source: issm/trunk-jpl/src/c/objects/Inputs/TransientInput.cpp@ 11197

Last change on this file since 11197 was 11197, checked in by Mathieu Morlighem, 13 years ago

Fixed many warnings on pleiades: missing return statement, uninitialized
values, wrong print argument,... etc

File size: 13.2 KB
Line 
1/*!\file TransientInput.c
2 * \brief: implementation of the TransientInput object
3 */
4/*Headers{{{1*/
5#ifdef HAVE_CONFIG_H
6 #include <config.h>
7#else
8#error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
9#endif
10
11#include <stdio.h>
12#include <string.h>
13#include "../objects.h"
14#include "../../EnumDefinitions/EnumDefinitions.h"
15#include "../../shared/shared.h"
16#include "../../Container/Container.h"
17#include "../../include/include.h"
18/*}}}*/
19
20/*TransientInput constructors and destructor*/
21/*FUNCTION TransientInput::TransientInput(){{{1*/
22TransientInput::TransientInput(){
23
24 enum_type=UNDEF;
25 inputs=NULL;
26 this->numtimesteps=0;
27 this->parameters=NULL;
28 this->timesteps=NULL;
29
30}
31/*}}}*/
32/*FUNCTION TransientInput::TransientInput(int in_enum_type){{{1*/
33TransientInput::TransientInput(int in_enum_type)
34{
35 /*Set Enum*/
36 enum_type=in_enum_type;
37
38 /*Allocate values and timesteps, and copy: */
39 this->numtimesteps=0;
40 this->timesteps=NULL;
41 inputs = new Inputs();
42 this->parameters=NULL;
43
44}
45/*}}}*/
46/*FUNCTION TransientInput::~TransientInput{{{1*/
47TransientInput::~TransientInput(){
48 xfree((void**)&this->timesteps);
49 this->timesteps=NULL;
50 this->numtimesteps=0;
51 parameters=NULL;
52 delete this->inputs;
53 return;
54}
55/*}}}*/
56/*FUNCTION void TransientInput::AddTimeInput(Input* input,double time){{{1*/
57void TransientInput::AddTimeInput(Input* input,double time){
58
59 /*insert values at time step: */
60 if (this->numtimesteps>0 && time<=this->timesteps[this->numtimesteps-1]) _assert_("timestep values must increase sequentially");
61
62 //copy timesteps, add the new time, delete previous timesteps, and add the new input: inputs->AddObject(input);
63 double* old_timesteps=NULL;
64
65 if (this->numtimesteps > 0){
66 old_timesteps=(double*)xmalloc(this->numtimesteps*sizeof(double));
67 memcpy(old_timesteps,this->timesteps,this->numtimesteps*sizeof(double));
68 xfree((void**)&this->timesteps);
69 }
70
71 this->numtimesteps=this->numtimesteps+1;
72 this->timesteps=(double*)xmalloc(this->numtimesteps*sizeof(double));
73
74 if (this->numtimesteps > 1){
75 memcpy(this->timesteps,old_timesteps,(this->numtimesteps-1)*sizeof(double));
76 xfree((void**)&old_timesteps);
77 }
78
79 /*go ahead and plug: */
80 this->timesteps[this->numtimesteps-1]=time;
81 inputs->AddObject(input);
82
83}
84/*}}}*/
85
86/*Object virtual functions definitions:*/
87/*FUNCTION TransientInput::Echo {{{1*/
88void TransientInput::Echo(void){
89 this->DeepEcho();
90}
91/*}}}*/
92/*FUNCTION TransientInput::DeepEcho{{{1*/
93void TransientInput::DeepEcho(void){
94
95 int i;
96
97 printf("TransientInput:\n");
98 printf(" enum: %i (%s)\n",this->enum_type,EnumToStringx(this->enum_type));
99 printf(" numtimesteps: %i\n",this->numtimesteps);
100 printf("---inputs: \n");
101 for(i=0;i<this->numtimesteps;i++){
102 printf(" time: %g \n",this->timesteps[i]);
103 ((Input*)this->inputs->GetObjectByOffset(i))->Echo();
104 }
105}
106/*}}}*/
107/*FUNCTION TransientInput::Id{{{1*/
108int TransientInput::Id(void){ return -1; }
109/*}}}*/
110/*FUNCTION TransientInput::MyRank{{{1*/
111int TransientInput::MyRank(void){
112 extern int my_rank;
113 return my_rank;
114}
115/*}}}*/
116#ifdef _SERIAL_
117/*FUNCTION TransientInput::Marshall{{{1*/
118void TransientInput::Marshall(char** pmarshalled_dataset){
119
120 char* marshalled_dataset=NULL;
121 char* marshalled_inputs=NULL;
122 int marshalled_inputs_size;
123 int enum_value=0;
124
125 /*recover marshalled_dataset: */
126 marshalled_dataset=*pmarshalled_dataset;
127
128 /*get enum value of TransientInput: */
129 enum_value=TransientInputEnum;
130
131 /*marshall enum: */
132 memcpy(marshalled_dataset,&enum_value,sizeof(enum_value));marshalled_dataset+=sizeof(enum_value);
133
134 /*marshall TransientInput data: */
135 memcpy(marshalled_dataset,&enum_type,sizeof(enum_type));marshalled_dataset+=sizeof(enum_type);
136 memcpy(marshalled_dataset,&numtimesteps,sizeof(numtimesteps));marshalled_dataset+=sizeof(numtimesteps);
137 memcpy(marshalled_dataset,timesteps,numtimesteps*sizeof(double));marshalled_dataset+=numtimesteps*sizeof(double);
138
139 /*marshal inputs*/
140 marshalled_inputs_size=inputs->MarshallSize();
141 marshalled_inputs=inputs->Marshall();
142 memcpy(marshalled_dataset,marshalled_inputs,marshalled_inputs_size*sizeof(char));
143 marshalled_dataset+=marshalled_inputs_size;
144
145 /*clean up and assign output pointer*/
146 xfree((void**)&marshalled_inputs);
147 *pmarshalled_dataset=marshalled_dataset;
148
149}
150/*}}}*
151/*FUNCTION TransientInput::MarshallSize{{{1*/
152int TransientInput::MarshallSize(){
153
154 return
155 +sizeof(enum_type)+
156 +sizeof(numtimesteps)+
157 +inputs->MarshallSize()
158 +numtimesteps*sizeof(double)+
159 +sizeof(int); //sizeof(int) for enum value
160}
161/*}}}*/
162/*FUNCTION TransientInput::Demarshall{{{1*/
163void TransientInput::Demarshall(char** pmarshalled_dataset){
164
165 char* marshalled_dataset=NULL;
166
167 /*recover marshalled_dataset: */
168 marshalled_dataset=*pmarshalled_dataset;
169
170 /*this time, no need to get enum type, the pointer directly points to the beginning of the
171 *object data (thanks to DataSet::Demarshall):*/
172 memcpy(&enum_type,marshalled_dataset,sizeof(enum_type));marshalled_dataset+=sizeof(enum_type);
173 memcpy(&numtimesteps,marshalled_dataset,sizeof(numtimesteps));marshalled_dataset+=sizeof(numtimesteps);
174
175 /*allocate: */
176 timesteps=(double*)xmalloc(numtimesteps*sizeof(double));
177 memcpy(timesteps,marshalled_dataset,numtimesteps*sizeof(double));marshalled_dataset+=numtimesteps*sizeof(double);
178
179 /*Demarshal values*/
180 inputs=(Inputs*)DataSetDemarshallRaw(&marshalled_dataset);
181
182 /*return: */
183 *pmarshalled_dataset=marshalled_dataset;
184 return;
185
186}
187/*}}}*/
188#endif
189/*FUNCTION TransientInput::ObjectEnum{{{1*/
190int TransientInput::ObjectEnum(void){
191
192 return TransientInputEnum;
193
194}
195/*}}}*/
196/*FUNCTION TransientInput::copy{{{1*/
197Object* TransientInput::copy() {
198
199 TransientInput* output=NULL;
200
201 output = new TransientInput();
202 output->enum_type=this->enum_type;
203 output->numtimesteps=this->numtimesteps;
204 output->timesteps=(double*)xmalloc(this->numtimesteps*sizeof(double));
205 memcpy(output->timesteps,this->timesteps,this->numtimesteps*sizeof(double));
206 output->inputs=(Inputs*)this->inputs->Copy();
207 output->parameters=this->parameters;
208
209 return output;
210
211}
212/*}}}*/
213
214/*TransientInput management*/
215/*FUNCTION TransientInput::InstanceEnum{{{1*/
216int TransientInput::InstanceEnum(void){
217
218 return this->enum_type;
219
220}
221/*}}}*/
222/*FUNCTION TransientInput::SpawnTriaInput{{{1*/
223Input* TransientInput::SpawnTriaInput(int* indices){
224
225 /*output*/
226 TransientInput* outinput=NULL;
227
228 /*Create new Transientinput (copy of current input)*/
229 outinput=new TransientInput();
230 outinput->enum_type=this->enum_type;
231 outinput->numtimesteps=this->numtimesteps;
232 outinput->timesteps=(double*)xmalloc(this->numtimesteps*sizeof(double));
233 memcpy(outinput->timesteps,this->timesteps,this->numtimesteps*sizeof(double));
234 outinput->inputs=(Inputs*)this->inputs->SpawnTriaInputs(indices);
235 outinput->parameters=this->parameters;
236
237 /*Assign output*/
238 return outinput;
239
240}
241/*}}}*/
242/*FUNCTION TransientInput::SpawnResult{{{1*/
243ElementResult* TransientInput::SpawnResult(int step, double time){
244
245 ElementResult* elementresult=NULL;
246
247 /*Ok, we want to spawn an ElementResult. We have the time, just get
248 *the correct values: */
249 Input* input=GetTimeInput(time);
250
251 elementresult=input->SpawnResult(step,time);
252
253 delete input;
254
255 return elementresult;
256}
257/*}}}*/
258
259/*Object functions*/
260/*FUNCTION TransientInput::GetInputValue(double* pvalue,GaussTria* gauss){{{1*/
261void TransientInput::GetInputValue(double* pvalue,GaussTria* gauss){
262
263 double time;
264
265 /*First, recover current time from parameters: */
266 this->parameters->FindParam(&time,TimeEnum);
267
268 /*Retrieve interpolated values for this time step: */
269 Input* input=GetTimeInput(time);
270
271 /*Call input function*/
272 input->GetInputValue(pvalue,gauss);
273
274 delete input;
275
276}
277/*}}}*/
278/*FUNCTION TransientInput::GetInputDerivativeValue(double* p, double* xyz_list, GaussTria* gauss){{{1*/
279void TransientInput::GetInputDerivativeValue(double* p, double* xyz_list, GaussTria* gauss){
280
281 double time;
282
283 /*First, recover current time from parameters: */
284 parameters->FindParam(&time,TimeEnum);
285
286 /*Retrieve interpolated values for this time step: */
287 Input* input=GetTimeInput(time);
288
289 /*Call input function*/
290 input->GetInputDerivativeValue(p,xyz_list,gauss);
291
292 delete input;
293
294}
295/*}}}*/
296/*FUNCTION TransientInput::ChangeEnum{{{1*/
297void TransientInput::ChangeEnum(int newenumtype){
298 this->enum_type=newenumtype;
299}
300/*}}}*/
301/*FUNCTION TransientInput::GetInputAverage{{{1*/
302void TransientInput::GetInputAverage(double* pvalue){
303
304 double time;
305
306 /*First, recover current time from parameters: */
307 parameters->FindParam(&time,TimeEnum);
308
309 /*Retrieve interpolated values for this time step: */
310 Input* input=GetTimeInput(time);
311
312 /*Call input function*/
313 input->GetInputAverage(pvalue);
314
315 delete input;
316
317}
318/*}}}*/
319
320/*Intermediary*/
321/*FUNCTION TransientInput::SquareMin{{{1*/
322void TransientInput::SquareMin(double* psquaremin, bool process_units,Parameters* parameters){
323
324 double time;
325
326 /*First, recover current time from parameters: */
327 parameters->FindParam(&time,TimeEnum);
328
329 /*Retrieve interpolated values for this time step: */
330 Input* input=GetTimeInput(time);
331
332 /*Call input function*/
333 input->SquareMin(psquaremin,process_units,parameters);
334
335 delete input;
336
337}
338/*}}}*/
339/*FUNCTION TransientInput::InfinityNorm{{{1*/
340double TransientInput::InfinityNorm(void){
341
342 double time;
343 double infnorm;
344
345 /*First, recover current time from parameters: */
346 parameters->FindParam(&time,TimeEnum);
347
348 /*Retrieve interpolated values for this time step: */
349 Input* input=GetTimeInput(time);
350
351 /*Call input function*/
352 infnorm=input->InfinityNorm();
353
354 /*Clean-up and return*/
355 delete input;
356 return infnorm;
357}
358/*}}}*/
359/*FUNCTION TransientInput::Max{{{1*/
360double TransientInput::Max(void){
361
362 double time;
363 double max;
364
365 /*First, recover current time from parameters: */
366 parameters->FindParam(&time,TimeEnum);
367
368 /*Retrieve interpolated values for this time step: */
369 Input* input=GetTimeInput(time);
370
371 /*Call input function*/
372 max=input->Max();
373
374 delete input;
375
376 return max;
377}
378/*}}}*/
379/*FUNCTION TransientInput::MaxAbs{{{1*/
380double TransientInput::MaxAbs(void){
381
382 double time;
383 double maxabs;
384
385 /*First, recover current time from parameters: */
386 parameters->FindParam(&time,TimeEnum);
387
388 /*Retrieve interpolated values for this time step: */
389 Input* input=GetTimeInput(time);
390
391 /*Call input function*/
392 maxabs=input->MaxAbs();
393
394 /*Clean-up and return*/
395 delete input;
396 return maxabs;
397
398}
399/*}}}*/
400/*FUNCTION TransientInput::Min{{{1*/
401double TransientInput::Min(void){
402
403 double time;
404 double min;
405
406 /*First, recover current time from parameters: */
407 parameters->FindParam(&time,TimeEnum);
408
409 /*Retrieve interpolated values for this time step: */
410 Input* input=GetTimeInput(time);
411
412 /*Call input function*/
413 min=input->Min();
414
415 /*Clean-up and return*/
416 delete input;
417 return min;
418
419}
420/*}}}*/
421/*FUNCTION TransientInput::MinAbs{{{1*/
422double TransientInput::MinAbs(void){
423
424 double time;
425 double minabs;
426
427 /*First, recover current time from parameters: */
428 parameters->FindParam(&time,TimeEnum);
429
430 /*Retrieve interpolated values for this time step: */
431 Input* input=GetTimeInput(time);
432
433 /*Call input function*/
434 minabs=input->MinAbs();
435
436 /*Clean-up and return*/
437 delete input;
438 return minabs;
439}
440/*}}}*/
441/*FUNCTION TransientInput::GetVectorFromInputs{{{1*/
442void TransientInput::GetVectorFromInputs(Vec vector,int* doflist){
443
444 double time;
445
446 /*First, recover current time from parameters: */
447 parameters->FindParam(&time,TimeEnum);
448
449 /*Retrieve interpolated values for this time step: */
450 Input* input=GetTimeInput(time);
451
452 /*Call input function*/
453 input->GetVectorFromInputs(vector,doflist);
454
455 delete input;
456
457} /*}}}*/
458/*FUNCTION TransientInput::GetTimeInput{{{1*/
459Input* TransientInput::GetTimeInput(double intime){
460
461 int i,j;
462 double deltat;
463 double alpha1,alpha2;
464 bool found=false;
465 Input* input=NULL;
466 Input* input1=NULL;
467 Input* input2=NULL;
468
469 /*Ok, we have the time, go through the timesteps, and figure out which interval we
470 *fall within. Then interpolate the values on this interval: */
471 if(intime<this->timesteps[0]){
472 /*get values for the first time: */
473 input=(Input*)((Input*)this->inputs->GetObjectByOffset(0))->copy();
474 found=true;
475 }
476 else if(intime>this->timesteps[this->numtimesteps-1]){
477 /*get values for the last time: */
478 input=(Input*)((Input*)this->inputs->GetObjectByOffset(numtimesteps-1))->copy();
479 found=true;
480 }
481 else{
482 /*Find which interval we fall within: */
483 for(i=0;i<this->numtimesteps;i++){
484 if(intime==this->timesteps[i]){
485 /*We are right on one step time: */
486 input=(Input*)((Input*)this->inputs->GetObjectByOffset(i))->copy();
487 found=true;
488 break; //we are done with the time interpolation.
489 }
490 else{
491 if(this->timesteps[i]<intime && intime<this->timesteps[i+1]){
492 /*ok, we have the interval ]i:i+1[. Interpolate linearly for now: */
493 deltat=this->timesteps[i+1]-this->timesteps[i];
494 alpha2=(intime-this->timesteps[i])/deltat;
495 alpha1=(1-alpha2);
496
497 input1=(Input*)this->inputs->GetObjectByOffset(i);
498 input2=(Input*)this->inputs->GetObjectByOffset(i+1);
499
500 input=(Input*)input1->copy();
501 input->Scale(alpha1);
502 input->AXPY(input2,alpha2);
503
504 found=true;
505 break;
506 }
507 else continue; //keep looking on the next interval
508 }
509 }
510 }
511 if(!found)_error_("did not find time interval on which to interpolate forcing values!");
512
513 /*Assign output pointer*/
514 return input;
515}
516/*}}}*/
517/*FUNCTION TransientInput::Configure{{{1*/
518void TransientInput::Configure(Parameters* parameters){
519 this->parameters=parameters;
520}
521/*}}}*/
Note: See TracBrowser for help on using the repository browser.