checpoint with allowing thomason also in the interface

NonLinearSolver/materialLaw/mlawNonLocalPorousCoupled.cpp

@@ -562,20 +562,58 @@ void mlawNonLocalPorousCoupledLaw::checkCoalescence(IPNonLocalPorosity* q1, cons
     else
     {
       // In the bulk:
-        // Coalescence never occurs
+      if (!q0Thom->getCoalescenceOnsetFlag()){
+        double yieldThomason = _mlawCoales->yieldFunction(kcorEq,pcor,R,yieldfV,q0,q1,T);
+
+        if (yieldThomason > _tol){
+          q1Thom->getRefToCoalescenceCriterion() = 0.;
+          // update offset
+          if (_CfTOffsetMethod == 1){
+            q1Thom->getRefToCrackOffsetOnCft() = 1. + yieldThomason*R0/(Cft*R);
+            Msg::Info("coalescence occurs, use Thomason yield surface in next step with offset = %e",q1Thom->getCrackOffsetOnCft());
+          }
+          else if(_CfTOffsetMethod == 0){
+            q1Thom->getRefToCrackOffsetOnCft() = 1.;
+          }
+
+          q1Thom->getRefToCoalescenceActiveFlag() = true;
+          q1Thom->getRefToCoalescenceOnsetFlag() = true;
+          //Msg::Info("coalescene yield surface is used");
+          // initial coalesence onset value
+          q1Thom->getRefToPorosityAtCoalescenceOnset() = yieldfV;
+          q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q1Thom->getLigamentRatio();
+          q1Thom->getRefToAspectRatioAtCoalescenceOnset() = q1Thom->getAspectRatio();
+          q1Thom->getRefToShapeFactorAtCoalescenceOnset() = q1Thom->getShapeFactor();
+          q1Thom->getRefToAccelerateRate() = 1.;
+        }
+        else{
+          q1Thom->getRefToCoalescenceCriterion() = yieldThomason; // save data
+          q1Thom->getRefToCoalescenceActiveFlag() = false;
           q1Thom->getRefToCoalescenceOnsetFlag() = false;
-      // Reset onset variables to their default value
-        // Porosity
+
+          // reset to defaut value
           q1Thom->getRefToPorosityAtCoalescenceOnset() = 0.;
-        // Geometrical parameters
           q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = 0.;
           q1Thom->getRefToAspectRatioAtCoalescenceOnset() = 0.;
           q1Thom->getRefToShapeFactorAtCoalescenceOnset() = 0.;
 
-      q1Thom->getRefToCoalescenceActiveFlag() = false;
           q1Thom->getRefToCrackOffsetOnCft() = 1.;
           q1Thom->getRefToAccelerateRate() = 1.;
-      q1Thom->getRefToYieldOffset() = 0.;
+        }
+      }
+      else{
+        q1Thom->getRefToCoalescenceCriterion() = 0.;
+        q1Thom->getRefToCoalescenceActiveFlag() = true;
+        q1Thom->getRefToCoalescenceOnsetFlag() = true;
+
+        q1Thom->getRefToLigamentRatioAtCoalescenceOnset() = q0Thom->getLigamentRatioAtCoalescenceOnset();
+        q1Thom->getRefToAspectRatioAtCoalescenceOnset() = q0Thom->getAspectRatioAtCoalescenceOnset();
+        q1Thom->getRefToShapeFactorAtCoalescenceOnset() = q0Thom->getShapeFactorAtCoalescenceOnset();
+        q1Thom->getRefToPorosityAtCoalescenceOnset() = q0Thom->getPorosityAtCoalescenceOnset();
+
+        q1Thom->getRefToCrackOffsetOnCft() = q0Thom->getCrackOffsetOnCft();
+        q1Thom->getRefToAccelerateRate() = q0Thom->getAccelerateRate();
+      }
     }
   }
   else{