diff --git a/NonLinearSolver/Domain/partDomain.cpp b/NonLinearSolver/Domain/partDomain.cpp
index 68cc85828453fa23fb5c71ec759628b30c6c73fe..a02b0b9bbcd275ea413d445d70b26b95290ce54e 100644
--- a/NonLinearSolver/Domain/partDomain.cpp
+++ b/NonLinearSolver/Domain/partDomain.cpp
@@ -338,34 +338,6 @@ void partDomain::createIPMap(){
}
}
};
-
-void partDomain::getBlockDissipationIps(std::vector<int>& blockDissipationIPs) const{
- blockDissipationIPs.reserve(blockDissipationIPs.size()+_blockedDissipationPs.size());
- for (std::set<int>::const_iterator it = _blockedDissipationPs.begin(); it != _blockedDissipationPs.end(); it++){
- blockDissipationIPs.push_back(*it);
- }
-};
-void partDomain::getBrokenIPs(std::vector<int>& brkIps) const{
- brkIps.reserve(brkIps.size() + _brokenIPs.size());
- for (std::set<int>::const_iterator it = _brokenIPs.begin(); it != _brokenIPs.end(); it++){
- brkIps.push_back(*it);
- }
-};
-
-void partDomain::blockDissipationIP(const int ip){
- _blockedDissipationPs.insert(ip);
-};
-
-void partDomain::clearBlockDissipationIPs(){
- _blockedDissipationPs.clear();
-}
-void partDomain::brokenIP(const int ip){
- _brokenIPs.insert(ip);
-};
-void partDomain::clearBrokenIPs(){
- _brokenIPs.clear();
-}
-
#endif //HAVE_MPI
void partDomain::initMicroMeshId(){
@@ -426,7 +398,7 @@ void partDomain::computeActiveDissipationAverageStressIncrement(const IPField* i
for (groupOfElements::elementContainer::iterator it= g->begin(); it!=g->end(); it++){
MElement* ele= *it;
- AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ele->getNum());
int npts= integBulk->getIntPoints(ele,&GP);
for (int i=0; i<npts; i++){
@@ -510,7 +482,7 @@ void partDomain::computeActiveDissipationCenter(const IPField* ipf, SPoint3& pos
// compute from damaging zone
for (groupOfElements::elementContainer::iterator it= g->begin(); it!=g->end(); it++){
MElement* ele= *it;
- AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ele->getNum());
int npts= integBulk->getIntPoints(ele,&GP);
for (int i=0; i<npts; i++){
const IPStateBase *ips = (*vips)[i];
@@ -539,7 +511,7 @@ void partDomain::computeActiveDissipationAverageStrainIncrement(const IPField* i
// compute from damaging zone
for (groupOfElements::elementContainer::iterator it= g->begin(); it!=g->end(); it++){
MElement* ele= *it;
- AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ele->getNum());
int npts= integBulk->getIntPoints(ele,&GP);
for (int i=0; i<npts; i++){
const IPStateBase *ips = (*vips)[i];
diff --git a/NonLinearSolver/Domain/partDomain.h b/NonLinearSolver/Domain/partDomain.h
index 40d1d3fc51452e933fb9d379f171f7e6dc23fb79..602f82ecc5b3d61bcf12009e7bef734b78d37e69 100644
--- a/NonLinearSolver/Domain/partDomain.h
+++ b/NonLinearSolver/Domain/partDomain.h
@@ -64,10 +64,6 @@ class partDomain{
std::set<int> _otherRanks;
std::set<int> _domainIPBulk; // bulk integration points lie on this domain
std::map<int,std::set<int> > _mapIPBulk; // map IP for bulk elements
-
- // for damage blockage
- std::set<int> _blockedDissipationPs;
- std::set<int> _brokenIPs;
#endif //HAVE_MPI
bool _distributedOtherRanks;
@@ -179,7 +175,7 @@ public:
int getPhysical() const{return _phys;}
//check cohesive element insertion, failure
- virtual void checkFailure(AllIPState* aips) = 0;
+ virtual void checkFailure(IPField* ipf) const = 0;
virtual void initializeTerms(unknownField *uf,IPField *ip)=0;
// alocate the object (new) but is not responsiblefor the delete (the NeumannBC is in charged of the delete)
@@ -294,13 +290,6 @@ public:
virtual void createIPMap();
virtual const std::map<int,std::set<int> >& getMapIPBulk() const {return _mapIPBulk;};
virtual const std::set<int>& getDomainIPBulk() const {return _domainIPBulk;};
-
- virtual void getBlockDissipationIps(std::vector<int>& blockDissipationIPs) const;
- virtual void getBrokenIPs(std::vector<int>& brkIps) const;
- virtual void blockDissipationIP(const int ip);
- virtual void clearBlockDissipationIPs();
- virtual void brokenIP(const int ip);
- virtual void clearBrokenIPs();
#endif //HAVE_MPI
virtual bool considerForTimeStep(AllIPState *aips,MElement *e) {return true;};
virtual double computeMeanValueDomain(const int num, const IPField* ipf, double& volume) const;
diff --git a/NonLinearSolver/internalPoints/ipField.cpp b/NonLinearSolver/internalPoints/ipField.cpp
index 693882a1896f1fc42568769f9be325a024a82b1e..520074ffe3e768b8bc459780e04cc52d39308c98 100644
--- a/NonLinearSolver/internalPoints/ipField.cpp
+++ b/NonLinearSolver/internalPoints/ipField.cpp
@@ -230,12 +230,21 @@ void IPField::compute1state(IPStateBase::whichState ws, bool stiff){
}
void IPField::checkFailure(){
+ #if defined(HAVE_MPI)
+ if (_solver->isMultiscale() and _solver->rankOnSolver(Msg::GetCommRank())){
+ // get data in order to send to other procs
+ // only new block and broken ips are filled in the group
+ this->clearBlockDissipationIPs();
+ this->clearBrokenIPs();
+ }
+ #endif //MPI
+
//Msg::Info("checking failure all IP"); //too many output in explicit
std::vector<partDomain*>* domainVector = _solver->getDomainVector();
// check failure, cohesive insertion at current time step
for (int idom=0; idom < domainVector->size(); idom++){
partDomain* dom = (*domainVector)[idom];
- dom->checkFailure(_AIPS);
+ dom->checkFailure(this);
}
#if defined(HAVE_MPI)
@@ -248,13 +257,9 @@ void IPField::checkFailure(){
std::set<int> otherRanks;
_solver->getOtherRanks(Msg::GetCommRank(),otherRanks);
if (otherRanks.size() > 0){
- std::vector<int> brokenIPs;
- std::vector<int> blockDissipationIPs;
- for (int idom=0; idom< domainVector->size(); idom++){
- (*domainVector)[idom]->getBrokenIPs(brokenIPs);
- (*domainVector)[idom]->getBlockDissipationIps(blockDissipationIPs);
- }
- // send data
+ std::vector<int> brokenIPs(_brokenIPs.begin(),_brokenIPs.end());
+ std::vector<int> blockDissipationIPs(_blockedDissipationPs.begin(),_blockedDissipationPs.end());
+ // send data
for (std::set<int>::const_iterator it = otherRanks.begin(); it != otherRanks.end(); it++){
int otherRank = *it;
int sizeIPData[2];
@@ -1580,3 +1585,47 @@ void IPField::resetFile(){
viparch[i].openFile();
}
};
+
+
+void IPField::blockDissipationBulkElement(const int elnum){
+ _blockedDissipationBulkElements.insert(elnum);
+};
+
+bool IPField::dissipationIsBlockedBulkElement(const int elnum) const{
+ return (_blockedDissipationBulkElements.find(elnum) != _blockedDissipationBulkElements.end());
+};
+
+void IPField::brokenBulkElement(const int elnum){
+ _brokenBulkElements.insert(elnum);
+};
+
+bool IPField::bulkElementIsBroken(const int elnum) const{
+ return (_brokenBulkElements.find(elnum) != _brokenBulkElements.end());
+};
+
+void IPField::brokenInterfaceElement(const int el1, const int el2){
+ TwoNum key(el1,el2);
+ _brokenInterfaceElements.insert(key);
+};
+
+bool IPField::interfaceElementIsBroken(const int el1, const int el2) const{
+ TwoNum key(el1,el2);
+ return (_brokenInterfaceElements.find(key) != _brokenInterfaceElements.end());
+};
+
+#if defined(HAVE_MPI)
+void IPField::blockDissipationIP(const int ip){
+ _blockedDissipationPs.insert(ip);
+};
+
+void IPField::clearBlockDissipationIPs(){
+ _blockedDissipationPs.clear();
+}
+void IPField::brokenIP(const int ip){
+ _brokenIPs.insert(ip);
+};
+void IPField::clearBrokenIPs(){
+ _brokenIPs.clear();
+}
+
+#endif //HAVE_MPI
\ No newline at end of file
diff --git a/NonLinearSolver/internalPoints/ipField.h b/NonLinearSolver/internalPoints/ipField.h
index 0453385d2fdadcf804742575193ea2663d7be469..88b06e28350769f2c7604544fea7958787b3cc2b 100644
--- a/NonLinearSolver/internalPoints/ipField.h
+++ b/NonLinearSolver/internalPoints/ipField.h
@@ -11,6 +11,8 @@
//
# ifndef _IPFIELD_H_
# define _IPFIELD_H_
+
+#ifndef SWIG
#include<vector>
#include"quadratureRules.h"
#include "unknownField.h"
@@ -18,7 +20,41 @@
#include "GModel.h"
#include "ipstate.h"
#include "partDomain.h"
-struct archiveIPVariable;
+#include "InterfaceKeys.h"
+
+class TwoNum{
+ public:
+ int small, large;
+ TwoNum(const int e1, const int e2){
+ if (e1 < e2) {
+ small = e1;
+ large = e2;
+ }
+ else{
+ small = e2;
+ large = e1;
+ }
+ }
+ TwoNum(const TwoNum& src): small(src.small),large(src.large){}
+ ~TwoNum(){};
+ TwoNum& operator = (const TwoNum& other) {
+ small = other.small;
+ large = other.large;
+ return *this;
+ };
+ bool operator < (const TwoNum &other) const{
+ if (small < other.small) return true;
+ if (small > other.small) return false;
+ if (large < other.large) return true;
+ return false;
+ }
+ bool operator == (const TwoNum &other) const{
+ return (small == other.small && large == other.large);
+ }
+
+};
+#endif //SWIG
+
class IPField : public elementsField {
public:
enum Output { SVM=0,SIG_XX, SIG_YY, SIG_ZZ, SIG_XY, SIG_YZ, SIG_XZ, DAMAGE, PLASTICSTRAIN, GL_EQUIVALENT_STRAIN,
@@ -183,6 +219,16 @@ class IPField : public elementsField {
nonLinearMechSolver* _solver; // solver to which ipField belongs
int _numOfActiveDissipationIPs;
int _numOfActiveDissipationIPsReached;
+
+ std::set<TwoNum> _brokenInterfaceElements; //all interface broken element is stored in order to avoid rechecking
+ std::set<int> _brokenBulkElements; // all broken element is stored in order to avoid rechecking
+ std::set<int> _blockedDissipationBulkElements; // all broken element is stored in order to avoid rechecking
+
+ #if defined(HAVE_MPI)
+ // for damage blockage
+ std::set<int> _blockedDissipationPs;
+ std::set<int> _brokenIPs;
+ #endif //HAVE_MPI
public :
IPField(nonLinearMechSolver* sl, std::vector<ip2archive> &vaip,
@@ -191,7 +237,9 @@ class IPField : public elementsField {
static std::string ToString(const int i);
- AllIPState* getAips() const {return _AIPS;}
+ AllIPState* getAips() {return _AIPS;}
+ const AllIPState* getAips() const {return _AIPS;}
+
const std::vector<partDomain*>& getDomainVector() const;
virtual ~IPField(){delete _AIPS;}
@@ -370,6 +418,23 @@ class IPField : public elementsField {
virtual void checkCritialInterfaceIP();
virtual void resetFile();
+
+ void blockDissipationBulkElement(const int elnum); // bulk element
+ bool dissipationIsBlockedBulkElement(const int elnum) const;
+
+ void brokenBulkElement(const int elnum);
+ bool bulkElementIsBroken(const int elnum) const;
+
+ void brokenInterfaceElement(const int elnum1, const int elnum2);
+ bool interfaceElementIsBroken(const int elnum1, const int elnum2) const;
+
+ #if defined(HAVE_MPI)
+ void blockDissipationIP(const int ip);
+ void clearBlockDissipationIPs();
+ void brokenIP(const int ip);
+ void clearBrokenIPs();
+ #endif //HAVE_MPI
+
#endif // swig
};
#endif // IPField
diff --git a/NonLinearSolver/internalPoints/ipstate.cpp b/NonLinearSolver/internalPoints/ipstate.cpp
index 5d66e692a4e88177c83f1dd2f2ae48b8335ea93a..f595d8ddd3f0cf114639852d142982eabfe2f54f 100644
--- a/NonLinearSolver/internalPoints/ipstate.cpp
+++ b/NonLinearSolver/internalPoints/ipstate.cpp
@@ -461,6 +461,27 @@ AllIPState::ipstateElementContainer* AllIPState::operator[](const long int num)
return &(_mapall.find(num)->second);
}
+const AllIPState::ipstateElementContainer* AllIPState::getIPstate(const long int num) const
+{
+
+ ipstateContainer::const_iterator it = _mapall.find(num);
+ #ifdef _DEBUG
+ if(it == _mapall.end()){
+ Msg::Error("Try to get an inexisting IPState on rank %d",Msg::GetCommRank());
+ return NULL;
+ }
+ else
+ #endif // _DEBUG
+ {
+ return &(it->second);
+ }
+}
+
+const AllIPState::ipstateElementContainer* AllIPState::operator[](const long int num) const
+{
+ return &(_mapall.find(num)->second);
+}
+
void AllIPState::restart()
{
for(ipstateContainer::iterator its=_mapall.begin(); its!=_mapall.end();++its)
diff --git a/NonLinearSolver/internalPoints/ipstate.h b/NonLinearSolver/internalPoints/ipstate.h
index 6e771cf120c9f4a747ff6e0b8fe7404d428d1861..3799fd7618883c41277bd728acad320e690c4537 100644
--- a/NonLinearSolver/internalPoints/ipstate.h
+++ b/NonLinearSolver/internalPoints/ipstate.h
@@ -118,6 +118,10 @@ class AllIPState{
~AllIPState();
ipstateElementContainer* getIPstate(const long int num);
ipstateElementContainer* operator[](const long int num);
+
+ const ipstateElementContainer* getIPstate(const long int num) const;
+ const ipstateElementContainer* operator[] (const long int num) const;
+
const std::vector<partDomain*>& getDomainVector() const {return _domainVector;};
void nextStep();
diff --git a/NonLinearSolver/materialLaw/mlaw.h b/NonLinearSolver/materialLaw/mlaw.h
index 037456df58faa1c5c65a24037fe7f22425ad054f..1712d2641d87424164bd7e0c333a1541debf2439 100644
--- a/NonLinearSolver/materialLaw/mlaw.h
+++ b/NonLinearSolver/materialLaw/mlaw.h
@@ -73,8 +73,8 @@ class materialLaw{
virtual const nonLinearMechSolver* getMacroSolver() const {return _macroSolver;};
virtual bool withEnergyDissipation() const = 0;
- // function to define when one IP is broken
- virtual bool brokenCheck(IPVariable* ipv) const { return false;}
+ // function to define when one IP is broken for block dissipation
+ virtual bool brokenCheck(const IPVariable* ipv) const { return false;}
virtual materialLaw* clone() const = 0;
#endif
};
diff --git a/dG3D/src/dG3DCohesiveMaterialLaw.cpp b/dG3D/src/dG3DCohesiveMaterialLaw.cpp
index 01f7c825c4ec3d12ea905ab9a4e8bd17ef774d2c..0fa3b172adaca3312f19c2c7ade43bfa9a295668 100644
--- a/dG3D/src/dG3DCohesiveMaterialLaw.cpp
+++ b/dG3D/src/dG3DCohesiveMaterialLaw.cpp
@@ -987,21 +987,6 @@ void BulkFollwedCohesive3DLaw::createIPVariable(IPVariable* &ipv,const MElement
ipv = new BulkFollowedCohesive3DIPVariable();
}
-bool BulkFollwedCohesive3DLaw::brokenCheck(IPVariable* ipv) const{
- dG3DIPVariableBase* dgip = dynamic_cast<dG3DIPVariableBase*>(ipv);
- if (dgip != NULL){
- const STensor43& L = dgip->getConstRefToTangentModuli();
- static SVector3 dir;
- dir*=0.;
- double cr;
- FailureDetection::getLostSolutionUniquenessCriterion(getMacroSolver()->getDim() ,L,dir,cr);
- if (cr <= 0) return true;
- else return false;
- }
- else
- return false;
-};
-
void BulkFollwedCohesive3DLaw::checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert) const{
// check failure onset from both negative and postive IPStates
FractureCohesive3DIPVariable* fMinusCur = dynamic_cast<FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::current));
diff --git a/dG3D/src/dG3DCohesiveMaterialLaw.h b/dG3D/src/dG3DCohesiveMaterialLaw.h
index 49aff4da87ae1414f057d19e5fdbac04121e33d5..86205a1cb43761348fd850818e34dc21923423be 100644
--- a/dG3D/src/dG3DCohesiveMaterialLaw.h
+++ b/dG3D/src/dG3DCohesiveMaterialLaw.h
@@ -57,6 +57,7 @@ class Cohesive3DLaw : public materialLaw{
virtual double getCharacteristicLength() const {return 1.;};
virtual materialLaw* clone() const = 0;
virtual bool withEnergyDissipation() const {return true;}
+ virtual bool brokenCheck(const IPVariable* ipv) const { return false;} // alway false
#endif
};
@@ -305,7 +306,6 @@ class GeneralBulkFollwedCohesive3DLaw : public Cohesive3DLaw{
virtual void stress(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff=true, const bool checkfrac=true) = 0;
virtual void checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert= false) const = 0;
virtual void transferInterfaceDataToBulk(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff) const = 0;
- virtual bool brokenCheck(IPVariable* ipv) const =0;
virtual materialLaw* clone() const = 0;
double getLostSolutionUniquenssTolerance() const {return _lostSolutionUniquenssTolerance;};
#endif // SWIG
@@ -326,7 +326,6 @@ class BulkFollwedCohesive3DLaw : public GeneralBulkFollwedCohesive3DLaw{
virtual void stress(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff=true, const bool checkfrac=true);
virtual void checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert= false) const;
virtual void transferInterfaceDataToBulk(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff) const;
- virtual bool brokenCheck(IPVariable* ipv) const;
virtual materialLaw* clone() const{return new BulkFollwedCohesive3DLaw(*this);};
#endif // SWIG
};
diff --git a/dG3D/src/dG3DDomain.cpp b/dG3D/src/dG3DDomain.cpp
index 9ffc082d6a0b92a5064d1947ac2570c23d481540..51b586aa6a55196c13e79fcd24f4cb00b3802ba3 100644
--- a/dG3D/src/dG3DDomain.cpp
+++ b/dG3D/src/dG3DDomain.cpp
@@ -153,10 +153,7 @@ void dG3DDomain::setImposedCrack(const bool flag){
void dG3DDomain::imposeCrack(const int em, const int ep){
printf("add crack between %d and %d \n",em,ep);
-
- int size = _imposedInterfaceCrack.size();
- std::pair<int,int> pa(em,ep);
- _imposedInterfaceCrack.insert(std::pair<int,std::pair<int,int> >(size,pa));
+ _imposedInterfaceCrack.insert(TwoNum(em,ep));
};
@@ -1353,20 +1350,12 @@ int dG3DDomain::getStressDimension() const{
return stresDim;
}
- void dG3DDomain::checkFailure(AllIPState* aips){
+ void dG3DDomain::checkFailure(IPField* ipf) const{
// check insertion of cohesive element
const dG3DMaterialLaw *mlawminus = static_cast<const dG3DMaterialLaw*>(this->getMaterialLawMinus());
const dG3DMaterialLaw *mlawplus = static_cast<const dG3DMaterialLaw*>(this->getMaterialLawPlus());
-
- #if defined(HAVE_MPI)
- // get data in order to send to other procs
- // only new block and broken ips are filled in the group
- if (mlawminus->isNumeric() or mlawplus->isNumeric()){
- this->clearBlockDissipationIPs();
- this->clearBrokenIPs();
- }
- #endif //MPI
-
+
+ // check failure
if (_fullDg and (mlawminus->getType() == materialLaw::fracture or mlawplus->getType() == materialLaw::fracture)){
const FractureByCohesive3DLaw * flawMinus = static_cast<const FractureByCohesive3DLaw*>(mlawminus);
const Cohesive3DLaw *mclaw = static_cast<const Cohesive3DLaw*>(flawMinus->getFractureLaw());
@@ -1379,182 +1368,186 @@ int dG3DDomain::getStressDimension() const{
// if only domain work only on rootranks
for(groupOfElements::elementContainer::iterator it=gi->begin(); it!=gi->end(); ++it){
MInterfaceElement *ie = dynamic_cast<MInterfaceElement*>(*it);
- bool interfaceElementBroken = false; // element broken flag (said element is broken if all Ip are already broken)
-
- std::map<int,std::pair<int,int> >::iterator itcrack = _imposedInterfaceCrack.end();
- bool checkHere = true;
- if (_imposeCrackFlag){
- checkHere =false;
- if (_imposedInterfaceCrack.size() > 0){
- for (itcrack = _imposedInterfaceCrack.begin(); itcrack!= _imposedInterfaceCrack.end(); itcrack++){
- std::pair<int,int>& twoNum = itcrack->second;
- if (((ie->getElem(0)->getNum() == twoNum.first) and (ie->getElem(1)->getNum() == twoNum.second)) or
- ((ie->getElem(1)->getNum() == twoNum.first) and (ie->getElem(0)->getNum() == twoNum.second))){
- checkHere = true;
- printf("checking failure at interface between ele %d and ele %d \n",ie->getElem(0)->getNum(), ie->getElem(1)->getNum());
- break;
- }
- }
- }
- }
-
- if (checkHere){
- MElement *ele = dynamic_cast<MElement*>(ie);
- IntPt* GP;
- int npts=integBound->getIntPoints(ele,&GP);
- int nbIpVBroken = 0; // Number of IPv broken at the interface element
- bool OneIpvBroken = false; // Crack insertion flag (true if one Ipv broken at this time step)
-
- /* Check for cohesive insertion on ip's of the interface element */
- AllIPState::ipstateElementContainer *vips = aips->getIPstate(ie->getNum());
- for(int j=0;j<npts;j++){
- IPStateBase* ipsm = (*vips)[j];
- IPStateBase* ipsp = (*vips)[j+npts];
-
- const FractureCohesive3DIPVariable *ipvfprevm = static_cast<const FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::previous));
- const FractureCohesive3DIPVariable *ipvfprevp = static_cast<const FractureCohesive3DIPVariable*>(ipsp->getState(IPStateBase::previous));
- FractureCohesive3DIPVariable *ipvfm = static_cast<FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::current));
- FractureCohesive3DIPVariable *ipvfp = static_cast<FractureCohesive3DIPVariable*>(ipsp->getState(IPStateBase::current));
-
- // check for current state only
- if (!ipvfm->isbroken()){
- mclaw->checkCohesiveInsertion(ipsm,ipsp);
- // Flag OneIpvBroken set to true if we have crack insertion at this time step at least at one gauss point
- if(ipvfm->isbroken()){
- OneIpvBroken = true;
- #if defined(HAVE_MPI)
- if (mlawminus->isNumeric() or mlawplus->isNumeric()){
- this->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j));
- this->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j+npts));
- }
- #endif //HAVE_MPI
- };
- }
- // Count the number of currently broken ipv
- if(ipvfm->isbroken()){nbIpVBroken++;};
- }
+ int numElMinus = ie->getElem(0)->getNum();
+ int numElPlus = ie->getElem(1)->getNum();
+ if (!(ipf->interfaceElementIsBroken(numElMinus,numElPlus))){
+ // if element is not yet broken
+ bool willCheckFalue = true;
+ if (_imposeCrackFlag){
+ willCheckFalue =false;
+ if (_imposedInterfaceCrack.size() > 0){
+ TwoNum key(numElMinus, numElPlus);
+ if (_imposedInterfaceCrack.find(key) != _imposedInterfaceCrack.end()){
+ willCheckFalue = true;
+ printf("checking failure at interface between ele %d and ele %d \n",numElMinus, numElPlus);
+ break;
+ }
+ }
+ }
+
+ if (willCheckFalue){
+ bool interfaceElementBroken = false; // element broken flag (said element is broken if all Ip are already broken)
+
+ MElement *ele = dynamic_cast<MElement*>(ie);
+
+ IntPt* GP;
+ int npts=integBound->getIntPoints(ele,&GP);
+ int nbIpVBroken = 0; // Number of IPv broken at the interface element
+
+ /* Check for cohesive insertion on ip's of the interface element */
+ AllIPState::ipstateElementContainer *vips = ipf->getAips()->getIPstate(ie->getNum());
+ for(int j=0;j<npts;j++){
+ IPStateBase* ipsm = (*vips)[j];
+ IPStateBase* ipsp = (*vips)[j+npts];
+
+ const FractureCohesive3DIPVariable *ipvfprevm = static_cast<const FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::previous));
+ const FractureCohesive3DIPVariable *ipvfprevp = static_cast<const FractureCohesive3DIPVariable*>(ipsp->getState(IPStateBase::previous));
+ FractureCohesive3DIPVariable *ipvfm = static_cast<FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::current));
+ FractureCohesive3DIPVariable *ipvfp = static_cast<FractureCohesive3DIPVariable*>(ipsp->getState(IPStateBase::current));
+
+ // check for current state only
+ if (!ipvfm->isbroken()){
+ // failure is check for both sides
+ mclaw->checkCohesiveInsertion(ipsm,ipsp);
+ if(ipvfm->isbroken()){
+ #if defined(HAVE_MPI)
+ if (mlawminus->isNumeric() or mlawplus->isNumeric()){
+ ipf->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j));
+ ipf->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j+npts));
+ }
+ #endif //HAVE_MPI
+ };
+ }
+ // Count the number of currently broken ipv
+ if(ipvfm->isbroken()){nbIpVBroken++;};
+ }
- /* Force cohesive insertion and block bulk damage evolution if needed */
- // only carried out during cohesive insertion process if we have crack insertion at this time step at least at one gauss point
- if (OneIpvBroken)
- {
- if (nbIpVBroken ==npts) // all GP are broken
- {
- interfaceElementBroken = true; // interfaceElement is naturally broken
- }
+
+ if (nbIpVBroken ==npts) // all GP are broken
+ {
+ interfaceElementBroken = true; // interfaceElement is naturally broken
+ }
- // Check if forcing is alowed and not already done
- if (_forceInterfaceElementBroken and !interfaceElementBroken)
- {
- // Determine the needed number of broken ipv before forcing
- int nptsToBreakInterfaceElement(0);
- if (_failureIPRatio <= 0.)
- {
- // break interface element when cohesive IP is inserted in first IP
- nptsToBreakInterfaceElement = 1;
- }
- else if (_failureIPRatio >= 1.)
- {
- // break interface element when cohesive IPs are inserted in all IP
- nptsToBreakInterfaceElement = npts;
- }
- else
- {
- // break interface element when cohesive IPs are inserted in agiven number of IPs
- nptsToBreakInterfaceElement = npts*_failureIPRatio;
- if (nptsToBreakInterfaceElement < 1) nptsToBreakInterfaceElement = 1;
- }
+ // Check if forcing is alowed and not already done
+ if (_forceInterfaceElementBroken and !interfaceElementBroken)
+ {
+ // Determine the needed number of broken ipv before forcing
+ int nptsToBreakInterfaceElement(0);
+ if (_failureIPRatio <= 0.)
+ {
+ // break interface element when cohesive IP is inserted in first IP
+ nptsToBreakInterfaceElement = 1;
+ }
+ else if (_failureIPRatio >= 1.)
+ {
+ // break interface element when cohesive IPs are inserted in all IP
+ nptsToBreakInterfaceElement = npts;
+ }
+ else
+ {
+ // break interface element when cohesive IPs are inserted in agiven number of IPs
+ nptsToBreakInterfaceElement = npts*_failureIPRatio;
+ if (nptsToBreakInterfaceElement < 1) nptsToBreakInterfaceElement = 1;
+ }
- // Check if a sufficient number of IP's are broken
- if(nbIpVBroken >= nptsToBreakInterfaceElement)
- {
- interfaceElementBroken = true; // element is force to break
- // Force cohesive insertion at all cohesive ipv if not already broken
- for(int j=0;j<npts;j++)
- {
- // Get Ipstate
- IPStateBase* ipsm = (*vips)[j];
- IPStateBase* ipsp = (*vips)[j+npts];
- // Check if failure already detected and force insertion if needed
- FractureCohesive3DIPVariable *ipvfm = static_cast<FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::current));
- if(!ipvfm->isbroken())
- {
- mclaw->checkCohesiveInsertion(ipsm,ipsp,true);
- #if defined(HAVE_MPI)
- if (mlawminus->isNumeric() or mlawplus->isNumeric()){
- this->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j));
- this->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j+npts));
- }
- #endif //HAVE_MPI
- }
- }
- }
- }
+ // Check if a sufficient number of IP's are broken
+ if(nbIpVBroken >= nptsToBreakInterfaceElement)
+ {
+ interfaceElementBroken = true; // element is force to break
+ // Force cohesive insertion at all cohesive ipv if not already broken
+ for(int j=0;j<npts;j++)
+ {
+ // Get Ipstate
+ IPStateBase* ipsm = (*vips)[j];
+ IPStateBase* ipsp = (*vips)[j+npts];
+ // Check if failure already detected and force insertion if needed
+ FractureCohesive3DIPVariable *ipvfm = static_cast<FractureCohesive3DIPVariable*>(ipsm->getState(IPStateBase::current));
+ if(!ipvfm->isbroken())
+ {
+ // force broken all remaining IPVs
+ mclaw->checkCohesiveInsertion(ipsm,ipsp,true);
+ #if defined(HAVE_MPI)
+ if (mlawminus->isNumeric() or mlawplus->isNumeric()){
+ ipf->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j));
+ ipf->brokenIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),j+npts));
+ }
+ #endif //HAVE_MPI
+ }
+ }
+ }
+ }
+
+ if (interfaceElementBroken){
+ // add to the list of failure elements
+ ipf->brokenInterfaceElement(numElMinus,numElPlus);
+ }
- // Control damage blocking inside bulk elements in terms of chosen method
- bool blockBulkDamage = false; // true if damage is blocked in bulk elements
- if (_bulkDamageBlockMethod == 0)
- {
- // block damage when cohesive crack is inserted in first IP
- if (nbIpVBroken > 0) blockBulkDamage = true;
- }
- else if (_bulkDamageBlockMethod == 1)
- {
- // block damage when cohesive crack is interface element is broken
- if (interfaceElementBroken) blockBulkDamage = true;
- }
- else{
- // no block damage in bulk ips
+ // Control damage blocking inside bulk elements in terms of chosen method
+ bool blockBulkDamage = false; // true if damage is blocked in bulk elements
+ if (_bulkDamageBlockMethod == 0)
+ {
+ // block damage when cohesive crack is inserted in first IP
+ if (nbIpVBroken > 0) blockBulkDamage = true;
+ }
+ else if (_bulkDamageBlockMethod == 1)
+ {
+ // block damage when cohesive crack is interface element is broken
+ if (interfaceElementBroken) blockBulkDamage = true;
+ }
+ else{
+ // no block damage in bulk ips
blockBulkDamage = false;
- }
+ }
- // Block damage in bulk in Plus and Minus bulk elements when using
- // bulk material with damage if needed
- if (blockBulkDamage){
- if (bulkLawMinus->withEnergyDissipation()){
- // Get bulk element and check if it exists
- const partDomain* dom = this->getMinusDomain();
- MElement* em = ie->getElem(0);
- if(dom->g_find(em)){
- // Loop on bulk ipv of the element
- int npts_bulk = dom->getBulkGaussIntegrationRule()->getIntPoints(em,&GP);
- AllIPState::ipstateElementContainer *vips_bulk = aips->getIPstate(em->getNum());
- for(int k=0;k<npts_bulk;k++)
- {
- // Update critical damage and flag in bulk ipv's
- IPStateBase* ips_bulk = (*vips_bulk)[k];
- dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
- const dG3DIPVariableBase* ipv_bulkPrev = static_cast<const dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::previous));
+ // Block damage in bulk in Plus and Minus bulk elements when using
+ // bulk material with damage if needed
+ if (blockBulkDamage){
+ if (bulkLawMinus->withEnergyDissipation()){
+ // Get bulk element and check if it exists
+ const partDomain* dom = this->getMinusDomain();
+ MElement* em = ie->getElem(0);
+ if(dom->g_find(em) and (!ipf->dissipationIsBlockedBulkElement(em->getNum()))){
+ ipf->blockDissipationBulkElement(em->getNum());
+ // Loop on bulk ipv of the element
+ int npts_bulk = dom->getBulkGaussIntegrationRule()->getIntPoints(em,&GP);
+ AllIPState::ipstateElementContainer *vips_bulk = ipf->getAips()->getIPstate(em->getNum());
+ for(int k=0;k<npts_bulk;k++)
+ {
+ // Update critical damage and flag in bulk ipv's
+ IPStateBase* ips_bulk = (*vips_bulk)[k];
+ dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
+ const dG3DIPVariableBase* ipv_bulkPrev = static_cast<const dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::previous));
if (!ipv_bulk->dissipationIsBlocked()){
ipv_bulk->setBulkCriticalDamage(ipv_bulk->getBulkCriticalDamage());
ipv_bulk->blockDissipation(true);
#if defined(HAVE_MPI)
if (bulkLawMinus->isNumeric()){
- this->blockDissipationIP(numericalMaterialBase::createTypeWithTwoInts(em->getNum(),k));
+ ipf->blockDissipationIP(numericalMaterialBase::createTypeWithTwoInts(em->getNum(),k));
}
#endif //HAVE_MPI
}
- }
- }
- }
+ }
+ }
+ }
- if (bulkLawPlus->withEnergyDissipation())
- {
- const partDomain* dom = this->getPlusDomain();
- MElement* ep = ie->getElem(1);
- if(dom->g_find(ep))
- {
- // Loop on bulk ipv of the element
- int npts_bulk = dom->getBulkGaussIntegrationRule()->getIntPoints(ep,&GP);
- AllIPState::ipstateElementContainer *vips_bulk = aips->getIPstate(ep->getNum());
- for(int k=0;k<npts_bulk;k++)
- {
- // Update critical damage and flag in bulk ipv's
- IPStateBase* ips_bulk = (*vips_bulk)[k];
- dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
- const dG3DIPVariableBase* ipv_bulkPrev = static_cast<const dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::previous));
+ if (bulkLawPlus->withEnergyDissipation())
+ {
+ const partDomain* dom = this->getPlusDomain();
+ MElement* ep = ie->getElem(1);
+ if(dom->g_find(ep) and (!ipf->dissipationIsBlockedBulkElement(ep->getNum())))
+ {
+ ipf->blockDissipationBulkElement(ep->getNum());
+ // Loop on bulk ipv of the element
+ int npts_bulk = dom->getBulkGaussIntegrationRule()->getIntPoints(ep,&GP);
+ AllIPState::ipstateElementContainer *vips_bulk = ipf->getAips()->getIPstate(ep->getNum());
+ for(int k=0;k<npts_bulk;k++)
+ {
+ // Update critical damage and flag in bulk ipv's
+ IPStateBase* ips_bulk = (*vips_bulk)[k];
+ dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
+ const dG3DIPVariableBase* ipv_bulkPrev = static_cast<const dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::previous));
if (!ipv_bulk->dissipationIsBlocked()){
ipv_bulk->setBulkCriticalDamage(ipv_bulk->getBulkCriticalDamage());
@@ -1562,98 +1555,63 @@ int dG3DDomain::getStressDimension() const{
#if defined(HAVE_MPI)
if (bulkLawPlus->isNumeric()){
- this->blockDissipationIP(numericalMaterialBase::createTypeWithTwoInts(ep->getNum(),k));
+ ipf->blockDissipationIP(numericalMaterialBase::createTypeWithTwoInts(ep->getNum(),k));
}
#endif //HAVE_MPI
}
- }
- }
- }
- }
- }
- }
-
- if (_imposeCrackFlag){
- if (interfaceElementBroken and itcrack != _imposedInterfaceCrack.end()){
- printf("element %d is fully broken\n", ie->getNum());
- _imposedInterfaceCrack.erase(itcrack);
+ }
+ }
+ }
+ }
}
+
}
- }
+ #ifdef _DEBUG
+ else{
+ printf("element %d with negative %d, positive %d is fully broken\n", ie->getNum(),numElMinus,numElPlus);
+ }
+ #endif //_DEBUG
+ }
}
// block damage in elements which are not related to crack mouth
if (this->groupOfElementsSize() > 0 and (_bulkDamageBlockMethod >=0)){
- materialLaw* mlaw = this->getMaterialLaw();
-
- static std::vector<materialLaw*> allFailureLaw;
- allFailureLaw.clear();
- if (mlaw->getType() == materialLaw::fracture){
- allFailureLaw.push_back(mlaw);
+ const materialLaw* mlaw = this->getMaterialLaw();
+ IntPt* GP;
+ for (groupOfElements::elementContainer::iterator ite = this->g_cbegin(); ite!= this->g_cend(); ite++){
+ MElement* ele = *ite;
+
+ if (!ipf->dissipationIsBlockedBulkElement(ele->getNum())){
+ int npts_bulk = this->getBulkGaussIntegrationRule()->getIntPoints(ele,&GP);
+ AllIPState::ipstateElementContainer *vips_bulk = ipf->getAips()->getIPstate(ele->getNum());
+
+ // if one IP is broken, block all element
+ bool willBlock = false;
+ for(int k=0;k<npts_bulk;k++){
+ IPStateBase* ips_bulk = (*vips_bulk)[k];
+ dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
+ if (mlaw->brokenCheck(ipv_bulk)){
+ willBlock = true;
+ break;
+ }
+ }
+ if (willBlock){
+ ipf->blockDissipationBulkElement(ele->getNum());
+ // block all IP
+ for(int k=0;k<npts_bulk;k++){
+ IPStateBase* ips_bulk = (*vips_bulk)[k];
+ dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
+ if (!ipv_bulk->dissipationIsBlocked()){
+ ipv_bulk->blockDissipation(true);
+ #if defined(HAVE_MPI)
+ ipf->blockDissipationIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),k));
+ #endif
+ }
+ }
+ }
+ }
+
}
-
- const std::vector<partDomain*>* domainVector = &(aips->getDomainVector());
- for (int jdom =0; jdom < domainVector->size(); jdom++){
- dgPartDomain* dgdom = dynamic_cast<dgPartDomain*>((*domainVector)[jdom]);
- if (dgdom != NULL){
- // if this dgdom is purely interfacial
- if (dgdom->getMinusDomain()->getPhysical() != dgdom->getPlusDomain()->getPhysical()){
- if ((dgdom->getMinusDomain()->getPhysical() == this->getPhysical()) and
- (dgdom->getMaterialLawMinus()->getType() == materialLaw::fracture)){
- allFailureLaw.push_back(dgdom->getMaterialLawMinus());
- }
- else if ((dgdom->getPlusDomain()->getPhysical() == this->getPhysical()) and
- (dgdom->getMaterialLawPlus()->getType() == materialLaw::fracture)){
- allFailureLaw.push_back(dgdom->getMaterialLawPlus());
- }
- }
- }
- }
-
- if (allFailureLaw.size() > 0){
- IntPt* GP;
- for (groupOfElements::elementContainer::iterator ite = this->g_cbegin(); ite!= this->g_cend(); ite++){
- MElement* ele = *ite;
- int npts_bulk = this->getBulkGaussIntegrationRule()->getIntPoints(ele,&GP);
- AllIPState::ipstateElementContainer *vips_bulk = aips->getIPstate(ele->getNum());
- int numFailureIP = 0;
- int numBlockedIP = 0;
- for(int k=0;k<npts_bulk;k++){
- IPStateBase* ips_bulk = (*vips_bulk)[k];
- dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
- if (ipv_bulk->dissipationIsBlocked()){
- numBlockedIP++;
- numFailureIP++;
- }
- else{
- bool broken = false;
- for (int ilaw = 0; ilaw < allFailureLaw.size(); ilaw++){
- if (allFailureLaw[ilaw]->brokenCheck(ipv_bulk)){
- broken = true;
- break;
- };
- }
- if (broken){
- numFailureIP++;
- }
- }
- }
- //printf("rank %d numBlockedIP %d numFailureIP %d \n",Msg::GetCommRank(),numBlockedIP,numFailureIP);
- if ((numFailureIP <= npts_bulk) and (numFailureIP>0) and (numBlockedIP < npts_bulk )){
- // block all IP
- for(int k=0;k<npts_bulk;k++){
- IPStateBase* ips_bulk = (*vips_bulk)[k];
- dG3DIPVariableBase *ipv_bulk = static_cast<dG3DIPVariableBase*>(ips_bulk->getState(IPStateBase::current));
- if (!ipv_bulk->dissipationIsBlocked()){
- ipv_bulk->blockDissipation(true);
- #if defined(HAVE_MPI)
- this->blockDissipationIP(numericalMaterialBase::createTypeWithTwoInts(ele->getNum(),k));
- #endif
- }
- }
- }
- }
- }
}
};
diff --git a/dG3D/src/dG3DDomain.h b/dG3D/src/dG3DDomain.h
index 1a48e563b7419d4554291d2d840d7638260358fc..0d2e1bfe9d09670bfc7343afa0ff776065efe85c 100644
--- a/dG3D/src/dG3DDomain.h
+++ b/dG3D/src/dG3DDomain.h
@@ -17,6 +17,7 @@
#include "dG3DIPVariable.h"
#include "partDomain.h"
#include "interfaceQuadrature.h"
+#include "ipField.h"
#endif
class dG3DDomain : public dgPartDomain{
public :
@@ -67,7 +68,7 @@ class dG3DDomain : public dgPartDomain{
// crack injection
bool _imposeCrackFlag;
- std::map<int, std::pair<int,int> > _imposedInterfaceCrack; // pair is positive and negative
+ std::set<TwoNum> _imposedInterfaceCrack; // pair is positive and negative
bool _accountPathFollowing;
@@ -97,7 +98,7 @@ class dG3DDomain : public dgPartDomain{
virtual int getStressDimension() const;
- virtual void checkFailure(AllIPState* aips);
+ virtual void checkFailure(IPField* aips) const;
virtual void initialIPVariable(AllIPState *aips,const unknownField *ufield,const IPStateBase::whichState ws, bool stiff);
diff --git a/dG3D/src/dG3DHomogenizedTangentTerms.cpp b/dG3D/src/dG3DHomogenizedTangentTerms.cpp
index af849fcbe46a9c06464d94049ba113354638af41..f20bcadef170e195c260d448548d7b680357c279 100644
--- a/dG3D/src/dG3DHomogenizedTangentTerms.cpp
+++ b/dG3D/src/dG3DHomogenizedTangentTerms.cpp
@@ -21,7 +21,7 @@ void dG3DHomogenizedTangent::get(MElement *ele,int npts,IntPt *GP,fullMatrix<dou
m.resize(stressDim,nbDof);
m.setAll(0.);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++){
const IPStateBase *ips = (*vips)[i];
@@ -104,7 +104,7 @@ void dG3DNonLocalHomogenizedTangent::get(MElement *ele,int npts,IntPt *GP,fullMa
int nbFF = ele->getNumShapeFunctions();
int stressDim = _dom->getStressDimension();
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
// get derivatives of average strain on active damaging zone
// over applied strain
@@ -200,7 +200,7 @@ void dG3DExtraDofHomogenizedTangent::get(MElement *ele,int npts,IntPt *GP,fullMa
int numMDof =nbFF*3;
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++){
const IPStateBase *ips = (*vips)[i];
diff --git a/dG3D/src/dG3DMultiscaleMaterialLaw.cpp b/dG3D/src/dG3DMultiscaleMaterialLaw.cpp
index 82a0b0f2b6fae25a2a4a45b91c161d8e4adc403c..084135079700f9317834bdb52b355eb7fe57ebd8 100644
--- a/dG3D/src/dG3DMultiscaleMaterialLaw.cpp
+++ b/dG3D/src/dG3DMultiscaleMaterialLaw.cpp
@@ -122,8 +122,8 @@ void dG3DMultiscaleMaterialLaw::initLaws(const std::map<int,materialLaw*> &mapla
}
};
-bool dG3DMultiscaleMaterialLaw::brokenCheck(IPVariable* ipv) const{
- dG3DMultiscaleIPVariable* mipv = dynamic_cast<dG3DMultiscaleIPVariable*>(ipv);
+bool dG3DMultiscaleMaterialLaw::brokenCheck(const IPVariable* ipv) const{
+ const dG3DMultiscaleIPVariable* mipv = dynamic_cast<const dG3DMultiscaleIPVariable*>(ipv);
if (mipv != NULL){
return mipv->solverIsBroken();
}
@@ -520,15 +520,6 @@ dG3DMultiscaleCohesiveLaw::dG3DMultiscaleCohesiveLaw(const int num, const double
dG3DMultiscaleCohesiveLaw::dG3DMultiscaleCohesiveLaw(const dG3DMultiscaleCohesiveLaw& src):
GeneralMultiscaleBulkFollwedCohesive3DLaw(src),_withNormalPart(src._withNormalPart),_fixedCrackFace(src._fixedCrackFace),_fixedCrackFaceFlag(src._fixedCrackFaceFlag){}
-bool dG3DMultiscaleCohesiveLaw::brokenCheck(IPVariable* ipv) const{
- dG3DMultiscaleIPVariable* mipv = dynamic_cast<dG3DMultiscaleIPVariable*>(ipv);
- if (mipv != NULL){
- return mipv->solverIsBroken();
- }
- else{
- return false;
- }
-};
void dG3DMultiscaleCohesiveLaw::checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert) const{
// check failure onset from both negative and postive IPStates
@@ -1138,15 +1129,6 @@ void TwoFieldMultiscaleCohesive3DLaw::transferInterfaceDataToBulk(IPVariable* ip
}
};
-bool TwoFieldMultiscaleCohesive3DLaw::brokenCheck(IPVariable* ipv) const{
- dG3DMultiscaleIPVariable* mipv = dynamic_cast<dG3DMultiscaleIPVariable*>(ipv);
- if (mipv != NULL){
- return mipv->solverIsBroken();
- }
- else{
- return false;
- }
-};
void TwoFieldMultiscaleCohesive3DLaw::stress(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff, const bool checkfrac){
FractureCohesive3DIPVariable *fipv = static_cast < FractureCohesive3DIPVariable *> (ipv);
diff --git a/dG3D/src/dG3DMultiscaleMaterialLaw.h b/dG3D/src/dG3DMultiscaleMaterialLaw.h
index 3f46732b06aa43146950f3900da9d45321cf2aa1..6caa74f9970577458c7b75e094aec67b64e157c2 100644
--- a/dG3D/src/dG3DMultiscaleMaterialLaw.h
+++ b/dG3D/src/dG3DMultiscaleMaterialLaw.h
@@ -110,7 +110,7 @@ class dG3DMultiscaleMaterialLaw : public dG3DMultiscaleMaterialLawBase{
virtual bool isNumeric() const{return true;}
virtual void initialIPVariable(IPVariable* ipv, bool stiff);
- virtual bool brokenCheck(IPVariable* ipv) const;
+ virtual bool brokenCheck(const IPVariable* ipv) const;
virtual materialLaw* clone() const {return new dG3DMultiscaleMaterialLaw(*this);};
#endif // SWIG
@@ -191,7 +191,6 @@ class GeneralMultiscaleBulkFollwedCohesive3DLaw : public GeneralBulkFollwedCohes
GeneralMultiscaleBulkFollwedCohesive3DLaw(const int num, const bool init = true);
GeneralMultiscaleBulkFollwedCohesive3DLaw(const GeneralMultiscaleBulkFollwedCohesive3DLaw& src);
virtual ~GeneralMultiscaleBulkFollwedCohesive3DLaw();
- virtual bool brokenCheck(IPVariable* ipv) const = 0;
virtual void initLaws(const std::map<int,materialLaw*> &maplaw){} // nothing to do for this law
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 = 0;
virtual void createIPVariable(IPVariable* &ipv,const MElement *ele, const int nbFF_, const IntPt *GP, const int gpt) const = 0;
@@ -232,7 +231,6 @@ class dG3DMultiscaleCohesiveLaw : public GeneralMultiscaleBulkFollwedCohesive3DL
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(IPVariable* &ipv,const MElement *ele, const int nbFF_, const IntPt *GP, const int gpt) const;
- virtual bool brokenCheck(IPVariable* ipv) const;
virtual void stress(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff=true, const bool checkfrac=true);
virtual void checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert= false) const;
virtual void transferInterfaceDataToBulk(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff) const;
@@ -269,7 +267,6 @@ class TwoFieldMultiscaleCohesive3DLaw : public GeneralMultiscaleBulkFollwedCohes
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(IPVariable* &ipv,const MElement *ele, const int nbFF_, const IntPt *GP, const int gpt) const;
- virtual bool brokenCheck(IPVariable* ipv) const;
virtual void stress(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff=true, const bool checkfrac=true);
virtual void checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert= false) const;
virtual void transferInterfaceDataToBulk(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff) const;
@@ -290,7 +287,6 @@ class TwoFieldMultiscaleAdhesive3DLaw : public GeneralMultiscaleBulkFollwedCohes
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(IPVariable* &ipv,const MElement *ele, const int nbFF_, const IntPt *GP, const int gpt) const;
- virtual bool brokenCheck(IPVariable* ipv) const {return false;};
virtual void stress(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff=true, const bool checkfrac=true);
virtual void checkCohesiveInsertion(IPStateBase* ipsm, IPStateBase* ipsp, const bool forcedInsert= false) const;
virtual void transferInterfaceDataToBulk(IPVariable* ipv, const IPVariable* ipvprev, const bool stiff) const;
diff --git a/dG3D/src/dG3DTerms.cpp b/dG3D/src/dG3DTerms.cpp
index cf567a68b7d90e8b41ac7b4851a05b39a18f9f7e..05233e7312f22184b59a615b3a19c7137cee4575 100644
--- a/dG3D/src/dG3DTerms.cpp
+++ b/dG3D/src/dG3DTerms.cpp
@@ -195,7 +195,7 @@ void dG3DForceInter::get(MElement *ele, int npts, IntPt *GP, fullVector<double>
const double b1hs = _beta1/hs;
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -375,7 +375,7 @@ void dG3DStiffnessInter::get(MElement *ele,int npts,IntPt *GP,fullMatrix<double>
const double b1hs = _beta1/hs;
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -981,7 +981,7 @@ void nonLocalDG3DForceInter::get(MElement *ele, int npts, IntPt *GP, fullVector<
if( _continuity)
{
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// value on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -1149,7 +1149,7 @@ void nonLocalDG3DStiffnessInter::get(MElement *ele,int npts,IntPt *GP,fullMatrix
if( _continuity)
{
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// values of Gauss points on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -1588,7 +1588,7 @@ void ExtraDofDiffusionDG3DForceInter::get(MElement *ele, int npts, IntPt *GP, fu
if( _continuity)
{
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// value on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -1798,7 +1798,7 @@ void ExtraDofDiffusionDG3DStiffnessInter::get(MElement *ele,int npts,IntPt *GP,f
if( _continuity)
{
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// values of Gauss points on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -2705,7 +2705,7 @@ void ElecTherMechDG3DForceInter::get(MElement *ele, int npts, IntPt *GP, fullVec
if( _continuity)
{
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// value on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
@@ -3225,7 +3225,7 @@ void ElecTherMechDG3DStiffnessInter::get(MElement *ele,int npts,IntPt *GP,fullMa
if( _continuity)
{
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// values of Gauss points on minus and plus elements
IntPt *GPm; IntPt *GPp;
_interQuad->getIntPoints(ie,GP,&GPm,&GPp);
diff --git a/dG3D/src/multipleFieldDG3DTerms.cpp b/dG3D/src/multipleFieldDG3DTerms.cpp
index 3c93dee8bb394fbed68ea0c27027dc7620bd71e3..f6f82f8b83bd42ac8a2bb91c103583f306fe3329 100644
--- a/dG3D/src/multipleFieldDG3DTerms.cpp
+++ b/dG3D/src/multipleFieldDG3DTerms.cpp
@@ -22,7 +22,7 @@ void twoFieldDG3DForceBulk::get(MElement *ele,int npts,IntPt *GP,fullVector<doub
vFor.resize(nbdof);
vFor.setAll(0.);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++){
const dG3DEnhancedStrainIPVariable *ipv = static_cast<const dG3DEnhancedStrainIPVariable*>((*vips)[i]->getState(IPStateBase::current));
@@ -66,7 +66,7 @@ void twoFieldDG3DStiffnessBulk::get(MElement *ele,int npts,IntPt *GP,fullMatrix<
mStiff.resize(nbdof,nbdof);
mStiff.setAll(0.);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++){
const dG3DEnhancedStrainIPVariable *ipv = static_cast<const dG3DEnhancedStrainIPVariable*>((*vips)[i]->getState(IPStateBase::current));
@@ -156,7 +156,7 @@ void twoFieldDG3DForceInter::get(MElement *ele, int npts, IntPt *GP, fullVector<
if(_dom->getFormulation()){
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_dom->getInterfaceQuadrature()->getIntPoints(ie,GP,&GPm,&GPp);
@@ -313,7 +313,7 @@ void twoFieldDG3DStiffnessInter::get(MElement *ele,int npts,IntPt *GP,fullMatrix
double hs = ie->getCharacteristicSize();
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_dom->getInterfaceQuadrature()->getIntPoints(ie,GP,&GPm,&GPp);
@@ -643,7 +643,7 @@ void dG3DTwoFieldForceInter::get(MElement *ele, int npts, IntPt *GP, fullVector<
if(_dom->getFormulation()){
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_dom->getInterfaceQuadrature()->getIntPoints(ie,GP,&GPm,&GPp);
@@ -854,7 +854,7 @@ void dG3DTwoFieldStiffnessInter::get(MElement *ele,int npts,IntPt *GP, fullMatri
if(_dom->getFormulation()){
// get value at gauss's point
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_dom->getInterfaceQuadrature()->getIntPoints(ie,GP,&GPm,&GPp);
diff --git a/dG3D/src/pathFollowingTerms.cpp b/dG3D/src/pathFollowingTerms.cpp
index 625ecc125814f5a7e91c44626cfac728ee41fce4..46b7e6f2314a31e4ff6f7a7f3f1ea9172267ab78 100644
--- a/dG3D/src/pathFollowingTerms.cpp
+++ b/dG3D/src/pathFollowingTerms.cpp
@@ -15,7 +15,7 @@
void dG3DDissipationPathFollowingBulkScalarTerm::get(MElement *ele, int npts, IntPt *GP, double &val) const{
val = 0.;
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++){
const dG3DIPVariableBase* ipv = dynamic_cast<const dG3DIPVariableBase*>((*vips)[i]->getState(IPStateBase::current));
const dG3DIPVariableBase* ipvprev = dynamic_cast<const dG3DIPVariableBase*>((*vips)[i]->getState(IPStateBase::previous));
@@ -35,7 +35,7 @@ void dG3DDissipationPathFollowingBulkLinearTerm::get(MElement *ele,int npts,IntP
m.resize(nbdof);
m.setAll(0.);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++)
{
@@ -70,7 +70,7 @@ void dG3DEnhancedStrainDissipationPathFollowingBulkLinearTerm::get(MElement *ele
FunctionSpaceBase* sp = _dom->getFunctionSpace();
int nbFF = ele->getNumShapeFunctions();
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
for (int i = 0; i < npts; i++)
{
@@ -104,7 +104,7 @@ void dG3DDissipationPathFollowingBoundScalarTerm::get(MElement *ele, int npts, I
}
if ((ie != NULL) and _dom->getFormulation()){
// get all ipvs
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
for (int i = 0; i < npts; i++){
const IPStateBase *ipsm = (*vips)[i];
@@ -152,7 +152,7 @@ void dG3DDissipationPathFollowingBoundLinearTerm::get(MElement *ele,int npts,Int
m.resize(nbdofm+nbdofp);
m.scale(0.);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_dom->getInterfaceQuadrature()->getIntPoints(ie,GP,&GPm,&GPp);
@@ -274,7 +274,7 @@ void twoFieldDG3DDissipationPathFollowingBoundLinearTerm::get(MElement *ele,int
m.scale(0.);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ie->getNum());
// Values on minus and plus elements
IntPt *GPm; IntPt *GPp;
_dom->getInterfaceQuadrature()->getIntPoints(ie,GP,&GPm,&GPp);
@@ -354,7 +354,7 @@ void dG3DNonLocalDissipationPathFollowingBulkScalarTerm::get(MElement *ele, int
void dG3DNonLocalDissipationPathFollowingBulkLinearTerm::get(MElement *ele,int npts,IntPt *GP,fullVector<double> &m) const{
dG3DDissipationPathFollowingBulkLinearTerm::get(ele,npts,GP,m);
- AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer *vips = _ipf->getAips()->getIPstate(ele->getNum());
FunctionSpaceBase* space = _dom->getFunctionSpace();
diff --git a/dgshell/src/dgC0ShellTerms.cpp b/dgshell/src/dgC0ShellTerms.cpp
index d60919a9cba190ca7da9a10b522e065e04627be7..2af8facbf9ab41c240aa04da58490db65379e877 100644
--- a/dgshell/src/dgC0ShellTerms.cpp
+++ b/dgshell/src/dgC0ShellTerms.cpp
@@ -259,7 +259,7 @@ void IsotropicForceBulkTermC0DgShell::get(MElement *ele,int npts,IntPt *GP,fullV
m2.resize(6*npts,nbdof,false); // set operation --> no need to setAll 0 here
// get gauss' point data
- AllIPState::ipstateElementContainer* vips = ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer* vips = ipf->getAips()->getIPstate(ele->getNum());
// set matrix m1 and m2 (one set by component x,y,z for m1)
// get the shape function (and derivative values)
@@ -273,7 +273,7 @@ void IsotropicForceBulkTermC0DgShell::get(MElement *ele,int npts,IntPt *GP,fullV
// Coordonate of Gauss' point i
// u = GP[i].pt[0]; v = GP[i].pt[1]; w = GP[i].pt[2];
// ipv data
- IPStateBase* ips = (*vips)[i];
+ const IPStateBase* ips = (*vips)[i];
const IPVariableShell *ipv = static_cast<const IPVariableShell*>(ips->getState(IPStateBase::current));
const shellLocalBasis *lb = ipv->getshellLocalBasis();
const Bvector *B = ipv->getBvector();
diff --git a/dgshell/src/dgNonLinearShellTerms.cpp b/dgshell/src/dgNonLinearShellTerms.cpp
index 2ccd6dc93da81fd266c44782a1ffe71bcc29c053..3c8b1b2314de529efc16ebf2fd132e135234714d 100644
--- a/dgshell/src/dgNonLinearShellTerms.cpp
+++ b/dgshell/src/dgNonLinearShellTerms.cpp
@@ -19,7 +19,7 @@ void dgNonLinearShellForceBulk::get(MElement *ele,int npts,IntPt *GP,fullVector<
m.resize(nbdof);
m.scale(0.);
// get IPVariable of element
- AllIPState::ipstateElementContainer* vips = _ipf->getAips()->getIPstate(ele->getNum());
+ const AllIPState::ipstateElementContainer* vips = _ipf->getAips()->getIPstate(ele->getNum());
// Shape functions values at Gauss points
nlsFunctionSpaceUVW<double>* sp1 = static_cast<nlsFunctionSpaceUVW<double>*>(&(this->space1));
std::vector<GaussPointSpaceValues<double>*> vgps;
diff --git a/dgshell/src/dgShellDomain.h b/dgshell/src/dgShellDomain.h
index 170a89faf739605b040f9d5f1631d04829679f79..629951225279161516e88fd0c83cdf2eee3adb61 100644
--- a/dgshell/src/dgShellDomain.h
+++ b/dgshell/src/dgShellDomain.h
@@ -55,7 +55,7 @@ class dgLinearShellDomain : public dgPartDomain{
fullVector<double> &dispp, fullVector<double> &dispExtra,
const bool virt, bool stiff, const bool checkfrac=true);
virtual void computeIPVariable(AllIPState *aips,const unknownField *ufield,const IPStateBase::whichState ws, bool stiff);
- virtual void checkFailure(AllIPState* aips){
+ virtual void checkFailure(IPField* ipf) const{
// Msg::Warning("This function checkFailure is not defined ");
}
virtual FunctionSpaceBase* getSpaceForBC(const nonLinearBoundaryCondition::type bc_type, const nonLinearBoundaryCondition::location bc_location,