16 #include "metaphysicl/dualsemidynamicsparsenumberarray.h" 17 #include "metaphysicl/parallel_dualnumber.h" 18 #include "metaphysicl/parallel_dynamic_std_array_wrapper.h" 19 #include "metaphysicl/parallel_semidynamicsparsenumberarray.h" 33 "The frictional Lagrange's multiplier for an addtional direction.");
36 "Coupled function to evaluate friction with values from contact pressure and relative " 37 "tangential velocities");
38 params.
addParam<
Real>(
"c_t", 1e0,
"Numerical parameter for tangential constraints");
42 "Minimum value of contact pressure that will trigger frictional enforcement");
44 "mu",
"mu > 0",
"The friction coefficient for the Coulomb friction law");
46 "The weighted tangential velocities user object.");
55 _c_t(getParam<
Real>(
"c_t")),
56 _secondary_x_dot(_secondary_var.adUDot()),
57 _primary_x_dot(_primary_var.adUDotNeighbor()),
58 _secondary_y_dot(adCoupledDot(
"disp_y")),
59 _primary_y_dot(adCoupledNeighborValueDot(
"disp_y")),
60 _secondary_z_dot(_has_disp_z ? &adCoupledDot(
"disp_z") : nullptr),
61 _primary_z_dot(_has_disp_z ? &adCoupledNeighborValueDot(
"disp_z") : nullptr),
62 _epsilon(getParam<
Real>(
"epsilon")),
63 _mu(isParamValid(
"function_friction") ?
std::numeric_limits<double>::quiet_NaN()
64 : getParam<
Real>(
"mu")),
65 _function_friction(isParamValid(
"function_friction") ? &getFunction(
"function_friction")
67 _has_friction_function(isParamValid(
"function_friction")),
74 "Please only provide friction either as a function or as a constant value, but not both.");
76 paramError(
"mu",
"Please provide a value or a function for the coefficient of friction.");
78 if (!getParam<bool>(
"use_displaced_mesh"))
80 "'use_displaced_mesh' must be true for the " 81 "ComputeFrictionalForceLMMechanicalContact object");
85 "Three-dimensional mortar frictional contact simulations require an additional " 86 "frictional Lagrange's multiplier to enforce a second tangential pressure");
97 "Frictional contact constraints only support elemental variables of CONSTANT order");
118 const auto & dof_to_weighted_tangential_velocity =
121 const std::unordered_map<const DofObject *, std::pair<ADReal, Real>> & dof_to_weighted_gap =
126 for (
const auto & [dof_object, weighted_velocities_pr] : dof_to_weighted_tangential_velocity)
131 const auto & [weighted_gap_pr, normalization] =
132 libmesh_map_find(dof_to_weighted_gap, dof_object);
149 const std::unordered_set<const Node *> & inactive_lm_nodes)
151 const auto & dof_to_weighted_tangential_velocity =
155 for (
const auto & [dof_object, weighted_velocities_pr] : dof_to_weighted_tangential_velocity)
158 if ((inactive_lm_nodes.find(static_cast<const Node *>(dof_object)) !=
159 inactive_lm_nodes.end()) ||
190 std::array<dof_id_type, 2> friction_dof_indices;
191 std::array<ADReal, 2> friction_lm_values;
193 const unsigned int num_tangents = 2;
216 dof_residual = friction_lm_values[0];
217 dof_residual_dir = friction_lm_values[1];
222 const Real epsilon_sqrt = 1.0e-48;
224 const auto lamdba_plus_cg = contact_pressure +
c * weighted_gap;
225 std::array<ADReal, 2> lambda_t_plus_ctu;
226 lambda_t_plus_ctu[0] = friction_lm_values[0] + c_t * *tangential_vel[0] *
_dt;
227 lambda_t_plus_ctu[1] = friction_lm_values[1] + c_t * *tangential_vel[1] *
_dt;
229 const auto term_1_x =
230 std::max(mu_ad * lamdba_plus_cg,
231 std::sqrt(lambda_t_plus_ctu[0] * lambda_t_plus_ctu[0] +
232 lambda_t_plus_ctu[1] * lambda_t_plus_ctu[1] + epsilon_sqrt)) *
233 friction_lm_values[0];
235 const auto term_1_y =
236 std::max(mu_ad * lamdba_plus_cg,
237 std::sqrt(lambda_t_plus_ctu[0] * lambda_t_plus_ctu[0] +
238 lambda_t_plus_ctu[1] * lambda_t_plus_ctu[1] + epsilon_sqrt)) *
239 friction_lm_values[1];
241 const auto term_2_x = mu_ad * std::max(0.0, lamdba_plus_cg) * lambda_t_plus_ctu[0];
243 const auto term_2_y = mu_ad * std::max(0.0, lamdba_plus_cg) * lambda_t_plus_ctu[1];
245 dof_residual = term_1_x - term_2_x;
246 dof_residual_dir = term_1_y - term_2_y;
250 std::array<ADReal, 1>{{dof_residual}},
251 std::array<dof_id_type, 1>{{friction_dof_indices[0]}},
254 std::array<ADReal, 1>{{dof_residual_dir}},
255 std::array<dof_id_type, 1>{{friction_dof_indices[1]}},
287 dof_residual = friction_lm_value;
290 const auto term_1 = std::max(mu_ad * (contact_pressure +
c * weighted_gap),
291 std::abs(friction_lm_value + c_t * tangential_vel *
_dt)) *
293 const auto term_2 = mu_ad * std::max(0.0, contact_pressure +
c * weighted_gap) *
294 (friction_lm_value + c_t * tangential_vel *
_dt);
296 dof_residual = term_1 - term_2;
300 std::array<ADReal, 1>{{dof_residual}},
301 std::array<dof_id_type, 1>{{friction_dof_index}},
307 const ADReal & tangential_vel,
308 const ADReal & tangential_vel_dir)
317 ADReal tangential_vel_magnitude = std::sqrt(tangential_vel * tangential_vel +
318 tangential_vel_dir * tangential_vel_dir + 1.0e-24);
virtual const NumericVector< Number > *const & currentSolution() const=0
void addResidualsAndJacobian(Assembly &assembly, const Residuals &residuals, const Indices &dof_indices, Real scaling_factor)
unsigned int number() const
MooseVariable * getVar(const std::string &var_name, unsigned int comp)
Creates dof object to weighted tangential velocities map.
DualNumber< Real, DNDerivativeType, true > ADReal
bool isParamValid(const std::string &name) const
const std::unordered_map< const DofObject *, std::array< ADReal, 2 > > & dofToWeightedVelocities() const
Get the degree of freedom to weighted velocities information.
void paramError(const std::string ¶m, Args... args) const
unsigned int number() const
const std::unordered_map< const DofObject *, std::pair< ADReal, Real > > & dofToWeightedGap() const
Get the degree of freedom to weighted gap information.
MooseVariable *const _var
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
IntRange< T > make_range(T beg, T end)
const InputParameters & parameters() const
void derivInsert(SemiDynamicSparseNumberArray< Real, libMesh::dof_id_type, NWrapper< N >> &derivs, libMesh::dof_id_type index, Real value)
virtual Real value(Real t, const Point &p) const
processor_id_type processor_id() const