1 /*---------------------------------------------------------------------------*\
3 \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
5 \\ / A nd | Copyright (C) 2011 OpenFOAM Foundation
7 -------------------------------------------------------------------------------
9 This file is part of OpenFOAM.
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.
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
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/>.
24 \*---------------------------------------------------------------------------*/
27 #include "PstreamReduceOps.H"
32 // * * * * * * * * * * * * * Static Member Data * * * * * * * * * * * * * * //
34 defineTypeNameAndDebug(Foam::Time, 0);
39 const char* Foam::NamedEnum
41 Foam::Time::stopAtControls,
52 const char* Foam::NamedEnum
54 Foam::Time::writeControls,
66 const Foam::NamedEnum<Foam::Time::stopAtControls, 4>
67 Foam::Time::stopAtControlNames_;
69 const Foam::NamedEnum<Foam::Time::writeControls, 5>
70 Foam::Time::writeControlNames_;
72 Foam::Time::fmtflags Foam::Time::format_(Foam::Time::general);
73 int Foam::Time::precision_(6);
75 Foam::word Foam::Time::controlDictName("controlDict");
78 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
80 void Foam::Time::adjustDeltaT()
82 if (writeControl_ == wcAdjustableRunTime)
84 scalar timeToNextWrite = max
87 (outputTimeIndex_ + 1)*writeInterval_ - (value() - startTime_)
90 scalar nSteps = timeToNextWrite/deltaT_ - SMALL;
92 // For tiny deltaT the label can overflow!
93 if (nSteps < labelMax)
95 label nStepsToNextWrite = label(nSteps) + 1;
97 scalar newDeltaT = timeToNextWrite/nStepsToNextWrite;
99 // Control the increase of the time step to within a factor of 2
100 // and the decrease within a factor of 5.
101 if (newDeltaT >= deltaT_)
103 deltaT_ = min(newDeltaT, 2.0*deltaT_);
107 deltaT_ = max(newDeltaT, 0.2*deltaT_);
114 void Foam::Time::setControls()
116 // default is to resume calculation from "latestTime"
117 const word startFrom = controlDict_.lookupOrDefault<word>
123 if (startFrom == "startTime")
125 controlDict_.lookup("startTime") >> startTime_;
129 // Search directory for valid time directories
130 instantList timeDirs = findTimes(path());
132 if (startFrom == "firstTime")
136 if (timeDirs[0].name() == constant() && timeDirs.size() >= 2)
138 startTime_ = timeDirs[1].value();
142 startTime_ = timeDirs[0].value();
146 else if (startFrom == "latestTime")
150 startTime_ = timeDirs.last().value();
155 FatalIOErrorIn("Time::setControls()", controlDict_)
156 << "expected startTime, firstTime or latestTime"
157 << " found '" << startFrom << "'"
158 << exit(FatalIOError);
162 setTime(startTime_, 0);
165 deltaTSave_ = deltaT_;
168 if (Pstream::parRun())
170 scalar sumStartTime = startTime_;
171 reduce(sumStartTime, sumOp<scalar>());
174 mag(Pstream::nProcs()*startTime_ - sumStartTime)
175 > Pstream::nProcs()*deltaT_/10.0
178 FatalIOErrorIn("Time::setControls()", controlDict_)
179 << "Start time is not the same for all processors" << nl
180 << "processor " << Pstream::myProcNo() << " has startTime "
181 << startTime_ << exit(FatalIOError);
185 IOdictionary timeDict
193 IOobject::READ_IF_PRESENT,
199 if (timeDict.readIfPresent("deltaT", deltaT_))
201 deltaTSave_ = deltaT_;
205 timeDict.readIfPresent("deltaT0", deltaT0_);
207 if (timeDict.readIfPresent("index", startTimeIndex_))
209 timeIndex_ = startTimeIndex_;
214 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
218 const word& controlDictName,
219 const fileName& rootPath,
220 const fileName& caseName,
221 const word& systemName,
222 const word& constantName,
223 const bool enableFunctionObjects
234 objectRegistry(*this),
245 IOobject::MUST_READ_IF_MODIFIED,
256 writeControl_(wcTimeStep),
257 writeInterval_(GREAT),
261 writeFormat_(IOstream::ASCII),
262 writeVersion_(IOstream::currentVersion),
263 writeCompression_(IOstream::UNCOMPRESSED),
265 runTimeModifiable_(true),
267 functionObjects_(*this, enableFunctionObjects)
269 libs_.open(controlDict_, "libs");
271 // Explicitly set read flags on objectRegistry so anything constructed
272 // from it reads as well (e.g. fvSolution).
273 readOpt() = IOobject::MUST_READ_IF_MODIFIED;
277 // Time objects not registered so do like objectRegistry::checkIn ourselves.
278 if (runTimeModifiable_)
284 regIOobject::fileModificationChecking == inotify
285 || regIOobject::fileModificationChecking == inotifyMaster
289 // File might not exist yet.
290 fileName f(controlDict_.filePath());
294 // We don't have this file but would like to re-read it.
295 // Possibly if in master-only reading mode. Use a non-existing
296 // file to keep fileMonitor synced.
297 f = controlDict_.objectPath();
300 controlDict_.watchIndex() = addWatch(f);
307 const word& controlDictName,
309 const word& systemName,
310 const word& constantName
321 objectRegistry(*this),
332 IOobject::MUST_READ_IF_MODIFIED,
343 writeControl_(wcTimeStep),
344 writeInterval_(GREAT),
348 writeFormat_(IOstream::ASCII),
349 writeVersion_(IOstream::currentVersion),
350 writeCompression_(IOstream::UNCOMPRESSED),
352 runTimeModifiable_(true),
354 functionObjects_(*this, !args.optionFound("noFunctionObjects"))
356 libs_.open(controlDict_, "libs");
358 // Explicitly set read flags on objectRegistry so anything constructed
359 // from it reads as well (e.g. fvSolution).
360 readOpt() = IOobject::MUST_READ_IF_MODIFIED;
364 // Time objects not registered so do like objectRegistry::checkIn ourselves.
365 if (runTimeModifiable_)
371 regIOobject::fileModificationChecking == inotify
372 || regIOobject::fileModificationChecking == inotifyMaster
376 // File might not exist yet.
377 fileName f(controlDict_.filePath());
381 // We don't have this file but would like to re-read it.
382 // Possibly if in master-only reading mode. Use a non-existing
383 // file to keep fileMonitor synced.
384 f = controlDict_.objectPath();
387 controlDict_.watchIndex() = addWatch(f);
394 const dictionary& dict,
395 const fileName& rootPath,
396 const fileName& caseName,
397 const word& systemName,
398 const word& constantName,
399 const bool enableFunctionObjects
410 objectRegistry(*this),
433 writeControl_(wcTimeStep),
434 writeInterval_(GREAT),
438 writeFormat_(IOstream::ASCII),
439 writeVersion_(IOstream::currentVersion),
440 writeCompression_(IOstream::UNCOMPRESSED),
442 runTimeModifiable_(true),
444 functionObjects_(*this, enableFunctionObjects)
446 libs_.open(controlDict_, "libs");
449 // Explicitly set read flags on objectRegistry so anything constructed
450 // from it reads as well (e.g. fvSolution).
451 readOpt() = IOobject::MUST_READ_IF_MODIFIED;
453 // Since could not construct regIOobject with setting:
454 controlDict_.readOpt() = IOobject::MUST_READ_IF_MODIFIED;
458 // Time objects not registered so do like objectRegistry::checkIn ourselves.
459 if (runTimeModifiable_)
465 regIOobject::fileModificationChecking == inotify
466 || regIOobject::fileModificationChecking == inotifyMaster
470 // File might not exist yet.
471 fileName f(controlDict_.filePath());
475 // We don't have this file but would like to re-read it.
476 // Possibly if in master-only reading mode. Use a non-existing
477 // file to keep fileMonitor synced.
478 f = controlDict_.objectPath();
481 controlDict_.watchIndex() = addWatch(f);
488 const fileName& rootPath,
489 const fileName& caseName,
490 const word& systemName,
491 const word& constantName,
492 const bool enableFunctionObjects
503 objectRegistry(*this),
525 writeControl_(wcTimeStep),
526 writeInterval_(GREAT),
530 writeFormat_(IOstream::ASCII),
531 writeVersion_(IOstream::currentVersion),
532 writeCompression_(IOstream::UNCOMPRESSED),
534 runTimeModifiable_(true),
536 functionObjects_(*this, enableFunctionObjects)
538 libs_.open(controlDict_, "libs");
542 // * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
546 if (controlDict_.watchIndex() != -1)
548 removeWatch(controlDict_.watchIndex());
551 // destroy function objects first
552 functionObjects_.clear();
556 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
558 Foam::label Foam::Time::addWatch(const fileName& fName) const
560 return monitorPtr_().addWatch(fName);
564 bool Foam::Time::removeWatch(const label watchIndex) const
566 return monitorPtr_().removeWatch(watchIndex);
569 const Foam::fileName& Foam::Time::getFile(const label watchIndex) const
571 return monitorPtr_().getFile(watchIndex);
575 Foam::fileMonitor::fileState Foam::Time::getState
580 return monitorPtr_().getState(watchFd);
584 void Foam::Time::setUnmodified(const label watchFd) const
586 monitorPtr_().setUnmodified(watchFd);
590 Foam::word Foam::Time::timeName(const scalar t)
592 std::ostringstream buf;
593 buf.setf(ios_base::fmtflags(format_), ios_base::floatfield);
594 buf.precision(precision_);
600 Foam::word Foam::Time::timeName() const
602 return dimensionedScalar::name();
606 // Search the construction path for times
607 Foam::instantList Foam::Time::times() const
609 return findTimes(path());
613 Foam::word Foam::Time::findInstancePath(const instant& t) const
615 instantList timeDirs = findTimes(path());
617 forAllReverse(timeDirs, timeI)
619 if (timeDirs[timeI] == t)
621 return timeDirs[timeI].name();
629 Foam::instant Foam::Time::findClosestTime(const scalar t) const
631 instantList timeDirs = findTimes(path());
633 // there is only one time (likely "constant") so return it
634 if (timeDirs.size() == 1)
639 if (t < timeDirs[1].value())
643 else if (t > timeDirs.last().value())
645 return timeDirs.last();
648 label nearestIndex = -1;
649 scalar deltaT = GREAT;
651 for (label timeI=1; timeI < timeDirs.size(); ++timeI)
653 scalar diff = mag(timeDirs[timeI].value() - t);
657 nearestIndex = timeI;
661 return timeDirs[nearestIndex];
665 // This should work too,
666 // if we don't worry about checking "constant" explicitly
668 // Foam::instant Foam::Time::findClosestTime(const scalar t) const
670 // instantList timeDirs = findTimes(path());
671 // label timeIndex = min(findClosestTimeIndex(timeDirs, t), 0);
672 // return timeDirs[timeIndex];
675 Foam::label Foam::Time::findClosestTimeIndex
677 const instantList& timeDirs,
681 label nearestIndex = -1;
682 scalar deltaT = GREAT;
684 forAll(timeDirs, timeI)
686 if (timeDirs[timeI].name() == "constant") continue;
688 scalar diff = mag(timeDirs[timeI].value() - t);
692 nearestIndex = timeI;
700 Foam::label Foam::Time::startTimeIndex() const
702 return startTimeIndex_;
706 Foam::dimensionedScalar Foam::Time::startTime() const
708 return dimensionedScalar("startTime", dimTime, startTime_);
712 Foam::dimensionedScalar Foam::Time::endTime() const
714 return dimensionedScalar("endTime", dimTime, endTime_);
718 bool Foam::Time::run() const
720 bool running = value() < (endTime_ - 0.5*deltaT_);
724 // only execute when the condition is no longer true
725 // ie, when exiting the control loop
726 if (!running && timeIndex_ != startTimeIndex_)
728 // Note, end() also calls an indirect start() as required
729 functionObjects_.end();
737 const_cast<Time&>(*this).readModifiedObjects();
739 if (timeIndex_ == startTimeIndex_)
741 functionObjects_.start();
745 functionObjects_.execute();
749 // Update the "running" status following the
750 // possible side-effects from functionObjects
751 running = value() < (endTime_ - 0.5*deltaT_);
758 bool Foam::Time::loop()
760 bool running = run();
771 bool Foam::Time::end() const
773 return value() > (endTime_ + 0.5*deltaT_);
777 bool Foam::Time::stopAt(const stopAtControls sa) const
779 const bool changed = (stopAt_ != sa);
785 controlDict_.lookup("endTime") >> endTime_;
795 void Foam::Time::setTime(const Time& t)
798 dimensionedScalar::name() = t.dimensionedScalar::name();
799 timeIndex_ = t.timeIndex_;
803 void Foam::Time::setTime(const instant& inst, const label newIndex)
805 value() = inst.value();
806 dimensionedScalar::name() = inst.name();
807 timeIndex_ = newIndex;
809 IOdictionary timeDict
817 IOobject::READ_IF_PRESENT,
823 timeDict.readIfPresent("deltaT", deltaT_);
824 timeDict.readIfPresent("deltaT0", deltaT0_);
825 timeDict.readIfPresent("index", timeIndex_);
829 void Foam::Time::setTime(const dimensionedScalar& newTime, const label newIndex)
831 setTime(newTime.value(), newIndex);
835 void Foam::Time::setTime(const scalar newTime, const label newIndex)
838 dimensionedScalar::name() = timeName(timeToUserTime(newTime));
839 timeIndex_ = newIndex;
843 void Foam::Time::setEndTime(const dimensionedScalar& endTime)
845 setEndTime(endTime.value());
849 void Foam::Time::setEndTime(const scalar endTime)
855 void Foam::Time::setDeltaT(const dimensionedScalar& deltaT)
857 setDeltaT(deltaT.value());
861 void Foam::Time::setDeltaT(const scalar deltaT)
864 deltaTchanged_ = true;
869 Foam::TimeState Foam::Time::subCycle(const label nSubCycles)
872 prevTimeState_.set(new TimeState(*this));
874 setTime(*this - deltaT(), (timeIndex() - 1)*nSubCycles);
875 deltaT_ /= nSubCycles;
876 deltaT0_ /= nSubCycles;
877 deltaTSave_ = deltaT0_;
879 return prevTimeState();
883 void Foam::Time::endSubCycle()
888 TimeState::operator=(prevTimeState());
889 prevTimeState_.clear();
894 // * * * * * * * * * * * * * * * Member Operators * * * * * * * * * * * * * //
896 Foam::Time& Foam::Time::operator+=(const dimensionedScalar& deltaT)
898 return operator+=(deltaT.value());
902 Foam::Time& Foam::Time::operator+=(const scalar deltaT)
909 Foam::Time& Foam::Time::operator++()
911 deltaT0_ = deltaTSave_;
912 deltaTSave_ = deltaT_;
914 // Save old time name
915 const word oldTimeName = dimensionedScalar::name();
917 setTime(value() + deltaT_, timeIndex_ + 1);
921 // If the time is very close to zero reset to zero
922 if (mag(value()) < 10*SMALL*deltaT_)
924 setTime(0.0, timeIndex_);
929 // Check that new time representation differs from old one
930 if (dimensionedScalar::name() == oldTimeName)
932 int oldPrecision = precision_;
936 setTime(value(), timeIndex());
938 while (precision_ < 100 && dimensionedScalar::name() == oldTimeName);
940 WarningIn("Time::operator++()")
941 << "Increased the timePrecision from " << oldPrecision
942 << " to " << precision_
943 << " to distinguish between timeNames at time " << value()
950 switch (writeControl_)
953 outputTime_ = !(timeIndex_ % label(writeInterval_));
957 case wcAdjustableRunTime:
959 label outputIndex = label
961 ((value() - startTime_) + 0.5*deltaT_)
965 if (outputIndex > outputTimeIndex_)
968 outputTimeIndex_ = outputIndex;
979 label outputIndex = label
981 returnReduce(elapsedCpuTime(), maxOp<double>())
984 if (outputIndex > outputTimeIndex_)
987 outputTimeIndex_ = outputIndex;
998 label outputIndex = label
1000 returnReduce(label(elapsedClockTime()), maxOp<label>())
1003 if (outputIndex > outputTimeIndex_)
1006 outputTimeIndex_ = outputIndex;
1010 outputTime_ = false;
1016 // see if endTime needs adjustment to stop at the next run()/end() check
1019 if (stopAt_ == saNoWriteNow)
1023 else if (stopAt_ == saWriteNow)
1028 else if (stopAt_ == saNextWrite && outputTime_ == true)
1039 Foam::Time& Foam::Time::operator++(int)
1041 return operator++();
1045 // ************************************************************************* //