idaholab/moose: contact coverage diff

Base 846a5c Head #33390 250e9c
Total Total +/- New
Rate 91.29% 91.35% +0.06% 100.00%
Hits 5021 5109 +88 117
Misses 479 484 +5 0
Filename Stmts Miss Cover
modules/contact/src/actions/ContactAction.C +3 -1 +0.14%
modules/contact/src/constraints/ComputeFrictionalForceLMMechanicalContact.C +1 0 +0.06%
modules/contact/src/constraints/ComputeWeightedGapLMMechanicalContact.C +5 0 +0.62%
modules/contact/src/constraints/NormalMortarMechanicalContact.C +7 +1 -2.94%
modules/contact/src/constraints/TangentialMortarMechanicalContact.C +4 +1 -1.85%
modules/contact/src/userobjects/CohesiveZoneModelBase.C +3 0 +0.09%
modules/contact/src/userobjects/LMWeightedGapUserObject.C +1 0 +0.37%
modules/contact/src/userobjects/LMWeightedVelocitiesUserObject.C +1 0 +0.15%
modules/contact/src/userobjects/PenaltyWeightedGapUserObject.C +19 0 +0.75%
modules/contact/src/userobjects/WeightedGapUserObject.C +42 +4 -2.80%
modules/contact/src/userobjects/WeightedVelocitiesUserObject.C +7 0 +0.41%
TOTAL +93 +5 +0.06%
code
coverage unchanged
code
coverage increased
code
coverage decreased
+
line added or modified

modules/contact/src/actions/ContactAction.C

868  
869  
870  
871 +
872 +
873 +
874  
875  
876  

    params.set<MooseEnum>("mortar_3d_subpatch_plane") =
        getParam<MooseEnum>("mortar_3d_subpatch_plane");
    const bool augmented_penalty =
        _formulation == ContactFormulation::MORTAR_PENALTY &&
        dynamic_cast<AugmentedLagrangianContactProblemInterface *>(_problem.get());
    params.set<bool>("ghost_point_neighbors") = !_mortar_dynamics && !augmented_penalty;
    addRelationshipManagers(input_rm_type, params);
  }
  else

1160  
1161  
1162  
1163  
1164  
1165  
1166  

      uo_params.set<bool>("correct_edge_dropping") = getParam<bool>("correct_edge_dropping");
      uo_params.set<std::vector<VariableName>>("disp_y") = {displacements[1]};
      if (ndisp > 2)
        uo_params.set<std::vector<VariableName>>("disp_z") = {displacements[2]};

      uo_params.set<VariableName>("secondary_variable") = displacements[0];
      uo_params.set<bool>("use_displaced_mesh") = true;

modules/contact/src/constraints/ComputeFrictionalForceLMMechanicalContact.C

73  
74  
75  
76 +
77  
78  
79  

    _3d(_has_disp_z)

{
  _weighted_velocities_uo.includeNodalNormalDerivatives();

  if (parameters.isParamSetByUser("mu") && _has_friction_function)
    paramError(

modules/contact/src/constraints/ComputeWeightedGapLMMechanicalContact.C

63  
64  
65  
66 +
67 +
68  
69  
70  

      "the value of c effectively depends on element size since in the constraint we compare nodal "
      "Lagrange Multiplier values to integrated gap values (LM nodal value is independent of "
      "element size, where integrated values are dependent on element size).");
  params.set<bool>("ghost_point_neighbors") = true;
  params.suppressParameter<bool>("ghost_point_neighbors");
  params.set<bool>("use_displaced_mesh") = true;
  params.set<bool>("interpolate_normals") = false;
  params.addRequiredParam<UserObjectName>("weighted_gap_uo", "The weighted gap user object");

89  
90  
91  
92 +
93  
94  
95  
96  
97  
98  
99 +
100 +
101  
102  
103  

    _disp_z_var(_has_disp_z ? getVar("disp_z", 0) : nullptr),
    _weighted_gap_uo(getUserObject<WeightedGapUserObject>("weighted_gap_uo"))
{
  _weighted_gap_uo.includeNodalNormalDerivatives();

  if (!getParam<bool>("use_displaced_mesh"))
    paramError(
        "use_displaced_mesh",
        "'use_displaced_mesh' must be true for the ComputeWeightedGapLMMechanicalContact object");

  if (getParam<bool>("interpolate_normals"))
    paramError("interpolate_normals",
               "Mechanical mortar contact uses normalized secondary nodal normals and cannot be "
               "combined with quadrature-point normal interpolation.");

modules/contact/src/constraints/NormalMortarMechanicalContact.C

31  
32  
33  
34 +
35  
36 +
37 +
38  
39  
40  

NormalMortarMechanicalContact::NormalMortarMechanicalContact(const InputParameters & parameters)
  : ADMortarLagrangeConstraint(parameters),
    _component(getParam<MooseEnum>("component")),
    _weighted_gap_uo(getUserObject<WeightedGapUserObject>("weighted_gap_uo"))
{
  if (getParam<bool>("interpolate_normals"))
    paramError("interpolate_normals",
               "Mechanical mortar contact uses normalized secondary nodal normals and cannot be "
               "combined with quadrature-point normal interpolation.");
}

57  
58  
59  
60 +
61 +
62 +
63  
64 +
65  
66  
67  

      {
        const auto normal_index = libmesh_map_find(_secondary_ip_lowerd_map, _i);
        const auto normal =
            _weighted_gap_uo.usesNodalNormalDerivatives()
                ? _weighted_gap_uo.contactNormal(*_lower_secondary_elem, normal_index)
                : ADRealVectorValue(_normals[normal_index]);
        return _test_secondary[_i][_qp] * _weighted_gap_uo.contactPressure()[_qp] *
               normal(_component);
      }

    case Moose::MortarType::Primary:

70  
71  
72  
73 +
74 +
75 +
76  
77 +
78  
79  
80  

      {
        const auto normal_index = libmesh_map_find(_primary_ip_lowerd_map, _i);
        const auto normal =
            _weighted_gap_uo.usesNodalNormalDerivatives()
                ? _weighted_gap_uo.contactNormal(*_lower_secondary_elem, normal_index)
                : ADRealVectorValue(_normals[normal_index]);
        return -_test_primary[_i][_qp] * _weighted_gap_uo.contactPressure()[_qp] *
               normal(_component);
      }
    default:
      return 0;

modules/contact/src/constraints/TangentialMortarMechanicalContact.C

41  
42  
43  
44 +
45  
46 +
47 +
48  
49  
50  

  : ADMortarLagrangeConstraint(parameters),
    _component(getParam<MooseEnum>("component")),
    _direction(getParam<MooseEnum>("direction")),
    _weighted_velocities_uo(getUserObject<WeightedVelocitiesUserObject>("weighted_velocities_uo"))
{
  if (getParam<bool>("interpolate_normals"))
    paramError("interpolate_normals",
               "Mechanical mortar contact uses tangents derived from normalized secondary nodal "
               "normals and cannot be combined with quadrature-point normal interpolation.");
}

52  
53  
54  
55 +
56  
57 +
58  
59  
60 +
61 +
62  
63  
64 +
65 +
66 +
67 +
68  
69  
70  

ADReal
TangentialMortarMechanicalContact::computeQpResidual(Moose::MortarType type)
{
  const auto tangent_component = [this](const unsigned int geometry_index)
  {
    if (_weighted_velocities_uo.usesNodalNormalDerivatives())
    {
      const auto tangents =
          _weighted_velocities_uo.contactTangents(*_lower_secondary_elem, geometry_index);
      return tangents[_direction](_component) / tangents[_direction].norm();
    }

    const auto & tangents = amg().getNodalTangents(*_lower_secondary_elem);
    return ADReal(tangents[_direction][geometry_index](_component) /
                  tangents[_direction][geometry_index].norm());
  };

  MooseEnum direction("direction_1 direction_2", "direction_1");


85  
86  
87  
88 +
89 +
90  
91  
92  
93 +
94 +
95  
96  
97  

      // means we want the residual to be negative in that case. So the sign of this residual should
      // be the same as the sign of lambda
      {
        const auto geometry_index = libmesh_map_find(_secondary_ip_lowerd_map, _i);
        return _test_secondary[_i][_qp] * tangential_pressure * tangent_component(geometry_index);
      }
    case Moose::MortarType::Primary:
    {
      const auto geometry_index = libmesh_map_find(_primary_ip_lowerd_map, _i);
      return -_test_primary[_i][_qp] * tangential_pressure * tangent_component(geometry_index);
    }
    default:
      return 0;

modules/contact/src/userobjects/CohesiveZoneModelBase.C

44  
45  
46  
47 +
48 +
49 +
50  
51  
52  

  params.suppressParameter<Real>("penalty_multiplier");
  params.suppressParameter<Real>("penetration_tolerance");
  params.setParameters("penalty", 0.0);
  params.set<bool>("allow_nodal_normal_derivatives") = false;
  params.set<bool>("use_nodal_normal_derivatives") = false;
  params.set<bool>("ghost_point_neighbors") = false;
  return params;
}

modules/contact/src/userobjects/LMWeightedGapUserObject.C

33  
34  
35  
36 +
37  
38  
39  

  InputParameters params = WeightedGapUserObject::validParams();
  params.addClassDescription(
      "Provides the mortar normal Lagrange multiplier for constraint enforcement.");
  params.set<bool>("allow_nodal_normal_derivatives") = true;
  params += LMWeightedGapUserObject::newParams();
  return params;
}

modules/contact/src/userobjects/LMWeightedVelocitiesUserObject.C

20  
21  
22  
23 +
24  
25  
26  

  params += LMWeightedGapUserObject::newParams();
  params.addClassDescription("Provides the mortar contact Lagrange multipliers (normal and "
                             "tangential) for constraint enforcement.");
  params.set<bool>("allow_nodal_normal_derivatives") = true;
  params.renameCoupledVar("lm_variable", "lm_variable_normal", "");
  params.addRequiredCoupledVar(
      "lm_variable_tangential_one",

modules/contact/src/userobjects/PenaltyWeightedGapUserObject.C

20  
21  
22  
23 +
24 +
25 +
26  
27  
28  
29 +
30  
31 +
32 +
33 +
34 +
35 +
36 +
37 +
38 +
39 +
40  
41  
42 +
43 +
44 +
45 +
46 +
47  
48  
49  
50 +
51  
52  
53  

PenaltyWeightedGapUserObject::validParams()
{
  InputParameters params = WeightedGapUserObject::validParams();
  params.set<bool>("allow_nodal_normal_derivatives") = true;
  params.set<bool>("use_nodal_normal_derivatives") = true;
  params.suppressParameter<bool>("ghost_point_neighbors");
  // Penalty contact has no LM constraint to contribute the nodal-normal one-ring to the matrix
  // graph, so its weighted-gap object supplies the coupling relationship manager directly.
  const auto configure_one_ring =
      [](const InputParameters & obj_params, InputParameters & rm_params)
  {
    rm_params.set<bool>("use_displaced_mesh") = obj_params.get<bool>("use_displaced_mesh");
    rm_params.set<BoundaryName>("secondary_boundary") =
        obj_params.get<BoundaryName>("secondary_boundary");
    rm_params.set<BoundaryName>("primary_boundary") =
        obj_params.get<BoundaryName>("primary_boundary");
    rm_params.set<SubdomainName>("secondary_subdomain") =
        obj_params.get<SubdomainName>("secondary_subdomain");
    rm_params.set<SubdomainName>("primary_subdomain") =
        obj_params.get<SubdomainName>("primary_subdomain");
    // penetration_tolerance is required exactly for augmented-Lagrange contact, which uses
    // frozen-normal geometry.
    rm_params.set<bool>("ghost_point_neighbors") =
        obj_params.get<bool>("use_nodal_normal_derivatives") &&
        !obj_params.isParamValid("penetration_tolerance");
  };
  params.addRelationshipManager("AugmentSparsityOnInterface",
                                Moose::RelationshipManagerType::GEOMETRIC |
                                    Moose::RelationshipManagerType::ALGEBRAIC,
                                configure_one_ring);
  params.addRelationshipManager(
      "AugmentSparsityOnInterface", Moose::RelationshipManagerType::COUPLING, configure_one_ring);
  params.addClassDescription("Computes the mortar normal contact force via a penalty approach.");
  params.addRequiredParam<Real>("penalty", "The penalty factor");

modules/contact/src/userobjects/WeightedGapUserObject.C

31  
32  
33  
34 +
35 +
36  
37  
38  

  params.addCoupledVar("disp_z", "The z displacement variable");
  params.set<bool>("use_displaced_mesh") = true;
  params.set<bool>("interpolate_normals") = false;
  params.addPrivateParam<bool>("allow_nodal_normal_derivatives", false);
  params.addPrivateParam<bool>("use_nodal_normal_derivatives", false);
  params.set<ExecFlagEnum>("execute_on") = {EXEC_LINEAR, EXEC_NONLINEAR};
  params.suppressParameter<ExecFlagEnum>("execute_on");
  return params;

52  
53  
54  
55 +
56 +
57 +
58  
59  
60  

    _primary_disp_y(_disp_y_var->adSlnNeighbor()),
    _secondary_disp_z(_has_disp_z ? &_disp_z_var->adSln() : nullptr),
    _primary_disp_z(_has_disp_z ? &_disp_z_var->adSlnNeighbor() : nullptr),
    _coord(_assembly.mortarCoordTransformation()),
    _allow_nodal_normal_derivatives(getParam<bool>("allow_nodal_normal_derivatives")),
    _use_nodal_normal_derivatives(false)
{
  if (!getParam<bool>("use_displaced_mesh"))
    paramError("use_displaced_mesh",

62  
63  
64  
65 +
66 +
67 +
68  
69  
70 +
71  
72 +
73 +
74  
75 +
76 +
77  
78 +
79 +
80  
81  
82  
83  
84  
85 +
86 +
87  
88 +
89 +
90  
91 +
92 +
93  
94  
95  
96  
97 +
98  
99  
100  
101 +
102  
103 +
104  
105  
106  


  // A valid penetration_tolerance identifies augmented-Lagrange contact, which retains frozen
  // normal directions.
  if (getParam<bool>("use_nodal_normal_derivatives") && !isParamValid("penetration_tolerance"))
    includeNodalNormalDerivatives();
}

void
WeightedGapUserObject::includeNodalNormalDerivatives() const
{
  if (_use_nodal_normal_derivatives)
    return;

  if (!_allow_nodal_normal_derivatives)
    mooseError("Nodal-normal derivatives are not supported by user object '", name(), "'.");

  if (getParam<bool>("interpolate_normals"))
    paramError("interpolate_normals",
               "Nodal-normal derivatives require normalized secondary nodal normals and cannot be "
               "combined with quadrature-point normal interpolation.");

  const std::array<std::pair<const MooseVariable *, const char *>, 3> displacement_variables{
      {{_disp_x_var, "disp_x"}, {_disp_y_var, "disp_y"}, {_disp_z_var, "disp_z"}}};
  for (const auto & [variable, parameter_name] : displacement_variables)
    if (variable)
    {
      if (!variable->isNodal())
        paramError(parameter_name,
                   "Nodal-normal derivatives require a nodal displacement variable.");
      if (&variable->sys() != &_sys)
        paramError(parameter_name,
                   "Nodal-normal derivatives require displacement variables in the nonlinear "
                   "system assembled by this contact object.");
    }

  _use_nodal_normal_derivatives = true;
}

bool
WeightedGapUserObject::usesNodalNormalDerivatives() const
{
  return _use_nodal_normal_derivatives && Moose::doDerivatives(_subproblem, _sys);
}

void

168  
169  
170  
171 +
172  
173  
174  

void
WeightedGapUserObject::computeQpIProperties()
{
  if (usesNodalNormalDerivatives())
    mooseAssert(_lower_secondary_elem->n_nodes() > _i,
                "Making sure that the nodal normal index is valid");
  else

179  
180  
181  
182 +
183 +
184  
185 +
186  
187  
188  


  auto & [weighted_gap, normalization] = _dof_to_weighted_gap[dof];

  if (usesNodalNormalDerivatives())
    weighted_gap += (*_test)[_i][_qp] * _qp_gap_nodal * contactNormal(*_lower_secondary_elem, _i);
  else
    weighted_gap += (*_test)[_i][_qp] * _qp_gap_nodal * _normals[_i];

  normalization += (*_test)[_i][_qp] * _qp_factor;


190  
191  
192  
193 +
194  
195  
196  

}

const ADRealVectorValue &
WeightedGapUserObject::contactNormal(const Elem & lower_secondary_elem,
                                     const unsigned int nodal_index) const
{
  mooseAssert(nodal_index < lower_secondary_elem.n_nodes(),

198  
199  
200  
201 +
202  
203 +
204 +
205  
206 +
207  
208 +
209 +
210 +
211  
212  
213  
214 +
215  
216  
217  
218  
219 +
220  
221  
222  
223 +
224 +
225  
226 +
227  
228 +
229 +
230 +
231  
232  
233 +
234  
235  
236  

  mooseAssert(usesNodalNormalDerivatives(),
              "AD contact normals should only be requested while recording nodal-normal "
              "derivatives.");
  const Node * const node = lower_secondary_elem.node_ptr(nodal_index);
  const auto normal_it = _ad_nodal_normals.find(node);
  if (normal_it != _ad_nodal_normals.end())
    return normal_it->second;

  amg().computeADNodalNormals(
      [this](const Node & coordinate_node, const Point & geometry_coordinate)
      { return nodalCoordinate(coordinate_node, geometry_coordinate); },
      _ad_nodal_normals);
  return libmesh_map_find(_ad_nodal_normals, node);
}

ADPoint
WeightedGapUserObject::nodalCoordinate(const Node & node, const Point & geometry_coordinate) const
{
  mooseAssert(usesNodalNormalDerivatives(),
              "AD nodal coordinates should only be requested while recording nodal-normal "
              "derivatives.");
  ADPoint point = geometry_coordinate;

  const std::array<std::pair<const MooseVariable *, unsigned int>, 3> displacement_variables{
      {{_disp_x_var, 0}, {_disp_y_var, 1}, {_disp_z_var, 2}}};
  for (const auto & [variable, component] : displacement_variables)
    if (variable)
    {
      const auto sys_num = variable->sys().number();
      const auto var_num = variable->number();
      if (node.n_dofs(sys_num, var_num))
        Moose::derivInsert(
            point(component).derivatives(), node.dof_number(sys_num, var_num, 0), 1.0);
    }

  return point;
}

void

modules/contact/src/userobjects/WeightedVelocitiesUserObject.C

55  
56  
57  
58 +
59  
60  
61  

}

std::array<ADRealVectorValue, 2>
WeightedVelocitiesUserObject::contactTangents(const Elem & lower_secondary_elem,
                                              const unsigned int nodal_index) const
{
  mooseAssert(nodal_index < lower_secondary_elem.n_nodes(),

63  
64  
65  
66 +
67  
68  
69  

  mooseAssert(usesNodalNormalDerivatives(),
              "AD contact tangents should only be requested while recording nodal-normal "
              "derivatives.");
  return Moose::Mortar::householderTangents(contactNormal(lower_secondary_elem, nodal_index));
}

void

122  
123  
124  
125 +
126  
127 +
128 +
129  
130  
131 +
132  
133 +
134  
135  
136  
137  
138 +
139  
140  
141 +
142  
143 +
144  
145  
146 +
147  
148  
149 +
150  
151  
152 +
153  
154  
155  

      _is_weighted_gap_nodal ? static_cast<const DofObject *>(_lower_secondary_elem->node_ptr(_i))
                             : static_cast<const DofObject *>(_lower_secondary_elem);

  if (usesNodalNormalDerivatives())
  {
    const auto tangents = contactTangents(*_lower_secondary_elem, _i);
    for (const auto direction : make_range(_3d ? 2 : 1))
    {
      _dof_to_weighted_tangential_velocity[dof][direction] +=
          (*_test)[_i][_qp] * _qp_tangential_velocity_nodal * tangents[direction];
      _dof_to_real_tangential_velocity[dof][direction] +=
          (*_test)[_i][_qp] * _qp_real_tangential_velocity_nodal * tangents[direction];
    }
  }
  else
  {
    const auto & nodal_tangents = amg().getNodalTangents(*_lower_secondary_elem);

    _dof_to_weighted_tangential_velocity[dof][0] +=
        (*_test)[_i][_qp] * _qp_tangential_velocity_nodal * nodal_tangents[0][_i];
    _dof_to_real_tangential_velocity[dof][0] +=
        (*_test)[_i][_qp] * _qp_real_tangential_velocity_nodal * nodal_tangents[0][_i];

    // Get the _dof_to_weighted_tangential_velocity map for a second direction
    if (_3d)
    {
      _dof_to_weighted_tangential_velocity[dof][1] +=
          (*_test)[_i][_qp] * _qp_tangential_velocity_nodal * nodal_tangents[1][_i];

      _dof_to_real_tangential_velocity[dof][1] +=
          (*_test)[_i][_qp] * _qp_real_tangential_velocity_nodal * nodal_tangents[1][_i];
    }
  }
}