Index: /issm/trunk/test/NightlyRun/test1201.m
===================================================================
--- /issm/trunk/test/NightlyRun/test1201.m	(revision 5374)
+++ /issm/trunk/test/NightlyRun/test1201.m	(revision 5375)
@@ -1,58 +1,74 @@
-%This test is a test from the EISMINT for Ice shelves 
-% Vincent Rommelaere 1996
+%This test is a test from the EISMINT for Ice shelves Vincent Rommelaere 1996
 
-%The goal is to test the prognostic model
-%md=bamg(model,'domain','DomainOutline.exp','hmax',4550);
-md=bamg(model,'domain','../Exp/SquareEISMINT.exp','hmax',3000);
-md=geography(md,'all','');
-md=parameterize(md,'../Par/SquareEISMINT.par');
-md=setelementstype(md,'macayeal','all');
+results={};
 
-disp('      initial velocity');
-md.vx=zeros(md.numberofgrids,1);
-md.vy=-400*ones(md.numberofgrids,1);
-md.vz=zeros(md.numberofgrids,1);
+for type=1:3;
+	%The goal is to test the prognostic model
+	md=bamg(model,'domain','../Exp/SquareEISMINT.exp','hmax',3000);
+	md=geography(md,'all','');
+	md=parameterize(md,'../Par/SquareEISMINT.par');
+	md=setelementstype(md,'macayeal','all');
 
-%analytical test
-md.spcvelocity(:,1:3)=1;
-md.spcvelocity(:,4)=0;
-md.spcvelocity(:,5)=-400;
-md.spcvelocity(1,1)=0;
-md.prognostic_DG=1;
+	disp('      initial velocity');
+	md.vx=zeros(md.numberofgrids,1);
+	md.vy=-400*ones(md.numberofgrids,1);
+	md.vz=zeros(md.numberofgrids,1);
 
-%Launch many prognostic as we don't really want a transient
-%FOR NOW:
-i=0;
-time=0;
-md.dummy=struct();
-pos=find(md.y>199999.9);
-connectivity=full(sparse(md.elements(:),1,1));
-while(time<500),
-	disp(['step ' num2str(i) '/' num2str(500/md.dt)]);
-	i=i+1;
-	time=time+md.dt;
-	md=solve(md,'analysis_type',PrognosticSolutionEnum);
-	thickness   =zeros(md.numberofgrids,1);
-	connectivity=sparse(md.results.PrognosticSolution.Thickness.index(:),1,1);
-	thickness   =sparse(md.results.PrognosticSolution.Thickness.index(:),1,md.results.PrognosticSolution.Thickness.value(:));
-	md.thickness=full(thickness./connectivity);
-	md.thickness(pos)=500+100*sin(2*pi*time/200)*ones(size(pos,1),1);
-	md.surface=md.bed+md.thickness;
-	md.dummy(i).thickness=md.thickness;
-	%plotmodel(md,'data',md.thickness,'sectionvalue','CrossLine.exp')
+	%analytical test
+	md.spcvelocity(:,1:3)=1;
+	md.spcvelocity(:,4)=0;
+	md.spcvelocity(:,5)=-400;
+	md.spcvelocity(1,1)=0;
+	if type==1,
+		md.artificial_diffusivity=1;
+		md.prognostic_DG=0;
+	elseif type==2,
+		md.artificial_diffusivity=0;
+		md.prognostic_DG=0;
+	elseif type==3,
+		md.prognostic_DG=1;
+	end
+
+	%Launch many prognostic as we don't really want a transient
+	%FOR NOW:
+	i=0;
+	time=0;
+	md.dummy=struct();
+	pos=find(md.y>199999.9);
+	connectivity=full(sparse(md.elements(:),1,1));
+	while(time<500),
+		disp(['step ' num2str(i) '/' num2str(500/md.dt)]);
+		i=i+1;
+		time=time+md.dt;
+		md=solve(md,'analysis_type',PrognosticSolutionEnum);
+		md.thickness=PatchToVec(md.results.PrognosticSolution.Thickness);
+		md.thickness(pos)=500+100*sin(2*pi*time/200)*ones(size(pos,1),1);
+		md.surface=md.bed+md.thickness;
+		md.dummy(i).thickness=md.thickness;
+	end
+	results{type}=md.thickness;
 end
 
 %plot results
-plotmodel(md,'data',md.thickness)
-%plotmodel(md,'data',[md.thickness 500+100*sin(2*pi/200*(500-md.y/400))],'sectionvalue','../Exp/CrossLineEISMINT.exp')
+[elements,x,y,z,s,h1]=SectionValues(md,results{1},'../Exp/CrossLineEISMINT.exp',100);
+[elements,x,y,z,s,h2]=SectionValues(md,results{2},'../Exp/CrossLineEISMINT.exp',100);
+[elements,x,y,z,s,h3]=SectionValues(md,results{3},'../Exp/CrossLineEISMINT.exp',100);
+[elements,x,y,z,s,hth]=SectionValues(md, 500+100*sin(2*pi/200*(500-md.y/400)),'../Exp/CrossLineEISMINT.exp',100);
+plot(s,h1,'r',s,h2,'b',s,h3,'g',s,hth,'k')
+legend('Art. diff.','No Art. diff.','Discontinue','Theorique')
+
+set(gcf,'Color','w')
+printmodel('eismintmasscthickness','png','margin','on','marginsize',25,'frame','off','resolution',2,'hardcopy','off');
 
 %Fields and tolerances to track changes
 field_names     ={ ...
-	'Thickness', ...
+	'ThicknessArtDigg','ThicknessNoArtDiff','ThicknessDG' ...
 };
 field_tolerances={...
-	1e-13, ...
+	1e-13, 1e-13, 1e-13...
 };
 field_values={
-	md.thickness, ...
+	results{1}, ...
+	results{2}, ...
+	results{3}, ...
 };
