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 "WeakPlaneStress.h"
      11             : 
      12             : #include "Material.h"
      13             : #include "MooseMesh.h"
      14             : #include "MooseVariable.h"
      15             : #include "RankTwoTensor.h"
      16             : #include "RankFourTensor.h"
      17             : 
      18             : registerMooseObject("SolidMechanicsApp", WeakPlaneStress);
      19             : 
      20             : InputParameters
      21          92 : WeakPlaneStress::validParams()
      22             : {
      23          92 :   InputParameters params = Kernel::validParams();
      24          92 :   params.addClassDescription("Plane stress kernel to provide out-of-plane strain contribution.");
      25         184 :   params.addCoupledVar("displacements",
      26             :                        "The string of displacements suitable for the problem statement");
      27         184 :   params.addCoupledVar("temperature",
      28             :                        "The name of the temperature variable used in the "
      29             :                        "ComputeThermalExpansionEigenstrain.  (Not required for "
      30             :                        "simulations without temperature coupling.)");
      31         184 :   params.addParam<std::vector<MaterialPropertyName>>(
      32             :       "eigenstrain_names",
      33             :       {},
      34             :       "List of eigenstrains used in the strain calculation. Used for computing their derivaties "
      35             :       "for off-diagonal Jacobian terms.");
      36         184 :   params.addParam<std::string>("base_name", "Material property base name");
      37             : 
      38         184 :   MooseEnum direction("x y z", "z");
      39         184 :   params.addParam<MooseEnum>("out_of_plane_strain_direction",
      40             :                              direction,
      41             :                              "The direction of the out-of-plane strain variable");
      42             : 
      43          92 :   params.set<bool>("use_displaced_mesh") = false;
      44             : 
      45          92 :   return params;
      46          92 : }
      47             : 
      48          46 : WeakPlaneStress::WeakPlaneStress(const InputParameters & parameters)
      49             :   : DerivativeMaterialInterface<Kernel>(parameters),
      50          46 :     _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
      51          92 :     _stress(getMaterialProperty<RankTwoTensor>(_base_name + "stress")),
      52          92 :     _Jacobian_mult(getMaterialProperty<RankFourTensor>(_base_name + "Jacobian_mult")),
      53          92 :     _direction(getParam<MooseEnum>("out_of_plane_strain_direction")),
      54          46 :     _disp_coupled(isCoupled("displacements")),
      55          46 :     _ndisp(_disp_coupled ? coupledComponents("displacements") : 0),
      56          46 :     _disp_var(_ndisp),
      57          46 :     _temp_coupled(isCoupled("temperature")),
      58         110 :     _temp_var(_temp_coupled ? coupled("temperature") : 0)
      59             : {
      60          46 :   if (_disp_coupled)
      61         138 :     for (unsigned int i = 0; i < _ndisp; ++i)
      62          92 :       _disp_var[i] = coupled("displacements", i);
      63             : 
      64          46 :   if (_temp_coupled)
      65             :   {
      66          54 :     for (auto eigenstrain_name : getParam<std::vector<MaterialPropertyName>>("eigenstrain_names"))
      67          54 :       _deigenstrain_dT.push_back(&getMaterialPropertyDerivative<RankTwoTensor>(
      68          36 :           eigenstrain_name, coupledName("temperature", 0)));
      69             :   }
      70             : 
      71             :   // Checking for consistency between mesh size and length of the provided displacements vector
      72          46 :   if (_disp_coupled && _ndisp != _mesh.dimension())
      73           0 :     mooseError("The number of displacement variables supplied must match the mesh dimension.");
      74          46 : }
      75             : 
      76             : Real
      77      144128 : WeakPlaneStress::computeQpResidual()
      78             : {
      79      144128 :   return _stress[_qp](_direction, _direction) * _test[_i][_qp];
      80             : }
      81             : 
      82             : Real
      83       62976 : WeakPlaneStress::computeQpJacobian()
      84             : {
      85       62976 :   return _Jacobian_mult[_qp](_direction, _direction, _direction, _direction) * _test[_i][_qp] *
      86       62976 :          _phi[_j][_qp];
      87             : }
      88             : 
      89             : Real
      90        2048 : WeakPlaneStress::computeQpOffDiagJacobian(unsigned int jvar)
      91             : {
      92             :   Real val = 0.0;
      93             : 
      94             :   // off-diagonal Jacobian with respect to a coupled displacement component
      95        2048 :   if (_disp_coupled)
      96             :   {
      97        6144 :     for (unsigned int coupled_direction = 0; coupled_direction < _ndisp; ++coupled_direction)
      98             :     {
      99        4096 :       if (jvar == _disp_var[coupled_direction])
     100             :       {
     101             :         unsigned int coupled_direction_index = 0;
     102        2048 :         switch (_direction)
     103             :         {
     104           0 :           case 0: // x
     105             :           {
     106           0 :             if (coupled_direction == 0)
     107             :               coupled_direction_index = 1;
     108             :             else
     109             :               coupled_direction_index = 2;
     110             :             break;
     111             :           }
     112           0 :           case 1: // y
     113             :           {
     114           0 :             if (coupled_direction == 0)
     115             :               coupled_direction_index = 0;
     116             :             else
     117             :               coupled_direction_index = 2;
     118             :             break;
     119             :           }
     120             :           default: // z
     121             :           {
     122             :             coupled_direction_index = coupled_direction;
     123             :             break;
     124             :           }
     125             :         }
     126             : 
     127        2048 :         val = _Jacobian_mult[_qp](
     128        2048 :                   _direction, _direction, coupled_direction_index, coupled_direction_index) *
     129        2048 :               _test[_i][_qp] * _grad_phi[_j][_qp](coupled_direction_index);
     130             :       }
     131             :     }
     132             :   }
     133             : 
     134             :   // off-diagonal Jacobian with respect to a coupled temperature variable
     135        2048 :   if (_temp_coupled && jvar == _temp_var)
     136             :   {
     137           0 :     RankTwoTensor total_deigenstrain_dT;
     138           0 :     for (const auto deigenstrain_dT : _deigenstrain_dT)
     139           0 :       total_deigenstrain_dT += (*deigenstrain_dT)[_qp];
     140             : 
     141             :     Real sum = 0.0;
     142           0 :     for (unsigned int i = 0; i < 3; ++i)
     143           0 :       sum += _Jacobian_mult[_qp](_direction, _direction, i, i) * total_deigenstrain_dT(i, i);
     144             : 
     145           0 :     val = -sum * _test[_i][_qp] * _phi[_j][_qp];
     146             :   }
     147             : 
     148        2048 :   return val;
     149             : }