https://mooseframework.inl.gov
MortarConstraintBase.C
Go to the documentation of this file.
1 //* This file is part of the MOOSE framework
2 //* https://mooseframework.inl.gov
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 "MortarConstraintBase.h"
11 #include "FEProblemBase.h"
12 #include "Assembly.h"
13 #include "MooseVariableFE.h"
14 
15 #include "libmesh/string_to_enum.h"
16 
19 {
23 
24  // Whether on a displaced or undisplaced mesh, coupling ghosting will only happen for
25  // cross-interface elements
26  params.addRelationshipManager("AugmentSparsityOnInterface",
28  [](const InputParameters & obj_params, InputParameters & rm_params)
29  {
30  rm_params.set<bool>("use_displaced_mesh") =
31  obj_params.get<bool>("use_displaced_mesh");
32  rm_params.set<BoundaryName>("secondary_boundary") =
33  obj_params.get<BoundaryName>("secondary_boundary");
34  rm_params.set<BoundaryName>("primary_boundary") =
35  obj_params.get<BoundaryName>("primary_boundary");
36  rm_params.set<SubdomainName>("secondary_subdomain") =
37  obj_params.get<SubdomainName>("secondary_subdomain");
38  rm_params.set<SubdomainName>("primary_subdomain") =
39  obj_params.get<SubdomainName>("primary_subdomain");
40  rm_params.set<bool>("ghost_higher_d_neighbors") =
41  obj_params.get<bool>("ghost_higher_d_neighbors");
42  });
43 
44  // If the LM is ever a Lagrange variable, we will attempt to obtain its dof indices from the
45  // process that owns the node. However, in order for the dof indices to get set for the Lagrange
46  // variable, the process that owns the node needs to have local copies of any lower-d elements
47  // that have the connected node. Note that the geometric ghosting done here is different than that
48  // done by the AugmentSparsityOnInterface RM, even when ghost_point_neighbors is true. The latter
49  // ghosts equal-manifold lower-dimensional secondary element point neighbors and their interface
50  // couplings. This ghosts lower-dimensional point neighbors of higher-dimensional elements.
51  // Neither is guaranteed to be a superset of the other. For instance ghosting of lower-d point
52  // neighbors (AugmentSparsityOnInterface with ghost_point_neighbors = true) is only guaranteed to
53  // ghost those lower-d point neighbors on *processes that own lower-d elements*. And you may have
54  // a process that only owns higher-dimensional elements
55  //
56  // Note that in my experience it is only important for the higher-d lower-d point neighbors to be
57  // ghosted when forming sparsity patterns and so I'm putting this here instead of at the
58  // MortarConsumerInterface level
59  params.addRelationshipManager("GhostHigherDLowerDPointNeighbors",
61 
62  params.addParam<VariableName>("secondary_variable", "Primal variable on secondary surface.");
63  params.addParam<VariableName>(
64  "primary_variable",
65  "Primal variable on primary surface. If this parameter is not provided then the primary "
66  "variable will be initialized to the secondary variable");
67  params.makeParamNotRequired<NonlinearVariableName>("variable");
68  params.setDocString(
69  "variable",
70  "The name of the lagrange multiplier variable that this constraint is applied to. This "
71  "parameter may not be supplied in the case of using penalty methods for example");
72  params.addParam<bool>(
73  "compute_primal_residuals", true, "Whether to compute residuals for the primal variable.");
74  params.addParam<bool>(
75  "compute_lm_residuals", true, "Whether to compute Lagrange Multiplier residuals");
76  params.addParam<MooseEnum>(
77  "quadrature",
78  MooseEnum("DEFAULT FIRST SECOND THIRD FOURTH", "DEFAULT"),
79  "Quadrature rule to use on mortar segments. For 2D mortar DEFAULT is recommended. "
80  "For 3D mortar, QUAD meshes are integrated using triangle mortar segments. "
81  "While DEFAULT quadrature order is typically sufficiently accurate, exact integration of "
82  "QUAD mortar faces requires SECOND order quadrature for FIRST variables and FOURTH order "
83  "quadrature for SECOND order variables.");
84  params.addParam<bool>(
85  "use_petrov_galerkin",
86  false,
87  "Whether to use the Petrov-Galerkin approach for the mortar-based constraints. If set to "
88  "true, we use the standard basis as the test function and dual basis as "
89  "the shape function for the interpolation of the Lagrange multiplier variable.");
90  params.addCoupledVar("aux_lm",
91  "Auxiliary Lagrange multiplier variable that is utilized together with the "
92  "Petrov-Galerkin approach.");
93  return params;
94 }
95 
97  : Constraint(parameters),
102  true,
103  isParamValid("variable") ? "variable" : "secondary_variable",
106  _fe_problem(*getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
107  _var(isParamValid("variable")
108  ? &_subproblem.getStandardVariable(_tid, parameters.getMooseType("variable"))
109  : nullptr),
110  _secondary_var(
111  isParamValid("secondary_variable")
112  ? _sys.getActualFieldVariable<Real>(_tid, parameters.getMooseType("secondary_variable"))
113  : _sys.getActualFieldVariable<Real>(_tid, parameters.getMooseType("primary_variable"))),
114  _primary_var(
115  isParamValid("primary_variable")
116  ? _sys.getActualFieldVariable<Real>(_tid, parameters.getMooseType("primary_variable"))
117  : _secondary_var),
118 
119  _compute_primal_residuals(getParam<bool>("compute_primal_residuals")),
120  _compute_lm_residuals(!_var ? false : getParam<bool>("compute_lm_residuals")),
121  _test_dummy(),
122  _use_dual(_var ? _var->useDual() : false),
123  _tangents(_assembly.tangents()),
124  _coord(_assembly.mortarCoordTransformation()),
125  _q_point(_assembly.qPointsMortar()),
126  _use_petrov_galerkin(getParam<bool>("use_petrov_galerkin")),
127  _aux_lm_var(isCoupled("aux_lm") ? getVar("aux_lm", 0) : nullptr),
128  _test(_var
129  ? ((_use_petrov_galerkin && _aux_lm_var) ? _aux_lm_var->phiLower() : _var->phiLower())
130  : _test_dummy),
131  _test_secondary(_secondary_var.phiFace()),
132  _test_primary(_primary_var.phiFaceNeighbor()),
133  _grad_test_secondary(_secondary_var.gradPhiFace()),
134  _grad_test_primary(_primary_var.gradPhiFaceNeighbor()),
135  _interior_secondary_elem(_assembly.elem()),
136  _interior_primary_elem(_assembly.neighbor()),
137  _displaced(getParam<bool>("use_displaced_mesh"))
138 {
139  if (_use_dual)
141 
143  paramError("use_petrov_galerkin",
144  "We need to set `use_dual = true` while using the Petrov-Galerkin approach");
145 
146  if (_use_petrov_galerkin && ((!isParamValid("aux_lm")) || _aux_lm_var == nullptr))
147  paramError("use_petrov_galerkin",
148  "We need to specify an auxiliary variable `aux_lm` while using the Petrov-Galerkin "
149  "approach");
150 
152  paramError("aux_lm",
153  "Auxiliary LM variable needs to use standard shape function, i.e., set `use_dual = "
154  "false`.");
155 
156  // Note parameter is discretization order, we then convert to quadrature order
157  const MooseEnum p_order = getParam<MooseEnum>("quadrature");
158  // If quadrature not DEFAULT, set mortar qrule
159  if (p_order != "DEFAULT")
160  {
161  Order q_order = static_cast<Order>(2 * Utility::string_to_enum<Order>(p_order) + 1);
162  _assembly.setMortarQRule(q_order);
163  }
164 
165  if (_var)
169 }
170 
171 void
173 {
175 
177  {
178  // Compute the residual for the secondary interior primal dofs
180 
181  // Compute the residual for the primary interior primal dofs.
183  }
184 
186  // Compute the residual for the lower dimensional LM dofs (if we even have an LM variable)
188 }
189 
190 void
192 {
194 
196  {
197  // Compute the jacobian for the secondary interior primal dofs
199 
200  // Compute the jacobian for the primary interior primal dofs.
202  }
203 
205  // Compute the jacobian for the lower dimensional LM dofs (if we even have an LM variable)
207 }
208 
209 void
210 MortarConstraintBase::zeroInactiveLMDofs(const std::unordered_set<const Node *> & inactive_lm_nodes,
211  const std::unordered_set<const Elem *> & inactive_lm_elems)
212 {
213  // If no LM variable has been defined, skip
214  if (!_var)
215  return;
216 
217  const auto sn = _sys.number();
218  const auto vn = _var->number();
219 
220  // If variable is nodal, zero DoFs based on inactive LM nodes
221  if (_var->isNodal())
222  {
223  for (const auto node : inactive_lm_nodes)
224  {
225  // Allow mixed Lagrange orders between primal and LM
226  if (!node->n_comp(sn, vn))
227  continue;
228 
229  const auto dof_index = node->dof_number(sn, vn, 0);
230  // No scaling; this is not physics
233  _assembly, /*element_value=*/1, dof_index, dof_index, /*scaling_factor=*/1);
235  {
236  const Real lm_value = _var->getNodalValue(*node);
238  std::array<Real, 1>{{lm_value}},
239  std::array<dof_id_type, 1>{{dof_index}},
240  /*scaling_factor=*/1);
241  }
242  }
243  }
244  // If variable is elemental, zero based on inactive LM elems
245  else
246  {
247  for (const auto el : inactive_lm_elems)
248  {
249  const auto n_comp = el->n_comp(sn, vn);
250 
251  for (const auto comp : make_range(n_comp))
252  {
253  const auto dof_index = el->dof_number(sn, vn, comp);
254  // No scaling; this is not physics
257  _assembly, /*element_value=*/1, dof_index, dof_index, /*scaling_factor=*/1);
259  {
260  const Real lm_value = _var->getElementalValue(el, comp);
262  std::array<Real, 1>{{lm_value}},
263  std::array<dof_id_type, 1>{{dof_index}},
264  /*scaling_factor=*/1);
265  }
266  }
267  }
268  }
269 }
VarFieldType
Definition: MooseTypes.h:722
MooseVariableField< Real > & _secondary_var
Reference to the secondary variable.
Order
void setMortarQRule(Order order)
Specifies a custom qrule for integration on mortar segment mesh.
Definition: Assembly.C:735
Intermediate base class that ties together all the interfaces for getting MooseVariables with the Moo...
void paramError(const std::string &param, Args... args) const
Emits an error prefixed with the file and line number of the given param (from the input file) along ...
Definition: MooseBase.h:435
void setDocString(const std::string &name, const std::string &doc)
Set the doc string of a parameter.
unsigned int number() const
Get variable number coming from libMesh.
std::vector< std::pair< R1, R2 > > get(const std::string &param1, const std::string &param2) const
Combine two vector parameters into a single vector of pairs.
void zeroInactiveLMDofs(const std::unordered_set< const Node *> &inactive_lm_nodes, const std::unordered_set< const Elem *> &inactive_lm_elems)
A post routine for zeroing all inactive LM DoFs.
Base class for all Constraint types.
Definition: Constraint.h:19
T & set(const std::string &name, bool quiet_mode=false)
Returns a writable reference to the named parameters.
void addResiduals(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
Add the provided incoming residuals corresponding to the provided dof indices.
virtual void precalculateResidual()
The main MOOSE class responsible for handling user-defined parameters in almost every MOOSE system...
const bool _compute_primal_residuals
Whether to compute primal residuals.
void addRelationshipManager(const std::string &name, Moose::RelationshipManagerType rm_type, Moose::RelationshipManagerInputParameterCallback input_parameter_callback=nullptr)
Tells MOOSE about a RelationshipManager that this object needs.
Specialization of SubProblem for solving nonlinear equations plus auxiliary equations.
MortarConstraintBase(const InputParameters &parameters)
virtual void computeResidual() override
Method for computing the residual.
An interface for accessing mortar mesh data.
static InputParameters validParams()
OutputData getElementalValue(const Elem *elem, unsigned int idx=0) const
Get the current value of this variable on an element.
SystemBase & _sys
Reference to the EquationSystem object.
bool useDual() const
Get dual mortar option.
VarKindType
Framework-wide stuff.
Definition: MooseTypes.h:715
This is a "smart" enum class intended to replace many of the shortcomings in the C++ enum type It sho...
Definition: MooseEnum.h:33
std::vector< BoundaryID > getBoundaryIDs(const libMesh::MeshBase &mesh, const std::vector< BoundaryName > &boundary_name, bool generate_unknown, const std::set< BoundaryID > &mesh_boundary_ids)
Gets the boundary IDs with their names.
unsigned int number() const
Gets the number of this system.
Definition: SystemBase.C:1149
const bool _use_dual
Whether to use the dual motar approach.
Assembly & _assembly
Reference to this Kernel&#39;s assembly object.
void addCoupledVar(const std::string &name, const std::string &doc_string)
This method adds a coupled variable name pair.
OutputData getNodalValue(const Node &node) const
Get the value of this variable at given node.
bool isNodal() const override
Is this variable nodal.
void makeParamNotRequired(const std::string &name)
Changes the parameter to not be required.
MooseVariable *const _var
Pointer to the lagrange multipler variable. nullptr if none.
void addMooseVariableDependency(MooseVariableFieldBase *var)
Call this function to add the passed in MooseVariableFieldBase as a variable that this object depends...
const bool _use_petrov_galerkin
Whether to use Petrov-Galerkin approach.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const std::set< SubdomainID > EMPTY_BLOCK_IDS
Definition: MooseTypes.h:684
bool computingResidual() const
Definition: Assembly.h:1901
const MooseVariable *const _aux_lm_var
The auxiliary Lagrange multiplier variable (used together whith the Petrov-Galerkin approach) ...
void activateDual()
Indicates that dual shape functions are used for mortar constraint.
Definition: Assembly.h:578
static InputParameters validParams()
virtual void computeJacobian() override
Method for computing the Jacobian.
Interface for objects that need to get values of MooseVariables.
IntRange< T > make_range(T beg, T end)
This interface is designed for DGKernel, InternalSideUserObject, InterfaceUserObject, where material properties on a side of both its primary side (face) and its secondary side (neighbor) all required.
void addParam(const std::string &name, const S &value, const std::string &doc_string)
These methods add an optional parameter and a documentation string to the InputParameters object...
bool isParamValid(const std::string &name) const
Test if the supplied parameter is valid.
Definition: MooseBase.h:195
MOOSE now contains C++17 code, so give a reasonable error message stating what the user can do to add...
bool computingJacobian() const
Definition: Assembly.h:1906
void addJacobianElement(Assembly &assembly, Real value, dof_id_type row_index, dof_id_type column_index, Real scaling_factor)
Add into a single Jacobian element.
static InputParameters validParams()
Definition: Constraint.C:15
const bool _compute_lm_residuals
Whether to compute lagrange multiplier residuals.
MooseVariableField< Real > & _primary_var
Reference to the primary variable.