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 "Compute2DFiniteStrain.h"
      11             : 
      12             : #include "libmesh/quadrature.h"
      13             : 
      14             : InputParameters
      15        1632 : Compute2DFiniteStrain::validParams()
      16             : {
      17        1632 :   InputParameters params = ComputeFiniteStrain::validParams();
      18        1632 :   params.addClassDescription(
      19             :       "Compute a strain increment and rotation increment for finite strains in 2D geometries.");
      20             : 
      21        3264 :   MooseEnum outOfPlaneDirection("x y z", "z");
      22        3264 :   params.addParam<MooseEnum>(
      23             :       "out_of_plane_direction", outOfPlaneDirection, "The direction of the out-of-plane strain.");
      24        1632 :   return params;
      25        1632 : }
      26             : 
      27        1224 : Compute2DFiniteStrain::Compute2DFiniteStrain(const InputParameters & parameters)
      28             :   : ComputeFiniteStrain(parameters),
      29        2448 :     _out_of_plane_direction(getParam<MooseEnum>("out_of_plane_direction"))
      30             : {
      31        1224 : }
      32             : 
      33             : void
      34        1158 : Compute2DFiniteStrain::initialSetup()
      35             : {
      36        4632 :   for (unsigned int i = 0; i < 3; ++i)
      37             :   {
      38        3474 :     if (_out_of_plane_direction == i)
      39             :     {
      40        1158 :       _disp[i] = &_zero;
      41        1158 :       _grad_disp[i] = &_grad_zero;
      42             :     }
      43             :     else
      44             :     {
      45        2316 :       _disp[i] = &coupledValue("displacements", i);
      46        2316 :       _grad_disp[i] = &coupledGradient("displacements", i);
      47             :     }
      48             : 
      49        3474 :     if (_fe_problem.isTransient() && i != _out_of_plane_direction)
      50        2316 :       _grad_disp_old[i] = &coupledGradientOld("displacements", i);
      51             :     else
      52        1158 :       _grad_disp_old[i] = &_grad_zero;
      53             :   }
      54        1158 : }
      55             : 
      56             : void
      57     2094826 : Compute2DFiniteStrain::computeProperties()
      58             : {
      59     2094826 :   RankTwoTensor ave_Fhat;
      60             :   Real ave_dfgrd_det = 0.0;
      61             : 
      62    10460690 :   for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
      63             :   {
      64             :     // Deformation gradient calculation for 2D problems
      65             :     auto A = RankTwoTensor::initializeFromRows(
      66     8365864 :         (*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
      67             : 
      68             :     // Old Deformation gradient
      69             :     auto Fbar = RankTwoTensor::initializeFromRows(
      70     8365864 :         (*_grad_disp_old[0])[_qp], (*_grad_disp_old[1])[_qp], (*_grad_disp_old[2])[_qp]);
      71             : 
      72             :     // Compute the displacement gradient for the out of plane direction for plane strain,
      73             :     // generalized plane strain, or axisymmetric problems
      74             : 
      75     8365864 :     A(_out_of_plane_direction, _out_of_plane_direction) = computeOutOfPlaneGradDisp();
      76     8365864 :     Fbar(_out_of_plane_direction, _out_of_plane_direction) = computeOutOfPlaneGradDispOld();
      77             : 
      78             :     // Gauss point deformation gradient
      79     8365864 :     _deformation_gradient[_qp] = A;
      80     8365864 :     _deformation_gradient[_qp].addIa(1.0);
      81             : 
      82             :     // deformation gradient midpoint
      83     8365864 :     if (_use_hw)
      84             :     {
      85           0 :       (*_def_grad_mid)[_qp].setToIdentity();
      86           0 :       (*_def_grad_mid)[_qp] += 0.5 * (A + Fbar);
      87             :     }
      88             : 
      89     8365864 :     A -= Fbar; // very nearly A = gradU - gradUold
      90             : 
      91             :     // _f_bar = dDU/dX_o
      92     8365864 :     if (_use_hw)
      93           0 :       (*_f_bar)[_qp] = A;
      94             : 
      95     8365864 :     Fbar.addIa(1.0); // Fbar = ( I + gradUold)
      96             : 
      97             :     // Incremental deformation gradient _Fhat = I + A Fbar^-1
      98     8365864 :     _Fhat[_qp] = A * Fbar.inverse();
      99     8365864 :     _Fhat[_qp].addIa(1.0);
     100             : 
     101     8365864 :     if (_volumetric_locking_correction)
     102             :     {
     103             :       // Calculate average _Fhat for volumetric locking correction
     104     6784144 :       ave_Fhat += _Fhat[_qp] * _JxW[_qp] * _coord[_qp];
     105             : 
     106             :       // Average deformation gradient
     107     6784144 :       ave_dfgrd_det += _deformation_gradient[_qp].det() * _JxW[_qp] * _coord[_qp];
     108             :     }
     109             :   }
     110     2094826 :   if (_volumetric_locking_correction)
     111             :   {
     112             :     // needed for volumetric locking correction
     113     1696036 :     ave_Fhat /= _current_elem_volume;
     114             :     // average deformation gradient
     115     1696036 :     ave_dfgrd_det /= _current_elem_volume;
     116             :   }
     117    10460690 :   for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
     118             :   {
     119     8365864 :     if (_volumetric_locking_correction)
     120             :     {
     121             :       // Finalize volumetric locking correction
     122     6784144 :       _Fhat[_qp] *= std::cbrt(ave_Fhat.det() / _Fhat[_qp].det());
     123             :       // Volumetric locking correction
     124     6784144 :       _deformation_gradient[_qp] *= std::cbrt(ave_dfgrd_det / _deformation_gradient[_qp].det());
     125             :     }
     126             : 
     127     8365864 :     computeQpStrain();
     128             :   }
     129     2094826 : }
     130             : 
     131             : void
     132           0 : Compute2DFiniteStrain::displacementIntegrityCheck()
     133             : {
     134           0 :   if (_out_of_plane_direction != 2 && _ndisp != 3)
     135           0 :     mooseError("For 2D simulations where the out-of-plane direction is x or y the number of "
     136             :                "supplied displacements must be three.");
     137           0 :   else if (_out_of_plane_direction == 2 && _ndisp != 2)
     138           0 :     mooseError("For 2D simulations where the out-of-plane direction is z the number of supplied "
     139             :                "displacements must be two.");
     140           0 : }