www.mooseframework.org
ComputeFullJacobianThread.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* DO NOT MODIFY THIS HEADER */
3 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
4 /* */
5 /* (c) 2010 Battelle Energy Alliance, LLC */
6 /* ALL RIGHTS RESERVED */
7 /* */
8 /* Prepared by Battelle Energy Alliance, LLC */
9 /* Under Contract No. DE-AC07-05ID14517 */
10 /* With the U. S. Department of Energy */
11 /* */
12 /* See COPYRIGHT for full restrictions */
13 /****************************************************************/
14 
16 #include "NonlinearSystem.h"
17 #include "FEProblem.h"
18 #include "KernelBase.h"
19 #include "IntegratedBC.h"
20 #include "DGKernel.h"
21 #include "InterfaceKernel.h"
22 #include "MooseVariable.h"
23 #include "MooseVariableScalar.h"
24 #include "NonlocalKernel.h"
25 #include "NonlocalIntegratedBC.h"
26 #include "libmesh/threads.h"
27 
29  SparseMatrix<Number> & jacobian,
30  Moose::KernelType kernel_type)
31  : ComputeJacobianThread(fe_problem, jacobian),
32  _nl(fe_problem.getNonlinearSystemBase()),
33  _integrated_bcs(_nl.getIntegratedBCWarehouse()),
34  _dg_kernels(_nl.getDGKernelWarehouse()),
35  _interface_kernels(_nl.getInterfaceKernelWarehouse()),
36  _kernel_type(kernel_type),
37  _warehouse(NULL)
38 {
39  switch (_kernel_type)
40  {
41  case Moose::KT_ALL:
43  break;
44 
45  case Moose::KT_TIME:
47  break;
48 
49  case Moose::KT_NONTIME:
51  break;
52 
53  case Moose::KT_EIGEN:
55  break;
56 
57  case Moose::KT_NONEIGEN:
59  break;
60 
61  default:
62  mooseError("Unknown kernel type \n");
63  }
64 }
65 
66 // Splitting Constructor
68  Threads::split split)
69  : ComputeJacobianThread(x, split),
70  _nl(x._nl),
76 {
77 }
78 
80 
81 void
83 {
84  std::vector<std::pair<MooseVariable *, MooseVariable *>> & ce = _fe_problem.couplingEntries(_tid);
85  for (const auto & it : ce)
86  {
87  MooseVariable & ivariable = *(it.first);
88  MooseVariable & jvariable = *(it.second);
89 
90  unsigned int ivar = ivariable.number();
91  unsigned int jvar = jvariable.number();
92 
93  if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) &&
95  {
96  // only if there are dofs for j-variable (if it is subdomain restricted var, there may not be
97  // any)
98  const std::vector<std::shared_ptr<KernelBase>> & kernels =
100  for (const auto & kernel : kernels)
101  if ((kernel->variable().number() == ivar) && kernel->isImplicit())
102  {
103  kernel->subProblem().prepareShapes(jvar, _tid);
104  kernel->computeOffDiagJacobian(jvar);
105  }
106  }
107  }
108 
111  {
112  std::vector<std::pair<MooseVariable *, MooseVariable *>> & cne =
114  for (const auto & it : cne)
115  {
116  MooseVariable & ivariable = *(it.first);
117  MooseVariable & jvariable = *(it.second);
118 
119  unsigned int ivar = ivariable.number();
120  unsigned int jvar = jvariable.number();
121 
122  if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) &&
124  {
125  const std::vector<std::shared_ptr<KernelBase>> & kernels =
127  for (const auto & kernel : kernels)
128  {
129  std::shared_ptr<NonlocalKernel> nonlocal_kernel =
130  std::dynamic_pointer_cast<NonlocalKernel>(kernel);
131  if (nonlocal_kernel)
132  if ((kernel->variable().number() == ivar) && kernel->isImplicit())
133  {
134  kernel->subProblem().prepareShapes(jvar, _tid);
135  kernel->computeNonlocalOffDiagJacobian(jvar);
136  }
137  }
138  }
139  }
140  }
141 
142  const std::vector<MooseVariableScalar *> & scalar_vars = _nl.getScalarVariables(_tid);
143  if (scalar_vars.size() > 0)
144  {
145  // go over nl-variables (non-scalar)
146  const std::vector<MooseVariable *> & vars = _nl.getVariables(_tid);
147  for (const auto & ivariable : vars)
148  if (ivariable->activeOnSubdomain(_subdomain) > 0 &&
150  {
151  // for each variable get the list of active kernels
152  const std::vector<std::shared_ptr<KernelBase>> & kernels =
154  for (const auto & kernel : kernels)
155  if (kernel->isImplicit())
156  {
157  // now, get the list of coupled scalar vars and compute their off-diag jacobians
158  const std::vector<MooseVariableScalar *> coupled_scalar_vars =
159  kernel->getCoupledMooseScalarVars();
160 
161  // Do: dvar / dscalar_var, only want to process only nl-variables (not aux ones)
162  for (const auto & jvariable : coupled_scalar_vars)
163  if (_nl.hasScalarVariable(jvariable->name()))
164  kernel->computeOffDiagJacobianScalar(jvariable->number());
165  }
166  }
167  }
168 }
169 
170 void
172 {
173  std::vector<std::pair<MooseVariable *, MooseVariable *>> & ce = _fe_problem.couplingEntries(_tid);
174  for (const auto & it : ce)
175  {
176  MooseVariable & ivar = *(it.first);
177  MooseVariable & jvar = *(it.second);
180  {
181  // only if there are dofs for j-variable (if it is subdomain restricted var, there may not be
182  // any)
183  const std::vector<std::shared_ptr<IntegratedBC>> & bcs =
185  for (const auto & bc : bcs)
186  if (bc->shouldApply() && bc->variable().number() == ivar.number() && bc->isImplicit())
187  {
188  bc->subProblem().prepareFaceShapes(jvar.number(), _tid);
189  bc->computeJacobianBlock(jvar.number());
190  }
191  }
192  }
193 
196  {
197  std::vector<std::pair<MooseVariable *, MooseVariable *>> & cne =
199  for (const auto & it : cne)
200  {
201  MooseVariable & ivariable = *(it.first);
202  MooseVariable & jvariable = *(it.second);
203 
204  unsigned int ivar = ivariable.number();
205  unsigned int jvar = jvariable.number();
206 
207  if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) &&
209  {
210  const std::vector<std::shared_ptr<IntegratedBC>> & integrated_bcs =
212  for (const auto & integrated_bc : integrated_bcs)
213  {
214  std::shared_ptr<NonlocalIntegratedBC> nonlocal_integrated_bc =
215  std::dynamic_pointer_cast<NonlocalIntegratedBC>(integrated_bc);
216  if (nonlocal_integrated_bc)
217  if ((integrated_bc->variable().number() == ivar) && integrated_bc->isImplicit())
218  {
219  integrated_bc->subProblem().prepareFaceShapes(jvar, _tid);
220  integrated_bc->computeNonlocalOffDiagJacobian(jvar);
221  }
222  }
223  }
224  }
225  }
226 
227  const std::vector<MooseVariableScalar *> & scalar_vars = _nl.getScalarVariables(_tid);
228  if (scalar_vars.size() > 0)
229  {
230  // go over nl-variables (non-scalar)
231  const std::vector<MooseVariable *> & vars = _nl.getVariables(_tid);
232  for (const auto & ivar : vars)
233  if (ivar->activeOnSubdomain(_subdomain) > 0 &&
235  {
236  // for each variable get the list of active kernels
237  const std::vector<std::shared_ptr<IntegratedBC>> & bcs =
239  for (const auto & bc : bcs)
240  if (bc->variable().number() == ivar->number() && bc->isImplicit())
241  {
242  // now, get the list of coupled scalar vars and compute their off-diag jacobians
243  const std::vector<MooseVariableScalar *> coupled_scalar_vars =
244  bc->getCoupledMooseScalarVars();
245 
246  // Do: dvar / dscalar_var, only want to process only nl-variables (not aux ones)
247  for (const auto & jvar : coupled_scalar_vars)
248  if (_nl.hasScalarVariable(jvar->name()))
249  bc->computeJacobianBlockScalar(jvar->number());
250  }
251  }
252  }
253 }
254 
255 void
257 {
259  {
260  const auto & ce = _fe_problem.couplingEntries(_tid);
261  for (const auto & it : ce)
262  {
263  const std::vector<std::shared_ptr<DGKernel>> & dgks =
265  for (const auto & dg : dgks)
266  {
267  MooseVariable & ivariable = *(it.first);
268  MooseVariable & jvariable = *(it.second);
269 
270  unsigned int ivar = ivariable.number();
271  unsigned int jvar = jvariable.number();
272 
273  if (dg->variable().number() == ivar && dg->isImplicit() &&
274  dg->hasBlocks(neighbor->subdomain_id()) && jvariable.activeOnSubdomain(_subdomain))
275  {
276  dg->subProblem().prepareFaceShapes(jvar, _tid);
277  dg->subProblem().prepareNeighborShapes(jvar, _tid);
278  dg->computeOffDiagJacobian(jvar);
279  }
280  }
281  }
282  }
283 }
284 
285 void
287 {
289  {
290  const auto & ce = _fe_problem.couplingEntries(_tid);
291  for (const auto & it : ce)
292  {
293  const std::vector<std::shared_ptr<InterfaceKernel>> & int_ks =
295  for (const auto & interface_kernel : int_ks)
296  {
297  if (!interface_kernel->isImplicit())
298  continue;
299 
300  unsigned int ivar = it.first->number();
301  unsigned int jvar = it.second->number();
302 
303  interface_kernel->subProblem().prepareFaceShapes(jvar, _tid);
304  interface_kernel->subProblem().prepareNeighborShapes(jvar, _tid);
305 
306  if (interface_kernel->variable().number() == ivar)
307  interface_kernel->computeElementOffDiagJacobian(jvar);
308 
309  if (interface_kernel->neighborVariable().number() == ivar)
310  interface_kernel->computeNeighborOffDiagJacobian(jvar);
311  }
312  }
313  }
314 }
SubProblem & subProblem()
Get a reference to the subproblem.
const std::vector< std::shared_ptr< KernelBase > > & getActiveVariableBlockObjects(unsigned int variable_id, SubdomainID block_id, THREAD_ID tid=0) const
const std::vector< MooseVariableScalar * > & getScalarVariables(THREAD_ID tid)
Definition: SystemBase.h:422
virtual void prepareFaceShapes(unsigned int var, THREAD_ID tid)=0
const MooseObjectWarehouse< DGKernel > & _dg_kernels
Class for stuff related to variables.
Definition: MooseVariable.h:43
virtual void prepareShapes(unsigned int var, THREAD_ID tid)=0
virtual bool hasScalarVariable(const std::string &var_name)
Definition: SystemBase.C:565
virtual void computeInternalFaceJacobian(const Elem *neighbor) override
const KernelWarehouse & getNonEigenKernelWarehouse()
SubProblem & subProblem()
Returns a reference to the SubProblem for which this Kernel is active.
Definition: KernelBase.C:162
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getActiveBoundaryObjects(THREAD_ID tid=0) const
NonlocalIntegratedBC is used for solving integral terms in integro-differential equations.
void mooseError(Args &&...args)
Emit an error message with the given stringified, concatenated args and terminate the application...
Definition: MooseError.h:182
const KernelWarehouse & getEigenKernelWarehouse()
static PetscErrorCode Vec x
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getBoundaryObjects(THREAD_ID tid=0) const
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
const KernelWarehouse * _warehouse
const MooseObjectWarehouse< InterfaceKernel > & _interface_kernels
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
std::vector< std::pair< MooseVariable *, MooseVariable * > > & nonlocalCouplingEntries(THREAD_ID tid)
bool activeOnSubdomain(SubdomainID subdomain) const
Is the variable active on the subdomain?
virtual bool checkNonlocalCouplingRequirement()
Definition: SubProblem.h:62
const std::vector< MooseVariable * > & getVariables(THREAD_ID tid)
Definition: SystemBase.h:418
virtual void computeInternalInterFaceJacobian(BoundaryID bnd_id) override
KernelType
Definition: MooseTypes.h:162
const std::map< SubdomainID, std::vector< std::shared_ptr< T > > > & getActiveBlockObjects(THREAD_ID tid=0) const
virtual void computeFaceJacobian(BoundaryID bnd_id) override
ComputeFullJacobianThread(FEProblemBase &fe_problem, SparseMatrix< Number > &jacobian, Moose::KernelType kernel_type=Moose::KT_ALL)
std::vector< std::pair< MooseVariable *, MooseVariable * > > & couplingEntries(THREAD_ID tid)
unsigned int number() const
Get variable number coming from libMesh.
NonlocalKernel is used for solving integral terms in integro-differential equations.
virtual void computeJacobian() override
bool hasActiveBlockObjects(THREAD_ID tid=0) const
SubdomainID _subdomain
The subdomain for the current element.
const KernelWarehouse & getTimeKernelWarehouse()
bool hasActiveVariableBlockObjects(unsigned int variable_id, SubdomainID block_id, THREAD_ID tid=0) const
Methods for checking/getting variable kernels for a variable and SubdomainID.
const KernelWarehouse & getNonTimeKernelWarehouse()
const MooseObjectWarehouse< IntegratedBC > & _integrated_bcs
boundary_id_type BoundaryID
Definition: MooseTypes.h:75
const KernelWarehouse & getKernelWarehouse()
Access functions to Warehouses from outside NonlinearSystemBase.