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