adapt coupled porous law to new selective update

NonLinearSolver/materialLaw/mlawNonLocalPorousCoalescence.cpp

@@ -130,10 +130,10 @@ double mlawNonLocalPorousThomasonLaw::yieldFunction(const double kcorEq, const d
   double Cft = q1->getConstRefToIPCoalescence().getConcentrationFactor();
   double CftOffset = q0->getConstRefToIPCoalescence().getCrackOffsetOnCft();
   if (_CfTOffsetMethod==2){
-    if (!q0->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
+    if (!q1->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
   }
   else if (_CfTOffsetMethod==3){
-    if (!q0->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
+    if (!q1->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
   }
 
   double R0 = _j2IH->getYield0();
@@ -176,10 +176,10 @@ void mlawNonLocalPorousThomasonLaw::computeYieldDerivatives( fullVector<double>&
   double Cft = q1Thom->getConcentrationFactor();
   double CftOffset = q0Thom->getCrackOffsetOnCft();
   if (_CfTOffsetMethod==2){
-    if (!q0->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
+    if (!q1->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
   }
   else if (_CfTOffsetMethod==3){
-    if (!q0->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
+    if (!q1->getConstRefToIPCoalescence().getCoalescenceOnsetFlag()){CftOffset = 1.;};
   }
   
   double DCftDChi = q1Thom->getDConcentrationFactorDLigamentRatio();

NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.cpp

@@ -373,7 +373,7 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
     if (q1->getLocation() == IPVariable::INTERFACE_MINUS or q1->getLocation() == IPVariable::INTERFACE_PLUS){
       // At the interface:
       // Check if transition has already occurs or not
-      if (!q0Thom->getCoalescenceOnsetFlag()){
+      if (!q1Thom->getCoalescenceOnsetFlag()){
         // If Coalescence has not yet occured at least once
         // Check if the normal traction force is sufficient to obtain coalescence in normal direction of the interface
 
@@ -464,7 +464,8 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
       }
       else{
         // Coalescence has occured at least once
-        q1Thom->getRefToCoalescenceOnsetFlag() = true;
+        q1Thom->getRefToCoalescenceOnsetFlag() = true; // already done....
+        
         // Update onset value (yield porosity and geometrical parameters)
         q1Thom->getRefToPorosityAtCoalescenceOnset() = q0Thom->getPorosityAtCoalescenceOnset();
         q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q0Thom->getLigamentRatioAtCoalescenceOnset();
@@ -476,6 +477,13 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
 
 
         // Determine which mode should be used for the next time step.
+        if(!q0Thom->getCoalescenceOnsetFlag()){
+          // First time that coalescence is detected : Thomason == ON
+          q1Thom->getRefToCoalescenceActiveFlag() = true;
+          q1Thom->getRefToAccelerateRate() = 1.0;
+        }
+        else{
+          // Not the first time that coalescence appears
           if (q0Thom->getCoalescenceActiveFlag()){
             // If Thomason is/was used for this time step : check if Gurson is more restrictive
             double yieldGurson = _mlawGrowth->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
@@ -498,11 +506,10 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
               q1Thom->getRefToCoalescenceActiveFlag() = true; // Thomason will be used at the next time step
               q1Thom->getRefToAccelerateRate() = 1.0;
               
-            if (_CfTOffsetMethod == 3){
-              q1Thom->getRefToCrackOffsetOnCft() = q0Thom->getCrackOffsetOnCft() + yieldThomason*R0/(Cft*R);
-              
-              if (q1Thom->getRefToCrackOffsetOnCft() > 2.){Msg::Error("High predicted value of Cft Offset = %e",q1Thom->getRefToCrackOffsetOnCft());};
-            }
+              //if (_CfTOffsetMethod == 3){
+              //  q1Thom->getRefToCrackOffsetOnCft() = q0Thom->getCrackOffsetOnCft() + yieldThomason*R0/(Cft*R);
+              //  if (q1Thom->getRefToCrackOffsetOnCft() > 2.){Msg::Error("High predicted value of Cft Offset = %e",q1Thom->getRefToCrackOffsetOnCft());};
+              //}
               
             }
             else{
@@ -512,6 +519,7 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
           }
         }
       }
+    }
     else
     {
       // In the bulk: