diff --git a/NonLinearSolver/internalPoints/ipHyperelastic.cpp b/NonLinearSolver/internalPoints/ipHyperelastic.cpp
index ca34ee5070f2ee933105f0566c1f9ca5f4dbdf40..87549c92c4316ecdc5bbeec0556cf4c406183d6d 100644
--- a/NonLinearSolver/internalPoints/ipHyperelastic.cpp
+++ b/NonLinearSolver/internalPoints/ipHyperelastic.cpp
@@ -12,7 +12,8 @@
#include "ipHyperelastic.h"
#include "restartManager.h"
-IPHyperViscoElastic::IPHyperViscoElastic(const int N):IPVariableMechanics(),_N(N),_elasticEnergy(0.),_Ee(0.),_kirchhoff(0.){
+IPHyperViscoElastic::IPHyperViscoElastic(const int N):IPVariableMechanics(),_N(N),_elasticEnergy(0.),_Ee(0.),_kirchhoff(0.),
+ _irreversibleEnergy(0.),_DirreversibleEnergyDF(0.){
_A.clear();
_B.clear();
for (int i=0; i< N; i++){
@@ -22,7 +23,8 @@ IPHyperViscoElastic::IPHyperViscoElastic(const int N):IPVariableMechanics(),_N(N
}
};
IPHyperViscoElastic::IPHyperViscoElastic(const IPHyperViscoElastic& src): IPVariableMechanics(src), _Ee(src._Ee),
- _N(src._N),_A(src._A),_B(src._B),_elasticEnergy(src._elasticEnergy), _kirchhoff(src._kirchhoff){};
+ _N(src._N),_A(src._A),_B(src._B),_elasticEnergy(src._elasticEnergy), _kirchhoff(src._kirchhoff),
+ _irreversibleEnergy(src._irreversibleEnergy),_DirreversibleEnergyDF(src._DirreversibleEnergyDF){};
IPHyperViscoElastic& IPHyperViscoElastic::operator =(const IPVariable& src){
IPVariableMechanics::operator =(src);
@@ -34,6 +36,8 @@ IPHyperViscoElastic& IPHyperViscoElastic::operator =(const IPVariable& src){
_A = psrc->_A;
_B = psrc->_B;
_elasticEnergy = psrc->_elasticEnergy;
+ _irreversibleEnergy = psrc->_irreversibleEnergy;
+ _DirreversibleEnergyDF = psrc->_DirreversibleEnergyDF;
}
return *this;
};
@@ -46,6 +50,8 @@ void IPHyperViscoElastic::restart() {
restartManager::restart(_B);
restartManager::restart(_Ee);
restartManager::restart(_kirchhoff);
+ restartManager::restart(_irreversibleEnergy);
+ restartManager::restart(_DirreversibleEnergyDF);
}
@@ -113,21 +119,42 @@ IPHyperViscoElastoPlastic& IPHyperViscoElastoPlastic::operator =(const IPVariabl
_gF = ps->_gF;
_dgFdF = ps->_dgFdF;
- if (_ipCompression != NULL) delete _ipCompression;
- _ipCompression = NULL;
- if ( ps->_ipCompression != NULL) _ipCompression = dynamic_cast<IPJ2IsotropicHardening*>(ps->_ipCompression->clone());
-
- if (_ipTraction != NULL) delete _ipTraction;
- _ipTraction = NULL;
- if (ps->_ipTraction != NULL) _ipTraction = dynamic_cast<IPJ2IsotropicHardening*>(ps->_ipTraction->clone());
-
- if (_ipShear != NULL) delete _ipShear;
- _ipShear = NULL;
- if (ps->_ipShear != NULL) _ipShear = dynamic_cast<IPJ2IsotropicHardening*>(ps->_ipShear->clone());
-
- if (_ipKinematic!= NULL) delete _ipKinematic;
- _ipKinematic = NULL;
- if (ps->_ipKinematic != NULL) _ipKinematic = dynamic_cast<IPKinematicHardening*>(ps->_ipKinematic->clone());
+ if ( ps->_ipCompression != NULL) {
+ if (_ipCompression == NULL){
+ _ipCompression = dynamic_cast<IPJ2IsotropicHardening*>(ps->_ipCompression->clone());
+ }
+ else{
+ _ipCompression->operator=(*dynamic_cast<const IPVariable*>(ps->_ipCompression));
+ }
+ }
+
+ if ( ps->_ipTraction != NULL) {
+ if (_ipTraction == NULL){
+ _ipTraction = dynamic_cast<IPJ2IsotropicHardening*>(ps->_ipTraction->clone());
+ }
+ else{
+ _ipTraction->operator=(*dynamic_cast<const IPVariable*>(ps->_ipTraction));
+ }
+ }
+
+ if ( ps->_ipShear != NULL) {
+ if (_ipShear == NULL){
+ _ipShear = dynamic_cast<IPJ2IsotropicHardening*>(ps->_ipShear->clone());
+ }
+ else{
+ _ipShear->operator=(*dynamic_cast<const IPVariable*>(ps->_ipShear));
+ }
+ }
+
+ if ( ps->_ipKinematic != NULL) {
+ if (_ipKinematic == NULL){
+ _ipKinematic = dynamic_cast<IPKinematicHardening*>(ps->_ipKinematic->clone());
+ }
+ else{
+ _ipKinematic->operator=(*dynamic_cast<const IPVariable*>(ps->_ipKinematic));
+ }
+ }
+
}
return *this;
};
diff --git a/NonLinearSolver/internalPoints/ipHyperelastic.h b/NonLinearSolver/internalPoints/ipHyperelastic.h
index 849442746b6f5bb7caf66394db3a1fa3919dbc64..0b665e40f386b620c2e4465c76c717c7163fc765 100644
--- a/NonLinearSolver/internalPoints/ipHyperelastic.h
+++ b/NonLinearSolver/internalPoints/ipHyperelastic.h
@@ -26,6 +26,9 @@ class IPHyperViscoElastic : public IPVariableMechanics{
STensor3 _Ee; // elastic strain
STensor3 _kirchhoff; // corotational Kirchhoff stress
+ double _irreversibleEnergy;
+ STensor3 _DirreversibleEnergyDF;
+
public:
IPHyperViscoElastic(const int N);
IPHyperViscoElastic(const IPHyperViscoElastic& src);
@@ -43,6 +46,13 @@ class IPHyperViscoElastic : public IPVariableMechanics{
virtual double defoEnergy() const{return _elasticEnergy;}
virtual double& getRefToElasticEnergy() {return _elasticEnergy;};
virtual double plasticEnergy() const{return 0.;}
+
+ virtual double irreversibleEnergy() const {return _irreversibleEnergy;};
+ virtual double & getRefToIrreversibleEnergy() {return _irreversibleEnergy;};
+
+ virtual STensor3& getRefToDIrreversibleEnergyDF() {return _DirreversibleEnergyDF;};
+ virtual const STensor3& getConstRefToDIrreversibleEnergyDF() const{return _DirreversibleEnergyDF;};
+
virtual void restart();
};
@@ -73,7 +83,7 @@ class IPHyperViscoElastoPlastic : public IPHyperViscoElastic{
double _r; // failure onset
double _gF; //
STensor3 _dgFdF;
-
+
public:
IPHyperViscoElastoPlastic(const J2IsotropicHardening* comp,
diff --git a/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.cpp b/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.cpp
index 97ee4fffc38c5f77dac797e090e50ca3a9ff6242..ae5d47515a852dc842798228340789821efe1e11 100644
--- a/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.cpp
+++ b/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.cpp
@@ -11,26 +11,18 @@
#include "ipNonLocalDamageHyperelastic.h"
#include "restartManager.h"
-IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage(const J2IsotropicHardening* comp,
+IPHyperViscoElastoPlasticLocalDamage::IPHyperViscoElastoPlasticLocalDamage(const J2IsotropicHardening* comp,
const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin, const int N,
- const CLengthLaw *cll, const DamageLaw *daml): IPHyperViscoElastoPlastic(comp,trac,shear,kin,N),
- _nonlocalPlasticStrain (0), _Dc(1.),_damageBlocked(false)
+ const DamageLaw *daml): IPHyperViscoElastoPlastic(comp,trac,shear,kin,N),
+ _Dc(1.),_damageBlocked(false)
{
- ipvCL=NULL;
- if(cll ==NULL) Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage has no cll");
- cll->createIPVariable(ipvCL);
-
ipvDam=NULL;
- if(daml ==NULL) Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage has no daml");
+ if(daml ==NULL) Msg::Error("IPHyperViscoElastoPlasticLocalDamage::IPHyperViscoElastoPlasticLocalDamage has no daml");
daml->createIPVariable(ipvDam);
};
-IPHyperViscoElastoPlasticNonLocalDamage::~IPHyperViscoElastoPlasticNonLocalDamage(){
- if(ipvCL != NULL)
- {
- delete ipvCL;
- }
+IPHyperViscoElastoPlasticLocalDamage::~IPHyperViscoElastoPlasticLocalDamage(){
if(ipvDam !=NULL)
{
delete ipvDam;
@@ -38,44 +30,24 @@ IPHyperViscoElastoPlasticNonLocalDamage::~IPHyperViscoElastoPlasticNonLocalDamag
};
-IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage(const IPHyperViscoElastoPlasticNonLocalDamage &source):
+IPHyperViscoElastoPlasticLocalDamage::IPHyperViscoElastoPlasticLocalDamage(const IPHyperViscoElastoPlasticLocalDamage &source):
IPHyperViscoElastoPlastic(source){
- _nonlocalPlasticStrain=source._nonlocalPlasticStrain;
_Dc = source._Dc;
_damageBlocked = source._damageBlocked;
- ipvCL = NULL;
- if(source.ipvCL != NULL)
- {
- ipvCL = dynamic_cast<IPCLength*>(source.ipvCL->clone());
- }
ipvDam = NULL;
if(source.ipvDam != NULL)
{
ipvDam = dynamic_cast<IPDamage*>(source.ipvDam->clone());
}
};
-IPHyperViscoElastoPlasticNonLocalDamage& IPHyperViscoElastoPlasticNonLocalDamage::operator=(const IPVariable &source){
+IPHyperViscoElastoPlasticLocalDamage& IPHyperViscoElastoPlasticLocalDamage::operator=(const IPVariable &source){
IPHyperViscoElastoPlastic::operator=(source);
- const IPHyperViscoElastoPlasticNonLocalDamage* src = dynamic_cast<const IPHyperViscoElastoPlasticNonLocalDamage*>(&source);
+ const IPHyperViscoElastoPlasticLocalDamage* src = dynamic_cast<const IPHyperViscoElastoPlasticLocalDamage*>(&source);
if(src != NULL)
{
- _nonlocalPlasticStrain=src->_nonlocalPlasticStrain;
_Dc = src->_Dc;
_damageBlocked = src->_damageBlocked;
-
- if(src->ipvCL != NULL)
- {
- if (ipvCL != NULL){
- ipvCL->operator=(*dynamic_cast<const IPVariable*>(src->ipvCL));
- }
- else
- ipvCL= dynamic_cast<IPCLength*>(src->ipvCL->clone());
- }
- else{
- if(ipvCL != NULL) delete ipvCL; ipvCL = NULL;
- }
-
if(src->ipvDam != NULL)
{
@@ -91,82 +63,177 @@ IPHyperViscoElastoPlasticNonLocalDamage& IPHyperViscoElastoPlasticNonLocalDamage
return *this;
};
-void IPHyperViscoElastoPlasticNonLocalDamage::restart()
+void IPHyperViscoElastoPlasticLocalDamage::restart()
{
IPHyperViscoElastoPlastic::restart();
- restartManager::restart(ipvCL);
restartManager::restart(ipvDam);
- restartManager::restart(_nonlocalPlasticStrain);
restartManager::restart(_Dc);
restartManager::restart(_damageBlocked);
return;
}
+IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage(const J2IsotropicHardening* comp,
+ const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin, const int N,
+ const CLengthLaw *cll, const DamageLaw *daml): IPHyperViscoElastoPlasticLocalDamage(comp,trac,shear,kin,N,daml),
+ _nonlocalPlasticStrain(0.), _DirreversibleEnergyDNonLocalVariable(0.)
+{
+ ipvCL=NULL;
+ if(cll ==NULL) Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage has no cll");
+ cll->createIPVariable(ipvCL);
+};
+
+IPHyperViscoElastoPlasticNonLocalDamage::~IPHyperViscoElastoPlasticNonLocalDamage(){
+ if(ipvCL != NULL)
+ {
+ delete ipvCL;
+ }
+};
+
+IPHyperViscoElastoPlasticNonLocalDamage::IPHyperViscoElastoPlasticNonLocalDamage(const IPHyperViscoElastoPlasticNonLocalDamage &source):
+ IPHyperViscoElastoPlasticLocalDamage(source), _nonlocalPlasticStrain(source._nonlocalPlasticStrain),
+ _DirreversibleEnergyDNonLocalVariable(source._DirreversibleEnergyDNonLocalVariable){
+ ipvCL = NULL;
+ if(source.ipvCL != NULL)
+ {
+ ipvCL = dynamic_cast<IPCLength*>(source.ipvCL->clone());
+ }
+};
+IPHyperViscoElastoPlasticNonLocalDamage& IPHyperViscoElastoPlasticNonLocalDamage::operator=(const IPVariable &source){
+ IPHyperViscoElastoPlasticLocalDamage::operator=(source);
+ const IPHyperViscoElastoPlasticNonLocalDamage* src = dynamic_cast<const IPHyperViscoElastoPlasticNonLocalDamage*>(&source);
+ if(src != NULL)
+ {
+ _nonlocalPlasticStrain=src->_nonlocalPlasticStrain;
+ _DirreversibleEnergyDNonLocalVariable = src->_DirreversibleEnergyDNonLocalVariable;
+ if(src->ipvCL != NULL)
+ {
+ if (ipvCL != NULL){
+ ipvCL->operator=(*dynamic_cast<const IPVariable*>(src->ipvCL));
+ }
+ else
+ ipvCL= dynamic_cast<IPCLength*>(src->ipvCL->clone());
+ }
+ else{
+ if(ipvCL != NULL) delete ipvCL; ipvCL = NULL;
+ }
+
+ }
+ return *this;
+};
+
+void IPHyperViscoElastoPlasticNonLocalDamage::restart()
+{
+ IPHyperViscoElastoPlasticLocalDamage::restart();
+ restartManager::restart(ipvCL);
+ restartManager::restart(_nonlocalPlasticStrain);
+ restartManager::restart(_DirreversibleEnergyDNonLocalVariable);
+};
+
+IPHyperViscoElastoPlasticMultipleLocalDamage::IPHyperViscoElastoPlasticMultipleLocalDamage(const J2IsotropicHardening* comp,
+ const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin,
+ const int N, const std::vector<DamageLaw*>& dl): IPHyperViscoElastoPlastic(comp,trac,shear,kin,N),
+ _Dc(1.){
+ numNonLocalVariable = dl.size();
+ ipvDam.clear();
+ for (int i=0; i< dl.size(); i++){
+ IPDamage* ipvD = NULL;
+ dl[i]->createIPVariable(ipvD);
+ ipvDam.push_back(ipvD);
+ }
+};
+IPHyperViscoElastoPlasticMultipleLocalDamage::~IPHyperViscoElastoPlasticMultipleLocalDamage(){
+ for (int i=0; i< numNonLocalVariable; i++){
+ delete ipvDam[i];
+ }
+ ipvDam.clear();
+};
+IPHyperViscoElastoPlasticMultipleLocalDamage::IPHyperViscoElastoPlasticMultipleLocalDamage(const IPHyperViscoElastoPlasticMultipleLocalDamage& source):
+ IPHyperViscoElastoPlastic(source){
+ numNonLocalVariable = source.numNonLocalVariable;
+ _Dc = source._Dc;
+
+ ipvDam.clear();
+ for (int i=0; i< numNonLocalVariable; i++){
+ if(source.ipvDam[i] != NULL){
+ ipvDam.push_back(dynamic_cast<IPDamage*>(source.ipvDam[i]->clone()));
+ }
+ }
+};
+IPHyperViscoElastoPlasticMultipleLocalDamage& IPHyperViscoElastoPlasticMultipleLocalDamage::operator=(const IPVariable& source){
+ IPHyperViscoElastoPlastic::operator=(source);
+ const IPHyperViscoElastoPlasticMultipleLocalDamage* ipsource = dynamic_cast<const IPHyperViscoElastoPlasticMultipleLocalDamage*>(&source);
+ if (ipsource != NULL){
+ numNonLocalVariable = ipsource->numNonLocalVariable;
+ _Dc = ipsource->_Dc;
+ for (int i=0; i< ipsource->numNonLocalVariable; i++){
+ if(ipsource->ipvDam[i] != NULL){
+ if (ipvDam[i] != NULL){
+ ipvDam[i]->operator=(*dynamic_cast<const IPVariable*>(ipsource->ipvDam[i]));
+ }
+ }
+ }
+ }
+ return *this;
+};
+
+
+void IPHyperViscoElastoPlasticMultipleLocalDamage::restart(){
+ IPHyperViscoElastoPlastic::restart();
+ restartManager::restart(numNonLocalVariable);
+ restartManager::restart(_Dc);
+ for (int i=0; i< numNonLocalVariable; i++){
+ ipvDam[i]->restart();
+ }
+};
+
+
IPHyperViscoElastoPlasticMultipleNonLocalDamage::IPHyperViscoElastoPlasticMultipleNonLocalDamage(const J2IsotropicHardening* comp,
const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin,
- const int N, const std::vector<CLengthLaw*>& cll, const std::vector<DamageLaw*>& dl): IPHyperViscoElastoPlastic(comp,trac,shear,kin,N),
- _nonlocalEqPlasticStrain(0.),_nonlocalFailurePlasticity(0.),_Dc(1.){
- numNonLocalVariable = cll.size();
+ const int N, const std::vector<CLengthLaw*>& cll, const std::vector<DamageLaw*>& dl):
+ IPHyperViscoElastoPlasticMultipleLocalDamage(comp,trac,shear,kin,N,dl),
+ _nonlocalEqPlasticStrain(0.),_nonlocalFailurePlasticity(0.){
ipvCL.clear();
- ipvDam.clear();
+ _DirreversibleEnergyDNonLocalVariable.clear();
for (int i=0; i< cll.size(); i++){
IPCLength* ipvLength=NULL;
cll[i]->createIPVariable(ipvLength);
ipvCL.push_back(ipvLength);
- IPDamage* ipvD = NULL;
- dl[i]->createIPVariable(ipvD);
- ipvDam.push_back(ipvD);
+ _DirreversibleEnergyDNonLocalVariable.push_back(0.);
}
};
IPHyperViscoElastoPlasticMultipleNonLocalDamage::~IPHyperViscoElastoPlasticMultipleNonLocalDamage(){
for (int i=0; i< numNonLocalVariable; i++){
delete ipvCL[i];
- delete ipvDam[i];
}
ipvCL.clear();
- ipvDam.clear();
};
IPHyperViscoElastoPlasticMultipleNonLocalDamage::IPHyperViscoElastoPlasticMultipleNonLocalDamage(const IPHyperViscoElastoPlasticMultipleNonLocalDamage& source):
- IPHyperViscoElastoPlastic(source){
- numNonLocalVariable = source.numNonLocalVariable;
+ IPHyperViscoElastoPlasticMultipleLocalDamage(source){
_nonlocalEqPlasticStrain= source._nonlocalEqPlasticStrain;
_nonlocalFailurePlasticity = source._nonlocalFailurePlasticity;
- _Dc = source._Dc;
-
+ _DirreversibleEnergyDNonLocalVariable = source._DirreversibleEnergyDNonLocalVariable;
ipvCL.clear();
- ipvDam.clear();
for (int i=0; i< numNonLocalVariable; i++){
if(source.ipvCL[i] != NULL){
ipvCL.push_back(dynamic_cast<IPCLength*>(source.ipvCL[i]->clone()));
}
- if(source.ipvDam[i] != NULL){
- ipvDam.push_back(dynamic_cast<IPDamage*>(source.ipvDam[i]->clone()));
- }
}
};
IPHyperViscoElastoPlasticMultipleNonLocalDamage& IPHyperViscoElastoPlasticMultipleNonLocalDamage::operator=(const IPVariable& source){
- IPHyperViscoElastoPlastic::operator=(source);
+ IPHyperViscoElastoPlasticMultipleLocalDamage::operator=(source);
const IPHyperViscoElastoPlasticMultipleNonLocalDamage* ipsource = dynamic_cast<const IPHyperViscoElastoPlasticMultipleNonLocalDamage*>(&source);
if (ipsource != NULL){
- for (int i=0; i< numNonLocalVariable; i++){
- if (ipvCL[i] != NULL) delete ipvCL[i];
- if (ipvDam[i] != NULL) delete ipvDam[i];
- }
- ipvCL.clear();
- ipvDam.clear();
-
- numNonLocalVariable = ipsource->numNonLocalVariable;
_nonlocalEqPlasticStrain = ipsource->_nonlocalEqPlasticStrain;
_nonlocalFailurePlasticity = ipsource->_nonlocalFailurePlasticity;
- _Dc = ipsource->_Dc;
+ _DirreversibleEnergyDNonLocalVariable = ipsource->_DirreversibleEnergyDNonLocalVariable;
+
for (int i=0; i< ipsource->numNonLocalVariable; i++){
if(ipsource->ipvCL[i] != NULL){
- ipvCL.push_back(dynamic_cast<IPCLength*>(ipsource->ipvCL[i]->clone()));
- }
- if(ipsource->ipvDam[i] != NULL){
- ipvDam.push_back(dynamic_cast<IPDamage*>(ipsource->ipvDam[i]->clone()));
+ if (ipvCL[i] != NULL){
+ ipvCL[i]->operator=(*dynamic_cast<const IPVariable*>(ipsource->ipvCL[i]));
+ }
}
}
}
@@ -175,14 +242,13 @@ IPHyperViscoElastoPlasticMultipleNonLocalDamage& IPHyperViscoElastoPlasticMultip
void IPHyperViscoElastoPlasticMultipleNonLocalDamage::restart(){
- IPHyperViscoElastoPlastic::restart();
- restartManager::restart(numNonLocalVariable);
+ IPHyperViscoElastoPlasticMultipleLocalDamage::restart();
restartManager::restart(_nonlocalEqPlasticStrain);
restartManager::restart(_nonlocalFailurePlasticity);
- restartManager::restart(_Dc);
+ for (int i=0; i< numNonLocalVariable; i++){
+ restartManager::restart(_DirreversibleEnergyDNonLocalVariable[i]);
+ }
for (int i=0; i< numNonLocalVariable; i++){
ipvCL[i]->restart();
- ipvDam[i]->restart();
}
-
};
diff --git a/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.h b/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.h
index a6cb6ccde928d3e44953d995fc5416cf8bc32e43..6ec8f8f6d89e005a3f36601714559252f63df560 100644
--- a/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.h
+++ b/NonLinearSolver/internalPoints/ipNonLocalDamageHyperelastic.h
@@ -18,27 +18,23 @@
#include "CLengthLaw.h"
#include "DamageLaw.h"
-class IPHyperViscoElastoPlasticNonLocalDamage : public IPHyperViscoElastoPlastic
+class IPHyperViscoElastoPlasticLocalDamage : public IPHyperViscoElastoPlastic
{
protected:
- IPCLength* ipvCL;
IPDamage* ipvDam;
-
- double _nonlocalPlasticStrain;
double _Dc; // critical damage
bool _damageBlocked;
public:
- IPHyperViscoElastoPlasticNonLocalDamage(const J2IsotropicHardening* comp,
+ IPHyperViscoElastoPlasticLocalDamage(const J2IsotropicHardening* comp,
const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin,
- const int N,
- const CLengthLaw *cll, const DamageLaw *dl);
- virtual ~IPHyperViscoElastoPlasticNonLocalDamage();
+ const int N, const DamageLaw *dl);
+ virtual ~IPHyperViscoElastoPlasticLocalDamage();
- IPHyperViscoElastoPlasticNonLocalDamage(const IPHyperViscoElastoPlasticNonLocalDamage &source);
- virtual IPHyperViscoElastoPlasticNonLocalDamage& operator=(const IPVariable &source);
+ IPHyperViscoElastoPlasticLocalDamage(const IPHyperViscoElastoPlasticLocalDamage &source);
+ virtual IPHyperViscoElastoPlasticLocalDamage& operator=(const IPVariable &source);
- virtual IPVariable* clone() const{return new IPHyperViscoElastoPlasticNonLocalDamage(*this);}
+ virtual IPVariable* clone() const{return new IPHyperViscoElastoPlasticLocalDamage(*this);}
virtual void blockDamage(const bool fl){_damageBlocked = fl;};
@@ -47,187 +43,234 @@ class IPHyperViscoElastoPlasticNonLocalDamage : public IPHyperViscoElastoPlastic
virtual void setCriticalDamage(const double DT){_Dc = DT;};
virtual double getCriticalDamage() const {return _Dc;};
-
virtual void restart();
- virtual double getCurrentPlasticStrain() const { return _epspbarre;}
- virtual double &getRefToCurrentPlasticStrain() { return _epspbarre;}
- virtual double getEffectivePlasticStrain() const { return _nonlocalPlasticStrain;}
- virtual double &getRefToEffectivePlasticStrain() { return _nonlocalPlasticStrain;}
-
- virtual const IPCLength &getConstRefToIPCLength() const
- {
- if(ipvCL==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
- return *ipvCL;
- }
- virtual IPCLength &getRefToIPCLength()
- {
- if(ipvCL==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
- return *ipvCL;
- }
- virtual const STensor3 &getConstRefToCharacteristicLength() const
- {
- if(ipvCL==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
- return ipvCL->getConstRefToCL();
- }
- virtual STensor3 &getRefToCharacteristicLength() const
- {
- if(ipvCL==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
- return ipvCL->getRefToCL();
- }
virtual const IPDamage &getConstRefToIPDamage() const
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return *ipvDam;
}
virtual IPDamage &getRefToIPDamage()
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return *ipvDam;
}
virtual double getDamage() const
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return ipvDam->getDamage();
}
virtual double getMaximalP() const
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return ipvDam->getMaximalP();
}
virtual double getDeltaDamage() const
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return ipvDam->getDeltaDamage();
}
virtual double getDDamageDp() const
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return ipvDam->getDDamageDp();
}
virtual const STensor3 &getConstRefToDDamageDFe() const
{
if(ipvDam==NULL)
- Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvDam not initialized");
+ Msg::Error("IPHyperViscoElastoPlasticLocalDamage: ipvDam not initialized");
return ipvDam->getConstRefToDDamageDFe();
}
};
-class IPHyperViscoElastoPlasticMultipleNonLocalDamage : public IPHyperViscoElastoPlastic{
+
+class IPHyperViscoElastoPlasticNonLocalDamage : public IPHyperViscoElastoPlasticLocalDamage
+{
+ protected:
+ IPCLength* ipvCL;
+ double _nonlocalPlasticStrain;
+ double _DirreversibleEnergyDNonLocalVariable;
+
+ public:
+ IPHyperViscoElastoPlasticNonLocalDamage(const J2IsotropicHardening* comp,
+ const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin,
+ const int N,
+ const CLengthLaw *cll, const DamageLaw *dl);
+ virtual ~IPHyperViscoElastoPlasticNonLocalDamage();
+
+ IPHyperViscoElastoPlasticNonLocalDamage(const IPHyperViscoElastoPlasticNonLocalDamage &source);
+ virtual IPHyperViscoElastoPlasticNonLocalDamage& operator=(const IPVariable &source);
+
+ virtual IPVariable* clone() const{return new IPHyperViscoElastoPlasticNonLocalDamage(*this);}
+
+ virtual const double& getDIrreversibleEnergyDNonLocalVariable() const {return _DirreversibleEnergyDNonLocalVariable;};
+ virtual double& getRefToDIrreversibleEnergyDNonLocalVariable() {return _DirreversibleEnergyDNonLocalVariable;};
+
+ virtual void restart();
+ virtual double getCurrentPlasticStrain() const { return _epspbarre;}
+ virtual double &getRefToCurrentPlasticStrain() { return _epspbarre;}
+ virtual double getEffectivePlasticStrain() const { return _nonlocalPlasticStrain;}
+ virtual double &getRefToEffectivePlasticStrain() { return _nonlocalPlasticStrain;}
+
+ virtual const IPCLength &getConstRefToIPCLength() const
+ {
+ if(ipvCL==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
+ return *ipvCL;
+ }
+ virtual IPCLength &getRefToIPCLength()
+ {
+ if(ipvCL==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
+ return *ipvCL;
+ }
+ virtual const STensor3 &getConstRefToCharacteristicLength() const
+ {
+ if(ipvCL==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
+ return ipvCL->getConstRefToCL();
+ }
+ virtual STensor3 &getRefToCharacteristicLength() const
+ {
+ if(ipvCL==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticNonLocalDamage: ipvCL not initialized");
+ return ipvCL->getRefToCL();
+ }
+};
+
+class IPHyperViscoElastoPlasticMultipleLocalDamage : public IPHyperViscoElastoPlastic{
protected:
int numNonLocalVariable;
- std::vector<IPCLength*> ipvCL;
std::vector<IPDamage*> ipvDam;
-
- // bareps
- double _nonlocalEqPlasticStrain;
- double _nonlocalFailurePlasticity;
-
+
double _Dc; // critical damage
public:
- IPHyperViscoElastoPlasticMultipleNonLocalDamage(const J2IsotropicHardening* comp,
+ IPHyperViscoElastoPlasticMultipleLocalDamage(const J2IsotropicHardening* comp,
const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin,
- const int N, const std::vector<CLengthLaw*>& cll, const std::vector<DamageLaw*>& dl);
- virtual ~IPHyperViscoElastoPlasticMultipleNonLocalDamage();
- IPHyperViscoElastoPlasticMultipleNonLocalDamage(const IPHyperViscoElastoPlasticMultipleNonLocalDamage& source);
- IPHyperViscoElastoPlasticMultipleNonLocalDamage& operator=(const IPVariable& source);
+ const int N, const std::vector<DamageLaw*>& dl);
+ virtual ~IPHyperViscoElastoPlasticMultipleLocalDamage();
+ IPHyperViscoElastoPlasticMultipleLocalDamage(const IPHyperViscoElastoPlasticMultipleLocalDamage& source);
+ virtual IPHyperViscoElastoPlasticMultipleLocalDamage& operator=(const IPVariable& source);
virtual void setCriticalDamage(const double DT){_Dc = DT;};
virtual double getCriticalDamage() const {return _Dc;};
virtual int getNumNonLocalVariable() const{return numNonLocalVariable;};
- virtual IPVariable* clone() const{return new IPHyperViscoElastoPlasticMultipleNonLocalDamage(*this);}
+ virtual IPVariable* clone() const{return new IPHyperViscoElastoPlasticMultipleLocalDamage(*this);}
virtual void restart();
- virtual double getCurrentPlasticStrain() const { return _epspbarre;}
- virtual double &getRefToCurrentPlasticStrain() { return _epspbarre;}
- virtual double getEffectivePlasticStrain() const { return _nonlocalEqPlasticStrain;}
- virtual double &getRefToEffectivePlasticStrain() { return _nonlocalEqPlasticStrain;}
-
- virtual double getNonLocalFailurePlasticity() const {return _nonlocalFailurePlasticity;};
- virtual double & getRefToNonLocalFailurePlasticity() {return _nonlocalFailurePlasticity;};
-
- virtual const IPCLength &getConstRefToIPCLength(const int i) const
- {
- if(ipvCL[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
- return *ipvCL[i];
- }
- virtual IPCLength &getRefToIPCLength(const int i)
- {
- if(ipvCL[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
- return *ipvCL[i];
- }
- virtual const STensor3 &getConstRefToCharacteristicLength(const int i) const
- {
- if(ipvCL[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
- return ipvCL[i]->getConstRefToCL();
- }
- virtual STensor3 &getRefToCharacteristicLength(const int i) const
- {
- if(ipvCL[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
- return ipvCL[i]->getRefToCL();
- }
virtual const IPDamage &getConstRefToIPDamage(const int i) const
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return *ipvDam[i];
}
virtual IPDamage &getRefToIPDamage(const int i)
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return *ipvDam[i];
}
virtual double getDamage(const int i) const
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return ipvDam[i]->getDamage();
}
virtual double getMaximalP(const int i) const
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return ipvDam[i]->getMaximalP();
}
virtual double getDeltaDamage(const int i) const
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return ipvDam[i]->getDeltaDamage();
}
virtual double getDDamageDp( const int i) const
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return ipvDam[i]->getDDamageDp();
}
virtual const STensor3 &getConstRefToDDamageDFe(const int i) const
{
if(ipvDam[i]==NULL)
- Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ Msg::Error("IPHyperViscoElastoPlasticMultipleLocalDamage: ipvCL[%d] not initialized",i);
return ipvDam[i]->getConstRefToDDamageDFe();
}
};
+
+
+class IPHyperViscoElastoPlasticMultipleNonLocalDamage : public IPHyperViscoElastoPlasticMultipleLocalDamage{
+ protected:
+ std::vector<IPCLength*> ipvCL;
+ // bareps
+ double _nonlocalEqPlasticStrain;
+ double _nonlocalFailurePlasticity;
+ std::vector<double> _DirreversibleEnergyDNonLocalVariable;
+
+ public:
+ IPHyperViscoElastoPlasticMultipleNonLocalDamage(const J2IsotropicHardening* comp,
+ const J2IsotropicHardening* trac,const J2IsotropicHardening* shear, const kinematicHardening* kin,
+ const int N, const std::vector<CLengthLaw*>& cll, const std::vector<DamageLaw*>& dl);
+ virtual ~IPHyperViscoElastoPlasticMultipleNonLocalDamage();
+ IPHyperViscoElastoPlasticMultipleNonLocalDamage(const IPHyperViscoElastoPlasticMultipleNonLocalDamage& source);
+ virtual IPHyperViscoElastoPlasticMultipleNonLocalDamage& operator=(const IPVariable& source);
+
+ virtual IPVariable* clone() const{return new IPHyperViscoElastoPlasticMultipleNonLocalDamage(*this);}
+
+ virtual const double& getDIrreversibleEnergyDNonLocalVariable(const int i) const {return _DirreversibleEnergyDNonLocalVariable[i];};
+ virtual double& getRefToDIrreversibleEnergyDNonLocalVariable(const int i) {return _DirreversibleEnergyDNonLocalVariable[i];};
+
+ virtual void restart();
+ virtual double getCurrentPlasticStrain() const { return _epspbarre;}
+ virtual double &getRefToCurrentPlasticStrain() { return _epspbarre;}
+
+ virtual double getEffectivePlasticStrain() const { return _nonlocalEqPlasticStrain;}
+ virtual double &getRefToEffectivePlasticStrain() { return _nonlocalEqPlasticStrain;}
+
+ virtual double getNonLocalFailurePlasticity() const {return _nonlocalFailurePlasticity;};
+ virtual double & getRefToNonLocalFailurePlasticity() {return _nonlocalFailurePlasticity;};
+
+ virtual const IPCLength &getConstRefToIPCLength(const int i) const
+ {
+ if(ipvCL[i]==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ return *ipvCL[i];
+ }
+ virtual IPCLength &getRefToIPCLength(const int i)
+ {
+ if(ipvCL[i]==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ return *ipvCL[i];
+ }
+ virtual const STensor3 &getConstRefToCharacteristicLength(const int i) const
+ {
+ if(ipvCL[i]==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ return ipvCL[i]->getConstRefToCL();
+ }
+ virtual STensor3 &getRefToCharacteristicLength(const int i) const
+ {
+ if(ipvCL[i]==NULL)
+ Msg::Error("IPHyperViscoElastoPlasticMultipleNonLocalDamage: ipvCL[%d] not initialized",i);
+ return ipvCL[i]->getRefToCL();
+ }
+};
+
#endif // IPNONLOCALDAMAGEHYPERELASTIC_H_
diff --git a/NonLinearSolver/materialLaw/mlaw.h b/NonLinearSolver/materialLaw/mlaw.h
index 8794466b052b2028a8eb8c8b4f9cddf5e60c26f1..d7076c8e0520163ab58995ac8da2977b7a0bbaf1 100644
--- a/NonLinearSolver/materialLaw/mlaw.h
+++ b/NonLinearSolver/materialLaw/mlaw.h
@@ -28,7 +28,7 @@ class materialLaw{
transverseIsoYarnB, Anisotropic, AnisotropicStoch, nonLocalDamage, vumat,
FSElastic, FSElastoPlastic, numeric, secondOrderElastic, j2smallstrain,nonLocalDamageGurson,nonLocalDamageJ2Hyper, nonLocalDamageIsotropicElasticity,LinearThermoMechanics,J2ThermoMechanics,SMP,LinearElecTherMech,
AnIsotropicElecTherMech, hyperelastic, quadYieldHyper, powerYieldLaw, powerYieldLawWithFailure, nonLocalDamageQuadYieldHyper,nonLocalDamagePowerYieldHyper,
- nonLocalDamagePowerYieldHyperWithFailure,ElecSMP,ThermalConducter,AnIsotropicTherMech, localDamageJ2Hyper,linearElastic};
+ localDamagePowerYieldHyperWithFailure,nonLocalDamagePowerYieldHyperWithFailure,ElecSMP,ThermalConducter,AnIsotropicTherMech, localDamageJ2Hyper,linearElastic};
protected :
diff --git a/NonLinearSolver/materialLaw/mlawHyperelastic.cpp b/NonLinearSolver/materialLaw/mlawHyperelastic.cpp
index ad7edaf4fdb27ca9fa71e520a8bf23cfd741bb41..05d19df95f8a51838b97b56a5e94000569cd9695 100644
--- a/NonLinearSolver/materialLaw/mlawHyperelastic.cpp
+++ b/NonLinearSolver/materialLaw/mlawHyperelastic.cpp
@@ -96,6 +96,23 @@ void mlawHyperViscoElastic::setViscoElasticData(const std::string filename){
fclose(file);
};
+double mlawHyperViscoElastic::deformationEnergy(const STensor3 &C) const
+{
+ static STensor3 logCdev;
+ STensorOperation::logSTensor3(C,_order,logCdev,NULL);
+ double trace = logCdev.trace();
+ double lnJ = 0.5*trace;
+ logCdev(0,0)-=trace/3.;
+ logCdev(1,1)-=trace/3.;
+ logCdev(2,2)-=trace/3.;
+
+ double Psy = _K*0.5*lnJ*lnJ+_mu*0.25*dot(logCdev,logCdev);
+ for (int i=0; i<_N; i++){
+ Psy += _Ki[i]*0.5*lnJ*lnJ+_Gi[i]*0.25*dot(logCdev,logCdev);
+ }
+ return Psy;;
+}
+
mlawHyperViscoElastic::mlawHyperViscoElastic(const int num,const double E,const double nu, const double rho,
const bool matrixbyPerturbation, const double pert):
@@ -380,6 +397,8 @@ void mlawHyperViscoElastic::predictorCorrector_ViscoElastic(const STensor3& F0,
STensorOperation::multSTensor3STensor43(corKir,dlnCdC,secondPK);
STensorOperation::multSTensor3(F,secondPK,P);
// first PK
+
+ q1->_elasticEnergy=deformationEnergy(C);
if (stiff){
static STensor43 DsecondPKdC;
@@ -479,27 +498,6 @@ void mlawPowerYieldHyper::setViscosityEffect(const viscosityLaw& v, const double
_viscosity = v.clone();
};
-double mlawPowerYieldHyper::deformationEnergy(const STensor3 &C) const
-{
- double Jac= sqrt((C(0,0) * (C(1,1) * C(2,2) - C(1,2) * C(2,1)) -
- C(0,1) * (C(1,0) * C(2,2) - C(1,2) * C(2,0)) +
- C(0,2) * (C(1,0) * C(2,1) - C(1,1) * C(2,0))));
- double lnJ = log(Jac);
- STensor3 logCdev;
- STensorOperation::logSTensor3(C,_order,logCdev,NULL);
- double trace = logCdev.trace();
- logCdev(0,0)-=trace/3.;
- logCdev(1,1)-=trace/3.;
- logCdev(2,2)-=trace/3.;
-
- double Psy = _K*0.5*lnJ*lnJ+_mu*0.25*dot(logCdev,logCdev);
- for (int i=0; i<_N; i++){
- Psy += _Ki[i]*0.5*lnJ*lnJ+_Gi[i]*0.25*dot(logCdev,logCdev);
- }
- return Psy;;
-}
-
-
void mlawPowerYieldHyper::hardening(IPHyperViscoElastoPlastic* q) const{
//Msg::Error("epspCompression = %e, epspTRaction = %e, epspShear = %e",q->_epspCompression,q->_epspTraction,q->_epspShear);
if (_compression != NULL && q->_ipCompression != NULL){
@@ -531,7 +529,7 @@ void mlawPowerYieldHyper::tangent_full_perturbation(
static STensor43 tmpSTensor43;
static STensor3 Fplus, Pplus;
- static IPHyperViscoElastoPlastic q11(*q1);
+ static IPHyperViscoElastoPlastic q11(*q0);
for (int i=0; i<3; i++){
for (int j=0; j<3; j++){
@@ -539,6 +537,7 @@ void mlawPowerYieldHyper::tangent_full_perturbation(
Fplus(i,j)+=_perturbationfactor;
this->predictorCorrector(Fplus,q0,&q11,Pplus,false,tmpSTensor43,tmpSTensor43,tmpSTensor43);
q1->_DgammaDF(i,j) = (q11._epspbarre - q1->_epspbarre)/_perturbationfactor;
+ q1->_DirreversibleEnergyDF(i,j) = (q11._irreversibleEnergy - q1->_irreversibleEnergy)/_perturbationfactor;
for (int k=0; k<3; k++){
for (int l=0; l<3; l++){
T_(k,l,i,j) = (Pplus(k,l)-P(k,l))/_perturbationfactor;
@@ -680,7 +679,7 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
q1->_nup = (9.-2.*_b)/(18.+2.*_b);
double kk = 1./sqrt(1.+2.*q1->_nup*q1->_nup);
- double Gamma = 0.; // q1->_Gamma; // flow rule parameter
+ double Gamma = 0.; // // flow rule parameter
double PhiEq = PhiEqpr;
@@ -712,12 +711,13 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
double dDgammaDGamma = 0.;
+ double Dgamma = 0.; // eqplastic strain
if (f>0){
// plasticity
int ite = 0;
int maxite = 100; // maximal number of iters
- double Dgamma = 0.; // eqplastic strain
+
//Msg::Error("plasticity occurs f = %e",f);
@@ -829,7 +829,10 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
q1->_kirchhoff(0,0) += (ptilde);
q1->_kirchhoff(1,1) += (ptilde);
q1->_kirchhoff(2,2) += (ptilde);
- };
+ }
+ else{
+ N *= 0.;
+ }
const STensor3& KS = q1->_kirchhoff;
@@ -845,10 +848,17 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
P(i,j) += Fe(i,k)*S(k,l)*Fpinv(j,l);
- // update current internal variables
+ // defo energy
q1->_elasticEnergy=deformationEnergy(Ce);
-
-
+
+ // dissipation
+ double& irrEnerg = q1->getRefToIrreversibleEnergy();
+ irrEnerg = q0->irreversibleEnergy();
+ if (Dgamma > 0){
+ double dotKSN = dot(KS,N);
+ irrEnerg += Gamma*dotKSN;
+ }
+
if (stiff){
static STensor3 DpprDCepr;
static STensor43 DdevKprDCepr;
@@ -864,6 +874,9 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
static STensor43 dFpDCepr;
static STensor3 DgamaDCepr;
+ static STensor3 DtrNDCepr;
+ static STensor43 DdevNDCepr;
+ static STensor3 DGDCepr;
if (Gamma >0){
// plastic
@@ -883,7 +896,7 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
}
}
- static STensor3 DGDCepr;
+
for (int i=0; i<3; i++){
for (int j=0; j<3; j++){
DGDCepr(i,j) = (-dfDCepr(i,j)-dfdDgamma*kk*Gamma*dAdCepr(i,j))/DfDGamma;
@@ -896,14 +909,11 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
}
}
-
- static STensor3 DtrNDCepr;
+
for (int i=0; i<3; i++)
for (int j=0; j<3; j++)
DtrNDCepr(i,j) = 2.*_b/v*DpprDCepr(i,j) - 2.*_b*ptildepr*(2.*_b*Kt)/(v*v)*DGDCepr(i,j);
-
- static STensor43 DdevNDCepr;
DdevNDCepr = (DdevKprDCepr);
DdevNDCepr *= (3./u);
for (int i=0; i<3; i++)
@@ -949,6 +959,9 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
// elastic
DgamaDCepr *= 0.;
dFpDCepr *= 0.;
+ DtrNDCepr *= 0.;
+ DdevNDCepr *= 0.;
+ DGDCepr *= 0.;
}
static STensor43 dKcorDcepr;
@@ -977,9 +990,15 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
STensor3& DgammaDF = q1->_DgammaDF;
static STensor43 DKcorDF;
- STensorOperation::multSTensor3STensor43(DgamaDCepr,CeprToF,DgammaDF);
STensorOperation::multSTensor43(dKcorDcepr,CeprToF,DKcorDF);
- STensorOperation::multSTensor43(dFpDCepr,CeprToF,dFpdF);
+ if (Gamma > 0){
+ STensorOperation::multSTensor3STensor43(DgamaDCepr,CeprToF,DgammaDF);
+ STensorOperation::multSTensor43(dFpDCepr,CeprToF,dFpdF);
+ }
+ else{
+ DgammaDF *= 0.;
+ dFpdF *= 0.;
+ }
static STensor43 DinvFpDF; //
for (int i=0; i<3; i++)
@@ -1055,7 +1074,29 @@ void mlawPowerYieldHyper::predictorCorrector_nonAssociatedFlow(const STensor3& F
}
}
}
-
+
+
+ // irreversible energy
+ STensor3& DirrEnergDF = q1->_DirreversibleEnergyDF;
+ if (Dgamma > 0){
+ static STensor3 DirrEnergDCepr;
+ double dotKSN = dot(KS,N);
+ DirrEnergDCepr = DGDCepr;
+ DirrEnergDCepr *= dotKSN;
+ for (int i=0; i<3; i++){
+ for (int j=0; j<3; j++){
+ for (int k=0; k<3; k++){
+ for (int l=0; l<3; l++){
+ DirrEnergDCepr(i,j) += (Gamma*dKcorDcepr(k,l,i,j)*N(k,l)+ Gamma*KS(k,l)*(DdevNDCepr(k,l,i,j)+_I(k,l)*DtrNDCepr(i,j)/3.));
+ }
+ }
+ }
+ }
+ STensorOperation::multSTensor3STensor43(DirrEnergDCepr,CeprToF,DirrEnergDF);
+ }
+ else{
+ DirrEnergDF *= 0.;
+ }
};
};
@@ -1288,7 +1329,16 @@ void mlawPowerYieldHyper::predictorCorrector_associatedFlow(const STensor3& F, c
for(int l=0; l<3; l++)
P(i,j) += Fe(i,k)*S(k,l)*Fpinv(j,l);
+ // defo energy
q1->_elasticEnergy=deformationEnergy(Ce);
+
+ // dissipation
+ double& irrEnerg = q1->getRefToIrreversibleEnergy();
+ irrEnerg = q0->irreversibleEnergy();
+ if (Dgamma > 0){
+ double dotKSN = dot(KS,N);
+ irrEnerg += Gamma*dotKSN;
+ }
if (stiff){
static STensor3 DpprDCepr;
@@ -1305,8 +1355,10 @@ void mlawPowerYieldHyper::predictorCorrector_associatedFlow(const STensor3& F, c
static STensor43 dFpDCepr;
static STensor3 DgamaDCepr;
-
- if (Dgamma > 0.){
+ static STensor3 DGDCepr;
+ static STensor3 DtrNDCepr;
+ static STensor43 DdevNDCepr;
+ if (Dgamma > 0.){
static STensor3 dg0dCepr;
for (int i=0; i<3; i++)
for (int j=0; j<3; j++){
@@ -1335,7 +1387,7 @@ void mlawPowerYieldHyper::predictorCorrector_associatedFlow(const STensor3& F, c
DsigVMDCepr(i,j) = -invJ(1,0)*DPhigammaDCepr(i,j) - invJ(1,1)*DPhiSigVMDCepr(i,j);
}
- static STensor3 DGDCepr;
+
for (int i=0; i<3; i++)
for (int j=0; j<3; j++){
DGDCepr(i,j) = (g0/(_n*a(1)))*(dzdDgamma*DgamaDCepr(i,j)+dzdsigVM*DsigVMDCepr(i,j)) +
@@ -1343,11 +1395,11 @@ void mlawPowerYieldHyper::predictorCorrector_associatedFlow(const STensor3& F, c
}
- static STensor3 DtrNDCepr;
+
DtrNDCepr = (DgamaDCepr);
DtrNDCepr *= -Da(1);
- static STensor43 DdevNDCepr;
+
DdevNDCepr = (DdevKprDCepr);
double B = a(2)*pow(sigVM,_n-2.);
double u = 1.+3.*Ge*Gamma*_n*B;
@@ -1406,6 +1458,9 @@ void mlawPowerYieldHyper::predictorCorrector_associatedFlow(const STensor3& F, c
else{
DgamaDCepr *= 0.;
dFpDCepr *= 0.;
+ DGDCepr *= 0.;
+ DdevNDCepr *= 0.;
+ DtrNDCepr *= 0.;
}
@@ -1518,6 +1573,27 @@ void mlawPowerYieldHyper::predictorCorrector_associatedFlow(const STensor3& F, c
}
}
}
+ // irreversible energy
+ STensor3& DirrEnergDF = q1->_DirreversibleEnergyDF;
+ if (Dgamma > 0){
+ static STensor3 DirrEnergDCepr;
+ double dotKSN = dot(KS,N);
+ DirrEnergDCepr = DGDCepr;
+ DirrEnergDCepr *= dotKSN;
+ for (int i=0; i<3; i++){
+ for (int j=0; j<3; j++){
+ for (int k=0; k<3; k++){
+ for (int l=0; l<3; l++){
+ DirrEnergDCepr(i,j) += (Gamma*dKcorDcepr(k,l,i,j)*N(k,l)+ Gamma*KS(k,l)*(DdevNDCepr(k,l,i,j)+_I(k,l)*DtrNDCepr(i,j)/3.));
+ }
+ }
+ }
+ }
+ STensorOperation::multSTensor3STensor43(DirrEnergDCepr,CeprToF,DirrEnergDF);
+ }
+ else{
+ DirrEnergDF *= 0.;
+ }
}
};
diff --git a/NonLinearSolver/materialLaw/mlawHyperelastic.h b/NonLinearSolver/materialLaw/mlawHyperelastic.h
index ebd6603f0ed39248e4e74bdeaf0b96ab701816a0..1974c61bcd9ae3ce7a7cbe57f6ea5ad34b5d6c8b 100644
--- a/NonLinearSolver/materialLaw/mlawHyperelastic.h
+++ b/NonLinearSolver/materialLaw/mlawHyperelastic.h
@@ -43,6 +43,7 @@ class mlawHyperViscoElastic : public materialLaw{
std::vector<double> _gi;
protected:
+ virtual double deformationEnergy(const STensor3 &C) const ;
void updateViscoElasticFlow(const IPHyperViscoElastic *q0, IPHyperViscoElastic *q1, double& Ke, double& Ge) const;
void viscoElasticPredictor(const STensor3& Ee, const STensor3& Ee0,
const IPHyperViscoElastic *q0, IPHyperViscoElastic *q1,
@@ -129,9 +130,6 @@ class mlawPowerYieldHyper : public mlawHyperViscoElastic{
protected:
-
- virtual double deformationEnergy(const STensor3 &C) const ;
-
virtual void hardening(IPHyperViscoElastoPlastic* q) const;
virtual void updateEqPlasticDeformation(IPHyperViscoElastoPlastic *q1, const IPHyperViscoElastoPlastic *q0,
diff --git a/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.cpp b/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.cpp
index a1addc68e97d0b80cae656dc21d0e7c1b5724b4c..100ab7cab8e6cac45acf777a1ab620c9de9cea2d 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.cpp
+++ b/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.cpp
@@ -152,6 +152,147 @@ void mlawNonLocalDamagePowerYieldHyper::constitutive(
}
+mlawLocalDamagePowerYieldHyperWithFailure::mlawLocalDamagePowerYieldHyperWithFailure(const int num,const double E,const double nu, const double rho,
+ const double tol, const bool matrixbyPerturbation , const double pert):
+ mlawPowerYieldHyperWithFailure(num,E,nu,rho,tol,matrixbyPerturbation,pert), _saturated(false){};
+
+void mlawLocalDamagePowerYieldHyperWithFailure::clearAllDamageLaw(){
+ for (int i=0; i< damLaw.size(); i++){
+ if (damLaw[i]!= NULL) delete damLaw[i];
+ }
+ damLaw.clear();
+};
+
+void mlawLocalDamagePowerYieldHyperWithFailure::setDamageLaw(const DamageLaw &_damLaw){
+ damLaw.push_back(_damLaw.clone());
+};
+
+mlawLocalDamagePowerYieldHyperWithFailure::mlawLocalDamagePowerYieldHyperWithFailure(const mlawLocalDamagePowerYieldHyperWithFailure &source):
+ mlawPowerYieldHyperWithFailure(source),_saturated(source._saturated){
+ damLaw.clear();
+ for (int i=0; i< source.damLaw.size(); i++){
+ if(source.damLaw[i] != NULL)
+ {
+ damLaw.push_back(source.damLaw[i]->clone());
+ };
+ }
+}
+
+mlawLocalDamagePowerYieldHyperWithFailure::~mlawLocalDamagePowerYieldHyperWithFailure(){
+ for (int i=0; i< damLaw.size(); i++){
+ if (damLaw[i]!= NULL) delete damLaw[i];
+ }
+ damLaw.clear();
+};
+
+void mlawLocalDamagePowerYieldHyperWithFailure::createIPState(IPStateBase* &ips,const bool* state_,const MElement *ele, const int nbFF_, const IntPt *GP, const int gpt) const
+{
+ IPVariable* ipvi = new IPHyperViscoElastoPlasticMultipleLocalDamage(_compression,_traction,_shear,_kinematic,_N,damLaw);
+ IPVariable* ipv1 = new IPHyperViscoElastoPlasticMultipleLocalDamage(_compression,_traction,_shear,_kinematic,_N,damLaw);
+ IPVariable* ipv2 = new IPHyperViscoElastoPlasticMultipleLocalDamage(_compression,_traction,_shear,_kinematic,_N,damLaw);
+ if(ips != NULL) delete ips;
+ ips = new IP3State(state_,ipvi,ipv1,ipv2);
+}
+void mlawLocalDamagePowerYieldHyperWithFailure::createIPState(IPHyperViscoElastoPlasticMultipleLocalDamage *ivi, IPHyperViscoElastoPlasticMultipleLocalDamage *iv1, IPHyperViscoElastoPlasticMultipleLocalDamage *iv2) const
+{
+
+}
+void mlawLocalDamagePowerYieldHyperWithFailure::createIPVariable(IPHyperViscoElastoPlasticMultipleLocalDamage *&ipv,const MElement *ele,const int nbFF,const IntPt *GP, const int gpt) const
+{
+
+}
+
+double mlawLocalDamagePowerYieldHyperWithFailure::soundSpeed() const
+{
+ return mlawPowerYieldHyperWithFailure::soundSpeed();
+}
+
+void mlawLocalDamagePowerYieldHyperWithFailure::constitutive(
+ const STensor3& F0, // initial deformation gradient (input @ time n)
+ const STensor3& Fn, // updated deformation gradient (input @ time n+1)
+ STensor3 &P, // updated 1st Piola-Kirchhoff stress tensor (output)
+ const IPHyperViscoElastoPlasticMultipleLocalDamage *ipvprev, // array of initial internal variable
+ IPHyperViscoElastoPlasticMultipleLocalDamage *ipvcur, // updated array of internal variable (in ipvcur on output),
+ STensor43 &Tangent, // constitutive tangents (output)
+ const bool stiff // if true compute the tangents
+ ) const
+{
+
+ if (isSaturatedHardening()){
+ // saturate when defined critical for damage law
+ if (ipvprev->getDamage(1) >= damLaw[1]->getCriticalDamage()){
+ ipvcur->saturateHardening(ipvprev);
+ }
+ }
+ static STensor43 dFedF, dFpdF;
+ static STensor3 Peff;
+ mlawPowerYieldHyperWithFailure::predictorCorrector(Fn,ipvprev,ipvcur,Peff,stiff,Tangent,dFedF,dFpdF);
+
+ // get result from effective law
+ const STensor3& Fe = ipvcur->_Fe;
+ const STensor3& dgammadF = ipvcur->_DgammaDF;
+ const STensor3 &Fp = ipvcur->getConstRefToFp();
+ const STensor3& dgFdF = ipvcur->getConstRefToDFailurePlasticityDF();
+
+ double ene = ipvcur->defoEnergy();
+ //Msg::Info("enery = %e",ene);
+
+ // saturation damage always develops
+ damLaw[0]->computeDamage(ipvcur->getConstRefToEqPlasticStrain(),
+ ipvprev->getConstRefToEqPlasticStrain(),
+ ene, Fe, Fp, Peff, Cel,
+ ipvprev->getConstRefToIPDamage(0),ipvcur->getRefToIPDamage(0));
+
+
+ if (ipvprev->getDamage(1) < ipvcur->getCriticalDamage()){
+ // set for current critical damage
+ ipvcur->setCriticalDamage(ipvcur->getCriticalDamage());
+ // damage evolution
+ damLaw[1]->computeDamage(ipvcur->getFailurePlasticity(),
+ ipvprev->getFailurePlasticity(),
+ ene, Fe, Fp, Peff, Cel,
+ ipvprev->getConstRefToIPDamage(1),ipvcur->getRefToIPDamage(1));
+ }
+ else{
+ // damage stop increasing
+ STensor3 dDDFe(0.);
+ IPDamage& curDama1 = ipvcur->getRefToIPDamage(1);
+ curDama1.setValues(ipvprev->getDamage(1),ipvcur->getMaximalP(1),0.,0.,dDDFe);
+ }
+
+
+ // computue total softening
+ double D1 = ipvcur->getDamage(0);
+ double D2 = ipvcur->getDamage(1);
+ double D = 1. - (1.- D1)*(1.-D2);
+ P = Peff;
+ P*= (1.- D);
+
+
+ if(stiff)
+ {
+ // we need to correct partial P/partial F: (1-D) partial P/partial F - Peff otimes partial D partial F
+ Tangent*=(1.-D);
+
+ for(int i=0;i<3;i++){
+ for(int j=0;j<3;j++){
+ for(int k=0;k<3;k++){
+ for(int l=0;l<3;l++){
+ Tangent(i,j,k,l)-=Peff(i,j)*(ipvcur->getDDamageDp(0)*dgammadF(k,l)*(1.- D2)+(1.-D1)*ipvcur->getDDamageDp(1)*dgFdF(k,l));
+ for(int m=0; m<3; m++){
+ for(int n=0; n<3; n++){
+ Tangent(i,j,k,l)-=Peff(i,j)*((1.-D2)*ipvcur->getConstRefToDDamageDFe(0)(m,n) +(1.-D1)*ipvcur->getConstRefToDDamageDFe(1)(m,n)) *dFedF(m,n,k,l);
+ }
+ }
+ }
+ }
+ }
+ }
+ }
+}
+
+
+
mlawNonLocalDamagePowerYieldHyperWithFailure::mlawNonLocalDamagePowerYieldHyperWithFailure(const int num,const double E,const double nu, const double rho,
@@ -301,6 +442,20 @@ void mlawNonLocalDamagePowerYieldHyperWithFailure::constitutive(
double D = 1. - (1.- D1)*(1.-D2);
P = Peff;
P*= (1.- D);
+
+ ipvcur->_elasticEnergy = (1-D)*ene;
+
+
+
+ double dIrrPlatic = ipvcur->irreversibleEnergy();
+ dIrrPlatic -= ipvprev->irreversibleEnergy();
+ double Dprev = 1- (1-ipvprev->getDamage(0))*(1-ipvprev->getDamage(1));
+
+ // dissipation
+ double& irrEnerg = ipvcur->getRefToIrreversibleEnergy();
+ irrEnerg = ipvprev->irreversibleEnergy();
+ irrEnerg += ((1-D)*dIrrPlatic + ene*(D - Dprev));
+
if(stiff)
@@ -308,18 +463,12 @@ void mlawNonLocalDamagePowerYieldHyperWithFailure::constitutive(
// we need to correct partial P/partial F: (1-D) partial P/partial F - Peff otimes partial D partial F
Tangent*=(1.-D);
- for(int i=0;i<3;i++)
- {
- for(int j=0;j<3;j++)
- {
- for(int k=0;k<3;k++)
- {
- for(int l=0;l<3;l++)
- {
- for(int m=0; m<3; m++)
- {
- for(int n=0; n<3; n++)
- {
+ for(int i=0;i<3;i++){
+ for(int j=0;j<3;j++){
+ for(int k=0;k<3;k++){
+ for(int l=0;l<3;l++){
+ for(int m=0; m<3; m++){
+ for(int n=0; n<3; n++){
Tangent(i,j,k,l)-=Peff(i,j)*((1.-D2)*ipvcur->getConstRefToDDamageDFe(0)(m,n) +(1.-D1)*ipvcur->getConstRefToDDamageDFe(1)(m,n)) *dFedF(m,n,k,l);
}
}
@@ -343,8 +492,28 @@ void mlawNonLocalDamagePowerYieldHyperWithFailure::constitutive(
// partial p partial tilde p (0 if no MFH)
dLocalVariableDNonLocalVariable.setAll(0.);
-
+
+
+ static STensor3 DIrrPlaticDF;
+ DIrrPlaticDF = ipvcur->getConstRefToDIrreversibleEnergyDF();
+
+ STensor3& DirrEnergDF = ipvcur->getRefToDIrreversibleEnergyDF();
+ for(int i=0;i<3;i++){
+ for(int j=0;j<3;j++){
+ DirrEnergDF(i,j) = (1-D)*DIrrPlaticDF(i,j);
+ for(int k=0;k<3;k++){
+ for(int l=0;l<3;l++){
+ for(int m=0; m<3; m++){
+ DirrEnergDF(i,j) += (D - Dprev)*Peff(k,m)*Fp(l,m)*dFedF(k,l,i,j);
+ };
+ }
+ }
+ }
+ }
+
+ ipvcur->getRefToDIrreversibleEnergyDNonLocalVariable(0) = (ene- dIrrPlatic)*(ipvcur->getDDamageDp(0)*(1.- D2));
+ ipvcur->getRefToDIrreversibleEnergyDNonLocalVariable(1) = (ene- dIrrPlatic)*(ipvcur->getDDamageDp(1)*(1.- D1));
+
}
-
}
diff --git a/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.h b/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.h
index b0d868fc54aae898eb0d617d6091e5e0da442878..bac77a429b5011aa31871ef506fc04cdf45fa868 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.h
+++ b/NonLinearSolver/materialLaw/mlawNonLocalDamageHyperelastic.h
@@ -77,6 +77,54 @@ class mlawNonLocalDamagePowerYieldHyper : public mlawPowerYieldHyper{
#endif // SWIG
};
+class mlawLocalDamagePowerYieldHyperWithFailure : public mlawPowerYieldHyperWithFailure{
+ protected:
+ std::vector<DamageLaw*> damLaw;
+ bool _saturated; // saturate after reaching critical state
+
+ public:
+ mlawLocalDamagePowerYieldHyperWithFailure(const int num,const double E,const double nu, const double rho,
+ const double tol=1.e-6, const bool matrixbyPerturbation = false, const double pert = 1e-8);
+
+ void clearAllDamageLaw();
+ void setDamageLaw(const DamageLaw &_damLaw);
+
+
+ #ifndef SWIG
+ mlawLocalDamagePowerYieldHyperWithFailure(const mlawLocalDamagePowerYieldHyperWithFailure &source);
+ virtual ~mlawLocalDamagePowerYieldHyperWithFailure();
+
+ virtual materialLaw* clone() const { return new mlawLocalDamagePowerYieldHyperWithFailure(*this);}
+
+ bool isSaturatedHardening() const {return _saturated;};
+ void saturatedHardening(const bool falg) {_saturated = falg;};
+
+ // function of materialLaw
+ virtual matname getType() const{return materialLaw::localDamagePowerYieldHyperWithFailure;}
+ virtual bool withDamage() const {return true;};
+ virtual void createIPState(IPHyperViscoElastoPlasticMultipleLocalDamage *ivi, IPHyperViscoElastoPlasticMultipleLocalDamage *iv1, IPHyperViscoElastoPlasticMultipleLocalDamage *iv2) const;
+ virtual void createIPState(IPStateBase* &ips,const bool* state_=NULL,const MElement *ele=NULL, const int nbFF_=0, const IntPt *GP=NULL, const int gpt = 0) const;
+ virtual void createIPVariable(IPHyperViscoElastoPlasticMultipleLocalDamage *&ipv,const MElement *ele,const int nbFF,const IntPt *GP, const int gpt) const;
+
+ virtual void initLaws(const std::map<int,materialLaw*> &maplaw){}; // this law is initialized so nothing to do
+ virtual double soundSpeed() const; // default but you can redefine it for your case
+ virtual const std::vector<DamageLaw*>& getDamageLaw() const {return damLaw; };
+ // specific function
+ public:
+
+ virtual void constitutive(
+ const STensor3& F0, // initial deformation gradient (input @ time n)
+ const STensor3& Fn, // updated deformation gradient (input @ time n+1)
+ STensor3 &P, // updated 1st Piola-Kirchhoff stress tensor (output)
+ const IPHyperViscoElastoPlasticMultipleLocalDamage *q0, // array of initial internal variable
+ IPHyperViscoElastoPlasticMultipleLocalDamage *q1, // updated array of internal variable (in ipvcur on output),
+ STensor43 &Tangent, // constitutive tangents (output)
+ const bool stiff // if true compute the tangents
+ ) const;
+ #endif // SWIG
+};
+
+
class mlawNonLocalDamagePowerYieldHyperWithFailure : public mlawPowerYieldHyperWithFailure{
protected:
std::vector<CLengthLaw*> cLLaw;
diff --git a/NonLinearSolver/materialLaw/mlawNonLocalDamageJ2Hyper.cpp b/NonLinearSolver/materialLaw/mlawNonLocalDamageJ2Hyper.cpp
index 267d29500584eb2abbd63f977d844492e53030df..346c5bd3d71141a8ff6e64a2a5f6c7e89d8c0190 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalDamageJ2Hyper.cpp
+++ b/NonLinearSolver/materialLaw/mlawNonLocalDamageJ2Hyper.cpp
@@ -282,7 +282,19 @@ void mlawNonLocalDamageJ2Hyper::constitutive(
}
else if (this->_pfApproxMethod == mlawJ2linear::Implicit){
double H = ipvcur->getConstRefToIPJ2IsotropicHardening().getDR();
- DdissEnergDF = Peff;
+ for(int i=0;i<3;i++){
+ for(int j=0;j<3;j++){
+ DdissEnergDF(i,j) = 0.;
+ for(int k=0;k<3;k++){
+ for(int l=0;l<3;l++){
+ for (int m=0; m<3; m++){
+ DdissEnergDF(i,j) += Peff(k,m)*Fp(l,m)*dFedF(k,l,i,j);
+ }
+ }
+ }
+ }
+ }
+
DdissEnergDF*= (D-D0);
double DdissDD = ene - Sy*dp;
DdissEnergDF.daxpy(DDamageDF,DdissDD);
diff --git a/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp b/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp
index 07ca904bc810a6d883e83f75eb2399790f5a0ba6..32e6938857eeb42a19b655ad7dcdc810de16da3c 100644
--- a/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp
+++ b/NonLinearSolver/nlsolver/nonLinearMechSolver.cpp
@@ -1058,7 +1058,10 @@ void nonLinearMechSolver::oneStepPreSolvePathFollowing(const double curtime, con
currentPathFollowingIncr = _loadStep;
}
else if (pathsys->getControlType() == pathFollowingSystemBase::LOCAL_CONTROL){
- double ff = double(_numNROptimalLocal)/double(_numNROptimal);
+ double ff = 1.;
+ if (_timeStepAdaptation){
+ ff = double(_numNROptimalLocal)/double(_numNROptimal);
+ }
_localStep = _localStepPrev*_timeStepAdaptationFractor*ff;
currentPathFollowingIncr = _localStep;
}
diff --git a/dG3D/src/dG3DIPVariable.h b/dG3D/src/dG3DIPVariable.h
index 1bfa746e1b782a31c183b11bafd92f77fd290c02..f31dad5de1339cd955c26d1bb46fdbc86d5ef575 100644
--- a/dG3D/src/dG3DIPVariable.h
+++ b/dG3D/src/dG3DIPVariable.h
@@ -362,6 +362,13 @@ public:
if (_ipViscoElastic)
_ipViscoElastic->setLocation(loc);
};
+
+ // for path following based on irreversibe energt
+ virtual double irreversibleEnergy() const {return _ipViscoElastic->irreversibleEnergy();};
+ virtual double& getRefToIrreversibleEnergy() {return _ipViscoElastic->getRefToIrreversibleEnergy();};
+
+ virtual const STensor3& getConstRefToDIrreversibleEnergyDDeformationGradient() const {return _ipViscoElastic->getConstRefToDIrreversibleEnergyDF();};
+ virtual STensor3& getRefToDIrreversibleEnergyDDeformationGradient() {return _ipViscoElastic->getRefToDIrreversibleEnergyDF();};
virtual double defoEnergy() const;
@@ -390,11 +397,16 @@ class HyperViscoElastoPlasticdG3DIPVariable : public dG3DIPVariable{
virtual double get(const int i) const;
+ // for path following based on irreversibe energt
+ virtual double irreversibleEnergy() const {return _ipViscoElastoPlastic->irreversibleEnergy();};
+ virtual double& getRefToIrreversibleEnergy() {return _ipViscoElastoPlastic->getRefToIrreversibleEnergy();};
+
+ virtual const STensor3& getConstRefToDIrreversibleEnergyDDeformationGradient() const {return _ipViscoElastoPlastic->getConstRefToDIrreversibleEnergyDF();};
+ virtual STensor3& getRefToDIrreversibleEnergyDDeformationGradient() {return _ipViscoElastoPlastic->getRefToDIrreversibleEnergyDF();};
+
virtual double defoEnergy() const {return _ipViscoElastoPlastic->defoEnergy();};
virtual double plasticEnergy() const {return _ipViscoElastoPlastic->plasticEnergy();};
- // for path following based on irreversibe energt
- virtual double irreversibleEnergy() const {return _ipViscoElastoPlastic->irreversibleEnergy();};
virtual IPVariable* clone() const {return new HyperViscoElastoPlasticdG3DIPVariable(*this);};
virtual void restart();
diff --git a/dG3D/src/nonLocalDamageHyperelasticDG3DIPVariable.h b/dG3D/src/nonLocalDamageHyperelasticDG3DIPVariable.h
index 162afa477ab0af71990a09dd288b80d95229ac0d..c3aae4e3513f277dd0f5d6bfedc210cbffe691b4 100644
--- a/dG3D/src/nonLocalDamageHyperelasticDG3DIPVariable.h
+++ b/dG3D/src/nonLocalDamageHyperelasticDG3DIPVariable.h
@@ -54,6 +54,32 @@ class nonLocalDamageHyperelasticDG3DIPVariable: public nonLocalDamageDG3DIPVaria
else
Msg::Fatal("the non-local damge id = %d is not defined",idex);
}
+
+ // for path following based on irreversibe energt
+ virtual double irreversibleEnergy() const {return _nldHyperipv->irreversibleEnergy();};
+ virtual double& getRefToIrreversibleEnergy() {return _nldHyperipv->getRefToIrreversibleEnergy();};
+
+ virtual const STensor3& getConstRefToDIrreversibleEnergyDDeformationGradient() const
+ {return _nldHyperipv->getConstRefToDIrreversibleEnergyDF();};
+ virtual STensor3& getRefToDIrreversibleEnergyDDeformationGradient()
+ {return _nldHyperipv->getRefToDIrreversibleEnergyDF();};
+
+ virtual const double& getConstRefToDIrreversibleEnergyDNonLocalVariable(const int index) const
+ {
+ if (index == 0)
+ {
+ return _nldHyperipv->getDIrreversibleEnergyDNonLocalVariable();
+ }
+ else{Msg::Fatal("the non-local variable %d is not defined",index);}
+ };
+ virtual double& getRefToDIrreversibleEnergyDDNonLocalVariable(const int index)
+ {
+ if (index == 0)
+ {
+ return _nldHyperipv->getRefToDIrreversibleEnergyDNonLocalVariable();
+ }
+ else{Msg::Fatal("the non-local variable %d is not defined",index);}
+ };
virtual IPVariable* clone() const {return new nonLocalDamageHyperelasticDG3DIPVariable(*this);};
virtual void restart();
@@ -97,6 +123,26 @@ class nonLocalDamageHyperelasticDG3DIPVariableWithFailure: public nonLocalDamage
{
return _nldHyperipv->getConstRefToCharacteristicLength(idex);
}
+
+ // for path following based on irreversibe energt
+ virtual double irreversibleEnergy() const {return _nldHyperipv->irreversibleEnergy();};
+ virtual double& getRefToIrreversibleEnergy() {return _nldHyperipv->getRefToIrreversibleEnergy();};
+
+ virtual const STensor3& getConstRefToDIrreversibleEnergyDDeformationGradient() const {
+ return _nldHyperipv->getConstRefToDIrreversibleEnergyDF();
+ };
+ virtual STensor3& getRefToDIrreversibleEnergyDDeformationGradient(){
+ return _nldHyperipv->getRefToDIrreversibleEnergyDF();
+ };
+
+ virtual const double& getConstRefToDIrreversibleEnergyDNonLocalVariable(const int index) const
+ {
+ return _nldHyperipv->getDIrreversibleEnergyDNonLocalVariable(index);
+ };
+ virtual double& getRefToDIrreversibleEnergyDDNonLocalVariable(const int index)
+ {
+ return _nldHyperipv->getRefToDIrreversibleEnergyDNonLocalVariable(index);
+ };
virtual IPVariable* clone() const {return new nonLocalDamageHyperelasticDG3DIPVariableWithFailure(*this);};
virtual void restart();
diff --git a/dG3D/src/pathFollowingTerms.cpp b/dG3D/src/pathFollowingTerms.cpp
index deeee14673719822ae7ca4aefc5636cdc8de59c6..17ee48375bdb705f7571c5d4d50c99f7fc3095fd 100644
--- a/dG3D/src/pathFollowingTerms.cpp
+++ b/dG3D/src/pathFollowingTerms.cpp
@@ -441,7 +441,6 @@ void dG3DNonLocalDissipationPathFollowingBulkLinearTerm::get(MElement *ele,int n
int numNonLocalVar = ipv->getNumberNonLocalVariable();
for (int inl =0; inl < numNonLocalVar; inl++){
const double & DdamEnergyDNonLocalVar = ipv->getConstRefToDIrreversibleEnergyDNonLocalVariable(inl);
-
for (int k=0; k< nbFF; k++){
m(k+inl*nbFF+3*nbFF) += DdamEnergyDNonLocalVar*Vals[k+ inl*nbFF+3*nbFF]*detJ*weight;
}