Index: /issm/trunk/src/m/solutions/jpl/ControlInitialization.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/ControlInitialization.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/ControlInitialization.m	(revision 3839)
@@ -1,3 +1,3 @@
-function [inputs models]=ControlInitialization(models,inputs);
+function models=ControlInitialization(models);
 
 %recover models
@@ -55,5 +55,5 @@
 % get pressure (reconditionned) and create 4d u_g
 displaystring(verbose,'\n%s',['computing pressure according to Pattyn...']);
-p_g=ComputePressure(m_dh.elements,m_dh.nodes,mdh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,inputs,DiagnosticAnalysisEnum(),HorizAnalysisEnum());
+p_g=ComputePressure(m_dh.elements,m_dh.nodes,mdh.vertices,m_dh.loads,m_dh.materials,m_dh.parameters,DiagnosticAnalysisEnum(),HorizAnalysisEnum());
 p_g=p_g/m_ds.parameters.stokesreconditioning;
 u_g_stokes=zeros(m_ds.nodesets.gsize,1);
Index: /issm/trunk/src/m/solutions/jpl/SpawnCore.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/SpawnCore.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/SpawnCore.m	(revision 3839)
@@ -1,7 +1,7 @@
-function responses=SpawnCore(models,inputs,variables,variabledescriptors,analysis_type,sub_analysis_type,counter);
+function responses=SpawnCore(models,variables,variabledescriptors,analysis_type,sub_analysis_type,counter);
 %SPAWNCORE - for Qmu analysis, using Dakota. Spawn the core solution.
 %
 %   Usage:
-%      responses=SpawnCore(models,inputs,variables,variabledescriptors)
+%      responses=SpawnCore(models,variables,variabledescriptors)
 %
 
@@ -47,5 +47,5 @@
 if analysis_type==DiagnosticAnalysisEnum(),
 
-	results=diagnostic_core(models,inputs);
+	results=diagnostic_core(models);
 
 else
Index: /issm/trunk/src/m/solutions/jpl/balancedthickness.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/balancedthickness.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/balancedthickness.m	(revision 3839)
@@ -25,5 +25,5 @@
 
 	displaystring(md.verbose,'\n%s',['call computational core:']);
-	results=balancedthickness_core(models,inputs,BalancedthicknessAnalysisEnum(),NoneAnalysisEnum());
+	results=balancedthickness_core(models,BalancedthicknessAnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(md.verbose,'\n%s',['load results...']);
Index: /issm/trunk/src/m/solutions/jpl/balancedthickness2.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/balancedthickness2.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/balancedthickness2.m	(revision 3839)
@@ -28,5 +28,5 @@
 
 	displaystring(md.verbose,'\n%s',['call computational core:']);
-	results=balancedthickness2_core(models,inputs,Balancedthickness2AnalysisEnum(),NoneAnalysisEnum());
+	results=balancedthickness2_core(models,Balancedthickness2AnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(md.verbose,'\n%s',['load results...']);
Index: /issm/trunk/src/m/solutions/jpl/balancedthickness2_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/balancedthickness2_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/balancedthickness2_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=balancedthickness2_core(models,inputs,analysis_type,sub_analysis_type)
+function results=balancedthickness2_core(models,analysis_type,sub_analysis_type)
 %BALANCEDTHICKNESS_CORE - linear solution sequence
 %
 %   Usage:
-%      h_g=balancedthickness2_core(m,inputs,analysis_type,sub_analysis_type)
+%      h_g=balancedthickness2_core(m,analysis_type,sub_analysis_type)
 
 	%get FE model
Index: /issm/trunk/src/m/solutions/jpl/balancedthickness_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/balancedthickness_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/balancedthickness_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=balancedthickness_core(models,inputs,analysis_type,sub_analysis_type)
+function results=balancedthickness_core(models,analysis_type,sub_analysis_type)
 %BALANCEDTHICKNESS_CORE - linear solution sequence
 %
 %   Usage:
-%      h_g=balancedthickness_core(m,inputs,analysis_type,sub_analysis_type)
+%      h_g=balancedthickness_core(m,analysis_type,sub_analysis_type)
 
 	%get FE model
Index: /issm/trunk/src/m/solutions/jpl/balancedvelocities.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/balancedvelocities.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/balancedvelocities.m	(revision 3839)
@@ -27,5 +27,5 @@
 
 	displaystring(md.verbose,'\n%s',['call computational core:']);
-	results=balancedvelocities_core(models,inputs,BalancedvelocitiesAnalysisEnum(),NoneAnalysisEnum());
+	results=balancedvelocities_core(models,BalancedvelocitiesAnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(md.verbose,'\n%s',['load results...']);
Index: /issm/trunk/src/m/solutions/jpl/balancedvelocities_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/balancedvelocities_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/balancedvelocities_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=balancedvelocities_core(models,inputs,analysis_type,sub_analysis_type)
+function results=balancedvelocities_core(models,analysis_type,sub_analysis_type)
 %BALANCEDVELOCITIES_CORE - linear solution sequence
 %
 %   Usage:
-%      v_g=balancedvelocities_core(m,inputs,analysis_type,sub_analysis_type)
+%      v_g=balancedvelocities_core(m,analysis_type,sub_analysis_type)
 
 	%get FE model
Index: /issm/trunk/src/m/solutions/jpl/control_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/control_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/control_core.m	(revision 3839)
@@ -1,11 +1,11 @@
-function results=control_core(models,inputs)
+function results=control_core(models)
 %CONTROL_CORE - compute the core inversion
 %
 %   Usage:
-%      results=control_core(models,inputs);
+%      results=control_core(models);
 %
 
 %Preprocess models
-[inputs models]=ControlInitialization(models,inputs);
+models=ControlInitialization(models);
 
 %recover active model.
@@ -37,5 +37,5 @@
 	%In case we are running a steady state control method, compute new temperature field using new parameter distribution: 
 	if model.parameters.ControlSteady;
-		steadystate_results=steadystate_core(models,inputs); t_g=steadystate_results.t_g; 
+		steadystate_results=steadystate_core(models); t_g=steadystate_results.t_g; 
 		inputs=add(inputs,'temperature',t_g,'doublevec',1,model.parameters.NumberOfNodes);
 	end
@@ -46,8 +46,8 @@
 
 	%Update inputs in datasets
-	[model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters]=UpdateFromInputs(model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters,inputs);
+	[model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters]=UpdateFromInputs(model.elements,model.nodes,model.vertices,model.loads,model.materials,model.parameters);
 
 	displaystring(verbose,'\n%s',['      computing gradJ...']);
-	results_grad=gradjcompute_core(models,inputs);
+	results_grad=gradjcompute_core(models);
 	u_g=results_grad.u_g; c(n).grad_g=results_grad.grad_g;
 	if dim==3,
@@ -78,5 +78,5 @@
 
 	displaystring(verbose,'\n%s',['      optimizing along gradient direction...']);
-	[search_scalar c(n).J]=ControlOptimization('objectivefunctionC',0,1,options,models,inputs,param_g,c(n).grad_g,n,model.parameters);
+	[search_scalar c(n).J]=ControlOptimization('objectivefunctionC',0,1,options,models,param_g,c(n).grad_g,n,model.parameters);
 
 	displaystring(verbose,'\n%s',['      updating parameter using optimized search scalar...']);
@@ -119,5 +119,5 @@
 if model.parameters.ControlSteady;
 	inputs=add(inputs,model.parameters.ControlType,param_g,'doublevec',1,model.parameters.NumberOfNodes);
-	steadystate_results=steadystate_core(models,inputs); t_g=steadystate_results.t_g; 
+	steadystate_results=steadystate_core(models); t_g=steadystate_results.t_g; 
 	u_g=steadystate_results.u_g;
 	t_g=steadystate_results.t_g;
@@ -125,5 +125,5 @@
 else
 	inputs=add(inputs,model.parameters.ControlType,param_g,'doublevec',1,model.parameters.NumberOfNodes);
-	results_diag=diagnostic_core(models,inputs);
+	results_diag=diagnostic_core(models);
 	u_g=results_diag.u_g;
 end
Index: /issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/diagnostic_core_linear.m	(revision 3839)
@@ -9,9 +9,9 @@
 
 	%Update inputs in datasets
-	[m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs);
+	[m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters);
 	
 	%system matrices
-	[K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
-	[K_gg, p_g,kmax]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+	[K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
+	[K_gg, p_g,kmax]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 	
 	%Reduce tangent matrix from g size to f size
Index: /issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/diagnostic_core_nonlinear.m	(revision 3839)
@@ -34,5 +34,5 @@
 	
 		%penalties
-		[K_gg , p_g, kmax]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+		[K_gg , p_g, kmax]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 
 		%Reduce tangent matrix from g size to f size
@@ -58,5 +58,5 @@
 		
 		%penalty constraints
-		[loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints( m.elements,m.nodes,m.vertices,loads, m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+		[loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints( m.elements,m.nodes,m.vertices,loads, m.materials,m.parameters,analysis_type,sub_analysis_type);
 		error;
 		
@@ -66,5 +66,5 @@
 		converged=convergence(K_ff,p_f,soln(count).u_f,soln(count-1).u_f,m.parameters);
 			
-		%add convergence status into inputs
+		%add convergence status into 
 		inputs=add(inputs,'converged',converged,'double');
 
@@ -85,5 +85,5 @@
 		inputs=add(inputs,'velocity',u_g,'doublevec',m.parameters.NumberOfDofsPerNode,m.parameters.NumberOfNodes);
 		m.parameters.Kflag=1; m.parameters.Pflag=0; 
-		[K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+		[K_gg, p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 		[K_ff, K_fs] = Reducematrixfromgtof( K_gg, m.Gmn, m.nodesets); 
 		varargout(1)={K_ff};
Index: /issm/trunk/src/m/solutions/jpl/gradjcompute_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/gradjcompute_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/gradjcompute_core.m	(revision 3839)
@@ -1,3 +1,3 @@
-function results=gradjcompute_core(models,inputs);
+function results=gradjcompute_core(models);
 
 %recover active models
@@ -21,5 +21,5 @@
 %Buid Du, difference between observed velocity and model velocity.
 displaystring(verbose,'%s','          computing Du...');
-[Du_g]=Du(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+[Du_g]=Du(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 
 %Reduce adjoint load from g-set to f-set
@@ -35,5 +35,5 @@
 
 %Compute gradJ 
-grad_g=Gradj(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+grad_g=Gradj(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 if (dim==3 & extrude_param),
 	displaystring(verbose,'%s','          extruding gradient...');
@@ -44,5 +44,5 @@
 	%compute initial velocity from diagnostic_core (horiz+vertical)
 	displaystring(verbose,'%s','          compute 3d initial velocity...');
-	results_diag=diagnostic_core(models,inputs);
+	results_diag=diagnostic_core(models);
 	u_g=results_diag.u_g;
 end
Index: /issm/trunk/src/m/solutions/jpl/objectivefunctionC.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/objectivefunctionC.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/objectivefunctionC.m	(revision 3839)
@@ -1,3 +1,3 @@
-function J =objectivefunctionC(search_scalar,models,inputs,p_g,grad_g,n,analysis_type,sub_analysis_type);
+function J =objectivefunctionC(search_scalar,models,p_g,grad_g,n,analysis_type,sub_analysis_type);
 
 %recover active model.
@@ -23,3 +23,3 @@
 
 %Compute misfit for this velocity field. 
-J=CostFunction(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+J=CostFunction(m.elements,m.nodes,m.vertices,m.loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
Index: /issm/trunk/src/m/solutions/jpl/prognostic.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/prognostic.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/prognostic.m	(revision 3839)
@@ -27,5 +27,5 @@
 
 	displaystring(md.verbose,'\n%s',['call computational core:']);
-	results=prognostic_core(models,inputs,PrognosticAnalysisEnum(),NoneAnalysisEnum());
+	results=prognostic_core(models,PrognosticAnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(md.verbose,'\n%s',['load results...']);
Index: /issm/trunk/src/m/solutions/jpl/prognostic2.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/prognostic2.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/prognostic2.m	(revision 3839)
@@ -28,5 +28,5 @@
 
 	displaystring(md.verbose,'\n%s',['call computational core:']);
-	results=prognostic2_core(models,inputs,Prognostic2AnalysisEnum(),NoneAnalysisEnum());
+	results=prognostic2_core(models,Prognostic2AnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(md.verbose,'\n%s',['load results...']);
Index: /issm/trunk/src/m/solutions/jpl/prognostic2_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/prognostic2_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/prognostic2_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=prognostic2_core(models,inputs,analysis_type,sub_analysis_type)
+function results=prognostic2_core(models,analysis_type,sub_analysis_type)
 %PROGNOSTIC2_CORE - linear solution sequence
 %
 %   Usage:
-%      h_g=prognostic2_core(m,inputs,analysis_type,sub_analysis_type)
+%      h_g=prognostic2_core(m,analysis_type,sub_analysis_type)
 
 	%get FE model
Index: /issm/trunk/src/m/solutions/jpl/prognostic_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/prognostic_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/prognostic_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=prognostic_core(models,inputs,analysis_type,sub_analysis_type)
+function results=prognostic_core(models,analysis_type,sub_analysis_type)
 %PROGNOSTIC_CORE - linear solution sequence
 %
 %   Usage:
-%      h_g=prognostic_core(m,inputs,analysis_type,sub_analysis_type)
+%      h_g=prognostic_core(m,analysis_type,sub_analysis_type)
 
 	%get FE model
Index: /issm/trunk/src/m/solutions/jpl/slopecompute.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/slopecompute.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/slopecompute.m	(revision 3839)
@@ -17,10 +17,7 @@
 	md.dof=modelsize(models);
 
-	%initialize inputs
-	inputs=inputlist;
-
 	%compute solution
 	displaystring(md.verbose,'\n%s',['call computational core:']);
-	results=slopecompute_core(models,inputs,SlopecomputeAnalysisEnum(),NoneAnalysisEnum());
+	results=slopecompute_core(models,SlopecomputeAnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(md.verbose,'\n%s',['load results...']);
Index: /issm/trunk/src/m/solutions/jpl/slopecompute_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/slopecompute_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/slopecompute_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=slopecompute_core(models,inputs,analysis_type,sub_analysis_type)
+function results=slopecompute_core(models,analysis_type,sub_analysis_type)
 %SLOPECOMPUTE_CORE - linear solution sequence
 %
 %   Usage:
-%      [sx_g sy_g]=slopecompute_core(m,inputs,analysis_type,sub_analysis_type)
+%      [sx_g sy_g]=slopecompute_core(m,analysis_type,sub_analysis_type)
 
 	%get FE model
Index: /issm/trunk/src/m/solutions/jpl/steadystate.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/steadystate.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/steadystate.m	(revision 3839)
@@ -55,5 +55,5 @@
 			models.ds.parameters.ControlSteady=1;
 			%launch core of control solution.
-			results=control_core(models,inputs);
+			results=control_core(models);
 
 			%process results
@@ -62,5 +62,5 @@
 		else,
 			%launch core of steadystate solution.
-			results=steadystate_core(models,inputs);
+			results=steadystate_core(models);
 		
 			%process results
@@ -70,5 +70,5 @@
 	else
 		%launch dakota driver for steadystate core solution
-		Qmu(models,inputs,models.dh.parameters);
+		Qmu(models,models.dh.parameters);
 	end
 
Index: /issm/trunk/src/m/solutions/jpl/steadystate_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/steadystate_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/steadystate_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=steadystate_core(models,inputs);
+function results=steadystate_core(models);
 %STEADYSTATE_CORE - compute the core temperature and velocity field  at thermal steady state.
 %
 %   Usage:
-%      results=steadystate_core(models,inputs);
+%      results=steadystate_core(models);
 %
 
@@ -42,5 +42,5 @@
 		inputs=add(inputs,'velocity',results_diagnostic.u_g,'doublevec',ndof,m_t.parameters.NumberOfNodes);
 	end
-	results_thermal=thermal_core(models,inputs);
+	results_thermal=thermal_core(models);
 	
 	%add temperature to inputs.
@@ -51,5 +51,5 @@
 	
 	%now compute diagnostic velocity using the steady state temperature.
-	results_diagnostic=diagnostic_core(models,inputs);
+	results_diagnostic=diagnostic_core(models);
 	
 	%convergence? 
Index: /issm/trunk/src/m/solutions/jpl/thermal.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/thermal.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/thermal.m	(revision 3839)
@@ -30,5 +30,5 @@
 	if ~models.t.parameters.QmuAnalysis,
 		%launch core of diagnostic solution.
-		results=thermal_core(models,inputs);
+		results=thermal_core(models);
 	
 		%process results
@@ -37,5 +37,5 @@
 	else
 		%launch dakota driver for diagnostic core solution
-		Qmu(models,inputs,models.t.parameters);
+		Qmu(models,models.t.parameters);
 	end
 
Index: /issm/trunk/src/m/solutions/jpl/thermal_core.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/thermal_core.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/thermal_core.m	(revision 3839)
@@ -1,7 +1,7 @@
-function results=thermal_core(models,inputs)
+function results=thermal_core(models)
 %THERMAL_CORE - core of thermal solution
 %
 %   Usage:
-%      solution=thermal_core(models,inputs)
+%      solution=thermal_core(models)
 
 %get FE model
@@ -15,5 +15,5 @@
 
 	displaystring(m_t.parameters.Verbose,'\n%s',['computing temperatures...']);
-	[results.t_g m_t.loads melting_offset]=thermal_core_nonlinear(m_t,inputs,ThermalAnalysisEnum(),NoneAnalysisEnum());
+	[results.t_g m_t.loads melting_offset]=thermal_core_nonlinear(m_t,,ThermalAnalysisEnum(),NoneAnalysisEnum());
 
 	displaystring(m_t.parameters.Verbose,'\n%s',['computing melting...']);
@@ -43,5 +43,5 @@
 		displaystring(m_t.parameters.Verbose,'\n%s',['    computing temperatures...']);
 		inputs=add(inputs,'temperature',results(n).t_g,'doublevec',1,m_t.parameters.NumberOfNodes);
-		[results(n+1).t_g m_t.loads melting_offset]=thermal_core_nonlinear(m_t,inputs,ThermalAnalysisEnum(),NoneAnalysisEnum());
+		[results(n+1).t_g m_t.loads melting_offset]=thermal_core_nonlinear(m_t,ThermalAnalysisEnum(),NoneAnalysisEnum());
 
 		displaystring(m_t.parameters.Verbose,'\n%s',['    computing melting...']);
Index: /issm/trunk/src/m/solutions/jpl/thermal_core_nonlinear.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/thermal_core_nonlinear.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/thermal_core_nonlinear.m	(revision 3839)
@@ -1,8 +1,8 @@
-function [t_g ,loads, melting_offset]=thermal_core_nonlinear(m,inputs,analysis_type,sub_analysis_type)
+function [t_g ,loads, melting_offset]=thermal_core_nonlinear(m,analysis_type,sub_analysis_type)
 %THERMAL_CORE_NONLINEAR - core of thermal solution sequence.
 %   model is return together with temperature
 %
 %   Usage:
-%      [t_g ,loads, melting_offset]=thermal_core_nonlinear(m,inputs,analysis_type,sub_analysis_type);
+%      [t_g ,loads, melting_offset]=thermal_core_nonlinear(m,analysis_type,sub_analysis_type);
 
 	count=1;
@@ -20,5 +20,5 @@
 
 		%Update inputs in datasets
-		[m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters,inputs);
+		[m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters);
 
 		%system matrices 
@@ -26,13 +26,13 @@
 			if count==1
 				displaystring(m.parameters.Verbose,'%s',['   system matrices']);
-				[K_gg_nopenalty, p_g_nopenalty]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+				[K_gg_nopenalty, p_g_nopenalty]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 			end
 			displaystring(m.parameters.Verbose,'%s',['   penalty system matrices']);
-			[K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+			[K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg_nopenalty,p_g_nopenalty,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 		else
 			displaystring(m.parameters.Verbose,'%s',['   system matrices']);
-			[K_gg , p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+			[K_gg , p_g]=SystemMatrices(m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 			displaystring(m.parameters.Verbose,'%s',['   penalty system matrices']);
-			[K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+			[K_gg , p_g, melting_offset]=PenaltySystemMatrices(K_gg,p_g,m.elements,m.nodes,m.vertices,loads,m.materials,m.parameters,analysis_type,sub_analysis_type);
 		end
 
@@ -54,9 +54,9 @@
 		inputs=add(inputs,'temperature',t_g,'doublevec',m.parameters.NumberOfDofsPerNode,m.parameters.NumberOfNodes);
 		displaystring(m.parameters.Verbose,'%s',['   update inputs']);
-		[m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters,inputs);
+		[m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters]=UpdateFromInputs(m.elements,m.nodes,m.vertices, loads,m.materials,m.parameters);
 	
 		%penalty constraints
 		displaystring(m.parameters.Verbose,'%s',['   penalty constraints']);
-		[loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints(m.elements,m.nodes,m.vertices,loads, m.materials,m.parameters,inputs,analysis_type,sub_analysis_type);
+		[loads,constraints_converged,num_unstable_constraints] =PenaltyConstraints(m.elements,m.nodes,m.vertices,loads, m.materials,m.parameters,analysis_type,sub_analysis_type);
 	
 		if ~converged,
Index: /issm/trunk/src/m/solutions/jpl/transient2d.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/transient2d.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/transient2d.m	(revision 3839)
@@ -71,5 +71,5 @@
 
 	%Get horizontal solution. 
-	rawresults=diagnostic_core(models,inputs);
+	rawresults=diagnostic_core(models);
 	solution(n+1).u_g=rawresults.u_g; solution(n+1).p_g=rawresults.p_g;
 
@@ -77,5 +77,5 @@
 	disp(sprintf('%s','   computing new thickness...'));
 	inputs=add(inputs,'velocity',solution(n+1).u_g,'doublevec',2,models.p.parameters.NumberOfNodes);
-	rawresults=prognostic_core(models,inputs,PrognosticAnalysisEnum(),NoneAnalysisEnum());
+	rawresults=prognostic_core(models,PrognosticAnalysisEnum(),NoneAnalysisEnum());
 	new_thickness=rawresults.h_g;
 
Index: /issm/trunk/src/m/solutions/jpl/transient3d.m
===================================================================
--- /issm/trunk/src/m/solutions/jpl/transient3d.m	(revision 3838)
+++ /issm/trunk/src/m/solutions/jpl/transient3d.m	(revision 3839)
@@ -84,5 +84,5 @@
 	%Deal with temperature first 
 	displaystring(md.verbose,'\n%s',['    computing temperatures...']);
-	[results(n+1).t_g models.t.loads melting_offset]=thermal_core_nonlinear(models.t,inputs,ThermalAnalysisEnum(),TransientAnalysisEnum());
+	[results(n+1).t_g models.t.loads melting_offset]=thermal_core_nonlinear(models.t,ThermalAnalysisEnum(),TransientAnalysisEnum());
 	inputs=add(inputs,'temperature',results(n+1).t_g,'doublevec',1,models.t.parameters.NumberOfNodes);
 	
@@ -96,5 +96,5 @@
 
 	%Deal with velocities.
-	rawresults=diagnostic_core(models,inputs);
+	rawresults=diagnostic_core(models);
 	results(n+1).u_g=rawresults.u_g; results(n+1).p_g=rawresults.p_g;
 
@@ -102,5 +102,5 @@
 	displaystring(md.verbose,'\n%s',['    computing new thickness...']);
 	inputs=add(inputs,'velocity',results(n+1).u_g,'doublevec',3,models.p.parameters.NumberOfNodes);
-	rawresults=prognostic_core(models,inputs,PrognosticAnalysisEnum(),NoneAnalysisEnum());
+	rawresults=prognostic_core(models,PrognosticAnalysisEnum(),NoneAnalysisEnum());
 	new_thickness=rawresults.h_g;