cpp_diff_for_cohesiveLaw

--- /home/lingran/YadeCodes/compare/old_CohesiveFrictionalContactLaw.cpp    2013-01-16 14:58:08.859546000 +0100
+++ /home/lingran/YadeCodes/compare/new_CohesiveFrictionalContactLaw.cpp    2013-08-26 10:40:56.821407000 +0200
@@ -41,6 +41,32 @@
    return shearEnergy;
 }

+Real Law2_ScGeom6D_CohFrictPhys_CohesionMoment::bendingElastEnergy()
+{
+   Real bendingEnergy=0;
+   FOREACH(const shared_ptr<Interaction>& I, *scene->interactions){
+       if(!I->isReal()) continue;
+       CohFrictPhys* phys = YADE_CAST<CohFrictPhys*>(I->phys.get());
+       if (phys) {
+           bendingEnergy += 0.5*(phys->moment_bending.squaredNorm()/phys->kr);
+       }
+   }
+   return bendingEnergy;
+}
+
+Real Law2_ScGeom6D_CohFrictPhys_CohesionMoment::twistElastEnergy()
+{
+   Real twistEnergy=0;
+   FOREACH(const shared_ptr<Interaction>& I, *scene->interactions){
+       if(!I->isReal()) continue;
+       CohFrictPhys* phys = YADE_CAST<CohFrictPhys*>(I->phys.get());
+       if (phys) {
+           twistEnergy += 0.5*(phys->moment_twist.squaredNorm()/phys->ktw);
+       }
+   }
+   return twistEnergy;
+}
+
 void CohesiveFrictionalContactLaw::action()
 {
    if(!functor) functor=shared_ptr<Law2_ScGeom6D_CohFrictPhys_CohesionMoment>(new Law2_ScGeom6D_CohFrictPhys_CohesionMoment);