Line data    Source code

       1             : //* This file is part of the MOOSE framework
       2             : //* https://www.mooseframework.org
       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 "StressDivergenceRSphericalTensors.h"
      11             : #include "ElasticityTensorTools.h"
      12             : #include "FEProblem.h"
      13             : #include "MooseMesh.h"
      14             : 
      15             : registerMooseObject("TensorMechanicsApp", StressDivergenceRSphericalTensors);
      16             : 
      17             : InputParameters
      18          60 : StressDivergenceRSphericalTensors::validParams()
      19             : {
      20          60 :   InputParameters params = StressDivergenceTensors::validParams();
      21          60 :   params.addClassDescription(
      22             :       "Calculate stress divergence for a spherically symmetric 1D problem in polar coordinates.");
      23          60 :   params.set<unsigned int>("component") = 0;
      24          60 :   params.set<bool>("use_displaced_mesh") = true;
      25          60 :   return params;
      26           0 : }
      27             : 
      28          30 : StressDivergenceRSphericalTensors::StressDivergenceRSphericalTensors(
      29          30 :     const InputParameters & parameters)
      30          30 :   : StressDivergenceTensors(parameters)
      31             : {
      32          30 :   if (_component != 0)
      33           0 :     mooseError("Invalid component for this 1D RSpherical problem.");
      34          30 : }
      35             : 
      36             : void
      37          30 : StressDivergenceRSphericalTensors::initialSetup()
      38             : {
      39          30 :   if (getBlockCoordSystem() != Moose::COORD_RSPHERICAL)
      40           0 :     mooseError("The coordinate system in the Problem block must be set to RSPHERICAL for 1D "
      41             :                "spherically symmetric geometries.");
      42          30 : }
      43             : 
      44             : Real
      45        4088 : StressDivergenceRSphericalTensors::computeQpResidual()
      46             : {
      47        4088 :   return _grad_test[_i][_qp](0) * _stress[_qp](0, 0) +               // stress_{rr} part 1
      48        4088 :          +(_test[_i][_qp] / _q_point[_qp](0)) * _stress[_qp](1, 1) + // stress_{\theta \theta}
      49        4088 :          +(_test[_i][_qp] / _q_point[_qp](0)) * _stress[_qp](2, 2);  // stress_{\phi \phi}
      50             : }
      51             : 
      52             : Real
      53        2956 : StressDivergenceRSphericalTensors::computeQpJacobian()
      54             : {
      55        2956 :   return calculateJacobian(_component, _component);
      56             : }
      57             : 
      58             : Real
      59           0 : StressDivergenceRSphericalTensors::computeQpOffDiagJacobian(unsigned int jvar)
      60             : {
      61           0 :   for (unsigned int i = 0; i < _ndisp; ++i)
      62           0 :     if (jvar == _disp_var[i])
      63           0 :       return calculateJacobian(_component, i);
      64             : 
      65             :   return 0.0;
      66             : }
      67             : 
      68             : Real
      69        2956 : StressDivergenceRSphericalTensors::calculateJacobian(unsigned int libmesh_dbg_var(ivar),
      70             :                                                      unsigned int libmesh_dbg_var(jvar))
      71             : {
      72             :   mooseAssert(ivar == 0 && jvar == 0,
      73             :               "Invalid component in Jacobian Calculation"); // Only nonzero case for a 1D simulation
      74             : 
      75        2956 :   const Real test = _grad_test[_i][_qp](0);
      76        2956 :   const Real test_r1 = _test[_i][_qp] / _q_point[_qp](0);
      77             :   // const Real test_r2 = test_r1;
      78             : 
      79        2956 :   const Real phi = _grad_phi[_j][_qp](0);
      80        2956 :   const Real phi_r1 = _phi[_j][_qp] / _q_point[_qp](0);
      81             :   const Real phi_r2 = phi_r1;
      82             : 
      83        2956 :   const Real term1 = test * _Jacobian_mult[_qp](0, 0, 0, 0) * phi;
      84        2956 :   const Real term2 = test * _Jacobian_mult[_qp](0, 0, 1, 1) * phi_r1; // same as term3
      85             :   // const Real term3 = test * _Jacobian_mult[_qp](0, 0, 2, 2) * phi_r2;
      86             : 
      87        2956 :   const Real term4 = test_r1 * _Jacobian_mult[_qp](1, 1, 0, 0) * phi;    // same as term7
      88        2956 :   const Real term5 = test_r1 * _Jacobian_mult[_qp](1, 1, 1, 1) * phi_r1; // same as term9
      89        2956 :   const Real term6 = test_r1 * _Jacobian_mult[_qp](1, 1, 2, 2) * phi_r2; // same as term8
      90             : 
      91             :   // const Real term7 = test_r2 * _Jacobian_mult[_qp](2, 2, 0, 0) * phi;
      92             :   // const Real term8 = test_r2 * _Jacobian_mult[_qp](2, 2, 1, 1) * phi_r1;
      93             :   // const Real term9 = test_r2 * _Jacobian_mult[_qp](2, 2, 2, 2) * phi_r2;
      94             : 
      95        2956 :   return term1 + 2 * (term2 + term4 + term5 + term6);
      96             : }