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