AdaptInstationary.cc 13.2 KB
Newer Older
1 2 3 4 5 6 7 8 9 10 11 12
//
// Software License for AMDiS
//
// Copyright (c) 2010 Dresden University of Technology 
// All rights reserved.
// Authors: Simon Vey, Thomas Witkowski et al.
//
// This file is part of AMDiS
//
// See also license.opensource.txt in the distribution.


13
#include "AdaptInstationary.h"
14
#include "Initfile.h"
15
#include "est/Estimator.h"
16 17 18 19
#include "ProblemIterationInterface.h"
#include "ProblemTimeInterface.h"
#include "Serializer.h"

20 21
#if HAVE_PARALLEL_DOMAIN_AMDIS
#include "parallel/MeshDistributor.h"
Naumann, Andreas's avatar
Naumann, Andreas committed
22
#ifndef HAVE_PARALLEL_MTL4
23 24
#include <petsc.h>
#endif
Naumann, Andreas's avatar
Naumann, Andreas committed
25
#endif
26

27 28
namespace AMDiS {

Thomas Witkowski's avatar
Thomas Witkowski committed
29
  AdaptInstationary::AdaptInstationary(std::string name,
30
				       ProblemIterationInterface *problemStat,  
31 32 33
				       AdaptInfo *info,
				       ProblemTimeInterface *problemInstat,
				       AdaptInfo *initialInfo,
34
				       time_t initialTimestampSet)
35
    : AdaptBase(name, problemStat, info, problemInstat, initialInfo),
36 37
      breakWhenStable(0),
      dbgMode(false)
38 39 40
  {
    FUNCNAME("AdaptInstationary::AdaptInstationary()");

41 42
    MSG("You make use of the obsolete constructor AdaptInstationary::AdaptInstationary(...)!\n");
    MSG("Please use the constructor that uses references instead of pointers!\n");
43

44 45
    initConstructor(problemStat, info, initialInfo, initialTimestampSet);
 }
46

Thomas Witkowski's avatar
Thomas Witkowski committed
47

48 49 50 51 52
  AdaptInstationary::AdaptInstationary(std::string name,
				       ProblemIterationInterface &problemStat,  
				       AdaptInfo &info,
				       ProblemTimeInterface &problemInstat,
				       AdaptInfo &initialInfo,
53
				       time_t initialTimestampSet)
54 55 56 57 58 59
    : AdaptBase(name, &problemStat, &info, &problemInstat, &initialInfo),
      breakWhenStable(0),
      dbgMode(false)
  {
    FUNCNAME("AdaptInstationary::AdaptInstationary()");

60
    initConstructor(&problemStat, &info, &initialInfo, initialTimestampSet);
61 62
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
63

64 65 66
  void AdaptInstationary::initConstructor(ProblemIterationInterface *problemStat,  
					  AdaptInfo *info,
					  AdaptInfo *initialInfo,
67
					  time_t initialTimestampSet)
68
  {
Thomas Witkowski's avatar
Thomas Witkowski committed
69
    initialize(name);
70

Thomas Witkowski's avatar
Thomas Witkowski committed
71
    fixedTimestep = (info->getMinTimestep() == info->getMaxTimestep());
72 73
 
    if (initialTimestampSet == 0)
Thomas Witkowski's avatar
Thomas Witkowski committed
74
      initialTimestamp = time(NULL);
75
    else
76 77
      initialTimestamp = initialTimestampSet;

78
    // Check if the problem should be deserialized because of the -rs parameter.
Thomas Witkowski's avatar
Thomas Witkowski committed
79
    std::string serializationFilename = "";
80
    Parameters::get("argv->rs", serializationFilename);
81 82 83 84

    if (serializationFilename.compare("")) {
      // The value of the -rs argument is ignored, because we want to use the 
      // serialization file mentioned in the used init file.
Thomas Witkowski's avatar
Thomas Witkowski committed
85
      MSG("Deserialization from file: %s\n", queueSerializationFilename.c_str());
86

87
      std::ifstream in(queueSerializationFilename.c_str() , ios::in);
88 89 90 91 92 93

      // Read the revision number of the AMDiS version which was used to create 
      // the serialization file.
      int revNumber = -1;
      SerUtil::deserialize(in, revNumber);
      
94 95 96
      deserialize(in);
      in.close();

97 98
      info->setIsDeserialized(true);
      initialInfo->setIsDeserialized(true);
99 100 101 102
    } else {
      int readSerialization = 0;
      int readSerializationWithAdaptInfo = 0;

103 104 105 106
      Parameters::get((*problemStat).getName() + "->input->read serialization",
		      readSerialization);
      Parameters::get((*problemStat).getName() + "->input->serialization with adaptinfo",
		      readSerializationWithAdaptInfo);
107 108

      if (readSerialization && readSerializationWithAdaptInfo) {
Thomas Witkowski's avatar
Thomas Witkowski committed
109
	std::string serializationFilename = "";
110

111 112
	Parameters::get((*problemStat).getName() + "->input->serialization filename", 
			serializationFilename);
113 114 115
	TEST_EXIT(serializationFilename != "")("no serialization file\n");

	MSG("Deserialization with AdaptInfo from file: %s\n", serializationFilename.c_str());
Thomas Witkowski's avatar
Thomas Witkowski committed
116
	std::ifstream in(serializationFilename.c_str());
117 118 119 120 121 122

	// Read the revision number of the AMDiS version which was used to create 
	// the serialization file.
	int revNumber = -1;
	SerUtil::deserialize(in, revNumber);

123 124 125 126 127 128
	deserialize(in);
	in.close();
      }
    }
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
129

130 131 132 133 134
  void AdaptInstationary::explicitTimeStrategy()
  {
    FUNCNAME("AdaptInstationary::explicitTimeStrategy()");

    // estimate before first adaption
Thomas Witkowski's avatar
Thomas Witkowski committed
135
    if (adaptInfo->getTime() <= adaptInfo->getStartTime())
Thomas Witkowski's avatar
Thomas Witkowski committed
136
      problemIteration->oneIteration(adaptInfo, ESTIMATE);
Thomas Witkowski's avatar
Thomas Witkowski committed
137 138


139
    // increment time
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
140
    adaptInfo->setTime(adaptInfo->getTime() + adaptInfo->getTimestep());
141

Thomas Witkowski's avatar
Thomas Witkowski committed
142
    problemTime->setTime(adaptInfo);
143

144
    INFO(info, 6)("time = %e, timestep = %e\n",
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
145
		  adaptInfo->getTime(), adaptInfo->getTimestep());
146

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
147
    adaptInfo->setSpaceIteration(0);
148 149
  
    // do the iteration
Thomas Witkowski's avatar
Thomas Witkowski committed
150 151 152
    problemIteration->beginIteration(adaptInfo);
    problemIteration->oneIteration(adaptInfo, FULL_ITERATION);
    problemIteration->endIteration(adaptInfo);
153
    adaptInfo->setLastProcessedTimestep(adaptInfo->getTimestep()); 
154 155
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
156

157 158 159 160 161
  void AdaptInstationary::implicitTimeStrategy()
  {
    FUNCNAME("AdaptInstationary::implicitTimeStrategy()");

    do {
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
162
      adaptInfo->setTime(adaptInfo->getTime() + adaptInfo->getTimestep());
Thomas Witkowski's avatar
Thomas Witkowski committed
163
      problemTime->setTime(adaptInfo);
164

Thomas Witkowski's avatar
Thomas Witkowski committed
165 166
      INFO(info,6)("time = %e, try timestep = %e\n",
		   adaptInfo->getTime(), adaptInfo->getTimestep());
167

Thomas Witkowski's avatar
Thomas Witkowski committed
168
      problemIteration->oneIteration(adaptInfo, NO_ADAPTION);
169

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
170
      adaptInfo->incTimestepIteration();
171

Thomas Witkowski's avatar
Thomas Witkowski committed
172
      if (!fixedTimestep && 
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
173
	  !adaptInfo->timeToleranceReached() &&
174
	  adaptInfo->getTimestepIteration() <= adaptInfo->getMaxTimestepIteration() &&
Thomas Witkowski's avatar
Thomas Witkowski committed
175
	  !(adaptInfo->getTimestep() <= adaptInfo->getMinTimestep())) {
176 177 178 179
  
	adaptInfo->setTime(adaptInfo->getTime() - adaptInfo->getTimestep());
	adaptInfo->setTimestep(adaptInfo->getTimestep() * timeDelta1);
	continue;
180
      }
181

182

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
183
      adaptInfo->setSpaceIteration(0);
184 185


186
      // === Do only space iterations only if the maximum is higher than 0. === 
187

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
188
      if (adaptInfo->getMaxSpaceIteration() > 0) {
189
    
190
	// === Space iterations. ===
191
	do {
Thomas Witkowski's avatar
Thomas Witkowski committed
192
	  problemIteration->beginIteration(adaptInfo);
193
	  
Thomas Witkowski's avatar
Thomas Witkowski committed
194 195
	  if (problemIteration->oneIteration(adaptInfo, FULL_ITERATION)) {
	    if (!fixedTimestep && 
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
196
		!adaptInfo->timeToleranceReached() &&
Thomas Witkowski's avatar
Thomas Witkowski committed
197 198
		!(adaptInfo->getTimestep() <= adaptInfo->getMinTimestep())) {
	      adaptInfo->setTime(adaptInfo->getTime() - adaptInfo->getTimestep());
Thomas Witkowski's avatar
Thomas Witkowski committed
199 200
	      adaptInfo->setTimestep(adaptInfo->getTimestep() * timeDelta2);
	      problemIteration->endIteration(adaptInfo);
Thomas Witkowski's avatar
Thomas Witkowski committed
201 202 203
	      adaptInfo->incSpaceIteration();
	      break;
	    }	
204 205
	  }

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
206
	  adaptInfo->incSpaceIteration();
Thomas Witkowski's avatar
Thomas Witkowski committed
207
	  problemIteration->endIteration(adaptInfo);
208
	  
Thomas Witkowski's avatar
Thomas Witkowski committed
209 210
	} while (!adaptInfo->spaceToleranceReached() && 
		 adaptInfo->getSpaceIteration() <= adaptInfo->getMaxSpaceIteration());
211 212

      } else {
Thomas Witkowski's avatar
Thomas Witkowski committed
213
	problemIteration->endIteration(adaptInfo);
214 215
      }

216

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
217
    } while(!adaptInfo->timeToleranceReached() &&
Thomas Witkowski's avatar
Thomas Witkowski committed
218
	    !(adaptInfo->getTimestep() <= adaptInfo->getMinTimestep()) && 
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
219
	    adaptInfo->getTimestepIteration() <= adaptInfo->getMaxTimestepIteration());  
220

221
    adaptInfo->setLastProcessedTimestep(adaptInfo->getTimestep()); 
222 223 224 225

    // After successful iteration/timestep the timestep will be changed according 
    // adaption rules for next timestep. 
    // First, check for increase of timestep
Thomas Witkowski's avatar
Thomas Witkowski committed
226 227
    if (!fixedTimestep && adaptInfo->timeErrorLow()) {
      adaptInfo->setTimestep(adaptInfo->getTimestep() * timeDelta2);
228 229 230 231 232 233 234
      if (dbgMode) {
	// print information about timestep increase
      }
    } else {
      if (dbgMode) {
	std::cout << "=== ADAPT INFO DEBUG MODE ===\n";
	std::cout << " Do not increase timestep: \n";
Thomas Witkowski's avatar
Thomas Witkowski committed
235
	if (fixedTimestep)
236
	  std::cout << "   fixedTimestep = true\n";	
Thomas Witkowski's avatar
Thomas Witkowski committed
237
	if (!adaptInfo->timeErrorLow())
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
238
	  adaptInfo->printTimeErrorLowInfo();
239
      }
240
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
241

242
    // Second, check for decrease of timestep
Thomas Witkowski's avatar
Thomas Witkowski committed
243
    if (!fixedTimestep &&
244
	!adaptInfo->timeToleranceReached() &&
Thomas Witkowski's avatar
Thomas Witkowski committed
245 246
	!(adaptInfo->getTimestep() <= adaptInfo->getMinTimestep()))
	adaptInfo->setTimestep(adaptInfo->getTimestep() * timeDelta1);    
247 248
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
249 250 251

  void AdaptInstationary::simpleAdaptiveTimeStrategy()
  {
252
    FUNCNAME("AdaptInstationary::simpleAdaptiveTimeStrategy()");
Thomas Witkowski's avatar
Thomas Witkowski committed
253 254 255 256 257 258 259 260 261 262 263 264 265 266

    // estimate before first adaption
    if (adaptInfo->getTime() <= adaptInfo->getStartTime())
      problemIteration->oneIteration(adaptInfo, ESTIMATE);

    adaptInfo->setTime(adaptInfo->getTime() + adaptInfo->getTimestep());
    problemTime->setTime(adaptInfo);
    
    INFO(info,6)("time = %e, timestep = %e\n",
		 adaptInfo->getTime(), adaptInfo->getTimestep());
    
    problemIteration->oneIteration(adaptInfo, FULL_ITERATION);

    adaptInfo->setLastProcessedTimestep(adaptInfo->getTimestep());        
267 268 269 270 271 272

    if (dbgMode) {
	std::cout << "=== ADAPT INFO DEBUG MODE           ===\n";
	std::cout << "=== in simpleAdaptiveTimeStrategy() ===\n";
	adaptInfo->printTimeErrorLowInfo();
    }
Thomas Witkowski's avatar
Thomas Witkowski committed
273 274 275 276 277 278 279 280 281 282 283 284 285
    
    // First, check for increase of timestep
    if (!fixedTimestep && adaptInfo->timeErrorLow())
      adaptInfo->setTimestep(adaptInfo->getTimestep() * timeDelta2);
    
    // Second, check for decrease of timestep
    if (!fixedTimestep &&
	!adaptInfo->timeToleranceReached() &&
	!(adaptInfo->getTimestep() <= adaptInfo->getMinTimestep()))
      adaptInfo->setTimestep(adaptInfo->getTimestep() * timeDelta1);    
  }


286 287
  void AdaptInstationary::oneTimestep()
  {
288
    FUNCNAME("AdaptInstationary::oneTimestep()");
289

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
290
    adaptInfo->setTimestepIteration(0);
291

292 293 294 295 296 297 298
    switch (strategy) {
    case 0:
      explicitTimeStrategy();
      break;
    case 1:
      implicitTimeStrategy();
      break;
Thomas Witkowski's avatar
Thomas Witkowski committed
299 300 301
    case 2:
      simpleAdaptiveTimeStrategy();
      break;
302
    default:
Thomas Witkowski's avatar
Thomas Witkowski committed
303
      ERROR_EXIT("Unknown strategy = %d!\n", strategy);
304
    }
305

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
306
    adaptInfo->incTimestepNumber();
307 308
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
309

310 311 312 313 314
  int AdaptInstationary::adapt()
  {
    FUNCNAME("AdaptInstationary::adapt()");
    int errorCode = 0;

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
315
    TEST_EXIT(adaptInfo->getTimestep() >= adaptInfo->getMinTimestep())
316
      ("timestep < min timestep\n");
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
317
    TEST_EXIT(adaptInfo->getTimestep() <= adaptInfo->getMaxTimestep())
318 319
      ("timestep > max timestep\n");

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
320
    TEST_EXIT(adaptInfo->getTimestep() > 0)("timestep <= 0!\n");
321

322 323 324 325
#if HAVE_PARALLEL_DOMAIN_AMDIS
    MeshDistributor::globalMeshDistributor->initParallelization(); 
#endif

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
326 327
    if (adaptInfo->getTimestepNumber() == 0) {
      adaptInfo->setTime(adaptInfo->getStartTime());
Thomas Witkowski's avatar
Thomas Witkowski committed
328 329
      initialAdaptInfo->setStartTime(adaptInfo->getStartTime());
      initialAdaptInfo->setTime(adaptInfo->getStartTime());
330

Thomas Witkowski's avatar
Thomas Witkowski committed
331
      problemTime->setTime(adaptInfo);
332 333

      // initial adaption
Thomas Witkowski's avatar
Thomas Witkowski committed
334 335
      problemTime->solveInitialProblem(initialAdaptInfo);
      problemTime->transferInitialSolution(adaptInfo);
336 337
    }

Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
338
    while (!adaptInfo->reachedEndTime()) {
Thomas Witkowski's avatar
Thomas Witkowski committed
339
      iterationTimestamp = time(NULL);
340

Thomas Witkowski's avatar
Thomas Witkowski committed
341
      problemTime->initTimestep(adaptInfo);
342
      oneTimestep();
Thomas Witkowski's avatar
Thomas Witkowski committed
343
      problemTime->closeTimestep(adaptInfo);
Thomas Witkowski's avatar
Thomas Witkowski committed
344

345
      if (breakWhenStable && (adaptInfo->getSolverIterations() == 0)) 
346
	break;
347

348 349 350 351 352 353 354 355 356 357 358 359
      // Check if there is a runtime limitation. If there is a runtime limitation
      // and there is no more time for a next adaption loop, than return the error
      // code for rescheduling the problem and break the adaption loop.
      if (checkQueueRuntime()) {
	errorCode = RescheduleErrorCode;
	break;
      }
    }

    return errorCode;
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
360

Thomas Witkowski's avatar
Thomas Witkowski committed
361
  void AdaptInstationary::initialize(std::string aName)
362 363 364 365
  {
    FUNCNAME("AdaptInstationary::initialize()");

    strategy = 0;
Thomas Witkowski's avatar
Thomas Witkowski committed
366 367 368 369
    timeDelta1 = 0.7071;
    timeDelta2 = 1.4142;
    queueRuntime = -1;
    queueSerializationFilename = "__serialized_problem.ser";
370

371 372 373 374 375
    Parameters::get(aName + "->strategy", strategy);
    Parameters::get(aName + "->time delta 1", timeDelta1);
    Parameters::get(aName + "->time delta 2", timeDelta2);
    Parameters::get(aName + "->info", info);
    Parameters::get(aName + "->break when stable", breakWhenStable);
376
    Parameters::get(aName + "->time adaptivity debug mode", dbgMode);
377 378 379
    Parameters::get(aName + "->queue->runtime", queueRuntime);
    Parameters::get(aName + "->queue->serialization filename", 
		    queueSerializationFilename);
380 381
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
382

Thomas Witkowski's avatar
Thomas Witkowski committed
383
  void AdaptInstationary::serialize(std::ostream &out)
384 385 386
  {
    FUNCNAME("AdaptInstationary::serialize()");

Thomas Witkowski's avatar
Thomas Witkowski committed
387
    problemIteration->serialize(out);
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
388
    adaptInfo->serialize(out);
Thomas Witkowski's avatar
Thomas Witkowski committed
389 390
    if (problemTime)
      problemTime->serialize(out);    
391 392
  }

Thomas Witkowski's avatar
Thomas Witkowski committed
393

Thomas Witkowski's avatar
Thomas Witkowski committed
394
  void AdaptInstationary::deserialize(std::istream &in)
395 396
  {
    FUNCNAME("AdaptInstationary::deserialize()");
397 398 399 400

    if (in.fail())
      ERROR_EXIT("File not found for deserialization!\n");

Thomas Witkowski's avatar
Thomas Witkowski committed
401
    problemIteration->deserialize(in);
Thomas Witkowski's avatar
* Bla  
Thomas Witkowski committed
402
    adaptInfo->deserialize(in);
Thomas Witkowski's avatar
Thomas Witkowski committed
403 404
    if (problemTime)
      problemTime->deserialize(in);    
405 406 407 408 409 410
  }


  bool AdaptInstationary::checkQueueRuntime()
  {
    // If there is no time limited runtime queue, there is also nothing to check.
Thomas Witkowski's avatar
Thomas Witkowski committed
411
    if (queueRuntime == -1) {
412 413 414 415 416 417
      return false;
    }
    // Get the current time.
    time_t currentTimestamp = time(NULL);

    // Update list with the last iteration runtimes.
Thomas Witkowski's avatar
Thomas Witkowski committed
418
    lastIterationsDuration.push(currentTimestamp - iterationTimestamp);
419
    // The list should not contain more than 5 elements. If so, delete the oldest one.
Thomas Witkowski's avatar
Thomas Witkowski committed
420 421
    if (lastIterationsDuration.size() > 5)
      lastIterationsDuration.pop();    
422 423

    // Calculate the avarage of the last iterations.
Thomas Witkowski's avatar
Thomas Witkowski committed
424
    std::queue<int> tmpQueue = lastIterationsDuration;
425 426
    int avrgLastIterations = 0;
    while (!tmpQueue.empty()) {
427 428
	avrgLastIterations += tmpQueue.front();
	tmpQueue.pop();
429
    } 
Thomas Witkowski's avatar
Thomas Witkowski committed
430
    avrgLastIterations /= lastIterationsDuration.size();
431

432
    // Check if there is enough time for a further iteration.
Thomas Witkowski's avatar
Thomas Witkowski committed
433 434
    if (initialTimestamp + queueRuntime - currentTimestamp < avrgLastIterations * 2) {
      std::ofstream out(queueSerializationFilename.c_str());
435 436 437 438 439 440 441 442 443 444
      serialize(out);
      out.close();

      return true;
    }

    return false;
  }

}