www.mooseframework.org
Public Member Functions | Static Public Member Functions | Protected Member Functions | Protected Attributes | Private Types | Private Attributes | List of all members
InteractionIntegral Class Reference

This postprocessor computes the Interaction Integral. More...

#include <InteractionIntegral.h>

Inheritance diagram for InteractionIntegral:
[legend]

Public Member Functions

 InteractionIntegral (const InputParameters &parameters)
 
virtual Real getValue ()
 

Static Public Member Functions

static MooseEnum qFunctionType ()
 
static MooseEnum sifModeType ()
 

Protected Member Functions

virtual void initialSetup ()
 
virtual Real computeQpIntegral ()
 
virtual Real computeIntegral ()
 
void computeAuxFields (RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
 
void computeTFields (RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
 

Protected Attributes

unsigned int _ndisp
 
const CrackFrontDefinition *const _crack_front_definition
 
bool _has_crack_front_point_index
 
const unsigned int _crack_front_point_index
 
bool _treat_as_2d
 
const MaterialProperty< RankTwoTensor > * _stress
 
const MaterialProperty< RankTwoTensor > * _strain
 
std::vector< const VariableGradient * > _grad_disp
 
const bool _has_temp
 
const VariableGradient & _grad_temp
 
Real _K_factor
 
bool _has_symmetry_plane
 
Real _poissons_ratio
 
Real _youngs_modulus
 
unsigned int _ring_index
 
const MaterialProperty< RankTwoTensor > * _total_deigenstrain_dT
 
std::vector< Real > _q_curr_elem
 
const std::vector< std::vector< Real > > * _phi_curr_elem
 
const std::vector< std::vector< RealGradient > > * _dphi_curr_elem
 
Real _kappa
 
Real _shear_modulus
 
Real _r
 
Real _theta
 

Private Types

enum  QMethod { QMethod::Geometry, QMethod::Topology }
 
enum  SifMethod { SifMethod::KI, SifMethod::KII, SifMethod::KIII, SifMethod::T }
 

Private Attributes

const QMethod _q_function_type
 
const SifMethod _sif_mode
 

Detailed Description

This postprocessor computes the Interaction Integral.

Definition at line 27 of file InteractionIntegral.h.

Member Enumeration Documentation

◆ QMethod

enum InteractionIntegral::QMethod
strongprivate
Enumerator
Geometry 
Topology 

Definition at line 68 of file InteractionIntegral.h.

69  {
70  Geometry,
71  Topology
72  };

◆ SifMethod

enum InteractionIntegral::SifMethod
strongprivate
Enumerator
KI 
KII 
KIII 

Definition at line 76 of file InteractionIntegral.h.

77  {
78  KI,
79  KII,
80  KIII,
81  T
82  };

Constructor & Destructor Documentation

◆ InteractionIntegral()

InteractionIntegral::InteractionIntegral ( const InputParameters &  parameters)

Definition at line 75 of file InteractionIntegral.C.

76  : ElementIntegralPostprocessor(parameters),
77  _ndisp(coupledComponents("displacements")),
78  _crack_front_definition(&getUserObject<CrackFrontDefinition>("crack_front_definition")),
79  _has_crack_front_point_index(isParamValid("crack_front_point_index")),
81  _has_crack_front_point_index ? getParam<unsigned int>("crack_front_point_index") : 0),
82  _treat_as_2d(false),
83  _stress(hasMaterialProperty<RankTwoTensor>("stress")
84  ? &getMaterialPropertyByName<RankTwoTensor>("stress")
85  : nullptr),
86  _strain(hasMaterialProperty<RankTwoTensor>("elastic_strain")
87  ? &getMaterialPropertyByName<RankTwoTensor>("elastic_strain")
88  : nullptr),
89  _grad_disp(3),
90  _has_temp(isCoupled("temperature")),
91  _grad_temp(_has_temp ? coupledGradient("temperature") : _grad_zero),
92  _K_factor(getParam<Real>("K_factor")),
93  _has_symmetry_plane(isParamValid("symmetry_plane")),
94  _poissons_ratio(getParam<Real>("poissons_ratio")),
95  _youngs_modulus(getParam<Real>("youngs_modulus")),
96  _ring_index(getParam<unsigned int>("ring_index")),
97  _total_deigenstrain_dT(hasMaterialProperty<RankTwoTensor>("total_deigenstrain_dT")
98  ? &getMaterialProperty<RankTwoTensor>("total_deigenstrain_dT")
99  : nullptr),
100  _q_function_type(getParam<MooseEnum>("q_function_type").getEnum<QMethod>()),
101  _sif_mode(getParam<MooseEnum>("sif_mode").getEnum<SifMethod>())
102 {
103  if (!hasMaterialProperty<RankTwoTensor>("stress"))
104  mooseError("InteractionIntegral Error: RankTwoTensor material property 'stress' not found. "
105  "This may be because solid mechanics system is being used to calculate a SymmTensor "
106  "'stress' material property. To use interaction integral calculation with solid "
107  "mechanics application, please set 'solid_mechanics = true' in the DomainIntegral "
108  "block.");
109 
110  if (!hasMaterialProperty<RankTwoTensor>("elastic_strain"))
111  mooseError("InteractionIntegral Error: RankTwoTensor material property 'elastic_strain' not "
112  "found. This may be because solid mechanics system is being used to calculate a "
113  "SymmTensor 'elastic_strain' material property. To use interaction integral "
114  "calculation with solid mechanics application, please set 'solid_mechanics = true' "
115  "in the DomainIntegral block.");
116 
118  mooseError("InteractionIntegral Error: To include thermal strain term in interaction integral, "
119  "must both couple temperature in DomainIntegral block and compute "
120  "total_deigenstrain_dT using ThermalFractureIntegral material model.");
121 
122  // plane strain
123  _kappa = 3.0 - 4.0 * _poissons_ratio;
124  _shear_modulus = _youngs_modulus / (2.0 * (1.0 + _poissons_ratio));
125 
126  // Checking for consistency between mesh size and length of the provided displacements vector
127  if (_ndisp != _mesh.dimension())
128  mooseError("InteractionIntegral Error: number of variables supplied in 'displacements' must "
129  "match the mesh dimension.");
130 
131  // fetch gradients of coupled variables
132  for (unsigned int i = 0; i < _ndisp; ++i)
133  _grad_disp[i] = &coupledGradient("displacements", i);
134 
135  // set unused dimensions to zero
136  for (unsigned i = _ndisp; i < 3; ++i)
137  _grad_disp[i] = &_grad_zero;
138 }
const unsigned int _crack_front_point_index
const QMethod _q_function_type
const MaterialProperty< RankTwoTensor > * _stress
const CrackFrontDefinition *const _crack_front_definition
std::vector< const VariableGradient * > _grad_disp
const MaterialProperty< RankTwoTensor > * _strain
const SifMethod _sif_mode
const MaterialProperty< RankTwoTensor > * _total_deigenstrain_dT
const VariableGradient & _grad_temp

Member Function Documentation

◆ computeAuxFields()

void InteractionIntegral::computeAuxFields ( RankTwoTensor &  aux_stress,
RankTwoTensor &  grad_disp 
)
protected

Definition at line 291 of file InteractionIntegral.C.

Referenced by computeQpIntegral().

292 {
293  RealVectorValue k(0.0);
294  if (_sif_mode == SifMethod::KI)
295  k(0) = 1.0;
296  else if (_sif_mode == SifMethod::KII)
297  k(1) = 1.0;
298  else if (_sif_mode == SifMethod::KIII)
299  k(2) = 1.0;
300 
301  Real t = _theta;
302  Real t2 = _theta / 2.0;
303  Real tt2 = 3.0 * _theta / 2.0;
304  Real st = std::sin(t);
305  Real ct = std::cos(t);
306  Real st2 = std::sin(t2);
307  Real ct2 = std::cos(t2);
308  Real stt2 = std::sin(tt2);
309  Real ctt2 = std::cos(tt2);
310  Real ct2sq = Utility::pow<2>(ct2);
311  Real ct2cu = Utility::pow<3>(ct2);
312  Real sqrt2PiR = std::sqrt(2.0 * libMesh::pi * _r);
313 
314  // Calculate auxiliary stress tensor
315  aux_stress.zero();
316 
317  aux_stress(0, 0) =
318  1.0 / sqrt2PiR * (k(0) * ct2 * (1.0 - st2 * stt2) - k(1) * st2 * (2.0 + ct2 * ctt2));
319  aux_stress(1, 1) = 1.0 / sqrt2PiR * (k(0) * ct2 * (1.0 + st2 * stt2) + k(1) * st2 * ct2 * ctt2);
320  aux_stress(0, 1) = 1.0 / sqrt2PiR * (k(0) * ct2 * st2 * ctt2 + k(1) * ct2 * (1.0 - st2 * stt2));
321  aux_stress(0, 2) = -1.0 / sqrt2PiR * k(2) * st2;
322  aux_stress(1, 2) = 1.0 / sqrt2PiR * k(2) * ct2;
323  // plane stress
324  // Real s33 = 0;
325  // plane strain
326  aux_stress(2, 2) = _poissons_ratio * (aux_stress(0, 0) + aux_stress(1, 1));
327 
328  aux_stress(1, 0) = aux_stress(0, 1);
329  aux_stress(2, 0) = aux_stress(0, 2);
330  aux_stress(2, 1) = aux_stress(1, 2);
331 
332  // Calculate x1 derivative of auxiliary displacements
333  grad_disp.zero();
334 
335  grad_disp(0, 0) = k(0) / (4.0 * _shear_modulus * sqrt2PiR) *
336  (ct * ct2 * _kappa + ct * ct2 - 2.0 * ct * ct2cu + st * st2 * _kappa +
337  st * st2 - 6.0 * st * st2 * ct2sq) +
338  k(1) / (4.0 * _shear_modulus * sqrt2PiR) *
339  (ct * st2 * _kappa + ct * st2 + 2.0 * ct * st2 * ct2sq - st * ct2 * _kappa +
340  3.0 * st * ct2 - 6.0 * st * ct2cu);
341 
342  grad_disp(0, 1) = k(0) / (4.0 * _shear_modulus * sqrt2PiR) *
343  (ct * st2 * _kappa + ct * st2 - 2.0 * ct * st2 * ct2sq - st * ct2 * _kappa -
344  5.0 * st * ct2 + 6.0 * st * ct2cu) +
345  k(1) / (4.0 * _shear_modulus * sqrt2PiR) *
346  (-ct * ct2 * _kappa + 3.0 * ct * ct2 - 2.0 * ct * ct2cu -
347  st * st2 * _kappa + 3.0 * st * st2 - 6.0 * st * st2 * ct2sq);
348 
349  grad_disp(0, 2) = k(2) / (_shear_modulus * sqrt2PiR) * (st2 * ct - ct2 * st);
350 }
const SifMethod _sif_mode

◆ computeIntegral()

Real InteractionIntegral::computeIntegral ( )
protectedvirtual

Definition at line 252 of file InteractionIntegral.C.

253 {
254  Real sum = 0.0;
255 
256  // calculate phi and dphi for this element
257  FEType fe_type(Utility::string_to_enum<Order>("first"),
258  Utility::string_to_enum<FEFamily>("lagrange"));
259  const unsigned int dim = _current_elem->dim();
260  std::unique_ptr<FEBase> fe(FEBase::build(dim, fe_type));
261  fe->attach_quadrature_rule(_qrule);
262  _phi_curr_elem = &fe->get_phi();
263  _dphi_curr_elem = &fe->get_dphi();
264  fe->reinit(_current_elem);
265 
266  // calculate q for all nodes in this element
267  _q_curr_elem.clear();
268  unsigned int ring_base = (_q_function_type == QMethod::Topology) ? 0 : 1;
269 
270  for (unsigned int i = 0; i < _current_elem->n_nodes(); ++i)
271  {
272  Node * this_node = _current_elem->get_node(i);
273  Real q_this_node;
274 
277  _crack_front_point_index, _ring_index - ring_base, this_node);
280  _crack_front_point_index, _ring_index - ring_base, this_node);
281 
282  _q_curr_elem.push_back(q_this_node);
283  }
284 
285  for (_qp = 0; _qp < _qrule->n_points(); _qp++)
286  sum += _JxW[_qp] * _coord[_qp] * computeQpIntegral();
287  return sum;
288 }
const std::vector< std::vector< Real > > * _phi_curr_elem
const unsigned int _crack_front_point_index
const QMethod _q_function_type
const CrackFrontDefinition *const _crack_front_definition
const std::vector< std::vector< RealGradient > > * _dphi_curr_elem
virtual Real computeQpIntegral()
Real DomainIntegralTopologicalQFunction(unsigned int crack_front_point_index, unsigned int ring_index, const Node *const current_node) const
std::vector< Real > _q_curr_elem
Real DomainIntegralQFunction(unsigned int crack_front_point_index, unsigned int ring_index, const Node *const current_node) const

◆ computeQpIntegral()

Real InteractionIntegral::computeQpIntegral ( )
protectedvirtual

Definition at line 160 of file InteractionIntegral.C.

Referenced by computeIntegral().

161 {
162  Real scalar_q = 0.0;
163  RealVectorValue grad_q(0.0, 0.0, 0.0);
164 
165  const std::vector<std::vector<Real>> & phi_curr_elem = *_phi_curr_elem;
166  const std::vector<std::vector<RealGradient>> & dphi_curr_elem = *_dphi_curr_elem;
167 
168  for (unsigned int i = 0; i < _current_elem->n_nodes(); ++i)
169  {
170  scalar_q += phi_curr_elem[i][_qp] * _q_curr_elem[i];
171 
172  for (unsigned int j = 0; j < _current_elem->dim(); ++j)
173  grad_q(j) += dphi_curr_elem[i][_qp](j) * _q_curr_elem[i];
174  }
175 
176  // In the crack front coordinate system, the crack direction is (1,0,0)
177  RealVectorValue crack_direction(0.0);
178  crack_direction(0) = 1.0;
179 
180  // Calculate (r,theta) position of qp relative to crack front
181  Point p(_q_point[_qp]);
183 
184  RankTwoTensor aux_stress;
185  RankTwoTensor aux_du;
186 
188  computeAuxFields(aux_stress, aux_du);
189  else if (_sif_mode == SifMethod::T)
190  computeTFields(aux_stress, aux_du);
191 
192  RankTwoTensor grad_disp((*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
193 
194  // Rotate stress, strain, displacement and temperature to crack front coordinate system
195  RealVectorValue grad_q_cf =
197  RankTwoTensor grad_disp_cf =
199  RankTwoTensor stress_cf =
201  RankTwoTensor strain_cf =
203  RealVectorValue grad_temp_cf =
205 
206  RankTwoTensor dq;
207  dq(0, 0) = crack_direction(0) * grad_q_cf(0);
208  dq(0, 1) = crack_direction(0) * grad_q_cf(1);
209  dq(0, 2) = crack_direction(0) * grad_q_cf(2);
210 
211  // Calculate interaction integral terms
212 
213  // Term1 = stress * x1-derivative of aux disp * dq
214  RankTwoTensor tmp1 = dq * stress_cf;
215  Real term1 = aux_du.doubleContraction(tmp1);
216 
217  // Term2 = aux stress * x1-derivative of disp * dq
218  RankTwoTensor tmp2 = dq * aux_stress;
219  Real term2 = grad_disp_cf(0, 0) * tmp2(0, 0) + grad_disp_cf(1, 0) * tmp2(0, 1) +
220  grad_disp_cf(2, 0) * tmp2(0, 2);
221 
222  // Term3 = aux stress * strain * dq_x (= stress * aux strain * dq_x)
223  Real term3 = dq(0, 0) * aux_stress.doubleContraction(strain_cf);
224 
225  // Term4 (thermal strain term) = q * aux_stress * alpha * dtheta_x
226  // - the term including the derivative of alpha is not implemented
227  Real term4 = 0.0;
228  if (_has_temp)
229  {
230  Real sigma_alpha = aux_stress.doubleContraction((*_total_deigenstrain_dT)[_qp]);
231  term4 = scalar_q * sigma_alpha * grad_temp_cf(0);
232  }
233 
234  Real q_avg_seg = 1.0;
236  {
237  q_avg_seg =
240  2.0;
241  }
242 
243  Real eq = term1 + term2 - term3 + term4;
244 
246  eq *= 2.0;
247 
248  return eq / q_avg_seg;
249 }
void computeAuxFields(RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
void calculateRThetaToCrackFront(const Point qp, const unsigned int point_index, Real &r, Real &theta) const
calculate r and theta in the crack front polar cooridnate
Real getCrackFrontBackwardSegmentLength(const unsigned int point_index) const
const std::vector< std::vector< Real > > * _phi_curr_elem
const unsigned int _crack_front_point_index
const MaterialProperty< RankTwoTensor > * _stress
const CrackFrontDefinition *const _crack_front_definition
std::vector< const VariableGradient * > _grad_disp
void computeTFields(RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
const std::vector< std::vector< RealGradient > > * _dphi_curr_elem
const MaterialProperty< RankTwoTensor > * _strain
RealVectorValue rotateToCrackFrontCoords(const RealVectorValue vector, const unsigned int point_index) const
const SifMethod _sif_mode
Real getCrackFrontForwardSegmentLength(const unsigned int point_index) const
const MaterialProperty< RankTwoTensor > * _total_deigenstrain_dT
const VariableGradient & _grad_temp
std::vector< Real > _q_curr_elem

◆ computeTFields()

void InteractionIntegral::computeTFields ( RankTwoTensor &  aux_stress,
RankTwoTensor &  grad_disp 
)
protected

Definition at line 353 of file InteractionIntegral.C.

Referenced by computeQpIntegral().

354 {
355 
356  Real t = _theta;
357  Real st = std::sin(t);
358  Real ct = std::cos(t);
359  Real stsq = Utility::pow<2>(st);
360  Real ctsq = Utility::pow<2>(ct);
361  Real ctcu = Utility::pow<3>(ct);
362  Real oneOverPiR = 1.0 / (libMesh::pi * _r);
363 
364  aux_stress.zero();
365  aux_stress(0, 0) = -oneOverPiR * ctcu;
366  aux_stress(0, 1) = -oneOverPiR * st * ctsq;
367  aux_stress(1, 0) = -oneOverPiR * st * ctsq;
368  aux_stress(1, 1) = -oneOverPiR * ct * stsq;
369  aux_stress(2, 2) = -oneOverPiR * _poissons_ratio * (ctcu + ct * stsq);
370 
371  grad_disp.zero();
372  grad_disp(0, 0) = oneOverPiR / (4.0 * _youngs_modulus) *
373  (ct * (4.0 * Utility::pow<2>(_poissons_ratio) - 3.0 + _poissons_ratio) -
374  std::cos(3.0 * t) * (1.0 + _poissons_ratio));
375  grad_disp(0, 1) = -oneOverPiR / (4.0 * _youngs_modulus) *
376  (st * (4.0 * Utility::pow<2>(_poissons_ratio) - 3.0 + _poissons_ratio) +
377  std::sin(3.0 * t) * (1.0 + _poissons_ratio));
378 }

◆ getValue()

Real InteractionIntegral::getValue ( )
virtual

Definition at line 147 of file InteractionIntegral.C.

148 {
149  gatherSum(_integral_value);
150 
152  _integral_value +=
155 
156  return _K_factor * _integral_value;
157 }
const unsigned int _crack_front_point_index
const CrackFrontDefinition *const _crack_front_definition
const SifMethod _sif_mode
Real getCrackFrontTangentialStrain(const unsigned int node_index) const

◆ initialSetup()

void InteractionIntegral::initialSetup ( )
protectedvirtual

Definition at line 141 of file InteractionIntegral.C.

142 {
144 }
const CrackFrontDefinition *const _crack_front_definition

◆ qFunctionType()

MooseEnum InteractionIntegral::qFunctionType ( )
static

Definition at line 21 of file InteractionIntegral.C.

Referenced by validParams< InteractionIntegral >().

22 {
23  return MooseEnum("Geometry Topology", "Geometry");
24 }

◆ sifModeType()

MooseEnum InteractionIntegral::sifModeType ( )
static

Definition at line 27 of file InteractionIntegral.C.

Referenced by validParams< InteractionIntegral >().

28 {
29  return MooseEnum("KI KII KIII T", "KI");
30 }

Member Data Documentation

◆ _crack_front_definition

const CrackFrontDefinition* const InteractionIntegral::_crack_front_definition
protected

Definition at line 44 of file InteractionIntegral.h.

Referenced by computeIntegral(), computeQpIntegral(), getValue(), and initialSetup().

◆ _crack_front_point_index

const unsigned int InteractionIntegral::_crack_front_point_index
protected

Definition at line 46 of file InteractionIntegral.h.

Referenced by computeIntegral(), computeQpIntegral(), and getValue().

◆ _dphi_curr_elem

const std::vector<std::vector<RealGradient> >* InteractionIntegral::_dphi_curr_elem
protected

Definition at line 61 of file InteractionIntegral.h.

Referenced by computeIntegral(), and computeQpIntegral().

◆ _grad_disp

std::vector<const VariableGradient *> InteractionIntegral::_grad_disp
protected

Definition at line 50 of file InteractionIntegral.h.

Referenced by computeQpIntegral(), and InteractionIntegral().

◆ _grad_temp

const VariableGradient& InteractionIntegral::_grad_temp
protected

Definition at line 52 of file InteractionIntegral.h.

Referenced by computeQpIntegral().

◆ _has_crack_front_point_index

bool InteractionIntegral::_has_crack_front_point_index
protected

Definition at line 45 of file InteractionIntegral.h.

◆ _has_symmetry_plane

bool InteractionIntegral::_has_symmetry_plane
protected

Definition at line 54 of file InteractionIntegral.h.

Referenced by computeQpIntegral().

◆ _has_temp

const bool InteractionIntegral::_has_temp
protected

Definition at line 51 of file InteractionIntegral.h.

Referenced by computeQpIntegral(), and InteractionIntegral().

◆ _K_factor

Real InteractionIntegral::_K_factor
protected

Definition at line 53 of file InteractionIntegral.h.

Referenced by getValue().

◆ _kappa

Real InteractionIntegral::_kappa
protected

Definition at line 62 of file InteractionIntegral.h.

Referenced by computeAuxFields(), and InteractionIntegral().

◆ _ndisp

unsigned int InteractionIntegral::_ndisp
protected

Definition at line 43 of file InteractionIntegral.h.

Referenced by InteractionIntegral().

◆ _phi_curr_elem

const std::vector<std::vector<Real> >* InteractionIntegral::_phi_curr_elem
protected

Definition at line 60 of file InteractionIntegral.h.

Referenced by computeIntegral(), and computeQpIntegral().

◆ _poissons_ratio

Real InteractionIntegral::_poissons_ratio
protected

◆ _q_curr_elem

std::vector<Real> InteractionIntegral::_q_curr_elem
protected

Definition at line 59 of file InteractionIntegral.h.

Referenced by computeIntegral(), and computeQpIntegral().

◆ _q_function_type

const QMethod InteractionIntegral::_q_function_type
private

Definition at line 74 of file InteractionIntegral.h.

Referenced by computeIntegral().

◆ _r

Real InteractionIntegral::_r
protected

Definition at line 64 of file InteractionIntegral.h.

Referenced by computeAuxFields(), computeQpIntegral(), and computeTFields().

◆ _ring_index

unsigned int InteractionIntegral::_ring_index
protected

Definition at line 57 of file InteractionIntegral.h.

Referenced by computeIntegral().

◆ _shear_modulus

Real InteractionIntegral::_shear_modulus
protected

Definition at line 63 of file InteractionIntegral.h.

Referenced by computeAuxFields(), and InteractionIntegral().

◆ _sif_mode

const SifMethod InteractionIntegral::_sif_mode
private

Definition at line 84 of file InteractionIntegral.h.

Referenced by computeAuxFields(), computeQpIntegral(), and getValue().

◆ _strain

const MaterialProperty<RankTwoTensor>* InteractionIntegral::_strain
protected

Definition at line 49 of file InteractionIntegral.h.

Referenced by computeQpIntegral().

◆ _stress

const MaterialProperty<RankTwoTensor>* InteractionIntegral::_stress
protected

Definition at line 48 of file InteractionIntegral.h.

Referenced by computeQpIntegral().

◆ _theta

Real InteractionIntegral::_theta
protected

Definition at line 65 of file InteractionIntegral.h.

Referenced by computeAuxFields(), computeQpIntegral(), and computeTFields().

◆ _total_deigenstrain_dT

const MaterialProperty<RankTwoTensor>* InteractionIntegral::_total_deigenstrain_dT
protected

Definition at line 58 of file InteractionIntegral.h.

Referenced by computeQpIntegral(), and InteractionIntegral().

◆ _treat_as_2d

bool InteractionIntegral::_treat_as_2d
protected

Definition at line 47 of file InteractionIntegral.h.

Referenced by getValue(), and initialSetup().

◆ _youngs_modulus

Real InteractionIntegral::_youngs_modulus
protected

Definition at line 56 of file InteractionIntegral.h.

Referenced by computeTFields(), and InteractionIntegral().


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