From 1c35e8a2f5c1cf93578b0c6c2a09208e4df8b050 Mon Sep 17 00:00:00 2001
From: Will Bainbridge <http://cfd.direct>
Date: Mon, 5 Nov 2018 09:28:03 +0000
Subject: [PATCH] reactingEulerFoam: Added update interval for population
 balance source terms.

This can be used to speedup simulations that converge to a steady state
by only updating the source terms once every few iterations. The number
of iterations between each update is set in fvSolution as follows:

    solvers
    {
        bubbles
        {
            sourceUpdateInterval 10;
        }

        // etc ...
    }

By default the interval is 1, and so the sources update every time.

Based on a patch contributed by Juho Peltola, VTT.
---
 .../populationBalanceModel.C                  | 38 +++++++++++--------
 .../populationBalanceModel.H                  | 11 ++++++
 .../populationBalanceModelI.H                 | 17 +++++++++
 3 files changed, 51 insertions(+), 15 deletions(-)

diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.C
index 5260a24433..9a333508fd 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.C
@@ -841,6 +841,15 @@ void Foam::diameterModels::populationBalanceModel::calcVelocity()
     }
 }
 
+bool Foam::diameterModels::populationBalanceModel::updateSources()
+{
+    const bool result = sourceUpdateCounter_ % sourceUpdateInterval() == 0;
+
+    ++ sourceUpdateCounter_;
+
+    return result;
+}
+
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
@@ -927,7 +936,11 @@ Foam::diameterModels::populationBalanceModel::populationBalanceModel
     nucleationRate_(),
     alphas_(),
     dsm_(),
-    U_()
+    U_(),
+    sourceUpdateCounter_
+    (
+        (mesh_.time().timeIndex()*nCorr())%sourceUpdateInterval()
+    )
 {
     this->registerVelocityAndSizeGroups();
 
@@ -1224,22 +1237,14 @@ Foam::diameterModels::populationBalanceModel::continuousTurbulence() const
 void Foam::diameterModels::populationBalanceModel::solve()
 {
     const dictionary& solutionControls = mesh_.solverDict(name_);
-    bool solveOnFinalIterOnly
-        (
-            solutionControls.lookupOrDefault<bool>
-            (
-                "solveOnFinalIterOnly",
-                false
-            )
-        );
+    bool solveOnFinalIterOnly =
+        solutionControls.lookupOrDefault<bool>("solveOnFinalIterOnly", false);
 
     if (!solveOnFinalIterOnly || pimple_.finalIter())
     {
-        label nCorr(readLabel(solutionControls.lookup("nCorr")));
-        scalar tolerance
-            (
-                readScalar(solutionControls.lookup("tolerance"))
-            );
+        const label nCorr = this->nCorr();
+        const scalar tolerance =
+            readScalar(solutionControls.lookup("tolerance"));
 
         if (nCorr > 0)
         {
@@ -1257,7 +1262,10 @@ void Foam::diameterModels::populationBalanceModel::solve()
                 << iCorr
                 << endl;
 
-            sources();
+            if (updateSources())
+            {
+                sources();
+            }
 
             dmdt();
 
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.H
index 347e332c73..89e30eb87f 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.H
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModel.H
@@ -277,6 +277,9 @@ class populationBalanceModel
         //- Average velocity
         autoPtr<volVectorField> U_;
 
+        //- Counter for interval between source term updates
+        label sourceUpdateCounter_;
+
 
     // Private member functions
 
@@ -316,6 +319,14 @@ class populationBalanceModel
 
         void calcVelocity();
 
+        //- Return whether the sources should be updated on this iteration
+        bool updateSources();
+
+        //- Return the number of corrections
+        inline label nCorr() const;
+
+        //- Return the interval at which the sources are updated
+        inline label sourceUpdateInterval() const;
 
 public:
 
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModelI.H b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModelI.H
index 5301f104cc..1de910ecf4 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModelI.H
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/populationBalanceModel/populationBalanceModel/populationBalanceModelI.H
@@ -23,6 +23,23 @@ License
 
 \*---------------------------------------------------------------------------*/
 
+// * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * * //
+
+inline Foam::label Foam::diameterModels::populationBalanceModel::nCorr() const
+{
+    return mesh_.solverDict(name_).lookupType<label>("nCorr");
+}
+
+
+inline Foam::label
+Foam::diameterModels::populationBalanceModel::sourceUpdateInterval() const
+{
+    return
+        mesh_.solverDict(name_)
+       .lookupOrDefault<label>("sourceUpdateInterval", 1);
+}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 inline const Foam::phaseSystem&