Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 22763) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 22764) @@ -2082,5 +2082,5 @@ /*Control Inputs*/ - if (control_analysis){ + if (control_analysis && !ad_analysis){ if(!ad_analysis)iomodel->FindConstant(&controls,NULL,"md.inversion.control_parameters"); if(ad_analysis)iomodel->FindConstant(&controls,NULL,"md.autodiff.independent_object_names"); Index: /issm/trunk-jpl/src/c/classes/ExternalResults/GenericExternalResult.h =================================================================== --- /issm/trunk-jpl/src/c/classes/ExternalResults/GenericExternalResult.h (revision 22763) +++ /issm/trunk-jpl/src/c/classes/ExternalResults/GenericExternalResult.h (revision 22764) @@ -179,4 +179,7 @@ } /*}}}*/ +void Transpose(void){ /*{{{*/ + _error_("not implemented yet"); +} /*}}}*/ char* GetResultName(void){ /*{{{*/ char* name = xNew(strlen(this->result_name)+1); @@ -594,183 +597,206 @@ } /*}}}*/ - -/*Specific instantiations for IssmDouble*: */ +template <> inline void GenericExternalResult::Transpose(void){/*{{{*/ + + + /*Perform transpose only if we have a matrix*/ + if(M>1 && N>1){ + IssmPDouble* temp=xNew(M*N); + for(int i=0;i(this->value); + this->value = temp; + } + + /*Switch dimensions*/ + int temp2 = this->N; + this->N = this->M; + this->M = temp2; + + + +} /*}}}*/ + + /*Specific instantiations for IssmDouble*: */ #if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_) //We hook off this specific specialization when not running ADOLC, otherwise we get a redeclaration with the next specialization. -template <> inline void GenericExternalResult::WriteData(FILE* fid,bool io_gather){ /*{{{*/ - - int i; - int my_rank; - int type; - int rows,cols; - char *name = NULL; - IssmPDouble passiveDouble; - IssmPDouble* passiveDoubles; - - /*recover my_rank:*/ - my_rank=IssmComm::GetRank(); - - if(io_gather){ - /*we are gathering the data on cpu 0, don't write on other cpus: */ - if(my_rank) return; - } - - /*First write enum: */ - int length=(strlen(this->result_name)+1)*sizeof(char); - fwrite(&length,sizeof(int),1,fid); - fwrite(this->result_name,length,1,fid); - - /*Now write time and step: */ - passiveDouble=reCast(time); - fwrite(&passiveDouble,sizeof(IssmPDouble),1,fid); - fwrite(&step,sizeof(int),1,fid); - - /*writing a IssmDouble array, type is 3:*/ - type=3; - fwrite(&type,sizeof(int),1,fid); - rows=this->M; - fwrite(&rows,sizeof(int),1,fid); - cols=this->N; - fwrite(&cols,sizeof(int),1,fid); - - passiveDoubles=xNew(this->M*this->N); - for (i=0;iM*this->N;i++)passiveDoubles[i]=reCast(value[i]); - fwrite(passiveDoubles,cols*rows*sizeof(IssmPDouble),1,fid); - xDelete(passiveDoubles); - -} -/*}}}*/ -#endif - -/*Specifics instantiations for Vector*/ -template <> inline GenericExternalResult*>::GenericExternalResult(int in_id, int in_enum_type,Vector* in_values, int in_M,int in_N,int in_step,IssmDouble in_time){/*{{{*/ - _error_("instanciation not correct"); -} -/*}}}*/ -template <> inline GenericExternalResult*>::GenericExternalResult(int in_id, int in_enum_type,Vector* in_value,int in_step, IssmDouble in_time){ /*{{{*/ - id = in_id; - M = 0; - N = 0; - - /*Convert enum to name*/ - EnumToStringx(&this->result_name,in_enum_type); - - step = in_step; - time = in_time; - - value = in_value; -} /*}}}*/ -template <> inline GenericExternalResult*>::~GenericExternalResult(){ /*{{{*/ - xDelete(this->result_name); - delete value; -} /*}}}*/ -template <> inline void GenericExternalResult*>::Echo(void){ /*{{{*/ - - _printf_("GenericExternalResult*>:\n"); - this->GenericEcho(); - this->value->Echo(); - -} /*}}}*/ -template <> inline void GenericExternalResult*>::DeepEcho(void){ /*{{{*/ - - this->Echo(); - -} /*}}}*/ -template <> inline Object* GenericExternalResult*>::copy(void){ /*{{{*/ - return new GenericExternalResult*>(this->id,StringToEnumx(this->result_name),this->value,this->step,this->time); -} /*}}}*/ -#if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_) //We hook off this specific specialization when not running ADOLC, otherwise we get a redeclaration with the next specialization. -template <> inline void GenericExternalResult*>::WriteData(FILE* fid,bool io_gather){ /*{{{*/ - - char *name = NULL; - int length,rows,cols=1; - - if(!io_gather){ - _error_("not supported yet"); - } - - /*Serialize vector*/ - IssmPDouble* serialvalues = this->value->ToMPISerial(); - this->value->GetSize(&rows); - - if(IssmComm::GetRank()==0){ - /*First write name: */ - length=(strlen(this->result_name)+1)*sizeof(char); + template <> inline void GenericExternalResult::WriteData(FILE* fid,bool io_gather){ /*{{{*/ + + int i; + int my_rank; + int type; + int rows,cols; + char *name = NULL; + IssmPDouble passiveDouble; + IssmPDouble* passiveDoubles; + + /*recover my_rank:*/ + my_rank=IssmComm::GetRank(); + + if(io_gather){ + /*we are gathering the data on cpu 0, don't write on other cpus: */ + if(my_rank) return; + } + + /*First write enum: */ + int length=(strlen(this->result_name)+1)*sizeof(char); fwrite(&length,sizeof(int),1,fid); fwrite(this->result_name,length,1,fid); /*Now write time and step: */ - IssmPDouble passiveDouble=reCast(time); + passiveDouble=reCast(time); fwrite(&passiveDouble,sizeof(IssmPDouble),1,fid); fwrite(&step,sizeof(int),1,fid); /*writing a IssmDouble array, type is 3:*/ - int type=3; + type=3; fwrite(&type,sizeof(int),1,fid); + rows=this->M; fwrite(&rows,sizeof(int),1,fid); + cols=this->N; fwrite(&cols,sizeof(int),1,fid); - fwrite(serialvalues,cols*rows*sizeof(IssmPDouble),1,fid); - } - - /*Clean up*/ - xDelete(serialvalues); - -} -/*}}}*/ + + passiveDoubles=xNew(this->M*this->N); + for (i=0;iM*this->N;i++)passiveDoubles[i]=reCast(value[i]); + fwrite(passiveDoubles,cols*rows*sizeof(IssmPDouble),1,fid); + xDelete(passiveDoubles); + + } + /*}}}*/ #endif -template <> inline int GenericExternalResult*>::ObjectEnum(void){ /*{{{*/ - return NoneEnum; - /*???? FIXME*/ -} /*}}}*/ - -/*Specifics instantiations for Vector*/ -template <> inline void GenericExternalResult*>::WriteData(FILE* fid,bool io_gather){ /*{{{*/ - - int i; - char *name = NULL; - int length,rows,cols=1; - IssmDouble* serialvalues = NULL; - IssmPDouble* pserialvalues = NULL; - - if(!io_gather){ - _error_("not supported yet"); - } - - /*Serialize vector*/ - serialvalues = this->value->ToMPISerial(); - this->value->GetSize(&rows); - - pserialvalues=xNew(rows); - for(i=0;i(serialvalues[i]); - - if(IssmComm::GetRank()==0){ - /*First write name: */ - length=(strlen(this->result_name)+1)*sizeof(char); - fwrite(&length,sizeof(int),1,fid); - fwrite(this->result_name,length,1,fid); - - /*Now write time and step: */ - IssmPDouble passiveDouble=reCast(time); - fwrite(&passiveDouble,sizeof(IssmPDouble),1,fid); - fwrite(&step,sizeof(int),1,fid); - - /*writing a IssmDouble array, type is 3:*/ - int type=3; - fwrite(&type,sizeof(int),1,fid); - fwrite(&rows,sizeof(int),1,fid); - fwrite(&cols,sizeof(int),1,fid); - fwrite(pserialvalues,cols*rows*sizeof(IssmPDouble),1,fid); - } - - /*Clean up*/ - xDelete(pserialvalues); - xDelete(serialvalues); - -} -/*}}}*/ -template <> inline void GenericExternalResult*>::Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){/*{{{*/ - - _error_("GenericExternalResult instantiated for type Vector* called " << result_name << " not implemented yet"); - -} /*}}}*/ + + /*Specifics instantiations for Vector*/ + template <> inline GenericExternalResult*>::GenericExternalResult(int in_id, int in_enum_type,Vector* in_values, int in_M,int in_N,int in_step,IssmDouble in_time){/*{{{*/ + _error_("instanciation not correct"); + } + /*}}}*/ + template <> inline GenericExternalResult*>::GenericExternalResult(int in_id, int in_enum_type,Vector* in_value,int in_step, IssmDouble in_time){ /*{{{*/ + id = in_id; + M = 0; + N = 0; + + /*Convert enum to name*/ + EnumToStringx(&this->result_name,in_enum_type); + + step = in_step; + time = in_time; + + value = in_value; + } /*}}}*/ + template <> inline GenericExternalResult*>::~GenericExternalResult(){ /*{{{*/ + xDelete(this->result_name); + delete value; + } /*}}}*/ + template <> inline void GenericExternalResult*>::Echo(void){ /*{{{*/ + + _printf_("GenericExternalResult*>:\n"); + this->GenericEcho(); + this->value->Echo(); + + } /*}}}*/ + template <> inline void GenericExternalResult*>::DeepEcho(void){ /*{{{*/ + + this->Echo(); + + } /*}}}*/ + template <> inline Object* GenericExternalResult*>::copy(void){ /*{{{*/ + return new GenericExternalResult*>(this->id,StringToEnumx(this->result_name),this->value,this->step,this->time); + } /*}}}*/ +#if defined(_HAVE_ADOLC_) && !defined(_WRAPPERS_) //We hook off this specific specialization when not running ADOLC, otherwise we get a redeclaration with the next specialization. + template <> inline void GenericExternalResult*>::WriteData(FILE* fid,bool io_gather){ /*{{{*/ + + char *name = NULL; + int length,rows,cols=1; + + if(!io_gather){ + _error_("not supported yet"); + } + + /*Serialize vector*/ + IssmPDouble* serialvalues = this->value->ToMPISerial(); + this->value->GetSize(&rows); + + if(IssmComm::GetRank()==0){ + /*First write name: */ + length=(strlen(this->result_name)+1)*sizeof(char); + fwrite(&length,sizeof(int),1,fid); + fwrite(this->result_name,length,1,fid); + + /*Now write time and step: */ + IssmPDouble passiveDouble=reCast(time); + fwrite(&passiveDouble,sizeof(IssmPDouble),1,fid); + fwrite(&step,sizeof(int),1,fid); + + /*writing a IssmDouble array, type is 3:*/ + int type=3; + fwrite(&type,sizeof(int),1,fid); + fwrite(&rows,sizeof(int),1,fid); + fwrite(&cols,sizeof(int),1,fid); + fwrite(serialvalues,cols*rows*sizeof(IssmPDouble),1,fid); + } + + /*Clean up*/ + xDelete(serialvalues); + + } + /*}}}*/ +#endif + template <> inline int GenericExternalResult*>::ObjectEnum(void){ /*{{{*/ + return NoneEnum; + /*???? FIXME*/ + } /*}}}*/ + + /*Specifics instantiations for Vector*/ + template <> inline void GenericExternalResult*>::WriteData(FILE* fid,bool io_gather){ /*{{{*/ + + int i; + char *name = NULL; + int length,rows,cols=1; + IssmDouble* serialvalues = NULL; + IssmPDouble* pserialvalues = NULL; + + if(!io_gather){ + _error_("not supported yet"); + } + + /*Serialize vector*/ + serialvalues = this->value->ToMPISerial(); + this->value->GetSize(&rows); + + pserialvalues=xNew(rows); + for(i=0;i(serialvalues[i]); + + if(IssmComm::GetRank()==0){ + /*First write name: */ + length=(strlen(this->result_name)+1)*sizeof(char); + fwrite(&length,sizeof(int),1,fid); + fwrite(this->result_name,length,1,fid); + + /*Now write time and step: */ + IssmPDouble passiveDouble=reCast(time); + fwrite(&passiveDouble,sizeof(IssmPDouble),1,fid); + fwrite(&step,sizeof(int),1,fid); + + /*writing a IssmDouble array, type is 3:*/ + int type=3; + fwrite(&type,sizeof(int),1,fid); + fwrite(&rows,sizeof(int),1,fid); + fwrite(&cols,sizeof(int),1,fid); + fwrite(pserialvalues,cols*rows*sizeof(IssmPDouble),1,fid); + } + + /*Clean up*/ + xDelete(pserialvalues); + xDelete(serialvalues); + + } + /*}}}*/ + template <> inline void GenericExternalResult*>::Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){/*{{{*/ + + _error_("GenericExternalResult instantiated for type Vector* called " << result_name << " not implemented yet"); + + } /*}}}*/ #endif /* _EXTERNAL_RESULTOBJECT_H */ Index: /issm/trunk-jpl/src/c/cores/controladm1qn3_core.cpp =================================================================== --- /issm/trunk-jpl/src/c/cores/controladm1qn3_core.cpp (revision 22763) +++ /issm/trunk-jpl/src/c/cores/controladm1qn3_core.cpp (revision 22764) @@ -71,4 +71,5 @@ int* N = NULL; int N_add = 0; + int* control_enum = NULL; if (solution_type == TransientSolutionEnum){ @@ -88,4 +89,5 @@ femmodel->parameters->FindParam(&scaling_factors,NULL,InversionControlScalingFactorsEnum); femmodel->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); + femmodel->parameters->FindParam(&control_enum,NULL,InversionControlParametersEnum); numberofvertices=femmodel->vertices->NumberOfVertices(); @@ -348,4 +350,5 @@ int offset = 0; int N_add; + int* control_enum; /*Recover some parameters*/ @@ -358,4 +361,5 @@ femmodel->parameters->FindParamAndMakePassive(>tol,InversionGttolEnum); femmodel->parameters->FindParam(&scaling_factors,NULL,InversionControlScalingFactorsEnum); + femmodel->parameters->FindParam(&control_enum,NULL,InversionControlParametersEnum); femmodel->parameters->SetParam(false,SaveResultsEnum); femmodel->parameters->FindParam(&N,NULL,ControlInputSizeNEnum); @@ -478,8 +482,31 @@ int step = 1; femmodel->parameters->SetParam(step,StepEnum); - femmodel->OutputControlsx(&femmodel->results); femmodel->results->AddObject(new GenericExternalResult(femmodel->results->Size()+1,JEnum,mystruct.Jlist,(*mystruct.i),mystruct.N,1,0)); - femmodel->results->AddObject(new GenericExternalResult(femmodel->results->Size()+1,Gradient1Enum,G,numberofvertices,intn/numberofvertices,1,0)); + + /*TEMP START*/ + int offset = 0; + for(int i=0;i* G_output = new GenericExternalResult(femmodel->results->Size()+1,Gradient1Enum+i,&G[offset],N[i],numberofvertices,1,0.); + GenericExternalResult* X_output = new GenericExternalResult(femmodel->results->Size()+1,control_enum[i],&X[offset],N[i],numberofvertices,1,0.); + + /*transpose for consistency with MATLAB's formating*/ + G_output->Transpose(); + X_output->Transpose(); + + /*Add to results*/ + femmodel->results->AddObject(G_output); + femmodel->results->AddObject(X_output); + + offset += N[i]*numberofvertices; + } + /*TEMP END*/ + + femmodel->results->AddObject(new GenericExternalResult(femmodel->results->Size()+1,GradientEnum,G,numberofvertices,intn/numberofvertices,1,0)); + femmodel->results->AddObject(new GenericExternalResult(femmodel->results->Size()+1,InversionControlParametersEnum,X,numberofvertices,intn/numberofvertices,1,0)); + + //femmodel->OutputControlsx(&femmodel->results); } else{ Index: /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateOutputDefinitions.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateOutputDefinitions.cpp (revision 22763) +++ /issm/trunk-jpl/src/c/modules/ModelProcessorx/CreateOutputDefinitions.cpp (revision 22764) @@ -341,5 +341,4 @@ iomodel->FetchMultipleData(&cfsurfacelogvel_vyobs_string,&num_cfsurfacelogvels, "md.cfsurfacelogvel.vyobs_string"); iomodel->FetchMultipleData(&cfsurfacelogvel_weights,&cfsurfacelogvel_weights_M,&cfsurfacelogvel_weights_N,&num_cfsurfacelogvels, "md.cfsurfacelogvel.weights"); iomodel->FetchMultipleData(&cfsurfacelogvel_weightstring,&num_cfsurfacelogvels, "md.cfsurfacelogvel.weights_string"); - _printf_("Num with weight string: "<FetchMultipleData(&cfsurfacelogvel_datatime,&num_cfsurfacelogvels, "md.cfsurfacelogvel.datatime");