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

StressDivergenceTensors mostly copies from StressDivergence. More...

#include <StressDivergenceTensors.h>

Inheritance diagram for StressDivergenceTensors:
[legend]

Public Member Functions

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

Protected Member Functions

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

Protected Attributes

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 unsigned int _out_of_plane_strain_var
 
const unsigned int _out_of_plane_direction
 
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

StressDivergenceTensors mostly copies from StressDivergence.

There are small changes to use RankFourTensor and RankTwoTensors instead of SymmElasticityTensors and SymmTensors. This is done to allow for more mathematical transparancy.

Definition at line 28 of file StressDivergenceTensors.h.

Constructor & Destructor Documentation

◆ StressDivergenceTensors()

StressDivergenceTensors::StressDivergenceTensors ( const InputParameters &  parameters)

Definition at line 62 of file StressDivergenceTensors.C.

63  : ALEKernel(parameters),
64  _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
65  _use_finite_deform_jacobian(getParam<bool>("use_finite_deform_jacobian")),
66  _stress(getMaterialPropertyByName<RankTwoTensor>(_base_name + "stress")),
67  _Jacobian_mult(getMaterialPropertyByName<RankFourTensor>(_base_name + "Jacobian_mult")),
68  _component(getParam<unsigned int>("component")),
69  _ndisp(coupledComponents("displacements")),
71  _temp_coupled(isCoupled("temperature")),
72  _temp_var(_temp_coupled ? coupled("temperature") : 0),
73  _deigenstrain_dT(_temp_coupled ? &getMaterialPropertyDerivative<RankTwoTensor>(
74  getParam<std::string>("thermal_eigenstrain_name"),
75  getVar("temperature", 0)->name())
76  : nullptr),
77  _out_of_plane_strain_coupled(isCoupled("out_of_plane_strain")),
78  _out_of_plane_strain_var(_out_of_plane_strain_coupled ? coupled("out_of_plane_strain") : 0),
79  _out_of_plane_direction(getParam<MooseEnum>("out_of_plane_direction")),
80  _avg_grad_test(_test.size(), std::vector<Real>(3, 0.0)),
81  _avg_grad_phi(_phi.size(), std::vector<Real>(3, 0.0)),
82  _volumetric_locking_correction(getParam<bool>("volumetric_locking_correction"))
83 {
84  for (unsigned int i = 0; i < _ndisp; ++i)
85  _disp_var[i] = coupled("displacements", i);
86 
87  // Checking for consistency between mesh size and length of the provided displacements vector
88  if (_out_of_plane_direction != 2 && _ndisp != 3)
89  mooseError("For 2D simulations where the out-of-plane direction is x or y coordinate "
90  "directions the number of supplied displacements must be three.");
91  else if (_out_of_plane_direction == 2 && _ndisp != _mesh.dimension())
92  mooseError("The number of displacement variables supplied must match the mesh dimension");
93 
95  {
97  &getMaterialProperty<RankTwoTensor>(_base_name + "deformation_gradient");
99  &getMaterialPropertyOld<RankTwoTensor>(_base_name + "deformation_gradient");
100  _rotation_increment = &getMaterialProperty<RankTwoTensor>(_base_name + "rotation_increment");
101  }
102 
103  // Error if volumetric locking correction is turned on for 1D problems
105  mooseError("Volumetric locking correction should be set to false for 1-D problems.");
106 
107  // Generate warning when volumetric locking correction is used with second order elements
108  if (_mesh.hasSecondOrderElements() && _volumetric_locking_correction)
109  mooseWarning("Volumteric locking correction is not required for second order elements. Using "
110  "volumetric locking with second order elements could cause zigzag patterns in "
111  "stresses and strains.");
112 }
std::vector< std::vector< Real > > _avg_grad_phi
Gradient of phi function averaged over the element. Used in volumetric locking correction calculation...
bool _volumetric_locking_correction
Flag for volumetric locking correction.
unsigned int _ndisp
Coupled displacement variables.
std::vector< std::vector< Real > > _avg_grad_test
Gradient of test function averaged over the element. Used in volumetric locking correction calculatio...
const unsigned int _out_of_plane_strain_var
const MaterialProperty< RankTwoTensor > * _deformation_gradient_old
const unsigned int _out_of_plane_direction
const MaterialProperty< RankTwoTensor > * _rotation_increment
std::vector< unsigned int > _disp_var
ALEKernel(const InputParameters &parameters)
Definition: ALEKernel.C:24
const std::string name
Definition: Setup.h:22
const MaterialProperty< RankTwoTensor > * _deformation_gradient
const MaterialProperty< RankFourTensor > & _Jacobian_mult
const MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< RankTwoTensor > *const _deigenstrain_dT
d(strain)/d(temperature), if computed by ComputeThermalExpansionEigenstrain

Member Function Documentation

◆ computeAverageGradientPhi()

void StressDivergenceTensors::computeAverageGradientPhi ( )
protectedvirtual

Reimplemented in StressDivergenceRZTensors.

Definition at line 408 of file StressDivergenceTensors.C.

Referenced by computeJacobian(), and computeOffDiagJacobian().

409 {
410  // Calculate volume average derivatives for phi
411  _avg_grad_phi.resize(_phi.size());
412  for (_i = 0; _i < _phi.size(); ++_i)
413  {
414  _avg_grad_phi[_i].resize(3);
415  for (unsigned int component = 0; component < _mesh.dimension(); ++component)
416  {
417  _avg_grad_phi[_i][component] = 0.0;
418  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
419  _avg_grad_phi[_i][component] += _grad_phi[_i][_qp](component) * _JxW[_qp] * _coord[_qp];
420 
421  _avg_grad_phi[_i][component] /= _current_elem_volume;
422  }
423  }
424 }
std::vector< std::vector< Real > > _avg_grad_phi
Gradient of phi function averaged over the element. Used in volumetric locking correction calculation...
Real component(const SymmTensor &symm_tensor, unsigned int index)

◆ computeAverageGradientTest()

void StressDivergenceTensors::computeAverageGradientTest ( )
protectedvirtual

Reimplemented in StressDivergenceRZTensors.

Definition at line 392 of file StressDivergenceTensors.C.

Referenced by computeJacobian(), computeOffDiagJacobian(), and computeResidual().

393 {
394  // Calculate volume averaged value of shape function derivative
395  _avg_grad_test.resize(_test.size());
396  for (_i = 0; _i < _test.size(); ++_i)
397  {
398  _avg_grad_test[_i].resize(3);
399  _avg_grad_test[_i][_component] = 0.0;
400  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
401  _avg_grad_test[_i][_component] += _grad_test[_i][_qp](_component) * _JxW[_qp] * _coord[_qp];
402 
403  _avg_grad_test[_i][_component] /= _current_elem_volume;
404  }
405 }
std::vector< std::vector< Real > > _avg_grad_test
Gradient of test function averaged over the element. Used in volumetric locking correction calculatio...

◆ computeFiniteDeformJacobian()

void StressDivergenceTensors::computeFiniteDeformJacobian ( )
protectedvirtual

Definition at line 335 of file StressDivergenceTensors.C.

Referenced by computeJacobian(), and computeOffDiagJacobian().

336 {
337  const RankTwoTensor I(RankTwoTensor::initIdentity);
338  const RankFourTensor II_ijkl = I.mixedProductIkJl(I);
339 
340  // Bring back to unrotated config
341  const RankTwoTensor unrotated_stress =
342  (*_rotation_increment)[_qp].transpose() * _stress[_qp] * (*_rotation_increment)[_qp];
343 
344  // Incremental deformation gradient Fhat
345  const RankTwoTensor Fhat =
346  (*_deformation_gradient)[_qp] * (*_deformation_gradient_old)[_qp].inverse();
347  const RankTwoTensor Fhatinv = Fhat.inverse();
348 
349  const RankTwoTensor rot_times_stress = (*_rotation_increment)[_qp] * unrotated_stress;
350  const RankFourTensor dstress_drot =
351  I.mixedProductIkJl(rot_times_stress) + I.mixedProductJkIl(rot_times_stress);
352  const RankFourTensor rot_rank_four =
353  (*_rotation_increment)[_qp].mixedProductIkJl((*_rotation_increment)[_qp]);
354  const RankFourTensor drot_dUhatinv = Fhat.mixedProductIkJl(I);
355 
356  const RankTwoTensor A = I - Fhatinv;
357 
358  // Ctilde = Chat^-1 - I
359  const RankTwoTensor Ctilde = A * A.transpose() - A - A.transpose();
360  const RankFourTensor dCtilde_dFhatinv =
361  -I.mixedProductIkJl(A) - I.mixedProductJkIl(A) + II_ijkl + I.mixedProductJkIl(I);
362 
363  // Second order approximation of Uhat - consistent with strain increment definition
364  // const RankTwoTensor Uhat = I - 0.5 * Ctilde - 3.0/8.0 * Ctilde * Ctilde;
365 
366  RankFourTensor dUhatinv_dCtilde =
367  0.5 * II_ijkl - 1.0 / 8.0 * (I.mixedProductIkJl(Ctilde) + Ctilde.mixedProductIkJl(I));
368  RankFourTensor drot_dFhatinv = drot_dUhatinv * dUhatinv_dCtilde * dCtilde_dFhatinv;
369 
370  drot_dFhatinv -= Fhat.mixedProductIkJl((*_rotation_increment)[_qp].transpose());
371  _finite_deform_Jacobian_mult[_qp] = dstress_drot * drot_dFhatinv;
372 
373  const RankFourTensor dstrain_increment_dCtilde =
374  -0.5 * II_ijkl + 0.25 * (I.mixedProductIkJl(Ctilde) + Ctilde.mixedProductIkJl(I));
376  rot_rank_four * _Jacobian_mult[_qp] * dstrain_increment_dCtilde * dCtilde_dFhatinv;
377  _finite_deform_Jacobian_mult[_qp] += Fhat.mixedProductJkIl(_stress[_qp]);
378 
379  const RankFourTensor dFhat_dFhatinv = -Fhat.mixedProductIkJl(Fhat.transpose());
380  const RankTwoTensor dJ_dFhatinv = dFhat_dFhatinv.innerProductTranspose(Fhat.ddet());
381 
382  // Component from Jacobian derivative
383  _finite_deform_Jacobian_mult[_qp] += _stress[_qp].outerProduct(dJ_dFhatinv);
384 
385  // Derivative of Fhatinv w.r.t. undisplaced coordinates
386  const RankTwoTensor Finv = (*_deformation_gradient)[_qp].inverse();
387  const RankFourTensor dFhatinv_dGradu = -Fhatinv.mixedProductIkJl(Finv.transpose());
388  _finite_deform_Jacobian_mult[_qp] = _finite_deform_Jacobian_mult[_qp] * dFhatinv_dGradu;
389 }
const MaterialProperty< RankTwoTensor > * _rotation_increment
std::vector< RankFourTensor > _finite_deform_Jacobian_mult
const MaterialProperty< RankFourTensor > & _Jacobian_mult
const MaterialProperty< RankTwoTensor > & _stress

◆ computeJacobian()

void StressDivergenceTensors::computeJacobian ( )
overridevirtual

Reimplemented from ALEKernel.

Definition at line 160 of file StressDivergenceTensors.C.

161 {
163  {
166  }
167 
169  {
170  _finite_deform_Jacobian_mult.resize(_qrule->n_points());
171 
172  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
174 
176  }
177  else
178  Kernel::computeJacobian();
179 }
bool _volumetric_locking_correction
Flag for volumetric locking correction.
std::vector< RankFourTensor > _finite_deform_Jacobian_mult
virtual void computeFiniteDeformJacobian()
virtual void computeJacobian() override
Definition: ALEKernel.C:35

◆ computeOffDiagJacobian()

void StressDivergenceTensors::computeOffDiagJacobian ( MooseVariableFEBase &  jvar)
overridevirtual

Reimplemented from ALEKernel.

Definition at line 182 of file StressDivergenceTensors.C.

183 {
185  {
188  }
189 
191  {
192  _finite_deform_Jacobian_mult.resize(_qrule->n_points());
193 
194  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
196 
198  }
199  else
200  Kernel::computeOffDiagJacobian(jvar);
201 }
bool _volumetric_locking_correction
Flag for volumetric locking correction.
virtual void computeOffDiagJacobian(MooseVariableFEBase &jvar) override
Definition: ALEKernel.C:42
std::vector< RankFourTensor > _finite_deform_Jacobian_mult
virtual void computeFiniteDeformJacobian()

◆ computeQpJacobian()

Real StressDivergenceTensors::computeQpJacobian ( )
overrideprotectedvirtual

Reimplemented in StressDivergenceRSphericalTensors, StressDivergenceRZTensors, and DynamicStressDivergenceTensors.

Definition at line 204 of file StressDivergenceTensors.C.

Referenced by DynamicStressDivergenceTensors::computeQpJacobian().

205 {
208  _component,
209  _component,
210  _grad_test[_i][_qp],
211  _grad_phi_undisplaced[_j][_qp]);
212 
213  Real sum_C3x3 = _Jacobian_mult[_qp].sum3x3();
214  RealGradient sum_C3x1 = _Jacobian_mult[_qp].sum3x1();
215 
216  Real jacobian = 0.0;
217  // B^T_i * C * B_j
219  _Jacobian_mult[_qp], _component, _component, _grad_test[_i][_qp], _grad_phi[_j][_qp]);
220 
222  {
223  // jacobian = Bbar^T_i * C * Bbar_j where Bbar = B + Bvol
224  // jacobian = B^T_i * C * B_j + Bvol^T_i * C * Bvol_j + Bvol^T_i * C * B_j + B^T_i * C * Bvol_j
225 
226  // Bvol^T_i * C * Bvol_j
227  jacobian += sum_C3x3 * (_avg_grad_test[_i][_component] - _grad_test[_i][_qp](_component)) *
228  (_avg_grad_phi[_j][_component] - _grad_phi[_j][_qp](_component)) / 9.0;
229 
230  // B^T_i * C * Bvol_j
231  jacobian += sum_C3x1(_component) * _grad_test[_i][_qp](_component) *
232  (_avg_grad_phi[_j][_component] - _grad_phi[_j][_qp](_component)) / 3.0;
233 
234  // Bvol^T_i * C * B_j
235  RankTwoTensor phi;
236  if (_component == 0)
237  {
238  phi(0, 0) = _grad_phi[_j][_qp](0);
239  phi(0, 1) = phi(1, 0) = _grad_phi[_j][_qp](1);
240  phi(0, 2) = phi(2, 0) = _grad_phi[_j][_qp](2);
241  }
242  else if (_component == 1)
243  {
244  phi(1, 1) = _grad_phi[_j][_qp](1);
245  phi(0, 1) = phi(1, 0) = _grad_phi[_j][_qp](0);
246  phi(1, 2) = phi(2, 1) = _grad_phi[_j][_qp](2);
247  }
248  else if (_component == 2)
249  {
250  phi(2, 2) = _grad_phi[_j][_qp](2);
251  phi(0, 2) = phi(2, 0) = _grad_phi[_j][_qp](0);
252  phi(1, 2) = phi(2, 1) = _grad_phi[_j][_qp](1);
253  }
254 
255  jacobian += (_Jacobian_mult[_qp] * phi).trace() *
256  (_avg_grad_test[_i][_component] - _grad_test[_i][_qp](_component)) / 3.0;
257  }
258  return jacobian;
259 }
std::vector< std::vector< Real > > _avg_grad_phi
Gradient of phi function averaged over the element. Used in volumetric locking correction calculation...
VectorValue< Real > RealGradient
bool _volumetric_locking_correction
Flag for volumetric locking correction.
std::vector< std::vector< Real > > _avg_grad_test
Gradient of test function averaged over the element. Used in volumetric locking correction calculatio...
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...
std::vector< RankFourTensor > _finite_deform_Jacobian_mult
const MaterialProperty< RankFourTensor > & _Jacobian_mult
const VariablePhiGradient & _grad_phi_undisplaced
Shape and test functions on the undisplaced mesh.
Definition: ALEKernel.h:39

◆ computeQpOffDiagJacobian()

Real StressDivergenceTensors::computeQpOffDiagJacobian ( unsigned int  jvar)
overrideprotectedvirtual

Reimplemented in StressDivergenceRSphericalTensors, StressDivergenceRZTensors, CosseratStressDivergenceTensors, and DynamicStressDivergenceTensors.

Definition at line 262 of file StressDivergenceTensors.C.

Referenced by CosseratStressDivergenceTensors::computeQpOffDiagJacobian(), and DynamicStressDivergenceTensors::computeQpOffDiagJacobian().

263 {
264  // off-diagonal Jacobian with respect to a coupled displacement component
265  for (unsigned int coupled_component = 0; coupled_component < _ndisp; ++coupled_component)
266  if (jvar == _disp_var[coupled_component])
267  {
268  if (_out_of_plane_direction != 2)
269  {
270  if (coupled_component == _out_of_plane_direction)
271  continue;
272  }
273 
276  _component,
277  coupled_component,
278  _grad_test[_i][_qp],
279  _grad_phi_undisplaced[_j][_qp]);
280 
281  const Real sum_C3x3 = _Jacobian_mult[_qp].sum3x3();
282  const RealGradient sum_C3x1 = _Jacobian_mult[_qp].sum3x1();
283  Real jacobian = 0.0;
284 
285  // B^T_i * C * B_j
287  _component,
288  coupled_component,
289  _grad_test[_i][_qp],
290  _grad_phi[_j][_qp]);
291 
293  {
294  // jacobian = Bbar^T_i * C * Bbar_j where Bbar = B + Bvol
295  // jacobian = B^T_i * C * B_j + Bvol^T_i * C * Bvol_j + Bvol^T_i * C * B_j + B^T_i * C *
296  // Bvol_j
297 
298  // Bvol^T_i * C * Bvol_j
299  jacobian += sum_C3x3 * (_avg_grad_test[_i][_component] - _grad_test[_i][_qp](_component)) *
300  (_avg_grad_phi[_j][coupled_component] - _grad_phi[_j][_qp](coupled_component)) /
301  9.0;
302 
303  // B^T_i * C * Bvol_j
304  jacobian += sum_C3x1(_component) * _grad_test[_i][_qp](_component) *
305  (_avg_grad_phi[_j][coupled_component] - _grad_phi[_j][_qp](coupled_component)) /
306  3.0;
307 
308  // Bvol^T_i * C * B_i
309  RankTwoTensor phi;
310  for (unsigned int i = 0; i < 3; ++i)
311  phi(coupled_component, i) = _grad_phi[_j][_qp](i);
312 
313  jacobian += (_Jacobian_mult[_qp] * phi).trace() *
314  (_avg_grad_test[_i][_component] - _grad_test[_i][_qp](_component)) / 3.0;
315  }
316 
317  return jacobian;
318  }
319 
320  // off-diagonal Jacobian with respect to a coupled out_of_plane_strain variable
322  return _Jacobian_mult[_qp](
324  _grad_test[_i][_qp](_component) * _phi[_j][_qp];
325 
326  // off-diagonal Jacobian with respect to a coupled temperature variable
327  if (_temp_coupled && jvar == _temp_var)
328  return -((_Jacobian_mult[_qp] * (*_deigenstrain_dT)[_qp]) *
329  _grad_test[_i][_qp])(_component)*_phi[_j][_qp];
330 
331  return 0.0;
332 }
std::vector< std::vector< Real > > _avg_grad_phi
Gradient of phi function averaged over the element. Used in volumetric locking correction calculation...
VectorValue< Real > RealGradient
bool _volumetric_locking_correction
Flag for volumetric locking correction.
unsigned int _ndisp
Coupled displacement variables.
std::vector< std::vector< Real > > _avg_grad_test
Gradient of test function averaged over the element. Used in volumetric locking correction calculatio...
const unsigned int _out_of_plane_strain_var
const unsigned int _out_of_plane_direction
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...
std::vector< unsigned int > _disp_var
std::vector< RankFourTensor > _finite_deform_Jacobian_mult
const MaterialProperty< RankFourTensor > & _Jacobian_mult
const VariablePhiGradient & _grad_phi_undisplaced
Shape and test functions on the undisplaced mesh.
Definition: ALEKernel.h:39

◆ computeQpResidual()

Real StressDivergenceTensors::computeQpResidual ( )
overrideprotectedvirtual

Reimplemented in StressDivergenceRSphericalTensors, StressDivergenceRZTensors, and DynamicStressDivergenceTensors.

Definition at line 148 of file StressDivergenceTensors.C.

Referenced by computeResidual().

149 {
150  Real residual = _stress[_qp].row(_component) * _grad_test[_i][_qp];
151  // volumetric locking correction
153  residual += _stress[_qp].trace() / 3.0 *
154  (_avg_grad_test[_i][_component] - _grad_test[_i][_qp](_component));
155 
156  return residual;
157 }
bool _volumetric_locking_correction
Flag for volumetric locking correction.
std::vector< std::vector< Real > > _avg_grad_test
Gradient of test function averaged over the element. Used in volumetric locking correction calculatio...
const MaterialProperty< RankTwoTensor > & _stress

◆ computeResidual()

void StressDivergenceTensors::computeResidual ( )
overrideprotectedvirtual

Definition at line 123 of file StressDivergenceTensors.C.

124 {
125  DenseVector<Number> & re = _assembly.residualBlock(_var.number());
126  _local_re.resize(re.size());
127  _local_re.zero();
128 
131 
132  precalculateResidual();
133  for (_i = 0; _i < _test.size(); ++_i)
134  for (_qp = 0; _qp < _qrule->n_points(); ++_qp)
135  _local_re(_i) += _JxW[_qp] * _coord[_qp] * computeQpResidual();
136 
137  re += _local_re;
138 
139  if (_has_save_in)
140  {
141  Threads::spin_mutex::scoped_lock lock(Threads::spin_mtx);
142  for (const auto & var : _save_in)
143  var->sys().solution().add_vector(_local_re, var->dofIndices());
144  }
145 }
bool _volumetric_locking_correction
Flag for volumetric locking correction.
virtual Real computeQpResidual() override

◆ initialSetup()

void StressDivergenceTensors::initialSetup ( )
overrideprotectedvirtual

Reimplemented in StressDivergenceRSphericalTensors, and StressDivergenceRZTensors.

Definition at line 115 of file StressDivergenceTensors.C.

116 {
117  if (getBlockCoordSystem() != Moose::COORD_XYZ)
118  mooseError(
119  "The coordinate system in the Problem block must be set to XYZ for cartesian geometries.");
120 }

Member Data Documentation

◆ _assembly_undisplaced

Assembly& ALEKernel::_assembly_undisplaced
protectedinherited

undisplaced problem

Definition at line 33 of file ALEKernel.h.

◆ _avg_grad_phi

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

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

Definition at line 81 of file StressDivergenceTensors.h.

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

◆ _avg_grad_test

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

◆ _base_name

std::string StressDivergenceTensors::_base_name
protected

Definition at line 49 of file StressDivergenceTensors.h.

Referenced by StressDivergenceTensors().

◆ _component

const unsigned int StressDivergenceTensors::_component
protected

◆ _deformation_gradient

const MaterialProperty<RankTwoTensor>* StressDivergenceTensors::_deformation_gradient
protected

Definition at line 56 of file StressDivergenceTensors.h.

Referenced by StressDivergenceTensors().

◆ _deformation_gradient_old

const MaterialProperty<RankTwoTensor>* StressDivergenceTensors::_deformation_gradient_old
protected

Definition at line 57 of file StressDivergenceTensors.h.

Referenced by StressDivergenceTensors().

◆ _deigenstrain_dT

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

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

Definition at line 71 of file StressDivergenceTensors.h.

◆ _disp_var

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

◆ _finite_deform_Jacobian_mult

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

◆ _grad_phi_undisplaced

const VariablePhiGradient& ALEKernel::_grad_phi_undisplaced
protectedinherited

Shape and test functions on the undisplaced mesh.

Definition at line 39 of file ALEKernel.h.

Referenced by computeQpJacobian(), and computeQpOffDiagJacobian().

◆ _grad_test_undisplaced

const VariableTestGradient& ALEKernel::_grad_test_undisplaced
protectedinherited

Definition at line 40 of file ALEKernel.h.

◆ _Jacobian_mult

const MaterialProperty<RankFourTensor>& StressDivergenceTensors::_Jacobian_mult
protected

◆ _ndisp

unsigned int StressDivergenceTensors::_ndisp
protected

◆ _out_of_plane_direction

const unsigned int StressDivergenceTensors::_out_of_plane_direction
protected

Definition at line 75 of file StressDivergenceTensors.h.

Referenced by computeQpOffDiagJacobian(), and StressDivergenceTensors().

◆ _out_of_plane_strain_coupled

const bool StressDivergenceTensors::_out_of_plane_strain_coupled
protected

Definition at line 73 of file StressDivergenceTensors.h.

Referenced by computeQpOffDiagJacobian().

◆ _out_of_plane_strain_var

const unsigned int StressDivergenceTensors::_out_of_plane_strain_var
protected

Definition at line 74 of file StressDivergenceTensors.h.

Referenced by computeQpOffDiagJacobian().

◆ _rotation_increment

const MaterialProperty<RankTwoTensor>* StressDivergenceTensors::_rotation_increment
protected

◆ _stress

const MaterialProperty<RankTwoTensor>& StressDivergenceTensors::_stress
protected

◆ _temp_coupled

const bool StressDivergenceTensors::_temp_coupled
protected

◆ _temp_var

const unsigned int StressDivergenceTensors::_temp_var
protected

◆ _use_finite_deform_jacobian

bool StressDivergenceTensors::_use_finite_deform_jacobian
protected

◆ _var_undisplaced

MooseVariable& ALEKernel::_var_undisplaced
protectedinherited

Reference to this Kernel's undisplaced MooseVariable object.

Definition at line 36 of file ALEKernel.h.

◆ _volumetric_locking_correction

bool StressDivergenceTensors::_volumetric_locking_correction
protected

The documentation for this class was generated from the following files: