www.mooseframework.org
EshelbyTensor.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 
8 #include "EshelbyTensor.h"
9 #include "RankTwoTensor.h"
10 #include "MooseMesh.h"
11 
12 template <>
13 InputParameters
15 {
16  InputParameters params = validParams<Material>();
17  params.addClassDescription("Stuff");
18  params.addRequiredCoupledVar(
19  "displacements",
20  "The displacements appropriate for the simulation geometry and coordinate system");
21  params.addParam<std::string>("base_name",
22  "Optional parameter that allows the user to define "
23  "multiple mechanics material systems on the same "
24  "block, i.e. for multiple phases");
25  params.addCoupledVar("temperature", "Coupled temperature");
26  return params;
27 }
28 
29 EshelbyTensor::EshelbyTensor(const InputParameters & parameters)
30  : DerivativeMaterialInterface<Material>(parameters),
31  _base_name(isParamValid("base_name") ? getParam<std::string>("base_name") + "_" : ""),
32  _sed(declareProperty<Real>(_base_name + "strain_energy_density")),
33  _sed_old(getMaterialPropertyOld<Real>(_base_name + "strain_energy_density")),
34  _eshelby_tensor(declareProperty<RankTwoTensor>(_base_name + "Eshelby_tensor")),
35  _stress(getMaterialProperty<RankTwoTensor>(_base_name + "stress")),
36  _stress_old(getMaterialPropertyOld<RankTwoTensor>(_base_name + "stress")),
37  _grad_disp(3),
38  _J_thermal_term_vec(declareProperty<RealVectorValue>("J_thermal_term_vec")),
39  _has_temp(isCoupled("temperature")),
40  _grad_temp(coupledGradient("temperature")),
41  _total_deigenstrain_dT(hasMaterialProperty<RankTwoTensor>("total_deigenstrain_dT")
42  ? &getMaterialProperty<RankTwoTensor>("total_deigenstrain_dT")
43  : nullptr)
44 {
45  unsigned int ndisp = coupledComponents("displacements");
46 
47  // Checking for consistency between mesh size and length of the provided displacements vector
48  if (ndisp != _mesh.dimension())
49  mooseError(
50  "The number of variables supplied in 'displacements' must match the mesh dimension.");
51 
52  // fetch coupled gradients
53  for (unsigned int i = 0; i < ndisp; ++i)
54  _grad_disp[i] = &coupledGradient("displacements", i);
55 
56  // set unused dimensions to zero
57  for (unsigned i = ndisp; i < 3; ++i)
58  _grad_disp[i] = &_grad_zero;
59 
61  mooseError("EshelbyTensor Error: To include thermal strain term in Fracture integral "
62  "calculation, must both couple temperature in DomainIntegral block and compute "
63  "total_deigenstrain_dT using ThermalFractureIntegral material model.");
64 
65  if (hasMaterialProperty<RankTwoTensor>(_base_name + "strain_increment"))
66  {
67  _strain_increment = &getMaterialProperty<RankTwoTensor>(_base_name + "strain_increment");
68  _mechanical_strain = nullptr;
69  }
70  else if (hasMaterialProperty<RankTwoTensor>(_base_name + "mechanical_strain"))
71  {
72  _mechanical_strain = &getMaterialProperty<RankTwoTensor>(_base_name + "mechanical_strain");
73  _strain_increment = nullptr;
74  }
75  else
76  mooseError("EshelbyTensor cannot find either mechanical_strain or strain_increment material "
77  "properties.");
78 }
79 
80 void
82 {
83  _sed[_qp] = 0.0;
84 }
85 
86 void
88 {
89  RankTwoTensor F((*_grad_disp[0])[_qp],
90  (*_grad_disp[1])[_qp],
91  (*_grad_disp[2])[_qp]); // Deformation gradient
92 
93  RankTwoTensor H(F);
94  F.addIa(1.0);
95  Real detF = F.det();
96  RankTwoTensor FinvT(F.inverse().transpose());
97 
98  if (_strain_increment != nullptr)
99  _sed[_qp] = _sed_old[_qp] + _stress[_qp].doubleContraction((*_strain_increment)[_qp]) / 2.0 +
100  _stress_old[_qp].doubleContraction((*_strain_increment)[_qp]) / 2.0;
101  else
102  _sed[_qp] = _stress[_qp].doubleContraction((*_mechanical_strain)[_qp]) / 2.0;
103 
104  // 1st Piola-Kirchoff Stress (P):
105  RankTwoTensor P = detF * _stress[_qp] * FinvT;
106 
107  // HTP = H^T * P = H^T * detF * sigma * FinvT;
108  RankTwoTensor HTP = H.transpose() * P;
109 
110  RankTwoTensor WI = RankTwoTensor(RankTwoTensor::initIdentity);
111  WI *= (_sed[_qp] * detF);
112 
113  _eshelby_tensor[_qp] = WI - HTP;
114 
115  if (_has_temp)
116  {
117  Real sigma_alpha = _stress[_qp].doubleContraction((*_total_deigenstrain_dT)[_qp]);
118  _J_thermal_term_vec[_qp] = sigma_alpha * _grad_temp[_qp];
119  }
120  else
121  _J_thermal_term_vec[_qp].zero();
122 }
const MaterialProperty< RankTwoTensor > * _total_deigenstrain_dT
Definition: EshelbyTensor.h:41
virtual void initQpStatefulProperties() override
Definition: EshelbyTensor.C:81
const MaterialProperty< Real > & _sed_old
Definition: EshelbyTensor.h:30
MaterialProperty< RealVectorValue > & _J_thermal_term_vec
Definition: EshelbyTensor.h:38
MaterialProperty< Real > & _sed
Definition: EshelbyTensor.h:29
std::vector< const VariableGradient * > _grad_disp
Definition: EshelbyTensor.h:36
const MaterialProperty< RankTwoTensor > * _mechanical_strain
Definition: EshelbyTensor.h:35
const VariableGradient & _grad_temp
Definition: EshelbyTensor.h:40
const MaterialProperty< RankTwoTensor > & _stress_old
Definition: EshelbyTensor.h:33
InputParameters validParams< EshelbyTensor >()
Definition: EshelbyTensor.C:14
const MaterialProperty< RankTwoTensor > & _stress
Definition: EshelbyTensor.h:32
const bool _has_temp
Definition: EshelbyTensor.h:39
EshelbyTensor(const InputParameters &parameters)
Definition: EshelbyTensor.C:29
virtual void computeQpProperties() override
Definition: EshelbyTensor.C:87
std::string _base_name
Definition: EshelbyTensor.h:27
const MaterialProperty< RankTwoTensor > * _strain_increment
Definition: EshelbyTensor.h:34
MaterialProperty< RankTwoTensor > & _eshelby_tensor
Definition: EshelbyTensor.h:31