www.mooseframework.org
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
StressDivergenceRZTensors Class Reference

StressDivergenceRZTensors is a modification of StressDivergenceTensors to accommodate the Axisymmetric material models that use cylindrical coordinates. More...

#include <StressDivergenceRZTensors.h>

Inheritance diagram for StressDivergenceRZTensors:
[legend]

Public Member Functions

 StressDivergenceRZTensors (const InputParameters &parameters)
 
virtual void computeJacobian () override
 
virtual void computeOffDiagJacobian (MooseVariableFEBase &jvar) override
 

Static Public Member Functions

static InputParameters validParams ()
 

Protected Member Functions

virtual void initialSetup () override
 
virtual Real computeQpResidual () override
 
virtual Real computeQpJacobian () override
 
virtual Real computeQpOffDiagJacobian (unsigned int jvar) override
 
virtual void computeAverageGradientTest () override
 
virtual void computeAverageGradientPhi () override
 
Real calculateJacobian (unsigned int ivar, unsigned int jvar)
 
virtual void computeResidual () override
 
virtual void computeFiniteDeformJacobian ()
 

Protected Attributes

const std::string _base_name
 
bool _use_finite_deform_jacobian
 
const MaterialProperty< RankTwoTensor > & _stress
 
const MaterialProperty< RankFourTensor > & _Jacobian_mult
 
std::vector< RankFourTensor_finite_deform_Jacobian_mult
 
const MaterialProperty< RankTwoTensor > * _deformation_gradient
 
const MaterialProperty< RankTwoTensor > * _deformation_gradient_old
 
const MaterialProperty< RankTwoTensor > * _rotation_increment
 
const unsigned int _component
 
unsigned int _ndisp
 Coupled displacement variables. More...
 
std::vector< unsigned int > _disp_var
 
const bool _temp_coupled
 
const unsigned int _temp_var
 
const MaterialProperty< RankTwoTensor > *const _deigenstrain_dT
 d(strain)/d(temperature), if computed by ComputeThermalExpansionEigenstrain More...
 
const bool _out_of_plane_strain_coupled
 
const VariableValue * _out_of_plane_strain
 
const unsigned int _out_of_plane_strain_var
 
const unsigned int _out_of_plane_direction
 
const bool _use_displaced_mesh
 Whether this object is acting on the displaced mesh. More...
 
std::vector< std::vector< Real > > _avg_grad_test
 Gradient of test function averaged over the element. Used in volumetric locking correction calculation. More...
 
std::vector< std::vector< Real > > _avg_grad_phi
 Gradient of phi function averaged over the element. Used in volumetric locking correction calculation. More...
 
bool _volumetric_locking_correction
 Flag for volumetric locking correction. More...
 
Assembly & _assembly_undisplaced
 undisplaced problem More...
 
MooseVariable & _var_undisplaced
 Reference to this Kernel's undisplaced MooseVariable object. More...
 
const VariablePhiGradient & _grad_phi_undisplaced
 Shape and test functions on the undisplaced mesh. More...
 
const VariableTestGradient & _grad_test_undisplaced
 

Detailed Description

StressDivergenceRZTensors is a modification of StressDivergenceTensors to accommodate the Axisymmetric material models that use cylindrical coordinates.

This kernel is for symmetrical loading only. The key modifications are a result of the circumferential stress' dependence on displacement in the axial direction. Reference: Cook et.al. Concepts and Applications of Finite Element Analysis, 4th Ed. 2002. p 510. Within this kernel, '_disp_x' refers to displacement in the radial direction, u_r, and '_disp_y' refers to displacement in the axial direction, u_z. The COORD_TYPE in the Problem block must be set to RZ.

Definition at line 31 of file StressDivergenceRZTensors.h.

Constructor & Destructor Documentation

◆ StressDivergenceRZTensors()

StressDivergenceRZTensors::StressDivergenceRZTensors ( const InputParameters &  parameters)

Definition at line 34 of file StressDivergenceRZTensors.C.

35  : StressDivergenceTensors(parameters)
36 {
37 }

Member Function Documentation

◆ calculateJacobian()

Real StressDivergenceRZTensors::calculateJacobian ( unsigned int  ivar,
unsigned int  jvar 
)
protected

Definition at line 123 of file StressDivergenceRZTensors.C.

124 {
125  // B^T_i * C * B_j
126  RealGradient test, test_z, phi, phi_z;
127  Real first_term = 0.0;
128  if (ivar == 0) // Case grad_test for x, requires contributions from stress_xx, stress_xy, and
129  // stress_zz
130  {
131  test(0) = _grad_test[_i][_qp](0);
132  test(1) = _grad_test[_i][_qp](1);
133  test_z(2) = _test[_i][_qp] / _q_point[_qp](0);
134  }
135  else // Case grad_test for y
136  {
137  test(0) = _grad_test[_i][_qp](0);
138  test(1) = _grad_test[_i][_qp](1);
139  }
140 
141  if (jvar == 0)
142  {
143  phi(0) = _grad_phi[_j][_qp](0);
144  phi(1) = _grad_phi[_j][_qp](1);
145  phi_z(2) = _phi[_j][_qp] / _q_point[_qp](0);
146  }
147  else
148  {
149  phi(0) = _grad_phi[_j][_qp](0);
150  phi(1) = _grad_phi[_j][_qp](1);
151  }
152 
153  if (ivar == 0 &&
154  jvar == 0) // Case when both phi and test are functions of x and z; requires four terms
155  {
156  const Real first_sum = ElasticityTensorTools::elasticJacobian(
157  _Jacobian_mult[_qp], ivar, jvar, test, phi); // test_x and phi_x
158  const Real second_sum = ElasticityTensorTools::elasticJacobian(
159  _Jacobian_mult[_qp], 2, 2, test_z, phi_z); // test_z and phi_z
160  const Real mixed_sum1 = ElasticityTensorTools::elasticJacobian(
161  _Jacobian_mult[_qp], ivar, 2, test, phi_z); // test_x and phi_z
162  const Real mixed_sum2 = ElasticityTensorTools::elasticJacobian(
163  _Jacobian_mult[_qp], 2, jvar, test_z, phi); // test_z and phi_x
164 
165  first_term = first_sum + second_sum + mixed_sum1 + mixed_sum2;
166  }
167  else if (ivar == 0 && jvar == 1)
168  {
169  const Real first_sum = ElasticityTensorTools::elasticJacobian(
170  _Jacobian_mult[_qp], ivar, jvar, test, phi); // test_x and phi_y
171  const Real mixed_sum2 = ElasticityTensorTools::elasticJacobian(
172  _Jacobian_mult[_qp], 2, jvar, test_z, phi); // test_z and phi_y
173 
174  first_term = first_sum + mixed_sum2;
175  }
176  else if (ivar == 1 && jvar == 0)
177  {
178  const Real second_sum = ElasticityTensorTools::elasticJacobian(
179  _Jacobian_mult[_qp], ivar, jvar, test, phi); // test_y and phi_x
180  const Real mixed_sum1 = ElasticityTensorTools::elasticJacobian(
181  _Jacobian_mult[_qp], ivar, 2, test, phi_z); // test_y and phi_z
182 
183  first_term = second_sum + mixed_sum1;
184  }
185  else if (ivar == 1 && jvar == 1)
187  _Jacobian_mult[_qp], ivar, jvar, test, phi); // test_y and phi_y
188  else
189  mooseError("Invalid component in Jacobian Calculation");
190 
191  Real val = 0.0;
192  // volumetric locking correction
193  // K = Bbar^T_i * C * Bbar^T_j where Bbar = B + Bvol
194  // K = B^T_i * C * B_j + Bvol^T_i * C * Bvol_j + B^T_i * C * Bvol_j + Bvol^T_i * C * B_j
196  {
197  RealGradient new_test(2, 0.0);
198  RealGradient new_phi(2, 0.0);
199 
200  new_test(0) = _grad_test[_i][_qp](0) + _test[_i][_qp] / _q_point[_qp](0);
201  new_test(1) = _grad_test[_i][_qp](1);
202  new_phi(0) = _grad_phi[_j][_qp](0) + _phi[_j][_qp] / _q_point[_qp](0);
203  new_phi(1) = _grad_phi[_j][_qp](1);
204 
205  // Bvol^T_i * C * Bvol_j
206  val += _Jacobian_mult[_qp].sum3x3() * (_avg_grad_test[_i][ivar] - new_test(ivar)) *
207  (_avg_grad_phi[_j][jvar] - new_phi(jvar)) / 3.0;
208 
209  // B^T_i * C * Bvol_j
210  RealGradient sum_3x1 = _Jacobian_mult[_qp].sum3x1();
211  if (ivar == 0 && jvar == 0)
212  val += (sum_3x1(0) * test(0) + sum_3x1(2) * test_z(2)) * (_avg_grad_phi[_j][0] - new_phi(0));
213  else if (ivar == 0 && jvar == 1)
214  val += (sum_3x1(0) * test(0) + sum_3x1(2) * test_z(2)) * (_avg_grad_phi[_j][1] - new_phi(1));
215  else if (ivar == 1 && jvar == 0)
216  val += sum_3x1(1) * test(1) * (_avg_grad_phi[_j][0] - new_phi(0));
217  else
218  val += sum_3x1(1) * test(1) * (_avg_grad_phi[_j][1] - new_phi(1));
219 
220  // Bvol^T_i * C * B_j
221  // val = trace (C * B_j) *(avg_grad_test[_i][ivar] - new_test(ivar))
222  if (jvar == 0)
223  for (unsigned int i = 0; i < 3; ++i)
224  val +=
225  (_Jacobian_mult[_qp](i, i, 0, 0) * phi(0) + _Jacobian_mult[_qp](i, i, 0, 1) * phi(1) +
226  _Jacobian_mult[_qp](i, i, 2, 2) * phi_z(2)) *
227  (_avg_grad_test[_i][ivar] - new_test(ivar));
228  else if (jvar == 1)
229  for (unsigned int i = 0; i < 3; ++i)
230  val +=
231  (_Jacobian_mult[_qp](i, i, 0, 1) * phi(0) + _Jacobian_mult[_qp](i, i, 1, 1) * phi(1)) *
232  (_avg_grad_test[_i][ivar] - new_test(ivar));
233  }
234 
235  return val / 3.0 + first_term;
236 }

Referenced by computeQpJacobian(), and computeQpOffDiagJacobian().

◆ computeAverageGradientPhi()

void StressDivergenceRZTensors::computeAverageGradientPhi ( )
overrideprotectedvirtual

Reimplemented from StressDivergenceTensors.

Definition at line 261 of file StressDivergenceRZTensors.C.

262 {
263  _avg_grad_phi.resize(_phi.size());
264  for (_i = 0; _i < _phi.size(); ++_i)
265  {
266  _avg_grad_phi[_i].resize(2);
267  for (unsigned int component = 0; component < 2; ++component)
268  {
269  _avg_grad_phi[_i][component] = 0.0;
270  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
271  {
272  if (component == 0)
273  _avg_grad_phi[_i][component] +=
274  (_grad_phi[_i][_qp](component) + _phi[_i][_qp] / _q_point[_qp](0)) * _JxW[_qp] *
275  _coord[_qp];
276  else
277  _avg_grad_phi[_i][component] += _grad_phi[_i][_qp](component) * _JxW[_qp] * _coord[_qp];
278  }
279  _avg_grad_phi[_i][component] /= _current_elem_volume;
280  }
281  }
282 }

◆ computeAverageGradientTest()

void StressDivergenceRZTensors::computeAverageGradientTest ( )
overrideprotectedvirtual

Reimplemented from StressDivergenceTensors.

Definition at line 239 of file StressDivergenceRZTensors.C.

240 {
241  // calculate volume averaged value of shape function derivative
242  _avg_grad_test.resize(_test.size());
243  for (_i = 0; _i < _test.size(); ++_i)
244  {
245  _avg_grad_test[_i].resize(2);
246  _avg_grad_test[_i][_component] = 0.0;
247  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
248  {
249  if (_component == 0)
250  _avg_grad_test[_i][_component] +=
251  (_grad_test[_i][_qp](_component) + _test[_i][_qp] / _q_point[_qp](0)) * _JxW[_qp] *
252  _coord[_qp];
253  else
254  _avg_grad_test[_i][_component] += _grad_test[_i][_qp](_component) * _JxW[_qp] * _coord[_qp];
255  }
256  _avg_grad_test[_i][_component] /= _current_elem_volume;
257  }
258 }

◆ computeFiniteDeformJacobian()

void StressDivergenceTensors::computeFiniteDeformJacobian ( )
protectedvirtualinherited

Definition at line 350 of file StressDivergenceTensors.C.

351 {
352  const RankTwoTensor I(RankTwoTensor::initIdentity);
353  const RankFourTensor II_ijkl = I.mixedProductIkJl(I);
354 
355  // Bring back to unrotated config
356  const RankTwoTensor unrotated_stress =
357  (*_rotation_increment)[_qp].transpose() * _stress[_qp] * (*_rotation_increment)[_qp];
358 
359  // Incremental deformation gradient Fhat
360  const RankTwoTensor Fhat =
361  (*_deformation_gradient)[_qp] * (*_deformation_gradient_old)[_qp].inverse();
362  const RankTwoTensor Fhatinv = Fhat.inverse();
363 
364  const RankTwoTensor rot_times_stress = (*_rotation_increment)[_qp] * unrotated_stress;
365  const RankFourTensor dstress_drot =
366  I.mixedProductIkJl(rot_times_stress) + I.mixedProductJkIl(rot_times_stress);
367  const RankFourTensor rot_rank_four =
368  (*_rotation_increment)[_qp].mixedProductIkJl((*_rotation_increment)[_qp]);
369  const RankFourTensor drot_dUhatinv = Fhat.mixedProductIkJl(I);
370 
371  const RankTwoTensor A = I - Fhatinv;
372 
373  // Ctilde = Chat^-1 - I
374  const RankTwoTensor Ctilde = A * A.transpose() - A - A.transpose();
375  const RankFourTensor dCtilde_dFhatinv =
376  -I.mixedProductIkJl(A) - I.mixedProductJkIl(A) + II_ijkl + I.mixedProductJkIl(I);
377 
378  // Second order approximation of Uhat - consistent with strain increment definition
379  // const RankTwoTensor Uhat = I - 0.5 * Ctilde - 3.0/8.0 * Ctilde * Ctilde;
380 
381  RankFourTensor dUhatinv_dCtilde =
382  0.5 * II_ijkl - 1.0 / 8.0 * (I.mixedProductIkJl(Ctilde) + Ctilde.mixedProductIkJl(I));
383  RankFourTensor drot_dFhatinv = drot_dUhatinv * dUhatinv_dCtilde * dCtilde_dFhatinv;
384 
385  drot_dFhatinv -= Fhat.mixedProductIkJl((*_rotation_increment)[_qp].transpose());
386  _finite_deform_Jacobian_mult[_qp] = dstress_drot * drot_dFhatinv;
387 
388  const RankFourTensor dstrain_increment_dCtilde =
389  -0.5 * II_ijkl + 0.25 * (I.mixedProductIkJl(Ctilde) + Ctilde.mixedProductIkJl(I));
391  rot_rank_four * _Jacobian_mult[_qp] * dstrain_increment_dCtilde * dCtilde_dFhatinv;
392  _finite_deform_Jacobian_mult[_qp] += Fhat.mixedProductJkIl(_stress[_qp]);
393 
394  const RankFourTensor dFhat_dFhatinv = -Fhat.mixedProductIkJl(Fhat.transpose());
395  const RankTwoTensor dJ_dFhatinv = dFhat_dFhatinv.innerProductTranspose(Fhat.ddet());
396 
397  // Component from Jacobian derivative
398  _finite_deform_Jacobian_mult[_qp] += _stress[_qp].outerProduct(dJ_dFhatinv);
399 
400  // Derivative of Fhatinv w.r.t. undisplaced coordinates
401  const RankTwoTensor Finv = (*_deformation_gradient)[_qp].inverse();
402  const RankFourTensor dFhatinv_dGradu = -Fhatinv.mixedProductIkJl(Finv.transpose());
403  _finite_deform_Jacobian_mult[_qp] = _finite_deform_Jacobian_mult[_qp] * dFhatinv_dGradu;
404 }

Referenced by StressDivergenceTensors::computeJacobian(), and StressDivergenceTensors::computeOffDiagJacobian().

◆ computeJacobian()

void StressDivergenceTensors::computeJacobian ( )
overridevirtualinherited

Reimplemented from ALEKernel.

Definition at line 168 of file StressDivergenceTensors.C.

169 {
171  {
174  }
175 
177  {
178  _finite_deform_Jacobian_mult.resize(_qrule->n_points());
179 
180  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
182 
184  }
185  else
186  Kernel::computeJacobian();
187 }

◆ computeOffDiagJacobian()

void StressDivergenceTensors::computeOffDiagJacobian ( MooseVariableFEBase &  jvar)
overridevirtualinherited

Reimplemented from ALEKernel.

Definition at line 190 of file StressDivergenceTensors.C.

191 {
193  {
196  }
197 
199  {
200  _finite_deform_Jacobian_mult.resize(_qrule->n_points());
201 
202  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
204 
206  }
207  else
208  Kernel::computeOffDiagJacobian(jvar);
209 }

◆ computeQpJacobian()

Real StressDivergenceRZTensors::computeQpJacobian ( )
overrideprotectedvirtual

Reimplemented from StressDivergenceTensors.

Definition at line 81 of file StressDivergenceRZTensors.C.

82 {
84 }

◆ computeQpOffDiagJacobian()

Real StressDivergenceRZTensors::computeQpOffDiagJacobian ( unsigned int  jvar)
overrideprotectedvirtual

Reimplemented from StressDivergenceTensors.

Definition at line 87 of file StressDivergenceRZTensors.C.

88 {
89  for (unsigned int i = 0; i < _ndisp; ++i)
90  {
91  if (jvar == _disp_var[i])
92  return calculateJacobian(_component, i);
93  }
94 
95  if (_temp_coupled && jvar == _temp_var)
96  {
97  Real jac = 0.0;
98  if (_component == 0)
99  {
100  for (unsigned k = 0; k < LIBMESH_DIM; ++k)
101  for (unsigned l = 0; l < LIBMESH_DIM; ++l)
102  jac -= (_grad_test[_i][_qp](0) * _Jacobian_mult[_qp](0, 0, k, l) +
103  _test[_i][_qp] / _q_point[_qp](0) * _Jacobian_mult[_qp](2, 2, k, l) +
104  _grad_test[_i][_qp](1) * _Jacobian_mult[_qp](0, 1, k, l)) *
105  (*_deigenstrain_dT)[_qp](k, l);
106  return jac * _phi[_j][_qp];
107  }
108  else if (_component == 1)
109  {
110  for (unsigned k = 0; k < LIBMESH_DIM; ++k)
111  for (unsigned l = 0; l < LIBMESH_DIM; ++l)
112  jac -= (_grad_test[_i][_qp](1) * _Jacobian_mult[_qp](1, 1, k, l) +
113  _grad_test[_i][_qp](0) * _Jacobian_mult[_qp](1, 0, k, l)) *
114  (*_deigenstrain_dT)[_qp](k, l);
115  return jac * _phi[_j][_qp];
116  }
117  }
118 
119  return 0.0;
120 }

◆ computeQpResidual()

Real StressDivergenceRZTensors::computeQpResidual ( )
overrideprotectedvirtual

Reimplemented from StressDivergenceTensors.

Definition at line 51 of file StressDivergenceRZTensors.C.

52 {
53  Real div = 0.0;
54  if (_component == 0)
55  {
56  div = _grad_test[_i][_qp](0) * _stress[_qp](0, 0) +
57  +(_test[_i][_qp] / _q_point[_qp](0)) * _stress[_qp](2, 2) +
58  +_grad_test[_i][_qp](1) * _stress[_qp](0, 1); // stress_{rz}
59 
60  // volumetric locking correction
62  div += (_avg_grad_test[_i][0] - _grad_test[_i][_qp](0) - _test[_i][_qp] / _q_point[_qp](0)) *
63  (_stress[_qp].trace()) / 3.0;
64  }
65  else if (_component == 1)
66  {
67  div = _grad_test[_i][_qp](1) * _stress[_qp](1, 1) +
68  +_grad_test[_i][_qp](0) * _stress[_qp](1, 0); // stress_{zr}
69 
70  // volumetric locking correction
72  div += (_avg_grad_test[_i][1] - _grad_test[_i][_qp](1)) * (_stress[_qp].trace()) / 3.0;
73  }
74  else
75  mooseError("Invalid component for this AxisymmetricRZ problem.");
76 
77  return div;
78 }

◆ computeResidual()

void StressDivergenceTensors::computeResidual ( )
overrideprotectedvirtualinherited

Definition at line 126 of file StressDivergenceTensors.C.

127 {
128  prepareVectorTag(_assembly, _var.number());
129 
132 
133  precalculateResidual();
134  for (_i = 0; _i < _test.size(); ++_i)
135  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
136  _local_re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual();
137 
138  accumulateTaggedLocalResidual();
139 
140  if (_has_save_in)
141  {
142  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
143  for (const auto & var : _save_in)
144  var->sys().solution().add_vector(_local_re, var->dofIndices());
145  }
146 }

◆ initialSetup()

void StressDivergenceRZTensors::initialSetup ( )
overrideprotectedvirtual

Reimplemented from StressDivergenceTensors.

Definition at line 40 of file StressDivergenceRZTensors.C.

41 {
42  if (getBlockCoordSystem() != Moose::COORD_RZ)
43  mooseError("The coordinate system in the Problem block must be set to RZ for axisymmetric "
44  "geometries.");
45 
46  if (getBlockCoordSystem() == Moose::COORD_RZ && _fe_problem.getAxisymmetricRadialCoord() != 0)
47  mooseError("rz_coord_axis=Y is the only supported option for StressDivergenceRZTensors");
48 }

◆ validParams()

InputParameters StressDivergenceRZTensors::validParams ( )
static

Definition at line 20 of file StressDivergenceRZTensors.C.

21 {
22  InputParameters params = StressDivergenceTensors::validParams();
23  params.addClassDescription(
24  "Calculate stress divergence for an axisymmetric problem in cylindrical coordinates.");
25  params.addRequiredRangeCheckedParam<unsigned int>(
26  "component",
27  "component < 2",
28  "An integer corresponding to the direction the variable this kernel acts in. (0 "
29  "refers to the radial and 1 to the axial displacement.)");
30  params.set<bool>("use_displaced_mesh") = true;
31  return params;
32 }

Member Data Documentation

◆ _assembly_undisplaced

Assembly& ALEKernel::_assembly_undisplaced
protectedinherited

undisplaced problem

Definition at line 34 of file ALEKernel.h.

◆ _avg_grad_phi

std::vector<std::vector<Real> > StressDivergenceTensors::_avg_grad_phi
protectedinherited

Gradient of phi function averaged over the element. Used in volumetric locking correction calculation.

Definition at line 86 of file StressDivergenceTensors.h.

Referenced by calculateJacobian(), computeAverageGradientPhi(), StressDivergenceTensors::computeAverageGradientPhi(), StressDivergenceTensors::computeQpJacobian(), and StressDivergenceTensors::computeQpOffDiagJacobian().

◆ _avg_grad_test

std::vector<std::vector<Real> > StressDivergenceTensors::_avg_grad_test
protectedinherited

◆ _base_name

const std::string StressDivergenceTensors::_base_name
protectedinherited

◆ _component

const unsigned int StressDivergenceTensors::_component
protectedinherited

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>* StressDivergenceTensors::_deformation_gradient
protectedinherited

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>* StressDivergenceTensors::_deformation_gradient_old
protectedinherited

◆ _deigenstrain_dT

const MaterialProperty<RankTwoTensor>* const StressDivergenceTensors::_deigenstrain_dT
protectedinherited

d(strain)/d(temperature), if computed by ComputeThermalExpansionEigenstrain

Definition at line 72 of file StressDivergenceTensors.h.

◆ _disp_var

std::vector<unsigned int> StressDivergenceTensors::_disp_var
protectedinherited

◆ _finite_deform_Jacobian_mult

std::vector<RankFourTensor> StressDivergenceTensors::_finite_deform_Jacobian_mult
protectedinherited

◆ _grad_phi_undisplaced

const VariablePhiGradient& ALEKernel::_grad_phi_undisplaced
protectedinherited

Shape and test functions on the undisplaced mesh.

Definition at line 40 of file ALEKernel.h.

Referenced by StressDivergenceTensors::computeQpJacobian(), and StressDivergenceTensors::computeQpOffDiagJacobian().

◆ _grad_test_undisplaced

const VariableTestGradient& ALEKernel::_grad_test_undisplaced
protectedinherited

Definition at line 41 of file ALEKernel.h.

◆ _Jacobian_mult

const MaterialProperty<RankFourTensor>& StressDivergenceTensors::_Jacobian_mult
protectedinherited

◆ _ndisp

unsigned int StressDivergenceTensors::_ndisp
protectedinherited

◆ _out_of_plane_direction

const unsigned int StressDivergenceTensors::_out_of_plane_direction
protectedinherited

◆ _out_of_plane_strain

const VariableValue* StressDivergenceTensors::_out_of_plane_strain
protectedinherited

◆ _out_of_plane_strain_coupled

const bool StressDivergenceTensors::_out_of_plane_strain_coupled
protectedinherited

◆ _out_of_plane_strain_var

const unsigned int StressDivergenceTensors::_out_of_plane_strain_var
protectedinherited

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>* StressDivergenceTensors::_rotation_increment
protectedinherited

◆ _stress

const MaterialProperty<RankTwoTensor>& StressDivergenceTensors::_stress
protectedinherited

◆ _temp_coupled

const bool StressDivergenceTensors::_temp_coupled
protectedinherited

◆ _temp_var

const unsigned int StressDivergenceTensors::_temp_var
protectedinherited

◆ _use_displaced_mesh

const bool StressDivergenceTensors::_use_displaced_mesh
protectedinherited

Whether this object is acting on the displaced mesh.

Definition at line 80 of file StressDivergenceTensors.h.

Referenced by StressDivergenceTensors::computeQpJacobian(), and StressDivergenceTensors::computeQpResidual().

◆ _use_finite_deform_jacobian

bool StressDivergenceTensors::_use_finite_deform_jacobian
protectedinherited

◆ _var_undisplaced

MooseVariable& ALEKernel::_var_undisplaced
protectedinherited

Reference to this Kernel's undisplaced MooseVariable object.

Definition at line 37 of file ALEKernel.h.

◆ _volumetric_locking_correction

bool StressDivergenceTensors::_volumetric_locking_correction
protectedinherited

The documentation for this class was generated from the following files:
ElasticityTensorTools::elasticJacobian
Real elasticJacobian(const RankFourTensor &r4t, unsigned int i, unsigned int k, const RealGradient &grad_test, const RealGradient &grad_phi)
This is used for the standard kernel stress_ij*d(test)/dx_j, when varied wrt u_k Jacobian entry: d(st...
Definition: ElasticityTensorTools.C:21
ALEKernel::computeOffDiagJacobian
virtual void computeOffDiagJacobian(MooseVariableFEBase &jvar) override
Definition: ALEKernel.C:43
StressDivergenceTensors::StressDivergenceTensors
StressDivergenceTensors(const InputParameters &parameters)
Definition: StressDivergenceTensors.C:62
StressDivergenceTensors::_component
const unsigned int _component
Definition: StressDivergenceTensors.h:62
StressDivergenceTensors::computeAverageGradientPhi
virtual void computeAverageGradientPhi()
Definition: StressDivergenceTensors.C:423
StressDivergenceTensors::_disp_var
std::vector< unsigned int > _disp_var
Definition: StressDivergenceTensors.h:66
libMesh::RealGradient
VectorValue< Real > RealGradient
Definition: GrainForceAndTorqueInterface.h:17
StressDivergenceTensors::_rotation_increment
const MaterialProperty< RankTwoTensor > * _rotation_increment
Definition: StressDivergenceTensors.h:59
StressDivergenceTensors::_temp_coupled
const bool _temp_coupled
Definition: StressDivergenceTensors.h:68
StressDivergenceTensors::_avg_grad_phi
std::vector< std::vector< Real > > _avg_grad_phi
Gradient of phi function averaged over the element. Used in volumetric locking correction calculation...
Definition: StressDivergenceTensors.h:86
ALEKernel::computeJacobian
virtual void computeJacobian() override
Definition: ALEKernel.C:36
StressDivergenceRZTensors::calculateJacobian
Real calculateJacobian(unsigned int ivar, unsigned int jvar)
Definition: StressDivergenceRZTensors.C:123
StressDivergenceTensors::computeFiniteDeformJacobian
virtual void computeFiniteDeformJacobian()
Definition: StressDivergenceTensors.C:350
StressDivergenceTensors::_volumetric_locking_correction
bool _volumetric_locking_correction
Flag for volumetric locking correction.
Definition: StressDivergenceTensors.h:89
StressDivergenceTensors::computeQpResidual
virtual Real computeQpResidual() override
Definition: StressDivergenceTensors.C:149
StressDivergenceTensors::validParams
static InputParameters validParams()
Definition: StressDivergenceTensors.C:26
MaterialTensorCalculatorTools::component
Real component(const SymmTensor &symm_tensor, unsigned int index)
Definition: MaterialTensorCalculatorTools.C:16
StressDivergenceTensors::_Jacobian_mult
const MaterialProperty< RankFourTensor > & _Jacobian_mult
Definition: StressDivergenceTensors.h:54
StressDivergenceTensors::_finite_deform_Jacobian_mult
std::vector< RankFourTensor > _finite_deform_Jacobian_mult
Definition: StressDivergenceTensors.h:56
StressDivergenceTensors::_stress
const MaterialProperty< RankTwoTensor > & _stress
Definition: StressDivergenceTensors.h:53
RankFourTensorTempl
Definition: ACGrGrElasticDrivingForce.h:20
StressDivergenceTensors::_avg_grad_test
std::vector< std::vector< Real > > _avg_grad_test
Gradient of test function averaged over the element. Used in volumetric locking correction calculatio...
Definition: StressDivergenceTensors.h:83
StressDivergenceTensors::computeAverageGradientTest
virtual void computeAverageGradientTest()
Definition: StressDivergenceTensors.C:407
StressDivergenceTensors::_ndisp
unsigned int _ndisp
Coupled displacement variables.
Definition: StressDivergenceTensors.h:65
StressDivergenceTensors::_use_finite_deform_jacobian
bool _use_finite_deform_jacobian
Definition: StressDivergenceTensors.h:51
RankTwoTensorTempl< Real >
StressDivergenceTensors::_temp_var
const unsigned int _temp_var
Definition: StressDivergenceTensors.h:69