diff --git a/src/OpenFOAM/primitives/functions/DataEntry/DataEntry/DataEntry.H b/src/OpenFOAM/primitives/functions/DataEntry/DataEntry/DataEntry.H
index 9cfdb66e25..a9b1686d58 100644
--- a/src/OpenFOAM/primitives/functions/DataEntry/DataEntry/DataEntry.H
+++ b/src/OpenFOAM/primitives/functions/DataEntry/DataEntry/DataEntry.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -142,62 +142,60 @@ public:
             virtual void convertTimeBase(const Time& t);
 
 
-public:
+        // Evaluation
+
+            //- Return value as a function of (scalar) independent variable
+            virtual Type value(const scalar x) const;
+
+            //- Return value as a function of (scalar) independent variable
+            virtual tmp<Field<Type> > value(const scalarField& x) const;
+
+            //- Integrate between two (scalar) values
+            virtual Type integrate(const scalar x1, const scalar x2) const;
+
+            //- Integrate between two (scalar) values
+            virtual tmp<Field<Type> > integrate
+            (
+                const scalarField& x1,
+                const scalarField& x2
+            ) const;
+
+            //- Return dimensioned type
+            virtual dimensioned<Type> dimValue(const scalar x) const;
+
+            //- Return dimensioned type as a function of (scalar)
+            virtual tmp<Field<dimensioned<Type> > > dimValue
+            (
+                const scalarField& x
+            ) const;
+
+            //- Integrate between two scalars and return a dimensioned type
+            virtual dimensioned<Type> dimIntegrate
+            (
+                const scalar x1,
+                const scalar x2
+            ) const;
+
+            //- Integrate between two scalar fields and return a field of
+            //  dimensioned type
+            virtual tmp<Field<dimensioned<Type> > > dimIntegrate
+            (
+                const scalarField& x1,
+                const scalarField& x2
+            ) const;
 
 
-    // Evaluation
+        // I/O
 
-        //- Return value as a function of (scalar) independent variable
-        virtual Type value(const scalar x) const;
+            //- Ostream Operator
+            friend Ostream& operator<< <Type>
+            (
+                Ostream& os,
+                const DataEntry<Type>& de
+            );
 
-        //- Return value as a function of (scalar) independent variable
-        virtual tmp<Field<Type> > value(const scalarField& x) const;
-
-        //- Integrate between two (scalar) values
-        virtual Type integrate(const scalar x1, const scalar x2) const;
-
-        //- Integrate between two (scalar) values
-        virtual tmp<Field<Type> > integrate
-        (
-            const scalarField& x1,
-            const scalarField& x2
-        ) const;
-
-        //- Return dimensioned type
-        virtual dimensioned<Type> dimValue(const scalar x) const;
-
-        //- Return dimensioned type as a function of (scalar)
-        virtual tmp<Field<dimensioned<Type> > > dimValue
-        (
-            const scalarField& x
-        ) const;
-
-        //- Integrate between two scalars and returns a dimensioned type
-        virtual dimensioned<Type> dimIntegrate
-        (
-            const scalar x1,
-            const scalar x2
-        ) const;
-
-        //- Integrate between two scalars and returns list of dimensioned type
-        virtual tmp<Field<dimensioned<Type> > > dimIntegrate
-        (
-            const scalarField& x1,
-            const scalarField& x2
-        ) const;
-
-
-    // I/O
-
-        //- Ostream Operator
-        friend Ostream& operator<< <Type>
-        (
-            Ostream& os,
-            const DataEntry<Type>& de
-        );
-
-        //- Write in dictionary format
-        virtual void writeData(Ostream& os) const;
+            //- Write in dictionary format
+            virtual void writeData(Ostream& os) const;
 };
 
 
diff --git a/src/finiteVolume/cfdTools/general/adjustPhi/adjustPhi.C b/src/finiteVolume/cfdTools/general/adjustPhi/adjustPhi.C
index 206e7c4e95..27622aeeb7 100644
--- a/src/finiteVolume/cfdTools/general/adjustPhi/adjustPhi.C
+++ b/src/finiteVolume/cfdTools/general/adjustPhi/adjustPhi.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -56,11 +56,7 @@ bool Foam::adjustPhi
 
             if (!isA<processorFvsPatchScalarField>(phip))
             {
-                if
-                (
-                    Up.fixesValue()
-                 && !isA<inletOutletFvPatchVectorField>(Up)
-                )
+                if (Up.fixesValue() && !isA<inletOutletFvPatchVectorField>(Up))
                 {
                     forAll(phip, i)
                     {
@@ -113,8 +109,12 @@ bool Foam::adjustPhi
         {
             FatalErrorIn
             (
-                "adjustPhi(surfaceScalarField& phi, const volVectorField& U,"
-                "const volScalarField& p"
+                "adjustPhi"
+                "("
+                    "surfaceScalarField&, "
+                    "const volVectorField&,"
+                    "volScalarField&"
+                ")"
             )   << "Continuity error cannot be removed by adjusting the"
                    " outflow.\nPlease check the velocity boundary conditions"
                    " and/or run potentialFoam to initialise the outflow." << nl
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C
index cebed76580..97c3c819a0 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/ConeNozzleInjection/ConeNozzleInjection.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -167,7 +167,8 @@ Foam::ConeNozzleInjection<CloudType>::ConeNozzleInjection
             "Foam::ConeNozzleInjection<CloudType>::ConeNozzleInjection"
             "("
                 "const dictionary&, "
-                "CloudType&"
+                "CloudType&, "
+                "const word&"
             ")"
         )<< "innerNozzleDiameter >= outerNozzleDiameter" << nl
          << exit(FatalError);
@@ -369,6 +370,7 @@ void Foam::ConeNozzleInjection<CloudType>::setPositionAndCell
                     "const scalar, "
                     "vector&, "
                     "label&, "
+                    "label&, "
                     "label&"
                 ")"
             )<< "Unknown injectionMethod type" << nl
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.C b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.C
index c40a95fbac..293fd7d740 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,7 @@ Foam::KinematicLookupTableInjection<CloudType>::KinematicLookupTableInjection
     (
         readScalar(this->coeffDict().lookup("parcelsPerSecond"))
     ),
+    randomise_(readBool(this->coeffDict().lookup("randomise"))),
     injectors_
     (
         IOobject
@@ -87,6 +88,7 @@ Foam::KinematicLookupTableInjection<CloudType>::KinematicLookupTableInjection
     inputFileName_(im.inputFileName_),
     duration_(im.duration_),
     parcelsPerSecond_(im.parcelsPerSecond_),
+    randomise_(im.randomise_),
     injectors_(im.injectors_),
     injectorCells_(im.injectorCells_),
     injectorTetFaces_(im.injectorTetFaces_),
@@ -177,7 +179,16 @@ void Foam::KinematicLookupTableInjection<CloudType>::setPositionAndCell
     label& tetPtI
 )
 {
-    label injectorI = parcelI*injectorCells_.size()/nParcels;
+    label injectorI = 0;
+    if (randomise_)
+    {
+        cachedRandom& rnd = this->owner().rndGen();
+        injectorI = rnd.position<label>(0, injectorCells_.size() - 1);
+    }
+    else
+    {
+        injectorI = parcelI*injectorCells_.size()/nParcels;
+    }
 
     position = injectors_[injectorI].x();
     cellOwner = injectorCells_[injectorI];
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.H b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.H
index d156f4793a..f17410bc1d 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.H
+++ b/src/lagrangian/intermediate/submodels/Kinematic/InjectionModel/KinematicLookupTableInjection/KinematicLookupTableInjection.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -78,6 +78,9 @@ class KinematicLookupTableInjection
         //- Number of parcels per injector - common to all injection sources
         const scalar parcelsPerSecond_;
 
+        //- Flag to indicate to randomise injection positions
+        bool randomise_;
+
         //- List of injectors
         kinematicParcelInjectionDataIOList injectors_;
 
diff --git a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Lift/LiftForce/LiftForce.C b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Lift/LiftForce/LiftForce.C
index 290d2f3c3d..a4faceeec7 100644
--- a/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Lift/LiftForce/LiftForce.C
+++ b/src/lagrangian/intermediate/submodels/Kinematic/ParticleForces/Lift/LiftForce/LiftForce.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -141,7 +141,7 @@ Foam::forceSuSp Foam::LiftForce<CloudType>::calcCoupled
 
     scalar Cl = this->Cl(p, curlUc, Re, muc);
 
-    value.Su() = mass/p.rho()*p.d()/2.0*p.rhoc()*Cl*((p.Uc() - p.U())^curlUc);
+    value.Su() = mass/p.rho()*p.rhoc()*Cl*((p.Uc() - p.U())^curlUc);
 
     return value;
 }
diff --git a/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.C b/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.C
index 711b3cec7a..1259da7e37 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.C
+++ b/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -42,6 +42,7 @@ Foam::ReactingLookupTableInjection<CloudType>::ReactingLookupTableInjection
     (
         readScalar(this->coeffDict().lookup("parcelsPerSecond"))
     ),
+    randomise_(readBool(this->coeffDict().lookup("randomise"))),
     injectors_
     (
         IOobject
@@ -86,6 +87,7 @@ Foam::ReactingLookupTableInjection<CloudType>::ReactingLookupTableInjection
     inputFileName_(im.inputFileName_),
     duration_(im.duration_),
     parcelsPerSecond_(im.parcelsPerSecond_),
+    randomise_(im.randomise_),
     injectors_(im.injectors_),
     injectorCells_(im.injectorCells_),
     injectorTetFaces_(im.injectorTetFaces_),
@@ -176,7 +178,16 @@ void Foam::ReactingLookupTableInjection<CloudType>::setPositionAndCell
     label& tetPtI
 )
 {
-    label injectorI = parcelI*injectorCells_.size()/nParcels;
+    label injectorI = 0;
+    if (randomise_)
+    {
+        cachedRandom& rnd = this->owner().rndGen();
+        injectorI = rnd.position<label>(0, injectorCells_.size() - 1);
+    }
+    else
+    {
+        injectorI = parcelI*injectorCells_.size()/nParcels;
+    }
 
     position = injectors_[injectorI].x();
     cellOwner = injectorCells_[injectorI];
diff --git a/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.H b/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.H
index ef062fd72b..654109c509 100644
--- a/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.H
+++ b/src/lagrangian/intermediate/submodels/Reacting/InjectionModel/ReactingLookupTableInjection/ReactingLookupTableInjection.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -81,6 +81,9 @@ class ReactingLookupTableInjection
         //- Number of parcels per injector - common to all injection sources
         const scalar parcelsPerSecond_;
 
+        //- Flag to indicate to randomise injection positions
+        bool randomise_;
+
         //- List of injectors
         reactingParcelInjectionDataIOList injectors_;
 
diff --git a/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.C b/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.C
index bd630d879d..179b846caf 100644
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.C
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,7 @@ ReactingMultiphaseLookupTableInjection
     (
         readScalar(this->coeffDict().lookup("parcelsPerSecond"))
     ),
+    randomise_(readBool(this->coeffDict().lookup("randomise"))),
     injectors_
     (
         IOobject
@@ -88,6 +89,7 @@ ReactingMultiphaseLookupTableInjection
     inputFileName_(im.inputFileName_),
     duration_(im.duration_),
     parcelsPerSecond_(im.parcelsPerSecond_),
+    randomise_(im.randomise_),
     injectors_(im.injectors_),
     injectorCells_(im.injectorCells_),
     injectorTetFaces_(im.injectorTetFaces_),
@@ -182,7 +184,16 @@ void Foam::ReactingMultiphaseLookupTableInjection<CloudType>::setPositionAndCell
     label& tetPtI
 )
 {
-    label injectorI = parcelI*injectorCells_.size()/nParcels;
+    label injectorI = 0;
+    if (randomise_)
+    {
+        cachedRandom& rnd = this->owner().rndGen();
+        injectorI = rnd.position<label>(0, injectorCells_.size() - 1);
+    }
+    else
+    {
+        injectorI = parcelI*injectorCells_.size()/nParcels;
+    }
 
     position = injectors_[injectorI].x();
     cellOwner = injectorCells_[injectorI];
diff --git a/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.H b/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.H
index 7e6c328558..81f41d42de 100644
--- a/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.H
+++ b/src/lagrangian/intermediate/submodels/ReactingMultiphase/InjectionModel/ReactingMultiphaseLookupTableInjection/ReactingMultiphaseLookupTableInjection.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -84,6 +84,9 @@ class ReactingMultiphaseLookupTableInjection
         //- Number of parcels per injector - common to all injection sources
         const scalar parcelsPerSecond_;
 
+        //- Flag to indicate to randomise injection positions
+        bool randomise_;
+
         //- List of injectors
         reactingMultiphaseParcelInjectionDataIOList injectors_;
 
diff --git a/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.C b/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.C
index fdc0a3b4a3..682fc9c3ea 100644
--- a/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.C
+++ b/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,7 @@ Foam::ThermoLookupTableInjection<CloudType>::ThermoLookupTableInjection
     (
         readScalar(this->coeffDict().lookup("parcelsPerSecond"))
     ),
+    randomise_(readBool(this->coeffDict().lookup("randomise"))),
     injectors_
     (
         IOobject
@@ -87,6 +88,7 @@ Foam::ThermoLookupTableInjection<CloudType>::ThermoLookupTableInjection
     inputFileName_(im.inputFileName_),
     duration_(im.duration_),
     parcelsPerSecond_(im.parcelsPerSecond_),
+    randomise_(im.randomise_),
     injectors_(im.injectors_),
     injectorCells_(im.injectorCells_),
     injectorTetFaces_(im.injectorTetFaces_),
@@ -177,7 +179,16 @@ void Foam::ThermoLookupTableInjection<CloudType>::setPositionAndCell
     label& tetPtI
 )
 {
-    label injectorI = parcelI*injectorCells_.size()/nParcels;
+    label injectorI = 0;
+    if (randomise_)
+    {
+        cachedRandom& rnd = this->owner().rndGen();
+        injectorI = rnd.position<label>(0, injectorCells_.size() - 1);
+    }
+    else
+    {
+        injectorI = parcelI*injectorCells_.size()/nParcels;
+    }
 
     position = injectors_[injectorI].x();
     cellOwner = injectorCells_[injectorI];
diff --git a/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.H b/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.H
index efb83c3614..42c48f290d 100644
--- a/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.H
+++ b/src/lagrangian/intermediate/submodels/Thermodynamic/InjectionModel/ThermoLookupTableInjection/ThermoLookupTableInjection.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -80,6 +80,9 @@ class ThermoLookupTableInjection
         //- Number of parcels per injector - common to all injection sources
         const scalar parcelsPerSecond_;
 
+        //- Flag to indicate to randomise injection positions
+        bool randomise_;
+
         //- List of injectors
         kinematicParcelInjectionDataIOList injectors_;
 
diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C
index d714747310..f19b7f69d1 100644
--- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C
+++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.C
@@ -311,6 +311,7 @@ Foam::fieldAverage::fieldAverage
     prevTimeIndex_(-1),
     resetOnRestart_(false),
     resetOnOutput_(false),
+    initialised_(false),
     faItems_(),
     meanScalarFields_(),
     meanVectorFields_(),
@@ -361,11 +362,7 @@ void Foam::fieldAverage::read(const dictionary& dict)
         dict.readIfPresent("resetOnOutput", resetOnOutput_);
         dict.lookup("fields") >> faItems_;
 
-        initialize();
         readAveragingProperties();
-
-        // ensure first averaging works unconditionally
-        prevTimeIndex_ = -1;
     }
 }
 
@@ -374,6 +371,16 @@ void Foam::fieldAverage::execute()
 {
     if (active_)
     {
+        if (!initialised_)
+        {
+            initialize();
+
+            // ensure first averaging works unconditionally
+            prevTimeIndex_ = -1;
+
+            initialised_ = true;
+        }
+
         calcAverages();
     }
 }
diff --git a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.H b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.H
index f8b689757d..4a6d75e06b 100644
--- a/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.H
+++ b/src/postProcessing/functionObjects/field/fieldAverage/fieldAverage/fieldAverage.H
@@ -171,6 +171,9 @@ protected:
         //- Reset the averaging process on output flag
         Switch resetOnOutput_;
 
+        //- Initialised flag
+        bool initialised_;
+
         //- List of field average items, describing what averages to be
         //  calculated and output
         List<fieldAverageItem> faItems_;
diff --git a/src/postProcessing/functionObjects/forces/forces/forces.C b/src/postProcessing/functionObjects/forces/forces/forces.C
index c159b652e8..dd164d99a4 100644
--- a/src/postProcessing/functionObjects/forces/forces/forces.C
+++ b/src/postProcessing/functionObjects/forces/forces/forces.C
@@ -28,16 +28,12 @@ License
 #include "dictionary.H"
 #include "Time.H"
 #include "wordReList.H"
-
-#include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H"
-#include "incompressible/RAS/RASModel/RASModel.H"
-#include "incompressible/LES/LESModel/LESModel.H"
-
-#include "fluidThermo.H"
-#include "compressible/RAS/RASModel/RASModel.H"
-#include "compressible/LES/LESModel/LESModel.H"
-
+#include "fvcGrad.H"
 #include "porosityModel.H"
+#include "fluidThermo.H"
+#include "incompressible/turbulenceModel/turbulenceModel.H"
+#include "compressible/turbulenceModel/turbulenceModel.H"
+#include "incompressible/singlePhaseTransportModel/singlePhaseTransportModel.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -70,38 +66,27 @@ void Foam::forces::writeFileHeader(const label i)
 
 Foam::tmp<Foam::volSymmTensorField> Foam::forces::devRhoReff() const
 {
-    if (obr_.foundObject<compressible::RASModel>("RASProperties"))
-    {
-        const compressible::RASModel& ras
-            = obr_.lookupObject<compressible::RASModel>("RASProperties");
+    typedef compressible::turbulenceModel cmpTurbModel;
+    typedef incompressible::turbulenceModel icoTurbModel;
 
-        return ras.devRhoReff();
-    }
-    else if (obr_.foundObject<incompressible::RASModel>("RASProperties"))
+    if (obr_.foundObject<cmpTurbModel>(cmpTurbModel::typeName))
     {
-        const incompressible::RASModel& ras
-            = obr_.lookupObject<incompressible::RASModel>("RASProperties");
+        const cmpTurbModel& turb =
+            obr_.lookupObject<cmpTurbModel>(cmpTurbModel::typeName);
 
-        return rho()*ras.devReff();
+        return turb.devRhoReff();
     }
-    else if (obr_.foundObject<compressible::LESModel>("LESProperties"))
+    else if (obr_.foundObject<icoTurbModel>(icoTurbModel::typeName))
     {
-        const compressible::LESModel& les =
-        obr_.lookupObject<compressible::LESModel>("LESProperties");
+        const incompressible::turbulenceModel& turb =
+            obr_.lookupObject<icoTurbModel>(icoTurbModel::typeName);
 
-        return les.devRhoReff();
+        return rho()*turb.devReff();
     }
-    else if (obr_.foundObject<incompressible::LESModel>("LESProperties"))
-    {
-        const incompressible::LESModel& les
-            = obr_.lookupObject<incompressible::LESModel>("LESProperties");
-
-        return rho()*les.devReff();
-    }
-    else if (obr_.foundObject<fluidThermo>("thermophysicalProperties"))
+    else if (obr_.foundObject<fluidThermo>(fluidThermo::typeName))
     {
         const fluidThermo& thermo =
-             obr_.lookupObject<fluidThermo>("thermophysicalProperties");
+            obr_.lookupObject<fluidThermo>(fluidThermo::typeName);
 
         const volVectorField& U = obr_.lookupObject<volVectorField>(UName_);
 
diff --git a/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.C b/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.C
index 607d52e503..181099c007 100644
--- a/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.C
+++ b/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.C
@@ -70,7 +70,7 @@ Foam::dimensionedScalar Foam::pressureTools::rhoScale
     }
     else
     {
-        return dimensionedScalar("rhoRef", dimDensity, rhoRef_);
+        return dimensionedScalar("rhoRef", dimDensity, rhoInf_);
     }
 }
 
@@ -100,7 +100,7 @@ Foam::tmp<Foam::volScalarField> Foam::pressureTools::rho
                         IOobject::NO_WRITE
                     ),
                     p.mesh(),
-                    dimensionedScalar("zero", dimDensity, rhoRef_)
+                    dimensionedScalar("zero", dimDensity, rhoInf_)
                 )
             );
     }
@@ -193,7 +193,6 @@ Foam::pressureTools::pressureTools
     pName_("p"),
     UName_("U"),
     rhoName_("rho"),
-    rhoRef_(1.0),
     calcTotal_(false),
     pRef_(0.0),
     calcCoeff_(false),
@@ -242,7 +241,7 @@ void Foam::pressureTools::read(const dictionary& dict)
 
         if (p.dimensions() != dimPressure)
         {
-            dict.lookup("rhoRef") >> rhoRef_;
+            dict.lookup("rhoRef") >> rhoInf_;
         } 
 
         dict.lookup("calcTotal") >> calcTotal_;
diff --git a/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.H b/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.H
index 11f490429d..82a44921ed 100644
--- a/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.H
+++ b/src/postProcessing/functionObjects/utilities/pressureTools/pressureTools.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -92,7 +92,6 @@ Description
     \table
         Property     | Description             | Required    | Default value
         type         | type name: pressureTools| yes         |
-        rhoRef       | Reference density for incompressible cases | no | 1
         calcTotal    | Calculate total coefficient | yes     |
         pRef         | Reference pressure for total pressure | no | 0.0
         calcCoeff    | Calculate pressure coefficient | yes  |
@@ -150,9 +149,6 @@ class pressureTools
         //- Name of density field, default is "rho"
         word rhoName_;
 
-        //- Reference density employed for incompressible cases
-        scalar rhoRef_;
-
 
         // Total pressure calculation
 
diff --git a/src/regionModels/surfaceFilmModels/Make/files b/src/regionModels/surfaceFilmModels/Make/files
index 1151e2f7fd..f14a789338 100644
--- a/src/regionModels/surfaceFilmModels/Make/files
+++ b/src/regionModels/surfaceFilmModels/Make/files
@@ -24,6 +24,10 @@ $(KINEMATICMODELS)/injectionModel/drippingInjection/drippingInjection.C
 $(KINEMATICMODELS)/injectionModel/removeInjection/removeInjection.C
 $(KINEMATICMODELS)/injectionModel/curvatureSeparation/curvatureSeparation.C
 
+$(KINEMATICMODELS)/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C
+$(KINEMATICMODELS)/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModelNew.C
+$(KINEMATICMODELS)/filmTurbulenceModel/laminar/laminar.C
+
 THERMOMODELS=submodels/thermo
 $(THERMOMODELS)/phaseChangeModel/phaseChangeModel/phaseChangeModel.C
 $(THERMOMODELS)/phaseChangeModel/phaseChangeModel/phaseChangeModelNew.C
diff --git a/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.C b/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.C
index fecefcd2e3..00a3ac39af 100644
--- a/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.C
+++ b/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.C
@@ -163,7 +163,7 @@ tmp<volScalarField> kinematicSingleLayer::pu()
         (
             IOobject
             (
-                "pu",
+                typeName + ":pu",
                 time_.timeName(),
                 regionMesh(),
                 IOobject::NO_READ,
@@ -185,7 +185,7 @@ tmp<volScalarField> kinematicSingleLayer::pp()
         (
             IOobject
             (
-                "pp",
+                typeName + ":pp",
                 time_.timeName(),
                 regionMesh(),
                 IOobject::NO_READ,
@@ -216,6 +216,8 @@ void kinematicSingleLayer::updateSubmodels()
     // Update source fields
     const dimensionedScalar deltaT = time().deltaT();
     rhoSp_ += cloudMassTrans_/magSf()/deltaT;
+
+    turbulence_->correct();
 }
 
 
@@ -282,9 +284,7 @@ void kinematicSingleLayer::updateSurfaceVelocities()
     Uw_ -= nHat()*(Uw_ & nHat());
     Uw_.correctBoundaryConditions();
 
-    // apply quadratic profile to surface velocity (scale by sqrt(2))
-    Us_ = 1.414*U_;
-    Us_.correctBoundaryConditions();
+    Us_ = turbulence_->Us();
 }
 
 
@@ -309,6 +309,7 @@ tmp<Foam::fvVectorMatrix> kinematicSingleLayer::solveMomentum
 //      - fvm::SuSp(rhoSp_, U_)
       - rhoSp_*U_
       + forces_.correct(U_)
+      + turbulence_->Su(U_)
     );
 
     fvVectorMatrix& UEqn = tUEqn();
@@ -779,6 +780,8 @@ kinematicSingleLayer::kinematicSingleLayer
 
     injection_(*this, coeffs_),
 
+    turbulence_(filmTurbulenceModel::New(*this, coeffs_)),
+
     forces_(*this, coeffs_),
 
     addedMassTotal_(0.0)
@@ -1020,7 +1023,7 @@ tmp<volScalarField> kinematicSingleLayer::primaryMassTrans() const
         (
             IOobject
             (
-                "kinematicSingleLayer::primaryMassTrans",
+                typeName + ":primaryMassTrans",
                 time().timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,
@@ -1075,7 +1078,7 @@ tmp<DimensionedField<scalar, volMesh> > kinematicSingleLayer::Srho() const
         (
             IOobject
             (
-                "kinematicSingleLayer::Srho",
+                typeName + ":Srho",
                 time().timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,
@@ -1100,7 +1103,7 @@ tmp<DimensionedField<scalar, volMesh> > kinematicSingleLayer::Srho
         (
             IOobject
             (
-                "kinematicSingleLayer::Srho(" + Foam::name(i) + ")",
+                typeName + ":Srho(" + Foam::name(i) + ")",
                 time().timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,
@@ -1122,7 +1125,7 @@ tmp<DimensionedField<scalar, volMesh> > kinematicSingleLayer::Sh() const
         (
             IOobject
             (
-                "kinematicSingleLayer::Sh",
+                typeName + ":Sh",
                 time().timeName(),
                 primaryMesh(),
                 IOobject::NO_READ,
diff --git a/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.H b/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.H
index 079999d0cd..2c9e1020a3 100644
--- a/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.H
+++ b/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayer.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,6 +43,7 @@ SourceFiles
 
 #include "injectionModelList.H"
 #include "forceList.H"
+#include "filmTurbulenceModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -197,6 +198,9 @@ protected:
             //- Cloud injection
             injectionModelList injection_;
 
+            //- Turbulence model
+            autoPtr<filmTurbulenceModel> turbulence_;
+
             //- List of film forces
             forceList forces_;
 
@@ -444,6 +448,9 @@ public:
             //- Injection
             inline injectionModelList& injection();
 
+            //- Turbulence
+            inline const filmTurbulenceModel& turbulence() const;
+
 
         // Helper functions
 
diff --git a/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayerI.H b/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayerI.H
index 6a77b614bc..15f6398ad8 100644
--- a/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayerI.H
+++ b/src/regionModels/surfaceFilmModels/kinematicSingleLayer/kinematicSingleLayerI.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -169,6 +169,12 @@ inline injectionModelList& kinematicSingleLayer::injection()
 }
 
 
+inline const filmTurbulenceModel& kinematicSingleLayer::turbulence() const
+{
+    return turbulence_();
+}
+
+
 inline tmp<volScalarField> kinematicSingleLayer::mass() const
 {
     return rho_*delta_*magSf();
diff --git a/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C
new file mode 100644
index 0000000000..a77307e881
--- /dev/null
+++ b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C
@@ -0,0 +1,73 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "filmTurbulenceModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace regionModels
+{
+namespace surfaceFilmModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(filmTurbulenceModel, 0);
+defineRunTimeSelectionTable(filmTurbulenceModel, dictionary);
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+filmTurbulenceModel::filmTurbulenceModel(const surfaceFilmModel& owner)
+:
+    subModelBase(owner)
+{}
+
+
+filmTurbulenceModel::filmTurbulenceModel
+(
+    const word& type,
+    const surfaceFilmModel& owner,
+    const dictionary& dict
+)
+:
+    subModelBase(type, owner, dict)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+filmTurbulenceModel::~filmTurbulenceModel()
+{}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace surfaceFilmModels
+} // End namespace regionModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H
new file mode 100644
index 0000000000..571bf83e8d
--- /dev/null
+++ b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H
@@ -0,0 +1,148 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::regionModels::surfaceFilmModels::filmTurbulenceModel
+
+Description
+    Base class for film turbulence models
+
+SourceFiles
+    filmTurbulenceModel.C
+    filmTurbulenceModelNew.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef filmTurbulenceModel_H
+#define filmTurbulenceModel_H
+
+#include "subModelBase.H"
+#include "runTimeSelectionTables.H"
+#include "fvMatricesFwd.H"
+#include "volFieldsFwd.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace regionModels
+{
+namespace surfaceFilmModels
+{
+
+/*---------------------------------------------------------------------------*\
+                       Class filmTurbulenceModel Declaration
+\*---------------------------------------------------------------------------*/
+
+class filmTurbulenceModel
+:
+    public subModelBase
+{
+private:
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        filmTurbulenceModel(const filmTurbulenceModel&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const filmTurbulenceModel&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("filmTurbulenceModel");
+
+
+    // Declare runtime constructor selection table
+
+         declareRunTimeSelectionTable
+         (
+             autoPtr,
+             filmTurbulenceModel,
+             dictionary,
+             (
+                const surfaceFilmModel& owner,
+                const dictionary& dict
+             ),
+             (owner, dict)
+         );
+
+    // Constructors
+
+        //- Construct null
+        filmTurbulenceModel(const surfaceFilmModel& owner);
+
+        //- Construct from type name, dictionary and surface film model
+        filmTurbulenceModel
+        (
+            const word& type,
+            const surfaceFilmModel& owner,
+            const dictionary& dict
+        );
+
+
+    // Selectors
+
+        //- Return a reference to the selected injection model
+        static autoPtr<filmTurbulenceModel> New
+        (
+            const surfaceFilmModel& owner,
+            const dictionary& dict
+        );
+
+
+    //- Destructor
+    virtual ~filmTurbulenceModel();
+
+
+    // Member Functions
+
+        // Evolution
+
+            //- Return the film surface velocity
+            virtual tmp<volVectorField> Us() const = 0;
+
+            //- Return the film turbulence viscosity
+            virtual tmp<volScalarField> mut() const = 0;
+
+            //- Correct/update the model
+            virtual void correct() = 0;
+
+            //- Return the source for the film momentum equation
+            virtual tmp<fvVectorMatrix> Su(volVectorField& U) const = 0;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace surfaceFilmModels
+} // End namespace regionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModelNew.C b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModelNew.C
new file mode 100644
index 0000000000..ba939ef798
--- /dev/null
+++ b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModelNew.C
@@ -0,0 +1,77 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "filmTurbulenceModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace regionModels
+{
+namespace surfaceFilmModels
+{
+
+// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
+
+autoPtr<filmTurbulenceModel> filmTurbulenceModel::New
+(
+    const surfaceFilmModel& model,
+    const dictionary& dict
+)
+{
+    const word modelType(dict.lookup("turbulence"));
+
+    Info<< "        " << modelType << endl;
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(modelType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "filmTurbulenceModel::New"
+            "("
+                "const surfaceFilmModel&, "
+                "const dictionary&"
+            ")"
+        )   << "Unknown filmTurbulenceModel type " << modelType
+            << nl << nl << "Valid filmTurbulenceModel types are:" << nl
+            << dictionaryConstructorTablePtr_->toc()
+            << exit(FatalError);
+    }
+
+    return autoPtr<filmTurbulenceModel>(cstrIter()(model, dict));
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace surfaceFilmModels
+} // End namespace regionModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/laminar/laminar.C b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/laminar/laminar.C
new file mode 100644
index 0000000000..108109702b
--- /dev/null
+++ b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/laminar/laminar.C
@@ -0,0 +1,160 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "laminar.H"
+#include "addToRunTimeSelectionTable.H"
+#include "fvMesh.H"
+#include "fvMatrices.H"
+#include "Time.H"
+#include "volFields.H"
+#include "fvmSup.H"
+#include "kinematicSingleLayer.H"
+#include "zeroGradientFvPatchFields.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace regionModels
+{
+namespace surfaceFilmModels
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(laminar, 0);
+addToRunTimeSelectionTable(filmTurbulenceModel, laminar, dictionary);
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+laminar::laminar
+(
+    const surfaceFilmModel& owner,
+    const dictionary& dict
+)
+:
+    filmTurbulenceModel(type(), owner, dict),
+    Cf_(readScalar(coeffs_.lookup("Cf")))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+laminar::~laminar()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+tmp<volVectorField> laminar::Us() const
+{
+    tmp<volVectorField> tUs
+    (
+        new volVectorField
+        (
+            IOobject
+            (
+                typeName + ".Us",
+                owner_.regionMesh().time().timeName(),
+                owner_.regionMesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            owner_.regionMesh(),
+            dimensionedVector("zero", dimVelocity, vector::zero),
+            zeroGradientFvPatchVectorField::typeName
+        )
+    );
+
+    // apply quadratic profile
+    tUs() = Foam::sqrt(2.0)*owner_.U();
+    tUs().correctBoundaryConditions();
+
+    return tUs;
+}
+
+
+tmp<volScalarField> laminar::mut() const
+{
+    return tmp<volScalarField>
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                typeName + ".mut",
+                owner_.regionMesh().time().timeName(),
+                owner_.regionMesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            owner_.regionMesh(),
+            dimensionedScalar("zero", dimMass/dimLength/dimTime, 0.0)
+        )
+    );
+}
+
+
+void laminar::correct()
+{
+    // do nothing
+}
+
+tmp<fvVectorMatrix> laminar::Su(volVectorField& U) const
+{
+    // local reference to film model
+    const kinematicSingleLayer& film =
+        static_cast<const kinematicSingleLayer&>(owner_);
+
+    // local references to film fields
+    const volScalarField& mu = film.mu();
+    const volVectorField& Uw = film.Uw();
+    const volVectorField& Us = film.Us();
+    const volScalarField& delta = film.delta();
+    const volVectorField& Up = film.UPrimary();
+    const volScalarField& rhop = film.rhoPrimary();
+
+    // employ simple coeff-based model
+    volScalarField Cs("Cs", Cf_*rhop*mag(Up - Us));
+
+    dimensionedScalar d0("SMALL", delta.dimensions(), SMALL);
+    volScalarField Cw("Cw", mu/(0.3333*(delta + d0)));
+    Cw.min(5000.0);
+
+    return
+    (
+       - fvm::Sp(Cs, U) + Cs*Us // surface contribution
+       - fvm::Sp(Cw, U) + Cw*Uw // wall contribution
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace surfaceFilmModels
+} // End namespace regionModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/laminar/laminar.H b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/laminar/laminar.H
new file mode 100644
index 0000000000..397030efff
--- /dev/null
+++ b/src/regionModels/surfaceFilmModels/submodels/kinematic/filmTurbulenceModel/laminar/laminar.H
@@ -0,0 +1,118 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::regionModels::surfaceFilmModels::laminar
+
+Description
+    Film laminar turbulence model.
+
+SourceFiles
+    laminar.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef laminar_H
+#define laminar_H
+
+#include "filmTurbulenceModel.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace regionModels
+{
+namespace surfaceFilmModels
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class laminar Declaration
+\*---------------------------------------------------------------------------*/
+
+class laminar
+:
+    public filmTurbulenceModel
+{
+private:
+
+    // Private Data
+
+        //- Surface roughness coefficient
+        scalar Cf_;
+
+
+    // Private member functions
+
+        //- Disallow default bitwise copy construct
+        laminar(const laminar&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const laminar&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("laminar");
+
+
+    // Constructors
+
+        //- Construct from surface film model
+        laminar(const surfaceFilmModel& owner, const dictionary& dict);
+
+
+    //- Destructor
+    virtual ~laminar();
+
+
+    // Member Functions
+
+        // Evolution
+
+            //- Return the film surface velocity
+            virtual tmp<volVectorField> Us() const;
+
+            //- Return the film turbulence viscosity
+            virtual tmp<volScalarField> mut() const;
+
+            //- Correct/update the model
+            virtual void correct();
+
+            //- Return the source for the film momentum equation
+            virtual tmp<fvVectorMatrix> Su(volVectorField& U) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace surfaceFilmModels
+} // End namespace regionModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.C b/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.C
index 6994f926e3..324d857a31 100644
--- a/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.C
+++ b/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -63,20 +63,51 @@ tmp<volScalarField> SpalartAllmaras::fv1() const
 
 tmp<volScalarField> SpalartAllmaras::fv2() const
 {
-    return 1.0/pow3(scalar(1) + rho()*nuTilda_/(mu()*Cv2_));
+    if (ashfordCorrection_)
+    {
+        return 1.0/pow3(scalar(1) + rho()*nuTilda_/(mu()*Cv2_));
+    }
+    else
+    {
+        const volScalarField chi("chi", rho()*nuTilda_/mu());
+        return 1.0 - chi/(1.0 + chi*fv1());
+    }
 }
 
 
 tmp<volScalarField> SpalartAllmaras::fv3() const
 {
-    volScalarField chi(rho()*nuTilda_/mu());
-    volScalarField chiByCv2((1/Cv2_)*chi);
+    if (ashfordCorrection_)
+    {
+        volScalarField chi(rho()*nuTilda_/mu());
+        volScalarField chiByCv2((1/Cv2_)*chi);
 
-    return
-        (scalar(1) + chi*fv1())
-       *(1/Cv2_)
-       *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
-       /pow3(scalar(1) + chiByCv2);
+        return
+            (scalar(1) + chi*fv1())
+           *(1/Cv2_)
+           *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
+           /pow3(scalar(1) + chiByCv2);
+    }
+    else
+    {
+        return tmp<volScalarField>
+        (
+            new volScalarField
+            (
+                IOobject
+                (
+                    "fv3",
+                    mesh_.time().timeName(),
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                ),
+                mesh_,
+                dimensionedScalar("fv3", dimless, 1),
+                zeroGradientFvPatchScalarField::typeName
+            )
+        );
+    }
 }
 
 
@@ -222,6 +253,8 @@ SpalartAllmaras::SpalartAllmaras
         )
     ),
 
+    ashfordCorrection_(coeffDict_.lookupOrDefault("ashfordCorrection", true)),
+
     nuTilda_
     (
         IOobject
@@ -265,6 +298,11 @@ SpalartAllmaras::SpalartAllmaras
     updateSubGridScaleFields();
 
     printCoeffs();
+
+    if (ashfordCorrection_)
+    {
+        Info<< "    Employing Ashford correction" << endl;
+    }
 }
 
 
@@ -344,16 +382,19 @@ bool SpalartAllmaras::read()
     {
         sigmaNut_.readIfPresent(coeffDict());
         Prt_.readIfPresent(coeffDict());
+
         Cb1_.readIfPresent(coeffDict());
         Cb2_.readIfPresent(coeffDict());
-        Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
-        Cw2_.readIfPresent(coeffDict());
-        Cw3_.readIfPresent(coeffDict());
         Cv1_.readIfPresent(coeffDict());
         Cv2_.readIfPresent(coeffDict());
         CDES_.readIfPresent(coeffDict());
         ck_.readIfPresent(coeffDict());
         kappa_.readIfPresent(*this);
+        Cw1_ = Cb1_/sqr(kappa_) + (1.0 + Cb2_)/sigmaNut_;
+        Cw2_.readIfPresent(coeffDict());
+        Cw3_.readIfPresent(coeffDict());
+
+        ashfordCorrection_.readIfPresent("ashfordCorrection", coeffDict());
 
         return true;
     }
diff --git a/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.H b/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.H
index 1897578942..11d24d3388 100644
--- a/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.H
+++ b/src/turbulenceModels/compressible/LES/SpalartAllmaras/SpalartAllmaras.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -30,6 +30,15 @@ Group
 Description
     SpalartAllmaras for compressible flows
 
+    Extended according to
+    \verbatim
+        "An Unstructured Grid Generation and Adaptive Solution Technique
+        for High Reynolds Number Compressible Flows"
+        G.A. Ashford,
+        Ph.D. thesis, University of Michigan, 1996.
+    \endverbatim
+    by using the optional flag \c ashfordCorrection
+
 SourceFiles
     SpalartAllmaras.C
 
@@ -77,12 +86,16 @@ class SpalartAllmaras
             dimensionedScalar Cw3_;
 
 
-    // Fields
+        //- Optional flag to activate the Ashford correction
+        Switch ashfordCorrection_;
 
-        volScalarField nuTilda_;
-        volScalarField dTilda_;
-        volScalarField muSgs_;
-        volScalarField alphaSgs_;
+
+        // Fields
+
+            volScalarField nuTilda_;
+            volScalarField dTilda_;
+            volScalarField muSgs_;
+            volScalarField alphaSgs_;
 
 
     // Private Member Functions
diff --git a/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.C b/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.C
index 6c5f3aae0c..9e5cf2d11f 100644
--- a/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.C
+++ b/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -63,7 +63,14 @@ tmp<volScalarField> SpalartAllmaras::fv2
     const volScalarField& fv1
 ) const
 {
-    return 1.0/pow3(scalar(1) + chi/Cv2_);
+    if (ashfordCorrection_)
+    {
+        return 1.0/pow3(scalar(1) + chi/Cv2_);
+    }
+    else
+    {
+        return 1.0 - chi/(1.0 + chi*fv1);
+    }
 }
 
 
@@ -73,13 +80,36 @@ tmp<volScalarField> SpalartAllmaras::fv3
     const volScalarField& fv1
 ) const
 {
-    const volScalarField chiByCv2((1/Cv2_)*chi);
+    if (ashfordCorrection_)
+    {
+        const volScalarField chiByCv2((1/Cv2_)*chi);
 
-    return
-        (scalar(1) + chi*fv1)
-       *(1/Cv2_)
-       *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
-       /pow3(scalar(1) + chiByCv2);
+        return
+            (scalar(1) + chi*fv1)
+           *(1/Cv2_)
+           *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
+           /pow3(scalar(1) + chiByCv2);
+    }
+    else
+    {
+        return tmp<volScalarField>
+        (
+            new volScalarField
+            (
+                IOobject
+                (
+                    "fv3",
+                    mesh_.time().timeName(),
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                ),
+                mesh_,
+                dimensionedScalar("fv3", dimless, 1),
+                zeroGradientFvPatchScalarField::typeName
+            )
+        );
+    }
 }
 
 
@@ -252,6 +282,11 @@ SpalartAllmaras::SpalartAllmaras
     alphat_.correctBoundaryConditions();
 
     printCoeffs();
+
+    if (ashfordCorrection_)
+    {
+        Info<< "    Employing Ashford correction" << endl;
+    }
 }
 
 
@@ -372,6 +407,8 @@ bool SpalartAllmaras::read()
         Cv1_.readIfPresent(coeffDict());
         Cv2_.readIfPresent(coeffDict());
 
+        ashfordCorrection_.readIfPresent("ashfordCorrection", coeffDict());
+
         return true;
     }
     else
diff --git a/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.H b/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.H
index ea829fc88d..ce1054c482 100644
--- a/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.H
+++ b/src/turbulenceModels/compressible/RAS/SpalartAllmaras/SpalartAllmaras.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -37,14 +37,16 @@ Description
         P.R. Spalart,
         S.R. Allmaras,
         La Recherche Aerospatiale, No. 1, 1994, pp. 5-21.
+    \endverbatim
 
-        Extended according to:
-
+    Extended according to:
+    \verbatim
         "An Unstructured Grid Generation and Adaptive Solution Technique
         for High Reynolds Number Compressible Flows"
         G.A. Ashford,
         Ph.D. thesis, University of Michigan, 1996.
     \endverbatim
+    by using the optional flag \c ashfordCorrection
 
     The default model coefficients correspond to the following:
     \verbatim
@@ -110,6 +112,10 @@ protected:
             dimensionedScalar Cv2_;
 
 
+        //- Optional flag to activate the Ashford correction
+        Switch ashfordCorrection_;
+
+
         // Fields
 
             volScalarField nuTilda_;
diff --git a/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.C b/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.C
index 4bbf67702a..39afb02f30 100644
--- a/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.C
+++ b/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -60,20 +60,51 @@ tmp<volScalarField> SpalartAllmaras::fv1() const
 
 tmp<volScalarField> SpalartAllmaras::fv2() const
 {
-    return 1/pow3(scalar(1) + nuTilda_/(Cv2_*nu()));
+    if (ashfordCorrection_)
+    {
+        return 1/pow3(scalar(1) + nuTilda_/(Cv2_*nu()));
+    }
+    else
+    {
+        const volScalarField chi("chi", nuTilda_/nu());
+        return 1.0 - chi/(1.0 + chi*fv1());
+    }
 }
 
 
 tmp<volScalarField> SpalartAllmaras::fv3() const
 {
-    const volScalarField chi("chi", nuTilda_/nu());
-    const volScalarField chiByCv2(chi/Cv2_);
+    if (ashfordCorrection_)
+    {
+        const volScalarField chi("chi", nuTilda_/nu());
+        const volScalarField chiByCv2(chi/Cv2_);
 
-    return
-        (scalar(1) + chi*fv1())
-       *(1/Cv2_)
-       *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
-       /pow3(scalar(1) + chiByCv2);
+        return
+            (scalar(1) + chi*fv1())
+           *(1/Cv2_)
+           *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
+           /pow3(scalar(1) + chiByCv2);
+    }
+    else
+    {
+        return tmp<volScalarField>
+        (
+            new volScalarField
+            (
+                IOobject
+                (
+                    "fv3",
+                    mesh_.time().timeName(),
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                ),
+                mesh_,
+                dimensionedScalar("fv3", dimless, 1),
+                zeroGradientFvPatchScalarField::typeName
+            )
+        );
+    }
 }
 
 
@@ -246,6 +277,8 @@ SpalartAllmaras::SpalartAllmaras
         )
     ),
 
+    ashfordCorrection_(coeffDict_.lookupOrDefault("ashfordCorrection", true)),
+
     y_(mesh_),
 
     nuTilda_
@@ -275,6 +308,11 @@ SpalartAllmaras::SpalartAllmaras
     )
 {
     updateSubGridScaleFields();
+
+    if (ashfordCorrection_)
+    {
+        Info<< "    Employing Ashford correction" << endl;
+    }
 }
 
 
@@ -376,6 +414,7 @@ bool SpalartAllmaras::read()
     {
         sigmaNut_.readIfPresent(coeffDict());
         kappa_.readIfPresent(*this);
+
         Cb1_.readIfPresent(coeffDict());
         Cb2_.readIfPresent(coeffDict());
         Cv1_.readIfPresent(coeffDict());
@@ -386,6 +425,8 @@ bool SpalartAllmaras::read()
         Cw2_.readIfPresent(coeffDict());
         Cw3_.readIfPresent(coeffDict());
 
+        ashfordCorrection_.readIfPresent("ashfordCorrection", coeffDict());
+
         return true;
     }
     else
diff --git a/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.H b/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.H
index fe4fe280f3..3f3b39e12b 100644
--- a/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.H
+++ b/src/turbulenceModels/incompressible/LES/SpalartAllmaras/SpalartAllmaras.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -30,6 +30,15 @@ Group
 Description
     SpalartAllmaras DES (SA + LES) turbulence model for incompressible flows
 
+    Extended according to
+    \verbatim
+        "An Unstructured Grid Generation and Adaptive Solution Technique
+        for High Reynolds Number Compressible Flows"
+        G.A. Ashford,
+        Ph.D. thesis, University of Michigan, 1996.
+    \endverbatim
+    by using the optional flag \c ashfordCorrection
+
 SourceFiles
     SpalartAllmaras.C
 
@@ -90,6 +99,10 @@ protected:
             dimensionedScalar Cw3_;
 
 
+        //- Optional flag to activate the Ashford correction
+        Switch ashfordCorrection_;
+
+
         // Fields
 
             wallDist y_;
diff --git a/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.C b/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.C
index 76046a80e2..999978f3c9 100644
--- a/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.C
+++ b/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -61,7 +61,14 @@ tmp<volScalarField> SpalartAllmaras::fv2
     const volScalarField& fv1
 ) const
 {
-    return 1.0/pow3(scalar(1) + chi/Cv2_);
+    if (ashfordCorrection_)
+    {
+        return 1.0/pow3(scalar(1) + chi/Cv2_);
+    }
+    else
+    {
+        return 1.0 - chi/(1.0 + chi*fv1);
+    }
 }
 
 
@@ -71,13 +78,36 @@ tmp<volScalarField> SpalartAllmaras::fv3
     const volScalarField& fv1
 ) const
 {
-    const volScalarField chiByCv2((1/Cv2_)*chi);
+    if (ashfordCorrection_)
+    {
+        const volScalarField chiByCv2((1/Cv2_)*chi);
 
-    return
-        (scalar(1) + chi*fv1)
-       *(1/Cv2_)
-       *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
-       /pow3(scalar(1) + chiByCv2);
+        return
+            (scalar(1) + chi*fv1)
+           *(1/Cv2_)
+           *(3*(scalar(1) + chiByCv2) + sqr(chiByCv2))
+           /pow3(scalar(1) + chiByCv2);
+    }
+    else
+    {
+        return tmp<volScalarField>
+        (
+            new volScalarField
+            (
+                IOobject
+                (
+                    "fv3",
+                    mesh_.time().timeName(),
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                ),
+                mesh_,
+                dimensionedScalar("fv3", dimless, 1),
+                zeroGradientFvPatchScalarField::typeName
+            )
+        );
+    }
 }
 
 
@@ -195,6 +225,8 @@ SpalartAllmaras::SpalartAllmaras
         )
     ),
 
+    ashfordCorrection_(coeffDict_.lookupOrDefault("ashfordCorrection", true)),
+
     nuTilda_
     (
         IOobject
@@ -224,6 +256,11 @@ SpalartAllmaras::SpalartAllmaras
     d_(mesh_)
 {
     printCoeffs();
+
+    if (ashfordCorrection_)
+    {
+        Info<< "    Employing Ashford correction" << endl;
+    }
 }
 
 
@@ -368,6 +405,8 @@ bool SpalartAllmaras::read()
         Cv1_.readIfPresent(coeffDict());
         Cv2_.readIfPresent(coeffDict());
 
+        ashfordCorrection_.readIfPresent("ashfordCorrection", coeffDict());
+
         return true;
     }
     else
diff --git a/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.H b/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.H
index d7479c2dd5..597b0859c4 100644
--- a/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.H
+++ b/src/turbulenceModels/incompressible/RAS/SpalartAllmaras/SpalartAllmaras.H
@@ -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-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -37,14 +37,16 @@ Description
         P.R. Spalart,
         S.R. Allmaras,
         La Recherche Aerospatiale, No. 1, 1994, pp. 5-21.
+    \endverbatim
 
-        Extended according to:
-
+    Extended according to
+    \verbatim
         "An Unstructured Grid Generation and Adaptive Solution Technique
         for High Reynolds Number Compressible Flows"
         G.A. Ashford,
         Ph.D. thesis, University of Michigan, 1996.
     \endverbatim
+    using the optional flag \c ashfordCorrection
 
     The default model coefficients correspond to the following:
     \verbatim
@@ -82,7 +84,7 @@ namespace RASModels
 {
 
 /*---------------------------------------------------------------------------*\
-                           Class SpalartAllmaras Declaration
+                       Class SpalartAllmaras Declaration
 \*---------------------------------------------------------------------------*/
 
 class SpalartAllmaras
@@ -108,6 +110,10 @@ protected:
             dimensionedScalar Cv2_;
 
 
+        //- Optional flag to activate the Ashford correction
+        Switch ashfordCorrection_;
+
+
         // Fields
 
             volScalarField nuTilda_;