www.mooseframework.org
ComputeAxisymmetric1DFiniteStrain.C
Go to the documentation of this file.
1 /****************************************************************/
2 /* MOOSE - Multiphysics Object Oriented Simulation Environment */
3 /* */
4 /* All contents are licensed under LGPL V2.1 */
5 /* See LICENSE for full restrictions */
6 /****************************************************************/
7 
9 
10 template <>
11 InputParameters
13 {
14  InputParameters params = validParams<Compute1DFiniteStrain>();
15  params.addClassDescription("Compute a strain increment and rotation increment for finite strains "
16  "in an axisymmetric 1D problem");
17  params.addParam<UserObjectName>("subblock_index_provider",
18  "SubblockIndexProvider user object name");
19  params.addCoupledVar("scalar_out_of_plane_strain", "Scalar variable for axisymmetric 1D problem");
20  params.addCoupledVar("out_of_plane_strain", "Nonlinear variable for axisymmetric 1D problem");
21 
22  return params;
23 }
24 
26  const InputParameters & parameters)
27  : Compute1DFiniteStrain(parameters),
28  _disp_old_0(coupledValueOld("displacements", 0)),
29  _subblock_id_provider(isParamValid("subblock_index_provider")
30  ? &getUserObject<SubblockIndexProvider>("subblock_index_provider")
31  : nullptr),
32  _has_out_of_plane_strain(isParamValid("out_of_plane_strain")),
33  _out_of_plane_strain(_has_out_of_plane_strain ? coupledValue("out_of_plane_strain") : _zero),
34  _out_of_plane_strain_old(_has_out_of_plane_strain ? coupledValueOld("out_of_plane_strain")
35  : _zero),
36  _has_scalar_out_of_plane_strain(isParamValid("scalar_out_of_plane_strain")),
37  _nscalar_strains(
38  _has_scalar_out_of_plane_strain ? coupledScalarComponents("scalar_out_of_plane_strain") : 0)
39 {
41  mooseError("Must define only one of out_of_plane_strain or scalar_out_of_plane_strain");
42 
44  mooseError("Must define either out_of_plane_strain or scalar_out_of_plane_strain");
45 
46  // in case when the provided scalar_out_of_plane_strain is not a coupled
47  // scalar variable, still set _nscalar_strains = 1 but return its default value 0
48  if (coupledScalarComponents("scalar_out_of_plane_strain") == 0)
49  _nscalar_strains = 1;
50 
52  {
55  for (unsigned int i = 0; i < _nscalar_strains; ++i)
56  {
57  _scalar_out_of_plane_strain[i] = &coupledScalarValue("scalar_out_of_plane_strain", i);
58  _scalar_out_of_plane_strain_old[i] = &coupledScalarValueOld("scalar_out_of_plane_strain", i);
59  }
60  }
61 }
62 
63 void
65 {
66  if (getBlockCoordSystem() != Moose::COORD_RZ)
67  mooseError("The coordinate system must be set to RZ for Axisymmetric geometries.");
68 }
69 
70 Real
72 {
74  return std::exp((*_scalar_out_of_plane_strain[getCurrentSubblockIndex()])[0]) - 1.0;
75  else
76  return std::exp(_out_of_plane_strain[_qp]) - 1.0;
77 }
78 
79 Real
81 {
83  return std::exp((*_scalar_out_of_plane_strain_old[getCurrentSubblockIndex()])[0]) - 1.0;
84  else
85  return std::exp(_out_of_plane_strain_old[_qp]) - 1.0;
86 }
87 
88 Real
90 {
91  if (!MooseUtils::absoluteFuzzyEqual(_q_point[_qp](0), 0.0))
92  return (*_disp[0])[_qp] / _q_point[_qp](0);
93  else
94  return 0.0;
95 }
96 
97 Real
99 {
100  if (!MooseUtils::absoluteFuzzyEqual(_q_point[_qp](0), 0.0))
101  return _disp_old_0[_qp] / _q_point[_qp](0);
102  else
103  return 0.0;
104 }
ComputeAxisymmetric1DFiniteStrain(const InputParameters &parameters)
Abstract base class for user objects that provide an index for a given element that is independent of...
Compute1DFiniteStrain defines a strain increment for finite strains in 1D problems, handling strains in other two directions.
std::vector< const VariableValue * > _disp
Real computeGradDispYY() override
Computes the current dUy/dy for axisymmetric problems.
Real computeGradDispZZ() override
Computes the current dUz/dz for axisymmetric problems, where .
const VariableValue & _disp_old_0
the old value of the first component of the displacements vector
std::vector< const VariableValue * > _scalar_out_of_plane_strain
Real computeGradDispZZOld() override
Computes the old dUz/dz for axisymmetric problems, where .
InputParameters validParams< Compute1DFiniteStrain >()
unsigned int getCurrentSubblockIndex() const
gets its subblock index for current element
std::vector< const VariableValue * > _scalar_out_of_plane_strain_old
Real computeGradDispYYOld() override
Computes the old dUy/dy for axisymmetric problems.
InputParameters validParams< ComputeAxisymmetric1DFiniteStrain >()