Changeset 22427


Ignore:
Timestamp:
02/16/18 10:52:23 (7 years ago)
Author:
erobo
Message:

CHG: adding cpp files for Nodalvalues, Misfit, and Regionaloutput

Location:
issm/trunk-jpl/src/c/classes
Files:
3 added
3 edited

Legend:

Unmodified
Added
Removed

  • issm/trunk-jpl/src/c/classes/Misfit.h

    r22137 r22427  
    77
    88/*Headers:*/
    9 /*{{{*/
    109#include "./Definition.h"
    1110#include "../datastructures/datastructures.h"
    … …  
    1716#include "../classes/Inputs/Input.h"
    1817#include "../classes/gauss/Gauss.h"
    19 /*}}}*/
     18
    2019IssmDouble OutputDefinitionsResponsex(FemModel* femmodel,int output_enum);
    21 void    GetVectorFromInputsx( IssmDouble** pvector, int* pvector_size, FemModel* femmodel,int name);
     20void  GetVectorFromInputsx( IssmDouble** pvector, int* pvector_size, FemModel* femmodel,int name);
    2221
    2322class Misfit: public Object, public Definition{
    … …  
    3736               
    3837                /*Misfit constructors, destructors :*/
    39                 Misfit(){/*{{{*/
    40 
    41                         this->definitionenum = -1;
    42                         this->name = NULL;
    43                         this->model_enum = UNDEF;
    44                         this->observation_enum = UNDEF;
    45                         this->weights_enum = UNDEF;
    46                         this->timeinterpolation=NULL;
    47                         this->local=1;
    48                         this->misfit=0;
    49                         this->lock=0;
    50 
    51                 }
    52                 /*}}}*/
    53                 Misfit(char* in_name, int in_definitionenum, int in_model_enum, int in_observation_enum, char* in_timeinterpolation, int in_local, int in_weights_enum){/*{{{*/
    54 
    55                         this->definitionenum=in_definitionenum;
    56                        
    57                         this->name              = xNew<char>(strlen(in_name)+1);
    58                         xMemCpy<char>(this->name,in_name,strlen(in_name)+1);
    59 
    60                         this->timeinterpolation = xNew<char>(strlen(in_timeinterpolation)+1);
    61                         xMemCpy<char>(this->timeinterpolation,in_timeinterpolation,strlen(in_timeinterpolation)+1);
    62                                                
    63                         this->model_enum=in_model_enum;
    64                         this->observation_enum=in_observation_enum;
    65                         this->weights_enum=in_weights_enum;
    66                         this->local=in_local;
    67                        
    68                         this->misfit=0;
    69                         this->lock=0;
    70                 }
    71                 /*}}}*/
    72                 ~Misfit(){/*{{{*/
    73                         if(this->name)xDelete(this->name);
    74                         if(this->timeinterpolation)xDelete(this->timeinterpolation);
    75                         this->misfit=0;
    76                         this->lock=0;
    77                 }
    78                 /*}}}*/
     38                Misfit();
     39                Misfit(char* in_name, int in_definitionenum, int in_model_enum, int in_observation_enum, char* in_timeinterpolation, int in_local, int in_weights_enum);
     40                ~Misfit();
     41               
    7942                /*Object virtual function resolutoin: */
    80                 Object* copy() {/*{{{*/
    81                         Misfit* mf = new Misfit(this->name,this->definitionenum, this->model_enum,this->observation_enum,this->timeinterpolation,this->local,this->weights_enum);
    82                         mf->misfit=this->misfit;
    83                         mf->lock=this->lock;
    84                         return (Object*) mf;
    85                 }
    86                 /*}}}*/
    87                 void DeepEcho(void){/*{{{*/
    88                         this->Echo();
    89                 }
    90                 /*}}}*/
    91                 void Echo(void){/*{{{*/
    92                         _printf_(" Misfit: " << name << " " << this->definitionenum << "\n");
    93                         _printf_("    model_enum: " << model_enum << " " << EnumToStringx(model_enum) << "\n");
    94                         _printf_("    observation_enum: " << observation_enum << " " << EnumToStringx(observation_enum) << "\n");
    95                         _printf_("    weights_enum: " << weights_enum << " " << EnumToStringx(weights_enum) << "\n");
    96                         _printf_("    timeinterpolation: " << timeinterpolation << "\n");
    97                         _printf_("    local: " << local << "\n");
    98                 }
    99                 /*}}}*/
    100                 int Id(void){/*{{{*/
    101                         return -1;
    102                 }
    103                 /*}}}*/
    104                 void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){/*{{{*/
    105                         _error_("not implemented yet!");
    106                 }
    107                 /*}}}*/
    108                 int ObjectEnum(void){/*{{{*/
    109                         return MisfitEnum;
    110                 }
    111                 /*}}}*/
     43                Object* copy();
     44                void DeepEcho(void);
     45                void Echo(void);
     46                int Id(void);
     47                void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction);
     48                int ObjectEnum(void);
     49               
    11250                /*Definition virtual function resolutoin: */
    113                 int DefinitionEnum(){/*{{{*/
    114                         return this->definitionenum;
    115                 }
    116                 /*}}}*/
    117                 char* Name(){/*{{{*/
    118                         char* name2=xNew<char>(strlen(this->name)+1);
    119                         xMemCpy(name2,this->name,strlen(this->name)+1);
    120 
    121                         return name2;
    122                 }
    123                 /*}}}*/
    124                  IssmDouble Response(FemModel* femmodel){/*{{{*/
    125                                  
    126                          /*diverse: */
    127                          IssmDouble time,starttime,finaltime;
    128                          IssmDouble dt;
    129                          
    130                          /*recover time parameters: */
    131                          femmodel->parameters->FindParam(&starttime,TimesteppingStartTimeEnum);
    132                          femmodel->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum);
    133                          femmodel->parameters->FindParam(&time,TimeEnum);
    134                          femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum);
    135 
    136                          if (this->local==1){ /*area integration using elements: {{{*/
    137 
    138                                  int i;
    139                                  IssmDouble misfit_t=0.;
    140                                  IssmDouble all_misfit_t=0.;
    141                                  IssmDouble area_t=0.;
    142                                  IssmDouble all_area_t;
    143 
    144                        
    145                                  /*If we are locked, return time average: */
    146                                  if(this->lock)return misfit/(time-starttime);
    147 
    148                                  for(i=0;i<femmodel->elements->Size();i++){
    149                                          Element* element=(Element*)femmodel->elements->GetObjectByOffset(i);
    150                                          misfit_t+=element->Misfit(model_enum,observation_enum,weights_enum);
    151                                          area_t+=element->MisfitArea(weights_enum);
    152                                  }
    153 
    154                                  ISSM_MPI_Allreduce ( (void*)&misfit_t,(void*)&all_misfit_t,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
    155                                  ISSM_MPI_Allreduce ( (void*)&area_t,(void*)&all_area_t,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
    156                                  area_t=all_area_t;
    157                                  misfit_t=all_misfit_t;
    158                                  
    159                                  /*Divide by surface area if not nill!: */
    160                                  if (area_t!=0) misfit_t=misfit_t/area_t;
    161 
    162                                  /*Add this time's contribution to curent misfit: */
    163                                  misfit+=dt*misfit_t;
    164 
    165                                  /*Do we lock? i.e. are we at final_time? :*/
    166                                  if(time==finaltime)this->lock=1;
    167 
    168                                  /*What we return is the value of misfit / time if transient*/
    169                                  if(time!=0.) return misfit/(time-starttime);
    170                                  return misfit;
    171                          } /*}}}*/
    172                          else if (this->local==2){ /*vertex by vertex computation: {{{*/
    173 
    174                                  IssmDouble* model = NULL;
    175                                  IssmDouble* observation= NULL;
    176                                  IssmDouble* weights= NULL;
    177                                  int msize,osize,wsize;
    178                                  
    179                                  /*Are we transient?:*/
    180                                  if (time==0){
    181                                          IssmDouble misfit_t=0.;
    182                                          
    183                                          /*get global vectors: */
    184                                          GetVectorFromInputsx(&model,&msize,femmodel,model_enum);
    185                                          GetVectorFromInputsx(&observation,&osize,femmodel,observation_enum);_assert_(msize==osize);
    186                                          GetVectorFromInputsx(&weights,&wsize,femmodel,weights_enum); _assert_(wsize==msize);
    187 
    188                                          int count=0;
    189                                          for (int i=0;i<msize;i++){
    190                                                  misfit_t += pow(model[i]-observation[i],2)*weights[i];
    191                                                  if (weights[i]!=0)count++;
    192                                          }
    193                                          misfit=sqrt(misfit_t/count);
    194 
    195                                          /*Free ressources:*/
    196                                          xDelete<IssmDouble>(model);
    197                                          xDelete<IssmDouble>(observation);
    198                                          xDelete<IssmDouble>(weights);
    199 
    200                                          /*return value: */
    201                                          return misfit;
    202                                  }
    203                                  else{
    204                                          
    205                                          IssmDouble misfit_t=0.;
    206                                          IssmDouble all_misfit_t=0.;
    207 
    208                                          /*If we are locked, return time average: */
    209                                          if(this->lock)return misfit/(time-starttime);
    210 
    211                                          /*get global vectors: */
    212                                          GetVectorFromInputsx(&model,&msize,femmodel,model_enum);
    213                                          GetVectorFromInputsx(&observation,&osize,femmodel,observation_enum);_assert_(msize==osize);
    214                                          GetVectorFromInputsx(&weights,&wsize,femmodel,weights_enum); _assert_(wsize==msize);
    215                                          
    216                                          int count=0;
    217                                          for (int i=0;i<msize;i++){
    218                                                  misfit_t += pow(model[i]-observation[i],2)*weights[i];
    219                                                  if (weights[i]!=0)count++;
    220                                          }
    221                                          misfit=sqrt(misfit_t/count);
    222 
    223                                          /*Add this time's contribution to curent misfit: */
    224                                          misfit=sqrt(misfit_t)/count;
    225                                          misfit+=dt*misfit_t;
    226 
    227                                          /*Do we lock? i.e. are we at final_time? :*/
    228                                          if(time==finaltime)this->lock=1;
    229                                          
    230                                          /*Free ressources:*/
    231                                          xDelete<IssmDouble>(model);
    232                                          xDelete<IssmDouble>(observation);
    233                                          xDelete<IssmDouble>(weights);
    234 
    235                                          /*What we return is the value of misfit / time: */
    236                                          return misfit/(time-starttime);
    237                                  }
    238 
    239                          } /*}}}*/
    240                          else{ /*global computation: {{{ */
    241                                  
    242                                  IssmDouble model, observation;
    243                                  
    244                                  /*If we are locked, return time average: */
    245                                  if(this->lock)return misfit/(time-starttime);
    246 
    247                                  /*First, the global  model response: */
    248                                  model=OutputDefinitionsResponsex(femmodel,this->model_enum);
    249                                  /*Now, the observation is buried inside the elements, go fish it in the first element (cludgy, needs fixing): */
    250                                  Element* element=(Element*)femmodel->elements->GetObjectByOffset(0); _assert_(element);
    251                                  Input* input = element->GetInput(observation_enum); _assert_(input);
    252                                  input->GetInputAverage(&observation);
    253                                  
    254                                  /*Add this time's contribution to curent misfit: */
    255                                  misfit+=dt*(model-observation);
    256                                  
    257                                  /*Do we lock? i.e. are we at final_time? :*/
    258                                  if(time==finaltime)this->lock=1;
    259                                  
    260                                  /*What we return is the value of misfit / time: */
    261                                  return misfit/(time-starttime);
    262                          } /*}}}*/
    263 
    264                  }
    265                         /*}}}*/
     51                int DefinitionEnum();
     52                char* Name();
     53                IssmDouble Response(FemModel* femmodel);
    26654};
    267 
    26855#endif  /* _MISFIT_H_ */

  • issm/trunk-jpl/src/c/classes/Nodalvalue.h

    r20810 r22427  
    3232               
    3333                /*Nodalvalue constructors, destructors :*/
    34                 Nodalvalue(){/*{{{*/
     34Nodalvalue();
     35Nodalvalue(char* in_name, int in_definitionenum, int in_model_enum, int in_node);
     36~Nodalvalue();
    3537
    36                         this->definitionenum = -1;
    37                         this->name = NULL;
    38                         this->model_enum = UNDEF;
    39                         this->node = -1;
     38/*Object virtual function resolutoin: */
     39Object* copy();
     40void DeepEcho(void);
     41void Echo(void);
     42int Id(void);
     43void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction);
     44int ObjectEnum(void);
    4045
    41                 }
    42                 /*}}}*/
    43                 Nodalvalue(char* in_name, int in_definitionenum, int in_model_enum, int in_node){/*{{{*/
    44 
    45                         this->definitionenum=in_definitionenum;
    46                         this->name   = xNew<char>(strlen(in_name)+1);
    47                         xMemCpy<char>(this->name,in_name,strlen(in_name)+1);
    48 
    49                         this->model_enum=in_model_enum;
    50                         this->node=in_node;
    51                 }
    52                 /*}}}*/
    53                 ~Nodalvalue(){/*{{{*/
    54                         if(this->name)xDelete(this->name);
    55                 }
    56                 /*}}}*/
    57                 /*Object virtual function resolutoin: */
    58                 Object* copy() {/*{{{*/
    59                         Nodalvalue* mf = new Nodalvalue(this->name,this->definitionenum, this->model_enum,this->node);
    60                         return (Object*) mf;
    61                 }
    62                 /*}}}*/
    63                 void DeepEcho(void){/*{{{*/
    64                         this->Echo();
    65                 }
    66                 /*}}}*/
    67                 void Echo(void){/*{{{*/
    68                         _printf_(" Nodalvalue: " << name << " " << this->definitionenum << "\n");
    69                         _printf_("    model_enum: " << model_enum << " " << EnumToStringx(model_enum) << "\n");
    70                         _printf_("    node: " << node << "\n");
    71                 }
    72                 /*}}}*/
    73                 int Id(void){/*{{{*/
    74                         return -1;
    75                 }
    76                 /*}}}*/
    77                 void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){/*{{{*/
    78                         _error_("not implemented yet!");
    79                 }
    80                 /*}}}*/
    81                 int ObjectEnum(void){/*{{{*/
    82                         return NodalvalueEnum;
    83                 }
    84                 /*}}}*/
    85                 /*Definition virtual function resolutoin: */
    86                 int DefinitionEnum(){/*{{{*/
    87 
    88                         return this->definitionenum;
    89                 }
    90                 /*}}}*/
    91                 char* Name(){/*{{{*/
    92 
    93                         char* name2=xNew<char>(strlen(this->name)+1);
    94                         xMemCpy(name2,this->name,strlen(this->name)+1);
    95 
    96                         return name2;
    97                 }
    98                 /*}}}*/
    99                  IssmDouble Response(FemModel* femmodel){/*{{{*/
    100                        
    101                          /*output:*/
    102                          IssmDouble value;
    103 
    104                          /*set index, which will be used by the NodalValue module: */
    105                          femmodel->parameters->SetParam(node,IndexEnum);
    106 
    107                          /*call Nodalvalue:*/
    108                          NodalValuex(&value, model_enum, femmodel->elements, femmodel->nodes, femmodel->vertices, femmodel->loads,
    109                                          femmodel->materials, femmodel->parameters);
    110 
    111                          /*done:*/
    112                          return value;
    113                  }
    114                  /*}}}*/
     46/*Definition virtual function resolutoin: */
     47int DefinitionEnum();
     48char* Name();
     49IssmDouble Response(FemModel* femmodel);
    11550};
    11651

  • issm/trunk-jpl/src/c/classes/Regionaloutput.h

    r22336 r22427  
    1919class Regionaloutput: public Object, public Definition{
    2020
    21         public:
     21public:
    2222
    23                 int         definitionenum;
    24                 char*       outputname;
    25                 char*       name;
    26                 IssmDouble* mask;
    27                 int         M;
    28                
    29                 /*Regionalicevolume: constructors, destructors :*/
    30                 Regionaloutput(){/*{{{*/
     23        int         definitionenum;
     24        char*       outputname;
     25        char*       name;
     26        IssmDouble* mask;
     27        int         M;
     28       
     29/*Regionalicevolume: constructors, destructors :*/
     30Regionaloutput();
     31Regionaloutput(char* in_name, int in_definitionenum, char* in_outputname, IssmDouble* maskin, int Min);
     32~Regionaloutput();
    3133
    32                         this->definitionenum = -1;
    33                         this->outputname = NULL;
    34                         this->name = NULL;
    35                         this->mask=NULL;
    36                         this->M=0;
     34/*Object virtual function resolutoin: */
     35Object* copy();
     36void DeepEcho(void);
     37void Echo(void);
     38int Id(void);
     39void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction);
     40int ObjectEnum(void);
    3741
    38                 }
    39                 /*}}}*/
    40                 Regionaloutput(char* in_name, int in_definitionenum, char* in_outputname, IssmDouble* maskin, int Min){ /*{{{*/
    41 
    42                         this->definitionenum=in_definitionenum;
    43                         this->outputname = xNew<char>(strlen(in_outputname)+1);
    44                         xMemCpy<char>(this->outputname,in_outputname,strlen(in_outputname)+1);
    45                         this->name = xNew<char>(strlen(in_name)+1);
    46                         xMemCpy<char>(this->name,in_name,strlen(in_name)+1);
    47 
    48                         this->mask   = xNew<IssmDouble>(Min);
    49                         xMemCpy<IssmDouble>(this->mask, maskin, Min);
    50 
    51                         this->M=Min;
    52 
    53                 }
    54                 /*}}}*/
    55                 ~Regionaloutput(){/*{{{*/
    56                         if(this->name)xDelete(this->name);
    57                         if(this->outputname)xDelete(this->outputname);
    58                         if(this->mask)xDelete(this->mask);
    59                 }
    60                 /*}}}*/
    61                 /*Object virtual function resolutoin: */
    62                 Object* copy() {/*{{{*/
    63                         Regionaloutput* mf = new Regionaloutput(this->name,this->definitionenum,this->outputname,this->mask,this->M);
    64                         return (Object*) mf;
    65                 }
    66                 /*}}}*/
    67                 void DeepEcho(void){/*{{{*/
    68                         this->Echo();
    69                 }
    70                 /*}}}*/
    71                 void Echo(void){/*{{{*/
    72                         _printf_(" Regionaloutput: " << this->name << " " << this->definitionenum << "\n");
    73                         _printf_("    outputname enum: " << this->outputname << "Enum\n");
    74                         _printf_("    mask: " << this->mask << "\n");
    75                         _printf_("    M: " << this->M << "\n");
    76                 }
    77                 /*}}}*/
    78                 int Id(void){/*{{{*/
    79                         return -1;
    80                 }
    81                 /*}}}*/
    82                 void Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){/*{{{*/
    83                         _error_("not implemented yet!");
    84                 }
    85                 /*}}}*/
    86                 int ObjectEnum(void){/*{{{*/
    87                         return RegionaloutputEnum;
    88                 }
    89                 /*}}}*/
    90                 /*Definition virtual function resolutoin: */
    91                 int DefinitionEnum(){/*{{{*/
    92 
    93                         return this->definitionenum;
    94                 }
    95                 /*}}}*/
    96                 char* Name(){/*{{{*/
    97 
    98                         char* name2=xNew<char>(strlen(this->name)+1);
    99                         xMemCpy(name2,this->name,strlen(this->name)+1);
    100 
    101                         return name2;
    102                 }
    103                 /*}}}*/
    104                 IssmDouble Response(FemModel* femmodel){/*{{{*/
    105 
    106                         int i;
    107                         IssmDouble val_t=0.;
    108                         IssmDouble all_val_t=0.;
    109                         int outputenum = StringToEnumx(this->outputname);
    110 
    111                         for(i=0;i<femmodel->elements->Size();i++){
    112                                 Element* element=(Element*)femmodel->elements->GetObjectByOffset(i);
    113                                 switch(outputenum){
    114                                         case GroundedAreaEnum:
    115                                                 val_t+=element->GroundedArea(this->mask,false);
    116                                                 break;
    117                                         case GroundedAreaScaledEnum:
    118                                                 val_t+=element->GroundedArea(this->mask,true);
    119                                                 break;
    120                                         case FloatingAreaEnum:
    121                                                 val_t+=element->FloatingArea(this->mask,false);
    122                                                 break;
    123                                         case FloatingAreaScaledEnum:
    124                                                 val_t+=element->FloatingArea(this->mask,true);
    125                                                 break;
    126                                         case IceMassEnum:
    127                                                 val_t+=element->IceMass(this->mask,false);
    128                                                 break;
    129                                         case IceMassScaledEnum:
    130                                                 val_t+=element->IceMass(this->mask,true);
    131                                                 break;
    132                                         case IceVolumeEnum:
    133                                                 val_t+=element->IceVolume(this->mask,false);
    134                                                 break;
    135                                         case IceVolumeScaledEnum:
    136                                                 val_t+=element->IceVolume(this->mask,true);
    137                                                 break;
    138                                         case IceVolumeAboveFloatationEnum:
    139                                                 val_t+=element->IceVolumeAboveFloatation(this->mask,false);
    140                                                 break;
    141                                         case IceVolumeAboveFloatationScaledEnum:
    142                                                 val_t+=element->IceVolumeAboveFloatation(this->mask,true);
    143                                                 break;
    144                                         case TotalFloatingBmbEnum:
    145                                                 val_t+=element->TotalFloatingBmb(this->mask,false);
    146                                                 break;
    147                                         case TotalFloatingBmbScaledEnum:
    148                                                 val_t+=element->TotalFloatingBmb(this->mask,true);
    149                                                 break;
    150                                         case TotalGroundedBmbEnum:
    151                                                 val_t+=element->TotalGroundedBmb(this->mask,false);
    152                                                 break;
    153                                         case TotalGroundedBmbScaledEnum:
    154                                                 val_t+=element->TotalGroundedBmb(this->mask,true);
    155                                                 break;
    156                                         case TotalSmbEnum:
    157                                                 val_t+=element->TotalSmb(this->mask,false);
    158                                                 break;
    159                                         case TotalSmbScaledEnum:
    160                                                 val_t+=element->TotalSmb(this->mask,true);
    161                                                 break;
    162                                         default:
    163                                                 _error_("Regional output type " << this->outputname << " not supported yet!");
    164                                 }
    165                         }
    166 
    167                         ISSM_MPI_Allreduce ( (void*)&val_t,(void*)&all_val_t,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
    168                         val_t=all_val_t;
    169 
    170                         return val_t;
    171                 }
    172                 /*}}}*/
     42/*Definition virtual function resolutoin: */
     43int DefinitionEnum();
     44char* Name();
     45IssmDouble Response(FemModel* femmodel);
    17346};
    174 
    17547#endif  /* _REGIONALOUTPUT_H_ */
Note: See TracChangeset for help on using the changeset viewer.