Index: /issm/trunk/src/m/model/ismodelselfconsistent.m
===================================================================
--- /issm/trunk/src/m/model/ismodelselfconsistent.m	(revision 7625)
+++ /issm/trunk/src/m/model/ismodelselfconsistent.m	(revision 7626)
@@ -11,83 +11,4 @@
 	help ismodelselfconsistent
 	error('ismodelselfconsistent error message: wrong usage');
-end
-%}}}
-
-%recursive call for TRANSIENT {{{1
-if (md.analysis_type==Transient2DSolutionEnum | md.analysis_type==Transient3DSolutionEnum),
-	if md.dt<=0,
-		error('model not consistent: field dt must be positive for a transient run')
-	end
-
-	%recursive call to ismodelselfconsistent
-	if (md.dim==2),
-		analysis=[DiagnosticSolutionEnum PrognosticSolutionEnum];
-	else
-		analysis=[DiagnosticSolutionEnum PrognosticSolutionEnum ThermalSolutionEnum];
-	end
-
-	for i=1:length(analysis),
-		md.analysis_type=analysis(i);
-		ismodelselfconsistent(md);
-	end
-
-	if(md.cfl_coefficient>1 | md.cfl_coefficient<0),
-		error(['model not consistent: model ' md.name ' cfl_coefficient field should between  0 and 1']);
-	end
-end
-%}}}
-%recursive call for STEADYSTATE{{{1
-if md.analysis_type==SteadystateSolutionEnum,
-
-	%NDT
-	if md.dt~=0,
-		error(['model not consistent: for a steadystate computation, dt must be zero.']);
-	end
-
-	%PRESSURE
-	if isnans(md.pressure),
-		error(['model not consistent: for a steadystate computation, the model must have an initial pressure, even lithostatic will do.']);
-	end
-
-	%eps: 
-	if isnan(md.eps_rel),
-		error(['model not consistent: for a steadystate computation, eps_rel (relative convergence criterion) must be defined!']);
-	end
-
-	%recursive call to ismodelselfconsistent
-	analysis=[DiagnosticSolutionEnum ThermalSolutionEnum];
-
-	for i=1:length(analysis),
-		md.analysis_type=analysis(i);
-		ismodelselfconsistent(md);
-	end
-
-end
-%}}}
-%recursive call for GROUNDINGLINEMIGRATION2D{{{1
-if md.analysis_type==GroundingLineMigration2DSolutionEnum,
-	if strcmpi(md.cluster.name,'none'),
-		error(['model not consistent: ' md.analysis_type ' is only implemented in parallel mode !'])
-	end
-
-	if md.dt<=0,
-		error('model not consistent: field dt must be positive for a transient run')
-	end
-
-	%recursive call to ismodelselfconsistent
-	if (md.dim==2),
-		analysis=[DiagnosticSolutionEnum PrognosticSolutionEnum];
-	else
-		error(['model not consistent: for a ' md.analysis_type ' computation, the grounding line module is only implemented in 2d !'])
-	end
-
-	for i=1:length(analysis),
-		md.analysis_type=analysis(i);
-		ismodelselfconsistent(md);
-	end
-
-	if(md.cfl_coefficient>1 | md.cfl_coefficient<0),
-		error(['model not consistent: model ' md.name ' cfl_coefficient field should between  0 and 1']);
-	end
 end
 %}}}
@@ -276,125 +197,4 @@
 IsConsistent(md.cluster);
 %}}}
-
-% Solution checks
-%DIAGNOSTIC{{{1
-if md.analysis_type==DiagnosticSolutionEnum,
-
-	%HUTTER ON ICESHELF WARNING
-	if any(md.elements_type==HutterApproximationEnum & md.elementoniceshelf),
-		disp(sprintf('\n !!! Warning: Hutter''s model is not consistent on ice shelves !!!\n'));
-	end
-
-	%SINGULAR
-	if ~any(sum(md.spcvelocity(:,1:2),2)==2),
-		error(['model not consistent: model ' md.name ' is not well posed (singular). You need at least one grid with fixed velocity!'])
-	end
-
-	%DIRICHLET IF THICKNESS <= 0
-	if any(md.thickness<=0),
-		pos=find(md.thickness<=0);
-		if any(find(md.spcthickness(pos,1)==0)),
-			error(['model not consistent: model ' md.name ' has some grids with 0 thickness']);
-		end
-	end
-
-	%INITIAL VELOCITY
-	if length(md.vx)==md.numberofgrids & length(md.vy)==md.numberofgrids,
-		fields={'vx','vy'};
-		checknan(md,fields);
-	end
-end
-%}}}
-%PROGNOSTIC{{{1
-if md.analysis_type==PrognosticSolutionEnum,
-
-	%INITIAL VELOCITIES
-	fields={'vx','vy'};
-	checksize(md,fields,[md.numberofgrids 1]);
-	checknan(md,fields);
-
-	%CHECK THAT WE ARE NOT FULLY CONSTRAINED
-	if (md.dim==2),
-		if isempty(find(~md.spcthickness(:,1))),
-			error(['model not consistent: model ' md.name ' is totally constrained for prognostic, no need to solve!']);
-		end
-	end
-
-end
-%}}}
-%THERMAL {{{1
-%THERMAL STEADY AND THERMAL TRANSIENT
-if md.analysis_type==ThermalSolutionEnum,
-
-	%EXTRUSION
-	if (md.dim==2),
-		error(['model not consistent: for a ' md.analysis_type ' computation, the model must be 3d, extrude it first!'])
-	end
-
-	%CHECK THAT WE ARE NOT FULLY CONSTRAINED
-	if isempty(find(~md.spctemperature(:,1))),
-		error(['model not consistent: model ' md.name ' is totally constrained for temperature, no need to solve!']);
-	end
-
-	%VELOCITIES AND PRESSURE
-	fields={'vx','vy','vz','pressure','geothermalflux'};
-	checksize(md,fields,[md.numberofgrids 1]);
-	checknan(md,fields);
-
-end
-
-%THERMAL TRANSIENT
-if md.analysis_type==ThermalSolutionEnum & md.dt~=0,
-
-	%DT and NDT
-	fields={'dt','ndt'};
-	checkgreaterstrict(md,fields,0);
-
-	%INITIAL TEMPERATURE, MELTING AND ACCUMULATION
-	fields={'temperature','accumulation_rate','melting_rate'};
-	checksize(md,fields,[md.numberofgrids 1]);
-	checknan(md,fields);
-
-	%INITIAL TEMPERATURE
-	fields={'temperature','spctemperature(:,2)','observed_temperature'};
-	checkgreater(md,fields,0)
-
-end
-%}}}
-%BALANCEDTHICKNESS{{{1
-if md.analysis_type==BalancedthicknessSolutionEnum
-
-	%VELOCITIES MELTING AND ACCUMULATION
-	fields={'vx','vy','accumulation_rate','melting_rate','dhdt'};
-	checksize(md,fields,[md.numberofgrids 1]);
-	checknan(md,fields);
-
-	%SPC				 
-%	if ~md.prognostic_DG,
-%		if any(md.spcthickness(find(md.gridonboundary))~=1),		 
-%			error(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcthickness']);			 
-%		end 
-%	end
-
-	%Triangle with zero velocity
-	if any(sum(abs(md.vx(md.elements)),2)==0 & sum(abs(md.vy(md.elements)),2)==0)
-		error('model not consistent: at least one triangle has all its vertices with a zero velocity');
-	end
-end
-%}}}
-%BALANCEDVELOCITIES{{{1
-if md.analysis_type==BalancedvelocitiesSolutionEnum
-
-	%VELOCITIES MELTING AND ACCUMULATION
-	fields={'vx','vy','accumulation_rate','melting_rate'};
-	checksize(md,fields,[md.numberofgrids 1]);
-	checknan(md,fields);
-
-	%SPC
-	if any(md.spcvelocity(find(md.gridonboundary),[1:2])~=1),
-		error(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcvelocity']);
-	end
-end
-%}}}
 %CONTROL{{{1
 if md.control_analysis,
@@ -488,5 +288,197 @@
 %}}}
 
-end
+%Solution specific check
+%TRANSIENT {{{1
+if (md.analysis_type==Transient2DSolutionEnum | md.analysis_type==Transient3DSolutionEnum),
+
+	if md.dt<=0,
+		error('model not consistent: field dt must be positive for a transient run')
+	end
+
+	if(md.cfl_coefficient>1 | md.cfl_coefficient<0),
+		error(['model not consistent: model ' md.name ' cfl_coefficient field should between  0 and 1']);
+	end
+end
+%}}}
+%STEADYSTATE{{{1
+if md.analysis_type==SteadystateSolutionEnum,
+
+	%NDT
+	if md.dt~=0,
+		error(['model not consistent: for a steadystate computation, dt must be zero.']);
+	end
+
+	%eps: 
+	if isnan(md.eps_rel),
+		error(['model not consistent: for a steadystate computation, eps_rel (relative convergence criterion) must be defined!']);
+	end
+end
+%}}}
+%GROUNDINGLINEMIGRATION2D{{{1
+if md.analysis_type==GroundingLineMigration2DSolutionEnum,
+	if strcmpi(md.cluster.name,'none'),
+		error(['model not consistent: ' md.analysis_type ' is only implemented in parallel mode !'])
+	end
+
+	if md.dt<=0,
+		error('model not consistent: field dt must be positive for a transient run')
+	end
+
+	%recursive call to ismodelselfconsistent
+	if (md.dim~=2),
+		error(['model not consistent: for a ' md.analysis_type ' computation, the grounding line module is only implemented in 2d !'])
+	end
+
+	if(md.cfl_coefficient>1 | md.cfl_coefficient<0),
+		error(['model not consistent: model ' md.name ' cfl_coefficient field should between  0 and 1']);
+	end
+end
+%}}}
+%ANALYSESCHECKS {{{1
+[analyses,numanalyses,solutioncore]=SolutionConfiguration(md.analysis_type);
+for i=1:numanalyses,
+	checkforanalysis(md,analyses(i));
+end
+%}}}
+end
+
+%analysis checks
+%checkforanalysis {{{1
+function checkforanalysis(md,analysis_enum)
+
+	switch(analysis_enum),
+		case DiagnosticHorizAnalysisEnum,
+			% {{{2
+			%SINGULAR
+			if ~any(sum(md.spcvelocity(:,1:2),2)==2),
+				error(['model not consistent: model ' md.name ' is not well posed (singular). You need at least one grid with fixed velocity!'])
+			end
+
+			%DIRICHLET IF THICKNESS <= 0
+			if any(md.thickness<=0),
+				pos=find(md.thickness<=0);
+				if any(find(md.spcthickness(pos,1)==0)),
+					error(['model not consistent: model ' md.name ' has some grids with 0 thickness']);
+				end
+			end
+
+			%INITIAL VELOCITY
+			if ~isempty(md.vx) & ~isempty(md.vy),
+				fields={'vx','vy'};
+				checknan(md,fields);
+				checksize(md,fields,[md.numberofgrids 1]);
+			end
+			%}}}
+		case DiagnosticVertAnalysisEnum,
+			% {{{2
+			if md.dim==3,
+				% No checks for now
+			end
+			%}}}
+		case DiagnosticHutterAnalysisEnum,
+			% {{{2
+			%HUTTER ON ICESHELF WARNING
+			if any(md.elements_type==HutterApproximationEnum & md.elementoniceshelf),
+				disp(sprintf('\n !!! Warning: Hutter''s model is not consistent on ice shelves !!!\n'));
+			end
+			%}}}
+		case BedSlopeAnalysisEnum,
+			% {{{2
+			% No checks for now
+			%}}}
+		case SurfaceSlopeAnalysisEnum,
+			% {{{2
+			% No checks for now
+			%}}}
+		case PrognosticAnalysisEnum,
+			% {{{2
+			%INITIAL VELOCITIES
+			fields={'vx','vy'};
+			checksize(md,fields,[md.numberofgrids 1]);
+			checknan(md,fields);
+
+			%CHECK THAT WE ARE NOT FULLY CONSTRAINED
+			if (md.dim==2),
+				if isempty(find(~md.spcthickness(:,1))),
+					error(['model not consistent: model ' md.name ' is totally constrained for prognostic, no need to solve!']);
+				end
+			end
+			%}}}
+		case ThermalAnalysisEnum,
+			% {{{2
+			%EXTRUSION
+			if (md.dim==2),
+				error(['model not consistent: for a ' md.analysis_type ' computation, the model must be 3d, extrude it first!'])
+			end
+
+			%CHECK THAT WE ARE NOT FULLY CONSTRAINED
+			if isempty(find(~md.spctemperature(:,1))),
+				error(['model not consistent: model ' md.name ' is totally constrained for temperature, no need to solve!']);
+			end
+
+			%VELOCITIES AND PRESSURE
+			fields={'vx','vy','vz','pressure','geothermalflux'};
+			checksize(md,fields,[md.numberofgrids 1]);
+			checknan(md,fields);
+
+			%THERMAL TRANSIENT
+			if md.dt~=0,
+
+				%DT and NDT
+				fields={'dt','ndt'};
+				checkgreaterstrict(md,fields,0);
+
+				%INITIAL TEMPERATURE, MELTING AND ACCUMULATION
+				fields={'temperature','accumulation_rate','melting_rate'};
+				checksize(md,fields,[md.numberofgrids 1]);
+				checknan(md,fields);
+
+				%INITIAL TEMPERATURE
+				fields={'temperature','spctemperature(:,2)','observed_temperature'};
+				checkgreater(md,fields,0)
+			end
+			%}}}
+		case MeltingAnalysisEnum,
+			% {{{2
+			% No checks for now
+			%}}}
+		case BalancedthicknessAnalysisEnum,
+			% {{{2
+			%VELOCITIES MELTING AND ACCUMULATION
+			fields={'vx','vy','accumulation_rate','melting_rate','dhdt'};
+			checksize(md,fields,[md.numberofgrids 1]);
+			checknan(md,fields);
+
+			%SPC				 
+			%	if ~md.prognostic_DG,
+			%		if any(md.spcthickness(find(md.gridonboundary))~=1),		 
+			%			error(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcthickness']);			 
+			%		end 
+			%	end
+
+			%Triangle with zero velocity
+			if any(sum(abs(md.vx(md.elements)),2)==0 & sum(abs(md.vy(md.elements)),2)==0)
+				error('model not consistent: at least one triangle has all its vertices with a zero velocity');
+			end
+			%}}}
+		case BalancedvelocitiesAnalysisEnum,
+			% {{{2
+			%VELOCITIES MELTING AND ACCUMULATION
+			fields={'vx','vy','accumulation_rate','melting_rate'};
+			checksize(md,fields,[md.numberofgrids 1]);
+			checknan(md,fields);
+
+			%SPC
+			if any(md.spcvelocity(find(md.gridonboundary),[1:2])~=1),
+				error(['model not consistent: model ' md.name ' should have all the nodes on boundary constrained in field spcvelocity']);
+			end
+			%}}}
+		otherwise
+			disp(['WARNING: no check implemented for analysis ' EnumToString(analysis_enum) '!']);
+			error('stop')
+	end
+
+end
+%}}}
 
 %checks additional functions