www.mooseframework.org
ComputeJacobianThread.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 
15 #include "ComputeJacobianThread.h"
16 
17 #include "DGKernel.h"
18 #include "FEProblem.h"
19 #include "IntegratedBC.h"
20 #include "InterfaceKernel.h"
21 #include "KernelWarehouse.h"
22 #include "MooseVariable.h"
23 #include "NonlinearSystem.h"
24 #include "NonlocalIntegratedBC.h"
25 #include "NonlocalKernel.h"
26 #include "SwapBackSentinel.h"
27 #include "TimeDerivative.h"
28 
29 #include "libmesh/threads.h"
30 
32  SparseMatrix<Number> & jacobian,
33  Moose::KernelType kernel_type)
34  : ThreadedElementLoop<ConstElemRange>(fe_problem),
35  _jacobian(jacobian),
36  _nl(fe_problem.getNonlinearSystemBase()),
37  _num_cached(0),
38  _integrated_bcs(_nl.getIntegratedBCWarehouse()),
39  _dg_kernels(_nl.getDGKernelWarehouse()),
40  _interface_kernels(_nl.getInterfaceKernelWarehouse()),
41  _kernels(_nl.getKernelWarehouse()),
42  _kernel_type(kernel_type)
43 {
44 }
45 
46 // Splitting Constructor
48  : ThreadedElementLoop<ConstElemRange>(x, split),
50  _nl(x._nl),
55  _kernels(x._kernels),
57 {
58 }
59 
61 
62 void
64 {
65  const MooseObjectWarehouse<KernelBase> * warehouse;
66  switch (_kernel_type)
67  {
68  case Moose::KT_ALL:
69  warehouse = &_nl.getKernelWarehouse();
70  break;
71 
72  case Moose::KT_TIME:
73  warehouse = &_nl.getTimeKernelWarehouse();
74  break;
75 
76  case Moose::KT_NONTIME:
77  warehouse = &_nl.getNonTimeKernelWarehouse();
78  break;
79 
80  case Moose::KT_EIGEN:
81  warehouse = &_nl.getEigenKernelWarehouse();
82  break;
83 
84  case Moose::KT_NONEIGEN:
85  warehouse = &_nl.getNonEigenKernelWarehouse();
86  break;
87 
88  default:
89  mooseError("Unknown kernel type \n");
90  }
91 
92  if (warehouse->hasActiveBlockObjects(_subdomain, _tid))
93  {
94  const std::vector<std::shared_ptr<KernelBase>> & kernels =
96  for (const auto & kernel : kernels)
97  if (kernel->isImplicit())
98  {
99  kernel->subProblem().prepareShapes(kernel->variable().number(), _tid);
100  kernel->computeJacobian();
103  {
104  std::shared_ptr<NonlocalKernel> nonlocal_kernel =
105  std::dynamic_pointer_cast<NonlocalKernel>(kernel);
106  if (nonlocal_kernel)
107  kernel->computeNonlocalJacobian();
108  }
109  }
110  }
111 }
112 
113 void
115 {
116  const std::vector<std::shared_ptr<IntegratedBC>> & bcs =
118  for (const auto & bc : bcs)
119  if (bc->shouldApply() && bc->isImplicit())
120  {
121  bc->subProblem().prepareFaceShapes(bc->variable().number(), _tid);
122  bc->computeJacobian();
125  {
126  std::shared_ptr<NonlocalIntegratedBC> nonlocal_integrated_bc =
127  std::dynamic_pointer_cast<NonlocalIntegratedBC>(bc);
128  if (nonlocal_integrated_bc)
130  }
131  }
132 }
133 
134 void
136 {
137  // No need to call hasActiveObjects, this is done in the calling method (see onInternalSide)
138  const std::vector<std::shared_ptr<DGKernel>> & dgks =
140  for (const auto & dg : dgks)
141  if (dg->isImplicit())
142  {
143  dg->subProblem().prepareFaceShapes(dg->variable().number(), _tid);
144  dg->subProblem().prepareNeighborShapes(dg->variable().number(), _tid);
145  if (dg->hasBlocks(neighbor->subdomain_id()))
146  dg->computeJacobian();
147  }
148 }
149 
150 void
152 {
153  // No need to call hasActiveObjects, this is done in the calling method (see onInterface)
154  const std::vector<std::shared_ptr<InterfaceKernel>> & intks =
156  for (const auto & intk : intks)
157  if (intk->isImplicit())
158  {
159  intk->subProblem().prepareFaceShapes(intk->variable().number(), _tid);
160  intk->subProblem().prepareNeighborShapes(intk->neighborVariable().number(), _tid);
161  intk->computeJacobian();
162  }
163 }
164 
165 void
167 {
169 
170  // Update variable Dependencies
171  std::set<MooseVariable *> needed_moose_vars;
176 
177  // Update material dependencies
178  std::set<unsigned int> needed_mat_props;
183 
185  _fe_problem.setActiveMaterialProperties(needed_mat_props, _tid);
187 }
188 
189 void
191 {
192  _fe_problem.prepare(elem, _tid);
193 
194  _fe_problem.reinitElem(elem, _tid);
195 
196  // Set up Sentinel class so that, even if reinitMaterials() throws, we
197  // still remember to swap back during stack unwinding.
200 
201  if (_nl.getScalarVariables(_tid).size() > 0)
203 
204  computeJacobian();
205 }
206 
207 void
208 ComputeJacobianThread::onBoundary(const Elem * elem, unsigned int side, BoundaryID bnd_id)
209 {
211  {
212  _fe_problem.reinitElemFace(elem, side, bnd_id, _tid);
213 
214  // Set up Sentinel class so that, even if reinitMaterials() throws, we
215  // still remember to swap back during stack unwinding.
217 
218  _fe_problem.reinitMaterialsFace(elem->subdomain_id(), _tid);
220 
221  computeFaceJacobian(bnd_id);
222  }
223 }
224 
225 void
226 ComputeJacobianThread::onInternalSide(const Elem * elem, unsigned int side)
227 {
229  {
230  // Pointer to the neighbor we are currently working on.
231  const Elem * neighbor = elem->neighbor_ptr(side);
232 
233  // Get the global id of the element and the neighbor
234  const dof_id_type elem_id = elem->id(), neighbor_id = neighbor->id();
235 
236  if ((neighbor->active() && (neighbor->level() == elem->level()) && (elem_id < neighbor_id)) ||
237  (neighbor->level() < elem->level()))
238  {
239  _fe_problem.reinitNeighbor(elem, side, _tid);
240 
241  // Set up Sentinels so that, even if one of the reinitMaterialsXXX() calls throws, we
242  // still remember to swap back during stack unwinding.
244  _fe_problem.reinitMaterialsFace(elem->subdomain_id(), _tid);
245 
247  _fe_problem.reinitMaterialsNeighbor(neighbor->subdomain_id(), _tid);
248 
249  computeInternalFaceJacobian(neighbor);
250 
251  {
252  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
254  }
255  }
256  }
257 }
258 
259 void
260 ComputeJacobianThread::onInterface(const Elem * elem, unsigned int side, BoundaryID bnd_id)
261 {
263  {
264  // Pointer to the neighbor we are currently working on.
265  const Elem * neighbor = elem->neighbor_ptr(side);
266 
267  if (neighbor->active())
268  {
269  _fe_problem.reinitNeighbor(elem, side, _tid);
270 
271  // Set up Sentinels so that, even if one of the reinitMaterialsXXX() calls throws, we
272  // still remember to swap back during stack unwinding.
274  _fe_problem.reinitMaterialsFace(elem->subdomain_id(), _tid);
275 
277  _fe_problem.reinitMaterialsNeighbor(neighbor->subdomain_id(), _tid);
278 
280 
281  {
282  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
284  }
285  }
286  }
287 }
288 
289 void
291 {
293  _num_cached++;
294 
295  if (_num_cached % 20 == 0)
296  {
297  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
299  }
300 }
301 
302 void
304 {
307 }
308 
309 void
311 {
312 }
virtual void post() override
Called after the element range loop.
void updateBlockVariableDependency(SubdomainID id, std::set< MooseVariable * > &needed_moose_vars, THREAD_ID tid=0) const
virtual void setActiveElementalMooseVariables(const std::set< MooseVariable * > &moose_vars, THREAD_ID tid) override
Set the MOOSE variables to be reinited on each element.
Moose::KernelType _kernel_type
virtual void postElement(const Elem *) override
Called after the element assembly is done (including surface assembling)
virtual void prepare(const Elem *elem, THREAD_ID tid) override
Base class for assembly-like calculations.
SparseMatrix< Number > & _jacobian
const std::vector< MooseVariableScalar * > & getScalarVariables(THREAD_ID tid)
Definition: SystemBase.h:422
virtual void onElement(const Elem *elem) override
Assembly of the element (not including surface assembly)
virtual void prepareMaterials(SubdomainID blk_id, THREAD_ID tid)
Add the MooseVariables that the current materials depend on to the dependency list.
void updateBlockMatPropDependency(SubdomainID id, std::set< unsigned int > &needed_mat_props, THREAD_ID tid=0) const
virtual void computeNonlocalJacobian()
computeNonlocalJacobian and computeNonlocalOffDiagJacobian methods are introduced for providing the j...
const KernelWarehouse & getNonEigenKernelWarehouse()
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.
virtual void reinitMaterialsNeighbor(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
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
virtual void reinitElem(const Elem *elem, THREAD_ID tid) override
virtual void reinitMaterials(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void reinitMaterialsBoundary(BoundaryID boundary_id, THREAD_ID tid, bool swap_stateful=true)
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
virtual void setActiveMaterialProperties(const std::set< unsigned int > &mat_prop_ids, THREAD_ID tid) override
Record and set the material properties required by the current computing thread.
virtual void subdomainSetup(SubdomainID subdomain, THREAD_ID tid)
void updateBoundaryVariableDependency(std::set< MooseVariable * > &needed_moose_vars, THREAD_ID tid=0) const
ComputeJacobianThread(FEProblemBase &fe_problem, SparseMatrix< Number > &jacobian, Moose::KernelType kernel_type=Moose::KT_ALL)
virtual void cacheJacobian(THREAD_ID tid) override
void updateBoundaryMatPropDependency(std::set< unsigned int > &needed_mat_props, THREAD_ID tid=0) const
virtual void computeInternalInterFaceJacobian(BoundaryID bnd_id)
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
virtual void clearActiveMaterialProperties(THREAD_ID tid) override
Clear the active material properties.
virtual bool checkNonlocalCouplingRequirement()
Definition: SubProblem.h:62
virtual void swapBackMaterials(THREAD_ID tid)
virtual void addCachedJacobian(SparseMatrix< Number > &jacobian, THREAD_ID tid) override
KernelType
Definition: MooseTypes.h:162
NonlinearSystemBase & _nl
virtual void swapBackMaterialsNeighbor(THREAD_ID tid)
const std::map< SubdomainID, std::vector< std::shared_ptr< T > > > & getActiveBlockObjects(THREAD_ID tid=0) const
virtual void subdomainChanged() override
Called every time the current subdomain changes (i.e.
const MooseObjectWarehouse< DGKernel > & _dg_kernels
virtual void reinitElemFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid) override
virtual void onBoundary(const Elem *elem, unsigned int side, BoundaryID bnd_id) override
Called when doing boundary assembling.
virtual void reinitMaterialsFace(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void computeFaceJacobian(BoundaryID bnd_id)
void join(const ComputeJacobianThread &)
NonlocalKernel is used for solving integral terms in integro-differential equations.
virtual void reinitNeighbor(const Elem *elem, unsigned int side, THREAD_ID tid) override
bool hasActiveBlockObjects(THREAD_ID tid=0) const
virtual void computeInternalFaceJacobian(const Elem *neighbor)
virtual void addJacobianNeighbor(SparseMatrix< Number > &jacobian, THREAD_ID tid) override
SubdomainID _subdomain
The subdomain for the current element.
const KernelWarehouse & getTimeKernelWarehouse()
virtual void computeNonlocalJacobian()
computeNonlocalJacobian and computeNonlocalOffDiagJacobian methods are introduced for providing the j...
const KernelWarehouse & _kernels
virtual void reinitOffDiagScalars(THREAD_ID tid) override
const KernelWarehouse & getNonTimeKernelWarehouse()
const MooseObjectWarehouse< InterfaceKernel > & _interface_kernels
virtual void onInterface(const Elem *elem, unsigned int side, BoundaryID bnd_id) override
Called when doing interface assembling.
const MooseObjectWarehouse< IntegratedBC > & _integrated_bcs
virtual void swapBackMaterialsFace(THREAD_ID tid)
virtual void onInternalSide(const Elem *elem, unsigned int side) override
Called when doing internal edge assembling.
virtual void clearActiveElementalMooseVariables(THREAD_ID tid) override
Clear the active elemental MooseVariable.
boundary_id_type BoundaryID
Definition: MooseTypes.h:75
const KernelWarehouse & getKernelWarehouse()
Access functions to Warehouses from outside NonlinearSystemBase.
The "SwapBackSentinel" class&#39;s destructor guarantees that FEProblemBase::swapBackMaterials{Face,Neighbor}() is called even when an exception is thrown from FEProblemBase::reinitMaterials{Face,Neighbor}.