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