TUT: added two unsteady adjoint optimisation tutorials

showcasing the use of compressedFullStorage and binomialCheckPointing
ENH: introduced unsteady adjoint functionality
2022-07-15 16:43:12 +03:00 · 2022-07-15 16:43:10 +03:00 · 2022-07-15 16:43:07 +03:00 · 2022-07-15 16:43:06 +03:00
382 changed files with 17244 additions and 4219 deletions
--- a/applications/solvers/combustion/fireFoam/createClouds.H
+++ b/applications/solvers/combustion/fireFoam/createClouds.H
@ -1,9 +1,9 @@
 Info<< "\nConstructing reacting cloud" << endl;
-reactingCloud parcels
+basicReactingCloud parcels
 (
    "reactingCloud1",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );
--- a/applications/solvers/combustion/fireFoam/fireFoam.C
+++ b/applications/solvers/combustion/fireFoam/fireFoam.C
@ -37,7 +37,7 @@ Description

 #include "fvCFD.H"
 #include "turbulentFluidThermoModel.H"
-#include "reactingCloud.H"
+#include "basicReactingCloud.H"
 #include "surfaceFilmModel.H"
 #include "pyrolysisModelCollection.H"
 #include "radiationModel.H"
--- a/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/DPMDyMFoam.C
+++ b/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/DPMDyMFoam.C
@ -42,11 +42,11 @@ Description
 #include "CorrectPhi.H"

 #ifdef MPPIC
-    #include "kinematicCloud.H"
-    #define kinematicTypeCloud kinematicCloud
+    #include "basicKinematicCloud.H"
+    #define basicKinematicTypeCloud basicKinematicCloud
 #else
-    #include "kinematicCollidingCloud.H"
-    #define kinematicTypeCloud kinematicCollidingCloud
+    #include "basicKinematicCollidingCloud.H"
+    #define basicKinematicTypeCloud basicKinematicCollidingCloud
 #endif

 int main(int argc, char *argv[])
@ -88,7 +88,7 @@ int main(int argc, char *argv[])
        Info<< "Time = " << runTime.timeName() << nl << endl;

        // Store the particle positions
-        kCloud.storeGlobalPositions();
+        kinematicCloud.storeGlobalPositions();

        mesh.update();

@ -111,15 +111,15 @@ int main(int argc, char *argv[])
        continuousPhaseTransport.correct();
        muc = rhoc*continuousPhaseTransport.nu();

-        kCloud.evolve();
+        kinematicCloud.evolve();

        // Update continuous phase volume fraction field
-        alphac = max(1.0 - kCloud.theta(), alphacMin);
+        alphac = max(1.0 - kinematicCloud.theta(), alphacMin);
        alphac.correctBoundaryConditions();
        alphacf = fvc::interpolate(alphac);
        alphaPhic = alphacf*phic;

-        fvVectorMatrix cloudSU(kCloud.SU(Uc));
+        fvVectorMatrix cloudSU(kinematicCloud.SU(Uc));
        volVectorField cloudVolSUSu
        (
            IOobject
--- a/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/MPPICDyMFoam/Make/options
+++ b/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/MPPICDyMFoam/Make/options
@ -12,8 +12,6 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
    -I$(LIB_SRC)/transportModels \
    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
--- a/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/Make/options
+++ b/applications/solvers/lagrangian/DPMFoam/DPMDyMFoam/Make/options
@ -11,11 +11,8 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
    -I$(LIB_SRC)/transportModels \
    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/phaseIncompressible/lnInclude \
--- a/applications/solvers/lagrangian/DPMFoam/DPMFoam.C
+++ b/applications/solvers/lagrangian/DPMFoam/DPMFoam.C
@ -43,11 +43,11 @@ Description
 #include "pimpleControl.H"

 #ifdef MPPIC
-    #include "kinematicCloud.H"
-    #define kinematicTypeCloud kinematicCloud
+    #include "basicKinematicCloud.H"
+    #define basicKinematicTypeCloud basicKinematicCloud
 #else
-    #include "kinematicCollidingCloud.H"
-    #define kinematicTypeCloud kinematicCollidingCloud
+    #include "basicKinematicCollidingCloud.H"
+    #define basicKinematicTypeCloud basicKinematicCollidingCloud
 #endif

 int main(int argc, char *argv[])
@ -91,16 +91,16 @@ int main(int argc, char *argv[])
        continuousPhaseTransport.correct();
        muc = rhoc*continuousPhaseTransport.nu();

-        Info<< "Evolving " << kCloud.name() << endl;
-        kCloud.evolve();
+        Info<< "Evolving " << kinematicCloud.name() << endl;
+        kinematicCloud.evolve();

        // Update continuous phase volume fraction field
-        alphac = max(1.0 - kCloud.theta(), alphacMin);
+        alphac = max(1.0 - kinematicCloud.theta(), alphacMin);
        alphac.correctBoundaryConditions();
        alphacf = fvc::interpolate(alphac);
        alphaPhic = alphacf*phic;

-        fvVectorMatrix cloudSU(kCloud.SU(Uc));
+        fvVectorMatrix cloudSU(kinematicCloud.SU(Uc));
        volVectorField cloudVolSUSu
        (
            IOobject
--- a/applications/solvers/lagrangian/DPMFoam/MPPICFoam/Make/options
+++ b/applications/solvers/lagrangian/DPMFoam/MPPICFoam/Make/options
@ -10,8 +10,6 @@ EXE_INC = \
    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
    -I$(LIB_SRC)/transportModels \
    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
--- a/applications/solvers/lagrangian/DPMFoam/Make/options
+++ b/applications/solvers/lagrangian/DPMFoam/Make/options
@ -10,11 +10,8 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
    -I$(LIB_SRC)/transportModels \
    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
    -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
--- a/applications/solvers/lagrangian/DPMFoam/createFields.H
+++ b/applications/solvers/lagrangian/DPMFoam/createFields.H
@ -125,13 +125,13 @@ const word kinematicCloudName
 );

 Info<< "Constructing kinematicCloud " << kinematicCloudName << endl;
-kinematicTypeCloud kCloud
+basicKinematicTypeCloud kinematicCloud
 (
    kinematicCloudName,
-    g,
    rhoc,
    Uc,
-    muc
+    muc,
+    g
 );

 // Particle fraction upper limit
@ -139,13 +139,13 @@ scalar alphacMin
 (
    1.0
  - (
-        kCloud.particleProperties().subDict("constantProperties")
+        kinematicCloud.particleProperties().subDict("constantProperties")
       .get<scalar>("alphaMax")
    )
 );

 // Update alphac from the particle locations
-alphac = max(1.0 - kCloud.theta(), alphacMin);
+alphac = max(1.0 - kinematicCloud.theta(), alphacMin);
 alphac.correctBoundaryConditions();

 surfaceScalarField alphacf("alphacf", fvc::interpolate(alphac));
--- a/applications/solvers/lagrangian/coalChemistryFoam/Make/options
+++ b/applications/solvers/lagrangian/coalChemistryFoam/Make/options
@ -6,6 +6,7 @@ EXE_INC = \
    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
    -I$(LIB_SRC)/lagrangian/basic/lnInclude \
    -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
+    -I$(LIB_SRC)/lagrangian/coalCombustion/lnInclude \
    -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
@ -34,6 +35,7 @@ EXE_LIBS = \
    -llagrangian \
    -llagrangianIntermediate \
    -llagrangianTurbulence \
+    -lcoalCombustion\
    -lspecie \
    -lcompressibleTransportModels \
    -lfluidThermophysicalModels \
--- a/applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C
+++ b/applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C
@ -37,8 +37,8 @@ Description

 #include "fvCFD.H"
 #include "turbulentFluidThermoModel.H"
-#include "thermoCloud.H"
-#include "reactingMultiphaseCloud.H"
+#include "basicThermoCloud.H"
+#include "coalCloud.H"
 #include "psiReactionThermo.H"
 #include "CombustionModel.H"
 #include "fvOptions.H"
--- a/applications/solvers/lagrangian/coalChemistryFoam/createClouds.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/createClouds.H
@ -1,19 +1,19 @@
 Info<< "\nConstructing coal cloud" << endl;
-reactingMultiphaseCloud coalParcels
+coalCloud coalParcels
 (
    "coalCloud1",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );

 Info<< "\nConstructing limestone cloud" << endl;
-thermoCloud limestoneParcels
+basicThermoCloud limestoneParcels
 (
    "limestoneCloud1",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );
--- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/Make/options
+++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/Make/options
@ -9,8 +9,6 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
    -I$(LIB_SRC)/transportModels \
--- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/createFields.H
@ -63,13 +63,13 @@ const word kinematicCloudName
 );

 Info<< "Constructing kinematicCloud " << kinematicCloudName << endl;
-kinematicCollidingCloud kCloud
+basicKinematicCollidingCloud kinematicCloud
 (
    kinematicCloudName,
-    g,
    rhoInf,
    U,
-    mu
+    mu,
+    g
 );

 IOobject Hheader
--- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/Make/options
+++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/Make/options
@ -10,8 +10,6 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
    -I$(LIB_SRC)/transportModels \
--- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/icoUncoupledKinematicParcelDyMFoam.C
+++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelDyMFoam/icoUncoupledKinematicParcelDyMFoam.C
@ -41,7 +41,7 @@ Description
 #include "dynamicFvMesh.H"
 #include "singlePhaseTransportModel.H"
 #include "turbulentTransportModel.H"
-#include "kinematicCollidingCloud.H"
+#include "basicKinematicCollidingCloud.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

@ -76,19 +76,19 @@ int main(int argc, char *argv[])
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;

-        kCloud.storeGlobalPositions();
+        kinematicCloud.storeGlobalPositions();

        mesh.update();

        U.correctBoundaryConditions();

-        Info<< "Evolving " << kCloud.name() << endl;
+        Info<< "Evolving " << kinematicCloud.name() << endl;

        laminarTransport.correct();

        mu = laminarTransport.nu()*rhoInfValue;

-        kCloud.evolve();
+        kinematicCloud.evolve();

        runTime.write();

--- a/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelFoam.C
+++ b/applications/solvers/lagrangian/icoUncoupledKinematicParcelFoam/icoUncoupledKinematicParcelFoam.C
@ -40,7 +40,7 @@ Description
 #include "fvCFD.H"
 #include "singlePhaseTransportModel.H"
 #include "turbulentTransportModel.H"
-#include "kinematicCollidingCloud.H"
+#include "basicKinematicCollidingCloud.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

@ -75,13 +75,13 @@ int main(int argc, char *argv[])
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;

-        Info<< "Evolving " << kCloud.name() << endl;
+        Info<< "Evolving " << kinematicCloud.name() << endl;

        laminarTransport.correct();

        mu = laminarTransport.nu()*rhoInfValue;

-        kCloud.evolve();
+        kinematicCloud.evolve();

        runTime.write();

--- a/applications/solvers/lagrangian/kinematicParcelFoam/Make/options
+++ b/applications/solvers/lagrangian/kinematicParcelFoam/Make/options
@ -8,12 +8,6 @@ EXE_INC = \
    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
    -I$(LIB_SRC)/transportModels \
    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
    -I$(LIB_SRC)/dynamicMesh/lnInclude \
    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
    -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
--- a/applications/solvers/lagrangian/kinematicParcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/kinematicParcelFoam/UEqn.H
@ -6,7 +6,7 @@
      + MRF.DDt(U)
      + turbulence->divDevReff(U)
     ==
-        scalar(1)/rhoInfValue*parcels.SU(U)
+        parcels.SU(U, true)
      + fvOptions(U)
    );

--- a/applications/solvers/lagrangian/kinematicParcelFoam/createClouds.H
+++ b/applications/solvers/lagrangian/kinematicParcelFoam/createClouds.H
@ -4,12 +4,13 @@ const word kinematicCloudName
 );

 Info<< "Constructing kinematicCloud " << kinematicCloudName << endl;
-kinematicCloud parcels
+
+basicKinematicCloud parcels
 (
    kinematicCloudName,
-    g,
    rhoInf,
    U,
-    muc
+    muc,
+    g
 );

--- a/applications/solvers/lagrangian/kinematicParcelFoam/kinematicParcelFoam.C
+++ b/applications/solvers/lagrangian/kinematicParcelFoam/kinematicParcelFoam.C
@ -40,7 +40,7 @@ Description
 #include "singlePhaseTransportModel.H"
 #include "turbulentTransportModel.H"
 #include "surfaceFilmModel.H"
-#include "kinematicCloud.H"
+#include "basicKinematicCloud.H"
 #include "fvOptions.H"
 #include "pimpleControl.H"
 #include "CorrectPhi.H"
@ -67,6 +67,7 @@ int main(int argc, char *argv[])
    #include "createDyMControls.H"
    #include "createFields.H"
    #include "createFieldRefs.H"
+    #include "createRegionControls.H"
    #include "createUfIfPresent.H"

    turbulence->validate();
@ -94,7 +95,7 @@ int main(int argc, char *argv[])
        // Do any mesh changes
        mesh.update();

-        if (pimple.solveFlow() && mesh.changing())
+        if (solvePrimaryRegion && mesh.changing())
        {
            MRF.update();

@ -119,7 +120,7 @@ int main(int argc, char *argv[])
        parcels.evolve();
        surfaceFilm.evolve();

-        if (pimple.solveFlow())
+        if (solvePrimaryRegion)
        {
            // --- PIMPLE loop
            while (pimple.loop())
--- a/applications/solvers/lagrangian/parcelFoam/Make/files
+++ b/applications/solvers/lagrangian/parcelFoam/Make/files
@ -1,3 +0,0 @@
-parcelFoam.C
-
-EXE = $(FOAM_USER_APPBIN)/parcelFoam
--- a/applications/solvers/lagrangian/parcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/parcelFoam/UEqn.H
@ -1,22 +0,0 @@
-    MRF.correctBoundaryVelocity(U);
-
-    fvVectorMatrix UEqn
-    (
-        fvm::ddt(U) + fvm::div(phi, U)
-      + MRF.DDt(U)
-      + turbulence->divDevReff(U)
-     ==
-        scalar(1)/rhoInfValue*parcels.SU(U)
-      + fvOptions(U)
-    );
-
-    UEqn.relax();
-
-    fvOptions.constrain(UEqn);
-
-    if (pimple.momentumPredictor())
-    {
-        solve(UEqn == -fvc::grad(p));
-
-        fvOptions.correct(U);
-    }
--- a/applications/solvers/lagrangian/parcelFoam/createClouds.H
+++ b/applications/solvers/lagrangian/parcelFoam/createClouds.H
@ -1,2 +0,0 @@
-auto parcels = parcelCloudModelList(g, rhoInf, U, muc);
-
--- a/applications/solvers/lagrangian/parcelFoam/createFieldRefs.H
+++ b/applications/solvers/lagrangian/parcelFoam/createFieldRefs.H
@ -1 +0,0 @@
-regionModels::surfaceFilmModel& surfaceFilm = tsurfaceFilm();
--- a/applications/solvers/lagrangian/parcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/parcelFoam/createFields.H
@ -1,85 +0,0 @@
-#include "readGravitationalAcceleration.H"
-
-Info<< "Reading field p\n" << endl;
-volScalarField p
-(
-    IOobject
-    (
-        "p",
-        runTime.timeName(),
-        mesh,
-        IOobject::MUST_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh
-);
-
-Info<< "\nReading field U\n" << endl;
-volVectorField U
-(
-    IOobject
-    (
-        "U",
-        runTime.timeName(),
-        mesh,
-        IOobject::MUST_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh
-);
-
-#include "createPhi.H"
-
-singlePhaseTransportModel laminarTransport(U, phi);
-
-dimensionedScalar rhoInfValue
-(
-    "rhoInf",
-    dimDensity,
-    laminarTransport
-);
-
-volScalarField rhoInf
-(
-    IOobject
-    (
-        "rho",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::NO_WRITE
-    ),
-    mesh,
-    rhoInfValue
-);
-
-volScalarField muc
-(
-    IOobject
-    (
-        "muc",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::NO_WRITE
-    ),
-    rhoInf*laminarTransport.nu()
-);
-
-Info<< "Creating turbulence model\n" << endl;
-autoPtr<incompressible::turbulenceModel> turbulence
-(
-    incompressible::turbulenceModel::New(U, phi, laminarTransport)
-);
-
-label pRefCell = 0;
-scalar pRefValue = 0.0;
-setRefCell(p, pimple.dict(), pRefCell, pRefValue);
-mesh.setFluxRequired(p.name());
-
-#include "createMRF.H"
-#include "createClouds.H"
-#include "createSurfaceFilmModel.H"
-#include "createFvOptions.H"
-
-
--- a/applications/solvers/lagrangian/parcelFoam/pEqn.H
+++ b/applications/solvers/lagrangian/parcelFoam/pEqn.H
@ -1,58 +0,0 @@
-
-volScalarField rAU(1.0/UEqn.A());
-volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
-
-surfaceScalarField phiHbyA("phiHbyA", fvc::flux(HbyA));
-
-if (pimple.ddtCorr())
-{
-    phiHbyA += MRF.zeroFilter(fvc::interpolate(rAU)*fvc::ddtCorr(U, phi, Uf));
-}
-
-MRF.makeRelative(phiHbyA);
-
-if (p.needReference())
-{
-    fvc::makeRelative(phiHbyA, U);
-    adjustPhi(phiHbyA, U, p);
-    fvc::makeAbsolute(phiHbyA, U);
-}
-
-
-// Update the pressure BCs to ensure flux consistency
-constrainPressure(p, U, phiHbyA, rAU, MRF);
-
-// Non-orthogonal pressure corrector loop
-while (pimple.correctNonOrthogonal())
-{
-    fvScalarMatrix pEqn
-    (
-        fvm::laplacian(rAU, p)
-     ==
-        fvc::div(phiHbyA)
-    );
-
-    pEqn.setReference(pRefCell, pRefValue);
-
-    pEqn.solve(mesh.solver(p.select(pimple.finalInnerIter())));
-
-    if (pimple.finalNonOrthogonalIter())
-    {
-        phi = phiHbyA - pEqn.flux();
-    }
-}
-
-#include "continuityErrs.H"
-
-p.relax();
-
-U = HbyA - rAU*fvc::grad(p);
-U.correctBoundaryConditions();
-fvOptions.correct(U);
-
-// Correct rhoUf if the mesh is moving
-fvc::correctUf(Uf, U, phi);
-
-// Make the fluxes relative to the mesh motion
-fvc::makeRelative(phi, U);
-
--- a/applications/solvers/lagrangian/parcelFoam/parcelFoam.C
+++ b/applications/solvers/lagrangian/parcelFoam/parcelFoam.C
@ -1,154 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | www.openfoam.com
-     \\/     M anipulation  |
-------------------------------------------------------------------------------
-    Copyright (C) 2021 OpenCFD Ltd.
-------------------------------------------------------------------------------
-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
-    parcelFoam
-
-Group
-    grpLagrangianSolvers
-
-Description
-    Transient solver for incompressible, turbulent flow with kinematic,
-    particle cloud, and surface film modelling.
-
-\*---------------------------------------------------------------------------*/
-
-#include "fvCFD.H"
-#include "dynamicFvMesh.H"
-#include "singlePhaseTransportModel.H"
-#include "turbulentTransportModel.H"
-#include "surfaceFilmModel.H"
-#include "parcelCloudModelList.H"
-#include "fvOptions.H"
-#include "pimpleControl.H"
-#include "CorrectPhi.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-int main(int argc, char *argv[])
-{
-    argList::addNote
-    (
-        "Transient solver for incompressible, turbulent flow"
-        " with kinematic particle clouds"
-        " and surface film modelling."
-    );
-
-    #define CREATE_MESH createMeshesPostProcess.H
-    #include "postProcess.H"
-
-    #include "addCheckCaseOptions.H"
-    #include "setRootCaseLists.H"
-    #include "createTime.H"
-    #include "createDynamicFvMesh.H"
-    #include "initContinuityErrs.H"
-    #include "createDyMControls.H"
-    #include "createFields.H"
-    #include "createFieldRefs.H"
-    #include "createUfIfPresent.H"
-
-    turbulence->validate();
-
-    #include "CourantNo.H"
-    #include "setInitialDeltaT.H"
-
-    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-    Info<< "\nStarting time loop\n" << endl;
-
-    while (runTime.run())
-    {
-        #include "readDyMControls.H"
-        #include "CourantNo.H"
-        #include "setMultiRegionDeltaT.H"
-
-        ++runTime;
-
-        Info<< "Time = " << runTime.timeName() << nl << endl;
-
-        // Store the particle positions
-        parcels.storeGlobalPositions();
-
-        // Do any mesh changes
-        mesh.update();
-
-        if (pimple.solveFlow() && mesh.changing())
-        {
-            MRF.update();
-
-            if (correctPhi)
-            {
-                // Calculate absolute flux
-                // from the mapped surface velocity
-                phi = mesh.Sf() & Uf();
-
-                #include "../../incompressible/pimpleFoam/correctPhi.H"
-
-                // Make the fluxes relative to the mesh-motion
-                fvc::makeRelative(phi, U);
-            }
-
-            if (checkMeshCourantNo)
-            {
-                #include "meshCourantNo.H"
-            }
-        }
-
-        parcels.evolve();
-        surfaceFilm.evolve();
-
-        if (pimple.solveFlow())
-        {
-            // --- PIMPLE loop
-            while (pimple.loop())
-            {
-                #include "UEqn.H"
-
-                // --- Pressure corrector loop
-                while (pimple.correct())
-                {
-                    #include "pEqn.H"
-                }
-
-                if (pimple.turbCorr())
-                {
-                    laminarTransport.correct();
-                    turbulence->correct();
-                }
-            }
-        }
-
-        runTime.write();
-
-        runTime.printExecutionTime(Info);
-    }
-
-    Info<< "End\n" << endl;
-
-    return 0;
-}
-
-
-// ************************************************************************* //
--- a/applications/solvers/lagrangian/reactingParcelFoam/createClouds.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/createClouds.H
@ -1,9 +1,9 @@
 Info<< "\nConstructing reacting cloud" << endl;
-reactingTypeCloud parcels
+basicReactingTypeCloud parcels
 (
    "reactingCloud1",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );
--- a/applications/solvers/lagrangian/reactingParcelFoam/reactingHeterogenousParcelFoam/reactingHeterogenousParcelFoam.C
+++ b/applications/solvers/lagrangian/reactingParcelFoam/reactingHeterogenousParcelFoam/reactingHeterogenousParcelFoam.C
@ -37,8 +37,8 @@ Description

 \*---------------------------------------------------------------------------*/

-#define CLOUD_BASE_TYPE heterogeneousReacting
-#define CLOUD_BASE_TYPE_NAME "heterogeneousReacting"
+#define CLOUD_BASE_TYPE HeterogeneousReacting
+#define CLOUD_BASE_TYPE_NAME "HeterogeneousReacting"

 #include "reactingParcelFoam.C"

--- a/applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C
+++ b/applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C
@ -53,12 +53,12 @@ Description
 #include "cloudMacros.H"

 #ifndef CLOUD_BASE_TYPE
-    #define CLOUD_BASE_TYPE reactingMultiphase
+    #define CLOUD_BASE_TYPE ReactingMultiphase
    #define CLOUD_BASE_TYPE_NAME "reacting"
 #endif

 #include CLOUD_INCLUDE_FILE(CLOUD_BASE_TYPE)
-#define reactingTypeCloud CLOUD_TYPE(CLOUD_BASE_TYPE)
+#define basicReactingTypeCloud CLOUD_TYPE(CLOUD_BASE_TYPE)

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

@ -81,6 +81,7 @@ int main(int argc, char *argv[])
    #include "createDyMControls.H"
    #include "createFields.H"
    #include "createFieldRefs.H"
+    #include "createRegionControls.H"
    #include "initContinuityErrs.H"
    #include "createRhoUfIfPresent.H"

@ -104,7 +105,7 @@ int main(int argc, char *argv[])
        // so that it can be mapped and used in correctPhi
        // to ensure the corrected phi has the same divergence
        autoPtr<volScalarField> divrhoU;
-        if (pimple.solveFlow() && correctPhi)
+        if (solvePrimaryRegion && correctPhi)
        {
            divrhoU.reset
            (
@ -132,7 +133,7 @@ int main(int argc, char *argv[])

        // Store momentum to set rhoUf for introduced faces.
        autoPtr<volVectorField> rhoU;
-        if (pimple.solveFlow() && rhoUf.valid())
+        if (solvePrimaryRegion && rhoUf.valid())
        {
            rhoU.reset(new volVectorField("rhoU", rho*U));
        }
@ -143,7 +144,7 @@ int main(int argc, char *argv[])
        // Do any mesh changes
        mesh.update();

-        if (pimple.solveFlow() && mesh.changing())
+        if (solvePrimaryRegion && mesh.changing())
        {
            gh = (g & mesh.C()) - ghRef;
            ghf = (g & mesh.Cf()) - ghRef;
@ -171,7 +172,7 @@ int main(int argc, char *argv[])
        parcels.evolve();
        surfaceFilm.evolve();

-        if (pimple.solveFlow())
+        if (solvePrimaryRegion)
        {
            if (pimple.nCorrPIMPLE() <= 1)
            {
--- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/Make/options
+++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/Make/options
@ -7,6 +7,7 @@ EXE_INC = \
    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
    -I$(LIB_SRC)/lagrangian/basic/lnInclude \
    -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
+    -I$(LIB_SRC)/lagrangian/coalCombustion/lnInclude \
    -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
--- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/createClouds.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/createClouds.H
@ -1,9 +1,9 @@
 Info<< "\nConstructing " << CLOUD_BASE_TYPE_NAME << " cloud" << endl;
-reactingTypeCloud parcels
+basicReactingTypeCloud parcels
 (
    word(CLOUD_BASE_TYPE_NAME) & "Cloud1",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );
--- a/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/simpleReactingParcelFoam.C
+++ b/applications/solvers/lagrangian/reactingParcelFoam/simpleReactingParcelFoam/simpleReactingParcelFoam.C
@ -45,13 +45,13 @@ Description
 #include "cloudMacros.H"

 #ifndef CLOUD_BASE_TYPE
-    #define CLOUD_BASE_TYPE reactingMultiphase
+    #define CLOUD_BASE_TYPE ReactingMultiphase
    //#define CLOUD_BASE_TYPE_NAME "reactingMultiphase" Backwards compat
    #define CLOUD_BASE_TYPE_NAME "reacting"
 #endif

 #include CLOUD_INCLUDE_FILE(CLOUD_BASE_TYPE)
-#define reactingTypeCloud CLOUD_TYPE(CLOUD_BASE_TYPE)
+#define basicReactingTypeCloud CLOUD_TYPE(CLOUD_BASE_TYPE)

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

--- a/applications/solvers/lagrangian/rhoParcelFoam/EEqn.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/EEqn.H
@ -1,41 +0,0 @@
-{
-    volScalarField& he = thermo.he();
-
-    fvScalarMatrix EEqn
-    (
-        fvm::ddt(rho, he) + mvConvection->fvmDiv(phi, he)
-      + fvc::ddt(rho, K) + fvc::div(phi, K)
-      + (
-            he.name() == "e"
-          ? fvc::div
-            (
-                fvc::absolute(phi/fvc::interpolate(rho), U),
-                p,
-                "div(phiv,p)"
-            )
-          : -dpdt
-        )
-      - fvm::laplacian(turbulence->alphaEff(), he)
-     ==
-        rho*(U&g)
-      + parcels.Sh(he)
-      + surfaceFilm.Sh()
-      + radiation->Sh(thermo, he)
-      + Qdot
-      + fvOptions(rho, he)
-    );
-
-    EEqn.relax();
-
-    fvOptions.constrain(EEqn);
-
-    EEqn.solve();
-
-    fvOptions.correct(he);
-
-    thermo.correct();
-    radiation->correct();
-
-    Info<< "T gas min/max   = " << min(T).value() << ", "
-        << max(T).value() << endl;
-}
--- a/applications/solvers/lagrangian/rhoParcelFoam/Make/files
+++ b/applications/solvers/lagrangian/rhoParcelFoam/Make/files
@ -1,3 +0,0 @@
-rhoParcelFoam.C
-
-EXE = $(FOAM_USER_APPBIN)/rhoParcelFoam
--- a/applications/solvers/lagrangian/rhoParcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/UEqn.H
@ -1,34 +0,0 @@
-    MRF.correctBoundaryVelocity(U);
-
-    fvVectorMatrix UEqn
-    (
-        fvm::ddt(rho, U) + fvm::div(phi, U)
-      + MRF.DDt(rho, U)
-      + turbulence->divDevRhoReff(U)
-     ==
-        parcels.SU(U)
-      + fvOptions(rho, U)
-    );
-
-    UEqn.relax();
-
-    fvOptions.constrain(UEqn);
-
-    if (pimple.momentumPredictor())
-    {
-        solve
-        (
-            UEqn
-          ==
-            fvc::reconstruct
-            (
-                (
-                  - ghf*fvc::snGrad(rho)
-                  - fvc::snGrad(p_rgh)
-                )*mesh.magSf()
-            )
-        );
-
-        fvOptions.correct(U);
-        K = 0.5*magSqr(U);
-    }
--- a/applications/solvers/lagrangian/rhoParcelFoam/YEqn.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/YEqn.H
@ -1,51 +0,0 @@
-tmp<fv::convectionScheme<scalar>> mvConvection
-(
-    fv::convectionScheme<scalar>::New
-    (
-        mesh,
-        fields,
-        phi,
-        mesh.divScheme("div(phi,Yi_h)")
-    )
-);
-
-
-{
-    combustion->correct();
-    Qdot = combustion->Qdot();
-    volScalarField Yt(0.0*Y[0]);
-
-    forAll(Y, i)
-    {
-        if (i != inertIndex && composition.active(i))
-        {
-            volScalarField& Yi = Y[i];
-
-            fvScalarMatrix YEqn
-            (
-                fvm::ddt(rho, Yi)
-              + mvConvection->fvmDiv(phi, Yi)
-              - fvm::laplacian(turbulence->muEff(), Yi)
-              ==
-                parcels.SYi(i, Yi)
-              + fvOptions(rho, Yi)
-              + combustion->R(Yi)
-              + surfaceFilm.Srho(i)
-            );
-
-            YEqn.relax();
-
-            fvOptions.constrain(YEqn);
-
-            YEqn.solve(mesh.solver("Yi"));
-
-            fvOptions.correct(Yi);
-
-            Yi.max(0.0);
-            Yt += Yi;
-        }
-    }
-
-    Y[inertIndex] = scalar(1) - Yt;
-    Y[inertIndex].max(0.0);
-}
--- a/applications/solvers/lagrangian/rhoParcelFoam/createClouds.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/createClouds.H
@ -1,3 +0,0 @@
-Info<< "\nConstructing cloud" << endl;
-
-auto parcels = parcelCloudModelList(g, rho, U, slgThermo);
--- a/applications/solvers/lagrangian/rhoParcelFoam/createFieldRefs.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/createFieldRefs.H
@ -1,5 +0,0 @@
-const volScalarField& T = thermo.T();
-const volScalarField& psi = thermo.psi();
-const label inertIndex(composition.species().find(inertSpecie));
-
-regionModels::surfaceFilmModel& surfaceFilm = tsurfaceFilm();
--- a/applications/solvers/lagrangian/rhoParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/createFields.H
@ -1,132 +0,0 @@
-#include "createRDeltaT.H"
-
-Info<< "Reading thermophysical properties\n" << endl;
-autoPtr<rhoReactionThermo> pThermo(rhoReactionThermo::New(mesh));
-rhoReactionThermo& thermo = pThermo();
-thermo.validate(args.executable(), "h", "e");
-
-SLGThermo slgThermo(mesh, thermo);
-
-basicSpecieMixture& composition = thermo.composition();
-PtrList<volScalarField>& Y = composition.Y();
-
-const word inertSpecie(thermo.get<word>("inertSpecie"));
-if
-(
-    !composition.species().found(inertSpecie)
-  && composition.species().size() > 0
-)
-{
-    FatalIOErrorIn(args.executable().c_str(), thermo)
-        << "Inert specie " << inertSpecie << " not found in available species "
-        << composition.species()
-        << exit(FatalIOError);
-}
-
-Info<< "Creating field rho\n" << endl;
-volScalarField rho
-(
-    IOobject
-    (
-        "rho",
-        runTime.timeName(),
-        mesh,
-        IOobject::NO_READ,
-        IOobject::AUTO_WRITE
-    ),
-    thermo.rho()
-);
-
-volScalarField& p = thermo.p();
-
-Info<< "\nReading field U\n" << endl;
-volVectorField U
-(
-    IOobject
-    (
-        "U",
-        runTime.timeName(),
-        mesh,
-        IOobject::MUST_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh
-);
-
-#include "compressibleCreatePhi.H"
-
-Info<< "Creating turbulence model\n" << endl;
-autoPtr<compressible::turbulenceModel> turbulence
-(
-    compressible::turbulenceModel::New
-    (
-        rho,
-        U,
-        phi,
-        thermo
-    )
-);
-
-Info<< "Creating combustion model\n" << endl;
-autoPtr<CombustionModel<rhoReactionThermo>> combustion
-(
-    CombustionModel<rhoReactionThermo>::New(thermo, turbulence())
-);
-
-#include "readGravitationalAcceleration.H"
-#include "readhRef.H"
-#include "gh.H"
-
-
-volScalarField p_rgh
-(
-    IOobject
-    (
-        "p_rgh",
-        runTime.timeName(),
-        mesh,
-        IOobject::MUST_READ,
-        IOobject::AUTO_WRITE
-    ),
-    mesh
-);
-
-// Force p_rgh to be consistent with p
-p_rgh = p - rho*gh;
-
-pressureControl pressureControl(p, rho, pimple.dict(), false);
-
-mesh.setFluxRequired(p_rgh.name());
-
-Info<< "Creating multi-variate interpolation scheme\n" << endl;
-multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
-
-forAll(Y, i)
-{
-    fields.add(Y[i]);
-}
-fields.add(thermo.he());
-
-volScalarField Qdot
-(
-    IOobject
-    (
-        "Qdot",
-        runTime.timeName(),
-        mesh,
-        IOobject::READ_IF_PRESENT,
-        IOobject::AUTO_WRITE
-    ),
-    mesh,
-    dimensionedScalar(dimEnergy/dimVolume/dimTime, Zero)
-);
-
-#include "createDpdt.H"
-
-#include "createK.H"
-
-#include "createMRF.H"
-#include "createRadiationModel.H"
-#include "createClouds.H"
-#include "createSurfaceFilmModel.H"
-#include "createFvOptions.H"
--- a/applications/solvers/lagrangian/rhoParcelFoam/createMeshesPostProcess.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/createMeshesPostProcess.H
@ -1,35 +0,0 @@
-#include "createMesh.H"
-
-dictionary filmDict;
-
-IOobject io
-(
-    "surfaceFilmProperties",
-    mesh.time().constant(),
-    mesh,
-    IOobject::READ_IF_PRESENT,
-    IOobject::NO_WRITE,
-    false
-);
-
-if (io.typeHeaderOk<IOdictionary>())
-{
-    IOdictionary propDict(io);
-
-    filmDict = std::move(propDict);
-
-    const word filmRegionName = filmDict.get<word>("region");
-
-    fvMesh filmMesh
-    (
-        IOobject
-        (
-            filmRegionName,
-            runTime.timeName(),
-            runTime,
-            IOobject::MUST_READ
-        )
-    );
-}
-
-
--- a/applications/solvers/lagrangian/rhoParcelFoam/createSurfaceFilmModel.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/createSurfaceFilmModel.H
@ -1,6 +0,0 @@
-Info<< "\nConstructing surface film model" << endl;
-
-autoPtr<regionModels::surfaceFilmModel> tsurfaceFilm
-(
-    regionModels::surfaceFilmModel::New(mesh, g)
-);
--- a/applications/solvers/lagrangian/rhoParcelFoam/pEqn.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/pEqn.H
@ -1,96 +0,0 @@
-if (!pimple.SIMPLErho())
-{
-    rho = thermo.rho();
-}
-
-// Thermodynamic density needs to be updated by psi*d(p) after the
-// pressure solution
-const volScalarField psip0(psi*p);
-
-volScalarField rAU(1.0/UEqn.A());
-surfaceScalarField rhorAUf("rhorAUf", fvc::interpolate(rho*rAU));
-volVectorField HbyA(constrainHbyA(rAU*UEqn.H(), U, p));
-
-surfaceScalarField phig(-rhorAUf*ghf*fvc::snGrad(rho)*mesh.magSf());
-
-surfaceScalarField phiHbyA
-(
-    "phiHbyA",
-    (
-        fvc::flux(rho*HbyA)
-      + MRF.zeroFilter(rhorAUf*fvc::ddtCorr(rho, U, phi))
-    )
-  + phig
-);
-
-fvc::makeRelative(phiHbyA, rho, U);
-MRF.makeRelative(fvc::interpolate(rho), phiHbyA);
-
-// Update the pressure BCs to ensure flux consistency
-constrainPressure(p_rgh, rho, U, phiHbyA, rhorAUf, MRF);
-
-fvScalarMatrix p_rghDDtEqn
-(
-    fvc::ddt(rho) + psi*correction(fvm::ddt(p_rgh))
-  + fvc::div(phiHbyA)
- ==
-    parcels.Srho()
-  + surfaceFilm.Srho()
-  + fvOptions(psi, p_rgh, rho.name())
-);
-
-while (pimple.correctNonOrthogonal())
-{
-    fvScalarMatrix p_rghEqn
-    (
-        p_rghDDtEqn
-      - fvm::laplacian(rhorAUf, p_rgh)
-    );
-
-    p_rghEqn.solve(mesh.solver(p_rgh.select(pimple.finalInnerIter())));
-
-    if (pimple.finalNonOrthogonalIter())
-    {
-        phi = phiHbyA + p_rghEqn.flux();
-
-        // Explicitly relax pressure for momentum corrector
-        p_rgh.relax();
-
-        U = HbyA + rAU*fvc::reconstruct((p_rghEqn.flux() + phig)/rhorAUf);
-        U.correctBoundaryConditions();
-        fvOptions.correct(U);
-        K = 0.5*magSqr(U);
-    }
-}
-
-p = p_rgh + rho*gh;
-
-// Thermodynamic density update
-thermo.correctRho(psi*p - psip0);
-
-#include "rhoEqn.H"
-#include "compressibleContinuityErrs.H"
-
-if (pressureControl.limit(p))
-{
-    p.correctBoundaryConditions();
-    rho = thermo.rho();
-    p_rgh = p - rho*gh;
-}
-else if (pimple.SIMPLErho())
-{
-    rho = thermo.rho();
-}
-
-// Correct rhoUf if the mesh is moving
-fvc::correctRhoUf(rhoUf, rho, U, phi);
-
-if (thermo.dpdt())
-{
-    dpdt = fvc::ddt(p);
-
-    if (mesh.moving())
-    {
-        dpdt -= fvc::div(fvc::meshPhi(rho, U), p);
-    }
-}
--- a/applications/solvers/lagrangian/rhoParcelFoam/rhoEqn.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/rhoEqn.H
@ -1,50 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | www.openfoam.com
-     \\/     M anipulation  |
-------------------------------------------------------------------------------
-    Copyright (C) 2011-2017 OpenFOAM Foundation
-------------------------------------------------------------------------------
-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
-    rhoEqn
-
-Description
-    Solve the continuity for density.
-
-\*---------------------------------------------------------------------------*/
-
-{
-    fvScalarMatrix rhoEqn
-    (
-        fvm::ddt(rho)
-      + fvc::div(phi)
-      ==
-        parcels.Srho(rho)
-      + surfaceFilm.Srho()
-      + fvOptions(rho)
-    );
-
-    rhoEqn.solve();
-
-    fvOptions.correct(rho);
-}
-
-// ************************************************************************* //
--- a/applications/solvers/lagrangian/rhoParcelFoam/rhoParcelFoam.C
+++ b/applications/solvers/lagrangian/rhoParcelFoam/rhoParcelFoam.C
@ -1,206 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | www.openfoam.com
-     \\/     M anipulation  |
-------------------------------------------------------------------------------
-    Copyright (C) 2011-2020 OpenFOAM Foundation
-    Copyright (C) 2018-2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
-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
-    rhoParcelFoam
-
-Group
-    grpLagrangianSolvers
-
-Description
-    Transient solver for compressible, turbulent flow with a reacting,
-    multiphase particle cloud, and surface film modelling.
-
-\*---------------------------------------------------------------------------*/
-
-#include "fvCFD.H"
-#include "dynamicFvMesh.H"
-#include "turbulentFluidThermoModel.H"
-#include "surfaceFilmModel.H"
-#include "rhoReactionThermo.H"
-#include "CombustionModel.H"
-#include "radiationModel.H"
-#include "SLGThermo.H"
-#include "fvOptions.H"
-#include "pimpleControl.H"
-#include "pressureControl.H"
-#include "CorrectPhi.H"
-#include "localEulerDdtScheme.H"
-#include "fvcSmooth.H"
-
-#include "parcelCloudModelList.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-int main(int argc, char *argv[])
-{
-    argList::addNote
-    (
-        "Transient solver for compressible, turbulent flow"
-        " with lagrangian parcel clouds and surface film modelling."
-    );
-
-    #define CREATE_MESH createMeshesPostProcess.H
-    #include "postProcess.H"
-
-    #include "addCheckCaseOptions.H"
-    #include "setRootCaseLists.H"
-    #include "createTime.H"
-    #include "createDynamicFvMesh.H"
-    #include "createDyMControls.H"
-    #include "createFields.H"
-    #include "createFieldRefs.H"
-    #include "initContinuityErrs.H"
-    #include "createRhoUfIfPresent.H"
-
-    turbulence->validate();
-
-    if (!LTS)
-    {
-        #include "compressibleCourantNo.H"
-        #include "setInitialDeltaT.H"
-    }
-
-    // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-    Info<< "\nStarting time loop\n" << endl;
-
-    while (runTime.run())
-    {
-        #include "readDyMControls.H"
-
-        // Store divrhoU from the previous mesh
-        // so that it can be mapped and used in correctPhi
-        // to ensure the corrected phi has the same divergence
-        autoPtr<volScalarField> divrhoU;
-        if (pimple.solveFlow() && correctPhi)
-        {
-            divrhoU.reset
-            (
-                new volScalarField
-                (
-                    "divrhoU",
-                    fvc::div(fvc::absolute(phi, rho, U))
-                )
-            );
-        }
-
-        if (LTS)
-        {
-            #include "setRDeltaT.H"
-        }
-        else
-        {
-            #include "compressibleCourantNo.H"
-            #include "setMultiRegionDeltaT.H"
-        }
-
-        ++runTime;
-
-        Info<< "Time = " << runTime.timeName() << nl << endl;
-
-        // Store momentum to set rhoUf for introduced faces.
-        autoPtr<volVectorField> rhoU;
-        if (pimple.solveFlow() && rhoUf.valid())
-        {
-            rhoU.reset(new volVectorField("rhoU", rho*U));
-        }
-
-        // Store the particle positions
-        parcels.storeGlobalPositions();
-
-        // Do any mesh changes
-        mesh.update();
-
-        if (pimple.solveFlow() && mesh.changing())
-        {
-            gh = (g & mesh.C()) - ghRef;
-            ghf = (g & mesh.Cf()) - ghRef;
-
-            MRF.update();
-
-            if (correctPhi)
-            {
-                // Calculate absolute flux
-                // from the mapped surface velocity
-                phi = mesh.Sf() & rhoUf();
-
-                #include "../../compressible/rhoPimpleFoam/correctPhi.H"
-
-                // Make the fluxes relative to the mesh-motion
-                fvc::makeRelative(phi, rho, U);
-            }
-
-            if (checkMeshCourantNo)
-            {
-                #include "meshCourantNo.H"
-            }
-        }
-
-        parcels.evolve();
-        surfaceFilm.evolve();
-
-        if (pimple.solveFlow())
-        {
-            if (pimple.nCorrPIMPLE() <= 1)
-            {
-                #include "rhoEqn.H"
-            }
-
-            // --- PIMPLE loop
-            while (pimple.loop())
-            {
-                #include "UEqn.H"
-                #include "YEqn.H"
-                #include "EEqn.H"
-
-                // --- Pressure corrector loop
-                while (pimple.correct())
-                {
-                    #include "pEqn.H"
-                }
-
-                if (pimple.turbCorr())
-                {
-                    turbulence->correct();
-                }
-            }
-
-            rho = thermo.rho();
-        }
-
-        runTime.write();
-
-        runTime.printExecutionTime(Info);
-    }
-
-    Info<< "End\n" << endl;
-
-    return 0;
-}
-
-
-// ************************************************************************* //
--- a/applications/solvers/lagrangian/rhoParcelFoam/setMultiRegionDeltaT.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/setMultiRegionDeltaT.H
@ -1,55 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | www.openfoam.com
-     \\/     M anipulation  |
-------------------------------------------------------------------------------
-    Copyright (C) 2011-2017 OpenFOAM Foundation
-------------------------------------------------------------------------------
-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
-    setMultiRegionDeltaT
-
-Description
-    Reset the timestep to maintain a constant maximum Courant numbers.
-    Reduction of time-step is immediate, but increase is damped to avoid
-    unstable oscillations.
-
-\*---------------------------------------------------------------------------*/
-
-if (adjustTimeStep)
-{
-    const scalar maxDeltaTFact =
-        min(maxCo/(CoNum + SMALL), maxCo/(surfaceFilm.CourantNumber() + SMALL));
-    const scalar deltaTFact =
-        min(min(maxDeltaTFact, 1.0 + 0.1*maxDeltaTFact), 1.2);
-
-    runTime.setDeltaT
-    (
-        min
-        (
-            deltaTFact*runTime.deltaTValue(),
-            maxDeltaT
-        )
-    );
-
-    Info<< "deltaT = " <<  runTime.deltaTValue() << endl;
-}
-
-// ************************************************************************* //
--- a/applications/solvers/lagrangian/rhoParcelFoam/setRDeltaT.H
+++ b/applications/solvers/lagrangian/rhoParcelFoam/setRDeltaT.H
@ -1,137 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | www.openfoam.com
-     \\/     M anipulation  |
-------------------------------------------------------------------------------
-    Copyright (C) 2011-2016 OpenFOAM Foundation
-    Copyright (C) 2020 OpenCFD Ltd.
-------------------------------------------------------------------------------
-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/>.
-
-\*---------------------------------------------------------------------------*/
-
-{
-    volScalarField& rDeltaT = trDeltaT.ref();
-
-    const dictionary& pimpleDict = pimple.dict();
-
-    // Maximum flow Courant number
-    scalar maxCo(pimpleDict.get<scalar>("maxCo"));
-
-    // Maximum time scale
-    scalar maxDeltaT(pimpleDict.getOrDefault<scalar>("maxDeltaT", GREAT));
-
-    // Smoothing parameter (0-1) when smoothing iterations > 0
-    scalar rDeltaTSmoothingCoeff
-    (
-        pimpleDict.getOrDefault<scalar>("rDeltaTSmoothingCoeff", 0.1)
-    );
-
-    // Damping coefficient (1-0)
-    scalar rDeltaTDampingCoeff
-    (
-        pimpleDict.getOrDefault<scalar>("rDeltaTDampingCoeff", 0.2)
-    );
-
-    // Maximum change in cell temperature per iteration
-    // (relative to previous value)
-    scalar alphaTemp(pimpleDict.getOrDefault("alphaTemp", 0.05));
-
-
-    Info<< "Time scales min/max:" << endl;
-
-    // Cache old reciprocal time scale field
-    volScalarField rDeltaT0("rDeltaT0", rDeltaT);
-
-    // Flow time scale
-    {
-        rDeltaT.ref() =
-        (
-            fvc::surfaceSum(mag(phi))()()
-           /((2*maxCo)*mesh.V()*rho())
-        );
-
-        // Limit the largest time scale
-        rDeltaT.max(1/maxDeltaT);
-
-        Info<< "    Flow        = "
-            << gMin(1/rDeltaT.primitiveField()) << ", "
-            << gMax(1/rDeltaT.primitiveField()) << endl;
-    }
-
-    // Reaction source time scale
-    {
-        volScalarField::Internal rDeltaTT
-        (
-            mag
-            (
-                parcels.hsTrans()/(mesh.V()*runTime.deltaT())
-              + Qdot()
-            )
-           /(
-               alphaTemp
-              *rho()
-              *thermo.Cp()()()
-              *T()
-           )
-        );
-
-        Info<< "    Temperature = "
-            << gMin(1/(rDeltaTT.field() + VSMALL)) << ", "
-            << gMax(1/(rDeltaTT.field() + VSMALL)) << endl;
-
-        rDeltaT.ref() = max
-        (
-            rDeltaT(),
-            rDeltaTT
-        );
-    }
-
-    // Update the boundary values of the reciprocal time-step
-    rDeltaT.correctBoundaryConditions();
-
-    // Spatially smooth the time scale field
-    if (rDeltaTSmoothingCoeff < 1.0)
-    {
-        fvc::smooth(rDeltaT, rDeltaTSmoothingCoeff);
-    }
-
-    // Limit rate of change of time scale
-    // - reduce as much as required
-    // - only increase at a fraction of old time scale
-    if
-    (
-        rDeltaTDampingCoeff < 1.0
-     && runTime.timeIndex() > runTime.startTimeIndex() + 1
-    )
-    {
-        rDeltaT = max
-        (
-            rDeltaT,
-            (scalar(1) - rDeltaTDampingCoeff)*rDeltaT0
-        );
-    }
-
-    Info<< "    Overall     = "
-        << gMin(1/rDeltaT.primitiveField())
-        << ", " << gMax(1/rDeltaT.primitiveField()) << endl;
-}
-
-
-// ************************************************************************* //
--- a/applications/solvers/lagrangian/simpleCoalParcelFoam/Make/options
+++ b/applications/solvers/lagrangian/simpleCoalParcelFoam/Make/options
@ -7,6 +7,7 @@ EXE_INC = \
    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
    -I$(LIB_SRC)/lagrangian/basic/lnInclude \
    -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
+    -I$(LIB_SRC)/lagrangian/coalCombustion/lnInclude \
    -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
@ -46,6 +47,7 @@ EXE_LIBS = \
    -lsurfaceFilmModels \
    -lcombustionModels \
    -lsampling \
+    -lcoalCombustion \
    -lregionFaModels \
    -lfiniteArea \
    -lfaOptions
--- a/applications/solvers/lagrangian/simpleCoalParcelFoam/createClouds.H
+++ b/applications/solvers/lagrangian/simpleCoalParcelFoam/createClouds.H
@ -1,9 +1,9 @@
 Info<< "\nConstructing coal cloud" << endl;
-reactingMultiphaseCloud parcels
+coalCloud parcels
 (
    "reactingCloud1",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );
--- a/applications/solvers/lagrangian/simpleCoalParcelFoam/simpleCoalParcelFoam.C
+++ b/applications/solvers/lagrangian/simpleCoalParcelFoam/simpleCoalParcelFoam.C
@ -37,7 +37,7 @@ Description

 #include "fvCFD.H"
 #include "turbulentFluidThermoModel.H"
-#include "reactingMultiphaseCloud.H"
+#include "coalCloud.H"
 #include "rhoReactionThermo.H"
 #include "CombustionModel.H"
 #include "radiationModel.H"
--- a/applications/solvers/lagrangian/sprayFoam/Make/options
+++ b/applications/solvers/lagrangian/sprayFoam/Make/options
@ -35,6 +35,7 @@ EXE_LIBS = \
    -llagrangian \
    -llagrangianIntermediate \
    -llagrangianTurbulence \
+    -llagrangianSpray \
    -lspecie \
    -lcompressibleTransportModels \
    -lfluidThermophysicalModels \
--- a/applications/solvers/lagrangian/sprayFoam/createClouds.H
+++ b/applications/solvers/lagrangian/sprayFoam/createClouds.H
@ -1,9 +1,9 @@
 Info<< "\nConstructing reacting cloud" << endl;
-sprayCloud parcels
+basicSprayCloud parcels
 (
    "sprayCloud",
-    g,
    rho,
    U,
+    g,
    slgThermo
 );
--- a/applications/solvers/lagrangian/sprayFoam/engineFoam/Make/options
+++ b/applications/solvers/lagrangian/sprayFoam/engineFoam/Make/options
@ -39,6 +39,7 @@ EXE_LIBS = \
    -llagrangian \
    -llagrangianIntermediate \
    -llagrangianTurbulence \
+    -llagrangianSpray \
    -lspecie \
    -lcompressibleTransportModels \
    -lfluidThermophysicalModels \
--- a/applications/solvers/lagrangian/sprayFoam/engineFoam/engineFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/engineFoam/engineFoam.C
@ -37,7 +37,7 @@ Description
 #include "engineTime.H"
 #include "engineMesh.H"
 #include "turbulentFluidThermoModel.H"
-#include "sprayCloud.H"
+#include "basicSprayCloud.H"
 #include "psiReactionThermo.H"
 #include "CombustionModel.H"
 #include "radiationModel.H"
--- a/applications/solvers/lagrangian/sprayFoam/simpleSprayFoam/Make/options
+++ b/applications/solvers/lagrangian/sprayFoam/simpleSprayFoam/Make/options
@ -34,6 +34,7 @@ EXE_LIBS = \
    -lthermoTools \
    -llagrangian \
    -llagrangianIntermediate \
+    -llagrangianSpray \
    -llagrangianTurbulence \
    -lspecie \
    -lcompressibleTransportModels \
--- a/applications/solvers/lagrangian/sprayFoam/simpleSprayFoam/simpleSprayFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/simpleSprayFoam/simpleSprayFoam.C
@ -35,7 +35,7 @@ Description

 \*---------------------------------------------------------------------------*/

-#define CLOUD_BASE_TYPE spray
+#define CLOUD_BASE_TYPE Spray
 #define CLOUD_BASE_TYPE_NAME "spray"

 #include "simpleReactingParcelFoam.C"
--- a/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/Make/options
+++ b/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/Make/options
@ -40,6 +40,7 @@ EXE_LIBS = \
    -llagrangian \
    -llagrangianIntermediate \
    -llagrangianTurbulence \
+    -llagrangianSpray \
    -lspecie \
    -lcompressibleTransportModels \
    -lfluidThermophysicalModels \
--- a/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/sprayDyMFoam/sprayDyMFoam.C
@ -38,7 +38,7 @@ Description
 #include "fvCFD.H"
 #include "dynamicFvMesh.H"
 #include "turbulenceModel.H"
-#include "sprayCloud.H"
+#include "basicSprayCloud.H"
 #include "psiReactionThermo.H"
 #include "CombustionModel.H"
 #include "radiationModel.H"
--- a/applications/solvers/lagrangian/sprayFoam/sprayFoam.C
+++ b/applications/solvers/lagrangian/sprayFoam/sprayFoam.C
@ -37,7 +37,7 @@ Description

 #include "fvCFD.H"
 #include "turbulentFluidThermoModel.H"
-#include "sprayCloud.H"
+#include "basicSprayCloud.H"
 #include "psiReactionThermo.H"
 #include "CombustionModel.H"
 #include "radiationModel.H"
--- a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/Make/options
+++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/Make/options
@ -9,8 +9,6 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
    -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
--- a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/createFields.H
@ -55,11 +55,11 @@ const word kinematicCloudName
 );

 Info<< "Constructing kinematicCloud " << kinematicCloudName << endl;
-kinematicCloud kCloud
+basicKinematicCloud kinematicCloud
 (
    kinematicCloudName,
-    g,
    rho,
    U,
-    thermo.mu()
+    thermo.mu(),
+    g
 );
--- a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/options
+++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/Make/options
@ -10,8 +10,6 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
    -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
--- a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/uncoupledKinematicParcelDyMFoam.C
+++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelDyMFoam/uncoupledKinematicParcelDyMFoam.C
@ -37,7 +37,7 @@ Description
 #include "dynamicFvMesh.H"
 #include "psiThermo.H"
 #include "turbulentFluidThermoModel.H"
-#include "kinematicCloud.H"
+#include "basicKinematicCloud.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

@ -73,14 +73,14 @@ int main(int argc, char *argv[])
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;

-        kCloud.storeGlobalPositions();
+        kinematicCloud.storeGlobalPositions();

        mesh.update();

        U.correctBoundaryConditions();

-        Info<< "Evolving " << kCloud.name() << endl;
-        kCloud.evolve();
+        Info<< "Evolving " << kinematicCloud.name() << endl;
+        kinematicCloud.evolve();

        runTime.write();

--- a/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelFoam.C
+++ b/applications/solvers/lagrangian/uncoupledKinematicParcelFoam/uncoupledKinematicParcelFoam.C
@ -39,7 +39,7 @@ Description
 #include "fvCFD.H"
 #include "psiThermo.H"
 #include "turbulentFluidThermoModel.H"
-#include "kinematicCloud.H"
+#include "basicKinematicCloud.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

@ -75,8 +75,8 @@ int main(int argc, char *argv[])
    {
        Info<< "Time = " << runTime.timeName() << nl << endl;

-        Info<< "Evolving " << kCloud.name() << endl;
-        kCloud.evolve();
+        Info<< "Evolving " << kinematicCloud.name() << endl;
+        kinematicCloud.evolve();

        runTime.write();

--- a/applications/solvers/multiphase/MPPICInterFoam/MPPICInterFoam.C
+++ b/applications/solvers/multiphase/MPPICInterFoam/MPPICInterFoam.C
@ -52,7 +52,7 @@ Description
 #include "CorrectPhi.H"
 #include "fvcSmooth.H"

-#include "kinematicCloud.H"
+#include "basicKinematicCloud.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

@ -111,12 +111,12 @@ int main(int argc, char *argv[])

        Info<< "Time = " << runTime.timeName() << nl << endl;

-        Info<< "Evolving " << kCloud.name() << endl;
+        Info<< "Evolving " << kinematicCloud.name() << endl;

-        kCloud.evolve();
+        kinematicCloud.evolve();

        // Update continuous phase volume fraction field
-        alphac = max(1.0 - kCloud.theta(), alphacMin);
+        alphac = max(1.0 - kinematicCloud.theta(), alphacMin);
        alphac.correctBoundaryConditions();

         Info<< "Continuous phase-1 volume fraction = "
@ -130,7 +130,7 @@ int main(int argc, char *argv[])
        alphaPhic = alphacf*phi;
        alphacRho = alphac*rho;

-        fvVectorMatrix cloudSU(kCloud.SU(U));
+        fvVectorMatrix cloudSU(kinematicCloud.SU(U));
        volVectorField cloudVolSUSu
        (
            IOobject
--- a/applications/solvers/multiphase/MPPICInterFoam/Make/options
+++ b/applications/solvers/multiphase/MPPICInterFoam/Make/options
@ -12,14 +12,11 @@ EXE_INC = \
    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
    -I$(LIB_SRC)/transportModels \
    -I$(LIB_SRC)/transportModels/immiscibleIncompressibleTwoPhaseMixture/lnInclude \
    -I$(LIB_SRC)/transportModels/twoPhaseMixture/lnInclude \
    -I$(LIB_SRC)/transportModels/incompressible/lnInclude \
    -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
-    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
    -I$(LIB_SRC)/TurbulenceModels/phaseCompressible/lnInclude \
--- a/applications/solvers/multiphase/MPPICInterFoam/createFields.H
+++ b/applications/solvers/multiphase/MPPICInterFoam/createFields.H
@ -147,13 +147,13 @@ volScalarField alphacRho(alphac*rho);
 alphacRho.oldTime();

 Info<< "Constructing kinematicCloud " << endl;
-kinematicCloud kCloud
+basicKinematicCloud kinematicCloud
 (
    "kinematicCloud",
-    g,
    rho,
    U,
-    mu
+    mu,
+    g
 );

 // Particle fraction upper limit
@ -161,13 +161,13 @@ scalar alphacMin
 (
    1.0
  - (
-        kCloud.particleProperties().subDict("constantProperties")
+        kinematicCloud.particleProperties().subDict("constantProperties")
       .get<scalar>("alphaMax")
    )
 );

 // Update alphac from the particle locations
-alphac = max(1.0 - kCloud.theta(), alphacMin);
+alphac = max(1.0 - kinematicCloud.theta(), alphacMin);
 alphac.correctBoundaryConditions();

 surfaceScalarField alphacf("alphacf", fvc::interpolate(alphac));
--- a/src/OpenFOAM/db/Time/Time.C
+++ b/src/OpenFOAM/db/Time/Time.C
@ -967,6 +967,34 @@ bool Foam::Time::end() const
 }


+bool Foam::Time::reverseRun() const
+{
+    // Must not solve for the 0 time step
+    bool isRunning = value() > (endTime_ + 0.5*deltaT_);
+
+    return isRunning;
+}
+
+
+bool Foam::Time::reverseLoop()
+{
+    const bool isRunning = reverseRun();
+
+    if (isRunning)
+    {
+        operator--();
+    }
+
+    return isRunning;
+}
+
+
+bool Foam::Time::reverseEnd() const
+{
+    return value() < (endTime_ - 0.5*deltaT_);
+}
+
+
 bool Foam::Time::stopAt(const stopAtControls stopCtrl) const
 {
    if (stopCtrl == stopAtControls::saUnknown)
@ -1365,4 +1393,242 @@ Foam::Time& Foam::Time::operator++(int)
 }


+Foam::Time&
+Foam::Time::operator-=(const dimensionedScalar& deltaT)
+{
+    return operator-=(deltaT.value());
+}
+
+
+Foam::Time& Foam::Time::operator-=(const scalar deltaT)
+{
+    setDeltaT(deltaT);
+    return operator--();
+}
+
+
+Foam::Time& Foam::Time::operator--()
+{
+    deltaT0_ = deltaTSave_;
+    deltaTSave_ = deltaT_;
+
+    // Save old time value and name
+    const scalar oldTimeValue = timeToUserTime(value());
+    const word oldTimeName = dimensionedScalar::name();
+
+    // Decrease time
+    setTime(value() - deltaT_, timeIndex_ - 1);
+
+    if (!subCycling_)
+    {
+        // If the time is very close to zero reset to zero
+        if (mag(value()) < 10*SMALL*deltaT_)
+        {
+            setTime(0.0, timeIndex_);
+        }
+
+        if (sigStopAtWriteNow_.active() || sigWriteNow_.active())
+        {
+            // A signal might have been sent on one processor only
+            // Reduce so all decide the same.
+
+            label flag = 0;
+            if (sigStopAtWriteNow_.active() && stopAt_ == saWriteNow)
+            {
+                flag += 1;
+            }
+            if (sigWriteNow_.active() && writeOnce_)
+            {
+                flag += 2;
+            }
+            reduce(flag, maxOp<label>());
+
+            if (flag & 1)
+            {
+                stopAt_ = saWriteNow;
+            }
+            if (flag & 2)
+            {
+                writeOnce_ = true;
+            }
+        }
+
+        writeTime_ = false;
+
+        switch (writeControl_)
+        {
+            case wcNone:
+            case wcUnknown:
+            break;
+
+            case wcTimeStep:
+                // Avoid writing just because timIndex_ is zero
+                writeTime_ = 
+                    !(timeIndex_ % label(writeInterval_)) && timeIndex_;
+            break;
+
+            case wcRunTime:
+            case wcAdjustableRunTime:
+            {
+                const label writeIndex = label
+                (
+                    ((value() - startTime_) - 0.5*deltaT_)
+                  / writeInterval_
+                );
+
+                if (writeIndex < writeTimeIndex_)
+                {
+                    writeTime_ = true;
+                    writeTimeIndex_ = writeIndex;
+                }
+            }
+            break;
+
+            case wcCpuTime:
+            {
+                const label writeIndex = label
+                (
+                    returnReduce(elapsedCpuTime(), maxOp<double>())
+                  / writeInterval_
+                );
+                if (writeIndex > writeTimeIndex_)
+                {
+                    writeTime_ = true;
+                    writeTimeIndex_ = writeIndex;
+                }
+            }
+            break;
+
+            case wcClockTime:
+            {
+                const label writeIndex = label
+                (
+                    returnReduce(elapsedClockTime(), maxOp<double>())
+                  / writeInterval_
+                );
+                if (writeIndex > writeTimeIndex_)
+                {
+                    writeTime_ = true;
+                    writeTimeIndex_ = writeIndex;
+                }
+            }
+            break;
+        }
+
+
+        // Check if endTime needs adjustment to stop at the next
+        // reverseRun()/reverseEnd()
+        if (!reverseEnd())
+        {
+            if (stopAt_ == saNoWriteNow)
+            {
+                endTime_ = value();
+            }
+            else if (stopAt_ == saWriteNow)
+            {
+                endTime_ = value();
+                writeTime_ = true;
+            }
+            else if (stopAt_ == saNextWrite && writeTime_ == true)
+            {
+                endTime_ = value();
+            }
+        }
+
+        // Override writeTime if one-shot writing
+        if (writeOnce_)
+        {
+            writeTime_ = true;
+            writeOnce_ = false;
+        }
+
+        // Adjust the precision of the time directory name if necessary
+        if (writeTime_)
+        {
+            // Tolerance used when testing time equivalence
+            const scalar timeTol =
+                max(min(pow(10.0, -precision_), 0.1*deltaT_), SMALL);
+
+            // User-time equivalent of deltaT
+            const scalar userDeltaT = timeToUserTime(deltaT_);
+
+            // Time value obtained by reading timeName
+            scalar timeNameValue = -VGREAT;
+
+            // Check that new time representation differs from old one
+            // reinterpretation of the word
+            if
+            (
+                readScalar(dimensionedScalar::name(), timeNameValue)
+             && (mag(-timeNameValue + oldTimeValue - userDeltaT) > timeTol)
+            )
+            {
+                int oldPrecision = precision_;
+                while
+                (
+                    precision_ < maxPrecision_
+                 && readScalar(dimensionedScalar::name(), timeNameValue)
+                 && (mag(-timeNameValue + oldTimeValue - userDeltaT) > timeTol)
+                )
+                {
+                    precision_++;
+                    setTime(value(), timeIndex());
+                }
+
+                if (precision_ != oldPrecision)
+                {
+                    WarningInFunction
+                        << "Increased the timePrecision from " << oldPrecision
+                        << " to " << precision_
+                        << " to distinguish between timeNames at time "
+                        << dimensionedScalar::name()
+                        << endl;
+
+                    if (precision_ == maxPrecision_)
+                    {
+                        // Reached maxPrecision limit
+                        WarningInFunction
+                            << "Current time name " << dimensionedScalar::name()
+                            << nl
+                            << "    The maximum time precision has been reached"
+                               " which might result in overwriting previous"
+                               " results."
+                            << endl;
+                    }
+
+                    // Check if round-off error caused time-reversal
+                    scalar oldTimeNameValue = -VGREAT;
+                    if
+                    (
+                        readScalar(oldTimeName, oldTimeNameValue)
+                     && (
+                            sign(timeNameValue - oldTimeNameValue)
+                         != sign(deltaT_)
+                        )
+                    )
+                    {
+                        WarningInFunction
+                            << "Current time name " << dimensionedScalar::name()
+                            << " is set to an instance prior to the "
+                               "previous one "
+                            << oldTimeName << nl
+                            << "    This might result in temporal "
+                               "discontinuities."
+                            << endl;
+                    }
+                }
+            }
+        }
+    }
+
+    return *this;
+}
+
+
+Foam::Time& Foam::Time::operator--(int)
+{
+    return operator--();
+}
+
+
 // ************************************************************************* //
--- a/src/OpenFOAM/db/Time/Time.H
+++ b/src/OpenFOAM/db/Time/Time.H
@ -7,6 +7,7 @@
 -------------------------------------------------------------------------------
    Copyright (C) 2011-2017 OpenFOAM Foundation
    Copyright (C) 2016-2020 OpenCFD Ltd.
+    Copyright (C) 2022      PCOpt/NTUA
 -------------------------------------------------------------------------------
 License
    This file is part of OpenFOAM.
@ -571,6 +572,43 @@ public:
            //      not yield the same result
            virtual bool end() const;

+            //- Return true if run should continue.
+            //  Does not invoke any functionObject methods
+            //  \note
+            //  For correct behaviour, the following style of time-loop
+            //  is recommended:
+            //  \code
+            //      while (runTime.reverseRun())
+            //      {
+            //          --runTime;
+            //          solve;
+            //          runTime.write();
+            //      }
+            //  \endcode
+            virtual bool reverseRun() const;
+
+            //- Return true if run should continue and if so decrease time
+            //  Does not invoke any functionObject methods
+            //  \note
+            //  For correct behaviour, the following style of time-loop
+            //  is recommended:
+            //  \code
+            //      while (runTime.reverseLoop())
+            //      {
+            //          solve;
+            //          runTime.write();
+            //      }
+            //  \endcode
+            virtual bool reverseLoop();
+
+            //- Return true if end of run,
+            //  does not invoke any functionObject methods
+            //  \note
+            //      The rounding heuristics near endTime mean that
+            //      \code reverseRun() \endcode and \code !reverseEnd() \endcode
+            //      may not yield the same result
+            virtual bool reverseEnd() const;
+

        // Edit

@ -653,6 +691,20 @@ public:

        //- Postfix increment, this is identical to the prefix increment
        virtual Time& operator++(int);
+
+        //- Set deltaT to that specified and decrease time via operator--()
+        virtual Time& operator-=(const dimensionedScalar& deltaT);
+
+        //- Set deltaT to that specified and decrease time via operator--()
+        virtual Time& operator-=(const scalar deltaT);
+
+        //- prefix decrease,
+        //  also invokes the functionobjectlist::start() or
+        //  functionobjectlist::execute() method, depending on the time-index
+        virtual Time& operator--();
+
+        //- Postfix decrease, this is identical to the decrease increment
+        virtual Time& operator--(int);
 };


--- a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.C
+++ b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.C
@ -545,7 +545,8 @@ Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
 (
    const IOobject& io,
    const Mesh& mesh,
-    const dictionary& dict
+    const dictionary& dict,
+    const bool readOldTime
 )
 :
    Internal(io, mesh, dimless, false),
@ -566,6 +567,11 @@ Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
            << exit(FatalIOError);
    }

+    if (readOldTime)
+    {
+        readOldTimeIfPresent();
+    }
+
    DebugInFunction
        << "Finishing dictionary-construct" << nl << this->info() << endl;
 }
@ -623,7 +629,8 @@ template<class Type, template<class> class PatchField, class GeoMesh>
 Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
 (
    const IOobject& io,
-    const GeometricField<Type, PatchField, GeoMesh>& gf
+    const GeometricField<Type, PatchField, GeoMesh>& gf,
+    const bool readOldTime
 )
 :
    Internal(io, gf),
@ -636,7 +643,7 @@ Foam::GeometricField<Type, PatchField, GeoMesh>::GeometricField
        << "Copy construct, resetting IO params" << nl
        << this->info() << endl;

-    if (!readIfPresent() && gf.field0Ptr_)
+    if (!readIfPresent() && gf.field0Ptr_ && readOldTime)
    {
        field0Ptr_ = new GeometricField<Type, PatchField, GeoMesh>
        (
--- a/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H
+++ b/src/OpenFOAM/fields/GeometricFields/GeometricField/GeometricField.H
@ -296,7 +296,8 @@ public:
        (
            const IOobject& io,
            const Mesh& mesh,
-            const dictionary& dict
+            const dictionary& dict,
+            const bool readOldTime = false
        );

        //- Copy construct
@ -315,7 +316,8 @@ public:
        GeometricField
        (
            const IOobject& io,
-            const GeometricField<Type, PatchField, GeoMesh>& gf
+            const GeometricField<Type, PatchField, GeoMesh>& gf,
+            const bool readOldTime = true
        );

        //- Construct from tmp\<GeometricField\> resetting IO parameters
--- a/src/functionObjects/lagrangian/Make/options
+++ b/src/functionObjects/lagrangian/Make/options
@ -12,14 +12,7 @@ EXE_INC = \
    -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
    -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
    -I$(LIB_SRC)/regionFaModels/lnInclude \
-    -I$(LIB_SRC)/faOptions/lnInclude \
-    \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
-    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude
+    -I$(LIB_SRC)/faOptions/lnInclude

 LIB_LIBS = \
    -lfiniteVolume \
--- a/src/functionObjects/lagrangian/cloudInfo/cloudInfo.C
+++ b/src/functionObjects/lagrangian/cloudInfo/cloudInfo.C
@ -28,7 +28,7 @@ License

 #include "cloudInfo.H"
 #include "cloud.H"
-#include "parcelCloudModel.H"
+#include "kinematicCloud.H"
 #include "dictionary.H"
 #include "mathematicalConstants.H"
 #include "PstreamReduceOps.H"
@ -147,7 +147,7 @@ bool Foam::functionObjects::cloudInfo::performAction(unsigned request)

        const word& cloudName = names()[cloudi];

-        const auto* kinCloudPtr = obr_.cfindObject<parcelCloudModel>(cloudName);
+        const auto* kinCloudPtr = obr_.cfindObject<kinematicCloud>(cloudName);

        if (!kinCloudPtr)
        {
--- a/src/functionObjects/lagrangian/icoUncoupledKinematicCloud/icoUncoupledKinematicCloud.C
+++ b/src/functionObjects/lagrangian/icoUncoupledKinematicCloud/icoUncoupledKinematicCloud.C
@ -92,10 +92,10 @@ Foam::functionObjects::icoUncoupledKinematicCloud::icoUncoupledKinematicCloud
    kinematicCloud_
    (
        kinematicCloudName_,
-        g_,
        rho_,
        U_,
-        mu_
+        mu_,
+        g_
    )
 {}

--- a/src/functionObjects/lagrangian/icoUncoupledKinematicCloud/icoUncoupledKinematicCloud.H
+++ b/src/functionObjects/lagrangian/icoUncoupledKinematicCloud/icoUncoupledKinematicCloud.H
@ -76,7 +76,7 @@ SourceFiles
 #define functionObjects_icoUncoupledKinematicCloud_H

 #include "fvMeshFunctionObject.H"
-#include "kinematicCollidingCloud.H"
+#include "basicKinematicCollidingCloud.H"
 #include "volFields.H"
 #include "uniformDimensionedFields.H"

@ -122,7 +122,7 @@ class icoUncoupledKinematicCloud
        word kinematicCloudName_;

        //- The kinematic cloud
-        kinematicCollidingCloud kinematicCloud_;
+        basicKinematicCollidingCloud kinematicCloud_;


    // Private member functions
--- a/src/lagrangian/Allwmake
+++ b/src/lagrangian/Allwmake
@ -9,7 +9,9 @@ wmake $targetType basic
 wmake $targetType solidParticle
 wmake $targetType intermediate
 wmake $targetType turbulence
+wmake $targetType spray
 wmake $targetType DSMC
+wmake $targetType coalCombustion

 molecularDynamics/Allwmake $targetType $*

--- a/src/lagrangian/coalCombustion/Make/files
+++ b/src/lagrangian/coalCombustion/Make/files
@ -0,0 +1,5 @@
+/* Coal parcel and sub-models */
+coalParcel/makeCoalParcelSubmodels.C
+coalCloudList/coalCloudList.C
+
+LIB = $(FOAM_LIBBIN)/libcoalCombustion
--- a/applications/solvers/lagrangian/parcelFoam/Make/options
+++ b/applications/solvers/lagrangian/parcelFoam/Make/options
@ -1,48 +1,49 @@
 EXE_INC = \
-    -I$(FOAM_SOLVERS)/lagrangian/reactingParcelFoam \
    -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/finiteArea/lnInclude \
    -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/sampling/lnInclude \
-    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
-    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
-    -I$(LIB_SRC)/transportModels \
-    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
-    -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
-    -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
-    -I$(LIB_SRC)/regionFaModels/lnInclude \
-    -I$(LIB_SRC)/faOptions/lnInclude \
    -I$(LIB_SRC)/lagrangian/basic/lnInclude \
    -I$(LIB_SRC)/lagrangian/intermediate/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
+    -I$(LIB_SRC)/lagrangian/distributionModels/lnInclude \
    -I$(LIB_SRC)/transportModels/compressible/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
    -I$(LIB_SRC)/thermophysicalModels/specie/lnInclude \
-    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude
+    -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/thermophysicalProperties/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/reactionThermo/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/SLGThermo/lnInclude \
+    -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/compressible/lnInclude \
+    -I$(LIB_SRC)/transportModels \
+    -I$(LIB_SRC)/regionModels/regionModel/lnInclude \
+    -I$(LIB_SRC)/regionModels/surfaceFilmModels/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/regionFaModels/lnInclude \
+    -I$(LIB_SRC)/finiteArea/lnInclude \
+    -I$(LIB_SRC)/faOptions/lnInclude

-
-EXE_LIBS = \
+LIB_LIBS = \
    -lfiniteVolume \
-    -lfvOptions \
    -lmeshTools \
-    -lsampling \
-    -lturbulenceModels \
-    -lincompressibleTurbulenceModels \
-    -lincompressibleTransportModels \
-    -ldynamicMesh \
-    -ldynamicFvMesh \
-    -ltopoChangerFvMesh \
-    -latmosphericModels \
-    -lregionModels \
-    -lsurfaceFilmModels \
-    -lsurfaceFilmDerivedFvPatchFields \
    -llagrangian \
    -llagrangianIntermediate \
    -llagrangianTurbulence \
+    -ldistributionModels \
+    -lspecie \
+    -lcompressibleTransportModels \
+    -lfluidThermophysicalModels \
+    -lthermophysicalProperties \
+    -lreactionThermophysicalModels \
+    -lSLGThermo \
+    -lradiationModels \
+    -lturbulenceModels \
+    -lincompressibleTurbulenceModels \
+    -lcompressibleTurbulenceModels \
+    -lincompressibleTransportModels \
+    -lregionModels \
+    -lsurfaceFilmModels \
+    -ldynamicMesh \
+    -ldynamicFvMesh \
    -lregionFaModels \
-    -lfiniteArea \
-    -lfaOptions
+    -lfiniteArea
--- a/src/lagrangian/intermediate/clouds/derived/reactingCloud/reactingCloud.H
+++ b/src/lagrangian/intermediate/clouds/derived/reactingCloud/reactingCloud.H
@ -5,7 +5,7 @@
    \\  /    A nd           | www.openfoam.com
     \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2011 OpenFOAM Foundation
+    Copyright (C) 2011-2015 OpenFOAM Foundation
 -------------------------------------------------------------------------------
 License
    This file is part of OpenFOAM.
@ -24,28 +24,28 @@ License
    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.

 Class
-    Foam::reactingCloud
+    Foam::coalCloud

 Description
-    Cloud class to introduce reacting parcels
+    Cloud class to introduce coal parcels

 \*---------------------------------------------------------------------------*/

-#ifndef Foam_reactingCloud_H
-#define Foam_reactingCloud_H
+#ifndef coalCloud_H
+#define coalCloud_H

-#include "BaseCloud.H"
+#include "Cloud.H"
 #include "KinematicCloud.H"
 #include "ThermoCloud.H"
 #include "ReactingCloud.H"
-#include "ParcelCloudModel.H"
-#include "reactingParcel.H"
+#include "ReactingMultiphaseCloud.H"
+#include "coalParcel.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

 namespace Foam
 {
-    typedef ParcelCloudModel
+    typedef ReactingMultiphaseCloud
    <
        ReactingCloud
        <
@ -53,14 +53,14 @@ namespace Foam
            <
                KinematicCloud
                <
-                    BaseCloud
+                    Cloud
                    <
-                        reactingParcel
+                        coalParcel
                    >
                >
            >
        >
-    > reactingCloud;
+    > coalCloud;
 }

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
--- a/src/lagrangian/intermediate/clouds/parcelCloudModelList/parcelCloudModelList.C
+++ b/src/lagrangian/intermediate/clouds/parcelCloudModelList/parcelCloudModelList.C
@ -26,102 +26,66 @@ License

 \*---------------------------------------------------------------------------*/

-#include "parcelCloudModelList.H"
+#include "coalCloudList.H"

-template<class ... Args>
-void Foam::parcelCloudModelList::initialise(const Args& ... args)
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::coalCloudList::coalCloudList
+(
+    const volScalarField& rho,
+    const volVectorField& U,
+    const dimensionedVector& g,
+    const SLGThermo& slgThermo
+)
+:
+    PtrList<coalCloud>(),
+    mesh_(rho.mesh())
 {
    IOdictionary props
    (
        IOobject
        (
-            "cloudProperties",
+            "coalCloudList",
            mesh_.time().constant(),
            mesh_,
            IOobject::MUST_READ
        )
    );

-    wordHashSet cloudTypes;
-    if (props.readIfPresent("cloudTypes", cloudTypes))
+    const wordHashSet cloudNames(props.get<wordList>("clouds"));
+
+    setSize(cloudNames.size());
+
+    label i = 0;
+    for (const word& name : cloudNames)
    {
-        setSize(cloudTypes.size());
+        Info<< "creating cloud: " << name << endl;

-        label i = 0;
-        for (const auto& name : cloudTypes)
-        {
-            Info<< "creating cloud: " << name << endl;
+        set
+        (
+            i++,
+            new coalCloud
+            (
+                name,
+                rho,
+                U,
+                g,
+                slgThermo
+            )
+        );

-            set(i++, parcelCloudModel::New(name, args ...));
-
-            Info<< endl;
-        }
+        Info<< endl;
    }
-    else
-    {
-        setSize(1);
-        set(0, parcelCloudModel::New("cloud", args ...));
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::parcelCloudModelList::parcelCloudModelList
-(
-    const dimensionedVector& g,
-    const volScalarField& rho,
-    const volVectorField& U,
-    const volScalarField& mu
-)
-:
-    PtrList<parcelCloudModel>(),
-    mesh_(rho.mesh())
-{
-    initialise(g, rho, U, mu);
-}
-
-
-Foam::parcelCloudModelList::parcelCloudModelList
-(
-    const dimensionedVector& g,
-    const volScalarField& rho,
-    const volVectorField& U,
-    const SLGThermo& slgThermo
-)
-:
-    PtrList<parcelCloudModel>(),
-    mesh_(rho.mesh())
-{
-    initialise(g, rho, U, slgThermo);
 }


 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

-void Foam::parcelCloudModelList::info()
+void Foam::coalCloudList::evolve()
 {
-    for (auto& c : *this)
+    forAll(*this, i)
    {
-        c.info();
-    }
-}
-
-
-void Foam::parcelCloudModelList::evolve()
-{
-    for (auto& c : *this)
-    {
-        c.evolveME();
-    }
-}
-
-
-void Foam::parcelCloudModelList::storeGlobalPositions()
-{
-    for (auto& c : *this)
-    {
-        c.storeGlobalPositionsME();
+        operator[](i).evolve();
    }
 }

--- a/src/lagrangian/intermediate/clouds/parcelCloudModelList/parcelCloudModelList.H
+++ b/src/lagrangian/intermediate/clouds/parcelCloudModelList/parcelCloudModelList.H
@ -5,8 +5,7 @@
    \\  /    A nd           | www.openfoam.com
     \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2012-2021 OpenFOAM Foundation
-    Copyright (C) 2022 OpenCFD Ltd.
+    Copyright (C) 2012-2016 OpenFOAM Foundation
 -------------------------------------------------------------------------------
 License
    This file is part of OpenFOAM.
@ -26,10 +25,10 @@ License

 \*---------------------------------------------------------------------------*/

-#ifndef Foam_parcelCloudModelList_H
-#define Foam_parcelCloudModelList_H
+#ifndef coalCloudList_H
+#define coalCloudList_H

-#include "parcelCloudModel.H"
+#include "coalCloud.H"
 #include "volFieldsFwd.H"
 #include "fvMatricesFwd.H"

@ -39,43 +38,28 @@ namespace Foam
 {

 /*---------------------------------------------------------------------------*\
-                   Class parcelCloudModelList Declaration
+                       Class coalCloudList Declaration
 \*---------------------------------------------------------------------------*/

-class parcelCloudModelList
+class coalCloudList
 :
-    public PtrList<parcelCloudModel>
+    public PtrList<coalCloud>
 {
 private:

    //- Reference to the mesh
    const fvMesh& mesh_;

-    // Private Member Functions
-
-        //- Initialise the clouds
-        template<class ... Args>
-        void initialise(const Args& ... args);
-

 public:

-    // Constructors
+    // Constructor

-        parcelCloudModelList
+        coalCloudList
        (
-            const dimensionedVector& g,
            const volScalarField& rho,
            const volVectorField& U,
-            const volScalarField& mu
-        );
-
-
-        parcelCloudModelList
-        (
            const dimensionedVector& g,
-            const volScalarField& rho,
-            const volVectorField& U,
            const SLGThermo& slgThermo
        );

@ -84,61 +68,46 @@ public:

        // Evolution

-            void info();
-
            //- Evolve the cloud collection
            void evolve();

-            void storeGlobalPositions();
-

        // Source terms

-            // Momentum
-
-            //- Return momentum transfer [kg m/s]
+            //- Return const reference to momentum source
            inline tmp<volVectorField::Internal> UTrans() const;

-            //- Momentum transfer coefficient [kg]
-            inline tmp<volScalarField::Internal> UCeoff() const;
-
-            //- Return momentum source term [kg m/s^2]
+            //- Return tmp momentum source term
            inline tmp<fvVectorMatrix> SU(volVectorField& U) const;

-
-            // Energy
-
-            //- Return sensible enthalpy transfer [J]
+            //- Sensible enthalpy transfer [J/kg]
            inline tmp<volScalarField::Internal> hsTrans() const;

-            //- Return sensible enthalpy transfer coefficient [J/K]
-            inline tmp<volScalarField::Internal> hsCoeff() const;
-
-            //- Return sensible enthalpy source term [J/s]
+            //- Return sensible enthalpy source term [J/kg/m3/s]
            inline tmp<fvScalarMatrix> Sh(volScalarField& hs) const;

-            //- Return equivalent particulate emission [kg/m/s^3]
-            inline tmp<volScalarField> Ep() const;

-            //- Return equivalent particulate absorption [1/m]
-            inline tmp<volScalarField> ap() const;
-
-            //- Return equivalent particulate scattering factor [1/m]
-            inline tmp<volScalarField> sigmap() const;
-
-
-            // Mass
-
-            //- Return mass source term for specie [kg/s]
+            //- Return mass source term for specie i - specie eqn
            inline tmp<fvScalarMatrix> SYi
            (
-                const label speciei,
+                const label i,
                volScalarField& Yi
            ) const;

-            //- Return total mass source term [kg/m^3/s]
+            //- Return total mass transfer [kg/m3]
+            inline tmp<volScalarField::Internal> rhoTrans() const;
+
+            //- Return tmp total mass source for carrier phase
+            //  - fully explicit
            inline tmp<volScalarField::Internal> Srho() const;

+            //- Return tmp total mass source for carrier phase specie i
+            //  - fully explicit
+            inline tmp<volScalarField::Internal> Srho
+            (
+                const label i
+            ) const;
+
            //- Return total mass source term [kg/m3/s]
            inline tmp<fvScalarMatrix> Srho(volScalarField& rho) const;
 };
@ -150,7 +119,7 @@ public:

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

-#include "parcelCloudModelListI.H"
+#include "coalCloudListI.H"

 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

--- a/src/lagrangian/coalCombustion/coalCloudList/coalCloudListI.H
+++ b/src/lagrangian/coalCombustion/coalCloudList/coalCloudListI.H
@ -0,0 +1,264 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2012-2016 OpenFOAM Foundation
+-------------------------------------------------------------------------------
+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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "fvMatrices.H"
+#include "volFields.H"
+#include "DimensionedField.H"
+
+Foam::tmp<Foam::DimensionedField<Foam::vector, Foam::volMesh>>
+Foam::coalCloudList::UTrans() const
+{
+    tmp<volVectorField::Internal> tfld
+    (
+        new volVectorField::Internal
+        (
+            IOobject
+            (
+                "UTransEff",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh_,
+            dimensionedVector(dimMass*dimVelocity, Zero)
+        )
+    );
+
+    volVectorField::Internal& fld = tfld.ref();
+
+    forAll(*this, i)
+    {
+        fld += operator[](i).UTrans();
+    }
+
+    return tfld;
+}
+
+
+Foam::tmp<Foam::fvVectorMatrix> Foam::coalCloudList::SU
+(
+    volVectorField& U
+) const
+{
+    tmp<fvVectorMatrix> tfvm(new fvVectorMatrix(U, dimForce));
+    fvVectorMatrix& fvm = tfvm.ref();
+
+    forAll(*this, i)
+    {
+        fvm += operator[](i).SU(U);
+    }
+
+    return tfvm;
+}
+
+
+Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh>>
+Foam::coalCloudList::hsTrans() const
+{
+    tmp<volScalarField::Internal> tfld
+    (
+        new volScalarField::Internal
+        (
+            IOobject
+            (
+                "hsTransEff",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh_,
+            dimensionedScalar(dimEnergy, Zero)
+        )
+    );
+
+    volScalarField::Internal& fld = tfld.ref();
+
+    forAll(*this, i)
+    {
+        fld += operator[](i).hsTrans();
+    }
+
+    return tfld;
+}
+
+
+Foam::tmp<Foam::fvScalarMatrix> Foam::coalCloudList::Sh
+(
+    volScalarField& hs
+) const
+{
+    tmp<fvScalarMatrix> tfvm(new fvScalarMatrix(hs, dimEnergy/dimTime));
+    fvScalarMatrix& fvm = tfvm.ref();
+
+    forAll(*this, i)
+    {
+        fvm += operator[](i).Sh(hs);
+    }
+
+    return tfvm;
+}
+
+
+Foam::tmp<Foam::fvScalarMatrix> Foam::coalCloudList::SYi
+(
+    const label ii,
+    volScalarField& Yi
+) const
+{
+    tmp<fvScalarMatrix> tfvm(new fvScalarMatrix(Yi, dimMass/dimTime));
+    fvScalarMatrix& fvm = tfvm.ref();
+
+    forAll(*this, i)
+    {
+        fvm += operator[](i).SYi(ii, Yi);
+    }
+
+    return tfvm;
+}
+
+
+Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh>>
+Foam::coalCloudList::rhoTrans() const
+{
+    tmp<volScalarField::Internal> tfld
+    (
+        new volScalarField::Internal
+        (
+            IOobject
+            (
+                "rhoTransEff",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh_,
+            dimensionedScalar(dimMass, Zero)
+        )
+    );
+
+    volScalarField::Internal& fld = tfld.ref();
+
+    forAll(*this, i)
+    {
+        forAll(operator[](i).rhoTrans(), j)
+        {
+            fld += operator[](i).rhoTrans()[j];
+        }
+    }
+
+    return tfld;
+}
+
+
+
+
+Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh>>
+Foam::coalCloudList::Srho() const
+{
+    tmp<volScalarField::Internal> tfld
+    (
+        new volScalarField::Internal
+        (
+            IOobject
+            (
+                "rhoTransEff",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh_,
+            dimensionedScalar(dimDensity/dimTime, Zero)
+        )
+    );
+
+    volScalarField::Internal& fld = tfld.ref();
+
+    forAll(*this, i)
+    {
+        fld += operator[](i).Srho();
+    }
+
+    return tfld;
+}
+
+
+Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh>>
+Foam::coalCloudList::Srho
+(
+    const label i
+) const
+{
+    tmp<volScalarField::Internal> tfld
+    (
+        new volScalarField::Internal
+        (
+            IOobject
+            (
+                "rhoTransEff",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            mesh_,
+            dimensionedScalar(dimDensity/dimTime, Zero)
+        )
+    );
+
+    volScalarField::Internal& fld = tfld.ref();
+
+    forAll(*this, j)
+    {
+        fld += operator[](j).Srho(i);
+    }
+
+    return tfld;
+}
+
+
+Foam::tmp<Foam::fvScalarMatrix> Foam::coalCloudList::Srho
+(
+    volScalarField& rho
+) const
+{
+    tmp<fvScalarMatrix> tfvm(new fvScalarMatrix(rho, dimMass/dimTime));
+    fvScalarMatrix& fvm = tfvm.ref();
+
+    forAll(*this, i)
+    {
+        fvm += operator[](i).Srho(rho);
+    }
+
+    return tfvm;
+}
+
+
+// ************************************************************************* //
--- a/src/lagrangian/coalCombustion/coalParcel/coalParcel.H
+++ b/src/lagrangian/coalCombustion/coalParcel/coalParcel.H
@ -0,0 +1,74 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011-2015 OpenFOAM Foundation
+    Copyright (C) 2019 OpenCFD Ltd.
+-------------------------------------------------------------------------------
+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/>.
+
+Class
+    Foam::coalParcel
+
+Description
+    Definition of coal parcel
+
+SourceFiles
+    CoalParcel.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef coalParcel_H
+#define coalParcel_H
+
+#include "contiguous.H"
+#include "particle.H"
+#include "KinematicParcel.H"
+#include "ThermoParcel.H"
+#include "ReactingParcel.H"
+#include "ReactingMultiphaseParcel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    typedef ReactingMultiphaseParcel
+    <
+        ReactingParcel
+        <
+            ThermoParcel
+            <
+                KinematicParcel
+                <
+                    particle
+                >
+            >
+        >
+    > coalParcel;
+
+    //- Non-contiguous data for coalParcel
+    template<> struct is_contiguous<coalParcel> : std::false_type {};
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
--- a/src/lagrangian/coalCombustion/coalParcel/makeCoalParcelSubmodels.C
+++ b/src/lagrangian/coalCombustion/coalParcel/makeCoalParcelSubmodels.C
@ -0,0 +1,39 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011 OpenFOAM Foundation
+-------------------------------------------------------------------------------
+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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "coalCloud.H"
+#include "makeCoalParcelSurfaceReactionModels.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    makeCoalParcelSurfaceReactionModels(coalCloud);
+}
+
+
+// ************************************************************************* //
--- a/src/lagrangian/coalCombustion/include/makeCoalParcelSurfaceReactionModels.H
+++ b/src/lagrangian/coalCombustion/include/makeCoalParcelSurfaceReactionModels.H
@ -0,0 +1,59 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | www.openfoam.com
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+    Copyright (C) 2011-2016 OpenFOAM Foundation
+-------------------------------------------------------------------------------
+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/>.
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef makeCoalParcelSurfaceReactionModels_H
+#define makeCoalParcelSurfaceReactionModels_H
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "NoSurfaceReaction.H"
+#include "COxidationDiffusionLimitedRate.H"
+#include "COxidationIntrinsicRate.H"
+#include "COxidationKineticDiffusionLimitedRate.H"
+#include "COxidationHurtMitchell.H"
+#include "COxidationMurphyShaddix.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#define makeCoalParcelSurfaceReactionModels(CloudType)                         \
+                                                                               \
+    makeSurfaceReactionModelType(COxidationDiffusionLimitedRate, CloudType);   \
+    makeSurfaceReactionModelType                                               \
+    (                                                                          \
+        COxidationKineticDiffusionLimitedRate,                                 \
+        CloudType                                                              \
+    );                                                                         \
+    makeSurfaceReactionModelType(COxidationIntrinsicRate, CloudType);          \
+    makeSurfaceReactionModelType(COxidationHurtMitchell, CloudType);           \
+    makeSurfaceReactionModelType(COxidationMurphyShaddix, CloudType);
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.C
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.H
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationDiffusionLimitedRate/COxidationDiffusionLimitedRate.H
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationHurtMitchell/COxidationHurtMitchell.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationHurtMitchell/COxidationHurtMitchell.C
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationHurtMitchell/COxidationHurtMitchell.H
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationHurtMitchell/COxidationHurtMitchell.H
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationIntrinsicRate/COxidationIntrinsicRate.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationIntrinsicRate/COxidationIntrinsicRate.C
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationIntrinsicRate/COxidationIntrinsicRate.H
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationIntrinsicRate/COxidationIntrinsicRate.H
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.C
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.H
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationKineticDiffusionLimitedRate/COxidationKineticDiffusionLimitedRate.H
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.C
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.H
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/SurfaceReactionModel/COxidationMurphyShaddix/COxidationMurphyShaddix.H
--- a/src/lagrangian/intermediate/Make/files
+++ b/src/lagrangian/intermediate/Make/files
@ -2,75 +2,59 @@ PARCELS=parcels
 BASEPARCELS=$(PARCELS)/baseClasses
 DERIVEDPARCELS=$(PARCELS)/derived

-
 CLOUDS=clouds
-TEMPLATECLOUDS=$(CLOUDS)/Templates
-$(TEMPLATECLOUDS)/KinematicCloud/KinematicCloudName.C
-$(TEMPLATECLOUDS)/ThermoCloud/ThermoCloudName.C
-$(TEMPLATECLOUDS)/ReactingCloud/ReactingCloudName.C
-$(TEMPLATECLOUDS)/ReactingMultiphaseCloud/ReactingMultiphaseCloudName.C
-$(TEMPLATECLOUDS)/ReactingHeterogeneousCloud/ReactingHeterogeneousCloudName.C
-$(TEMPLATECLOUDS)/CollidingCloud/CollidingCloudName.C
-$(TEMPLATECLOUDS)/MPPICCloud/MPPICCloudName.C
-$(TEMPLATECLOUDS)/SprayCloud/SprayCloudName.C
-
-
+BASECLOUDS=$(CLOUDS)/baseClasses
 DERIVEDCLOUDS=$(CLOUDS)/derived
-$(DERIVEDCLOUDS)/kinematicCloud/kinematicCloud.C
-$(DERIVEDCLOUDS)/thermoCloud/thermoCloud.C
-$(DERIVEDCLOUDS)/reactingCloud/reactingCloud.C
-$(DERIVEDCLOUDS)/reactingMultiphaseCloud/reactingMultiphaseCloud.C
-$(DERIVEDCLOUDS)/heterogeneousReactingCloud/heterogeneousReactingCloud.C
-$(DERIVEDCLOUDS)/kinematicCollidingCloud/kinematicCollidingCloud.C
-$(DERIVEDCLOUDS)/kinematicMPPICCloud/kinematicMPPICCloud.C
-$(DERIVEDCLOUDS)/sprayCloud/sprayCloud.C
+
+
+/* Cloud base classes */
+$(BASECLOUDS)/kinematicCloud/kinematicCloud.C
+$(BASECLOUDS)/thermoCloud/thermoCloud.C
+$(BASECLOUDS)/reactingCloud/reactingCloud.C
+$(BASECLOUDS)/reactingMultiphaseCloud/reactingMultiphaseCloud.C
+$(BASECLOUDS)/reactingHeterogeneousCloud/reactingHeterogeneousCloud.C


 /* kinematic parcel sub-models */
-KINEMATICPARCEL=$(DERIVEDPARCELS)/kinematicParcel
-$(KINEMATICPARCEL)/defineKinematicParcel.C
-$(KINEMATICPARCEL)/makeKinematicParcelSubmodels.C
+KINEMATICPARCEL=$(DERIVEDPARCELS)/basicKinematicParcel
+$(KINEMATICPARCEL)/defineBasicKinematicParcel.C
+$(KINEMATICPARCEL)/makeBasicKinematicParcelSubmodels.C


 /* kinematic colliding parcel sub-models */
-KINEMATICCOLLIDINGPARCEL=$(DERIVEDPARCELS)/kinematicCollidingParcel
-$(KINEMATICCOLLIDINGPARCEL)/defineKinematicCollidingParcel.C
-$(KINEMATICCOLLIDINGPARCEL)/makeKinematicCollidingParcelSubmodels.C
+KINEMATICCOLLIDINGPARCEL=$(DERIVEDPARCELS)/basicKinematicCollidingParcel
+$(KINEMATICCOLLIDINGPARCEL)/defineBasicKinematicCollidingParcel.C
+$(KINEMATICCOLLIDINGPARCEL)/makeBasicKinematicCollidingParcelSubmodels.C


 /* thermo parcel sub-models */
-THERMOPARCEL=$(DERIVEDPARCELS)/thermoParcel
-$(THERMOPARCEL)/defineThermoParcel.C
-$(THERMOPARCEL)/makeThermoParcelSubmodels.C
+THERMOPARCEL=$(DERIVEDPARCELS)/basicThermoParcel
+$(THERMOPARCEL)/defineBasicThermoParcel.C
+$(THERMOPARCEL)/makeBasicThermoParcelSubmodels.C


 /* reacting parcel sub-models */
-REACTINGPARCEL=$(DERIVEDPARCELS)/reactingParcel
-$(REACTINGPARCEL)/defineReactingParcel.C
-$(REACTINGPARCEL)/makeReactingParcelSubmodels.C
+REACTINGPARCEL=$(DERIVEDPARCELS)/basicReactingParcel
+$(REACTINGPARCEL)/defineBasicReactingParcel.C
+$(REACTINGPARCEL)/makeBasicReactingParcelSubmodels.C


 /* reacting multiphase parcel sub-models */
-REACTINGMPPARCEL=$(DERIVEDPARCELS)/reactingMultiphaseParcel
-$(REACTINGMPPARCEL)/defineReactingMultiphaseParcel.C
-$(REACTINGMPPARCEL)/makeReactingMultiphaseParcelSubmodels.C
+REACTINGMPPARCEL=$(DERIVEDPARCELS)/basicReactingMultiphaseParcel
+$(REACTINGMPPARCEL)/defineBasicReactingMultiphaseParcel.C
+$(REACTINGMPPARCEL)/makeBasicReactingMultiphaseParcelSubmodels.C


-/* heterogeneous reacting parcel sub-models */
-REACTINGHETERMPPARCEL=$(DERIVEDPARCELS)/heterogeneousReactingParcel
-$(REACTINGHETERMPPARCEL)/defineHeterogeneousReactingParcel.C
-$(REACTINGHETERMPPARCEL)/makeHeterogeneousReactingParcelSubmodels.C
+/* heterogeous reacting parcel sub-models */
+REACTINGHETERMPPARCEL=$(DERIVEDPARCELS)/basicHeterogeneousReactingParcel
+$(REACTINGHETERMPPARCEL)/defineBasicHeterogeneousReactingParcel.C
+$(REACTINGHETERMPPARCEL)/makeBasicHeterogeneousReactingParcelSubmodels.C


 /* kinematic MPPIC parcel sub-models */
-KINEMATICMPPICPARCEL=$(DERIVEDPARCELS)/kinematicMPPICParcel
-$(KINEMATICMPPICPARCEL)/defineKinematicMPPICParcel.C
-$(KINEMATICMPPICPARCEL)/makeKinematicMPPICParcelSubmodels.C
-
-/* spray parcel sub-models */
-SPRAYPARCEL=$(DERIVEDPARCELS)/sprayParcel
-$(SPRAYPARCEL)/defineSprayParcel.C
-$(SPRAYPARCEL)/makeSprayParcelSubmodels.C
+KINEMATICMPPICPARCEL=$(DERIVEDPARCELS)/basicKinematicMPPICParcel
+$(KINEMATICMPPICPARCEL)/defineBasicKinematicMPPICParcel.C
+$(KINEMATICMPPICPARCEL)/makeBasicKinematicMPPICParcelSubmodels.C


 /* bolt-on models */
@ -146,10 +130,4 @@ conversion/ensight/ensightOutputCloud.C
 conversion/vtk/foamVtkLagrangianWriter.C


-clouds/baseCloudInterface.C
-clouds/parcelCloudModel.C
-clouds/parcelCloudModelNew.C
-clouds/parcelCloudModelList/parcelCloudModelList.C
-
-
 LIB = $(FOAM_LIBBIN)/liblagrangianIntermediate
--- a/src/lagrangian/intermediate/clouds/BaseCloud.H
+++ b/src/lagrangian/intermediate/clouds/BaseCloud.H
@ -1,275 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | www.openfoam.com
-     \\/     M anipulation  |
-------------------------------------------------------------------------------
-    Copyright (C) 2020-2021 OpenFOAM Foundation
-    Copyright (C) 2022 OpenCFD Ltd.
-------------------------------------------------------------------------------
-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/>.
-
-Class
-    Foam::BaseCloud
-
-Description
-    Base template for parcel clouds. Inserts the baseCloudInterface
-    virtualisation layer into the class. Also defines default zero-return
-    source methods to enable the less functional clouds to be used in more
-    complex situations.
-
-\*---------------------------------------------------------------------------*/
-
-#ifndef Foam_BaseCloud_H
-#define Foam_BaseCloud_H
-
-#include "Cloud.H"
-#include "baseCloudInterface.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-class SLGThermo;
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-template<class ParticleType>
-class BaseCloud
-:
-    public Cloud<ParticleType>,
-    virtual public baseCloudInterface
-{
-    // Private data
-
-        //- Reference to the mesh
-        const fvMesh& mesh_;
-
-
-public:
-
-    // Constructors
-
-        BaseCloud
-        (
-            const word& cloudName,
-            const dimensionedVector& g,
-            const volScalarField& rho,
-            const volVectorField& U,
-            const volScalarField& mu,
-            const bool readFields = true
-        )
-        :
-            Cloud<ParticleType>(rho.mesh(), cloudName, false),
-            mesh_(rho.mesh())
-        {}
-
-
-        BaseCloud
-        (
-            const word& cloudName,
-            const dimensionedVector& g,
-            const volScalarField& rho,
-            const volVectorField& U,
-            const SLGThermo& thermo,
-            const bool readFields = true
-        )
-        :
-            Cloud<ParticleType>(rho.mesh(), cloudName, false),
-            mesh_(rho.mesh())
-        {}
-
-        //- Copy constructor with new name
-        BaseCloud
-        (
-            BaseCloud<ParticleType>& c,
-            const word& name
-        )
-        :
-            Cloud<ParticleType>(c.mesh(), name, c),
-            mesh_(c.mesh_)
-        {}
-
-        //- Copy constructor with new name - creates bare cloud
-        BaseCloud
-        (
-            const fvMesh& mesh,
-            const word& name,
-            const BaseCloud<ParticleType>& c
-        )
-        :
-            Cloud<ParticleType>(mesh, name, IDLList<ParticleType>()),
-            mesh_(mesh)
-        {}
-
-    //- Destructor
-    virtual ~BaseCloud() = default;
-
-
-    // Member Functions
-
-        // Access
-
-            //- Return reference to the mesh
-            const fvMesh& mesh() const
-            {
-                return mesh_;
-            }
-
-
-        // Sources
-
-            // Momentum
-
-                //- Momentum transfer [kg m/s]
-                virtual tmp<volVectorField::Internal> UTrans() const
-                {
-                    return volVectorField::Internal::New
-                    (
-                        IOobject::scopedName(this->name(), ":UTrans"),
-                        mesh_,
-                        dimensionedVector(dimMass*dimVelocity, Zero)
-                    );
-                }
-
-                //- Momentum transfer coefficient [kg]
-                virtual tmp<volScalarField::Internal> UCoeff() const
-                {
-                    return volScalarField::Internal::New
-                    (
-                        IOobject::scopedName(this->name(), "UCoeff"),
-                        mesh_,
-                        dimensionedScalar(dimMass, Zero)
-                    );
-                }
-
-                //- Return momentum source term [kg m/s^2]
-                virtual tmp<fvVectorMatrix> SU(volVectorField& U) const
-                {
-                    return tmp<fvVectorMatrix>::New(U, dimMass*dimAcceleration);
-                }
-
-
-            // Energy
-
-                //- Return sensible enthalpy transfer [J]
-                virtual tmp<volScalarField::Internal> hsTrans() const
-                {
-                    return volScalarField::Internal::New
-                    (
-                        IOobject::scopedName(this->name(), "hsTrans"),
-                        mesh_,
-                        dimensionedScalar(dimEnergy, Zero)
-                    );
-                }
-
-                //- Return sensible enthalpy transfer coefficient [J/K]
-                virtual tmp<volScalarField::Internal> hsCoeff() const
-                {
-                    return volScalarField::Internal::New
-                    (
-                        IOobject::scopedName(this->name(), "hsCoeff"),
-                        mesh_,
-                        dimensionedScalar(dimEnergy/dimTemperature, Zero)
-                    );
-                }
-
-                //- Return sensible enthalpy source term [J/s]
-                virtual tmp<fvScalarMatrix> Sh(volScalarField& hs) const
-                {
-                    return tmp<fvScalarMatrix>::New(hs, dimEnergy/dimTime);
-                }
-
-
-            // Radiation - overrides thermoCloud virtual abstract members
-
-                //- Return equivalent particulate emission [kg/m/s^3]
-                virtual tmp<volScalarField> Ep() const
-                {
-                    return volScalarField::New
-                    (
-                        IOobject::scopedName(this->name(), "radiation_Ep"),
-                        mesh_,
-                        dimensionedScalar(dimMass/dimLength/pow3(dimTime), Zero)
-                    );
-                }
-
-                //- Return equivalent particulate absorption [1/m]
-                virtual tmp<volScalarField> ap() const
-                {
-                    return volScalarField::New
-                    (
-                        IOobject::scopedName(this->name(), "radiation_ap"),
-                        mesh_,
-                        dimensionedScalar(dimless/dimLength, Zero)
-                    );
-                }
-
-                //- Return equivalent particulate scattering factor [1/m]
-                virtual tmp<volScalarField> sigmap() const
-                {
-                    return volScalarField::New
-                    (
-                        IOobject::scopedName(this->name(), "radiation_sigmap"),
-                        mesh_,
-                        dimensionedScalar(dimless/dimLength, Zero)
-                    );
-                }
-
-
-            // Mass
-
-                //- Return mass source term for specie [kg/s]
-                virtual tmp<fvScalarMatrix> SYi
-                (
-                    const label speciei,
-                    volScalarField& Yi
-                ) const
-                {
-                    return tmp<fvScalarMatrix>::New(Yi, dimMass/dimTime);
-                }
-
-                //- Return total mass source term [kg/m^3/s]
-                virtual tmp<volScalarField::Internal> Srho() const
-                {
-                    return volScalarField::Internal::New
-                    (
-                        IOobject::scopedName(this->name(), "Srho"),
-                        mesh_,
-                        dimensionedScalar(dimDensity/dimTime, Zero)
-                    );
-                }
-
-                //- Return total mass source term [kg/s]
-                virtual tmp<fvScalarMatrix> Srho(volScalarField& rho) const
-                {
-                    return tmp<fvScalarMatrix>::New(rho, dimMass/dimTime);
-                }
-};
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-#endif
-
-// ************************************************************************* //
--- a/src/lagrangian/intermediate/clouds/Templates/CollidingCloud/CollidingCloud.C
+++ b/src/lagrangian/intermediate/clouds/Templates/CollidingCloud/CollidingCloud.C
@ -28,7 +28,6 @@ License
 #include "CollidingCloud.H"
 #include "CollisionModel.H"
 #include "NoCollision.H"
-#include "subCycleTime.H"

 // * * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * //

@ -89,14 +88,14 @@ template<class CloudType>
 Foam::CollidingCloud<CloudType>::CollidingCloud
 (
    const word& cloudName,
-    const dimensionedVector& g,
    const volScalarField& rho,
    const volVectorField& U,
    const volScalarField& mu,
+    const dimensionedVector& g,
    bool readFields
 )
 :
-    CloudType(cloudName, g, rho, U, mu, false),
+    CloudType(cloudName, rho, U, mu, g, false),
    constProps_(this->particleProperties()),
    collisionModel_(nullptr)
 {
@ -149,6 +148,13 @@ Foam::CollidingCloud<CloudType>::CollidingCloud
 {}


+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template<class CloudType>
+Foam::CollidingCloud<CloudType>::~CollidingCloud()
+{}
+
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //

 template<class CloudType>
--- a/src/lagrangian/intermediate/clouds/Templates/CollidingCloud/CollidingCloud.H
+++ b/src/lagrangian/intermediate/clouds/Templates/CollidingCloud/CollidingCloud.H
@ -58,8 +58,6 @@ namespace Foam
 template<class CloudType>
 class CollisionModel;

-TemplateName(CollidingCloud);
-
 /*---------------------------------------------------------------------------*\
                       Class CollidingCloud Declaration
 \*---------------------------------------------------------------------------*/
@ -67,8 +65,7 @@ TemplateName(CollidingCloud);
 template<class CloudType>
 class CollidingCloud
 :
-    public CloudType,
-    public CollidingCloudName
+    public CloudType
 {
 public:

@ -112,7 +109,8 @@ protected:
        // References to the cloud sub-models

            //- Collision model
-            autoPtr<CollisionModel<CollidingCloud<CloudType>>> collisionModel_;
+            autoPtr<CollisionModel<CollidingCloud<CloudType>>>
+                collisionModel_;


        // Initialisation
@ -144,10 +142,10 @@ public:
        CollidingCloud
        (
            const word& cloudName,
-            const dimensionedVector& g,
            const volScalarField& rho,
            const volVectorField& U,
            const volScalarField& mu,
+            const dimensionedVector& g,
            bool readFields = true
        );

@ -186,7 +184,7 @@ public:


    //- Destructor
-    virtual ~CollidingCloud() = default;
+    virtual ~CollidingCloud();


    // Member Functions
--- a/Show More
+++ b/Show More
Author	SHA1	Message	Date
Andreas Margetis	07664c0de8	TUT: added two unsteady adjoint optimisation tutorials showcasing the use of compressedFullStorage and binomialCheckPointing	2022-07-15 16:43:12 +03:00
Andreas Margetis	32c8fb57b7	ENH: introduced unsteady adjoint functionality - The unsteady adjoint equations are integrated backwards in time. Since each adjoint time-step requires the primal solution of that time-step to be known, schemes for managing the storage/retrieval of the entire flow series are necessary. These are implemented through the primalStorage class and its derived ones. The latter manipulate a new class of fields, called compressedGeometricFields, which provide hooks for compressing/decompressing a field during the time integration of the primal/adjoint equations. The method used for compressing/decompressing is run-time selectable. - The current commit provides the shortGeometricField implementation which avoids the storage of patchFields that can be retrieved from the internalField (e.g. coupled, zeroGradient, symmetry, etc) , to cut on the storage requirements. More elaborate compression approaches will be included in the future, during the exaFoam project. - Two primalStorage options are included: compressedFullStorage and binomialCheckPointing. - compressedFullStorage stores the entire flow time-series, potentially by compressing each time-step (only the above-mentioned short approach is available for the moment). - binomialCheckPointing is based on the homonymous algorithm found in \verbatim Wang, Q., Moin, P., & Iaccarino, G.. Minimal Repetition Dynamic Checkpointing Algorithm for Unsteady Adjoint Calculation (2009). SIAM Journal on Scientific Computing, 31(4), 2549-2567. 10.1137/080727890, \endverbatim which stores the solution of the flow equations in a predefined number of time-steps, named checkpoints. During the backwards-in-time integration of the adjoint equations, if the primal solution at a certain time-step is not available, it is retrieved by re-computing the primal flow field starting from the closest checkpoint. Checkpoints are optimally distributed throughout the time-series to invoke the least number of flow recomputations during the backwards-in-time solution of the adjoint equations. Binomial checkpointing is the current state of the art though its re-computation cost frequently amounts for an extra solution of the flow equations in medium-to-large cases. - The adjoint to the PISO and PIMPLE solvers, along with their solverControl variants, are additionally included. - Objective functions are integrated in time, through appropriate entries in the dictionaries defining them. Authored by Andreas Margetis and reviewed by Vaggelis Papoutsis, with earlier contributions from Dr. Ioannis Kavvadias.	2022-07-15 16:43:10 +03:00
Andreas Margetis	4ee7dd50ee	ENH: optionally read oldTimes in two GeometricField constructors These are used by the adjoint code and are necessary for unsteady adjoint simulations. Both additions are implemented through optional arguments with default values, to maintain backwards compatibility for the rest of the code base.	2022-07-15 16:43:07 +03:00
Andreas Margetis	162f5f29ec	ENH: modifications in Time to support unsteady adjoint Since the unsteady adjoint equations are integrated backwards in time, the -- operator and the reverseEnd and reverseLoop methods were added to control the flow of time and the ending criteria.	2022-07-15 16:43:06 +03:00
				`@ -1,2 +0,0 @@`
				`auto parcels = parcelCloudModelList(g, rhoInf, U, muc);`
				`@ -1 +0,0 @@`
				`regionModels::surfaceFilmModel& surfaceFilm = tsurfaceFilm();`