diff --git a/NonLinearSolver/materialLaw/mlawNonLinearTVE.cpp b/NonLinearSolver/materialLaw/mlawNonLinearTVE.cpp
index 0df08f388cb7828ddb739681823d4a16761a5f9b..a814172f69a60454e4fd5a0d20b2e3abf379045f 100644
--- a/NonLinearSolver/materialLaw/mlawNonLinearTVE.cpp
+++ b/NonLinearSolver/materialLaw/mlawNonLinearTVE.cpp
@@ -582,6 +582,28 @@ double mlawNonLinearTVE::freeEnergyMechanical(const IPNonLinearTVE& q, const dou
return Psy_mech;
}
+void mlawNonLinearTVE::calculateMullinsEffectScaler(const IPNonLinearTVE* q0, IPNonLinearTVE *q1, const double T) const{
+
+ // update eta - Mullins Effect Scalar
+ if (_mullinsEffect != NULL && q1->_ipMullinsEffect != NULL){
+ _mullinsEffect->mullinsEffectScaling(q0->_elasticEnergy, *q0->_ipMullinsEffect, q1->_elasticEnergy,*q1->_ipMullinsEffect);
+ }
+
+ // get eta - Mullins Effect Scalar (at unloading intitiation-> eta = 1.; unloading end-> eta = eta_min; reloading -> eta = eta_min;)
+ double eta = 1., DetaDpsi = 0., DDetaDpsipsi = 0., DetaDT = 0., DDetaDTT = 0.;
+ if (q1->_ipMullinsEffect != NULL){
+ eta = q1->_ipMullinsEffect->getEta();
+ DetaDpsi = q1->_ipMullinsEffect->getDetaDpsi();
+ DetaDT = q1->_ipMullinsEffect->getDetaDT();
+ DDetaDTT = q1->_ipMullinsEffect->getDDetaDTT();
+ }
+
+ // update Bulk/Shear PropertyScaleFactor
+ q1->getRefToElasticBulkPropertyScaleFactor() = eta;
+ q1->getRefToElasticShearPropertyScaleFactor() = eta;
+
+};
+
void mlawNonLinearTVE::corKirInfinity(const IPNonLinearTVE *q1, const STensor3& devEe, const double& trEe, const double T1, STensor3& CorKirDevInf, double& CorKirTrInf) const{
double KT, GT, AlphaT;
getK(q1,KT,T1); getG(q1,GT,T1); getAlpha(AlphaT,T1);
@@ -685,10 +707,10 @@ void mlawNonLinearTVE::getTVEdCorKirDT(const IPNonLinearTVE *q0, IPNonLinearTVE
if (_TemFuncOpt == 0){ // Use Shift Factor
getKi(q1,KiT,T); getGi(q1,GiT,T);
- DKiDT[i] = -_Ki[i]/(_ki[i]) * Ddt_shiftDT*expdtkby2; // Removed _ki/2 and _gi/2
- DDKiDT[i] = -_Ki[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2))*expdtkby2;
- DGiDT[i] = -_Gi[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
- DDGiDT[i] = -_Gi[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
+ DKiDT[i] = -KiT[i]/(_ki[i]) * Ddt_shiftDT*expdtkby2; // Removed _ki/2 and _gi/2
+ DDKiDT[i] = -KiT[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2))*expdtkby2;
+ DGiDT[i] = -GiT[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
+ DDGiDT[i] = -GiT[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
}
for (int j=0; j<3; j++){
@@ -838,10 +860,10 @@ void mlawNonLinearTVE::getMechSourceTVE(const IPNonLinearTVE *q0, IPNonLinearTVE
if (_TemFuncOpt == 0){ // Use Shift Factor
getKi(q1,KiT,T); getGi(q1,GiT,T);
- DKiDT[i] = -_Ki[i]/(_ki[i]) * Ddt_shiftDT*expdtkby2; // Removed _ki/2 and _gi/2
- DDKiDT[i] = -_Ki[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2))*expdtkby2;
- DGiDT[i] = -_Gi[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
- DDGiDT[i] = -_Gi[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
+ DKiDT[i] = -KiT[i]/(_ki[i]) * Ddt_shiftDT*expdtkby2; // Removed _ki/2 and _gi/2
+ DDKiDT[i] = -KiT[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2))*expdtkby2;
+ DGiDT[i] = -GiT[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
+ DDGiDT[i] = -GiT[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
}
for (int j=0; j<3; j++){
@@ -891,9 +913,9 @@ void mlawNonLinearTVE::getMechSourceTVE(const IPNonLinearTVE *q0, IPNonLinearTVE
STensorOperation::zero(DGammai); STensorOperation::zero(DDGammaiDT);
DQiDT -= DOiDT; DDQiDT -= DDOiDT; DDGammaiDT -= DDOGammaiDT;
for (int j=0; j<3; j++){
- Qi(j,j) += (_Ki[i]*trEe - trOi[i]);
+ Qi(j,j) += (KiT[i]*trEe - trOi[i]);
for (int k=0; k<3; k++){
- Qi(j,k) += (2*_Gi[i]*devEe(j,k) - devOi[i](j,k));
+ Qi(j,k) += (2*GiT[i]*devEe(j,k) - devOi[i](j,k));
}
}
DGammai += DEe;
@@ -906,8 +928,8 @@ void mlawNonLinearTVE::getMechSourceTVE(const IPNonLinearTVE *q0, IPNonLinearTVE
static STensor43 temp_DDQiDTDE;
static STensor43 temp_DDgammiDtDE;
- double Ki_term = (_Ki[i]*(1-expdtkby2)+T*DKiDT[i]);
- double Gi_term = (_Gi[i]*(1-expdtgby2)+T*DGiDT[i]);
+ double Ki_term = (KiT[i]*(1-expdtkby2)+T*DKiDT[i]);
+ double Gi_term = (GiT[i]*(1-expdtgby2)+T*DGiDT[i]);
isotropicHookTensor( Ki_term, Gi_term , temp_DDQiDTDE );
double Ki_term_2 = 1/3*expdtkby2;
@@ -1007,6 +1029,11 @@ void mlawNonLinearTVE::ThermoViscoElasticPredictor(const STensor3& Ee, const STe
double& DKe, double& DGe, double& DKde, double& DGde,
double& KTsum, double& GTsum, double& DKDTsum, double& DGDTsum,
const double T0, const double T1, const bool calculateStress) const{
+
+ // Mullins Effect Scaler -> should be before any calls to get functions
+ if (_mullinsEffect != NULL && q1->_ipMullinsEffect != NULL){
+ mlawNonLinearTVE::calculateMullinsEffectScaler(q0, q1, T1);
+ }
// Set Temperature to calculate temperature updated properties
double T_set = setTemp(T0,T1,1); // 1->T1, 2->T_mid, 3->T_midshift
@@ -1078,8 +1105,8 @@ void mlawNonLinearTVE::ThermoViscoElasticPredictor(const STensor3& Ee, const STe
// Maxwell
if (_TemFuncOpt == 0){ // Use Shift Factor
getGi(q1,GiT,T_set); // Removed 2*_gi[i]
- DGiDT[i] = -_Gi[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
- DDGiDT[i] = -_Gi[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
+ DGiDT[i] = -GiT[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
+ DDGiDT[i] = -GiT[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
}
// For Single Convolution
@@ -1165,8 +1192,8 @@ void mlawNonLinearTVE::ThermoViscoElasticPredictor(const STensor3& Ee, const STe
// Maxwell
if (_TemFuncOpt == 0){ // Use Shift Factor
getKi(q1,KiT,T_set); getAlphai(AlphaiT,T_set); // Removed 2*_ki[i]
- DKiDT[i] = -_Ki[i]/(_ki[i]) * Ddt_shiftDT *expdtkby2;
- DDKiDT[i] = -_Ki[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2)) *expdtkby2;
+ DKiDT[i] = -KiT[i]/(_ki[i]) * Ddt_shiftDT *expdtkby2;
+ DDKiDT[i] = -KiT[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2)) *expdtkby2;
}
// Single Convolution
@@ -1614,10 +1641,10 @@ void mlawNonLinearTVE::predictorCorrector_ThermoViscoElastic(
if (_TemFuncOpt == 0){ // Use Shift Factor
getKi(q1,KiT,T_set); getGi(q1,GiT,T_set); getAlphai(AlphaiT,T_set);
- DKiDT[i] = -_Ki[i]/(_ki[i]) * Ddt_shiftDT*expdtkby2; // Removed _ki/2 and _gi/2
- DDKiDT[i] = -_Ki[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2))*expdtkby2;
- DGiDT[i] = -_Gi[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
- DDGiDT[i] = -_Gi[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
+ DKiDT[i] = -KiT[i]/(_ki[i]) * Ddt_shiftDT*expdtkby2; // Removed _ki/2 and _gi/2
+ DDKiDT[i] = -KiT[i]/(_ki[i]) * (DDdt_shiftDT - 1/(_ki[i])*pow(Ddt_shiftDT,2))*expdtkby2;
+ DGiDT[i] = -GiT[i]/(_gi[i]) * Ddt_shiftDT*expdtgby2;
+ DDGiDT[i] = -GiT[i]/(_gi[i]) * (DDdt_shiftDT - 1/(_gi[i])*pow(Ddt_shiftDT,2))*expdtgby2;
}
for (int j=0; j<3; j++){
@@ -1673,9 +1700,9 @@ void mlawNonLinearTVE::predictorCorrector_ThermoViscoElastic(
STensorOperation::zero(DGammai); STensorOperation::zero(DDGammaiDT);
DQiDT -= DOiDT; DDQiDT -= DDOiDT; DDGammaiDT -= DDOGammaiDT;
for (int j=0; j<3; j++){
- Qi(j,j) += (_Ki[i]*trEe - trOi[i]);
+ Qi(j,j) += (KiT[i]*trEe - trOi[i]);
for (int k=0; k<3; k++){
- Qi(j,k) += (2*_Gi[i]*devEe(j,k) - devOi[i](j,k));
+ Qi(j,k) += (2*GiT[i]*devEe(j,k) - devOi[i](j,k));
}
}
DGammai += DEe;
@@ -1691,8 +1718,8 @@ void mlawNonLinearTVE::predictorCorrector_ThermoViscoElastic(
static STensor43 temp_DDQiDTDE;
static STensor43 temp_DDgammiDtDE;
- double Ki_term = (_Ki[i]*(1-expdtkby2)+T*DKiDT[i]);
- double Gi_term = (_Gi[i]*(1-expdtgby2)+T*DGiDT[i]);
+ double Ki_term = (KiT[i]*(1-expdtkby2)+T*DKiDT[i]);
+ double Gi_term = (GiT[i]*(1-expdtgby2)+T*DGiDT[i]);
isotropicHookTensor( Ki_term, Gi_term , temp_DDQiDTDE );
double Ki_term_2 = 1/3*expdtkby2;
diff --git a/NonLinearSolver/materialLaw/mlawNonLinearTVE.h b/NonLinearSolver/materialLaw/mlawNonLinearTVE.h
index 0e9b11cca833ccf6686443f0cee1372adcac2bde..90ae75b93e749f22388e8b8fce874349673cd7a4 100644
--- a/NonLinearSolver/materialLaw/mlawNonLinearTVE.h
+++ b/NonLinearSolver/materialLaw/mlawNonLinearTVE.h
@@ -79,6 +79,8 @@ class mlawNonLinearTVE : public mlawPowerYieldHyper{
virtual double linearInterpTemp(const double t_interp, const double t0, const double t1, const double T0, const double T1) const;
virtual double freeEnergyMechanical(const IPNonLinearTVE& q, const double Tc) const ;
virtual double setTemp(const double T0, const double T1, const int para) const ;
+
+ virtual void calculateMullinsEffectScaler(const IPNonLinearTVE* q0, IPNonLinearTVE *q1, const double T) const;
void corKirInfinity(const IPNonLinearTVE *q1, const STensor3& devEe, const double& trEe, const double T, STensor3& CorKirDevInf, double& CorKirTrInf) const;
void ThermoViscoElasticPredictor(const STensor3& Ee, const STensor3& Ee0,
@@ -226,21 +228,23 @@ class mlawNonLinearTVE : public mlawPowerYieldHyper{
}
virtual void getK(const IPNonLinearTVE *q1, double& KT, const double Tc, double *dK = NULL, double *ddK = NULL) const{
- KT = _K*_temFunc_K->getVal(Tc);
+ double K = _K*q1->getConstRefToElasticBulkPropertyScaleFactor();
+ KT = K*_temFunc_K->getVal(Tc);
if(dK!=NULL){
- *dK=_K*_temFunc_K->getDiff(Tc); // Here, *dK can be imagined as *(dK+0) or dK[0] if dK were a vector. See dKi below.
+ *dK = K*_temFunc_K->getDiff(Tc); // Here, *dK can be imagined as *(dK+0) or dK[0] if dK were a vector. See dKi below.
}
if(ddK!=NULL){
- *ddK = _K*_temFunc_K->getDoubleDiff(Tc);
+ *ddK = K*_temFunc_K->getDoubleDiff(Tc);
}
}
virtual void getG(const IPNonLinearTVE *q1, double& GT, const double Tc, double *dG = NULL, double *ddG = NULL) const{
- GT = _G*_temFunc_G->getVal(Tc);
+ double G = _G*q1->getConstRefToElasticShearPropertyScaleFactor();
+ GT = G*_temFunc_G->getVal(Tc);
if(dG!=NULL){
- *dG = _G*_temFunc_G->getDiff(Tc);
+ *dG = G*_temFunc_G->getDiff(Tc);
}
if(ddG!=NULL){
- *ddG = _G*_temFunc_G->getDoubleDiff(Tc);
+ *ddG = G*_temFunc_G->getDoubleDiff(Tc);
}
}
virtual void getAlpha(double& alphaT, const double Tc, double *dAlpha = NULL, double *ddAlpha = NULL) const{
@@ -270,28 +274,42 @@ class mlawNonLinearTVE : public mlawPowerYieldHyper{
// The 3rd argument is a pointer to a vector.
virtual void getKi(const IPNonLinearTVE *q1, std::vector<double>& Ki_T, const double Tc, std::vector<double>* dKi=NULL, std::vector<double>* ddKi=NULL) const{
+
+ std::vector<double> Ki;
+ Ki.clear(); Ki.resize(_N,0.);
+ for (int i=0; i<_Ki.size();i++){
+ Ki[i] = _Ki[i]*q1->getConstRefToElasticBulkPropertyScaleFactor();
+ }
+
Ki_T.resize(_N,0.);
for (int i=0; i<_Ki.size();i++){
- Ki_T[i] = _Ki[i]*_temFunc_Ki[i]->getVal(Tc);
+ Ki_T[i] = Ki[i]*_temFunc_Ki[i]->getVal(Tc);
if(dKi!=NULL){
- dKi->at(i) = _Ki[i]*_temFunc_Ki[i]->getDiff(Tc);
+ dKi->at(i) = Ki[i]*_temFunc_Ki[i]->getDiff(Tc);
}
if(ddKi!=NULL){
- ddKi->at(i) = _Ki[i]*_temFunc_Ki[i]->getDoubleDiff(Tc);
+ ddKi->at(i) = Ki[i]*_temFunc_Ki[i]->getDoubleDiff(Tc);
}
}
}
// The 3rd argument is a pointer to a vector.
virtual void getGi(const IPNonLinearTVE *q1, std::vector<double>& Gi_T, const double Tc, std::vector<double>* dGi=NULL, std::vector<double>* ddGi=NULL) const{
+
+ std::vector<double> Gi;
+ Gi.clear(); Gi.resize(_N,0.);
+ for (int i=0; i<_Gi.size();i++){
+ Gi[i] = _Gi[i]*q1->getConstRefToElasticShearPropertyScaleFactor();
+ }
+
Gi_T.resize(_N,0.);
for (int i=0; i<_Gi.size();i++){
- Gi_T[i] = _Gi[i]*_temFunc_Gi[i]->getVal(Tc);
+ Gi_T[i] = Gi[i]*_temFunc_Gi[i]->getVal(Tc);
if(dGi!=NULL){
- dGi->at(i) = _Gi[i]*_temFunc_Gi[i]->getDiff(Tc); // *(dGi+i) is equivalent to dGi[i]
+ dGi->at(i) = Gi[i]*_temFunc_Gi[i]->getDiff(Tc); // *(dGi+i) is equivalent to dGi[i]
}
if(ddGi!=NULL){
- ddGi->at(i) = _Gi[i]*_temFunc_Gi[i]->getDoubleDiff(Tc);
+ ddGi->at(i) = Gi[i]*_temFunc_Gi[i]->getDoubleDiff(Tc);
}
}
}
diff --git a/NonLinearSolver/materialLaw/mlawNonLinearTVP.cpp b/NonLinearSolver/materialLaw/mlawNonLinearTVP.cpp
index 8469054b8af14fabf738d25d917a1e5c2909ca2f..873300d9660130da66ba4632e51691ed142aaf97 100644
--- a/NonLinearSolver/materialLaw/mlawNonLinearTVP.cpp
+++ b/NonLinearSolver/materialLaw/mlawNonLinearTVP.cpp
@@ -380,12 +380,12 @@ void mlawNonLinearTVP::getFreeEnergyTVM(const IPNonLinearTVP *q0, IPNonLinearTVP
std::vector<double> KiT, GiT; KiT.resize(_N,0.); GiT.resize(_N,0.);
getKi(q1,KiT,T); getGi(q1,GiT,T);
for (int k=0; k<_Gi.size(); k++){
- DpsiVEdT += 1/9 * 1/_Ki[k] * trOi[k] * trDOiDT[k];
- DDpsiVEdTT += 1/9 * 1/_Ki[k] * (trOi[k]*trDDOiDTT[k] + pow(trDOiDT[k],2));
+ DpsiVEdT += 1/9 * 1/KiT[k] * trOi[k] * trDOiDT[k];
+ DDpsiVEdTT += 1/9 * 1/KiT[k] * (trOi[k]*trDDOiDTT[k] + pow(trDOiDT[k],2));
for (int i=0; i<3; i++)
for (int j=0; j<3; j++){
- DpsiVEdT += 2/_Gi[k] * devOi[k](i,j) * devDOiDT[k](i,j);
- DDpsiVEdTT += 2/_Gi[k] * (devOi[k](i,j) * devDDOiDTT[k](i,j) + pow(devDOiDT[k](i,j),2));
+ DpsiVEdT += 2/GiT[k] * devOi[k](i,j) * devDOiDT[k](i,j);
+ DDpsiVEdTT += 2/GiT[k] * (devOi[k](i,j) * devDDOiDTT[k](i,j) + pow(devDOiDT[k](i,j),2));
}
}
}
diff --git a/NonLinearSolver/materialLaw/mullinsEffect.cpp b/NonLinearSolver/materialLaw/mullinsEffect.cpp
index 022b622a808ce228e9b79149ce1691940ad1238d..47f32a8b9178899e13a5d6e790d7eb1cd15432a7 100644
--- a/NonLinearSolver/materialLaw/mullinsEffect.cpp
+++ b/NonLinearSolver/materialLaw/mullinsEffect.cpp
@@ -91,24 +91,30 @@ void negativeExponentialScaler::createIPVariable(IPMullinsEffect* &ipv) const
ipv = new IPMullinsEffect();
}
-void negativeExponentialScaler::mullinsEffectScaling(double psi0, double _psi, IPMullinsEffect &ipv) const{
+void negativeExponentialScaler::mullinsEffectScaling(double psi0, const IPMullinsEffect &ipvprev, double _psi, IPMullinsEffect &ipv) const{
+ // at unloading intitiation-> eta = 1.; unloading end-> eta = eta_min; reloading -> eta = eta_min;
+
+ double eta0 = ipvprev.getEta();
+
double psi = _psi;
- double eta = 0., DetaDpsi = 0., DDetaDpsipsi = 0., DetaDT = 0., DDetaDTT = 0.;
+ double eta = 1., DetaDpsi = 0., DDetaDpsipsi = 0., DetaDT = 0., DDetaDTT = 0.;
if (psi<psi0){
eta = 1. - _r*( 1. - exp( - sqrt(_m*(psi0-psi)) ) );
}
else{
- eta = 1.;
+ eta = eta0;
}
ipv.set(eta,DetaDpsi,DDetaDpsipsi,DetaDT,DDetaDTT);
}
-void negativeExponentialScaler::mullinsEffectScaling(double psi0, double _psi, double T, IPMullinsEffect &ipv) const{
+void negativeExponentialScaler::mullinsEffectScaling(double psi0, const IPMullinsEffect &ipvprev, double _psi, double T, IPMullinsEffect &ipv) const{
+ double eta0 = ipvprev.getEta();
+
double psi = _psi;
double r = _r*_temFunc_r->getVal(T);
double m = _m*_temFunc_m->getVal(T);
@@ -119,7 +125,7 @@ void negativeExponentialScaler::mullinsEffectScaling(double psi0, double _psi, d
eta = 1. - r*( 1. - exp( - sqrt(m*(psi0-psi)) ) );
}
else{
- eta = 1.;
+ eta = eta0;
}
ipv.set(eta,DetaDpsi,DDetaDpsipsi,DetaDT,DDetaDTT);
diff --git a/NonLinearSolver/materialLaw/mullinsEffect.h b/NonLinearSolver/materialLaw/mullinsEffect.h
index bbca308bc5dda7e811531773d623c4fae7ddffe3..9e06ee651f89844e61f4e9ad7a13eed49e76e590 100644
--- a/NonLinearSolver/materialLaw/mullinsEffect.h
+++ b/NonLinearSolver/materialLaw/mullinsEffect.h
@@ -33,8 +33,8 @@ class mullinsEffect {
virtual void createIPVariable(IPMullinsEffect* &ipv) const=0;
virtual void initLaws(const std::map<int,mullinsEffect*> &maplaw)=0;
virtual const bool isInitialized() {return _initialized;}
- virtual void mullinsEffectScaling(double psi0, double _psi, IPMullinsEffect &ipv) const=0;
- virtual void mullinsEffectScaling(double psi0, double _psi, double T, IPMullinsEffect &ipv) const=0;
+ virtual void mullinsEffectScaling(double psi0, const IPMullinsEffect &ipvprev, double _psi, IPMullinsEffect &ipv) const=0;
+ virtual void mullinsEffectScaling(double psi0, const IPMullinsEffect &ipvprev, double _psi, double T, IPMullinsEffect &ipv) const=0;
virtual mullinsEffect * clone() const=0;
#endif
@@ -61,8 +61,8 @@ class negativeExponentialScaler : public mullinsEffect{
virtual scalingFunction getType() const{return mullinsEffect::negativeExponentialScaler; };
virtual void createIPVariable(IPMullinsEffect* &ipv) const;
virtual void initLaws(const std::map<int,mullinsEffect*> &maplaw) {}; //nothing now, we will see if we use the mapping
- virtual void mullinsEffectScaling(double psi0, double _psi, IPMullinsEffect &ipv) const;
- virtual void mullinsEffectScaling(double psi0, double _psi, double T, IPMullinsEffect &ipv) const;
+ virtual void mullinsEffectScaling(double psi0, const IPMullinsEffect &ipvprev, double _psi, IPMullinsEffect &ipv) const;
+ virtual void mullinsEffectScaling(double psi0, const IPMullinsEffect &ipvprev, double _psi, double T, IPMullinsEffect &ipv) const;
virtual mullinsEffect * clone() const;
#endif // SWIG
};
diff --git a/NonLinearSolver/nlmechsolpy.i b/NonLinearSolver/nlmechsolpy.i
index f030bdaa3e77ba0604619e63e6b81eaef0123e25..c929327301d993acaf92445f2ebcd5371481d0e2 100644
--- a/NonLinearSolver/nlmechsolpy.i
+++ b/NonLinearSolver/nlmechsolpy.i
@@ -15,6 +15,7 @@
#include "j2IsotropicHardening.h"
#include "kinematicHardening.h"
#include "viscosityLaw.h"
+ #include "mullinsEffect.h"
#include "CLengthLaw.h"
#include "NucleationLaw.h"
#include "nucleationCriterionLaw.h"
@@ -81,6 +82,7 @@
%nodefaultctor J2IsotropicHardening;
%nodefaultctor kinematicHardening;
%nodefaultctor viscosityLaw;
+%nodefaultctor mullinsEffect;
%nodefaultctor RateFailureLaw;
%nodefaultctor CLengthLaw;
%nodefaultctor DamageLaw;
@@ -160,6 +162,7 @@
%include "j2IsotropicHardening.h"
%include "kinematicHardening.h"
%include "viscosityLaw.h"
+%include "mullinsEffect.h"
%include "CLengthLaw.h"
%include "NucleationLaw.h"
%include "nucleationCriterionLaw.h"
diff --git a/dG3D/benchmarks/nonLinearTVP_cube/cubeTVP.py b/dG3D/benchmarks/nonLinearTVP_cube/cubeTVP.py
index c38ec3d627645b798615c7d9526e4b6e3abd32eb..d4acb99e9ca9ae13a1e20a824c9b37337fc0e19a 100644
--- a/dG3D/benchmarks/nonLinearTVP_cube/cubeTVP.py
+++ b/dG3D/benchmarks/nonLinearTVP_cube/cubeTVP.py
@@ -30,7 +30,7 @@ E = relSpec[0] #MPa # 2700
nu = 0.33 #0.33
K = E/3./(1.-2.*nu) # Bulk mudulus
mu = E/2./(1.+nu) # Shear mudulus
-rho = 7850e-9 #905e-12 # Bulk mass 905 kg/m3 -> 905e-12 tonne/mm3
+rho = 7850e-12 #905e-12 # Bulk mass 905 kg/m3 -> 905e-12 tonne/mm3
Cp = rho*1900e+6 # 1900 J/kg/K -> 1900e+6 Nmm/tonne/K
KThCon = 0.14 # 0.14 W/m/K -> 0.14 Nmm/s/mm/K
Alpha = 0.6e-6 # 1/K
diff --git a/dG3D/src/dG3DMaterialLaw.h b/dG3D/src/dG3DMaterialLaw.h
index 7248a08d22c735ad97653c8722c0303774a8d83a..3d1f07b6536fff28d4fa7839c5e44bb96a31ccc4 100644
--- a/dG3D/src/dG3DMaterialLaw.h
+++ b/dG3D/src/dG3DMaterialLaw.h
@@ -52,6 +52,7 @@
#include "mlawNonLinearTVM.h"
#include "mlawNonLinearTVE.h"
#include "mlawNonLinearTVP.h"
+#include "mullinsEffect.h"
#endif
#if defined(HAVE_TORCH)
#include <torch/torch.h>