www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ComputeJacobianBlocksThread Class Reference

Specialization for filling multiple "small" preconditioning matrices simulatenously. More...

#include <ComputeJacobianBlocksThread.h>

Inheritance diagram for ComputeJacobianBlocksThread:
[legend]

Public Member Functions

 ComputeJacobianBlocksThread (FEProblemBase &fe_problem, std::vector< JacobianBlock * > &blocks)
 
 ComputeJacobianBlocksThread (ComputeJacobianBlocksThread &x, Threads::split split)
 
virtual ~ComputeJacobianBlocksThread ()
 
void join (const ComputeJacobianThread &)
 
virtual void subdomainChanged () override
 Called every time the current subdomain changes (i.e. More...
 
virtual void onElement (const Elem *elem) override
 Assembly of the element (not including surface assembly) More...
 
virtual void onBoundary (const Elem *elem, unsigned int side, BoundaryID bnd_id) override
 Called when doing boundary assembling. More...
 
virtual void onInternalSide (const Elem *elem, unsigned int side) override
 Called when doing internal edge assembling. More...
 
virtual void onInterface (const Elem *elem, unsigned int side, BoundaryID bnd_id) override
 Called when doing interface assembling. More...
 
virtual void post () override
 Called after the element range loop. More...
 
virtual void caughtMooseException (MooseException &e) override
 Called if a MooseException is caught anywhere during the computation. More...
 
virtual bool keepGoing () override
 Whether or not the loop should continue. More...
 
virtual void preElement (const Elem *elem) override
 Called before the element assembly. More...
 
virtual void preInternalSide (const Elem *elem, unsigned int side) override
 Called before evaluations on an element internal side. More...
 
virtual void neighborSubdomainChanged () override
 Called every time the neighbor subdomain changes (i.e. More...
 
void operator() (const ConstElemRange &range, bool bypass_threading=false)
 
virtual void pre ()
 Called before the element range loop. More...
 
virtual void postInternalSide (const Elem *elem, unsigned int side)
 Called after evaluations on an element internal side. More...
 

Protected Member Functions

virtual void postElement (const Elem *elem) override
 Called after the element assembly is done (including surface assembling) More...
 
virtual void computeJacobian () override
 
virtual void computeFaceJacobian (BoundaryID bnd_id) override
 
virtual void computeInternalFaceJacobian (const Elem *neighbor) override
 
virtual void computeInternalInterFaceJacobian (BoundaryID bnd_id) override
 

Protected Attributes

std::vector< JacobianBlock * > _blocks
 
NonlinearSystemBase_nl
 
const MooseObjectWarehouse< IntegratedBC > & _integrated_bcs
 
const MooseObjectWarehouse< DGKernel > & _dg_kernels
 
const MooseObjectWarehouse< InterfaceKernel > & _interface_kernels
 
Moose::KernelType _kernel_type
 
const KernelWarehouse_warehouse
 
SparseMatrix< Number > & _jacobian
 
unsigned int _num_cached
 
const KernelWarehouse_kernels
 
FEProblemBase_fe_problem
 
MooseMesh_mesh
 
THREAD_ID _tid
 
SubdomainID _subdomain
 The subdomain for the current element. More...
 
SubdomainID _old_subdomain
 The subdomain for the last element. More...
 
SubdomainID _neighbor_subdomain
 The subdomain for the current neighbor. More...
 
SubdomainID _old_neighbor_subdomain
 The subdomain for the last neighbor. More...
 

Detailed Description

Specialization for filling multiple "small" preconditioning matrices simulatenously.

Definition at line 45 of file ComputeJacobianBlocksThread.h.

Constructor & Destructor Documentation

ComputeJacobianBlocksThread::ComputeJacobianBlocksThread ( FEProblemBase fe_problem,
std::vector< JacobianBlock * > &  blocks 
)

Definition at line 26 of file ComputeJacobianBlocksThread.C.

28  : ComputeFullJacobianThread(fe_problem, blocks[0]->_jacobian /* have to pass something */),
29  _blocks(blocks)
30 {
31 }
SparseMatrix< Number > & _jacobian
std::vector< JacobianBlock * > _blocks
ComputeFullJacobianThread(FEProblemBase &fe_problem, SparseMatrix< Number > &jacobian, Moose::KernelType kernel_type=Moose::KT_ALL)
ComputeJacobianBlocksThread::ComputeJacobianBlocksThread ( ComputeJacobianBlocksThread x,
Threads::split  split 
)

Definition at line 34 of file ComputeJacobianBlocksThread.C.

37 {
38 }
std::vector< JacobianBlock * > _blocks
ComputeFullJacobianThread(FEProblemBase &fe_problem, SparseMatrix< Number > &jacobian, Moose::KernelType kernel_type=Moose::KT_ALL)
ComputeJacobianBlocksThread::~ComputeJacobianBlocksThread ( )
virtual

Definition at line 40 of file ComputeJacobianBlocksThread.C.

40 {}

Member Function Documentation

virtual void ThreadedElementLoop< ConstElemRange >::caughtMooseException ( MooseException e)
overridevirtualinherited

Called if a MooseException is caught anywhere during the computation.

The single input parameter taken is a MooseException object.

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

void ComputeFullJacobianThread::computeFaceJacobian ( BoundaryID  bnd_id)
overrideprotectedvirtualinherited

done only when nonlocal integrated_bcs exist in the system

Reimplemented from ComputeJacobianThread.

Definition at line 171 of file ComputeFullJacobianThread.C.

Referenced by ComputeFullJacobianThread::join().

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 }
SubProblem & subProblem()
Get a reference to the subproblem.
const std::vector< MooseVariableScalar * > & getScalarVariables(THREAD_ID tid)
Definition: SystemBase.h:422
virtual void prepareFaceShapes(unsigned int var, THREAD_ID tid)=0
Class for stuff related to variables.
Definition: MooseVariable.h:43
virtual bool hasScalarVariable(const std::string &var_name)
Definition: SystemBase.C:565
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.
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getBoundaryObjects(THREAD_ID tid=0) const
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
std::vector< std::pair< MooseVariable *, MooseVariable * > > & couplingEntries(THREAD_ID tid)
unsigned int number() const
Get variable number coming from libMesh.
SubdomainID _subdomain
The subdomain for the current element.
const MooseObjectWarehouse< IntegratedBC > & _integrated_bcs
void ComputeFullJacobianThread::computeInternalFaceJacobian ( const Elem *  neighbor)
overrideprotectedvirtualinherited

Reimplemented from ComputeJacobianThread.

Definition at line 256 of file ComputeFullJacobianThread.C.

Referenced by ComputeFullJacobianThread::join().

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 }
const MooseObjectWarehouse< DGKernel > & _dg_kernels
Class for stuff related to variables.
Definition: MooseVariable.h:43
bool activeOnSubdomain(SubdomainID subdomain) const
Is the variable active on the subdomain?
const std::map< SubdomainID, std::vector< std::shared_ptr< T > > > & getActiveBlockObjects(THREAD_ID tid=0) const
std::vector< std::pair< MooseVariable *, MooseVariable * > > & couplingEntries(THREAD_ID tid)
unsigned int number() const
Get variable number coming from libMesh.
bool hasActiveBlockObjects(THREAD_ID tid=0) const
SubdomainID _subdomain
The subdomain for the current element.
void ComputeFullJacobianThread::computeInternalInterFaceJacobian ( BoundaryID  bnd_id)
overrideprotectedvirtualinherited

Reimplemented from ComputeJacobianThread.

Definition at line 286 of file ComputeFullJacobianThread.C.

Referenced by ComputeFullJacobianThread::join().

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 }
const std::map< BoundaryID, std::vector< std::shared_ptr< T > > > & getActiveBoundaryObjects(THREAD_ID tid=0) const
const MooseObjectWarehouse< InterfaceKernel > & _interface_kernels
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
std::vector< std::pair< MooseVariable *, MooseVariable * > > & couplingEntries(THREAD_ID tid)
void ComputeFullJacobianThread::computeJacobian ( )
overrideprotectedvirtualinherited

done only when nonlocal kernels exist in the system

Reimplemented from ComputeJacobianThread.

Definition at line 82 of file ComputeFullJacobianThread.C.

Referenced by ComputeFullJacobianThread::join().

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 }
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
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
SubProblem & subProblem()
Returns a reference to the SubProblem for which this Kernel is active.
Definition: KernelBase.C:162
const std::string & name() const
Get the variable name.
const KernelWarehouse * _warehouse
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
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.
SubdomainID _subdomain
The subdomain for the current element.
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.
void ComputeJacobianBlocksThread::join ( const ComputeJacobianThread )
inline

Definition at line 55 of file ComputeJacobianBlocksThread.h.

55 {}
virtual bool ThreadedElementLoop< ConstElemRange >::keepGoing ( )
inlineoverridevirtualinherited

Whether or not the loop should continue.

Returns
true to keep going, false to stop.

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 50 of file ThreadedElementLoop.h.

50 { return !_fe_problem.hasException(); }
virtual bool hasException()
Whether or not an exception has occurred.
virtual void ThreadedElementLoop< ConstElemRange >::neighborSubdomainChanged ( )
overridevirtualinherited

Called every time the neighbor subdomain changes (i.e.

the subdomain of this neighbor is not the same as the subdomain of the last neighbor). Beware of over-using this! You might think that you can do some expensive stuff in here and get away with it... but there are applications that have TONS of subdomains....

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

void ComputeJacobianThread::onBoundary ( const Elem *  elem,
unsigned int  side,
BoundaryID  bnd_id 
)
overridevirtualinherited

Called when doing boundary assembling.

Parameters
elem- The element we are checking is on the boundary.
side- The side of the element in question.
bnd_id- ID of the boundary we are at

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 208 of file ComputeJacobianThread.C.

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 }
virtual void reinitMaterialsBoundary(BoundaryID boundary_id, THREAD_ID tid, bool swap_stateful=true)
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
virtual void reinitElemFace(const Elem *elem, unsigned int side, BoundaryID bnd_id, THREAD_ID tid) override
virtual void reinitMaterialsFace(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void computeFaceJacobian(BoundaryID bnd_id)
const MooseObjectWarehouse< IntegratedBC > & _integrated_bcs
virtual void swapBackMaterialsFace(THREAD_ID tid)
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}.
void ComputeJacobianThread::onElement ( const Elem *  elem)
overridevirtualinherited

Assembly of the element (not including surface assembly)

Parameters
elem- active element

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 190 of file ComputeJacobianThread.C.

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 }
virtual void prepare(const Elem *elem, THREAD_ID tid) override
const std::vector< MooseVariableScalar * > & getScalarVariables(THREAD_ID tid)
Definition: SystemBase.h:422
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 swapBackMaterials(THREAD_ID tid)
NonlinearSystemBase & _nl
SubdomainID _subdomain
The subdomain for the current element.
virtual void reinitOffDiagScalars(THREAD_ID tid) override
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}.
void ComputeJacobianThread::onInterface ( const Elem *  elem,
unsigned int  side,
BoundaryID  bnd_id 
)
overridevirtualinherited

Called when doing interface assembling.

Parameters
elem- Element we are on
side- local side number of the element 'elem'
bnd_id- ID of the interface we are at

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 260 of file ComputeJacobianThread.C.

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 }
SparseMatrix< Number > & _jacobian
virtual void reinitMaterialsNeighbor(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void computeInternalInterFaceJacobian(BoundaryID bnd_id)
bool hasActiveBoundaryObjects(THREAD_ID tid=0) const
virtual void swapBackMaterialsNeighbor(THREAD_ID tid)
virtual void reinitMaterialsFace(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void reinitNeighbor(const Elem *elem, unsigned int side, THREAD_ID tid) override
virtual void addJacobianNeighbor(SparseMatrix< Number > &jacobian, THREAD_ID tid) override
const MooseObjectWarehouse< InterfaceKernel > & _interface_kernels
virtual void swapBackMaterialsFace(THREAD_ID tid)
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}.
void ComputeJacobianThread::onInternalSide ( const Elem *  elem,
unsigned int  side 
)
overridevirtualinherited

Called when doing internal edge assembling.

Parameters
elem- Element we are on
side- local side number of the element 'elem'

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 226 of file ComputeJacobianThread.C.

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 }
SparseMatrix< Number > & _jacobian
virtual void reinitMaterialsNeighbor(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
virtual void swapBackMaterialsNeighbor(THREAD_ID tid)
const MooseObjectWarehouse< DGKernel > & _dg_kernels
virtual void reinitMaterialsFace(SubdomainID blk_id, THREAD_ID tid, bool swap_stateful=true)
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.
virtual void swapBackMaterialsFace(THREAD_ID tid)
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}.
void ThreadedElementLoopBase< ConstElemRange >::operator() ( const ConstElemRange &  range,
bool  bypass_threading = false 
)
inherited
void ComputeJacobianThread::post ( )
overridevirtualinherited

Called after the element range loop.

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 303 of file ComputeJacobianThread.C.

304 {
307 }
virtual void clearActiveMaterialProperties(THREAD_ID tid) override
Clear the active material properties.
virtual void clearActiveElementalMooseVariables(THREAD_ID tid) override
Clear the active elemental MooseVariable.
void ComputeJacobianBlocksThread::postElement ( const Elem *  elem)
overrideprotectedvirtual

Called after the element assembly is done (including surface assembling)

Parameters
elem- active element

Reimplemented from ComputeJacobianThread.

Definition at line 43 of file ComputeJacobianBlocksThread.C.

44 {
45  std::vector<dof_id_type>
46  dof_indices; // Do this out here to avoid creating and destroying it thousands of times
47 
48  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
49 
50  for (const auto & block : _blocks)
51  {
52  const DofMap & dof_map = block->_precond_system.get_dof_map();
53  dof_map.dof_indices(elem, dof_indices);
54 
56  block->_jacobian, block->_ivar, block->_jvar, dof_map, dof_indices, _tid);
57  }
58 }
DofMap & dof_map
std::vector< JacobianBlock * > _blocks
virtual void addJacobianBlock(SparseMatrix< Number > &jacobian, unsigned int ivar, unsigned int jvar, const DofMap &dof_map, std::vector< dof_id_type > &dof_indices, THREAD_ID tid) override
virtual void ThreadedElementLoopBase< ConstElemRange >::postInternalSide ( const Elem *  elem,
unsigned int  side 
)
virtualinherited

Called after evaluations on an element internal side.

Parameters
elem- Element we are on
side- local side number of the element 'elem'
virtual void ThreadedElementLoopBase< ConstElemRange >::pre ( )
virtualinherited

Called before the element range loop.

virtual void ThreadedElementLoop< ConstElemRange >::preElement ( const Elem *  elem)
overridevirtualinherited

Called before the element assembly.

Parameters
elem- active element

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

virtual void ThreadedElementLoop< ConstElemRange >::preInternalSide ( const Elem *  elem,
unsigned int  side 
)
overridevirtualinherited

Called before evaluations on an element internal side.

Parameters
elem- Element we are on
side- local side number of the element 'elem'

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

void ComputeJacobianThread::subdomainChanged ( )
overridevirtualinherited

Called every time the current subdomain changes (i.e.

the subdomain of this element is not the same as the subdomain of the last element). Beware of over-using this! You might think that you can do some expensive stuff in here and get away with it... but there are applications that have TONS of subdomains....

Reimplemented from ThreadedElementLoopBase< ConstElemRange >.

Definition at line 166 of file ComputeJacobianThread.C.

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 }
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.
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 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
void updateBoundaryMatPropDependency(std::set< unsigned int > &needed_mat_props, THREAD_ID tid=0) const
const MooseObjectWarehouse< DGKernel > & _dg_kernels
SubdomainID _subdomain
The subdomain for the current element.
const KernelWarehouse & _kernels
const MooseObjectWarehouse< InterfaceKernel > & _interface_kernels
const MooseObjectWarehouse< IntegratedBC > & _integrated_bcs

Member Data Documentation

std::vector<JacobianBlock *> ComputeJacobianBlocksThread::_blocks
protected

Definition at line 60 of file ComputeJacobianBlocksThread.h.

Referenced by postElement().

const MooseObjectWarehouse<DGKernel>& ComputeFullJacobianThread::_dg_kernels
protectedinherited
FEProblemBase& ThreadedElementLoop< ConstElemRange >::_fe_problem
protectedinherited
const MooseObjectWarehouse<IntegratedBC>& ComputeFullJacobianThread::_integrated_bcs
protectedinherited
const MooseObjectWarehouse<InterfaceKernel>& ComputeFullJacobianThread::_interface_kernels
protectedinherited
SparseMatrix<Number>& ComputeJacobianThread::_jacobian
protectedinherited
Moose::KernelType ComputeFullJacobianThread::_kernel_type
protectedinherited
const KernelWarehouse& ComputeJacobianThread::_kernels
protectedinherited

Definition at line 68 of file ComputeJacobianThread.h.

Referenced by ComputeJacobianThread::subdomainChanged().

MooseMesh& ThreadedElementLoopBase< ConstElemRange >::_mesh
protectedinherited

Definition at line 141 of file ThreadedElementLoopBase.h.

SubdomainID ThreadedElementLoopBase< ConstElemRange >::_neighbor_subdomain
protectedinherited

The subdomain for the current neighbor.

Definition at line 151 of file ThreadedElementLoopBase.h.

NonlinearSystemBase& ComputeFullJacobianThread::_nl
protectedinherited
unsigned int ComputeJacobianThread::_num_cached
protectedinherited

Definition at line 56 of file ComputeJacobianThread.h.

Referenced by ComputeJacobianThread::postElement().

SubdomainID ThreadedElementLoopBase< ConstElemRange >::_old_neighbor_subdomain
protectedinherited

The subdomain for the last neighbor.

Definition at line 154 of file ThreadedElementLoopBase.h.

SubdomainID ThreadedElementLoopBase< ConstElemRange >::_old_subdomain
protectedinherited

The subdomain for the last element.

Definition at line 148 of file ThreadedElementLoopBase.h.

SubdomainID ThreadedElementLoopBase< ConstElemRange >::_subdomain
protectedinherited
THREAD_ID ThreadedElementLoopBase< ConstElemRange >::_tid
protectedinherited

Definition at line 142 of file ThreadedElementLoopBase.h.

Referenced by ComputeFullJacobianThread::computeFaceJacobian(), ComputeJacobianThread::computeFaceJacobian(), ComputeFullJacobianThread::computeInternalFaceJacobian(), ComputeJacobianThread::computeInternalFaceJacobian(), ComputeFullJacobianThread::computeInternalInterFaceJacobian(), ComputeJacobianThread::computeInternalInterFaceJacobian(), ComputeFullJacobianThread::computeJacobian(), ComputeJacobianThread::computeJacobian(), ComputeResidualThread::onBoundary(), ComputeJacobianThread::onBoundary(), ComputeUserObjectsThread::onBoundary(), ComputeMaterialsObjectThread::onBoundary(), ComputeMarkerThread::onElement(), ComputeElemDampingThread::onElement(), ComputeElemAuxVarsThread::onElement(), ComputeResidualThread::onElement(), ComputeJacobianThread::onElement(), ComputeIndicatorThread::onElement(), ComputeUserObjectsThread::onElement(), ComputeMaterialsObjectThread::onElement(), ComputeResidualThread::onInterface(), ComputeJacobianThread::onInterface(), ComputeJacobianThread::onInternalSide(), ComputeResidualThread::onInternalSide(), ComputeIndicatorThread::onInternalSide(), ComputeMaterialsObjectThread::onInternalSide(), ComputeUserObjectsThread::onInternalSide(), ComputeMarkerThread::post(), ComputeElemAuxVarsThread::post(), ComputeMaterialsObjectThread::post(), ComputeIndicatorThread::post(), ComputeResidualThread::post(), ComputeJacobianThread::post(), ComputeUserObjectsThread::post(), ComputeResidualThread::postElement(), ComputeJacobianThread::postElement(), postElement(), ComputeMarkerThread::subdomainChanged(), ComputeElemAuxVarsThread::subdomainChanged(), ComputeResidualThread::subdomainChanged(), ComputeJacobianThread::subdomainChanged(), ComputeIndicatorThread::subdomainChanged(), ComputeMaterialsObjectThread::subdomainChanged(), and ComputeUserObjectsThread::subdomainChanged().

const KernelWarehouse* ComputeFullJacobianThread::_warehouse
protectedinherited

The documentation for this class was generated from the following files: