Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://www.mooseframework.org
       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 "ComputeLagrangianStrain.h"
      11             : 
      12             : template <class G>
      13             : InputParameters
      14        2589 : ComputeLagrangianStrainBase<G>::baseParams()
      15             : {
      16        2589 :   InputParameters params = Material::validParams();
      17             : 
      18        5178 :   params.addRequiredCoupledVar("displacements", "Displacement variables");
      19        5178 :   params.addParam<bool>(
      20        5178 :       "large_kinematics", false, "Use large displacement kinematics in the kernel.");
      21        5178 :   params.addParam<bool>("stabilize_strain", false, "Average the volumetric strains");
      22        5178 :   params.addParam<std::vector<MaterialPropertyName>>(
      23             :       "eigenstrain_names", {}, "List of eigenstrains to account for");
      24        5178 :   params.addParam<std::vector<MaterialPropertyName>>(
      25             :       "homogenization_gradient_names",
      26             :       {},
      27             :       "List of homogenization gradients to add to the displacement gradient");
      28             : 
      29        5178 :   params.addParam<std::string>("base_name", "Material property base name");
      30             : 
      31             :   // We rely on this *not* having use_displaced mesh on
      32        2589 :   params.suppressParameter<bool>("use_displaced_mesh");
      33             : 
      34        2589 :   return params;
      35           0 : }
      36             : 
      37             : template <class G>
      38        1941 : ComputeLagrangianStrainBase<G>::ComputeLagrangianStrainBase(const InputParameters & parameters)
      39             :   : Material(parameters),
      40        1941 :     _ndisp(coupledComponents("displacements")),
      41        1941 :     _disp(coupledValues("displacements")),
      42        1941 :     _grad_disp(coupledGradients("displacements")),
      43        3918 :     _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
      44        3882 :     _large_kinematics(getParam<bool>("large_kinematics")),
      45        3882 :     _stabilize_strain(getParam<bool>("stabilize_strain")),
      46        3882 :     _eigenstrain_names(getParam<std::vector<MaterialPropertyName>>("eigenstrain_names")),
      47        1941 :     _eigenstrains(_eigenstrain_names.size()),
      48        1941 :     _eigenstrains_old(_eigenstrain_names.size()),
      49        1941 :     _total_strain(declareProperty<RankTwoTensor>(_base_name + "total_strain")),
      50        3882 :     _total_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "total_strain")),
      51        1941 :     _mechanical_strain(declareProperty<RankTwoTensor>(_base_name + "mechanical_strain")),
      52        3882 :     _mechanical_strain_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "mechanical_strain")),
      53        1941 :     _strain_increment(declareProperty<RankTwoTensor>(_base_name + "strain_increment")),
      54        1941 :     _vorticity_increment(declareProperty<RankTwoTensor>(_base_name + "vorticity_increment")),
      55        1941 :     _F_ust(declareProperty<RankTwoTensor>(_base_name + "unstabilized_deformation_gradient")),
      56        1941 :     _F_avg(declareProperty<RankTwoTensor>(_base_name + "average_deformation_gradient")),
      57        1941 :     _F(declareProperty<RankTwoTensor>(_base_name + "deformation_gradient")),
      58        3882 :     _F_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "deformation_gradient")),
      59        1941 :     _F_inv(declareProperty<RankTwoTensor>(_base_name + "inverse_deformation_gradient")),
      60        1941 :     _f_inv(declareProperty<RankTwoTensor>(_base_name + "inverse_incremental_deformation_gradient")),
      61        3882 :     _homogenization_gradient_names(
      62             :         getParam<std::vector<MaterialPropertyName>>("homogenization_gradient_names")),
      63        1941 :     _homogenization_contributions(_homogenization_gradient_names.size()),
      64        3882 :     _rotation_increment(declareProperty<RankTwoTensor>(_base_name + "rotation_increment"))
      65             : {
      66             :   // Setup eigenstrains
      67        2115 :   for (auto i : make_range(_eigenstrain_names.size()))
      68             :   {
      69         174 :     _eigenstrains[i] = &getMaterialProperty<RankTwoTensor>(_eigenstrain_names[i]);
      70         174 :     _eigenstrains_old[i] = &getMaterialPropertyOld<RankTwoTensor>(_eigenstrain_names[i]);
      71             :   }
      72             : 
      73             :   // In the future maybe there is a reason to have more than one, but for now
      74        1941 :   if (_homogenization_gradient_names.size() > 1)
      75           0 :     mooseError("ComputeLagrangianStrainBase cannot accommodate more than one "
      76             :                "homogenization gradient");
      77             : 
      78             :   // Setup homogenization contributions
      79        2208 :   for (unsigned int i = 0; i < _homogenization_gradient_names.size(); i++)
      80         267 :     _homogenization_contributions[i] =
      81         267 :         &getMaterialProperty<RankTwoTensor>(_homogenization_gradient_names[i]);
      82        1941 : }
      83             : 
      84             : template <class G>
      85             : void
      86      465124 : ComputeLagrangianStrainBase<G>::initQpStatefulProperties()
      87             : {
      88      465124 :   _total_strain[_qp].zero();
      89      465124 :   _mechanical_strain[_qp].zero();
      90      465124 :   _F[_qp].setToIdentity();
      91      465124 : }
      92             : 
      93             : template <class G>
      94             : void
      95     6012171 : ComputeLagrangianStrainBase<G>::computeProperties()
      96             : {
      97             :   // Average the volumetric terms, if required
      98     6012171 :   computeDeformationGradient();
      99             : 
     100    49070671 :   for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
     101    43058500 :     computeQpProperties();
     102     6012171 : }
     103             : 
     104             : template <class G>
     105             : void
     106    43058500 : ComputeLagrangianStrainBase<G>::computeQpProperties()
     107             : {
     108             :   // Add in the macroscale gradient contribution
     109    61858172 :   for (auto contribution : _homogenization_contributions)
     110    18799672 :     _F[_qp] += (*contribution)[_qp];
     111             : 
     112             :   // If the kernel is large deformation then we need the "actual"
     113             :   // kinematic quantities
     114    43058500 :   RankTwoTensor dL;
     115    43058500 :   if (_large_kinematics)
     116             :   {
     117    20241820 :     _F_inv[_qp] = _F[_qp].inverse();
     118    20241820 :     _f_inv[_qp] = _F_old[_qp] * _F_inv[_qp];
     119    20241820 :     dL = RankTwoTensor::Identity() - _f_inv[_qp];
     120             :   }
     121             :   // For small deformations we just provide the identity
     122             :   else
     123             :   {
     124    22816680 :     _F_inv[_qp] = RankTwoTensor::Identity();
     125    22816680 :     _f_inv[_qp] = RankTwoTensor::Identity();
     126    22816680 :     dL = _F[_qp] - _F_old[_qp];
     127             :   }
     128             : 
     129    43058500 :   computeQpIncrementalStrains(dL);
     130    43058500 : }
     131             : 
     132             : template <class G>
     133             : void
     134    43058500 : ComputeLagrangianStrainBase<G>::computeQpIncrementalStrains(const RankTwoTensor & dL)
     135             : {
     136             :   // Get the deformation increments
     137    43058500 :   _strain_increment[_qp] = (dL + dL.transpose()) / 2.0;
     138    43058500 :   _vorticity_increment[_qp] = (dL - dL.transpose()) / 2.0;
     139             : 
     140             :   // Increment the total strain
     141    43058500 :   _total_strain[_qp] = _total_strain_old[_qp] + _strain_increment[_qp];
     142             : 
     143             :   // Get rid of the eigenstrains
     144             :   // Note we currently do not alter the deformation gradient -- this will be
     145             :   // needed in the future for a "complete" system
     146    43058500 :   subtractQpEigenstrainIncrement(_strain_increment[_qp]);
     147             : 
     148             :   // Increment the mechanical strain
     149    43058500 :   _mechanical_strain[_qp] = _mechanical_strain_old[_qp] + _strain_increment[_qp];
     150             : 
     151             :   // Faked rotation increment for ComputeStressBase materials
     152    43058500 :   _rotation_increment[_qp] = RankTwoTensor::Identity();
     153    43058500 : }
     154             : 
     155             : template <class G>
     156             : void
     157    43058500 : ComputeLagrangianStrainBase<G>::subtractQpEigenstrainIncrement(RankTwoTensor & strain)
     158             : {
     159    43471396 :   for (auto i : make_range(_eigenstrain_names.size()))
     160      412896 :     strain -= (*_eigenstrains[i])[_qp] - (*_eigenstrains_old[i])[_qp];
     161    43058500 : }
     162             : 
     163             : template <class G>
     164             : void
     165    43058500 : ComputeLagrangianStrainBase<G>::computeQpUnstabilizedDeformationGradient()
     166             : {
     167    43058500 :   _F_ust[_qp].setToIdentity();
     168   171099640 :   for (auto component : make_range(_ndisp))
     169   128041140 :     G::addGradOp(_F_ust[_qp],
     170             :                  component,
     171             :                  (*_grad_disp[component])[_qp],
     172   128041140 :                  (*_disp[component])[_qp],
     173      336648 :                  _q_point[_qp]);
     174    43058500 : }
     175             : 
     176             : template <class G>
     177             : void
     178     6012171 : ComputeLagrangianStrainBase<G>::computeDeformationGradient()
     179             : {
     180             :   // First calculate the unstabilized deformation gradient at each qp
     181    49070671 :   for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
     182             :   {
     183    43058500 :     computeQpUnstabilizedDeformationGradient();
     184    43058500 :     _F[_qp] = _F_ust[_qp];
     185             :   }
     186             : 
     187             :   // If stabilization is on do the volumetric correction
     188     6012171 :   if (_stabilize_strain)
     189             :   {
     190     1850354 :     const auto F_avg = StabilizationUtils::elementAverage(
     191    15564754 :         [this](unsigned int qp) { return _F_ust[qp]; }, _JxW, _coord);
     192             :     // All quadrature points have the same F_avg
     193     1850354 :     _F_avg.set().setAllValues(F_avg);
     194             :     // Make the appropriate modification, depending on small or large
     195             :     // deformations
     196    15564754 :     for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
     197             :     {
     198    13714400 :       if (_large_kinematics)
     199     4766700 :         _F[_qp] *= std::pow(F_avg.det() / _F[_qp].det(), 1.0 / 3.0);
     200             :       else
     201     8947700 :         _F[_qp] += (F_avg.trace() - _F[_qp].trace()) * RankTwoTensor::Identity() / 3.0;
     202             :     }
     203             :   }
     204     6012171 : }
     205             : 
     206             : template class ComputeLagrangianStrainBase<GradientOperatorCartesian>;
     207             : template class ComputeLagrangianStrainBase<GradientOperatorAxisymmetricCylindrical>;
     208             : template class ComputeLagrangianStrainBase<GradientOperatorCentrosymmetricSpherical>;