From c9fea670b6e3180b492bd54e7201ecfda21f716e Mon Sep 17 00:00:00 2001
From: Vanja Skuric <vanja.skuric@fsb.hr>
Date: Fri, 27 May 2016 22:07:30 +0200
Subject: [PATCH] Backported lagrangian solvers and tutorials (vanilla OF
 3.0.1)

---
 .../lagrangian/coalChemistryFoam/UEqn.H        |  2 +-
 .../coalChemistryFoam/coalChemistryFoam.C      |  9 ++++++---
 .../coalChemistryFoam/createFields.H           |  2 ++
 .../lagrangian/coalChemistryFoam/pEqn.H        | 10 +++++-----
 .../UEqn.H                                     |  2 +-
 .../createFields.H                             |  2 ++
 .../pEqn.H                                     | 13 +++----------
 .../porousExplicitSourceReactingParcelFoam.C   |  7 +++++--
 .../lagrangian/reactingParcelFoam/UEqn.H       |  2 +-
 .../reactingParcelFoam/createFields.H          |  2 ++
 .../lagrangian/reactingParcelFoam/pEqn.H       | 10 +++++-----
 .../reactingParcelFoam/reactingParcelFoam.C    |  9 ++++++---
 .../simplifiedSiwek/constant/polyMesh/boundary |  4 ++--
 .../simplifiedSiwek/system/fvSchemes           |  8 +-------
 .../simplifiedSiwek/system/fvSolution          |  2 +-
 .../constant/polyMesh/boundary                 |  5 ++---
 .../channelParticles/system/fvSchemes          |  6 ------
 .../icoLagrangianFoam/createFields.H           |  2 ++
 .../icoLagrangianFoam/icoLagrangianFoam.C      | 11 +++++++----
 .../filter/Allclean                            |  2 +-
 .../filter/constant/polyMesh/boundary          | 18 +++++++++---------
 .../filter/system/fvSchemes                    |  8 +-------
 .../parcelInBox/0/G                            |  4 ++++
 .../parcelInBox/0/H2O                          |  4 ++++
 .../parcelInBox/0/T                            |  4 ++++
 .../parcelInBox/0/U                            |  4 ++++
 .../parcelInBox/0/air                          |  4 ++++
 .../parcelInBox/0/p                            |  4 ++++
 .../parcelInBox/constant/polyMesh/boundary     |  4 ++--
 .../parcelInBox/system/fvSchemes               |  8 +-------
 .../verticalChannel/system/fvSchemes           |  8 +-------
 .../evaporationTest/constant/polyMesh/boundary |  4 ++--
 .../evaporationTest/system/fvSchemes           |  8 +-------
 .../evaporationTest/system/fvSolution          |  2 +-
 .../rhoPisoTwinParcelFoam/UEqn.H               |  2 +-
 .../rhoPisoTwinParcelFoam/createFields.H       |  2 ++
 .../rhoPisoTwinParcelFoam/pEqn.H               | 10 +++++-----
 .../rhoPisoTwinParcelFoam.C                    | 12 +++++++-----
 .../simplifiedSiwek/constant/polyMesh/boundary |  4 ++--
 .../simplifiedSiwek/system/fvSchemes           |  8 +-------
 .../simplifiedSiwek/system/fvSolution          |  2 +-
 41 files changed, 116 insertions(+), 118 deletions(-)

diff --git a/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H
index b02a7c61e..deb710cce 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/UEqn.H
@@ -11,7 +11,7 @@
 
     UEqn.relax();
 
-    if (momentumPredictor)
+    if (pimple.momentumPredictor())
     {
         solve(UEqn == -fvc::grad(p));
     }
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C b/applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C
index 4d031e521..570327101 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C
+++ b/applications/solvers/lagrangian/coalChemistryFoam/coalChemistryFoam.C
@@ -44,6 +44,7 @@ Description
 #include "chemistrySolver.H"
 #include "timeActivatedExplicitSource.H"
 #include "radiationModel.H"
+#include "pimpleControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -53,6 +54,9 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
+
+    pimpleControl pimple(mesh);
+
     #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
@@ -70,7 +74,6 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readTimeControls.H"
-        #include "readPISOControls.H"
         #include "compressibleCourantNo.H"
         #include "setDeltaT.H"
 
@@ -89,14 +92,14 @@ int main(int argc, char *argv[])
         #include "rhoEqn.H"
 
         // --- PIMPLE loop
-        for (int ocorr=1; ocorr<=nOuterCorr; ocorr++)
+        while (pimple.loop())
         {
             #include "UEqn.H"
             #include "YEqn.H"
             #include "hsEqn.H"
 
             // --- PISO loop
-            for (int corr=1; corr<=nCorr; corr++)
+            while (pimple.correct())
             {
                 #include "pEqn.H"
             }
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/createFields.H b/applications/solvers/lagrangian/coalChemistryFoam/createFields.H
index 65be09775..12a6e9548 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/createFields.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/createFields.H
@@ -149,3 +149,5 @@
         mesh,
         dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
     );
+
+    mesh.schemesDict().setFluxRequired(p.name());
diff --git a/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H b/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H
index b31ad45ed..9878707eb 100644
--- a/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H
+++ b/applications/solvers/lagrangian/coalChemistryFoam/pEqn.H
@@ -3,7 +3,7 @@ rho = thermo.rho();
 volScalarField rUA = 1.0/UEqn.A();
 U = rUA*UEqn.H();
 
-if (transonic)
+if (pimple.transonic())
 {
     surfaceScalarField phid
     (
@@ -15,7 +15,7 @@ if (transonic)
         )
     );
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (pimple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -28,7 +28,7 @@ if (transonic)
 
         pEqn.solve();
 
-        if (nonOrth == nNonOrthCorr)
+        if (pimple.finalNonOrthogonalIter())
         {
             phi == pEqn.flux();
         }
@@ -43,7 +43,7 @@ else
           + fvc::ddtPhiCorr(rUA, rho, U, phi)
         );
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (pimple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -56,7 +56,7 @@ else
 
         pEqn.solve();
 
-        if (nonOrth == nNonOrthCorr)
+        if (pimple.finalNonOrthogonalIter())
         {
             phi += pEqn.flux();
         }
diff --git a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H
index e77fe75da..f46d1d902 100644
--- a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/UEqn.H
@@ -12,7 +12,7 @@
 
     pZones.addResistance(UEqn);
 
-    if (momentumPredictor)
+    if (piso.momentumPredictor())
     {
         solve(UEqn == -fvc::grad(p));
     }
diff --git a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/createFields.H b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/createFields.H
index cf3f484bd..7a5b48ac9 100644
--- a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/createFields.H
@@ -103,3 +103,5 @@
         mesh,
         dimensionedScalar("chemistrySh", dimEnergy/dimTime/dimVolume, 0.0)
     );
+
+    mesh.schemesDict().setFluxRequired(p.name());
diff --git a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H
index 43d677609..ff590e59e 100644
--- a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H
+++ b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/pEqn.H
@@ -23,7 +23,7 @@
             );
     }
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (piso.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -35,16 +35,9 @@
           + massSource.SuTot()
         );
 
-        if (corr == nCorr-1 && nonOrth == nNonOrthCorr)
-        {
-            pEqn.solve(mesh.solutionDict().solver("pFinal"));
-        }
-        else
-        {
-            pEqn.solve();
-        }
+        pEqn.solve(mesh.solutionDict().solver(p.select(piso.finalInnerIter())));
 
-        if (nonOrth == nNonOrthCorr)
+        if (piso.finalNonOrthogonalIter())
         {
             phi += pEqn.flux();
         }
diff --git a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/porousExplicitSourceReactingParcelFoam.C b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/porousExplicitSourceReactingParcelFoam.C
index de63f2fbb..68ef32448 100644
--- a/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/porousExplicitSourceReactingParcelFoam.C
+++ b/applications/solvers/lagrangian/porousExplicitSourceReactingParcelFoam/porousExplicitSourceReactingParcelFoam.C
@@ -49,6 +49,7 @@ Description
 #include "radiationModel.H"
 #include "porousZones.H"
 #include "timeActivatedExplicitSource.H"
+#include "pisoControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -58,6 +59,9 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
+
+    pisoControl piso(mesh);
+
     #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
@@ -77,7 +81,6 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readTimeControls.H"
-        #include "readPISOControls.H"
         #include "readAdditionalSolutionControls.H"
         #include "compressibleCourantNo.H"
         #include "setDeltaT.H"
@@ -95,7 +98,7 @@ int main(int argc, char *argv[])
         #include "hsEqn.H"
 
         // --- PISO loop
-        for (int corr = 0; corr < nCorr; corr++)
+        while (piso.correct())
         {
             #include "pEqn.H"
         }
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/UEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/UEqn.H
index 3c4a92709..1c1c1fd0c 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/UEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/UEqn.H
@@ -10,7 +10,7 @@
 
     UEqn.relax();
 
-    if (momentumPredictor)
+    if (pimple.momentumPredictor())
     {
         solve(UEqn == -fvc::grad(p));
     }
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/createFields.H b/applications/solvers/lagrangian/reactingParcelFoam/createFields.H
index e071ac8b3..9579be78c 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/createFields.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/createFields.H
@@ -109,3 +109,5 @@
         mesh,
         dimensionedScalar("chemistry::Sh", dimEnergy/dimTime/dimVolume, 0.0)
     );
+
+    mesh.schemesDict().setFluxRequired(p.name());
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/pEqn.H b/applications/solvers/lagrangian/reactingParcelFoam/pEqn.H
index d27afcce8..50185f5fc 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/pEqn.H
+++ b/applications/solvers/lagrangian/reactingParcelFoam/pEqn.H
@@ -3,7 +3,7 @@ rho = thermo.rho();
 volScalarField rUA = 1.0/UEqn.A();
 U = rUA*UEqn.H();
 
-if (transonic)
+if (pimple.transonic())
 {
     surfaceScalarField phid
     (
@@ -15,7 +15,7 @@ if (transonic)
         )
     );
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (pimple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -28,7 +28,7 @@ if (transonic)
 
         pEqn.solve();
 
-        if (nonOrth == nNonOrthCorr)
+        if (pimple.finalNonOrthogonalIter())
         {
             phi == pEqn.flux();
         }
@@ -43,7 +43,7 @@ else
           + fvc::ddtPhiCorr(rUA, rho, U, phi)
         );
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (pimple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -56,7 +56,7 @@ else
 
         pEqn.solve();
 
-        if (nonOrth == nNonOrthCorr)
+        if (pimple.finalNonOrthogonalIter())
         {
             phi += pEqn.flux();
         }
diff --git a/applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C b/applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C
index 154774cca..93971a15b 100644
--- a/applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C
+++ b/applications/solvers/lagrangian/reactingParcelFoam/reactingParcelFoam.C
@@ -37,6 +37,7 @@ Description
 #include "psiChemistryModel.H"
 #include "chemistrySolver.H"
 #include "radiationModel.H"
+#include "pimpleControl.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -46,6 +47,9 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
+
+    pimpleControl pimple(mesh);
+
     #include "readChemistryProperties.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
@@ -63,7 +67,6 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readTimeControls.H"
-        #include "readPISOControls.H"
         #include "compressibleCourantNo.H"
         #include "setDeltaT.H"
 
@@ -77,13 +80,13 @@ int main(int argc, char *argv[])
         #include "rhoEqn.H"
 
         // --- PIMPLE loop
-        for (int ocorr=1; ocorr<=nOuterCorr; ocorr++)
+        while (pimple.loop())
         {
             #include "UEqn.H"
             #include "YEqn.H"
 
             // --- PISO loop
-            for (int corr=1; corr<=nCorr; corr++)
+            while (pimple.correct())
             {
                 #include "hsEqn.H"
                 #include "pEqn.H"
diff --git a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/polyMesh/boundary b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/polyMesh/boundary
index dff7e8226..eedf0fdb7 100644
--- a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/polyMesh/boundary
+++ b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/constant/polyMesh/boundary
@@ -1,9 +1,9 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | foam-extend: Open Source CFD                    |
-|  \\    /   O peration     | Version:     3.2                                |
+|  \\    /   O peration     | Version:     4.0                                |
 |   \\  /    A nd           | Web:         http://www.foam-extend.org         |
-|    \\/     M anipulation  |                                                 |
+|    \\/     M anipulation  | For copyright notice see file Copyright         |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
diff --git a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSchemes b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSchemes
index c538134fd..d98d660a2 100644
--- a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSchemes
+++ b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSchemes
@@ -36,7 +36,7 @@ divSchemes
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
     div(U)          Gauss linear;
-    div((muEff*dev2(grad(U).T()))) Gauss linear;
+    div((muEff*dev2(T(grad(U))))) Gauss linear;
     div(phi,Yi_h)   Gauss upwind;
 }
 
@@ -62,10 +62,4 @@ snGradSchemes
     default         corrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p               ;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSolution b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSolution
index e48e0f6b4..d459df935 100644
--- a/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSolution
+++ b/tutorials/lagrangian/coalChemistryFoam/simplifiedSiwek/system/fvSolution
@@ -50,7 +50,7 @@ solvers
     }
 }
 
-PISO
+PIMPLE
 {
     transonic       yes;
     nCorrectors     2;
diff --git a/tutorials/lagrangian/icoLagrangianFoam/channelParticles/constant/polyMesh/boundary b/tutorials/lagrangian/icoLagrangianFoam/channelParticles/constant/polyMesh/boundary
index 0beee5a65..8e0826950 100644
--- a/tutorials/lagrangian/icoLagrangianFoam/channelParticles/constant/polyMesh/boundary
+++ b/tutorials/lagrangian/icoLagrangianFoam/channelParticles/constant/polyMesh/boundary
@@ -1,9 +1,9 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | foam-extend: Open Source CFD                    |
-|  \\    /   O peration     | Version:     3.2                                |
+|  \\    /   O peration     | Version:     4.0                                |
 |   \\  /    A nd           | Web:         http://www.foam-extend.org         |
-|    \\/     M anipulation  |                                                 |
+|    \\/     M anipulation  | For copyright notice see file Copyright         |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
@@ -38,7 +38,6 @@ FoamFile
     frontAndBack
     {
         type            empty;
-        physicalType    empty;
         nFaces          2000;
         startFace       2110;
     }
diff --git a/tutorials/lagrangian/icoLagrangianFoam/channelParticles/system/fvSchemes b/tutorials/lagrangian/icoLagrangianFoam/channelParticles/system/fvSchemes
index 277df7d18..e807e470c 100644
--- a/tutorials/lagrangian/icoLagrangianFoam/channelParticles/system/fvSchemes
+++ b/tutorials/lagrangian/icoLagrangianFoam/channelParticles/system/fvSchemes
@@ -49,10 +49,4 @@ snGradSchemes
     default         corrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/createFields.H b/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/createFields.H
index 794eb514b..54c52cb09 100644
--- a/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/createFields.H
+++ b/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/createFields.H
@@ -53,6 +53,8 @@
     label pRefCell = 0;
     scalar pRefValue = 0.0;
     setRefCell(p, mesh.solutionDict().subDict("PISO"), pRefCell, pRefValue);
+    mesh.schemesDict().setFluxRequired(p.name());
+
     dimensionedScalar rhoVal
     (
         transportProperties.lookup("rho")
diff --git a/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/icoLagrangianFoam.C b/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/icoLagrangianFoam.C
index 34efa82f2..36f4462a0 100644
--- a/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/icoLagrangianFoam.C
+++ b/tutorials/lagrangian/icoLagrangianFoam/icoLagrangianFoam/icoLagrangianFoam.C
@@ -31,6 +31,7 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
+#include "pisoControl.H"
 
 #include "basicKinematicCloud.H"
 
@@ -42,6 +43,9 @@ int main(int argc, char *argv[])
 
 #   include "createTime.H"
 #   include "createMesh.H"
+
+    pisoControl piso(mesh);
+
 #   include "createFields.H"
 #   include "initContinuityErrs.H"
 
@@ -53,7 +57,6 @@ int main(int argc, char *argv[])
     {
         Info<< "Time = " << runTime.timeName() << nl << endl;
 
-#       include "readPISOControls.H"
 #       include "CourantNo.H"
 
         Info<< "Evolving kinematicCloud" << endl;
@@ -73,7 +76,7 @@ int main(int argc, char *argv[])
 
         // --- PISO loop
 
-        for (int corr=0; corr<nCorr; corr++)
+        while (piso.correct())
         {
             volScalarField rUA = 1.0/UEqn.A();
 
@@ -83,7 +86,7 @@ int main(int argc, char *argv[])
 
             adjustPhi(phi, U, p);
 
-            for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+            while (piso.correctNonOrthogonal())
             {
                 fvScalarMatrix pEqn
                 (
@@ -93,7 +96,7 @@ int main(int argc, char *argv[])
                 pEqn.setReference(pRefCell, pRefValue);
                 pEqn.solve();
 
-                if (nonOrth == nNonOrthCorr)
+                if (piso.finalNonOrthogonalIter())
                 {
                     phi -= pEqn.flux();
                 }
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/Allclean b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/Allclean
index cb015f101..80fcf5310 100755
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/Allclean
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/Allclean
@@ -6,7 +6,7 @@
 cleanCase
 
 # remove 0 directory
-rm -rf 0
+rm -rf 0 postProcessing
 
 # copy 0.org to 0
 cp -r 0.org 0
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/constant/polyMesh/boundary b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/constant/polyMesh/boundary
index 9be6801c5..211bb44ae 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/constant/polyMesh/boundary
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/constant/polyMesh/boundary
@@ -1,9 +1,9 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | foam-extend: Open Source CFD                    |
-|  \\    /   O peration     | Version:     3.2                                |
+|  \\    /   O peration     | Version:     4.0                                |
 |   \\  /    A nd           | Web:         http://www.foam-extend.org         |
-|    \\/     M anipulation  |                                                 |
+|    \\/     M anipulation  | For copyright notice see file Copyright         |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
@@ -21,39 +21,39 @@ FoamFile
     {
         type            wall;
         nFaces          172;
-        startFace       3334;
+        startFace       3314;
     }
     inlet
     {
         type            patch;
         nFaces          20;
-        startFace       3506;
+        startFace       3486;
     }
     outlet
     {
         type            patch;
         nFaces          20;
-        startFace       3526;
+        startFace       3506;
     }
     cycLeft
     {
         type            cyclic;
-        nFaces          0;
-        startFace       3546;
+        nFaces          40;
+        startFace       3526;
         featureCos      0.9;
     }
     cycRight
     {
         type            cyclic;
         nFaces          0;
-        startFace       3546;
+        startFace       3566;
         featureCos      0.9;
     }
     frontAndBack
     {
         type            empty;
         nFaces          3440;
-        startFace       3546;
+        startFace       3566;
     }
 )
 
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/system/fvSchemes b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/system/fvSchemes
index c538134fd..d98d660a2 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/system/fvSchemes
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/filter/system/fvSchemes
@@ -36,7 +36,7 @@ divSchemes
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
     div(U)          Gauss linear;
-    div((muEff*dev2(grad(U).T()))) Gauss linear;
+    div((muEff*dev2(T(grad(U))))) Gauss linear;
     div(phi,Yi_h)   Gauss upwind;
 }
 
@@ -62,10 +62,4 @@ snGradSchemes
     default         corrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p               ;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/G b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/G
index c5b215263..0ef856417 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/G
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/G
@@ -47,6 +47,10 @@ boundaryField
         emissivity      1;
         value           uniform 0;
     }
+    defaultFaces
+    {
+        type            empty;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/H2O b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/H2O
index c9b2b3a93..8b893a122 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/H2O
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/H2O
@@ -44,6 +44,10 @@ boundaryField
         type            fixedValue;
         value           uniform 0.01;
     }
+    defaultFaces
+    {
+        type            empty;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/T b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/T
index b2aab8347..5e1a51629 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/T
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/T
@@ -43,6 +43,10 @@ boundaryField
         type            fixedValue;
         value           uniform 473.0;
     }
+    defaultFaces
+    {
+        type            empty;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/U b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/U
index 84938a483..8aa7ead8e 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/U
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/U
@@ -43,6 +43,10 @@ boundaryField
         type            fixedValue;
         value           uniform (0 1e-8 0);
     }
+    defaultFaces
+    {
+        type            empty;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/air b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/air
index 5b0d480d9..9fa907810 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/air
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/air
@@ -44,6 +44,10 @@ boundaryField
         type            fixedValue;
         value           uniform 0.99;
     }
+    defaultFaces
+    {
+        type            empty;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/p b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/p
index cae73b9ca..155b08669 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/p
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/0/p
@@ -42,6 +42,10 @@ boundaryField
     {
         type            zeroGradient;
     }
+    defaultFaces
+    {
+        type            empty;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/constant/polyMesh/boundary b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/constant/polyMesh/boundary
index d960c9845..299824a75 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/constant/polyMesh/boundary
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/constant/polyMesh/boundary
@@ -1,9 +1,9 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | foam-extend: Open Source CFD                    |
-|  \\    /   O peration     | Version:     3.2                                |
+|  \\    /   O peration     | Version:     4.0                                |
 |   \\  /    A nd           | Web:         http://www.foam-extend.org         |
-|    \\/     M anipulation  |                                                 |
+|    \\/     M anipulation  | For copyright notice see file Copyright         |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/system/fvSchemes b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/system/fvSchemes
index 0aff01154..233870cb1 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/system/fvSchemes
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/parcelInBox/system/fvSchemes
@@ -37,7 +37,7 @@ divSchemes
     div(phi,epsilon) Gauss upwind;
     div(phi,omega) Gauss upwind;
 //    div(U)          Gauss linear;
-    div((muEff*dev2(grad(U).T()))) Gauss linear;
+    div((muEff*dev2(T(grad(U))))) Gauss linear;
     div(phi,Yi_h)   Gauss upwind;
 }
 
@@ -56,10 +56,4 @@ snGradSchemes
     default         uncorrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/verticalChannel/system/fvSchemes b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/verticalChannel/system/fvSchemes
index e61a42743..a4b2cfc04 100644
--- a/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/verticalChannel/system/fvSchemes
+++ b/tutorials/lagrangian/porousExplicitSourceReactingParcelFoam/verticalChannel/system/fvSchemes
@@ -38,7 +38,7 @@ divSchemes
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
     div(phi,omega) Gauss upwind;
-    div((muEff*dev2(grad(U).T()))) Gauss linear;
+    div((muEff*dev2(T(grad(U))))) Gauss linear;
     div(phi,Yi_h)   Gauss upwind;
 }
 
@@ -57,10 +57,4 @@ snGradSchemes
     default         uncorrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/polyMesh/boundary b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/polyMesh/boundary
index 8451880b5..9803a53e2 100644
--- a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/polyMesh/boundary
+++ b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/constant/polyMesh/boundary
@@ -1,9 +1,9 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | foam-extend: Open Source CFD                    |
-|  \\    /   O peration     | Version:     3.2                                |
+|  \\    /   O peration     | Version:     4.0                                |
 |   \\  /    A nd           | Web:         http://www.foam-extend.org         |
-|    \\/     M anipulation  |                                                 |
+|    \\/     M anipulation  | For copyright notice see file Copyright         |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
diff --git a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSchemes b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSchemes
index c0c891c23..d98d660a2 100644
--- a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSchemes
+++ b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSchemes
@@ -36,7 +36,7 @@ divSchemes
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
     div(U)          Gauss linear;
-    div((muEff*dev2(grad(U).T()))) Gauss linear;
+    div((muEff*dev2(T(grad(U))))) Gauss linear;
     div(phi,Yi_h)   Gauss upwind;
 }
 
@@ -62,10 +62,4 @@ snGradSchemes
     default         corrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSolution b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSolution
index 11c0ed404..fbc3a5ee3 100644
--- a/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSolution
+++ b/tutorials/lagrangian/reactingParcelFoam/evaporationTest/system/fvSolution
@@ -96,7 +96,7 @@ solvers
     }
 }
 
-PISO
+PIMPLE
 {
     transonic       no;
     nCorrectors     2;
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/UEqn.H b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/UEqn.H
index 5f2e597a1..d913e3fa1 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/UEqn.H
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/UEqn.H
@@ -11,7 +11,7 @@
 
     UEqn.relax();
 
-    if (momentumPredictor)
+    if (pimple.momentumPredictor())
     {
         solve(UEqn == -fvc::grad(p));
     }
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/createFields.H b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/createFields.H
index 46bbe2d91..39a93262b 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/createFields.H
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/createFields.H
@@ -57,3 +57,5 @@
     Info<< "Creating field DpDt\n" << endl;
     volScalarField DpDt =
         fvc::DDt(surfaceScalarField("phiU", phi/fvc::interpolate(rho)), p);
+
+    mesh.schemesDict().setFluxRequired(p.name());
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/pEqn.H b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/pEqn.H
index 9443f909a..2881cfbf8 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/pEqn.H
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/pEqn.H
@@ -3,7 +3,7 @@ rho = thermo.rho();
 volScalarField rUA = 1.0/UEqn.A();
 U = rUA*UEqn.H();
 
-if (transonic)
+if (pimple.transonic())
 {
     surfaceScalarField phid
     (
@@ -15,7 +15,7 @@ if (transonic)
         )
     );
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (pimple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -26,7 +26,7 @@ if (transonic)
 
         pEqn.solve();
 
-        if (nonOrth == nNonOrthCorr)
+        if (pimple.finalNonOrthogonalIter())
         {
             phi == pEqn.flux();
         }
@@ -41,7 +41,7 @@ else
           + fvc::ddtPhiCorr(rUA, rho, U, phi)
         );
 
-    for (int nonOrth=0; nonOrth<=nNonOrthCorr; nonOrth++)
+    while (pimple.correctNonOrthogonal())
     {
         fvScalarMatrix pEqn
         (
@@ -52,7 +52,7 @@ else
 
         pEqn.solve();
 
-        if (nonOrth == nNonOrthCorr)
+        if (pimple.finalNonOrthogonalIter())
         {
             phi += pEqn.flux();
         }
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam.C b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam.C
index 4ca4709fe..6161d1fee 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam.C
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam/rhoPisoTwinParcelFoam.C
@@ -32,6 +32,7 @@ Description
 #include "fvCFD.H"
 #include "basicPsiThermo.H"
 #include "turbulenceModel.H"
+#include "pimpleControl.H"
 
 #include "basicThermoCloud.H"
 #include "basicKinematicCloud.H"
@@ -44,12 +45,14 @@ int main(int argc, char *argv[])
 
     #include "createTime.H"
     #include "createMesh.H"
+
+    pimpleControl pimple(mesh);
+
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
     #include "createClouds.H"
-    #include "readPISOControls.H"
     #include "initContinuityErrs.H"
-    #include "readTimeControls.H"
+    #include "createTimeControls.H"
     #include "compressibleCourantNo.H"
     #include "setInitialDeltaT.H"
 
@@ -60,7 +63,6 @@ int main(int argc, char *argv[])
     while (runTime.run())
     {
         #include "readTimeControls.H"
-        #include "readPISOControls.H"
         #include "compressibleCourantNo.H"
         #include "setDeltaT.H"
 
@@ -78,12 +80,12 @@ int main(int argc, char *argv[])
         #include "rhoEqn.H"
 
         // --- PIMPLE loop
-        for (int ocorr=1; ocorr<=nOuterCorr; ocorr++)
+        while (pimple.loop())
         {
             #include "UEqn.H"
 
             // --- PISO loop
-            for (int corr=1; corr<=nCorr; corr++)
+            while (pimple.correct())
             {
                 #include "hsEqn.H"
                 #include "pEqn.H"
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/constant/polyMesh/boundary b/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/constant/polyMesh/boundary
index dff7e8226..eedf0fdb7 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/constant/polyMesh/boundary
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/constant/polyMesh/boundary
@@ -1,9 +1,9 @@
 /*--------------------------------*- C++ -*----------------------------------*\
 | =========                 |                                                 |
 | \\      /  F ield         | foam-extend: Open Source CFD                    |
-|  \\    /   O peration     | Version:     3.2                                |
+|  \\    /   O peration     | Version:     4.0                                |
 |   \\  /    A nd           | Web:         http://www.foam-extend.org         |
-|    \\/     M anipulation  |                                                 |
+|    \\/     M anipulation  | For copyright notice see file Copyright         |
 \*---------------------------------------------------------------------------*/
 FoamFile
 {
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSchemes b/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSchemes
index c538134fd..d98d660a2 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSchemes
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSchemes
@@ -36,7 +36,7 @@ divSchemes
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
     div(U)          Gauss linear;
-    div((muEff*dev2(grad(U).T()))) Gauss linear;
+    div((muEff*dev2(T(grad(U))))) Gauss linear;
     div(phi,Yi_h)   Gauss upwind;
 }
 
@@ -62,10 +62,4 @@ snGradSchemes
     default         corrected;
 }
 
-fluxRequired
-{
-    default         no;
-    p               ;
-}
-
 // ************************************************************************* //
diff --git a/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSolution b/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSolution
index 1fbde64c4..480f43c40 100644
--- a/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSolution
+++ b/tutorials/lagrangian/rhoPisoTwinParcelFoam/simplifiedSiwek/system/fvSolution
@@ -74,7 +74,7 @@ solvers
     }
 }
 
-PISO
+PIMPLE
 {
     transonic       yes;
     nCorrectors     2;