From 6c336965c60a0c1c94c7b536093a7d2b77a42020 Mon Sep 17 00:00:00 2001 From: Antoine Jacquey Date: Tue, 5 Mar 2019 09:46:13 +0100 Subject: [PATCH] Capped yield working. --- include/materials/LynxDamageDeformation.h | 4 ++ include/materials/LynxDeformationBase.h | 1 + src/materials/LynxDamageDeformation.C | 61 ++++++++++++++++++++++- src/materials/LynxDeformationBase.C | 6 +++ 4 files changed, 70 insertions(+), 2 deletions(-) diff --git a/include/materials/LynxDamageDeformation.h b/include/materials/LynxDamageDeformation.h index 374aa46..61d0e0c 100644 --- a/include/materials/LynxDamageDeformation.h +++ b/include/materials/LynxDamageDeformation.h @@ -36,8 +36,10 @@ public: protected: virtual void initQpStatefulProperties() override; + virtual void initializeQpDeformation() override; virtual void elasticModuli() override; virtual void plasticCorrection(Real & pressure, RankTwoTensor & stress_dev) override; + virtual void damageCorrection() override; virtual Real computePlasticityYield(const Real & pressure, const Real & eqv_stress); virtual Real plasticIncrement(const Real & /*pressure*/, const Real & eqv_stress); virtual Real computeConvexPlasticityYield2(const Real & pressure, const Real & eqv_stress); @@ -45,6 +47,7 @@ protected: virtual Real convexPlasticIncrement(Real & vol_plastic_incr, Real & eqv_plastic_incr); virtual void computeDamageProperties(const Real & pressure, const Real & eqv_stress); virtual void updateDamageParameters(); + virtual void initializeDamageParameters(); virtual void updateDamageConvexParameters(const Real & pressure, const Real & eqv_stress); virtual Real convexReferencePressure(); virtual Real dConvexPlasticYield2(const Real & rho); @@ -52,6 +55,7 @@ protected: virtual Real dConvexPlasticYield2_dq(const Real & pressure, const Real & eqv_stress); virtual Real getConvexProjection(const Real & x1, const Real & x2); virtual Real strainRatio(const RankTwoTensor & elastic_strain); + virtual RankTwoTensor rotatedElasticStrain(const RankTwoTensor & elastic_strain); // Coupled variables bool _coupled_dam; diff --git a/include/materials/LynxDeformationBase.h b/include/materials/LynxDeformationBase.h index 9cda22d..38e9f74 100644 --- a/include/materials/LynxDeformationBase.h +++ b/include/materials/LynxDeformationBase.h @@ -55,6 +55,7 @@ protected: virtual void deviatoricDeformation(const Real & pressure, RankTwoTensor & stress_dev); virtual void reformStressTensor(const Real & pressure, const RankTwoTensor & stress_dev); virtual void plasticCorrection(Real & pressure, RankTwoTensor & stress_dev) = 0; + virtual void damageCorrection(); virtual void tangentOperator(); virtual void plasticTangentOperator(const RankTwoTensor & flow_direction, const RankFourTensor & flow_direction_dyad); diff --git a/src/materials/LynxDamageDeformation.C b/src/materials/LynxDamageDeformation.C index e59521c..ee8fb60 100644 --- a/src/materials/LynxDamageDeformation.C +++ b/src/materials/LynxDamageDeformation.C @@ -109,14 +109,32 @@ LynxDamageDeformation::initQpStatefulProperties() LynxDeformationBase::initQpStatefulProperties(); } +void +LynxDamageDeformation::initializeQpDeformation() +{ + // Need to get _xi0 + initializeDamageParameters(); + + LynxDeformationBase::initializeQpDeformation(); +} + void LynxDamageDeformation::elasticModuli() { + // Strain ratio + const Real xi = strainRatio(_elastic_strain_old[_qp]); + // Bulk modulus - _K[_qp] = _has_bulk_modulus ? averageProperty(_bulk_modulus) : 0.0; + _K[_qp] = _has_bulk_modulus ? averageProperty(_bulk_modulus) - + 2.0 / 3.0 * _damage_old[_qp] * _damage_plasticity->_gamma * + (xi / 2.0 - _damage_plasticity->_xi0) + : 0.0; // Shear modulus - _G[_qp] = _has_shear_modulus ? averageProperty(_shear_modulus) : 0.0; + _G[_qp] = _has_shear_modulus + ? averageProperty(_shear_modulus) - _damage_old[_qp] * _damage_plasticity->_gamma * + (xi / 2.0 - _damage_plasticity->_xi0) + : 0.0; } void @@ -169,6 +187,18 @@ LynxDamageDeformation::plasticCorrection(Real & pressure, RankTwoTensor & stress rho_0; } } + + // Update elastic strain + _elastic_strain[_qp] = spinRotation(_elastic_strain_old[_qp].deviatoric()); + _elastic_strain[_qp].addIa(_elastic_strain_old[_qp].trace() / 3.0); + + _elastic_strain[_qp] += + _strain_increment[_qp] - _viscous_strain_incr[_qp] - _plastic_strain_incr[_qp]; +} + +void +LynxDamageDeformation::damageCorrection() +{ } Real @@ -330,6 +360,23 @@ LynxDamageDeformation::updateDamageParameters() averageProperty(_one_on_damage_eta)); } +void +LynxDamageDeformation::initializeDamageParameters() +{ + Real sin_phi = std::sin(averageProperty(_friction_angle) * libMesh::pi / 180.0); + Real xi0 = -std::sqrt(3.0) / + std::sqrt(1.0 + 1.5 * Utility::pow<2>(averageProperty(_bulk_modulus) / + averageProperty(_shear_modulus) * sin_phi)); + Real porous_coeff = + averageProperty(_porous_coeff) + averageProperty(_porous_coeff_linear) * _porosity[_qp]; + Real p_cr = _K[_qp] * (std::sqrt(3.0) + xi0) / (3.0 * porous_coeff); + _damage_plasticity->fill(xi0, + averageProperty(_damage_modulus), + p_cr, + averageProperty(_one_on_plastic_eta), + averageProperty(_one_on_damage_eta)); +} + void LynxDamageDeformation::updateDamageConvexParameters(const Real & pressure, const Real & eqv_stress) { @@ -519,3 +566,13 @@ LynxDamageDeformation::strainRatio(const RankTwoTensor & elastic_strain) else return -std::sqrt(3.0); } + +RankTwoTensor +LynxDamageDeformation::rotatedElasticStrain(const RankTwoTensor & elastic_strain) +{ + const Real strain_v = elastic_strain.trace(); + RankTwoTensor strain_rot = spinRotation(elastic_strain.deviatoric()); + strain_rot.addIa(strain_v / 3.0); + + return strain_rot; +} diff --git a/src/materials/LynxDeformationBase.C b/src/materials/LynxDeformationBase.C index b68d60c..d33fdb2 100644 --- a/src/materials/LynxDeformationBase.C +++ b/src/materials/LynxDeformationBase.C @@ -340,6 +340,7 @@ LynxDeformationBase::computeQpDeformation() reformStressTensor(pressure, stress_dev); // Additional correction + damageCorrection(); // Update tangent operator modulus if (_fe_problem.currentlyComputingJacobian()) @@ -426,6 +427,11 @@ LynxDeformationBase::reformStressTensor(const Real & pressure, const RankTwoTens _stress[_qp].addIa(-pressure); } +void +LynxDeformationBase::damageCorrection() +{ +} + void LynxDeformationBase::tangentOperator() { -- GitLab