www.mooseframework.org
Public Member Functions | Protected Member Functions | Protected Attributes | List of all members
ComputeIsotropicLinearElasticPFFractureStress Class Reference

Phase-field fracture This class computes the stress and energy contribution for the small strain Isotropic Elastic formulation of phase field fracture. More...

#include <ComputeIsotropicLinearElasticPFFractureStress.h>

Inheritance diagram for ComputeIsotropicLinearElasticPFFractureStress:
[legend]

Public Member Functions

 ComputeIsotropicLinearElasticPFFractureStress (const InputParameters &parameters)
 

Protected Member Functions

virtual void initQpStatefulProperties ()
 Function required to initialize statefull material properties. More...
 
virtual void computeQpStress ()
 Primary funciton of this material, computes stress and elastic energy. More...
 
virtual void computeQpProperties () override
 

Protected Attributes

const VariableValue & _c
 Coupled order parameter defining the crack. More...
 
const Real _kdamage
 Small number to avoid non-positive definiteness at or near complete damage. More...
 
bool _use_current_hist
 Use current value of history variable. More...
 
const MaterialProperty< Real > & _l
 Material property defining crack width, declared elsewhere. More...
 
const MaterialProperty< Real > & _gc
 Material property defining gc parameter, declared elsewhere. More...
 
MaterialProperty< Real > & _F
 Elastic energy and derivatives, declared in this material. More...
 
MaterialProperty< Real > & _dFdc
 
MaterialProperty< Real > & _d2Fdc2
 
MaterialProperty< RankTwoTensor > & _dstress_dc
 
MaterialProperty< RankTwoTensor > & _d2Fdcdstrain
 
MaterialProperty< Real > & _hist
 History variable that prevents crack healing, declared in this material. More...
 
const MaterialProperty< Real > & _hist_old
 Old value of history variable. More...
 
const std::string _base_name
 
const std::string _elasticity_tensor_name
 
const MaterialProperty< RankTwoTensor > & _mechanical_strain
 
MaterialProperty< RankTwoTensor > & _stress
 
MaterialProperty< RankTwoTensor > & _elastic_strain
 
const MaterialProperty< RankFourTensor > & _elasticity_tensor
 
const MaterialProperty< RankTwoTensor > & _extra_stress
 Extra stress tensor. More...
 
std::vector< Function * > _initial_stress_fcn
 initial stress components More...
 
MaterialProperty< RankFourTensor > & _Jacobian_mult
 derivative of stress w.r.t. strain (_dstress_dstrain) More...
 

Detailed Description

Phase-field fracture This class computes the stress and energy contribution for the small strain Isotropic Elastic formulation of phase field fracture.

Definition at line 25 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Constructor & Destructor Documentation

◆ ComputeIsotropicLinearElasticPFFractureStress()

ComputeIsotropicLinearElasticPFFractureStress::ComputeIsotropicLinearElasticPFFractureStress ( const InputParameters &  parameters)

Definition at line 31 of file ComputeIsotropicLinearElasticPFFractureStress.C.

33  : ComputeStressBase(parameters),
34  _c(coupledValue("c")),
35  _kdamage(getParam<Real>("kdamage")),
36  _use_current_hist(getParam<bool>("use_current_history_variable")),
37  _l(getMaterialProperty<Real>("l")),
38  _gc(getMaterialProperty<Real>("gc_prop")),
39  _F(declareProperty<Real>(getParam<MaterialPropertyName>("F_name"))),
40  _dFdc(declarePropertyDerivative<Real>(getParam<MaterialPropertyName>("F_name"),
41  getVar("c", 0)->name())),
42  _d2Fdc2(declarePropertyDerivative<Real>(
43  getParam<MaterialPropertyName>("F_name"), getVar("c", 0)->name(), getVar("c", 0)->name())),
45  declarePropertyDerivative<RankTwoTensor>(_base_name + "stress", getVar("c", 0)->name())),
46  _d2Fdcdstrain(declareProperty<RankTwoTensor>("d2Fdcdstrain")),
47  _hist(declareProperty<Real>("hist")),
48  _hist_old(getMaterialPropertyOld<Real>("hist"))
49 {
50 }
const VariableValue & _c
Coupled order parameter defining the crack.
const Real _kdamage
Small number to avoid non-positive definiteness at or near complete damage.
const MaterialProperty< Real > & _gc
Material property defining gc parameter, declared elsewhere.
const std::string name
Definition: Setup.h:22
const MaterialProperty< Real > & _l
Material property defining crack width, declared elsewhere.
const std::string _base_name
ComputeStressBase(const InputParameters &parameters)
MaterialProperty< Real > & _hist
History variable that prevents crack healing, declared in this material.
MaterialProperty< Real > & _F
Elastic energy and derivatives, declared in this material.
const MaterialProperty< Real > & _hist_old
Old value of history variable.

Member Function Documentation

◆ computeQpProperties()

void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Definition at line 51 of file ComputeStressBase.C.

52 {
54 
55  // Add in extra stress
56  _stress[_qp] += _extra_stress[_qp];
57 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.

◆ computeQpStress()

void ComputeIsotropicLinearElasticPFFractureStress::computeQpStress ( )
protectedvirtual

Primary funciton of this material, computes stress and elastic energy.

Implements ComputeStressBase.

Definition at line 59 of file ComputeIsotropicLinearElasticPFFractureStress.C.

60 {
61  const Real c = _c[_qp];
62 
63  // Isotropic elasticity is assumed and should be enforced
64  const Real lambda = _elasticity_tensor[_qp](0, 0, 1, 1);
65  const Real mu = _elasticity_tensor[_qp](0, 1, 0, 1);
66 
67  // Compute eigenvectors and eigenvalues of mechanical strain and projection tensor
68  RankTwoTensor eigvec;
69  std::vector<Real> eigval(LIBMESH_DIM);
70  RankFourTensor proj_pos =
71  _mechanical_strain[_qp].positiveProjectionEigenDecomposition(eigval, eigvec);
72  RankFourTensor I4sym(RankFourTensor::initIdentitySymmetricFour);
73  RankFourTensor proj_neg = I4sym - proj_pos;
74 
75  // Calculate tensors of outerproduct of eigen vectors
76  std::vector<RankTwoTensor> etens(LIBMESH_DIM);
77 
78  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
79  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
80  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
81  etens[i](j, k) = eigvec(j, i) * eigvec(k, i);
82 
83  // Separate out positive and negative eigen values
84  std::vector<Real> epos(LIBMESH_DIM), eneg(LIBMESH_DIM);
85  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
86  {
87  epos[i] = (std::abs(eigval[i]) + eigval[i]) / 2.0;
88  eneg[i] = (std::abs(eigval[i]) - eigval[i]) / 2.0;
89  }
90 
91  // Seprate positive and negative sums of all eigenvalues
92  Real etr = 0.0;
93  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
94  etr += eigval[i];
95 
96  const Real etrpos = (std::abs(etr) + etr) / 2.0;
97  const Real etrneg = (std::abs(etr) - etr) / 2.0;
98 
99  // Calculate the tensile (postive) and compressive (negative) parts of stress
100  RankTwoTensor stress0pos, stress0neg;
101  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
102  {
103  stress0pos += etens[i] * (lambda * etrpos + 2.0 * mu * epos[i]);
104  stress0neg += etens[i] * (lambda * etrneg + 2.0 * mu * eneg[i]);
105  }
106 
107  // sum squares of epos and eneg
108  Real pval(0.0), nval(0.0);
109  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
110  {
111  pval += epos[i] * epos[i];
112  nval += eneg[i] * eneg[i];
113  }
114 
115  // Energy with positive principal strains
116  const Real G0_pos = lambda * etrpos * etrpos / 2.0 + mu * pval;
117  const Real G0_neg = lambda * etrneg * etrneg / 2.0 + mu * nval;
118 
119  // Assign history variable and derivative
120  if (G0_pos > _hist_old[_qp])
121  _hist[_qp] = G0_pos;
122  else
123  _hist[_qp] = _hist_old[_qp];
124 
125  Real hist_variable = _hist_old[_qp];
126  if (_use_current_hist)
127  hist_variable = _hist[_qp];
128 
129  // Damage associated with positive component of stress
130  _stress[_qp] = stress0pos * ((1.0 - c) * (1.0 - c) * (1 - _kdamage) + _kdamage) - stress0neg;
131 
132  // Elastic free energy density
133  _F[_qp] = hist_variable * ((1.0 - c) * (1.0 - c) * (1 - _kdamage) + _kdamage) - G0_neg +
134  _gc[_qp] / (2 * _l[_qp]) * c * c;
135 
136  // derivative of elastic free energy density wrt c
137  _dFdc[_qp] = -hist_variable * 2.0 * (1.0 - c) * (1 - _kdamage) + _gc[_qp] / _l[_qp] * c;
138 
139  // 2nd derivative of elastic free energy density wrt c
140  _d2Fdc2[_qp] = hist_variable * 2.0 * (1 - _kdamage) + _gc[_qp] / _l[_qp];
141 
142  // 2nd derivative wrt c and strain = 0.0 if we used the previous step's history varible
143  if (_use_current_hist)
144  _d2Fdcdstrain[_qp] = -stress0pos * 2.0 * (1.0 - c) * (1 - _kdamage);
145 
146  // Used in StressDivergencePFFracTensors off-diagonal Jacobian
147  _dstress_dc[_qp] = -stress0pos * 2.0 * (1.0 - c) * (1 - _kdamage);
148 
149  // Calculate dstress/dstrain
150  // dstress/dstrain = dstress_pos/dstrain_pos * dstrain_pos/dstrain
151  // + dstress_neg/dstrain_neg * dstrain_neg/dstrain
152  // dstrain_pos/dstrain = proj_pos (fourth order projection tensor)
153  // dstrain_neg/dstrain = proj_neg (fourth order projection tensor)
154  // Note that proj_pos + proj_neg = symmetric identify fourth order tensor
155 
156  // d trace(A) / dA
157  RankTwoTensor I(RankTwoTensor::initIdentity);
158  RankFourTensor dtraceAdA = I.outerProduct(I);
159 
160  _Jacobian_mult[_qp] = ((1.0 - c) * (1.0 - c) * (1 - _kdamage) + _kdamage) *
161  (lambda * dtraceAdA * MathUtils::heavyside(etr) + 2 * mu * proj_pos) +
162  (lambda * dtraceAdA * MathUtils::heavyside(-etr) + 2 * mu * proj_neg);
163 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const VariableValue & _c
Coupled order parameter defining the crack.
const Real _kdamage
Small number to avoid non-positive definiteness at or near complete damage.
MaterialProperty< RankTwoTensor > & _stress
const MaterialProperty< Real > & _gc
Material property defining gc parameter, declared elsewhere.
const MaterialProperty< RankTwoTensor > & _mechanical_strain
const MaterialProperty< Real > & _l
Material property defining crack width, declared elsewhere.
MaterialProperty< Real > & _hist
History variable that prevents crack healing, declared in this material.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
MaterialProperty< Real > & _F
Elastic energy and derivatives, declared in this material.
const MaterialProperty< Real > & _hist_old
Old value of history variable.

◆ initQpStatefulProperties()

void ComputeIsotropicLinearElasticPFFractureStress::initQpStatefulProperties ( )
protectedvirtual

Function required to initialize statefull material properties.

Reimplemented from ComputeStressBase.

Definition at line 53 of file ComputeIsotropicLinearElasticPFFractureStress.C.

54 {
55  _hist[_qp] = 0.0;
56 }
MaterialProperty< Real > & _hist
History variable that prevents crack healing, declared in this material.

Member Data Documentation

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

◆ _c

const VariableValue& ComputeIsotropicLinearElasticPFFractureStress::_c
protected

Coupled order parameter defining the crack.

Definition at line 38 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _d2Fdc2

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_d2Fdc2
protected

Definition at line 55 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _d2Fdcdstrain

MaterialProperty<RankTwoTensor>& ComputeIsotropicLinearElasticPFFractureStress::_d2Fdcdstrain
protected

Definition at line 57 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _dFdc

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_dFdc
protected

Definition at line 54 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _dstress_dc

MaterialProperty<RankTwoTensor>& ComputeIsotropicLinearElasticPFFractureStress::_dstress_dc
protected

Definition at line 56 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _elastic_strain

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited

◆ _elasticity_tensor

const MaterialProperty<RankFourTensor>& ComputeStressBase::_elasticity_tensor
protectedinherited

◆ _elasticity_tensor_name

const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited

◆ _extra_stress

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited

Extra stress tensor.

Definition at line 47 of file ComputeStressBase.h.

Referenced by ComputeStressBase::computeQpProperties().

◆ _F

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_F
protected

Elastic energy and derivatives, declared in this material.

Definition at line 53 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _gc

const MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_gc
protected

Material property defining gc parameter, declared elsewhere.

Definition at line 50 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _hist

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_hist
protected

History variable that prevents crack healing, declared in this material.

Definition at line 60 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

◆ _hist_old

const MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_hist_old
protected

Old value of history variable.

Definition at line 63 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _initial_stress_fcn

std::vector<Function *> ComputeStressBase::_initial_stress_fcn
protectedinherited

initial stress components

Definition at line 50 of file ComputeStressBase.h.

◆ _Jacobian_mult

MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited

◆ _kdamage

const Real ComputeIsotropicLinearElasticPFFractureStress::_kdamage
protected

Small number to avoid non-positive definiteness at or near complete damage.

Definition at line 41 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _l

const MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_l
protected

Material property defining crack width, declared elsewhere.

Definition at line 47 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

◆ _use_current_hist

bool ComputeIsotropicLinearElasticPFFractureStress::_use_current_hist
protected

Use current value of history variable.

Definition at line 44 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().


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