Index: /issm/trunk/src/c/parallel/control_core.cpp
===================================================================
--- /issm/trunk/src/c/parallel/control_core.cpp	(revision 2218)
+++ /issm/trunk/src/c/parallel/control_core.cpp	(revision 2219)
@@ -22,4 +22,5 @@
 	DataSet* diagnostic_results=NULL;
 	DataSet* gradjcompute_results=NULL;
+	DataSet* steadystate_results=NULL;
 	Vec     u_g=NULL;
 	double  search_scalar;
@@ -48,4 +49,5 @@
 	int analysis_type;
 	int sub_analysis_type;
+	int     control_steady;
 	int debug=0;
 	int convergence=0;
@@ -71,4 +73,5 @@
 	model->FindParam(&sub_analysis_type,"sub_analysis_type");
 	model->FindParam(&numberofnodes,"numberofnodes");
+	model->FindParam(&control_steady,"control_steady");
 	gsize=fem_model->nodes->NumberOfDofs();
 
@@ -121,4 +124,12 @@
 		ControlConstrainx(param_g,numberofnodes,mincontrolconstraint,maxcontrolconstraint,control_type);
 		_printf_("%s\n","      done.");
+
+		/*In case we are running a steady state control method, compute new temperature field using new parameter 
+		 * distribution: */
+		steadystate_results=new DataSet(ResultsEnum());
+		steadystate_core(steadystate_results,model,inputs);
+		VecFree(&t_g); steadystate_results->FindResult(&t_g,"t_g");
+		delete results_steadystate;
+		inputs->Add("temperature",t_g,1,numberofnodes);
 
 		_printf_("%s%i%s%g\n","      value of misfit J after optimization #",n+1,": ",J[n]);
Index: /issm/trunk/src/m/solutions/cielo/control_core.m
===================================================================
--- /issm/trunk/src/m/solutions/cielo/control_core.m	(revision 2218)
+++ /issm/trunk/src/m/solutions/cielo/control_core.m	(revision 2219)
@@ -75,4 +75,10 @@
 	param_g=ControlConstrain(param_g,model.parameters);
 
+	%In case we are running a steady state control method, compute new temperature field using new parameter distribution: 
+	if model.parameters.control_steady;
+		steadystate_results=steadystate_core(models,inputs); t_g=results.t_g; 
+		inputs=add(inputs,'temperature',t_g,'doublevec',1,,model.parameters.numberofnodes);
+	end
+
 	disp(['      value of misfit J after optimization #' num2str(n) ':' num2str(c(n).J)]);