timeControlFunctionObject.C
Go to the documentation of this file.
1 /*---------------------------------------------------------------------------*\
2  ========= |
3  \\ / F ield | OpenFOAM: The Open Source CFD Toolbox
4  \\ / O peration |
5  \\ / A nd | www.openfoam.com
6  \\/ M anipulation |
7 -------------------------------------------------------------------------------
8  Copyright (C) 2016 OpenFOAM Foundation
9  Copyright (C) 2017-2019 OpenCFD Ltd.
10 -------------------------------------------------------------------------------
11 License
12  This file is part of OpenFOAM.
13 
14  OpenFOAM is free software: you can redistribute it and/or modify it
15  under the terms of the GNU General Public License as published by
16  the Free Software Foundation, either version 3 of the License, or
17  (at your option) any later version.
18 
19  OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
20  ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
21  FITNESS FOR A PARTICULAR PURPOSE. See the GNU General Public License
22  for more details.
23 
24  You should have received a copy of the GNU General Public License
25  along with OpenFOAM. If not, see <http://www.gnu.org/licenses/>.
26 
27 \*---------------------------------------------------------------------------*/
28 
30 #include "polyMesh.H"
31 #include "mapPolyMesh.H"
33 
34 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
35 
36 namespace Foam
37 {
38 namespace functionObjects
39 {
40  defineTypeNameAndDebug(timeControl, 0);
41 }
42 }
43 
46 ({
47  { controlMode::TIME, "time" },
48  { controlMode::TRIGGER, "trigger" },
49  { controlMode::TIME_OR_TRIGGER, "timeOrTrigger" },
50  { controlMode::TIME_AND_TRIGGER, "timeAndTrigger" }
51 });
52 
53 
54 // * * * * * * * * * * * * * * * Private Members * * * * * * * * * * * * * * //
55 
56 void Foam::functionObjects::timeControl::readControls()
57 {
58  if (dict_.readIfPresent("timeStart", timeStart_))
59  {
60  timeStart_ = time_.userTimeToTime(timeStart_);
61  }
62  if (dict_.readIfPresent("timeEnd", timeEnd_))
63  {
64  timeEnd_ = time_.userTimeToTime(timeEnd_);
65  }
66 
67  if (dict_.readIfPresent("triggerStart", triggerStart_))
68  {
69  dict_.readIfPresent("triggerEnd", triggerEnd_);
70  controlMode_ = controlModeNames_.get("controlMode", dict_);
71  }
72 
73  deltaTCoeff_ = GREAT;
74  if (dict_.readIfPresent("deltaTCoeff", deltaTCoeff_))
75  {
76  nStepsToStartTimeChange_ = labelMax;
77  }
78  else
79  {
80  nStepsToStartTimeChange_ = 3;
81  dict_.readIfPresent
82  (
83  "nStepsToStartTimeChange",
84  nStepsToStartTimeChange_
85  );
86  }
87 }
88 
89 
90 bool Foam::functionObjects::timeControl::active() const
91 {
92  label triggeri = time_.functionObjects().triggerIndex();
93 
94  bool inTime =
95  time_.value() >= (timeStart_ - 0.5*time_.deltaTValue())
96  && time_.value() <= (timeEnd_ + 0.5*time_.deltaTValue());
97 
98  bool inTrigger = triggeri >= triggerStart_ && triggeri <= triggerEnd_;
99 
100  switch (controlMode_)
101  {
102  case controlMode::TIME:
103  {
104  return inTime;
105  }
106  case controlMode::TRIGGER:
107  {
108  return inTrigger;
109  }
110  case controlMode::TIME_OR_TRIGGER:
111  {
112  return inTime || inTrigger;
113  }
114  case controlMode::TIME_AND_TRIGGER:
115  {
116  return inTime && inTrigger;
117  }
118  default:
119  {
121  << "Unhandled enumeration: " << controlModeNames_[controlMode_]
122  << abort(FatalError);
123  }
124  }
125 
126  return false;
127 }
128 
129 
130 Foam::scalar Foam::functionObjects::timeControl::calcExpansion
131 (
132  const scalar startRatio,
133  const scalar y,
134  const label n
135 )
136 {
137  scalar ratio = startRatio;
138  // This function is used to calculate the 'expansion ratio' or
139  // 'time step growth factor'. Inputs:
140  // - ratio: wanted ratio
141  // - y: overall time required (divided by deltaT) to get from 1 to ratio
142  // - n: number of steps in which to get to wanted ratio
143 
144  // Let a0 = first term in geometric progression (GP), an = last term,
145  // n = number of terms, ratio = ratio of GP, Sn = sum of series for
146  // first n terms.
147  // We have an=a0*ratio^(n-1) Sn = a0(ratio^n -1)/(ratio -1)
148  // -> Sn=an/ratio^(n-1)*(ratio^n -1)/(ratio -1), assume y=Sn/an
149  // -> y*ratio^(n-1)=(ratio^n -1)/(ratio -1)
150  // -> y*ratio^n - y*ratio^(n-1) -ratio^n +1 = 0
151 
152  // Solve using Newton-Raphson iteration to determine the expansion
153  // ratio giving the desired overall timestep change.
154  // Note max iteration count to avoid hanging; function generally
155  // converges in 5 iterations or so.
156  for (label iter = 0; iter < 100; iter++)
157  {
158  // Dimensionless equation
159  scalar f = (y-1)*pow(ratio, n)+1-y*pow(ratio, n-1);
160  scalar dfdratio = (y-1)*n*pow(ratio, n-1)-y*(n-1)*pow(ratio, n-2);
161  scalar dratio = f/(dfdratio + SMALL);
162  ratio -= dratio;
163  // Exit if function is satisfied
164  if (mag(f) < 1e-6)
165  {
166  return ratio;
167  }
168  }
169 
170  if (debug)
171  {
173  << "Did not converge to find new timestep growth factor given "
174  << "overall factor " << y << " and " << n << " steps to do it in."
175  << endl << " Returning current best guess " << ratio << endl;
176  }
177  return ratio;
178 }
179 
180 
181 void Foam::functionObjects::timeControl::calcDeltaTCoeff
182 (
183  scalar& requiredDeltaTCoeff,
184  const scalar wantedDT,
185  const label nSteps,
186  const scalar presentTime,
187  const scalar timeToNextWrite,
188  const bool rampDirectionUp
189 )
190 {
191  const scalar writeInterval = writeControl_.interval();
192 
193 
194  // Set new deltaT approximated to last step value
195  // adjust to include geometric series sum of nSteps terms
196  scalar newDeltaT = deltaT0_;
197 
198  if (seriesDTCoeff_ != GREAT)
199  {
200  newDeltaT *= seriesDTCoeff_;
201  }
202 
203  // Recalculate new deltaTCoeff_ based on rounded steps
204  requiredDeltaTCoeff = Foam::exp(Foam::log(wantedDT/newDeltaT)/nSteps);
205 
206  // Calculate total time required with given dT increment
207  // to reach specified dT value
208  scalar requiredTimeInterval = newDeltaT;
209 
210  // For nSteps = 1 during we get coeff = 1.0. Adding SMALL in denominator
211  // might lead to SMALL requiredTimeInterval, and further unnecessary
212  // calculations
213  if (requiredDeltaTCoeff != 1.0)
214  {
215  requiredTimeInterval *=
216  (pow(requiredDeltaTCoeff, nSteps) - 1)
217  /(requiredDeltaTCoeff - 1);
218  }
219 
220  // Nearest time which is multiple of wantedDT, after
221  // requiredTimeInterval
222  scalar timeToNextMultiple = -presentTime;
223 
224  if (rampDirectionUp)
225  {
226  timeToNextMultiple +=
227  ceil
228  (
229  (presentTime + requiredTimeInterval)
230  /wantedDT
231  )
232  *wantedDT;
233  }
234  else
235  {
236  timeToNextMultiple +=
237  floor
238  (
239  (presentTime - requiredTimeInterval)
240  /wantedDT
241  )
242  *wantedDT;
243  }
244 
245  // Check nearest time and decide parameters
246  if (timeToNextWrite <= timeToNextMultiple)
247  {
248  // As the ramping can't be achieved, use current, unadapted deltaT
249  // (from e.g. maxCo calculation in solver) and adjust when
250  // writeInterval is closer than this deltaT
251  newDeltaT = deltaT0_;
252  if (timeToNextWrite < newDeltaT)
253  {
254  newDeltaT = timeToNextWrite;
255  }
256 
257  // Corner case when required time span is less than writeInterval
258  if (requiredTimeInterval > writeInterval)
259  {
261  << "With given ratio needed time span "
262  << requiredTimeInterval
263  << " exceeds available writeInterval "
264  << writeInterval << nl
265  << "Disabling all future time step ramping"
266  << endl;
267  deltaTCoeff_ = GREAT;
268  newDeltaT = wantedDT;
269  }
270  else
271  {
272  // Reset the series ratio, as this could change later
273  seriesDTCoeff_ = GREAT;
274  if (debug)
275  {
276  Info<< "Disabling ramping until time "
277  << presentTime+timeToNextWrite << endl;
278  }
279  }
280  requiredDeltaTCoeff = newDeltaT/deltaT0_;
281  }
282  else
283  {
284  // Find the minimum number of time steps (with constant deltaT
285  // variation) to get us to the writeInterval and multiples of wantedDT.
286  // Increases the number of nSteps to do the adjustment in until it
287  // has reached one that is possible. Returns the new expansionratio.
288 
289  scalar y = timeToNextMultiple/wantedDT;
290  label requiredSteps = nSteps;
291 
292  scalar ratioEstimate = deltaTCoeff_;
293  scalar ratioMax = deltaTCoeff_;
294 
295  if (seriesDTCoeff_ != GREAT)
296  {
297  ratioEstimate = seriesDTCoeff_;
298  }
299 
300  if (!rampDirectionUp)
301  {
302  ratioEstimate = 1/ratioEstimate;
303  ratioMax = 1/ratioMax;
304  requiredSteps *= -1;
305  }
306  // Provision for fall-back value if we can't satisfy requirements
307  bool searchConverged = false;
308  for (label iter = 0; iter < 100; iter++)
309  {
310  // Calculate the new expansion and the overall time step changes
311  const scalar newRatio = calcExpansion
312  (
313  ratioEstimate,
314  y,
315  requiredSteps
316  );
317  const scalar deltaTLoop =
318  wantedDT
319  / pow(newRatio, mag(requiredSteps)-1);
320  scalar firstDeltaRatio = deltaTLoop/deltaT0_;
321  // Avoid division by zero for ratio = 1.0
322  scalar Sn =
323  deltaTLoop
324  *(pow(newRatio, mag(requiredSteps))-1)
325  /(newRatio-1+SMALL);
326 
327  if (debug)
328  {
329  Info<< " nSteps " << requiredSteps
330  << " ratioEstimate " << ratioEstimate
331  << " a0 " << deltaTLoop
332  << " firstDeltaRatio " << firstDeltaRatio
333  << " Sn " << Sn << " Sn+Time " << Sn+presentTime
334  << " seriesRatio "<< newRatio << endl;
335  }
336 
337  // If series ratio becomes unity check next deltaT multiple
338  // Do the same if first ratio is decreasing when ramping up!
339  if
340  (
341  (firstDeltaRatio < 1.0 && rampDirectionUp)
342  || (firstDeltaRatio > 1.0 && !rampDirectionUp)
343  || newRatio == 1.0
344  )
345  {
346  y += 1;
347  requiredSteps = mag(nSteps);
348  if (debug)
349  {
350  Info<< "firstDeltaRatio " << firstDeltaRatio << " rampDir"
351  << rampDirectionUp << " newRatio " << newRatio
352  << " y " << y << " steps " << requiredSteps <<endl;
353  }
354  continue;
355  }
356 
357  if (firstDeltaRatio > ratioMax && rampDirectionUp)
358  {
359  requiredSteps++;
360  if (debug)
361  {
362  Info<< "First ratio " << firstDeltaRatio
363  << " exceeds threshold " << ratioMax << nl
364  << "Increasing required steps " << requiredSteps
365  << endl;
366  }
367  }
368  else if (firstDeltaRatio < ratioMax && !rampDirectionUp)
369  {
370  requiredSteps--;
371  if (debug)
372  {
373  Info<< "First ratio " << firstDeltaRatio
374  << " exceeds threshold " << ratioMax << nl
375  << "Decreasing required steps " << requiredSteps
376  << endl;
377  }
378  }
379  else if
380  (
381  (newRatio > ratioMax && rampDirectionUp)
382  || (newRatio < ratioMax && !rampDirectionUp)
383  )
384  {
385  y += 1;
386  requiredSteps = nSteps;
387  if (debug)
388  {
389  Info<< "Series ratio " << newRatio << "exceeds threshold "
390  << ratioMax << nl << "Consider next deltaT multiple "
391  << y*wantedDT + presentTime << endl;
392  }
393  }
394  else
395  {
396  seriesDTCoeff_ = newRatio;
397  // Reset future series expansion at last step
398  if (requiredSteps == 1)
399  {
400  Sn = deltaT0_*firstDeltaRatio;
401  seriesDTCoeff_ = GREAT;
402  }
403  // Situations where we achieve wantedDT but fail to achieve
404  // multiple of writeInterval
405  scalar jumpError =
406  remainder(Sn + presentTime, wantedDT) / wantedDT;
407 
408  if (mag(jumpError) > ROOTSMALL)
409  {
410  requiredSteps = label(timeToNextWrite/wantedDT);
411  firstDeltaRatio = timeToNextWrite/requiredSteps/deltaT0_;
412  }
413  if (debug)
414  {
415  Info<< "All conditions satisfied deltaT0_:" << deltaT0_
416  << " calculated deltaTCoeff:" << firstDeltaRatio
417  << " series ratio set to:" << seriesDTCoeff_ << endl;
418  }
419  searchConverged = true;
420  requiredDeltaTCoeff = firstDeltaRatio;
421  break;
422  }
423  }
424 
425  if (!searchConverged)
426  {
427  if (debug)
428  {
429  // Not converged. Return fall-back value
431  << "Did not converge to find new timestep growth factor"
432  << " given overall factor " << y
433  << " and " << requiredSteps << " steps to do it in."
434  << nl << " Falling back to non-adjusted deltatT "
435  << deltaT0_ << endl;
436  }
437  requiredDeltaTCoeff = 1.0;
438  }
439  }
440 }
441 
442 
443 // * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
444 
445 Foam::functionObjects::timeControl::timeControl
446 (
447  const word& name,
448  const Time& runTime,
449  const dictionary& dict
450 )
451 :
453  dict_(dict),
454  controlMode_(controlMode::TIME),
455  timeStart_(-VGREAT),
456  timeEnd_(VGREAT),
457  triggerStart_(labelMax),
458  triggerEnd_(labelMax),
459  nStepsToStartTimeChange_(labelMax),
460  executeControl_(runTime, dict, "execute"),
461  writeControl_(runTime, dict, "write"),
462  foPtr_(functionObject::New(name, runTime, dict_)),
463  executeTimeIndex_(-1),
464  deltaT0_(0),
465  seriesDTCoeff_(GREAT)
466 {
467  readControls();
468 }
469 
470 
471 // * * * * * * * * * * * * * * * Member Functions * * * * * * * * * * * * * //
472 
474 {
475  if
476  (
478  || Foam::timeControl::entriesPresent(dict, "output") // backwards compat
480  || dict.found("timeStart")
481  || dict.found("timeEnd")
482  || dict.found("triggerStart")
483  || dict.found("triggerEnd")
484  )
485  {
486  return true;
487  }
488 
489  return false;
490 }
491 
492 
494 {
495  // Store old deltaT for reference
496  deltaT0_ = time_.deltaTValue();
497 
498  if (active() && (postProcess || executeControl_.execute()))
499  {
500  executeTimeIndex_ = time_.timeIndex();
501  foPtr_->execute();
502  }
503 
504  return true;
505 }
506 
507 
509 {
510  if (active())
511  {
512  // Call underlying function object directly
513  foPtr_->execute(subIndex);
514  }
515 
516  return true;
517 }
518 
519 
521 {
522  if (active() && (postProcess || writeControl_.execute()))
523  {
524  // Ensure written results reflect the current state
525  if (executeTimeIndex_ != time_.timeIndex())
526  {
527  executeTimeIndex_ = time_.timeIndex();
528  foPtr_->execute();
529  }
530 
531  foPtr_->write();
532  }
533 
534  return true;
535 }
536 
537 
539 {
540  if (active() && (executeControl_.execute() || writeControl_.execute()))
541  {
542  foPtr_->end();
543  }
544 
545  return true;
546 }
547 
548 
550 {
551  if (active())
552  {
553  if
554  (
555  writeControl_.control()
557  )
558  {
559  const scalar presentTime = time_.value();
560 
561  //const scalar oldDeltaT = time_.deltaTValue();
562 
563  // Call underlying fo to do time step adjusting
564  foPtr_->adjustTimeStep();
565  const scalar wantedDT = time_.deltaTValue();
566 
567  //Pout<< "adjustTimeStep by " << foPtr_->name()
568  // << " from " << oldDeltaT << " to " << wantedDT << endl;
569 
570  const label writeTimeIndex = writeControl_.executionIndex();
571  const scalar writeInterval = writeControl_.interval();
572 
573  // Get time to next writeInterval
574  scalar timeToNextWrite =
575  (writeTimeIndex+1)*writeInterval
576  - (presentTime - time_.startTime().value());
577  if (timeToNextWrite <= 0.0)
578  {
579  timeToNextWrite = writeInterval;
580  }
581 
582  // The target DT (coming from fixed dT or maxCo etc) may not be
583  // an integral multiple of writeInterval. If so ignore any step
584  // calculation and give priority to writeInterval.
585  // Note looser tolerance on the relative intervalError - SMALL is
586  // too strict.
587  scalar intervalError =
588  remainder(writeInterval, wantedDT)/writeInterval;
589  if
590  (
591  mag(intervalError) > ROOTSMALL
592  || deltaTCoeff_ == GREAT
593  )
594  {
595  scalar deltaT = time_.deltaTValue();
596  scalar nSteps = timeToNextWrite/deltaT;
597 
598  // For tiny deltaT the label can overflow!
599  if
600  (
601  nSteps < scalar(labelMax)
602  && (
603  deltaTCoeff_ != GREAT
604  || nSteps < nStepsToStartTimeChange_
605  )
606  )
607  {
608  // nSteps can be < 1 so make sure at least 1
609 
610  label nStepsToNextWrite = max(1, round(nSteps));
611  scalar newDeltaT = timeToNextWrite/nStepsToNextWrite;
612 
613  // Backwards compatibility: clip deltaT to 0.2 .. 2
614  scalar clipThreshold = 2;
615  if (deltaTCoeff_ != GREAT)
616  {
617  clipThreshold = deltaTCoeff_;
618  }
619  // Adjust time step
620  if (newDeltaT >= deltaT)
621  {
622  deltaT = min(newDeltaT, clipThreshold*deltaT);
623  }
624  else
625  {
626  clipThreshold = 1/clipThreshold;
627  deltaT = max(newDeltaT, clipThreshold*deltaT);
628  }
629 
630  const_cast<Time&>(time_).setDeltaT(deltaT, false);
631  }
632  }
633  else
634  {
635  // Set initial approximation to user defined ratio.
636  scalar requiredDeltaTCoeff = deltaTCoeff_;
637 
638  // Use if we have series expansion ratio from last attempt.
639  // Starting from user defined coefficient, might get different
640  // steps and convergence time (could be recursive - worst case)
641  // This is more likely to happen when coefficient limit is high
642  if (seriesDTCoeff_ != GREAT)
643  {
644  requiredDeltaTCoeff = seriesDTCoeff_;
645  }
646  // Avoid divide by zero if we need ratio = 1
647  if (1/Foam::log(requiredDeltaTCoeff) > scalar(labelMax))
648  {
649  requiredDeltaTCoeff = deltaTCoeff_;
650  }
651 
652  // Try and observe the deltaTCoeff and the writeInterval
653  // and the wanted deltaT. Below code calculates the
654  // minimum number of time steps in which to end up on
655  // writeInterval+n*deltaT with the minimum possible deltaTCoeff.
656  // This is non-trivial!
657 
658  // Approx steps required from present dT to required dT
659  label nSteps = 0;
660  if (wantedDT < deltaT0_)
661  {
662  nSteps = label
663  (
664  floor
665  (
666  Foam::log(wantedDT/deltaT0_)
667  /Foam::log(requiredDeltaTCoeff + SMALL)
668  )
669  );
670  }
671  else
672  {
673  nSteps = label
674  (
675  ceil
676  (
677  Foam::log(wantedDT/deltaT0_)
678  /Foam::log(requiredDeltaTCoeff + SMALL)
679  )
680  );
681  }
682  // Negative steps -> ramp down needed,
683  // zero steps -> we are at writeInterval-multiple time!
684  if (nSteps < 1)
685  {
686  // Not enough steps to do gradual time step adjustment
687  // if earlier deltaT larger than wantedDT
688  if (wantedDT < deltaT0_)
689  {
690  requiredDeltaTCoeff = 1/requiredDeltaTCoeff;
691  calcDeltaTCoeff
692  (
693  requiredDeltaTCoeff,
694  wantedDT,
695  nSteps,
696  presentTime,
697  timeToNextWrite,
698  false
699  );
700  }
701  else
702  {
703  if (timeToNextWrite > wantedDT)
704  {
705  requiredDeltaTCoeff = wantedDT/deltaT0_;
706  }
707  else
708  {
709  requiredDeltaTCoeff = timeToNextWrite/deltaT0_;
710  }
711  }
712  // When allowed deltaTCoeff can't give required ramp
713  if (requiredDeltaTCoeff > deltaTCoeff_ && debug)
714  {
716  << "Required deltaTCoeff "<< requiredDeltaTCoeff
717  << " is larger than allowed value "<< deltaTCoeff_
718  << endl;
719  }
720  }
721  else
722  {
723  calcDeltaTCoeff
724  (
725  requiredDeltaTCoeff,
726  wantedDT,
727  nSteps,
728  presentTime,
729  timeToNextWrite,
730  true
731  );
732  }
733 
734  // Calculate deltaT to be used based on ramp
735  scalar newDeltaT = deltaT0_*requiredDeltaTCoeff;
736 
737  // Set time, deltaT adjustments for writeInterval purposes
738  // are already taken care. Hence disable auto-update
739  const_cast<Time&>(time_).setDeltaT(newDeltaT, false);
740 
741  if (seriesDTCoeff_ < 1.0)
742  {
743  requiredDeltaTCoeff = 1/requiredDeltaTCoeff;
744  seriesDTCoeff_ = 1/seriesDTCoeff_;
745  }
746  }
747  }
748  else
749  {
750  foPtr_->adjustTimeStep();
751  const scalar wantedDT = time_.deltaTValue();
752 
753  if (deltaTCoeff_ != GREAT)
754  {
755  // Clip time step change to deltaTCoeff
756 
757  scalar requiredDeltaTCoeff =
758  (
759  max
760  (
761  1.0/deltaTCoeff_,
762  min(deltaTCoeff_, wantedDT/deltaT0_)
763  )
764  );
765 
766  scalar newDeltaT = deltaT0_*requiredDeltaTCoeff;
767 
768  // Set time, deltaT adjustments for writeInterval purposes
769  // are already taken care. Hence disable auto-update
770  const_cast<Time&>( time_).setDeltaT(newDeltaT, false);
771  }
772  else
773  {
774  // Set the DT without any ramping
775  const_cast<Time&>( time_).setDeltaT(wantedDT, false);
776  }
777  }
778  }
779 
780  return true;
781 }
782 
783 
785 {
786  if (dict != dict_)
787  {
788  dict_ = dict;
789 
790  writeControl_.read(dict);
791  executeControl_.read(dict);
792  readControls();
793 
794  // Forward to underlying function object
795  return foPtr_->read(dict);
796  }
797 
798  return false;
799 }
800 
801 
803 {
804  return active() ? foPtr_->filesModified() : false;
805 }
806 
807 
809 {
810  if (active())
811  {
812  foPtr_->updateMesh(mpm);
813  }
814 }
815 
816 
818 {
819  if (active())
820  {
821  foPtr_->movePoints(mesh);
822  }
823 }
824 
825 
826 // ************************************************************************* //
Foam::expressions::patchExpr::debug
int debug
Static debugging option.
Foam::functionObjects::timeControl::movePoints
virtual void movePoints(const polyMesh &mesh)
Update for changes of mesh.
Definition: timeControlFunctionObject.C:817
runTime
engineTime & runTime
Definition: createEngineTime.H:13
Foam::Enum
Enum is a wrapper around a list of names/values that represent particular enumeration (or int) values...
Definition: IOstreamOption.H:57
Foam::Time
Class to control time during OpenFOAM simulations that is also the top-level objectRegistry.
Definition: Time.H:73
Foam::labelMax
constexpr label labelMax
Definition: label.H:61
Foam::word
A class for handling words, derived from Foam::string.
Definition: word.H:65
Foam::functionObject::New
static autoPtr< functionObject > New(const word &name, const Time &runTime, const dictionary &dict)
Select from dictionary, based on its "type" entry.
Definition: functionObject.C:79
Foam::functionObjects::timeFunctionObject::time_
const Time & time_
Reference to the time database.
Definition: timeFunctionObject.H:65
mapPolyMesh.H
Foam::functionObjects::timeControl::end
virtual bool end()
Called when Time::run() determines that the time-loop exits.
Definition: timeControlFunctionObject.C:538
Foam::TimeState::userTimeToTime
virtual scalar userTimeToTime(const scalar theta) const
Convert the user-time (e.g. CA deg) to real-time (s).
Definition: TimeState.C:49
Foam::functionObjects::timeFunctionObject
Virtual base class for function objects with a reference to Time.
Definition: timeFunctionObject.H:56
Foam::timeControl::entriesPresent
static bool entriesPresent(const dictionary &dict, const word &prefix)
Identify if a timeControl object is present in the dictionary.
Definition: timeControl.C:87
Foam::endl
Ostream & endl(Ostream &os)
Add newline and flush stream.
Definition: Ostream.H:369
Foam::exp
dimensionedScalar exp(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:261
Foam::Enum::get
EnumType get(const word &enumName) const
The enumeration corresponding to the given name.
Definition: Enum.C:75
polyMesh.H
Foam::min
label min(const labelHashSet &set, label minValue=labelMax)
Find the min value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:33
Foam::polyMesh
Mesh consisting of general polyhedral cells.
Definition: polyMesh.H:77
Foam::TimeState::deltaTValue
scalar deltaTValue() const noexcept
Return time step value.
Definition: TimeStateI.H:43
Foam::functionObjects::timeControl::read
virtual bool read(const dictionary &)
Read and set the function object if its data have changed.
Definition: timeControlFunctionObject.C:784
n
label n
Definition: TABSMDCalcMethod2.H:31
Foam::functionObjects::timeControl::controlModeNames_
static const Enum< controlMode > controlModeNames_
Definition: timeControlFunctionObject.H:93
Foam::Info
messageStream Info
Information stream (stdout output on master, null elsewhere)
Foam::functionObjects::timeControl::adjustTimeStep
virtual bool adjustTimeStep()
Called at the end of Time::adjustDeltaT() if adjustTime is true.
Definition: timeControlFunctionObject.C:549
Foam::pow
dimensionedScalar pow(const dimensionedScalar &ds, const dimensionedScalar &expt)
Definition: dimensionedScalar.C:75
Foam::max
label max(const labelHashSet &set, label maxValue=labelMin)
Find the max value in labelHashSet, optionally limited by second argument.
Definition: hashSets.C:47
dict
dictionary dict
Definition: searchingEngine.H:14
Foam::FatalError
error FatalError
Foam::dictionary
A list of keyword definitions, which are a keyword followed by a number of values (eg,...
Definition: dictionary.H:123
Foam::log
dimensionedScalar log(const dimensionedScalar &ds)
Definition: dimensionedScalar.C:262
mesh
dynamicFvMesh & mesh
Definition: createDynamicFvMesh.H:6
addToRunTimeSelectionTable.H
Macros for easy insertion into run-time selection tables.
Foam::functionObjects::timeControl::filesModified
virtual bool filesModified() const
Did any file get changed during execution?
Definition: timeControlFunctionObject.C:802
Foam
Namespace for OpenFOAM.
Definition: atmBoundaryLayer.C:33
Foam::abort
errorManip< error > abort(error &err)
Definition: errorManip.H:144
Foam::functionObjects::timeControl::write
virtual bool write()
Called at each ++ or += of the time-loop.
Definition: timeControlFunctionObject.C:520
Foam::timeControl::ocAdjustableRunTime
Currently identical to "runTime".
Definition: timeControl.H:73
FatalErrorInFunction
#define FatalErrorInFunction
Report an error message using Foam::FatalError.
Definition: error.H:453
Foam::nl
constexpr char nl
Definition: Ostream.H:404
Foam::Sn
static scalar Sn(const scalar a, const scalar x)
Definition: invIncGamma.C:68
f
labelList f(nPoints)
Foam::functionObjects::timeControl::updateMesh
virtual void updateMesh(const mapPolyMesh &mpm)
Update for changes of mesh.
Definition: timeControlFunctionObject.C:808
Foam::mag
dimensioned< typename typeOfMag< Type >::type > mag(const dimensioned< Type > &dt)
Foam::constant::electromagnetic::e
const dimensionedScalar e
Elementary charge.
Definition: createFields.H:11
timeControlFunctionObject.H
Foam::functionObjects::defineTypeNameAndDebug
defineTypeNameAndDebug(ObukhovLength, 0)
Foam::mapPolyMesh
Class containing mesh-to-mesh mapping information after a change in polyMesh topology.
Definition: mapPolyMesh.H:161
Foam::name
word name(const expressions::valueTypeCode typeCode)
A word representation of a valueTypeCode. Empty for INVALID.
Definition: exprTraits.C:59
Foam::functionObjects::timeControl::entriesPresent
static bool entriesPresent(const dictionary &dict)
Helper function to identify if a timeControl object is present.
Definition: timeControlFunctionObject.C:473
Foam::functionObjects::timeControl::execute
virtual bool execute()
Called at each ++ or += of the time-loop.
Definition: timeControlFunctionObject.C:493
WarningInFunction
#define WarningInFunction
Report a warning using Foam::Warning.
Definition: messageStream.H:328
Foam::dictionary::readIfPresent
bool readIfPresent(const word &keyword, T &val, enum keyType::option matchOpt=keyType::REGEX) const
Definition: dictionaryTemplates.C:405
y
scalar y
Definition: LISASMDCalcMethod1.H:14