www.mooseframework.org
NonlocalKernel.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 "NonlocalKernel.h"
16 #include "Assembly.h"
17 #include "MooseVariable.h"
18 #include "Problem.h"
19 #include "SubProblem.h"
20 #include "SystemBase.h"
21 #include "MooseMesh.h"
22 
23 #include "libmesh/threads.h"
24 #include "libmesh/quadrature.h"
25 
26 template <>
29 {
31  return params;
32 }
33 
34 NonlocalKernel::NonlocalKernel(const InputParameters & parameters) : Kernel(parameters)
35 {
36  _mesh.errorIfDistributedMesh("NonlocalKernel");
37  mooseWarning("NonlocalKernel is a computationally expensive experimental capability used only "
38  "for integral terms.");
39 }
40 
41 void
43 {
44  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
45  _local_ke.resize(ke.m(), ke.n());
46  _local_ke.zero();
47 
49  for (_j = 0; _j < _phi.size();
50  _j++) // looping order for _i & _j are reversed for performance improvement
51  {
53  for (_i = 0; _i < _test.size(); _i++)
54  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
56  }
57 
58  ke += _local_ke;
59 
61  {
62  unsigned int rows = ke.m();
63  DenseVector<Number> diag(rows);
64  for (unsigned int i = 0; i < rows; i++)
65  diag(i) = _local_ke(i, i);
66 
67  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
68  for (const auto & var : _diag_save_in)
69  var->sys().solution().add_vector(diag, var->dofIndices());
70  }
71 }
72 
73 void
75 {
76  if (jvar == _var.number())
78  else
79  {
80  MooseVariable & jv = _sys.getVariable(_tid, jvar);
81  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar);
82 
84  for (_j = 0; _j < _phi.size();
85  _j++) // looping order for _i & _j are reversed for performance improvement
86  {
87  getUserObjectJacobian(jvar, jv.dofIndices()[_j]);
88  for (_i = 0; _i < _test.size(); _i++)
89  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
90  ke(_i, _j) += _JxW[_qp] * _coord[_qp] * computeQpOffDiagJacobian(jvar);
91  }
92  }
93 }
94 
95 void
97 {
98  DenseMatrix<Number> & keg = _assembly.jacobianBlockNonlocal(_var.number(), _var.number());
99  // compiling set of global IDs for the local DOFs on the element
100  std::set<dof_id_type> local_dofindices(_var.dofIndices().begin(), _var.dofIndices().end());
101  // storing the global IDs for all the DOFs of the variable
102  const std::vector<dof_id_type> & var_alldofindices = _var.allDofIndices();
103  unsigned int n_total_dofs = var_alldofindices.size();
104 
106  for (_k = 0; _k < n_total_dofs;
107  _k++) // looping order for _i & _k are reversed for performance improvement
108  {
109  // eliminating the local components
110  auto it = local_dofindices.find(var_alldofindices[_k]);
111  if (it == local_dofindices.end())
112  {
113  getUserObjectJacobian(_var.number(), var_alldofindices[_k]);
114  // skip global DOFs that do not contribute to the jacobian
115  if (!globalDoFEnabled(_var, var_alldofindices[_k]))
116  continue;
117 
118  for (_i = 0; _i < _test.size(); _i++)
119  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
120  keg(_i, _k) += _JxW[_qp] * _coord[_qp] * computeQpNonlocalJacobian(var_alldofindices[_k]);
121  }
122  }
123 }
124 
125 void
127 {
128  if (jvar == _var.number())
130  else
131  {
132  MooseVariable & jv = _sys.getVariable(_tid, jvar);
133  DenseMatrix<Number> & keg = _assembly.jacobianBlockNonlocal(_var.number(), jvar);
134  // compiling set of global IDs for the local DOFs on the element
135  std::set<dof_id_type> local_dofindices(jv.dofIndices().begin(), jv.dofIndices().end());
136  // storing the global IDs for all the DOFs of the variable
137  const std::vector<dof_id_type> & jv_alldofindices = jv.allDofIndices();
138  unsigned int n_total_dofs = jv_alldofindices.size();
139 
141  for (_k = 0; _k < n_total_dofs;
142  _k++) // looping order for _i & _k are reversed for performance improvement
143  {
144  // eliminating the local components
145  auto it = local_dofindices.find(jv_alldofindices[_k]);
146  if (it == local_dofindices.end())
147  {
148  getUserObjectJacobian(jvar, jv_alldofindices[_k]);
149  // skip global DOFs that do not contribute to the jacobian
150  if (!globalDoFEnabled(jv, jv_alldofindices[_k]))
151  continue;
152 
153  for (_i = 0; _i < _test.size(); _i++)
154  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
155  keg(_i, _k) += _JxW[_qp] * _coord[_qp] *
156  computeQpNonlocalOffDiagJacobian(jvar, jv_alldofindices[_k]);
157  }
158  }
159  }
160 }
const std::vector< dof_id_type > & allDofIndices() const
Get all global dofindices for the variable.
void mooseWarning(Args &&...args) const
Definition: MooseObject.h:89
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: Kernel.C:126
std::vector< MooseVariable * > _diag_save_in
Definition: KernelBase.h:176
virtual void getUserObjectJacobian(unsigned int, dof_id_type)
Optimization option for getting jocobinas from userobject once per dof.
Class for stuff related to variables.
Definition: MooseVariable.h:43
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:140
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
Definition: KernelBase.h:167
const VariableTestValue & _test
the current test function
Definition: KernelBase.h:152
SystemBase & _sys
Reference to the EquationSystem object.
Definition: KernelBase.h:111
THREAD_ID _tid
The thread ID for this kernel.
Definition: KernelBase.h:114
QBase *& _qrule
active quadrature rule
Definition: KernelBase.h:137
virtual void computeOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar... storing the result in Ke.
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:2473
virtual bool globalDoFEnabled(MooseVariable &, dof_id_type)
optimization option for executing nonlocal jacobian calculation only for nonzero elements ...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:175
NonlocalKernel(const InputParameters &parameters)
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
Definition: KernelBase.h:117
virtual MooseVariable & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:103
InputParameters validParams< NonlocalKernel >()
unsigned int _i
current index for the test function
Definition: KernelBase.h:146
virtual Real computeQpOffDiagJacobian(unsigned int jvar)
This is the virtual that derived classes should override for computing an off-diagonal Jacobian compo...
Definition: Kernel.C:132
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:143
virtual Real computeQpNonlocalOffDiagJacobian(unsigned int, dof_id_type)
const VariablePhiValue & _phi
the current shape functions
Definition: KernelBase.h:158
MooseMesh & _mesh
Reference to this Kernel&#39;s mesh object.
Definition: KernelBase.h:123
unsigned int _j
current index for the shape function
Definition: KernelBase.h:149
InputParameters validParams< Kernel >()
Definition: Kernel.C:30
virtual void computeNonlocalOffDiagJacobian(unsigned int jvar)
Computes d-residual / d-jvar...
unsigned int number() const
Get variable number coming from libMesh.
std::vector< dof_id_type > & dofIndices()
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:250
virtual void precalculateOffDiagJacobian(unsigned int)
Definition: Kernel.h:48
Definition: Kernel.h:25
virtual void computeNonlocalJacobian()
computeNonlocalJacobian and computeNonlocalOffDiagJacobian methods are introduced for providing the j...
unsigned int _k
DenseMatrix< Number > & jacobianBlock(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:887
virtual Real computeQpNonlocalJacobian(dof_id_type)
Compute this Kernel&#39;s contribution to the Jacobian corresponding to nolocal dof at the current quadra...
MooseVariable & _var
Reference to this Kernel&#39;s MooseVariable object.
Definition: KernelBase.h:120
DenseMatrix< Number > & jacobianBlockNonlocal(unsigned int ivar, unsigned int jvar)
Definition: Assembly.C:894
virtual void precalculateJacobian()
Definition: Kernel.h:47
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:131
virtual void computeJacobian()
computeJacobian and computeQpOffDiagJacobian methods are almost same as Kernel except for few additio...