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

This postprocessor computes the Interaction Integral. More...

#include <InteractionIntegralSM.h>

Inheritance diagram for InteractionIntegralSM:
[legend]

Public Member Functions

 InteractionIntegralSM (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< SymmTensor > & _stress
 
const MaterialProperty< SymmTensor > & _strain
 
std::vector< const VariableGradient * > _grad_disp
 
const bool _has_temp
 
const VariableGradient & _grad_temp
 
const MaterialProperty< Real > * _current_instantaneous_thermal_expansion_coef
 
Real _K_factor
 
bool _has_symmetry_plane
 
Real _poissons_ratio
 
Real _youngs_modulus
 
unsigned int _ring_index
 
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 25 of file InteractionIntegralSM.h.

Member Enumeration Documentation

enum InteractionIntegralSM::QMethod
strongprivate
Enumerator
Geometry 
Topology 

Definition at line 67 of file InteractionIntegralSM.h.

68  {
69  Geometry,
70  Topology
71  };
Enumerator
KI 
KII 
KIII 

Definition at line 75 of file InteractionIntegralSM.h.

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

Constructor & Destructor Documentation

InteractionIntegralSM::InteractionIntegralSM ( const InputParameters &  parameters)

Definition at line 66 of file InteractionIntegralSM.C.

67  : ElementIntegralPostprocessor(parameters),
68  _ndisp(coupledComponents("displacements")),
69  _crack_front_definition(&getUserObject<CrackFrontDefinition>("crack_front_definition")),
70  _has_crack_front_point_index(isParamValid("crack_front_point_index")),
72  _has_crack_front_point_index ? getParam<unsigned int>("crack_front_point_index") : 0),
73  _treat_as_2d(false),
74  _stress(getMaterialPropertyByName<SymmTensor>("stress")),
75  _strain(getMaterialPropertyByName<SymmTensor>("elastic_strain")),
76  _grad_disp(3),
77  _has_temp(isCoupled("temp")),
78  _grad_temp(_has_temp ? coupledGradient("temp") : _grad_zero),
80  hasMaterialProperty<Real>("current_instantaneous_thermal_expansion_coef")
81  ? &getMaterialProperty<Real>("current_instantaneous_thermal_expansion_coef")
82  : NULL),
83  _K_factor(getParam<Real>("K_factor")),
84  _has_symmetry_plane(isParamValid("symmetry_plane")),
85  _poissons_ratio(getParam<Real>("poissons_ratio")),
86  _youngs_modulus(getParam<Real>("youngs_modulus")),
87  _ring_index(getParam<unsigned int>("ring_index")),
88  _q_function_type(getParam<MooseEnum>("q_function_type").getEnum<QMethod>()),
89  _sif_mode(getParam<MooseEnum>("sif_mode").getEnum<SifMethod>())
90 {
92  mooseError("To include thermal strain term in interaction integral, must both couple "
93  "temperature in DomainIntegral block and compute thermal expansion property in "
94  "material model using compute_InteractionIntegral = true.");
95 
96  // plane strain
97  _kappa = 3.0 - 4.0 * _poissons_ratio;
99 
100  // Checking for consistency between mesh size and length of the provided displacements vector
101  if (_ndisp != _mesh.dimension())
102  mooseError(
103  "The number of variables supplied in 'displacements' must match the mesh dimension.");
104  // fetch gradients of coupled variables
105  for (unsigned int i = 0; i < _ndisp; ++i)
106  _grad_disp[i] = &coupledGradient("displacements", i);
107 
108  // set unused dimensions to zero
109  for (unsigned i = _ndisp; i < 3; ++i)
110  _grad_disp[i] = &_grad_zero;
111 }
const unsigned int _crack_front_point_index
const CrackFrontDefinition *const _crack_front_definition
const MaterialProperty< SymmTensor > & _stress
const MaterialProperty< Real > * _current_instantaneous_thermal_expansion_coef
std::vector< const VariableGradient * > _grad_disp
const VariableGradient & _grad_temp
const MaterialProperty< SymmTensor > & _strain

Member Function Documentation

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

Definition at line 284 of file InteractionIntegralSM.C.

Referenced by computeQpIntegral().

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

Definition at line 245 of file InteractionIntegralSM.C.

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

Definition at line 133 of file InteractionIntegralSM.C.

Referenced by computeIntegral().

134 {
135  Real scalar_q = 0.0;
136  RealVectorValue grad_q(0.0, 0.0, 0.0);
137 
138  const std::vector<std::vector<Real>> & phi_curr_elem = *_phi_curr_elem;
139  const std::vector<std::vector<RealGradient>> & dphi_curr_elem = *_dphi_curr_elem;
140 
141  for (unsigned int i = 0; i < _current_elem->n_nodes(); ++i)
142  {
143  scalar_q += phi_curr_elem[i][_qp] * _q_curr_elem[i];
144 
145  for (unsigned int j = 0; j < _current_elem->dim(); ++j)
146  grad_q(j) += dphi_curr_elem[i][_qp](j) * _q_curr_elem[i];
147  }
148 
149  // In the crack front coordinate system, the crack direction is (1,0,0)
150  RealVectorValue crack_direction(0.0);
151  crack_direction(0) = 1.0;
152 
153  // Calculate (r,theta) position of qp relative to crack front
154  Point p(_q_point[_qp]);
156 
157  RankTwoTensor aux_stress;
158  RankTwoTensor aux_du;
159 
161  computeAuxFields(aux_stress, aux_du);
162  else if (_sif_mode == SifMethod::T)
163  computeTFields(aux_stress, aux_du);
164 
165  RankTwoTensor stress;
166  stress(0, 0) = _stress[_qp].xx();
167  stress(0, 1) = _stress[_qp].xy();
168  stress(0, 2) = _stress[_qp].xz();
169  stress(1, 0) = _stress[_qp].xy();
170  stress(1, 1) = _stress[_qp].yy();
171  stress(1, 2) = _stress[_qp].yz();
172  stress(2, 0) = _stress[_qp].xz();
173  stress(2, 1) = _stress[_qp].yz();
174  stress(2, 2) = _stress[_qp].zz();
175 
176  RankTwoTensor strain;
177  strain(0, 0) = _strain[_qp].xx();
178  strain(0, 1) = _strain[_qp].xy();
179  strain(0, 2) = _strain[_qp].xz();
180  strain(1, 0) = _strain[_qp].xy();
181  strain(1, 1) = _strain[_qp].yy();
182  strain(1, 2) = _strain[_qp].yz();
183  strain(2, 0) = _strain[_qp].xz();
184  strain(2, 1) = _strain[_qp].yz();
185  strain(2, 2) = _strain[_qp].zz();
186 
187  RankTwoTensor grad_disp((*_grad_disp[0])[_qp], (*_grad_disp[1])[_qp], (*_grad_disp[2])[_qp]);
188 
189  // Rotate stress, strain, displacement and temperature to crack front coordinate system
190  RealVectorValue grad_q_cf =
192  RankTwoTensor grad_disp_cf =
194  RankTwoTensor stress_cf =
196  RankTwoTensor strain_cf =
198  RealVectorValue grad_temp_cf =
200 
201  RankTwoTensor dq;
202  dq(0, 0) = crack_direction(0) * grad_q_cf(0);
203  dq(0, 1) = crack_direction(0) * grad_q_cf(1);
204  dq(0, 2) = crack_direction(0) * grad_q_cf(2);
205 
206  // Calculate interaction integral terms
207 
208  // Term1 = stress * x1-derivative of aux disp * dq
209  RankTwoTensor tmp1 = dq * stress_cf;
210  Real term1 = aux_du.doubleContraction(tmp1);
211 
212  // Term2 = aux stress * x1-derivative of disp * dq
213  RankTwoTensor tmp2 = dq * aux_stress;
214  Real term2 = grad_disp_cf(0, 0) * tmp2(0, 0) + grad_disp_cf(1, 0) * tmp2(0, 1) +
215  grad_disp_cf(2, 0) * tmp2(0, 2);
216 
217  // Term3 = aux stress * strain * dq_x (= stress * aux strain * dq_x)
218  Real term3 = dq(0, 0) * aux_stress.doubleContraction(strain_cf);
219 
220  // Term4 (thermal strain term) = q * aux_stress * alpha * dtheta_x
221  // - the term including the derivative of alpha is not implemented
222  Real term4 = 0.0;
223  if (_has_temp)
224  term4 = scalar_q * aux_stress.trace() * (*_current_instantaneous_thermal_expansion_coef)[_qp] *
225  grad_temp_cf(0);
226 
227  Real q_avg_seg = 1.0;
229  {
230  q_avg_seg =
233  2.0;
234  }
235 
236  Real eq = term1 + term2 - term3 + term4;
237 
239  eq *= 2.0;
240 
241  return eq / q_avg_seg;
242 }
const unsigned int _crack_front_point_index
void computeAuxFields(RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
const CrackFrontDefinition *const _crack_front_definition
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
void computeTFields(RankTwoTensor &aux_stress, RankTwoTensor &grad_disp)
const MaterialProperty< SymmTensor > & _stress
const std::vector< std::vector< Real > > * _phi_curr_elem
Real getCrackFrontBackwardSegmentLength(const unsigned int point_index) const
RealVectorValue rotateToCrackFrontCoords(const RealVectorValue vector, const unsigned int point_index) const
std::vector< const VariableGradient * > _grad_disp
Real getCrackFrontForwardSegmentLength(const unsigned int point_index) const
std::vector< Real > _q_curr_elem
const std::vector< std::vector< RealGradient > > * _dphi_curr_elem
const VariableGradient & _grad_temp
const MaterialProperty< SymmTensor > & _strain
void InteractionIntegralSM::computeTFields ( RankTwoTensor &  aux_stress,
RankTwoTensor &  grad_disp 
)
protected

Definition at line 346 of file InteractionIntegralSM.C.

Referenced by computeQpIntegral().

347 {
348 
349  Real t = _theta;
350  Real st = std::sin(t);
351  Real ct = std::cos(t);
352  Real stsq = Utility::pow<2>(st);
353  Real ctsq = Utility::pow<2>(ct);
354  Real ctcu = Utility::pow<3>(ct);
355  Real oneOverPiR = 1.0 / (libMesh::pi * _r);
356 
357  aux_stress.zero();
358  aux_stress(0, 0) = -oneOverPiR * ctcu;
359  aux_stress(0, 1) = -oneOverPiR * st * ctsq;
360  aux_stress(1, 0) = -oneOverPiR * st * ctsq;
361  aux_stress(1, 1) = -oneOverPiR * ct * stsq;
362  aux_stress(2, 2) = -oneOverPiR * _poissons_ratio * (ctcu + ct * stsq);
363 
364  grad_disp.zero();
365  grad_disp(0, 0) = oneOverPiR / (4.0 * _youngs_modulus) *
366  (ct * (4.0 * Utility::pow<2>(_poissons_ratio) - 3.0 + _poissons_ratio) -
367  std::cos(3.0 * t) * (1.0 + _poissons_ratio));
368  grad_disp(0, 1) = -oneOverPiR / (4.0 * _youngs_modulus) *
369  (st * (4.0 * Utility::pow<2>(_poissons_ratio) - 3.0 + _poissons_ratio) +
370  std::sin(3.0 * t) * (1.0 + _poissons_ratio));
371 }
Real InteractionIntegralSM::getValue ( )
virtual

Definition at line 120 of file InteractionIntegralSM.C.

121 {
122  gatherSum(_integral_value);
123 
125  _integral_value +=
128 
129  return _K_factor * _integral_value;
130 }
const unsigned int _crack_front_point_index
Real getCrackFrontTangentialStrain(const unsigned int node_index) const
const CrackFrontDefinition *const _crack_front_definition
void InteractionIntegralSM::initialSetup ( )
protectedvirtual

Definition at line 114 of file InteractionIntegralSM.C.

115 {
117 }
const CrackFrontDefinition *const _crack_front_definition
MooseEnum InteractionIntegralSM::qFunctionType ( )
static

Definition at line 17 of file InteractionIntegralSM.C.

Referenced by validParams< InteractionIntegralSM >().

18 {
19  return MooseEnum("Geometry Topology", "Geometry");
20 }
MooseEnum InteractionIntegralSM::sifModeType ( )
static

Definition at line 23 of file InteractionIntegralSM.C.

Referenced by validParams< InteractionIntegralSM >().

24 {
25  return MooseEnum("KI KII KIII T", "KI");
26 }

Member Data Documentation

const CrackFrontDefinition* const InteractionIntegralSM::_crack_front_definition
protected

Definition at line 43 of file InteractionIntegralSM.h.

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

const unsigned int InteractionIntegralSM::_crack_front_point_index
protected

Definition at line 45 of file InteractionIntegralSM.h.

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

const MaterialProperty<Real>* InteractionIntegralSM::_current_instantaneous_thermal_expansion_coef
protected

Definition at line 52 of file InteractionIntegralSM.h.

Referenced by InteractionIntegralSM().

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

Definition at line 60 of file InteractionIntegralSM.h.

Referenced by computeIntegral(), and computeQpIntegral().

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

Definition at line 49 of file InteractionIntegralSM.h.

Referenced by computeQpIntegral(), and InteractionIntegralSM().

const VariableGradient& InteractionIntegralSM::_grad_temp
protected

Definition at line 51 of file InteractionIntegralSM.h.

Referenced by computeQpIntegral().

bool InteractionIntegralSM::_has_crack_front_point_index
protected

Definition at line 44 of file InteractionIntegralSM.h.

bool InteractionIntegralSM::_has_symmetry_plane
protected

Definition at line 54 of file InteractionIntegralSM.h.

Referenced by computeQpIntegral().

const bool InteractionIntegralSM::_has_temp
protected

Definition at line 50 of file InteractionIntegralSM.h.

Referenced by computeQpIntegral(), and InteractionIntegralSM().

Real InteractionIntegralSM::_K_factor
protected

Definition at line 53 of file InteractionIntegralSM.h.

Referenced by getValue().

Real InteractionIntegralSM::_kappa
protected

Definition at line 61 of file InteractionIntegralSM.h.

Referenced by computeAuxFields(), and InteractionIntegralSM().

unsigned int InteractionIntegralSM::_ndisp
protected

Definition at line 42 of file InteractionIntegralSM.h.

Referenced by InteractionIntegralSM().

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

Definition at line 59 of file InteractionIntegralSM.h.

Referenced by computeIntegral(), and computeQpIntegral().

Real InteractionIntegralSM::_poissons_ratio
protected
std::vector<Real> InteractionIntegralSM::_q_curr_elem
protected

Definition at line 58 of file InteractionIntegralSM.h.

Referenced by computeIntegral(), and computeQpIntegral().

const QMethod InteractionIntegralSM::_q_function_type
private

Definition at line 73 of file InteractionIntegralSM.h.

Referenced by computeIntegral().

Real InteractionIntegralSM::_r
protected

Definition at line 63 of file InteractionIntegralSM.h.

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

unsigned int InteractionIntegralSM::_ring_index
protected

Definition at line 57 of file InteractionIntegralSM.h.

Referenced by computeIntegral().

Real InteractionIntegralSM::_shear_modulus
protected

Definition at line 62 of file InteractionIntegralSM.h.

Referenced by computeAuxFields(), and InteractionIntegralSM().

const SifMethod InteractionIntegralSM::_sif_mode
private

Definition at line 83 of file InteractionIntegralSM.h.

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

const MaterialProperty<SymmTensor>& InteractionIntegralSM::_strain
protected

Definition at line 48 of file InteractionIntegralSM.h.

Referenced by computeQpIntegral().

const MaterialProperty<SymmTensor>& InteractionIntegralSM::_stress
protected

Definition at line 47 of file InteractionIntegralSM.h.

Referenced by computeQpIntegral().

Real InteractionIntegralSM::_theta
protected

Definition at line 64 of file InteractionIntegralSM.h.

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

bool InteractionIntegralSM::_treat_as_2d
protected

Definition at line 46 of file InteractionIntegralSM.h.

Referenced by getValue(), and initialSetup().

Real InteractionIntegralSM::_youngs_modulus
protected

Definition at line 56 of file InteractionIntegralSM.h.

Referenced by computeTFields(), and InteractionIntegralSM().


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