interFoam family: Added run-time selectable LTS support

LTS is selected by the ddt scheme e.g. in the
tutorials/multiphase/interFoam/ras/DTCHull case:

ddtSchemes
{
    default         localEuler rDeltaT;
}

LTSInterFoam is no longer needed now that interFoam includes LTS
support.

applications/solvers/compressible/rhoCentralFoam/setRDeltaT.H

@@ -15,7 +15,7 @@
     (
         1/dimensionedScalar("maxDeltaT", dimTime, maxDeltaT),
         fvc::surfaceSum(amaxSf)().dimensionedInternalField()
-       /(2*maxCo*mesh.V())
+       /((2*maxCo)*mesh.V())
     );
 
     // Update tho boundary values of the reciprocal time-step

applications/solvers/multiphase/compressibleMultiphaseInterFoam/multiphaseMixtureThermo/Make/options

 EXE_INC = \
     -I$(LIB_SRC)/transportModels/compressible/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-    -I$(LIB_SRC)/finiteVolume/lnInclude
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude
 
 LIB_LIBS = \
     -lcompressibleTransportModels \
     -lfluidThermophysicalModels \
     -lspecie \
-    -lfiniteVolume
+    -lfiniteVolume \
+    -lmeshTools

applications/solvers/multiphase/driftFluxFoam/alphaEqn.H

@@ -6,11 +6,7 @@
     {
         fvScalarMatrix alpha1Eqn
         (
-            #ifdef LTSSOLVE
-            fv::localEulerDdtScheme<scalar>(mesh, rDeltaT.name()).fvmDdt(alpha1)
-            #else
             fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
-            #endif
           + fv::gaussConvectionScheme<scalar>
             (
                 mesh,
@@ -34,16 +30,6 @@
         {
             Info<< "Applying the previous iteration correction flux" << endl;
 
-            #ifdef LTSSOLVE
-            MULES::LTScorrect
-            (
-                alpha1,
-                phiAlpha,
-                tphiAlphaCorr0(),
-                mixture.alphaMax(),
-                0
-            );
-            #else
             MULES::correct
             (
                 alpha1,
@@ -52,7 +38,6 @@
                 mixture.alphaMax(),
                 0
             );
-            #endif
 
             phiAlpha += tphiAlphaCorr0();
         }
@@ -84,16 +69,6 @@
             tmp<surfaceScalarField> tphiAlphaCorr(tphiAlphaUn() - phiAlpha);
             volScalarField alpha10("alpha10", alpha1);
 
-            #ifdef LTSSOLVE
-            MULES::LTScorrect
-            (
-                alpha1,
-                tphiAlphaUn(),
-                tphiAlphaCorr(),
-                mixture.alphaMax(),
-                0
-            );
-            #else
             MULES::correct
             (
                 alpha1,
@@ -102,7 +77,6 @@
                 mixture.alphaMax(),
                 0
             );
-            #endif
 
             // Under-relax the correction for all but the 1st corrector
             if (aCorr == 0)
@@ -119,16 +93,6 @@
         {
             phiAlpha = tphiAlphaUn;
 
-            #ifdef LTSSOLVE
-            MULES::explicitLTSSolve
-            (
-                alpha1,
-                phi,
-                phiAlpha,
-                mixture.alphaMax(),
-                0
-            );
-            #else
             MULES::explicitSolve
             (
                 alpha1,
@@ -137,7 +101,6 @@
                 mixture.alphaMax(),
                 0
             );
-            #endif
         }
     }
 

applications/solvers/multiphase/interFoam/LTSInterFoam/LTSInterFoam.C

-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-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/>.
-
-Application
-    interFoam
-
-Description
-    Solver for 2 incompressible, isothermal immiscible fluids using a VOF
-    (volume of fluid) phase-fraction based interface capturing approach.
-
-    The momentum and other fluid properties are of the "mixture" and a single
-    momentum equation is solved.
-
-    Turbulence modelling is generic, i.e.  laminar, RAS or LES may be selected.
-
-    For a two-fluid approach see twoPhaseEulerFoam.
-
-\*---------------------------------------------------------------------------*/
-
-#include "fvCFD.H"
-#include "CMULES.H"
-#include "EulerDdtScheme.H"
-#include "localEulerDdtScheme.H"
-#include "CrankNicolsonDdtScheme.H"
-#include "subCycle.H"
-#include "immiscibleIncompressibleTwoPhaseMixture.H"
-#include "turbulentTransportModel.H"
-#include "fvcSmooth.H"
-#include "pimpleControl.H"
-#include "fvIOoptionList.H"
-#include "CorrectPhi.H"
-#include "fixedFluxPressureFvPatchScalarField.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-int main(int argc, char *argv[])
-{
-    #include "setRootCase.H"
-    #include "createTime.H"
-    #include "createMesh.H"
-
-    pimpleControl pimple(mesh);
-
-    #include "initContinuityErrs.H"
-    #include "createFields.H"
-    #include "createMRF.H"
-    #include "createFvOptions.H"
-    #include "correctPhi.H"
-    #include "CourantNo.H"
-    #include "createRDeltaT.H"
-
-    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-    Info<< "\nStarting time loop\n" << endl;
-
-    while (runTime.run())
-    {
-        runTime++;
-
-        Info<< "Time = " << runTime.timeName() << nl << endl;
-
-        #include "setrDeltaT.H"
-
-        // --- Pressure-velocity PIMPLE corrector loop
-        while (pimple.loop())
-        {
-            #include "alphaControls.H"
-
-            #define LTSSOLVE
-            #include "alphaEqnSubCycle.H"
-            #undef LTSSOLVE
-
-            mixture.correct();
-
-            turbulence->correct();
-
-            #include "UEqn.H"
-
-            // --- Pressure corrector loop
-            while (pimple.correct())
-            {
-                #include "pEqn.H"
-            }
-        }
-
-        runTime.write();
-
-        Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
-            << "  ClockTime = " << runTime.elapsedClockTime() << " s"
-            << nl << endl;
-    }
-
-    Info<< "End\n" << endl;
-
-    return 0;
-}
-
-
-// ************************************************************************* //

applications/solvers/multiphase/interFoam/LTSInterFoam/Make/files

-LTSInterFoam.C
-
-EXE = $(FOAM_APPBIN)/LTSInterFoam

applications/solvers/multiphase/interFoam/LTSInterFoam/Make/options

-EXE_INC = \
-    -I.. \
-    -I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
-    -I$(LIB_SRC)/transportModels \
-    -I$(LIB_SRC)/transportModels/incompressible/lnInclude \
-    -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-    -I$(LIB_SRC)/transportModels/immiscibleIncompressibleTwoPhaseMixture/lnInclude \
-    -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/fvOptions/lnInclude \
-    -I$(LIB_SRC)/sampling/lnInclude
-
-EXE_LIBS = \
-    -limmiscibleIncompressibleTwoPhaseMixture \
-    -lturbulenceModels \
-    -lincompressibleTurbulenceModels \
-    -lfiniteVolume \
-    -lmeshTools \
-    -lfvOptions \
-    -lsampling

applications/solvers/multiphase/interFoam/LTSInterFoam/createRDeltaT.H

-scalar maxDeltaT
-(
-    pimple.dict().lookupOrDefault<scalar>("maxDeltaT", GREAT)
-);
-
-volScalarField rDeltaT
-(
-    IOobject
-    (
-        "rDeltaT",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh,
-    1/dimensionedScalar("maxDeltaT", dimTime, maxDeltaT),
-    zeroGradientFvPatchScalarField::typeName
-);
-
-volScalarField rSubDeltaT
-(
-    IOobject
-    (
-        "rSubDeltaT",
-        runTime.timeName(),
-        mesh
-    ),
-    mesh,
-    1/dimensionedScalar("maxDeltaT", dimTime, maxDeltaT)
-);

applications/solvers/multiphase/interFoam/alphaEqn.H

@@ -77,11 +77,15 @@
     {
         fvScalarMatrix alpha1Eqn
         (
-            #ifdef LTSSOLVE
-            fv::localEulerDdtScheme<scalar>(mesh, rDeltaT.name()).fvmDdt(alpha1)
-            #else
-            fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
-            #endif
+            (
+                LTS
+              ? fv::localEulerDdtScheme<scalar>
+                (
+                    mesh,
+                    trSubDeltaT().name()
+                ).fvmDdt(alpha1)
+              : fv::EulerDdtScheme<scalar>(mesh).fvmDdt(alpha1)
+            )
           + fv::gaussConvectionScheme<scalar>
             (
                 mesh,
@@ -104,11 +108,7 @@
         if (alphaApplyPrevCorr && tphiAlphaCorr0.valid())
         {
             Info<< "Applying the previous iteration compression flux" << endl;
-            #ifdef LTSSOLVE
-            MULES::LTScorrect(alpha1, phiAlpha, tphiAlphaCorr0(), 1, 0);
-            #else
             MULES::correct(alpha1, phiAlpha, tphiAlphaCorr0(), 1, 0);
-            #endif
 
             phiAlpha += tphiAlphaCorr0();
         }
@@ -154,11 +154,7 @@
             tmp<surfaceScalarField> tphiAlphaCorr(tphiAlphaUn() - phiAlpha);
             volScalarField alpha10("alpha10", alpha1);
 
-            #ifdef LTSSOLVE
-            MULES::LTScorrect(alpha1, tphiAlphaUn(), tphiAlphaCorr(), 1, 0);
-            #else
             MULES::correct(alpha1, tphiAlphaUn(), tphiAlphaCorr(), 1, 0);
-            #endif
 
             // Under-relax the correction for all but the 1st corrector
             if (aCorr == 0)
@@ -175,11 +171,7 @@
         {
             phiAlpha = tphiAlphaUn;
 
-            #ifdef LTSSOLVE
-            MULES::explicitLTSSolve(alpha1, phi, phiAlpha, 1, 0);
-            #else
             MULES::explicitSolve(alpha1, phiCN, phiAlpha, 1, 0);
-            #endif
         }
 
         alpha2 = 1.0 - alpha1;

applications/solvers/multiphase/interFoam/createRDeltaT.H

+bool LTS =
+    word(mesh.ddtScheme("default"))
+ == fv::localEulerDdtScheme<scalar>::typeName;
+
+tmp<volScalarField> trDeltaT;
+tmp<volScalarField> trSubDeltaT;
+
+if (LTS)
+{
+    scalar maxDeltaT
+    (
+        pimple.dict().lookupOrDefault<scalar>("maxDeltaT", GREAT)
+    );
+
+    trDeltaT = tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "rDeltaT",
+                runTime.timeName(),
+                mesh,
+                IOobject::NO_READ,
+                IOobject::AUTO_WRITE
+            ),
+            mesh,
+            1/dimensionedScalar("maxDeltaT", dimTime, maxDeltaT),
+            zeroGradientFvPatchScalarField::typeName
+        )
+    );
+
+    trSubDeltaT = tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "rSubDeltaT",
+                runTime.timeName(),
+                mesh
+            ),
+            mesh,
+            1/dimensionedScalar("maxDeltaT", dimTime, maxDeltaT)
+        )
+    );
+}

applications/solvers/multiphase/interFoam/interDyMFoam/interDyMFoam.C

@@ -45,6 +45,8 @@ Description
 #include "fvIOoptionList.H"
 #include "CorrectPhi.H"
 #include "fixedFluxPressureFvPatchScalarField.H"
+#include "localEulerDdtScheme.H"
+#include "fvcSmooth.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -60,7 +62,6 @@ int main(int argc, char *argv[])
     #include "createFields.H"
     #include "createMRF.H"
     #include "createFvOptions.H"
-    #include "readTimeControls.H"
 
     volScalarField rAU
     (
@@ -78,8 +79,14 @@ int main(int argc, char *argv[])
 
     #include "correctPhi.H"
     #include "createUf.H"
-    #include "CourantNo.H"
-    #include "setInitialDeltaT.H"
+    #include "createRDeltaT.H"
+
+    if (!LTS)
+    {
+        #include "readTimeControls.H"
+        #include "CourantNo.H"
+        #include "setInitialDeltaT.H"
+    }
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
     Info<< "\nStarting time loop\n" << endl;
@@ -87,10 +94,17 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readControls.H"
-        #include "alphaCourantNo.H"
-        #include "CourantNo.H"
 
-        #include "setDeltaT.H"
+        if (LTS)
+        {
+            #include "setRDeltaT.H"
+        }
+        else
+        {
+            #include "CourantNo.H"
+            #include "alphaCourantNo.H"
+            #include "setDeltaT.H"
+        }
 
         runTime++;
 

applications/solvers/multiphase/interFoam/interFoam.C

@@ -49,6 +49,8 @@ Description
 #include "fvIOoptionList.H"
 #include "CorrectPhi.H"
 #include "fixedFluxPressureFvPatchScalarField.H"
+#include "localEulerDdtScheme.H"
+#include "fvcSmooth.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -64,10 +66,15 @@ int main(int argc, char *argv[])
     #include "createFields.H"
     #include "createMRF.H"
     #include "createFvOptions.H"
-    #include "readTimeControls.H"
     #include "correctPhi.H"
-    #include "CourantNo.H"
-    #include "setInitialDeltaT.H"
+    #include "createRDeltaT.H"
+
+    if (!LTS)
+    {
+        #include "readTimeControls.H"
+        #include "CourantNo.H"
+        #include "setInitialDeltaT.H"
+    }
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -76,9 +83,17 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readTimeControls.H"
-        #include "CourantNo.H"
-        #include "alphaCourantNo.H"
-        #include "setDeltaT.H"
+
+        if (LTS)
+        {
+            #include "setRDeltaT.H"
+        }
+        else
+        {
+            #include "CourantNo.H"
+            #include "alphaCourantNo.H"
+            #include "setDeltaT.H"
+        }
 
         runTime++;
 

applications/solvers/multiphase/interFoam/interMixingFoam/alphaEqns.H

@@ -72,19 +72,41 @@
         phiAlpha1 -= phiAlpha1BD;
 
         // Calculate the limiter for alpha1
-        MULES::limiter
-        (
-            allLambda,
-            1.0/runTime.deltaT().value(),
-            geometricOneField(),
-            alpha1,
-            phiAlpha1BD,
-            phiAlpha1,
-            zeroField(),
-            zeroField(),
-            1,
-            0
-        );
+        if (LTS)
+        {
+            const volScalarField& rDeltaT =
+                mesh.objectRegistry::lookupObject<volScalarField>("rSubDeltaT");
+
+            MULES::limiter
+            (
+                allLambda,
+                rDeltaT,
+                geometricOneField(),
+                alpha1,
+                phiAlpha1BD,
+                phiAlpha1,
+                zeroField(),
+                zeroField(),
+                1,
+                0
+            );
+        }
+        else
+        {
+            MULES::limiter
+            (
+                allLambda,
+                1.0/runTime.deltaT().value(),
+                geometricOneField(),
+                alpha1,
+                phiAlpha1BD,
+                phiAlpha1,
+                zeroField(),
+                zeroField(),
+                1,
+                0
+            );
+        }
 
         // Create the complete flux for alpha2
         surfaceScalarField phiAlpha2
@@ -113,19 +135,41 @@
         phiAlpha2 -= phiAlpha2BD;
 
         // Further limit the limiter for alpha2
-        MULES::limiter
-        (
-            allLambda,
-            1.0/runTime.deltaT().value(),
-            geometricOneField(),
-            alpha2,
-            phiAlpha2BD,
-            phiAlpha2,
-            zeroField(),
-            zeroField(),
-            1,
-            0
-        );
+        if (LTS)
+        {
+            const volScalarField& rDeltaT =
+                mesh.objectRegistry::lookupObject<volScalarField>("rSubDeltaT");
+
+            MULES::limiter
+            (
+                allLambda,
+                rDeltaT,
+                geometricOneField(),
+                alpha2,
+                phiAlpha2BD,
+                phiAlpha2,
+                zeroField(),
+                zeroField(),
+                1,
+                0
+            );
+        }
+        else
+        {
+            MULES::limiter
+            (
+                allLambda,
+                1.0/runTime.deltaT().value(),
+                geometricOneField(),
+                alpha2,
+                phiAlpha2BD,
+                phiAlpha2,
+                zeroField(),
+                zeroField(),
+                1,
+                0
+            );
+        }
 
         // Construct the limited fluxes
         phiAlpha1 = phiAlpha1BD + lambda*phiAlpha1;