[8757] | 1 | /*!\file TransientInput.c
|
---|
| 2 | * \brief: implementation of the TransientInput object
|
---|
| 3 | */
|
---|
| 4 | /*Headers{{{1*/
|
---|
| 5 | #ifdef HAVE_CONFIG_H
|
---|
[9320] | 6 | #include <config.h>
|
---|
[8757] | 7 | #else
|
---|
| 8 | #error "Cannot compile with HAVE_CONFIG_H symbol! run configure first!"
|
---|
| 9 | #endif
|
---|
| 10 |
|
---|
[9320] | 11 | #include <stdio.h>
|
---|
[8757] | 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*/
|
---|
| 22 | TransientInput::TransientInput(){
|
---|
[8791] | 23 |
|
---|
[8967] | 24 | enum_type=UNDEF;
|
---|
[8757] | 25 | inputs=NULL;
|
---|
[8791] | 26 | this->numtimesteps=0;
|
---|
| 27 | this->parameters=NULL;
|
---|
| 28 | this->timesteps=NULL;
|
---|
| 29 |
|
---|
[8757] | 30 | }
|
---|
| 31 | /*}}}*/
|
---|
| 32 | /*FUNCTION TransientInput::TransientInput(int in_enum_type){{{1*/
|
---|
| 33 | TransientInput::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*/
|
---|
| 47 | TransientInput::~TransientInput(){
|
---|
| 48 | xfree((void**)&this->timesteps);
|
---|
[8791] | 49 | this->timesteps=NULL;
|
---|
[8757] | 50 | this->numtimesteps=0;
|
---|
| 51 | parameters=NULL;
|
---|
| 52 | delete this->inputs;
|
---|
| 53 | return;
|
---|
| 54 | }
|
---|
| 55 | /*}}}*/
|
---|
[8967] | 56 | /*FUNCTION void TransientInput::AddTimeInput(Input* input,double time){{{1*/
|
---|
[8757] | 57 | void 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*/
|
---|
| 88 | void TransientInput::Echo(void){
|
---|
| 89 | this->DeepEcho();
|
---|
| 90 | }
|
---|
| 91 | /*}}}*/
|
---|
| 92 | /*FUNCTION TransientInput::DeepEcho{{{1*/
|
---|
| 93 | void 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*/
|
---|
| 108 | int TransientInput::Id(void){ return -1; }
|
---|
| 109 | /*}}}*/
|
---|
| 110 | /*FUNCTION TransientInput::MyRank{{{1*/
|
---|
| 111 | int TransientInput::MyRank(void){
|
---|
| 112 | extern int my_rank;
|
---|
| 113 | return my_rank;
|
---|
| 114 | }
|
---|
| 115 | /*}}}*/
|
---|
[9883] | 116 | /*FUNCTION TransientInput::ObjectEnum{{{1*/
|
---|
| 117 | int TransientInput::ObjectEnum(void){
|
---|
[8757] | 118 |
|
---|
| 119 | return TransientInputEnum;
|
---|
| 120 |
|
---|
| 121 | }
|
---|
| 122 | /*}}}*/
|
---|
| 123 | /*FUNCTION TransientInput::copy{{{1*/
|
---|
| 124 | Object* TransientInput::copy() {
|
---|
| 125 |
|
---|
| 126 | TransientInput* output=NULL;
|
---|
| 127 |
|
---|
| 128 | output = new TransientInput();
|
---|
| 129 | output->enum_type=this->enum_type;
|
---|
| 130 | output->numtimesteps=this->numtimesteps;
|
---|
[8791] | 131 | output->timesteps=(double*)xmalloc(this->numtimesteps*sizeof(double));
|
---|
[8757] | 132 | memcpy(output->timesteps,this->timesteps,this->numtimesteps*sizeof(double));
|
---|
| 133 | output->inputs=(Inputs*)this->inputs->Copy();
|
---|
| 134 | output->parameters=this->parameters;
|
---|
| 135 |
|
---|
| 136 | return output;
|
---|
| 137 |
|
---|
| 138 | }
|
---|
| 139 | /*}}}*/
|
---|
| 140 |
|
---|
| 141 | /*TransientInput management*/
|
---|
[9883] | 142 | /*FUNCTION TransientInput::InstanceEnum{{{1*/
|
---|
| 143 | int TransientInput::InstanceEnum(void){
|
---|
[8757] | 144 |
|
---|
| 145 | return this->enum_type;
|
---|
| 146 |
|
---|
| 147 | }
|
---|
| 148 | /*}}}*/
|
---|
| 149 | /*FUNCTION TransientInput::SpawnTriaInput{{{1*/
|
---|
| 150 | Input* TransientInput::SpawnTriaInput(int* indices){
|
---|
| 151 |
|
---|
| 152 | /*output*/
|
---|
| 153 | TransientInput* outinput=NULL;
|
---|
| 154 |
|
---|
| 155 | /*Create new Transientinput (copy of current input)*/
|
---|
| 156 | outinput=new TransientInput();
|
---|
| 157 | outinput->enum_type=this->enum_type;
|
---|
[8791] | 158 | outinput->numtimesteps=this->numtimesteps;
|
---|
| 159 | outinput->timesteps=(double*)xmalloc(this->numtimesteps*sizeof(double));
|
---|
| 160 | memcpy(outinput->timesteps,this->timesteps,this->numtimesteps*sizeof(double));
|
---|
[8757] | 161 | outinput->inputs=(Inputs*)this->inputs->SpawnTriaInputs(indices);
|
---|
[8791] | 162 | outinput->parameters=this->parameters;
|
---|
[8757] | 163 |
|
---|
| 164 | /*Assign output*/
|
---|
| 165 | return outinput;
|
---|
| 166 |
|
---|
| 167 | }
|
---|
| 168 | /*}}}*/
|
---|
| 169 | /*FUNCTION TransientInput::SpawnResult{{{1*/
|
---|
| 170 | ElementResult* TransientInput::SpawnResult(int step, double time){
|
---|
| 171 |
|
---|
| 172 | ElementResult* elementresult=NULL;
|
---|
| 173 |
|
---|
| 174 | /*Ok, we want to spawn an ElementResult. We have the time, just get
|
---|
| 175 | *the correct values: */
|
---|
| 176 | Input* input=GetTimeInput(time);
|
---|
| 177 |
|
---|
| 178 | elementresult=input->SpawnResult(step,time);
|
---|
| 179 |
|
---|
[8791] | 180 | delete input;
|
---|
[8757] | 181 |
|
---|
| 182 | return elementresult;
|
---|
| 183 | }
|
---|
| 184 | /*}}}*/
|
---|
| 185 |
|
---|
| 186 | /*Object functions*/
|
---|
[10135] | 187 | /*FUNCTION TransientInput::GetInputValue(double* pvalue,GaussTria* gauss){{{1*/
|
---|
| 188 | void TransientInput::GetInputValue(double* pvalue,GaussTria* gauss){
|
---|
[8757] | 189 |
|
---|
| 190 | double time;
|
---|
| 191 |
|
---|
| 192 | /*First, recover current time from parameters: */
|
---|
[8791] | 193 | this->parameters->FindParam(&time,TimeEnum);
|
---|
[8757] | 194 |
|
---|
| 195 | /*Retrieve interpolated values for this time step: */
|
---|
| 196 | Input* input=GetTimeInput(time);
|
---|
| 197 |
|
---|
| 198 | /*Call input function*/
|
---|
[10135] | 199 | input->GetInputValue(pvalue,gauss);
|
---|
[8757] | 200 |
|
---|
[8791] | 201 | delete input;
|
---|
[8757] | 202 |
|
---|
| 203 | }
|
---|
| 204 | /*}}}*/
|
---|
[10135] | 205 | /*FUNCTION TransientInput::GetInputDerivativeValue(double* p, double* xyz_list, GaussTria* gauss){{{1*/
|
---|
| 206 | void TransientInput::GetInputDerivativeValue(double* p, double* xyz_list, GaussTria* gauss){
|
---|
[8757] | 207 |
|
---|
| 208 | double time;
|
---|
| 209 |
|
---|
| 210 | /*First, recover current time from parameters: */
|
---|
| 211 | parameters->FindParam(&time,TimeEnum);
|
---|
| 212 |
|
---|
| 213 | /*Retrieve interpolated values for this time step: */
|
---|
| 214 | Input* input=GetTimeInput(time);
|
---|
| 215 |
|
---|
| 216 | /*Call input function*/
|
---|
[10135] | 217 | input->GetInputDerivativeValue(p,xyz_list,gauss);
|
---|
[8757] | 218 |
|
---|
[8791] | 219 | delete input;
|
---|
| 220 |
|
---|
[8757] | 221 | }
|
---|
| 222 | /*}}}*/
|
---|
| 223 | /*FUNCTION TransientInput::ChangeEnum{{{1*/
|
---|
| 224 | void TransientInput::ChangeEnum(int newenumtype){
|
---|
| 225 | this->enum_type=newenumtype;
|
---|
| 226 | }
|
---|
| 227 | /*}}}*/
|
---|
[10135] | 228 | /*FUNCTION TransientInput::GetInputAverage{{{1*/
|
---|
| 229 | void TransientInput::GetInputAverage(double* pvalue){
|
---|
[8757] | 230 |
|
---|
| 231 | double time;
|
---|
| 232 |
|
---|
| 233 | /*First, recover current time from parameters: */
|
---|
| 234 | parameters->FindParam(&time,TimeEnum);
|
---|
| 235 |
|
---|
| 236 | /*Retrieve interpolated values for this time step: */
|
---|
| 237 | Input* input=GetTimeInput(time);
|
---|
| 238 |
|
---|
| 239 | /*Call input function*/
|
---|
[10135] | 240 | input->GetInputAverage(pvalue);
|
---|
[8757] | 241 |
|
---|
[8791] | 242 | delete input;
|
---|
[8757] | 243 |
|
---|
| 244 | }
|
---|
| 245 | /*}}}*/
|
---|
| 246 |
|
---|
| 247 | /*Intermediary*/
|
---|
| 248 | /*FUNCTION TransientInput::SquareMin{{{1*/
|
---|
| 249 | void TransientInput::SquareMin(double* psquaremin, bool process_units,Parameters* parameters){
|
---|
| 250 |
|
---|
| 251 | double time;
|
---|
| 252 |
|
---|
| 253 | /*First, recover current time from parameters: */
|
---|
| 254 | parameters->FindParam(&time,TimeEnum);
|
---|
| 255 |
|
---|
| 256 | /*Retrieve interpolated values for this time step: */
|
---|
| 257 | Input* input=GetTimeInput(time);
|
---|
| 258 |
|
---|
| 259 | /*Call input function*/
|
---|
| 260 | input->SquareMin(psquaremin,process_units,parameters);
|
---|
| 261 |
|
---|
[8791] | 262 | delete input;
|
---|
[8757] | 263 |
|
---|
| 264 | }
|
---|
| 265 | /*}}}*/
|
---|
| 266 | /*FUNCTION TransientInput::InfinityNorm{{{1*/
|
---|
| 267 | double TransientInput::InfinityNorm(void){
|
---|
| 268 |
|
---|
| 269 | double time;
|
---|
| 270 | double infnorm;
|
---|
| 271 |
|
---|
| 272 | /*First, recover current time from parameters: */
|
---|
| 273 | parameters->FindParam(&time,TimeEnum);
|
---|
| 274 |
|
---|
| 275 | /*Retrieve interpolated values for this time step: */
|
---|
| 276 | Input* input=GetTimeInput(time);
|
---|
| 277 |
|
---|
| 278 | /*Call input function*/
|
---|
| 279 | infnorm=input->InfinityNorm();
|
---|
| 280 |
|
---|
[11197] | 281 | /*Clean-up and return*/
|
---|
[8791] | 282 | delete input;
|
---|
[11197] | 283 | return infnorm;
|
---|
[8757] | 284 | }
|
---|
| 285 | /*}}}*/
|
---|
| 286 | /*FUNCTION TransientInput::Max{{{1*/
|
---|
| 287 | double TransientInput::Max(void){
|
---|
| 288 |
|
---|
| 289 | double time;
|
---|
| 290 | double max;
|
---|
| 291 |
|
---|
| 292 | /*First, recover current time from parameters: */
|
---|
| 293 | parameters->FindParam(&time,TimeEnum);
|
---|
| 294 |
|
---|
| 295 | /*Retrieve interpolated values for this time step: */
|
---|
| 296 | Input* input=GetTimeInput(time);
|
---|
| 297 |
|
---|
| 298 | /*Call input function*/
|
---|
| 299 | max=input->Max();
|
---|
| 300 |
|
---|
[8791] | 301 | delete input;
|
---|
[8757] | 302 |
|
---|
| 303 | return max;
|
---|
| 304 | }
|
---|
| 305 | /*}}}*/
|
---|
| 306 | /*FUNCTION TransientInput::MaxAbs{{{1*/
|
---|
| 307 | double TransientInput::MaxAbs(void){
|
---|
| 308 |
|
---|
| 309 | double time;
|
---|
| 310 | double maxabs;
|
---|
| 311 |
|
---|
| 312 | /*First, recover current time from parameters: */
|
---|
| 313 | parameters->FindParam(&time,TimeEnum);
|
---|
| 314 |
|
---|
| 315 | /*Retrieve interpolated values for this time step: */
|
---|
| 316 | Input* input=GetTimeInput(time);
|
---|
| 317 |
|
---|
| 318 | /*Call input function*/
|
---|
| 319 | maxabs=input->MaxAbs();
|
---|
| 320 |
|
---|
[11197] | 321 | /*Clean-up and return*/
|
---|
[8791] | 322 | delete input;
|
---|
[8757] | 323 | return maxabs;
|
---|
[8791] | 324 |
|
---|
[8757] | 325 | }
|
---|
| 326 | /*}}}*/
|
---|
| 327 | /*FUNCTION TransientInput::Min{{{1*/
|
---|
| 328 | double TransientInput::Min(void){
|
---|
| 329 |
|
---|
| 330 | double time;
|
---|
| 331 | double min;
|
---|
| 332 |
|
---|
| 333 | /*First, recover current time from parameters: */
|
---|
| 334 | parameters->FindParam(&time,TimeEnum);
|
---|
| 335 |
|
---|
| 336 | /*Retrieve interpolated values for this time step: */
|
---|
| 337 | Input* input=GetTimeInput(time);
|
---|
| 338 |
|
---|
| 339 | /*Call input function*/
|
---|
| 340 | min=input->Min();
|
---|
| 341 |
|
---|
[11197] | 342 | /*Clean-up and return*/
|
---|
[8791] | 343 | delete input;
|
---|
[8757] | 344 | return min;
|
---|
[8791] | 345 |
|
---|
[8757] | 346 | }
|
---|
| 347 | /*}}}*/
|
---|
| 348 | /*FUNCTION TransientInput::MinAbs{{{1*/
|
---|
| 349 | double TransientInput::MinAbs(void){
|
---|
| 350 |
|
---|
| 351 | double time;
|
---|
| 352 | double minabs;
|
---|
| 353 |
|
---|
| 354 | /*First, recover current time from parameters: */
|
---|
| 355 | parameters->FindParam(&time,TimeEnum);
|
---|
| 356 |
|
---|
| 357 | /*Retrieve interpolated values for this time step: */
|
---|
| 358 | Input* input=GetTimeInput(time);
|
---|
| 359 |
|
---|
| 360 | /*Call input function*/
|
---|
| 361 | minabs=input->MinAbs();
|
---|
| 362 |
|
---|
[11197] | 363 | /*Clean-up and return*/
|
---|
[8791] | 364 | delete input;
|
---|
[8757] | 365 | return minabs;
|
---|
| 366 | }
|
---|
| 367 | /*}}}*/
|
---|
| 368 | /*FUNCTION TransientInput::GetVectorFromInputs{{{1*/
|
---|
[11695] | 369 | void TransientInput::GetVectorFromInputs(Vector* vector,int* doflist){
|
---|
[8757] | 370 |
|
---|
| 371 | double time;
|
---|
| 372 |
|
---|
| 373 | /*First, recover current time from parameters: */
|
---|
| 374 | parameters->FindParam(&time,TimeEnum);
|
---|
| 375 |
|
---|
| 376 | /*Retrieve interpolated values for this time step: */
|
---|
| 377 | Input* input=GetTimeInput(time);
|
---|
| 378 |
|
---|
| 379 | /*Call input function*/
|
---|
| 380 | input->GetVectorFromInputs(vector,doflist);
|
---|
| 381 |
|
---|
[8791] | 382 | delete input;
|
---|
[8757] | 383 |
|
---|
| 384 | } /*}}}*/
|
---|
| 385 | /*FUNCTION TransientInput::GetTimeInput{{{1*/
|
---|
| 386 | Input* TransientInput::GetTimeInput(double intime){
|
---|
| 387 |
|
---|
| 388 | int i,j;
|
---|
| 389 | double deltat;
|
---|
| 390 | double alpha1,alpha2;
|
---|
| 391 | bool found=false;
|
---|
| 392 | Input* input=NULL;
|
---|
| 393 | Input* input1=NULL;
|
---|
| 394 | Input* input2=NULL;
|
---|
| 395 |
|
---|
| 396 | /*Ok, we have the time, go through the timesteps, and figure out which interval we
|
---|
| 397 | *fall within. Then interpolate the values on this interval: */
|
---|
| 398 | if(intime<this->timesteps[0]){
|
---|
| 399 | /*get values for the first time: */
|
---|
| 400 | input=(Input*)((Input*)this->inputs->GetObjectByOffset(0))->copy();
|
---|
| 401 | found=true;
|
---|
| 402 | }
|
---|
| 403 | else if(intime>this->timesteps[this->numtimesteps-1]){
|
---|
| 404 | /*get values for the last time: */
|
---|
| 405 | input=(Input*)((Input*)this->inputs->GetObjectByOffset(numtimesteps-1))->copy();
|
---|
| 406 | found=true;
|
---|
| 407 | }
|
---|
| 408 | else{
|
---|
| 409 | /*Find which interval we fall within: */
|
---|
| 410 | for(i=0;i<this->numtimesteps;i++){
|
---|
| 411 | if(intime==this->timesteps[i]){
|
---|
| 412 | /*We are right on one step time: */
|
---|
| 413 | input=(Input*)((Input*)this->inputs->GetObjectByOffset(i))->copy();
|
---|
| 414 | found=true;
|
---|
| 415 | break; //we are done with the time interpolation.
|
---|
| 416 | }
|
---|
| 417 | else{
|
---|
| 418 | if(this->timesteps[i]<intime && intime<this->timesteps[i+1]){
|
---|
| 419 | /*ok, we have the interval ]i:i+1[. Interpolate linearly for now: */
|
---|
| 420 | deltat=this->timesteps[i+1]-this->timesteps[i];
|
---|
| 421 | alpha2=(intime-this->timesteps[i])/deltat;
|
---|
| 422 | alpha1=(1-alpha2);
|
---|
| 423 |
|
---|
| 424 | input1=(Input*)this->inputs->GetObjectByOffset(i);
|
---|
| 425 | input2=(Input*)this->inputs->GetObjectByOffset(i+1);
|
---|
| 426 |
|
---|
| 427 | input=(Input*)input1->copy();
|
---|
| 428 | input->Scale(alpha1);
|
---|
| 429 | input->AXPY(input2,alpha2);
|
---|
| 430 |
|
---|
| 431 | found=true;
|
---|
| 432 | break;
|
---|
| 433 | }
|
---|
| 434 | else continue; //keep looking on the next interval
|
---|
| 435 | }
|
---|
| 436 | }
|
---|
| 437 | }
|
---|
| 438 | if(!found)_error_("did not find time interval on which to interpolate forcing values!");
|
---|
| 439 |
|
---|
| 440 | /*Assign output pointer*/
|
---|
| 441 | return input;
|
---|
| 442 | }
|
---|
| 443 | /*}}}*/
|
---|
| 444 | /*FUNCTION TransientInput::Configure{{{1*/
|
---|
| 445 | void TransientInput::Configure(Parameters* parameters){
|
---|
| 446 | this->parameters=parameters;
|
---|
| 447 | }
|
---|
| 448 | /*}}}*/
|
---|