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-2018 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  libs_(),
335 
336  controlDict_
337  (
338  IOobject
339  (
340  controlDictName,
341  system(),
342  *this,
345  false
346  )
347  ),
348 
349  startTimeIndex_(0),
350  startTime_(0),
351  endTime_(0),
352 
353  stopAt_(stopAtControl::endTime),
354  writeControl_(writeControl::timeStep),
355  writeInterval_(great),
356  purgeWrite_(0),
357  writeOnce_(false),
358  subCycling_(false),
359  sigWriteNow_(true, *this),
360  sigStopAtWriteNow_(true, *this),
361 
362  writeFormat_(IOstream::ASCII),
363  writeVersion_(IOstream::currentVersion),
364  writeCompression_(IOstream::UNCOMPRESSED),
365  graphFormat_("raw"),
366  runTimeModifiable_(false),
367 
368  functionObjects_(*this, enableFunctionObjects)
369 {
370  libs_.open(controlDict_, "libs");
371 
372  // Explicitly set read flags on objectRegistry so anything constructed
373  // from it reads as well (e.g. fvSolution).
375 
376  setControls();
377 
378  // Time objects not registered so do like objectRegistry::checkIn ourselves.
379  if (runTimeModifiable_)
380  {
381  // Monitor all files that controlDict depends on
382  fileHandler().addWatches(controlDict_, controlDict_.files());
383  }
384 
385  // Clear dependent files
386  controlDict_.files().clear();
387 }
388 
389 
391 (
392  const word& controlDictName,
393  const argList& args,
394  const word& systemName,
395  const word& constantName
396 )
397 :
398  TimePaths
399  (
400  args.parRunControl().parRun(),
401  args.rootPath(),
402  args.globalCaseName(),
403  args.caseName(),
404  systemName,
405  constantName
406  ),
407 
408  objectRegistry(*this),
409 
410  libs_(),
411 
412  controlDict_
413  (
414  IOobject
415  (
416  controlDictName,
417  system(),
418  *this,
421  false
422  )
423  ),
424 
425  startTimeIndex_(0),
426  startTime_(0),
427  endTime_(0),
428 
429  stopAt_(stopAtControl::endTime),
430  writeControl_(writeControl::timeStep),
431  writeInterval_(great),
432  purgeWrite_(0),
433  writeOnce_(false),
434  subCycling_(false),
435  sigWriteNow_(true, *this),
436  sigStopAtWriteNow_(true, *this),
437 
438  writeFormat_(IOstream::ASCII),
439  writeVersion_(IOstream::currentVersion),
440  writeCompression_(IOstream::UNCOMPRESSED),
441  graphFormat_("raw"),
442  runTimeModifiable_(false),
443 
444  functionObjects_
445  (
446  *this,
447  argList::validOptions.found("withFunctionObjects")
448  ? args.optionFound("withFunctionObjects")
449  : !args.optionFound("noFunctionObjects")
450  )
451 {
452  libs_.open(controlDict_, "libs");
453 
454  // Explicitly set read flags on objectRegistry so anything constructed
455  // from it reads as well (e.g. fvSolution).
457 
458  setControls();
459 
460  // Time objects not registered so do like objectRegistry::checkIn ourselves.
461  if (runTimeModifiable_)
462  {
463  // Monitor all files that controlDict depends on
464  fileHandler().addWatches(controlDict_, controlDict_.files());
465  }
466 
467  // Clear dependent files since not needed
468  controlDict_.files().clear();
469 }
470 
471 
473 (
474  const dictionary& dict,
475  const fileName& rootPath,
476  const fileName& caseName,
477  const word& systemName,
478  const word& constantName,
479  const bool enableFunctionObjects
480 )
481 :
482  TimePaths
483  (
484  rootPath,
485  caseName,
486  systemName,
487  constantName
488  ),
489 
490  objectRegistry(*this),
491 
492  libs_(),
493 
494  controlDict_
495  (
496  IOobject
497  (
498  controlDictName,
499  system(),
500  *this,
503  false
504  ),
505  dict
506  ),
507 
508  startTimeIndex_(0),
509  startTime_(0),
510  endTime_(0),
511 
512  stopAt_(stopAtControl::endTime),
513  writeControl_(writeControl::timeStep),
514  writeInterval_(great),
515  purgeWrite_(0),
516  writeOnce_(false),
517  subCycling_(false),
518  sigWriteNow_(true, *this),
519  sigStopAtWriteNow_(true, *this),
520 
521  writeFormat_(IOstream::ASCII),
522  writeVersion_(IOstream::currentVersion),
523  writeCompression_(IOstream::UNCOMPRESSED),
524  graphFormat_("raw"),
525  runTimeModifiable_(false),
526 
527  functionObjects_(*this, enableFunctionObjects)
528 {
529  libs_.open(controlDict_, "libs");
530 
531 
532  // Explicitly set read flags on objectRegistry so anything constructed
533  // from it reads as well (e.g. fvSolution).
535 
536  // Since could not construct regIOobject with setting:
537  controlDict_.readOpt() = IOobject::MUST_READ_IF_MODIFIED;
538 
539  setControls();
540 
541  // Time objects not registered so do like objectRegistry::checkIn ourselves.
542  if (runTimeModifiable_)
543  {
544  // Monitor all files that controlDict depends on
545  fileHandler().addWatches(controlDict_, controlDict_.files());
546  }
547 
548  // Clear dependent files since not needed
549  controlDict_.files().clear();
550 }
551 
552 
554 (
555  const fileName& rootPath,
556  const fileName& caseName,
557  const word& systemName,
558  const word& constantName,
559  const bool enableFunctionObjects
560 )
561 :
562  TimePaths
563  (
564  rootPath,
565  caseName,
566  systemName,
567  constantName
568  ),
569 
570  objectRegistry(*this),
571 
572  libs_(),
573 
574  controlDict_
575  (
576  IOobject
577  (
578  controlDictName,
579  system(),
580  *this,
583  false
584  )
585  ),
586 
587  startTimeIndex_(0),
588  startTime_(0),
589  endTime_(0),
590 
591  stopAt_(stopAtControl::endTime),
592  writeControl_(writeControl::timeStep),
593  writeInterval_(great),
594  purgeWrite_(0),
595  writeOnce_(false),
596  subCycling_(false),
597 
598  writeFormat_(IOstream::ASCII),
599  writeVersion_(IOstream::currentVersion),
600  writeCompression_(IOstream::UNCOMPRESSED),
601  graphFormat_("raw"),
602  runTimeModifiable_(false),
603 
604  functionObjects_(*this, enableFunctionObjects)
605 {
606  libs_.open(controlDict_, "libs");
607 }
608 
609 
610 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
611 
613 {
614  forAllReverse(controlDict_.watchIndices(), i)
615  {
616  fileHandler().removeWatch(controlDict_.watchIndices()[i]);
617  }
618 
619  // Destroy function objects first
620  functionObjects_.clear();
621 }
622 
623 
624 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
625 
626 Foam::word Foam::Time::timeName(const scalar t, const int precision)
627 {
628  std::ostringstream buf;
629  buf.setf(ios_base::fmtflags(format_), ios_base::floatfield);
630  buf.precision(precision);
631  buf << t;
632  return buf.str();
633 }
634 
635 
637 {
638  return dimensionedScalar::name();
639 }
640 
641 
643 {
644  return findTimes(path(), constant());
645 }
646 
647 
649 (
650  const fileName& dir,
651  const word& name,
652  const IOobject::readOption rOpt,
653  const word& stopInstance
654 ) const
655 {
656  IOobject startIO
657  (
658  name, // name might be empty!
659  timeName(),
660  dir,
661  *this,
662  rOpt
663  );
664 
665  IOobject io
666  (
667  fileHandler().findInstance
668  (
669  startIO,
670  timeOutputValue(),
671  stopInstance
672  )
673  );
674  return io.instance();
675 }
676 
677 
679 (
680  const fileName& directory,
681  const instant& t
682 ) const
683 {
684  // Simplified version: use findTimes (readDir + sort). The expensive
685  // bit is the readDir, not the sorting. Tbd: avoid calling findInstancePath
686  // from filePath.
687 
688  instantList timeDirs = findTimes(path(), constant());
689  // Note:
690  // - times will include constant (with value 0) as first element.
691  // For backwards compatibility make sure to find 0 in preference
692  // to constant.
693  // - list is sorted so could use binary search
694 
695  forAllReverse(timeDirs, i)
696  {
697  if (t.equal(timeDirs[i].value()))
698  {
699  return timeDirs[i].name();
700  }
701  }
702 
703  return word::null;
704 }
705 
706 
708 {
709  return findInstancePath(path(), t);
710 }
711 
712 
714 {
715  instantList timeDirs = findTimes(path(), constant());
716 
717  // There is only one time (likely "constant") so return it
718  if (timeDirs.size() == 1)
719  {
720  return timeDirs[0];
721  }
722 
723  if (t < timeDirs[1].value())
724  {
725  return timeDirs[1];
726  }
727  else if (t > timeDirs.last().value())
728  {
729  return timeDirs.last();
730  }
731 
732  label nearestIndex = -1;
733  scalar deltaT = great;
734 
735  for (label timei=1; timei < timeDirs.size(); ++timei)
736  {
737  scalar diff = mag(timeDirs[timei].value() - t);
738  if (diff < deltaT)
739  {
740  deltaT = diff;
741  nearestIndex = timei;
742  }
743  }
744 
745  return timeDirs[nearestIndex];
746 }
747 
748 
750 (
751  const instantList& timeDirs,
752  const scalar t,
753  const word& constantName
754 )
755 {
756  label nearestIndex = -1;
757  scalar deltaT = great;
758 
759  forAll(timeDirs, timei)
760  {
761  if (timeDirs[timei].name() == constantName) continue;
762 
763  scalar diff = mag(timeDirs[timei].value() - t);
764  if (diff < deltaT)
765  {
766  deltaT = diff;
767  nearestIndex = timei;
768  }
769  }
770 
771  return nearestIndex;
772 }
773 
774 
776 {
777  return startTimeIndex_;
778 }
779 
780 
782 {
783  return dimensionedScalar("startTime", dimTime, startTime_);
784 }
785 
786 
788 {
789  return dimensionedScalar("endTime", dimTime, endTime_);
790 }
791 
792 
794 {
795  return value() < (endTime_ - 0.5*deltaT_);
796 }
797 
798 
799 bool Foam::Time::run() const
800 {
801  bool running = this->running();
802 
803  if (!subCycling_)
804  {
805  if (!running && timeIndex_ != startTimeIndex_)
806  {
807  functionObjects_.execute();
808  functionObjects_.end();
809  }
810  }
811 
812  if (running)
813  {
814  if (!subCycling_)
815  {
816  const_cast<Time&>(*this).readModifiedObjects();
817 
818  if (timeIndex_ == startTimeIndex_)
819  {
820  functionObjects_.start();
821  }
822  else
823  {
824  functionObjects_.execute();
825  }
826  }
827 
828  // Re-evaluate if running in case a function object has changed things
829  running = this->running();
830  }
831 
832  return running;
833 }
834 
835 
837 {
838  bool running = run();
839 
840  if (running)
841  {
842  operator++();
843  }
844 
845  return running;
846 }
847 
848 
849 bool Foam::Time::end() const
850 {
851  return value() > (endTime_ + 0.5*deltaT_);
852 }
853 
854 
855 bool Foam::Time::stopAt(const stopAtControl sa) const
856 {
857  const bool changed = (stopAt_ != sa);
858  stopAt_ = sa;
859 
860  // adjust endTime
861  if (sa == stopAtControl::endTime)
862  {
863  controlDict_.lookup("endTime") >> endTime_;
864  }
865  else
866  {
867  endTime_ = great;
868  }
869  return changed;
870 }
871 
872 
873 void Foam::Time::setTime(const Time& t)
874 {
875  value() = t.value();
876  dimensionedScalar::name() = t.dimensionedScalar::name();
877  timeIndex_ = t.timeIndex_;
878  fileHandler().setTime(*this);
879 }
880 
881 
882 void Foam::Time::setTime(const instant& inst, const label newIndex)
883 {
884  value() = inst.value();
885  dimensionedScalar::name() = inst.name();
886  timeIndex_ = newIndex;
887 
888  IOdictionary timeDict
889  (
890  IOobject
891  (
892  "time",
893  timeName(),
894  "uniform",
895  *this,
898  false
899  )
900  );
901 
902  timeDict.readIfPresent("deltaT", deltaT_);
903  timeDict.readIfPresent("deltaT0", deltaT0_);
904  timeDict.readIfPresent("index", timeIndex_);
905  fileHandler().setTime(*this);
906 }
907 
908 
909 void Foam::Time::setTime(const dimensionedScalar& newTime, const label newIndex)
910 {
911  setTime(newTime.value(), newIndex);
912 }
913 
914 
915 void Foam::Time::setTime(const scalar newTime, const label newIndex)
916 {
917  value() = newTime;
918  dimensionedScalar::name() = timeName(timeToUserTime(newTime));
919  timeIndex_ = newIndex;
920  fileHandler().setTime(*this);
921 }
922 
923 
925 {
926  setEndTime(endTime.value());
927 }
928 
929 
930 void Foam::Time::setEndTime(const scalar endTime)
931 {
932  endTime_ = endTime;
933 }
934 
935 
937 {
938  setDeltaT(deltaT.value());
939 }
940 
941 
942 void Foam::Time::setDeltaT(const scalar deltaT)
943 {
944  setDeltaTNoAdjust(deltaT);
945 
946  functionObjects_.setTimeStep();
947 
948  if (writeControl_ == writeControl::adjustableRunTime)
949  {
950  adjustDeltaT();
951  }
952 }
953 
954 
955 void Foam::Time::setDeltaTNoAdjust(const scalar deltaT)
956 {
957  deltaT_ = deltaT;
958  deltaTchanged_ = true;
959 }
960 
961 
963 {
964  subCycling_ = true;
965  prevTimeState_.set(new TimeState(*this));
966 
967  setTime(*this - deltaT(), (timeIndex() - 1)*nSubCycles);
968  deltaT_ /= nSubCycles;
969  deltaT0_ /= nSubCycles;
970  deltaTSave_ = deltaT0_;
971 
972  return prevTimeState();
973 }
974 
975 
977 {
978  if (subCycling_)
979  {
980  subCycling_ = false;
981  TimeState::operator=(prevTimeState());
982  prevTimeState_.clear();
983  }
984 }
985 
986 
987 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
988 
990 {
991  return operator+=(deltaT.value());
992 }
993 
994 
995 Foam::Time& Foam::Time::operator+=(const scalar deltaT)
996 {
997  setDeltaT(deltaT);
998  return operator++();
999 }
1000 
1001 
1003 {
1004  deltaT0_ = deltaTSave_;
1005  deltaTSave_ = deltaT_;
1006 
1007  // Save old time value and name
1008  const scalar oldTimeValue = timeToUserTime(value());
1009  const word oldTimeName = dimensionedScalar::name();
1010 
1011  // Increment time
1012  setTime(value() + deltaT_, timeIndex_ + 1);
1013 
1014  if (!subCycling_)
1015  {
1016  // If the time is very close to zero reset to zero
1017  if (mag(value()) < 10*small*deltaT_)
1018  {
1019  setTime(0, timeIndex_);
1020  }
1021 
1022  if (sigStopAtWriteNow_.active() || sigWriteNow_.active())
1023  {
1024  // A signal might have been sent on one processor only
1025  // Reduce so all decide the same.
1026 
1027  label flag = 0;
1028  if
1029  (
1030  sigStopAtWriteNow_.active()
1031  && stopAt_ == stopAtControl::writeNow
1032  )
1033  {
1034  flag += 1;
1035  }
1036  if (sigWriteNow_.active() && writeOnce_)
1037  {
1038  flag += 2;
1039  }
1040  reduce(flag, maxOp<label>());
1041 
1042  if (flag & 1)
1043  {
1044  stopAt_ = stopAtControl::writeNow;
1045  }
1046  if (flag & 2)
1047  {
1048  writeOnce_ = true;
1049  }
1050  }
1051 
1052  writeTime_ = false;
1053 
1054  switch (writeControl_)
1055  {
1056  case writeControl::timeStep:
1057  writeTime_ = !(timeIndex_ % label(writeInterval_));
1058  break;
1059 
1060  case writeControl::runTime:
1061  case writeControl::adjustableRunTime:
1062  {
1063  label writeIndex = label
1064  (
1065  ((value() - startTime_) + 0.5*deltaT_)
1066  / writeInterval_
1067  );
1068 
1069  if (writeIndex > writeTimeIndex_)
1070  {
1071  writeTime_ = true;
1072  writeTimeIndex_ = writeIndex;
1073  }
1074  }
1075  break;
1076 
1077  case writeControl::cpuTime:
1078  {
1079  label writeIndex = label
1080  (
1081  returnReduce(elapsedCpuTime(), maxOp<double>())
1082  / writeInterval_
1083  );
1084  if (writeIndex > writeTimeIndex_)
1085  {
1086  writeTime_ = true;
1087  writeTimeIndex_ = writeIndex;
1088  }
1089  }
1090  break;
1091 
1092  case writeControl::clockTime:
1093  {
1094  label writeIndex = label
1095  (
1096  returnReduce(label(elapsedClockTime()), maxOp<label>())
1097  / writeInterval_
1098  );
1099  if (writeIndex > writeTimeIndex_)
1100  {
1101  writeTime_ = true;
1102  writeTimeIndex_ = writeIndex;
1103  }
1104  }
1105  break;
1106  }
1107 
1108 
1109  // Check if endTime needs adjustment to stop at the next run()/end()
1110  if (!end())
1111  {
1112  if (stopAt_ == stopAtControl::noWriteNow)
1113  {
1114  endTime_ = value();
1115  }
1116  else if (stopAt_ == stopAtControl::writeNow)
1117  {
1118  endTime_ = value();
1119  writeTime_ = true;
1120  }
1121  else if (stopAt_ == stopAtControl::nextWrite && writeTime_ == true)
1122  {
1123  endTime_ = value();
1124  }
1125  }
1126 
1127  // Override writeTime if one-shot writing
1128  if (writeOnce_)
1129  {
1130  writeTime_ = true;
1131  writeOnce_ = false;
1132  }
1133 
1134  // Adjust the precision of the time directory name if necessary
1135  if (writeTime_)
1136  {
1137  // Tolerance used when testing time equivalence
1138  const scalar timeTol =
1139  max(min(pow(10.0, -precision_), 0.1*deltaT_), small);
1140 
1141  // User-time equivalent of deltaT
1142  const scalar userDeltaT = timeToUserTime(deltaT_);
1143 
1144  // Time value obtained by reading timeName
1145  scalar timeNameValue = -vGreat;
1146 
1147  // Check that new time representation differs from old one
1148  // reinterpretation of the word
1149  if
1150  (
1151  readScalar(dimensionedScalar::name().c_str(), timeNameValue)
1152  && (mag(timeNameValue - oldTimeValue - userDeltaT) > timeTol)
1153  )
1154  {
1155  int oldPrecision = precision_;
1156  while
1157  (
1158  precision_ < maxPrecision_
1159  && readScalar(dimensionedScalar::name().c_str(), timeNameValue)
1160  && (mag(timeNameValue - oldTimeValue - userDeltaT) > timeTol)
1161  )
1162  {
1163  precision_++;
1164  setTime(value(), timeIndex());
1165  }
1166 
1167  if (precision_ != oldPrecision)
1168  {
1170  << "Increased the timePrecision from " << oldPrecision
1171  << " to " << precision_
1172  << " to distinguish between timeNames at time "
1174  << endl;
1175 
1176  if (precision_ == maxPrecision_)
1177  {
1178  // Reached maxPrecision limit
1180  << "Current time name " << dimensionedScalar::name()
1181  << nl
1182  << " The maximum time precision has been reached"
1183  " which might result in overwriting previous"
1184  " results."
1185  << endl;
1186  }
1187 
1188  // Check if round-off error caused time-reversal
1189  scalar oldTimeNameValue = -vGreat;
1190  if
1191  (
1192  readScalar(oldTimeName.c_str(), oldTimeNameValue)
1193  && (
1194  sign(timeNameValue - oldTimeNameValue)
1195  != sign(deltaT_)
1196  )
1197  )
1198  {
1200  << "Current time name " << dimensionedScalar::name()
1201  << " is set to an instance prior to the "
1202  "previous one "
1203  << oldTimeName << nl
1204  << " This might result in temporal "
1205  "discontinuities."
1206  << endl;
1207  }
1208  }
1209  }
1210  }
1211  }
1212 
1213  return *this;
1214 }
1215 
1216 
1218 {
1219  return operator++();
1220 }
1221 
1222 
1223 // ************************************************************************* //
virtual Time & operator++()
Prefix increment,.
Definition: Time.C:1002
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:48
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:612
const fileName & rootPath() const
Return root path.
Definition: argListI.H:36
void size(const label)
Override size to be inconsistent with allocated storage.
Definition: ListI.H:163
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:836
virtual dimensionedScalar startTime() const
Return start time.
Definition: Time.C:781
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:256
bool optionFound(const word &opt) const
Return true if the named option is found.
Definition: argListI.H:108
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:649
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:787
stopAtControl
Stop-run control options.
Definition: Time.H:98
readOption
Enumeration defining the read options.
Definition: IOobject.H:107
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:636
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:68
virtual bool run() const
Return true if run should continue,.
Definition: Time.C:799
const ParRunControl & parRunControl() const
Return parRunControl.
Definition: argListI.H:54
static const int maxPrecision_
Maximum time directory name precision.
Definition: Time.H:160
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:962
virtual bool stopAt(const stopAtControl) const
Adjust the current stopAtControl. Note that this value.
Definition: Time.C:855
static int precision_
Time directory name precision.
Definition: Time.H:157
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:793
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:873
format
Supported time directory name formats.
Definition: Time.H:107
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:68
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:391
static word controlDictName
The default control dictionary name (normally "controlDict")
Definition: Time.H:198
const fileName & caseName() const
Return case name (parallel run) or global case (serial run)
Definition: argListI.H:42
static format format_
Time directory name format.
Definition: Time.H:154
static const char nl
Definition: Ostream.H:265
virtual void setDeltaTNoAdjust(const scalar)
Reset time step without additional adjustment or modification.
Definition: Time.C:955
defineTypeNameAndDebug(combustionModel, 0)
scalar value() const
Value (const access)
Definition: instant.H:114
static const NamedEnum< stopAtControl, 4 > stopAtControlNames_
Definition: Time.H:123
const word & name() const
Return const reference to name.
dimensioned< Type > min(const dimensioned< Type > &, const dimensioned< Type > &)
void readModifiedObjects()
Read the objects that have been modified.
Definition: TimeIO.C:460
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:924
static const NamedEnum< writeControl, 5 > writeControlNames_
Definition: Time.H:126
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:378
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:849
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:642
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:713
virtual void setDeltaT(const dimensionedScalar &)
Reset time step.
Definition: Time.C:936
writeControl
Write control options.
Definition: Time.H:88
dimensioned< scalar > mag(const dimensioned< Type > &)
virtual void endSubCycle()
Reset time after sub-cycling back to previous TimeState.
Definition: Time.C:976
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:679
virtual label startTimeIndex() const
Return start time index.
Definition: Time.C:775
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:989
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:750
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