Index: /issm/trunk/src/m/solutions/NewFemModel.m
===================================================================
--- /issm/trunk/src/m/solutions/NewFemModel.m	(revision 4187)
+++ /issm/trunk/src/m/solutions/NewFemModel.m	(revision 4188)
@@ -19,4 +19,11 @@
    [femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.constraints,femmodel.loads,femmodel.materials,femmodel.parameters]=ModelProcessor(md,solution_type,femmodel.analysis_type_list);
 
+	%Initialize some fiels with empty celils
+	femmodel.m_yg=cell(nummodels,1);
+	femmodel.m_Rmg=cell(nummodels,1);
+	femmodel.m_nodesets=cell(nummodels,1);
+	femmodel.m_ys=cell(nummodels,1);
+	femmodel.m_Gmn=cell(nummodels,1);
+
    %now, go through all analyses types and post-process datasets
    for i=1:length(analysis_types),
@@ -26,25 +33,30 @@
 
 	   displaystring(md.verbose,'%s','      generating degrees of freedofemmodel...');
-	   if ~isfield(femmodel,'part') [femmodel.vertices,femmodel.part,femmodel.tpart]=VerticesDof(femmodel.vertices, femmodel.parameters); %do not create partition vector twice! we only have one set of vertices!
+		if ~isfield(femmodel,'part'),
+			[femmodel.vertices,femmodel.part,femmodel.tpart]=VerticesDof(femmodel.vertices, femmodel.parameters); %do not create partition vector twice! we only have one set of vertices!
+		end
 		   
-		   [femmodel.nodes]=NodesDof(femmodel.nodes,femmodel.parameters);
+		[femmodel.nodes]=NodesDof(femmodel.nodes,femmodel.parameters);
 
-		   displaystring(md.verbose,'%s','      generating single point constraints...');
-		   [femmodel.nodes,femmodel.m_yg(i)]=SpcNodes(femmodel.nodes,femmodel.constraints,analysis_type);
+		displaystring(md.verbose,'%s','      generating single point constraints...');
+		[femmodel.nodes,femmodel.m_yg{i}]=SpcNodes(femmodel.nodes,femmodel.constraints,analysis_type);
 
-		   displaystring(md.verbose,'%s','      generating rigid body constraints...');
-		   [femmodel.m_Rmg(i),femmodel.nodes]=MpcNodes(femmodel.nodes,femmodel.constraints,analysis_types);
+		displaystring(md.verbose,'%s','      generating rigid body constraints...');
+		[femmodel.m_Rmg{i},femmodel.nodes]=MpcNodes(femmodel.nodes,femmodel.constraints,analysis_types);
 
-		   displaystring(md.verbose,'%s','      generating node sets...');
-		   femmodel.m_nodesets(i)=BuildNodeSets(femmodel.nodes,analysis_type);
+		displaystring(md.verbose,'%s','      generating node sets...');
+		femmodel.m_nodesets{i}=BuildNodeSets(femmodel.nodes,analysis_type);
 
-		   displaystring(md.verbose,'%s','      reducing single point constraints vector...');
-		   femmodel.m_ys(i)=Reducevectorgtos(femmodel.m_yg(i).vector,femmodel.m_nodesets(i));
+		displaystring(md.verbose,'%s','      reducing single point constraints vector...');
+		femmodel.m_ys{i}=Reducevectorgtos(femmodel.m_yg{i},femmodel.m_nodesets{i});
 
-		   displaystring(md.verbose,'%s','      normalizing rigid body constraints matrix...');
-		   femmodel.m_Gmn(i)= NormalizeConstraints(femmodel.m_Rmg(i),femmodel.m_nodesets(i));
+		displaystring(md.verbose,'%s','      normalizing rigid body constraints matrix...');
+		femmodel.m_Gmn{i}= NormalizeConstraints(femmodel.m_Rmg{i},femmodel.m_nodesets{i});
 
-		   displaystring(md.verbose,'%s','      configuring element and loads...');
-		   [femmodel.elements,femmodel.loads,femmodel.nodes,femmodel.parameters] = ConfigureObjects( femmodel.elements, femmodel.loads, femmodel.nodes, femmodel.vertices,femmodel.materials,femmodel.parameters);
-	   end
+		displaystring(md.verbose,'%s','      configuring element and loads...');
+		[femmodel.elements,femmodel.loads,femmodel.nodes,femmodel.parameters] = ConfigureObjects( femmodel.elements, femmodel.loads, femmodel.nodes, femmodel.vertices,femmodel.materials,femmodel.parameters);
+
    end
+
+	displaystring(md.verbose,'%s','   process parameters...');
+	femmodel.parameters=ProcessParams(femmodel.parameters);
Index: /issm/trunk/src/m/solutions/SetCurrentAnalysis.m
===================================================================
--- /issm/trunk/src/m/solutions/SetCurrentAnalysis.m	(revision 4187)
+++ /issm/trunk/src/m/solutions/SetCurrentAnalysis.m	(revision 4188)
@@ -11,12 +11,13 @@
 
 	%first, look for analysis: 
-	int found=0;
+	found=0;
 	for i=1:length(femmodel.analysis_type_list),
-		if analysis_type_list(i)==analysis_enum,
+		if femmodel.analysis_type_list(i)==analysis_enum,
 			found=i;
 			break;
+		end
 	end
 
-	if (found!=1),
+	if ~found==1,
 		error('SetCurrentAnalysis error message: could not find analysis_type in list of FemModel analyses');
 	end
Index: /issm/trunk/src/m/solutions/diagnostic_core.m
===================================================================
--- /issm/trunk/src/m/solutions/diagnostic_core.m	(revision 4187)
+++ /issm/trunk/src/m/solutions/diagnostic_core.m	(revision 4188)
@@ -16,5 +16,5 @@
 	ismacayealpattyn=femmodel.parameters.IsMacAyealPattyn;
 	isstokes=femmodel.parameters.IsStokes;
-	stokesreconditioning=femmodel.parameters.stokesreconditioning;
+	if(isstokes) stokesreconditioning=femmodel.parameters.stokesreconditioning; end
 	qmu_analysis=femmodel.parameters.QmuAnalysis;
 
@@ -58,6 +58,6 @@
 
 			%"recondition" pressure computed previously:
-			femmodel.elements=InputDuplicate(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureEnum,PressureStokesEnum);
-			femmodel.elements=InputScale(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureStokesEnum,1.0/stokesreconditioning);
+			femmodel.elements=InputDuplicate(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,PressureEnum,PressureStokesEnum);
+			femmodel.elements=InputScale(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,PressureStokesEnum,1.0/stokesreconditioning);
 
 			displaystring(verbose,'\n%s',['update boundary conditions for stokes using velocities previously computed...']);
@@ -71,6 +71,6 @@
 
 	displaystring(verbose,'\n%s',['saving results...']);
-	femmodel.elements=InputToResult(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,VxEnum);
-	femmodel.elements=InputToResult(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,VyEnum);
-	femmodel.elements=InputToResult(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,PressureEnum);
-	if(dim==3) femmodel.elements=InputToResult(femmodel->elements,femmodel->nodes,femmodel->vertices,femmodel->loads,femmodel->materials,femmodel->parameters,VxEnum);
+	femmodel.elements=InputToResult(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,VxEnum);
+	femmodel.elements=InputToResult(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,VyEnum);
+	femmodel.elements=InputToResult(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,PressureEnum);
+	if(dim==3) femmodel.elements=InputToResult(femmodel.elements,femmodel.nodes,femmodel.vertices,femmodel.loads,femmodel.materials,femmodel.parameters,VxEnum); end