Index: /issm/trunk/src/c/parallel/diagnostic_core.cpp
===================================================================
--- /issm/trunk/src/c/parallel/diagnostic_core.cpp	(revision 3751)
+++ /issm/trunk/src/c/parallel/diagnostic_core.cpp	(revision 3752)
@@ -91,7 +91,7 @@
 	//for qmu analysis, be sure the velocity input we are starting from  is the one in the parameters: */
 	if(qmu_analysis){
-		model->FindParamByAnalysisAndSub(&vx,&dummy,VxEnum,DiagnosticAnalysisEnum,HorizAnalysisEnum); model->UpdateInputsFromVector(vx,VxEnum,VertexEnum);
-		model->FindParamByAnalysisAndSub(&vy,&dummy,VyEnum,DiagnosticAnalysisEnum,HorizAnalysisEnum); model->UpdateInputsFromVector(vy,VyEnum,VertexEnum);
-		model->FindParamByAnalysisAndSub(&vz,&dummy,VzEnum,DiagnosticAnalysisEnum,HorizAnalysisEnum); model->UpdateInputsFromVector(vz,VzEnum,VertexEnum);
+		fem_dh->FindParam(&vx,&dummy,VxEnum); model->UpdateInputsFromVector(vx,VxEnum,VertexEnum);
+		fem_dh->FindParam(&vy,&dummy,VyEnum); model->UpdateInputsFromVector(vy,VyEnum,VertexEnum);
+		fem_dh->FindParam(&vz,&dummy,VzEnum); model->UpdateInputsFromVector(vz,VzEnum,VertexEnum);
 	}
 
Index: /issm/trunk/src/m/solutions/jpl/diagnostic.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/diagnostic.m	(revision 3751)
+++ /issm/trunk/src/m/solutions/jpl/diagnostic.m	(revision 3752)
@@ -29,12 +29,4 @@
 	md.dof=modelsize(models);
 
-	%initialize inputs
-	inputs=inputlist;
-	inputs=add(inputs,'velocity',models.dh.parameters.u_g,'doublevec',3,models.dh.parameters.numberofnodes);
-	if md.control_analysis,
-		inputs=add(inputs,'velocity_obs',models.dh.parameters.u_g_obs,'doublevec',2,models.dh.parameters.numberofnodes);
-		inputs=add(inputs,'weights',models.dh.parameters.weights,'doublevec',1,models.dh.parameters.numberofnodes);
-	end
-
 	%compute solution
 	if ~models.dh.parameters.qmu_analysis,
Index: /issm/trunk/src/m/solutions/jpl/diagnostic_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/diagnostic_core.m	(revision 3751)
+++ /issm/trunk/src/m/solutions/jpl/diagnostic_core.m	(revision 3752)
@@ -1,7 +1,7 @@
-function results=diagnostic_core(models,inputs);
+function results=diagnostic_core(model);
 %DIAGNOSTIC_CORE - compute the core velocity field 
 %
 %   Usage:
-%      results=diagnostic_core(models,inputs);
+%      results=diagnostic_core(models);
 %
 
@@ -20,28 +20,33 @@
 isstokes=m_ds.parameters.isstokes;
 numrifts=m_dhu.parameters.numrifts;
+qmu_analysis=m_dh.parameters.qmu_analysis;
 
-debug=1;
+	
+%for qmu analysis, be sure the velocity input we are starting from  is the one in the parameters: 
+if qmu_analysis,
+	UpdateInputsFromVector(models,m_dh.vx,VxEnum,VertexEnum);
+	UpdateInputsFromVector(models,m_dh.vy,VyEnum,VertexEnum);
+	UpdateInputsFromVector(models,m_dh.vz,VzEnum,VertexEnum);
+end
+
+%Compute slopes: 
+[surfaceslopex,surfaceslopey]=slope_core(m_sl,SurfaceAnalysisEnum);
+[bedslopex,bedslopey]=slope_core(m_sl,BedAnalysisEnum);
+
+%Update:
+UpdateInputsFromVector(models,surfaceslopex,SurfaceSlopexEnum,VertexEnum);
+UpdateInputsFromVector(models,surfaceslopey,SurfaceSlopeyEnum,VertexEnum);
+UpdateInputsFromVector(models,bedslopex,BedSlopexEnum,VertexEnum);
+UpdateInputsFromVector(models,bedslopey,BedSlopeyEnum,VertexEnum);
+
+
 
 if ishutter,
 
-	displaystring(verbose,'\n%s',['computing surface slope (x and y derivatives)...']);
-	slopex=diagnostic_core_linear(m_sl,inputs,SlopecomputeAnalysisEnum(),SurfaceXAnalysisEnum());
-	slopey=diagnostic_core_linear(m_sl,inputs,SlopecomputeAnalysisEnum(),SurfaceYAnalysisEnum());
-
-	if dim==3,
-		displaystring(verbose,'\n%s',['extruding slopes in 3d...']);
-		slopex=FieldExtrude(m_sl.elements,m_sl.nodes,m_sl.vertices,m_sl.loads,m_sl.materials,m_sl.parameters,slopex,'slopex',0);
-		slopey=FieldExtrude(m_sl.elements,m_sl.nodes,m_sl.vertices,m_sl.loads,m_sl.materials,m_sl.parameters,slopey,'slopey',0);
-	end
-
-	displaystring(verbose,'\n%s',['computing slopes...']);
-	inputs=add(inputs,'surfaceslopex',slopex,'doublevec',m_sl.parameters.numberofdofspernode,m_sl.parameters.numberofnodes);
-	inputs=add(inputs,'surfaceslopey',slopey,'doublevec',m_sl.parameters.numberofdofspernode,m_sl.parameters.numberofnodes);
-
 	displaystring(verbose,'\n%s',['computing hutter velocities...']);
-	u_g=diagnostic_core_linear(m_dhu,inputs,DiagnosticAnalysisEnum(),HutterAnalysisEnum());
+	u_g=diagnostic_core_linear(m_dhu,DiagnosticAnalysisEnum(),HutterAnalysisEnum());
 
 	displaystring(verbose,'\n%s',['computing pressure according to MacAyeal...']);
-	p_g=ComputePressure(m_dhu.elements,m_dhu.nodes,m_dhu.vertices,m_dhu.loads,m_dhu.materials,m_dhu.parameters,inputs,DiagnosticAnalysisEnum(),HutterAnalysisEnum());
+	p_g=ComputePressure(m_dhu.elements,m_dhu.nodes,m_dhu.vertices,m_dhu.loads,m_dhu.materials,m_dhu.parameters,DiagnosticAnalysisEnum(),HutterAnalysisEnum());
 
 	displaystring(verbose,'\n%s',['update boundary conditions for macyeal pattyn using hutter results...']);
@@ -56,8 +61,10 @@
 
 	displaystring(verbose,'\n%s',['computing horizontal velocities...']);
-	[u_g m_dh.loads]=diagnostic_core_nonlinear(m_dh,inputs,DiagnosticAnalysisEnum(),HorizAnalysisEnum());
+	[u_g m_dh.loads]=diagnostic_core_nonlinear(m_dh,DiagnosticAnalysisEnum(),HorizAnalysisEnum());
 
-	displaystring(verbose,'\n%s',['computing pressure according to MacAyeal...']);
-	p_g=ComputePressure(m_dh.elements,m_dh.nodes,m_dh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,inputs,DiagnosticAnalysisEnum(),HorizAnalysisEnum());
+	if dim==2,
+		displaystring(verbose,'\n%s',['computing pressure according to MacAyeal...']);
+		p_g=ComputePressure(m_dh.elements,m_dh.nodes,m_dh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,DiagnosticAnalysisEnum(),HorizAnalysisEnum());
+	end
 end