diff --git a/src/faOptions/Make/files b/src/faOptions/Make/files
index cd8c59745df5cbd23f79099fe327b62161876a73..af26eff51b25f37b90bad9d5a4d62ce25918ba80 100644
--- a/src/faOptions/Make/files
+++ b/src/faOptions/Make/files
@@ -13,4 +13,9 @@ $(derivedSources)/jouleHeatingSource/jouleHeatingSource.C
 $(derivedSources)/contactHeatFluxSource/contactHeatFluxSource.C
 $(derivedSources)/externalFileSource/externalFileSource.C
 
+/* Corrections */
+corrections=corrections
+
+$(corrections)/limitVelocity/limitVelocity.C
+
 LIB = $(FOAM_LIBBIN)/libfaOptions
diff --git a/src/faOptions/corrections/limitVelocity/limitVelocity.C b/src/faOptions/corrections/limitVelocity/limitVelocity.C
index a7440549a6096317429eda309fd391201e00723d..8fbb5756c0f4cabf3daa4344404628b821ac875c 100644
--- a/src/faOptions/corrections/limitVelocity/limitVelocity.C
+++ b/src/faOptions/corrections/limitVelocity/limitVelocity.C
@@ -96,7 +96,7 @@ void Foam::fa::limitVelocity::correct(areaVectorField& U)
     }
 
     // Handle boundaries in the case of 'all'
-    if (!cellSetOption::useSubMesh())
+    if (!faceSetOption::useSubMesh())
     {
         for (faPatchVectorField& Up : U.boundaryFieldRef())
         {
diff --git a/src/regionFaModels/liquidFilm/liquidFilmBase.C b/src/regionFaModels/liquidFilm/liquidFilmBase.C
index 5ae39554858c9408d6e33e6e98e1db7145ed9722..d5a49af1255a787bc4b1393dbd42e399e63e6df9 100644
--- a/src/regionFaModels/liquidFilm/liquidFilmBase.C
+++ b/src/regionFaModels/liquidFilm/liquidFilmBase.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -551,6 +551,13 @@ scalar liquidFilmBase::pRef()
     return pRef_;
 }
 
+
+word liquidFilmBase::UName() const
+{
+    return UName_;
+}
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace areaSurfaceFilmModels
diff --git a/src/regionFaModels/liquidFilm/liquidFilmBase.H b/src/regionFaModels/liquidFilm/liquidFilmBase.H
index 8bbce9f9b5213acb2d24a4c7b13fe89613ba4eb1..b51d8d6f457bc3ef8143d7e154e182704e72db51 100644
--- a/src/regionFaModels/liquidFilm/liquidFilmBase.H
+++ b/src/regionFaModels/liquidFilm/liquidFilmBase.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -261,6 +261,9 @@ public:
             //- Access to pRef
             scalar pRef();
 
+            //- Name of the U field
+            word UName() const;
+
 
             // Transfer fields - to the primary region (lagragian injection)
 
diff --git a/src/regionFaModels/liquidFilm/liquidFilmModel/liquidFilmModel.C b/src/regionFaModels/liquidFilm/liquidFilmModel/liquidFilmModel.C
index 15a4149e013fe98d169d17f967e905f08d685eb6..14d462c640e5444092b91190fa538beb39b73d43 100644
--- a/src/regionFaModels/liquidFilm/liquidFilmModel/liquidFilmModel.C
+++ b/src/regionFaModels/liquidFilm/liquidFilmModel/liquidFilmModel.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -180,7 +180,7 @@ liquidFilmModel::liquidFilmModel
             primaryMesh().time().timeName(),
             primaryMesh(),
             IOobject::NO_READ,
-            IOobject::AUTO_WRITE
+            IOobject::NO_WRITE
         ),
         regionMesh(),
         dimensionedScalar(dimVelocity, Zero)
@@ -215,7 +215,7 @@ liquidFilmModel::liquidFilmModel
             primaryMesh().time().timeName(),
             primaryMesh(),
             IOobject::NO_READ,
-            IOobject::AUTO_WRITE
+            IOobject::NO_WRITE
         ),
         primaryMesh(),
         dimensionedScalar(dimMass, Zero),
@@ -230,7 +230,7 @@ liquidFilmModel::liquidFilmModel
             primaryMesh().time().timeName(),
             primaryMesh(),
             IOobject::NO_READ,
-            IOobject::AUTO_WRITE
+            IOobject::NO_WRITE
         ),
         primaryMesh(),
         dimensionedScalar(dimLength, Zero),
diff --git a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C
index cec1b7acb9a4301e045f25185efe9cf63f1348be..44a078987865b5b34f3962762802d89113441442 100644
--- a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C
+++ b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -27,6 +27,8 @@ License
 
 #include "filmTurbulenceModel.H"
 #include "gravityMeshObject.H"
+#include "turbulentTransportModel.H"
+#include "turbulentFluidThermoModel.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -54,6 +56,18 @@ filmTurbulenceModel::frictionMethodTypeNames_
     { frictionMethodType::mManningStrickler, "ManningStrickler" }
 };
 
+
+const Enum
+<
+    filmTurbulenceModel::shearMethodType
+>
+filmTurbulenceModel::shearMethodTypeNames_
+{
+    { shearMethodType::msimple, "simple" },
+    { shearMethodType::mwallFunction, "wallFunction" }
+};
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 filmTurbulenceModel::filmTurbulenceModel
@@ -65,14 +79,16 @@ filmTurbulenceModel::filmTurbulenceModel
 :
     film_(film),
     dict_(dict.subDict(modelType + "Coeffs")),
-    method_(frictionMethodTypeNames_.get("friction", dict_))
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-filmTurbulenceModel::~filmTurbulenceModel()
-{}
+    method_(frictionMethodTypeNames_.get("friction", dict_)),
+    shearMethod_(shearMethodTypeNames_.get("shearStress", dict_)),
+    rhoName_(dict_.getOrDefault<word>("rho", "rho")),
+    rhoRef_(VGREAT)
+{
+    if (rhoName_ == "rhoInf")
+    {
+        rhoRef_ = dict_.get<scalar>("rhoInf");
+    }
+}
 
 
 // * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
@@ -99,34 +115,45 @@ tmp<areaScalarField> filmTurbulenceModel::Cw() const
         );
     auto& Cw = tCw.ref();
 
+
     switch (method_)
     {
         case mquadraticProfile:
         {
             const scalarField& mu = film_.mu().primitiveField();
             const scalarField& h = film_.h().primitiveField();
+            const scalarField& rho = film_.rho().primitiveField();
+
             const scalar h0 = film_.h0().value();
-            Cw.primitiveFieldRef() = 3*mu/(h + h0);
+
+            Cw.primitiveFieldRef() = 3*mu/((h + h0)*rho);
             Cw.min(5000.0);
+
             break;
         }
         case mlinearProfile:
         {
             const scalarField& mu = film_.mu().primitiveField();
             const scalarField& h = film_.h().primitiveField();
+            const scalarField& rho = film_.rho().primitiveField();
+
             const scalar h0 = film_.h0().value();
-            Cw.primitiveFieldRef() = 2*mu/(h + h0);
+
+            Cw.primitiveFieldRef() = 2*mu/((h + h0)*rho);
+
             break;
         }
         case mDarcyWeisbach:
         {
             const uniformDimensionedVectorField& g =
                 meshObjects::gravity::New(film_.primaryMesh().time());
-
             const vectorField& Uf = film_.Uf().primitiveField();
+            const scalarField& rho = film_.rho().primitiveField();
+
+            const scalar Cf = dict_.get<scalar>("DarcyWeisbach");
+
+            Cw.primitiveFieldRef() = Cf*mag(g.value())*mag(Uf)/rho;
 
-            scalar Cf = dict_.get<scalar>("Cf");
-            Cw.primitiveFieldRef() = Cf*mag(g.value())*mag(Uf);
             break;
         }
         case mManningStrickler:
@@ -134,14 +161,15 @@ tmp<areaScalarField> filmTurbulenceModel::Cw() const
             const uniformDimensionedVectorField& g =
                 meshObjects::gravity::New(film_.primaryMesh().time());
 
-            const scalar n = dict_.get<scalar>("n");
-
             const vectorField& Uf = film_.Uf().primitiveField();
             const scalarField& h = film_.h().primitiveField();
             const scalar h0 = film_.h0().value();
 
+            const scalar n = dict_.get<scalar>("n");
+
             Cw.primitiveFieldRef() =
-                sqr(n)*mag(g.value())*mag(Uf)/cbrt(h+h0);
+                sqr(n)*mag(g.value())*mag(Uf)/cbrt(h + h0);
+
             break;
         }
         default:
@@ -161,6 +189,167 @@ tmp<areaScalarField> filmTurbulenceModel::Cw() const
 }
 
 
+tmp<faVectorMatrix> filmTurbulenceModel::primaryRegionFriction
+(
+    areaVectorField& U
+) const
+{
+    tmp<faVectorMatrix> tshearStress
+    (
+        new faVectorMatrix(U, sqr(U.dimensions())*sqr(dimLength))
+    );
+
+    switch (shearMethod_)
+    {
+        case msimple:
+        {
+            tmp<areaVectorField> Up = film_.Up();
+
+            const dimensionedScalar Cf
+            (
+                "Cf",
+                dimVelocity,
+                dict_.get<scalar>("Cf")
+            );
+
+            tshearStress.ref() += - fam::Sp(Cf, U) + Cf*Up();
+
+            break;
+        }
+        case mwallFunction:
+        {
+            tmp<volSymmTensorField> tdevRhoReff = devRhoReff();
+
+            const volSymmTensorField::Boundary& devRhoReffb
+                = tdevRhoReff().boundaryField();
+
+            const label patchi = film_.patchID();
+
+            const surfaceVectorField::Boundary& Sfb =
+                film_.primaryMesh().Sf().boundaryField();
+
+            vectorField fT(Sfb[patchi] & devRhoReffb[patchi]);
+
+            const vectorField& nHat =
+                film_.regionMesh().faceAreaNormals().internalField();
+
+            // Substract normal component
+            fT -= nHat*(fT & nHat);
+
+            auto taForce = tmp<areaVectorField>::New
+            (
+                IOobject
+                (
+                    "taForce",
+                    film_.primaryMesh().time().timeName(),
+                    film_.primaryMesh()
+                ),
+                film_.regionMesh(),
+                dimensionedVector(sqr(dimVelocity), Zero)
+            );
+            vectorField& aForce = taForce.ref().primitiveFieldRef();
+
+            // Map ft to surface
+            const vectorField afT(film_.vsm().mapToSurface(fT));
+
+            const DimensionedField<scalar, areaMesh>& magSf =
+                film_.regionMesh().S();
+
+            aForce = afT/(film_.rho().primitiveField()*magSf);
+
+            tshearStress.ref() += taForce();
+
+            if (film_.regionMesh().time().writeTime())
+            {
+                taForce().write();
+            }
+
+            break;
+        }
+    }
+
+    return tshearStress;
+}
+
+
+tmp<Foam::volSymmTensorField> filmTurbulenceModel::devRhoReff() const
+{
+    typedef compressible::turbulenceModel cmpTurbModel;
+    typedef incompressible::turbulenceModel icoTurbModel;
+
+    const fvMesh& m = film_.primaryMesh();
+
+    const auto& U = m.lookupObject<volVectorField>(film_.UName());
+
+    if (m.foundObject<cmpTurbModel>(cmpTurbModel::propertiesName))
+    {
+        const auto& turb =
+            m.lookupObject<cmpTurbModel>(cmpTurbModel::propertiesName);
+
+        return turb.devRhoReff();
+    }
+    else if (m.foundObject<icoTurbModel>(icoTurbModel::propertiesName))
+    {
+        const auto& turb =
+            m.lookupObject<icoTurbModel>(icoTurbModel::propertiesName);
+
+        return rho()*turb.devReff();
+    }
+    else if (m.foundObject<fluidThermo>(fluidThermo::dictName))
+    {
+        const auto& thermo =
+            m.lookupObject<fluidThermo>(fluidThermo::dictName);
+
+        return -thermo.mu()*dev(twoSymm(fvc::grad(U)));
+    }
+    else if (m.foundObject<transportModel>("transportProperties"))
+    {
+        const auto& laminarT =
+            m.lookupObject<transportModel>("transportProperties");
+
+        return -rho()*laminarT.nu()*dev(twoSymm(fvc::grad(U)));
+    }
+    else if (m.foundObject<dictionary>("transportProperties"))
+    {
+        const auto& transportProperties =
+            m.lookupObject<dictionary>("transportProperties");
+
+        const dimensionedScalar nu("nu", dimViscosity, transportProperties);
+
+        return -rho()*nu*dev(twoSymm(fvc::grad(U)));
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "No valid model for viscous stress calculation"
+            << exit(FatalError);
+
+        return volSymmTensorField::null();
+    }
+}
+
+
+tmp<Foam::volScalarField> filmTurbulenceModel::rho() const
+{
+    const fvMesh& m = film_.primaryMesh();
+    if (rhoName_ == "rhoInf")
+    {
+        return tmp<volScalarField>::New
+        (
+            IOobject
+            (
+                "rho",
+                m.time().timeName(),
+                m
+            ),
+            m,
+            dimensionedScalar(dimDensity, rhoRef_)
+        );
+    }
+
+    return m.lookupObject<volScalarField>(rhoName_);
+}
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace areaSurfaceFilmModels
diff --git a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H
index 5277983eecc508f2507272905a8c714c7cd8020e..69ce7cea2674e5a96fd2231ece9c5d345c44081f 100644
--- a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H
+++ b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/filmTurbulenceModel/filmTurbulenceModel.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -82,6 +82,13 @@ public:
             mManningStrickler
         };
 
+        //- Options for the shear stress models
+        enum shearMethodType
+        {
+            msimple,
+            mwallFunction
+        };
+
 
 protected:
 
@@ -93,12 +100,24 @@ protected:
         //- Names for friction models
         static const Enum<frictionMethodType> frictionMethodTypeNames_;
 
+        //- Names for shear stress models
+        static const Enum<shearMethodType> shearMethodTypeNames_;
+
         //- Model dictionary
         const dictionary dict_;
 
         //- Method used
         const frictionMethodType method_;
 
+        //- Shear method used
+        const shearMethodType shearMethod_;
+
+        //- Name of density field (optional)
+        word rhoName_;
+
+        //- Reference density needed for incompressible calculations
+        scalar rhoRef_;
+
 
 public:
 
@@ -143,7 +162,7 @@ public:
 
 
     //- Destructor
-    virtual ~filmTurbulenceModel();
+    virtual ~filmTurbulenceModel() = default;
 
 
     // Member Functions
@@ -154,7 +173,19 @@ public:
             const liquidFilmBase& film() const;
 
 
-        // Evolution
+        // Turbulence
+
+            //- Return the effective viscous stress (laminar + turbulent)
+            tmp<volSymmTensorField> devRhoReff() const;
+
+            //- Return primary region friction
+            tmp<faVectorMatrix> primaryRegionFriction
+            (
+                areaVectorField& U
+            ) const;
+
+            //- Return rho if specified otherwise rhoRef
+            tmp<volScalarField> rho() const;
 
             //- Return the wall film surface friction
             virtual tmp<areaScalarField> Cw() const;
@@ -162,6 +193,9 @@ public:
             //- Return the film turbulence viscosity
             virtual tmp<areaScalarField> mut() const = 0;
 
+
+        // Evaluation
+
             //- Correct/update the model
             virtual void correct() = 0;
 
diff --git a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.C b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.C
index 7630e254661109ce28688d7d664dd581ad1a3489..7859246570decd4140b48369f740a9c7e5d1ca34 100644
--- a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.C
+++ b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.C
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -51,14 +51,7 @@ laminar::laminar
     const dictionary& dict
 )
 :
-    filmTurbulenceModel(type(), film, dict),
-    Cf_(dict_.get<scalar>("Cf"))
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-laminar::~laminar()
+    filmTurbulenceModel(type(), film, dict)
 {}
 
 
@@ -66,20 +59,19 @@ laminar::~laminar()
 
 tmp<areaScalarField> laminar::mut() const
 {
-    auto tmut =
-        tmp<areaScalarField>::New
+    auto tmut = tmp<areaScalarField>::New
+    (
+        IOobject
         (
-            IOobject
-            (
-                "mut",
-                film().primaryMesh().time().timeName(),
-                film().primaryMesh(),
-                IOobject::NO_READ,
-                IOobject::NO_WRITE
-            ),
-            film().regionMesh(),
-            dimensionedScalar(dimMass/dimLength/dimTime)
-        );
+            "mut",
+            film().primaryMesh().time().timeName(),
+            film().primaryMesh(),
+            IOobject::NO_READ,
+            IOobject::NO_WRITE
+        ),
+        film().regionMesh(),
+        dimensionedScalar(dimMass/dimLength/dimTime)
+    );
 
     return tmut;
 }
@@ -91,18 +83,18 @@ void laminar::correct()
 
 tmp<faVectorMatrix> laminar::Su(areaVectorField& U) const
 {
-    // local references to film fields
-    tmp<areaVectorField> Uw = film_.Uw();
-    tmp<areaVectorField> Up = film_.Up();
+    return primaryRegionFriction(U) + wallFriction(U);
+}
 
-    // employ simple coeff-based model
-    const dimensionedScalar Cf("Cf", dimVelocity, Cf_);
 
+tmp<faVectorMatrix> laminar::wallFriction(areaVectorField& U) const
+{
+    // local references to film fields
+    tmp<areaVectorField> Uw = film_.Uw();
     tmp<areaScalarField> wf = Cw();
 
     return
     (
-       - fam::Sp(Cf, U) + Cf*Up()     // surface contribution
        - fam::Sp(wf(), U) + wf()*Uw() // wall contribution
     );
 }
diff --git a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.H b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.H
index b0e6b7054136866f3adf18ce0953ac63d5f09663..64488aa23a1b64089a1c2a6f65cccfc770131106 100644
--- a/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.H
+++ b/src/regionFaModels/liquidFilm/subModels/kinematic/filmTurbulenceModel/laminar/laminar.H
@@ -5,7 +5,7 @@
     \\  /    A nd           | www.openfoam.com
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
-    Copyright (C) 2020 OpenCFD Ltd.
+    Copyright (C) 2020-2021 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -57,12 +57,6 @@ class laminar
 :
     public filmTurbulenceModel
 {
-    // Private Data
-
-        //- Surface roughness coefficient
-        scalar Cf_;
-
-
     // Private Member Functions
 
         //- No copy construct
@@ -85,21 +79,27 @@ public:
 
 
     //- Destructor
-    virtual ~laminar();
+    virtual ~laminar() = default;
 
 
     // Member Functions
 
-        // Evolution
+        // Turbulence
+
+            //- Wall friction
+            tmp<faVectorMatrix> wallFriction(areaVectorField& U) const;
 
             //- Return the film turbulence viscosity
             virtual tmp<areaScalarField> mut() const;
 
-            //- Correct/update the model
-            virtual void correct();
-
             //- Return the source for the film momentum equation
             virtual tmp<faVectorMatrix> Su(areaVectorField& U) const;
+
+
+        // Evaluation
+
+            //- Correct/update the model
+            virtual void correct();
 };
 
 
diff --git a/tutorials/incompressible/pimpleFoam/laminar/inclinedPlaneFilm/0/U b/tutorials/incompressible/pimpleFoam/laminar/inclinedPlaneFilm/0/U
index a40b782ffe7c02258a403134f9f56346c91cafff..b4f2534dd51900be005f74031dd464eb7e7a3190 100644
--- a/tutorials/incompressible/pimpleFoam/laminar/inclinedPlaneFilm/0/U
+++ b/tutorials/incompressible/pimpleFoam/laminar/inclinedPlaneFilm/0/U
@@ -51,6 +51,7 @@ boundaryField
 
         laminarCoeffs
         {
+            shearStress     simple;
             friction        ManningStrickler;
             n               0.1; // Manning number
             Cf              0.9; // Gas friction
diff --git a/tutorials/lagrangian/kinematicParcelFoam/pitzDailyWithSprinklers/0/U b/tutorials/lagrangian/kinematicParcelFoam/pitzDailyWithSprinklers/0/U
index fc2264f16cbde48e2f1564f3cb7910cc2798bbba..e250ac7e0552b4307ec9a419b1a9901f9479c755 100644
--- a/tutorials/lagrangian/kinematicParcelFoam/pitzDailyWithSprinklers/0/U
+++ b/tutorials/lagrangian/kinematicParcelFoam/pitzDailyWithSprinklers/0/U
@@ -57,6 +57,7 @@ boundaryField
 
         laminarCoeffs
         {
+            shearStress     simple;
             friction        ManningStrickler; // Wall friction model
             n               0.005;  // Manning number
             Cf              0;      // Gas friction