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 "SolidMechanicsPlasticOrthotropic.h"
      11             : #include "RankFourTensor.h"
      12             : #include "libmesh/utility.h"
      13             : 
      14             : registerMooseObject("SolidMechanicsApp", SolidMechanicsPlasticOrthotropic);
      15             : registerMooseObjectRenamed("SolidMechanicsApp",
      16             :                            TensorMechanicsPlasticOrthotropic,
      17             :                            "01/01/2025 00:00",
      18             :                            SolidMechanicsPlasticOrthotropic);
      19             : 
      20             : InputParameters
      21          28 : SolidMechanicsPlasticOrthotropic::validParams()
      22             : {
      23          28 :   InputParameters params = SolidMechanicsPlasticIsotropicSD::validParams();
      24          56 :   params.addRequiredParam<std::vector<Real>>("c1", "The six coefficients of L prime");
      25          56 :   params.addRequiredParam<std::vector<Real>>("c2", "The six coefficients of L prime prime");
      26          28 :   params.addClassDescription("Orthotropic plasticity for pressure sensitive materials and also "
      27             :                              "models the strength differential effect");
      28          28 :   return params;
      29           0 : }
      30             : 
      31          14 : SolidMechanicsPlasticOrthotropic::SolidMechanicsPlasticOrthotropic(
      32          14 :     const InputParameters & parameters)
      33             :   : SolidMechanicsPlasticIsotropicSD(parameters),
      34          14 :     _c1(getParam<std::vector<Real>>("c1")),
      35          56 :     _c2(getParam<std::vector<Real>>("c2"))
      36             : {
      37          14 :   _c = 1.0;
      38          14 :   _l1(0, 0, 0, 0) = (_c1[1] + _c1[2]) / 3.0;
      39          14 :   _l1(0, 0, 1, 1) = -_c1[2] / 3.0;
      40          14 :   _l1(0, 0, 2, 2) = -_c1[1] / 3.0;
      41          14 :   _l1(1, 1, 0, 0) = -_c1[2] / 3.0;
      42          14 :   _l1(1, 1, 1, 1) = (_c1[0] + _c1[2]) / 3.0;
      43          14 :   _l1(1, 1, 2, 2) = -_c1[0] / 3.0;
      44          14 :   _l1(2, 2, 0, 0) = -_c1[1] / 3.0;
      45          14 :   _l1(2, 2, 1, 1) = -_c1[0] / 3.0;
      46          14 :   _l1(2, 2, 2, 2) = (_c1[0] + _c1[1]) / 3.0;
      47          14 :   _l1(0, 1, 1, 0) = _c1[5] / 2.0;
      48          14 :   _l1(0, 1, 0, 1) = _c1[5] / 2.0;
      49          14 :   _l1(1, 0, 1, 0) = _c1[5] / 2.0;
      50          14 :   _l1(1, 0, 0, 1) = _c1[5] / 2.0;
      51          14 :   _l1(0, 2, 0, 2) = _c1[4] / 2.0;
      52          14 :   _l1(0, 2, 2, 0) = _c1[4] / 2.0;
      53          14 :   _l1(2, 0, 2, 0) = _c1[4] / 2.0;
      54          14 :   _l1(2, 0, 0, 2) = _c1[4] / 2.0;
      55          14 :   _l1(1, 2, 2, 1) = _c1[3] / 2.0;
      56          14 :   _l1(1, 2, 1, 2) = _c1[3] / 2.0;
      57          14 :   _l1(2, 1, 1, 2) = _c1[3] / 2.0;
      58          14 :   _l1(2, 1, 2, 1) = _c1[3] / 2.0;
      59             : 
      60          14 :   _l2(0, 0, 0, 0) = (_c2[1] + _c2[2]) / 3.0;
      61          14 :   _l2(0, 0, 1, 1) = -_c2[2] / 3.0;
      62          14 :   _l2(0, 0, 2, 2) = -_c2[1] / 3.0;
      63          14 :   _l2(1, 1, 0, 0) = -_c2[2] / 3.0;
      64          14 :   _l2(1, 1, 1, 1) = (_c2[0] + _c2[2]) / 3.0;
      65          14 :   _l2(1, 1, 2, 2) = -_c2[0] / 3.0;
      66          14 :   _l2(2, 2, 0, 0) = -_c2[1] / 3.0;
      67          14 :   _l2(2, 2, 1, 1) = -_c2[0] / 3.0;
      68          14 :   _l2(2, 2, 2, 2) = (_c2[0] + _c2[1]) / 3.0;
      69          14 :   _l2(0, 1, 1, 0) = _c2[5] / 2.0;
      70          14 :   _l2(0, 1, 0, 1) = _c2[5] / 2.0;
      71          14 :   _l2(1, 0, 1, 0) = _c2[5] / 2.0;
      72          14 :   _l2(1, 0, 0, 1) = _c2[5] / 2.0;
      73          14 :   _l2(0, 2, 0, 2) = _c2[4] / 2.0;
      74          14 :   _l2(0, 2, 2, 0) = _c2[4] / 2.0;
      75          14 :   _l2(2, 0, 2, 0) = _c2[4] / 2.0;
      76          14 :   _l2(2, 0, 0, 2) = _c2[4] / 2.0;
      77          14 :   _l2(1, 2, 2, 1) = _c2[3] / 2.0;
      78          14 :   _l2(1, 2, 1, 2) = _c2[3] / 2.0;
      79          14 :   _l2(2, 1, 1, 2) = _c2[3] / 2.0;
      80          14 :   _l2(2, 1, 2, 1) = _c2[3] / 2.0;
      81          14 : }
      82             : 
      83             : Real
      84       11456 : SolidMechanicsPlasticOrthotropic::yieldFunction(const RankTwoTensor & stress, Real intnl) const
      85             : {
      86       11456 :   const RankTwoTensor j2prime = _l1 * stress;
      87       11456 :   const RankTwoTensor j3prime = _l2 * stress;
      88       11456 :   return _b * stress.trace() +
      89       11456 :          std::pow(std::pow(-j2prime.generalSecondInvariant(), 3.0 / 2.0) - j3prime.det(),
      90             :                   1.0 / 3.0) -
      91       11456 :          yieldStrength(intnl);
      92             : }
      93             : 
      94             : RankTwoTensor
      95        9184 : SolidMechanicsPlasticOrthotropic::dyieldFunction_dstress(const RankTwoTensor & stress,
      96             :                                                          Real /*intnl*/) const
      97             : {
      98        9184 :   const RankTwoTensor j2prime = _l1 * stress;
      99        9184 :   const RankTwoTensor j3prime = _l2 * stress;
     100        9184 :   const Real j2 = -j2prime.generalSecondInvariant();
     101        9184 :   const Real j3 = j3prime.det();
     102        9184 :   return _b * dI_sigma() + dphi_dj2(j2, j3) * _l1.innerProductTranspose(dj2_dSkl(j2prime)) +
     103        9184 :          dphi_dj3(j2, j3) * _l2.innerProductTranspose(j3prime.ddet());
     104             : }
     105             : 
     106             : RankFourTensor
     107        2296 : SolidMechanicsPlasticOrthotropic::dflowPotential_dstress(const RankTwoTensor & stress,
     108             :                                                          Real /*intnl*/) const
     109             : {
     110        2296 :   if (_associative)
     111             :   {
     112        2296 :     const RankTwoTensor j2prime = _l1 * stress;
     113        2296 :     const RankTwoTensor j3prime = _l2 * stress;
     114        2296 :     const RankTwoTensor dj2 = dj2_dSkl(j2prime);
     115        2296 :     const RankTwoTensor dj3 = j3prime.ddet();
     116        2296 :     const Real j2 = -j2prime.generalSecondInvariant();
     117        2296 :     const Real j3 = j3prime.det();
     118             :     const RankFourTensor dr =
     119        2296 :         dfj2_dj2(j2, j3) *
     120        2296 :             _l1.innerProductTranspose(dj2).outerProduct(_l1.innerProductTranspose(dj2)) +
     121        2296 :         dfj2_dj3(j2, j3) *
     122        4592 :             _l1.innerProductTranspose(dj2).outerProduct(_l2.innerProductTranspose(dj3)) +
     123        2296 :         dfj3_dj2(j2, j3) *
     124        4592 :             _l2.innerProductTranspose(dj3).outerProduct(_l1.innerProductTranspose(dj2)) +
     125        2296 :         dfj3_dj3(j2, j3) *
     126        4592 :             _l2.innerProductTranspose(dj3).outerProduct(_l2.innerProductTranspose(dj3));
     127        2296 :     const RankTwoTensor r = _b * dI_sigma() +
     128        2296 :                             dphi_dj2(j2, j3) * _l1.innerProductTranspose(dj2_dSkl(j2prime)) +
     129        2296 :                             dphi_dj3(j2, j3) * _l2.innerProductTranspose(j3prime.ddet());
     130        2296 :     const Real norm = r.L2norm();
     131        4592 :     return dr / norm - (r / Utility::pow<3>(norm)).outerProduct(dr.innerProductTranspose(r));
     132             :   }
     133             :   else
     134           0 :     return SolidMechanicsPlasticJ2::dflowPotential_dstress(stress, 0);
     135             : }
     136             : 
     137             : RankTwoTensor
     138        6888 : SolidMechanicsPlasticOrthotropic::flowPotential(const RankTwoTensor & stress, Real intnl) const
     139             : {
     140        6888 :   if (_associative)
     141             :   {
     142        6888 :     const RankTwoTensor a = dyieldFunction_dstress(stress, intnl);
     143        6888 :     return a / a.L2norm();
     144             :   }
     145             :   else
     146           0 :     return SolidMechanicsPlasticJ2::dyieldFunction_dstress(stress, intnl);
     147             : }