Index: /issm/trunk/src/c/modules/ModelProcessorx/Control/CreateParametersControl.cpp =================================================================== --- /issm/trunk/src/c/modules/ModelProcessorx/Control/CreateParametersControl.cpp (revision 8600) +++ /issm/trunk/src/c/modules/ModelProcessorx/Control/CreateParametersControl.cpp (revision 8601) @@ -23,6 +23,7 @@ if(iomodel->control_analysis){ - /*How many controls?*/ + /*How many controls and how many responses?*/ parameters->AddObject(new IntParam(NumControlsEnum,iomodel->num_control_type)); + parameters->AddObject(new IntParam(NumResponsesEnum,iomodel->num_cm_responses)); /*What control type?*/ @@ -52,7 +53,7 @@ IoModelFetchData(&iomodel->maxiter,NULL,NULL,iomodel_handle,"maxiter"); - parameters->AddObject(new DoubleVecParam(CmResponsesEnum,iomodel->cm_responses,iomodel->nsteps)); + parameters->AddObject(new DoubleMatParam(OptScalEnum,iomodel->optscal,iomodel->nsteps,iomodel->num_control_type)); + parameters->AddObject(new DoubleMatParam(CmResponsesEnum,iomodel->cm_responses,iomodel->nsteps,iomodel->num_cm_responses)); parameters->AddObject(new DoubleVecParam(CmJumpEnum,iomodel->cm_jump,iomodel->nsteps)); - parameters->AddObject(new DoubleMatParam(OptScalEnum,iomodel->optscal,iomodel->nsteps,iomodel->num_control_type)); parameters->AddObject(new DoubleVecParam(MaxIterEnum,iomodel->maxiter,iomodel->nsteps)); Index: /issm/trunk/src/c/objects/Elements/Penta.cpp =================================================================== --- /issm/trunk/src/c/objects/Elements/Penta.cpp (revision 8600) +++ /issm/trunk/src/c/objects/Elements/Penta.cpp (revision 8601) @@ -6295,5 +6295,4 @@ name==VzObsEnum || name==TemperatureEnum || - name==CmResponseEnum || name==DragCoefficientEnum || name==GradientEnum || Index: /issm/trunk/src/c/objects/Elements/Tria.cpp =================================================================== --- /issm/trunk/src/c/objects/Elements/Tria.cpp (revision 8600) +++ /issm/trunk/src/c/objects/Elements/Tria.cpp (revision 8601) @@ -1650,5 +1650,7 @@ /*Intermediaries */ - int i,ig,response; + int i,ig; + int *responses=NULL; + int response; double Jdet; double obs_velocity_mag,velocity_mag; @@ -1667,4 +1669,5 @@ this->parameters->FindParam(&meanvel,MeanVelEnum); this->parameters->FindParam(&epsvel,EpsVelEnum); + this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum); response=responses[0]; Input* weights_input=inputs->GetInput(WeightsEnum); _assert_(weights_input); Input* vx_input =inputs->GetInput(VxEnum); _assert_(vx_input); @@ -1673,5 +1676,4 @@ Input* vyobs_input =inputs->GetInput(VyObsEnum); _assert_(vyobs_input); - inputs->GetParameterValue(&response,CmResponseEnum); if(response==SurfaceAverageVelMisfitEnum) inputs->GetParameterValue(&S,SurfaceAreaEnum); @@ -1796,4 +1798,5 @@ /*Clean up and return*/ delete gauss; + xfree((void**)&responses); return pe; } @@ -1803,5 +1806,7 @@ /*Intermediaries */ - int i,ig,response; + int i,ig; + int *responses=NULL; + int response; double Jdet; double obs_velocity_mag,velocity_mag; @@ -1820,4 +1825,5 @@ this->parameters->FindParam(&meanvel,MeanVelEnum); this->parameters->FindParam(&epsvel,EpsVelEnum); + this->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum); response=responses[0]; Input* weights_input=inputs->GetInput(WeightsEnum); _assert_(weights_input); Input* vx_input =inputs->GetInput(VxEnum); _assert_(vx_input); @@ -1826,5 +1832,4 @@ Input* vyobs_input =inputs->GetInput(VyObsEnum); _assert_(vyobs_input); - inputs->GetParameterValue(&response,CmResponseEnum); if(response==SurfaceAverageVelMisfitEnum) inputs->GetParameterValue(&S,SurfaceAreaEnum); @@ -1949,4 +1954,5 @@ /*Clean up and return*/ delete gauss; + xfree((void**)&responses); return pe; } @@ -4103,5 +4109,4 @@ name==VxObsEnum || name==VyObsEnum || - name==CmResponseEnum || name==DragCoefficientEnum || name==GradientEnum || Index: /issm/trunk/src/c/objects/OptArgs.h =================================================================== --- /issm/trunk/src/c/objects/OptArgs.h (revision 8600) +++ /issm/trunk/src/c/objects/OptArgs.h (revision 8601) @@ -9,17 +9,14 @@ #include "mex.h" - struct OptArgs{ - char* function_name; mxArray* femmodel; - mxArray* n; +}; -}; #else + class Model; struct OptArgs{ FemModel* femmodel; - int n; }; #endif Index: /issm/trunk/src/c/objects/objects.h =================================================================== --- /issm/trunk/src/c/objects/objects.h (revision 8600) +++ /issm/trunk/src/c/objects/objects.h (revision 8601) @@ -125,4 +125,5 @@ #include "./Params/IntParam.h" #include "./Params/IntVecParam.h" +#include "./Params/IntMatParam.h" #include "./Params/FileParam.h" #include "./Params/Param.h" Index: /issm/trunk/src/c/shared/Numerics/OptFunc.cpp =================================================================== --- /issm/trunk/src/c/shared/Numerics/OptFunc.cpp (revision 8600) +++ /issm/trunk/src/c/shared/Numerics/OptFunc.cpp (revision 8601) @@ -20,5 +20,5 @@ double J; - mxArray* inputs[3]; + mxArray* inputs[2]; mxArray* psearch_scalar=NULL; mxArray* mxJ=NULL; @@ -27,7 +27,6 @@ inputs[0]=psearch_scalar; inputs[1]=optargs->femmodel; - inputs[2]=optargs->n; - mexCallMATLAB( 1, &mxJ, 3, (mxArray**)inputs, optargs->function_name); + mexCallMATLAB(1,&mxJ,2,(mxArray**)inputs, optargs->function_name); /*extract misfit from mxArray*/ Index: /issm/trunk/src/c/solutions/control_core.cpp =================================================================== --- /issm/trunk/src/c/solutions/control_core.cpp (revision 8600) +++ /issm/trunk/src/c/solutions/control_core.cpp (revision 8601) @@ -17,5 +17,5 @@ /*parameters: */ - int num_controls; + int num_controls,num_responses; int nsteps; double eps_cm; @@ -29,7 +29,7 @@ int* control_type = NULL; double* responses=NULL; + int* step_responses=NULL; double* maxiter=NULL; double* cm_jump=NULL; - /*intermediary: */ @@ -46,8 +46,9 @@ /*Recover parameters used throughout the solution:{{{1*/ + femmodel->parameters->FindParam(&num_controls,NumControlsEnum); + femmodel->parameters->FindParam(&num_responses,NumResponsesEnum); + femmodel->parameters->FindParam(&control_type,NULL,ControlTypeEnum); + femmodel->parameters->FindParam(&responses,NULL,NULL,CmResponsesEnum); femmodel->parameters->FindParam(&nsteps,NStepsEnum); - femmodel->parameters->FindParam(&num_controls,NumControlsEnum); - femmodel->parameters->FindParam(&control_type,NULL,ControlTypeEnum); - femmodel->parameters->FindParam(&responses,NULL,CmResponsesEnum); femmodel->parameters->FindParam(&maxiter,NULL,MaxIterEnum); femmodel->parameters->FindParam(&cm_jump,NULL,CmJumpEnum); @@ -71,4 +72,5 @@ /*Initialize responses: */ J=(double*)xmalloc(nsteps*sizeof(double)); + step_responses=(int*)xmalloc(num_responses*sizeof(double)); /*Initialize some of the BrentSearch arguments: */ @@ -81,5 +83,6 @@ /*Display info*/ _printf_(VerboseControl(),"\n%s%i%s%i\n"," control method step ",n+1,"/",nsteps); - InputUpdateFromConstantx(femmodel->elements,femmodel->nodes, femmodel->vertices, femmodel->loads, femmodel->materials, femmodel->parameters,(int)responses[n],CmResponseEnum); + for(i=0;iparameters->SetParam(step_responses,1,num_responses,StepResponsesEnum); /*In case we are running a steady state control method, compute new temperature field using new parameter distribution: */ @@ -98,5 +101,5 @@ _printf_(VerboseControl(),"%s\n"," optimizing along gradient direction"); - optargs.n=n; optpars.maxiter=(int)maxiter[n]; optpars.cm_jump=cm_jump[n]; + optpars.maxiter=(int)maxiter[n]; optpars.cm_jump=cm_jump[n]; BrentSearch(&search_scalar,J+n,&optpars,&objectivefunctionC,&optargs); //OptimalSearch(&search_scalar,J+n,&optpars,&objectivefunctionC,&optargs); @@ -129,4 +132,5 @@ xfree((void**)&control_type); xfree((void**)&responses); + xfree((void**)&step_responses); xfree((void**)&maxiter); xfree((void**)&cm_jump); Index: /issm/trunk/src/c/solutions/objectivefunctionC.cpp =================================================================== --- /issm/trunk/src/c/solutions/objectivefunctionC.cpp (revision 8600) +++ /issm/trunk/src/c/solutions/objectivefunctionC.cpp (revision 8601) @@ -25,14 +25,13 @@ /*output: */ - double J; + double J=0,Jplus; /*parameters: */ - FemModel *femmodel = NULL; - int n; - double *responses = NULL; - int solution_type; - int analysis_type; + int num_responses; + int solution_type,analysis_type; bool isstokes = false; bool conserve_loads = true; + int *responses = NULL; + FemModel *femmodel = NULL; /*Recover finite element model: */ @@ -40,7 +39,6 @@ /*Recover parameters: */ - n=optargs->n; - - femmodel->parameters->FindParam(&responses,NULL,CmResponsesEnum); + femmodel->parameters->FindParam(&num_responses,NumResponsesEnum); + femmodel->parameters->FindParam(&responses,NULL,NULL,StepResponsesEnum); femmodel->parameters->FindParam(&isstokes,IsStokesEnum); femmodel->parameters->FindParam(&analysis_type,AnalysisTypeEnum); @@ -76,5 +74,8 @@ /*Compute misfit for this velocity field.*/ - CostFunctionx( &J, femmodel->elements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters,(int)responses[n]); + for(int i=0;ielements,femmodel->nodes, femmodel->vertices,femmodel->loads, femmodel->materials, femmodel->parameters,responses[i]); + J+=Jplus; + } /*Free ressources:*/