/*!\file Misfit.h * \brief: header file for Misfit object */ #ifndef _MISFIT_H_ #define _MISFIT_H_ /*Headers:*/ /*{{{*/ #include "./Definition.h" #include "../datastructures/datastructures.h" #include "./Elements/Element.h" #include "./Elements/Elements.h" #include "./FemModel.h" #include "../modules/SurfaceAreax/SurfaceAreax.h" #include "../classes/Params/Parameters.h" /*}}}*/ class Misfit: public Object, public Definition{ public: char* name; int model_enum; int observation_enum; int weights_enum; char* timeinterpolation; IssmDouble misfit; //value carried over in time. int lock; // if lock is on, we just return the value stored in "misfit". this is used so we don't compute misfit past the final_time /*Misfit constructors, destructors :*/ /*FUNCTION Misfit() {{{*/ Misfit(){ this->name = NULL; this->model_enum = UNDEF; this->observation_enum = UNDEF; this->weights_enum = UNDEF; this->timeinterpolation=NULL; this->misfit=0; this->lock=0; } /*}}}*/ /*FUNCTION Misfit(char* in_name, int in_model_enum, int in_observation_enum char* in_timeinterpolation, int in_weights_enum) {{{*/ Misfit(char* in_name, int in_model_enum, int in_observation_enum, char* in_timeinterpolation, int in_weights_enum){ this->name = xNew(strlen(in_name)+1); xMemCpy(this->name,in_name,strlen(in_name)+1); this->timeinterpolation = xNew(strlen(in_timeinterpolation)+1); xMemCpy(this->timeinterpolation,in_timeinterpolation,strlen(in_timeinterpolation)+1); this->model_enum=in_model_enum; this->observation_enum=in_observation_enum; this->weights_enum=in_weights_enum; this->misfit=0; this->lock=0; } /*}}}*/ /*FUNCTION ~Misfit() {{{*/ ~Misfit(){ if(this->name)xDelete(this->name); if(this->timeinterpolation)xDelete(this->timeinterpolation); this->misfit=0; this->lock=0; } /*}}}*/ /*Object virtual function resolutoin: */ /*FUNCTION Echo(){{{*/ void Echo(void){ _printf_(" Misfit: " << name << "\n"); _printf_(" model_enum: " << model_enum << " " << EnumToStringx(model_enum) << "\n"); _printf_(" observation_enum: " << observation_enum << " " << EnumToStringx(observation_enum) << "\n"); _printf_(" weights_enum: " << weights_enum << " " << EnumToStringx(weights_enum) << "\n"); _printf_(" timeinterpolation: " << timeinterpolation << "\n"); } /*}}}*/ /*FUNCTION DeepEcho(){{{*/ void DeepEcho(void){ this->Echo(); } /*}}}*/ /*FUNCTION Id(){{{*/ int Id(void){ return -1; } /*}}}*/ /*FUNCTION ObjectEnum{{{*/ int ObjectEnum(void){ return MisfitEnum; } /*}}}*/ /*FUNCTION copy {{{*/ Object* copy() { return new Misfit(this->name,this->model_enum,this->observation_enum,this->timeinterpolation,this->weights_enum); } /*}}}*/ /*Definition virtual function resolutoin: */ /*FUNCTION char* Name() {{{*/ char* Name(){ char* name2=xNew(strlen(this->name)+1); xMemCpy(name2,this->name,strlen(this->name)+1); return name2; } /*}}}*/ /*FUNCTION IssmDouble Response(FemModel* femmodel) {{{*/ IssmDouble Response(FemModel* femmodel){ int i; IssmDouble misfit_t=0.; IssmDouble all_misfit_t=0.; IssmDouble dt; IssmDouble area_t=0.; IssmDouble all_area_t; IssmDouble time,starttime,finaltime; femmodel->parameters->FindParam(&starttime,TimesteppingStartTimeEnum); femmodel->parameters->FindParam(&finaltime,TimesteppingFinalTimeEnum); femmodel->parameters->FindParam(&time,TimeEnum); /*If we are locked, return time average: */ if(this->lock)return misfit/(time-starttime); for(i=0;ielements->Size();i++){ Element* element=(Element*)femmodel->elements->GetObjectByOffset(i); misfit_t+=element->Misfit(model_enum,observation_enum,weights_enum); area_t+=element->MisfitArea(weights_enum); } ISSM_MPI_Allreduce ( (void*)&misfit_t,(void*)&all_misfit_t,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm()); ISSM_MPI_Allreduce ( (void*)&area_t,(void*)&all_area_t,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm()); area_t=all_area_t; misfit_t=all_misfit_t; /*Divide by surface area if not nill!: */ if (area_t!=0) misfit_t=misfit_t/area_t; /*Recover delta_t: */ femmodel->parameters->FindParam(&dt,TimesteppingTimeStepEnum); /*Add this time's contribution to curent misfit: */ misfit+=dt*misfit_t; /*Do we lock? i.e. are we at final_time? :*/ if(time==finaltime)this->lock=1; /*What we return is the value of misfit / time: */ return misfit/(time-starttime); } /*}}}*/ }; #endif /* _MISFIT_H_ */