Merge branch 'feature-interIsoFoam-update' into 'develop'

Update of interIsoFoam and isoAdvection extending to work with porous media

See merge request Development/openfoam!511

applications/solvers/multiphase/interIsoFoam/UEqn.H

@@ -0,0 +1,35 @@
+    MRF.correctBoundaryVelocity(U);
+
+    fvVectorMatrix UEqn
+    (
+        fvm::ddt(rho, U) + fvm::div(rhoPhi, U)
+      + MRF.DDt(rho, U)
+      + turbulence->divDevRhoReff(rho, U)
+     ==
+        fvOptions(rho, U)
+    );
+
+    #include "UEqnAddPorosity.H"
+
+    UEqn.relax();
+
+    fvOptions.constrain(UEqn);
+
+    if (pimple.momentumPredictor())
+    {
+        solve
+        (
+            UEqn
+         ==
+            fvc::reconstruct
+            (
+                (
+                    mixture.surfaceTensionForce()
+                  - ghf*fvc::snGrad(rho)
+                  - fvc::snGrad(p_rgh)
+                ) * mesh.magSf()
+            )
+        );
+
+        fvOptions.correct(U);
+    }

applications/solvers/multiphase/interIsoFoam/UEqnAddPorosity.H

@@ -0,0 +1,48 @@
+// Including porosity effects in UEqn following:
+// Jensen, B., Jacobsen, N. G., & Christensen, E. D. (2014).
+// Investigations on the porous media equations and resistance
+// coefficients for coastal structures. Coastal Engineering, 84, 56-72.
+
+if (porosityEnabled)
+{
+    const volScalarField& porosity = tporosity.cref();
+
+    const word porosityModel("JensenEtAl2014");
+    const dictionary& dict =
+        porosityProperties.subDict(porosityModel + "Coeffs");
+    const dimensionedScalar alpha(dimless/dimArea, dict.get<scalar>("alpha"));
+    const dimensionedScalar beta(dimless/dimLength, dict.get<scalar>("beta"));
+    const dimensionedScalar d50(dimless, dict.get<scalar>("d50"));
+    const dimensionedScalar KC(dimless, dict.get<scalar>("KC"));
+
+    // Generating Darcy-Forchheimer coefficient: F = rho*U*(a + b*|U|)
+    // Shoud it be mu or muEff in the equation below?
+    {
+        // Darcy term
+        volScalarField DarcyForchheimerCoeff
+        (
+            alpha*sqr(1 - porosity)*mixture.mu()/sqr(porosity)/sqr(d50)
+        );
+
+        // Adding Forchheimer term
+        DarcyForchheimerCoeff += rho*mag(U)
+            *beta*(1 + pos(KC)*7.5/KC)*(1 - porosity)/sqr(porosity)/d50;
+
+        // Adding Darcy-Forchheimer term as implicit source term
+        UEqn += fvm::Sp(DarcyForchheimerCoeff, U);
+    }
+
+    {
+        // Generating added mass force coefficient
+        const dimensionedScalar gamma_p(dimless, dict.get<scalar>("gamma_p"));
+        const volScalarField Cm(gamma_p*(1 - porosity));
+
+        UEqn += Cm*fvm::ddt(rho, U);
+        UEqn += Cm*MRF.DDt(rho, U);
+    }
+
+    // Dividing both matrix entries and source term by porosity to compensate
+    // for the fact that the FVM cell volume averages use division by cell
+    // volume V whereas only the cell pore volume, porosity*V, is accessible.
+    UEqn *= scalar(1)/porosity;
+}

applications/solvers/multiphase/interIsoFoam/alphaEqn.H

@@ -22,8 +22,20 @@
     mixture.correct();
 }
 
+scalar domainFraction = 0;
+if (porosityEnabled)
+{
+    const volScalarField& porosity = tporosity.cref();
+    rhoPhi *= scalar(1)/fvc::interpolate(porosity);
+    domainFraction = alpha1.weightedAverage(mesh.Vsc()*porosity).value();
+}
+else
+{
+    domainFraction = alpha1.weightedAverage(mesh.Vsc()).value();
+}
+
 Info<< "Phase-1 volume fraction = "
-    << alpha1.weightedAverage(mesh.Vsc()).value()
+    << domainFraction
     << "  Min(" << alpha1.name() << ") = " << min(alpha1).value()
     << "  Max(" << alpha1.name() << ") = " << max(alpha1).value()
     << endl;

applications/solvers/multiphase/interIsoFoam/createFields.H

@@ -130,5 +130,6 @@ mesh.setFluxRequired(alpha1.name());
 
 #include "createMRF.H"
 #include "createFvOptions.H"
+#include "createPorosity.H"
 
 isoAdvection advector(alpha1, phi, U);

applications/solvers/multiphase/interIsoFoam/createPorosity.H

@@ -0,0 +1,34 @@
+// Reading porosity properties from constant directory
+IOdictionary porosityProperties
+(
+    IOobject
+    (
+        "porosityProperties",
+        runTime.constant(),
+        runTime,
+        IOobject::READ_IF_PRESENT,
+        IOobject::NO_WRITE
+    )
+);
+
+const bool porosityEnabled
+(
+    porosityProperties.getOrDefault<bool>("porosityEnabled", false)
+);
+
+tmp<volScalarField> tporosity;
+if (porosityEnabled)
+{
+    tporosity = tmp<volScalarField>::New
+    (
+        IOobject
+        (
+            "porosity",
+            runTime.timeName(),
+            mesh,
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE
+        ),
+        mesh
+    );
+}

applications/solvers/multiphase/interIsoFoam/interIsoFoam.C

@@ -92,7 +92,7 @@ int main(int argc, char *argv[])
     #include "initCorrectPhi.H"
     #include "createUfIfPresent.H"
 
-    #include "CourantNo.H"
+    #include "porousCourantNo.H"
     #include "setInitialDeltaT.H"
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -101,8 +101,8 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readDyMControls.H"
-        #include "CourantNo.H"
-        #include "alphaCourantNo.H"
+        #include "porousCourantNo.H"
+        #include "porousAlphaCourantNo.H"
         #include "setDeltaT.H"
 
         ++runTime;

applications/solvers/multiphase/interIsoFoam/porousAlphaCourantNo.H

@@ -0,0 +1,69 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 Johan Roenby
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Global
+    porousAlphaCourantNo
+
+Description
+    Calculates and outputs the mean and maximum Courant Numbers.
+
+\*---------------------------------------------------------------------------*/
+
+scalar maxAlphaCo = runTime.controlDict().get<scalar>("maxAlphaCo");
+scalar alphaCoNum = 0;
+scalar meanAlphaCoNum = 0;
+
+if (mesh.nInternalFaces())
+{
+    scalarField sumPhi
+    (
+        mixture.nearInterface()().primitiveField()
+       *fvc::surfaceSum(mag(phi))().primitiveField()
+    );
+
+    if (porosityEnabled)
+    {
+        const volScalarField& porosity = tporosity.cref();
+
+        alphaCoNum =
+            0.5*gMax(sumPhi/mesh.V().field()/porosity)*runTime.deltaTValue();
+
+        meanAlphaCoNum =
+            0.5*(gSum(sumPhi)/gSum(mesh.V().field()*porosity))
+           *runTime.deltaTValue();
+    }
+    else
+    {
+        alphaCoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue();
+
+        meanAlphaCoNum =
+            0.5*(gSum(sumPhi)/gSum(mesh.V().field()))*runTime.deltaTValue();
+    }
+}
+
+Info<< "Interface Courant Number mean: " << meanAlphaCoNum
+    << " max: " << alphaCoNum << endl;
+
+// ************************************************************************* //

applications/solvers/multiphase/interIsoFoam/porousCourantNo.H

@@ -0,0 +1,67 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2021 Johan Roenby
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Global
+    porousCourantNo
+
+Description
+    Calculates and outputs the mean and maximum Courant Numbers.
+
+\*---------------------------------------------------------------------------*/
+
+scalar CoNum = 0;
+scalar meanCoNum = 0;
+
+{
+    scalarField sumPhi
+    (
+        fvc::surfaceSum(mag(phi))().primitiveField()
+    );
+
+    if (porosityEnabled)
+    {
+        const volScalarField& porosity = tporosity.cref();
+
+        CoNum =
+            0.5*gMax(sumPhi/mesh.V().field()/porosity)
+           *runTime.deltaTValue();
+
+        meanCoNum =
+            0.5*(gSum(sumPhi)/gSum(mesh.V().field()*porosity))
+           *runTime.deltaTValue();
+    }
+    else
+    {
+        CoNum = 0.5*gMax(sumPhi/mesh.V().field())*runTime.deltaTValue();
+
+        meanCoNum =
+            0.5*(gSum(sumPhi)/gSum(mesh.V().field()))*runTime.deltaTValue();
+    }
+}
+
+Info<< "Courant Number mean: " << meanCoNum
+    << " max: " << CoNum << endl;
+
+// ************************************************************************* //

src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvection.C

@@ -7,8 +7,8 @@
 -------------------------------------------------------------------------------
     Copyright (C) 2016-2017 DHI
     Modified code Copyright (C) 2016-2017 OpenCFD Ltd.
-    Modified code Copyright (C) 2019 Johan Roenby
     Modified code Copyright (C) 2019-2020 DLR
+    Modified code Copyright (C) 2018, 2021 Johan Roenby
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -39,6 +39,7 @@ License
 #include "meshTools.H"
 #include "OBJstream.H"
 #include "syncTools.H"
+#include "profiling.H"
 
 #include "addToRunTimeSelectionTable.H"
 
@@ -108,6 +109,10 @@ Foam::isoAdvection::isoAdvection
     bsn0_(bsFaces_.size()),
     bsUn0_(bsFaces_.size()),
 
+    // Porosity
+    porosityEnabled_(dict_.getOrDefault<bool>("porosityEnabled", false)),
+    porosityPtr_(nullptr),
+
     // Parallel run data
     procPatchLabels_(mesh_.boundary().size()),
     surfaceCellFacesOnProcPatches_(0)
@@ -132,6 +137,48 @@ Foam::isoAdvection::isoAdvection
         // Get boundary mesh and resize the list for parallel comms
         setProcessorPatches();
     }
+
+    // Reading porosity properties from constant directory
+    IOdictionary porosityProperties
+    (
+        IOobject
+        (
+            "porosityProperties",
+            mesh_.time().constant(),
+            mesh_,
+            IOobject::READ_IF_PRESENT,
+            IOobject::NO_WRITE
+        )
+    );
+
+    porosityEnabled_ =
+        porosityProperties.getOrDefault<bool>("porosityEnabled", false);
+
+    if (porosityEnabled_)
+    {
+        if (mesh_.foundObject<volScalarField>("porosity"))
+        {
+            porosityPtr_ = mesh_.getObjectPtr<volScalarField>("porosity");
+
+            if
+            (
+                gMin(porosityPtr_->primitiveField()) <= 0
+             || gMax(porosityPtr_->primitiveField()) > 1 + SMALL
+            )
+            {
+                FatalErrorInFunction
+                    << "Porosity field has values <= 0 or > 1"
+                    << exit(FatalError);
+            }
+        }
+        else
+        {
+            FatalErrorInFunction
+                << "Porosity enabled in constant/porosityProperties "
+                << "but no porosity field is found in object registry."
+                << exit(FatalError);
+        }
+    }
 }
 
 
@@ -195,6 +242,7 @@ void Foam::isoAdvection::extendMarkedCells
 
 void Foam::isoAdvection::timeIntegratedFlux()
 {
+    addProfilingInFunction(geometricVoF);
     // Get time step
     const scalar dt = mesh_.time().deltaTValue();
 
@@ -222,6 +270,26 @@ void Foam::isoAdvection::timeIntegratedFlux()
     DynamicList<List<point>> isoFacePts;
     const DynamicField<label>& interfaceLabels = surf_->interfaceLabels();
 
+    // Calculating isoface normal velocity
+    scalarField Un0(interfaceLabels.size());
+    forAll(Un0, i)
+    {
+        const label celli = interfaceLabels[i];
+        const point x0(surf_->centre()[celli]);
+        const vector n0(normalised(-surf_->normal()[celli]));
+        Un0[i] = UInterp.interpolate(x0, celli) & n0;
+    }
+
+    // Taking acount of porosity if enabled
+    if (porosityEnabled_)
+    {
+        forAll(Un0, i)
+        {
+            const label celli = interfaceLabels[i];
+            Un0[i] /= porosityPtr_->primitiveField()[celli];
+        }
+    }
+
     // Loop through cells
     forAll(interfaceLabels, i)
     {
@@ -232,6 +300,8 @@ void Foam::isoAdvection::timeIntegratedFlux()
             // This is a surface cell, increment counter, append and mark cell
             nSurfaceCells++;
             surfCells_.append(celli);
+            const point x0(surf_->centre()[celli]);
+            const vector n0(normalised(-surf_->normal()[celli]));
 
             DebugInfo
                 << "\n------------ Cell " << celli << " with alpha1 = "
@@ -239,19 +309,6 @@ void Foam::isoAdvection::timeIntegratedFlux()
                 << 1.0 - alpha1In_[celli] << " ------------"
                 << endl;
 
-            // Cell is cut
-            const point x0 = surf_->centre()[celli];
-            const vector n0(normalised(-surf_->normal()[celli]));
-
-            // Get the speed of the isoface by interpolating velocity and
-            // dotting it with isoface unit normal
-            const scalar Un0 = UInterp.interpolate(x0, celli) & n0;
-
-            DebugInfo
-                << "calcIsoFace gives initial surface: \nx0 = " << x0
-                << ", \nn0 = " << n0 << ", \nUn0 = "
-                << Un0 << endl;
-
             // Estimate time integrated flux through each downwind face
             // Note: looping over all cell faces - in reduced-D, some of
             //       these faces will be on empty patches
@@ -284,7 +341,7 @@ void Foam::isoAdvection::timeIntegratedFlux()
                             facei,
                             x0,
                             n0,
-                            Un0,
+                            Un0[i],
                             dt,
                             phiIn[facei],
                             magSfIn[facei]
@@ -297,7 +354,7 @@ void Foam::isoAdvection::timeIntegratedFlux()
                     bsFaces_.append(facei);
                     bsx0_.append(x0);
                     bsn0_.append(n0);
-                    bsUn0_.append(Un0);
+                    bsUn0_.append(Un0[i]);
 
                     // Note: we must not check if the face is on the
                     // processor patch here.
@@ -552,6 +609,7 @@ void Foam::isoAdvection::checkIfOnProcPatch(const label facei)
 
 void Foam::isoAdvection::applyBruteForceBounding()
 {
+    addProfilingInFunction(geometricVoF);
     bool alpha1Changed = false;
 
     const scalar snapAlphaTol = dict_.getOrDefault<scalar>("snapTol", 0);

src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvection.H

@@ -7,8 +7,8 @@
 -------------------------------------------------------------------------------
     Copyright (C) 2016-2017 DHI
     Modified code Copyright (C) 2016-2019 OpenCFD Ltd.
-    Modified code Copyright (C) 2018 Johan Roenby
     Modified code Copyright (C) 2019 DLR
+    Modified code Copyright (C) 2018, 2021 Johan Roenby
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -30,14 +30,13 @@ Class
     Foam::isoAdvection
 
 Description
-    Calculates the new VOF (alpha) field after time step dt given the initial
-    VOF field and a velocity field U and face fluxes phi. The fluid transport
-    calculation is based on an idea of using isosurfaces to estimate the
-    internal distribution of fluid in cells and advecting such isosurfaces
-    across the mesh faces with the velocity field interpolated to the
-    isosurfaces.
+    An implementation of the isoAdvector geometric Volume-of-Fluid method
+    advancing the provided volume fraction field (alpha1) in time using the
+    given velocity field, U, and corresponding face fluxes, phi.
 
-    Reference:
+    References:
+
+    Main isoAdvector idea:
     \verbatim
         Roenby, J., Bredmose, H. and Jasak, H. (2016).
         A computational method for sharp interface advection
@@ -45,6 +44,20 @@ Description
         doi 10.1098/rsos.160405
     \endverbatim
 
+    Calculation of rhoPhi:
+    \verbatim
+        Roenby, J., Bredmose, H., & Jasak, H. (2019).
+        IsoAdvector: Geometric VOF on general meshes.
+        OpenFOAM® (pp. 281-296). Springer, Cham.
+    \endverbatim
+
+    Extension to porous media flows:
+    \verbatim
+        Missios, K., Jacobsen, N. G., Moeller, K., & Roenby, J.
+        Using the isoAdvector Geometric VOF Method for Interfacial Flows 
+        Through Porous Media. MARINE 2021.
+    \endverbatim
+
     Original code supplied by Johan Roenby, DHI (2016)
     Modified Henning Scheufler, DLR
 
@@ -157,7 +170,11 @@ class isoAdvection
             //- Storage for boundary surface iso face speed
             DynamicScalarList bsUn0_;
 
+            //- Switch telling if porosity is enabled
+            bool porosityEnabled_;
 
+            //- Pointer to porosity field (read from objectRegistry if exists)
+            volScalarField* porosityPtr_;
 
         // Additional data for parallel runs
 

src/transportModels/geometricVoF/advectionSchemes/isoAdvection/isoAdvectionTemplates.C

@@ -8,6 +8,7 @@
     Copyright (C) 2016-2017 DHI
     Modified code Copyright (C) 2016-2017 OpenCFD Ltd.
     Modified code Copyright (C) 2019-2020 DLR
+    Modified code Copyright (C) 2021 Johan Roenby
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -118,6 +119,7 @@ void Foam::isoAdvection::limitFluxes
     const SuType& Su
 )
 {
+    addProfilingInFunction(geometricVoF);
     DebugInFunction << endl;
 
     const scalar aTol = 100*SMALL;          // Note: tolerances
@@ -162,14 +164,23 @@ void Foam::isoAdvection::limitFluxes
                 {
                     checkIfOnProcPatch(facei);
                     const label own = owner[facei];
+                    scalar Vown = mesh_.V()[own];
+                    if (porosityEnabled_)
+                    {
+                        Vown *= porosityPtr_->primitiveField()[own];
+                    }
+                    alpha1_[own] -= faceValue(dVfcorrectionValues, facei)/Vown;
 
-                    alpha1_[own] -=
-                        faceValue(dVfcorrectionValues, facei)/mesh_.V()[own];
                     if (mesh_.isInternalFace(facei))
                     {
                         const label nei = neighbour[facei];
+                        scalar Vnei = mesh_.V()[nei];
+                        if (porosityEnabled_)
+                        {
+                            Vnei *= porosityPtr_->primitiveField()[nei];
+                        }
                         alpha1_[nei] +=
-                            faceValue(dVfcorrectionValues, facei)/mesh_.V()[nei];
+                            faceValue(dVfcorrectionValues, facei)/Vnei;
                     }
 
                     // Change to treat boundaries consistently
@@ -222,6 +233,7 @@ void Foam::isoAdvection::boundFlux
     const SuType& Su
 )
 {
+    addProfilingInFunction(geometricVoF);
     DebugInFunction << endl;
     scalar rDeltaT = 1/mesh_.time().deltaTValue();
 
@@ -242,7 +254,12 @@ void Foam::isoAdvection::boundFlux
     {
         if (alpha1_[celli] < -aTol || alpha1_[celli] > 1 + aTol)
         {
-            const scalar Vi = meshV[celli];
+            scalar Vi = meshV[celli];
+            if (porosityEnabled_)
+            {
+                Vi *= porosityPtr_->primitiveField()[celli];
+            }
+
             scalar alphaOvershoot =
                 pos0(alpha1_[celli] - 1)*(alpha1_[celli] - 1)
               + neg0(alpha1_[celli])*alpha1_[celli];
@@ -374,6 +391,7 @@ void Foam::isoAdvection::boundFlux
 template<class SpType, class SuType>
 void Foam::isoAdvection::advect(const SpType& Sp, const SuType& Su)
 {
+    addProfilingInFunction(geometricVoF);
     DebugInFunction << endl;
 
     if (mesh_.topoChanging())
@@ -409,21 +427,30 @@ void Foam::isoAdvection::advect(const SpType& Sp, const SuType& Su)
     }
 
     // Advect the free surface
-    alpha1_.primitiveFieldRef() =
-    (
-        alpha1_.oldTime().primitiveField()*rDeltaT
-      + Su.field()
-      - fvc::surfaceIntegrate(dVf_)().primitiveField()*rDeltaT
-    )/(rDeltaT - Sp.field());
+    if (porosityEnabled_)
+    {
+        // Should Su and Sp also be divided by porosity?
+        alpha1_.primitiveFieldRef() =
+        (
+            alpha1_.oldTime().primitiveField()*rDeltaT + Su.field()
+          - fvc::surfaceIntegrate(dVf_)().primitiveField()*rDeltaT
+            /porosityPtr_->primitiveField()
+        )/(rDeltaT - Sp.field());
+    }
+    else
+    {
+        alpha1_.primitiveFieldRef() =
+        (
+            alpha1_.oldTime().primitiveField()*rDeltaT
+          + Su.field()
+          - fvc::surfaceIntegrate(dVf_)().primitiveField()*rDeltaT
+        )/(rDeltaT - Sp.field());
+    }
 
     alpha1_.correctBoundaryConditions();
 
     // Adjust dVf for unbounded cells
-    limitFluxes
-    (
-        Sp,
-        Su
-    );
+    limitFluxes(Sp, Su);
 
     scalar maxAlphaMinus1 = gMax(alpha1In_) - 1;
     scalar minAlpha = gMin(alpha1In_);

src/transportModels/geometricVoF/cellCuts/cutFace/cutFace.C

@@ -275,10 +275,10 @@ void Foam::cutFace::calcSubFaceCentreAndArea
     if (nPoints == 3)
     {
         subFaceCentre =
-                (1.0 / 3.0) * (subFacePoints[0] + subFacePoints[1] + subFacePoints[2]);
+            (1.0/3.0)*(subFacePoints[0] + subFacePoints[1] + subFacePoints[2]);
 
-        subFaceArea = 0.5 * ((subFacePoints[1] - subFacePoints[0]) ^
-                             (subFacePoints[2] - subFacePoints[0]));
+        subFaceArea = 0.5*((subFacePoints[1] - subFacePoints[0]) ^
+            (subFacePoints[2] - subFacePoints[0]));
     }
     else if (nPoints > 0)
     {

src/transportModels/geometricVoF/reconstructionSchemes/isoSchemes/isoAlpha/isoAlpha.C

@@ -28,6 +28,7 @@ License
 #include "isoAlpha.H"
 #include "addToRunTimeSelectionTable.H"
 #include "cutCellPLIC.H"
+#include "profiling.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -76,6 +77,7 @@ Foam::reconstruction::isoAlpha::isoAlpha
 
 void Foam::reconstruction::isoAlpha::reconstruct(bool forceUpdate)
 {
+    addProfilingInFunction(geometricVoF);
     const bool uptodate = alreadyReconstructed(forceUpdate);
 
     if (uptodate && !forceUpdate)

src/transportModels/geometricVoF/reconstructionSchemes/plicSchemes/gradAlpha/gradAlpha.C

@@ -29,7 +29,7 @@ License
 #include "fvc.H"
 #include "leastSquareGrad.H"
 #include "addToRunTimeSelectionTable.H"
-
+#include "profiling.H"
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
 namespace Foam
@@ -46,6 +46,7 @@ namespace reconstruction
 
 void Foam::reconstruction::gradAlpha::gradSurf(const volScalarField& phi)
 {
+    addProfilingInFunction(geometricVoF);
     leastSquareGrad<scalar> lsGrad("polyDegree1",mesh_.geometricD());
 
     zoneDistribute& exchangeFields = zoneDistribute::New(mesh_);
@@ -123,6 +124,7 @@ Foam::reconstruction::gradAlpha::gradAlpha
 
 void Foam::reconstruction::gradAlpha::reconstruct(bool forceUpdate)
 {
+    addProfilingInFunction(geometricVoF);
     const bool uptodate = alreadyReconstructed(forceUpdate);
 
     if (uptodate && !forceUpdate)

src/transportModels/geometricVoF/reconstructionSchemes/plicSchemes/plicRDF/plicRDF.C

@@ -30,6 +30,7 @@ License
 #include "fvc.H"
 #include "leastSquareGrad.H"
 #include "addToRunTimeSelectionTable.H"
+#include "profiling.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -47,6 +48,7 @@ namespace reconstruction
 
 void Foam::reconstruction::plicRDF::interpolateNormal()
 {
+    addProfilingInFunction(geometricVoF);
     scalar dt = mesh_.time().deltaTValue();
     zoneDistribute& exchangeFields = zoneDistribute::New(mesh_);
 
@@ -159,6 +161,7 @@ void Foam::reconstruction::plicRDF::interpolateNormal()
 
 void Foam::reconstruction::plicRDF::gradSurf(const volScalarField& phi)
 {
+    addProfilingInFunction(geometricVoF);
     leastSquareGrad<scalar> lsGrad("polyDegree1", mesh_.geometricD());
     zoneDistribute& exchangeFields = zoneDistribute::New(mesh_);
 
@@ -205,6 +208,7 @@ void Foam::reconstruction::plicRDF::gradSurf(const volScalarField& phi)
 
 void Foam::reconstruction::plicRDF::setInitNormals(bool interpolate)
 {
+    addProfilingInFunction(geometricVoF);
     zoneDistribute& exchangeFields = zoneDistribute::New(mesh_);
 
     interfaceLabels_.clear();
@@ -239,6 +243,7 @@ void Foam::reconstruction::plicRDF::calcResidual
     List<normalRes>& normalResidual
 )
 {
+    addProfilingInFunction(geometricVoF);
     zoneDistribute& exchangeFields = zoneDistribute::New(mesh_);
     exchangeFields.setUpCommforZone(interfaceCell_,false);
 
@@ -377,6 +382,7 @@ Foam::reconstruction::plicRDF::plicRDF
 
 void Foam::reconstruction::plicRDF::reconstruct(bool forceUpdate)
 {
+    addProfilingInFunction(geometricVoF);
     zoneDistribute& exchangeFields = zoneDistribute::New(mesh_);
     const bool uptodate = alreadyReconstructed(forceUpdate);
 

tutorials/multiphase/interIsoFoam/damBreak/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss limitedLinearV 1;
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss interfaceCompression;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
 }
 

tutorials/multiphase/interIsoFoam/damBreakWithObstacle/constant/dynamicMeshDict

@@ -47,7 +47,7 @@ correctFluxes
     (phi none)
     (nHatf none)
     (rhoPhi none)
-    (alphaPhi none)
+    (alphaPhi_ none)
     (ghf none)
     (phi0 none)
     (dVf_ none)

tutorials/multiphase/interIsoFoam/damBreakWithObstacle/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss upwind;
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
     div(phi,k)       Gauss upwind;
     div(phi,omega)   Gauss upwind;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;

tutorials/multiphase/interIsoFoam/damBreakWithObstacle/system/fvSolution

@@ -29,9 +29,8 @@ solvers
         writeIsoFaces   false;
         reconstructionScheme plicRDF; // isoAlpha
 
-        nAlphaCorr      1;
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read by interIsoFoam but not used
     }
 
     p_rgh

tutorials/multiphase/interIsoFoam/discInConstantFlow/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rho*phi,U)   Gauss limitedLinearV 1;
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss interfaceCompression;
     div((muEff*dev(T(grad(U))))) Gauss linear;
 }
 

tutorials/multiphase/interIsoFoam/discInConstantFlow/system/fvSolution

@@ -19,12 +19,12 @@ solvers
 {
     "alpha.*"
     {
-        isoFaceTol      1e-8;
-        surfCellTol     1e-8;
-        nAlphaBounds    3;
+        isoFaceTol      1e-10;
+        surfCellTol     1e-10;
+        nAlphaBounds    5;
         snapTol         0;
         clip            false;
-        cAlpha          1;
+        cAlpha          1; //Read by interIsoFoam but not used
         nAlphaSubCycles 1;
         reconstructionScheme plicRDF; // isoAlpha
     }

tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss linearUpwind grad(U);
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
 }
 

tutorials/multiphase/interIsoFoam/discInConstantFlowCyclicBCs/system/fvSolution

@@ -27,17 +27,8 @@ solvers
         clip            true;
         reconstructionScheme isoAlpha; // isoAlpha
 
-        nAlphaCorr      2;
         nAlphaSubCycles 1;
-        cAlpha          1;
-
-        MULESCorr       yes;
-        nLimiterIter    3;
-
-        solver          smoothSolver;
-        smoother        symGaussSeidel;
-        tolerance       1e-8;
-        relTol          0;
+        cAlpha          1; // Read but not used by interIsoFoam
     }
 
     "pcorr.*"

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/0.orig/U

@@ -0,0 +1,40 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    object      U;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (1 0 0.5);
+
+boundaryField
+{
+    rim
+    {
+        type            zeroGradient;
+    }
+
+    front
+    {
+        type            empty;
+    }
+
+    back
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/0.orig/alpha.water

@@ -0,0 +1,40 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      alpha.water;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    rim
+    {
+        type    zeroGradient;
+    }
+
+    front
+    {
+        type    empty;
+    }
+
+    back
+    {
+        type    empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/0.orig/p_rgh

@@ -0,0 +1,40 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      p_rgh;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [1 -1 -2 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    rim
+    {
+        type    zeroGradient;
+    }
+
+    front
+    {
+        type    empty;
+    }
+
+    back
+    {
+        type    empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/0.orig/porosity

@@ -0,0 +1,40 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      porosity;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 1;
+
+boundaryField
+{
+    rim
+    {
+        type    zeroGradient;
+    }
+
+    front
+    {
+        type    empty;
+    }
+
+    back
+    {
+        type    empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/Allclean

@@ -0,0 +1,8 @@
+#!/bin/sh
+cd "${0%/*}" || exit                                # Run from this directory
+. ${WM_PROJECT_DIR:?}/bin/tools/CleanFunctions      # Tutorial clean functions
+#------------------------------------------------------------------------------
+
+cleanCase0
+
+#------------------------------------------------------------------------------

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/Allrun

@@ -0,0 +1,34 @@
+#!/bin/sh
+cd "${0%/*}" || exit                                # Run from this directory
+. ${WM_PROJECT_DIR:?}/bin/tools/RunFunctions        # Tutorial run functions
+#------------------------------------------------------------------------------
+
+restore0Dir
+
+runApplication blockMesh
+
+runApplication setAlphaField
+
+runApplication setFields
+
+runApplication decomposePar
+
+runParallel $(getApplication)
+
+runApplication reconstructPar
+
+# restart
+
+latestTime=$(foamListTimes -latestTime)
+
+mv -f "$latestTime" "$latestTime".bak
+
+rm -rf processor*
+
+runParallel -s decompose redistributePar -decompose -latestTime
+
+runParallel -s 2 $(getApplication)
+
+runParallel -s reconstruct redistributePar -reconstruct -latestTime
+
+#------------------------------------------------------------------------------

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/constant/g

@@ -0,0 +1,21 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       uniformDimensionedVectorField;
+    object      g;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -2 0 0 0 0];
+value           ( 0 -9.81 0 );
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/constant/porosityProperties

@@ -0,0 +1,20 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      porosityProperties;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+porosityEnabled true;
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/constant/transportProperties

@@ -0,0 +1,36 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      transportProperties;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+phases          (water air);
+
+water
+{
+    transportModel  Newtonian;
+    nu              0;
+    rho             1000;
+}
+
+air
+{
+    transportModel  Newtonian;
+    nu              0;
+    rho             1;
+}
+
+sigma           0;
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/blockMeshDict

@@ -0,0 +1,82 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      blockMeshDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+scale   1;
+
+L       5;
+nx      100;
+
+y1      -0.05;
+y2      0.05;
+ny      1;
+
+H       3;
+nz      60;
+
+vertices
+(
+    (0 $y1 0)
+    ($L $y1 0)
+    ($L $y2 0)
+    (0 $y2 0)
+    (0 $y1 $H)
+    ($L $y1 $H)
+    ($L $y2 $H)
+    (0 $y2 $H)
+);
+
+blocks
+(
+    hex (0 1 2 3 4 5 6 7) ($nx $ny $nz) simpleGrading (1 1 1)
+);
+
+edges
+(
+);
+
+boundary
+(
+    rim
+    {
+        type patch;
+        faces
+        (
+            (4 5 6 7)
+            (0 4 7 3)
+            (1 2 6 5)
+            (0 3 2 1)
+        );
+    }
+    front
+    {
+        type empty;
+        faces
+        (
+            (0 1 5 4)
+        );
+    }
+    back
+    {
+        type empty;
+        faces
+        (
+            (2 3 7 6)
+        );
+    }
+);
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/controlDict

@@ -0,0 +1,58 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version         2.0;
+    format          ascii;
+    class           dictionary;
+    object          controlDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+application     interIsoFoam;
+
+startFrom       latestTime;
+
+startTime       0.0;
+
+stopAt          endTime;
+
+endTime         2.5;
+
+writeControl    adjustable;
+
+writeInterval   0.1;
+
+deltaT          1e-3;
+
+purgeWrite      0;
+
+writeFormat     ascii;
+
+writePrecision  14;
+
+writeCompression off;
+
+timeFormat      general;
+
+timePrecision   6;
+
+graphFormat     raw;
+
+runTimeModifiable no;
+
+adjustTimeStep  yes;
+
+maxCo           1e6;
+
+maxAlphaCo      0.5;
+
+maxDeltaT       1e6;
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/decomposeParDict

@@ -0,0 +1,22 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      decomposeParDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+numberOfSubdomains 4;
+
+method          scotch;
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/fvSchemes

@@ -0,0 +1,57 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      fvSchemes;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+ddtSchemes
+{
+    default          Euler;
+}
+
+gradSchemes
+{
+    default          Gauss linear;
+}
+
+divSchemes
+{
+    div(rhoPhi,U)    Gauss limitedLinearV 1;
+    div((muEff*dev(T(grad(U))))) Gauss linear;
+}
+
+laplacianSchemes
+{
+    default          Gauss linear corrected;
+}
+
+interpolationSchemes
+{
+    default          linear;
+}
+
+snGradSchemes
+{
+    default          corrected;
+}
+
+fluxRequired
+{
+    default          no;
+    p_rgh;
+    pcorr;
+    alpha1;
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/fvSolution

@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      fvSolution;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+solvers
+{
+    "alpha.*"
+    {
+        isoFaceTol      1e-8;
+        surfCellTol     1e-8;
+        nAlphaBounds    3;
+        snapTol         0;
+        clip            false;
+        cAlpha          1; // Read but not used by interIsoFoam
+        nAlphaSubCycles 1;
+        reconstructionScheme plicRDF; // isoAlpha
+    }
+}
+
+PIMPLE
+{
+    frozenFlow      yes;
+    momentumPredictor no;
+    nCorrectors     -1;
+    nNonOrthogonalCorrectors -1;
+    pRefCell        0;
+    pRefValue       0;
+}
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/setAlphaFieldDict

@@ -0,0 +1,24 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      setAlphaFieldDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+field           alpha.water;
+type            cylinder;
+radius          0.25;
+direction       (0 1 0);
+origin          (0.5 0 0.5);
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInConstantPorousFlow/system/setFieldsDict

@@ -0,0 +1,36 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      setFieldsDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+defaultFieldValues
+(
+    volScalarFieldValue porosity 1
+);
+
+regions
+(
+   boxToCell
+    {
+        box (1 -0.05 0) (4 0.05 3);
+
+        fieldValues
+        (
+            volScalarFieldValue porosity 0.5
+        );
+    }
+);
+
+
+// ************************************************************************* //

tutorials/multiphase/interIsoFoam/discInReversedVortexFlow/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)   Gauss limitedLinearV 1;
-    div(phi,alpha)  Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
 }
 
 laplacianSchemes

tutorials/multiphase/interIsoFoam/discInReversedVortexFlow/system/fvSolution

@@ -26,7 +26,7 @@ solvers
         clip            false;
 
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read by interIsoFoam but not used
         reconstructionScheme plicRDF; // isoAlpha
     }
 }

tutorials/multiphase/interIsoFoam/iobasin/system/fvSchemes

@@ -29,8 +29,6 @@ divSchemes
 {
     default             none;
     div(rhoPhi,U)       Gauss limitedLinearV 1;
-    div(phi,alpha)      Gauss vanLeer;
-    div(phirb,alpha)    Gauss linear;
 
     "div\(phi,(k|omega)\)"      Gauss upwind;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;

tutorials/multiphase/interIsoFoam/iobasin/system/fvSolution

@@ -20,7 +20,7 @@ solvers
     "alpha.water.*"
     {
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read by interIsoFoam but not used
         reconstructionScheme isoAlpha; // plicRDF
         isoFaceTol      1e-8;
         surfCellTol     1e-6;

tutorials/multiphase/interIsoFoam/notchedDiscInSolidBodyRotation/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)   Gauss limitedLinearV 1;
-    div(phi,alpha)  Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
 }
 

tutorials/multiphase/interIsoFoam/notchedDiscInSolidBodyRotation/system/fvSolution

@@ -25,7 +25,7 @@ solvers
         nAlphaBounds    3;
         clip            false;
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read by interIsoFoam but not used
         reconstructionScheme isoAlpha; // plicRDF
     }
 }

tutorials/multiphase/interIsoFoam/sloshingTank2D/system/fvSchemes

@@ -27,8 +27,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss vanLeerV;
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss vanLeer;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
 }
 

tutorials/multiphase/interIsoFoam/sloshingTank2D/system/fvSolution

@@ -18,9 +18,8 @@ solvers
 {
     alpha.water
     {
-        nAlphaCorr      1;
         nAlphaSubCycles 3;
-        cAlpha          1.5;
+        cAlpha          1; // Read by interIsoFoam but not used
 
         isoFaceTol      1e-10;
         surfCellTol     1e-6;

tutorials/multiphase/interIsoFoam/sphereInReversedVortexFlow/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss limitedLinearV 1;
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
 }
 
 laplacianSchemes

tutorials/multiphase/interIsoFoam/sphereInReversedVortexFlow/system/fvSolution

@@ -25,7 +25,7 @@ solvers
         nAlphaBounds    8;
         clip            false;
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read but not used by interIsoFoam
         reconstructionScheme isoAlpha;
     }
 }

tutorials/multiphase/interIsoFoam/standingWave/system/fvSolution

@@ -26,7 +26,7 @@ solvers
         clip            true;
 
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read by interIsoFoam but not used
         reconstructionScheme isoAlpha; // plicRDF
     }
 

tutorials/multiphase/interIsoFoam/waveExampleStreamFunction/system/fvSchemes

@@ -27,8 +27,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss linearUpwind grad(U);
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
     div(((rho*nuEff)*dev2(T(grad(U))))) Gauss linear;
 }
 

tutorials/multiphase/interIsoFoam/waveExampleStreamFunction/system/fvSolution

@@ -24,9 +24,8 @@ solvers
         snapTol         1e-12;
         clip            true;
 
-        nAlphaCorr      1;
         nAlphaSubCycles 1;
-        cAlpha          1;
+        cAlpha          1; // Read by interIsoFoam but not used
 
         reconstructionScheme plicRDF; // plicRDF
     }

tutorials/multiphase/interIsoFoam/weirOverflow/system/fvSchemes

@@ -28,8 +28,6 @@ gradSchemes
 divSchemes
 {
     div(rhoPhi,U)    Gauss linear;
-    div(phi,alpha)   Gauss vanLeer;
-    div(phirb,alpha) Gauss linear;
     div(phi,k)       Gauss upwind;
     div(phi,epsilon) Gauss upwind;
     div(phi,R)       Gauss upwind;

tutorials/verificationAndValidation/multiphase/StefanProblem/setups.orig/common/constant/turbulenceProperties

@@ -0,0 +1,20 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      turbulenceProperties;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+simulationType  laminar;
+
+
+// ************************************************************************* //

tutorials/verificationAndValidation/multiphase/interIsoFoam/porousDamBreak/0.orig/U

@@ -0,0 +1,51 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volVectorField;
+    object      U;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 1 -1 0 0 0 0];
+
+internalField   uniform (0 0 0);
+
+boundaryField
+{
+    leftWall
+    {
+        type            slip;
+    }
+
+    rightWall
+    {
+        type            slip;
+    }
+
+    lowerWall
+    {
+        type            slip;
+    }
+
+    atmosphere
+    {
+        type            pressureInletOutletVelocity;
+        value           uniform (0 0 0);
+    }
+
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //

tutorials/verificationAndValidation/multiphase/interIsoFoam/porousDamBreak/0.orig/alpha.water

@@ -0,0 +1,52 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  v2106                                 |
+|   \\  /    A nd           | Website:  www.openfoam.com                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       volScalarField;
+    object      alpha.water;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dimensions      [0 0 0 0 0 0 0];
+
+internalField   uniform 0;
+
+boundaryField
+{
+    leftWall
+    {
+        type            zeroGradient;
+    }
+
+    rightWall
+    {
+        type            zeroGradient;
+    }
+
+    lowerWall
+    {
+        type            zeroGradient;
+    }
+
+    atmosphere
+    {
+        type            inletOutlet;
+        inletValue      $internalField;
+        value           $internalField;
+    }
+
+    defaultFaces
+    {
+        type            empty;
+    }
+}
+
+
+// ************************************************************************* //