Index: /issm/trunk-jpl/src/c/classes/Cfsurfacelogvel.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Cfsurfacelogvel.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Cfsurfacelogvel.cpp (revision 25406) @@ -101,32 +101,27 @@ IssmDouble Cfsurfacelogvel::Response(FemModel* femmodel){/*{{{*/ - /*diverse: */ - IssmDouble time; - - /*recover time parameters: */ - femmodel->parameters->FindParam(&time,TimeEnum); - if(time < last_time) timepassedflag = false; - last_time = time; - - int i; - IssmDouble J=0.; - IssmDouble J_sum=0.; - - if(datatime<=time && !timepassedflag){ - for(i=0;ielements->Size();i++){ - Element* element=(Element*)femmodel->elements->GetObjectByOffset(i); - J+=this->Cfsurfacelogvel_Calculation(element,definitionenum); - } - - ISSM_MPI_Allreduce ( (void*)&J,(void*)&J_sum,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm()); - ISSM_MPI_Bcast(&J_sum,1,ISSM_MPI_DOUBLE,0,IssmComm::GetComm()); - J=J_sum; - - timepassedflag = true; - return J; - } - else return J; - } - /*}}}*/ + + /*recover time parameters: */ + IssmDouble time; + femmodel->parameters->FindParam(&time,TimeEnum); + if(time < last_time) this->timepassedflag = false; //FIXME: do we really need this? looks like an old piece of code before copy + this->last_time = time; + + IssmDouble J=0.; + IssmDouble J_sum=0.; + if(this->datatime<=time && !timepassedflag){ + for(int i=0;ielements->Size();i++){ + Element* element=(Element*)femmodel->elements->GetObjectByOffset(i); + J+=this->Cfsurfacelogvel_Calculation(element,definitionenum); + } + + ISSM_MPI_Allreduce ( (void*)&J,(void*)&J_sum,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm()); + ISSM_MPI_Bcast(&J_sum,1,ISSM_MPI_DOUBLE,0,IssmComm::GetComm()); + J=J_sum; + + this->timepassedflag = true; + } + return J; +}/*}}}*/ IssmDouble Cfsurfacelogvel::Cfsurfacelogvel_Calculation(Element* element, int definitionenum){/*{{{*/ @@ -165,5 +160,5 @@ Input *vy_input = NULL; if(numcomponents==2){ - vy_input = topelement->GetInput(VyEnum); _assert_(vy_input); + vy_input = topelement->GetInput(VyEnum); _assert_(vy_input); } @@ -192,7 +187,7 @@ * * [ vel + eps ] 2 * * J = 4 \bar{v}^2 | log ( ----------- ) | - * * [ vel + eps ] - * * obs - * */ + * [ vel + eps ] + * obs + */ if(numcomponents==1){ velocity_mag =fabs(vx)+epsvel; @@ -205,4 +200,6 @@ misfit=4*pow(meanvel,2)*pow(log(velocity_mag/obs_velocity_mag),2); + _printf_(velocity_mag<<" "<GetNumberOfVertices(); + IssmDouble value,value_min,value_max; int *vertexids = xNew(numvertices); @@ -1779,5 +1781,35 @@ this->AddControlInput(input_enum,inputs,iomodel,values,values_min,values_max,P0Enum,id); } - + else if(M==iomodel->numberofelements+1 && N>1){ + + /*create transient input: */ + IssmDouble* times = xNew(N); + for(int t=0;tSetTransientControlInput(input_enum,id,times,N); + + /*Get input*/ + ControlInput* control_input = inputs->GetControlInput(input_enum); _assert_(control_input); + TransientInput* transientinput_val = control_input->GetTransientInput("value"); + TransientInput* transientinput_min = control_input->GetTransientInput("lowerbound"); + TransientInput* transientinput_max = control_input->GetTransientInput("upperbound"); + for(int t=0;tSid()+t]; + value_min = scale*min_vector[N*this->Sid()+t]; + value_max = scale*max_vector[N*this->Sid()+t]; + switch(this->ObjectEnum()){ + case TriaEnum: + transientinput_val->AddTriaTimeInput( t,1,&(this->lid),&value,P0Enum); + transientinput_min->AddTriaTimeInput( t,1,&(this->lid),&value_min,P0Enum); + transientinput_max->AddTriaTimeInput( t,1,&(this->lid),&value_max,P0Enum); + break; + case PentaEnum: + //transientinput->AddPentaTimeInput(t,1,&(this->lid),&value,P0Enum); + _error_("to be implemented"); + break; + default: _error_("Not implemented yet"); + } + } + xDelete(times); + } else if(M==iomodel->numberofvertices+1 && N>1){ _error_("not supported tet"); Index: /issm/trunk-jpl/src/c/classes/Elements/Element.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Element.h (revision 25406) @@ -259,6 +259,5 @@ virtual int GetNumberOfNodes(int enum_type)=0; virtual int GetNumberOfVertices(void)=0; - virtual void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data,int offset)=0; - virtual void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data)=0; + virtual void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,int N,const char* data,int offset)=0; virtual void GetVerticesCoordinatesBase(IssmDouble** xyz_list)=0; virtual void GetVerticesCoordinatesTop(IssmDouble** xyz_list)=0; Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.cpp (revision 25406) @@ -215,8 +215,12 @@ } + /*Create Control Input*/ + inputs->SetControlInput(input_enum,PentaInputEnum,interpolation_enum,id); + ControlInput* control_input = inputs->GetControlInput(input_enum); _assert_(input_enum); + /*Call inputs method*/ switch(interpolation_enum){ case P1Enum: - inputs->SetPentaControlInput(input_enum,PentaInputEnum,interpolation_enum,id,NUMVERTICES,vertexlids,values,values_min,values_max); + control_input->SetControl(interpolation_enum,NUMVERTICES,&vertexlids[0],values,values_min,values_max); break; default: @@ -1830,7 +1834,6 @@ } /*}}}*/ -void Penta::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data){/*{{{*/ - - _error_("NOT NEEDED ANYMORE"); +void Penta::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,int N,const char* data,int offset){/*{{{*/ + /*Get input*/ if(control_enum==MaterialsRheologyBbarEnum) control_enum=MaterialsRheologyBEnum; @@ -1878,72 +1881,5 @@ for(int t=0;tnumtimesteps;t++) { IssmDouble time = transientinput->GetTimeByOffset(t); - _error_("not implemented"); - //PentaInput* timeinput = xDynamicCast(transientinput->GetTimeInput(time)); - //if(timeinput->interpolation_type!=P1Enum) _error_("not supported yet"); - //input->Serve(NUMVERTICES,&lidlist[0]); - ///*Create list of indices and values for global vector*/ - //for(int i=0;ivertices[i]->Sid()+t*M[control_index]; - // values[N*i+t] = timeinput->values[i]; - //} - } - vector->SetValues(NUMVERTICES*transientinput->numtimesteps,idlist,values,INS_VAL); - xDelete(M); - xDelete(idlist); - xDelete(values); - break; - } - default: _error_("input "<ObjectEnum())<<" not supported yet"); - } -} -/*}}}*/ -void Penta::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data,int offset){/*{{{*/ - - /*Get input*/ - if(control_enum==MaterialsRheologyBbarEnum) control_enum=MaterialsRheologyBEnum; - ElementInput* input=this->inputs->GetControlInputData(control_enum,data); _assert_(input); - - /*Lid list once for all*/ - int lidlist[NUMVERTICES]; - for(int i=0;ilid; - - /*Check what input we are dealing with*/ - switch(input->ObjectEnum()){ - case PentaInputEnum: - { - IssmDouble values[NUMVERTICES]; - int idlist[NUMVERTICES]; - - PentaInput* triainput = xDynamicCast(input); - - /*Create list of indices and values for global vector*/ - GradientIndexing(&idlist[0],control_index); - - if(triainput->GetInputInterpolationType()==P1Enum){ - input->Serve(NUMVERTICES,&lidlist[0]); - for(int i=0;ielement_values[i]; - vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); - } - else if(triainput->GetInputInterpolationType()==P0Enum){ - input->Serve(1,&this->lid); - vector->SetValue(idlist[0],triainput->element_values[0],INS_VAL); - } - else{ - _error_("not supported yet"); - } - break; - } - - case TransientInputEnum: - { - TransientInput* transientinput = xDynamicCast(input); - int N = transientinput->numtimesteps; - int* M = NULL; - parameters->FindParam(&M,NULL,ControlInputSizeMEnum); - int* idlist = xNew(NUMVERTICES*N); - IssmDouble* values = xNew(NUMVERTICES*N); - for(int t=0;tnumtimesteps;t++) { - IssmDouble time = transientinput->GetTimeByOffset(t); - _error_("not implemented"); + _error_("not implemented, SEE TRIA!"); //PentaInput* timeinput = xDynamicCast(transientinput->GetTimeInput(time)); //if(timeinput->interpolation_type!=P1Enum) _error_("not supported yet"); Index: /issm/trunk-jpl/src/c/classes/Elements/Penta.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Penta.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Penta.h (revision 25406) @@ -99,6 +99,5 @@ Penta* GetLowerPenta(void); Penta* GetUpperPenta(void); - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data); - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data,int offset); + void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,int N,const char* data,int offset); int GetVertexIndex(Vertex* vertex); void GetVerticesCoordinatesBase(IssmDouble** pxyz_list); Index: /issm/trunk-jpl/src/c/classes/Elements/Seg.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Seg.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Seg.h (revision 25406) @@ -73,6 +73,5 @@ int GetNumberOfNodes(int enum_type){_error_("not implemented yet");}; int GetNumberOfVertices(void); - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data,int offset){_error_("not implemented yet");}; - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data){_error_("not implemented yet");}; + void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,int N,const char* data,int offset){_error_("not implemented yet");}; void GetVerticesCoordinates(IssmDouble** pxyz_list); void GetVerticesCoordinatesBase(IssmDouble** pxyz_list){_error_("not implemented yet");}; Index: /issm/trunk-jpl/src/c/classes/Elements/Tetra.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tetra.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Tetra.h (revision 25406) @@ -77,6 +77,5 @@ int GetNumberOfNodes(int enum_type){_error_("not implemented yet");}; int GetNumberOfVertices(void); - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data,int offset){_error_("not implemented yet");}; - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data){_error_("not implemented yet");}; + void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,int N,const char* data,int offset){_error_("not implemented yet");}; void GetVerticesCoordinatesBase(IssmDouble** pxyz_list); void GetVerticesCoordinatesTop(IssmDouble** pxyz_list); Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.cpp (revision 25406) @@ -220,11 +220,15 @@ } + /*Create Control Input*/ + inputs->SetControlInput(input_enum,TriaInputEnum,interpolation_enum,id); + ControlInput* control_input = inputs->GetControlInput(input_enum); _assert_(input_enum); + /*Call inputs method*/ switch(interpolation_enum){ case P0Enum: - inputs->SetTriaControlInput(input_enum,TriaInputEnum,interpolation_enum,id,1,&this->lid,values,values_min,values_max); + control_input->SetControl(interpolation_enum,1,&this->lid,values,values_min,values_max); break; case P1Enum: - inputs->SetTriaControlInput(input_enum,TriaInputEnum,interpolation_enum,id,NUMVERTICES,vertexlids,values,values_min,values_max); + control_input->SetControl(interpolation_enum,NUMVERTICES,&vertexlids[0],values,values_min,values_max); break; default: @@ -1104,31 +1108,51 @@ IssmDouble values[NUMVERTICES]; int lidlist[NUMVERTICES]; - int idlist[NUMVERTICES]; - - ElementInput* input=this->inputs->GetControlInputData(control_enum,"gradient"); _assert_(input); + + /*Get list of ids for this element and this control*/ + int* idlist = xNew(NUMVERTICES*N); + GradientIndexing(&idlist[0],control_index); + + ControlInput* control_input=this->inputs->GetControlInput(control_enum); _assert_(control_input); this->GetVerticesLidList(&lidlist[0]); - GradientIndexing(&idlist[0],control_index); /*Get values on vertices*/ - if(input->ObjectEnum()==TriaInputEnum && input->GetInputInterpolationType()==P1Enum){ - _assert_(N==1); - for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); - } - else if(input->ObjectEnum()==TriaInputEnum && input->GetInputInterpolationType()==P0Enum){ - _assert_(N==1); - input->SetInput(P0Enum,this->lid,gradient[idlist[0]]); - } - else if(input->ObjectEnum()==TransientInputEnum){ + if(control_input->layout_enum==TriaInputEnum){ + ElementInput* gradient_input = control_input->GetInput("gradient"); _assert_(gradient_input); + if(gradient_input->GetInputInterpolationType()==P1Enum){ + _assert_(N==1); + for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); + } + else if(gradient_input->GetInputInterpolationType()==P0Enum){ + _assert_(N==1); + gradient_input->SetInput(P0Enum,this->lid,gradient[idlist[0]]); + } + else{ + _error_("not implemented yet"); + } + } + else if(control_input->layout_enum==TransientInputEnum){ + _assert_(N>1); + TransientInput* gradient_input = control_input->GetTransientInput("gradient"); _assert_(gradient_input); + for(int n=0;nSetInput(new_input,n); - //controlinput->Configure(parameters); + TriaInput* input_n = gradient_input->GetTriaInput(n); _assert_(input_n); + if(input_n->GetInputInterpolationType()==P1Enum){ + _error_("not implemented"); + } + else if(input_n->GetInputInterpolationType()==P0Enum){ + input_n->SetInput(P0Enum,this->lid,gradient[idlist[n]]); + } + else{ + _error_("not implemented yet"); + } } } else _error_("Type not supported"); + + /*Clean up*/ + xDelete(idlist); }/*}}}*/ @@ -2407,117 +2431,61 @@ } /*}}}*/ -void Tria::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data){/*{{{*/ - - _error_("DO WE NEED THIS FUNCTION?? REALLY?? see below"); - - /*Get out if this is not an element input*/ - if(!IsInputEnum(control_enum)) _error_("Enum "<inputs->GetControlInputData(control_enum,data); _assert_(input); - - /*Intermediaries*/ - int numindices; - int indices[NUMVERTICES]; - - /*Check interpolation*/ - int interpolation = input->GetInterpolation(); - switch(interpolation){ - case P0Enum: +void Tria::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,int N,const char* data,int offset){/*{{{*/ + + /*Get list of ids for this element and this control*/ + int* idlist = xNew(NUMVERTICES*N); + int lidlist[NUMVERTICES]; + this->GradientIndexing(&idlist[0],control_index); + this->GetVerticesLidList(&lidlist[0]); + + /*Get Control*/ + ControlInput* control_input=this->inputs->GetControlInput(control_enum); _assert_(control_input); + + /*Get values on vertices*/ + if(control_input->layout_enum==TriaInputEnum){ + _assert_(N==1); + TriaInput* input = xDynamicCast(control_input->GetInput("value")); _assert_(input); + + if(input->GetInputInterpolationType()==P1Enum){ + IssmDouble values[NUMVERTICES]; + input->Serve(NUMVERTICES,&lidlist[0]); + for(int i=0;ielement_values[i]; + vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); + } + else if(input->GetInputInterpolationType()==P0Enum){ input->Serve(1,&this->lid); - break; - case P1Enum: - numindices = NUMVERTICES; - for(int i=0;ilid; - input->Serve(numindices,&indices[0]); - break; - default: _error_("interpolation "<iscollapsed){ - xDynamicCast(input)->SetServeCollapsed(true); - } - - /* Start looping on the number of vertices: */ - IssmDouble values[NUMVERTICES]; - Gauss*gauss=this->NewGauss(); - for(int iv=0;ivGaussVertex(iv); - input->GetInputValue(&values[iv],gauss); - } - delete gauss; - - vector->SetValues(NUMVERTICES,idlist,&values[0],INS_VAL); -} -/*}}}*/ -void Tria::GetVectorFromControlInputs(Vector* vector,int control_enum,int control_index,const char* data,int offset){/*{{{*/ - - /*Get input*/ - ElementInput* input=this->inputs->GetControlInputData(control_enum,data); _assert_(input); - - /*Lid list once for all*/ - int lidlist[NUMVERTICES]; - for(int i=0;ilid; - - /*Check what input we are dealing with*/ - switch(input->ObjectEnum()){ - case TriaInputEnum: - { - IssmDouble values[NUMVERTICES]; - int idlist[NUMVERTICES]; - - TriaInput* triainput = xDynamicCast(input); - - /*Create list of indices and values for global vector*/ - GradientIndexing(&idlist[0],control_index); - - if(triainput->GetInputInterpolationType()==P1Enum){ - input->Serve(NUMVERTICES,&lidlist[0]); - for(int i=0;ielement_values[i]; - vector->SetValues(NUMVERTICES,idlist,values,INS_VAL); - } - else if(triainput->GetInputInterpolationType()==P0Enum){ - input->Serve(1,&this->lid); - vector->SetValue(idlist[0],triainput->element_values[0],INS_VAL); - } - else{ - _error_("not supported yet"); - } - break; - } - - case TransientInputEnum: - { - TransientInput* transientinput = xDynamicCast(input); - int N = transientinput->numtimesteps; - int* M = NULL; - parameters->FindParam(&M,NULL,ControlInputSizeMEnum); - int* idlist = xNew(NUMVERTICES*N); - IssmDouble* values = xNew(NUMVERTICES*N); - for(int t=0;tnumtimesteps;t++) { - IssmDouble time = transientinput->GetTimeByOffset(t); - _error_("not implemented"); - //TriaInput* timeinput = xDynamicCast(transientinput->GetTimeInput(time)); - //if(timeinput->interpolation_type!=P1Enum) _error_("not supported yet"); - //input->Serve(NUMVERTICES,&lidlist[0]); - ///*Create list of indices and values for global vector*/ - //for(int i=0;ivertices[i]->Sid()+t*M[control_index]; - // values[N*i+t] = timeinput->values[i]; - //} - } - vector->SetValues(NUMVERTICES*transientinput->numtimesteps,idlist,values,INS_VAL); - xDelete(M); - xDelete(idlist); - xDelete(values); - break; - } - default: _error_("input "<ObjectEnum())<<" not supported yet"); - } + vector->SetValue(idlist[0],input->element_values[0],INS_VAL); + } + else{ + _error_("not supported yet"); + } + } + else if(control_input->layout_enum==TransientInputEnum){ + TransientInput* input = control_input->GetTransientInput("value"); _assert_(input); + _assert_(N==input->numtimesteps); + IssmDouble* values = xNew(NUMVERTICES*N); + + int count = 0; + for(int n=0;nGetTriaInput(n); _assert_(input_n); + if(input_n->GetInputInterpolationType()==P1Enum){ + _error_("not implemented"); + } + else if(input_n->GetInputInterpolationType()==P0Enum){ + input_n->Serve(1,&this->lid); + values[n] = input_n->element_values[0]; + count++; + } + else{ + _error_("not implemented yet"); + } + } + vector->SetValues(count,idlist,values,INS_VAL); + xDelete(values); + } + else _error_("Type not supported"); + + /*Clean up*/ + xDelete(idlist); } /*}}}*/ @@ -3944,5 +3912,4 @@ IssmDouble values[NUMVERTICES]; int lidlist[NUMVERTICES]; - int idlist[NUMVERTICES]; /*Get Domain type*/ @@ -3966,29 +3933,49 @@ if(!IsInputEnum(control_enum)) return; - ElementInput* input=this->inputs->GetControlInputData(control_enum,"value"); _assert_(input); + /*Get list of ids for this element and this control*/ + int* idlist = xNew(NUMVERTICES*N); + GradientIndexing(&idlist[0],control_index); + + ControlInput* control_input=this->inputs->GetControlInput(control_enum); _assert_(control_input); this->GetVerticesLidList(&lidlist[0]); - GradientIndexing(&idlist[0],control_index); - - /*Get values on vertices*/ - if(input->ObjectEnum()==TriaInputEnum && input->GetInputInterpolationType()==P1Enum){ - _assert_(N==1); - for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); - } - else if(input->ObjectEnum()==TriaInputEnum && input->GetInputInterpolationType()==P0Enum){ - _assert_(N==1); - input->SetInput(P0Enum,this->lid,vector[idlist[0]]); - } - else if(input->ObjectEnum()==TransientInputEnum){ + + if(control_input->layout_enum==TriaInputEnum){ + ElementInput* input = control_input->GetInput("value"); _assert_(input); + if(input->GetInputInterpolationType()==P1Enum){ + _assert_(N==1); + for(int i=0;iSetInput(P1Enum,NUMVERTICES,&lidlist[0],&values[0]); + } + else if(input->GetInputInterpolationType()==P0Enum){ + _assert_(N==1); + input->SetInput(P0Enum,this->lid,vector[idlist[0]]); + } + else{ + _error_("not implemented yet"); + } + } + else if(control_input->layout_enum==TransientInputEnum){ + _assert_(N>1); + TransientInput* input = control_input->GetTransientInput("value"); _assert_(input); + for(int n=0;nSetInput(new_input,n); - //controlinput->Configure(parameters); + TriaInput* input_n = input->GetTriaInput(n); _assert_(input_n); + if(input_n->GetInputInterpolationType()==P1Enum){ + _error_("not implemented"); + } + else if(input_n->GetInputInterpolationType()==P0Enum){ + input_n->SetInput(P0Enum,this->lid,vector[idlist[n]]); + } + else{ + _error_("not implemented yet"); + } } } else _error_("Type not supported"); + + /*Clean up*/ + xDelete(idlist); } /*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Elements/Tria.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Elements/Tria.h (revision 25406) @@ -91,6 +91,5 @@ int GetNumberOfNodes(int enum_type); int GetNumberOfVertices(void); - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data,int offset); - void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,const char* data); + void GetVectorFromControlInputs(Vector* gradient,int control_enum,int control_index,int N,const char* data,int offset); void GetVerticesCoordinatesBase(IssmDouble** pxyz_list); void GetVerticesCoordinatesTop(IssmDouble** pxyz_list); Index: /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.cpp (revision 25406) @@ -13,4 +13,5 @@ #include "./TriaInput.h" #include "./PentaInput.h" +#include "./TransientInput.h" //#include "../../toolkits/objects/Vector.h" @@ -50,4 +51,17 @@ } /*}}}*/ +ControlInput::ControlInput(int enum_in,int nbe, int nbv,int id,IssmDouble* times, int numtimes){/*{{{*/ + + this->enum_type = enum_in; + this->control_id = id; + this->layout_enum = TransientInputEnum; /*Tria or Penta?*/ + + this->values =new TransientInput(enum_in,nbe,nbv,times,numtimes); + this->savedvalues=new TransientInput(enum_in,nbe,nbv,times,numtimes); + this->minvalues =new TransientInput(enum_in,nbe,nbv,times,numtimes); + this->maxvalues =new TransientInput(enum_in,nbe,nbv,times,numtimes); + this->gradient =new TransientInput(enum_in,nbe,nbv,times,numtimes); +} +/*}}}*/ ControlInput::~ControlInput(){/*{{{*/ if(values) delete values; @@ -67,13 +81,21 @@ output->enum_type=this->enum_type; output->control_id=this->control_id; - output->layout_enum = this->control_id; - - if(values) output->values = xDynamicCast(this->values->copy()); - if(savedvalues) output->savedvalues = xDynamicCast(this->savedvalues->copy()); - if(minvalues) output->minvalues = xDynamicCast(this->minvalues->copy()); - if(maxvalues) output->maxvalues = xDynamicCast(this->maxvalues->copy()); - if(gradient) output->gradient = xDynamicCast(this->gradient->copy()); + output->layout_enum = this->layout_enum; + + if(values) output->values = this->values->copy(); + if(savedvalues) output->savedvalues = this->savedvalues->copy(); + if(minvalues) output->minvalues = this->minvalues->copy(); + if(maxvalues) output->maxvalues = this->maxvalues->copy(); + if(gradient) output->gradient = this->gradient->copy(); return output; +} +/*}}}*/ +void ControlInput::Configure(Parameters* params){/*{{{*/ + this->values->Configure(params); + this->savedvalues->Configure(params); + this->minvalues->Configure(params); + this->maxvalues->Configure(params); + this->gradient->Configure(params); } /*}}}*/ @@ -82,5 +104,5 @@ _printf_("ControlInput:\n"); _printf_(setw(15)<<" ControlInput "<enum_type)<<"\n"); - _printf_(setw(15)<<" ControlInput "<layout_enum)<<"\n"); + _printf_(setw(15)<<" Layout "<layout_enum)<<"\n"); _printf_("---values: \n"); if (values) values->Echo(); _printf_("---savedvalues: \n");if (savedvalues) savedvalues->Echo(); @@ -114,8 +136,12 @@ /*Set input*/ - //TriaInput* input = xDynamicCast(this->inputs[id]); - this->values->SetInput(interp,numindices,indices,values_in); - this->minvalues->SetInput(interp,numindices,indices,values_min); - this->maxvalues->SetInput(interp,numindices,indices,values_max); + if(this->values->ObjectEnum()==TriaInputEnum || this->values->ObjectEnum()==PentaInputEnum){ + xDynamicCast(this->values)->SetInput(interp,numindices,indices,values_in); + xDynamicCast(this->minvalues)->SetInput(interp,numindices,indices,values_min); + xDynamicCast(this->maxvalues)->SetInput(interp,numindices,indices,values_max); + } + else{ + _error_("not supported"); + } } /*}}}*/ @@ -124,5 +150,10 @@ _assert_(this); _assert_(this->gradient); - this->gradient->SetInput(interp,numindices,indices,values_in); + if(this->gradient->ObjectEnum()==TriaInputEnum || this->gradient->ObjectEnum()==PentaInputEnum){ + xDynamicCast(this->gradient)->SetInput(interp,numindices,indices,values_in); + } + else{ + _error_("not supported"); + } } /*}}}*/ @@ -161,8 +192,13 @@ /*Cast and return*/ - if(this->values->ObjectEnum()!=TriaInputEnum){ + if(this->values->ObjectEnum()==TriaInputEnum){ + return xDynamicCast(this->values); + } + else if(this->values->ObjectEnum()==TransientInputEnum){ + return xDynamicCast(this->values)->GetTriaInput(); + } + else{ _error_("Cannot return a TriaInput"); } - return xDynamicCast(this->values); } @@ -180,28 +216,67 @@ ElementInput* ControlInput::GetInput(const char* data){/*{{{*/ - if(strcmp(data,"value")==0){ - _assert_(values); - return values; - } - else if(strcmp(data,"savedvalues")==0){ - _assert_(savedvalues); - return values; - } - else if (strcmp(data,"lowerbound")==0){ - _assert_(minvalues); - return minvalues; - } - else if (strcmp(data,"upperbound")==0){ - _assert_(maxvalues); - return maxvalues; - } - else if (strcmp(data,"gradient")==0){ - _assert_(gradient); - return gradient; - } - else{ - _error_("Data " << data << " not supported yet"); - } - -} -/*}}}*/ + if(this->gradient->ObjectEnum()==TriaInputEnum || this->gradient->ObjectEnum()==PentaInputEnum){ + if(strcmp(data,"value")==0){ + _assert_(values); + return xDynamicCast(values); + } + else if(strcmp(data,"savedvalues")==0){ + _assert_(savedvalues); + return xDynamicCast(values); + } + else if (strcmp(data,"lowerbound")==0){ + _assert_(minvalues); + return xDynamicCast(minvalues); + } + else if (strcmp(data,"upperbound")==0){ + _assert_(maxvalues); + return xDynamicCast(maxvalues); + } + else if (strcmp(data,"gradient")==0){ + _assert_(gradient); + return xDynamicCast(gradient); + } + else{ + _error_("Data " << data << " not supported yet"); + } + } + else{ + int* temp = xNew(3); + _error_("not supported"); + } + +} +/*}}}*/ +TransientInput* ControlInput::GetTransientInput(const char* data){/*{{{*/ + + if(this->values->ObjectEnum()==TransientInputEnum){ + if(strcmp(data,"value")==0){ + _assert_(values); + return xDynamicCast(values); + } + else if(strcmp(data,"savedvalues")==0){ + _assert_(savedvalues); + return xDynamicCast(values); + } + else if (strcmp(data,"lowerbound")==0){ + _assert_(minvalues); + return xDynamicCast(minvalues); + } + else if (strcmp(data,"upperbound")==0){ + _assert_(maxvalues); + return xDynamicCast(maxvalues); + } + else if (strcmp(data,"gradient")==0){ + _assert_(gradient); + return xDynamicCast(gradient); + } + else{ + _error_("Data " << data << " not supported yet"); + } + } + else{ + _error_("not supported"); + } + +} +/*}}}*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Inputs/ControlInput.h (revision 25406) @@ -10,33 +10,37 @@ class Gauss; class ElementInput; +class TransientInput; class ControlInput: public Input{ public: - int control_id; - int enum_type; - int layout_enum; - ElementInput *gradient; - ElementInput *maxvalues; - ElementInput *minvalues; - ElementInput *savedvalues; - ElementInput *values; + int control_id; + int enum_type; + int layout_enum; + Input *gradient; + Input *maxvalues; + Input *minvalues; + Input *savedvalues; + Input *values; /*ControlInput constructors, destructors: {{{*/ ControlInput(); ControlInput(int nbe, int nbv,int input_layout_enum,int interp,int id); + ControlInput(int enum_in,int nbe, int nbv,int id,IssmDouble* times, int numtimes); ~ControlInput(); /*}}}*/ /*Object virtual functions definitions:{{{ */ Input* copy(); - void DeepEcho(); - void Echo(); - int Id(); - void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction); - int ObjectEnum(); + void Configure(Parameters* params); + void DeepEcho(); + void Echo(); + int Id(); + void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction); + int ObjectEnum(); /*}}}*/ void SetInput(Input* in_input){_error_("not impelemented");}; void SetInput(Input* in_input,int timeoffset){_error_("not impelemented");}; ElementInput* GetInput(const char* data); + TransientInput* GetTransientInput(const char* data); void SetControl(int interp,int numindices,int* indices,IssmDouble* values_in,IssmDouble* values_min,IssmDouble* values_max); void SetGradient(int interp,int numindices,int* indices,IssmDouble* values_in); Index: /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Inputs/Inputs.cpp (revision 25406) @@ -690,7 +690,8 @@ } }/*}}}*/ -void Inputs::SetTriaControlInput(int enum_in,int layout,int interpolation,int control_id,int numindices,int* indices,IssmDouble* values,IssmDouble* values_min,IssmDouble* values_max){/*{{{*/ +void Inputs::SetControlInput(int enum_in,int layout,int interpolation,int control_id){/*{{{*/ bool recreate = false; + /*Get input id*/ int id = EnumToIndex(enum_in); @@ -711,7 +712,27 @@ } - /*Set input*/ - ControlInput* input = xDynamicCast(this->inputs[id]); - input->SetControl(interpolation,numindices,indices,values,values_min,values_max); +}/*}}}*/ +void Inputs::SetTransientControlInput(int enum_in,int control_id,IssmDouble* times,int numtimes){/*{{{*/ + + bool recreate = false; + + /*Get input id*/ + int id = EnumToIndex(enum_in); + + /*Create it if necessary*/ + if(this->inputs[id]){ + if(this->inputs[id]->ObjectEnum()!=ControlInputEnum){ + delete this->inputs[id]; + recreate = true; + } + } + else{ + recreate = true; + } + + if(recreate){ + this->inputs[id] = new ControlInput(enum_in,this->numberofelements_local,this->numberofvertices_local,control_id,times,numtimes); + } + }/*}}}*/ void Inputs::SetTriaControlInputGradient(int enum_in,int interpolation,int numindices,int* indices,IssmDouble* values){/*{{{*/ @@ -822,29 +843,4 @@ input->SetInput(interpolation,row,numindices,values); }/*}}}*/ -void Inputs::SetPentaControlInput(int enum_in,int layout,int interpolation,int control_id,int numindices,int* indices,IssmDouble* values,IssmDouble* values_min,IssmDouble* values_max){/*{{{*/ - - bool recreate = false; - /*Get input id*/ - int id = EnumToIndex(enum_in); - - /*Create it if necessary*/ - if(this->inputs[id]){ - if(this->inputs[id]->ObjectEnum()!=ControlInputEnum){ - delete this->inputs[id]; - recreate = true; - } - } - else{ - recreate = true; - } - - if(recreate){ - this->inputs[id] = new ControlInput(this->numberofelements_local,this->numberofvertices_local,layout,interpolation,control_id); - } - - /*Set input*/ - ControlInput* input = xDynamicCast(this->inputs[id]); - input->SetControl(interpolation,numindices,indices,values,values_min,values_max); -}/*}}}*/ void Inputs::SetPentaControlInputGradient(int enum_in,int interpolation,int numindices,int* indices,IssmDouble* values){/*{{{*/ Index: /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Inputs/Inputs.h (revision 25406) @@ -73,7 +73,8 @@ void SetDoubleInput(int enum_in,int index,IssmDouble value); void SetTransientInput(int enum_in,IssmDouble* times,int numtimes); + void SetControlInput(int enum_in,int layout,int interpolation,int control_id); + void SetTransientControlInput(int enum_in,int control_id,IssmDouble* times,int numtimes); TransientInput* SetDatasetTransientInput(int enum_in,int id,IssmDouble* times,int numtimes); void SetArrayInput(int enum_in,int row,IssmDouble* layers,int numlayers); - void SetTriaControlInput(int enum_in,int layout,int interpolation,int id,int numindices,int* indices,IssmDouble* values,IssmDouble* values_min,IssmDouble* values_max); void SetTriaControlInputGradient(int enum_in,int interpolation,int numindices,int* indices,IssmDouble* values); void SetTriaControlInputGradient(int enum_in,int interpolation,int numindices,int* indices,IssmDouble* values,int n); @@ -82,5 +83,4 @@ void SetTriaInput(int enum_in,int interpolation,int numindices,int* indices,IssmDouble* values); void SetTriaInput(int enum_in,int interpolation,int row,int numindices,IssmDouble* values); - void SetPentaControlInput(int enum_in,int layout,int interpolation,int id,int numindices,int* indices,IssmDouble* values,IssmDouble* values_min,IssmDouble* values_max); void SetPentaControlInputGradient(int enum_in,int interpolation,int numindices,int* indices,IssmDouble* values); void SetPentaDatasetInput(int enum_in,int id,int interpolation,int numindices,int* indices,IssmDouble* values); Index: /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.cpp (revision 25406) @@ -80,5 +80,10 @@ output->inputs = xNew(this->numtimesteps); for(int i=0;inumtimesteps;i++){ - output->inputs[i] = this->inputs[i]->copy(); + if(this->inputs[i]){ + output->inputs[i] = this->inputs[i]->copy(); + } + else{ + output->inputs[i] = NULL; + } } } Index: /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.h =================================================================== --- /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.h (revision 25405) +++ /issm/trunk-jpl/src/c/classes/Inputs/TransientInput.h (revision 25406) @@ -38,5 +38,5 @@ /*}}}*/ /*Object virtual functions definitions:{{{*/ - Input* copy(); + Input* copy(); void Configure(Parameters* params); void DeepEcho(); Index: /issm/trunk-jpl/src/c/modules/GetVectorFromControlInputsx/GetVectorFromControlInputsx.cpp =================================================================== --- /issm/trunk-jpl/src/c/modules/GetVectorFromControlInputsx/GetVectorFromControlInputsx.cpp (revision 25405) +++ /issm/trunk-jpl/src/c/modules/GetVectorFromControlInputsx/GetVectorFromControlInputsx.cpp (revision 25406) @@ -32,5 +32,5 @@ for(int j=0;jSize();j++){ Element* element=(Element*)elements->GetObjectByOffset(j); - element->GetVectorFromControlInputs(vector,control_type[i],i,data,offset); + element->GetVectorFromControlInputs(vector,control_type[i],i,N[i],data,offset); } offset += M[i]*N[i]; @@ -47,5 +47,5 @@ /*intermediary: */ - int N; + int N; Vector* vec_vector=NULL;