diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options
index 373084db4b8e80e55f8514db4f2888eb05ca7959..e3fb4e28e5ed6c6e652daa709df7f6017522fa99 100644
--- a/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options
+++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/Make/options
@@ -1,6 +1,7 @@
 interPhaseChangePath = $(FOAM_SOLVERS)/multiphase/interPhaseChangeFoam
 
 EXE_INC = \
+    -I. \
     -ItemperaturePhaseChangeTwoPhaseMixtures/lnInclude \
     -I$(interPhaseChangePath) \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H b/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H
index c0e2adcb51a99b341fb7ccaa4e4903d6ca477d70..4d43cf3d5da3ad3a666b2759787c1b21aac6a687 100644
--- a/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H
+++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/createFields.H
@@ -1,143 +1,142 @@
-    Info<< "Reading field p_rgh\n" << endl;
-    volScalarField p_rgh
+Info<< "Reading field p_rgh\n" << endl;
+volScalarField p_rgh
+(
+    IOobject
     (
-        IOobject
-        (
-            "p_rgh",
-            runTime.timeName(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh
-    );
-
-    Info<< "Reading field U\n" << endl;
-    volVectorField U
+        "p_rgh",
+        runTime.timeName(),
+        mesh,
+        IOobject::MUST_READ,
+        IOobject::AUTO_WRITE
+    ),
+    mesh
+);
+
+Info<< "Reading field U\n" << endl;
+volVectorField U
+(
+    IOobject
     (
-        IOobject
-        (
-            "U",
-            runTime.timeName(),
-            mesh,
-            IOobject::MUST_READ,
-            IOobject::AUTO_WRITE
-        ),
-        mesh
-    );
-
-    #include "createPhi.H"
-
-    // Creating e based thermo
-    autoPtr<twoPhaseMixtureEThermo> thermo;
-    thermo.set(new twoPhaseMixtureEThermo(U, phi));
+        "U",
+        runTime.timeName(),
+        mesh,
+        IOobject::MUST_READ,
+        IOobject::AUTO_WRITE
+    ),
+    mesh
+);
 
-    // Create mixture and
-    Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n" << endl;
-    autoPtr<temperaturePhaseChangeTwoPhaseMixture> mixture =
-        temperaturePhaseChangeTwoPhaseMixture::New(thermo(), mesh);
+#include "createPhi.H"
 
+// Creating e based thermo
+autoPtr<twoPhaseMixtureEThermo> thermo;
+thermo.set(new twoPhaseMixtureEThermo(U, phi));
 
-    volScalarField& T = thermo->T();
-    volScalarField& e = thermo->he();
-    e.oldTime();
+// Create mixture and
+Info<< "Creating temperaturePhaseChangeTwoPhaseMixture\n" << endl;
+autoPtr<temperaturePhaseChangeTwoPhaseMixture> mixture =
+    temperaturePhaseChangeTwoPhaseMixture::New(thermo(), mesh);
 
-    // Correct e from T and alpha
-    //thermo->correct();
 
-    volScalarField& alpha1(thermo->alpha1());
-    volScalarField& alpha2(thermo->alpha2());
+// Correct e from T and alpha
+//thermo->correct();
 
-    const dimensionedScalar& rho1 = thermo->rho1();
-    const dimensionedScalar& rho2 = thermo->rho2();
+volScalarField& alpha1(thermo->alpha1());
+volScalarField& alpha2(thermo->alpha2());
 
-    // Need to store rho for ddt(rho, U)
-    volScalarField rho
-    (
-        IOobject
-        (
-            "rho",
-            runTime.timeName(),
-            mesh,
-            IOobject::READ_IF_PRESENT,
-            IOobject::AUTO_WRITE
-        ),
-        alpha1*rho1 + alpha2*rho2,
-        alpha1.boundaryField().types()
-    );
-    rho.oldTime();
+const dimensionedScalar& rho1 = thermo->rho1();
+const dimensionedScalar& rho2 = thermo->rho2();
 
-    // Construct interface from alpha1 distribution
-    interfaceProperties interface
+// Need to store rho for ddt(rho, U)
+volScalarField rho
+(
+    IOobject
     (
-        alpha1,
-        U,
-        thermo->transportPropertiesDict()
-    );
-
-    // Construct incompressible turbulence model
-    autoPtr<incompressible::turbulenceModel> turbulence
+        "rho",
+        runTime.timeName(),
+        mesh,
+        IOobject::READ_IF_PRESENT,
+        IOobject::AUTO_WRITE
+    ),
+    alpha1*rho1 + alpha2*rho2,
+    alpha1.boundaryField().types()
+);
+rho.oldTime();
+
+// Construct interface from alpha1 distribution
+interfaceProperties interface
+(
+    alpha1,
+    U,
+    thermo->transportPropertiesDict()
+);
+
+// Construct incompressible turbulence model
+autoPtr<incompressible::turbulenceModel> turbulence
+(
+    incompressible::turbulenceModel::New(U, phi, thermo())
+);
+
+#include "readGravitationalAcceleration.H"
+#include "readhRef.H"
+#include "gh.H"
+
+volScalarField& p = thermo->p();
+
+label pRefCell = 0;
+scalar pRefValue = 0.0;
+setRefCell
+(
+    p,
+    p_rgh,
+    pimple.dict(),
+    pRefCell,
+    pRefValue
+);
+
+if (p_rgh.needReference())
+{
+    p += dimensionedScalar
     (
-        incompressible::turbulenceModel::New(U, phi, thermo())
+        "p",
+        p.dimensions(),
+        pRefValue - getRefCellValue(p, pRefCell)
     );
+    p_rgh = p - rho*gh;
+}
 
+mesh.setFluxRequired(p_rgh.name());
+mesh.setFluxRequired(alpha1.name());
 
-    Info<< "Calculating field g.h\n" << endl;
-    volScalarField gh("gh", g & mesh.C());
-    surfaceScalarField ghf("ghf", g & mesh.Cf());
-
-    volScalarField& p = thermo->p();
+// Turbulent Prandtl number
+dimensionedScalar Prt("Prt", dimless, thermo->transportPropertiesDict());
 
-    label pRefCell = 0;
-    scalar pRefValue = 0.0;
-    setRefCell
+volScalarField kappaEff
+(
+    IOobject
     (
-        p,
-        p_rgh,
-        pimple.dict(),
-        pRefCell,
-        pRefValue
-    );
-
-    if (p_rgh.needReference())
-    {
-        p += dimensionedScalar
-        (
-            "p",
-            p.dimensions(),
-            pRefValue - getRefCellValue(p, pRefCell)
-        );
-        p_rgh = p - rho*gh;
-    }
-
-    // Turbulent Prandtl number
-    dimensionedScalar Prt("Prt", dimless, thermo->transportPropertiesDict());
-
-    volScalarField kappaEff
+        "kappaEff",
+        runTime.timeName(),
+        mesh,
+        IOobject::NO_READ,
+        IOobject::NO_WRITE
+    ),
+    thermo->kappa()
+);
+
+Info<< "Creating field kinetic energy K\n" << endl;
+volScalarField K("K", 0.5*magSqr(U));
+
+Info<< "Creating field pDivU\n" << endl;
+volScalarField pDivU
+(
+    IOobject
     (
-        IOobject
-        (
-            "kappaEff",
-            runTime.timeName(),
-            mesh,
-            IOobject::NO_READ,
-            IOobject::NO_WRITE
-        ),
-        thermo->kappa()
-    );
-
-    Info<< "Creating field kinetic energy K\n" << endl;
-    volScalarField K("K", 0.5*magSqr(U));
+        "pDivU",
+        runTime.timeName(),
+        mesh
+    ),
+    mesh,
+    dimensionedScalar("pDivU", p.dimensions()/dimTime, 0)
+);
 
-    Info<< "Creating field pDivU\n" << endl;
-    volScalarField pDivU
-    (
-        IOobject
-        (
-            "pDivU",
-            runTime.timeName(),
-            mesh
-        ),
-        mesh,
-        dimensionedScalar("pDivU", p.dimensions()/dimTime, 0)
-    );
diff --git a/applications/solvers/multiphase/interCondensingEvaporatingFoam/interCondensatingEvaporatingFoam.C b/applications/solvers/multiphase/interCondensingEvaporatingFoam/interCondensatingEvaporatingFoam.C
index 362e4d013254778bbfc6f7439dffecd3de1ab37c..d9981659b8db1fe93eef0fb9646e3623c5dfee67 100644
--- a/applications/solvers/multiphase/interCondensingEvaporatingFoam/interCondensatingEvaporatingFoam.C
+++ b/applications/solvers/multiphase/interCondensingEvaporatingFoam/interCondensatingEvaporatingFoam.C
@@ -56,19 +56,22 @@ Description
 
 int main(int argc, char *argv[])
 {
+    #include "postProcess.H"
+
     #include "setRootCase.H"
     #include "createTime.H"
     #include "createMesh.H"
-
-    pimpleControl pimple(mesh);
-
-    #include "readGravitationalAcceleration.H"
+    #include "createControl.H"
     #include "createFields.H"
     #include "createFvOptions.H"
     #include "createTimeControls.H"
     #include "CourantNo.H"
     #include "setInitialDeltaT.H"
 
+    volScalarField& T = thermo->T();
+    volScalarField& e = thermo->he();
+    e.oldTime();
+
     turbulence->validate();
 
     // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //