ractingEulerFoam: populationBalance: Refactored scaling

Refactored the function for scaling the size group volume fractions to
better handle situations in which their sum drifts away from unity.
Scaling is now turned on by default, and can be turned off in the
solution dictionary for the population balance.

Additional revision and renaming of *Polydisperse tutorials

Patch contributed by Institute of Fluid Dynamics,
Helmholtz-Zentrum Dresden - Rossendorf (HZDR)

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/diameterModels/velocityGroup/velocityGroup.C

@@ -98,29 +98,28 @@ Foam::diameterModels::velocityGroup::fSum() const
 }
 
 
-void Foam::diameterModels::velocityGroup::renormalize()
+void Foam::diameterModels::velocityGroup::scale()
 {
-        Info<< "Renormalizing sizeGroups for velocityGroup "
+    Info<< "Scaling sizeGroups for velocityGroup " << phase_.name() << endl;
-            << phase_.name()
-            << endl;
 
-        // Set negative values to zero
+    forAll(sizeGroups_, i)
-        forAll(sizeGroups_, i)
+    {
-        {
+        sizeGroups_[i].max(0);
-            sizeGroups_[i] *= pos(sizeGroups_[i]);
+    };
-        };
 
-        forAll(sizeGroups_, i)
+    f_ = fSum();
-        {
+
-            sizeGroups_[i] /= fSum_;
+    forAll(sizeGroups_, i)
-        };
+    {
+        sizeGroups_[i] /= f_;
+    };
 }
 
 
 Foam::tmp<Foam::fv::convectionScheme<Foam::scalar>>
 Foam::diameterModels::velocityGroup::mvconvection() const
 {
-    tmp<fv::convectionScheme<Foam::scalar> > mvConvection
+    tmp<fv::convectionScheme<Foam::scalar>> mvConvection
     (
         fv::convectionScheme<Foam::scalar>::New
         (
@@ -174,24 +173,6 @@ Foam::diameterModels::velocityGroup::velocityGroup
         diameterProperties.lookup("sizeGroups"),
         sizeGroup::iNew(phase, *this)
     ),
-    fSum_
-    (
-        IOobject
-        (
-            IOobject::groupName
-            (
-                "fsum",
-                IOobject::groupName
-                (
-                    phase.name(),
-                    popBalName_
-                )
-            ),
-            phase.time().timeName(),
-            phase.mesh()
-        ),
-        fSum()
-    ),
     d_
     (
         IOobject
@@ -217,49 +198,6 @@ Foam::diameterModels::velocityGroup::velocityGroup
         dimensionedScalar(dimDensity/dimTime, Zero)
     )
 {
-    if
-    (
-        phase_.mesh().solverDict(popBalName_).lookupOrDefault<Switch>
-        (
-            "renormalizeAtRestart",
-            false
-        )
-     ||
-        phase_.mesh().solverDict(popBalName_).lookupOrDefault<Switch>
-        (
-            "renormalize",
-            false
-        )
-    )
-    {
-        renormalize();
-    }
-
-    fSum_ = fSum();
-
-    if
-    (
-        mag(1 - fSum_.weightedAverage(fSum_.mesh().V()).value()) >= 1e-5
-     || mag(1 - max(fSum_).value()) >= 1e-5
-     || mag(1 - min(fSum_).value()) >= 1e-5
-    )
-    {
-        FatalErrorInFunction
-            << " Initial values of the sizeGroups belonging to velocityGroup "
-            << this->phase().name()
-            << " must add to" << nl << " unity. This condition might be"
-            << " violated due to wrong entries in the" << nl
-            << " velocityGroupCoeffs subdictionary or bad initial conditions in"
-            << " the startTime" << nl
-            << " directory. The sizeGroups can be renormalized at every"
-            << " timestep or at restart" << nl
-            << " only by setting the corresponding switch renormalize or"
-            << " renormalizeAtRestart" << nl
-            << " in the fvSolution subdictionary " << popBalName_ << "."
-            << " Note that boundary conditions are not" << nl << "renormalized."
-            << exit(FatalError);
-    }
-
     forAll(sizeGroups_, i)
     {
         fields_.add(sizeGroups_[i]);
@@ -286,6 +224,28 @@ void Foam::diameterModels::velocityGroup::preSolve()
 
 void Foam::diameterModels::velocityGroup::postSolve()
 {
+    if
+    (
+        phase_.mesh().solverDict(popBalName_).lookupOrDefault<Switch>
+        (
+            "scale",
+            true
+        )
+    )
+    {
+        scale();
+    }
+
+    f_ = fSum();
+
+    f_.correctBoundaryConditions();
+
+    Info<< phase_.name() << " sizeGroups-sum volume fraction, min, max = "
+        << f_.weightedAverage(phase_.mesh().V()).value()
+        << ' ' << min(f_).value()
+        << ' ' << max(f_).value()
+        << endl;
+
     d_ = dsm();
 
     Info<< this->phase().name() << " Sauter mean diameter, min, max = "
@@ -293,26 +253,6 @@ void Foam::diameterModels::velocityGroup::postSolve()
         << ' ' << min(d_).value()
         << ' ' << max(d_).value()
         << endl;
-
-    fSum_ = fSum();
-
-    Info<< phase_.name() << " sizeGroups-sum volume fraction, min, max = "
-        << fSum_.weightedAverage(phase_.mesh().V()).value()
-        << ' ' << min(fSum_).value()
-        << ' ' << max(fSum_).value()
-        << endl;
-
-    if
-    (
-        phase_.mesh().solverDict(popBalName_).lookupOrDefault<Switch>
-        (
-            "renormalize",
-            false
-        )
-    )
-    {
-        renormalize();
-    }
 }
 
 

applications/solvers/multiphase/reactingEulerFoam/phaseSystems/diameterModels/velocityGroup/velocityGroup.H

@@ -100,7 +100,8 @@ class velocityGroup
         //- Name of the populationBalance this velocityGroup belongs to
         word popBalName_;
 
-        //- Reference field from which the sizeGroup fields are derived
+        //- Sum of the sizeGroup volume fractions and reference field from which
+        //  the sizeGroup fields are derived
         volScalarField f_;
 
         //- Form factor relating diameter and volume
@@ -109,9 +110,6 @@ class velocityGroup
         //- sizeGroups belonging to this velocityGroup
         PtrList<sizeGroup> sizeGroups_;
 
-        //- Sum of sizeGroup volume fractions
-        volScalarField fSum_;
-
         //- Number-based Sauter-mean diameter of the phase
         volScalarField d_;
 
@@ -131,7 +129,7 @@ class velocityGroup
 
         tmp<volScalarField> fSum() const;
 
-        void renormalize();
+        void scale();
 
         tmp<Foam::fv::convectionScheme<Foam::scalar>> mvconvection() const;
 

test/multiphase/reactingMultiphaseEulerFoam/populationBalanceModeling/binaryBreakup/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr           1;
         tolerance       1e-4;
-        renormalize     false;
+        scale           false;
     }
 
     p_rgh

test/multiphase/reactingMultiphaseEulerFoam/populationBalanceModeling/breakup/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr           1;
         tolerance       1e-4;
-        renormalize     false;
+        scale           false;
     }
 
     p_rgh

test/multiphase/reactingMultiphaseEulerFoam/populationBalanceModeling/coalescence/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr           1;
         tolerance       1e-4;
-        renormalize     false;
+        scale           false;
     }
 
     p_rgh

test/multiphase/reactingMultiphaseEulerFoam/populationBalanceModeling/simultaneousCoalescenceAndBreakup/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr           1;
         tolerance       1e-4;
-        renormalize     true;
+        scale           true;
     }
 
     p_rgh

test/multiphase/reactingTwoPhaseEulerFoam/populationBalanceModeling/drift/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr           1;
         tolerance       1e-4;
-        renormalize     false;
+        scale           false;
     }
 
     p_rgh

test/multiphase/reactingTwoPhaseEulerFoam/populationBalanceModeling/negativeDrift/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr           1;
         tolerance       1e-4;
-        renormalize     false;
+        scale           false;
     }
 
     p_rgh

tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumnFixedPolydisperse/Allrun

@@ -1,14 +1,11 @@
 #!/bin/sh
-cd ${0%/*} || exit 1    # run from this directory
+cd ${0%/*} || exit 1    # Run from this directory
 
 # Source tutorial run functions
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
-# Set application name
-application=$(getApplication)
-
 runApplication blockMesh
 runApplication setFields
-runApplication $application
+runApplication $(getApplication)
 
-# ----------------------------------------------------------------- end-of-file
+#------------------------------------------------------------------------------

tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumnPolydisperse/0/f.air1.bubbles

@@ -17,7 +17,7 @@ FoamFile
 
 dimensions      [0 0 0 0 0 0 0];
 
-internalField   uniform 0;
+internalField   uniform 1;
 
 boundaryField
 {

tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumnPolydisperse/0/f.air2.bubbles

@@ -17,7 +17,7 @@ FoamFile
 
 dimensions      [0 0 0 0 0 0 0];
 
-internalField   uniform 0;
+internalField   uniform 1;
 
 boundaryField
 {

tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumnPolydisperse/Allrun

@@ -1,15 +1,12 @@
 #!/bin/sh
-cd ${0%/*} || exit 1    # run from this directory
+cd ${0%/*} || exit 1    # Run from this directory
 
 # Source tutorial run functions
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
-# Set application name
-application=$(getApplication)
-
 runApplication blockMesh
 runApplication topoSet
 runApplication setFields
-runApplication $application
+runApplication $(getApplication)
 
-# ----------------------------------------------------------------- end-of-file
+#------------------------------------------------------------------------------

tutorials/multiphase/reactingMultiphaseEulerFoam/laminar/bubbleColumnPolydisperse/system/fvSolution

@@ -25,9 +25,10 @@ solvers
 
     bubbles
     {
-        nCorr           1;
+        nCorr                1;
-        tolerance       1e-4;
+        tolerance            1e-4;
-        renormalize     true;
+        scale                true;
+        solveOnFinalIterOnly true;
     }
 
     p_rgh

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/bubbleColumnPolydisperse/0/f.air.bubbles

@@ -17,7 +17,7 @@ FoamFile
 
 dimensions      [0 0 0 0 0 0 0];
 
-internalField   uniform 0;
+internalField   uniform 1;
 
 boundaryField
 {

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/bubbleColumnPolydisperse/Allrun

@@ -1,15 +1,12 @@
 #!/bin/sh
-cd ${0%/*} || exit 1    # run from this directory
+cd ${0%/*} || exit 1    # Run from this directory
 
 # Source tutorial run functions
 . $WM_PROJECT_DIR/bin/tools/RunFunctions
 
-# Set application name
-application=$(getApplication)
-
 runApplication blockMesh
 runApplication topoSet
 runApplication setFields
-runApplication $application
+runApplication $(getApplication)
 
 #------------------------------------------------------------------------------

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/bubbleColumnPolydisperse/constant/phaseProperties

@@ -81,21 +81,16 @@ populationBalanceCoeffs
 
         coalescenceModels
         (
-            CoulaloglouTavlarides{}
+            LehrMilliesMewes{}
         );
 
         binaryBreakupModels
-        ();
+        (
+            LehrMilliesMewes{}
+        );
 
         breakupModels
-        (
+        ();
-            exponential
-            {
-                C 10.0;
-                exponent 1.0;
-                daughterSizeDistributionModel uniform;
-            }
-        );
 
         driftModels
         (

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/bubbleColumnPolydisperse/constant/turbulenceProperties.air

@@ -24,11 +24,11 @@ RAS
     turbulence      on;
     printCoeffs     on;
 
-    // mixtureKEpsilonCoeffs
+    mixtureKEpsilonCoeffs
-    // {
+    {
-    //     Cp  1;
+        Cp  1;
-    //     C3  1;
+        C3  1;
-    // }
+    }
 }
 
 // ************************************************************************* //

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/bubbleColumnPolydisperse/system/fvSolution

@@ -27,7 +27,7 @@ solvers
     {
         nCorr                1;
         tolerance            1e-4;
-        renormalize          false;
+        scale                true;
         solveOnFinalIterOnly true;
     }
 

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/wallBoilingPolydisperse/constant/phaseProperties

@@ -76,24 +76,16 @@ populationBalanceCoeffs
 
         coalescenceModels
         (
-            hydrodynamic
+            LehrMilliesMewes{}
-            {
-                C 0.25;
-            }
         );
 
         binaryBreakupModels
-        ();
+        (
+            LehrMilliesMewes{}
+        );
 
         breakupModels
-        (
+        ();
-            exponential
-            {
-                C 0.5;
-                exponent 0.01;
-                daughterSizeDistributionModel uniform;
-            }
-        );
 
         driftModels
         (
@@ -125,7 +117,13 @@ blending
 }
 
 surfaceTension
-();
+(
+    (gas and liquid)
+    {
+        type            constant;
+        sigma           0.07;
+    }
+);
 
 saturationModel
 {

tutorials/multiphase/reactingTwoPhaseEulerFoam/RAS/wallBoilingPolydisperse/system/fvSolution

@@ -25,11 +25,10 @@ solvers
 
     bubbles
     {
-        nCorr                   1;
+        nCorr                1;
-        tolerance               1e-4;
+        tolerance            1e-4;
-        renormalizeOnRestart    true;
+        scale                true;
-        renormalize             true;
+        solveOnFinalIterOnly true;
-        solveOnFinalIterOnly    true;
     }
 
     p_rgh