diff --git a/NonLinearSolver/materialLaw/mlawNonLocalDamage.cpp b/NonLinearSolver/materialLaw/mlawNonLocalDamage.cpp
index 1bf302b69fdfa84393481300ab0cf2316aab5501..12088c4c56ba576a8dad3a76571bf7f352850ae1 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalDamage.cpp
+++ b/NonLinearSolver/materialLaw/mlawNonLocalDamage.cpp
@@ -1087,7 +1087,7 @@ void mlawNonLocalDamage::constitutive(const STensor3& F0,const STensor3& Fn,STen
     DirrEnergDF(2,1) = dEnePathFollowingdStrain[5]/sq2;
     
     STensor43 Cel, Sel;
-    ElasticStiffness(&Cel);
+    ElasticStiffness(statev_1,&Cel);
     STensorOperation::inverseSTensor43(Cel,Sel);     
     static STensor3 I(1.);
     for (int i=0; i<3; i++){
@@ -2109,6 +2109,34 @@ const double mlawNonLocalDamage::poissonRatio() const
 {
   return _nu0;
 }
+void mlawNonLocalDamage::ElasticStiffness(double *statev, STensor43* elasticTensor) const
+{
+  fullMatrix<double> Cel_fullMatrix(6,6);
+  static double **Cel_matrix;
+  static bool init=false;
+  if(init==false)
+  {
+    mallocmatrix(&Cel_matrix,6,6);
+    init=true;
+  }
+  mat->get_elOp(statev,Cel_matrix);
+  for(int i=0;i<6;i++)
+  {
+    for(int j=0;j<6;j++){
+      Cel_fullMatrix(i,j)=Cel_matrix[i][j];
+    }
+  }
+  for(int i=3;i<6;i++)
+  {
+    for(int j=0;j<3;j++)
+    {
+      Cel_fullMatrix(i,j)*=(1./sq2);
+      Cel_fullMatrix(j,i)*=2.;
+      Cel_fullMatrix(j,i)*=(1./sq2);
+    }
+  }
+  STensorOperation::fromFullMatrixToSTensor43(Cel_fullMatrix,(*elasticTensor));
+};
 
 void mlawNonLocalDamage::ElasticStiffness(STensor43* elasticTensor) const
 {
diff --git a/NonLinearSolver/materialLaw/mlawNonLocalDamage.h b/NonLinearSolver/materialLaw/mlawNonLocalDamage.h
index 2ded6164f8ccdf67f7a070217b723f0bd29e2bc4..971bd9caf2e8cce570c54c8927ea9d13e640bff9 100644
--- a/NonLinearSolver/materialLaw/mlawNonLocalDamage.h
+++ b/NonLinearSolver/materialLaw/mlawNonLocalDamage.h
@@ -97,6 +97,7 @@ class mlawNonLocalDamage : public materialLaw
   virtual int getNbNonLocalVariables() const {return nlVar;};
   virtual void updateCharacteristicLengthTensor(double** chara_Length, double *statev) const {};
   virtual void ElasticStiffness(STensor43* elasticTensor) const;
+  virtual void ElasticStiffness(double *statev, STensor43* elasticTensor) const;
   virtual void ElasticStiffnessMtx(double *statev, STensor43* elasticTensor) const;
   virtual void ElasticStiffnessInc(double *statev, STensor43* elasticTensor) const;
   virtual void ElasticStiffness_incOri(double e1, double e2, double e3,STensor43* elasticTensor) const;
diff --git a/dG3D/benchmarks/shiftMultiSystem/fullPlateHoleRes.py b/dG3D/benchmarks/shiftMultiSystem/fullPlateHoleRes.py
index 2679591a0e10f4efdd9f144f07c3acbcca492d87..1946f9eb0f4a2a6138f6c74feca3407321a258af 100644
--- a/dG3D/benchmarks/shiftMultiSystem/fullPlateHoleRes.py
+++ b/dG3D/benchmarks/shiftMultiSystem/fullPlateHoleRes.py
@@ -268,7 +268,7 @@ check.equal(6.859167e+03,mysolver.getArchivedForceOnPhysicalGroup("Face", 102, 0
 
 data21 = csv.reader(open('force102comp0_part2.csv'), delimiter=';')
 force2 = list(data21)
-check.equal(1.236485e+03,float(force2[-1][1]),1e-5)
+check.equal(1.235198e+03,float(force2[-1][1]),1e-5)
 
 
 #mysolver.disableResetRestart() #only battery: the second solve mimics a restart
diff --git a/dG3D/benchmarks/shiftMultiSystem/run.py b/dG3D/benchmarks/shiftMultiSystem/run.py
index f8fa40f7bbe04b10f2cbbae1cf9e04114024e553..9bc6a8f77a2aeda75030b8b60dc14076357431bc 100644
--- a/dG3D/benchmarks/shiftMultiSystem/run.py
+++ b/dG3D/benchmarks/shiftMultiSystem/run.py
@@ -13,7 +13,7 @@ else:
 
 data1 = csv.reader(open('force102comp0_part2.csv'), delimiter=';')
 force = list(data1)
-checkEqual(1.236485e+03,float(force[-1][1]),1e-5)
+checkEqual(1.235198e+03,float(force[-1][1]),1e-5)
 
 if sys.version_info[0] < 3:
   os.system('mpiexec -np 4 python fullPlateHoleRes.py')
@@ -22,5 +22,5 @@ else:
 
 data21 = csv.reader(open('force102comp0_part2.csv'), delimiter=';')
 force2 = list(data21)
-checkEqual(1.236485e+03,float(force2[-1][1]),1e-5)
+checkEqual(1.235198e+03,float(force2[-1][1]),1e-5)
 
diff --git a/dG3D/benchmarks/smpPhenomenological_Mech/SMP_mech.py b/dG3D/benchmarks/smpPhenomenological_Mech/SMP_mech.py
index 4913b895c9bb887be43ca06ade8cf624c71db358..c186e3f64f7cedd2e555131ad5d704f1df1e10dd 100644
--- a/dG3D/benchmarks/smpPhenomenological_Mech/SMP_mech.py
+++ b/dG3D/benchmarks/smpPhenomenological_Mech/SMP_mech.py
@@ -491,7 +491,7 @@ mysolver.solve()
 
 # ===Test==========================================================================================
 check= TestCheck()		
-check.equal(-7.604566e-05,mysolver.getArchivedNodalValue(7,0,mysolver.displacement),1e-5)
+check.equal(-7.612266e-05,mysolver.getArchivedNodalValue(7,0,mysolver.displacement),1e-5)
 
 
 
diff --git a/dG3D/benchmarks/smpPhenomenological_Thermo/SMP_thermo.py b/dG3D/benchmarks/smpPhenomenological_Thermo/SMP_thermo.py
index bd135e6b30c6517956531014203bbdf363b5bdf9..2bc87c102b309288c5a211435916108927d012c1 100644
--- a/dG3D/benchmarks/smpPhenomenological_Thermo/SMP_thermo.py
+++ b/dG3D/benchmarks/smpPhenomenological_Thermo/SMP_thermo.py
@@ -484,7 +484,7 @@ mysolver.solve()
 
 
 check = TestCheck()
-check.equal(-6.772246e-01,mysolver.getArchivedForceOnPhysicalGroup("Face", 5678, 3),1.e-6)
-check.equal(8.404647e+02,mysolver.getArchivedNodalValue(1,3,mysolver.displacement),1.e-6)
-check.equal(1.637999e-04,mysolver.getArchivedNodalValue(7,0,mysolver.displacement),1.e-6)
+check.equal(-6.735810e-01,mysolver.getArchivedForceOnPhysicalGroup("Face", 5678, 3),1.e-6)
+check.equal(8.380272e+02,mysolver.getArchivedNodalValue(1,3,mysolver.displacement),1.e-6)
+check.equal(1.606705e-04,mysolver.getArchivedNodalValue(7,0,mysolver.displacement),1.e-6)
 
diff --git a/dG3D/benchmarks/voidedRVE/voidedRVE.py b/dG3D/benchmarks/voidedRVE/voidedRVE.py
index 188b80ae1ec624ab485339ebf4caa2f2759fc27d..d4d3856afb867952d14ca7915587f34e35b43dd7 100644
--- a/dG3D/benchmarks/voidedRVE/voidedRVE.py
+++ b/dG3D/benchmarks/voidedRVE/voidedRVE.py
@@ -206,5 +206,5 @@ except:
     os._exit(1)
 else:
     st = float(lstress.split(';')[2])
-    check.equal(4.640682e+01,st,1.e-4)
+    check.equal(3.951369e+01,st,1.e-4)