www.mooseframework.org
NonlocalIntegratedBC.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 "NonlocalIntegratedBC.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 
30  : IntegratedBC(parameters)
31 {
32  _mesh.errorIfDistributedMesh("NonlocalIntegratedBC");
33  mooseWarning("NonlocalIntegratedBC is a computationally expensive experimental capability used "
34  "only for integral terms.");
35 }
36 
37 void
39 {
40  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), _var.number());
41  _local_ke.resize(ke.m(), ke.n());
42  _local_ke.zero();
43 
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  MooseVariableFEBase & jv = _sys.getVariable(_tid, jvar_num);
77  DenseMatrix<Number> & ke = _assembly.jacobianBlock(_var.number(), jvar_num);
78 
79  for (_j = 0; _j < jvar.phiFaceSize();
80  _j++) // looping order for _i & _j are reversed for performance improvement
81  {
82  getUserObjectJacobian(jvar_num, jv.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 
100  for (_k = 0; _k < n_total_dofs;
101  _k++) // looping order for _i & _k are reversed for performance improvement
102  {
103  // eliminating the local components
104  auto it = local_dofindices.find(var_alldofindices[_k]);
105  if (it == local_dofindices.end())
106  {
107  getUserObjectJacobian(_var.number(), var_alldofindices[_k]);
108  // skip global DOFs that do not contribute to the jacobian
109  if (!globalDoFEnabled(_var, var_alldofindices[_k]))
110  continue;
111 
112  for (_i = 0; _i < _test.size(); _i++)
113  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
114  keg(_i, _k) += _JxW[_qp] * _coord[_qp] * computeQpNonlocalJacobian(var_alldofindices[_k]);
115  }
116  }
117 }
118 
119 void
121 {
122  if (jvar == _var.number())
124  else
125  {
127  DenseMatrix<Number> & keg = _assembly.jacobianBlockNonlocal(_var.number(), jvar);
128  // compiling set of global IDs for the local DOFs on the element
129  std::set<dof_id_type> local_dofindices(jv.dofIndices().begin(), jv.dofIndices().end());
130  // storing the global IDs for all the DOFs of the variable
131  const std::vector<dof_id_type> & jv_alldofindices = jv.allDofIndices();
132  unsigned int n_total_dofs = jv_alldofindices.size();
133 
134  for (_k = 0; _k < n_total_dofs;
135  _k++) // looping order for _i & _k are reversed for performance improvement
136  {
137  // eliminating the local components
138  auto it = local_dofindices.find(jv_alldofindices[_k]);
139  if (it == local_dofindices.end())
140  {
141  getUserObjectJacobian(jvar, jv_alldofindices[_k]);
142  // skip global DOFs that do not contribute to the jacobian
143  if (!globalDoFEnabled(jv, jv_alldofindices[_k]))
144  continue;
145 
146  for (_i = 0; _i < _test.size(); _i++)
147  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
148  keg(_i, _k) += _JxW[_qp] * _coord[_qp] *
149  computeQpNonlocalOffDiagJacobian(jvar, jv_alldofindices[_k]);
150  }
151  }
152  }
153 }
const VariableTestValue & _test
test function values (in QPs)
Definition: IntegratedBC.h:68
InputParameters validParams< NonlocalIntegratedBC >()
virtual Real computeQpNonlocalOffDiagJacobian(unsigned int, dof_id_type)
InputParameters validParams< IntegratedBC >()
Definition: IntegratedBC.C:23
std::vector< MooseVariableFEBase * > _diag_save_in
void mooseWarning(Args &&... args) const
Definition: MooseObject.h:155
virtual Real computeQpJacobian()
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
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
unsigned int _i
i-th, j-th index for enumerating test and shape functions
virtual void getUserObjectJacobian(unsigned int, dof_id_type)
Optimization option for getting jocobinas from userobject once per dof.
const VariablePhiValue & _phi
shape function values (in QPs)
Definition: IntegratedBC.h:61
DenseMatrix< Number > _local_ke
Holds residual entries as they are accumulated by this Kernel.
unsigned int _qp
quadrature point index
NonlocalIntegratedBC(const InputParameters &parameters)
MooseVariableFEBase & getVariable(THREAD_ID tid, const std::string &var_name)
Gets a reference to a variable of with specified name.
Definition: SystemBase.C:105
void errorIfDistributedMesh(std::string name) const
Generate a unified error message if the underlying libMesh mesh is a DistributedMesh.
Definition: MooseMesh.C:2685
virtual void computeNonlocalJacobian() override
computeNonlocalJacobian and computeNonlocalOffDiagJacobian methods are introduced for providing the j...
DenseMatrix< Number > & jacobianBlockNonlocal(unsigned int ivar, unsigned int jvar, TagID tag=0)
Definition: Assembly.C:2026
unsigned int size() const
The number of elements that can currently be stored in the array.
Definition: MooseArray.h:259
virtual void computeNonlocalOffDiagJacobian(unsigned int jvar) override
Computes d-residual / d-jvar...
const std::vector< dof_id_type > & dofIndices() const final
Get local DoF indices.
virtual bool globalDoFEnabled(MooseVariableFEBase &, dof_id_type)
optimization option for executing nonlocal jacobian calculation only for nonzero elements ...
virtual void computeJacobianBlock(MooseVariableFEBase &jvar) override
Computes d-ivar-residual / d-jvar...
virtual const std::vector< dof_id_type > & dofIndices() const
Get local DoF indices.
Assembly & _assembly
Reference to assembly.
Base class for deriving any boundary condition of a integrated type.
Definition: IntegratedBC.h:24
const MooseArray< Real > & _coord
coordinate transformation
virtual size_t phiFaceSize() const =0
Return phiFace size.
virtual void computeJacobian() override
computeJacobian and computeQpOffDiagJacobian methods are almost same as IntegratedBC except for few a...
MooseVariable & _var
Definition: IntegratedBC.h:53
const QBase *const & _qrule
active quadrature rule
SystemBase & _sys
Reference to SystemBase.
const std::vector< dof_id_type > & allDofIndices() const
Get all global dofindices for the variable.
virtual Real computeQpOffDiagJacobian(unsigned int)
This is the virtual that derived classes should override for computing an off-diagonal jacobian compo...
virtual Real computeQpNonlocalJacobian(dof_id_type)
Compute this IntegratedBC&#39;s contribution to the Jacobian corresponding to nolocal dof at the current ...
bool _has_diag_save_in
The aux variables to save the diagonal Jacobian contributions to.
THREAD_ID _tid
Thread id.
const MooseArray< Real > & _JxW
transformed Jacobian weights
MooseMesh & _mesh
Mesh this BC is defined on.