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 "ADComputeShellStress.h"
      11             : #include "RankTwoTensor.h"
      12             : #include "RankFourTensor.h"
      13             : #include "ADComputeIsotropicElasticityTensorShell.h"
      14             : 
      15             : #include "libmesh/quadrature.h"
      16             : #include "libmesh/utility.h"
      17             : #include "libmesh/enum_quadrature_type.h"
      18             : #include "libmesh/fe_type.h"
      19             : #include "libmesh/string_to_enum.h"
      20             : #include "libmesh/quadrature_gauss.h"
      21             : 
      22             : registerMooseObject("SolidMechanicsApp", ADComputeShellStress);
      23             : 
      24             : InputParameters
      25         644 : ADComputeShellStress::validParams()
      26             : {
      27         644 :   InputParameters params = Material::validParams();
      28         644 :   params.addClassDescription("Compute in-plane stress using elasticity for shell");
      29        1288 :   params.addRequiredParam<std::string>("through_thickness_order",
      30             :                                        "Quadrature order in out of plane direction");
      31         644 :   return params;
      32           0 : }
      33             : 
      34         483 : ADComputeShellStress::ADComputeShellStress(const InputParameters & parameters)
      35         483 :   : Material(parameters)
      36             : {
      37             :   // get number of quadrature points along thickness based on order
      38             :   std::unique_ptr<libMesh::QGauss> t_qrule = std::make_unique<libMesh::QGauss>(
      39         966 :       1, Utility::string_to_enum<Order>(getParam<std::string>("through_thickness_order")));
      40         483 :   _t_points = t_qrule->get_points();
      41         483 :   _elasticity_tensor.resize(_t_points.size());
      42         483 :   _stress.resize(_t_points.size());
      43         483 :   _stress_old.resize(_t_points.size());
      44         483 :   _strain_increment.resize(_t_points.size());
      45         483 :   _local_transformation_matrix.resize(_t_points.size());
      46         483 :   _covariant_transformation_matrix.resize(_t_points.size());
      47         483 :   _global_stress.resize(_t_points.size());
      48         483 :   _local_stress.resize(_t_points.size());
      49        1449 :   for (unsigned int t = 0; t < _t_points.size(); ++t)
      50             :   {
      51         966 :     _elasticity_tensor[t] =
      52         966 :         &getADMaterialProperty<RankFourTensor>("elasticity_tensor_t_points_" + std::to_string(t));
      53         966 :     _stress[t] = &declareADProperty<RankTwoTensor>("stress_t_points_" + std::to_string(t));
      54         966 :     _stress_old[t] =
      55         966 :         &getMaterialPropertyOldByName<RankTwoTensor>("stress_t_points_" + std::to_string(t));
      56         966 :     _strain_increment[t] =
      57         966 :         &getADMaterialProperty<RankTwoTensor>("strain_increment_t_points_" + std::to_string(t));
      58             :     // rotation matrix and stress for output purposes only
      59         966 :     _local_transformation_matrix[t] =
      60         966 :         &getMaterialProperty<RankTwoTensor>("local_transformation_t_points_" + std::to_string(t));
      61         966 :     _covariant_transformation_matrix[t] = &getMaterialProperty<RankTwoTensor>(
      62         966 :         "covariant_transformation_t_points_" + std::to_string(t));
      63         966 :     _global_stress[t] =
      64         966 :         &declareProperty<RankTwoTensor>("global_stress_t_points_" + std::to_string(t));
      65         966 :     _local_stress[t] =
      66        1932 :         &declareProperty<RankTwoTensor>("local_stress_t_points_" + std::to_string(t));
      67             :   }
      68         483 : }
      69             : 
      70             : void
      71       58896 : ADComputeShellStress::initQpStatefulProperties()
      72             : {
      73             :   // initialize stress tensor to zero
      74      176688 :   for (unsigned int i = 0; i < _t_points.size(); ++i)
      75      117792 :     (*_stress[i])[_qp].zero();
      76       58896 : }
      77             : 
      78             : void
      79      159480 : ADComputeShellStress::computeQpProperties()
      80             : {
      81      478440 :   for (unsigned int i = 0; i < _t_points.size(); ++i)
      82             :   {
      83      318960 :     (*_stress[i])[_qp] =
      84      318960 :         (*_stress_old[i])[_qp] + (*_elasticity_tensor[i])[_qp] * (*_strain_increment[i])[_qp];
      85             : 
      86     1275840 :     for (unsigned int ii = 0; ii < 3; ++ii)
      87     3827520 :       for (unsigned int jj = 0; jj < 3; ++jj)
      88     2870640 :         _unrotated_stress(ii, jj) = MetaPhysicL::raw_value((*_stress[i])[_qp](ii, jj));
      89      318960 :     (*_global_stress[i])[_qp] = (*_covariant_transformation_matrix[i])[_qp].transpose() *
      90      318960 :                                 _unrotated_stress * (*_covariant_transformation_matrix[i])[_qp];
      91      318960 :     (*_local_stress[i])[_qp] = (*_local_transformation_matrix[i])[_qp] * (*_global_stress[i])[_qp] *
      92      318960 :                                (*_local_transformation_matrix[i])[_qp].transpose();
      93             :   }
      94      159480 : }