Line data Source code
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 :
17 : InputParameters
18 216431 : MortarConstraintBase::validParams()
19 : {
20 216431 : InputParameters params = Constraint::validParams();
21 216431 : params += MortarConsumerInterface::validParams();
22 216431 : params += TwoMaterialPropertyInterface::validParams();
23 :
24 : // Whether on a displaced or undisplaced mesh, coupling ghosting will only happen for
25 : // cross-interface elements
26 216431 : params.addRelationshipManager("AugmentSparsityOnInterface",
27 : Moose::RelationshipManagerType::COUPLING,
28 0 : [](const InputParameters & obj_params, InputParameters & rm_params)
29 : {
30 3688 : rm_params.set<bool>("use_displaced_mesh") =
31 3688 : obj_params.get<bool>("use_displaced_mesh");
32 7376 : rm_params.set<BoundaryName>("secondary_boundary") =
33 7376 : obj_params.get<BoundaryName>("secondary_boundary");
34 7376 : rm_params.set<BoundaryName>("primary_boundary") =
35 7376 : obj_params.get<BoundaryName>("primary_boundary");
36 7376 : rm_params.set<SubdomainName>("secondary_subdomain") =
37 7376 : obj_params.get<SubdomainName>("secondary_subdomain");
38 7376 : rm_params.set<SubdomainName>("primary_subdomain") =
39 7376 : obj_params.get<SubdomainName>("primary_subdomain");
40 3688 : rm_params.set<bool>("ghost_higher_d_neighbors") =
41 3688 : obj_params.get<bool>("ghost_higher_d_neighbors");
42 3688 : });
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 216431 : params.addRelationshipManager("GhostHigherDLowerDPointNeighbors",
60 : Moose::RelationshipManagerType::GEOMETRIC);
61 :
62 216431 : params.addParam<VariableName>("secondary_variable", "Primal variable on secondary surface.");
63 216431 : 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 216431 : params.makeParamNotRequired<NonlinearVariableName>("variable");
68 216431 : 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 649293 : params.addParam<bool>(
73 432862 : "compute_primal_residuals", true, "Whether to compute residuals for the primal variable.");
74 649293 : params.addParam<bool>(
75 432862 : "compute_lm_residuals", true, "Whether to compute Lagrange Multiplier residuals");
76 649293 : params.addParam<MooseEnum>(
77 : "quadrature",
78 432862 : 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 649293 : params.addParam<bool>(
85 : "use_petrov_galerkin",
86 432862 : 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 216431 : params.addCoupledVar("aux_lm",
91 : "Auxiliary Lagrange multiplier variable that is utilized together with the "
92 : "Petrov-Galerkin approach.");
93 216431 : return params;
94 0 : }
95 :
96 1232 : MortarConstraintBase::MortarConstraintBase(const InputParameters & parameters)
97 : : Constraint(parameters),
98 : NeighborCoupleableMooseVariableDependencyIntermediateInterface(this, false, false),
99 : MortarConsumerInterface(this),
100 : TwoMaterialPropertyInterface(this, Moose::EMPTY_BLOCK_IDS, getBoundaryIDs()),
101 : MooseVariableInterface<Real>(this,
102 : true,
103 2464 : isParamValid("variable") ? "variable" : "secondary_variable",
104 : Moose::VarKindType::VAR_SOLVER,
105 : Moose::VarFieldType::VAR_FIELD_STANDARD),
106 1232 : _fe_problem(*getCheckedPointerParam<FEProblemBase *>("_fe_problem_base")),
107 2464 : _var(isParamValid("variable")
108 1232 : ? &_subproblem.getStandardVariable(_tid, parameters.getMooseType("variable"))
109 : : nullptr),
110 1232 : _secondary_var(
111 2464 : isParamValid("secondary_variable")
112 2464 : ? _sys.getActualFieldVariable<Real>(_tid, parameters.getMooseType("secondary_variable"))
113 1232 : : _sys.getActualFieldVariable<Real>(_tid, parameters.getMooseType("primary_variable"))),
114 1232 : _primary_var(
115 2464 : isParamValid("primary_variable")
116 1232 : ? _sys.getActualFieldVariable<Real>(_tid, parameters.getMooseType("primary_variable"))
117 : : _secondary_var),
118 :
119 1232 : _compute_primal_residuals(getParam<bool>("compute_primal_residuals")),
120 1232 : _compute_lm_residuals(!_var ? false : getParam<bool>("compute_lm_residuals")),
121 1232 : _test_dummy(),
122 1232 : _use_dual(_var ? _var->useDual() : false),
123 1232 : _tangents(_assembly.tangents()),
124 1232 : _coord(_assembly.mortarCoordTransformation()),
125 1232 : _q_point(_assembly.qPointsMortar()),
126 1232 : _use_petrov_galerkin(getParam<bool>("use_petrov_galerkin")),
127 1232 : _aux_lm_var(isCoupled("aux_lm") ? getVar("aux_lm", 0) : nullptr),
128 2464 : _test(_var
129 1232 : ? ((_use_petrov_galerkin && _aux_lm_var) ? _aux_lm_var->phiLower() : _var->phiLower())
130 : : _test_dummy),
131 1232 : _test_secondary(_secondary_var.phiFace()),
132 1232 : _test_primary(_primary_var.phiFaceNeighbor()),
133 1232 : _grad_test_secondary(_secondary_var.gradPhiFace()),
134 1232 : _grad_test_primary(_primary_var.gradPhiFaceNeighbor()),
135 1232 : _interior_secondary_elem(_assembly.elem()),
136 1232 : _interior_primary_elem(_assembly.neighbor()),
137 4928 : _displaced(getParam<bool>("use_displaced_mesh"))
138 : {
139 1232 : if (_use_dual)
140 124 : _assembly.activateDual();
141 :
142 1232 : if (_use_petrov_galerkin && (!_use_dual))
143 0 : paramError("use_petrov_galerkin",
144 : "We need to set `use_dual = true` while using the Petrov-Galerkin approach");
145 :
146 1232 : if (_use_petrov_galerkin && ((!isParamValid("aux_lm")) || _aux_lm_var == nullptr))
147 0 : paramError("use_petrov_galerkin",
148 : "We need to specify an auxiliary variable `aux_lm` while using the Petrov-Galerkin "
149 : "approach");
150 :
151 1232 : if (_use_petrov_galerkin && _aux_lm_var->useDual())
152 0 : 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 1232 : const MooseEnum p_order = getParam<MooseEnum>("quadrature");
158 : // If quadrature not DEFAULT, set mortar qrule
159 1232 : if (p_order != "DEFAULT")
160 : {
161 0 : Order q_order = static_cast<Order>(2 * Utility::string_to_enum<Order>(p_order) + 1);
162 0 : _assembly.setMortarQRule(q_order);
163 : }
164 :
165 1232 : if (_var)
166 668 : addMooseVariableDependency(_var);
167 1232 : addMooseVariableDependency(&_secondary_var);
168 1232 : addMooseVariableDependency(&_primary_var);
169 1232 : }
170 :
171 : void
172 188598 : MortarConstraintBase::computeResidual()
173 : {
174 188598 : precalculateResidual();
175 :
176 188598 : if (_compute_primal_residuals)
177 : {
178 : // Compute the residual for the secondary interior primal dofs
179 188598 : computeResidual(Moose::MortarType::Secondary);
180 :
181 : // Compute the residual for the primary interior primal dofs.
182 188598 : computeResidual(Moose::MortarType::Primary);
183 : }
184 :
185 188598 : if (_compute_lm_residuals)
186 : // Compute the residual for the lower dimensional LM dofs (if we even have an LM variable)
187 89334 : computeResidual(Moose::MortarType::Lower);
188 188598 : }
189 :
190 : void
191 89628 : MortarConstraintBase::computeJacobian()
192 : {
193 89628 : precalculateResidual();
194 :
195 89628 : if (_compute_primal_residuals)
196 : {
197 : // Compute the jacobian for the secondary interior primal dofs
198 89628 : computeJacobian(Moose::MortarType::Secondary);
199 :
200 : // Compute the jacobian for the primary interior primal dofs.
201 89628 : computeJacobian(Moose::MortarType::Primary);
202 : }
203 :
204 89628 : if (_compute_lm_residuals)
205 : // Compute the jacobian for the lower dimensional LM dofs (if we even have an LM variable)
206 48972 : computeJacobian(Moose::MortarType::Lower);
207 89628 : }
208 :
209 : void
210 15741 : 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 15741 : if (!_var)
215 6503 : return;
216 :
217 9238 : const auto sn = _sys.number();
218 9238 : const auto vn = _var->number();
219 :
220 : // If variable is nodal, zero DoFs based on inactive LM nodes
221 9238 : if (_var->isNodal())
222 : {
223 6058 : for (const auto node : inactive_lm_nodes)
224 : {
225 : // Allow mixed Lagrange orders between primal and LM
226 0 : if (!node->n_comp(sn, vn))
227 0 : continue;
228 :
229 0 : const auto dof_index = node->dof_number(sn, vn, 0);
230 : // No scaling; this is not physics
231 0 : if (_assembly.computingJacobian())
232 0 : addJacobianElement(
233 : _assembly, /*element_value=*/1, dof_index, dof_index, /*scaling_factor=*/1);
234 0 : if (_assembly.computingResidual())
235 : {
236 0 : const Real lm_value = _var->getNodalValue(*node);
237 0 : addResiduals(_assembly,
238 0 : std::array<Real, 1>{{lm_value}},
239 0 : 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 3540 : for (const auto el : inactive_lm_elems)
248 : {
249 360 : const auto n_comp = el->n_comp(sn, vn);
250 :
251 720 : for (const auto comp : make_range(n_comp))
252 : {
253 360 : const auto dof_index = el->dof_number(sn, vn, comp);
254 : // No scaling; this is not physics
255 360 : if (_assembly.computingJacobian())
256 120 : addJacobianElement(
257 : _assembly, /*element_value=*/1, dof_index, dof_index, /*scaling_factor=*/1);
258 360 : if (_assembly.computingResidual())
259 : {
260 240 : const Real lm_value = _var->getElementalValue(el, comp);
261 240 : addResiduals(_assembly,
262 0 : std::array<Real, 1>{{lm_value}},
263 240 : std::array<dof_id_type, 1>{{dof_index}},
264 : /*scaling_factor=*/1);
265 : }
266 : }
267 : }
268 : }
269 : }
|