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 "CZMComputeDisplacementJumpTotalLagrangian.h"
      11             : #include "CohesiveZoneModelTools.h"
      12             : #include "FactorizedRankTwoTensor.h"
      13             : 
      14             : registerMooseObject("SolidMechanicsApp", CZMComputeDisplacementJumpTotalLagrangian);
      15             : registerMooseObject("SolidMechanicsApp", ADCZMComputeDisplacementJumpTotalLagrangian);
      16             : 
      17             : template <bool is_ad>
      18             : InputParameters
      19          88 : CZMComputeDisplacementJumpTotalLagrangianTempl<is_ad>::validParams()
      20             : {
      21          88 :   InputParameters params = CZMComputeDisplacementJumpBase<is_ad>::validParams();
      22          88 :   params.addClassDescription(
      23             :       "Compute the displacement jump increment across a czm interface in local "
      24             :       "coordinates for the Total Lagrangian kinematic formulation");
      25             : 
      26          88 :   return params;
      27           0 : }
      28             : 
      29             : template <bool is_ad>
      30          44 : CZMComputeDisplacementJumpTotalLagrangianTempl<
      31             :     is_ad>::CZMComputeDisplacementJumpTotalLagrangianTempl(const InputParameters & parameters)
      32             :   : CZMComputeDisplacementJumpBase<is_ad>(parameters),
      33          44 :     _F(this->template declareGenericPropertyByName<RankTwoTensor, is_ad>(_base_name + "F_czm")),
      34          44 :     _R(this->template declareGenericPropertyByName<RankTwoTensor, is_ad>(_base_name +
      35             :                                                                          "czm_rotation")),
      36          44 :     _czm_reference_rotation(this->template declareGenericPropertyByName<RankTwoTensor, is_ad>(
      37          44 :         _base_name + "czm_reference_rotation"))
      38             : {
      39             :   // initializing the displacement gradient vectors
      40         164 :   for (unsigned int i = 0; i < _ndisp; ++i)
      41             :   {
      42         120 :     _grad_disp.push_back(&this->template coupledGenericGradient<is_ad>("displacements", i));
      43         120 :     _grad_disp_neighbor.push_back(
      44         360 :         &this->template coupledGenericNeighborGradient<is_ad>("displacements", i));
      45             :   }
      46             : 
      47             :   // All others zero (so this will work naturally for 2D and 1D problems)
      48          56 :   for (unsigned int i = _ndisp; i < 3; i++)
      49             :   {
      50             :     if constexpr (is_ad)
      51             :     {
      52           0 :       _grad_disp.push_back(&_ad_grad_zero);
      53           0 :       _grad_disp_neighbor.push_back(&_ad_grad_zero);
      54             :     }
      55             :     else
      56             :     {
      57          12 :       _grad_disp.push_back(&_grad_zero);
      58          12 :       _grad_disp_neighbor.push_back(&_grad_zero);
      59             :     }
      60             :   }
      61          44 : }
      62             : 
      63             : template <bool is_ad>
      64             : void
      65         480 : CZMComputeDisplacementJumpTotalLagrangianTempl<is_ad>::initQpStatefulProperties()
      66             : {
      67         480 :   CZMComputeDisplacementJumpBase<is_ad>::initQpStatefulProperties();
      68         480 :   _displacement_jump_global[_qp] = 0;
      69         480 :   _F[_qp].setToIdentity();
      70         480 :   _R[_qp].setToIdentity();
      71         480 : }
      72             : 
      73             : template <bool is_ad>
      74             : void
      75       32912 : CZMComputeDisplacementJumpTotalLagrangianTempl<is_ad>::computeLocalDisplacementJump()
      76             : {
      77       39392 :   _interface_displacement_jump[_qp] =
      78       32912 :       _czm_total_rotation[_qp].transpose() * _displacement_jump_global[_qp];
      79       32912 : }
      80             : 
      81             : template <bool is_ad>
      82             : void
      83       32912 : CZMComputeDisplacementJumpTotalLagrangianTempl<is_ad>::computeRotationMatrices()
      84             : {
      85       32912 :   _czm_reference_rotation[_qp] = CohesiveZoneModelTools::computeReferenceRotation<is_ad>(
      86       32912 :       _normals[_qp], this->_mesh.dimension());
      87        6480 :   computeFandR();
      88       32912 :   _czm_total_rotation[_qp] = _R[_qp] * _czm_reference_rotation[_qp];
      89       32912 : }
      90             : 
      91             : template <bool is_ad>
      92             : void
      93       32912 : CZMComputeDisplacementJumpTotalLagrangianTempl<is_ad>::computeFandR()
      94             : {
      95       65824 :   const auto F = (GenericRankTwoTensor<is_ad>::Identity() +
      96             :                   GenericRankTwoTensor<is_ad>::initializeFromRows(
      97       32912 :                       (*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]));
      98       32912 :   const auto F_neighbor =
      99             :       (GenericRankTwoTensor<is_ad>::Identity() +
     100       32912 :        GenericRankTwoTensor<is_ad>::initializeFromRows((*_grad_disp_neighbor[0])[_qp],
     101       32912 :                                                        (*_grad_disp_neighbor[1])[_qp],
     102       32912 :                                                        (*_grad_disp_neighbor[2])[_qp]));
     103             : 
     104       39392 :   _F[_qp] = 0.5 * (F + F_neighbor);
     105             :   // According to Cody mooseError are always fatal, so nothing we can do about
     106             :   // them. The norm of the tensor might work, but there is the risk of an
     107             :   // unwanted overflow (I tried and it happens). So checking component by
     108             :   // component remains the only reasonable strategy. If someone finds a  better
     109             :   // way this could be changed in the future
     110      131648 :   for (unsigned int i = 0; i < 3; i++)
     111      394944 :     for (unsigned int j = 0; j < 3; j++)
     112      296208 :       if (!std::isfinite(MetaPhysicL::raw_value(_F[_qp](i, j))))
     113           0 :         throw MooseException("CZMMaterialBaseIncremental _F is not finite, reducing time step");
     114             : 
     115             :   using FR2T = std::conditional_t<is_ad, ADFactorizedRankTwoTensor, FactorizedRankTwoTensor>;
     116       32912 :   const FR2T C = _F[_qp].transpose() * _F[_qp];
     117       32912 :   const auto Uinv = MathUtils::sqrt(C).inverse().get();
     118       32912 :   _R[_qp] = _F[_qp] * Uinv;
     119       32912 : }
     120             : 
     121             : template class CZMComputeDisplacementJumpTotalLagrangianTempl<false>;
     122             : template class CZMComputeDisplacementJumpTotalLagrangianTempl<true>;