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
 
void addQpInitialStress ()
 InitialStress Deprecation: remove this method. More...
 

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...
 
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 > & _d2Fdcdstrain
 
MaterialProperty< RankTwoTensor > & _dstress_dc
 
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...
 
const bool _store_stress_old
 Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping. More...
 
const bool _initial_stress_provided
 Whether initial stress was provided. InitialStress Deprecation: remove this. More...
 
MaterialProperty< RankTwoTensor > * _initial_stress
 Initial stress, if provided. InitialStress Deprecation: remove this. More...
 
const MaterialProperty< RankTwoTensor > * _initial_stress_old
 Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this. 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 17 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Constructor & Destructor Documentation

ComputeIsotropicLinearElasticPFFractureStress::ComputeIsotropicLinearElasticPFFractureStress ( const InputParameters &  parameters)

Definition at line 23 of file ComputeIsotropicLinearElasticPFFractureStress.C.

25  : ComputeStressBase(parameters),
26  _c(coupledValue("c")),
27  _kdamage(getParam<Real>("kdamage")),
28  _l(getMaterialProperty<Real>("l")),
29  _gc(getMaterialProperty<Real>("gc_prop")),
30  _F(declareProperty<Real>(getParam<MaterialPropertyName>("F_name"))),
31  _dFdc(declarePropertyDerivative<Real>(getParam<MaterialPropertyName>("F_name"),
32  getVar("c", 0)->name())),
33  _d2Fdc2(declarePropertyDerivative<Real>(
34  getParam<MaterialPropertyName>("F_name"), getVar("c", 0)->name(), getVar("c", 0)->name())),
35  _d2Fdcdstrain(declareProperty<RankTwoTensor>("d2Fdcdstrain")),
37  declarePropertyDerivative<RankTwoTensor>(_base_name + "stress", getVar("c", 0)->name())),
38  _hist(declareProperty<Real>("hist")),
39  _hist_old(getMaterialPropertyOld<Real>("hist"))
40 {
41 }
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 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

void ComputeStressBase::addQpInitialStress ( )
protectedinherited

InitialStress Deprecation: remove this method.

Adds initial stress, if it is provided, to _stress[_qp]. This function should NOT be used if you calculate stress using

stress = stress_old + elasticity * strain_increment

because stress_old will already include initial stress. However this function SHOULD be used if your code uses

stress = elasticity * (elastic_strain_old + strain_increment) or stress = elasticity * elastic_strain

since in these cases the elastic_strain and elastic_strain_old will not include any contribution from initial stress.

Definition at line 112 of file ComputeStressBase.C.

Referenced by ComputeLinearElasticStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeStressBase::initQpStatefulProperties(), ComputeMultipleInelasticStress::updateQpState(), and ComputeMultipleInelasticStress::updateQpStateSingleModel().

113 {
115  _stress[_qp] += (*_initial_stress)[_qp];
116 }
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
void ComputeStressBase::computeQpProperties ( )
overrideprotectedvirtualinherited

Reimplemented in ComputeBirchMurnaghanEquationOfStress, and ComputeStressEosBase.

Definition at line 99 of file ComputeStressBase.C.

100 {
101  // InitialStress Deprecation: remove the following 2 lines
103  (*_initial_stress)[_qp] = (*_initial_stress_old)[_qp];
104 
105  computeQpStress();
106 
107  // Add in extra stress
108  _stress[_qp] += _extra_stress[_qp];
109 }
virtual void computeQpStress()=0
MaterialProperty< RankTwoTensor > & _stress
const bool _initial_stress_provided
Whether initial stress was provided. InitialStress Deprecation: remove this.
const MaterialProperty< RankTwoTensor > & _extra_stress
Extra stress tensor.
void ComputeIsotropicLinearElasticPFFractureStress::computeQpStress ( )
protectedvirtual

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

Implements ComputeStressBase.

Definition at line 50 of file ComputeIsotropicLinearElasticPFFractureStress.C.

51 {
52  const Real c = _c[_qp];
53 
54  // Compute cfactor to have c > 1 behave the same as c = 1
55  Real cfactor = 0.0;
56  if (c < 1.0)
57  cfactor = 1.0;
58 
59  // Isotropic elasticity is assumed and should be enforced
60  const Real lambda = _elasticity_tensor[_qp](0, 0, 1, 1);
61  const Real mu = _elasticity_tensor[_qp](0, 1, 0, 1);
62 
63  // Compute eigenvectors and eigenvalues of mechanical strain
64  RankTwoTensor eigvec;
65  std::vector<Real> eigval(LIBMESH_DIM);
66  _mechanical_strain[_qp].symmetricEigenvaluesEigenvectors(eigval, eigvec);
67 
68  // Calculate tensors of outerproduct of eigen vectors
69  std::vector<RankTwoTensor> etens(LIBMESH_DIM);
70 
71  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
72  for (unsigned int j = 0; j < LIBMESH_DIM; ++j)
73  for (unsigned int k = 0; k < LIBMESH_DIM; ++k)
74  etens[i](j, k) = eigvec(j, i) * eigvec(k, i);
75 
76  // Separate out positive and negative eigen values
77  std::vector<Real> epos(LIBMESH_DIM), eneg(LIBMESH_DIM);
78  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
79  {
80  epos[i] = (std::abs(eigval[i]) + eigval[i]) / 2.0;
81  eneg[i] = (std::abs(eigval[i]) - eigval[i]) / 2.0;
82  }
83 
84  // Seprate positive and negative sums of all eigenvalues
85  Real etr = 0.0;
86  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
87  etr += eigval[i];
88 
89  const Real etrpos = (std::abs(etr) + etr) / 2.0;
90  const Real etrneg = (std::abs(etr) - etr) / 2.0;
91 
92  // Calculate the tensile (postive) and compressive (negative) parts of stress
93  RankTwoTensor stress0pos, stress0neg;
94  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
95  {
96  stress0pos += etens[i] * (lambda * etrpos + 2.0 * mu * epos[i]);
97  stress0neg += etens[i] * (lambda * etrneg + 2.0 * mu * eneg[i]);
98  }
99 
100  // sum squares of epos and eneg
101  Real pval(0.0), nval(0.0);
102  for (unsigned int i = 0; i < LIBMESH_DIM; ++i)
103  {
104  pval += epos[i] * epos[i];
105  nval += eneg[i] * eneg[i];
106  }
107 
108  // Energy with positive principal strains
109  const Real G0_pos = lambda * etrpos * etrpos / 2.0 + mu * pval;
110  const Real G0_neg = lambda * etrneg * etrneg / 2.0 + mu * nval;
111 
112  // Assign history variable and derivative
113  if (G0_pos > _hist_old[_qp])
114  _hist[_qp] = G0_pos;
115  else
116  _hist[_qp] = _hist_old[_qp];
117 
118  // Damage associated with positive component of stress
119  _stress[_qp] =
120  stress0pos * cfactor * ((1.0 - c) * (1.0 - c) * (1 - _kdamage) + _kdamage) - stress0neg;
121 
122  // Elastic free energy density
123  _F[_qp] = _hist[_qp] * cfactor * ((1.0 - c) * (1.0 - c) * (1 - _kdamage) + _kdamage) - G0_neg +
124  _gc[_qp] / (2 * _l[_qp]) * c * c;
125 
126  // derivative of elastic free energy density wrt c
127  _dFdc[_qp] = -_hist[_qp] * 2.0 * cfactor * (1.0 - c) * (1 - _kdamage) + _gc[_qp] / _l[_qp] * c;
128 
129  // 2nd derivative of elastic free energy density wrt c
130  _d2Fdc2[_qp] = _hist[_qp] * 2.0 * cfactor * (1 - _kdamage) + _gc[_qp] / _l[_qp];
131 
132  // 2nd derivative wrt c and strain. Note that I am ignoring the history variable here, but this
133  // approach gets the fastest convergence
134  _d2Fdcdstrain[_qp] = -stress0pos * cfactor * (1.0 - c) * (1 - _kdamage);
135 
136  // Used in StressDivergencePFFracTensors off-diagonal Jacobian
137  _dstress_dc[_qp] = -stress0pos * 2.0 * cfactor * (1.0 - c) * (1 - _kdamage);
138 
139  // Estimate Jacobian mult, exact when c = 0 and VERY APPROXIMATE when c = 1
141 }
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.
void ComputeIsotropicLinearElasticPFFractureStress::initQpStatefulProperties ( )
protectedvirtual

Function required to initialize statefull material properties.

Reimplemented from ComputeStressBase.

Definition at line 44 of file ComputeIsotropicLinearElasticPFFractureStress.C.

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

Member Data Documentation

const std::string ComputeStressBase::_base_name
protectedinherited
const VariableValue& ComputeIsotropicLinearElasticPFFractureStress::_c
protected

Coupled order parameter defining the crack.

Definition at line 30 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_d2Fdc2
protected

Definition at line 44 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<RankTwoTensor>& ComputeIsotropicLinearElasticPFFractureStress::_d2Fdcdstrain
protected

Definition at line 45 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_dFdc
protected

Definition at line 43 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<RankTwoTensor>& ComputeIsotropicLinearElasticPFFractureStress::_dstress_dc
protected

Definition at line 46 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<RankTwoTensor>& ComputeStressBase::_elastic_strain
protectedinherited
const MaterialProperty<RankFourTensor>& ComputeStressBase::_elasticity_tensor
protectedinherited
const std::string ComputeStressBase::_elasticity_tensor_name
protectedinherited
const MaterialProperty<RankTwoTensor>& ComputeStressBase::_extra_stress
protectedinherited
MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_F
protected

Elastic energy and derivatives, declared in this material.

Definition at line 42 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

const MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_gc
protected

Material property defining gc parameter, declared elsewhere.

Definition at line 39 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_hist
protected

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

Definition at line 49 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress(), and initQpStatefulProperties().

const MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_hist_old
protected

Old value of history variable.

Definition at line 52 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress
protectedinherited

Initial stress, if provided. InitialStress Deprecation: remove this.

Definition at line 74 of file ComputeStressBase.h.

Referenced by ComputeStressBase::initQpStatefulProperties().

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

initial stress components

Definition at line 62 of file ComputeStressBase.h.

Referenced by ComputeStressBase::ComputeStressBase(), and ComputeStressBase::initQpStatefulProperties().

const MaterialProperty<RankTwoTensor>* ComputeStressBase::_initial_stress_old
protectedinherited

Old value of initial stress, which is needed to correctly implement finite-strain rotations. InitialStress Deprecation: remove this.

Definition at line 77 of file ComputeStressBase.h.

const bool ComputeStressBase::_initial_stress_provided
protectedinherited
MaterialProperty<RankFourTensor>& ComputeStressBase::_Jacobian_mult
protectedinherited
const Real ComputeIsotropicLinearElasticPFFractureStress::_kdamage
protected

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

Definition at line 33 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

const MaterialProperty<Real>& ComputeIsotropicLinearElasticPFFractureStress::_l
protected

Material property defining crack width, declared elsewhere.

Definition at line 36 of file ComputeIsotropicLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited
const bool ComputeStressBase::_store_stress_old
protectedinherited

Parameter which decides whether to store old stress. This is required for HHT time integration and Rayleigh damping.

Definition at line 68 of file ComputeStressBase.h.

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

Definition at line 53 of file ComputeStressBase.h.

Referenced by ComputeStressBase::addQpInitialStress(), ComputeMultipleInelasticCosseratStress::computeAdmissibleState(), ComputeMultipleInelasticStress::computeAdmissibleState(), ComputeStressBase::computeQpProperties(), ComputeStressEosBase::computeQpProperties(), ComputeBirchMurnaghanEquationOfStress::computeQpProperties(), ComputeLinearElasticStress::computeQpStress(), ComputeStrainIncrementBasedStress::computeQpStress(), ComputeLinearElasticPFFractureStress::computeQpStress(), ComputeCosseratLinearElasticStress::computeQpStress(), ComputeFiniteStrainElasticStress::computeQpStress(), computeQpStress(), ComputeFiniteStrainElasticStressBirchMurnaghan::computeQpStress(), FiniteStrainPlasticMaterial::computeQpStress(), ComputeMultiPlasticityStress::computeQpStress(), ComputeSmearedCrackingStress::computeQpStress(), ComputeLinearViscoelasticStress::computeQpStress(), ComputeMultipleInelasticStress::computeQpStressIntermediateConfiguration(), ComputeSmearedCrackingStress::crackingStressRotation(), ComputeMultipleInelasticStress::finiteStrainRotation(), ComputeStressBase::initQpStatefulProperties(), FiniteStrainCrystalPlasticity::initQpStatefulProperties(), FiniteStrainUObasedCP::initQpStatefulProperties(), FiniteStrainHyperElasticViscoPlastic::initQpStatefulProperties(), ComputeMultiPlasticityStress::postReturnMap(), FiniteStrainUObasedCP::postSolveQp(), FiniteStrainHyperElasticViscoPlastic::postSolveQp(), FiniteStrainCrystalPlasticity::postSolveQp(), ComputeMultipleInelasticStress::updateQpState(), ComputeMultipleInelasticStress::updateQpStateSingleModel(), and LinearIsoElasticPFDamage::updateVar().


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