diff --git a/applications/solvers/combustion/SLFMFoam/SLFMFoam.C b/applications/solvers/combustion/SLFMFoam/SLFMFoam.C
index 0dea1ccf0..8d1d48931 100644
--- a/applications/solvers/combustion/SLFMFoam/SLFMFoam.C
+++ b/applications/solvers/combustion/SLFMFoam/SLFMFoam.C
@@ -99,10 +99,9 @@ int main(int argc, char *argv[])
 
         turbulence->correct();
 
-        if(runTime.write() == true)
+        if (runTime.write())
         {
             #include "updateYi.H" //update species
-            rho.write();
         }
 
         Info<< "ExecutionTime = " << runTime.elapsedCpuTime() << " s"
diff --git a/applications/solvers/combustion/SLFMFoam/createFields.H b/applications/solvers/combustion/SLFMFoam/createFields.H
index a5383e293..deb56052c 100644
--- a/applications/solvers/combustion/SLFMFoam/createFields.H
+++ b/applications/solvers/combustion/SLFMFoam/createFields.H
@@ -26,15 +26,31 @@ if (!composition.contains(inertSpecie))
 volScalarField& p = thermo.p();
 volScalarField& T = thermo.T();
 
+wordList rhoBoundaryTypes
+(
+	T.boundaryField().size(),
+	zeroGradientFvPatchScalarField::typeName
+);
+
+forAll(composition.Y(inertSpecie).boundaryField(), patchi)
+{
+	rhoBoundaryTypes[patchi]
+        = composition.Y(inertSpecie).boundaryField()[patchi].type();
+}
+
+
 volScalarField rho
 (
     IOobject
     (
         "rho",
         runTime.timeName(),
-        mesh
+        mesh,
+        IOobject::NO_READ,
+        IOobject::AUTO_WRITE
     ),
-    thermo.rho()
+    thermo.rho(),
+    rhoBoundaryTypes
 );
 
 Info<< "\nReading field U\n" << endl;
@@ -92,11 +108,6 @@ autoPtr<compressible::turbulenceModel> turbulence
 // Set the turbulence into the combustion model
 combustion->setTurbulence(turbulence());
 
-
-Info<< "Creating field kinetic energy K\n" << endl;
-volScalarField K("K", 0.5*magSqr(U));
-
-
 Info<< "Creating multi-variate interpolation scheme\n" << endl;
 multivariateSurfaceInterpolationScheme<scalar>::fieldTable fields;
 
diff --git a/applications/solvers/combustion/SLFMFoam/updateT_RHO.H b/applications/solvers/combustion/SLFMFoam/updateT_RHO.H
index 9bf4cfb73..4e35bedcd 100644
--- a/applications/solvers/combustion/SLFMFoam/updateT_RHO.H
+++ b/applications/solvers/combustion/SLFMFoam/updateT_RHO.H
@@ -58,18 +58,3 @@ forAll(rho, cellI)
 T.correctBoundaryConditions();
 rho.correctBoundaryConditions();
 
-//update outlet boundary value of rho
-const fvPatchList& patches = mesh.boundary();
-forAll(patches, patchI)
-{
-    const fvPatch& curPatch = patches[patchI];
-    if(curPatch.name() == outletName)
-    {
-        forAll(curPatch, faceI)
-        {
-            label cellI = curPatch.faceCells()[faceI];
-            rho.boundaryFieldRef()[patchI][faceI] = rho[cellI];
-        }
-    }
-}
-