Index: ../trunk-jpl/src/c/classes/Elements/PentaRef.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/Elements/PentaRef.cpp	(revision 18522)
+++ ../trunk-jpl/src/c/classes/Elements/PentaRef.cpp	(revision 18523)
@@ -1,4 +1,4 @@
-/*!\file PentaRef.c
+/*!\file PentaRef.cpp
  * \brief: implementation of the PentaRef object
  */
 
@@ -25,6 +25,7 @@
 #define NUMNODESP2    18
 #define NUMNODESP2b   19
 #define NUMNODESP2xP4 30
+#define NUMNODESMAX   30
 
 /*Object constructors and destructor*/
 PentaRef::PentaRef(){/*{{{*/
@@ -885,24 +886,25 @@
 }
 /*}}}*/
 void PentaRef::GetInputValue(IssmDouble* pvalue,IssmDouble* plist,Gauss* gauss,int finiteelement){/*{{{*/
+	/* WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.*/
 
-	/*Output*/
-	IssmDouble value =0.;
+	/*Allocate basis functions*/
+	IssmDouble  basis[NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions*/
-	IssmDouble* basis=xNew<IssmDouble>(numnodes);
-	GetNodalFunctions(basis, gauss,finiteelement);
+	/*Get basis functions at this point*/
+	GetNodalFunctions(&basis[0],gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble value =0.;
 	for(int i=0;i<numnodes;i++) value += basis[i]*plist[i];
 
 	/*Assign output pointer*/
-	xDelete<IssmDouble>(basis);
 	*pvalue = value;
-
 }
 /*}}}*/
 void PentaRef::GetInputDerivativeValue(IssmDouble* p, IssmDouble* plist,IssmDouble* xyz_list, Gauss* gauss,int finiteelement){/*{{{*/
@@ -914,31 +916,33 @@
 	 *   dp/dz=p_list[0]*dh1/dz+p_list[1]*dh2/dz+p_list[2]*dh3/dz+p_list[3]*dh4/dz+p_list[4]*dh5/dz+p_list[5]*dh6/dz;
 	 *
 	 *   p is a vector of size 3x1 already allocated.
+	 *
+	 * WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.
 	 */
 
-	/*Output*/
-	IssmDouble dpx=0.;
-	IssmDouble dpy=0.;
-	IssmDouble dpz=0.;
+	/*Allocate derivatives of basis functions*/
+	IssmDouble  dbasis[3*NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions derivatives*/
-	IssmDouble* dbasis=xNew<IssmDouble>(3*numnodes);
-	GetNodalFunctionsDerivatives(dbasis,xyz_list,gauss,finiteelement);
+	/*Get basis functions derivatives at this point*/
+	GetNodalFunctionsDerivatives(&dbasis[0],xyz_list,gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble dpx=0.;
+	IssmDouble dpy=0.;
+	IssmDouble dpz=0.;
 	for(int i=0;i<numnodes;i++) dpx += dbasis[0*numnodes+i]*plist[i];
 	for(int i=0;i<numnodes;i++) dpy += dbasis[1*numnodes+i]*plist[i];
 	for(int i=0;i<numnodes;i++) dpz += dbasis[2*numnodes+i]*plist[i];
 
 	/*Assign values*/
-	xDelete<IssmDouble>(dbasis);
 	p[0]=dpx;
 	p[1]=dpy;
 	p[2]=dpz;
-
 }
 /*}}}*/
 int  PentaRef::NumberofNodes(int finiteelement){/*{{{*/
Index: ../trunk-jpl/src/c/classes/Elements/SegRef.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/Elements/SegRef.cpp	(revision 18522)
+++ ../trunk-jpl/src/c/classes/Elements/SegRef.cpp	(revision 18523)
@@ -17,6 +17,7 @@
 /*Element macros*/
 #define NUMNODESP0  1
 #define NUMNODESP1  2
+#define NUMNODESMAX 2
 
 /*Object constructors and destructor*/
 SegRef::SegRef(){/*{{{*/
@@ -114,44 +115,48 @@
 	 *   dp/dx=plist[0]*dh1/dx+plist[1]*dh2/dx
 	 *
 	 * p is a vector already allocated.
+	 *
+	 * WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.
 	 */
 
-	/*Output*/
-	IssmDouble dpx = 0.;
+	/*Allocate derivatives of basis functions*/
+	IssmDouble  dbasis[1*NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions derivatives*/
-	IssmDouble* dbasis=xNew<IssmDouble>(1*numnodes);
-	GetNodalFunctionsDerivatives(dbasis,xyz_list,gauss,finiteelement);
+	/*Get basis functions derivatives at this point*/
+	GetNodalFunctionsDerivatives(&dbasis[0],xyz_list,gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
-	for(int i=0;i<numnodes;i++) dpx += dbasis[i]*plist[i];
+	IssmDouble dpx=0.;
+	for(int i=0;i<numnodes;i++) dpx += dbasis[0*numnodes+i]*plist[i];
 
 	/*Assign values*/
-	xDelete<IssmDouble>(dbasis);
-	*p=dpx;
-
+	p[0]=dpx;
 }
 /*}}}*/
 void SegRef::GetInputValue(IssmDouble* p, IssmDouble* plist, GaussSeg* gauss,int finiteelement){/*{{{*/
+	/* WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.*/
 
-	/*Output*/
-	IssmDouble value =0.;
+	/*Allocate basis functions*/
+	IssmDouble  basis[NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions*/
-	IssmDouble* basis=xNew<IssmDouble>(numnodes);
-	GetNodalFunctions(basis, gauss,finiteelement);
+	/*Get basis functions at this point*/
+	GetNodalFunctions(&basis[0],gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble value =0.;
 	for(int i=0;i<numnodes;i++) value += basis[i]*plist[i];
 
 	/*Assign output pointer*/
-	xDelete<IssmDouble>(basis);
 	*p = value;
 }
 /*}}}*/
Index: ../trunk-jpl/src/c/classes/Elements/TetraRef.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/Elements/TetraRef.cpp	(revision 18522)
+++ ../trunk-jpl/src/c/classes/Elements/TetraRef.cpp	(revision 18523)
@@ -19,6 +19,7 @@
 #define NUMNODESP1  4
 #define NUMNODESP1b 5
 #define NUMNODESP2  10
+#define NUMNODESMAX 10
 
 /*Object constructors and destructor*/
 TetraRef::TetraRef(){/*{{{*/
@@ -216,49 +217,54 @@
 	 *   dp/dz=p_list[0]*dh1/dz+p_list[1]*dh2/dz+p_list[2]*dh3/dz+p_list[3]*dh4/dz+p_list[4]*dh5/dz+p_list[5]*dh6/dz;
 	 *
 	 *   p is a vector of size 3x1 already allocated.
+	 *
+	 * WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.
 	 */
 
-	/*Output*/
-	IssmDouble dpx=0.;
-	IssmDouble dpy=0.;
-	IssmDouble dpz=0.;
+	/*Allocate derivatives of basis functions*/
+	IssmDouble  dbasis[3*NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions derivatives*/
-	IssmDouble* dbasis=xNew<IssmDouble>(3*numnodes);
-	GetNodalFunctionsDerivatives(dbasis,xyz_list,gauss,finiteelement);
+	/*Get basis functions derivatives at this point*/
+	GetNodalFunctionsDerivatives(&dbasis[0],xyz_list,gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble dpx=0.;
+	IssmDouble dpy=0.;
+	IssmDouble dpz=0.;
 	for(int i=0;i<numnodes;i++) dpx += dbasis[0*numnodes+i]*plist[i];
 	for(int i=0;i<numnodes;i++) dpy += dbasis[1*numnodes+i]*plist[i];
 	for(int i=0;i<numnodes;i++) dpz += dbasis[2*numnodes+i]*plist[i];
 
 	/*Assign values*/
-	xDelete<IssmDouble>(dbasis);
 	p[0]=dpx;
 	p[1]=dpy;
 	p[2]=dpz;
 }
 /*}}}*/
 void TetraRef::GetInputValue(IssmDouble* p, IssmDouble* plist, Gauss* gauss,int finiteelement){/*{{{*/
+	/* WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.*/
 
-	/*Output*/
-	IssmDouble value =0.;
+	/*Allocate basis functions*/
+	IssmDouble  basis[NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions*/
-	IssmDouble* basis=xNew<IssmDouble>(numnodes);
-	GetNodalFunctions(basis, gauss,finiteelement);
+	/*Get basis functions at this point*/
+	GetNodalFunctions(&basis[0],gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble value =0.;
 	for(int i=0;i<numnodes;i++) value += basis[i]*plist[i];
 
 	/*Assign output pointer*/
-	xDelete<IssmDouble>(basis);
 	*p = value;
 }
 /*}}}*/
Index: ../trunk-jpl/src/c/classes/Elements/TriaRef.cpp
===================================================================
--- ../trunk-jpl/src/c/classes/Elements/TriaRef.cpp	(revision 18522)
+++ ../trunk-jpl/src/c/classes/Elements/TriaRef.cpp	(revision 18523)
@@ -20,6 +20,7 @@
 #define NUMNODESP1b 4
 #define NUMNODESP2  6
 #define NUMNODESP2b 7
+#define NUMNODESMAX 7
 
 /*Object constructors and destructor*/
 TriaRef::TriaRef(){/*{{{*/
@@ -353,54 +354,58 @@
 }
 /*}}}*/
 void TriaRef::GetInputDerivativeValue(IssmDouble* p, IssmDouble* plist,IssmDouble* xyz_list, Gauss* gauss,int finiteelement){/*{{{*/
-
 	/*From node values of parameter p (plist[0],plist[1],plist[2]), return parameter derivative value at gaussian 
 	 * point specified by gauss_basis:
 	 *   dp/dx=plist[0]*dh1/dx+plist[1]*dh2/dx+plist[2]*dh3/dx
 	 *   dp/dx=plist[0]*dh1/dx+plist[1]*dh2/dx+plist[2]*dh3/dx
 	 *
 	 * p is a vector already allocated.
+	 *
+	 * WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.
 	 */
 
-	/*Output*/
-	IssmDouble dpx=0.;
-	IssmDouble dpy=0.;
+	/*Allocate derivatives of basis functions*/
+	IssmDouble  dbasis[2*NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions derivatives*/
-	IssmDouble* dbasis=xNew<IssmDouble>(2*numnodes);
-	GetNodalFunctionsDerivatives(dbasis,xyz_list,gauss,finiteelement);
+	/*Get basis functions derivatives at this point*/
+	GetNodalFunctionsDerivatives(&dbasis[0],xyz_list,gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble dpx=0.;
+	IssmDouble dpy=0.;
 	for(int i=0;i<numnodes;i++) dpx += dbasis[0*numnodes+i]*plist[i];
 	for(int i=0;i<numnodes;i++) dpy += dbasis[1*numnodes+i]*plist[i];
 
 	/*Assign values*/
-	xDelete<IssmDouble>(dbasis);
 	p[0]=dpx;
 	p[1]=dpy;
 
 }
 /*}}}*/
 void TriaRef::GetInputValue(IssmDouble* p, IssmDouble* plist, Gauss* gauss,int finiteelement){/*{{{*/
+	/* WARNING: For a significant gain in performance, it is better to use
+	 * static memory allocation instead of dynamic.*/
 
-	/*Output*/
-	IssmDouble value =0.;
+	/*Allocate basis functions*/
+	IssmDouble  basis[NUMNODESMAX];
 
 	/*Fetch number of nodes for this finite element*/
 	int numnodes = this->NumberofNodes(finiteelement);
+	_assert_(numnodes<=NUMNODESMAX);
 
-	/*Get nodal functions*/
-	IssmDouble* basis=xNew<IssmDouble>(numnodes);
-	GetNodalFunctions(basis, gauss,finiteelement);
+	/*Get basis functions at this point*/
+	GetNodalFunctions(&basis[0],gauss,finiteelement);
 
 	/*Calculate parameter for this Gauss point*/
+	IssmDouble value =0.;
 	for(int i=0;i<numnodes;i++) value += basis[i]*plist[i];
 
 	/*Assign output pointer*/
-	xDelete<IssmDouble>(basis);
 	*p = value;
 }
 /*}}}*/
