Ice Sheet System Model  4.18
Code documentation
Cfsurfacesquare.cpp
Go to the documentation of this file.
1 
5 /*Headers:*/
6 /*{{{*/
7 #ifdef HAVE_CONFIG_H
8  #include <config.h>
9 #else
10 #error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
11 #endif
12 
13 #include "./classes.h"
16 #include "../datastructures/datastructures.h"
17 #include "./Elements/Element.h"
18 #include "./Elements/Elements.h"
19 #include "./FemModel.h"
20 #include "../modules/SurfaceAreax/SurfaceAreax.h"
21 #include "../classes/Params/Parameters.h"
22 #include "../classes/gauss/Gauss.h"
24 /*}}}*/
25 
26 /*Cfsurfacesquare constructors, destructors :*/
28 
29  this->definitionenum = -1;
30  this->name = NULL;
31  this->model_enum = UNDEF;
32  this->observation_enum = UNDEF;
33  this->weights_enum = UNDEF;
34  this->misfit=0;
35  this->lock=0;
36  this->datatime=0.;
37  this->timepassedflag = false;
38  this->last_time = 0.;
39 
40 }
41 /*}}}*/
42 Cfsurfacesquare::Cfsurfacesquare(char* in_name, int in_definitionenum, int in_model_enum, int in_observation_enum, int in_weights_enum, IssmDouble in_datatime, bool in_timepassedflag){/*{{{*/
43 
44  this->definitionenum=in_definitionenum;
45 
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->observation_enum=in_observation_enum;
51  this->weights_enum=in_weights_enum;
52  this->datatime=in_datatime;
53  this->timepassedflag=in_timepassedflag;
54 
55  this->misfit=0;
56  this->lock=0;
57 }
58 /*}}}*/
60  if(this->name)xDelete(this->name);
61  this->misfit=0;
62  this->lock=0;
63 }
64 /*}}}*/
65 /*Object virtual function resolutoin: */
68  mf->misfit=this->misfit;
69  mf->lock=this->lock;
70  return (Object*) mf;
71 }
72 /*}}}*/
73 void Cfsurfacesquare::DeepEcho(void){/*{{{*/
74  this->Echo();
75 }
76 /*}}}*/
77 void Cfsurfacesquare::Echo(void){/*{{{*/
78  _printf_(" Cfsurfacesquare: " << name << " " << this->definitionenum << "\n");
79  _printf_(" model_enum: " << model_enum << " " << EnumToStringx(model_enum) << "\n");
80  _printf_(" observation_enum: " << observation_enum << " " << EnumToStringx(observation_enum) << "\n");
81  _printf_(" weights_enum: " << weights_enum << " " << EnumToStringx(weights_enum) << "\n");
82  _printf_(" datatime: " << datatime << "\n");
83  _printf_(" timepassedflag: "<<timepassedflag<<"\n");
84 }
85 /*}}}*/
86 int Cfsurfacesquare::Id(void){/*{{{*/
87  return -1;
88 }
89 /*}}}*/
90 void Cfsurfacesquare::Marshall(char** pmarshalled_data,int* pmarshalled_data_size, int marshall_direction){/*{{{*/
91  _error_("not implemented yet!");
92 }
93 /*}}}*/
94 int Cfsurfacesquare::ObjectEnum(void){/*{{{*/
95  return CfsurfacesquareEnum;
96 }
97 /*}}}*/
98 /*Definition virtual function resolutoin: */
100  return this->definitionenum;
101 }
102 /*}}}*/
103 char* Cfsurfacesquare::Name(){/*{{{*/
104  char* name2=xNew<char>(strlen(this->name)+1);
105  xMemCpy(name2,this->name,strlen(this->name)+1);
106 
107  return name2;
108 }
109 /*}}}*/
111  /*diverse: */
112  IssmDouble time;
113 
114  /*recover time parameters: */
116  if(time < last_time) timepassedflag = false;
117  last_time = time;
118 
119  int i;
120  IssmDouble J=0.;
121  IssmDouble J_sum=0.;
122 
123  if(datatime<=time && !timepassedflag){
124  for(i=0;i<femmodel->elements->Size();i++){
127  }
128 
129  ISSM_MPI_Allreduce ( (void*)&J,(void*)&J_sum,1,ISSM_MPI_DOUBLE,ISSM_MPI_SUM,IssmComm::GetComm());
131  J=J_sum;
132 
133  timepassedflag = true;
134  return J;
135  }
136  else return J;
137  }
138  /*}}}*/
139 IssmDouble Cfsurfacesquare::Cfsurfacesquare_Calculation(Element* element, int model_enum, int observation_enum, int weights_enum){/*{{{*/
140 
141  int domaintype,numcomponents;
142  IssmDouble Jelem=0.;
143  IssmDouble misfit,Jdet;
144  IssmDouble model,obs,weight;
145  IssmDouble* xyz_list = NULL;
146 
147  /*Get basal element*/
148  if(!element->IsOnSurface()) return 0.;
149 
150  /*If on water, return 0: */
151  if(!element->IsIceInElement()) return 0.;
152 
153  /*Get problem dimension*/
154  element->FindParam(&domaintype,DomainTypeEnum);
155  switch(domaintype){
156  case Domain2DverticalEnum: numcomponents = 1; break;
157  case Domain3DEnum: numcomponents = 2; break;
158  case Domain2DhorizontalEnum: numcomponents = 2; break;
159  default: _error_("not supported yet");
160  }
161 
162  /*Spawn surface element*/
163  Element* topelement = element->SpawnTopElement();
164 
165  /* Get node coordinates*/
166  topelement->GetVerticesCoordinates(&xyz_list);
167 
168  /*Retrieve all inputs we will be needing: */
169  DatasetInput2 *datasetinput = topelement->GetDatasetInput2(definitionenum); _assert_(datasetinput);
170  Input2 *model_input = topelement->GetInput2(model_enum); _assert_(model_input);
171 
172  /* Start looping on the number of gaussian points: */
173  Gauss* gauss=topelement->NewGauss(2);
174  for(int ig=gauss->begin();ig<gauss->end();ig++){
175 
176  gauss->GaussPoint(ig);
177 
178  /* Get Jacobian determinant: */
179  topelement->JacobianDeterminant(&Jdet,xyz_list,gauss);
180 
181  /*Get all parameters at gaussian point*/
182  datasetinput->GetInputValue(&weight,gauss,WeightsSurfaceObservationEnum);
183  model_input->GetInputValue(&model,gauss);
184  datasetinput->GetInputValue(&obs,gauss,SurfaceObservationEnum);
185 
186  /*Compute SurfaceAbsVelMisfitEnum:
187  * *
188  * * 1 [ 2 2 ]
189  * * J = --- | (u - u ) + (v - v ) |
190  * * 2 [ obs obs ]
191  * *
192  * */
193  misfit=0.5*(model-obs)*(model-obs);
194 
195  /*Add to cost function*/
196  Jelem+=misfit*weight*Jdet*gauss->weight;
197  }
198 
199  /*clean up and Return: */
200  if(domaintype!=Domain2DhorizontalEnum){topelement->DeleteMaterials(); delete topelement;};
201  xDelete<IssmDouble>(xyz_list);
202  delete gauss;
203  return Jelem;
204 }/*}}}*/
DataSet::Size
int Size()
Definition: DataSet.cpp:399
SurfaceObservationEnum
@ SurfaceObservationEnum
Definition: EnumDefinitions.h:827
Cfsurfacesquare::DefinitionEnum
int DefinitionEnum()
Definition: Cfsurfacesquare.cpp:99
_assert_
#define _assert_(ignore)
Definition: exceptions.h:37
IssmDouble
double IssmDouble
Definition: types.h:37
DatasetInput2
Definition: DatasetInput2.h:14
ISSM_MPI_Allreduce
int ISSM_MPI_Allreduce(void *sendbuf, void *recvbuf, int count, ISSM_MPI_Datatype datatype, ISSM_MPI_Op op, ISSM_MPI_Comm comm)
Definition: issmmpi.cpp:116
Element::FindParam
void FindParam(bool *pvalue, int paramenum)
Definition: Element.cpp:933
Cfsurfacesquare::Marshall
void Marshall(char **pmarshalled_data, int *pmarshalled_data_size, int marshall_direction)
Definition: Cfsurfacesquare.cpp:90
_printf_
#define _printf_(StreamArgs)
Definition: Print.h:22
Cfsurfacesquare::observation_enum
int observation_enum
Definition: Cfsurfacesquare.h:23
ISSM_MPI_SUM
#define ISSM_MPI_SUM
Definition: issmmpi.h:134
FemModel::parameters
Parameters * parameters
Definition: FemModel.h:46
TimeEnum
@ TimeEnum
Definition: EnumDefinitions.h:427
Cfsurfacesquare::Cfsurfacesquare_Calculation
IssmDouble Cfsurfacesquare_Calculation(Element *element, int model_enum, int observation_enum, int weights_enum)
Definition: Cfsurfacesquare.cpp:139
IssmComm::GetComm
static ISSM_MPI_Comm GetComm(void)
Definition: IssmComm.cpp:30
Cfsurfacesquare::Id
int Id(void)
Definition: Cfsurfacesquare.cpp:86
Cfsurfacesquare::datatime
IssmDouble datatime
Definition: Cfsurfacesquare.h:26
Cfsurfacesquare::model_enum
int model_enum
Definition: Cfsurfacesquare.h:21
DatasetInput2::GetInputValue
void GetInputValue(IssmDouble *pvalue, Gauss *gauss, int index)
Definition: DatasetInput2.cpp:199
Cfsurfacesquare::~Cfsurfacesquare
~Cfsurfacesquare()
Definition: Cfsurfacesquare.cpp:59
Element::GetInput2
virtual Input2 * GetInput2(int inputenum)=0
Element::DeleteMaterials
void DeleteMaterials(void)
Definition: Element.cpp:429
Element
Definition: Element.h:41
Elements.h
Element::GetDatasetInput2
virtual DatasetInput2 * GetDatasetInput2(int inputenum)
Definition: Element.h:250
Domain2DhorizontalEnum
@ Domain2DhorizontalEnum
Definition: EnumDefinitions.h:534
ISSM_MPI_DOUBLE
#define ISSM_MPI_DOUBLE
Definition: issmmpi.h:125
Element.h
abstract class for Element object This class is a place holder for the Tria and the Penta elements....
Cfsurfacesquare::last_time
IssmDouble last_time
Definition: Cfsurfacesquare.h:28
Object
Definition: Object.h:13
DomainTypeEnum
@ DomainTypeEnum
Definition: EnumDefinitions.h:124
xDelete
void xDelete(T **&aT_pp)
Definition: MemOps.h:97
Element::NewGauss
virtual Gauss * NewGauss(void)=0
Cfsurfacesquare::timepassedflag
bool timepassedflag
Definition: Cfsurfacesquare.h:27
EnumToStringx
const char * EnumToStringx(int enum_in)
Definition: EnumToStringx.cpp:15
Cfsurfacesquare
Definition: Cfsurfacesquare.h:15
ExternalResult.h
abstract class for ExternalResult object
UNDEF
#define UNDEF
Definition: constants.h:8
Element::GetVerticesCoordinates
void GetVerticesCoordinates(IssmDouble **xyz_list)
Definition: Element.cpp:1446
Domain3DEnum
@ Domain3DEnum
Definition: EnumDefinitions.h:536
FemModel::elements
Elements * elements
Definition: FemModel.h:44
ISSM_MPI_Bcast
int ISSM_MPI_Bcast(void *buffer, int count, ISSM_MPI_Datatype datatype, int root, ISSM_MPI_Comm comm)
Definition: issmmpi.cpp:162
Cfsurfacesquare::Cfsurfacesquare
Cfsurfacesquare()
Definition: Cfsurfacesquare.cpp:27
Input2
Definition: Input2.h:18
FemModel
Definition: FemModel.h:31
Cfsurfacesquare::name
char * name
Definition: Cfsurfacesquare.h:22
Element::IsIceInElement
bool IsIceInElement()
Definition: Element.cpp:2021
_error_
#define _error_(StreamArgs)
Definition: exceptions.h:49
CfsurfacesquareEnum
@ CfsurfacesquareEnum
Definition: EnumDefinitions.h:1007
Element::SpawnTopElement
virtual Element * SpawnTopElement(void)=0
Cfsurfacesquare::lock
int lock
Definition: Cfsurfacesquare.h:30
Cfsurfacesquare::weights_enum
int weights_enum
Definition: Cfsurfacesquare.h:25
Gauss::begin
virtual int begin(void)=0
Cfsurfacesquare::DeepEcho
void DeepEcho(void)
Definition: Cfsurfacesquare.cpp:73
DataSet::GetObjectByOffset
Object * GetObjectByOffset(int offset)
Definition: DataSet.cpp:334
FemModel.h
WeightsSurfaceObservationEnum
@ WeightsSurfaceObservationEnum
Definition: EnumDefinitions.h:864
Cfsurfacesquare::Echo
void Echo(void)
Definition: Cfsurfacesquare.cpp:77
Element::JacobianDeterminant
virtual void JacobianDeterminant(IssmDouble *Jdet, IssmDouble *xyz_list, Gauss *gauss)=0
Gauss::GaussPoint
virtual void GaussPoint(int ig)=0
Cfsurfacesquare::copy
Object * copy()
Definition: Cfsurfacesquare.cpp:66
Parameters::FindParam
void FindParam(bool *pinteger, int enum_type)
Definition: Parameters.cpp:262
Cfsurfacesquare::Name
char * Name()
Definition: Cfsurfacesquare.cpp:103
Cfsurfacesquare::Response
IssmDouble Response(FemModel *femmodel)
Definition: Cfsurfacesquare.cpp:110
Cfsurfacesquare::ObjectEnum
int ObjectEnum(void)
Definition: Cfsurfacesquare.cpp:94
Input2::GetInputValue
virtual void GetInputValue(IssmDouble *pvalue, Gauss *gauss)
Definition: Input2.h:38
Gauss::weight
IssmDouble weight
Definition: Gauss.h:11
Element::IsOnSurface
bool IsOnSurface()
Definition: Element.cpp:1981
xMemCpy
T * xMemCpy(T *dest, const T *src, unsigned int size)
Definition: MemOps.h:152
Cfsurfacesquare::misfit
IssmDouble misfit
Definition: Cfsurfacesquare.h:31
Cfsurfacesquare::definitionenum
int definitionenum
Definition: Cfsurfacesquare.h:19
Domain2DverticalEnum
@ Domain2DverticalEnum
Definition: EnumDefinitions.h:535
Results.h
Gauss
Definition: Gauss.h:8
classes.h
DatasetInput2.h
: header file for datasetinput object
femmodel
FemModel * femmodel
Definition: esmfbinders.cpp:16