https://mooseframework.inl.gov
QuasiStaticSolidMechanicsPhysics.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 "Conversion.h"
11 #include "FEProblem.h"
12 #include "Factory.h"
13 #include "MooseMesh.h"
14 #include "MooseObjectAction.h"
16 #include "Material.h"
18 
19 #include "BlockRestrictable.h"
20 
22 #include "AddVariableAction.h"
23 
24 #include "libmesh/string_to_enum.h"
25 #include <algorithm>
26 
27 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "meta_action");
28 
29 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "setup_mesh_complete");
30 
31 registerMooseAction("SolidMechanicsApp",
33  "validate_coordinate_systems");
34 
35 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_variable");
36 
37 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_aux_variable");
38 
39 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_kernel");
40 
41 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_aux_kernel");
42 
43 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_material");
44 
45 registerMooseAction("SolidMechanicsApp",
47  "add_master_action_material");
48 
49 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_scalar_kernel");
50 
51 registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "add_user_object");
52 
55 {
57  params.addClassDescription("Set up stress divergence kernels with coordinate system aware logic");
58 
59  // parameters specified here only appear in the input file sub-blocks of the solid mechanics
60  // action, not in the common parameters area
61  params.addParam<std::vector<SubdomainName>>("block",
62  {},
63  "The list of ids of the blocks (subdomain) "
64  "that the stress divergence kernels will be "
65  "applied to");
66  params.addParamNamesToGroup("block", "Advanced");
67 
68  params.addParam<MultiMooseEnum>("additional_generate_output",
70  "Add scalar quantity output for stress and/or strain (will be "
71  "appended to the list in `generate_output`)");
72  params.addParam<MultiMooseEnum>(
73  "additional_material_output_order",
75  "Specifies the order of the FE shape function to use for this variable.");
76 
77  params.addParam<MultiMooseEnum>(
78  "additional_material_output_family",
80  "Specifies the family of FE shape functions to use for this variable.");
81 
82  params.addParamNamesToGroup("additional_generate_output additional_material_output_order "
83  "additional_material_output_family",
84  "Output");
85  params.addParam<std::string>(
86  "strain_base_name",
87  "The base name used for the strain. If not provided, it will be set equal to base_name");
88  params.addParam<std::vector<TagName>>(
89  "extra_vector_tags",
90  "The tag names for extra vectors that residual data should be saved into");
91  params.addParam<std::vector<TagName>>("absolute_value_vector_tags",
92  "The tag names for extra vectors that the absolute value "
93  "of the residual should be accumulated into");
94  params.addParam<Real>("scaling", "The scaling to apply to the displacement variables");
95  params.addParam<Point>(
96  "cylindrical_axis_point1",
97  "Starting point for direction of axis of rotation for cylindrical stress/strain.");
98  params.addParam<Point>(
99  "cylindrical_axis_point2",
100  "Ending point for direction of axis of rotation for cylindrical stress/strain.");
101  params.addParam<Point>("spherical_center_point",
102  "Center point of the spherical coordinate system.");
103  params.addParam<Point>("direction", "Direction stress/strain is calculated in");
104  params.addParam<bool>("automatic_eigenstrain_names",
105  false,
106  "Collects all material eigenstrains and passes to required strain "
107  "calculator within TMA internally.");
108 
109  // Homogenization system input
110  params.addParam<MultiMooseEnum>(
111  "constraint_types",
113  "Type of each constraint: strain, stress, or none. The types are specified in the "
114  "column-major order, and there must be 9 entries in total.");
115  params.addParam<std::vector<FunctionName>>(
116  "targets", {}, "Functions giving the targets to hit for constraint types that are not none.");
117 
118  params.addParamNamesToGroup("scaling", "Variables");
119  params.addParamNamesToGroup("strain_base_name automatic_eigenstrain_names", "Strain");
120  params.addParamNamesToGroup(
121  "cylindrical_axis_point1 cylindrical_axis_point2 spherical_center_point direction",
122  "Coordinate system");
123 
124  return params;
125 }
126 
129  _displacements(getParam<std::vector<VariableName>>("displacements")),
130  _ndisp(_displacements.size()),
131  _coupled_displacements(_ndisp),
132  _save_in(getParam<std::vector<AuxVariableName>>("save_in")),
133  _diag_save_in(getParam<std::vector<AuxVariableName>>("diag_save_in")),
134  _subdomain_names(getParam<std::vector<SubdomainName>>("block")),
135  _subdomain_ids(),
136  _strain(getParam<MooseEnum>("strain").getEnum<Strain>()),
137  _planar_formulation(getParam<MooseEnum>("planar_formulation").getEnum<PlanarFormulation>()),
138  _out_of_plane_direction(
139  getParam<MooseEnum>("out_of_plane_direction").getEnum<OutOfPlaneDirection>()),
140  _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
141  _material_output_order(getParam<MultiMooseEnum>("material_output_order")),
142  _material_output_family(getParam<MultiMooseEnum>("material_output_family")),
143  _cylindrical_axis_point1_valid(params.isParamSetByUser("cylindrical_axis_point1")),
144  _cylindrical_axis_point2_valid(params.isParamSetByUser("cylindrical_axis_point2")),
145  _direction_valid(params.isParamSetByUser("direction")),
146  _verbose(getParam<bool>("verbose")),
147  _spherical_center_point_valid(params.isParamSetByUser("spherical_center_point")),
148  _auto_eigenstrain(getParam<bool>("automatic_eigenstrain_names")),
149  _lagrangian_kernels(getParam<bool>("new_system")),
150  _lk_large_kinematics(_strain == Strain::Finite),
151  _lk_formulation(getParam<MooseEnum>("formulation").getEnum<LKFormulation>()),
152  _lk_locking(getParam<bool>("volumetric_locking_correction")),
153  _lk_homogenization(false),
154  _constraint_types(getParam<MultiMooseEnum>("constraint_types")),
155  _targets(getParam<std::vector<FunctionName>>("targets"))
156 {
157  // determine if incremental strains are to be used
158  if (isParamValid("incremental"))
159  {
160  const bool incremental = getParam<bool>("incremental");
161  if (!incremental && _strain == Strain::Small)
163  else if (!incremental && _strain == Strain::Finite)
165  else if (incremental && _strain == Strain::Small)
167  else if (incremental && _strain == Strain::Finite)
169  else
170  mooseError("Internal error");
171  }
172  else
173  {
174  if (_strain == Strain::Small)
175  {
177  mooseInfo("SolidMechanics Action: selecting 'total small strain' formulation. Use "
178  "`incremental = true` to select 'incremental small strain' instead.");
179  }
180  else if (_strain == Strain::Finite)
181  {
183  mooseInfo("SolidMechanics Action: selecting 'incremental finite strain' formulation.");
184  }
185  else
186  mooseError("Internal error");
187  }
188 
189  // determine if displaced mesh is to be used
191  if (params.isParamSetByUser("use_displaced_mesh") && !_lagrangian_kernels)
192  {
193  bool use_displaced_mesh_param = getParam<bool>("use_displaced_mesh");
194  if (use_displaced_mesh_param != _use_displaced_mesh && params.isParamSetByUser("strain"))
195  mooseError("Wrong combination of use displaced mesh and strain model");
196  _use_displaced_mesh = use_displaced_mesh_param;
197  }
198 
199  // convert vector of VariableName to vector of VariableName
200  for (unsigned int i = 0; i < _ndisp; ++i)
202 
203  if (_save_in.size() != 0 && _save_in.size() != _ndisp)
204  mooseError("Number of save_in variables should equal to the number of displacement variables ",
205  _ndisp);
206 
207  if (_diag_save_in.size() != 0 && _diag_save_in.size() != _ndisp)
208  mooseError(
209  "Number of diag_save_in variables should equal to the number of displacement variables ",
210  _ndisp);
211 
212  // plane strain consistency check
214  {
217  "moose",
218  27340,
219  "Planar formulations are not yet available through the Physics syntax with new_system = "
220  "true. They can be enabled by manually setting up the appropriate objects. Please refer "
221  "to the documentation and regression tests for examples.");
222  if (params.isParamSetByUser("out_of_plane_strain") &&
224  mooseError(
225  "out_of_plane_strain should only be specified with planar_formulation=WEAK_PLANE_STRESS");
226  else if (!params.isParamSetByUser("out_of_plane_strain") &&
228  mooseError("out_of_plane_strain must be specified with planar_formulation=WEAK_PLANE_STRESS");
229  }
230 
231  // convert output variable names to lower case
232  for (const auto & out : getParam<MultiMooseEnum>("generate_output"))
233  {
234  std::string lower(out);
235  std::transform(lower.begin(), lower.end(), lower.begin(), ::tolower);
236  _generate_output.push_back(lower);
237  }
238 
239  if (!_generate_output.empty())
241 
242  // Error if volumetric locking correction is true for 1D problems
243  if (_ndisp == 1 && getParam<bool>("volumetric_locking_correction"))
244  mooseError("Volumetric locking correction should be set to false for 1D problems.");
245 
246  if (!getParam<bool>("add_variables") && params.isParamSetByUser("scaling"))
247  paramError("scaling",
248  "The scaling parameter has no effect unless add_variables is set to true. Did you "
249  "mean to set 'add_variables = true'?");
250 
251  // Get cylindrical axis points if set by user
253  {
254  _cylindrical_axis_point1 = getParam<Point>("cylindrical_axis_point1");
255  _cylindrical_axis_point2 = getParam<Point>("cylindrical_axis_point2");
256  }
257 
258  // Get spherical center point if set by user
260  _spherical_center_point = getParam<Point>("spherical_center_point");
261 
262  // Get direction for tensor component if set by user
263  if (_direction_valid)
264  _direction = getParam<Point>("direction");
265 
266  // Get eigenstrain names if passed by user
267  _eigenstrain_names = getParam<std::vector<MaterialPropertyName>>("eigenstrain_names");
268 
269  // Determine if we're going to use the homogenization system for the new
270  // lagrangian kernels
271  bool ctype_set = params.isParamSetByUser("constraint_types");
272  bool targets_set = params.isParamSetByUser("targets");
273  if (ctype_set || targets_set)
274  {
275  if (!(ctype_set && targets_set))
276  mooseError("To use the Lagrangian kernel homogenization system you "
277  "most provide both the constraint_types and the targets "
278  "parameters");
279  _lk_homogenization = true;
280  // Do consistency checking on the kernel options
282  mooseError("The Lagrangian kernel homogenization system requires the "
283  "use of formulation = TOTAL");
284  }
285 
286  // Cross check options to weed out incompatible choices for the new lagrangian
287  // kernel system
289  {
290  if (_use_ad)
291  mooseError("The Lagrangian kernel system is not yet compatible with AD. "
292  "Do not set the use_automatic_differentiation flag.");
293 
294  if (params.isParamSetByUser("use_finite_deform_jacobian"))
295  mooseError("The Lagrangian kernel system always produces the exact "
296  "Jacobian. use_finite_deform_jacobian is redundant and "
297  " should not be set");
298  if (params.isParamSetByUser("global_strain"))
299  mooseError("The Lagrangian kernel system is not compatible with "
300  "the global_strain option. Use the homogenization "
301  " system instead");
302  if (params.isParamSetByUser("decomposition_method"))
303  mooseError("The decomposition_method parameter should not be used "
304  " with the Lagrangian kernel system. Similar options "
305  " for native small deformation material models are "
306  " available as part of the ComputeLagrangianStress "
307  " material system.");
308  }
309  else
310  {
311  if (params.isParamSetByUser("formulation"))
312  mooseError("The StressDiveregenceTensor system always uses an "
313  " updated Lagrangian formulation. Do not set the "
314  " formulation parameter, it is only used with the "
315  " new Lagrangian kernel system.");
316  if (_lk_homogenization)
317  mooseError("The homogenization system can only be used with the "
318  "new Lagrangian kernels");
319  }
320 }
321 
322 void
324 {
325  std::string ad_prepend = "";
326  if (_use_ad)
327  ad_prepend = "AD";
328 
329  // Consistency checks across subdomains
331 
332  // Gather info from all other QuasiStaticSolidMechanicsPhysics
334 
335  // Deal with the optional AuxVariable based tensor quantity output
337 
338  // Meta action which optionally spawns other actions
339  if (_current_task == "meta_action")
340  {
342  {
343  if (_use_ad)
344  paramError("use_automatic_differentiation", "AD not setup for use with PlaneStrain");
345  // Set the action parameters
346  const std::string type = "GeneralizedPlaneStrainAction";
347  auto action_params = _action_factory.getValidParams(type);
348  action_params.set<bool>("_built_by_moose") = true;
349  action_params.set<std::string>("registered_identifier") = "(AutoBuilt)";
350 
351  // Skipping selected parameters in applyParameters() and then manually setting them only if
352  // they are set by the user is just to prevent both the current and deprecated variants of
353  // these parameters from both getting passed to the UserObject. Once we get rid of the
354  // deprecated versions, we can just set them all with applyParameters().
355  action_params.applyParameters(parameters(),
356  {"use_displaced_mesh",
357  "out_of_plane_pressure",
358  "out_of_plane_pressure_function",
359  "factor",
360  "pressure_factor"});
361  action_params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
362 
363  if (parameters().isParamSetByUser("out_of_plane_pressure"))
364  action_params.set<FunctionName>("out_of_plane_pressure") =
365  getParam<FunctionName>("out_of_plane_pressure");
366  if (parameters().isParamSetByUser("out_of_plane_pressure_function"))
367  action_params.set<FunctionName>("out_of_plane_pressure_function") =
368  getParam<FunctionName>("out_of_plane_pressure_function");
369  if (parameters().isParamSetByUser("factor"))
370  action_params.set<Real>("factor") = getParam<Real>("factor");
371  if (parameters().isParamSetByUser("pressure_factor"))
372  action_params.set<Real>("pressure_factor") = getParam<Real>("pressure_factor");
373 
374  // Create and add the action to the warehouse
375  auto action = MooseSharedNamespace::static_pointer_cast<MooseObjectAction>(
376  _action_factory.create(type, name() + "_gps", action_params));
377  _awh.addActionBlock(action);
378  }
379  }
380 
381  // Add variables
382  else if (_current_task == "add_variable")
383  {
384  // Add variables here only if the CommonSolidMechanicsAction does not exist.
385  // This happens notably if the QuasiStaticSolidMechanics was created by a meta_action
386  const auto common_actions = _awh.getActions<CommonSolidMechanicsAction>();
387  if (common_actions.empty() && getParam<bool>("add_variables"))
388  {
389  auto params = _factory.getValidParams("MooseVariable");
390  // determine necessary order
391  const bool second = _problem->mesh().hasSecondOrderElements();
392 
393  params.set<MooseEnum>("order") = second ? "SECOND" : "FIRST";
394  params.set<MooseEnum>("family") = "LAGRANGE";
395  if (isParamValid("scaling"))
396  params.set<std::vector<Real>>("scaling") = {getParam<Real>("scaling")};
397 
398  // Note how we do not add the block restriction because BISON's meta-actions
399  // currently rely on them not being added.
400 
401  // Loop through the displacement variables
402  for (const auto & disp : _displacements)
403  {
404  // Create displacement variables
405  _problem->addVariable("MooseVariable", disp, params);
406  }
407  }
408 
409  // Homogenization scalar
410  if (_lk_homogenization)
411  {
412  InputParameters params = _factory.getValidParams("MooseVariable");
413  const std::map<bool, std::vector<unsigned int>> mg_order_max{{true, {1, 4, 9}},
414  {false, {1, 3, 6}}};
415  std::size_t mg_order = 0;
416  for (auto i : index_range(_constraint_types))
417  {
418  const auto ctype = static_cast<Homogenization::ConstraintType>(_constraint_types.get(i));
420  mg_order++;
421  }
422  if (mg_order > mg_order_max.at(_lk_large_kinematics)[_ndisp - 1])
423  paramError("constraint_types",
424  "Number of non-none constraint types must not be greater than ",
425  mg_order_max.at(_lk_large_kinematics)[_ndisp - 1],
426  ", but ",
427  mg_order,
428  " are provided.");
429  params.set<MooseEnum>("family") = "SCALAR";
430  params.set<MooseEnum>("order") = mg_order;
431  auto fe_type = AddVariableAction::feType(params);
432  auto var_type = AddVariableAction::variableType(fe_type);
433  _problem->addVariable(var_type, _hname, params);
434  }
435  }
436  // Add Materials
437  else if (_current_task == "add_master_action_material")
438  {
439  // Automatic eigenstrain names
440  if (_auto_eigenstrain)
442 
443  // Easiest just to branch on type here, as the strain systems are completely
444  // different
447  else
449  }
450 
451  // Add Stress Divergence (and optionally WeakPlaneStress) Kernels
452  else if (_current_task == "add_kernel")
453  {
454  for (unsigned int i = 0; i < _ndisp; ++i)
455  {
456  auto tensor_kernel_type = getKernelType();
457  auto params = getKernelParameters(ad_prepend + tensor_kernel_type);
458 
459  std::string kernel_name = "TM_" + name() + Moose::stringify(i);
460 
461  // Set appropriate components for kernels, including in the cases where a planar model is
462  // running in planes other than the x-y plane (defined by _out_of_plane_strain_direction).
464  continue;
465  else if (_out_of_plane_direction == OutOfPlaneDirection::y && i == 1)
466  continue;
467 
468  params.set<unsigned int>("component") = i;
469 
470  params.set<NonlinearVariableName>("variable") = _displacements[i];
471 
472  if (_save_in.size() == _ndisp)
473  params.set<std::vector<AuxVariableName>>("save_in") = {_save_in[i]};
474  if (_diag_save_in.size() == _ndisp)
475  params.set<std::vector<AuxVariableName>>("diag_save_in") = {_diag_save_in[i]};
476  if (isParamValid("out_of_plane_strain") && !_lagrangian_kernels)
477  params.set<std::vector<VariableName>>("out_of_plane_strain") = {
478  getParam<VariableName>("out_of_plane_strain")};
479 
480  if (_lk_homogenization)
481  {
482  params.set<std::vector<VariableName>>("scalar_variable") = {_hname};
483  params.set<MultiMooseEnum>("constraint_types") = _constraint_types;
484  params.set<std::vector<FunctionName>>("targets") = _targets;
485  }
486 
487  _problem->addKernel(ad_prepend + tensor_kernel_type, kernel_name, params);
488  }
489 
491  {
492  auto params = getKernelParameters(ad_prepend + "WeakPlaneStress");
493  std::string wps_kernel_name = "TM_WPS_" + name();
494  params.set<NonlinearVariableName>("variable") = getParam<VariableName>("out_of_plane_strain");
495 
496  _problem->addKernel(ad_prepend + "WeakPlaneStress", wps_kernel_name, params);
497  }
498  }
499 }
500 
501 void
503 {
504  // Do the coordinate system check only once the problem is created
505  if (_current_task == "setup_mesh_complete")
506  {
507  // get subdomain IDs
508  for (auto & name : _subdomain_names)
509  {
510  auto id = _mesh->getSubdomainID(name);
511  if (id == Moose::INVALID_BLOCK_ID)
512  paramError("block", "Subdomain \"" + name + "\" not found in mesh.");
513  else
514  _subdomain_ids.insert(id);
515  }
516  }
517 
518  if (_current_task == "validate_coordinate_systems")
519  {
520  // use either block restriction list or list of all subdomains in the mesh
521  const auto & check_subdomains =
522  _subdomain_ids.empty() ? _problem->mesh().meshSubdomains() : _subdomain_ids;
523  if (check_subdomains.empty())
524  mooseError("No subdomains found");
525 
526  // make sure all subdomains are using the same coordinate system
527  _coord_system = _problem->getCoordSystem(*check_subdomains.begin());
528  for (auto subdomain : check_subdomains)
529  if (_problem->getCoordSystem(subdomain) != _coord_system)
530  mooseError("The SolidMechanics action requires all subdomains to have the same coordinate "
531  "system.");
532 
534  {
536  mooseError("'out_of_plane_direction' must be 'z' for axisymmetric simulations");
537  }
539  {
541  mooseError(
542  "Must specify two displacements for plane strain when the out of plane direction is z");
544  mooseError("Must specify three displacements for plane strain when the out of plane "
545  "direction is x or y");
546  }
547  }
548 }
549 
550 void
552 {
553  if (_current_task == "add_material")
555 
556  // Add variables (optional)
557  if (_current_task == "add_aux_variable")
558  {
559  unsigned int index = 0;
560  for (auto out : _generate_output)
561  {
562  const auto & order = _material_output_order[index];
563  const auto & family = _material_output_family[index];
564 
565  std::string type = (order == "CONSTANT" && family == "MONOMIAL")
566  ? "MooseVariableConstMonomial"
567  : "MooseVariable";
568 
569  // Create output helper aux variables
570  auto params = _factory.getValidParams(type);
571  params.set<MooseEnum>("order") = order;
572  params.set<MooseEnum>("family") = family;
573 
574  if (family == "MONOMIAL")
575  _problem->addAuxVariable(type, _base_name + out, params);
576  else
577  _problem->addVariable(type, _base_name + out, params);
578 
579  index++;
580  }
581  }
582 
583  // Add output AuxKernels
584  else if (_current_task == "add_aux_kernel")
585  {
586  std::string ad_prepend = _use_ad ? "AD" : "";
587  // Loop through output aux variables
588  unsigned int index = 0;
589  for (auto out : _generate_output)
590  {
591  if (_material_output_family[index] == "MONOMIAL")
592  {
594 
595  params = _factory.getValidParams("MaterialRealAux");
596  params.applyParameters(parameters());
597  params.set<MaterialPropertyName>("property") = _base_name + out;
598  params.set<AuxVariableName>("variable") = _base_name + out;
599  params.set<ExecFlagEnum>("execute_on") = EXEC_TIMESTEP_END;
600 
601  _problem->addAuxKernel(
602  ad_prepend + "MaterialRealAux", _base_name + out + '_' + name(), params);
603  }
604  index++;
605  }
606  }
607  else if (_current_task == "add_kernel")
608  {
609  std::string ad_prepend = _use_ad ? "AD" : "";
610  // Loop through output aux variables
611  unsigned int index = 0;
612  for (auto out : _generate_output)
613  {
614  if (_material_output_family[index] != "MONOMIAL")
615  {
617 
618  params = _factory.getValidParams("MaterialPropertyValue");
619  params.applyParameters(parameters());
620  params.set<MaterialPropertyName>("prop_name") = _base_name + out;
621  params.set<NonlinearVariableName>("variable") = _base_name + out;
622 
623  _problem->addKernel(
624  ad_prepend + "MaterialPropertyValue", _base_name + out + '_' + name(), params);
625  }
626  index++;
627  }
628  }
629 }
630 
631 void
633 {
634  // Create containers for collecting blockIDs and eigenstrain names from materials
635  std::map<std::string, std::set<SubdomainID>> material_eigenstrain_map;
636  std::set<std::string> eigenstrain_set;
637 
638  std::set<MaterialPropertyName> verified_eigenstrain_names;
639 
640  std::map<std::string, std::string> remove_add_map;
641  std::set<std::string> remove_reduced_set;
642 
643  // Loop over all the materials(eigenstrains) already created
644  auto materials = _problem->getMaterialWarehouse().getObjects();
645  for (auto & mat : materials)
646  {
647  std::shared_ptr<BlockRestrictable> blk = std::dynamic_pointer_cast<BlockRestrictable>(mat);
648  const InputParameters & mat_params = mat->parameters();
649  auto & mat_name = mat->type();
650 
651  // Check for eigenstrain names, only deal with those materials
652  if (mat_params.isParamValid("eigenstrain_name"))
653  {
654  std::shared_ptr<MaterialData> mat_dat;
655  auto name = mat_params.get<std::string>("eigenstrain_name");
656 
657  // Check for base_name prefix
658  if (mat_params.isParamValid("base_name"))
659  name = mat_params.get<std::string>("base_name") + '_' + name;
660 
661  // Check block restrictions
662  if (!blk)
663  mooseError("Internal error, Material object that does not inherit form BlockRestricted");
664  const std::set<SubdomainID> & blocks =
665  blk->blockRestricted() ? blk->blockIDs() : blk->meshBlockIDs();
666 
667  if (std::includes(blocks.begin(), blocks.end(), _subdomain_ids.begin(), _subdomain_ids.end()))
668  {
669  material_eigenstrain_map[name].insert(blocks.begin(), blocks.end());
670  eigenstrain_set.insert(name);
671  }
672  }
673 
674  // Account for reduced eigenstrains and CompositeEigenstrains
675  if (mat_name == "ComputeReducedOrderEigenstrain")
676  {
677  auto input_eigenstrain_names =
678  mat_params.get<std::vector<MaterialPropertyName>>("input_eigenstrain_names");
679  remove_reduced_set.insert(input_eigenstrain_names.begin(), input_eigenstrain_names.end());
680  }
681  // Account for CompositeEigenstrains
682  if (mat_name == "CompositeEigenstrain")
683  {
684  auto remove_list = mat_params.get<std::vector<MaterialPropertyName>>("tensors");
685  for (auto i : remove_list)
686  remove_reduced_set.insert(i);
687  }
688 
689  // Account for MaterialADConverter , add or remove later
690  if (mat_name == "RankTwoTensorMaterialADConverter")
691  {
692  std::vector<MaterialPropertyName> remove_list;
693  std::vector<MaterialPropertyName> add_list;
694 
695  if (mat_params.isParamValid("ad_props_out") && mat_params.isParamValid("reg_props_in") &&
696  _use_ad)
697  {
698  remove_list = mat_params.get<std::vector<MaterialPropertyName>>("reg_props_in");
699  add_list = mat_params.get<std::vector<MaterialPropertyName>>("ad_props_out");
700  }
701  if (mat_params.isParamValid("ad_props_in") && mat_params.isParamValid("reg_props_out") &&
702  !_use_ad)
703  {
704  remove_list = mat_params.get<std::vector<MaterialPropertyName>>("ad_props_in");
705  add_list = mat_params.get<std::vector<MaterialPropertyName>>("reg_props_out");
706  }
707 
708  // These vectors are the same size as checked in MaterialADConverter
709  for (unsigned int index = 0; index < remove_list.size(); index++)
710  remove_add_map.emplace(remove_list[index], add_list[index]);
711  }
712  }
713  // All the materials have been accounted for, now remove or add parts
714 
715  // Remove names which aren't eigenstrains (converter properties)
716  for (auto remove_add_index : remove_add_map)
717  {
718  const bool is_in = eigenstrain_set.find(remove_add_index.first) != eigenstrain_set.end();
719  if (is_in)
720  {
721  eigenstrain_set.erase(remove_add_index.first);
722  eigenstrain_set.insert(remove_add_index.second);
723  }
724  }
725  for (auto index : remove_reduced_set)
726  eigenstrain_set.erase(index);
727 
728  // Compare the blockIDs set of eigenstrain names with the vector of _eigenstrain_names for the
729  // current subdomainID
730  std::set_union(eigenstrain_set.begin(),
731  eigenstrain_set.end(),
732  _eigenstrain_names.begin(),
733  _eigenstrain_names.end(),
734  std::inserter(verified_eigenstrain_names, verified_eigenstrain_names.begin()));
735 
736  // Ensure the eigenstrain names previously passed include any missing names
737  _eigenstrain_names.resize(verified_eigenstrain_names.size());
738  std::copy(verified_eigenstrain_names.begin(),
739  verified_eigenstrain_names.end(),
740  _eigenstrain_names.begin());
741 
742  Moose::out << COLOR_CYAN << "*** Automatic Eigenstrain Names ***"
743  << "\n"
744  << _name << ": " << Moose::stringify(_eigenstrain_names) << "\n"
745  << COLOR_DEFAULT << std::flush;
746 }
747 
748 void
750 {
751  // Ensure material output order and family vectors are same size as generate output
752 
753  // check number of supplied orders and families
755  paramError("material_output_order",
756  "The number of orders assigned to material outputs must be: 0 to be assigned "
757  "CONSTANT; 1 to assign all outputs the same value, or the same size as the number "
758  "of generate outputs listed.");
759 
760  if (_material_output_family.size() > 1 &&
762  paramError("material_output_family",
763  "The number of families assigned to material outputs must be: 0 to be assigned "
764  "MONOMIAL; 1 to assign all outputs the same value, or the same size as the number "
765  "of generate outputs listed.");
766 
767  // if no value was provided, chose the default CONSTANT
768  if (_material_output_order.size() == 0)
770 
771  // For only one order, make all orders the same magnitude
772  if (_material_output_order.size() == 1)
774  std::vector<std::string>(_generate_output.size(), _material_output_order[0]);
775 
776  if (_verbose)
777  Moose::out << COLOR_CYAN << "*** Automatic applied material output orders ***"
778  << "\n"
779  << _name << ": " << Moose::stringify(_material_output_order) << "\n"
780  << COLOR_DEFAULT << std::flush;
781 
782  // if no value was provided, chose the default MONOMIAL
783  if (_material_output_family.size() == 0)
785 
786  // For only one family, make all families that value
787  if (_material_output_family.size() == 1)
789  std::vector<std::string>(_generate_output.size(), _material_output_family[0]);
790 
791  if (_verbose)
792  Moose::out << COLOR_CYAN << "*** Automatic applied material output families ***"
793  << "\n"
794  << _name << ": " << Moose::stringify(_material_output_family) << "\n"
795  << COLOR_DEFAULT << std::flush;
796 }
797 
798 void
800 {
801  std::string ad_prepend = _use_ad ? "AD" : "";
802 
803  if (_current_task == "add_material")
804  {
805  // Add output Materials
806  for (auto out : _generate_output)
807  {
809 
810  // RankTwoCartesianComponent
811  if (
812  [&]()
813  {
814  for (const auto & r2q : _rank_two_cartesian_component_table)
815  for (unsigned int a = 0; a < 3; ++a)
816  for (unsigned int b = 0; b < 3; ++b)
817  if (r2q.first + '_' + _component_table[a] + _component_table[b] == out)
818  {
819  auto type = ad_prepend + "RankTwoCartesianComponent";
820  params = _factory.getValidParams(type);
821  params.set<MaterialPropertyName>("rank_two_tensor") = _base_name + r2q.second;
822  params.set<unsigned int>("index_i") = a;
823  params.set<unsigned int>("index_j") = b;
824 
825  params.applyParameters(parameters());
826  params.set<MaterialPropertyName>("property_name") = _base_name + out;
827  _problem->addMaterial(type, _base_name + out + '_' + name(), params);
828  return true;
829  }
830  return false;
831  }())
832  continue;
833 
834  // RankTwoDirectionalComponent
835  if (setupOutput(out,
837  [&](std::string prop_name, std::string invariant)
838  {
839  auto type = ad_prepend + "RankTwoDirectionalComponent";
840  params = _factory.getValidParams(type);
841  params.set<MaterialPropertyName>("rank_two_tensor") =
842  _base_name + prop_name;
843  params.set<MooseEnum>("invariant") = invariant;
844  params.applyParameters(parameters());
845  params.set<MaterialPropertyName>("property_name") = _base_name + out;
846  _problem->addMaterial(type, _base_name + out + '_' + name(), params);
847  }))
848  continue;
849 
850  // RankTwoInvariant
851  if (setupOutput(out,
853  [&](std::string prop_name, std::string invariant)
854  {
855  auto type = ad_prepend + "RankTwoInvariant";
856  params = _factory.getValidParams(type);
857  params.set<MaterialPropertyName>("rank_two_tensor") =
858  _base_name + prop_name;
859  params.set<MooseEnum>("invariant") = invariant;
860  params.applyParameters(parameters());
861  params.set<MaterialPropertyName>("property_name") = _base_name + out;
862  _problem->addMaterial(type, _base_name + out + '_' + name(), params);
863  }))
864  continue;
865 
866  // RankTwoCylindricalComponent
867  if (setupOutput(
868  out,
870  [&](std::string prop_name, std::string component)
871  {
873  mooseError(
874  "Cannot use cylindrical component output in a spherical coordinate system.");
875  auto type = ad_prepend + "RankTwoCylindricalComponent";
876  params = _factory.getValidParams(type);
877  params.set<MaterialPropertyName>("rank_two_tensor") = _base_name + prop_name;
878  params.set<MooseEnum>("cylindrical_component") = component;
879  params.applyParameters(parameters());
880  params.set<MaterialPropertyName>("property_name") = _base_name + out;
881  _problem->addMaterial(type, _base_name + out + '_' + name(), params);
882  }))
883  continue;
884 
885  // RankTwoSphericalComponent
886  if (setupOutput(out,
888  [&](std::string prop_name, std::string component)
889  {
890  auto type = ad_prepend + "RankTwoSphericalComponent";
891  params = _factory.getValidParams(type);
892  params.set<MaterialPropertyName>("rank_two_tensor") =
893  _base_name + prop_name;
894  params.set<MooseEnum>("spherical_component") = component;
895  params.applyParameters(parameters());
896  params.set<MaterialPropertyName>("property_name") = _base_name + out;
897  _problem->addMaterial(type, _base_name + out + '_' + name(), params);
898  }))
899  continue;
900 
901  paramError("generate_output", "Unable to add output Material for '", out, "'");
902  }
903  }
904 }
905 
906 void
908 {
909  // Gather info about all other solid mechanics physics when we add variables
910  if (_current_task == "validate_coordinate_systems" && getParam<bool>("add_variables"))
911  {
913  for (const auto & action : actions)
914  {
915  const auto size_before = _subdomain_id_union.size();
916  const auto added_size = action->_subdomain_ids.size();
917  _subdomain_id_union.insert(action->_subdomain_ids.begin(), action->_subdomain_ids.end());
918  const auto size_after = _subdomain_id_union.size();
919 
920  if (size_after != size_before + added_size)
921  mooseError("The block restrictions in the SolidMechanics/QuasiStatic actions must be "
922  "non-overlapping.");
923 
924  if (added_size == 0 && actions.size() > 1)
925  mooseError(
926  "No SolidMechanics/QuasiStatic action can be block unrestricted if more than one "
927  "SolidMechanics/QuasiStatic action is specified.");
928  }
929  }
930 }
931 
932 void
934 {
935  std::string type;
937  type = "ComputeLagrangianStrain";
938  else if (_coord_system == Moose::COORD_RZ)
939  type = "ComputeLagrangianStrainAxisymmetricCylindrical";
941  type = "ComputeLagrangianStrainCentrosymmetricSpherical";
942  else
943  mooseError("Unsupported coordinate system");
944 
945  auto params = _factory.getValidParams(type);
946 
947  if (isParamValid("strain_base_name"))
948  params.set<std::string>("base_name") = getParam<std::string>("strain_base_name");
949 
950  params.set<std::vector<VariableName>>("displacements") = _coupled_displacements;
951  params.set<std::vector<MaterialPropertyName>>("eigenstrain_names") = _eigenstrain_names;
952  params.set<bool>("large_kinematics") = _lk_large_kinematics;
953  params.set<std::vector<SubdomainName>>("block") = _subdomain_names;
954 
955  // Error if volumetric locking correction is on for higher-order elements
956  if (_problem->mesh().hasSecondOrderElements() && _lk_locking)
957  mooseError("Volumetric locking correction should not be used for "
958  "higher-order elements.");
959 
960  params.set<bool>("stabilize_strain") = _lk_locking;
961 
962  if (_lk_homogenization)
963  {
964  params.set<std::vector<MaterialPropertyName>>("homogenization_gradient_names") = {
966  }
967 
968  _problem->addMaterial(type, name() + "_strain", params);
969 
970  // Add the homogenization strain calculator
971  if (_lk_homogenization)
972  {
973  std::string type = "ComputeHomogenizedLagrangianStrain";
974  auto params = _factory.getValidParams(type);
975 
976  params.set<MaterialPropertyName>("homogenization_gradient_name") = _homogenization_strain_name;
977  params.set<std::vector<VariableName>>("macro_gradient") = {_hname};
978  params.set<MultiMooseEnum>("constraint_types") = _constraint_types;
979  params.set<std::vector<FunctionName>>("targets") = _targets;
980 
981  _problem->addMaterial(type, name() + "_compute_" + _homogenization_strain_name, params);
982  }
983 }
984 
985 void
987 {
988  std::string ad_prepend = _use_ad ? "AD" : "";
989 
990  std::string type;
991 
992  // no plane strain
994  {
995  std::map<std::pair<Moose::CoordinateSystemType, StrainAndIncrement>, std::string> type_map = {
996  {{Moose::COORD_XYZ, StrainAndIncrement::SmallTotal}, "ComputeSmallStrain"},
997  {{Moose::COORD_XYZ, StrainAndIncrement::SmallIncremental}, "ComputeIncrementalStrain"},
998  {{Moose::COORD_XYZ, StrainAndIncrement::FiniteIncremental}, "ComputeFiniteStrain"},
999  {{Moose::COORD_RZ, StrainAndIncrement::SmallTotal}, "ComputeAxisymmetricRZSmallStrain"},
1001  "ComputeAxisymmetricRZIncrementalStrain"},
1003  "ComputeAxisymmetricRZFiniteStrain"},
1004  {{Moose::COORD_RSPHERICAL, StrainAndIncrement::SmallTotal}, "ComputeRSphericalSmallStrain"},
1006  "ComputeRSphericalIncrementalStrain"},
1008  "ComputeRSphericalFiniteStrain"}};
1009 
1010  auto type_it = type_map.find(std::make_pair(_coord_system, _strain_and_increment));
1011  if (type_it != type_map.end())
1012  type = type_it->second;
1013  else
1014  mooseError("Unsupported strain formulation");
1015  }
1019  {
1022  paramError("use_automatic_differentiation",
1023  "AD not setup for use with PlaneStrain or GeneralizedPlaneStrain");
1024 
1025  std::map<std::pair<Moose::CoordinateSystemType, StrainAndIncrement>, std::string> type_map = {
1026  {{Moose::COORD_XYZ, StrainAndIncrement::SmallTotal}, "ComputePlaneSmallStrain"},
1027  {{Moose::COORD_XYZ, StrainAndIncrement::SmallIncremental}, "ComputePlaneIncrementalStrain"},
1028  {{Moose::COORD_XYZ, StrainAndIncrement::FiniteIncremental}, "ComputePlaneFiniteStrain"},
1029  {{Moose::COORD_RZ, StrainAndIncrement::SmallTotal}, "ComputeAxisymmetric1DSmallStrain"},
1031  "ComputeAxisymmetric1DIncrementalStrain"},
1033  "ComputeAxisymmetric1DFiniteStrain"}};
1034 
1035  // choose kernel type based on coordinate system
1036  auto type_it = type_map.find(std::make_pair(_coord_system, _strain_and_increment));
1037  if (type_it != type_map.end())
1038  type = type_it->second;
1039  else
1040  mooseError("Unsupported coordinate system for plane strain.");
1041  }
1042  else
1043  mooseError("Unsupported planar formulation");
1044 
1045  // set material parameters
1046  auto params = _factory.getValidParams(ad_prepend + type);
1047  params.applyParameters(
1048  parameters(),
1049  {"displacements", "use_displaced_mesh", "out_of_plane_strain", "scalar_out_of_plane_strain"});
1050 
1051  if (isParamValid("strain_base_name"))
1052  params.set<std::string>("base_name") = getParam<std::string>("strain_base_name");
1053 
1054  params.set<std::vector<VariableName>>("displacements") = _coupled_displacements;
1055  params.set<bool>("use_displaced_mesh") = false;
1056 
1057  if (isParamValid("scalar_out_of_plane_strain"))
1058  params.set<std::vector<VariableName>>("scalar_out_of_plane_strain") = {
1059  getParam<VariableName>("scalar_out_of_plane_strain")};
1060 
1061  if (isParamValid("out_of_plane_strain"))
1062  params.set<std::vector<VariableName>>("out_of_plane_strain") = {
1063  getParam<VariableName>("out_of_plane_strain")};
1064 
1065  params.set<std::vector<MaterialPropertyName>>("eigenstrain_names") = _eigenstrain_names;
1066 
1067  _problem->addMaterial(ad_prepend + type, name() + "_strain", params);
1068 }
1069 
1070 std::string
1072 {
1073  if (_lagrangian_kernels)
1074  {
1075  std::string type;
1077  {
1078  if (_lk_homogenization)
1079  type = "HomogenizedTotalLagrangianStressDivergence";
1081  type = "TotalLagrangianStressDivergence";
1083  type = "UpdatedLagrangianStressDivergence";
1084  else
1085  mooseError("Unknown formulation type");
1086  }
1087  else if (_coord_system == Moose::COORD_RZ)
1088  {
1089  if (_lk_homogenization)
1090  mooseError("The Lagrangian mechanics kernels do not yet support homogenization in "
1091  "coordinate systems other than Cartesian.");
1093  type = "TotalLagrangianStressDivergenceAxisymmetricCylindrical";
1095  mooseError("The Lagrangian mechanics kernels do not yet support the updated Lagrangian "
1096  "formulation in RZ coordinates.");
1097  }
1099  {
1100  if (_lk_homogenization)
1101  mooseError("The Lagrangian mechanics kernels do not yet support homogenization in "
1102  "coordinate systems other than Cartesian.");
1104  type = "TotalLagrangianStressDivergenceCentrosymmetricSpherical";
1106  mooseError("The Lagrangian mechanics kernels do not yet support the updated Lagrangian "
1107  "formulation in RZ coordinates.");
1108  }
1109  else
1110  mooseError("Unsupported coordinate system");
1111  return type;
1112  }
1113  else
1114  {
1115  std::map<Moose::CoordinateSystemType, std::string> type_map = {
1116  {Moose::COORD_XYZ, "StressDivergenceTensors"},
1117  {Moose::COORD_RZ, "StressDivergenceRZTensors"},
1118  {Moose::COORD_RSPHERICAL, "StressDivergenceRSphericalTensors"}};
1119 
1120  // choose kernel type based on coordinate system
1121  auto type_it = type_map.find(_coord_system);
1122  if (type_it != type_map.end())
1123  return type_it->second;
1124  else
1125  mooseError("Unsupported coordinate system");
1126  }
1127 }
1128 
1131 {
1133  params.applyParameters(
1134  parameters(),
1135  {"displacements", "use_displaced_mesh", "save_in", "diag_save_in", "out_of_plane_strain"});
1136 
1137  params.set<std::vector<VariableName>>("displacements") = _coupled_displacements;
1138 
1139  if (_lagrangian_kernels)
1140  {
1141  params.set<bool>("use_displaced_mesh") =
1143  params.set<bool>("large_kinematics") = _lk_large_kinematics;
1144  params.set<bool>("stabilize_strain") = _lk_locking;
1145  }
1146  else
1147  params.set<bool>("use_displaced_mesh") = _use_displaced_mesh;
1148 
1149  return params;
1150 }
const bool _lk_large_kinematics
Simplified flag for small/large deformations, Lagrangian kernel system.
void mooseInfo(Args &&... args) const
std::set< SubdomainID > _subdomain_id_union
set generated from the combined block restrictions of all SolidMechanics/Master action blocks ...
enum QuasiStaticSolidMechanicsPhysics::StrainAndIncrement _strain_and_increment
std::vector< AuxVariableName > _save_in
residual debugging
virtual InputParameters getKernelParameters(std::string type)
const std::string & _name
bool setupOutput(std::string out, T table, T2 setup)
Helper function to decode generate_outputs options using a "table" of scalar output quantities and a ...
ActionWarehouse & _awh
void paramError(const std::string &param, Args... args) const
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
std::vector< AuxVariableName > _diag_save_in
static const std::map< std::string, std::pair< std::string, std::vector< std::string > > > _rank_two_invariant_table
char ** blocks
std::set< SubdomainID > _subdomain_ids
set generated from the passed in vector of subdomain names
InputParameters getValidParams(const std::string &name)
void setAdditionalValue(const std::string &names)
Factory & _factory
unsigned int _ndisp
Number of displacement variables.
static const std::string component
Definition: NS.h:153
const InputParameters & parameters() const
ActionFactory & _action_factory
LKFormulation
New kernel system kinematics types.
static const std::map< std::string, std::pair< std::string, std::vector< std::string > > > _rank_two_cylindrical_component_table
T & set(const std::string &name, bool quiet_mode=false)
InputParameters getValidParams(const std::string &name) const
const bool _auto_eigenstrain
automatically gather names of eigenstrain tensors provided by simulation objects
registerMooseAction("SolidMechanicsApp", QuasiStaticSolidMechanicsPhysics, "meta_action")
unsigned int size() const
void applyParameters(const InputParameters &common, const std::vector< std::string > &exclude={}, const bool allow_private=false)
void mooseDocumentedError(const std::string &repo_name, const unsigned int issue_num, Args &&... args) const
void addActionBlock(std::shared_ptr< Action > blk)
Store common tensor mechanics parameters.
COORD_RSPHERICAL
const ExecFlagType EXEC_TIMESTEP_END
std::vector< VariableName > _displacements
displacement variables
std::vector< std::string > _generate_output
output materials to generate scalar stress/strain tensor quantities
static const std::vector< char > _component_table
table data for output generation
std::shared_ptr< Action > create(const std::string &action, const std::string &action_name, InputParameters &parameters)
InputParameters emptyInputParameters()
enum QuasiStaticSolidMechanicsPhysics::PlanarFormulation _planar_formulation
QuasiStaticSolidMechanicsPhysics(const InputParameters &params)
const SubdomainID INVALID_BLOCK_ID
std::vector< SubdomainName > _subdomain_names
if this vector is not empty the variables, kernels and materials are restricted to these subdomains ...
bool _use_displaced_mesh
use displaced mesh (true unless _strain is SMALL)
const std::string & name() const
enum QuasiStaticSolidMechanicsPhysics::Strain _strain
std::vector< VariableName > _coupled_displacements
Coupled displacement variables.
bool _cylindrical_axis_point1_valid
booleans used to determine if cylindrical axis points are passed
const std::string & type() const
bool _lk_locking
Simplified volumetric locking correction flag for new kernels.
const std::string & _current_task
static std::string variableType(const libMesh::FEType &fe_type, const bool is_fv=false, const bool is_array=false)
std::string stringify(const T &t)
static const std::map< std::string, std::pair< std::string, std::vector< std::string > > > _rank_two_directional_component_table
static std::map< std::string, std::string > _rank_two_cartesian_component_table
std::vector< MaterialPropertyName > _eigenstrain_names
bool isParamSetByUser(const std::string &name) const
std::shared_ptr< MooseMesh > & _mesh
static libMesh::FEType feType(const InputParameters &params)
unsigned int get(unsigned int i) const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
ConstraintType
Constraint type: stress/PK stress or strain/deformation gradient.
OStreamProxy out
const bool _lagrangian_kernels
New or old kernel system.
void mooseError(Args &&... args) const
void addClassDescription(const std::string &doc_string)
std::shared_ptr< FEProblemBase > & _problem
Point _spherical_center_point
center point for spherical stress/strain quantities
const MultiMooseEnum constraintType("strain stress none")
bool isParamValid(const std::string &name) const
bool _spherical_center_point_valid
booleans used to determine if spherical center point is passed
const std::string _base_name
base name for the current master action block
std::vector< const T *> getActions()
auto index_range(const T &sizable)
PlanarFormulation
use an out of plane stress/strain formulation
static const std::map< std::string, std::pair< std::string, std::vector< std::string > > > _rank_two_spherical_component_table
Point _cylindrical_axis_point1
points used to determine axis of rotation for cylindrical stress/strain quantities ...
void addParamNamesToGroup(const std::string &space_delim_names, const std::string group_name)
bool _lk_homogenization
Flag indicating if the homogenization system is present for new kernels.