www.mooseframework.org
ComputeFullJacobianThread.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://www.mooseframework.org
3 //*
4 //* All rights reserved, see COPYRIGHT for full restrictions
5 //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 //*
7 //* Licensed under LGPL 2.1, please see LICENSE for details
8 //* https://www.gnu.org/licenses/lgpl-2.1.html
9 
11 #include "NonlinearSystem.h"
12 #include "FEProblem.h"
13 #include "KernelBase.h"
14 #include "IntegratedBCBase.h"
15 #include "DGKernel.h"
16 #include "InterfaceKernel.h"
17 #include "MooseVariableFE.h"
18 #include "MooseVariableScalar.h"
19 #include "NonlocalKernel.h"
20 #include "NonlocalIntegratedBC.h"
21 #include "libmesh/threads.h"
22 
24  const std::set<TagID> & tags)
25  : ComputeJacobianThread(fe_problem, tags)
26 {
27 }
28 
29 // Splitting Constructor
33 {
34 }
35 
37 
38 void
40 {
41  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & ce =
43  for (const auto & it : ce)
44  {
45  MooseVariableFEBase & ivariable = *(it.first);
46  MooseVariableFEBase & jvariable = *(it.second);
47 
48  unsigned int ivar = ivariable.number();
49  unsigned int jvar = jvariable.number();
50 
51  if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) &&
53  {
54  // only if there are dofs for j-variable (if it is subdomain restricted var, there may not be
55  // any)
56  const auto & kernels = _warehouse->getActiveVariableBlockObjects(ivar, _subdomain, _tid);
57  for (const auto & kernel : kernels)
58  if ((kernel->variable().number() == ivar) && kernel->isImplicit())
59  {
60  kernel->subProblem().prepareShapes(jvar, _tid);
61  kernel->computeOffDiagJacobian(jvariable);
62  }
63  }
64  }
65 
67  {
69  for (const auto & kernel : kernels)
70  if (kernel->isImplicit())
71  kernel->computeADOffDiagJacobian();
72  }
73 
76  {
77  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & cne =
79  for (const auto & it : cne)
80  {
81  MooseVariableFEBase & ivariable = *(it.first);
82  MooseVariableFEBase & jvariable = *(it.second);
83 
84  unsigned int ivar = ivariable.number();
85  unsigned int jvar = jvariable.number();
86 
87  if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) &&
89  {
90  const auto & kernels = _warehouse->getActiveVariableBlockObjects(ivar, _subdomain, _tid);
91  for (const auto & kernel : kernels)
92  {
93  std::shared_ptr<NonlocalKernel> nonlocal_kernel =
95  if (nonlocal_kernel)
96  if ((kernel->variable().number() == ivar) && kernel->isImplicit())
97  {
98  kernel->subProblem().prepareShapes(jvar, _tid);
99  kernel->computeNonlocalOffDiagJacobian(jvar);
100  }
101  }
102  }
103  }
104  }
105 
106  const std::vector<MooseVariableScalar *> & scalar_vars = _nl.getScalarVariables(_tid);
107  if (scalar_vars.size() > 0)
108  {
109  // go over nl-variables (non-scalar)
110  const std::vector<MooseVariableFEBase *> & vars = _nl.getVariables(_tid);
111  for (const auto & ivariable : vars)
112  if (ivariable->activeOnSubdomain(_subdomain) > 0 &&
114  {
115  // for each variable get the list of active kernels
116  const auto & kernels =
118  for (const auto & kernel : kernels)
119  if (kernel->isImplicit())
120  {
121  // now, get the list of coupled scalar vars and compute their off-diag jacobians
122  const std::vector<MooseVariableScalar *> coupled_scalar_vars =
123  kernel->getCoupledMooseScalarVars();
124 
125  // Do: dvar / dscalar_var, only want to process only nl-variables (not aux ones)
126  for (const auto & jvariable : coupled_scalar_vars)
127  if (_nl.hasScalarVariable(jvariable->name()))
128  kernel->computeOffDiagJacobianScalar(jvariable->number());
129  }
130  }
131  }
132 }
133 
134 void
136 {
137  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & ce =
139  for (const auto & it : ce)
140  {
141  MooseVariableFEBase & ivar = *(it.first);
142  MooseVariableFEBase & jvar = *(it.second);
145  {
146  // only if there are dofs for j-variable (if it is subdomain restricted var, there may not be
147  // any)
148  const auto & bcs = _ibc_warehouse->getActiveBoundaryObjects(bnd_id, _tid);
149  for (const auto & bc : bcs)
150  if (bc->shouldApply() && bc->variable().number() == ivar.number() && bc->isImplicit())
151  {
152  bc->subProblem().prepareFaceShapes(jvar.number(), _tid);
153  bc->computeJacobianBlock(jvar);
154  }
155  }
156  }
157 
160  {
161  std::vector<std::pair<MooseVariableFEBase *, MooseVariableFEBase *>> & cne =
163  for (const auto & it : cne)
164  {
165  MooseVariableFEBase & ivariable = *(it.first);
166  MooseVariableFEBase & jvariable = *(it.second);
167 
168  unsigned int ivar = ivariable.number();
169  unsigned int jvar = jvariable.number();
170 
171  if (ivariable.activeOnSubdomain(_subdomain) && jvariable.activeOnSubdomain(_subdomain) &&
173  {
174  const std::vector<std::shared_ptr<IntegratedBCBase>> & integrated_bcs =
176  for (const auto & integrated_bc : integrated_bcs)
177  {
178  std::shared_ptr<NonlocalIntegratedBC> nonlocal_integrated_bc =
180  if (nonlocal_integrated_bc)
181  if ((integrated_bc->variable().number() == ivar) && integrated_bc->isImplicit())
182  {
183  integrated_bc->subProblem().prepareFaceShapes(jvar, _tid);
184  integrated_bc->computeNonlocalOffDiagJacobian(jvar);
185  }
186  }
187  }
188  }
189  }
190 
191  const std::vector<MooseVariableScalar *> & scalar_vars = _nl.getScalarVariables(_tid);
192  if (scalar_vars.size() > 0)
193  {
194  // go over nl-variables (non-scalar)
195  const std::vector<MooseVariableFEBase *> & vars = _nl.getVariables(_tid);
196  for (const auto & ivar : vars)
197  if (ivar->activeOnSubdomain(_subdomain) > 0 &&
199  {
200  // for each variable get the list of active kernels
201  const auto & bcs = _ibc_warehouse->getActiveBoundaryObjects(bnd_id, _tid);
202  for (const auto & bc : bcs)
203  if (bc->variable().number() == ivar->number() && bc->isImplicit())
204  {
205  // now, get the list of coupled scalar vars and compute their off-diag jacobians
206  const std::vector<MooseVariableScalar *> coupled_scalar_vars =
207  bc->getCoupledMooseScalarVars();
208 
209  // Do: dvar / dscalar_var, only want to process only nl-variables (not aux ones)
210  for (const auto & jvar : coupled_scalar_vars)
211  if (_nl.hasScalarVariable(jvar->name()))
212  bc->computeJacobianBlockScalar(jvar->number());
213  }
214  }
215  }
216 }
217 
218 void
220 {
222  {
223  const auto & ce = _fe_problem.couplingEntries(_tid);
224  for (const auto & it : ce)
225  {
226  const auto & dgks = _dg_warehouse->getActiveBlockObjects(_subdomain, _tid);
227  for (const auto & dg : dgks)
228  {
229  MooseVariableFEBase & ivariable = *(it.first);
230  MooseVariableFEBase & jvariable = *(it.second);
231 
232  unsigned int ivar = ivariable.number();
233  unsigned int jvar = jvariable.number();
234 
235  if (dg->variable().number() == ivar && dg->isImplicit() &&
236  dg->hasBlocks(neighbor->subdomain_id()) && jvariable.activeOnSubdomain(_subdomain))
237  {
238  dg->subProblem().prepareFaceShapes(jvar, _tid);
239  dg->subProblem().prepareNeighborShapes(jvar, _tid);
240  dg->computeOffDiagJacobian(jvar);
241  }
242  }
243  }
244  }
245 }
246 
247 void
249 {
251  {
252  const auto & ce = _fe_problem.couplingEntries(_tid);
253  for (const auto & it : ce)
254  {
255  const std::vector<std::shared_ptr<InterfaceKernel>> & int_ks =
257  for (const auto & interface_kernel : int_ks)
258  {
259  if (!interface_kernel->isImplicit())
260  continue;
261 
262  unsigned int ivar = it.first->number();
263  unsigned int jvar = it.second->number();
264 
265  interface_kernel->subProblem().prepareFaceShapes(jvar, _tid);
266  interface_kernel->subProblem().prepareNeighborShapes(jvar, _tid);
267 
268  if (interface_kernel->variable().number() == ivar)
269  interface_kernel->computeElementOffDiagJacobian(jvar);
270 
271  if (interface_kernel->neighborVariable().number() == ivar)
272  interface_kernel->computeNeighborOffDiagJacobian(jvar);
273  }
274  }
275  }
276 }
const std::vector< std::shared_ptr< T > > & getActiveVariableBlockObjects(unsigned int variable_id, SubdomainID block_id, THREAD_ID tid=0) const
SubProblem & subProblem()
Get a reference to the subproblem.
const std::vector< MooseVariableScalar * > & getScalarVariables(THREAD_ID tid)
Definition: SystemBase.h:614
bool hasActiveBlockObjects(THREAD_ID tid=0) const
virtual void prepareFaceShapes(unsigned int var, THREAD_ID tid)=0
const std::map< SubdomainID, std::vector< std::shared_ptr< T > > > & getActiveBlockObjects(THREAD_ID tid=0) const
unsigned int number() const
Get variable number coming from libMesh.
virtual void prepareShapes(unsigned int var, THREAD_ID tid)=0
ComputeFullJacobianThread(FEProblemBase &fe_problem, const std::set< TagID > &tags)
virtual void computeInternalFaceJacobian(const Elem *neighbor) override
NonlocalIntegratedBC is used for solving integral terms in integro-differential equations.
virtual bool activeOnSubdomain(SubdomainID subdomain) const =0
Is the variable active on the subdomain?
const std::vector< MooseVariableFEBase * > & getVariables(THREAD_ID tid)
Definition: SystemBase.h:609
std::unique_ptr< T_DEST, T_DELETER > dynamic_pointer_cast(std::unique_ptr< T_SRC, T_DELETER > &src)
These are reworked from https://stackoverflow.com/a/11003103.
static PetscErrorCode Vec x
MooseObjectWarehouse< DGKernelBase > * _dg_warehouse
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & couplingEntries(THREAD_ID tid)
MooseObjectWarehouse< InterfaceKernel > * _ik_warehouse
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
std::vector< std::string > split(const std::string &str, const std::string &delimiter)
Python like split function for strings.
Definition: MooseUtils.C:736
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.
MooseObjectWarehouse< KernelBase > * _adjk_warehouse
Pointer to tags.
boundary_id_type BoundaryID
virtual bool checkNonlocalCouplingRequirement()
Definition: SubProblem.h:69
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getBoundaryObjects(THREAD_ID tid=0) const
virtual void computeInternalInterFaceJacobian(BoundaryID bnd_id) override
MooseObjectWarehouse< KernelBase > * _warehouse
SubProblem & subProblem()
Returns a reference to the SubProblem for which this Kernel is active.
Definition: KernelBase.h:105
MooseObjectWarehouse< IntegratedBCBase > * _ibc_warehouse
NonlinearSystemBase & _nl
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getActiveBoundaryObjects(THREAD_ID tid=0) const
virtual void computeFaceJacobian(BoundaryID bnd_id) override
std::vector< std::pair< MooseVariableFEBase *, MooseVariableFEBase * > > & nonlocalCouplingEntries(THREAD_ID tid)
NonlocalKernel is used for solving integral terms in integro-differential equations.
virtual void computeJacobian() override
SubdomainID _subdomain
The subdomain for the current element.
virtual bool hasScalarVariable(const std::string &var_name) const
Definition: SystemBase.C:677