diff --git a/Solver/dofManager.h b/Solver/dofManager.h
index e65223a469ac7da27dcd5a6b6e4f2e2ecbadac09..be34e9706df35e18da0795d13e4e219bdc57a2e2 100644
--- a/Solver/dofManager.h
+++ b/Solver/dofManager.h
@@ -323,10 +323,10 @@ class dofManager{
         itConstraint = constraints.find(R[i]);
         if (itConstraint != constraints.end())
         {
-          for (unsigned i=0;i<(itConstraint->second).linear.size();i++)
+          for (unsigned j=0;j<(itConstraint->second).linear.size();j++)
           {
-                  std::map<Dof, int>::iterator itC = unknown.find((itConstraint->second).linear[i].first); // lin dep in unknown ?!
-                  _current->addToRightHandSide(itC->second, m(i)*(itConstraint->second).linear[i].second);
+                  std::map<Dof, int>::iterator itC = unknown.find((itConstraint->second).linear[j].first); // lin dep in unknown ?!
+                  _current->addToRightHandSide(itC->second, m(i)*(itConstraint->second).linear[j].second);
           }
         }
       }
diff --git a/Solver/solverAlgorithms.h b/Solver/solverAlgorithms.h
index 33f5a2b88665b3f20f998fdbc9b85378207089d9..1c76468a9f22925d7d23865489e709d8b5a309af 100644
--- a/Solver/solverAlgorithms.h
+++ b/Solver/solverAlgorithms.h
@@ -70,7 +70,25 @@ template<class Iterator,class Assembler> void Assemble(BilinearTermBase &term,
     term.get(e, npts, GP, localMatrix);
     shapeFcts.getKeys(e, R);
     testFcts.getKeys(e, C);
+//    std::cout << "assembling normal test function ; lagrange trial function : " << std::endl;
+//    for (int i = 0 ; i < R.size() ; i++)
+//    {
+//      std::cout << "tests : " << R[i].getEntity() << ":" << R[i].getType() << std::endl ;
+//    }
+//    for (int i = 0 ; i < R.size() ; i++)
+//    {
+//      std::cout << "trial : " << C[i].getEntity() << ":" << C[i].getType() << std::endl ;
+//    }
     assembler.assemble(R, C, localMatrix);
+//    std::cout << "assembling lagrange test function ; normal trial function : " << std::endl;
+//    for (int i = 0 ; i < R.size() ; i++)
+//    {
+//      std::cout << "tests : " << C[i].getEntity() << ":" << C[i].getType() << std::endl ;
+//    }
+//    for (int i = 0 ; i < R.size() ; i++)
+//    {
+//      std::cout << "trial : " << R[i].getEntity() << ":" << R[i].getType() << std::endl ;
+//    }
     assembler.assemble(C, R, localMatrix.transpose());
   }
 }
diff --git a/Solver/terms.h b/Solver/terms.h
index 7a4ba1bdfd0fa8f912472c5d2d3fb86098e9dd30..36cf2f69deae0cb6140cec381f395ee5fae2f49c 100644
--- a/Solver/terms.h
+++ b/Solver/terms.h
@@ -201,7 +201,7 @@ class IsotropicElasticTerm : public BilinearTerm<SVector3,SVector3>
   virtual ~IsotropicElasticTerm() {}
   virtual void get(MElement *ele,int npts,IntPt *GP,fullMatrix<double> &m)
   {
-    if (ele->getParent()) ele=ele->getParent();    
+    if (ele->getParent()) ele=ele->getParent();
     if (sym)
     {
       int nbFF = BilinearTerm<SVector3,SVector3>::space1.getNumKeys(ele);
@@ -314,6 +314,7 @@ template<class T1> class LoadTerm : public LinearTerm<T1>
   }
 };
 
+
 class LagrangeMultiplierTerm : public BilinearTerm<SVector3,double>
 {
   SVector3 _d;
@@ -344,4 +345,70 @@ class LagrangeMultiplierTerm : public BilinearTerm<SVector3,double>
   }
 };
 
+
+class LagMultTerm : public BilinearTerm<SVector3, SVector3>
+{
+  SVector3 _d;
+ public :
+  LagMultTerm(FunctionSpace<SVector3>& space1_, FunctionSpace<SVector3>& space2_, const SVector3 &d) :
+    BilinearTerm<SVector3,SVector3>(space1_, space2_) {for(int i=0; i < 3; i++) _d(i) = d(i);}
+  virtual ~LagMultTerm() {}
+  virtual void get(MElement *ele, int npts, IntPt *GP, fullMatrix<double> &m)
+  {
+    int nbFF1 = BilinearTerm<SVector3,SVector3>::space1.getNumKeys(ele); //nbVertices*nbcomp of parent
+    int nbFF2 = BilinearTerm<SVector3,SVector3>::space2.getNumKeys(ele); //nbVertices of boundary
+    double jac[3][3];
+    m.resize(nbFF1, nbFF2);
+    m.setAll(0.);
+    for (int i = 0; i < npts; i++) {
+      double u = GP[i].pt[0]; double v = GP[i].pt[1]; double w = GP[i].pt[2];
+      const double weight = GP[i].weight; const double detJ = ele->getJacobian(u, v, w, jac);
+      std::vector<TensorialTraits<SVector3>::ValType> Vals;
+      std::vector<TensorialTraits<SVector3>::ValType> ValsT;
+      BilinearTerm<SVector3,SVector3>::space1.f(ele, u, v, w, Vals);
+      BilinearTerm<SVector3,SVector3>::space2.f(ele, u, v, w, ValsT);
+      for (int j = 0; j < nbFF1; j++) {
+        for (int k = 0; k < nbFF2; k++) {
+          m(j, k) += dot(Vals[j], ValsT[k]) * weight * detJ;
+        }
+      }
+    }
+    m.print();
+  }
+};
+
+
+template<class T1> class LoadTermOnBorder : public LinearTerm<T1>
+{
+  simpleFunction<typename TensorialTraits<T1>::ValType> &Load;
+ public :
+  LoadTerm(FunctionSpace<T1>& space1_,simpleFunction<typename TensorialTraits<T1>::ValType> &Load_) :LinearTerm<T1>(space1_),Load(Load_) {}
+  virtual ~LoadTerm() {}
+
+  virtual void get(MElement *ele,int npts,IntPt *GP,fullVector<double> &m)
+  {
+    MElement * elep;
+    if (ele->getParent()) elep=ele->getParent();
+    int nbFF=LinearTerm<T1>::space1.getNumKeys(ele);
+    double jac[3][3];
+    m.resize(nbFF);
+    m.scale(0.);
+    for (int i = 0; i < npts; i++)
+    {
+      const double u = GP[i].pt[0]; const double v = GP[i].pt[1]; const double w = GP[i].pt[2];
+      const double weight = GP[i].weight; const double detJ = ele->getJacobian(u, v, w, jac);
+      std::vector<typename TensorialTraits<T1>::ValType> Vals;
+      LinearTerm<T1>::space1.f(ele, u, v, w, Vals);
+      SPoint3 p;
+      ele->pnt(u, v, w, p);
+      typename TensorialTraits<T1>::ValType load=Load(p.x(),p.y(),p.z());
+      for (int j = 0; j < nbFF ; ++j)
+      {
+        m(j) += dot(Vals[j], load) * weight * detJ;
+      }
+    }
+  }
+};
+
+
 #endif// _TERMS_H_