diff --git a/Graphics/IsoSimplex.cpp b/Graphics/IsoSimplex.cpp
new file mode 100644
index 0000000000000000000000000000000000000000..bef930f9c95f28335435659c8b57b50e7f584233
--- /dev/null
+++ b/Graphics/IsoSimplex.cpp
@@ -0,0 +1,186 @@
+#include "Gmsh.h"
+#include "GmshUI.h"
+#include "Geo.h"
+#include "Mesh.h"
+#include "Draw.h"
+#include "Iso.h"
+#include "Context.h"
+#include "Numeric.h"
+#include "Views.h"
+
+extern Context_T CTX;
+
+/* ------------------------------------------------------------------------ */
+/* S i m p l e x */
+/* ------------------------------------------------------------------------ */
+
+/*
+ compute the gradient of a linear interpolation in a tetrahedron
+*/
+void gradSimplex (double *x, double *y, double *z, double *v, double *grad)
+{
+ /*
+ p = p1 * (1-u-v-w) + p2 u + p3 v + p4 w
+ */
+
+ double mat[3][3];
+ double det,b[3];
+ mat[0][0] = x[1]-x[0];
+ mat[1][0] = x[2]-x[0];
+ mat[2][0] = x[3]-x[0];
+ mat[0][1] = y[1]-y[0];
+ mat[1][1] = y[2]-y[0];
+ mat[2][1] = y[3]-y[0];
+ mat[0][2] = z[1]-z[0];
+ mat[1][2] = z[2]-z[0];
+ mat[2][2] = z[3]-z[0];
+ b[0] = v[1]-v[0];
+ b[1] = v[2]-v[0];
+ b[2] = v[3]-v[0];
+ sys3x3 (mat, b, grad, &det);
+}
+
+void IsoSimplex( Post_View *View,
+ int preproNormals,
+ double *X, double *Y, double *Z, double *Val,
+ double V, double Vmin, double Vmax,
+ double *Offset, double Raise[3][5], int shade){
+ int nb,i;
+ double Xp[6],Yp[6],Zp[6];
+ double Xpi[6],Ypi[6],Zpi[6];
+ double norms[12];
+
+ if(V != Vmax){
+ nb = 0;
+ if((Val[0] > V && Val[1] <= V) || (Val[1] > V && Val[0] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,0,1,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[0] > V && Val[2] <= V) || (Val[2] > V && Val[0] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,0,2,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[0] > V && Val[3] <= V) || (Val[3] > V && Val[0] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,0,3,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[1] > V && Val[2] <= V) || (Val[2] > V && Val[1] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,1,2,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[1] > V && Val[3] <= V) || (Val[3] > V && Val[1] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,1,3,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[2] > V && Val[3] <= V) || (Val[3] > V && Val[2] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,2,3,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ }
+ else{
+ nb=0;
+ if((Val[0] < V && Val[1] <= V) || (Val[1] < V && Val[0] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,0,1,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[0] < V && Val[2] <= V) || (Val[2] < V && Val[0] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,0,2,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[0] < V && Val[3] <= V) || (Val[3] < V && Val[0] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,0,3,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[1] < V && Val[2] <= V) || (Val[2] < V && Val[1] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,1,2,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[1] < V && Val[3] <= V) || (Val[3] < V && Val[1] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,1,3,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ if((Val[2] < V && Val[3] <= V) || (Val[3] < V && Val[2] <= V)){
+ InterpolateIso(X,Y,Z,Val,V,2,3,&Xp[nb],&Yp[nb],&Zp[nb]); nb++;
+ }
+ }
+
+ /*
+ 3 possibilities for quads
+ -) 0,2,5,3
+ -) 0,1,5,4
+ -) 1,2,4,3
+ in all cases, simply invert the 2 last ones
+ for having the quads ordered
+ */
+
+ if(nb == 4)
+ {
+ double xx = Xp[3];
+ double yy = Yp[3];
+ double zz = Zp[3];
+ Xp[3] = Xp[2];
+ Yp[3] = Yp[2];
+ Zp[3] = Zp[2];
+ Xp[2] = xx;
+ Yp[2] = yy;
+ Zp[2] = zz;
+ }
+
+ /*
+ for having a nice isosurface, we should have n . grad v > 0
+ n = normal to the polygon
+ v = unknown field we wanna draw
+ */
+
+ if(nb > 2)
+ {
+ double v1[3] = {Xp[2]-Xp[0],Yp[2]-Yp[0],Zp[2]-Zp[0]};
+ double v2[3] = {Xp[1]-Xp[0],Yp[1]-Yp[0],Zp[1]-Zp[0]};
+ double gr[3];
+ double n[3],xx;
+ prodve(v1,v2,n);
+ norme(n);
+ gradSimplex(X,Y,Z,Val,gr);
+ prosca(gr,n,&xx);
+
+ if(xx > 0)
+ {
+ for(i=0;i<nb;i++)
+ {
+ Xpi[i] = Xp[i];
+ Ypi[i] = Yp[i];
+ Zpi[i] = Zp[i];
+ }
+ for(i=0;i<nb;i++)
+ {
+ Xp[i] = Xpi[nb-i-1];
+ Yp[i] = Ypi[nb-i-1];
+ Zp[i] = Zpi[nb-i-1];
+ }
+ }
+ else
+ {
+ n[0] = -n[0];
+ n[1] = -n[1];
+ n[2] = -n[2];
+ }
+ if(preproNormals)
+ {
+ for(i=0;i<nb;i++)
+ {
+ View->add_normal(Xp[i],Yp[i],Zp[i],n[0],n[1],n[2]);
+ }
+ return;
+ }
+ else
+ {
+ for(i=0;i<nb;i++)
+ {
+ if(!View->get_normal(Xp[i],Yp[i],Zp[i],norms[3*i],norms[3*i+1],norms[3*i+2]))
+ {
+ //printf("coucou\n");
+ norms[3*i] = n[0];
+ norms[3*i+1] = n[1];
+ norms[3*i+2] = n[2];
+ }
+ }
+ }
+ }
+
+
+ if(nb == 3)
+ Draw_Triangle(Xp,Yp,Zp,norms,Offset,Raise,shade);
+ else if(nb == 4)
+ Draw_Quadrangle(Xp,Yp,Zp,norms,Offset,Raise,shade);
+
+}
+
diff --git a/Graphics/IsoSimplex.h b/Graphics/IsoSimplex.h
new file mode 100644
index 0000000000000000000000000000000000000000..bc76112bc9bb5369c35324add006f7193c23ec1e
--- /dev/null
+++ b/Graphics/IsoSimplex.h
@@ -0,0 +1,10 @@
+#ifndef _ISOSIMPLEX_H_
+#define _ISOSIMPLEX_H_
+
+struct Post_View;
+void IsoSimplex (Post_View *View,
+ int preproNormals,
+ double *X, double *Y, double *Z, double *Val,
+ double V, double Vmin, double Vmax,
+ double *Offset, double Raise[3][5], int shade);
+#endif