www.mooseframework.org
EshelbyTensor.C
Go to the documentation of this file.
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 "EshelbyTensor.h"
11 #include "RankTwoTensor.h"
12 #include "MooseMesh.h"
13 
14 registerMooseObject("SolidMechanicsApp", EshelbyTensor);
15 registerMooseObject("SolidMechanicsApp", ADEshelbyTensor);
16 
17 template <bool is_ad>
20 {
22  params.addClassDescription("Computes the Eshelby tensor as a function of "
23  "strain energy density and the first "
24  "Piola-Kirchhoff stress");
25  params.addRequiredCoupledVar(
26  "displacements",
27  "The displacements appropriate for the simulation geometry and coordinate system");
28  params.addParam<std::string>("base_name",
29  "Optional parameter that allows the user to define "
30  "multiple mechanics material systems on the same "
31  "block, i.e. for multiple phases");
32  params.addParam<bool>(
33  "compute_dissipation",
34  false,
35  "Whether to compute Eshelby tensor's dissipation (or rate of change). This tensor"
36  "yields the increase in dissipation per unit crack advanced");
37  params.addCoupledVar("temperature", "Coupled temperature");
38  return params;
39 }
40 
41 template <bool is_ad>
44  _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
45  _compute_dissipation(getParam<bool>("compute_dissipation")),
46  _sed(getMaterialPropertyByName<Real>(_base_name + "strain_energy_density")),
47  _serd(_compute_dissipation
48  ? &getMaterialPropertyByName<Real>(_base_name + "strain_energy_rate_density")
49  : nullptr),
50  _eshelby_tensor(declareProperty<RankTwoTensor>(_base_name + "Eshelby_tensor")),
51  _eshelby_tensor_dissipation(
52  _compute_dissipation
53  ? &declareProperty<RankTwoTensor>(_base_name + "Eshelby_tensor_dissipation")
54  : nullptr),
55  _stress(getGenericMaterialProperty<RankTwoTensor, is_ad>(_base_name + "stress")),
56  _stress_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "stress")),
57  _grad_disp(3),
58  _grad_disp_old(3),
59  _J_thermal_term_vec(declareProperty<RealVectorValue>("J_thermal_term_vec")),
60  _grad_temp(coupledGradient("temperature")),
61  _has_temp(isCoupled("temperature")),
62  _total_deigenstrain_dT(getOptionalMaterialProperty<RankTwoTensor>("total_deigenstrain_dT"))
63 {
64  unsigned int ndisp = coupledComponents("displacements");
65 
66  // Checking for consistency between mesh size and length of the provided displacements vector
67  if (ndisp != _mesh.dimension())
68  mooseError(
69  "The number of variables supplied in 'displacements' must match the mesh dimension.");
70 
71  // fetch coupled gradients
72  for (unsigned int i = 0; i < ndisp; ++i)
73  _grad_disp[i] = &coupledGradient("displacements", i);
74 
75  // set unused dimensions to zero
76  for (unsigned i = ndisp; i < 3; ++i)
77  _grad_disp[i] = &_grad_zero;
78 
79  // Need previous step's displacements to compute deformation gradient time rate
81  {
82  // fetch coupled gradients previous step
83  for (unsigned int i = 0; i < ndisp; ++i)
84  _grad_disp_old[i] = &coupledGradientOld("displacements", i);
85 
86  // set unused dimensions to zero
87  for (unsigned i = ndisp; i < 3; ++i)
89  }
90 }
91 
92 template <bool is_ad>
93 void
95 {
96  if (_has_temp && !_total_deigenstrain_dT)
97  mooseError("EshelbyTensor Error: To include thermal strain term in Fracture integral "
98  "calculation, must both couple temperature in DomainIntegral block and compute "
99  "total_deigenstrain_dT using ThermalFractureIntegral material model.");
100 }
101 
102 template <bool is_ad>
103 void
105 {
106 }
107 
108 template <bool is_ad>
109 void
111 {
112  // Deformation gradient
114  (*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
115 
116  RankTwoTensor H(F);
117  F.addIa(1.0);
118  Real detF = F.det();
119  RankTwoTensor FinvT(F.inverse().transpose());
120 
121  // 1st Piola-Kirchoff Stress (P):
122  RankTwoTensor P = detF * MetaPhysicL::raw_value(_stress[_qp]) * FinvT;
123 
124  // HTP = H^T * P = H^T * detF * sigma * FinvT;
125  RankTwoTensor HTP = H.transpose() * P;
126 
128  WI *= (_sed[_qp] * detF);
129 
130  _eshelby_tensor[_qp] = WI - HTP;
131 
132  // Compute deformation gradient rate
133  if (_compute_dissipation == true)
134  {
136  (*_grad_disp_old[0])[_qp], (*_grad_disp_old[1])[_qp], (*_grad_disp_old[2])[_qp]);
137 
139  Wdot *= ((*_serd)[_qp] * detF);
140 
141  // F_dot = (F - F_old)/dt
142  RankTwoTensor F_dot = (H - H_old) / _dt;
143 
144  // FdotTP = Fdot^T * P = Fdot^T * detF * sigma * FinvT;
145  RankTwoTensor FdotTP = F_dot.transpose() * P;
146 
147  (*_eshelby_tensor_dissipation)[_qp] = Wdot - FdotTP;
148  }
149 
150  if (_has_temp)
151  {
152  const Real sigma_alpha =
153  MetaPhysicL::raw_value(_stress[_qp]).doubleContraction(_total_deigenstrain_dT[_qp]);
154  _J_thermal_term_vec[_qp] = sigma_alpha * _grad_temp[_qp];
155  }
156  else
157  _J_thermal_term_vec[_qp].zero();
158 }
registerMooseObject("SolidMechanicsApp", EshelbyTensor)
void addParam(const std::string &name, const std::initializer_list< typename T::value_type > &value, const std::string &doc_string)
virtual void initQpStatefulProperties() override
const VariableGradient & _grad_zero
EshelbyTensorTempl(const InputParameters &parameters)
Definition: EshelbyTensor.C:42
void mooseError(Args &&... args)
EshelbyTensor defines a strain increment and rotation increment, for finite strains.
Definition: EshelbyTensor.h:21
virtual void computeQpProperties() override
auto raw_value(const Eigen::Map< T > &in)
virtual const VariableGradient & coupledGradient(const std::string &var_name, unsigned int comp=0) const
std::vector< const VariableGradient * > _grad_disp_old
Definition: EshelbyTensor.h:52
static const std::string F
Definition: NS.h:150
static RankTwoTensorTempl initializeFromRows(const libMesh::TypeVector< Real > &row0, const libMesh::TypeVector< Real > &row1, const libMesh::TypeVector< Real > &row2)
static InputParameters validParams()
Definition: EshelbyTensor.C:19
virtual const VariableGradient & coupledGradientOld(const std::string &var_name, unsigned int comp=0) const
static InputParameters validParams()
virtual unsigned int dimension() const
const bool _compute_dissipation
Whether to also compute Eshelby tensor&#39;s dissipation for C(t) integral.
Definition: EshelbyTensor.h:37
std::vector< const VariableGradient * > _grad_disp
Definition: EshelbyTensor.h:51
void addCoupledVar(const std::string &name, const std::string &doc_string)
void addRequiredCoupledVar(const std::string &name, const std::string &doc_string)
unsigned int coupledComponents(const std::string &var_name) const
RankTwoTensorTempl< Real > transpose() const
DIE A HORRIBLE DEATH HERE typedef LIBMESH_DEFAULT_SCALAR_TYPE Real
void mooseError(Args &&... args) const
void addClassDescription(const std::string &doc_string)
virtual void initialSetup() override
Definition: EshelbyTensor.C:94