diff --git a/applications/solvers/compressible/rhoPimpleFoam/pEqn.H b/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
index a4e5fa20da..e0619a513b 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/pEqn.H
@@ -21,9 +21,9 @@ if (pimple.transonic())
         "phid",
         fvc::interpolate(psi)
        *(
-            (fvc::interpolate(rho*HbyA) & mesh.Sf())
-          + rhorAUf*fvc::ddtCorr(rho, U, phi)
-        )/fvc::interpolate(rho)
+            (fvc::interpolate(HbyA) & mesh.Sf())
+          + rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
+        )
     );
 
     fvOptions.makeRelative(fvc::interpolate(psi), phid);
diff --git a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H
index 1c58a02c33..e9fab27bc7 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/rhoPimpleDyMFoam/pEqn.H
@@ -21,9 +21,9 @@ if (pimple.transonic())
         "phid",
         fvc::interpolate(psi)
        *(
-            (fvc::interpolate(rho*HbyA) & mesh.Sf())
-          + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)
-        )/fvc::interpolate(rho)
+            (fvc::interpolate(HbyA) & mesh.Sf())
+          + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)/fvc::interpolate(rho)
+        )
     );
 
     fvc::makeRelative(phid, psi, U);
diff --git a/applications/solvers/compressible/rhoPimpleFoam/rhoPimplecFoam/pEqn.H b/applications/solvers/compressible/rhoPimpleFoam/rhoPimplecFoam/pEqn.H
index 0503db4207..7363add2b9 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/rhoPimplecFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoPimpleFoam/rhoPimplecFoam/pEqn.H
@@ -20,9 +20,10 @@ if (pimple.transonic())
         "phid",
         fvc::interpolate(psi)
        *(
-            (fvc::interpolate(rho*HbyA) & mesh.Sf())
+            (fvc::interpolate(HbyA) & mesh.Sf())
           + fvc::interpolate(rho*rAU)*fvc::ddtCorr(rho, U, phi)
-        )/fvc::interpolate(rho)
+           /fvc::interpolate(rho)
+        )
     );
 
     fvOptions.makeRelative(fvc::interpolate(psi), phid);
diff --git a/applications/solvers/compressible/rhoSimpleFoam/pEqn.H b/applications/solvers/compressible/rhoSimpleFoam/pEqn.H
index 444f5e996a..b1fc930be9 100644
--- a/applications/solvers/compressible/rhoSimpleFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoSimpleFoam/pEqn.H
@@ -13,8 +13,7 @@
         (
             "phid",
             fvc::interpolate(psi)
-           *(fvc::interpolate(rho*HbyA) & mesh.Sf())
-           /fvc::interpolate(rho)
+           *(fvc::interpolate(HbyA) & mesh.Sf())
         );
 
         fvOptions.makeRelative(fvc::interpolate(psi), phid);
diff --git a/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H b/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H
index c656b3c3d7..78878dcfaf 100644
--- a/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H
+++ b/applications/solvers/compressible/rhoSimpleFoam/rhoSimplecFoam/pEqn.H
@@ -14,8 +14,7 @@ if (simple.transonic())
     (
         "phid",
         fvc::interpolate(psi)
-       *(fvc::interpolate(rho*HbyA) & mesh.Sf())
-       /fvc::interpolate(rho)
+       *(fvc::interpolate(HbyA) & mesh.Sf())
     );
 
     surfaceScalarField phic
diff --git a/applications/solvers/compressible/sonicFoam/pEqn.H b/applications/solvers/compressible/sonicFoam/pEqn.H
index edde4b198e..022abc466a 100644
--- a/applications/solvers/compressible/sonicFoam/pEqn.H
+++ b/applications/solvers/compressible/sonicFoam/pEqn.H
@@ -11,9 +11,9 @@ surfaceScalarField phid
     "phid",
     fvc::interpolate(psi)*
     (
-        (mesh.Sf() & fvc::interpolate(rho*HbyA))
-      + rhorAUf*fvc::ddtCorr(rho, U, phi)
-    )/fvc::interpolate(rho)
+        (mesh.Sf() & fvc::interpolate(HbyA))
+      + rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
+    )
 );
 
 // Non-orthogonal pressure corrector loop
diff --git a/applications/solvers/compressible/sonicFoam/sonicDyMFoam/pEqn.H b/applications/solvers/compressible/sonicFoam/sonicDyMFoam/pEqn.H
index 4799f64567..b375f0e82a 100644
--- a/applications/solvers/compressible/sonicFoam/sonicDyMFoam/pEqn.H
+++ b/applications/solvers/compressible/sonicFoam/sonicDyMFoam/pEqn.H
@@ -12,9 +12,9 @@ surfaceScalarField phid
     fvc::interpolate(psi)
    *(
         (
-            (fvc::interpolate(rho*HbyA) & mesh.Sf())
-          + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)
-        )/fvc::interpolate(rho)
+            (fvc::interpolate(HbyA) & mesh.Sf())
+          + rhorAUf*fvc::ddtCorr(rho, U, rhoUf)/fvc::interpolate(rho)
+        )
     )
 );
 
diff --git a/applications/solvers/compressible/sonicFoam/sonicLiquidFoam/sonicLiquidFoam.C b/applications/solvers/compressible/sonicFoam/sonicLiquidFoam/sonicLiquidFoam.C
index df3f3cdc61..eb0a4a0d9a 100644
--- a/applications/solvers/compressible/sonicFoam/sonicLiquidFoam/sonicLiquidFoam.C
+++ b/applications/solvers/compressible/sonicFoam/sonicLiquidFoam/sonicLiquidFoam.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -89,9 +89,9 @@ int main(int argc, char *argv[])
                     "phid",
                     psi
                    *(
-                       (fvc::interpolate(rho*U) & mesh.Sf())
-                     + rhorAUf*fvc::ddtCorr(rho, U, phi)
-                    )/fvc::interpolate(rho)
+                       (fvc::interpolate(U) & mesh.Sf())
+                     + rhorAUf*fvc::ddtCorr(rho, U, phi)/fvc::interpolate(rho)
+                    )
                 );
 
                 phi = (rhoO/psi)*phid;
diff --git a/applications/solvers/lagrangian/DPMFoam/DPMFoam.C b/applications/solvers/lagrangian/DPMFoam/DPMFoam.C
index 513c123ddf..ee8ff02742 100644
--- a/applications/solvers/lagrangian/DPMFoam/DPMFoam.C
+++ b/applications/solvers/lagrangian/DPMFoam/DPMFoam.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -26,7 +26,7 @@ Application
 
 Description
     Transient solver for the coupled transport of a single kinematic particle
-    could including the effect of the volume fraction of particles on the
+    cloud including the effect of the volume fraction of particles on the
     continuous phase.
 
 \*---------------------------------------------------------------------------*/
diff --git a/applications/solvers/lagrangian/DPMFoam/MPPICFoam/MPPICFoam.C b/applications/solvers/lagrangian/DPMFoam/MPPICFoam/MPPICFoam.C
index 5a036f8f50..69048cf5a6 100644
--- a/applications/solvers/lagrangian/DPMFoam/MPPICFoam/MPPICFoam.C
+++ b/applications/solvers/lagrangian/DPMFoam/MPPICFoam/MPPICFoam.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -26,7 +26,7 @@ Application
 
 Description
     Transient solver for the coupled transport of a single kinematic particle
-    could including the effect of the volume fraction of particles on the
+    cloud including the effect of the volume fraction of particles on the
     continuous phase. Multi-Phase Particle In Cell (MPPIC) modeling is used to
     represent collisions without resolving particle-particle interactions.
 
diff --git a/applications/solvers/multiphase/interFoam/createFields.H b/applications/solvers/multiphase/interFoam/createFields.H
index 00a7d12db6..1c91f280bc 100644
--- a/applications/solvers/multiphase/interFoam/createFields.H
+++ b/applications/solvers/multiphase/interFoam/createFields.H
@@ -56,8 +56,6 @@
 
 
     // Mass flux
-    // Initialisation does not matter because rhoPhi is reset after the
-    // alpha1 solution before it is used in the U equation.
     surfaceScalarField rhoPhi
     (
         IOobject
@@ -68,8 +66,7 @@
             IOobject::NO_READ,
             IOobject::NO_WRITE
         ),
-        mesh,
-        dimensionedScalar("0", dimMass/dimTime, 0)
+        fvc::interpolate(rho)*phi
     );
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C
index 6c0fcbecec..2e7d971679 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.C
@@ -60,20 +60,16 @@ Foam::dragModels::Ergun::~Ergun()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::Ergun::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::Ergun::CdRe() const
 {
     return
         (4/3)
        *(
             150
-           *max(pair_.dispersed(), residualAlpha_)
-           *pair_.continuous().nu()
-           /(
-                max(scalar(1) - pair_.dispersed(), residualAlpha_)
-               *pair_.dispersed().d()
-               *max(pair_.magUr(), residualSlip_)
-            )
+           *max(scalar(1) - pair_.continuous(), residualAlpha_)
+           /max(pair_.continuous(), residualAlpha_)
           + 1.75
+           *pair_.Re()
         );
 }
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.H
index f8f29a8482..89e47dbb63 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Ergun/Ergun.H
@@ -82,7 +82,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C
index daaaf0bdc4..51e28ffdfb 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.C
@@ -48,8 +48,7 @@ Foam::dragModels::Gibilaro::Gibilaro
     const bool registerObject
 )
 :
-    dragModel(dict, pair, registerObject),
-    residualRe_("residualRe", dimless, dict.lookup("residualRe"))
+    dragModel(dict, pair, registerObject)
 {}
 
 
@@ -61,13 +60,13 @@ Foam::dragModels::Gibilaro::~Gibilaro()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::Gibilaro::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::Gibilaro::CdRe() const
 {
     volScalarField alpha2(max(scalar(1) - pair_.dispersed(), residualAlpha_));
 
     return
         (4/3)
-       *(17.3/(alpha2*max(pair_.Re(), residualRe_)) + scalar(0.336))
+       *(17.3/alpha2 + 0.336*pair_.Re())
        *max(pair_.continuous(), residualAlpha_)
        *pow(alpha2, -2.8);
 }
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.H
index c2e103e911..6b1d3ede03 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/Gibilaro/Gibilaro.H
@@ -58,14 +58,6 @@ class Gibilaro
 :
     public dragModel
 {
-private:
-
-    // Private data
-
-        //- Residual Reynolds number
-        const dimensionedScalar residualRe_;
-
-
 public:
 
     //- Runtime type information
@@ -90,7 +82,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C
index 9eb46e5060..ed8a871416 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.C
@@ -68,8 +68,7 @@ Foam::dragModels::GidaspowErgunWenYu::GidaspowErgunWenYu
             pair,
             false
         )
-    ),
-    residualRe_("residualRe", dimless, dict.lookup("residualRe"))
+    )
 {}
 
 
@@ -82,11 +81,11 @@ Foam::dragModels::GidaspowErgunWenYu::~GidaspowErgunWenYu()
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 Foam::tmp<Foam::volScalarField>
-Foam::dragModels::GidaspowErgunWenYu::Cd() const
+Foam::dragModels::GidaspowErgunWenYu::CdRe() const
 {
     return
-        pos(pair_.continuous() - 0.8)*WenYu_->Cd()
-      + neg(pair_.continuous() - 0.8)*Ergun_->Cd();
+        pos(pair_.continuous() - 0.8)*WenYu_->CdRe()
+      + neg(pair_.continuous() - 0.8)*Ergun_->CdRe();
 }
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.H
index e29bc28dcb..4bc2666d9d 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowErgunWenYu/GidaspowErgunWenYu.H
@@ -69,9 +69,6 @@ private:
         //- Wen Yu drag model
         autoPtr<WenYu> WenYu_;
 
-        //- Residual Reynolds number
-        const dimensionedScalar residualRe_;
-
 
 public:
 
@@ -97,7 +94,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C
index ce35a4ac25..0571242d55 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.C
@@ -62,18 +62,18 @@ Foam::dragModels::GidaspowSchillerNaumann::~GidaspowSchillerNaumann()
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 Foam::tmp<Foam::volScalarField>
-Foam::dragModels::GidaspowSchillerNaumann::Cd() const
+Foam::dragModels::GidaspowSchillerNaumann::CdRe() const
 {
     volScalarField alpha2(max(scalar(1) - pair_.dispersed(), residualAlpha_));
-    volScalarField Re(max(alpha2*pair_.Re(), residualRe_));
-    volScalarField Cds
+    volScalarField Re(alpha2*pair_.Re());
+    volScalarField CdsRe
     (
-        neg(Re - 1000)*(24.0*(1.0 + 0.15*pow(Re, 0.687))/Re)
-      + pos(Re - 1000)*0.44
+        neg(Re - 1000)*24.0*(1.0 + 0.15*pow(Re, 0.687))/alpha2
+      + pos(Re - 1000)*0.44*max(Re, residualRe_)
     );
 
     return
-        Cds
+        CdsRe
        *pow(alpha2, -2.65)
        *max(pair_.continuous(), residualAlpha_);
 }
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.H
index b14ccf06b6..f89d15d076 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/GidaspowSchillerNaumann/GidaspowSchillerNaumann.H
@@ -69,7 +69,7 @@ private:
 
     // Private data
 
-        //- Residual Reynolds number
+        //- Residual Reynolds Number
         const dimensionedScalar residualRe_;
 
 
@@ -97,7 +97,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C
index 5504d1a523..2d4a654030 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.C
@@ -61,13 +61,13 @@ Foam::dragModels::SchillerNaumann::~SchillerNaumann()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::SchillerNaumann::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::SchillerNaumann::CdRe() const
 {
-    volScalarField Re(pair_.Re() + residualRe_);
+    volScalarField Re(pair_.Re());
 
     return
-        neg(Re - 1000)*(24.0*(1.0 + 0.15*pow(Re, 0.687))/Re)
-      + pos(Re - 1000)*0.44;
+        neg(Re - 1000)*24.0*(1.0 + 0.15*pow(Re, 0.687))
+      + pos(Re - 1000)*0.44*max(Re, residualRe_);
 }
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.H
index 21d13358ca..e3c6175b63 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SchillerNaumann/SchillerNaumann.H
@@ -58,7 +58,7 @@ private:
 
     // Private data
 
-        //- Residual Reynolds number
+        //- Residual Reynolds Number
         const dimensionedScalar residualRe_;
 
 
@@ -86,7 +86,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C
index fbf839178f..86f38db155 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.C
@@ -48,8 +48,7 @@ Foam::dragModels::SyamlalOBrien::SyamlalOBrien
     const bool registerObject
 )
 :
-    dragModel(dict, pair, registerObject),
-    residualRe_("residualRe", dimless, dict.lookup("residualRe"))
+    dragModel(dict, pair, registerObject)
 {}
 
 
@@ -61,17 +60,16 @@ Foam::dragModels::SyamlalOBrien::~SyamlalOBrien()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::SyamlalOBrien::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::SyamlalOBrien::CdRe() const
 {
     volScalarField alpha2(max(scalar(1) - pair_.dispersed(), residualAlpha_));
-    volScalarField magUr(max(pair_.magUr(), residualSlip_));
     volScalarField A(pow(alpha2, 4.14));
     volScalarField B
     (
         neg(alpha2 - 0.85)*(0.8*pow(alpha2, 1.28))
       + pos(alpha2 - 0.85)*(pow(alpha2, 2.65))
     );
-    volScalarField Re(max(pair_.Re(), residualRe_));
+    volScalarField Re(pair_.Re());
     volScalarField Vr
     (
         0.5
@@ -79,10 +77,10 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::SyamlalOBrien::Cd() const
             A - 0.06*Re + sqrt(sqr(0.06*Re) + 0.12*Re*(2.0*B - A) + sqr(A))
         )
     );
-    volScalarField Cds(sqr(0.63 + 4.8*sqrt(Vr/Re)));
+    volScalarField CdsRe(sqr(0.63*sqrt(Re) + 4.8*sqrt(Vr)));
 
     return
-        Cds
+        CdsRe
        *max(pair_.continuous(), residualAlpha_)
        /sqr(Vr);
 }
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.H
index 4902997245..0b45dc9421 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/SyamlalOBrien/SyamlalOBrien.H
@@ -57,14 +57,6 @@ class SyamlalOBrien
 :
     public dragModel
 {
-private:
-
-    // Private data
-
-        //- Residual Reynolds number
-        const dimensionedScalar residualRe_;
-
-
 public:
 
     //- Runtime type information
@@ -89,7 +81,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.C
index 21a48c6673..0531e30bf0 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.C
@@ -49,6 +49,7 @@ Foam::dragModels::TomiyamaAnalytic::TomiyamaAnalytic
 )
 :
     dragModel(dict, pair, registerObject),
+    residualRe_("residualRe", dimless, dict.lookup("residualRe")),
     residualEo_("residualEo", dimless, dict.lookup("residualEo")),
     residualE_("residualE", dimless, dict.lookup("residualE"))
 {}
@@ -62,7 +63,8 @@ Foam::dragModels::TomiyamaAnalytic::~TomiyamaAnalytic()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::TomiyamaAnalytic::Cd() const
+Foam::tmp<Foam::volScalarField>
+Foam::dragModels::TomiyamaAnalytic::CdRe() const
 {
     volScalarField Eo(max(pair_.Eo(), residualEo_));
     volScalarField E(max(pair_.E(), residualE_));
@@ -73,12 +75,14 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::TomiyamaAnalytic::Cd() const
     volScalarField F((asin(rtOmEsq) - E*rtOmEsq)/OmEsq);
 
     return
-        (8.0/3.0)*Eo
+        (8.0/3.0)
+       *Eo
        /(
             Eo*pow(E, 2.0/3.0)/OmEsq
           + 16*pow(E, 0.75)
         )
-       /sqr(F);
+       /sqr(F)
+       *max(pair_.Re(), residualRe_);
 }
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.H
index 6c5b72fdca..237f6e9447 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaAnalytic/TomiyamaAnalytic.H
@@ -68,6 +68,9 @@ private:
 
     // Private data
 
+        //- Residual Reynolds Number
+        const dimensionedScalar residualRe_;
+
         //- Residual Eotvos number
         const dimensionedScalar residualEo_;
 
@@ -99,7 +102,7 @@ public:
     // Member Functions
 
         // Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.C
index d01af7ac9d..50fdc667f0 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.C
@@ -49,7 +49,6 @@ Foam::dragModels::TomiyamaCorrelated::TomiyamaCorrelated
 )
 :
     dragModel(dict, pair, registerObject),
-    residualRe_("residualRe", dimless, dict.lookup("residualRe")),
     A_("A", dimless, dict.lookup("A"))
 {}
 
@@ -63,21 +62,21 @@ Foam::dragModels::TomiyamaCorrelated::~TomiyamaCorrelated()
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 Foam::tmp<Foam::volScalarField>
-Foam::dragModels::TomiyamaCorrelated::Cd() const
+Foam::dragModels::TomiyamaCorrelated::CdRe() const
 {
-    volScalarField Re(pair_.Re() + residualRe_);
+    volScalarField Re(pair_.Re());
     volScalarField Eo(pair_.Eo());
 
     return
         max
         (
-            A_/Re
+            A_
            *min
             (
                 (1 + 0.15*pow(Re, 0.687)),
                 scalar(3)
             ),
-            8*Eo/(3*Eo + 12)
+            8*Eo*Re/(3*Eo + 12)
         );
 
 }
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.H
index c779d3262e..f3b37e121d 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/TomiyamaCorrelated/TomiyamaCorrelated.H
@@ -68,9 +68,6 @@ private:
 
     // Private data
 
-        //- Residual Reynolds number
-        const dimensionedScalar residualRe_;
-
         //- Coefficient
         const dimensionedScalar A_;
 
@@ -99,7 +96,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
index 41aeb231d6..34eab7db37 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.C
@@ -61,18 +61,18 @@ Foam::dragModels::WenYu::~WenYu()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::WenYu::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::WenYu::CdRe() const
 {
     volScalarField alpha2(max(scalar(1) - pair_.dispersed(), residualAlpha_));
-    volScalarField Re(max(pair_.Re(), residualRe_));
-    volScalarField Cds
+    volScalarField Re(pair_.Re());
+    volScalarField CdsRe
     (
-        neg(Re - 1000)*(24.0*(1.0 + 0.15*pow(Re, 0.687))/Re)
-      + pos(Re - 1000)*0.44
+        neg(Re - 1000)*24.0*(1.0 + 0.15*pow(Re, 0.687))
+      + pos(Re - 1000)*0.44*max(Re, residualRe_)
     );
 
     return
-        Cds
+        CdsRe
        *pow(alpha2, -2.65)
        *max(pair_.continuous(), residualAlpha_);
 }
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.H
index 48df2ae526..da3cac2cfc 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/WenYu/WenYu.H
@@ -72,7 +72,7 @@ private:
 
     // Private data
 
-        //- Residual Reynolds number
+        //- Residual Reynolds Number
         const dimensionedScalar residualRe_;
 
 
@@ -100,7 +100,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 };
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.C
index 754931df12..ce2a53a0bf 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.C
@@ -90,8 +90,7 @@ Foam::dragModel::dragModel
             pair
         )
     ),
-    residualAlpha_("residualAlpha", dimless, dict.lookup("residualAlpha")),
-    residualSlip_("residualSlip", dimVelocity, dict.lookup("residualSlip"))
+    residualAlpha_("residualAlpha", dimless, dict.lookup("residualAlpha"))
 {}
 
 
@@ -107,14 +106,14 @@ Foam::tmp<Foam::volScalarField> Foam::dragModel::K() const
 {
     return
         0.75
-       *Cd()
+       *CdRe()
        *swarmCorrection_->Cs()
        *pair_.continuous().rho()
+       *pair_.continuous().nu()
        /(
             max(pair_.continuous(), residualAlpha_)
-           *pair_.dispersed().d()
-        )
-       *max(pair_.magUr(), residualSlip_);
+           *sqr(pair_.dispersed().d())
+        );
 }
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.H
index ef0462e735..954d60c396 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/dragModel/dragModel.H
@@ -68,9 +68,6 @@ protected:
         //- Residual phase fraction
         const dimensionedScalar residualAlpha_;
 
-        //- Residual slip velocity
-        const dimensionedScalar residualSlip_;
-
 
 public:
 
@@ -134,7 +131,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const = 0;
+        virtual tmp<volScalarField> CdRe() const = 0;
 
         //- The drag function K used in the momentum equation
         //    ddt(alpha1*rho1*U1) + ... = ... alpha1*alpha2*K*(U1-U2)
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.C
index f57979c349..d955c9e796 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.C
@@ -60,7 +60,7 @@ Foam::dragModels::noDrag::~noDrag()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::noDrag::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::noDrag::CdRe() const
 {
     const fvMesh& mesh(this->pair_.phase1().mesh());
 
@@ -84,7 +84,7 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::noDrag::Cd() const
 
 Foam::tmp<Foam::volScalarField> Foam::dragModels::noDrag::K() const
 {
-    return Cd()*dimensionedScalar("zero", dimensionSet(1, -3, -1, 0, 0), 0);
+    return CdRe()*dimensionedScalar("zero", dimensionSet(1, -3, -1, 0, 0), 0);
 }
 
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.H
index 74c30e6745..b24267ee1b 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/noDrag/noDrag.H
@@ -78,7 +78,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 
         //- The drag function used in the momentum equation
         virtual tmp<volScalarField> K() const;
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.C
index 08b59fbf6b..92b82e537d 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.C
@@ -50,7 +50,6 @@ Foam::dragModels::segregated::segregated
 )
 :
     dragModel(dict, pair, registerObject),
-    residualRe_("residualRe", dimless, dict.lookup("residualRe")),
     m_("m", dimless, dict.lookup("m")),
     n_("n", dimless, dict.lookup("n"))
 {}
@@ -64,9 +63,9 @@ Foam::dragModels::segregated::~segregated()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-Foam::tmp<Foam::volScalarField> Foam::dragModels::segregated::Cd() const
+Foam::tmp<Foam::volScalarField> Foam::dragModels::segregated::CdRe() const
 {
-    FatalErrorIn("Foam::dragModels::segregated::Cd() const")
+    FatalErrorIn("Foam::dragModels::segregated::CdRe() const")
         << "Not implemented."
         << "Drag coefficient not defined for the segregated model."
         << exit(FatalError);
@@ -137,16 +136,12 @@ Foam::tmp<Foam::volScalarField> Foam::dragModels::segregated::K() const
 
     volScalarField ReI
     (
-        max
-        (
-            pair_.rho()
-           *pair_.magUr()
-           /(
-               magGradI
-              *max(alpha1*alpha2, sqr(residualAlpha_))
-              *muI
-            ),
-            residualRe_
+        pair_.rho()
+       *pair_.magUr()
+       /(
+           magGradI
+          *max(alpha1*alpha2, sqr(residualAlpha_))
+          *muI
         )
     );
 
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.H
index 1b68d94d25..13bf2a4d10 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/dragModels/segregated/segregated.H
@@ -66,9 +66,6 @@ private:
 
     // Private data
 
-        //- Residual reynolds number
-        const dimensionedScalar residualRe_;
-
         //- M coefficient
         const dimensionedScalar m_;
 
@@ -100,7 +97,7 @@ public:
     // Member Functions
 
         //- Drag coefficient
-        virtual tmp<volScalarField> Cd() const;
+        virtual tmp<volScalarField> CdRe() const;
 
         //- The drag function used in the momentum equation
         virtual tmp<volScalarField> K() const;
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.C
index 019926b34b..1dd2127338 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.C
@@ -47,8 +47,7 @@ Foam::heatTransferModels::RanzMarshall::RanzMarshall
     const phasePair& pair
 )
 :
-    heatTransferModel(dict, pair),
-    residualRe_("residualRe", dimless, dict.lookup("residualRe"))
+    heatTransferModel(dict, pair)
 {}
 
 
@@ -63,8 +62,7 @@ Foam::heatTransferModels::RanzMarshall::~RanzMarshall()
 Foam::tmp<Foam::volScalarField>
 Foam::heatTransferModels::RanzMarshall::K() const
 {
-    volScalarField Re(pair_.Re() + residualRe_);
-    volScalarField Nu(scalar(2) + 0.6*sqrt(Re)*cbrt(pair_.Pr()));
+    volScalarField Nu(scalar(2) + 0.6*pair_.Re()*cbrt(pair_.Pr()));
 
     return 6.0*pair_.continuous().kappa()*Nu/sqr(pair_.dispersed().d());
 }
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.H b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.H
index 0188d2e210..a52ff75b79 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/heatTransferModels/RanzMarshall/RanzMarshall.H
@@ -54,14 +54,6 @@ class RanzMarshall
 :
     public heatTransferModel
 {
-private:
-
-    // Private data
-        
-        //- Residual Reynolds number
-        const dimensionedScalar residualRe_;
-
-
 public:
 
     //- Runtime type information
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/turbulentDispersionModels/Gosman/Gosman.C b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/turbulentDispersionModels/Gosman/Gosman.C
index 33a144df68..eccbf4e3f0 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/turbulentDispersionModels/Gosman/Gosman.C
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/interfacialModels/turbulentDispersionModels/Gosman/Gosman.C
@@ -84,12 +84,12 @@ Foam::turbulentDispersionModels::Gosman::F() const
 
     return
       - 0.75
-       *drag.Cd()
+       *drag.CdRe()
+       *pair_.continuous().nu()
        *pair_.continuous().turbulence().nut()
-       *pair_.magUr()
        /(
             sigma_
-           *pair_.dispersed().d()
+           *sqr(pair_.dispersed().d())
         )
        *pair_.continuous().rho()
        *fvc::grad(pair_.dispersed());
diff --git a/src/OpenFOAM/matrices/solution/solution.C b/src/OpenFOAM/matrices/solution/solution.C
index ace5564e97..f3251d3c50 100644
--- a/src/OpenFOAM/matrices/solution/solution.C
+++ b/src/OpenFOAM/matrices/solution/solution.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,7 +35,7 @@ License
 
 namespace Foam
 {
-defineDebugSwitchWithName(solution, "solution", 0);
+    defineDebugSwitchWithName(solution, "solution", 0);
 }
 
 // List of sub-dictionaries to rewrite
@@ -60,10 +60,17 @@ void Foam::solution::read(const dictionary& dict)
     if (dict.found("relaxationFactors"))
     {
         const dictionary& relaxDict(dict.subDict("relaxationFactors"));
-        if (relaxDict.found("fields"))
+        if (relaxDict.found("fields") || relaxDict.found("equations"))
         {
-            fieldRelaxDict_ = relaxDict.subDict("fields");
-            eqnRelaxDict_ = relaxDict.subDict("equations");
+            if (relaxDict.found("fields"))
+            {
+                fieldRelaxDict_ = relaxDict.subDict("fields");
+            }
+
+            if (relaxDict.found("equations"))
+            {
+                eqnRelaxDict_ = relaxDict.subDict("equations");
+            }
         }
         else
         {
diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.C b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.C
index 137339e631..927c0e80b7 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -98,6 +98,22 @@ Foam::processorCyclicPolyPatch::processorCyclicPolyPatch
 {}
 
 
+Foam::processorCyclicPolyPatch::processorCyclicPolyPatch
+(
+    const processorCyclicPolyPatch& pp,
+    const polyBoundaryMesh& bm,
+    const label index,
+    const label newSize,
+    const label newStart
+)
+:
+    processorPolyPatch(pp, bm, index, newSize, newStart),
+    referPatchName_(pp.referPatchName_),
+    tag_(pp.tag()),
+    referPatchID_(-1)
+{}
+
+
 Foam::processorCyclicPolyPatch::processorCyclicPolyPatch
 (
     const processorCyclicPolyPatch& pp,
diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H
index 76a833ba61..2bfc2c9830 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H
+++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -145,6 +145,17 @@ public:
         //- Construct as given the original patch and resetting the
         //  face list and boundary mesh information
         processorCyclicPolyPatch
+        (
+            const processorCyclicPolyPatch& pp,
+            const polyBoundaryMesh& bm,
+            const label index,
+            const label newSize,
+            const label newStart
+        );
+
+        //- Construct as given the original patch and resetting the
+        //  face list, boundary mesh information and referPatch
+        processorCyclicPolyPatch
         (
             const processorCyclicPolyPatch& pp,
             const polyBoundaryMesh& bm,
@@ -171,6 +182,29 @@ public:
             return autoPtr<polyPatch>(new processorCyclicPolyPatch(*this, bm));
         }
 
+        //- Construct and return a clone, resetting the face list
+        //  and boundary mesh
+        virtual autoPtr<polyPatch> clone
+        (
+            const polyBoundaryMesh& bm,
+            const label index,
+            const label newSize,
+            const label newStart
+        ) const
+        {
+            return autoPtr<polyPatch>
+            (
+                new processorCyclicPolyPatch
+                (
+                    *this,
+                    bm,
+                    index,
+                    newSize,
+                    newStart
+                )
+            );
+        }
+
         //- Construct and return a clone, resetting the face list
         //  and boundary mesh
         virtual autoPtr<polyPatch> clone
diff --git a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/symmetryPlane/symmetryPlanePolyPatch.C b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/symmetryPlane/symmetryPlanePolyPatch.C
index d3c78904b9..63268e6832 100644
--- a/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/symmetryPlane/symmetryPlanePolyPatch.C
+++ b/src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/symmetryPlane/symmetryPlanePolyPatch.C
@@ -42,29 +42,38 @@ namespace Foam
 
 void Foam::symmetryPlanePolyPatch::calcGeometry(PstreamBuffers&)
 {
-    if (returnReduce(size(), sumOp<label>()) == 0)
+    if (n_ == vector::rootMax)
     {
-        // No faces in patch. Avoid gAverage complaining and set
-        // normal to nonsense value to catch any use
-        n_ = vector::rootMax;
-    }
-    else
-    {
-        const vectorField& nf(faceNormals());
-        n_ = gAverage(nf);
-
-        // Check the symmetry plane is planar
-        forAll(nf, facei)
+        if (returnReduce(size(), sumOp<label>()))
         {
-            if (magSqr(n_ - nf[facei]) > SMALL)
+            const vectorField& nf(faceNormals());
+            n_ = gAverage(nf);
+
+            if (debug)
             {
-                FatalErrorIn("symmetryPlanePolyPatch::n()")
-                    << "Symmetry plane '" << name() << "' is not planar."
-                    << endl
-                    << " Either split the patch into planar parts"
-                    << " or use the " << symmetryPolyPatch::typeName
-                    << " patch type"
-                    << exit(FatalError);
+                Info<< "Patch " << name() << " calculated average normal "
+                    << n_ << endl;
+            }
+
+
+            // Check the symmetry plane is planar
+            forAll(nf, facei)
+            {
+                if (magSqr(n_ - nf[facei]) > SMALL)
+                {
+                    FatalErrorIn("symmetryPlanePolyPatch::n()")
+                        << "Symmetry plane '" << name() << "' is not planar."
+                        << endl
+                        << "At local face at "
+                        << primitivePatch::faceCentres()[facei]
+                        << " the normal " << nf[facei]
+                        << " differs from the average normal " << n_
+                        << " by " << magSqr(n_ - nf[facei]) << endl
+                        << "Either split the patch into planar parts"
+                        << " or use the " << symmetryPolyPatch::typeName
+                        << " patch type"
+                        << exit(FatalError);
+                }
             }
         }
     }
@@ -109,7 +118,7 @@ Foam::symmetryPlanePolyPatch::symmetryPlanePolyPatch
 )
 :
     polyPatch(pp, bm),
-    n_(vector::rootMax)
+    n_(pp.n_)
 {}
 
 
@@ -123,7 +132,7 @@ Foam::symmetryPlanePolyPatch::symmetryPlanePolyPatch
 )
 :
     polyPatch(pp, bm, index, newSize, newStart),
-    n_(vector::rootMax)
+    n_(pp.n_)
 {}
 
 
@@ -137,7 +146,7 @@ Foam::symmetryPlanePolyPatch::symmetryPlanePolyPatch
 )
 :
     polyPatch(pp, bm, index, mapAddressing, newStart),
-    n_(vector::rootMax)
+    n_(pp.n_)
 {}
 
 
diff --git a/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C b/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C
index bfc6b567f4..92b7117169 100644
--- a/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C
+++ b/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -945,6 +945,56 @@ void Foam::dynamicRefineFvMesh::extendMarkedCells
 }
 
 
+void Foam::dynamicRefineFvMesh::checkEightAnchorPoints
+(
+    PackedBoolList& protectedCell,
+    label& nProtected
+) const
+{
+    const labelList& cellLevel = meshCutter_.cellLevel();
+    const labelList& pointLevel = meshCutter_.pointLevel();
+
+    labelList nAnchorPoints(nCells(), 0);
+
+    forAll(pointLevel, pointI)
+    {
+        const labelList& pCells = pointCells(pointI);
+
+        forAll(pCells, pCellI)
+        {
+            label cellI = pCells[pCellI];
+
+            if (pointLevel[pointI] <= cellLevel[cellI])
+            {
+                // Check if cell has already 8 anchor points -> protect cell
+                if (nAnchorPoints[cellI] == 8)
+                {
+                    if (protectedCell.set(cellI, true))
+                    {
+                        nProtected++;
+                    }
+                }
+
+                if (!protectedCell[cellI])
+                {
+                    nAnchorPoints[cellI]++;
+                }
+            }
+        }
+    }
+
+
+    forAll(protectedCell, cellI)
+    {
+        if (!protectedCell[cellI] && nAnchorPoints[cellI] != 8)
+        {
+            protectedCell.set(cellI, true);
+            nProtected++;
+        }
+    }
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::dynamicRefineFvMesh::dynamicRefineFvMesh(const IOobject& io)
@@ -1067,6 +1117,37 @@ Foam::dynamicRefineFvMesh::dynamicRefineFvMesh(const IOobject& io)
                 nProtected++;
             }
         }
+
+        // Also protect any cells that are less than hex
+        forAll(cells(), cellI)
+        {
+            const cell& cFaces = cells()[cellI];
+
+            if (cFaces.size() < 6)
+            {
+                if (protectedCell_.set(cellI, 1))
+                {
+                    nProtected++;
+                }
+            }
+            else
+            {
+                forAll(cFaces, cFaceI)
+                {
+                    if (faces()[cFaces[cFaceI]].size() < 4)
+                    {
+                        if (protectedCell_.set(cellI, 1))
+                        {
+                            nProtected++;
+                        }
+                        break;
+                    }
+                }
+            }
+        }
+
+        // Check cells for 8 corner points
+        checkEightAnchorPoints(protectedCell_, nProtected);
     }
 
     if (returnReduce(nProtected, sumOp<label>()) == 0)
diff --git a/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H b/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H
index 287e26ca7b..7611b86c09 100644
--- a/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H
+++ b/src/dynamicFvMesh/dynamicRefineFvMesh/dynamicRefineFvMesh.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -179,6 +179,12 @@ protected:
             //- Extend markedCell with cell-face-cell.
             void extendMarkedCells(PackedBoolList& markedCell) const;
 
+            //- Check all cells have 8 anchor points
+            void checkEightAnchorPoints
+            (
+                PackedBoolList& protectedCell,
+                label& nProtected
+            ) const;
 
 private:
 
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C
index 1f07b1bb14..20fb99a0b3 100644
--- a/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C
+++ b/src/finiteVolume/fields/fvPatchFields/derived/uniformFixedValue/uniformFixedValueFvPatchField.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2014 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -66,10 +66,10 @@ uniformFixedValueFvPatchField<Type>::uniformFixedValueFvPatchField
     const fvPatchFieldMapper& mapper
 )
 :
-    fixedValueFvPatchField<Type>(ptf, p, iF, mapper),
+    fixedValueFvPatchField<Type>(p, iF),  // bypass mapper
     uniformValue_(ptf.uniformValue_().clone().ptr())
 {
-    // For safety re-evaluate
+    // Evaluate since value not mapped
     const scalar t = this->db().time().timeOutputValue();
     fvPatchField<Type>::operator==(uniformValue_->value(t));
 }
diff --git a/tutorials/DNS/dnsFoam/boxTurb16/system/fvSolution b/tutorials/DNS/dnsFoam/boxTurb16/system/fvSolution
index f35d60d17c..06aba48225 100644
--- a/tutorials/DNS/dnsFoam/boxTurb16/system/fvSolution
+++ b/tutorials/DNS/dnsFoam/boxTurb16/system/fvSolution
@@ -27,8 +27,8 @@ solvers
 
     U
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
diff --git a/tutorials/compressible/rhoLTSPimpleFoam/angledDuct/system/fvSolution b/tutorials/compressible/rhoLTSPimpleFoam/angledDuct/system/fvSolution
index e4c1f816b5..db290d427b 100644
--- a/tutorials/compressible/rhoLTSPimpleFoam/angledDuct/system/fvSolution
+++ b/tutorials/compressible/rhoLTSPimpleFoam/angledDuct/system/fvSolution
@@ -34,8 +34,8 @@ solvers
 
     "(rho|U|h|k|epsilon|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoPimpleFoam/ras/angledDuct/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/ras/angledDuct/system/fvSolution
index d87df4b492..c89dc18278 100644
--- a/tutorials/compressible/rhoPimpleFoam/ras/angledDuct/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/ras/angledDuct/system/fvSolution
@@ -34,8 +34,8 @@ solvers
 
     "(rho|U|h|k|epsilon|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoPimpleFoam/ras/cavity/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/ras/cavity/system/fvSolution
index 9f65b0ec80..77cea8e603 100644
--- a/tutorials/compressible/rhoPimpleFoam/ras/cavity/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/ras/cavity/system/fvSolution
@@ -40,8 +40,8 @@ solvers
 
     "(U|h|R|k|epsilon|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoPimpleFoam/ras/mixerVessel2D/system/fvSolution b/tutorials/compressible/rhoPimpleFoam/ras/mixerVessel2D/system/fvSolution
index 5bd627be1b..825a3a30b2 100644
--- a/tutorials/compressible/rhoPimpleFoam/ras/mixerVessel2D/system/fvSolution
+++ b/tutorials/compressible/rhoPimpleFoam/ras/mixerVessel2D/system/fvSolution
@@ -34,16 +34,16 @@ solvers
 
     U
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-5;
         relTol          0.1;
     }
 
     UFinal
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-5;
         relTol          0;
     }
@@ -74,8 +74,8 @@ solvers
 
     "(U|h|k|epsilon)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-5;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoPimplecFoam/angledDuct/system/fvSolution b/tutorials/compressible/rhoPimplecFoam/angledDuct/system/fvSolution
index 3f77479c31..a89d48787c 100644
--- a/tutorials/compressible/rhoPimplecFoam/angledDuct/system/fvSolution
+++ b/tutorials/compressible/rhoPimplecFoam/angledDuct/system/fvSolution
@@ -34,8 +34,8 @@ solvers
 
     "(rho|U|h|k|epsilon|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoPorousSimpleFoam/angledDuctExplicit/system/fvSolution b/tutorials/compressible/rhoPorousSimpleFoam/angledDuctExplicit/system/fvSolution
index a61ae9a85d..961f865211 100644
--- a/tutorials/compressible/rhoPorousSimpleFoam/angledDuctExplicit/system/fvSolution
+++ b/tutorials/compressible/rhoPorousSimpleFoam/angledDuctExplicit/system/fvSolution
@@ -40,8 +40,8 @@ solvers
 
     e
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/system/fvSolution b/tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/system/fvSolution
index 13a08e5e1b..b22f827c12 100644
--- a/tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/system/fvSolution
+++ b/tutorials/compressible/rhoPorousSimpleFoam/angledDuctImplicit/system/fvSolution
@@ -31,8 +31,8 @@ solvers
 
     e
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvSolution b/tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvSolution
index c259abc052..290de2d85d 100644
--- a/tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvSolution
+++ b/tutorials/compressible/rhoSimpleFoam/angledDuctExplicitFixedCoeff/system/fvSolution
@@ -40,8 +40,8 @@ solvers
 
     e
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0.1;
     }
diff --git a/tutorials/compressible/sonicFoam/laminar/forwardStep/system/fvSolution b/tutorials/compressible/sonicFoam/laminar/forwardStep/system/fvSolution
index e07aa4bafc..2875f8c7aa 100644
--- a/tutorials/compressible/sonicFoam/laminar/forwardStep/system/fvSolution
+++ b/tutorials/compressible/sonicFoam/laminar/forwardStep/system/fvSolution
@@ -19,8 +19,8 @@ solvers
 {
     "p.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }
diff --git a/tutorials/compressible/sonicFoam/laminar/shockTube/system/fvSolution b/tutorials/compressible/sonicFoam/laminar/shockTube/system/fvSolution
index 9ed6e1a982..11b41801ae 100644
--- a/tutorials/compressible/sonicFoam/laminar/shockTube/system/fvSolution
+++ b/tutorials/compressible/sonicFoam/laminar/shockTube/system/fvSolution
@@ -19,8 +19,8 @@ solvers
 {
     "(p|U|e).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-15;
         relTol          0;
     }
diff --git a/tutorials/compressible/sonicFoam/ras/prism/system/controlDict b/tutorials/compressible/sonicFoam/ras/prism/system/controlDict
index 4c4191a48c..f5dcfb9881 100644
--- a/tutorials/compressible/sonicFoam/ras/prism/system/controlDict
+++ b/tutorials/compressible/sonicFoam/ras/prism/system/controlDict
@@ -17,7 +17,7 @@ FoamFile
 
 application     sonicFoam;
 
-startFrom       latestTime;
+startFrom       startTime;
 
 startTime       0;
 
diff --git a/tutorials/compressible/sonicFoam/ras/prism/system/fvSolution b/tutorials/compressible/sonicFoam/ras/prism/system/fvSolution
index 56886cc6e1..b1a21b1f9a 100644
--- a/tutorials/compressible/sonicFoam/ras/prism/system/fvSolution
+++ b/tutorials/compressible/sonicFoam/ras/prism/system/fvSolution
@@ -24,8 +24,8 @@ solvers
 
     "p.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/compressible/sonicLiquidFoam/decompressionTank/system/fvSolution b/tutorials/compressible/sonicLiquidFoam/decompressionTank/system/fvSolution
index 5f012bcea7..438b2dd3c4 100644
--- a/tutorials/compressible/sonicLiquidFoam/decompressionTank/system/fvSolution
+++ b/tutorials/compressible/sonicLiquidFoam/decompressionTank/system/fvSolution
@@ -19,8 +19,8 @@ solvers
 {
     "p.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }
diff --git a/tutorials/electromagnetics/electrostaticFoam/chargedWire/system/fvSolution b/tutorials/electromagnetics/electrostaticFoam/chargedWire/system/fvSolution
index de29f27df6..bb403e3cf2 100644
--- a/tutorials/electromagnetics/electrostaticFoam/chargedWire/system/fvSolution
+++ b/tutorials/electromagnetics/electrostaticFoam/chargedWire/system/fvSolution
@@ -27,8 +27,8 @@ solvers
 
     rho
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.2;
     }
diff --git a/tutorials/electromagnetics/mhdFoam/hartmann/system/fvSolution b/tutorials/electromagnetics/mhdFoam/hartmann/system/fvSolution
index 8e277b7c41..0f31523126 100644
--- a/tutorials/electromagnetics/mhdFoam/hartmann/system/fvSolution
+++ b/tutorials/electromagnetics/mhdFoam/hartmann/system/fvSolution
@@ -34,8 +34,8 @@ solvers
 
     "(U|B)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
diff --git a/tutorials/heatTransfer/buoyantSimpleFoam/externalCoupledCavity/externalSolver b/tutorials/heatTransfer/buoyantSimpleFoam/externalCoupledCavity/externalSolver
index 7b511f5dd9..03c25c644a 100755
--- a/tutorials/heatTransfer/buoyantSimpleFoam/externalCoupledCavity/externalSolver
+++ b/tutorials/heatTransfer/buoyantSimpleFoam/externalCoupledCavity/externalSolver
@@ -18,12 +18,12 @@ timeOut=10
 refGrad=0
 valueFraction=1
 
-function log()
+log()
 {
     echo "External: $@"
 }
 
-function init()
+init()
 {
     log "initialisation: creating ${dataFile}.in"
 
@@ -54,7 +54,7 @@ step=0
 while [ 1 ]; do
     if [ -f $lockFile ]; then
         log "found lock file ${lockFile} - waiting"
-        ((totalWait += waitSec))
+        totalWait=$(expr $totalWait + $waitSec)
         if [ $totalWait -gt $timeOut ]; then
             log "timeout"
             break
@@ -63,7 +63,7 @@ while [ 1 ]; do
         fi
     else
         totalWait=0
-        ((step++))
+        step=$(expr $step + 1)
         log "step $step"
         log "lock not present - taking control"
 
diff --git a/tutorials/mesh/foamyHexMesh/mixerVessel/Alltest b/tutorials/mesh/foamyHexMesh/mixerVessel/Alltest
new file mode 100755
index 0000000000..8fb0604162
--- /dev/null
+++ b/tutorials/mesh/foamyHexMesh/mixerVessel/Alltest
@@ -0,0 +1,9 @@
+#!/bin/sh
+cd ${0%/*} || exit 1    # run from this directory
+
+# Source tutorial run functions
+. $WM_PROJECT_DIR/bin/tools/RunFunctions
+
+# Do nothing for now.
+
+# ----------------------------------------------------------------- end-of-file
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/U b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/U
similarity index 85%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/0/U
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/U
index 88c7022844..507cfa947d 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/U
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/U
@@ -15,14 +15,18 @@ FoamFile
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-Umean 1;
+Umean 1.668;
+mUmean -1.668;
 
 dimensions      [0 1 -1 0 0 0 0];
 
-internalField   uniform ($Umean 0 0);
+internalField   uniform ($mUmean 0 0);
 
 boundaryField
 {
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
+
     inlet
     {
         type            fixedValue;
@@ -37,21 +41,15 @@ boundaryField
         value           $internalField;
     }
 
-    sides
-    {
-        type            symmetry;
-    }
-
     atmosphere
     {
         type            pressureInletOutletVelocity;
-        tangentialVelocity $internalField;
         value           uniform (0 0 0);
     }
 
     hull
     {
-        type            fixedValue;
+        type            movingWallVelocity;
         value           uniform (0 0 0);
     }
 }
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/alpha.water.org b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/alpha.water
similarity index 79%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/0/alpha.water.org
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/alpha.water
index 6f60c040ca..06ee368c9c 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/alpha.water.org
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/alpha.water
@@ -11,7 +11,7 @@ FoamFile
     format      ascii;
     class       volScalarField;
     location    "0";
-    object      alpha.water;
+    object      alpha;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -21,6 +21,9 @@ internalField   uniform 0;
 
 boundaryField
 {
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
+
     inlet
     {
         type            fixedValue;
@@ -29,18 +32,16 @@ boundaryField
 
     outlet
     {
-        type            zeroGradient;
-    }
-
-    sides
-    {
-        type            symmetry;
+        type            variableHeightFlowRate;
+        lowerBound      0;
+        upperBound      1;
+        value           $internalField;
     }
 
     atmosphere
     {
         type            inletOutlet;
-        inletValue      uniform 0;
+        inletValue      $internalField;
         value           $internalField;
     }
 
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/k b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/k
similarity index 85%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/0/k
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/k
index 853af9b0ef..392241cd6f 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/k
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/k
@@ -21,6 +21,9 @@ internalField   uniform 0.00015;
 
 boundaryField
 {
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
+
     inlet
     {
         type            fixedValue;
@@ -29,12 +32,9 @@ boundaryField
 
     outlet
     {
-        type            zeroGradient;
-    }
-
-    sides
-    {
-        type            symmetry;
+       type            inletOutlet;
+       inletValue      $internalField;
+       value           $internalField;
     }
 
     atmosphere
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/nut b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/nut
similarity index 83%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/0/nut
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/nut
index 87b4bf3c28..4abe529d46 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/nut
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/nut
@@ -21,6 +21,9 @@ internalField   uniform 5e-07;
 
 boundaryField
 {
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
+
     inlet
     {
         type            fixedValue;
@@ -32,11 +35,6 @@ boundaryField
         type            zeroGradient;
     }
 
-    sides
-    {
-        type            symmetry;
-    }
-
     atmosphere
     {
         type            zeroGradient;
@@ -44,7 +42,9 @@ boundaryField
 
     hull
     {
-        type            nutkWallFunction;
+        type            nutkRoughWallFunction;
+        Ks              uniform 100e-6;
+        Cs              uniform 0.5;
         value           $internalField;
     }
 }
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/omega b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/omega
similarity index 91%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/0/omega
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/omega
index a515138656..0dde05c3ac 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/omega
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/omega
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      omega;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -20,6 +21,9 @@ internalField   uniform 2;
 
 boundaryField
 {
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
+
     inlet
     {
         type            fixedValue;
@@ -33,16 +37,6 @@ boundaryField
        value           $internalField;
     }
 
-    outlet
-    {
-        type            zeroGradient;
-    }
-
-    sides
-    {
-        type            symmetry;
-    }
-
     atmosphere
     {
         type            inletOutlet;
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/p_rgh b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/p_rgh
similarity index 91%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/0/p_rgh
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/p_rgh
index 027c48884c..0c0eda7c6f 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/0/p_rgh
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/p_rgh
@@ -10,6 +10,7 @@ FoamFile
     version     2.0;
     format      ascii;
     class       volScalarField;
+    location    "0";
     object      p_rgh;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -20,6 +21,9 @@ internalField   uniform 0;
 
 boundaryField
 {
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
+
     inlet
     {
         type            fixedFluxPressure;
@@ -31,11 +35,6 @@ boundaryField
         type            zeroGradient;
     }
 
-    sides
-    {
-        type            symmetry;
-    }
-
     atmosphere
     {
         type            totalPressure;
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/polyMesh/blockMeshDict b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/pointDisplacement
similarity index 53%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/constant/polyMesh/blockMeshDict
rename to tutorials/multiphase/LTSInterFoam/DTCHull/0.org/pointDisplacement
index 989ecac8b0..350676e605 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/polyMesh/blockMeshDict
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/0.org/pointDisplacement
@@ -9,74 +9,44 @@ FoamFile
 {
     version     2.0;
     format      ascii;
-    class       dictionary;
-    object      blockMeshDict;
+    class       pointVectorField;
+    location    "0";
+    object      pointDisplacement;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-convertToMeters 1;
+dimensions      [0 1 0 0 0 0 0];
 
-vertices
-(
-    (-0.75 -0.75 -0.5)
-    ( 1.25 -0.75 -0.5)
-    ( 1.25  0.75 -0.5)
-    (-0.75  0.75 -0.5)
-    (-0.75 -0.75  0.0399)
-    ( 1.25 -0.75 0.0399)
-    ( 1.25  0.75  0.0399)
-    (-0.75  0.75  0.0399)
-);
+internalField   uniform (0 0 0);
 
-blocks
-(
-    hex (0 1 2 3 4 5 6 7) (40 30 10) simpleGrading (1 1 1)
-);
+boundaryField
+{
+    //- Set patchGroups for constraint patches
+    #include "${WM_PROJECT_DIR}/etc/caseDicts/setConstraintTypes"
 
-edges
-(
-);
-
-boundary
-(
     inlet
     {
-        type patch;
-        faces
-        (
-            (0 4 7 3)
-        );
+        type            fixedValue;
+        value           uniform (0 0 0);
     }
 
     outlet
     {
-        type patch;
-        faces
-        (
-            (2 6 5 1)
-        );
+        type            fixedValue;
+        value           uniform (0 0 0);
     }
 
     atmosphere
     {
-        type patch;
-        faces
-        (
-            (4 5 6 7)
-        );
+        type            fixedValue;
+        value           uniform (0 0 0);
     }
 
-    sides
+    hull
     {
-        type symmetry;
-        faces
-        (
-            (1 5 4 0)
-            (3 7 6 2)
-            (0 3 2 1)
-        );
+        type            calculated;
     }
-);
+}
 
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/Allclean b/tutorials/multiphase/LTSInterFoam/DTCHull/Allclean
new file mode 100755
index 0000000000..7ae6e656d7
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/Allclean
@@ -0,0 +1,17 @@
+#!/bin/sh
+cd ${0%/*} || exit 1    # run from this directory
+
+# Source tutorial clean functions
+. $WM_PROJECT_DIR/bin/tools/CleanFunctions
+
+# remove surface
+\rm -f constant/triSurface/DTC-scaled.stl.gz > /dev/null 2>&1
+\rm -rf constant/extendedFeatureEdgeMesh > /dev/null 2>&1
+\rm -f constant/triSurface/DTC-scaled.eMesh > /dev/null 2>&1
+
+cleanCase
+\rm constant/polyMesh/boundary > /dev/null 2>&1
+\rm -rf 0 > /dev/null 2>&1
+\rm system/topoSetDict > /dev/null 2>&1
+
+# ----------------------------------------------------------------- end-of-file
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/Allrun b/tutorials/multiphase/LTSInterFoam/DTCHull/Allrun
new file mode 100755
index 0000000000..b4f5560174
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/Allrun
@@ -0,0 +1,48 @@
+#!/bin/sh
+cd ${0%/*} || exit 1    # run from this directory
+
+# Source tutorial run functions
+. $WM_PROJECT_DIR/bin/tools/RunFunctions
+
+runRefineMesh ()
+{
+    echo "Running refineMesh on $PWD"
+    refineMesh -dict system/refineMeshDict -overwrite > log.refineMesh 2>&1
+}
+
+# copy DTC hull surface from resources folder
+cp $FOAM_TUTORIALS/resources/geometry/DTC-scaled.stl.gz constant/triSurface/
+
+runApplication surfaceFeatureExtract
+
+runApplication blockMesh
+
+i=1
+if [ -f log.topoSet.6 ] ; then
+   i=7
+fi
+while [ "$i" -lt 7 ] ; do
+    cp system/topoSetDict.${i} system/topoSetDict
+    runApplication topoSet
+    runRefineMesh
+    mv log.topoSet log.topoSet.$i
+    mv log.refineMesh log.refineMesh.$i
+    i=$(expr $i + 1)
+done
+
+runApplication snappyHexMesh -overwrite
+
+\rm -rf 0
+\cp -r 0.org 0
+
+runApplication setFields
+
+runApplication decomposePar
+
+runParallel renumberMesh 8 -overwrite
+
+runParallel $(getApplication) 8
+
+runApplication reconstructPar
+
+# ----------------------------------------------------------------- end-of-file
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/RASProperties b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/RASProperties
similarity index 100%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/constant/RASProperties
rename to tutorials/multiphase/LTSInterFoam/DTCHull/constant/RASProperties
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/constant/dynamicMeshDict b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/dynamicMeshDict
new file mode 100644
index 0000000000..7e7f4969af
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/dynamicMeshDict
@@ -0,0 +1,67 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      dynamicMeshDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+dynamicFvMesh       dynamicMotionSolverFvMesh;
+
+motionSolverLibs    ("libsixDoFRigidBodyMotion.so");
+
+solver              sixDoFRigidBodyMotion;
+
+sixDoFRigidBodyMotionCoeffs
+{
+    patches         (hull);
+    innerDistance   0.3;
+    outerDistance   1;
+
+    centreOfMass    (2.929541 0 0.2);
+    mass            412.73;
+    momentOfInertia (40 921 921);
+    rhoInf          1;
+    report          on;
+    accelerationRelaxation 0.3;
+    value           uniform (0 0 0);
+
+    constraints
+    {
+        zAxis
+        {
+            sixDoFRigidBodyMotionConstraint line;
+            direction     (0 0 1);
+        }
+        yPlane
+        {
+            sixDoFRigidBodyMotionConstraint axis;
+            axis          (0 1 0);
+        }
+    }
+
+    restraints
+    {
+        translationDamper
+        {
+            sixDoFRigidBodyMotionRestraint linearDamper;
+            coeff         8596;
+        }
+        rotationDamper
+        {
+            sixDoFRigidBodyMotionRestraint sphericalAngularDamper;
+            coeff         11586;
+        }
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/g b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/g
similarity index 100%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/constant/g
rename to tutorials/multiphase/LTSInterFoam/DTCHull/constant/g
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/constant/polyMesh/blockMeshDict b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/polyMesh/blockMeshDict
new file mode 100644
index 0000000000..8ac42309c3
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/polyMesh/blockMeshDict
@@ -0,0 +1,147 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      blockMeshDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+convertToMeters 1;
+
+vertices
+(
+    (-26 -19 -16)
+    (16 -19 -16)
+    (16 0 -16)
+    (-26 0 -16)
+
+    (-26 -19 -1)
+    (16 -19 -1)
+    (16 0 -1)
+    (-26 0 -1)
+
+    (-26 -19 0.185)
+    (16 -19 0.185)
+    (16 0 0.185)
+    (-26 0 0.185)
+
+    (-26 -19 0.244)
+    (16 -19 0.244)
+    (16 0 0.244)
+    (-26 0 0.244)
+
+    (-26 -19 0.3)
+    (16 -19 0.3)
+    (16 0 0.3)
+    (-26 0 0.3)
+
+    (-26 -19 1.6)
+    (16 -19 1.6)
+    (16 0 1.6)
+    (-26 0 1.6)
+
+    (-26 -19 4)
+    (16 -19 4)
+    (16 0 4)
+    (-26 0 4)
+);
+
+blocks
+(
+    hex (0 1 2 3 4 5 6 7) (42 19 50) simpleGrading (1 1 0.05)
+    hex (4 5 6 7 8 9 10 11) (42 19 50) simpleGrading (1 1 1)
+    hex (8 9 10 11 12 13 14 15) (42 19 4) simpleGrading (1 1 1)
+    hex (12 13 14 15 16 17 18 19) (42 19 4) simpleGrading (1 1 1)
+    hex (16 17 18 19 20 21 22 23) (42 19 40) simpleGrading (1 1 1)
+    hex (20 21 22 23 24 25 26 27) (42 19 20) simpleGrading (1 1 5)
+);
+
+edges
+(
+);
+
+boundary
+(
+    atmosphere
+    {
+        type patch;
+        faces
+        (
+            (24 25 26 27)
+        );
+    }
+    inlet
+    {
+        type patch;
+        faces
+        (
+            (1 2 6 5)
+            (5 6 10 9)
+            (9 10 14 13)
+            (13 14 18 17)
+            (17 18 22 21)
+            (21 22 26 25)
+        );
+    }
+    outlet
+    {
+        type patch;
+        faces
+        (
+            (0 4 7 3)
+            (4 8 11 7)
+            (8 12 15 11)
+            (12 16 19 15)
+            (16 20 23 19)
+            (20 24 27 23)
+        );
+    }
+    bottom
+    {
+        type symmetryPlane;
+        faces
+        (
+            (0 3 2 1)
+        );
+    }
+    side
+    {
+        type symmetryPlane;
+        faces
+        (
+            (0 1 5 4)
+            (4 5 9 8)
+            (8 9 13 12)
+            (12 13 17 16)
+            (16 17 21 20)
+            (20 21 25 24)
+        );
+    }
+    midPlane
+    {
+        type symmetryPlane;
+        faces
+        (
+            (3 7 6 2)
+            (7 11 10 6)
+            (11 15 14 10)
+            (15 19 18 14)
+            (19 23 22 18)
+            (23 27 26 22)
+        );
+    }
+);
+
+mergePatchPairs
+(
+);
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/transportProperties b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/transportProperties
similarity index 91%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/constant/transportProperties
rename to tutorials/multiphase/LTSInterFoam/DTCHull/constant/transportProperties
index 8c6e7a2008..0bfa2f19b9 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/transportProperties
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/transportProperties
@@ -20,8 +20,8 @@ phases (water air);
 water
 {
     transportModel  Newtonian;
-    nu              nu [ 0 2 -1 0 0 0 0 ]  1e-06;
-    rho             rho [ 1 -3 0 0 0 0 0 ] 1000;
+    nu              nu [ 0 2 -1 0 0 0 0 ]  1.09e-06;
+    rho             rho [ 1 -3 0 0 0 0 0 ] 998.8;
 }
 
 air
@@ -33,5 +33,4 @@ air
 
 sigma           sigma [ 1 0 -2 0 0 0 0 ] 0;
 
-
 // ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/triSurface/README b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/triSurface/README
similarity index 66%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/constant/triSurface/README
rename to tutorials/multiphase/LTSInterFoam/DTCHull/constant/triSurface/README
index de3f29cc01..dbf3b73d5f 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/triSurface/README
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/triSurface/README
@@ -1,4 +1,4 @@
-Directory to house tri-surfaces
+Folder to house tri-surfaces
 
 The Allrun script copies the surface from the $FOAM_TUTORIALS/resources/geometry
-directory
+folder
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/constant/turbulenceProperties b/tutorials/multiphase/LTSInterFoam/DTCHull/constant/turbulenceProperties
similarity index 100%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/constant/turbulenceProperties
rename to tutorials/multiphase/LTSInterFoam/DTCHull/constant/turbulenceProperties
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/controlDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/controlDict
similarity index 75%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/system/controlDict
rename to tutorials/multiphase/LTSInterFoam/DTCHull/system/controlDict
index d2aa8b3c71..c2713242ce 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/controlDict
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/controlDict
@@ -23,7 +23,7 @@ startTime       0;
 
 stopAt          endTime;
 
-endTime         1000;
+endTime         4000;
 
 deltaT          1;
 
@@ -33,11 +33,11 @@ writeInterval   100;
 
 purgeWrite      0;
 
-writeFormat     ascii;
+writeFormat     binary;
 
 writePrecision  6;
 
-writeCompression compressed;
+writeCompression uncompressed;
 
 timeFormat      general;
 
@@ -45,5 +45,22 @@ timePrecision   6;
 
 runTimeModifiable yes;
 
+functions
+{
+    forces
+    {
+        type forces;
+        functionObjectLibs ( "libforces.so" );
+        patches (hull);
+        rhoInf  998.8;
+        rhoName rho;
+        UName   U;
+        log     on;
+        outputControl   timeStep;
+        outputInterval  1;
+        CofR (2.929541 0 0.2);
+    }
+}
+
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/decomposeParDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/decomposeParDict
new file mode 100644
index 0000000000..a932227ac1
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/decomposeParDict
@@ -0,0 +1,58 @@
+/*---------------------------------------------------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      http://www.openfoam.org               |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+
+FoamFile
+{
+    version         2.0;
+    format          ascii;
+
+    root            "";
+    case            "";
+    instance        "";
+    local           "";
+
+    class           dictionary;
+    object          decomposeParDict;
+}
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+numberOfSubdomains 8;
+
+method          scotch;
+
+simpleCoeffs
+{
+    n               (2 1 1);
+    delta           0.001;
+}
+
+hierarchicalCoeffs
+{
+    n               (1 1 1);
+    delta           0.001;
+    order           xyz;
+}
+
+metisCoeffs
+{
+}
+
+manualCoeffs
+{
+    dataFile        "";
+}
+
+distributed     no;
+
+roots
+(
+);
+
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/fvSchemes b/tutorials/multiphase/LTSInterFoam/DTCHull/system/fvSchemes
similarity index 90%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/system/fvSchemes
rename to tutorials/multiphase/LTSInterFoam/DTCHull/system/fvSchemes
index aaf9168e32..0806a381d5 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/fvSchemes
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/fvSchemes
@@ -23,6 +23,7 @@ ddtSchemes
 gradSchemes
 {
     default         Gauss linear;
+    limitedGrad     cellLimited Gauss linear 1;
 }
 
 divSchemes
@@ -30,8 +31,8 @@ divSchemes
     div(rhoPhi,U)   Gauss linearUpwind grad(U);
     div(phi,alpha)  Gauss vanLeer;
     div(phirb,alpha) Gauss linear;
-    div(phi,k)      Gauss upwind;
-    div(phi,omega)  Gauss upwind;
+    div(phi,k)      Gauss linearUpwind limitedGrad;
+    div(phi,omega)  Gauss linearUpwind limitedGrad;
     div((muEff*dev(T(grad(U))))) Gauss linear;
 }
 
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/fvSolution b/tutorials/multiphase/LTSInterFoam/DTCHull/system/fvSolution
similarity index 80%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/system/fvSolution
rename to tutorials/multiphase/LTSInterFoam/DTCHull/system/fvSolution
index df8102dc06..7480f95218 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/fvSolution
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/fvSolution
@@ -19,19 +19,20 @@ solvers
 {
     "alpha.water.*"
     {
-        nAlphaCorr      1;
+        nAlphaCorr      2;
         nAlphaSubCycles 1;
         cAlpha          1;
-        icAlpha         0.25;
+        icAlpha         0;
 
         MULESCorr       yes;
-        nLimiterIter    8;
-        alphaApplyPrevCorr  true;
+        nLimiterIter    10;
+        alphaApplyPrevCorr  yes;
 
         solver          smoothSolver;
         smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
+        minIter         1;
     }
 
     pcorr
@@ -42,7 +43,7 @@ solvers
         {
             preconditioner  GAMG;
 
-            smoother        DICGaussSeidel;
+            smoother        GaussSeidel;
             agglomerator    faceAreaPair;
             mergeLevels     1;
             nCellsInCoarsestLevel 10;
@@ -60,7 +61,7 @@ solvers
     {
         solver          GAMG;
 
-        smoother        GaussSeidel;
+        smoother        DIC;
         agglomerator    faceAreaPair;
         mergeLevels     1;
         nCellsInCoarsestLevel 10;
@@ -73,7 +74,6 @@ solvers
     p_rghFinal
     {
         $p_rgh;
-        tolerance       1e-7;
         relTol          0;
     }
 
@@ -81,29 +81,31 @@ solvers
     {
         solver          smoothSolver;
 
-        smoother        GaussSeidel;
+        smoother        symGaussSeidel;
         nSweeps         1;
 
         tolerance       1e-7;
         relTol          0.1;
+        minIter         1;
     };
 }
 
 PIMPLE
 {
-    momentumPredictor no;
+    momentumPredictor   no;
 
-    nCorrectors     1;
+    nOuterCorrectors    1;
+    nCorrectors         2;
     nNonOrthogonalCorrectors 0;
 
-    maxCo           5;
-    maxAlphaCo      5;
+    maxCo               10;
+    maxAlphaCo          5;
 
     rDeltaTSmoothingCoeff 0.05;
     rDeltaTDampingCoeff 0.5;
-    nAlphaSweepIter  0;
-    nAlphaSpreadIter 0;
-    maxDeltaT       1;
+    nAlphaSpreadIter    0;
+    nAlphaSweepIter     0;
+    maxDeltaT           1;
 }
 
 relaxationFactors
@@ -117,5 +119,10 @@ relaxationFactors
     }
 }
 
+cache
+{
+    grad(U);
+}
+
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/meshQualityDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/meshQualityDict
new file mode 100644
index 0000000000..09d5f3928e
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/meshQualityDict
@@ -0,0 +1,86 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      meshQualityDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+//- Maximum non-orthogonality allowed. Set to 180 to disable.
+maxNonOrtho 70;
+
+//- Max skewness allowed. Set to <0 to disable.
+maxBoundarySkewness 20;
+maxInternalSkewness 4;
+
+//- Max concaveness allowed. Is angle (in degrees) below which concavity
+//  is allowed. 0 is straight face, <0 would be convex face.
+//  Set to 180 to disable.
+maxConcave 80;
+
+//- Minimum pyramid volume. Is absolute volume of cell pyramid.
+//  Set to a sensible fraction of the smallest cell volume expected.
+//  Set to very negative number (e.g. -1E30) to disable.
+minVol 1e-13;
+
+//- Minimum quality of the tet formed by the face-centre
+//  and variable base point minimum decomposition triangles and
+//  the cell centre.  Set to very negative number (e.g. -1E30) to
+//  disable.
+//     <0 = inside out tet,
+//      0 = flat tet
+//      1 = regular tet
+minTetQuality 1e-30;
+
+//- Minimum face area. Set to <0 to disable.
+minArea -1;
+
+//- Minimum face twist. Set to <-1 to disable. dot product of face normal
+//- and face centre triangles normal
+minTwist 0.05;
+
+//- minimum normalised cell determinant
+//- 1 = hex, <= 0 = folded or flattened illegal cell
+minDeterminant 0.001;
+
+//- minFaceWeight (0 -> 0.5)
+minFaceWeight 0.05;
+
+//- minVolRatio (0 -> 1)
+minVolRatio 0.01;
+
+//must be >0 for Fluent compatibility
+minTriangleTwist -1;
+
+//- if >0 : preserve single cells with all points on the surface if the
+//  resulting volume after snapping (by approximation) is larger than
+//  minVolCollapseRatio times old volume (i.e. not collapsed to flat cell).
+//  If <0 : delete always.
+//minVolCollapseRatio 0.5;
+
+
+// Advanced
+
+//- Number of error distribution iterations
+nSmoothScale 4;
+//- amount to scale back displacement at error points
+errorReduction 0.75;
+
+
+
+// Optional : some meshing phases allow usage of relaxed rules.
+// See e.g. addLayersControls::nRelaxedIter.
+relaxed
+{
+    //- Maximum non-orthogonality allowed. Set to 180 to disable.
+    maxNonOrtho 75;
+}
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/refineMeshDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/refineMeshDict
new file mode 100644
index 0000000000..ce5b7fe9cc
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/refineMeshDict
@@ -0,0 +1,43 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      refineMeshDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+set             c0;
+
+coordinateSystem global;
+
+globalCoeffs
+{
+    tan1            ( 1 0 0 );
+    tan2            ( 0 1 0 );
+}
+
+patchLocalCoeffs
+{
+    patch           outside;
+    tan1            ( 1 0 0 );
+}
+
+directions      ( tan1 tan2 );
+
+useHexTopology  no;
+
+geometricCut    yes;
+
+writeMesh       no;
+
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/setFieldsDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/setFieldsDict
similarity index 88%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/system/setFieldsDict
rename to tutorials/multiphase/LTSInterFoam/DTCHull/system/setFieldsDict
index 18f4b921a8..c19c3b991d 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/setFieldsDict
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/setFieldsDict
@@ -22,9 +22,12 @@ defaultFieldValues
 
 regions
 (
+    // Set cell values
+    // (does zerogradient on boundaries)
     boxToCell
     {
-        box (-100 -100 -100) (100 100 0);
+        box (-999 -999 -999) (999 999 0.244);
+
         fieldValues
         (
             volScalarFieldValue alpha.water 1
@@ -34,7 +37,7 @@ regions
     // Set patch values (using ==)
     boxToFace
     {
-        box (-100 -100 -100) (100 100 0);
+        box (-999 -999 -999) (999 999 0.244);
 
         fieldValues
         (
@@ -43,5 +46,4 @@ regions
     }
 );
 
-
 // ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/snappyHexMeshDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/snappyHexMeshDict
similarity index 72%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/system/snappyHexMeshDict
rename to tutorials/multiphase/LTSInterFoam/DTCHull/system/snappyHexMeshDict
index 36bd09b946..fbb5e12522 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/snappyHexMeshDict
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/snappyHexMeshDict
@@ -18,7 +18,7 @@ FoamFile
 // Which of the steps to run
 castellatedMesh true;
 snap            true;
-addLayers       false;
+addLayers       true;
 
 
 // Geometry. Definition of all surfaces. All surfaces are of class
@@ -29,7 +29,7 @@ addLayers       false;
 // - to 'snap' the mesh boundary to the surface
 geometry
 {
-    wigley-scaled-oriented.stl
+    DTC-scaled.stl
     {
         type triSurfaceMesh;
         name hull;
@@ -39,27 +39,6 @@ geometry
             type wall;
         }
     }
-
-    box1
-    {
-        type searchableBox;
-        min (-1 -0.5 -0.25);
-        max ( 1  0.5  1 );
-    }
-
-    box2
-    {
-        type searchableBox;
-        min (-1    -0.2  -0.15);
-        max ( 0.75  0.2   1);
-    }
-
-    box3
-    {
-        type searchableBox;
-        min (-10    -10  -0.03);
-        max ( 10  10   1);
-    }
 };
 
 
@@ -103,10 +82,10 @@ castellatedMeshControls
     // This is a featureEdgeMesh, read from constant/triSurface for now.
     features
     (
-        // {
-        //     file "wigley-scaled.eMesh";
-        //     level 3;
-        // }
+         {
+             file "DTC-scaled.eMesh";
+             level 0;
+         }
     );
 
 
@@ -125,7 +104,7 @@ castellatedMeshControls
         hull
         {
             // Surface-wise min and max refinement level
-            level (3 3);
+            level (0 0);
         }
     }
 
@@ -147,21 +126,6 @@ castellatedMeshControls
 
     refinementRegions
     {
-        box1
-        {
-            mode inside;
-            levels ((1e15 1));
-        }
-        box2
-        {
-            mode inside;
-            levels ((1e15 2));
-        }
-        box3
-        {
-            mode inside;
-            levels ((1e15 2));
-        }
     }
 
 
@@ -194,15 +158,17 @@ snapControls
     //- Relative distance for points to be attracted by surface feature point
     //  or edge. True distance is this factor times local
     //  maximum edge length.
-//    tolerance 4.0;
-tolerance 1.0;
+    //    tolerance 4.0;
+    tolerance 1.0;
 
     //- Number of mesh displacement relaxation iterations.
-    nSolveIter 300;
+    nSolveIter 100;
 
     //- Maximum number of snapping relaxation iterations. Should stop
     //  before upon reaching a correct mesh.
     nRelaxIter 5;
+
+    nFeatureSnapIter 10;
 }
 
 
@@ -219,19 +185,19 @@ addLayersControls
     {
         hull
         {
-            nSurfaceLayers 1;
+            nSurfaceLayers 3;
 
         }
     }
 
     // Expansion factor for layer mesh
-    expansionRatio 1.0;
+    expansionRatio 1.5;
 
 
     // Wanted thickness of final added cell layer. If multiple layers
     // is the thickness of the layer furthest away from the wall.
     // See relativeSizes parameter.
-    finalLayerThickness 0.3;
+    finalLayerThickness 0.7;
 
     // Minimum thickness of cell layer. If for any reason layer
     // cannot be above minThickness do not add layer.
@@ -296,80 +262,20 @@ addLayersControls
 // where to undo.
 meshQualityControls
 {
-    //- Maximum non-orthogonality allowed. Set to 180 to disable.
-    maxNonOrtho 65;
-
-    //- Max skewness allowed. Set to <0 to disable.
-    maxBoundarySkewness 20;
-    maxInternalSkewness 4;
-
-    //- Max concaveness allowed. Is angle (in degrees) below which concavity
-    //  is allowed. 0 is straight face, <0 would be convex face.
-    //  Set to 180 to disable.
-    maxConcave 80;
-
-    //- Minimum pyramid volume. Is absolute volume of cell pyramid.
-    //  Set to a sensible fraction of the smallest cell volume expected.
-    //  Set to very negative number (e.g. -1E30) to disable.
-    minVol 1e-13;
-
-    //- Minimum quality of the tet formed by the face-centre
-    //  and variable base point minimum decomposition triangles and
-    //  the cell centre.  Set to very negative number (e.g. -1E30) to
-    //  disable.
-    //     <0 = inside out tet,
-    //      0 = flat tet
-    //      1 = regular tet
-    minTetQuality 1e-30;
-
-    //- Minimum face area. Set to <0 to disable.
-    minArea -1;
-
-    //- Minimum face twist. Set to <-1 to disable. dot product of face normal
-    //- and face centre triangles normal
-    minTwist 0.05;
-
-    //- minimum normalised cell determinant
-    //- 1 = hex, <= 0 = folded or flattened illegal cell
-    minDeterminant 0.001;
-
-    //- minFaceWeight (0 -> 0.5)
-    minFaceWeight 0.05;
-
-    //- minVolRatio (0 -> 1)
-    minVolRatio 0.01;
-
-    //must be >0 for Fluent compatibility
-    minTriangleTwist -1;
-
-    //- if >0 : preserve single cells with all points on the surface if the
-    //  resulting volume after snapping (by approximation) is larger than
-    //  minVolCollapseRatio times old volume (i.e. not collapsed to flat cell).
-    //  If <0 : delete always.
-    //minVolCollapseRatio 0.5;
-
-
-    // Advanced
-
-    //- Number of error distribution iterations
-    nSmoothScale 4;
-    //- amount to scale back displacement at error points
-    errorReduction 0.75;
-
-
-
-    // Optional : some meshing phases allow usage of relaxed rules.
-    // See e.g. addLayersControls::nRelaxedIter.
-    relaxed
-    {
-        //- Maximum non-orthogonality allowed. Set to 180 to disable.
-        maxNonOrtho 75;
-    }
+    #include "meshQualityDict"
 }
 
 
 // Advanced
 
+// Flags for optional output
+// 0 : only write final meshes
+// 1 : write intermediate meshes
+// 2 : write volScalarField with cellLevel for postprocessing
+// 4 : write current intersections as .obj files
+debug 0;
+
+
 // Merge tolerance. Is fraction of overall bounding box of initial mesh.
 // Note: the write tolerance needs to be higher than this.
 mergeTolerance 1E-6;
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/surfaceFeatureExtractDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/surfaceFeatureExtractDict
new file mode 100644
index 0000000000..64ef9b9205
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/surfaceFeatureExtractDict
@@ -0,0 +1,48 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    object      surfaceFeatureExtractDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+DTC-scaled.stl
+{
+    // How to obtain raw features (extractFromFile || extractFromSurface)
+    extractionMethod    extractFromSurface;
+
+    extractFromSurfaceCoeffs
+    {
+        // Mark edges whose adjacent surface normals are at an angle less
+        // than includedAngle as features
+        // - 0  : selects no edges
+        // - 180: selects all edges
+        includedAngle   150;
+    }
+
+    subsetFeatures
+    {
+        // Keep nonManifold edges (edges with >2 connected faces)
+        nonManifoldEdges       yes;
+
+        // Keep open edges (edges with 1 connected face)
+        openEdges       yes;
+    }
+
+
+    // Write options
+
+        // Write features to obj format for postprocessing
+        writeObj                yes;
+}
+
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/decomposeParDict b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.1
similarity index 77%
rename from tutorials/multiphase/LTSInterFoam/wigleyHull/system/decomposeParDict
rename to tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.1
index 148aecf784..92ecf2d5bf 100644
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/system/decomposeParDict
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.1
@@ -11,23 +11,22 @@ FoamFile
     format      ascii;
     class       dictionary;
     location    "system";
-    object      decomposeParDict;
+    object      topoSetDict;
 }
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-numberOfSubdomains 4;
-
-method             simple; // Seriously don't actually use this method
-
-simpleCoeffs
-{
-    n               ( 2 2 1 );
-    delta           0.001;
-}
-
-distributed     no;
-
-roots           ( );
-
+actions
+(
+    {
+        name    c0;
+        type    cellSet;
+        action  new;
+        source  boxToCell;
+        sourceInfo
+        {
+            box (-10 -6 -3) (10 0 3);
+        }
+    }
+);
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.2 b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.2
new file mode 100644
index 0000000000..06f0616542
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.2
@@ -0,0 +1,32 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      topoSetDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+actions
+(
+    {
+        name    c0;
+        type    cellSet;
+        action  new;
+        source  boxToCell;
+        sourceInfo
+        {
+            box (-5 -3 -2.5) (9 0 2);
+        }
+    }
+);
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.3 b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.3
new file mode 100644
index 0000000000..04eaa8b97f
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.3
@@ -0,0 +1,32 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      topoSetDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+actions
+(
+    {
+        name    c0;
+        type    cellSet;
+        action  new;
+        source  boxToCell;
+        sourceInfo
+        {
+            box (-3 -1.5 -1) (8 0 1.5);
+        }
+    }
+);
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.4 b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.4
new file mode 100644
index 0000000000..23a54c9438
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.4
@@ -0,0 +1,32 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      topoSetDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+actions
+(
+    {
+        name    c0;
+        type    cellSet;
+        action  new;
+        source  boxToCell;
+        sourceInfo
+        {
+            box (-2 -1 -0.6) (7 0 1);
+        }
+    }
+);
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.5 b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.5
new file mode 100644
index 0000000000..a2e19ca8c1
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.5
@@ -0,0 +1,32 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      topoSetDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+actions
+(
+    {
+        name    c0;
+        type    cellSet;
+        action  new;
+        source  boxToCell;
+        sourceInfo
+        {
+            box (-1 -0.6 -0.3) (6.5 0 0.8);
+        }
+    }
+);
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.6 b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.6
new file mode 100644
index 0000000000..6af2ed95fb
--- /dev/null
+++ b/tutorials/multiphase/LTSInterFoam/DTCHull/system/topoSetDict.6
@@ -0,0 +1,32 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  dev                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "system";
+    object      topoSetDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+actions
+(
+    {
+        name    c0;
+        type    cellSet;
+        action  new;
+        source  boxToCell;
+        sourceInfo
+        {
+            box (-0.5 -0.55 -0.15) (6.25 0 0.65);
+        }
+    }
+);
+
+// ************************************************************************* //
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/Allclean b/tutorials/multiphase/LTSInterFoam/wigleyHull/Allclean
deleted file mode 100755
index 99318580f3..0000000000
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/Allclean
+++ /dev/null
@@ -1,14 +0,0 @@
-#!/bin/sh
-cd ${0%/*} || exit 1    # run from this directory
-
-# Source tutorial clean functions
-. $WM_PROJECT_DIR/bin/tools/CleanFunctions
-
-# remove surface
-\rm -f constant/triSurface/wigley-scaled-oriented.stl.gz
-
-cleanCase
-rm 0/alpha.water.gz 0/cellLevel.gz 0/pointLevel.gz constant/polyMesh/boundary \
-    > /dev/null 2>&1
-
-# ----------------------------------------------------------------- end-of-file
diff --git a/tutorials/multiphase/LTSInterFoam/wigleyHull/Allrun b/tutorials/multiphase/LTSInterFoam/wigleyHull/Allrun
deleted file mode 100755
index f57a6542d2..0000000000
--- a/tutorials/multiphase/LTSInterFoam/wigleyHull/Allrun
+++ /dev/null
@@ -1,20 +0,0 @@
-#!/bin/sh
-cd ${0%/*} || exit 1    # run from this directory
-
-# Source tutorial run functions
-. $WM_PROJECT_DIR/bin/tools/RunFunctions
-
-# copy wigley surface from resources directory
-cp $FOAM_TUTORIALS/resources/geometry/wigley-scaled-oriented.stl.gz constant/triSurface/
-
-runApplication blockMesh
-runApplication snappyHexMesh -overwrite
-if [ ! -e log.transformPoints ]
-then
-    transformPoints -scale '(2 1 1)' > log.transformPoints 2>&1
-fi
-cp 0/alpha.water.org 0/alpha.water
-runApplication setFields
-runApplication `getApplication`
-
-# ----------------------------------------------------------------- end-of-file
diff --git a/tutorials/multiphase/MRFInterFoam/mixerVessel2D/system/fvSolution b/tutorials/multiphase/MRFInterFoam/mixerVessel2D/system/fvSolution
index 4d0d891558..26571a081c 100644
--- a/tutorials/multiphase/MRFInterFoam/mixerVessel2D/system/fvSolution
+++ b/tutorials/multiphase/MRFInterFoam/mixerVessel2D/system/fvSolution
@@ -50,8 +50,8 @@ solvers
 
     U
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }
diff --git a/tutorials/multiphase/MRFMultiphaseInterFoam/mixerVessel2D/system/fvSolution b/tutorials/multiphase/MRFMultiphaseInterFoam/mixerVessel2D/system/fvSolution
index c4f3932da6..a19b09ca2a 100644
--- a/tutorials/multiphase/MRFMultiphaseInterFoam/mixerVessel2D/system/fvSolution
+++ b/tutorials/multiphase/MRFMultiphaseInterFoam/mixerVessel2D/system/fvSolution
@@ -50,8 +50,8 @@ solvers
 
     "(U|alpha)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }
diff --git a/tutorials/multiphase/cavitatingFoam/les/throttle/system/fvSolution b/tutorials/multiphase/cavitatingFoam/les/throttle/system/fvSolution
index f66b2caecf..1785af6f4e 100644
--- a/tutorials/multiphase/cavitatingFoam/les/throttle/system/fvSolution
+++ b/tutorials/multiphase/cavitatingFoam/les/throttle/system/fvSolution
@@ -19,16 +19,16 @@ solvers
 {
     rho
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.1;
     }
 
     rhoFinal
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
@@ -57,16 +57,16 @@ solvers
 
     "(U|k|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.1;
     }
 
     "(U|k|omega)Final"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/cavitatingFoam/les/throttle3D/system/fvSolution b/tutorials/multiphase/cavitatingFoam/les/throttle3D/system/fvSolution
index c29330bd28..7133194be0 100644
--- a/tutorials/multiphase/cavitatingFoam/les/throttle3D/system/fvSolution
+++ b/tutorials/multiphase/cavitatingFoam/les/throttle3D/system/fvSolution
@@ -19,8 +19,8 @@ solvers
 {
     rho
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.1;
     }
@@ -55,8 +55,8 @@ solvers
 
     "(U|k|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.1;
     }
diff --git a/tutorials/multiphase/cavitatingFoam/ras/throttle/system/fvSolution b/tutorials/multiphase/cavitatingFoam/ras/throttle/system/fvSolution
index f66b2caecf..1785af6f4e 100644
--- a/tutorials/multiphase/cavitatingFoam/ras/throttle/system/fvSolution
+++ b/tutorials/multiphase/cavitatingFoam/ras/throttle/system/fvSolution
@@ -19,16 +19,16 @@ solvers
 {
     rho
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.1;
     }
 
     rhoFinal
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
@@ -57,16 +57,16 @@ solvers
 
     "(U|k|omega)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0.1;
     }
 
     "(U|k|omega)Final"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSchemes b/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSchemes
index 0e66203335..045d871d38 100644
--- a/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSchemes
+++ b/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSchemes
@@ -60,7 +60,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSolution b/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSolution
index 2bec428fa5..5230b25c41 100644
--- a/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSolution
+++ b/tutorials/multiphase/compressibleInterDyMFoam/ras/sloshingTank2D/system/fvSolution
@@ -102,8 +102,8 @@ solvers
 
     "(T|k|B|nuTilda).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSolution b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSolution
index 62218a50fb..592be355bb 100644
--- a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSolution
+++ b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge2D/system/fvSolution
@@ -100,8 +100,8 @@ solvers
 
     "(T|k|B|nuTilda).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSolution b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSolution
index 563ce52459..fef915a850 100644
--- a/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSolution
+++ b/tutorials/multiphase/compressibleInterFoam/laminar/depthCharge3D/system/fvSolution
@@ -100,8 +100,8 @@ solvers
 
     "(T|k|B|nuTilda).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSchemes
index 9fe758f3d1..2ddd6a894a 100644
--- a/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSolution b/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSolution
index 1a1b04728a..e590cb4a8a 100644
--- a/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSolution
+++ b/tutorials/multiphase/interDyMFoam/ras/damBreakWithObstacle/system/fvSolution
@@ -97,8 +97,8 @@ solvers
 
     "(k|B|nuTilda)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution b/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution
index e920dcf8d9..ee177ba593 100644
--- a/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution
+++ b/tutorials/multiphase/interDyMFoam/ras/floatingObject/system/fvSolution
@@ -29,8 +29,8 @@ solvers
         nLimiterIter    5;
         alphaApplyPrevCorr  true;
 
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D/system/fvSchemes
index 90379ecc54..e9bac3f710 100644
--- a/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D3DoF/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D3DoF/system/fvSchemes
index 90379ecc54..e9bac3f710 100644
--- a/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D3DoF/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/sloshingTank2D3DoF/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D/system/fvSchemes
index 90379ecc54..e9bac3f710 100644
--- a/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D3DoF/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D3DoF/system/fvSchemes
index 90379ecc54..e9bac3f710 100644
--- a/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D3DoF/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D3DoF/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D6DoF/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D6DoF/system/fvSchemes
index 90379ecc54..e9bac3f710 100644
--- a/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D6DoF/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/sloshingTank3D6DoF/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interDyMFoam/ras/testTubeMixer/system/fvSchemes b/tutorials/multiphase/interDyMFoam/ras/testTubeMixer/system/fvSchemes
index 90379ecc54..e9bac3f710 100644
--- a/tutorials/multiphase/interDyMFoam/ras/testTubeMixer/system/fvSchemes
+++ b/tutorials/multiphase/interDyMFoam/ras/testTubeMixer/system/fvSchemes
@@ -53,7 +53,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interFoam/laminar/capillaryRise/system/fvSolution b/tutorials/multiphase/interFoam/laminar/capillaryRise/system/fvSolution
index 36156e8614..40102d5cbe 100644
--- a/tutorials/multiphase/interFoam/laminar/capillaryRise/system/fvSolution
+++ b/tutorials/multiphase/interFoam/laminar/capillaryRise/system/fvSolution
@@ -49,8 +49,8 @@ solvers
 
     U
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interFoam/laminar/damBreak/system/fvSolution b/tutorials/multiphase/interFoam/laminar/damBreak/system/fvSolution
index 0836eb66e5..bf1202d6ea 100644
--- a/tutorials/multiphase/interFoam/laminar/damBreak/system/fvSolution
+++ b/tutorials/multiphase/interFoam/laminar/damBreak/system/fvSolution
@@ -27,8 +27,8 @@ solvers
         MULESCorr       yes;
         nLimiterIter    3;
 
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSchemes b/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSchemes
index ddc0a7b6e4..07bd316e0e 100644
--- a/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSchemes
+++ b/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSchemes
@@ -57,7 +57,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSolution b/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSolution
index cb839db99a..47a4540f6c 100644
--- a/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSolution
+++ b/tutorials/multiphase/interFoam/les/nozzleFlow2D/system/fvSolution
@@ -95,8 +95,8 @@ solvers
 
     "(U|k)Final"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-08;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interFoam/ras/damBreak/system/controlDict b/tutorials/multiphase/interFoam/ras/damBreak/system/controlDict
index 8d4bf8a10b..b507f78b95 100644
--- a/tutorials/multiphase/interFoam/ras/damBreak/system/controlDict
+++ b/tutorials/multiphase/interFoam/ras/damBreak/system/controlDict
@@ -47,8 +47,8 @@ runTimeModifiable yes;
 
 adjustTimeStep  on;
 
-maxCo           0.2;
-maxAlphaCo      0.2;
+maxCo           1;
+maxAlphaCo      1;
 
 maxDeltaT       1;
 
diff --git a/tutorials/multiphase/interFoam/ras/damBreak/system/fvSchemes b/tutorials/multiphase/interFoam/ras/damBreak/system/fvSchemes
index 86fad259ec..3bc5347598 100644
--- a/tutorials/multiphase/interFoam/ras/damBreak/system/fvSchemes
+++ b/tutorials/multiphase/interFoam/ras/damBreak/system/fvSchemes
@@ -27,14 +27,11 @@ gradSchemes
 
 divSchemes
 {
-    div(rhoPhi,U)   Gauss linear;
+    div(rhoPhi,U)   Gauss linearUpwind grad(U);
     div(phi,alpha)  Gauss vanLeer;
     div(phirb,alpha) Gauss linear;
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
-    div(phi,R)      Gauss upwind;
-    div(R)          Gauss linear;
-    div(phi,nuTilda) Gauss upwind;
     div((muEff*dev(T(grad(U))))) Gauss linear;
 }
 
@@ -58,7 +55,7 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
+    alpha.water;
 }
 
 
diff --git a/tutorials/multiphase/interFoam/ras/damBreak/system/fvSolution b/tutorials/multiphase/interFoam/ras/damBreak/system/fvSolution
index f6725a947b..8cc95c0df2 100644
--- a/tutorials/multiphase/interFoam/ras/damBreak/system/fvSolution
+++ b/tutorials/multiphase/interFoam/ras/damBreak/system/fvSolution
@@ -17,18 +17,27 @@ FoamFile
 
 solvers
 {
-    alpha.water
+    "alpha.water.*"
     {
-        nAlphaCorr      1;
-        nAlphaSubCycles 4;
-        cAlpha          2;
+        nAlphaCorr      2;
+        nAlphaSubCycles 1;
+        alphaOuterCorrectors yes;
+        cAlpha          1;
+
+        MULESCorr       yes;
+        nLimiterIter    3;
+
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
+        tolerance       1e-8;
+        relTol          0;
     }
 
     pcorr
     {
         solver          PCG;
         preconditioner  DIC;
-        tolerance       1e-10;
+        tolerance       1e-5;
         relTol          0;
     }
 
@@ -42,35 +51,34 @@ solvers
 
     p_rghFinal
     {
-        solver          PCG;
-        preconditioner  DIC;
-        tolerance       1e-07;
+        $p_rgh;
         relTol          0;
     }
 
-    "(U|k|epsilon)"
+    "(U|k|epsilon).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
-    }
-
-    "(U|k|epsilon)Final"
-    {
-        solver          PBiCG;
-        preconditioner  DILU;
-        tolerance       1e-08;
-        relTol          0;
+        minIter         1;
     }
 }
 
 PIMPLE
 {
-    momentumPredictor no;
-    nCorrectors     3;
+    momentumPredictor   no;
+    nOuterCorrectors    1;
+    nCorrectors         3;
     nNonOrthogonalCorrectors 0;
 }
 
+relaxationFactors
+{
+    equations
+    {
+        ".*"                1;
+    }
+}
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/0/p_rgh b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/0/p_rgh
index 0b3ae667f0..f24e6e1ede 100644
--- a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/0/p_rgh
+++ b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/0/p_rgh
@@ -58,18 +58,18 @@ boundaryField
         patchType       cyclic;
         jump            uniform 0;
         value           uniform 0;
-        D               700;
+        D               1000;
         I               500;
-        length          1.05;
+        length          0.15;
     }
     porous_half1
     {
         type            porousBafflePressure;
         patchType       cyclic;
         value           uniform 0;
-        D               700;
+        D               1000;
         I               500;
-        length          1.05;
+        length          0.15;
     }
 }
 
diff --git a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/controlDict b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/controlDict
index bf6c452b46..53bb12ff93 100644
--- a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/controlDict
+++ b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/controlDict
@@ -29,7 +29,7 @@ deltaT          0.001;
 
 writeControl    adjustableRunTime;
 
-writeInterval   0.1;
+writeInterval   0.05;
 
 purgeWrite      0;
 
diff --git a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/createBafflesDict b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/createBafflesDict
index 46e0c66bf8..789f61dbc4 100644
--- a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/createBafflesDict
+++ b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/createBafflesDict
@@ -49,9 +49,9 @@ baffles
                     {
                         type            porousBafflePressure;
                         patchType       cyclic;
-                        D               700;
+                        D               1000;
                         I               500;
-                        length          1.05;
+                        length          0.15;
                         jump            uniform 0;
                         value           uniform 0;
                     }
diff --git a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSchemes b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSchemes
index 86fad259ec..3bc5347598 100644
--- a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSchemes
+++ b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSchemes
@@ -27,14 +27,11 @@ gradSchemes
 
 divSchemes
 {
-    div(rhoPhi,U)   Gauss linear;
+    div(rhoPhi,U)   Gauss linearUpwind grad(U);
     div(phi,alpha)  Gauss vanLeer;
     div(phirb,alpha) Gauss linear;
     div(phi,k)      Gauss upwind;
     div(phi,epsilon) Gauss upwind;
-    div(phi,R)      Gauss upwind;
-    div(R)          Gauss linear;
-    div(phi,nuTilda) Gauss upwind;
     div((muEff*dev(T(grad(U))))) Gauss linear;
 }
 
@@ -58,7 +55,7 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
+    alpha.water;
 }
 
 
diff --git a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSolution b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSolution
index a7103f6d62..8cc95c0df2 100644
--- a/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSolution
+++ b/tutorials/multiphase/interFoam/ras/damBreakPorousBaffle/system/fvSolution
@@ -17,18 +17,27 @@ FoamFile
 
 solvers
 {
-    alpha.water
+    "alpha.water.*"
     {
-        nAlphaCorr      1;
-        nAlphaSubCycles 4;
-        cAlpha          2;
+        nAlphaCorr      2;
+        nAlphaSubCycles 1;
+        alphaOuterCorrectors yes;
+        cAlpha          1;
+
+        MULESCorr       yes;
+        nLimiterIter    3;
+
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
+        tolerance       1e-8;
+        relTol          0;
     }
 
     pcorr
     {
         solver          PCG;
         preconditioner  DIC;
-        tolerance       1e-10;
+        tolerance       1e-5;
         relTol          0;
     }
 
@@ -42,41 +51,34 @@ solvers
 
     p_rghFinal
     {
-        solver          PCG;
-        preconditioner  DIC;
-        tolerance       1e-08;
+        $p_rgh;
         relTol          0;
     }
 
-    "(U|k|epsilon)"
+    "(U|k|epsilon).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
-    }
-
-    "(U|k|epsilon)Final"
-    {
-        solver          PBiCG;
-        preconditioner  DILU;
-        tolerance       1e-07;
-        relTol          0;
+        minIter         1;
     }
 }
 
 PIMPLE
 {
-    momentumPredictor yes;
-    nCorrectors     3;
+    momentumPredictor   no;
+    nOuterCorrectors    1;
+    nCorrectors         3;
     nNonOrthogonalCorrectors 0;
 }
 
 relaxationFactors
 {
-    "U.*"           1;
-    "k.*"           1;
-    "epsilon.*"     1;
+    equations
+    {
+        ".*"                1;
+    }
 }
 
 // ************************************************************************* //
diff --git a/tutorials/multiphase/interFoam/ras/waterChannel/system/fvSolution b/tutorials/multiphase/interFoam/ras/waterChannel/system/fvSolution
index ba6bfb0754..622dafc86c 100644
--- a/tutorials/multiphase/interFoam/ras/waterChannel/system/fvSolution
+++ b/tutorials/multiphase/interFoam/ras/waterChannel/system/fvSolution
@@ -26,8 +26,8 @@ solvers
         MULESCorr       yes;
         nLimiterIter    3;
 
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSchemes b/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSchemes
index a23cf9d48b..40392f6f90 100644
--- a/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSchemes
+++ b/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSchemes
@@ -58,7 +58,6 @@ fluxRequired
     default         no;
     p_rgh;
     pcorr;
-    alpha;
 }
 
 
diff --git a/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSolution b/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSolution
index dc48e429db..046c25c42c 100644
--- a/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSolution
+++ b/tutorials/multiphase/interFoam/ras/weirOverflow/system/fvSolution
@@ -47,8 +47,8 @@ solvers
 
     "(U|k|epsilon)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0.1;
     }
diff --git a/tutorials/multiphase/interMixingFoam/laminar/damBreak/system/fvSolution b/tutorials/multiphase/interMixingFoam/laminar/damBreak/system/fvSolution
index 3403701035..a1518045c8 100644
--- a/tutorials/multiphase/interMixingFoam/laminar/damBreak/system/fvSolution
+++ b/tutorials/multiphase/interMixingFoam/laminar/damBreak/system/fvSolution
@@ -56,8 +56,8 @@ solvers
 
     U
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }
diff --git a/tutorials/multiphase/interPhaseChangeFoam/cavitatingBullet/system/fvSolution b/tutorials/multiphase/interPhaseChangeFoam/cavitatingBullet/system/fvSolution
index 6e09921b2e..40709cfb30 100644
--- a/tutorials/multiphase/interPhaseChangeFoam/cavitatingBullet/system/fvSolution
+++ b/tutorials/multiphase/interPhaseChangeFoam/cavitatingBullet/system/fvSolution
@@ -25,16 +25,16 @@ solvers
         MULESCorr       yes;
         nLimiterIter    5;
 
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     };
 
     "U.*"
     {
-        solver           PBiCG;
-        preconditioner   DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance        1e-6;
         relTol           0;
     };
diff --git a/tutorials/multiphase/multiphaseEulerFoam/bubbleColumn/system/fvSolution b/tutorials/multiphase/multiphaseEulerFoam/bubbleColumn/system/fvSolution
index 73cb43356b..74edd57bbf 100644
--- a/tutorials/multiphase/multiphaseEulerFoam/bubbleColumn/system/fvSolution
+++ b/tutorials/multiphase/multiphaseEulerFoam/bubbleColumn/system/fvSolution
@@ -53,40 +53,40 @@ solvers
 
     "U.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "T.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     }
 
     "Theta.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "k.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "epsilon.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
diff --git a/tutorials/multiphase/settlingFoam/ras/dahl/system/fvSolution b/tutorials/multiphase/settlingFoam/ras/dahl/system/fvSolution
index 892618eaee..9194cde3c6 100644
--- a/tutorials/multiphase/settlingFoam/ras/dahl/system/fvSolution
+++ b/tutorials/multiphase/settlingFoam/ras/dahl/system/fvSolution
@@ -34,8 +34,8 @@ solvers
 
     "(U|Alpha|k|epsilon)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0.1;
     }
diff --git a/tutorials/multiphase/settlingFoam/ras/tank3D/system/fvSolution b/tutorials/multiphase/settlingFoam/ras/tank3D/system/fvSolution
index a399939cc3..8ad569733c 100644
--- a/tutorials/multiphase/settlingFoam/ras/tank3D/system/fvSolution
+++ b/tutorials/multiphase/settlingFoam/ras/tank3D/system/fvSolution
@@ -34,8 +34,8 @@ solvers
 
     "(U|Alpha|k|epsilon)"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-07;
         relTol          0.1;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/constant/phaseProperties
index 48d8185c6a..ce92781035 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/constant/phaseProperties
@@ -75,7 +75,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -87,7 +86,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -99,8 +97,6 @@ drag
     {
         type            segregated;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
-        residualRe      1e-3;
         m               0.5;
         n               8;
         swarmCorrection
@@ -130,13 +126,11 @@ heatTransfer
     (air in water)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 
     (water in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/system/fvSolution
index 31c2347993..64be5668a5 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/LES/bubbleColumn/system/fvSolution
@@ -47,24 +47,24 @@ solvers
 
     "U.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-7;
         relTol          0;
     }
 
     "h.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-7;
         relTol          0;
     }
 
     "(k|epsilon|Theta).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-7;
         relTol          0;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/constant/phaseProperties
index 48d8185c6a..ce92781035 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/constant/phaseProperties
@@ -75,7 +75,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -87,7 +86,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -99,8 +97,6 @@ drag
     {
         type            segregated;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
-        residualRe      1e-3;
         m               0.5;
         n               8;
         swarmCorrection
@@ -130,13 +126,11 @@ heatTransfer
     (air in water)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 
     (water in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/system/fvSolution
index 31c2347993..64be5668a5 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/RAS/bubbleColumn/system/fvSolution
@@ -47,24 +47,24 @@ solvers
 
     "U.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-7;
         relTol          0;
     }
 
     "h.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-7;
         relTol          0;
     }
 
     "(k|epsilon|Theta).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-7;
         relTol          0;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/constant/phaseProperties
index d79bf8cefa..07252b178c 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/constant/phaseProperties
@@ -40,7 +40,7 @@ blending
     default
     {
         type            none;
-        residualAlpha   1e-6;
+        residualAlpha   1e-3;
         continuousPhase air;
     }
 }
@@ -60,7 +60,6 @@ drag
     {
         type            GidaspowErgunWenYu;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -74,7 +73,7 @@ virtualMass
     (particles in air)
     {
         type            constantCoefficient;
-        Cvm             0.5;
+        Cvm             0;
     }
 );
 
@@ -83,7 +82,6 @@ heatTransfer
     (particles in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/system/fvSolution
index f5519a78a9..7e637fcfb4 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/RAS/fluidisedBed/system/fvSolution
@@ -53,32 +53,32 @@ solvers
 
     "U.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "(h|e).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-6;
         relTol          0;
     }
 
     "Theta.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-6;
         relTol          0;
     }
 
     "(k|epsilon).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties
index 48d8185c6a..ce92781035 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/constant/phaseProperties
@@ -75,7 +75,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -87,7 +86,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -99,8 +97,6 @@ drag
     {
         type            segregated;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
-        residualRe      1e-3;
         m               0.5;
         n               8;
         swarmCorrection
@@ -130,13 +126,11 @@ heatTransfer
     (air in water)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 
     (water in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/system/fvSolution
index a25ccff8c6..128aa1fa38 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumn/system/fvSolution
@@ -47,16 +47,16 @@ solvers
 
     "U.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "h.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/constant/phaseProperties
index 48d8185c6a..ce92781035 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/constant/phaseProperties
@@ -75,7 +75,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -87,7 +86,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -99,8 +97,6 @@ drag
     {
         type            segregated;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
-        residualRe      1e-3;
         m               0.5;
         n               8;
         swarmCorrection
@@ -130,13 +126,11 @@ heatTransfer
     (air in water)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 
     (water in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/system/fvSolution
index 0a266db563..61905c606a 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/bubbleColumnIATE/system/fvSolution
@@ -47,16 +47,16 @@ solvers
 
     "(U|kappai).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "h.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-8;
         relTol          0;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/constant/phaseProperties
index d79bf8cefa..a1aec908d9 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/constant/phaseProperties
@@ -40,7 +40,7 @@ blending
     default
     {
         type            none;
-        residualAlpha   1e-6;
+        residualAlpha   1e-3;
         continuousPhase air;
     }
 }
@@ -60,7 +60,6 @@ drag
     {
         type            GidaspowErgunWenYu;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -83,7 +82,6 @@ heatTransfer
     (particles in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/system/fvSolution
index f5519a78a9..7e637fcfb4 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/fluidisedBed/system/fvSolution
@@ -53,32 +53,32 @@ solvers
 
     "U.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
 
     "(h|e).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-6;
         relTol          0;
     }
 
     "Theta.*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-6;
         relTol          0;
     }
 
     "(k|epsilon).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-05;
         relTol          0;
     }
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/constant/phaseProperties b/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/constant/phaseProperties
index 48d8185c6a..ce92781035 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/constant/phaseProperties
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/constant/phaseProperties
@@ -75,7 +75,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -87,7 +86,6 @@ drag
     {
         type            SchillerNaumann;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
         residualRe      1e-3;
         swarmCorrection
         {
@@ -99,8 +97,6 @@ drag
     {
         type            segregated;
         residualAlpha   1e-6;
-        residualSlip    1e-3;
-        residualRe      1e-3;
         m               0.5;
         n               8;
         swarmCorrection
@@ -130,13 +126,11 @@ heatTransfer
     (air in water)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 
     (water in air)
     {
         type            RanzMarshall;
-        residualRe      1e-3;
     }
 );
 
diff --git a/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/system/fvSolution b/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/system/fvSolution
index 1646b9187f..380ea509a2 100644
--- a/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/system/fvSolution
+++ b/tutorials/multiphase/twoPhaseEulerFoam/laminar/mixerVessel2D/system/fvSolution
@@ -53,8 +53,8 @@ solvers
 
     "(k|epsilon|Theta|h).*"
     {
-        solver          PBiCG;
-        preconditioner  DILU;
+        solver          smoothSolver;
+        smoother        symGaussSeidel;
         tolerance       1e-06;
         relTol          0;
     }