diff --git a/NonLinearSolver/internalPoints/ipHardening.cpp b/NonLinearSolver/internalPoints/ipHardening.cpp
index e4c86895857edcb14cc98cf9da9ab35951558cb1..c4c24885750d061bcc9ce264816d59c8488a0617 100644
--- a/NonLinearSolver/internalPoints/ipHardening.cpp
+++ b/NonLinearSolver/internalPoints/ipHardening.cpp
@@ -249,6 +249,38 @@ IPVariable * IPLinearFollowedByExponentialJ2IsotropicHardening::clone() const
+/* IPLinearFollowedByPowerLawJ2IsotropicHardening */
+IPLinearFollowedByPowerLawJ2IsotropicHardening::IPLinearFollowedByPowerLawJ2IsotropicHardening(): IPJ2IsotropicHardening()
+{
+}
+
+IPLinearFollowedByPowerLawJ2IsotropicHardening::IPLinearFollowedByPowerLawJ2IsotropicHardening(const IPLinearFollowedByPowerLawJ2IsotropicHardening &source) : IPJ2IsotropicHardening(source)
+{
+}
+
+IPLinearFollowedByPowerLawJ2IsotropicHardening &IPLinearFollowedByPowerLawJ2IsotropicHardening::operator=(const IPVariable &source)
+{
+ IPJ2IsotropicHardening::operator=(source);
+ const IPLinearFollowedByPowerLawJ2IsotropicHardening* src = dynamic_cast<const IPLinearFollowedByPowerLawJ2IsotropicHardening*>(&source);
+ if(src != NULL)
+ {
+ }
+ return *this;
+}
+
+void IPLinearFollowedByPowerLawJ2IsotropicHardening::restart()
+{
+ IPJ2IsotropicHardening::restart();
+}
+
+IPVariable * IPLinearFollowedByPowerLawJ2IsotropicHardening::clone() const
+{
+ return new IPLinearFollowedByPowerLawJ2IsotropicHardening(*this);
+}
+
+
+
+
IPPolynomialJ2IsotropicHardening::IPPolynomialJ2IsotropicHardening(): IPJ2IsotropicHardening()
{
diff --git a/NonLinearSolver/internalPoints/ipHardening.h b/NonLinearSolver/internalPoints/ipHardening.h
index d4f41c209bbf608fd86b8721b786f013757c120a..c586ae0f4c8dacc6ed61a1cb7c28a9a2a7daaeda 100644
--- a/NonLinearSolver/internalPoints/ipHardening.h
+++ b/NonLinearSolver/internalPoints/ipHardening.h
@@ -160,6 +160,7 @@ class IPLinearExponentialJ2IsotropicHardening : public IPJ2IsotropicHardening
};
+
/* IPLinearFollowedByExponentialJ2IsotropicHardening */
class IPLinearFollowedByExponentialJ2IsotropicHardening : public IPJ2IsotropicHardening
{
@@ -177,6 +178,23 @@ class IPLinearFollowedByExponentialJ2IsotropicHardening : public IPJ2IsotropicHa
+/* IPLinearFollowedByPowerLawJ2IsotropicHardening */
+class IPLinearFollowedByPowerLawJ2IsotropicHardening : public IPJ2IsotropicHardening
+{
+
+ protected:
+ public:
+ IPLinearFollowedByPowerLawJ2IsotropicHardening();
+ IPLinearFollowedByPowerLawJ2IsotropicHardening(const IPLinearFollowedByPowerLawJ2IsotropicHardening &source);
+ virtual IPLinearFollowedByPowerLawJ2IsotropicHardening &operator=(const IPVariable &source);
+ virtual ~IPLinearFollowedByPowerLawJ2IsotropicHardening(){}
+ virtual void restart();
+ virtual IPVariable *clone() const;
+
+};
+
+
+
class IPPolynomialJ2IsotropicHardening : public IPJ2IsotropicHardening{
public:
IPPolynomialJ2IsotropicHardening();
diff --git a/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.cpp b/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.cpp
index 2a2ff8905dee412b9cbf55f9dcbec55f12fa1a50..1a0dd8c8ecd752353326ccf311bab2f85b6d4414 100644
--- a/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.cpp
+++ b/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.cpp
@@ -24,7 +24,7 @@ IPNonLocalDamageGurson::IPNonLocalDamageGurson() : IPVariableMechanics(), _elas
_Fp(0,0)=1.;
_Fp(1,1)=1.;
_Fp(2,2)=1.;
- _DFvStarDStrain=0;
+ _DFvStarDStrain=0.;
_logsqrtCe = 0.;
};
@@ -32,8 +32,11 @@ IPNonLocalDamageGurson::IPNonLocalDamageGurson(double fVinitial, const J2Isotrop
const CLengthLaw *cll,
const std::vector<GursonDamageNucleation*> *gdnLawContainer) :
IPVariableMechanics(), _elasticEnergy(0.),
- _eplmatrix(0.), _fV(fVinitial), _fVstar(fVinitial),
- _nldfVstar(fVinitial), _DFvStarDNldfVstar(1.),_dissipationBlocked(false)
+ _eplmatrix(0.), _fV(fVinitial),
+ _fVstar(fVinitial),_nldfVstar(fVinitial),
+ _DFvStarDNldfVstar(1.),
+ _dissipationBlocked(false),
+ _NonLocalToLocal(false)
{
ipvJ2IsotropicHardening=NULL;
if(j2IH ==NULL) Msg::Error("IPNonLocalDamageGurson::IPNonLocalDamageGurson has no j2IH");
@@ -54,14 +57,14 @@ IPNonLocalDamageGurson::IPNonLocalDamageGurson(double fVinitial, const J2Isotrop
if(fVinitial <0. or fVinitial >1.) Msg::Error("IPNonLocalDamageGurson::IPNonLocalDamageGurson wrong initial porosity");
- _Fp=0.;
+ _Fp = 0.;
_Fp(0,0)=1.;
_Fp(1,1)=1.;
_Fp(2,2)=1.;
- _logsqrtCe=0.;
+ _logsqrtCe = 0.;
- _DFvStarDStrain=0;
+ _DFvStarDStrain = 0.;
};
@@ -93,6 +96,7 @@ IPNonLocalDamageGurson::IPNonLocalDamageGurson(const IPNonLocalDamageGurson &sou
_DFvStarDStrain=source._DFvStarDStrain;
_DFvStarDNldfVstar=source._DFvStarDNldfVstar;
_dissipationBlocked = source._dissipationBlocked;
+ _NonLocalToLocal = source._NonLocalToLocal;
}
IPNonLocalDamageGurson& IPNonLocalDamageGurson::operator=(const IPVariable &source)
@@ -154,26 +158,12 @@ IPNonLocalDamageGurson& IPNonLocalDamageGurson::operator=(const IPVariable &sour
_DFvStarDStrain=src->_DFvStarDStrain;
_DFvStarDNldfVstar=src->_DFvStarDNldfVstar;
_dissipationBlocked = src->_dissipationBlocked;
+ _NonLocalToLocal = src->_NonLocalToLocal;
}
return *this;
}
-double IPNonLocalDamageGurson::defoEnergy() const
-{
- return _elasticEnergy;
-}
-
-
-double IPNonLocalDamageGurson::plasticEnergy() const
-{
- if(ipvJ2IsotropicHardening != NULL)
- return ipvJ2IsotropicHardening->getIntegR();
- else
- return 0.;
-}
-
-
void IPNonLocalDamageGurson::restart()
{
IPVariableMechanics::restart();
@@ -190,6 +180,7 @@ void IPNonLocalDamageGurson::restart()
restartManager::restart(_logsqrtCe);
restartManager::restart(_DFvStarDStrain);
restartManager::restart(_dissipationBlocked);
+ restartManager::restart(_NonLocalToLocal);
return;
}
diff --git a/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.h b/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.h
index c1d1798dfcb7b67f2ead49924ad59fe66513e369..21d8f239230f0222b31893bd4db9747d1dc1e792 100644
--- a/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.h
+++ b/NonLinearSolver/internalPoints/ipNonLocalDamageGurson.h
@@ -25,9 +25,6 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
IPJ2IsotropicHardening* ipvJ2IsotropicHardening; // Hardening law
std::vector<IPGursonNucleation*> ipvgdnContainer; // Nucleation law vector
- // Energies
- double _elasticEnergy; // elastic energy stored
-
// Internal variables
double _eplmatrix; // equivalent plastic strain in the matrix (hatp)
double _fV; // local porosity
@@ -42,7 +39,14 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
double _DFvStarDNldfVstar; //derivative of the local corrected porosity with the non local corrected porosity
// Damage and transition managing
- bool _dissipationBlocked;
+ bool _dissipationBlocked; // True if dissipation is blocked at the IPv
+ bool _NonLocalToLocal; // True if transition from local to non-local allowed
+
+ // Energies
+ double _elasticEnergy; // elastic energy stored
+ // plastic energy : stored in hardening law
+
+
public:
@@ -65,10 +69,7 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
}
}
-
// General functions
- virtual double defoEnergy() const;
- virtual double plasticEnergy() const;
virtual void restart();
virtual IPVariable* clone() const {return new IPNonLocalDamageGurson(*this);};
@@ -79,23 +80,31 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
virtual void getValuesFromMPI(const double *arrayMPI){_fV=arrayMPI[0];};
#endif // HAVE_MPI
-
+
+
// Damage managing
virtual void blockDissipation(const bool fl){_dissipationBlocked = fl;};
virtual bool dissipationIsBlocked() const {return _dissipationBlocked;};
+
- // Access functions - members values
+ // Access functions - Energy
virtual double &getRefToElasticEnergy() {return _elasticEnergy;};
-
+ virtual double defoEnergy() const { return _elasticEnergy;};
+ virtual double plasticEnergy() const
+ {
+ if(ipvJ2IsotropicHardening != NULL) {return ipvJ2IsotropicHardening->getIntegR();}
+ else {return 0.;}
+ };
+ // Those need to be implemented ???
+ virtual double damageEnergy() const {return 0.;}; // dissipation by damage
+ virtual double irreversibleEnergy() const {return 0.;} // irreversible energy for path following
+
+
+
+ // Access functions - Internal variables
virtual double getDamage() const{return getLocalCorrectedPorosity();};
- virtual const STensor3& getConstRefToFp() const {return _Fp;}
- virtual STensor3& getRefToFp(){return _Fp;}
-
- virtual const STensor3& getConstRefToLogsqrtCe() const{return _logsqrtCe;}
- virtual STensor3& getRefToLogsqrtCe(){return _logsqrtCe;}
-
virtual double getMatrixPlasticStrain() const { return _eplmatrix;}
virtual double &getRefToMatrixPlasticStrain() { return _eplmatrix;}
@@ -107,7 +116,16 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
virtual double getNonLocalCorrectedPorosity() const { return _nldfVstar;}
virtual double &getRefToNonLocalCorrectedPorosity() { return _nldfVstar;}
-
+
+ virtual const STensor3& getConstRefToFp() const {return _Fp;}
+ virtual STensor3& getRefToFp(){return _Fp;}
+
+ virtual const STensor3& getConstRefToLogsqrtCe() const{return _logsqrtCe;}
+ virtual STensor3& getRefToLogsqrtCe(){return _logsqrtCe;}
+
+
+
+ // Access functions - Derivatives
virtual double &getRefToDLocalCorrectedPorosityDNonLocalCorrectedPorosity(){return _DFvStarDNldfVstar;}
virtual double getDLocalCorrectedPorosityDNonLocalCorrectedPorosity() const{return _DFvStarDNldfVstar;}
@@ -115,6 +133,7 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
virtual STensor3 & getRefToDLocalCorrectedPorosityDStrain(){return _DFvStarDStrain;}
+
// Access functions - Included IPVs
virtual const IPJ2IsotropicHardening &getConstRefToIPJ2IsotropicHardening() const
{
@@ -176,7 +195,6 @@ class IPNonLocalDamageGurson : public IPVariableMechanics
{
return ipvgdnContainer;
}
-
};
#endif // IPNONLOCALDAMAGEGURSON_H_
diff --git a/NonLinearSolver/materialLaw/j2IsotropicHardening.cpp b/NonLinearSolver/materialLaw/j2IsotropicHardening.cpp
index b8b8b3e8f5aa04acc7f3201a25dd2bdbe096cddb..281ee14933c31049df9c72b2e35b188a62081b4b 100644
--- a/NonLinearSolver/materialLaw/j2IsotropicHardening.cpp
+++ b/NonLinearSolver/materialLaw/j2IsotropicHardening.cpp
@@ -315,7 +315,7 @@ J2IsotropicHardening * LinearExponentialJ2IsotropicHardening::clone() const
-/* LinearExponentialJ2IsotropicHardening */
+/* LinearFollowedByExponentialJ2IsotropicHardening */
LinearFollowedByExponentialJ2IsotropicHardening::LinearFollowedByExponentialJ2IsotropicHardening(const int num, double yield0,
double h1, double pexp, double h2, double hexp2, bool init) :
J2IsotropicHardening(num,yield0,init), _h1(h1), _pexp(pexp), _h2(h2), _hexp2(hexp2)
@@ -361,7 +361,7 @@ void LinearFollowedByExponentialJ2IsotropicHardening::hardening(double p, IPJ2Is
dR = _h1;
ddR = 0.;
}
- else if (p < _pexp) // Plastic case: linear part
+ else if (p <= _pexp) // Plastic case: linear part
{
R += _h1*p;
dR = _h1;
@@ -371,13 +371,12 @@ void LinearFollowedByExponentialJ2IsotropicHardening::hardening(double p, IPJ2Is
else // Plastic case: exponentional (saturation) part
{
double tmp = exp(-(p-_pexp)/_hexp2);
- R += _h1*_pexp + _h2*(1.-tmp);
- dR = _h2 * tmp /_hexp2;
- ddR = - dR /_hexp2;
- // ddR = _h2*exp(-(p-_pexp)/_hexp) * (p-_pexp)/_hexp * -(p-_pexp)/_hexp + _h2*exp(-(p-_pexp)/_hexp) * p/_hexp;
+ R += _h1*p + _h2*(1.-tmp);
+ dR = _h1 + _h2 * tmp /_hexp2;
+ ddR = -_h2 * tmp /_hexp2 /_hexp2;
intR = getYield0()*p;
- intR += 0.5*_h1*_pexp*_pexp + _h1*_pexp*(p-_pexp);
+ intR += 0.5*_h1*p*p;
intR += _h2*(p-_pexp) + _h2*(tmp-exp(0.)) *_hexp2;
}
ipv.set(R,dR,ddR,intR);
@@ -391,6 +390,89 @@ J2IsotropicHardening* LinearFollowedByExponentialJ2IsotropicHardening::clone() c
+
+/* LinearFollowedByPowerLawJ2IsotropicHardening */
+LinearFollowedByPowerLawJ2IsotropicHardening::LinearFollowedByPowerLawJ2IsotropicHardening(const int num, double yield0, double h, double pexp, double hexp, bool init):
+ J2IsotropicHardening(num,yield0,init), _h1(h), _pexp(pexp), _hexp(hexp)
+{
+ if(_h1 < 0. or _hexp < 0.) Msg::Error("LinearFollowedByPowerLawJ2IsotropicHardening:: negative hardening parameters");
+ if(_pexp <= 1.e-9 ) Msg::Error("LinearFollowedByPowerLawJ2IsotropicHardening:: negative or too small hardening parameters p_exp, use PowerLawJ2IsotropicHardening instead");
+}
+
+LinearFollowedByPowerLawJ2IsotropicHardening::LinearFollowedByPowerLawJ2IsotropicHardening(const LinearFollowedByPowerLawJ2IsotropicHardening &source) :
+ J2IsotropicHardening(source)
+{
+ _h1 = source._h1;
+ _pexp = source._pexp;
+ _hexp = source._hexp;
+}
+
+LinearFollowedByPowerLawJ2IsotropicHardening& LinearFollowedByPowerLawJ2IsotropicHardening::operator=(const J2IsotropicHardening &source)
+{
+ J2IsotropicHardening::operator=(source);
+ const LinearFollowedByPowerLawJ2IsotropicHardening* src = dynamic_cast<const LinearFollowedByPowerLawJ2IsotropicHardening*>(&source);
+ if(src != NULL)
+ {
+ _h1 = src->_h1;
+ _pexp = src->_pexp;
+ _hexp = src->_hexp;
+ }
+ return *this;
+}
+void LinearFollowedByPowerLawJ2IsotropicHardening::createIPVariable(IPJ2IsotropicHardening* &ipv) const
+{
+ if(ipv != NULL) delete ipv;
+ ipv = new IPLinearFollowedByPowerLawJ2IsotropicHardening();
+}
+
+void LinearFollowedByPowerLawJ2IsotropicHardening::hardening(double p, IPJ2IsotropicHardening &ipv) const
+{
+ double tol=1.e-16;
+ double dR=0., ddR=0., intR=0.;
+ double R = getYield0();
+ // elastic case
+ if(p < tol)
+ {
+ // why ????
+ if(_h1<1.)
+ {
+ dR = 1e20;
+ ddR = -1e20;
+ }
+ else
+ {
+ dR = _h1;
+ ddR = 0.;
+ }
+ }
+ // linear part
+ else if (p <= _pexp)
+ {
+ R = getYield0() + _h1*p;
+ dR = _h1;
+ ddR = 0.;
+ intR = getYield0()*p + 0.5*_h1*p*p;
+ }
+ // power part
+ else
+ {
+ R = (getYield0()+_h1*_pexp) * pow(p/_pexp,_hexp);
+ dR = R*_hexp/p;
+ ddR = dR*(_hexp-1.)/p;
+
+ intR = getYield0()*_pexp + 0.5*_h1*_pexp*_pexp;
+ intR += R*p /(_hexp+1.) -((getYield0()+_h1*_pexp)*_pexp/(_hexp+1.));
+ }
+ ipv.set(R,dR,ddR,intR);
+}
+J2IsotropicHardening * LinearFollowedByPowerLawJ2IsotropicHardening::clone() const
+{
+ return new LinearFollowedByPowerLawJ2IsotropicHardening(*this);
+}
+
+
+
+/* PolynomialJ2IsotropicHardening*/
PolynomialJ2IsotropicHardening::PolynomialJ2IsotropicHardening(const int num, double yield0, int order, bool init):
J2IsotropicHardening(num,yield0,init),_order(order){
_coefficients.resize(_order+1);
diff --git a/NonLinearSolver/materialLaw/j2IsotropicHardening.h b/NonLinearSolver/materialLaw/j2IsotropicHardening.h
index df27ae6da533ac284cb2715dd7febd06863edf5a..1b209d5c65b466fbf3548a3927b003d038103c41 100644
--- a/NonLinearSolver/materialLaw/j2IsotropicHardening.h
+++ b/NonLinearSolver/materialLaw/j2IsotropicHardening.h
@@ -20,7 +20,8 @@ class J2IsotropicHardening{
public :
enum hardeningname{perfectlyPlasticJ2IsotropicHardening,powerLawJ2IsotropicHardening, exponentialJ2IsotropicHardening,
swiftJ2IsotropicHardening, linearExponentialJ2IsotropicHardening,
- linearFollowedByExponentialJ2IsotropicHardening, polynomialJ2IsotropicHardening,
+ linearFollowedByExponentialJ2IsotropicHardening, linearFollowedByPowerLawJ2IsotropicHardening,
+ polynomialJ2IsotropicHardening,
twoExpJ2IsotropicHaderning, tanhJ2IsotropicHaderning};
protected :
int _num; // number of law (must be unique !)
@@ -155,7 +156,7 @@ class LinearFollowedByExponentialJ2IsotropicHardening : public J2IsotropicHarden
{
// R = yield0 + h1 * p
// Then when p > p_exp
- // R = yield0 + h1 * pexp + h2* (1- exp(-hexp2*(p-pexp)))
+ // R = yield0 + h1 * p + h2* (1- exp(-hexp2*(p-pexp)))
protected :
double _h1, _pexp, _h2, _hexp2;
@@ -176,6 +177,31 @@ class LinearFollowedByExponentialJ2IsotropicHardening : public J2IsotropicHarden
};
+class LinearFollowedByPowerLawJ2IsotropicHardening : public J2IsotropicHardening
+{
+ // R = yield0 + h1 * p
+ // Then when p > pexp
+ // R = yield0 * (1 + h1 * pexp) * (p/pexp) ^ h_exp
+ protected :
+ double _h1, _pexp, _hexp;
+
+ public:
+ // constructor
+ LinearFollowedByPowerLawJ2IsotropicHardening(const int num, double yield0, double h, double pexp, double hexp, bool init=true);
+#ifndef SWIG
+ virtual ~LinearFollowedByPowerLawJ2IsotropicHardening(){}
+ LinearFollowedByPowerLawJ2IsotropicHardening(const LinearFollowedByPowerLawJ2IsotropicHardening &source);
+ LinearFollowedByPowerLawJ2IsotropicHardening& operator=(const J2IsotropicHardening &source);
+ virtual int getNum() const{return _num;}
+ virtual hardeningname getType() const{return J2IsotropicHardening::linearFollowedByPowerLawJ2IsotropicHardening; };
+ virtual void createIPVariable(IPJ2IsotropicHardening* &ipv) const;
+ virtual void initLaws(const std::map<int,J2IsotropicHardening*> &maplaw) {}; //nothing now, we will see if we use the mapping
+ virtual void hardening(double p, IPJ2IsotropicHardening &ipv) const;
+ virtual J2IsotropicHardening * clone() const;
+#endif
+};
+
+
class PolynomialJ2IsotropicHardening : public J2IsotropicHardening{
#ifndef SWIG
protected: