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

Phase-field fracture This class computes the stress and energy contribution to fracture Small strain Anisotropic Elastic formulation Stiffness matrix scaled for heterogeneous elasticity property. More...

#include <ComputeLinearElasticPFFractureStress.h>

Inheritance diagram for ComputeLinearElasticPFFractureStress:
[legend]

Public Member Functions

 ComputeLinearElasticPFFractureStress (const InputParameters &parameters)
 

Protected Member Functions

virtual void computeQpStress ()
 
virtual void initQpStatefulProperties () override
 
virtual void computeQpProperties () override
 

Protected Attributes

bool _use_current_hist
 Use current value of history variable. More...
 
const std::string _uncracked_base_name
 Base name of the stress and strain modified to include cracks. More...
 
const VariableValue & _c
 Variable defining the phase field damage parameter. More...
 
const MaterialProperty< Real > & _gc_prop
 Critical energy release rate for fracture. More...
 
const MaterialProperty< Real > & _l
 Characteristic length, controls damage zone thickness. More...
 
const MaterialProperty< Real > & _visco
 
Real _kdamage
 Small number to avoid non-positive definiteness at or near complete damage. More...
 
MaterialProperty< Real > & _F
 
MaterialProperty< Real > & _dFdc
 
MaterialProperty< Real > & _d2Fdc2
 
MaterialProperty< RankTwoTensor > & _d2Fdcdstrain
 
MaterialProperty< RankTwoTensor > & _dstress_dc
 
MaterialProperty< Real > & _hist
 
const MaterialProperty< Real > & _hist_old
 
MaterialProperty< Real > & _kappa
 Property where the value for kappa will be defined. More...
 
MaterialProperty< Real > & _L
 Property where the value for L will be defined. 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 to fracture Small strain Anisotropic Elastic formulation Stiffness matrix scaled for heterogeneous elasticity property.

Definition at line 26 of file ComputeLinearElasticPFFractureStress.h.

Constructor & Destructor Documentation

◆ ComputeLinearElasticPFFractureStress()

ComputeLinearElasticPFFractureStress::ComputeLinearElasticPFFractureStress ( const InputParameters &  parameters)

Definition at line 36 of file ComputeLinearElasticPFFractureStress.C.

38  : ComputeStressBase(parameters),
39  _use_current_hist(getParam<bool>("use_current_history_variable")),
40  _c(coupledValue("c")),
41  _gc_prop(getMaterialProperty<Real>("gc_prop")),
42  _l(getMaterialProperty<Real>("l")),
43  _visco(getMaterialProperty<Real>("visco")),
44  _kdamage(getParam<Real>("kdamage")),
45  _F(declareProperty<Real>(getParam<MaterialPropertyName>("F_name"))),
46  _dFdc(declarePropertyDerivative<Real>(getParam<MaterialPropertyName>("F_name"),
47  getVar("c", 0)->name())),
48  _d2Fdc2(declarePropertyDerivative<Real>(
49  getParam<MaterialPropertyName>("F_name"), getVar("c", 0)->name(), getVar("c", 0)->name())),
50  _d2Fdcdstrain(declareProperty<RankTwoTensor>("d2Fdcdstrain")),
51  _dstress_dc(declarePropertyDerivative<RankTwoTensor>("stress", getVar("c", 0)->name())),
52  _hist(declareProperty<Real>("hist")),
53  _hist_old(getMaterialPropertyOld<Real>("hist")),
54  _kappa(declareProperty<Real>(getParam<MaterialPropertyName>("kappa_name"))),
55  _L(declareProperty<Real>(getParam<MaterialPropertyName>("mobility_name")))
56 {
57 }
const VariableValue & _c
Variable defining the phase field damage parameter.
bool _use_current_hist
Use current value of history variable.
MaterialProperty< Real > & _kappa
Property where the value for kappa will be defined.
const MaterialProperty< Real > & _l
Characteristic length, controls damage zone thickness.
const std::string name
Definition: Setup.h:22
MaterialProperty< Real > & _L
Property where the value for L will be defined.
ComputeStressBase(const InputParameters &parameters)
Real _kdamage
Small number to avoid non-positive definiteness at or near complete damage.
const MaterialProperty< Real > & _gc_prop
Critical energy release rate for fracture.

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 ComputeLinearElasticPFFractureStress::computeQpStress ( )
protectedvirtual

Implements ComputeStressBase.

Definition at line 60 of file ComputeLinearElasticPFFractureStress.C.

61 {
62  const Real c = _c[_qp];
63 
64  // Zero out values when c > 1
65  Real cfactor = 1.0;
66  if (c > 1.0)
67  cfactor = 0.0;
68 
69  // Compute Uncracked stress
70  RankTwoTensor uncracked_stress = _elasticity_tensor[_qp] * _mechanical_strain[_qp];
71 
72  // Create the positive and negative projection tensors
73  RankFourTensor I4sym(RankFourTensor::initIdentitySymmetricFour);
74  std::vector<Real> eigval;
75  RankTwoTensor eigvec;
76  RankFourTensor Ppos = uncracked_stress.positiveProjectionEigenDecomposition(eigval, eigvec);
77  RankFourTensor Pneg = I4sym - Ppos;
78 
79  // Project the positive and negative stresses
80  RankTwoTensor stress0pos = Ppos * uncracked_stress;
81  RankTwoTensor stress0neg = Pneg * uncracked_stress;
82 
83  // Compute the positive and negative elastic energies
84  Real G0_pos = (stress0pos).doubleContraction(_mechanical_strain[_qp]) / 2.0;
85  Real G0_neg = (stress0neg).doubleContraction(_mechanical_strain[_qp]) / 2.0;
86 
87  // Update the history variable
88  if (G0_pos > _hist_old[_qp])
89  _hist[_qp] = G0_pos;
90  else
91  _hist[_qp] = _hist_old[_qp];
92 
93  Real hist_variable = _hist_old[_qp];
95  hist_variable = _hist[_qp];
96 
97  // Compute degradation function and derivatives
98  Real h = cfactor * (1.0 - c) * (1.0 - c) * (1.0 - _kdamage) + _kdamage;
99  Real dhdc = -2.0 * cfactor * (1.0 - c) * (1.0 - _kdamage);
100  Real d2hdc2 = 2.0 * cfactor * (1.0 - _kdamage);
101 
102  // Compute stress and its derivatives
103  _stress[_qp] = stress0pos * h + stress0neg; // equivalent to (Ppos * h + Pneg) * uncracked_stress;
104  _dstress_dc[_qp] = stress0pos * dhdc;
105  _Jacobian_mult[_qp] = (Ppos * h + Pneg) * _elasticity_tensor[_qp];
106 
107  // Compute energy and its derivatives
108  _F[_qp] = hist_variable * h - G0_neg + _gc_prop[_qp] * c * c / (2 * _l[_qp]);
109  _dFdc[_qp] = hist_variable * dhdc + _gc_prop[_qp] * c / _l[_qp];
110  _d2Fdc2[_qp] = hist_variable * d2hdc2 + _gc_prop[_qp] / _l[_qp];
111 
112  // 2nd derivative wrt c and strain = 0.0 if we used the previous step's history varible
113  if (_use_current_hist)
114  _d2Fdcdstrain[_qp] = _dstress_dc[_qp];
115 
116  // Assign L and kappa
117  _kappa[_qp] = _gc_prop[_qp] * _l[_qp];
118  _L[_qp] = 1.0 / (_gc_prop[_qp] * _visco[_qp]);
119 }
MaterialProperty< RankFourTensor > & _Jacobian_mult
derivative of stress w.r.t. strain (_dstress_dstrain)
const VariableValue & _c
Variable defining the phase field damage parameter.
bool _use_current_hist
Use current value of history variable.
MaterialProperty< RankTwoTensor > & _stress
MaterialProperty< Real > & _kappa
Property where the value for kappa will be defined.
const MaterialProperty< Real > & _l
Characteristic length, controls damage zone thickness.
const MaterialProperty< RankTwoTensor > & _mechanical_strain
MaterialProperty< Real > & _L
Property where the value for L will be defined.
Real _kdamage
Small number to avoid non-positive definiteness at or near complete damage.
const MaterialProperty< RankFourTensor > & _elasticity_tensor
const MaterialProperty< Real > & _gc_prop
Critical energy release rate for fracture.

◆ initQpStatefulProperties()

void ComputeStressBase::initQpStatefulProperties ( )
overrideprotectedvirtualinherited

Member Data Documentation

◆ _base_name

const std::string ComputeStressBase::_base_name
protectedinherited

◆ _c

const VariableValue& ComputeLinearElasticPFFractureStress::_c
protected

Variable defining the phase field damage parameter.

Definition at line 41 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _d2Fdc2

MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_d2Fdc2
protected

Definition at line 57 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _d2Fdcdstrain

MaterialProperty<RankTwoTensor>& ComputeLinearElasticPFFractureStress::_d2Fdcdstrain
protected

Definition at line 58 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _dFdc

MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_dFdc
protected

Definition at line 56 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _dstress_dc

MaterialProperty<RankTwoTensor>& ComputeLinearElasticPFFractureStress::_dstress_dc
protected

Definition at line 59 of file ComputeLinearElasticPFFractureStress.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>& ComputeLinearElasticPFFractureStress::_F
protected

Definition at line 55 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _gc_prop

const MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_gc_prop
protected

Critical energy release rate for fracture.

Definition at line 44 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _hist

MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_hist
protected

Definition at line 61 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _hist_old

const MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_hist_old
protected

Definition at line 62 of file ComputeLinearElasticPFFractureStress.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

◆ _kappa

MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_kappa
protected

Property where the value for kappa will be defined.

Definition at line 65 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _kdamage

Real ComputeLinearElasticPFFractureStress::_kdamage
protected

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

Definition at line 53 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _l

const MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_l
protected

Characteristic length, controls damage zone thickness.

Definition at line 47 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _L

MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_L
protected

Property where the value for L will be defined.

Definition at line 68 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _mechanical_strain

const MaterialProperty<RankTwoTensor>& ComputeStressBase::_mechanical_strain
protectedinherited

◆ _stress

MaterialProperty<RankTwoTensor>& ComputeStressBase::_stress
protectedinherited

◆ _uncracked_base_name

const std::string ComputeLinearElasticPFFractureStress::_uncracked_base_name
protected

Base name of the stress and strain modified to include cracks.

Definition at line 38 of file ComputeLinearElasticPFFractureStress.h.

◆ _use_current_hist

bool ComputeLinearElasticPFFractureStress::_use_current_hist
protected

Use current value of history variable.

Definition at line 35 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().

◆ _visco

const MaterialProperty<Real>& ComputeLinearElasticPFFractureStress::_visco
protected

Definition at line 50 of file ComputeLinearElasticPFFractureStress.h.

Referenced by computeQpStress().


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