Time.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration | Website: https://openfoam.org
5  \\ / A nd | Copyright (C) 2011-2023 OpenFOAM Foundation
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"
27 #include "timeIOdictionary.H"
28 #include "argList.H"
29 
30 // * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * * * //
31 
32 namespace Foam
33 {
35 
36  template<>
37  const char* Foam::NamedEnum
38  <
40  4
41  >::names[] =
42  {
43  "endTime",
44  "noWriteNow",
45  "writeNow",
46  "nextWrite"
47  };
48 
49  template<>
50  const char* Foam::NamedEnum
51  <
53  5
54  >::names[] =
55  {
56  "timeStep",
57  "runTime",
58  "adjustableRunTime",
59  "clockTime",
60  "cpuTime"
61  };
62 }
63 
66 
69 
71 
73 
75 
76 const int Foam::Time::maxPrecision_(3 - log10(small));
77 
79 
80 
81 // * * * * * * * * * * * * Protected Member Functions * * * * * * * * * * * //
82 
84 {
85  const scalar timeToNextAction = min
86  (
87  max
88  (
89  0,
90  (writeTimeIndex_ + 1)*writeInterval_ - (value() - beginTime_)
91  ),
92  functionObjects_.timeToNextAction()
93  );
94 
95  const scalar nSteps = timeToNextAction/deltaT_;
96 
97  // Ensure nStepsToNextWrite does not overflow
98  if (nSteps < labelMax)
99  {
100  // Allow the time-step to increase by up to 1%
101  // to accommodate the next write time before splitting
102  const label nStepsToNextWrite = label(max(nSteps, 1) + 0.99);
103 
104  deltaT_ = timeToNextAction/nStepsToNextWrite;
105  }
106 }
107 
108 
110 {
111  // default is to resume calculation from "latestTime"
112  const word startFrom = controlDict_.lookupOrDefault<word>
113  (
114  "startFrom",
115  "latestTime"
116  );
117 
118  if (startFrom == "startTime")
119  {
120  controlDict_.lookup("startTime") >> startTime_;
121  startTime_ = userTimeToTime(startTime_);
122  }
123  else
124  {
125  // Search directory for valid time directories
126  instantList timeDirs = findTimes(path(), constant());
127 
128  if (startFrom == "firstTime")
129  {
130  if (timeDirs.size())
131  {
132  if (timeDirs[0].name() == constant() && timeDirs.size() >= 2)
133  {
134  startTime_ = userTimeToTime(timeDirs[1].value());
135  }
136  else
137  {
138  startTime_ = userTimeToTime(timeDirs[0].value());
139  }
140  }
141  }
142  else if (startFrom == "latestTime")
143  {
144  if (timeDirs.size())
145  {
146  startTime_ = userTimeToTime(timeDirs.last().value());
147  }
148  }
149  else
150  {
151  FatalIOErrorInFunction(controlDict_)
152  << "expected startTime, firstTime or latestTime"
153  << " found '" << startFrom << "'"
154  << exit(FatalIOError);
155  }
156  }
157 
158  setTime(startTime_, 0);
159  readDict();
160  deltaTSave_ = deltaT_;
161  deltaT0_ = deltaT_;
162 
163  // Check if time directory exists
164  // If not increase time precision to see if it is formatted differently.
165  if (!fileHandler().exists(timePath(), false, false))
166  {
167  int oldPrecision = curPrecision_;
168  int requiredPrecision = -1;
169 
170  for
171  (
172  curPrecision_ = maxPrecision_;
173  curPrecision_ > oldPrecision;
174  curPrecision_--
175  )
176  {
177  // Update the time formatting
178  setTime(startTime_, 0);
179 
180  // Check the existence of the time directory with the new format
181  if (fileHandler().exists(timePath(), false, false))
182  {
183  requiredPrecision = curPrecision_;
184  }
185  }
186 
187  if (requiredPrecision > 0)
188  {
189  // Update the time precision
190  curPrecision_ = requiredPrecision;
191 
192  // Update the time formatting
193  setTime(startTime_, 0);
194 
196  << "Increasing the timePrecision from " << oldPrecision
197  << " to " << curPrecision_
198  << " to support the formatting of the current time directory "
199  << name() << nl << endl;
200  }
201  else
202  {
203  // Could not find time directory so assume it is not present
204  curPrecision_ = oldPrecision;
205 
206  // Revert the time formatting
207  setTime(startTime_, 0);
208  }
209  }
210 
211  if (Pstream::parRun())
212  {
213  scalar sumStartTime = startTime_;
214  reduce(sumStartTime, sumOp<scalar>());
215  if
216  (
217  mag(Pstream::nProcs()*startTime_ - sumStartTime)
218  > Pstream::nProcs()*deltaT_/10.0
219  )
220  {
221  FatalIOErrorInFunction(controlDict_)
222  << "Start time is not the same for all processors" << nl
223  << "processor " << Pstream::myProcNo() << " has startTime "
224  << startTime_ << exit(FatalIOError);
225  }
226  }
227 
228  timeIOdictionary timeDict
229  (
230  IOobject
231  (
232  "time",
233  name(),
234  "uniform",
235  *this,
238  false
239  )
240  );
241 
242  if (controlDict_.found("beginTime"))
243  {
244  beginTime_ = userTimeToTime(controlDict_.lookup<scalar>("beginTime"));
245  }
246  else if (timeDict.found("beginTime"))
247  {
248  beginTime_ = userTimeToTime(timeDict.lookup<scalar>("beginTime"));
249  }
250  else
251  {
252  beginTime_ = startTime_;
253  }
254 
255  // Read and set the deltaT only if time-step adjustment is active
256  // otherwise use the deltaT from the controlDict
257  if (controlDict_.lookupOrDefault<Switch>("adjustTimeStep", false))
258  {
259  if (timeDict.found("deltaT"))
260  {
261  deltaT_ = userTimeToTime(timeDict.lookup<scalar>("deltaT"));
262  deltaTSave_ = deltaT_;
263  deltaT0_ = deltaT_;
264  }
265  }
266 
267  if (timeDict.found("deltaT0"))
268  {
269  deltaT0_ = userTimeToTime(timeDict.lookup<scalar>("deltaT0"));
270  }
271 
272  if (timeDict.readIfPresent("index", startTimeIndex_))
273  {
274  timeIndex_ = startTimeIndex_;
275  }
276 
277  // Set writeTimeIndex_ to correspond to beginTime_
278  if
279  (
280  (
281  writeControl_ == writeControl::runTime
282  || writeControl_ == writeControl::adjustableRunTime
283  )
284  )
285  {
286  writeTimeIndex_ = label
287  (
288  ((value() - beginTime_) + 0.5*deltaT_)/writeInterval_
289  );
290  }
291 
292 
293  // Check if values stored in time dictionary are consistent
294 
295  // 1. Based on time name
296  bool checkValue = true;
297 
298  string storedTimeName;
299  if (timeDict.readIfPresent("name", storedTimeName))
300  {
301  if (storedTimeName == name())
302  {
303  // Same time. No need to check stored value
304  checkValue = false;
305  }
306  }
307 
308  // 2. Based on time value
309  // (consistent up to the current time writing precision so it won't
310  // trigger if we just change the write precision)
311  if (checkValue)
312  {
313  scalar storedTimeValue;
314  if (timeDict.readIfPresent("value", storedTimeValue))
315  {
316  word storedTimeName(timeName(storedTimeValue));
317 
318  if (storedTimeName != name())
319  {
320  IOWarningInFunction(timeDict)
321  << "Time read from time dictionary " << storedTimeName
322  << " differs from actual time " << name() << '.' << nl
323  << " This may cause unexpected database behaviour."
324  << " If you are not interested" << nl
325  << " in preserving time state delete"
326  << " the time dictionary."
327  << endl;
328  }
329  }
330  }
331 }
332 
333 
334 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
335 
337 (
338  const word& controlDictName,
339  const argList& args,
340  const bool enableFunctionObjects
341 )
342 :
343  TimePaths
344  (
345  args.parRunControl().parRun(),
346  args.rootPath(),
347  args.globalCaseName(),
348  args.caseName()
349  ),
350 
351  objectRegistry(*this),
352 
353  runTimeModifiable_(false),
354 
355  controlDict_
356  (
357  IOobject
358  (
359  controlDictName,
360  system(),
361  *this,
362  IOobject::MUST_READ_IF_MODIFIED,
363  IOobject::NO_WRITE
364  )
365  ),
366 
367  startTimeIndex_(0),
368  startTime_(0),
369  endTime_(0),
370  beginTime_(startTime_),
371 
372  userTime_(userTimes::userTime::New(controlDict_)),
373 
374  stopAt_(stopAtControl::endTime),
375  writeControl_(writeControl::timeStep),
376  writeInterval_(great),
377  purgeWrite_(0),
378  writeOnce_(false),
379 
380  subCycling_(false),
381 
382  sigWriteNow_(writeInfoHeader, *this),
383  sigStopAtWriteNow_(writeInfoHeader, *this),
384 
385  writeFormat_(IOstream::ASCII),
386  writeVersion_(IOstream::currentVersion),
387  writeCompression_(IOstream::UNCOMPRESSED),
388  cacheTemporaryObjects_(true),
389 
390  functionObjects_
391  (
392  *this,
393  enableFunctionObjects
394  ? argList::validOptions.found("withFunctionObjects")
395  ? args.optionFound("withFunctionObjects")
396  : !args.optionFound("noFunctionObjects")
397  : false
398  )
399 {
400  libs.open(controlDict_, "libs");
401 
402  // Explicitly set read flags on objectRegistry so anything constructed
403  // from it reads as well (e.g. fvSolution).
405 
406  if (args.options().found("case"))
407  {
408  const wordList switchSets
409  (
410  {
411  "InfoSwitches",
412  "OptimisationSwitches",
413  "DebugSwitches",
414  "DimensionedConstants",
415  "DimensionSets"
416  }
417  );
418 
419  forAll(switchSets, i)
420  {
421  if (controlDict_.found(switchSets[i]))
422  {
423  IOWarningInFunction(controlDict_)
424  << switchSets[i]
425  << " in system/controlDict are only processed if "
426  << args.executable() << " is run in the "
427  << args.path() << " directory" << endl;
428  }
429  }
430  }
431 
432  setControls();
433 
434  // Add a watch on the controlDict file after runTimeModifiable_ is set
435  controlDict_.addWatch();
436 }
437 
438 
440 (
441  const word& controlDictName,
442  const fileName& rootPath,
443  const fileName& caseName,
444  const bool enableFunctionObjects
445 )
446 :
447  TimePaths(rootPath, caseName),
448 
449  objectRegistry(*this),
450 
451  runTimeModifiable_(false),
452 
453  controlDict_
454  (
455  IOobject
456  (
457  controlDictName,
458  system(),
459  *this,
460  IOobject::MUST_READ_IF_MODIFIED,
461  IOobject::NO_WRITE
462  )
463  ),
464 
465  startTimeIndex_(0),
466  startTime_(0),
467  endTime_(0),
468  beginTime_(startTime_),
469 
470  userTime_(userTimes::userTime::New(controlDict_)),
471 
472  stopAt_(stopAtControl::endTime),
473  writeControl_(writeControl::timeStep),
474  writeInterval_(great),
475  purgeWrite_(0),
476  writeOnce_(false),
477 
478  subCycling_(false),
479 
480  sigWriteNow_(writeInfoHeader, *this),
481  sigStopAtWriteNow_(writeInfoHeader, *this),
482 
483  writeFormat_(IOstream::ASCII),
484  writeVersion_(IOstream::currentVersion),
485  writeCompression_(IOstream::UNCOMPRESSED),
486  cacheTemporaryObjects_(true),
487 
488  functionObjects_(*this, enableFunctionObjects)
489 {
490  libs.open(controlDict_, "libs");
491 
492  // Explicitly set read flags on objectRegistry so anything constructed
493  // from it reads as well (e.g. fvSolution).
495 
496  setControls();
497 
498  // Add a watch on the controlDict file after runTimeModifiable_ is set
499  controlDict_.addWatch();
500 }
501 
502 
504 (
505  const dictionary& dict,
506  const fileName& rootPath,
507  const fileName& caseName,
508  const bool enableFunctionObjects
509 )
510 :
511  TimePaths(rootPath, caseName),
512 
513  objectRegistry(*this),
514 
515  runTimeModifiable_(false),
516 
517  controlDict_
518  (
519  IOobject
520  (
521  controlDictName,
522  system(),
523  *this,
524  IOobject::MUST_READ_IF_MODIFIED,
525  IOobject::NO_WRITE
526  ),
527  dict
528  ),
529 
530  startTimeIndex_(0),
531  startTime_(0),
532  endTime_(0),
533  beginTime_(startTime_),
534 
535  userTime_(userTimes::userTime::New(controlDict_)),
536 
537  stopAt_(stopAtControl::endTime),
538  writeControl_(writeControl::timeStep),
539  writeInterval_(great),
540  purgeWrite_(0),
541  writeOnce_(false),
542 
543  subCycling_(false),
544 
545  sigWriteNow_(writeInfoHeader, *this),
546  sigStopAtWriteNow_(writeInfoHeader, *this),
547 
548  writeFormat_(IOstream::ASCII),
549  writeVersion_(IOstream::currentVersion),
550  writeCompression_(IOstream::UNCOMPRESSED),
551  cacheTemporaryObjects_(true),
552 
553  functionObjects_(*this, enableFunctionObjects)
554 {
555  libs.open(controlDict_, "libs");
556 
557  // Explicitly set read flags on objectRegistry so anything constructed
558  // from it reads as well (e.g. fvSolution).
560 
561  setControls();
562 
563  // Add a watch on the controlDict file after runTimeModifiable_ is set
564  controlDict_.addWatch();
565 }
566 
567 
569 (
570  const fileName& rootPath,
571  const fileName& caseName,
572  const bool enableFunctionObjects
573 )
574 :
575  TimePaths(rootPath, caseName),
576 
577  objectRegistry(*this),
578 
579  runTimeModifiable_(false),
580 
581  controlDict_
582  (
583  IOobject
584  (
585  controlDictName,
586  system(),
587  *this,
588  IOobject::NO_READ,
589  IOobject::NO_WRITE
590  )
591  ),
592 
593  startTimeIndex_(0),
594  startTime_(0),
595  endTime_(0),
596  beginTime_(startTime_),
597 
598  userTime_(userTimes::userTime::New(controlDict_)),
599 
600  stopAt_(stopAtControl::endTime),
601  writeControl_(writeControl::timeStep),
602  writeInterval_(great),
603  purgeWrite_(0),
604  writeOnce_(false),
605 
606  subCycling_(false),
607 
608  writeFormat_(IOstream::ASCII),
609  writeVersion_(IOstream::currentVersion),
610  writeCompression_(IOstream::UNCOMPRESSED),
611  cacheTemporaryObjects_(true),
612 
613  functionObjects_(*this, enableFunctionObjects)
614 {
615  libs.open(controlDict_, "libs");
616 }
617 
618 
619 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
620 
622 {
623  // Destroy function objects first
624  functionObjects_.clear();
625 }
626 
627 
628 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
629 
630 Foam::word Foam::Time::timeName(const scalar t, const int precision)
631 {
632  std::ostringstream buf;
633  buf.setf(ios_base::fmtflags(format_), ios_base::floatfield);
634  buf.precision(precision);
635  buf << t;
636  return buf.str();
637 }
638 
639 
641 {
642  return findTimes(path(), constant());
643 }
644 
645 
647 (
648  const fileName& dir,
649  const word& name,
650  const IOobject::readOption rOpt,
651  const word& stopInstance
652 ) const
653 {
654  IOobject startIO
655  (
656  name, // name might be empty!
658  dir,
659  *this,
660  rOpt
661  );
662 
663  IOobject io
664  (
665  fileHandler().findInstance
666  (
667  startIO,
668  userTimeValue(),
669  stopInstance
670  )
671  );
672  return io.instance();
673 }
674 
675 
677 (
678  const fileName& directory,
679  const instant& t
680 ) const
681 {
682  // Simplified version: use findTimes (readDir + sort). The expensive
683  // bit is the readDir, not the sorting. Tbd: avoid calling findInstancePath
684  // from filePath.
685 
686  instantList timeDirs = findTimes(path(), constant());
687  // Note:
688  // - times will include constant (with value 0) as first element.
689  // For backwards compatibility make sure to find 0 in preference
690  // to constant.
691  // - list is sorted so could use binary search
692 
694  {
695  if (t.equal(timeDirs[i].value()))
696  {
697  return timeDirs[i].name();
698  }
699  }
700 
701  return word::null;
702 }
703 
704 
706 {
707  return findInstancePath(path(), t);
708 }
709 
710 
712 {
713  instantList timeDirs = findTimes(path(), constant());
714 
715  // There is only one time (likely "constant") so return it
716  if (timeDirs.size() == 1)
717  {
718  return timeDirs[0];
719  }
720 
721  if (t < timeDirs[1].value())
722  {
723  return timeDirs[1];
724  }
725  else if (t > timeDirs.last().value())
726  {
727  return timeDirs.last();
728  }
729 
730  label nearestIndex = -1;
731  scalar deltaT = great;
732 
733  for (label timei=1; timei < timeDirs.size(); ++timei)
734  {
735  scalar diff = mag(timeDirs[timei].value() - t);
736  if (diff < deltaT)
737  {
738  deltaT = diff;
739  nearestIndex = timei;
740  }
741  }
742 
743  return timeDirs[nearestIndex];
744 }
745 
746 
748 (
749  const instantList& timeDirs,
750  const scalar t,
751  const word& constantName
752 )
753 {
754  label nearestIndex = -1;
755  scalar deltaT = great;
756 
757  forAll(timeDirs, timei)
758  {
759  if (timeDirs[timei].name() == constantName) continue;
760 
761  scalar diff = mag(timeDirs[timei].value() - t);
762  if (diff < deltaT)
763  {
764  deltaT = diff;
765  nearestIndex = timei;
766  }
767  }
768 
769  return nearestIndex;
770 }
771 
772 
774 {
775  return startTimeIndex_;
776 }
777 
778 
780 {
781  return dimensionedScalar
782  (
783  timeName(timeToUserTime(beginTime_)),
784  dimTime,
785  beginTime_
786  );
787 }
788 
789 
791 {
792  return dimensionedScalar
793  (
794  timeName(timeToUserTime(startTime_)),
795  dimTime,
796  startTime_
797  );
798 }
799 
800 
802 {
803  return dimensionedScalar
804  (
805  timeName(timeToUserTime(endTime_)),
806  dimTime,
807  endTime_
808  );
809 }
810 
811 
813 {
814  return *userTime_;
815 }
816 
817 
818 Foam::scalar Foam::Time::userTimeValue() const
819 {
820  return userTime_->timeToUserTime(value());
821 }
822 
823 
824 Foam::scalar Foam::Time::userTimeToTime(const scalar tau) const
825 {
826  return userTime_->userTimeToTime(tau);
827 }
828 
829 
830 Foam::scalar Foam::Time::timeToUserTime(const scalar t) const
831 {
832  return userTime_->timeToUserTime(t);
833 }
834 
835 
837 {
838  if (name() == constant())
839  {
840  return constant();
841  }
842  else
843  {
844  return timeName(userTimeValue()) + userTime_->unit();
845  }
846 }
847 
848 
850 {
851  return value() < (endTime_ - 0.5*deltaT_);
852 }
853 
854 
855 bool Foam::Time::run() const
856 {
857  bool running = this->running();
858 
859  if (!subCycling_)
860  {
861  if (!running && timeIndex_ != startTimeIndex_)
862  {
863  functionObjects_.execute();
864  functionObjects_.end();
865 
866  if (cacheTemporaryObjects_)
867  {
868  cacheTemporaryObjects_ = checkCacheTemporaryObjects();
869  }
870  }
871  }
872 
873  if (running)
874  {
875  if (!subCycling_)
876  {
877  const_cast<Time&>(*this).readModifiedObjects();
878 
879  if (timeIndex_ == startTimeIndex_)
880  {
881  functionObjects_.start();
882  }
883  else
884  {
885  functionObjects_.execute();
886 
887  if (cacheTemporaryObjects_)
888  {
889  cacheTemporaryObjects_ = checkCacheTemporaryObjects();
890  }
891  }
892  }
893 
894  // Re-evaluate if running in case a function object has changed things
895  running = this->running();
896  }
897 
898  return running;
899 }
900 
901 
903 {
904  bool running = run();
905 
906  if (running)
907  {
908  operator++();
909  }
910 
911  return running;
912 }
913 
914 
915 bool Foam::Time::end() const
916 {
917  return value() > (endTime_ + 0.5*deltaT_);
918 }
919 
920 
921 bool Foam::Time::stopAt(const stopAtControl sa) const
922 {
923  const bool changed = (stopAt_ != sa);
924  stopAt_ = sa;
925 
926  // adjust endTime
927  if (sa == stopAtControl::endTime)
928  {
929  controlDict_.lookup("endTime") >> endTime_;
930  }
931  else
932  {
933  endTime_ = great;
934  }
935  return changed;
936 }
937 
938 
939 void Foam::Time::setTime(const Time& t)
940 {
941  value() = t.value();
942  dimensionedScalar::name() = t.dimensionedScalar::name();
943  timeIndex_ = t.timeIndex_;
944  fileHandler().setTime(*this);
945 }
946 
947 
948 void Foam::Time::setTime(const instant& inst, const label newIndex)
949 {
950  value() = userTimeToTime(inst.value());
951  dimensionedScalar::name() = inst.name();
952  timeIndex_ = newIndex;
953 
954  timeIOdictionary timeDict
955  (
956  IOobject
957  (
958  "time",
959  name(),
960  "uniform",
961  *this,
964  false
965  )
966  );
967 
968  if (timeDict.found("deltaT"))
969  {
970  deltaT_ = userTimeToTime(timeDict.lookup<scalar>("deltaT"));
971  }
972 
973  if (timeDict.found("deltaT0"))
974  {
975  deltaT0_ = userTimeToTime(timeDict.lookup<scalar>("deltaT0"));
976  }
977 
978  timeDict.readIfPresent("index", timeIndex_);
979 
980  fileHandler().setTime(*this);
981 }
982 
983 
984 void Foam::Time::setTime(const dimensionedScalar& newTime, const label newIndex)
985 {
986  setTime(newTime.value(), newIndex);
987 }
988 
989 
990 void Foam::Time::setTime(const scalar newTime, const label newIndex)
991 {
992  value() = newTime;
993  dimensionedScalar::name() = timeName(timeToUserTime(newTime));
994  timeIndex_ = newIndex;
995  fileHandler().setTime(*this);
996 }
997 
998 
1000 {
1001  setEndTime(endTime.value());
1002 }
1003 
1004 
1005 void Foam::Time::setEndTime(const scalar endTime)
1006 {
1007  endTime_ = endTime;
1008 }
1009 
1010 
1012 {
1013  setDeltaT(deltaT.value());
1014 }
1015 
1016 
1017 void Foam::Time::setDeltaT(const scalar deltaT)
1018 {
1019  setDeltaTNoAdjust(deltaT);
1020 
1021  if (writeControl_ == writeControl::adjustableRunTime)
1022  {
1023  adjustDeltaT();
1024  }
1025 }
1026 
1027 
1028 void Foam::Time::setDeltaTNoAdjust(const scalar deltaT)
1029 {
1030  deltaT_ = deltaT;
1031  deltaTchanged_ = true;
1032 }
1033 
1034 
1035 void Foam::Time::setWriteInterval(const scalar writeInterval)
1036 {
1037  if (writeInterval_ == great || !equal(writeInterval, writeInterval_))
1038  {
1039  writeInterval_ = writeInterval;
1040 
1041  if
1042  (
1043  writeControl_ == writeControl::runTime
1044  || writeControl_ == writeControl::adjustableRunTime
1045  )
1046  {
1047  // Recalculate writeTimeIndex_ for consistency with the new
1048  // writeInterval
1049  writeTimeIndex_ = label
1050  (
1051  ((value() - beginTime_) + 0.5*deltaT_)/writeInterval_
1052  );
1053  }
1054  else if (writeControl_ == writeControl::timeStep)
1055  {
1056  // Set to the nearest integer
1057  writeInterval_ = label(writeInterval + 0.5);
1058  }
1059  }
1060 }
1061 
1062 
1064 {
1065  subCycling_ = true;
1066  prevTimeState_.set(new TimeState(*this));
1067 
1068  setTime(*this - deltaT(), (timeIndex() - 1)*nSubCycles);
1069  deltaT_ /= nSubCycles;
1070  deltaT0_ /= nSubCycles;
1071  deltaTSave_ = deltaT0_;
1072 
1073  return prevTimeState();
1074 }
1075 
1076 
1078 {
1079  if (subCycling_)
1080  {
1081  subCycling_ = false;
1082  TimeState::operator=(prevTimeState());
1083  prevTimeState_.clear();
1084  }
1085 }
1086 
1087 
1088 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
1089 
1091 {
1092  return operator+=(deltaT.value());
1093 }
1094 
1095 
1096 Foam::Time& Foam::Time::operator+=(const scalar deltaT)
1097 {
1098  setDeltaT(deltaT);
1099  return operator++();
1100 }
1101 
1102 
1104 {
1105  deltaT0_ = deltaTSave_;
1106  deltaTSave_ = deltaT_;
1107 
1108  // Save old time value and name
1109  const scalar oldTimeValue = timeToUserTime(value());
1110  const word oldTimeName = dimensionedScalar::name();
1111 
1112  // Increment time
1113  setTime(value() + deltaT_, timeIndex_ + 1);
1114 
1115  if (!subCycling_)
1116  {
1117  // If the time is very close to zero reset to zero
1118  if (mag(value()) < 10*small*deltaT_)
1119  {
1120  setTime(0, timeIndex_);
1121  }
1122 
1123  if (sigStopAtWriteNow_.active() || sigWriteNow_.active())
1124  {
1125  // A signal might have been sent on one processor only
1126  // Reduce so all decide the same.
1127 
1128  label flag = 0;
1129  if
1130  (
1131  sigStopAtWriteNow_.active()
1132  && stopAt_ == stopAtControl::writeNow
1133  )
1134  {
1135  flag += 1;
1136  }
1137  if (sigWriteNow_.active() && writeOnce_)
1138  {
1139  flag += 2;
1140  }
1141  reduce(flag, maxOp<label>());
1142 
1143  if (flag & 1)
1144  {
1145  stopAt_ = stopAtControl::writeNow;
1146  }
1147  if (flag & 2)
1148  {
1149  writeOnce_ = true;
1150  }
1151  }
1152 
1153  writeTime_ = false;
1154 
1155  switch (writeControl_)
1156  {
1157  case writeControl::timeStep:
1158  writeTime_ = !(timeIndex_ % label(writeInterval_));
1159  break;
1160 
1161  case writeControl::runTime:
1162  case writeControl::adjustableRunTime:
1163  {
1164  label writeIndex = label
1165  (
1166  ((value() - beginTime_) + 0.5*deltaT_)
1167  / writeInterval_
1168  );
1169 
1170  if (writeIndex > writeTimeIndex_)
1171  {
1172  writeTime_ = true;
1173  writeTimeIndex_ = writeIndex;
1174  }
1175  }
1176  break;
1177 
1178  case writeControl::cpuTime:
1179  {
1180  label writeIndex = label
1181  (
1182  returnReduce(elapsedCpuTime(), maxOp<double>())
1183  / writeInterval_
1184  );
1185  if (writeIndex > writeTimeIndex_)
1186  {
1187  writeTime_ = true;
1188  writeTimeIndex_ = writeIndex;
1189  }
1190  }
1191  break;
1192 
1193  case writeControl::clockTime:
1194  {
1195  label writeIndex = label
1196  (
1197  returnReduce(label(elapsedClockTime()), maxOp<label>())
1198  / writeInterval_
1199  );
1200  if (writeIndex > writeTimeIndex_)
1201  {
1202  writeTime_ = true;
1203  writeTimeIndex_ = writeIndex;
1204  }
1205  }
1206  break;
1207  }
1208 
1209 
1210  // Check if endTime needs adjustment to stop at the next run()/end()
1211  if (!end())
1212  {
1213  if (stopAt_ == stopAtControl::noWriteNow)
1214  {
1215  endTime_ = value();
1216  }
1217  else if (stopAt_ == stopAtControl::writeNow)
1218  {
1219  endTime_ = value();
1220  writeTime_ = true;
1221  }
1222  else if (stopAt_ == stopAtControl::nextWrite && writeTime_ == true)
1223  {
1224  endTime_ = value();
1225  }
1226  }
1227 
1228  // Override writeTime if one-shot writing
1229  if (writeOnce_)
1230  {
1231  writeTime_ = true;
1232  writeOnce_ = false;
1233  }
1234 
1235  // Adjust the precision of the time directory name if necessary
1236  if (writeTime_)
1237  {
1238  // User-time equivalent of deltaT
1239  const scalar userDeltaT =
1240  timeToUserTime(value()) - timeToUserTime(value() - deltaT_);
1241 
1242  // Tolerance used when testing time equivalence
1243  const scalar timeTol =
1244  max(min(pow(10.0, -precision_), 0.1*userDeltaT), small);
1245 
1246  // Time value obtained by reading timeName
1247  scalar timeNameValue = -vGreat;
1248 
1249  // Check that new time representation differs from old one
1250  // reinterpretation of the word
1251  if
1252  (
1253  readScalar(dimensionedScalar::name().c_str(), timeNameValue)
1254  && (mag(timeNameValue - oldTimeValue - userDeltaT) > timeTol)
1255  )
1256  {
1257  int oldPrecision = curPrecision_;
1258  while
1259  (
1260  curPrecision_ < maxPrecision_
1261  && readScalar(dimensionedScalar::name().c_str(), timeNameValue)
1262  && (mag(timeNameValue - oldTimeValue - userDeltaT) > timeTol)
1263  )
1264  {
1265  curPrecision_++;
1266  setTime(value(), timeIndex());
1267  }
1268 
1269  if (curPrecision_ != oldPrecision)
1270  {
1272  << "Increased the timePrecision from " << oldPrecision
1273  << " to " << curPrecision_
1274  << " to distinguish between timeNames at time "
1276  << endl;
1277 
1278  if (curPrecision_ == maxPrecision_)
1279  {
1280  // Reached maxPrecision limit
1282  << "Current time name " << dimensionedScalar::name()
1283  << nl
1284  << " The maximum time precision has been reached"
1285  " which might result in overwriting previous"
1286  " results."
1287  << endl;
1288  }
1289 
1290  // Check if round-off error caused time-reversal
1291  scalar oldTimeNameValue = -vGreat;
1292  if
1293  (
1294  readScalar(oldTimeName.c_str(), oldTimeNameValue)
1295  && (
1296  sign(timeNameValue - oldTimeNameValue)
1297  != sign(deltaT_)
1298  )
1299  )
1300  {
1302  << "Current time name " << dimensionedScalar::name()
1303  << " is set to an instance prior to the "
1304  "previous one "
1305  << oldTimeName << nl
1306  << " This might result in temporal "
1307  "discontinuities."
1308  << endl;
1309  }
1310  }
1311  }
1312  }
1313  }
1314 
1315  return *this;
1316 }
1317 
1318 
1320 {
1321  return operator++();
1322 }
1323 
1324 
1325 // ************************************************************************* //
bool found
#define forAll(list, i)
Loop across all elements in list.
Definition: UList.H:434
#define forAllReverse(list, i)
Reverse loop across all elements in list.
Definition: UList.H:446
IOobject defines the attributes of an object for which implicit objectRegistry management is supporte...
Definition: IOobject.H:99
readOption
Enumeration defining the read options.
Definition: IOobject.H:117
@ MUST_READ_IF_MODIFIED
Definition: IOobject.H:119
fileName & instance() const
Return the instance directory, constant, system, <time> etc.
Definition: IOobject.C:355
readOption readOpt() const
Definition: IOobject.H:360
An IOstream is an abstract base class for all input/output systems; be they streams,...
Definition: IOstream.H:72
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:164
Initialise the NamedEnum HashTable from the static list of names.
Definition: NamedEnum.H:54
A simple wrapper around bool so that it can be read as a word: true/false, on/off,...
Definition: Switch.H:61
A class for addressing time paths without using the Time class.
Definition: TimePaths.H:49
The time value with time-stepping information, user-defined remapping, etc.
Definition: TimeState.H:55
friend class Time
Declare friendship with the Time class.
Definition: TimeState.H:70
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:76
virtual bool run() const
Return true if run should continue,.
Definition: Time.C:855
static int precision_
Time directory name precision.
Definition: Time.H:161
scalar timeToUserTime(const scalar t) const
Convert the real-time (s) into user-time (e.g. CA deg)
Definition: Time.C:830
static const NamedEnum< stopAtControl, 4 > stopAtControlNames_
Definition: Time.H:128
instantList times() const
Search the case for valid time directories.
Definition: Time.C:640
virtual void setTime(const Time &)
Reset the time and time-index to those of the given time.
Definition: Time.C:939
virtual dimensionedScalar startTime() const
Return start time.
Definition: Time.C:790
const word & timeName() const
Return current time name (for backwards-compatibility)
Definition: Time.H:382
static int curPrecision_
Current time directory name precision adjusted as necessary.
Definition: Time.H:166
void adjustDeltaT()
Adjust the time step so that writing occurs at the specified time.
Definition: Time.C:83
format
Supported time directory name formats.
Definition: Time.H:109
scalar userTimeToTime(const scalar tau) const
Convert the user-time (e.g. CA deg) to real-time (s).
Definition: Time.C:824
scalar writeInterval_
Definition: Time.H:134
virtual dimensionedScalar beginTime() const
Return begin time (initial start time)
Definition: Time.C:779
static const NamedEnum< writeControl, 5 > writeControlNames_
Definition: Time.H:131
scalar userTimeValue() const
Return current user time value.
Definition: Time.C:818
static const int maxPrecision_
Maximum time directory name precision.
Definition: Time.H:169
word findInstance(const fileName &dir, const word &name=word::null, const IOobject::readOption rOpt=IOobject::MUST_READ, const word &stopInstance=word::null) const
Return the location of "dir" containing the file "name".
Definition: Time.C:647
virtual void setWriteInterval(const scalar writeInterval)
Reset the write interval.
Definition: Time.C:1035
virtual void setEndTime(const dimensionedScalar &)
Reset end time.
Definition: Time.C:999
word userTimeName() const
Return current user time name with units.
Definition: Time.C:836
virtual bool running() const
Return true if run should continue without any side effects.
Definition: Time.C:849
void setControls()
Set the controls from the current controlDict.
Definition: Time.C:109
word findInstancePath(const fileName &path, const instant &) const
Search the case for the time directory path.
Definition: Time.C:677
virtual void setDeltaT(const dimensionedScalar &)
Reset time step.
Definition: Time.C:1011
static word controlDictName
The default control dictionary name (normally "controlDict")
Definition: Time.H:204
stopAtControl
Stop-run control options.
Definition: Time.H:100
static format format_
Time directory name format.
Definition: Time.H:158
virtual bool stopAt(const stopAtControl) const
Adjust the current stopAtControl. Note that this value.
Definition: Time.C:921
const userTimes::userTime & userTime() const
Return the userTime.
Definition: Time.C:812
virtual label startTimeIndex() const
Return start time index.
Definition: Time.C:773
virtual void endSubCycle()
Reset time after sub-cycling back to previous TimeState.
Definition: Time.C:1077
virtual void setDeltaTNoAdjust(const scalar)
Reset time step without additional adjustment or modification.
Definition: Time.C:1028
virtual Time & operator++()
Prefix increment,.
Definition: Time.C:1103
void readModifiedObjects()
Read the objects that have been modified.
Definition: TimeIO.C:226
scalar beginTime_
Definition: Time.H:123
virtual TimeState subCycle(const label nSubCycles)
Set time to sub-cycle for the given number of steps.
Definition: Time.C:1063
writeControl
Write control options.
Definition: Time.H:90
virtual bool loop()
Return true if run should continue and if so increment time.
Definition: Time.C:902
virtual dimensionedScalar endTime() const
Return end time.
Definition: Time.C:801
virtual ~Time()
Destructor.
Definition: Time.C:621
static label findClosestTimeIndex(const instantList &, const scalar, const word &constantName="constant")
Search instantList for the time index closest to the given time.
Definition: Time.C:748
virtual Time & operator+=(const dimensionedScalar &)
Set deltaT to that specified and increment time via operator++()
Definition: Time.C:1090
virtual bool end() const
Return true if end of run,.
Definition: Time.C:915
instant findClosestTime(const scalar) const
Search the case for the time closest to the given time.
Definition: Time.C:711
T & last()
Return the last element of the list.
Definition: UListI.H:128
static label nProcs(const label communicator=0)
Number of processes in parallel run.
Definition: UPstream.H:411
static bool & parRun()
Is this a parallel run?
Definition: UPstream.H:399
static int myProcNo(const label communicator=0)
Number of this process (starting from masterNo() = 0)
Definition: UPstream.H:429
Extract command arguments and options from the supplied argc and argv parameters.
Definition: argList.H:103
const word & executable() const
Name of executable without the path.
Definition: argListI.H:36
const Foam::HashTable< string > & options() const
Return options.
Definition: argListI.H:96
fileName path() const
Return the path to the caseName.
Definition: argListI.H:66
A list of keyword definitions, which are a keyword followed by any number of values (e....
Definition: dictionary.H:160
ITstream & lookup(const word &, bool recursive=false, bool patternMatch=true) const
Find and return an entry data stream.
Definition: dictionary.C:860
bool readIfPresent(const word &, T &, bool recursive=false, bool patternMatch=true) const
Find an entry if present, and assign to T.
bool found(const word &, bool recursive=false, bool patternMatch=true) const
Search dictionary for given keyword.
Definition: dictionary.C:659
const scalar & value() const
Return const reference to value.
const word & name() const
Return const reference to name.
bool open(const fileName &libName, const bool verbose=true)
Open the named library, optionally with warnings if problems occur.
A class for handling file names.
Definition: fileName.H:82
virtual void setTime(const Time &) const
Callback for time change.
scalar timeToNextAction()
Return the time to the next write.
An instant of time. Contains the time value and name.
Definition: instant.H:67
scalar value() const
Value (const access)
Definition: instant.H:114
const word & name() const
Name (const access)
Definition: instant.H:126
bool equal(const scalar) const
Comparison used for instants to be equal.
Definition: instant.C:59
Registry of regIOobjects.
void addWatch()
Add file watch on object (if registered and READ_IF_MODIFIED)
Definition: regIOobject.C:259
timeIOdictionary derived from IOdictionary with globalFile set false to enable writing to processor t...
An abstract class for the user time description.
A class for handling words, derived from string.
Definition: word.H:62
static const word null
An empty word.
Definition: word.H:77
#define FatalIOErrorInFunction(ios)
Report an error message using Foam::FatalIOError.
Definition: error.H:318
static instantList timeDirs
Definition: globalFoam.H:44
runTimeSource setTime(sourceTimes[sourceTimeIndex], sourceTimeIndex)
word timeName
Definition: getTimeIndex.H:3
#define IOWarningInFunction(ios)
Report an IO warning using Foam::Warning.
#define WarningInFunction
Report a warning using Foam::Warning.
Namespace for OpenFOAM.
const fileOperation & fileHandler()
Get current file handler.
errorManipArg< error, int > exit(error &err, const int errNo=1)
Definition: errorManip.H:124
dlLibraryTable libs
Table of loaded dynamic libraries.
dimensioned< scalar > dimensionedScalar
Dimensioned scalar obtained from generic dimensioned type.
void adjustDeltaT(Time &runTime, const PtrList< solver > &solvers)
Adjust the time-step according to the solver maxDeltaT.
int system(const std::string &command)
Execute the specified command.
Definition: POSIX.C:1230
bool equal(const T &s1, const T &s2)
Definition: doubleFloat.H:62
void setDeltaT(Time &runTime, const PtrList< solver > &solvers)
Set the initial time-step according to the solver maxDeltaT.
dimensionedScalar sign(const dimensionedScalar &ds)
intWM_LABEL_SIZE_t label
A label is an int32_t or int64_t as specified by the pre-processor macro WM_LABEL_SIZE.
Definition: label.H:59
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:251
bool writeInfoHeader
bool exists(const fileName &, const bool checkVariants=true, const bool followLink=true)
Does the name exist (as directory or file) in the file system?
Definition: POSIX.C:520
dimensionedScalar log10(const dimensionedScalar &ds)
layerAndWeight min(const layerAndWeight &a, const layerAndWeight &b)
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
const dimensionSet dimTime
void reduce(const List< UPstream::commsStruct > &comms, T &Value, const BinaryOp &bop, const int tag, const label comm)
bool readScalar(const char *buf, doubleScalar &s)
Read whole of buf as a scalar. Return true if successful.
Definition: doubleScalar.H:75
T returnReduce(const T &Value, const BinaryOp &bop, const int tag=Pstream::msgType(), const label comm=UPstream::worldComm)
scalar diff(const triad &A, const triad &B)
Return a quantity of the difference between two triads.
Definition: triad.C:403
dimensioned< scalar > mag(const dimensioned< Type > &)
defineTypeNameAndDebug(combustionModel, 0)
tmp< DimensionedField< TypeR, GeoMesh > > New(const tmp< DimensionedField< TypeR, GeoMesh >> &tdf1, const word &name, const dimensionSet &dimensions)
layerAndWeight max(const layerAndWeight &a, const layerAndWeight &b)
IOerror FatalIOError
static const label labelMax
Definition: label.H:62
word name(const complex &)
Return a string representation of a complex.
Definition: complex.C:47
static const char nl
Definition: Ostream.H:260
label timeIndex
Definition: getTimeIndex.H:4
dictionary dict
Foam::argList args(argc, argv)