Merge pull request #48 from mattfalcone1997/cfl_time_stepper

Add time stepper that allows OpenFOAM to control the time step so that CFL adaptive stepping can be used

Author: mattfalcone1997

include/base/FoamSolver.h

@@ -1,9 +1,38 @@
 #pragma once
 
-#include "fvCFD_moose.h"
+#include "solver.H"
+#include "functionObject.H"
 
 #include <vector>
 
+namespace Foam
+{
+namespace functionObjects
+{
+// Function object to tell OpenFOAM what MOOSE's dt is
+class mooseDeltaT : public functionObject
+{
+private:
+  const scalar & dt_;
+
+public:
+  TypeName("mooseDeltaT");
+  mooseDeltaT(const word & name, const Time & runTime, const scalar & dt)
+    : functionObject(name, runTime), dt_(dt)
+  {
+  }
+
+  wordList fields() const { return wordList::null(); }
+
+  bool executeAtStart() const { return false; }
+
+  bool execute() { return true; }
+  bool write() { return true; }
+  scalar maxDeltaT() const { return dt_; }
+};
+}
+}
+
 namespace Hippo
 {
 
@@ -32,6 +61,19 @@ class FoamSolver
   void setCurrentTime(double time) { runTime().setTime(time, runTime().timeIndex()); }
   // Set the time at which the solver should terminate.
   void setEndTime(double time) { runTime().setEndTime(time); }
+  // Run the presolve from MOOSE objects.
+  void preSolve();
+  // Provide access to the openfoam solver.
+  Foam::solver & solver() { return *_solver; };
+  // Calculate OpenFOAM's time step.
+  Foam::scalar computeDeltaT();
+  // check whether OpenFOAM has variable time step.
+  bool isDeltaTAdjustable() const;
+  // creates function object that tells OpenFOAM what MOOSE's
+  // time step is.
+  void appendDeltaTFunctionObject(const Foam::scalar & dt);
+  // get the current deltaT.
+  Foam::scalar getTimeDelta() { return runTime().deltaTValue(); };
 
 private:
   Foam::solver * _solver = nullptr;

include/solvers/timesteppers/FoamControlledTimeStepper.h

@@ -0,0 +1,42 @@
+#pragma once
+
+#include "FoamProblem.h"
+#include "FoamSolver.h"
+
+#include <InputParameters.h>
+#include <TimeStepper.h>
+
+/*
+Time stepper that allows OpenFOAM to control the time step enabling features such as CFL
+ daptive time steps. The intention is to allows the current time step in OpenFOAM
+ to be exposed to MOOSE
+ */
+
+class FoamControlledTimeStepper : public TimeStepper
+{
+public:
+  FoamControlledTimeStepper(InputParameters const & params);
+  static InputParameters validParams();
+
+  // Get initial time step from OpenFOAM input file
+  virtual Real computeInitialDT() { return computeDT(); };
+
+  /* Read time step from OpenFOAM
+    - Make sure the time step duration is computed in the current step
+  */
+  virtual Real computeDT();
+
+  // Consider how to communicate starting time for this stepper
+  // e.g. after a restart this would need to be executed after the
+  // OpenFOAM restart.
+  virtual void init();
+
+private:
+  // These two variables are needed depending on how the time-stepper is initialised
+  Hippo::FoamSolver & solver() { return problem()->solver(); }
+  FoamProblem * problem();
+  // Variables to determine whether an adjustable time step is used in OF and
+  // what it is.
+  bool _dt_adjustable = false;
+  Real _foam_initial_dt = 0.;
+};

src/base/FoamSolver.C

@@ -11,6 +11,15 @@
 #include <cassert>
 #include <iterator>
 
+namespace Foam
+{
+namespace functionObjects
+{
+defineTypeNameAndDebug(mooseDeltaT, 0);
+
+}
+}
+
 namespace Hippo
 {
 namespace
@@ -91,7 +100,6 @@ FoamSolver::run()
 
   // Adjust the time-step according to the solver maxDeltaT
   adjustDeltaT(time, solver);
-
   time++;
 
   // TODO: replace std::cout with MOOSE output or a dependency-injected stream.
@@ -206,4 +214,36 @@ FoamSolver::wallHeatFlux(const int patch_id, std::vector<double> & fill_vector)
   return whf_boundary.size();
 }
 
+void
+FoamSolver::preSolve()
+{
+  _solver->pimple.readIfModified();
+  _solver->preSolve();
+}
+
+Foam::scalar
+FoamSolver::computeDeltaT()
+{
+  // This code has been adapted from OpenFOAM's adjustDeltaT to determine the time-step that
+  // OpenFOAM will use on the next time step so MOOSE can predict it.
+  Foam::scalar deltaT =
+      std::min(_solver->maxDeltaT(), _solver->runTime.functionObjects().maxDeltaT());
+
+  if (deltaT < Foam::rootVGreat)
+    return std::min(Foam::solver::deltaTFactor * _solver->runTime.deltaTValue(), deltaT);
+  return _solver->runTime.deltaTValue();
+}
+
+bool
+FoamSolver::isDeltaTAdjustable() const
+{
+  return _solver->runTime.controlDict().lookupOrDefault("adjustTimeStep", false);
+}
+
+void
+FoamSolver::appendDeltaTFunctionObject(const Foam::scalar & dt)
+{
+  runTime().functionObjects().append(
+      new Foam::functionObjects::mooseDeltaT("Moose time step", runTime(), dt));
+}
 } // namespace Hippo

src/problems/FoamProblem.C

@@ -112,6 +112,7 @@ FoamProblem::externalSolve()
 {
   if (parameters().get<bool>("solve"))
   {
+    _solver.setTimeDelta(_dt); // Needed for constant deltaT cases
     _solver.run();
   }
 }

src/timesteppers/FoamControlledTimeStepper.C

@@ -0,0 +1,82 @@
+#include "FoamProblem.h"
+#include "FoamControlledTimeStepper.h"
+#include <pimpleSingleRegionControl.H>
+
+#include <scalar.H>
+#include <TimeStepper.h>
+#include <Transient.h>
+
+registerMooseObject("hippoApp", FoamControlledTimeStepper);
+
+InputParameters
+FoamControlledTimeStepper::validParams()
+{
+  auto params = TimeStepper::validParams();
+  return params;
+}
+
+FoamControlledTimeStepper::FoamControlledTimeStepper(InputParameters const & params)
+  : TimeStepper(params)
+{
+  auto problem = dynamic_cast<FoamProblem *>(&_app.feProblem());
+  if (!problem)
+  {
+    mooseError("FoamTimeStepper expects to be used with FoamProblem");
+  }
+}
+
+Real
+FoamControlledTimeStepper::computeDT()
+{
+
+  if (!_dt_adjustable)
+    return _foam_initial_dt;
+
+  // Not ideal, but for MOOSE to get an accurate deltaT
+  // preSolve must be called as this updates the BCs.
+  solver().preSolve();
+
+  // Ensure MOOSE gets OpenFOAM's time step unaffected by the mooseDeltaT
+  // functionObject.
+  Real dt_tmp = _dt;
+  _dt = Foam::rootVGreat;
+
+  // compute OpenFOAM's desired time step
+  Real deltaT = solver().computeDeltaT();
+
+  // reset MOOSE's time step and return
+  _dt = dt_tmp;
+  return deltaT;
+}
+
+FoamProblem *
+FoamControlledTimeStepper::problem()
+{
+  auto problem = dynamic_cast<FoamProblem *>(&_app.feProblem());
+  if (!problem)
+  {
+    mooseError("FoamTimeStepper expects to be used with FoamProblem");
+  }
+  return problem;
+}
+
+void
+FoamControlledTimeStepper::init()
+{
+  TimeStepper::init();
+
+  // determine if OpenFOAM's time-step is adjustable in controlDict
+  _dt_adjustable = solver().isDeltaTAdjustable();
+
+  // The key idea is that runTime.functionObjects().maxDeltaT() in adjustDeltaT
+  // loops over the function objects and chooses the minimum, so by having
+  // a function Object that returns what MOOSE wants, OpenFOAM will use the
+  // MOOSE time step if it is smaller than what OpenFOAM wants. As a result,
+  // if MOOSE wants to add a synchronisation step OpenFOAM will also use it too.
+
+  // create function object and append it to the solver's function object list
+  solver().appendDeltaTFunctionObject(_dt);
+
+  if (!_dt_adjustable)
+    _foam_initial_dt = solver().getTimeDelta();
+}

test/tests/timesteppers/foam_controlled_tstep_cfl_insert/fluid-openfoam/0/T

@@ -0,0 +1,56 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  10
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       volScalarField;
+    location    "0";
+    object      T;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 1 0 0 0];
+
+internalField   uniform 300;
+
+boundaryField
+{
+    inlet
+    {
+        type            fixedValue;
+        value           uniform 300;
+    }
+    outlet
+    {
+        type            zeroGradient;
+    }
+    top
+    {
+        type            zeroGradient;
+    }
+    slip-bottom
+    {
+        type            zeroGradient;
+    }
+    bottom
+    {
+        type            zeroGradient;
+    }
+    interface
+    {
+        type            fixedGradient;
+        gradient        uniform 0;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

test/tests/timesteppers/foam_controlled_tstep_cfl_insert/fluid-openfoam/0/U

@@ -0,0 +1,55 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     | Website:  https://openfoam.org
+    \\  /    A nd           | Version:  10
+     \\/     M anipulation  |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    format      ascii;
+    class       volVectorField;
+    location    "0";
+    object      U;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0.1 0 0);
+
+boundaryField
+{
+    inlet
+    {
+        type            fixedValue;
+        value           uniform (0.1 0 0);
+    }
+    outlet
+    {
+        type            zeroGradient;
+    }
+    top
+    {
+        type            slip;
+    }
+    slip-bottom
+    {
+        type            slip;
+    }
+    bottom
+    {
+        type            noSlip;
+    }
+    interface
+    {
+        type            noSlip;
+    }
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //