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 
23 {
25  return params;
26 }
27 
28 NonlocalKernel::NonlocalKernel(const InputParameters & parameters) : Kernel(parameters)
29 {
30  _mesh.errorIfDistributedMesh("NonlocalKernel");
31  mooseWarning("NonlocalKernel is a computationally expensive experimental capability used only "
32  "for integral terms.");
33 }
34 
35 void
37 {
40  for (_j = 0; _j < _phi.size();
41  _j++) // looping order for _i & _j are reversed for performance improvement
42  {
44  for (_i = 0; _i < _test.size(); _i++)
45  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
47  }
49 
51  {
52  unsigned int rows = _local_ke.m();
53  DenseVector<Number> diag(rows);
54  for (unsigned int i = 0; i < rows; i++)
55  diag(i) = _local_ke(i, i);
56 
57  for (const auto & var : _diag_save_in)
58  var->sys().solution().add_vector(diag, var->dofIndices());
59  }
60 }
61 
62 void
63 NonlocalKernel::computeOffDiagJacobian(const unsigned int jvar_num)
64 {
65  if (jvar_num == _var.number())
67  else
68  {
69  const auto & jvar = getVariable(jvar_num);
70 
71  prepareMatrixTag(_assembly, _var.number(), jvar_num);
72 
73  // This (undisplaced) jvar could potentially yield the wrong phi size if this object is acting
74  // on the displaced mesh
75  const auto phi_size = jvar.dofIndices().size();
76 
78  for (_j = 0; _j < phi_size;
79  _j++) // looping order for _i & _j are reversed for performance improvement
80  {
81  getUserObjectJacobian(jvar_num, jvar.dofIndices()[_j]);
82  for (_i = 0; _i < _test.size(); _i++)
83  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
85  }
87  }
88 }
89 
90 void
92 {
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  _nonlocal_ke(_i, _k) +=
116  _JxW[_qp] * _coord[_qp] * computeQpNonlocalJacobian(var_alldofindices[_k]);
117  }
118  }
120 }
121 
122 void
124 {
125  if (jvar == _var.number())
127  else
128  {
131  // compiling set of global IDs for the local DOFs on the element
132  std::set<dof_id_type> local_dofindices(jv.dofIndices().begin(), jv.dofIndices().end());
133  // storing the global IDs for all the DOFs of the variable
134  const std::vector<dof_id_type> & jv_alldofindices = jv.allDofIndices();
135  unsigned int n_total_dofs = jv_alldofindices.size();
136 
138  for (_k = 0; _k < n_total_dofs;
139  _k++) // looping order for _i & _k are reversed for performance improvement
140  {
141  // eliminating the local components
142  auto it = local_dofindices.find(jv_alldofindices[_k]);
143  if (it == local_dofindices.end())
144  {
145  getUserObjectJacobian(jvar, jv_alldofindices[_k]);
146  // skip global DOFs that do not contribute to the jacobian
147  if (!globalDoFEnabled(jv, jv_alldofindices[_k]))
148  continue;
149 
150  for (_i = 0; _i < _test.size(); _i++)
151  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
152  _nonlocal_ke(_i, _k) += _JxW[_qp] * _coord[_qp] *
153  computeQpNonlocalOffDiagJacobian(jvar, jv_alldofindices[_k]);
154  }
155  }
157  }
158 }
MooseMesh & _mesh
Reference to this Kernel&#39;s mesh object.
void accumulateTaggedNonlocalMatrix()
Nonlocal Jacobian blocks will be appended by adding the current nonlocal kernel Jacobian.
static InputParameters validParams()
Definition: Kernel.C:23
MooseVariable & _var
This is a regular kernel so we cast to a regular MooseVariable.
Definition: Kernel.h:72
virtual void getUserObjectJacobian(unsigned int, dof_id_type)
Optimization option for getting jocobinas from userobject once per dof.
const MooseArray< Real > & _JxW
The current quadrature point weight value.
Definition: KernelBase.h:51
virtual void precalculateOffDiagJacobian(unsigned int)
unsigned int number() const
Get variable number coming from libMesh.
virtual Real computeQpOffDiagJacobian(unsigned int)
For coupling standard variables.
Definition: Kernel.h:56
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int m() const
This class provides an interface for common operations on field variables of both FE and FV types wit...
virtual void computeJacobian() override
computeJacobian and computeQpOffDiagJacobian methods are almost same as Kernel except for few additio...
THREAD_ID _tid
The thread ID for this kernel.
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
Definition: KernelBase.h:68
DenseMatrix< Number > _local_ke
Holds local Jacobian entries as they are accumulated by this Kernel.
const MooseVariableFieldBase & getVariable(unsigned int jvar_num) const
Retrieve the variable object from our system associated with jvar_num.
void mooseWarning(Args &&... args) const
Emits a warning prefixed with object name and type.
NonlocalKernel(const InputParameters &parameters)
SystemBase & _sys
Reference to the EquationSystem object.
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:3367
void prepareMatrixTagNonlocal(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing nonlocal element jacobian according to the active tags.
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:256
const VariableTestValue & _test
the current test function
Definition: Kernel.h:75
const std::vector< dof_id_type > & dofIndices() const final
Get local DoF indices.
static InputParameters validParams()
virtual void precalculateJacobian()
std::vector< MooseVariableFEBase * > _diag_save_in
Definition: KernelBase.h:69
const QBase *const & _qrule
active quadrature rule
Definition: KernelBase.h:48
virtual bool globalDoFEnabled(MooseVariableFEBase &, dof_id_type)
optimization option for executing nonlocal jacobian calculation only for nonzero elements ...
virtual Real computeQpJacobian()
Compute this Kernel&#39;s contribution to the Jacobian at the current quadrature point.
Definition: Kernel.h:51
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:57
void accumulateTaggedLocalMatrix()
Local Jacobian blocks will be appended by adding the current local kernel Jacobian.
const MooseArray< Real > & _coord
The scaling factor to convert from cartesian to another coordinate system (e.g rz, spherical, etc.)
Definition: KernelBase.h:54
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
virtual void computeOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar... storing the result in Ke.
virtual Real computeQpNonlocalOffDiagJacobian(unsigned int, dof_id_type)
unsigned int _j
current index for the shape function
Definition: KernelBase.h:60
virtual void computeNonlocalOffDiagJacobian(unsigned int jvar) override
Computes Jacobian entries corresponding to nonlocal dofs of the jvar.
const std::vector< dof_id_type > & allDofIndices() const
Get all global dofindices for the variable.
Definition: Kernel.h:15
unsigned int _k
MooseVariableFieldBase & getVariable(THREAD_ID tid, const std::string &var_name) const
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:79
virtual void computeNonlocalJacobian() override
computeNonlocalJacobian and computeNonlocalOffDiagJacobian methods are introduced for providing the j...
void prepareMatrixTag(Assembly &assembly, unsigned int ivar, unsigned int jvar)
Prepare data for computing element jacobian according to the active tags.
DenseMatrix< Number > _nonlocal_ke
Holds nonlocal Jacobian entries as they are accumulated by this Kernel.
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:81
unsigned int _qp
The current quadrature point index.
Definition: KernelBase.h:42