LCOV - code coverage report
Current view: top level - src/physics - QuasiStaticSolidMechanicsPhysics.C (source / functions) Hit Total Coverage
Test: idaholab/moose solid_mechanics: #33416 (b10b36) with base 9fbd27 Lines: 473 555 85.2 %
Date: 2026-07-23 16:20:48 Functions: 17 18 94.4 %
Legend: Lines: hit not hit

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

Generated by: LCOV version 1.14