diff --git a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/UEqns.H b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/UEqns.H
index ed150abb88..fe40ce2f49 100644
--- a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/UEqns.H
+++ b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/UEqns.H
@@ -18,7 +18,7 @@ fvVectorMatrix U2Eqn(U2, U2.dimensions()*dimVol/dimTime);
         volTensorField Rc1
         (
             "Rc1",
-            ((2.0/3.0)*I)*(sqr(Ct)*k + nuEff1*tr(gradU1T)) - nuEff1*gradU1T
+            (((2.0/3.0)*I)*nuEff1)*tr(gradU1T) - nuEff1*gradU1T
         );
 
         if (kineticTheory.on())
@@ -28,19 +28,24 @@ fvVectorMatrix U2Eqn(U2, U2.dimensions()*dimVol/dimTime);
 
         U1Eqn =
         (
-            (scalar(1) + Cvm*rho2*alpha2/rho1)*
+            fvm::ddt(alpha1, U1)
+          + fvm::div(alphaPhi1, U1)
+          - fvm::Sp(fvc::ddt(alpha1) + fvc::div(alphaPhi1), U1)
+
+          + Cvm*rho2*alpha1*alpha2/rho1*
             (
-                fvm::ddt(alpha1, U1)
-              + fvm::div(alphaPhi1, U1)
-              - fvm::Sp(fvc::ddt(alpha1) + fvc::div(alphaPhi1), U1)
+                fvm::ddt(U1)
+              + fvm::div(phi1, U1)
+              - fvm::Sp(fvc::div(phi1), U1)
             )
+
           - fvm::laplacian(alpha1*nuEff1, U1)
           + fvc::div(alpha1*Rc1)
          ==
           - fvm::Sp(dragCoeff/rho1, U1)
           - alpha1*alpha2/rho1*(liftForce - Cvm*rho2*DDtU2)
         );
-
+        mrfZones.addCoriolis(alpha1, U1Eqn);
         U1Eqn.relax();
     }
 
@@ -49,24 +54,29 @@ fvVectorMatrix U2Eqn(U2, U2.dimensions()*dimVol/dimTime);
         volTensorField Rc2
         (
             "Rc2",
-            ((2.0/3.0)*I)*(k + nuEff2*tr(gradU2T)) - nuEff2*gradU2T
+            (((2.0/3.0)*I)*nuEff2)*tr(gradU2T) - nuEff2*gradU2T
         );
 
         U2Eqn =
         (
-            (scalar(1) + Cvm*rho2*alpha1/rho2)*
+            fvm::ddt(alpha2, U2)
+          + fvm::div(alphaPhi2, U2)
+          - fvm::Sp(fvc::ddt(alpha2) + fvc::div(alphaPhi2), U2)
+
+          + Cvm*rho2*alpha1*alpha2/rho2*
             (
-                fvm::ddt(alpha2, U2)
-              + fvm::div(alphaPhi2, U2)
-              - fvm::Sp(fvc::ddt(alpha2) + fvc::div(alphaPhi2), U2)
+                fvm::ddt(U2)
+              + fvm::div(phi2, U2)
+              - fvm::Sp(fvc::div(phi2), U2)
             )
+
           - fvm::laplacian(alpha2*nuEff2, U2)
           + fvc::div(alpha2*Rc2)
          ==
           - fvm::Sp(dragCoeff/rho2, U2)
           + alpha1*alpha2/rho2*(liftForce + Cvm*rho2*DDtU1)
         );
-
+        mrfZones.addCoriolis(alpha2, U2Eqn);
         U2Eqn.relax();
     }
 }
diff --git a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/compressibleTwoPhaseEulerFoam.C b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/compressibleTwoPhaseEulerFoam.C
index d29da482ba..78d70a1537 100644
--- a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/compressibleTwoPhaseEulerFoam.C
+++ b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/compressibleTwoPhaseEulerFoam.C
@@ -45,6 +45,8 @@ Description
 
 #include "pimpleControl.H"
 
+#include "MRFZones.H"
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 int main(int argc, char *argv[])
@@ -55,6 +57,7 @@ int main(int argc, char *argv[])
     #include "createMesh.H"
     #include "readGravitationalAcceleration.H"
     #include "createFields.H"
+    #include "createMRFZones.H"
     #include "readPPProperties.H"
     #include "initContinuityErrs.H"
     #include "readTimeControls.H"
diff --git a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/createMRFZones.H b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/createMRFZones.H
new file mode 100644
index 0000000000..4d5c2bab72
--- /dev/null
+++ b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/createMRFZones.H
@@ -0,0 +1,4 @@
+    MRFZones mrfZones(mesh);
+    mrfZones.correctBoundaryVelocity(U1);
+    mrfZones.correctBoundaryVelocity(U2);
+    mrfZones.correctBoundaryVelocity(U);
diff --git a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/pEqn.H b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/pEqn.H
index d2ef15776d..aabbc6d4c2 100644
--- a/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/pEqn.H
+++ b/applications/solvers/multiphase/compressibleTwoPhaseEulerFoam/pEqn.H
@@ -2,6 +2,11 @@
     rho1 = rho10 + psi1*p;
     rho2 = rho20 + psi2*p;
 
+    mrfZones.absoluteFlux(phi1.oldTime());
+    mrfZones.absoluteFlux(phi1);
+    mrfZones.absoluteFlux(phi2.oldTime());
+    mrfZones.absoluteFlux(phi2);
+
     surfaceScalarField alpha1f(fvc::interpolate(alpha1));
     surfaceScalarField alpha2f(scalar(1) - alpha1f);
 
@@ -25,6 +30,7 @@
       + fvc::interpolate((1.0/rho1)*rAU1*dragCoeff)*phi2
       + rAlphaAU1f*(g & mesh.Sf())
     );
+    mrfZones.relativeFlux(phiHbyA1);
 
     surfaceScalarField phiHbyA2
     (
@@ -34,9 +40,13 @@
       + fvc::interpolate((1.0/rho2)*rAU2*dragCoeff)*phi1
       + rAlphaAU2f*(g & mesh.Sf())
     );
+    mrfZones.relativeFlux(phiHbyA2);
 
     surfaceScalarField phiHbyA("phiHbyA", alpha1f*phiHbyA1 + alpha2f*phiHbyA2);
 
+    HbyA1 += (1.0/rho1)*rAU1*dragCoeff*U2;
+    HbyA2 += (1.0/rho2)*rAU2*dragCoeff*U1;
+
     surfaceScalarField Dp
     (
         "Dp",
@@ -94,10 +104,12 @@
             phi1.boundaryField() ==
                 (fvc::interpolate(U1) & mesh.Sf())().boundaryField();
             phi1 = phiHbyA1 + rAlphaAU1f*mSfGradp/fvc::interpolate(rho1);
+            mrfZones.relativeFlux(phi1.oldTime());
 
             phi2.boundaryField() ==
                 (fvc::interpolate(U2) & mesh.Sf())().boundaryField();
             phi2 = phiHbyA2 + rAlphaAU2f*mSfGradp/fvc::interpolate(rho2);
+            mrfZones.relativeFlux(phi2.oldTime());
 
             phi = alpha1f*phi1 + alpha2f*phi2;
 
@@ -110,10 +122,7 @@
             p.relax();
             mSfGradp = pEqnIncomp.flux()/Dp;
 
-            volVectorField U10 = U1;
-
             U1 = HbyA1
-               + (1.0/rho1)*rAU1*dragCoeff*U2
                + fvc::reconstruct
                  (
                      rAlphaAU1f*(g & mesh.Sf())
@@ -122,7 +131,6 @@
             U1.correctBoundaryConditions();
 
             U2 = HbyA2
-               + (1.0/rho2)*rAU2*dragCoeff*U10
                + fvc::reconstruct
                  (
                      rAlphaAU2f*(g & mesh.Sf())
diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/UEqns.H b/applications/solvers/multiphase/multiphaseEulerFoam/UEqns.H
index bce95f7d5a..b612e6b5d9 100644
--- a/applications/solvers/multiphase/multiphaseEulerFoam/UEqns.H
+++ b/applications/solvers/multiphase/multiphaseEulerFoam/UEqns.H
@@ -15,12 +15,17 @@ forAllIter(PtrDictionary<phaseModel>, fluid.phases(), iter)
         phasei,
         new fvVectorMatrix
         (
-            (scalar(1) + fluid.Cvm(phase)/phase.rho())*
+            fvm::ddt(alpha, U)
+          + fvm::div(phase.phiAlpha(), U)
+          - fvm::Sp(fvc::ddt(alpha) + fvc::div(phase.phiAlpha()), U)
+
+          + (alpha/phase.rho())*fluid.Cvm(phase)*
             (
-                fvm::ddt(alpha, U)
-              + fvm::div(phase.phiAlpha(), U)
-              - fvm::Sp(fvc::ddt(alpha) + fvc::div(phase.phiAlpha()), U)
+                fvm::ddt(U)
+              + fvm::div(phase.phi(), U)
+              - fvm::Sp(fvc::div(phase.phi()), U)
             )
+
           - fvm::laplacian(alpha*nuEff, U)
           - fvc::div
             (
diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H b/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H
index 4a817a1001..1c42c7b91a 100644
--- a/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H
+++ b/applications/solvers/multiphase/multiphaseEulerFoam/pEqn.H
@@ -76,6 +76,7 @@
           + fvc::ddtPhiCorr(rAUs[phasei], alpha, phase.U(), phase.phi())
           + rAlphaAUfs[phasei]*(g & mesh.Sf())
         );
+
         mrfZones.relativeFlux(phiHbyAs[phasei]);
 
         multiphaseSystem::dragModelTable::const_iterator dmIter =
@@ -89,11 +90,20 @@
             ++dmIter, ++dcIter
         )
         {
-            const phaseModel *phase2Ptr = &dmIter()->phase2();
-            if (phase2Ptr == &phase)
+            const phaseModel *phase2Ptr = NULL;
+
+            if (&phase == &dmIter()->phase1())
+            {
+                phase2Ptr = &dmIter()->phase2();
+            }
+            else if (&phase == &dmIter()->phase2())
             {
                 phase2Ptr = &dmIter()->phase1();
             }
+            else
+            {
+                continue;
+            }
 
             phiHbyAs[phasei] +=
                 fvc::interpolate
@@ -105,6 +115,9 @@
                 (1.0/phase.rho())*rAUs[phasei]*(*dcIter())
                *phase2Ptr->U();
 
+            // Alternative flux-pressure consistent drag
+            // but not momentum conservative
+            //
             // HbyAs[phasei] += fvc::reconstruct
             // (
             //     fvc::interpolate
diff --git a/applications/solvers/multiphase/twoPhaseEulerFoam/MRFtwoPhaseEulerFoam/createMRFZones.H b/applications/solvers/multiphase/twoPhaseEulerFoam/MRFtwoPhaseEulerFoam/createMRFZones.H
index 2fc0eff9eb..4d5c2bab72 100644
--- a/applications/solvers/multiphase/twoPhaseEulerFoam/MRFtwoPhaseEulerFoam/createMRFZones.H
+++ b/applications/solvers/multiphase/twoPhaseEulerFoam/MRFtwoPhaseEulerFoam/createMRFZones.H
@@ -1,3 +1,4 @@
     MRFZones mrfZones(mesh);
     mrfZones.correctBoundaryVelocity(U1);
     mrfZones.correctBoundaryVelocity(U2);
+    mrfZones.correctBoundaryVelocity(U);
diff --git a/applications/utilities/mesh/manipulation/createPatch/createPatch.C b/applications/utilities/mesh/manipulation/createPatch/createPatch.C
index 50c9f30839..7a969e15d1 100644
--- a/applications/utilities/mesh/manipulation/createPatch/createPatch.C
+++ b/applications/utilities/mesh/manipulation/createPatch/createPatch.C
@@ -111,10 +111,18 @@ void filterPatches(polyMesh& mesh, const HashSet<word>& addedPatchNames)
         // Note: reduce possible since non-proc patches guaranteed in same order
         if (!isA<processorPolyPatch>(pp))
         {
+
+            // Add if
+            // - non zero size
+            // - or added from the createPatchDict
+            // - or cyclic (since referred to by other cyclic half or
+            //   proccyclic)
+
             if
             (
                 addedPatchNames.found(pp.name())
              || returnReduce(pp.size(), sumOp<label>()) > 0
+             || isA<coupledPolyPatch>(pp)
             )
             {
                 allPatches.append
@@ -131,6 +139,7 @@ void filterPatches(polyMesh& mesh, const HashSet<word>& addedPatchNames)
             else
             {
                 Info<< "Removing zero-sized patch " << pp.name()
+                    << " type " << pp.type()
                     << " at position " << patchI << endl;
             }
         }
diff --git a/applications/utilities/preProcessing/setFields/setFields.C b/applications/utilities/preProcessing/setFields/setFields.C
index c28df31a2a..7b690d46f2 100644
--- a/applications/utilities/preProcessing/setFields/setFields.C
+++ b/applications/utilities/preProcessing/setFields/setFields.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -216,14 +216,20 @@ bool setFaceFieldType
         }
 
         // Override
+        bool hasWarned = false;
         labelList nChanged(field.boundaryField().size(), 0);
         forAll(selectedFaces, i)
         {
             label facei = selectedFaces[i];
             if (mesh.isInternalFace(facei))
             {
-                WarningIn("setFaceFieldType(..)")
-                    << "Ignoring internal face " << facei << endl;
+                if (!hasWarned)
+                {
+                    hasWarned = true;
+                    WarningIn("setFaceFieldType(..)")
+                        << "Ignoring internal face " << facei
+                        << ". Suppressing further warnings." << endl;
+                }
             }
             else
             {
diff --git a/applications/utilities/preProcessing/setFields/setFieldsDict b/applications/utilities/preProcessing/setFields/setFieldsDict
index 6f4b9cd490..695ab38383 100644
--- a/applications/utilities/preProcessing/setFields/setFieldsDict
+++ b/applications/utilities/preProcessing/setFields/setFieldsDict
@@ -16,7 +16,7 @@ FoamFile
 
 defaultFieldValues
 (
-    volScalarFieldValue gamma 0
+    volScalarFieldValue alpha1 0
     volVectorFieldValue U (0 0 0)
 );
 
@@ -30,7 +30,7 @@ regions
 
         fieldValues
         (
-            volScalarFieldValue gamma 1
+            volScalarFieldValue alpha1 1
         );
     }
 
@@ -41,7 +41,7 @@ regions
 
         fieldValues
         (
-            volScalarFieldValue gamma 1
+            volScalarFieldValue alpha1 1
         );
     }
 );
diff --git a/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/axisRotationMotion/axisRotationMotion.C b/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/axisRotationMotion/axisRotationMotion.C
index 456441026f..2de3fab7e0 100644
--- a/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/axisRotationMotion/axisRotationMotion.C
+++ b/src/dynamicFvMesh/solidBodyMotionFvMesh/solidBodyMotionFunctions/axisRotationMotion/axisRotationMotion.C
@@ -73,8 +73,6 @@ Foam::solidBodyMotionFunctions::axisRotationMotion::transformation() const
 {
     scalar t = time_.value();
 
-    // Motion around a centre of gravity
-
     // Rotation around centre of gravity (in radians)
     vector omega
     (
@@ -84,10 +82,7 @@ Foam::solidBodyMotionFunctions::axisRotationMotion::transformation() const
     );
 
     scalar magOmega = mag(omega);
-
-    scalar cosHalfTheta = cos(0.5*magOmega);
-
-    quaternion R(cosHalfTheta, omega/magOmega);
+    quaternion R(omega/magOmega, magOmega);
     septernion TR(septernion(CofG_)*R*septernion(-CofG_));
 
     Info<< "solidBodyMotionFunctions::axisRotationMotion::transformation(): "
diff --git a/src/thermophysicalModels/basic/basicThermo/basicThermo.C b/src/thermophysicalModels/basic/basicThermo/basicThermo.C
index ad40123d6e..59c0056040 100644
--- a/src/thermophysicalModels/basic/basicThermo/basicThermo.C
+++ b/src/thermophysicalModels/basic/basicThermo/basicThermo.C
@@ -350,6 +350,38 @@ Foam::tmp<Foam::volScalarField> Foam::basicThermo::hc() const
 }
 
 
+Foam::tmp<Foam::scalarField> Foam::basicThermo::TH
+(
+    const scalarField& h,
+    const scalarField& T0,      // starting temperature
+    const labelList& cells
+) const
+{
+    notImplemented
+    (
+        "basicThermo::TH"
+        "(const scalarField&, const scalarField&, const labelList&) const"
+    );
+    return tmp<scalarField>(NULL);
+}
+
+
+Foam::tmp<Foam::scalarField> Foam::basicThermo::TH
+(
+    const scalarField& h,
+    const scalarField& T0,      // starting temperature
+    const label patchi
+) const
+{
+    notImplemented
+    (
+        "basicThermo::TH"
+        "(const scalarField&, const scalarField&, const label) const"
+    );
+    return tmp<scalarField>(NULL);
+}
+
+
 Foam::volScalarField& Foam::basicThermo::e()
 {
     notImplemented("basicThermo::e()");
diff --git a/src/thermophysicalModels/basic/basicThermo/basicThermo.H b/src/thermophysicalModels/basic/basicThermo/basicThermo.H
index 0bb67f8b49..92c4b01258 100644
--- a/src/thermophysicalModels/basic/basicThermo/basicThermo.H
+++ b/src/thermophysicalModels/basic/basicThermo/basicThermo.H
@@ -198,6 +198,22 @@ public:
             //- Chemical enthalpy [J/kg]
             virtual tmp<volScalarField> hc() const;
 
+            //- Temperature from enthalpy for cell-set
+            virtual tmp<scalarField> TH
+            (
+                const scalarField& h,
+                const scalarField& T0,      // starting temperature
+                const labelList& cells
+            ) const;
+
+            //- Temperature from enthalpy for patch
+            virtual tmp<scalarField> TH
+            (
+                const scalarField& h,
+                const scalarField& T0,      // starting temperature
+                const label patchi
+            ) const;
+
             //- Internal energy [J/kg]
             //  Non-const access allowed for transport equations
             virtual volScalarField& e();
diff --git a/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.C b/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.C
index 9d26c95015..9d3f3c2ba7 100644
--- a/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.C
+++ b/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -205,6 +205,50 @@ Foam::tmp<Foam::scalarField> Foam::hRhoThermo<MixtureType>::h
 }
 
 
+template<class MixtureType>
+Foam::tmp<Foam::scalarField> Foam::hRhoThermo<MixtureType>::TH
+(
+    const scalarField& h,
+    const scalarField& T0,
+    const labelList& cells
+) const
+{
+    tmp<scalarField> tT(new scalarField(h.size()));
+    scalarField& T = tT();
+
+    forAll(h, celli)
+    {
+        T[celli] = this->cellMixture(cells[celli]).TH(h[celli], T0[celli]);
+    }
+
+    return tT;
+}
+
+
+template<class MixtureType>
+Foam::tmp<Foam::scalarField> Foam::hRhoThermo<MixtureType>::TH
+(
+    const scalarField& h,
+    const scalarField& T0,
+    const label patchi
+) const
+{
+    tmp<scalarField> tT(new scalarField(h.size()));
+    scalarField& T = tT();
+
+    forAll(h, facei)
+    {
+        T[facei] = this->patchFaceMixture
+        (
+            patchi,
+            facei
+        ).TH(h[facei], T0[facei]);
+    }
+
+    return tT;
+}
+
+
 template<class MixtureType>
 Foam::tmp<Foam::scalarField> Foam::hRhoThermo<MixtureType>::Cp
 (
diff --git a/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.H b/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.H
index ee3cd4ba7f..1f85793f6b 100644
--- a/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.H
+++ b/src/thermophysicalModels/basic/rhoThermo/hRhoThermo/hRhoThermo.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -134,6 +134,22 @@ public:
                 const label patchi
             ) const;
 
+            //- Temperature from enthalpy for cell-set
+            virtual tmp<scalarField> TH
+            (
+                const scalarField& h,
+                const scalarField& T0,      // starting temperature
+                const labelList& cells
+            ) const;
+
+            //- Temperature from enthalpy for patch
+            virtual tmp<scalarField> TH
+            (
+                const scalarField& h,
+                const scalarField& T0,      // starting temperature
+                const label patchi
+            ) const;
+
             //- Heat capacity at constant pressure for patch [J/kg/K]
             virtual tmp<scalarField> Cp
             (
diff --git a/src/triSurface/triSurface/interfaces/VTK/readVTK.C b/src/triSurface/triSurface/interfaces/VTK/readVTK.C
index 71b3e791aa..a5ffc2e247 100644
--- a/src/triSurface/triSurface/interfaces/VTK/readVTK.C
+++ b/src/triSurface/triSurface/interfaces/VTK/readVTK.C
@@ -35,7 +35,7 @@ namespace Foam
 
 bool Foam::triSurface::readVTK(const fileName& fName)
 {
-    // Read (and triagulate) point, faces, zone info
+    // Read (and triangulate) point, faces, zone info
     fileFormats::VTKsurfaceFormat<triFace> surf(fName);
 
     List<labelledTri> tris(surf.faces().size());
diff --git a/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/0/alpha1 b/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/0/alpha1
index 4d25044844..9a514b7807 100644
--- a/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/0/alpha1
+++ b/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/0/alpha1
@@ -970,381 +970,381 @@ internalField   nonuniform List<scalar>
 0
 0
 0
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
-1
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
+0
 1
 1
 1
diff --git a/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/system/setFieldsDict b/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/system/setFieldsDict
index 84545091c6..5b6d65bbec 100644
--- a/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/system/setFieldsDict
+++ b/tutorials/multiphase/compressibleTwoPhaseEulerFoam/bubbleColumn/system/setFieldsDict
@@ -24,7 +24,7 @@ regions
 (
     boxToCell
     {
-        box (0 0 -0.1) (0.15 0.501 0.1);
+        box (0 0 -0.1) (0.15 0.701 0.1);
         fieldValues
         (
             volScalarFieldValue alpha1 0