FreeFOAM The Cross-Platform CFD Toolkit
Time.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | Copyright (C) 1991-2011 OpenCFD Ltd.
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8 License
9  This file is part of OpenFOAM.
10 
11  OpenFOAM is free software: you can redistribute it and/or modify it
12  under the terms of the GNU General Public License as published by
13  the Free Software Foundation, either version 3 of the License, or
14  (at your option) any later version.
15 
16  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
17  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
18  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
19  for more details.
20 
21  You should have received a copy of the GNU General Public License
22  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
23 
24 \*---------------------------------------------------------------------------*/
25 
26 #include "Time.H"
28 
29 #include <sstream>
30 
31 // * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * * * //
32 
34 
35 template<>
37 {
38  "endTime",
39  "noWriteNow",
40  "writeNow",
41  "nextWrite"
42 };
43 
46 
47 template<>
49 {
50  "timeStep",
51  "runTime",
52  "adjustableRunTime",
53  "clockTime",
54  "cpuTime"
55 };
56 
59 
62 
64 
65 
66 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
67 
69 {
71  {
72  scalar timeToNextWrite = max
73  (
74  0.0,
76  );
77 
78  scalar nSteps = timeToNextWrite/deltaT_ - SMALL;
79 
80  // For tiny deltaT the label can overflow!
81  if (nSteps < labelMax)
82  {
83  label nStepsToNextWrite = label(nSteps) + 1;
84 
85  scalar newDeltaT = timeToNextWrite/nStepsToNextWrite;
86 
87  // Control the increase of the time step to within a factor of 2
88  // and the decrease within a factor of 5.
89  if (newDeltaT >= deltaT_)
90  {
91  deltaT_ = min(newDeltaT, 2.0*deltaT_);
92  }
93  else
94  {
95  deltaT_ = max(newDeltaT, 0.2*deltaT_);
96  }
97  }
98  }
99 }
100 
101 
103 {
104  // default is to resume calculation from "latestTime"
105  word startFrom = controlDict_.lookupOrDefault<word>
106  (
107  "startFrom",
108  "latestTime"
109  );
110 
111  if (startFrom == "startTime")
112  {
113  controlDict_.lookup("startTime") >> startTime_;
114  }
115  else
116  {
117  // Search directory for valid time directories
118  instantList timeDirs = findTimes(path());
119 
120  if (startFrom == "firstTime")
121  {
122  if (timeDirs.size())
123  {
124  startTime_ = timeDirs[0].value();
125  }
126  }
127  else if (startFrom == "latestTime")
128  {
129  if (timeDirs.size())
130  {
131  startTime_ = timeDirs[timeDirs.size()-1].value();
132  }
133  }
134  else
135  {
136  FatalIOErrorIn("Time::setControls()", controlDict_)
137  << "expected startTime, firstTime or latestTime"
138  << " found '" << startFrom << "'"
139  << exit(FatalIOError);
140  }
141  }
142 
143  setTime(startTime_, 0);
144 
145  readDict();
146  deltaTSave_ = deltaT_;
147  deltaT0_ = deltaTSave_;
148 
149  if (Pstream::parRun())
150  {
151  scalar sumStartTime = startTime_;
152  reduce(sumStartTime, sumOp<scalar>());
153  if
154  (
155  mag(Pstream::nProcs()*startTime_ - sumStartTime)
156  > Pstream::nProcs()*deltaT_/10.0
157  )
158  {
159  FatalIOErrorIn("Time::setControls()", controlDict_)
160  << "Start time is not the same for all processors" << nl
161  << "processor " << Pstream::myProcNo() << " has startTime "
162  << startTime_ << exit(FatalIOError);
163  }
164  }
165 
166  IOdictionary timeDict
167  (
168  IOobject
169  (
170  "time",
171  timeName(),
172  "uniform",
173  *this,
176  false
177  )
178  );
179 
180  if (timeDict.readIfPresent("deltaT", deltaTSave_))
181  {
182  deltaT0_ = deltaTSave_;
183  }
184 
185  if (timeDict.readIfPresent("index", startTimeIndex_))
186  {
187  timeIndex_ = startTimeIndex_;
188  }
189 }
190 
191 
192 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
193 
195 (
196  const word& controlDictName,
197  const fileName& rootPath,
198  const fileName& caseName,
199  const word& systemName,
200  const word& constantName
201 )
202 :
203  TimePaths
204  (
205  rootPath,
206  caseName,
207  systemName,
208  constantName
209  ),
210 
211  objectRegistry(*this),
212 
213  controlDict_
214  (
215  IOobject
216  (
217  controlDictName,
218  system(),
219  *this,
222  false
223  )
224  ),
225 
226  startTimeIndex_(0),
227  startTime_(0),
228  endTime_(0),
229 
230  stopAt_(saEndTime),
231  writeControl_(wcTimeStep),
232  writeInterval_(GREAT),
233  purgeWrite_(0),
234  subCycling_(false),
235 
236  writeFormat_(IOstream::ASCII),
237  writeVersion_(IOstream::currentVersion),
238  writeCompression_(IOstream::UNCOMPRESSED),
239  graphFormat_("raw"),
240  runTimeModifiable_(true),
241 
242  readLibs_(controlDict_, "libs"),
243  functionObjects_(*this)
244 {
245  setControls();
246 }
247 
248 
250 (
251  const dictionary& dict,
252  const fileName& rootPath,
253  const fileName& caseName,
254  const word& systemName,
255  const word& constantName
256 )
257 :
258  TimePaths
259  (
260  rootPath,
261  caseName,
262  systemName,
263  constantName
264  ),
265 
266  objectRegistry(*this),
267 
268  controlDict_
269  (
270  IOobject
271  (
272  controlDictName,
273  system(),
274  *this,
277  false
278  ),
279  dict
280  ),
281 
282  startTimeIndex_(0),
283  startTime_(0),
284  endTime_(0),
285 
286  stopAt_(saEndTime),
287  writeControl_(wcTimeStep),
288  writeInterval_(GREAT),
289  purgeWrite_(0),
290  subCycling_(false),
291 
292  writeFormat_(IOstream::ASCII),
293  writeVersion_(IOstream::currentVersion),
294  writeCompression_(IOstream::UNCOMPRESSED),
295  graphFormat_("raw"),
296  runTimeModifiable_(true),
297 
298  readLibs_(controlDict_, "libs"),
299  functionObjects_(*this)
300 {
301  setControls();
302 }
303 
304 
306 (
307  const fileName& rootPath,
308  const fileName& caseName,
309  const word& systemName,
310  const word& constantName
311 )
312 :
313  TimePaths
314  (
315  rootPath,
316  caseName,
317  systemName,
318  constantName
319  ),
320 
321  objectRegistry(*this),
322 
323  controlDict_
324  (
325  IOobject
326  (
327  controlDictName,
328  system(),
329  *this,
332  false
333  )
334  ),
335 
336  startTimeIndex_(0),
337  startTime_(0),
338  endTime_(0),
339 
340  stopAt_(saEndTime),
341  writeControl_(wcTimeStep),
342  writeInterval_(GREAT),
343  purgeWrite_(0),
344  subCycling_(false),
345 
346  writeFormat_(IOstream::ASCII),
347  writeVersion_(IOstream::currentVersion),
348  writeCompression_(IOstream::UNCOMPRESSED),
349  graphFormat_("raw"),
350  runTimeModifiable_(true),
351 
352  readLibs_(controlDict_, "libs"),
353  functionObjects_(*this)
354 {}
355 
356 
357 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
358 
360 {
361  // destroy function objects first
362  functionObjects_.clear();
363 }
364 
365 
366 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
367 
369 {
370  std::ostringstream buf;
371  buf.setf(ios_base::fmtflags(format_), ios_base::floatfield);
372  buf.precision(precision_);
373  buf << t;
374  return buf.str();
375 }
376 
377 
379 {
380  return dimensionedScalar::name();
381 }
382 
383 
384 // Search the construction path for times
386 {
387  return findTimes(path());
388 }
389 
390 
392 {
393  instantList timeDirs = findTimes(path());
394 
395  forAllReverse(timeDirs, timeI)
396  {
397  if (timeDirs[timeI] == t)
398  {
399  return timeDirs[timeI].name();
400  }
401  }
402 
403  return word::null;
404 }
405 
406 
408 {
409  instantList timeDirs = findTimes(path());
410 
411  // there is only one time (likely "constant") so return it
412  if (timeDirs.size() == 1)
413  {
414  return timeDirs[0];
415  }
416 
417  if (t < timeDirs[1].value())
418  {
419  return timeDirs[1];
420  }
421  else if (t > timeDirs[timeDirs.size()-1].value())
422  {
423  return timeDirs[timeDirs.size()-1];
424  }
425 
426  label nearestIndex = -1;
427  scalar deltaT = GREAT;
428 
429  for (label timeI=1; timeI < timeDirs.size(); ++timeI)
430  {
431  scalar diff = mag(timeDirs[timeI].value() - t);
432  if (diff < deltaT)
433  {
434  deltaT = diff;
435  nearestIndex = timeI;
436  }
437  }
438 
439  return timeDirs[nearestIndex];
440 }
441 
442 
443 // This should work too,
444 // if we don't worry about checking "constant" explicitly
445 //
446 // Foam::instant Foam::Time::findClosestTime(const scalar t) const
447 // {
448 // instantList timeDirs = findTimes(path());
449 // label timeIndex = min(findClosestTimeIndex(timeDirs, t), 0);
450 // return timeDirs[timeIndex];
451 // }
452 
454 (
455  const instantList& timeDirs,
456  const scalar t
457 )
458 {
459  label nearestIndex = -1;
460  scalar deltaT = GREAT;
461 
462  forAll(timeDirs, timeI)
463  {
464  if (timeDirs[timeI].name() == "constant") continue;
465 
466  scalar diff = mag(timeDirs[timeI].value() - t);
467  if (diff < deltaT)
468  {
469  deltaT = diff;
470  nearestIndex = timeI;
471  }
472  }
473 
474  return nearestIndex;
475 }
476 
477 
478 Foam::label Foam::Time::startTimeIndex() const
479 {
480  return startTimeIndex_;
481 }
482 
483 
485 {
486  return dimensionedScalar("startTime", dimTime, startTime_);
487 }
488 
489 
491 {
492  return dimensionedScalar("endTime", dimTime, endTime_);
493 }
494 
495 
496 bool Foam::Time::run() const
497 {
498  bool running = value() < (endTime_ - 0.5*deltaT_);
499 
500  if (!subCycling_)
501  {
502  // only execute when the condition is no longer true
503  // ie, when exiting the control loop
504  if (!running && timeIndex_ != startTimeIndex_)
505  {
506  // Note, end() also calls an indirect start() as required
507  functionObjects_.end();
508  }
509  }
510 
511  return running;
512 }
513 
514 
516 {
517  bool running = run();
518 
519  if (running)
520  {
521  operator++();
522  }
523 
524  return running;
525 }
526 
527 
528 bool Foam::Time::end() const
529 {
530  return value() > (endTime_ + 0.5*deltaT_);
531 }
532 
533 
534 void Foam::Time::setTime(const Time& t)
535 {
536  value() = t.value();
537  dimensionedScalar::name() = t.dimensionedScalar::name();
538  timeIndex_ = t.timeIndex_;
539 }
540 
541 
542 void Foam::Time::setTime(const instant& inst, const label newIndex)
543 {
544  value() = inst.value();
545  dimensionedScalar::name() = inst.name();
546  timeIndex_ = newIndex;
547 
548  IOdictionary timeDict
549  (
550  IOobject
551  (
552  "time",
553  timeName(),
554  "uniform",
555  *this,
558  false
559  )
560  );
561 
562  timeDict.readIfPresent("deltaT", deltaT_);
563  timeDict.readIfPresent("deltaT0", deltaT0_);
564  timeDict.readIfPresent("index", timeIndex_);
565 }
566 
567 
568 void Foam::Time::setTime(const dimensionedScalar& newTime, const label newIndex)
569 {
570  setTime(newTime.value(), newIndex);
571 }
572 
573 
574 void Foam::Time::setTime(const scalar newTime, const label newIndex)
575 {
576  value() = newTime;
577  dimensionedScalar::name() = timeName(timeToUserTime(newTime));
578  timeIndex_ = newIndex;
579 }
580 
581 
583 {
584  setEndTime(endTime.value());
585 }
586 
587 
588 void Foam::Time::setEndTime(const scalar endTime)
589 {
590  endTime_ = endTime;
591 }
592 
593 
595 {
596  setDeltaT(deltaT.value());
597 }
598 
599 
600 void Foam::Time::setDeltaT(const scalar deltaT)
601 {
602  deltaT_ = deltaT;
603  deltaTchanged_ = true;
604  adjustDeltaT();
605 }
606 
607 
608 Foam::TimeState Foam::Time::subCycle(const label nSubCycles)
609 {
610  subCycling_ = true;
611  prevTimeState_.set(new TimeState(*this));
612 
613  setTime(*this - deltaT(), (timeIndex() - 1)*nSubCycles);
614  deltaT_ /= nSubCycles;
615  deltaT0_ /= nSubCycles;
616  deltaTSave_ = deltaT0_;
617 
618  return prevTimeState();
619 }
620 
621 
623 {
624  if (subCycling_)
625  {
626  subCycling_ = false;
627  TimeState::operator=(prevTimeState());
628  prevTimeState_.clear();
629  }
630 }
631 
632 
633 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
634 
636 {
637  return operator+=(deltaT.value());
638 }
639 
640 
641 Foam::Time& Foam::Time::operator+=(const scalar deltaT)
642 {
643  setDeltaT(deltaT);
644  return operator++();
645 }
646 
647 
649 {
650  readModifiedObjects();
651 
652  if (!subCycling_)
653  {
654  if (timeIndex_ == startTimeIndex_)
655  {
656  functionObjects_.start();
657  }
658  else
659  {
660  functionObjects_.execute();
661  }
662  }
663 
664  deltaT0_ = deltaTSave_;
665  deltaTSave_ = deltaT_;
666 
667  const word oldTimeName = dimensionedScalar::name();
668 
669  setTime(value() + deltaT_, timeIndex_ + 1);
670 
671  // If the time is very close to zero reset to zero
672  if (mag(value()) < 10*SMALL*deltaT_)
673  {
674  setTime(0.0, timeIndex_);
675  }
676 
677  // Check that new time representation differs from old one
678  if (dimensionedScalar::name() == oldTimeName)
679  {
680  int oldPrecision = precision_;
681  do
682  {
683  precision_++;
684  setTime(value(), timeIndex());
685  }
686  while (precision_ < 100 && dimensionedScalar::name() == oldTimeName);
687 
688  WarningIn("Time::operator++()")
689  << "Increased the timePrecision from " << oldPrecision
690  << " to " << precision_
691  << " to distinguish between timeNames at time " << value()
692  << endl;
693  }
694 
695  switch (writeControl_)
696  {
697  case wcTimeStep:
698  outputTime_ = !(timeIndex_ % label(writeInterval_));
699  break;
700 
701  case wcRunTime:
702  case wcAdjustableRunTime:
703  {
704  label outputIndex =
705  label(((value() - startTime_) + 0.5*deltaT_)/writeInterval_);
706 
707  if (outputIndex > outputTimeIndex_)
708  {
709  outputTime_ = true;
710  outputTimeIndex_ = outputIndex;
711  }
712  else
713  {
714  outputTime_ = false;
715  }
716  }
717  break;
718 
719  case wcCpuTime:
720  {
721  label outputIndex = label
722  (
723  returnReduce(elapsedCpuTime(), maxOp<double>())
724  / writeInterval_
725  );
726  if (outputIndex > outputTimeIndex_)
727  {
728  outputTime_ = true;
729  outputTimeIndex_ = outputIndex;
730  }
731  else
732  {
733  outputTime_ = false;
734  }
735  }
736  break;
737 
738  case wcClockTime:
739  {
740  label outputIndex = label
741  (
742  returnReduce(label(elapsedClockTime()), maxOp<label>())
743  / writeInterval_
744  );
745  if (outputIndex > outputTimeIndex_)
746  {
747  outputTime_ = true;
748  outputTimeIndex_ = outputIndex;
749  }
750  else
751  {
752  outputTime_ = false;
753  }
754  }
755  break;
756  }
757 
758  // see if endTime needs adjustment to stop at the next run()/end() check
759  if (!end())
760  {
761  if (stopAt_ == saNoWriteNow)
762  {
763  endTime_ = value();
764  }
765  else if (stopAt_ == saWriteNow)
766  {
767  endTime_ = value();
768  outputTime_ = true;
769  }
770  else if (stopAt_ == saNextWrite && outputTime_ == true)
771  {
772  endTime_ = value();
773  }
774  }
775 
776  return *this;
777 }
778 
779 
781 {
782  return operator++();
783 }
784 
785 
786 // ************************ vim: set sw=4 sts=4 et: ************************ //