Loading [MathJax]/extensions/tex2jax.js
https://mooseframework.inl.gov
All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends
Public Types | Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
HomogenizedTotalLagrangianStressDivergenceS Class Reference

Total Lagrangian formulation with all homogenization terms (one disp_xyz field and macro_gradient scalar) More...

#include <HomogenizedTotalLagrangianStressDivergenceS.h>

Inheritance diagram for HomogenizedTotalLagrangianStressDivergenceS:
[legend]

Public Types

typedef std::vector< intJvarMap
 

Public Member Functions

 HomogenizedTotalLagrangianStressDivergenceS (const InputParameters &parameters)
 
virtual void initialSetup () override
 
template<>
void initialSetup ()
 
template<>
void initialSetup ()
 
template<>
void initialSetup ()
 
virtual void computeOffDiagJacobian (unsigned int jvar) override
 
unsigned int mapJvarToCvar (unsigned int jvar)
 
int mapJvarToCvar (unsigned int jvar, const JvarMap &jvar_map)
 
bool mapJvarToCvar (unsigned int jvar, unsigned int &cvar)
 
const JvarMapgetJvarMap ()
 
const JvarMapgetParameterJvarMap (std::string parameter_name)
 

Static Public Member Functions

static InputParameters validParams ()
 
static InputParameters baseParams ()
 

Protected Member Functions

virtual void computeScalarResidual () override
 Method for computing the scalar part of residual for _kappa. More...
 
virtual void computeScalarJacobian () override
 Method for computing the scalar variable part of Jacobian for d-_kappa-residual / d-_kappa. More...
 
virtual void computeScalarOffDiagJacobian (const unsigned int jvar_num) override
 Method for computing an off-diagonal jacobian component d-_kappa-residual / d-jvar. More...
 
virtual Real computeScalarQpOffDiagJacobian (const unsigned int jvar_num) override
 Method for computing an off-diagonal jacobian component at quadrature points. More...
 
virtual void computeOffDiagJacobianScalarLocal (const unsigned int svar_num) override
 Method for computing an off-diagonal jacobian component d-_var-residual / d-svar. More...
 
virtual Real computeQpOffDiagJacobianScalar (const unsigned int svar_num) override
 Method for computing d-_var-residual / d-svar at quadrature points. More...
 
virtual RankTwoTensor gradTest (unsigned int component) override
 
virtual RankTwoTensor gradTrial (unsigned int component) override
 
virtual void precalculateJacobianDisplacement (unsigned int component) override
 Prepare the average shape function gradients for stabilization. More...
 
virtual Real computeQpResidual () override
 
virtual Real computeQpJacobianDisplacement (unsigned int alpha, unsigned int beta) override
 
virtual Real computeQpJacobianTemperature (unsigned int cvar) override
 
virtual Real computeQpJacobianOutOfPlaneStrain () override
 
virtual void precalculateJacobian () override
 
virtual void precalculateOffDiagJacobian (unsigned int jvar) override
 
virtual Real computeQpJacobian () override
 
virtual Real computeQpOffDiagJacobian (unsigned int jvar) override
 

Protected Attributes

HomogenizationS::ConstraintMap _cmap
 Type of each constraint (stress or strain) for each component. More...
 
HomogenizationS::ConstraintType _ctype = HomogenizationS::ConstraintType::None
 The constraint type; initialize with 'none'. More...
 
unsigned int _m
 Used internally to iterate over each scalar component. More...
 
unsigned int _n
 
unsigned int _a
 
unsigned int _b
 
const MaterialProperty< RankTwoTensor > & _pk1
 The 1st Piola-Kirchhoff stress. More...
 
const MaterialProperty< RankFourTensor > & _dpk1
 The derivative of the PK1 stress with respect to the deformation gradient. More...
 
const bool _large_kinematics
 If true use large deformation kinematics. More...
 
const bool _stabilize_strain
 If true calculate the deformation gradient derivatives for F_bar. More...
 
const std::string _base_name
 Prepend to the material properties. More...
 
const unsigned int _alpha
 Which component of the vector residual this kernel is responsible for. More...
 
const unsigned int _ndisp
 Total number of displacements/size of residual vector. More...
 
std::vector< unsigned int_disp_nums
 The displacement numbers. More...
 
std::vector< std::vector< RankTwoTensor > > _avg_grad_trial
 
const MaterialProperty< RankTwoTensor > & _F_ust
 The unmodified deformation gradient. More...
 
const MaterialProperty< RankTwoTensor > & _F_avg
 The element-average deformation gradient. More...
 
const MaterialProperty< RankTwoTensor > & _f_inv
 The inverse increment deformation gradient. More...
 
const MaterialProperty< RankTwoTensor > & _F_inv
 The inverse deformation gradient. More...
 
const MaterialProperty< RankTwoTensor > & _F
 The actual (stabilized) deformation gradient. More...
 
const MooseVariable_temperature
 Temperature, if provided. This is used only to get the trial functions. More...
 
const MooseVariable_out_of_plane_strain
 Out-of-plane strain, if provided. More...
 
std::vector< std::vector< const MaterialProperty< RankTwoTensor > * > > _deigenstrain_dargs
 Eigenstrain derivatives wrt generate coupleds. More...
 
const unsigned int _n_args
 

Detailed Description

Total Lagrangian formulation with all homogenization terms (one disp_xyz field and macro_gradient scalar)

Definition at line 32 of file HomogenizedTotalLagrangianStressDivergenceS.h.

Constructor & Destructor Documentation

◆ HomogenizedTotalLagrangianStressDivergenceS()

HomogenizedTotalLagrangianStressDivergenceS::HomogenizedTotalLagrangianStressDivergenceS ( const InputParameters parameters)

Definition at line 37 of file HomogenizedTotalLagrangianStressDivergenceS.C.

40 {
41  // Constraint types
42  auto types = getParam<MultiMooseEnum>("constraint_types");
43  if (types.size() != Moose::dim * Moose::dim)
44  mooseError("Number of constraint types must equal dim * dim. ", types.size(), " are provided.");
45 
46  // Targets to hit
47  const std::vector<FunctionName> & fnames = getParam<std::vector<FunctionName>>("targets");
48 
49  // Prepare the constraint map
50  unsigned int fcount = 0;
51  for (const auto j : make_range(Moose::dim))
52  for (const auto i : make_range(Moose::dim))
53  {
54  const auto idx = i + Moose::dim * j;
55  const auto ctype = static_cast<HomogenizationS::ConstraintType>(types.get(idx));
57  {
58  const Function * const f = &getFunctionByName(fnames[fcount++]);
59  _cmap[{i, j}] = {ctype, f};
60  }
61  }
62 }
void mooseError(Args &&... args)
static constexpr std::size_t dim
Real f(Real x)
Test function for Brents method.
HomogenizationS::ConstraintMap _cmap
Type of each constraint (stress or strain) for each component.
ConstraintType
Constraint type: stress/PK stress or strain/deformation gradient.
TotalLagrangianStressDivergenceBaseS< GradientOperatorCartesian > TotalLagrangianStressDivergenceS
IntRange< T > make_range(T beg, T end)
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
unsigned int idx(const ElemType type, const unsigned int nx, const unsigned int i, const unsigned int j)

Member Function Documentation

◆ baseParams()

template<class G >
InputParameters TotalLagrangianStressDivergenceBaseS< G >::baseParams ( )
staticinherited

Definition at line 14 of file TotalLagrangianStressDivergenceBaseS.C.

Referenced by TotalLagrangianStressDivergenceBaseS< G >::validParams().

15 {
17  // This kernel requires use_displaced_mesh to be off
18  params.suppressParameter<bool>("use_displaced_mesh");
19  return params;
20 }
void suppressParameter(const std::string &name)

◆ computeOffDiagJacobianScalarLocal()

void HomogenizedTotalLagrangianStressDivergenceS::computeOffDiagJacobianScalarLocal ( const unsigned int  svar_num)
overrideprotectedvirtual

Method for computing an off-diagonal jacobian component d-_var-residual / d-svar.

svar is looped over all scalar variables, which herein is just _kappa

Definition at line 192 of file HomogenizedTotalLagrangianStressDivergenceS.C.

194 {
195 
196  // Get dofs and order of this scalar; at least one will be _kappa_var
197  const auto & svar = _sys.getScalarVariable(_tid, svar_num);
198  const unsigned int s_order = svar.order();
199  _local_ke.resize(_test.size(), s_order);
200 
201  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
202  {
203  unsigned int l = 0;
204  Real dV = _JxW[_qp] * _coord[_qp];
205  for (auto && [indices, constraint] : _cmap)
206  {
207  // copy constraint indices to protected variables to pass to Qp routine
208  _m = indices.first;
209  _n = indices.second;
210  _ctype = constraint.first;
211  initScalarQpJacobian(svar_num);
212  for (_i = 0; _i < _test.size(); _i++)
213  {
214  _local_ke(_i, l) += dV * computeQpOffDiagJacobianScalar(svar_num);
215  }
216  l++;
217  }
218  }
219 
220  addJacobian(_assembly, _local_ke, _var.dofIndices(), svar.dofIndices(), _var.scalingFactor());
221 }
unsigned int _m
Used internally to iterate over each scalar component.
HomogenizationS::ConstraintType _ctype
The constraint type; initialize with &#39;none&#39;.
HomogenizationS::ConstraintMap _cmap
Type of each constraint (stress or strain) for each component.
virtual Real computeQpOffDiagJacobianScalar(const unsigned int svar_num) override
Method for computing d-_var-residual / d-svar at quadrature points.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real

◆ computeQpJacobian()

Real LagrangianStressDivergenceBaseS::computeQpJacobian ( )
overrideprotectedvirtualinherited

Definition at line 115 of file LagrangianStressDivergenceBaseS.C.

116 {
118 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
virtual Real computeQpJacobianDisplacement(unsigned int alpha, unsigned int beta)=0

◆ computeQpJacobianDisplacement()

template<class G >
Real TotalLagrangianStressDivergenceBaseS< G >::computeQpJacobianDisplacement ( unsigned int  alpha,
unsigned int  beta 
)
overrideprotectedvirtualinherited

Implements LagrangianStressDivergenceBaseS.

Reimplemented in HomogenizedTotalLagrangianStressDivergenceA, and HomogenizedTotalLagrangianStressDivergenceR.

Definition at line 101 of file TotalLagrangianStressDivergenceBaseS.C.

103 {
104  // J_{alpha beta} = phi^alpha_{i, J} T_{iJkL} G^beta_{kL}
105  return gradTest(alpha).doubleContraction(_dpk1[_qp] * gradTrial(beta));
106 }
virtual RankTwoTensor gradTest(unsigned int component) override
Real doubleContraction(const RankTwoTensorTempl< Real > &a) const
static const std::string alpha
Definition: NS.h:134
const MaterialProperty< RankFourTensor > & _dpk1
The derivative of the PK1 stress with respect to the deformation gradient.
virtual RankTwoTensor gradTrial(unsigned int component) override

◆ computeQpJacobianOutOfPlaneStrain()

template<class G >
Real TotalLagrangianStressDivergenceBaseS< G >::computeQpJacobianOutOfPlaneStrain ( )
overrideprotectedvirtualinherited

Implements LagrangianStressDivergenceBaseS.

Definition at line 128 of file TotalLagrangianStressDivergenceBaseS.C.

129 {
130  return _dpk1[_qp].contractionKl(2, 2, gradTest(_alpha)) * _out_of_plane_strain->phi()[_j][_qp];
131 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
const MooseVariable * _out_of_plane_strain
Out-of-plane strain, if provided.
const FieldVariablePhiValue & phi() const override
virtual RankTwoTensor gradTest(unsigned int component) override
const MaterialProperty< RankFourTensor > & _dpk1
The derivative of the PK1 stress with respect to the deformation gradient.

◆ computeQpJacobianTemperature()

template<class G >
Real TotalLagrangianStressDivergenceBaseS< G >::computeQpJacobianTemperature ( unsigned int  cvar)
overrideprotectedvirtualinherited

Implements LagrangianStressDivergenceBaseS.

Definition at line 110 of file TotalLagrangianStressDivergenceBaseS.C.

111 {
112  usingTensorIndices(i_, j_, k_, l_);
113  // Multiple eigenstrains may depend on the same coupled var
114  RankTwoTensor total_deigen;
115  for (const auto deigen_darg : _deigenstrain_dargs[cvar])
116  total_deigen += (*deigen_darg)[_qp];
117 
118  const auto A = _f_inv[_qp].inverse();
119  const auto B = _F_inv[_qp].inverse();
120  const auto U = 0.5 * (A.template times<i_, k_, l_, j_>(B) + A.template times<i_, l_, k_, j_>(B));
121 
122  return -(_dpk1[_qp] * U * total_deigen).doubleContraction(gradTest(_alpha)) *
123  _temperature->phi()[_j][_qp];
124 }
const MooseVariable * _temperature
Temperature, if provided. This is used only to get the trial functions.
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
const MaterialProperty< RankTwoTensor > & _f_inv
The inverse increment deformation gradient.
const FieldVariablePhiValue & phi() const override
virtual RankTwoTensor gradTest(unsigned int component) override
std::vector< std::vector< const MaterialProperty< RankTwoTensor > * > > _deigenstrain_dargs
Eigenstrain derivatives wrt generate coupleds.
const MaterialProperty< RankTwoTensor > & _F_inv
The inverse deformation gradient.
const MaterialProperty< RankFourTensor > & _dpk1
The derivative of the PK1 stress with respect to the deformation gradient.

◆ computeQpOffDiagJacobian()

Real LagrangianStressDivergenceBaseS::computeQpOffDiagJacobian ( unsigned int  jvar)
overrideprotectedvirtualinherited

Definition at line 121 of file LagrangianStressDivergenceBaseS.C.

122 {
123  // Bail if jvar not coupled
124  if (getJvarMap()[jvar] < 0)
125  return 0.0;
126 
127  // Off diagonal terms for other displacements
128  for (auto beta : make_range(_ndisp))
129  if (jvar == _disp_nums[beta])
131 
132  // Off diagonal temperature term due to eigenstrain
133  if (_temperature && jvar == _temperature->number())
135 
136  // Off diagonal term due to weak plane stress
139 
140  return 0;
141 }
const MooseVariable * _temperature
Temperature, if provided. This is used only to get the trial functions.
std::vector< unsigned int > _disp_nums
The displacement numbers.
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
unsigned int number() const
const MooseVariable * _out_of_plane_strain
Out-of-plane strain, if provided.
const unsigned int _ndisp
Total number of displacements/size of residual vector.
virtual Real computeQpJacobianOutOfPlaneStrain()=0
virtual Real computeQpJacobianDisplacement(unsigned int alpha, unsigned int beta)=0
virtual Real computeQpJacobianTemperature(unsigned int cvar)=0
IntRange< T > make_range(T beg, T end)

◆ computeQpOffDiagJacobianScalar()

Real HomogenizedTotalLagrangianStressDivergenceS::computeQpOffDiagJacobianScalar ( const unsigned int  svar_num)
overrideprotectedvirtual

Method for computing d-_var-residual / d-svar at quadrature points.

Definition at line 224 of file HomogenizedTotalLagrangianStressDivergenceS.C.

Referenced by computeOffDiagJacobianScalarLocal().

225 {
226  // Just in case, skip any other scalar variables
227  if (svar_num == _kappa_var)
228  return _dpk1[_qp].contractionKl(_m, _n, gradTest(_alpha));
229  else
230  return 0.;
231 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
unsigned int _m
Used internally to iterate over each scalar component.
virtual RankTwoTensor gradTest(unsigned int component) override
const MaterialProperty< RankFourTensor > & _dpk1
The derivative of the PK1 stress with respect to the deformation gradient.

◆ computeQpResidual()

template<class G >
Real TotalLagrangianStressDivergenceBaseS< G >::computeQpResidual ( )
overrideprotectedvirtualinherited

Reimplemented in HomogenizedTotalLagrangianStressDivergenceA, and HomogenizedTotalLagrangianStressDivergenceR.

Definition at line 94 of file TotalLagrangianStressDivergenceBaseS.C.

95 {
96  return gradTest(_alpha).doubleContraction(_pk1[_qp]);
97 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
virtual RankTwoTensor gradTest(unsigned int component) override
Real doubleContraction(const RankTwoTensorTempl< Real > &a) const
const MaterialProperty< RankTwoTensor > & _pk1
The 1st Piola-Kirchhoff stress.

◆ computeScalarJacobian()

void HomogenizedTotalLagrangianStressDivergenceS::computeScalarJacobian ( )
overrideprotectedvirtual

Method for computing the scalar variable part of Jacobian for d-_kappa-residual / d-_kappa.

Definition at line 111 of file HomogenizedTotalLagrangianStressDivergenceS.C.

112 {
113  _local_ke.resize(_k_order, _k_order);
114 
115  // only assemble scalar residual once; i.e. when handling the first displacement component
116  if (_alpha == 0)
117  {
118  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
119  {
120  initScalarQpJacobian(_kappa_var);
121  Real dV = _JxW[_qp] * _coord[_qp];
122 
123  _h = 0;
124  for (auto && [indices1, constraint1] : _cmap)
125  {
126  auto && [i, j] = indices1;
127  auto && ctype = constraint1.first;
128  _l = 0;
129  for (auto && indices2 : _cmap)
130  {
131  auto && [a, b] = indices2.first;
133  _local_ke(_h, _l++) += dV * (_dpk1[_qp](i, j, a, b));
134  else if (ctype == HomogenizationS::ConstraintType::Strain)
135  {
136  if (_large_kinematics)
137  _local_ke(_h, _l++) += dV * (Real(i == a && j == b));
138  else
139  _local_ke(_h, _l++) +=
140  dV * (0.5 * Real(i == a && j == b) + 0.5 * Real(i == b && j == a));
141  }
142  else
143  mooseError("Unknown constraint type in Jacobian calculator!");
144  }
145  _h++;
146  }
147  }
148  }
149 
150  addJacobian(_assembly,
151  _local_ke,
152  _kappa_var_ptr->dofIndices(),
153  _kappa_var_ptr->dofIndices(),
154  _kappa_var_ptr->scalingFactor());
155 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
const bool _large_kinematics
If true use large deformation kinematics.
void mooseError(Args &&... args)
HomogenizationS::ConstraintMap _cmap
Type of each constraint (stress or strain) for each component.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
const MaterialProperty< RankFourTensor > & _dpk1
The derivative of the PK1 stress with respect to the deformation gradient.

◆ computeScalarOffDiagJacobian()

void HomogenizedTotalLagrangianStressDivergenceS::computeScalarOffDiagJacobian ( const unsigned int  jvar_num)
overrideprotectedvirtual

Method for computing an off-diagonal jacobian component d-_kappa-residual / d-jvar.

jvar is looped over all field variables, which herein is just disp_x and disp_y

Definition at line 158 of file HomogenizedTotalLagrangianStressDivergenceS.C.

160 {
161  const auto & jvar = getVariable(jvar_num);
162  // Get dofs and order of this variable; at least one will be _var
163  const auto jvar_size = jvar.phiSize();
164  _local_ke.resize(_k_order, jvar_size);
165 
166  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
167  {
168  // single index for Jacobian column; double indices for constraint tensor component
169  unsigned int h = 0;
170  Real dV = _JxW[_qp] * _coord[_qp];
171  for (auto && [indices, constraint] : _cmap)
172  {
173  // copy constraint indices to protected variables to pass to Qp routine
174  _m = indices.first;
175  _n = indices.second;
176  _ctype = constraint.first;
177  initScalarQpOffDiagJacobian(jvar);
178  for (_j = 0; _j < jvar_size; _j++)
179  _local_ke(h, _j) += dV * computeScalarQpOffDiagJacobian(jvar_num);
180  h++;
181  }
182  }
183 
184  addJacobian(_assembly,
185  _local_ke,
186  _kappa_var_ptr->dofIndices(),
187  jvar.dofIndices(),
188  _kappa_var_ptr->scalingFactor());
189 }
unsigned int _m
Used internally to iterate over each scalar component.
HomogenizationS::ConstraintType _ctype
The constraint type; initialize with &#39;none&#39;.
HomogenizationS::ConstraintMap _cmap
Type of each constraint (stress or strain) for each component.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
virtual Real computeScalarQpOffDiagJacobian(const unsigned int jvar_num) override
Method for computing an off-diagonal jacobian component at quadrature points.

◆ computeScalarQpOffDiagJacobian()

Real HomogenizedTotalLagrangianStressDivergenceS::computeScalarQpOffDiagJacobian ( const unsigned int  jvar_num)
overrideprotectedvirtual

Method for computing an off-diagonal jacobian component at quadrature points.

Definition at line 234 of file HomogenizedTotalLagrangianStressDivergenceS.C.

Referenced by computeScalarOffDiagJacobian().

235 {
236  // ONLY assemble the contribution from _alpha component, which is connected with _var
237  // The other components are handled by other kernel instances with other _alpha
238  if (jvar_num == _var.number())
239  {
241  return _dpk1[_qp].contractionIj(_m, _n, gradTrial(_alpha));
243  if (_large_kinematics)
244  return Real(_m == _alpha) * gradTrial(_alpha)(_m, _n);
245  else
246  return 0.5 * (Real(_m == _alpha) * gradTrial(_alpha)(_m, _n) +
247  Real(_n == _alpha) * gradTrial(_alpha)(_n, _m));
248  else
249  mooseError("Unknown constraint type in kernel calculation!");
250  }
251  else
252  return 0.;
253 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
const bool _large_kinematics
If true use large deformation kinematics.
void mooseError(Args &&... args)
unsigned int _m
Used internally to iterate over each scalar component.
HomogenizationS::ConstraintType _ctype
The constraint type; initialize with &#39;none&#39;.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
const MaterialProperty< RankFourTensor > & _dpk1
The derivative of the PK1 stress with respect to the deformation gradient.
virtual RankTwoTensor gradTrial(unsigned int component) override

◆ computeScalarResidual()

void HomogenizedTotalLagrangianStressDivergenceS::computeScalarResidual ( )
overrideprotectedvirtual

Method for computing the scalar part of residual for _kappa.

Definition at line 65 of file HomogenizedTotalLagrangianStressDivergenceS.C.

66 {
67  std::vector<Real> scalar_residuals(_k_order);
68 
69  // only assemble scalar residual once; i.e. when handling the first displacement component
70  if (_alpha == 0)
71  {
72  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
73  {
74  initScalarQpResidual();
75  Real dV = _JxW[_qp] * _coord[_qp];
76  _h = 0; // single index for residual vector; double indices for constraint tensor component
77  for (auto && [indices, constraint] : _cmap)
78  {
79  auto && [i, j] = indices;
80  auto && [ctype, ctarget] = constraint;
81 
83  {
85  scalar_residuals[_h++] += dV * (_pk1[_qp](i, j) - ctarget->value(_t, _q_point[_qp]));
87  scalar_residuals[_h++] +=
88  dV * (_F[_qp](i, j) - (Real(i == j) + ctarget->value(_t, _q_point[_qp])));
89  else
90  mooseError("Unknown constraint type in the integral!");
91  }
92  else
93  {
95  scalar_residuals[_h++] += dV * (_pk1[_qp](i, j) - ctarget->value(_t, _q_point[_qp]));
97  scalar_residuals[_h++] += dV * (0.5 * (_F[_qp](i, j) + _F[_qp](j, i)) -
98  (Real(i == j) + ctarget->value(_t, _q_point[_qp])));
99  else
100  mooseError("Unknown constraint type in the integral!");
101  }
102  }
103  }
104  }
105 
106  addResiduals(
107  _assembly, scalar_residuals, _kappa_var_ptr->dofIndices(), _kappa_var_ptr->scalingFactor());
108 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
const bool _large_kinematics
If true use large deformation kinematics.
void mooseError(Args &&... args)
HomogenizationS::ConstraintMap _cmap
Type of each constraint (stress or strain) for each component.
const MaterialProperty< RankTwoTensor > & _F
The actual (stabilized) deformation gradient.
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")
const MaterialProperty< RankTwoTensor > & _pk1
The 1st Piola-Kirchhoff stress.

◆ gradTest()

template<class G >
RankTwoTensor TotalLagrangianStressDivergenceBaseS< G >::gradTest ( unsigned int  component)
overrideprotectedvirtualinherited

◆ gradTrial()

template<class G >
RankTwoTensor TotalLagrangianStressDivergenceBaseS< G >::gradTrial ( unsigned int  component)
overrideprotectedvirtualinherited

◆ initialSetup() [1/4]

template<>
void TotalLagrangianStressDivergenceBaseS< GradientOperatorCartesian >::initialSetup ( )
inlineinherited

Definition at line 26 of file TotalLagrangianStressDivergenceS.h.

27 {
28  if (getBlockCoordSystem() != Moose::COORD_XYZ)
29  mooseError("This kernel should only act in Cartesian coordinates.");
30 }
void mooseError(Args &&... args)

◆ initialSetup() [2/4]

Definition at line 26 of file TotalLagrangianStressDivergenceAxisymmetricCylindricalS.h.

27 {
28  if (getBlockCoordSystem() != Moose::COORD_RZ)
29  mooseError("This kernel should only act in axisymmetric cylindrical coordinates.");
30 }
void mooseError(Args &&... args)

◆ initialSetup() [3/4]

Definition at line 26 of file TotalLagrangianStressDivergenceCentrosymmetricSphericalS.h.

27 {
28  if (getBlockCoordSystem() != Moose::COORD_RSPHERICAL)
29  mooseError("This kernel should only act in centrosymmetric spherical coordinates.");
30 }
void mooseError(Args &&... args)
COORD_RSPHERICAL

◆ initialSetup() [4/4]

template<class G >
virtual void TotalLagrangianStressDivergenceBaseS< G >::initialSetup ( )
overridevirtualinherited

◆ precalculateJacobian()

void LagrangianStressDivergenceBaseS::precalculateJacobian ( )
overrideprotectedvirtualinherited

Definition at line 85 of file LagrangianStressDivergenceBaseS.C.

86 {
87  // Skip if we are not doing stabilization
88  if (!_stabilize_strain)
89  return;
90 
91  // We need the gradients of shape functions in the reference frame
92  _fe_problem.prepareShapes(_var.number(), _tid);
93  _avg_grad_trial[_alpha].resize(_phi.size());
95 }
const unsigned int _alpha
Which component of the vector residual this kernel is responsible for.
const bool _stabilize_strain
If true calculate the deformation gradient derivatives for F_bar.
std::vector< std::vector< RankTwoTensor > > _avg_grad_trial
virtual void precalculateJacobianDisplacement(unsigned int component)=0
Prepare the average shape function gradients for stabilization.

◆ precalculateJacobianDisplacement()

template<class G >
void TotalLagrangianStressDivergenceBaseS< G >::precalculateJacobianDisplacement ( unsigned int  component)
overrideprotectedvirtualinherited

Prepare the average shape function gradients for stabilization.

Implements LagrangianStressDivergenceBaseS.

Definition at line 80 of file TotalLagrangianStressDivergenceBaseS.C.

81 {
82  // For total Lagrangian, the averaging is taken on the reference frame regardless of geometric
83  // nonlinearity. Convenient!
84  for (auto j : make_range(_phi.size()))
86  [this, component, j](unsigned int qp)
87  { return G::gradOp(component, _grad_phi[j][qp], _phi[j][qp], _q_point[qp]); },
88  _JxW,
89  _coord);
90 }
static const std::string component
Definition: NS.h:153
std::vector< std::vector< RankTwoTensor > > _avg_grad_trial
IntRange< T > make_range(T beg, T end)
auto elementAverage(const Functor &f, const MooseArray< Real > &JxW, const MooseArray< Real > &coord)
static const std::complex< double > j(0, 1)
Complex number "j" (also known as "i")

◆ precalculateOffDiagJacobian()

void LagrangianStressDivergenceBaseS::precalculateOffDiagJacobian ( unsigned int  jvar)
overrideprotectedvirtualinherited

Definition at line 98 of file LagrangianStressDivergenceBaseS.C.

99 {
100  // Skip if we are not doing stabilization
101  if (!_stabilize_strain)
102  return;
103 
104  for (auto beta : make_range(_ndisp))
105  if (jvar == _disp_nums[beta])
106  {
107  // We need the gradients of shape functions in the reference frame
108  _fe_problem.prepareShapes(jvar, _tid);
109  _avg_grad_trial[beta].resize(_phi.size());
111  }
112 }
std::vector< unsigned int > _disp_nums
The displacement numbers.
const unsigned int _ndisp
Total number of displacements/size of residual vector.
const bool _stabilize_strain
If true calculate the deformation gradient derivatives for F_bar.
std::vector< std::vector< RankTwoTensor > > _avg_grad_trial
IntRange< T > make_range(T beg, T end)
virtual void precalculateJacobianDisplacement(unsigned int component)=0
Prepare the average shape function gradients for stabilization.

◆ validParams()

InputParameters HomogenizedTotalLagrangianStressDivergenceS::validParams ( )
static

Definition at line 19 of file HomogenizedTotalLagrangianStressDivergenceS.C.

20 {
22  params.addClassDescription("Total Lagrangian stress equilibrium kernel with "
23  "homogenization constraint Jacobian terms");
24  params.renameCoupledVar(
25  "scalar_variable", "macro_var", "Optional scalar field with the macro gradient");
27  "constraint_types",
29  "Type of each constraint: strain, stress, or none. The types are specified in the "
30  "column-major order, and there must be 9 entries in total.");
31  params.addRequiredParam<std::vector<FunctionName>>(
32  "targets", "Functions giving the targets to hit for constraint types that are not none.");
33 
34  return params;
35 }
static InputParameters validParams()
void renameCoupledVar(const std::string &old_name, const std::string &new_name, const std::string &new_docstring)
const MultiMooseEnum constraintType("strain stress none")
Moose constraint type, for input.
void addRequiredParam(const std::string &name, const std::string &doc_string)
void addClassDescription(const std::string &doc_string)

Member Data Documentation

◆ _a

unsigned int HomogenizedTotalLagrangianStressDivergenceS::_a
protected

◆ _alpha

const unsigned int LagrangianStressDivergenceBaseS::_alpha
protectedinherited

◆ _avg_grad_trial

std::vector<std::vector<RankTwoTensor> > LagrangianStressDivergenceBaseS::_avg_grad_trial
protectedinherited

◆ _b

unsigned int HomogenizedTotalLagrangianStressDivergenceS::_b
protected

◆ _base_name

const std::string LagrangianStressDivergenceBaseS::_base_name
protectedinherited

Prepend to the material properties.

Definition at line 69 of file LagrangianStressDivergenceBaseS.h.

◆ _cmap

HomogenizationS::ConstraintMap HomogenizedTotalLagrangianStressDivergenceS::_cmap
protected

◆ _ctype

HomogenizationS::ConstraintType HomogenizedTotalLagrangianStressDivergenceS::_ctype = HomogenizationS::ConstraintType::None
protected

The constraint type; initialize with 'none'.

Definition at line 76 of file HomogenizedTotalLagrangianStressDivergenceS.h.

Referenced by computeOffDiagJacobianScalarLocal(), computeScalarOffDiagJacobian(), and computeScalarQpOffDiagJacobian().

◆ _deigenstrain_dargs

std::vector<std::vector<const MaterialProperty<RankTwoTensor> *> > LagrangianStressDivergenceBaseS::_deigenstrain_dargs
protectedinherited

Eigenstrain derivatives wrt generate coupleds.

Definition at line 107 of file LagrangianStressDivergenceBaseS.h.

Referenced by LagrangianStressDivergenceBaseS::LagrangianStressDivergenceBaseS().

◆ _disp_nums

std::vector<unsigned int> LagrangianStressDivergenceBaseS::_disp_nums
protectedinherited

◆ _dpk1

template<class G >
const MaterialProperty<RankFourTensor>& TotalLagrangianStressDivergenceBaseS< G >::_dpk1
protectedinherited

◆ _F

const MaterialProperty<RankTwoTensor>& LagrangianStressDivergenceBaseS::_F
protectedinherited

◆ _F_avg

const MaterialProperty<RankTwoTensor>& LagrangianStressDivergenceBaseS::_F_avg
protectedinherited

The element-average deformation gradient.

Definition at line 89 of file LagrangianStressDivergenceBaseS.h.

◆ _f_inv

const MaterialProperty<RankTwoTensor>& LagrangianStressDivergenceBaseS::_f_inv
protectedinherited

The inverse increment deformation gradient.

Definition at line 92 of file LagrangianStressDivergenceBaseS.h.

◆ _F_inv

const MaterialProperty<RankTwoTensor>& LagrangianStressDivergenceBaseS::_F_inv
protectedinherited

The inverse deformation gradient.

Definition at line 95 of file LagrangianStressDivergenceBaseS.h.

◆ _F_ust

const MaterialProperty<RankTwoTensor>& LagrangianStressDivergenceBaseS::_F_ust
protectedinherited

The unmodified deformation gradient.

Definition at line 86 of file LagrangianStressDivergenceBaseS.h.

◆ _large_kinematics

const bool LagrangianStressDivergenceBaseS::_large_kinematics
protectedinherited

◆ _m

unsigned int HomogenizedTotalLagrangianStressDivergenceS::_m
protected

◆ _n

unsigned int HomogenizedTotalLagrangianStressDivergenceS::_n
protected

◆ _ndisp

const unsigned int LagrangianStressDivergenceBaseS::_ndisp
protectedinherited

◆ _out_of_plane_strain

const MooseVariable* LagrangianStressDivergenceBaseS::_out_of_plane_strain
protectedinherited

Out-of-plane strain, if provided.

Definition at line 104 of file LagrangianStressDivergenceBaseS.h.

Referenced by LagrangianStressDivergenceBaseS::computeQpOffDiagJacobian().

◆ _pk1

template<class G >
const MaterialProperty<RankTwoTensor>& TotalLagrangianStressDivergenceBaseS< G >::_pk1
protectedinherited

◆ _stabilize_strain

const bool LagrangianStressDivergenceBaseS::_stabilize_strain
protectedinherited

If true calculate the deformation gradient derivatives for F_bar.

Definition at line 66 of file LagrangianStressDivergenceBaseS.h.

Referenced by LagrangianStressDivergenceBaseS::precalculateJacobian(), and LagrangianStressDivergenceBaseS::precalculateOffDiagJacobian().

◆ _temperature

const MooseVariable* LagrangianStressDivergenceBaseS::_temperature
protectedinherited

Temperature, if provided. This is used only to get the trial functions.

Definition at line 101 of file LagrangianStressDivergenceBaseS.h.

Referenced by LagrangianStressDivergenceBaseS::computeQpOffDiagJacobian().


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