diff --git a/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.cpp b/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.cpp
index ef3dfaa923cb46ba117bf088dd89676adea7ca53..9b907611461f86bb7b6378727f7126c2d818f253 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.cpp
+++ b/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.cpp
@@ -18,7 +18,7 @@ mlawNonLocalPorousCoupledLaw::mlawNonLocalPorousCoupledLaw(const int num,const d
const double fVinitial, const double lambda0, const J2IsotropicHardening &j2IH, const CLengthLaw &cLLaw,
const double tol, const bool matrixbyPerturbation, const double pert) :
mlawNonLocalPorosity(num, E, nu, rho, fVinitial, j2IH, cLLaw, tol, matrixbyPerturbation, pert),
- _withCrackTransition(false), _withBothYieldSurfaces(false)
+ _withCrackTransition(false), _withBothYieldSurfacesInBulk(false), _withBothYieldSurfacesInInterface(false), _beta(0.)
{
_mlawGrowth = new mlawNonLocalDamageGurson(num, E, nu, rho, q1, q2, q3, fVinitial, j2IH, cLLaw, tol, matrixbyPerturbation, pert);
@@ -34,7 +34,7 @@ mlawNonLocalPorousCoupledLaw::mlawNonLocalPorousCoupledLaw(const int num,const d
const double fVinitial, const double lambda0, const double kappa, const J2IsotropicHardening &j2IH, const CLengthLaw &cLLaw,
const double tol, const bool matrixbyPerturbation, const double pert) :
mlawNonLocalPorosity(num, E, nu, rho, fVinitial, j2IH, cLLaw, tol, matrixbyPerturbation, pert),
- _withCrackTransition(false), _withBothYieldSurfaces(false)
+ _withCrackTransition(false), _withBothYieldSurfacesInBulk(false), _withBothYieldSurfacesInInterface(false), _beta(0.)
{
_mlawGrowth = new mlawNonLocalDamageGurson(num, E, nu, rho, q1, q2, q3, fVinitial, j2IH, cLLaw, tol, matrixbyPerturbation, pert);
@@ -49,7 +49,10 @@ mlawNonLocalPorousCoupledLaw::mlawNonLocalPorousCoupledLaw(const int num,const d
mlawNonLocalPorousCoupledLaw::mlawNonLocalPorousCoupledLaw(const mlawNonLocalPorousCoupledLaw &source) :
mlawNonLocalPorosity(source),
- _withCrackTransition(source._withCrackTransition), _withBothYieldSurfaces(source._withBothYieldSurfaces)
+ _withCrackTransition(source._withCrackTransition), _withBothYieldSurfacesInBulk(source._withBothYieldSurfacesInBulk),
+ _withBothYieldSurfacesInInterface(source._withBothYieldSurfacesInInterface),
+ _beta(source._beta)
+
{
_mlawGrowth = NULL;
if(source._mlawGrowth != NULL) _mlawGrowth = dynamic_cast<mlawNonLocalDamageGurson*>(source._mlawGrowth->clone());
@@ -132,14 +135,15 @@ void mlawNonLocalPorousCoupledLaw::setCfTOffsetMethod(const int MethodNumber){
}
};
-void mlawNonLocalPorousCoupledLaw::setCouplingBehaviour(const bool flag){
- _withBothYieldSurfaces = flag;
- if(_withBothYieldSurfaces){Msg::Info("mlawNonLocalPorousCoupledLaw:: Gurson is also allowed during transition");};
-}
+void mlawNonLocalPorousCoupledLaw::setCouplingBehaviour(const bool flagBulk, const bool flagInterface){
+ _withBothYieldSurfacesInBulk = flagBulk;
+ _withBothYieldSurfacesInInterface = flagInterface;
+};
-void mlawNonLocalPorousCoupledLaw::setCrackTransition(const bool fl){
+void mlawNonLocalPorousCoupledLaw::setCrackTransition(const bool fl, const double beta){
_withCrackTransition = fl;
+ _beta = beta;
};
void mlawNonLocalPorousCoupledLaw::setBalanceEquationType(const int type){
@@ -413,7 +417,6 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
// Get effective traction stress on the interface from tau stress
static SVector3 tractionForce; // traction force vector
double tauEff = 0.; double tauN = 0.; double tauTanSq = 0.; // tranction force components
- double beta = 1.; // ratio => should be moved into law creation as a parameter
static SVector3 Norm;
Norm = q1->getConstRefToReferenceOutwardNormal();
if (Norm.norm() > 0){ // is it usefulll ???
@@ -433,10 +436,10 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
}
// Compute effective stress
if (tauN > 0.){
- tauEff = sqrt(tauN*tauN + beta*tauTanSq);
+ tauEff = sqrt(tauN*tauN + _beta*tauTanSq);
}
else{
- tauEff = sqrt(beta*tauTanSq);
+ tauEff = sqrt(_beta*tauTanSq);
}
}
else{
@@ -456,7 +459,7 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
double yieldThomason = _mlawCoales->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
// Offset management
- if(_CfTOffsetMethod == 2 or _CfTOffsetMethod == 3){
+ if(_CfTOffsetMethod == 2){
// Method 2: Offset is always computed to ensure contuinity if coalescence is forced (offset can be <1)
q1Thom->getRefToCrackOffsetOnCft() = 1. + yieldThomason*R0/(Cft*R);
}
@@ -487,8 +490,6 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
q1Thom->getRefToCrackOffsetOnCft() = 1.;
}
Msg::Info("Criterion satified with a traction force angle =%e rad", atan(sqrt(tauTanSq)/tauN));
-
-
}
else{
// Coalescence never occurs and is not occuring at the next step
@@ -525,7 +526,7 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
// Determine which mode should be used for the next time step.
if(!q0Thom->getCoalescenceOnsetFlag()){
// First time that coalescence is detected: Thomason == ON
- // This case is only appearing with with selective update
+ // This case is only appearing with selective update
q1Thom->getRefToCoalescenceActiveFlag() = true;
q1Thom->getRefToAccelerateRate() = 1.0;
q1Thom->getRefToYieldOffset() = _mlawGrowth->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
@@ -538,8 +539,8 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
q1Thom->getRefToYieldOffset() = yieldGurson;
if(yieldGurson > _tol){
Msg::Info("Change yield surface to Gurson fG =%e",yieldGurson);
- if(_withBothYieldSurfaces){
- q1Thom->getRefToCoalescenceActiveFlag() = false; // Gurson will be used at the next time step
+ if(_withBothYieldSurfacesInInterface){
+ q1Thom->getRefToCoalescenceActiveFlag() = false; // Gurson is forced to be used at the next time step
}
else{
q1Thom->getRefToCoalescenceActiveFlag() = true;
@@ -559,12 +560,6 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
Msg::Info("Change yield surface to Thomason fT = %e",yieldThomason);
q1Thom->getRefToCoalescenceActiveFlag() = true; // Thomason will be used at the next time step
q1Thom->getRefToAccelerateRate() = 1.0;
-
- if (_CfTOffsetMethod == 3){
- q1Thom->getRefToCrackOffsetOnCft() = q0Thom->getCrackOffsetOnCft() + q1Thom->getYieldOffset()*R0/(Cft*R);
- if (q1Thom->getRefToCrackOffsetOnCft() > 2.){Msg::Error("High predicted value of Cft Offset = %e",q1Thom->getRefToCrackOffsetOnCft());};
- }
-
}
else{
q1Thom->getRefToCoalescenceActiveFlag() = false; // Gurson will be used at the next time step
@@ -577,81 +572,114 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
else
{
// In the bulk:
-/*
- // Only growth model is allowed
- q1Thom->getRefToCoalescenceOnsetFlag() = false;
- // Reset onset variables to their default value
- q1Thom->getRefToPorosityAtCoalescenceOnset() = 0.;
- q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = 0.;
- q1Thom->getRefToAspectRatioAtCoalescenceOnset() = 0.;
- q1Thom->getRefToShapeFactorAtCoalescenceOnset() = 0.;
-
- q1Thom->getRefToCoalescenceActiveFlag() = false;
- q1Thom->getRefToCrackOffsetOnCft() = 1.;
- q1Thom->getRefToAccelerateRate() = 1.;
- q1Thom->getRefToYieldOffset() = 0.;
-
-*/
- if (!q0Thom->getCoalescenceOnsetFlag()){
- double yieldThomason = _mlawCoales->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
-
- if (yieldThomason > _tol){
- q1Thom->getRefToCoalescenceCriterion() = 0.;
- // update offset
- if (_CfTOffsetMethod == 1){
- q1Thom->getRefToCrackOffsetOnCft() = 1. + yieldThomason*R0/(Cft*R);
- Msg::Info("coalescence occurs, use Thomason yield surface in next step with offset = %e",q1Thom->getCrackOffsetOnCft());
+ if(_withBothYieldSurfacesInBulk){
+ if (!q0Thom->getCoalescenceOnsetFlag()){
+
+ double yieldThomason = _mlawCoales->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
+ if (yieldThomason > _tol){
+ q1Thom->getRefToCoalescenceCriterion() = 0.;
+
+ q1Thom->getRefToCoalescenceActiveFlag() = true;
+ q1Thom->getRefToCoalescenceOnsetFlag() = true;
+ //Msg::Info("coalescene yield surface is used");
+ // initial coalesence onset value
+ q1Thom->getRefToPorosityAtCoalescenceOnset() = yieldfV;
+ q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q1Thom->getLigamentRatio();
+ q1Thom->getRefToAspectRatioAtCoalescenceOnset() = q1Thom->getAspectRatio();
+ q1Thom->getRefToShapeFactorAtCoalescenceOnset() = q1Thom->getShapeFactor();
+ q1Thom->getRefToAccelerateRate() = 1.;
}
- else if(_CfTOffsetMethod == 0){
+ else{
+ q1Thom->getRefToCoalescenceCriterion() = yieldThomason; // save data
+ q1Thom->getRefToCoalescenceActiveFlag() = false;
+ q1Thom->getRefToCoalescenceOnsetFlag() = false;
+
+ // reset to defaut value
+ q1Thom->getRefToPorosityAtCoalescenceOnset() = 0.;
+ q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = 0.;
+ q1Thom->getRefToAspectRatioAtCoalescenceOnset() = 0.;
+ q1Thom->getRefToShapeFactorAtCoalescenceOnset() = 0.;
+
q1Thom->getRefToCrackOffsetOnCft() = 1.;
+ q1Thom->getRefToAccelerateRate() = 1.;
}
-
- q1Thom->getRefToCoalescenceActiveFlag() = true;
- q1Thom->getRefToCoalescenceOnsetFlag() = true;
- //Msg::Info("coalescene yield surface is used");
- // initial coalesence onset value
- q1Thom->getRefToPorosityAtCoalescenceOnset() = yieldfV;
- q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q1Thom->getLigamentRatio();
- q1Thom->getRefToAspectRatioAtCoalescenceOnset() = q1Thom->getAspectRatio();
- q1Thom->getRefToShapeFactorAtCoalescenceOnset() = q1Thom->getShapeFactor();
- q1Thom->getRefToAccelerateRate() = 1.;
}
else{
- q1Thom->getRefToCoalescenceCriterion() = yieldThomason; // save data
- q1Thom->getRefToCoalescenceActiveFlag() = false;
- q1Thom->getRefToCoalescenceOnsetFlag() = false;
-
- // reset to defaut value
- q1Thom->getRefToPorosityAtCoalescenceOnset() = 0.;
- q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = 0.;
- q1Thom->getRefToAspectRatioAtCoalescenceOnset() = 0.;
- q1Thom->getRefToShapeFactorAtCoalescenceOnset() = 0.;
-
- q1Thom->getRefToCrackOffsetOnCft() = 1.;
- q1Thom->getRefToAccelerateRate() = 1.;
+ // Coalescence has occured at least once
+ q1Thom->getRefToCoalescenceOnsetFlag() = true; // already done....
+ // Update onset value (yield porosity and geometrical parameters)
+ q1Thom->getRefToPorosityAtCoalescenceOnset() = q0Thom->getPorosityAtCoalescenceOnset();
+ q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q0Thom->getLigamentRatioAtCoalescenceOnset();
+ q1Thom->getRefToAspectRatioAtCoalescenceOnset() = q0Thom->getAspectRatioAtCoalescenceOnset();
+ q1Thom->getRefToShapeFactorAtCoalescenceOnset() = q0Thom->getShapeFactorAtCoalescenceOnset();
+
+ // Cft onset - Method 0-1-2: keep it always constant once the coalescence occurs
+ q1Thom->getRefToCrackOffsetOnCft() = q0Thom->getCrackOffsetOnCft();
+ // Determine which mode should be used for the next time step.
+ if(!q0Thom->getCoalescenceOnsetFlag()){
+ // First time that coalescence is detected: Thomason == ON
+ // This case is only appearing with selective update
+ q1Thom->getRefToCoalescenceActiveFlag() = true;
+ q1Thom->getRefToAccelerateRate() = 1.0;
+ q1Thom->getRefToYieldOffset() = _mlawGrowth->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
+ }
+ else{
+ // Not the first time that coalescence appears
+ if (q0Thom->getCoalescenceActiveFlag()){
+ // If Thomason is/was used for this time step : check if Gurson is more restrictive
+ double yieldGurson = _mlawGrowth->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
+ q1Thom->getRefToYieldOffset() = yieldGurson;
+ if(yieldGurson > _tol){
+ Msg::Info("Change yield surface to Gurson fG =%e",yieldGurson);
+ if(_withBothYieldSurfacesInInterface){
+ q1Thom->getRefToCoalescenceActiveFlag() = false; // Gurson is forced to be used at the next time step
+ }
+ else{
+ q1Thom->getRefToCoalescenceActiveFlag() = true;
+ }
+ q1Thom->getRefToAccelerateRate() = 1.0;
+ }
+ else{
+ q1Thom->getRefToCoalescenceActiveFlag() = true; // Thomason will be used at the next time step
+ q1Thom->getRefToAccelerateRate() = 1.0;
+ }
+ }
+ else{
+ // If Gurson was/is used for this time step : check if Thomason is more restrictive
+ double yieldThomason = _mlawCoales->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
+ q1Thom->getRefToYieldOffset() = yieldThomason;
+ if(yieldThomason > _tol){
+ Msg::Info("Change yield surface to Thomason fT = %e",yieldThomason);
+ q1Thom->getRefToCoalescenceActiveFlag() = true; // Thomason will be used at the next time step
+ q1Thom->getRefToAccelerateRate() = 1.0;
+ }
+ else{
+ q1Thom->getRefToCoalescenceActiveFlag() = false; // Gurson will be used at the next time step
+ q1Thom->getRefToAccelerateRate() = 1.0;
+ }
+ }
+ }
}
}
else{
- q1Thom->getRefToCoalescenceCriterion() = 0.;
- q1Thom->getRefToCoalescenceActiveFlag() = true;
- q1Thom->getRefToCoalescenceOnsetFlag() = true;
-
- q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q0Thom->getLigamentRatioAtCoalescenceOnset();
- q1Thom->getRefToAspectRatioAtCoalescenceOnset() = q0Thom->getAspectRatioAtCoalescenceOnset();
- q1Thom->getRefToShapeFactorAtCoalescenceOnset() = q0Thom->getShapeFactorAtCoalescenceOnset();
- q1Thom->getRefToPorosityAtCoalescenceOnset() = q0Thom->getPorosityAtCoalescenceOnset();
+ // Only growth model is allowed
+ q1Thom->getRefToCoalescenceOnsetFlag() = false;
+ // Reset onset variables to their default value
+ q1Thom->getRefToPorosityAtCoalescenceOnset() = 0.;
+ q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = 0.;
+ q1Thom->getRefToAspectRatioAtCoalescenceOnset() = 0.;
+ q1Thom->getRefToShapeFactorAtCoalescenceOnset() = 0.;
- q1Thom->getRefToCrackOffsetOnCft() = q0Thom->getCrackOffsetOnCft();
- q1Thom->getRefToAccelerateRate() = q0Thom->getAccelerateRate();
+ q1Thom->getRefToCoalescenceActiveFlag() = false;
+ q1Thom->getRefToCrackOffsetOnCft() = 1.;
+ q1Thom->getRefToAccelerateRate() = 1.;
+ q1Thom->getRefToYieldOffset() = 0.;
}
}
-
-
-
- // In the bulk: end
-
}
else{
+ // Without crack transition:
+
// If crack transition is not used: no distinction between bulk and interface ipvs
// if Thomason yield surface is first used, it always is used
// onset and active are the same
diff --git a/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.h b/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.h
index 4508914c6f27f1e0a3b39d5d8fb473039099494d..1c14be44f9813cf19f9d913a9fd4f6226fac32ce 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.h
+++ b/NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.h
@@ -26,7 +26,9 @@ class mlawNonLocalPorousCoupledLaw : public mlawNonLocalPorosity
// options
bool _withCrackTransition;
- bool _withBothYieldSurfaces;
+ bool _withBothYieldSurfacesInBulk;
+ bool _withBothYieldSurfacesInInterface;
+ double _beta;
public:
@@ -48,12 +50,12 @@ class mlawNonLocalPorousCoupledLaw : public mlawNonLocalPorosity
virtual void setScatterredInitialPorosity(double f0min, double f0max);
virtual void setCoalescenceLaw(const CoalescenceLaw& added_coalsLaw);
virtual void setYieldSurfaceExponent(const double newN);
- virtual void setCrackTransition(const bool fl);
+ virtual void setCrackTransition(const bool fl, const double beta=0.);
virtual void setBalanceEquationType(const int type);
virtual void setRoundedYieldSurfaceMethod(const int method);
virtual void setOnsetTriaxialityForRoundedYieldSurface(const double newTc);
virtual void setCfTOffsetMethod(const int MethodNumber);
- virtual void setCouplingBehaviour(const bool flag);
+ virtual void setCouplingBehaviour(const bool flagBulk, const bool flagInterface=false);
virtual void setShearPorosityGrowthFactor(const double k);
virtual void setLocalRegularizedFunction(const scalarFunction& fct);
diff --git a/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp b/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp
index 229deb0f3cc7037b57e48b2b010418b05b6189fe..ff2f60aff687428a8aca84277df1e5fbd7c0bc51 100644
--- a/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp
+++ b/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp
@@ -1332,12 +1332,12 @@ void NonLocalPorousCoupledDG3DMaterialLaw::setCfTOffsetMethod(const int method){
static_cast<mlawNonLocalPorousCoupledLaw*>(_nldPorous)->setCfTOffsetMethod(method);
};
-void NonLocalPorousCoupledDG3DMaterialLaw::setCouplingBehaviour(const bool flag){
- static_cast<mlawNonLocalPorousCoupledLaw*>(_nldPorous)->setCouplingBehaviour(flag);
+void NonLocalPorousCoupledDG3DMaterialLaw::setCouplingBehaviour(const bool flagBulk, const bool flagInterface){
+ static_cast<mlawNonLocalPorousCoupledLaw*>(_nldPorous)->setCouplingBehaviour(flagBulk, flagInterface);
};
-void NonLocalPorousCoupledDG3DMaterialLaw::setCrackTransition(const bool fl){
- static_cast<mlawNonLocalPorousCoupledLaw*>(_nldPorous)->setCrackTransition(fl);
+void NonLocalPorousCoupledDG3DMaterialLaw::setCrackTransition(const bool fl, const double beta){
+ static_cast<mlawNonLocalPorousCoupledLaw*>(_nldPorous)->setCrackTransition(fl,beta);
}
void NonLocalPorousCoupledDG3DMaterialLaw::setBalanceEquationType(const int type){
diff --git a/dG3D/src/nonLocalDamageDG3DMaterialLaw.h b/dG3D/src/nonLocalDamageDG3DMaterialLaw.h
index 320f12748f86fc62a73b199de71e1becf0ceac5c..4eeb91d885b5a00b1514e337a2f113bda80bf831 100644
--- a/dG3D/src/nonLocalDamageDG3DMaterialLaw.h
+++ b/dG3D/src/nonLocalDamageDG3DMaterialLaw.h
@@ -443,12 +443,12 @@ class NonLocalPorousCoupledDG3DMaterialLaw : public NonLocalPorousDuctileDamageD
const double tol=1.e-8, const bool matrixbyPerturbation = false, const double pert = 1e-8);
void setYieldSurfaceExponent(const double newN);
- void setCrackTransition(const bool fl);
+ void setCrackTransition(const bool fl, const double beta=0.);
void setBalanceEquationType(const int type);
void setRoundedYieldSurfaceMethod(const int method);
void setOnsetTriaxialityForRoundedYieldSurface(const double newTc);
void setCfTOffsetMethod(const int method);
- void setCouplingBehaviour(const bool flag);
+ void setCouplingBehaviour(const bool flagBulk, const bool flagInterface=false);
#ifndef SWIG
NonLocalPorousCoupledDG3DMaterialLaw(const NonLocalPorousCoupledDG3DMaterialLaw &source);
virtual ~NonLocalPorousCoupledDG3DMaterialLaw();