diff --git a/Numeric/Numeric.cpp b/Numeric/Numeric.cpp
index ad675906f89e03511aacd2c30403b5e2aa3c782a..04fa3ecb80bc75f3bbe019404d7375866cb0b0da 100644
--- a/Numeric/Numeric.cpp
+++ b/Numeric/Numeric.cpp
@@ -1,4 +1,4 @@
-// $Id: Numeric.cpp,v 1.21 2005-01-08 20:15:13 geuzaine Exp $
+// $Id: Numeric.cpp,v 1.22 2005-01-12 19:10:41 geuzaine Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -101,6 +101,13 @@ void prodve(double a[3], double b[3], double c[3])
   c[0] = a[1] * b[2] - a[2] * b[1];
 }
 
+void matvec(double mat[3][3], double vec[3], double res[3])
+{
+  res[0] = mat[0][0] * vec[0] + mat[0][1] * vec[1] + mat[0][2] * vec[2];
+  res[1] = mat[1][0] * vec[0] + mat[1][1] * vec[1] + mat[1][2] * vec[2];
+  res[2] = mat[2][0] * vec[0] + mat[2][1] * vec[1] + mat[2][2] * vec[2];
+}
+
 void prosca(double a[3], double b[3], double *c)
 {
   *c = a[0] * b[0] + a[1] * b[1] + a[2] * b[2];
@@ -232,27 +239,28 @@ int sys3x3_with_tol(double mat[3][3], double b[3], double res[3], double *det)
   return out;
 }
 
-int inv3x3(double mat[3][3], double inv[3][3], double *det)
+double inv3x3(double mat[3][3], double inv[3][3])
 {
-  double ud;
-
-  *det = det3x3(mat);
-
-  if(*det == 0.0)
-    return (0);
-
-  ud = 1. / (*det);
-
-  inv[0][0] = ud * (mat[1][1] * mat[2][2] - mat[1][2] * mat[2][1]);
-  inv[0][1] = -ud * (mat[1][0] * mat[2][2] - mat[1][2] * mat[2][0]);
-  inv[0][2] = ud * (mat[1][0] * mat[2][1] - mat[1][1] * mat[2][0]);
-  inv[1][0] = -ud * (mat[0][1] * mat[2][2] - mat[0][2] * mat[2][1]);
-  inv[1][1] = ud * (mat[0][0] * mat[2][2] - mat[0][2] * mat[2][0]);
-  inv[1][2] = -ud * (mat[0][0] * mat[2][1] - mat[0][1] * mat[2][0]);
-  inv[2][0] = ud * (mat[0][1] * mat[1][2] - mat[0][2] * mat[1][1]);
-  inv[2][1] = -ud * (mat[0][0] * mat[1][2] - mat[0][2] * mat[1][0]);
-  inv[2][2] = ud * (mat[0][0] * mat[1][1] - mat[0][1] * mat[1][0]);
-  return (1);
+  double det = det3x3(mat);
+  if(det){
+    double ud = 1. / det;
+    inv[0][0] =  ( mat[1][1] * mat[2][2] - mat[1][2] * mat[2][1] ) * ud ;
+    inv[1][0] = -( mat[1][0] * mat[2][2] - mat[1][2] * mat[2][0] ) * ud ;
+    inv[2][0] =  ( mat[1][0] * mat[2][1] - mat[1][1] * mat[2][0] ) * ud ;
+    inv[0][1] = -( mat[0][1] * mat[2][2] - mat[0][2] * mat[2][1] ) * ud ;
+    inv[1][1] =  ( mat[0][0] * mat[2][2] - mat[0][2] * mat[2][0] ) * ud ;
+    inv[2][1] = -( mat[0][0] * mat[2][1] - mat[0][1] * mat[2][0] ) * ud ;
+    inv[0][2] =  ( mat[0][1] * mat[1][2] - mat[0][2] * mat[1][1] ) * ud ;
+    inv[1][2] = -( mat[0][0] * mat[1][2] - mat[0][2] * mat[1][0] ) * ud ;
+    inv[2][2] =  ( mat[0][0] * mat[1][1] - mat[0][1] * mat[1][0] ) * ud ;
+  }
+  else{
+    Msg(GERROR, "Singular matrix");
+    for(int i = 0; i < 3; i++)
+      for(int j = 0; j < 3; j++)
+	inv[i][j] = 0.;
+  }
+  return det;
 }
 
 double angle_02pi(double A3)
diff --git a/Numeric/Numeric.h b/Numeric/Numeric.h
index 109c42c9c04c82069e2e905650e723538adccb8e..a04f8b790928e6b2993410492ba18bfdbe2c577a 100644
--- a/Numeric/Numeric.h
+++ b/Numeric/Numeric.h
@@ -61,6 +61,7 @@ double myasin(double a);
 double myacos(double a);
 void prodve(double a[3], double b[3], double c[3]);
 void prosca(double a[3], double b[3], double *c);
+void matvec(double mat[3][3], double vec[3], double res[3]);
 double norme(double a[3]);
 void normal3points(double x0, double y0, double z0,
 		   double x1, double y1, double z1,
@@ -72,7 +73,7 @@ int sys3x3_with_tol(double mat[3][3], double b[3], double res[3], double *det);
 double det3x3(double mat[3][3]);
 double trace3x3(double mat[3][3]);
 double trace2 (double mat[3][3]);
-int inv3x3(double mat[3][3], double inv[3][3], double *det);
+double inv3x3(double mat[3][3], double inv[3][3]);
 double angle_02pi(double A3);
 void eigenvalue(double mat[3][3], double re[3]);
 void FindCubicRoots(const double coeff[4], double re[3], double im[3]);
diff --git a/Plugin/Curl.cpp b/Plugin/Curl.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..6ccb7da6c48c235d27dabd7987cc7f9c67e10f26
--- /dev/null
+++ b/Plugin/Curl.cpp
@@ -0,0 +1,146 @@
+// $Id: Curl.cpp,v 1.1 2005-01-12 19:10:41 geuzaine Exp $
+//
+// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
+// USA.
+// 
+// Please report all bugs and problems to <gmsh@geuz.org>.
+
+#include "Plugin.h"
+#include "Curl.h"
+#include "List.h"
+#include "Views.h"
+#include "Context.h"
+#include "Numeric.h"
+#include "ShapeFunctions.h"
+
+extern Context_T CTX;
+
+StringXNumber CurlOptions_Number[] = {
+  {GMSH_FULLRC, "iView", NULL, -1.}
+};
+
+extern "C"
+{
+  GMSH_Plugin *GMSH_RegisterCurlPlugin()
+  {
+    return new GMSH_CurlPlugin();
+  }
+}
+
+GMSH_CurlPlugin::GMSH_CurlPlugin()
+{
+  ;
+}
+
+void GMSH_CurlPlugin::getName(char *name) const
+{
+  strcpy(name, "Curl");
+}
+
+void GMSH_CurlPlugin::getInfos(char *author, char *copyright,
+				    char *help_text) const
+{
+  strcpy(author, "C. Geuzaine (geuz@geuz.org)");
+  strcpy(copyright, "DGR (www.multiphysics.com)");
+  strcpy(help_text,
+	 "Plugin(Curl) computes the curl of the field\n"
+	 "in the view `iView'. If `iView' < 0, the plugin\n"
+	 "is run on the current view.\n"
+	 "\n"
+	 "Plugin(Curl) creates one new view.\n");
+}
+
+int GMSH_CurlPlugin::getNbOptions() const
+{
+  return sizeof(CurlOptions_Number) / sizeof(StringXNumber);
+}
+
+StringXNumber *GMSH_CurlPlugin::getOption(int iopt)
+{
+  return &CurlOptions_Number[iopt];
+}
+
+void GMSH_CurlPlugin::catchErrorMessage(char *errorMessage) const
+{
+  strcpy(errorMessage, "Curl failed...");
+}
+
+static void curl(int inNb, List_T *inList, int *outNb, List_T *outList, 
+		     int dim, int nbNod, int nbTime)
+{
+  if(!inNb) return;
+  
+  int nb = List_Nbr(inList) / inNb;
+  for(int i = 0; i < List_Nbr(inList); i += nb) {
+    double *x = (double *)List_Pointer_Fast(inList, i);
+    double *y = (double *)List_Pointer_Fast(inList, i + nbNod);
+    double *z = (double *)List_Pointer_Fast(inList, i + 2 * nbNod);
+    elementFactory factory;
+    element *element = factory.create(nbNod, dim, x, y, z);
+    if(!element) return;
+    for(int j = 0; j < 3 * nbNod; j++)
+      List_Add(outList, &x[j]);
+    for(int j = 0; j < nbTime; j++){
+      double *val = (double *)List_Pointer(inList, i + 3 * nbNod + nbNod * 3 * j);
+      for(int k = 0; k < nbNod; k++){
+	double u, v, w, f[3];
+	element->getNode(k, u, v, w);
+	element->interpolateCurl(val, u, v, w, f, 3);
+	List_Add(outList, &f[0]);
+	List_Add(outList, &f[1]);
+	List_Add(outList, &f[2]);
+      }
+    }
+    delete element;
+    (*outNb)++;
+  }
+}
+
+Post_View *GMSH_CurlPlugin::execute(Post_View * v)
+{
+  int iView = (int)CurlOptions_Number[0].def;
+  
+  if(iView < 0)
+    iView = v ? v->Index : 0;
+  
+  if(!List_Pointer_Test(CTX.post.list, iView)) {
+    Msg(GERROR, "View[%d] does not exist", iView);
+    return v;
+  }
+  
+  Post_View *v1 = *(Post_View **)List_Pointer(CTX.post.list, iView);
+  Post_View *v2 = BeginView(1);
+
+  curl(v1->NbVL, v1->VL, &v2->NbVL, v2->VL, 1, 2, v1->NbTimeStep);
+  curl(v1->NbVT, v1->VT, &v2->NbVT, v2->VT, 2, 3, v1->NbTimeStep);
+  curl(v1->NbVQ, v1->VQ, &v2->NbVQ, v2->VQ, 2, 4, v1->NbTimeStep);
+  curl(v1->NbVS, v1->VS, &v2->NbVS, v2->VS, 3, 4, v1->NbTimeStep);
+  curl(v1->NbVH, v1->VH, &v2->NbVH, v2->VH, 3, 8, v1->NbTimeStep);
+  curl(v1->NbVI, v1->VI, &v2->NbVI, v2->VI, 3, 6, v1->NbTimeStep);
+  curl(v1->NbVY, v1->VY, &v2->NbVY, v2->VY, 3, 5, v1->NbTimeStep);
+
+  // copy time data
+  for(int i = 0; i < List_Nbr(v1->Time); i++)
+    List_Add(v2->Time, List_Pointer(v1->Time, i));
+  // finalize
+  char name[1024], filename[1024];
+  sprintf(name, "%s_Curl", v1->Name);
+  sprintf(filename, "%s_Curl.pos", v1->Name);
+  EndView(v2, 1, filename, name);
+  
+  return v2;
+}
diff --git a/Plugin/Curl.h b/Plugin/Curl.h
new file mode 100644
index 0000000000000000000000000000000000000000..1742e9a91adc640e539b01512e48def2bcf51da8
--- /dev/null
+++ b/Plugin/Curl.h
@@ -0,0 +1,42 @@
+#ifndef _CURL_H_
+#define _CURL_H_
+
+// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
+// USA.
+// 
+// Please report all bugs and problems to <gmsh@geuz.org>.
+
+#include "Plugin.h"
+
+extern "C"
+{
+  GMSH_Plugin *GMSH_RegisterCurlPlugin();
+}
+
+class GMSH_CurlPlugin : public GMSH_Post_Plugin
+{
+ public:
+  GMSH_CurlPlugin();
+  void getName(char *name) const;
+  void getInfos(char *author, char *copyright, char *help_text) const;
+  void catchErrorMessage(char *errorMessage) const;
+  int getNbOptions() const;
+  StringXNumber* getOption(int iopt);  
+  Post_View *execute(Post_View *);
+};
+
+#endif
diff --git a/Plugin/Divergence.cpp b/Plugin/Divergence.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..61f160e743498edc6110ebb76cb36bb731b1e9cd
--- /dev/null
+++ b/Plugin/Divergence.cpp
@@ -0,0 +1,144 @@
+// $Id: Divergence.cpp,v 1.1 2005-01-12 19:10:41 geuzaine Exp $
+//
+// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
+// USA.
+// 
+// Please report all bugs and problems to <gmsh@geuz.org>.
+
+#include "Plugin.h"
+#include "Divergence.h"
+#include "List.h"
+#include "Views.h"
+#include "Context.h"
+#include "Numeric.h"
+#include "ShapeFunctions.h"
+
+extern Context_T CTX;
+
+StringXNumber DivergenceOptions_Number[] = {
+  {GMSH_FULLRC, "iView", NULL, -1.}
+};
+
+extern "C"
+{
+  GMSH_Plugin *GMSH_RegisterDivergencePlugin()
+  {
+    return new GMSH_DivergencePlugin();
+  }
+}
+
+GMSH_DivergencePlugin::GMSH_DivergencePlugin()
+{
+  ;
+}
+
+void GMSH_DivergencePlugin::getName(char *name) const
+{
+  strcpy(name, "Divergence");
+}
+
+void GMSH_DivergencePlugin::getInfos(char *author, char *copyright,
+				    char *help_text) const
+{
+  strcpy(author, "C. Geuzaine (geuz@geuz.org)");
+  strcpy(copyright, "DGR (www.multiphysics.com)");
+  strcpy(help_text,
+	 "Plugin(Divergence) computes the divergence of the\n"
+	 "field in the view `iView'. If `iView' < 0, the plugin\n"
+	 "is run on the current view.\n"
+	 "\n"
+	 "Plugin(Divergence) creates one new view.\n");
+}
+
+int GMSH_DivergencePlugin::getNbOptions() const
+{
+  return sizeof(DivergenceOptions_Number) / sizeof(StringXNumber);
+}
+
+StringXNumber *GMSH_DivergencePlugin::getOption(int iopt)
+{
+  return &DivergenceOptions_Number[iopt];
+}
+
+void GMSH_DivergencePlugin::catchErrorMessage(char *errorMessage) const
+{
+  strcpy(errorMessage, "Divergence failed...");
+}
+
+static void divergence(int inNb, List_T *inList, int *outNb, List_T *outList, 
+		     int dim, int nbNod, int nbTime)
+{
+  if(!inNb) return;
+  
+  int nb = List_Nbr(inList) / inNb;
+  for(int i = 0; i < List_Nbr(inList); i += nb) {
+    double *x = (double *)List_Pointer_Fast(inList, i);
+    double *y = (double *)List_Pointer_Fast(inList, i + nbNod);
+    double *z = (double *)List_Pointer_Fast(inList, i + 2 * nbNod);
+    elementFactory factory;
+    element *element = factory.create(nbNod, dim, x, y, z);
+    if(!element) return;
+    for(int j = 0; j < 3 * nbNod; j++)
+      List_Add(outList, &x[j]);
+    for(int j = 0; j < nbTime; j++){
+      double *val = (double *)List_Pointer(inList, i + 3 * nbNod + nbNod * 3 * j);
+      for(int k = 0; k < nbNod; k++){
+	double u, v, w;
+	element->getNode(k, u, v, w);
+	double f = element->interpolateDiv(val, u, v, w, 3);
+	List_Add(outList, &f);
+      }
+    }
+    delete element;
+    (*outNb)++;
+  }
+}
+
+Post_View *GMSH_DivergencePlugin::execute(Post_View * v)
+{
+  int iView = (int)DivergenceOptions_Number[0].def;
+  
+  if(iView < 0)
+    iView = v ? v->Index : 0;
+  
+  if(!List_Pointer_Test(CTX.post.list, iView)) {
+    Msg(GERROR, "View[%d] does not exist", iView);
+    return v;
+  }
+  
+  Post_View *v1 = *(Post_View **)List_Pointer(CTX.post.list, iView);
+  Post_View *v2 = BeginView(1);
+
+  divergence(v1->NbVL, v1->VL, &v2->NbSL, v2->SL, 1, 2, v1->NbTimeStep);
+  divergence(v1->NbVT, v1->VT, &v2->NbST, v2->ST, 2, 3, v1->NbTimeStep);
+  divergence(v1->NbVQ, v1->VQ, &v2->NbSQ, v2->SQ, 2, 4, v1->NbTimeStep);
+  divergence(v1->NbVS, v1->VS, &v2->NbSS, v2->SS, 3, 4, v1->NbTimeStep);
+  divergence(v1->NbVH, v1->VH, &v2->NbSH, v2->SH, 3, 8, v1->NbTimeStep);
+  divergence(v1->NbVI, v1->VI, &v2->NbSI, v2->SI, 3, 6, v1->NbTimeStep);
+  divergence(v1->NbVY, v1->VY, &v2->NbSY, v2->SY, 3, 5, v1->NbTimeStep);
+
+  // copy time data
+  for(int i = 0; i < List_Nbr(v1->Time); i++)
+    List_Add(v2->Time, List_Pointer(v1->Time, i));
+  // finalize
+  char name[1024], filename[1024];
+  sprintf(name, "%s_Divergence", v1->Name);
+  sprintf(filename, "%s_Divergence.pos", v1->Name);
+  EndView(v2, 1, filename, name);
+  
+  return v2;
+}
diff --git a/Plugin/Divergence.h b/Plugin/Divergence.h
new file mode 100644
index 0000000000000000000000000000000000000000..bc7e7e63fc1d81789746614fc36c843ec34e97c0
--- /dev/null
+++ b/Plugin/Divergence.h
@@ -0,0 +1,42 @@
+#ifndef _DIVERGENCE_H_
+#define _DIVERGENCE_H_
+
+// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
+// USA.
+// 
+// Please report all bugs and problems to <gmsh@geuz.org>.
+
+#include "Plugin.h"
+
+extern "C"
+{
+  GMSH_Plugin *GMSH_RegisterDivergencePlugin();
+}
+
+class GMSH_DivergencePlugin : public GMSH_Post_Plugin
+{
+ public:
+  GMSH_DivergencePlugin();
+  void getName(char *name) const;
+  void getInfos(char *author, char *copyright, char *help_text) const;
+  void catchErrorMessage(char *errorMessage) const;
+  int getNbOptions() const;
+  StringXNumber* getOption(int iopt);  
+  Post_View *execute(Post_View *);
+};
+
+#endif
diff --git a/Plugin/Gradient.cpp b/Plugin/Gradient.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..d2e94fc4233c0f8c5397b81a6c0eda498851b50c
--- /dev/null
+++ b/Plugin/Gradient.cpp
@@ -0,0 +1,156 @@
+// $Id: Gradient.cpp,v 1.5 2005-01-12 19:10:41 geuzaine Exp $
+//
+// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
+// USA.
+// 
+// Please report all bugs and problems to <gmsh@geuz.org>.
+
+#include "Plugin.h"
+#include "Gradient.h"
+#include "List.h"
+#include "Views.h"
+#include "Context.h"
+#include "Numeric.h"
+#include "ShapeFunctions.h"
+
+extern Context_T CTX;
+
+StringXNumber GradientOptions_Number[] = {
+  {GMSH_FULLRC, "iView", NULL, -1.}
+};
+
+extern "C"
+{
+  GMSH_Plugin *GMSH_RegisterGradientPlugin()
+  {
+    return new GMSH_GradientPlugin();
+  }
+}
+
+GMSH_GradientPlugin::GMSH_GradientPlugin()
+{
+  ;
+}
+
+void GMSH_GradientPlugin::getName(char *name) const
+{
+  strcpy(name, "Gradient");
+}
+
+void GMSH_GradientPlugin::getInfos(char *author, char *copyright,
+				    char *help_text) const
+{
+  strcpy(author, "C. Geuzaine (geuz@geuz.org)");
+  strcpy(copyright, "DGR (www.multiphysics.com)");
+  strcpy(help_text,
+	 "Plugin(Gradient) computes the gradient of the\n"
+	 "field in the view `iView'. If `iView' < 0, the\n"
+	 "plugin is run on the current view.\n"
+	 "\n"
+	 "Plugin(Gradient) creates one new view.\n");
+}
+
+int GMSH_GradientPlugin::getNbOptions() const
+{
+  return sizeof(GradientOptions_Number) / sizeof(StringXNumber);
+}
+
+StringXNumber *GMSH_GradientPlugin::getOption(int iopt)
+{
+  return &GradientOptions_Number[iopt];
+}
+
+void GMSH_GradientPlugin::catchErrorMessage(char *errorMessage) const
+{
+  strcpy(errorMessage, "Gradient failed...");
+}
+
+static void gradient(int inNb, List_T *inList, int *outNb, List_T *outList, 
+		     int dim, int nbNod, int nbComp, int nbTime)
+{
+  if(!inNb) return;
+  
+  int nb = List_Nbr(inList) / inNb;
+  for(int i = 0; i < List_Nbr(inList); i += nb) {
+    double *x = (double *)List_Pointer_Fast(inList, i);
+    double *y = (double *)List_Pointer_Fast(inList, i + nbNod);
+    double *z = (double *)List_Pointer_Fast(inList, i + 2 * nbNod);
+    elementFactory factory;
+    element *element = factory.create(nbNod, dim, x, y, z);
+    if(!element) return;
+    for(int j = 0; j < 3 * nbNod; j++)
+      List_Add(outList, &x[j]);
+    for(int j = 0; j < nbTime; j++){
+      double *val = (double *)List_Pointer(inList, i + 3 * nbNod + nbNod * nbComp * j);
+      for(int k = 0; k < nbNod; k++){
+	double u, v, w, f[3];
+	element->getNode(k, u, v, w);
+	for(int l = 0; l < nbComp; l++){
+	  element->interpolateGrad(val+l, u, v, w, f, nbComp);
+	  List_Add(outList, &f[0]);
+	  List_Add(outList, &f[1]);
+	  List_Add(outList, &f[2]);
+	}
+      }
+    }
+    delete element;
+    (*outNb)++;
+  }
+}
+
+Post_View *GMSH_GradientPlugin::execute(Post_View * v)
+{
+  int iView = (int)GradientOptions_Number[0].def;
+  
+  if(iView < 0)
+    iView = v ? v->Index : 0;
+  
+  if(!List_Pointer_Test(CTX.post.list, iView)) {
+    Msg(GERROR, "View[%d] does not exist", iView);
+    return v;
+  }
+  
+  Post_View *v1 = *(Post_View **)List_Pointer(CTX.post.list, iView);
+  Post_View *v2 = BeginView(1);
+  
+  gradient(v1->NbSL, v1->SL, &v2->NbVL, v2->VL, 1, 2, 1, v1->NbTimeStep);
+  gradient(v1->NbST, v1->ST, &v2->NbVT, v2->VT, 2, 3, 1, v1->NbTimeStep);
+  gradient(v1->NbSQ, v1->SQ, &v2->NbVQ, v2->VQ, 2, 4, 1, v1->NbTimeStep);
+  gradient(v1->NbSS, v1->SS, &v2->NbVS, v2->VS, 3, 4, 1, v1->NbTimeStep);
+  gradient(v1->NbSH, v1->SH, &v2->NbVH, v2->VH, 3, 8, 1, v1->NbTimeStep);
+  gradient(v1->NbSI, v1->SI, &v2->NbVI, v2->VI, 3, 6, 1, v1->NbTimeStep);
+  gradient(v1->NbSY, v1->SY, &v2->NbVY, v2->VY, 3, 5, 1, v1->NbTimeStep);
+
+  gradient(v1->NbVL, v1->VL, &v2->NbTL, v2->TL, 1, 2, 3, v1->NbTimeStep);
+  gradient(v1->NbVT, v1->VT, &v2->NbTT, v2->TT, 2, 3, 3, v1->NbTimeStep);
+  gradient(v1->NbVQ, v1->VQ, &v2->NbTQ, v2->TQ, 2, 4, 3, v1->NbTimeStep);
+  gradient(v1->NbVS, v1->VS, &v2->NbTS, v2->TS, 3, 4, 3, v1->NbTimeStep);
+  gradient(v1->NbVH, v1->VH, &v2->NbTH, v2->TH, 3, 8, 3, v1->NbTimeStep);
+  gradient(v1->NbVI, v1->VI, &v2->NbTI, v2->TI, 3, 6, 3, v1->NbTimeStep);
+  gradient(v1->NbVY, v1->VY, &v2->NbTY, v2->TY, 3, 5, 3, v1->NbTimeStep);
+
+  // copy time data
+  for(int i = 0; i < List_Nbr(v1->Time); i++)
+    List_Add(v2->Time, List_Pointer(v1->Time, i));
+  // finalize
+  char name[1024], filename[1024];
+  sprintf(name, "%s_Gradient", v1->Name);
+  sprintf(filename, "%s_Gradient.pos", v1->Name);
+  EndView(v2, 1, filename, name);
+  
+  return v2;
+}
diff --git a/Plugin/Gradient.h b/Plugin/Gradient.h
new file mode 100644
index 0000000000000000000000000000000000000000..3e38b587e34c940b4b613d4d1abb6891fb84088c
--- /dev/null
+++ b/Plugin/Gradient.h
@@ -0,0 +1,42 @@
+#ifndef _GRADIENT_H_
+#define _GRADIENT_H_
+
+// Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
+//
+// This program is free software; you can redistribute it and/or modify
+// it under the terms of the GNU General Public License as published by
+// the Free Software Foundation; either version 2 of the License, or
+// (at your option) any later version.
+//
+// This program is distributed in the hope that it will be useful,
+// but WITHOUT ANY WARRANTY; without even the implied warranty of
+// MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
+// GNU General Public License for more details.
+//
+// You should have received a copy of the GNU General Public License
+// along with this program; if not, write to the Free Software
+// Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA 02111-1307
+// USA.
+// 
+// Please report all bugs and problems to <gmsh@geuz.org>.
+
+#include "Plugin.h"
+
+extern "C"
+{
+  GMSH_Plugin *GMSH_RegisterGradientPlugin();
+}
+
+class GMSH_GradientPlugin : public GMSH_Post_Plugin
+{
+ public:
+  GMSH_GradientPlugin();
+  void getName(char *name) const;
+  void getInfos(char *author, char *copyright, char *help_text) const;
+  void catchErrorMessage(char *errorMessage) const;
+  int getNbOptions() const;
+  StringXNumber* getOption(int iopt);  
+  Post_View *execute(Post_View *);
+};
+
+#endif
diff --git a/Plugin/Integrate.cpp b/Plugin/Integrate.cpp
index b99d7aae8463c6bb6cded22beecc1490c7bea987..f961d480e533f08a89a1afe81288121751a1a776 100644
--- a/Plugin/Integrate.cpp
+++ b/Plugin/Integrate.cpp
@@ -1,4 +1,4 @@
-// $Id: Integrate.cpp,v 1.13 2005-01-08 20:15:19 geuzaine Exp $
+// $Id: Integrate.cpp,v 1.14 2005-01-12 19:10:41 geuzaine Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -101,10 +101,7 @@ static double integrate(int nbList, List_T *list, int dim,
 					    nbNod * nbComp * step);
     elementFactory factory;
     element *element = factory.create(nbNod, dim, x, y, z);
-    if(!element){
-      Msg(GERROR, "Unknown type of element (dim=%d, nbNod=%d)", dim, nbNod);
-      return 0.;
-    }
+    if(!element) return 0.;
     if(nbComp == 1){
       if(!levelsetPositive) 
 	res += element->integrate(v);
diff --git a/Plugin/Lambda2.cpp b/Plugin/Lambda2.cpp
index 4c270cdc37f68c0a838c73ab6797de33b1d7f199..dbe3d98ca2f8a745e441f34e1988bfe544f014f4 100644
--- a/Plugin/Lambda2.cpp
+++ b/Plugin/Lambda2.cpp
@@ -1,4 +1,4 @@
-// $Id: Lambda2.cpp,v 1.8 2005-01-08 20:15:19 geuzaine Exp $
+// $Id: Lambda2.cpp,v 1.9 2005-01-12 19:10:41 geuzaine Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -93,6 +93,29 @@ void GMSH_Lambda2Plugin::catchErrorMessage(char *errorMessage) const
   strcpy(errorMessage, "Lambda2 failed...");
 }
 
+static int inv3x3tran(double mat[3][3], double inv[3][3], double *det)
+{
+  double ud;
+
+  *det = det3x3(mat);
+
+  if(*det == 0.0)
+    return (0);
+
+  ud = 1. / (*det);
+
+  inv[0][0] = ud * (mat[1][1] * mat[2][2] - mat[1][2] * mat[2][1]);
+  inv[0][1] = -ud * (mat[1][0] * mat[2][2] - mat[1][2] * mat[2][0]);
+  inv[0][2] = ud * (mat[1][0] * mat[2][1] - mat[1][1] * mat[2][0]);
+  inv[1][0] = -ud * (mat[0][1] * mat[2][2] - mat[0][2] * mat[2][1]);
+  inv[1][1] = ud * (mat[0][0] * mat[2][2] - mat[0][2] * mat[2][0]);
+  inv[1][2] = -ud * (mat[0][0] * mat[2][1] - mat[0][1] * mat[2][0]);
+  inv[2][0] = ud * (mat[0][1] * mat[1][2] - mat[0][2] * mat[1][1]);
+  inv[2][1] = -ud * (mat[0][0] * mat[1][2] - mat[0][2] * mat[1][0]);
+  inv[2][2] = ud * (mat[0][0] * mat[1][1] - mat[0][1] * mat[1][0]);
+  return 1;
+}
+
 static void eigen(List_T *inList, int inNb, 
 		  List_T *outList, int *outNb,
 		  int nbTime, int nbNod, int nbComp, int lam)
@@ -159,7 +182,7 @@ static void eigen(List_T *inList, int inNb,
 	  dx_dksi[0][0] = x[1] - x[0]; dx_dksi[0][1] = x[2]-x[0]; dx_dksi[0][2] = cross[0];
 	  dx_dksi[1][0] = y[1] - y[0]; dx_dksi[1][1] = y[2]-y[0]; dx_dksi[1][2] = cross[1];
 	  dx_dksi[2][0] = z[1] - z[0]; dx_dksi[2][1] = z[2]-z[0]; dx_dksi[2][2] = cross[2];
-	  inv3x3(dx_dksi, dksi_dx, &det);
+	  inv3x3tran(dx_dksi, dksi_dx, &det);
 	  GradPhi_ksi[0][0]= -1;  GradPhi_ksi[0][1]= -1;  GradPhi_ksi[0][2]= 0;  
 	  GradPhi_ksi[1][0]=  1;  GradPhi_ksi[1][1]=  0;  GradPhi_ksi[1][2]= 0;
 	  GradPhi_ksi[2][0]=  0;  GradPhi_ksi[2][1]=  1;  GradPhi_ksi[2][2]= 0;
@@ -168,7 +191,7 @@ static void eigen(List_T *inList, int inNb,
 	  dx_dksi[0][0] = x[1] - x[0]; dx_dksi[0][1] = x[2]-x[0]; dx_dksi[0][2] = x[3]-x[0];
 	  dx_dksi[1][0] = y[1] - y[0]; dx_dksi[1][1] = y[2]-y[0]; dx_dksi[1][2] = y[3]-y[0];
 	  dx_dksi[2][0] = z[1] - z[0]; dx_dksi[2][1] = z[2]-z[0]; dx_dksi[2][2] = z[3]-z[0];   
-	  inv3x3(dx_dksi, dksi_dx, &det);
+	  inv3x3tran(dx_dksi, dksi_dx, &det);
 	  GradPhi_ksi[0][0]= -1;  GradPhi_ksi[0][1]= -1; GradPhi_ksi[0][2]= -1; 
 	  GradPhi_ksi[1][0]=  1;  GradPhi_ksi[1][1]=  0; GradPhi_ksi[1][2]=  0; 
 	  GradPhi_ksi[2][0]=  0;  GradPhi_ksi[2][1]=  1; GradPhi_ksi[2][2]=  0; 
diff --git a/Plugin/Makefile b/Plugin/Makefile
index 752746ad1c6d3404539c5a18c0c7d68d23abc8f4..2b3a4804eb7af48e83cc7186f49e26c3061823d7 100644
--- a/Plugin/Makefile
+++ b/Plugin/Makefile
@@ -1,4 +1,4 @@
-# $Id: Makefile,v 1.78 2005-01-09 02:18:59 geuzaine Exp $
+# $Id: Makefile,v 1.79 2005-01-12 19:10:41 geuzaine Exp $
 #
 # Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 #
@@ -43,7 +43,7 @@ SRC = Plugin.cpp\
         Extract.cpp\
         DecomposeInSimplex.cpp\
         Evaluate.cpp\
-        Integrate.cpp\
+        Integrate.cpp Gradient.cpp Curl.cpp Divergence.cpp\
         Annotate.cpp Remove.cpp\
         Probe.cpp\
         HarmonicToTime.cpp
@@ -77,14 +77,15 @@ Plugin.o: Plugin.cpp Plugin.h ../Common/Options.h ../Common/Message.h \
   ../Common/GmshMatrix.h ../Common/AdaptiveViews.h PluginManager.h \
   CutMap.h Levelset.h CutGrid.h StreamLines.h CutPlane.h CutParametric.h \
   CutSphere.h Skin.h ../DataStr/Tree.h ../DataStr/avl.h Extract.h \
-  HarmonicToTime.h Integrate.h Annotate.h Remove.h DecomposeInSimplex.h \
-  Smooth.h Transform.h Triangulate.h SphericalRaise.h DisplacementRaise.h \
-  StructuralSolver.h ../Geo/Geo.h ../Mesh/Mesh.h ../Mesh/Vertex.h \
-  ../Mesh/Element.h ../Mesh/Simplex.h ../Mesh/Face.h ../Mesh/Edge.h \
-  ../Geo/ExtrudeParams.h ../Mesh/DiscreteSurface.h ../Mesh/Metric.h \
-  ../Mesh/Matrix.h ../Common/GmshUI.h Eigenvectors.h Eigenvalues.h \
-  Lambda2.h ../Numeric/Numeric.h Evaluate.h OctreePost.h Octree.h \
-  OctreeInternals.h Probe.h ../Common/Context.h
+  HarmonicToTime.h Integrate.h Gradient.h Curl.h Divergence.h Annotate.h \
+  Remove.h DecomposeInSimplex.h Smooth.h Transform.h Triangulate.h \
+  SphericalRaise.h DisplacementRaise.h StructuralSolver.h ../Geo/Geo.h \
+  ../Mesh/Mesh.h ../Mesh/Vertex.h ../Mesh/Element.h ../Mesh/Simplex.h \
+  ../Mesh/Face.h ../Mesh/Edge.h ../Geo/ExtrudeParams.h \
+  ../Mesh/DiscreteSurface.h ../Mesh/Metric.h ../Mesh/Matrix.h \
+  ../Common/GmshUI.h Eigenvectors.h Eigenvalues.h Lambda2.h \
+  ../Numeric/Numeric.h Evaluate.h OctreePost.h Octree.h OctreeInternals.h \
+  Probe.h ../Common/Context.h
 Levelset.o: Levelset.cpp Levelset.h Plugin.h ../Common/Options.h \
   ../Common/Message.h ../Common/Views.h ../Common/ColorTable.h \
   ../DataStr/List.h ../Common/VertexArray.h ../Common/SmoothNormals.h \
@@ -228,6 +229,21 @@ Integrate.o: Integrate.cpp Plugin.h ../Common/Options.h \
   ../DataStr/List.h ../Common/VertexArray.h ../Common/SmoothNormals.h \
   ../Common/GmshMatrix.h ../Common/AdaptiveViews.h Integrate.h \
   ../Common/Context.h ../Numeric/Numeric.h ShapeFunctions.h
+Gradient.o: Gradient.cpp Plugin.h ../Common/Options.h ../Common/Message.h \
+  ../Common/Views.h ../Common/ColorTable.h ../DataStr/List.h \
+  ../Common/VertexArray.h ../Common/SmoothNormals.h \
+  ../Common/GmshMatrix.h ../Common/AdaptiveViews.h Gradient.h \
+  ../Common/Context.h ../Numeric/Numeric.h ShapeFunctions.h
+Curl.o: Curl.cpp Plugin.h ../Common/Options.h ../Common/Message.h \
+  ../Common/Views.h ../Common/ColorTable.h ../DataStr/List.h \
+  ../Common/VertexArray.h ../Common/SmoothNormals.h \
+  ../Common/GmshMatrix.h ../Common/AdaptiveViews.h Curl.h \
+  ../Common/Context.h ../Numeric/Numeric.h ShapeFunctions.h
+Divergence.o: Divergence.cpp Plugin.h ../Common/Options.h \
+  ../Common/Message.h ../Common/Views.h ../Common/ColorTable.h \
+  ../DataStr/List.h ../Common/VertexArray.h ../Common/SmoothNormals.h \
+  ../Common/GmshMatrix.h ../Common/AdaptiveViews.h Divergence.h \
+  ../Common/Context.h ../Numeric/Numeric.h ShapeFunctions.h
 Annotate.o: Annotate.cpp Plugin.h ../Common/Options.h ../Common/Message.h \
   ../Common/Views.h ../Common/ColorTable.h ../DataStr/List.h \
   ../Common/VertexArray.h ../Common/SmoothNormals.h \
diff --git a/Plugin/Plugin.cpp b/Plugin/Plugin.cpp
index 05d62797158870b1a3be62b110461a747092a7a4..5be5e88abdf77a3f3db7c398883c81df07e4beab 100644
--- a/Plugin/Plugin.cpp
+++ b/Plugin/Plugin.cpp
@@ -1,4 +1,4 @@
-// $Id: Plugin.cpp,v 1.72 2005-01-09 02:19:00 geuzaine Exp $
+// $Id: Plugin.cpp,v 1.73 2005-01-12 19:10:41 geuzaine Exp $
 //
 // Copyright (C) 1997-2005 C. Geuzaine, J.-F. Remacle
 //
@@ -42,6 +42,9 @@
 #include "Extract.h"
 #include "HarmonicToTime.h"
 #include "Integrate.h"
+#include "Gradient.h"
+#include "Curl.h"
+#include "Divergence.h"
 #include "Annotate.h"
 #include "Remove.h"
 #include "DecomposeInSimplex.h"
@@ -200,6 +203,12 @@ void GMSH_PluginManager::registerDefaultPlugins()
 		      ("HarmonicToTime", GMSH_RegisterHarmonicToTimePlugin()));
     allPlugins.insert(std::pair < char *, GMSH_Plugin * >
 		      ("Integrate", GMSH_RegisterIntegratePlugin()));
+    allPlugins.insert(std::pair < char *, GMSH_Plugin * >
+		      ("Gradient", GMSH_RegisterGradientPlugin()));
+    allPlugins.insert(std::pair < char *, GMSH_Plugin * >
+		      ("Curl", GMSH_RegisterCurlPlugin()));
+    allPlugins.insert(std::pair < char *, GMSH_Plugin * >
+		      ("Divergence", GMSH_RegisterDivergencePlugin()));
     allPlugins.insert(std::pair < char *, GMSH_Plugin * >
 		      ("Annotate", GMSH_RegisterAnnotatePlugin()));
     allPlugins.insert(std::pair < char *, GMSH_Plugin * >
diff --git a/Plugin/ShapeFunctions.h b/Plugin/ShapeFunctions.h
index 66258a86eae04355664e2b6f399a1bcefbc987f4..3b4585b72ad443ce62f64e9740c9e0545d8982fd 100644
--- a/Plugin/ShapeFunctions.h
+++ b/Plugin/ShapeFunctions.h
@@ -59,6 +59,13 @@ public:
 	jac[0][0] += _x[i] * s[0]; jac[0][1] += _y[i] * s[0]; jac[0][2] += _z[i] * s[0];
 	jac[1][0] += _x[i] * s[1]; jac[1][1] += _y[i] * s[1]; jac[1][2] += _z[i] * s[1];
       }
+      {
+	double a[3], b[3], c[3];
+	a[0]= _x[1] - _x[0]; a[1]= _y[1] - _y[0]; a[2]= _z[1] - _z[0];
+	b[0]= _x[2] - _x[0]; b[1]= _y[2] - _y[0]; b[2]= _z[2] - _z[0];
+	prodve(a, b, c);
+	jac[2][0] = c[0]; jac[2][1] = c[1]; jac[2][2] = c[2]; 
+      }
       return sqrt(DSQR(jac[0][0] * jac[1][1] - jac[0][1] * jac[1][0]) +
 		  DSQR(jac[0][2] * jac[1][0] - jac[0][0] * jac[1][2]) +
 		  DSQR(jac[0][1] * jac[1][2] - jac[0][2] * jac[1][1]));
@@ -67,19 +74,25 @@ public:
 	getGradShapeFunction(i, u, v, w, s);
 	jac[0][0] += _x[i] * s[0]; jac[0][1] += _y[i] * s[0]; jac[0][2] += _z[i] * s[0];
       }
+      {
+	double a[3], b[3], c[3];
+	a[0]= _x[1] - _x[0]; a[1]= _y[1] - _y[0]; a[2]= _z[1] - _z[0];	
+	if((fabs(a[0]) >= fabs(a[1]) && fabs(a[0]) >= fabs(a[2])) ||
+	   (fabs(a[1]) >= fabs(a[0]) && fabs(a[1]) >= fabs(a[2]))) {
+	  b[0] = a[1]; b[1] = -a[0]; b[2] = 0.;
+	}
+	else {
+	  b[0] = 0.; b[1] = a[2]; b[2] = -a[1];
+	}
+	prodve(a, b, c);
+	jac[1][0] = b[0]; jac[1][1] = b[1]; jac[1][2] = b[2]; 
+	jac[2][0] = c[0]; jac[2][1] = c[1]; jac[2][2] = c[2]; 
+      }
       return sqrt(DSQR(jac[0][0])+DSQR(jac[0][1])+DSQR(jac[0][2]));
     default:
       return 1.;
     }
   }
-  double getInverseJacobian(double jac[3][3], double invjac[3][3])
-  {
-    // maybe we should use the classic approach to define always
-    // non-singular jacobians, so that we can simply invert them with
-    // inv3x3
-    Msg(GERROR, "getInverseJacibian not done yet");
-    return 0.;
-  }
   double interpolate(double val[], double u, double v, double w, int stride=1)
   {
     double sum = 0;
@@ -87,29 +100,55 @@ public:
     for(int i = 0; i < getNumNodes(); i++){
       double s;
       getShapeFunction(i, u, v, w, s);
-      sum += val[i] * s;
+      sum += val[j] * s;
       j += stride;
     }
     return sum;
   }
-  void interpolateGrad(double val[], double u, double v, double w, double f[3], int stride=1)
+  void interpolateGrad(double val[], double u, double v, double w, double f[3], int stride=1,
+		       double invjac[3][3]=NULL)
   {
-    double fu[3] = {0., 0., 0.};
+    double dfdu[3] = {0., 0., 0.};
     int j = 0;
     for(int i = 0; i < getNumNodes(); i++){
       double s[3];
       getGradShapeFunction(i, u, v, w, s);
-      fu[0] += val[j] * s[0];
-      fu[1] += val[j] * s[1];
-      fu[2] += val[j] * s[2];
+      dfdu[0] += val[j] * s[0];
+      dfdu[1] += val[j] * s[1];
+      dfdu[2] += val[j] * s[2];
       j += stride;
     }
-    double jac[3][3], invjac[3][3];
+    if(invjac){
+      matvec(invjac, dfdu, f);
+    }
+    else{
+      double jac[3][3], inv[3][3];
+      getJacobian(u, v, w, jac);
+      inv3x3(jac, inv);
+      matvec(inv, dfdu, f);
+    }
+  }
+  void interpolateCurl(double val[], double u, double v, double w, double f[3], int stride=3)
+  {
+    double fx[3], fy[3], fz[3], jac[3][3], inv[3][3];
+    getJacobian(u, v, w, jac);
+    inv3x3(jac, inv);
+    interpolateGrad(&val[0], u, v, w, fx, stride, inv);
+    interpolateGrad(&val[1], u, v, w, fy, stride, inv);
+    interpolateGrad(&val[2], u, v, w, fz, stride, inv);
+    f[0] = fz[1] - fy[2];
+    f[1] = -(fz[0] - fx[2]);
+    f[2] = fy[0] - fx[1];
+  }
+  double interpolateDiv(double val[], double u, double v, double w, int stride=3)
+  {
+    double fx[3], fy[3], fz[3], jac[3][3], inv[3][3];
     getJacobian(u, v, w, jac);
-    getInverseJacobian(jac, invjac);
-    f[0] = fu[0] * invjac[0][0] + fu[1] * invjac[0][1] + fu[2] * invjac[0][2];
-    f[1] = fu[0] * invjac[1][0] + fu[1] * invjac[1][1] + fu[2] * invjac[1][2];
-    f[2] = fu[0] * invjac[2][0] + fu[1] * invjac[2][1] + fu[2] * invjac[2][2];
+    inv3x3(jac, inv);
+    interpolateGrad(&val[0], u, v, w, fx, stride, inv);
+    interpolateGrad(&val[1], u, v, w, fy, stride, inv);
+    interpolateGrad(&val[2], u, v, w, fz, stride, inv);
+    return fx[0] + fy[1] + fz[2];
   }
   double integrate(double val[], int stride=1)
   {
@@ -219,12 +258,8 @@ public:
     double t[3] = {_x[1]-_x[0], _y[1]-_y[0], _z[1]-_z[0]};
     norme(t);
     double v[3];
-    for(int i = 0; i < 3; i++){
-      double tmp[2]; 
-      tmp[0] = val[i];
-      tmp[1] = val[3+i];
-      v[i] = integrate(tmp);
-    }
+    for(int i = 0; i < 3; i++)
+      v[i] = integrate(&val[i], 3);
     double d;
     prosca(t, v, &d);
     return d;
@@ -284,13 +319,8 @@ public:
     prodve(t1, t2, n);
     norme(n);
     double v[3];
-    for(int i = 0; i < 3; i++){
-      double tmp[3]; 
-      tmp[0] = val[i];
-      tmp[1] = val[3+i];
-      tmp[2] = val[6+i];
-      v[i] = integrate(tmp);
-    }
+    for(int i = 0; i < 3; i++)
+      v[i] = integrate(&val[i], 3);
     double d;
     prosca(n, v, &d);
     return d;
@@ -353,14 +383,8 @@ public:
     prodve(t1, t2, n);
     norme(n);
     double v[3];
-    for(int i = 0; i < 3; i++){
-      double tmp[4]; 
-      tmp[0] = val[i];
-      tmp[1] = val[3+i];
-      tmp[2] = val[6+i];
-      tmp[3] = val[9+i];
-      v[i] = integrate(tmp);
-    }
+    for(int i = 0; i < 3; i++)
+      v[i] = integrate(&val[i], 3);
     double d;
     prosca(n, v, &d);
     return d;
@@ -670,7 +694,9 @@ class elementFactory{
       else if(numNodes == 6) return new prism(x, y, z);
       else if(numNodes == 5) return new pyramid(x, y, z);
       else return new tetrahedron(x, y, z);
-    default: return NULL;
+    default: 
+      Msg(GERROR, "Unknown type of element in factory");
+      return NULL;
     }
   }
 };
diff --git a/TODO b/TODO
index 016f94d6668f7e2e344be2c2ecb0191ea9be879b..51db9db1bf6b8007279642a873e22cddb3914800 100644
--- a/TODO
+++ b/TODO
@@ -1,14 +1,19 @@
-$Id: TODO,v 1.75 2005-01-09 02:18:59 geuzaine Exp $
+$Id: TODO,v 1.76 2005-01-12 19:10:40 geuzaine Exp $
 
-create Gradient, Curl and Divergence plugins (based on
-Plugins/ShapeFunctions.h)
+********************************************************************
+
+change the Tranform plugin so that it takes expression arguments?
 
 ********************************************************************
 
-Add dynamic variables in labels? E.g., if a string contains
-%string.string, replace it dynamically in Draw_Text2D/3D by calling
+Labels:
+- add dynamic variables? E.g., if a string contains %string.string,
+replace it dynamically in Draw_Text2D/3D by calling
 Get_OptionsCategory/Get_OptionValue. (provide a global option to
 disable all replacements)
+- position with the mouse
+- select/move once positioned?
+- add ability to add arrows
 
 ********************************************************************
 
diff --git a/doc/VERSIONS b/doc/VERSIONS
index 7d894bccbbf76b0febd78ef8c469a93714fa971f..11fc743347e42e4040903bbc82e32607fb58f01c 100644
--- a/doc/VERSIONS
+++ b/doc/VERSIONS
@@ -1,4 +1,4 @@
-$Id: VERSIONS,v 1.299 2005-01-10 04:24:55 geuzaine Exp $
+$Id: VERSIONS,v 1.300 2005-01-12 19:10:41 geuzaine Exp $
 
 New since 1.58: added support for discrete (triangulated) surfaces,
 either in STL format or with the new "Discrete Surface" command; added
@@ -6,8 +6,9 @@ STL and Text output format for post-processing views and STL output
 format for surface meshes; all levelset-based plugins can now also
 compute isovolumes; generalized Plugin(Evaluate) to handle external
 view data (based on the same or on a different mesh); generalized
-Plugin(CutGrid); new plugin (Eigenvalues); changed default colormap to
-match Matlab's "Jet" colormap; fixed small bugs.
+Plugin(CutGrid); new plugins (Eigenvalues, Gradient, Curl,
+Divergence); changed default colormap to match Matlab's "Jet"
+colormap; fixed small bugs.
 
 New in 1.58: fixed UNIX socket interface on Windows (broken by the TCP
 solver patch in 1.57); bumped version number of default
diff --git a/doc/texinfo/opt_plugin.texi b/doc/texinfo/opt_plugin.texi
index 1e0c8585c528e46cd2159a0e0e4b2e57f32e1097..1d42ac0cc2a3c9b22911e99ede34c9d3c9019d8e 100644
--- a/doc/texinfo/opt_plugin.texi
+++ b/doc/texinfo/opt_plugin.texi
@@ -36,6 +36,19 @@ Default value: @code{14}
 Default value: @code{-1}
 @end table
 
+@item Plugin(Curl)
+Plugin(Curl) computes the curl of the field
+in the view `iView'. If `iView' < 0, the plugin
+is run on the current view.
+
+Plugin(Curl) creates one new view.
+
+Numeric options:
+@table @code
+@item iView
+Default value: @code{-1}
+@end table
+
 @item Plugin(CutGrid)
 Plugin(CutGrid) cuts a triangle/tetrahedron view
 with a rectangular grid defined by the 3 points
@@ -249,6 +262,19 @@ Default value: @code{-1}
 Default value: @code{-1}
 @end table
 
+@item Plugin(Divergence)
+Plugin(Divergence) computes the divergence of the
+field in the view `iView'. If `iView' < 0, the plugin
+is run on the current view.
+
+Plugin(Divergence) creates one new view.
+
+Numeric options:
+@table @code
+@item iView
+Default value: @code{-1}
+@end table
+
 @item Plugin(Eigenvalues)
 Plugin(Eigenvalues) computes the three real
 eigenvalues of each tensor in the view `iView'.
@@ -381,6 +407,19 @@ Numeric options:
 Default value: @code{-1}
 @end table
 
+@item Plugin(Gradient)
+Plugin(Gradient) computes the gradient of the
+field in the view `iView'. If `iView' < 0, the
+plugin is run on the current view.
+
+Plugin(Gradient) creates one new view.
+
+Numeric options:
+@table @code
+@item iView
+Default value: @code{-1}
+@end table
+
 @item Plugin(HarmonicToTime)
 Plugin(HarmonicToTime) takes the values in the
 time steps `RealPart' and `ImaginaryPart' of