Index: /issm/trunk-jpl/test/NightlyRun/test1205.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test1205.py	(revision 14157)
+++ /issm/trunk-jpl/test/NightlyRun/test1205.py	(revision 14157)
@@ -0,0 +1,109 @@
+import numpy
+from model import *
+from roundmesh import *
+from setmask import *
+from parameterize import *
+from setflowequation import *
+from EnumDefinitions import *
+from solve import *
+from MatlabFuncs import *
+
+"""
+The aim of this program is to compare a model with an analytical solution given in MacAyeal EISMINT : Lessons in Ice-Sheet Modeling.
+"""
+
+printingflag=False
+
+numlayers=10
+resolution=30000.
+
+#To begin with the numerical model
+md=model()
+md=roundmesh(md,750000.,resolution)
+md=setmask(md,'','')    #We can not test iceshelves nor ice rises with this analytical solution
+md=parameterize(md,'../Par/RoundSheetStaticEISMINT.py')
+
+#Calculation of the analytical 2d velocity field
+constant=0.3
+vx_obs=constant/2.*md.mesh.x*(md.geometry.thickness)**-1
+vy_obs=constant/2.*md.mesh.y*(md.geometry.thickness)**-1
+vel_obs=numpy.sqrt((md.inversion.vx_obs)**2+(md.inversion.vy_obs)**2)
+
+#We extrude the model to have a 3d model
+md.extrude(numlayers,1.);
+md=setflowequation(md,'hutter','all')
+
+#Spc the nodes on the bed
+pos=numpy.nonzero(md.mesh.vertexonbed)
+md.diagnostic.spcvx[pos]=0.
+md.diagnostic.spcvy[pos]=0.
+md.diagnostic.spcvz[pos]=0.
+
+#Now we can solve the problem 
+md.cluster=generic('name',oshostname(),'np',8)
+md=solve(md,DiagnosticSolutionEnum())
+
+#Calculate the depth averaged velocity field (2d):
+vx=md.results.DiagnosticSolution.Vx
+vy=md.results.DiagnosticSolution.Vy
+vel=numpy.zeros((md.mesh.numberofvertices2d,1))
+
+for i in xrange(0,md.mesh.numberofvertices2d):
+	node_vel=0.
+	for j in xrange(0,md.mesh.numberoflayers-1):
+		node_vel=node_vel+1./(2.*(md.mesh.numberoflayers-1))*\
+			(numpy.sqrt(vx[i+(j+1)*md.mesh.numberofvertices2d,0]**2+vy[i+(j+1)*md.mesh.numberofvertices2d,0]**2)+\
+			numpy.sqrt(vx[i+j*md.mesh.numberofvertices2d,0]**2+vy[i+j*md.mesh.numberofvertices2d,0]**2))
+	vel[i,0]=node_vel
+
+#Plot of the velocity from the exact and calculated solutions
+#figure(1)
+#set(gcf,'Position',[1 1 1580 1150])
+#subplot(2,2,1)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#vel,'FaceColor','interp','EdgeColor','none');
+#title('Modelled velocity','FontSize',14,'FontWeight','bold')
+#colorbar; 
+#caxis([0 200]);
+   
+#subplot(2,2,2)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#vel_obs,'FaceColor','interp','EdgeColor','none');
+#title('Analytical velocity','FontSize',14,'FontWeight','bold')
+#colorbar; 
+#caxis([0 200]);
+
+#subplot(2,2,3)
+#hold on;
+#plot(sqrt((md.mesh.x(1:md.mesh.numberofvertices2d)).^2+(md.mesh.y(1:md.mesh.numberofvertices2d)).^2),vel,'r.');
+#plot(sqrt((md.mesh.x2d).^2+(md.mesh.y2d).^2),vel_obs,'b.');
+#title('Analytical vs calculated velocity','FontSize',14,'FontWeight','bold');
+#xlabel('distance to the center of the icesheet [m]','FontSize',14,'FontWeight','bold');
+#ylabel('velocity [m/yr]','FontSize',14,'FontWeight','bold');
+#legend('calculated velocity','exact velocity');
+#axis([0 750000 0 200]);
+#hold off;
+
+#subplot(2,2,4)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#abs(vel-vel_obs)./vel_obs*100,'FaceColor','interp','EdgeColor','none');
+#title('Relative misfit [%]','FontSize',14,'FontWeight','bold')
+#colorbar;
+#caxis([0 100]);
+
+if printingflag:
+	pass
+#	set(gcf,'Color','w')
+#	printmodel('hutterstatic','png','margin','on','marginsize',25,'frame','off','resolution',0.7,'hardcopy','off');
+#	system(['mv hutterstatic.png ' ISSM_DIR '/website/doc_pdf/validation/Images/EISMINT/IceSheet']);
+
+#Fields and tolerances to track changes
+field_names     =[ \
+	'Vx','Vy','Vel', \
+]
+field_tolerances=[ \
+	1e-13,1e-13,1e-13, \
+]
+field_values=[ \
+	vx,vy,vel, \
+]
Index: /issm/trunk-jpl/test/NightlyRun/test1206.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test1206.py	(revision 14157)
+++ /issm/trunk-jpl/test/NightlyRun/test1206.py	(revision 14157)
@@ -0,0 +1,108 @@
+import numpy
+from model import *
+from roundmesh import *
+from setmask import *
+from parameterize import *
+from setflowequation import *
+from EnumDefinitions import *
+from solve import *
+from MatlabFuncs import *
+
+"""
+The aim of this program is to compare a model with an analytical solution given in MacAyeal EISMINT : Lessons in Ice-Sheet Modeling.
+"""
+
+printingflag=False
+
+numlayers=10
+resolution=30000.
+
+#To begin with the numerical model
+md=model()
+md=roundmesh(md,750000.,resolution)
+md=setmask(md,'','')    #We can not test iceshelves nor ice rises with this analytical solution
+md=parameterize(md,'../Par/RoundSheetStaticEISMINT.py')
+
+#Calculation of the analytical 2d velocity field
+constant=0.3
+vx_obs=constant/2.*md.mesh.x.reshape(-1,1)*(md.geometry.thickness)**-1
+vy_obs=constant/2.*md.mesh.y.reshape(-1,1)*(md.geometry.thickness)**-1
+vel_obs=numpy.sqrt((md.inversion.vx_obs)**2+(md.inversion.vy_obs)**2)
+
+#We extrude the model to have a 3d model
+md.extrude(numlayers,1.)
+md=setflowequation(md,'pattyn','all')
+
+#Spc the nodes on the bed
+pos=numpy.nonzero(md.mesh.vertexonbed)
+md.diagnostic.spcvx[pos]=0.
+md.diagnostic.spcvy[pos]=0.
+md.diagnostic.spcvz[pos]=0.
+
+#Now we can solve the problem 
+md.cluster=generic('name',oshostname(),'np',8)
+md=solve(md,DiagnosticSolutionEnum())
+
+#Calculate the depth averaged velocity field (2d):
+vx=md.results.DiagnosticSolution.Vx
+vy=md.results.DiagnosticSolution.Vy
+vel=numpy.zeros((md.mesh.numberofvertices2d,1))
+
+for i in xrange(0,md.mesh.numberofvertices2d):
+	node_vel=0.
+	for j in xrange(0,md.mesh.numberoflayers-1):
+		node_vel=node_vel+1./(2.*(md.mesh.numberoflayers-1))*\
+			(numpy.sqrt(vx[i+(j+1)*md.mesh.numberofvertices2d,0]**2+vy[i+(j+1)*md.mesh.numberofvertices2d,0]**2)+\
+			numpy.sqrt(vx[i+j*md.mesh.numberofvertices2d,0]**2+vy[i+j*md.mesh.numberofvertices2d,0]**2))
+	vel[i,0]=node_vel
+
+#Plot of the velocity from the exact and calculated solutions
+#figure(1)
+#subplot(2,2,1)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#vel,'FaceColor','interp','EdgeColor','none');
+#title('Modelled velocity','FontSize',14,'FontWeight','bold')
+#colorbar; 
+#caxis([0 200]);
+   
+#subplot(2,2,2)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#vel_obs,'FaceColor','interp','EdgeColor','none');
+#title('Analytical velocity','FontSize',14,'FontWeight','bold')
+#colorbar; 
+#caxis([0 200]);
+
+#subplot(2,2,3)
+#hold on;
+#plot(sqrt((md.mesh.x(1:md.mesh.numberofvertices2d)).^2+(md.mesh.y(1:md.mesh.numberofvertices2d)).^2),vel,'r.');
+#plot(sqrt((md.mesh.x2d).^2+(md.mesh.y2d).^2),vel_obs,'b.');
+#title('Analytical vs calculated velocity','FontSize',14,'FontWeight','bold');
+#xlabel('distance to the center of the icesheet [m]','FontSize',14,'FontWeight','bold');
+#ylabel('velocity [m/yr]','FontSize',14,'FontWeight','bold');
+#legend('calculated velocity','exact velocity');
+#axis([0 750000 0 200]);
+#hold off;
+
+#subplot(2,2,4)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#abs(vel-vel_obs)./vel_obs*100,'FaceColor','interp','EdgeColor','none');
+#title('Relative misfit [%]','FontSize',14,'FontWeight','bold')
+#colorbar;
+#caxis([0 100]);
+
+if printingflag:
+	pass
+#	set(gcf,'Color','w')
+#	printmodel('pattynstatic','png','margin','on','marginsize',25,'frame','off','resolution',0.7,'hardcopy','off');
+#	system(['mv pattynstatic.png ' ISSM_DIR '/website/doc_pdf/validation/Images/EISMINT/IceSheet']);
+
+#Fields and tolerances to track changes
+field_names     =[ \
+	'Vx','Vy','Vel', \
+]
+field_tolerances=[ \
+	1e-12,1e-12,1e-12, \
+]
+field_values=[ \
+	vx,vy,vel, \
+]
Index: /issm/trunk-jpl/test/NightlyRun/test1207.py
===================================================================
--- /issm/trunk-jpl/test/NightlyRun/test1207.py	(revision 14157)
+++ /issm/trunk-jpl/test/NightlyRun/test1207.py	(revision 14157)
@@ -0,0 +1,108 @@
+import numpy
+from model import *
+from roundmesh import *
+from setmask import *
+from parameterize import *
+from setflowequation import *
+from EnumDefinitions import *
+from solve import *
+from MatlabFuncs import *
+
+"""
+The aim of this program is to compare a model with an analytical solution given in MacAyeal EISMINT : Lessons in Ice-Sheet Modeling.
+"""
+
+printingflag=False
+
+numlayers=10
+resolution=30000.
+
+#To begin with the numerical model
+md=model()
+md=roundmesh(md,750000.,resolution)
+md=setmask(md,'','')    #We can not test iceshelves nor ice rises with this analytical solution
+md=parameterize(md,'../Par/RoundSheetStaticEISMINT.py')
+
+#Calculation of the analytical 2d velocity field
+constant=0.3
+vx_obs=constant/2.*md.mesh.x.reshape(-1,1)*(md.geometry.thickness)**-1
+vy_obs=constant/2.*md.mesh.y.reshape(-1,1)*(md.geometry.thickness)**-1
+vel_obs=numpy.sqrt((md.inversion.vx_obs)**2+(md.inversion.vy_obs)**2)
+
+#We extrude the model to have a 3d model
+md.extrude(numlayers,1.)
+md=setflowequation(md,'stokes','all')
+
+#Spc the nodes on the bed
+pos=numpy.nonzero(md.mesh.vertexonbed)
+md.diagnostic.spcvx[pos]=0.
+md.diagnostic.spcvy[pos]=0.
+md.diagnostic.spcvz[pos]=0.
+
+#Now we can solve the problem 
+md.cluster=generic('name',oshostname(),'np',8)
+md=solve(md,DiagnosticSolutionEnum())
+
+#Calculate the depth averaged velocity field (2d):
+vx=md.results.DiagnosticSolution.Vx
+vy=md.results.DiagnosticSolution.Vy
+vel=numpy.zeros((md.mesh.numberofvertices2d,1))
+
+for i in xrange(0,md.mesh.numberofvertices2d):
+	node_vel=0.
+	for j in xrange(0,md.mesh.numberoflayers-1):
+		node_vel=node_vel+1./(2.*(md.mesh.numberoflayers-1))*\
+			(numpy.sqrt(vx[i+(j+1)*md.mesh.numberofvertices2d,0]**2+vy[i+(j+1)*md.mesh.numberofvertices2d,0]**2)+\
+			numpy.sqrt(vx[i+j*md.mesh.numberofvertices2d,0]**2+vy[i+j*md.mesh.numberofvertices2d,0]**2))
+	vel[i,0]=node_vel
+
+#Plot of the velocity from the exact and calculated solutions
+#figure(1)
+#subplot(2,2,1)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#vel,'FaceColor','interp','EdgeColor','none');
+#title('Modelled velocity','FontSize',14,'FontWeight','bold')
+#colorbar; 
+#caxis([0 200]);
+   
+#subplot(2,2,2)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#vel_obs,'FaceColor','interp','EdgeColor','none');
+#title('Analytical velocity','FontSize',14,'FontWeight','bold')
+#colorbar; 
+#caxis([0 200]);
+
+#subplot(2,2,3)
+#hold on;
+#plot(sqrt((md.mesh.x(1:md.mesh.numberofvertices2d)).^2+(md.mesh.y(1:md.mesh.numberofvertices2d)).^2),vel,'r.');
+#plot(sqrt((md.mesh.x2d).^2+(md.mesh.y2d).^2),vel_obs,'b.');
+#title('Analytical vs calculated velocity','FontSize',14,'FontWeight','bold');
+#xlabel('distance to the center of the icesheet [m]','FontSize',14,'FontWeight','bold');
+#ylabel('velocity [m/yr]','FontSize',14,'FontWeight','bold');
+#legend('calculated velocity','exact velocity');
+#axis([0 750000 0 200]);
+#hold off;
+
+#subplot(2,2,4)
+#p=patch('Faces',md.mesh.elements2d,'Vertices',[md.mesh.x2d md.mesh.y2d],'FaceVertexCData',...
+#abs(vel-vel_obs)./vel_obs*100,'FaceColor','interp','EdgeColor','none');
+#title('Relative misfit [%]','FontSize',14,'FontWeight','bold')
+#colorbar;
+#caxis([0 100]);
+
+if printingflag:
+	pass
+#	set(gcf,'Color','w')
+#	printmodel('stokesstatic','png','margin','on','marginsize',25,'frame','off','resolution',0.7,'hardcopy','off');
+#	system(['mv stokesstatic.png ' ISSM_DIR '/website/doc_pdf/validation/Images/EISMINT/IceSheet']);
+
+#Fields and tolerances to track changes
+field_names     =[ \
+	'Vx','Vy','Vel', \
+]
+field_tolerances=[ \
+	1e-12,1e-12,1e-12, \
+]
+field_values=[ \
+	vx,vy,vel, \
+]
