diff --git a/dG3D/src/nonLocalDamageDG3DIPVariable.cpp b/dG3D/src/nonLocalDamageDG3DIPVariable.cpp index e8a6b4c02baf2567b3e06e03db93efe8a3f40680..8caf6ff9bc1fbfa552e82a446ba59234d0d2c520 100644 --- a/dG3D/src/nonLocalDamageDG3DIPVariable.cpp +++ b/dG3D/src/nonLocalDamageDG3DIPVariable.cpp @@ -263,7 +263,7 @@ nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable::nonLocalDamageTorchANNDG3DIP _D36comp(source._D36comp), _DTensor43(source._DTensor43), _nlDTensor43(source._nlDTensor43), _internalVars_stress(source._internalVars_stress), _internalVars_elTang(source._internalVars_elTang), _internalVars_NonDamTang(source._internalVars_NonDamTang), _kinematicVariables(source._kinematicVariables),_irreversibleEnergy(source._irreversibleEnergy), _DirrEnergyDF(source._DirrEnergyDF),_DirrEnergyDVar_nl0(source._DirrEnergyDVar_nl0), _DirrEnergyDVar_nl1(source._DirrEnergyDVar_nl1), _plasticEnergy(source._plasticEnergy), _dExx(source._dExx), _dEyy(source._dEyy), _dExy(source._dExy), - _defoEnergy(source._defoEnergy), _damageEnergy(source._damageEnergy), _ipMeca(NULL), _DirrEnergyDVar_nl(source._DirrEnergyDVar_nl) + _defoEnergy(source._defoEnergy), _damageEnergy(source._damageEnergy), _ipMeca(NULL), _DirrEnergyDVar_nl(source._DirrEnergyDVar_nl), _S(source._S), _CelD(source._CelD), _Cel(source._Cel) { for( std::vector< IPCLength * > ::const_iterator it=source.ipvCL.begin(); it!=source.ipvCL.end(); it++) @@ -330,6 +330,8 @@ nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable& nonLocalDamageTorchANNDG3DIP _Eplast = psrc->_Eplast; _Eplastrec = psrc->_Eplastrec; //plastic strain reconstructed _Cel = psrc->_Cel; //elastic stiffness without damage + _CelD = psrc->_CelD; //elastic stiffness with damage + _S = psrc->_S; //elastic stiffness without damage _DS_nlDE = psrc->_DS_nlDE; //DStressreconstructedDStrain _DEplastDE = psrc->_DEplastDE; //DEplastDE _dExx = psrc->_dExx; //strain increment following xx (added on 21/05/25) @@ -428,6 +430,8 @@ void nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable::restart() //plastic strain reconstructed restartManager::restart(_Eplastrec); //elastic stiffness without damage + restartManager::restart(_S); + restartManager::restart(_CelD); restartManager::restart(_Cel); restartManager::restart(_DTensor43); restartManager::restart(_nlDTensor43); diff --git a/dG3D/src/nonLocalDamageDG3DIPVariable.h b/dG3D/src/nonLocalDamageDG3DIPVariable.h index 95f784df80b3404c978d13df990e3243ab09a6a0..5f9f405635b422abf2bb845b7fde5875cc68319c 100644 --- a/dG3D/src/nonLocalDamageDG3DIPVariable.h +++ b/dG3D/src/nonLocalDamageDG3DIPVariable.h @@ -1582,6 +1582,12 @@ class nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable : public dG3DIPVariable //plastic strain reconstructed STensor3 _Eplastrec; + + //local stress (added on 26/05/25) + STensor3 _S; + + //elastic stiffness with damage (added on 26/05/25) + STensor43 _CelD; //elastic stiffness without damage STensor43 _Cel; @@ -1694,6 +1700,14 @@ class nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable : public dG3DIPVariable //plastic strain reconstructed virtual const STensor3 &getConstRefToEplastrec() const {return _Eplastrec;}; virtual STensor3 &getRefToEplastrec() {return _Eplastrec;}; + + //local stress (added on 26/05/25) + virtual const STensor3 &getConstRefToLocalStress() const {return _S;}; + virtual STensor3 &getRefToLocalStress() {return _S;}; + + //elastic stiffness with damage (added on 26/05/25) + virtual const STensor43 &getConstRefToCelD() const {return _CelD;}; + virtual STensor43 &getRefToCelD() {return _CelD;}; //elastic stiffness without damage virtual const STensor43 &getConstRefToCel() const {return _Cel;}; diff --git a/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp b/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp index 8f1c9da44fe1d0caf9e3cc72154d5f71f4f7dd6c..47b151347c2d2d72a0c7339d676a979f9ae8152f 100644 --- a/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp +++ b/dG3D/src/nonLocalDamageDG3DMaterialLaw.cpp @@ -2823,12 +2823,12 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::initLaws(const std::map<int,mat E1Vector_plus_pert1(s) = E1_plus(0,s); } STensorOperation::fromFullVectorToSTensor3(E1Vector_plus_pert1, E1Tensor_plus_pert1); - compute_local_values(E1Tensor, E1Tensor_plus_pert1,h_stress_init, h_stress_tmp, h_elTang_init, h_elTang_tmp, h_NonDamTang_init, h_NonDamTang_tmp, S_plus, NonDamTangPred_plus,elTangPred_plus, D_plus, Eplast_plus,false, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); + //compute_local_values(E1Tensor, E1Tensor_plus_pert1,h_stress_init, h_stress_tmp, h_elTang_init, h_elTang_tmp, h_NonDamTang_init, h_NonDamTang_tmp, S_plus, NonDamTangPred_plus,elTangPred_plus, D_plus, Eplast_plus,false, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); (commented on 26/05/25 because wanted to add the ipvp) } } else { - compute_local_values(E1Tensor, E1Tensor, h_stress_init, h_stress_tmp, h_elTang_init, h_elTang_tmp, h_NonDamTang_init, h_NonDamTang_tmp, S,NonDamTangPred,elTangPred,D, Eplast, true, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); + //compute_local_values(E1Tensor, E1Tensor, h_stress_init, h_stress_tmp, h_elTang_init, h_elTang_tmp, h_NonDamTang_init, h_NonDamTang_tmp, S,NonDamTangPred,elTangPred,D, Eplast, true, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); (commented on 26/05/25 because wanted to add the ipvp) } _initialized = true; }; @@ -3114,7 +3114,11 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::stress(IPVariable* ipv, const I if(_tangentByPerturbation){ //beginning step 3.1: we compute the RNN predictions at this time step and the values deduced by the RNN predictions (local damage D and plastic strain Eplast) - compute_local_values(E0Tensor, E1Tensor, h_stress_0, h_stress_1, h_elTang_0, h_elTang_1, h_NonDamTang_0, h_NonDamTang_1, S,NonDamTangPred,elTangPred, D, Eplast, false, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); + compute_local_values(ipvprev,E0Tensor, E1Tensor, h_stress_0, h_stress_1, h_elTang_0, h_elTang_1, h_NonDamTang_0, h_NonDamTang_1, S,NonDamTangPred,elTangPred, D, Eplast, false, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); + + ipvcur->getRefToLocalStress()=S; //we retrieve the values predicted by RNN (added on 26/05/25) + ipvcur->getRefToCelD()=elTangPred; //we retrieve the values predicted by RNN (added on 26/05/25) + ipvcur->getRefToCel()=NonDamTangPred; //we retrieve the values predicted by RNN (added on 26/05/25) //to have D always increasing for (int i=0;i<3;i++){ @@ -3176,7 +3180,7 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::stress(IPVariable* ipv, const I E1Tensor_plus(i,j)=0.5*(F_plus(i,j)+Ft_plus(i,j))-I(i,j); } } - compute_local_values(E0Tensor, E1Tensor_plus, h_stress_0, h_stress_tmp, h_elTang_0, h_elTang_tmp, h_NonDamTang_0, h_NonDamTang_tmp, S_plus,NonDamTangPred_plus,elTangPred_plus, D_plus, Eplast_plus, false, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l_plus); + compute_local_values(ipvprev,E0Tensor, E1Tensor_plus, h_stress_0, h_stress_tmp, h_elTang_0, h_elTang_tmp, h_NonDamTang_0, h_NonDamTang_tmp, S_plus,NonDamTangPred_plus,elTangPred_plus, D_plus, Eplast_plus, false, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l_plus); //beginning step 3.3.1: computation of DEplastDE and DEplast_lDE @@ -3269,7 +3273,7 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::stress(IPVariable* ipv, const I //end step 3.3 } }else{ //be careful: it is the differentiation automatic case that does not work for now - compute_local_values(E0Tensor, E1Tensor, h_stress_0, h_stress_1, h_elTang_0, h_elTang_1, h_NonDamTang_0, h_NonDamTang_1, S,NonDamTangPred,elTangPred, D, Eplast, true, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); + compute_local_values(ipvprev,E0Tensor, E1Tensor, h_stress_0, h_stress_1, h_elTang_0, h_elTang_1, h_NonDamTang_0, h_NonDamTang_1, S,NonDamTangPred,elTangPred, D, Eplast, true, DSDE, B0, B1, B2, B3, B4, B5, Eplast_l); } } //end step 3 @@ -5420,7 +5424,7 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::stress(IPVariable* ipv, const I #if defined(HAVE_TORCH) -void NonlocalDamageTorchANNBasedDG3DMaterialLaw::compute_local_values(const STensor3& E0Tensor, const STensor3& E1Tensor, const torch::Tensor& h_stress_0, torch::Tensor& h_stress_1, const torch::Tensor& h_elTang_0, torch::Tensor& h_elTang_1, const torch::Tensor& h_NonDamTang_0, torch::Tensor& h_NonDamTang_1, STensor3& S, +void NonlocalDamageTorchANNBasedDG3DMaterialLaw::compute_local_values(const nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable* ipvprev,const STensor3& E0Tensor, const STensor3& E1Tensor, const torch::Tensor& h_stress_0, torch::Tensor& h_stress_1, const torch::Tensor& h_elTang_0, torch::Tensor& h_elTang_1, const torch::Tensor& h_NonDamTang_0, torch::Tensor& h_NonDamTang_1, STensor3& S, STensor43& NonDamTangPred,STensor43& elTangPred, STensor43& D, STensor3& Eplast, const bool stiff, fullMatrix<double>& DSDE, const STensor3& B0, const STensor3& B1, const STensor3& B2, const STensor3& B3, const STensor3& B4, const STensor3& B5, STensor3& Eplast_l) { torch::Tensor h_stress_tmp_0 = h_stress_0; @@ -5477,6 +5481,13 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::compute_local_values(const STen h_elTang_1 = h_elTang_tmp_1; h_NonDamTang_1 = h_NonDamTang_tmp_1; + //try to handle the case in which we are out of good RNN predictions, but it is not well defined, and maybe the code can be simplified (added on 26/05/25) + if (E1Tensor(0,0)>0.09){ + S=ipvprev->getConstRefToLocalStress(); + elTangPred=ipvprev->getConstRefToCelD(); + NonDamTangPred=ipvprev->getConstRefToCel(); + } + computation_D(NonDamTangPred, elTangPred, D); //computation of D computation_Eplast(E1Tensor, S, elTangPred, Eplast); //computation of Eplast computation_Eplast_l(B0, B1, B2, B3, B4, B5, Eplast, Eplast_l); //computation of Eplast_l @@ -5643,6 +5654,7 @@ void NonlocalDamageTorchANNBasedDG3DMaterialLaw::RNNstress_stiff(const STensor3& h_NonDamTang_1 = outputs_NonDamTang->elements()[1].toTensor(); InverseNormalize_NonDamTang(NonDamTang_norm, NonDamTangPred); + /*if (_RNN) diff --git a/dG3D/src/nonLocalDamageDG3DMaterialLaw.h b/dG3D/src/nonLocalDamageDG3DMaterialLaw.h index 06a8d460f0c7e4550896e5241b95b19df50237ad..89b80597507a5dd77a4bd7ab338b246ed4e5167f 100644 --- a/dG3D/src/nonLocalDamageDG3DMaterialLaw.h +++ b/dG3D/src/nonLocalDamageDG3DMaterialLaw.h @@ -1183,7 +1183,7 @@ class NonlocalDamageTorchANNBasedDG3DMaterialLaw : public dG3DMaterialLaw{ void computation_Eplast_l(const STensor3& B0, const STensor3& B1, const STensor3& B2, const STensor3& B3, const STensor3& B4, const STensor3& B5, const STensor3& Eplast, STensor3& Eplast_l); void computation_D_l(const STensor43& A0, const STensor43& A1, const STensor43& A2, const STensor43& A3, const STensor43& A4, const STensor43& A5, const STensor43& A6, const STensor43& A7, const STensor43& A8, const STensor43& A9, const STensor43& A10, const STensor43& A11, const STensor43& A12, const STensor43& A13, const STensor43& A14, const STensor43& A15, const STensor43& A16, const STensor43& A17, const STensor43& A18, const STensor43& A19, const STensor43& A20, const STensor43& D, STensor43& D_l); void computation_D_nl__Eplast_nl(const STensor43& A0, const STensor43& A1, const STensor43& A2, const STensor43& A3, const STensor43& A4, const STensor43& A5, const STensor43& A6, const STensor43& A7, const STensor43& A8, const STensor43& A9, const STensor43& A10, const STensor43& A11, const STensor43& A12, const STensor43& A13, const STensor43& A14, const STensor43& A15, const STensor43& A16, const STensor43& A17, const STensor43& A18, const STensor43& A19, const STensor43& A20, const STensor3& B0, const STensor3& B1, const STensor3& B2, const STensor3& B3, const STensor3& B4, const STensor3& B5, const STensor43& D, const STensor3& Eplast, STensor3& Eplast_nl, STensor43& D_nl, nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable* ipvcur); - void compute_local_values(const STensor3& E0Tensor, const STensor3& E1Tensor, const torch::Tensor& h_stress_0, torch::Tensor& h_stress_1, const torch::Tensor& h_elTang_0, torch::Tensor& h_elTang_1, const torch::Tensor& h_NonDamTang_0, torch::Tensor& h_NonDamTang_1, STensor3& S, + void compute_local_values(const nonLocalDamageTorchANNDG3DIPVariableDG3DIPVariable* ipvprev,const STensor3& E0Tensor, const STensor3& E1Tensor, const torch::Tensor& h_stress_0, torch::Tensor& h_stress_1, const torch::Tensor& h_elTang_0, torch::Tensor& h_elTang_1, const torch::Tensor& h_NonDamTang_0, torch::Tensor& h_NonDamTang_1, STensor3& S, STensor43& NonDamTangPred,STensor43& elTangPred, STensor43& D, STensor3& Eplast, const bool stiff, fullMatrix<double>& DSDE, const STensor3& B0, const STensor3& B1, const STensor3& B2, const STensor3& B3, const STensor3& B4, const STensor3& B5, STensor3& Eplast_l); void RNNstress_stiff(const STensor3& E0Tensor, const STensor3& E1Tensor, const torch::Tensor& h_stress_0, torch::Tensor& h_stress_1, const torch::Tensor& h_elTang_0, torch::Tensor& h_elTang_1, const torch::Tensor& h_NonDamTang_0, torch::Tensor& h_NonDamTang_1, STensor3& S, STensor43& NonDamTangPred,STensor43& elTangPred, const bool stiff, fullMatrix<double>& DSDE); //the function in whichn we compute the values predicted by RNN