From e0912a058cb68b0538ef6af6342a98ab1a331b81 Mon Sep 17 00:00:00 2001
From: sergio <sergio>
Date: Fri, 18 May 2018 01:27:13 -0700
Subject: [PATCH 1/8] ENH: Adding pressure work terms for compressibleInterFoam
 Removing temporary fix from compressibleInterDyFOAM to re-calculte Uf with
 morphing mesh

---
 .../solvers/multiphase/compressibleInterFoam/TEqn.H   | 11 +++++------
 .../compressibleInterDyMFoam.C                        |  6 +-----
 2 files changed, 6 insertions(+), 11 deletions(-)

diff --git a/applications/solvers/multiphase/compressibleInterFoam/TEqn.H b/applications/solvers/multiphase/compressibleInterFoam/TEqn.H
index 26f046545e..1561ae35cc 100644
--- a/applications/solvers/multiphase/compressibleInterFoam/TEqn.H
+++ b/applications/solvers/multiphase/compressibleInterFoam/TEqn.H
@@ -1,16 +1,15 @@
 {
     fvScalarMatrix TEqn
     (
-        fvm::ddt(rho, T) + fvm::div(rhoPhi, T)
-      - fvm::Sp(contErr, T)
+        fvm::ddt(rho, T) + fvm::div(rhoPhi, T) - fvm::Sp(contErr, T)
       - fvm::laplacian(turbulence.alphaEff(), T)
       + (
-            divU*p
-          + fvc::ddt(rho, K) + fvc::div(rhoPhi, K)
+             (divU*p)() - contErr/rho*p
+          +  (fvc::ddt(rho, K) + fvc::div(rhoPhi, K))() - contErr*K
         )
        *(
-           alpha1/mixture.thermo1().Cv()
-         + alpha2/mixture.thermo2().Cv()
+           alpha1()/mixture.thermo1().Cv()()
+         + alpha2()/mixture.thermo2().Cv()()
         )
      ==
         fvOptions(rho, T)
diff --git a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C
index 8d3d5a2702..816c495443 100644
--- a/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C
+++ b/applications/solvers/multiphase/compressibleInterFoam/compressibleInterDyMFoam/compressibleInterDyMFoam.C
@@ -119,12 +119,8 @@ int main(int argc, char *argv[])
                     ghf = (g & mesh.Cf()) - ghRef;
                 }
 
-                if ((mesh.changing() && correctPhi) || mesh.topoChanging())
+                if ((mesh.changing() && correctPhi))
                 {
-                    // Calculate absolute flux from the mapped surface velocity
-                    // Note: temporary fix until mapped Uf is assessed
-                    Uf = fvc::interpolate(U);
-
                     // Calculate absolute flux from the mapped surface velocity
                     phi = mesh.Sf() & Uf;
 

From f246a81f8af2a897e1fcd1db1ac5c3f772a9c4d4 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://openfoam.org>
Date: Wed, 21 Mar 2018 18:03:29 +0000
Subject: [PATCH 2/8] BUG: wideBandAbsorptionEmission: Corrected errors

Resolves bug-reports
https://bugs.openfoam.org/view.php?id=2881
https://bugs.openfoam.org/view.php?id=2882

Patches contributed by Kevin Nordin-Bates
---
 .../fvDOM/absorptionCoeffs/absorptionCoeffs.C |   6 +-
 .../blackBodyEmission/blackBodyEmission.C     |  56 +++++--
 .../blackBodyEmission/blackBodyEmission.H     |  12 +-
 .../radiationModels/fvDOM/fvDOM/fvDOM.C       |  97 +++++++++---
 .../wideBandAbsorptionEmission.C              | 140 ++++++++++++------
 .../wideBandAbsorptionEmission.H              |  41 ++---
 6 files changed, 240 insertions(+), 112 deletions(-)

diff --git a/src/thermophysicalModels/radiation/radiationModels/fvDOM/absorptionCoeffs/absorptionCoeffs.C b/src/thermophysicalModels/radiation/radiationModels/fvDOM/absorptionCoeffs/absorptionCoeffs.C
index 1b50a08cfd..3edf6b2c7a 100644
--- a/src/thermophysicalModels/radiation/radiationModels/fvDOM/absorptionCoeffs/absorptionCoeffs.C
+++ b/src/thermophysicalModels/radiation/radiationModels/fvDOM/absorptionCoeffs/absorptionCoeffs.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -26,7 +26,6 @@ License
 #include "absorptionCoeffs.H"
 #include "IOstreams.H"
 
-
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * //
 
 Foam::radiation::absorptionCoeffs::~absorptionCoeffs()
@@ -40,7 +39,7 @@ void Foam::radiation::absorptionCoeffs::checkT(const scalar T) const
     if (T < Tlow_ || T > Thigh_)
     {
         WarningInFunction
-            << "usinf absCoeff out of temperature range:" << nl
+            << "using absorptionCoeffs out of temperature range:" << nl
             << "    " << Tlow_ << " -> " << Thigh_ << ";  T = " << T
             << nl << endl;
     }
@@ -72,7 +71,6 @@ void Foam::radiation::absorptionCoeffs::initialise(const dictionary& dict)
     dict.lookup("Tlow") >> Tlow_;
     dict.lookup("Thigh") >> Thigh_;
     dict.lookup("invTemp") >> invTemp_;
-
     dict.lookup("loTcoeffs") >> lowACoeffs_;
     dict.lookup("hiTcoeffs") >> highACoeffs_;
 }
diff --git a/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.C b/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.C
index dcc1029a8a..57e4f0a16d 100644
--- a/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.C
+++ b/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -192,6 +192,44 @@ Foam::scalar Foam::radiation::blackBodyEmission::fLambdaT
 }
 
 
+Foam::tmp<Foam::volScalarField>
+Foam::radiation::blackBodyEmission::deltaLambdaT
+(
+    const volScalarField& T,
+    const Vector2D<scalar>& band
+) const
+{
+    tmp<volScalarField> deltaLambdaT
+    (
+        new volScalarField
+        (
+            IOobject
+            (
+                "deltaLambdaT",
+                T.mesh().time().timeName(),
+                T.mesh(),
+                IOobject::NO_READ,
+                IOobject::NO_WRITE
+            ),
+            T.mesh(),
+            dimensionedScalar("deltaLambdaT", dimless, 1.0)
+        )
+    );
+
+    if (band != Vector2D<scalar>::one)
+    {
+        scalarField& deltaLambdaTf = deltaLambdaT.ref();
+
+        forAll(T, i)
+        {
+            deltaLambdaTf[i] = fLambdaT(band[1]*T[i]) - fLambdaT(band[0]*T[i]);
+        }
+    }
+
+    return deltaLambdaT;
+}
+
+
 Foam::tmp<Foam::volScalarField>
 Foam::radiation::blackBodyEmission::EbDeltaLambdaT
 (
@@ -215,21 +253,13 @@ Foam::radiation::blackBodyEmission::EbDeltaLambdaT
         )
     );
 
-
-    if (band == Vector2D<scalar>::one)
-    {
-        return Eb;
-    }
-    else
+    if (band != Vector2D<scalar>::one)
     {
         scalarField& Ebif = Eb.ref();
 
         forAll(T, i)
         {
-            const scalar T1 = fLambdaT(band[1]*T[i]);
-            const scalar T2 = fLambdaT(band[0]*T[i]);
-
-            Ebif[i] *= T1 - T2;
+            Ebif[i] *= fLambdaT(band[1]*T[i]) - fLambdaT(band[0]*T[i]);
         }
 
         volScalarField::Boundary& EbBf = Eb.ref().boundaryFieldRef();
@@ -251,9 +281,9 @@ Foam::radiation::blackBodyEmission::EbDeltaLambdaT
                 }
             }
         }
-
-        return Eb;
     }
+
+    return Eb;
 }
 
 
diff --git a/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.H b/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.H
index c22319fc24..56fa3c441b 100644
--- a/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.H
+++ b/src/thermophysicalModels/radiation/radiationModels/fvDOM/blackBodyEmission/blackBodyEmission.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -27,7 +27,7 @@ Class
 Description
     Class black body emission
 
-    Table of black body emissive power taken from:
+    Table of black body emissive power from:
         Modest, "Radiative Heat Transfer", pp.775-777, 1993
 
 SourceFiles
@@ -123,6 +123,13 @@ public:
                 return (C1_/(pow5(lambda)*(exp(C2_/(lambda*T)) - 1.0)));
             }
 
+            //- Proportion of total energy at T from lambda1 to lambda2
+            tmp<Foam::volScalarField> deltaLambdaT
+            (
+                const volScalarField& T,
+                const Vector2D<scalar>& band
+            ) const;
+
             //- Integral energy at T from lambda1 to lambda2
             tmp<Foam::volScalarField> EbDeltaLambdaT
             (
@@ -130,7 +137,6 @@ public:
                 const Vector2D<scalar>& band
             ) const;
 
-
     // Edit
 
         // Update black body emission
diff --git a/src/thermophysicalModels/radiation/radiationModels/fvDOM/fvDOM/fvDOM.C b/src/thermophysicalModels/radiation/radiationModels/fvDOM/fvDOM/fvDOM.C
index f6506064c1..78402206e4 100644
--- a/src/thermophysicalModels/radiation/radiationModels/fvDOM/fvDOM/fvDOM.C
+++ b/src/thermophysicalModels/radiation/radiationModels/fvDOM/fvDOM/fvDOM.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -54,15 +54,15 @@ void Foam::radiation::fvDOM::initialise()
     {
         nRay_ = 4*nPhi_*nTheta_;
         IRay_.setSize(nRay_);
-        scalar deltaPhi = pi/(2.0*nPhi_);
-        scalar deltaTheta = pi/nTheta_;
+        const scalar deltaPhi = pi/(2.0*nPhi_);
+        const scalar deltaTheta = pi/nTheta_;
         label i = 0;
         for (label n = 1; n <= nTheta_; n++)
         {
             for (label m = 1; m <= 4*nPhi_; m++)
             {
-                scalar thetai = (2.0*n - 1.0)*deltaTheta/2.0;
-                scalar phii = (2.0*m - 1.0)*deltaPhi/2.0;
+                const scalar thetai = (2*n - 1)*deltaTheta/2.0;
+                const scalar phii = (2*m - 1)*deltaPhi/2.0;
                 IRay_.set
                 (
                     i,
@@ -87,15 +87,15 @@ void Foam::radiation::fvDOM::initialise()
     // 2D
     else if (mesh_.nSolutionD() == 2)
     {
-        scalar thetai = piByTwo;
-        scalar deltaTheta = pi;
+        const scalar thetai = piByTwo;
+        const scalar deltaTheta = pi;
         nRay_ = 4*nPhi_;
         IRay_.setSize(nRay_);
-        scalar deltaPhi = pi/(2.0*nPhi_);
+        const scalar deltaPhi = pi/(2.0*nPhi_);
         label i = 0;
         for (label m = 1; m <= 4*nPhi_; m++)
         {
-            scalar phii = (2.0*m - 1.0)*deltaPhi/2.0;
+            const scalar phii = (2*m - 1)*deltaPhi/2.0;
             IRay_.set
             (
                 i,
@@ -119,15 +119,15 @@ void Foam::radiation::fvDOM::initialise()
     // 1D
     else
     {
-        scalar thetai = piByTwo;
-        scalar deltaTheta = pi;
+        const scalar thetai = piByTwo;
+        const scalar deltaTheta = pi;
         nRay_ = 2;
         IRay_.setSize(nRay_);
-        scalar deltaPhi = pi;
+        const scalar deltaPhi = pi;
         label i = 0;
         for (label m = 1; m <= 2; m++)
         {
-            scalar phii = (2.0*m - 1.0)*deltaPhi/2.0;
+            const scalar phii = (2*m - 1)*deltaPhi/2.0;
             IRay_.set
             (
                 i,
@@ -444,14 +444,14 @@ void Foam::radiation::fvDOM::calculate()
     // Set rays convergence false
     List<bool> rayIdConv(nRay_, false);
 
-    scalar maxResidual = 0.0;
+    scalar maxResidual = 0;
     label radIter = 0;
     do
     {
         Info<< "Radiation solver iter: " << radIter << endl;
 
         radIter++;
-        maxResidual = 0.0;
+        maxResidual = 0;
         forAll(IRay_, rayI)
         {
             if (!rayIdConv[rayI])
@@ -474,7 +474,8 @@ void Foam::radiation::fvDOM::calculate()
 
 Foam::tmp<Foam::volScalarField> Foam::radiation::fvDOM::Rp() const
 {
-    return tmp<volScalarField>
+    // Construct using contribution from first frequency band
+    tmp<volScalarField> tRp
     (
         new volScalarField
         (
@@ -487,21 +488,75 @@ Foam::tmp<Foam::volScalarField> Foam::radiation::fvDOM::Rp() const
                 IOobject::NO_WRITE,
                 false
             ),
-            4.0*a_*physicoChemical::sigma //absorptionEmission_->a()
+            (
+                4
+               *physicoChemical::sigma
+               *(aLambda_[0] - absorptionEmission_->aDisp(0)())
+               *blackBody_.deltaLambdaT(T_, absorptionEmission_->bands(0))
+            )
         )
     );
+
+    volScalarField& Rp=tRp.ref();
+
+    // Add contributions over remaining frequency bands
+    for (label j=1; j < nLambda_; j++)
+    {
+        Rp +=
+        (
+            4
+           *physicoChemical::sigma
+           *(aLambda_[j] - absorptionEmission_->aDisp(j)())
+           *blackBody_.deltaLambdaT(T_, absorptionEmission_->bands(j))
+        );
+    }
+
+    return tRp;
 }
 
 
 Foam::tmp<Foam::DimensionedField<Foam::scalar, Foam::volMesh>>
 Foam::radiation::fvDOM::Ru() const
 {
+    tmp<DimensionedField<scalar, volMesh>> tRu
+    (
+        new DimensionedField<scalar, volMesh>
+        (
+            IOobject
+            (
+                "Ru",
+                mesh_.time().timeName(),
+                mesh_,
+                IOobject::NO_READ,
+                IOobject::NO_WRITE,
+                false
+            ),
+            mesh_,
+            dimensionedScalar("zero", dimensionSet(1, -1, -3, 0, 0), 0)
+        )
+    );
 
-    const volScalarField::Internal& G = G_();
-    const volScalarField::Internal E = absorptionEmission_->ECont()()();
-    const volScalarField::Internal a = a_.internalField();
+    DimensionedField<scalar, volMesh>& Ru=tRu.ref();
 
-    return a*G - E;
+    // Sum contributions over all frequency bands
+    for (label j=0; j < nLambda_; j++)
+    {
+        // Compute total incident radiation within frequency band
+        tmp<DimensionedField<scalar, volMesh>> Gj
+        (
+            IRay_[0].ILambda(j)()*IRay_[0].omega()
+        );
+
+        for (label rayI=1; rayI < nRay_; rayI++)
+        {
+            Gj.ref() += IRay_[rayI].ILambda(j)()*IRay_[rayI].omega();
+        }
+
+        Ru += (aLambda_[j]() - absorptionEmission_->aDisp(j)()())*Gj
+             - absorptionEmission_->ECont(j)()();
+    }
+
+    return tRu;
 }
 
 
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
index f37c3fc010..cbb35e3a99 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -25,6 +25,8 @@ License
 
 #include "wideBandAbsorptionEmission.H"
 #include "addToRunTimeSelectionTable.H"
+#include "basicSpecieMixture.H"
+#include "unitConversion.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -56,12 +58,7 @@ Foam::radiation::wideBandAbsorptionEmission::wideBandAbsorptionEmission
     coeffsDict_((dict.optionalSubDict(typeName + "Coeffs"))),
     speciesNames_(0),
     specieIndex_(label(0)),
-    lookUpTable_
-    (
-        fileName(coeffsDict_.lookup("lookUpTableFileName")),
-        mesh.time().constant(),
-        mesh
-    ),
+    lookUpTablePtr_(),
     thermo_(mesh.lookupObject<fluidThermo>(basicThermo::dictName)),
     Yj_(nSpecies_),
     totalWaveLength_(0)
@@ -95,7 +92,7 @@ Foam::radiation::wideBandAbsorptionEmission::wideBandAbsorptionEmission
                 if (!speciesNames_.found(key))
                 {
                     FatalErrorInFunction
-                        << "specie: " << key << "is not in all the bands"
+                        << "specie: " << key << " is not in all the bands"
                         << nl << exit(FatalError);
                 }
             }
@@ -106,37 +103,80 @@ Foam::radiation::wideBandAbsorptionEmission::wideBandAbsorptionEmission
     }
     nBands_ = nBand;
 
+    if (coeffsDict_.found("lookUpTableFileName"))
+    {
+        const word name = coeffsDict_.lookup("lookUpTableFileName");
+        if (name != "none")
+        {
+            lookUpTablePtr_.set
+            (
+                new interpolationLookUpTable<scalar>
+                (
+                    fileName(coeffsDict_.lookup("lookUpTableFileName")),
+                    mesh.time().constant(),
+                    mesh
+                )
+            );
+
+            if (!mesh.foundObject<volScalarField>("ft"))
+            {
+                FatalErrorInFunction
+                    << "specie ft is not present to use with "
+                    << "lookUpTableFileName " << nl
+                    << exit(FatalError);
+            }
+        }
+    }
+
     // Check that all the species on the dictionary are present in the
-    // look-up table  and save the corresponding indices of the look-up table
+    // look-up table and save the corresponding indices of the look-up table
 
     label j = 0;
     forAllConstIter(HashTable<label>, speciesNames_, iter)
     {
-        if (lookUpTable_.found(iter.key()))
+        if (!lookUpTablePtr_.empty())
         {
-            label index = lookUpTable_.findFieldIndex(iter.key());
-            Info<< "specie: " << iter.key() << " found in look-up table "
-                << " with index: " << index << endl;
-            specieIndex_[iter()] = index;
+            if (lookUpTablePtr_().found(iter.key()))
+            {
+                const label index =
+                    lookUpTablePtr_().findFieldIndex(iter.key());
+
+                Info<< "specie: " << iter.key() << " found on look-up table "
+                    << " with index: " << index << endl;
+
+                specieIndex_[iter()] = index;
+            }
+            else if (mesh.foundObject<volScalarField>(iter.key()))
+            {
+                Yj_.set(j, &mesh.lookupObjectRef<volScalarField>(iter.key()));
+                specieIndex_[iter()] = 0;
+                j++;
+                Info<< "specie: " << iter.key() << " is being solved" << endl;
+            }
+            else
+            {
+                FatalErrorInFunction
+                    << "specie: " << iter.key()
+                    << " is neither in look-up table: "
+                    << lookUpTablePtr_().tableName()
+                    << " nor is being solved" << nl
+                    << exit(FatalError);
+            }
         }
         else if (mesh.foundObject<volScalarField>(iter.key()))
         {
-            volScalarField& Y = const_cast<volScalarField&>
-                (mesh.lookupObject<volScalarField>(iter.key()));
-
-            Yj_.set(j, &Y);
-
-            specieIndex_[iter()] = 0.0;
+            Yj_.set(j, &mesh.lookupObjectRef<volScalarField>(iter.key()));
+            specieIndex_[iter()] = 0;
             j++;
-            Info<< "species: " << iter.key() << " is being solved" << endl;
         }
         else
         {
             FatalErrorInFunction
-                << "specie: " << iter.key()
-                << " is neither in look-up table : "
-                << lookUpTable_.tableName() << " nor is being solved"
+                << " there is no lookup table and the specie" << nl
+                << iter.key() << nl
+                << " is not found " << nl
                 << exit(FatalError);
+
         }
     }
 }
@@ -153,11 +193,11 @@ Foam::radiation::wideBandAbsorptionEmission::~wideBandAbsorptionEmission()
 Foam::tmp<Foam::volScalarField>
 Foam::radiation::wideBandAbsorptionEmission::aCont(const label bandi) const
 {
+    const basicSpecieMixture& mixture =
+        dynamic_cast<const basicSpecieMixture&>(thermo_);
+
     const volScalarField& T = thermo_.T();
     const volScalarField& p = thermo_.p();
-    const volScalarField& ft = mesh_.lookupObject<volScalarField>("ft");
-
-    label nSpecies = speciesNames_.size();
 
     tmp<volScalarField> ta
     (
@@ -178,38 +218,50 @@ Foam::radiation::wideBandAbsorptionEmission::aCont(const label bandi) const
 
     scalarField& a = ta.ref().primitiveFieldRef();
 
-    forAll(a, i)
+    forAll(a, celli)
     {
-        const List<scalar>& species = lookUpTable_.lookUp(ft[i]);
-
-        for (label n=0; n<nSpecies; n++)
+        forAllConstIter(HashTable<label>, speciesNames_, iter)
         {
-            label l = 0;
-            scalar Yipi = 0.0;
+            const label n = iter();
+            scalar Xipi = 0;
             if (specieIndex_[n] != 0)
             {
-                // moles x pressure [atm]
-                Yipi = species[specieIndex_[n]]*p[i]*9.869231e-6;
+                const volScalarField& ft =
+                    mesh_.lookupObject<volScalarField>("ft");
+
+                const List<scalar>& Ynft = lookUpTablePtr_().lookUp(ft[celli]);
+
+                // moles*pressure [atm]
+                Xipi = Ynft[specieIndex_[n]]*paToAtm(p[celli]);
             }
             else
             {
-                // mass fraction from species being solved
-                Yipi = Yj_[l][i];
-                l++;
+                scalar invWt = 0;
+                forAll(mixture.Y(), s)
+                {
+                    invWt += mixture.Y(s)[celli]/mixture.W(s);
+                }
+
+                const label index = mixture.species()[iter.key()];
+
+                const scalar Xk =
+                    mixture.Y(index)[celli]/(mixture.W(index)*invWt);
+
+                Xipi = Xk*paToAtm(p[celli]);
             }
 
-            scalar Ti = T[i];
+            scalar Ti = T[celli];
 
             const absorptionCoeffs::coeffArray& b =
-                coeffs_[n][bandi].coeffs(T[i]);
+                coeffs_[bandi][n].coeffs(T[celli]);
 
-            if (coeffs_[n][bandi].invTemp())
+            if (coeffs_[bandi][n].invTemp())
             {
-                Ti = 1.0/T[i];
+                Ti = 1.0/T[celli];
             }
 
-            a[i]+=
-                Yipi
+            a[celli]+=
+                Xipi
                *(
                     ((((b[5]*Ti + b[4])*Ti + b[3])*Ti + b[2])*Ti + b[1])*Ti
                   + b[0]
diff --git a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H
index 4149f25f58..ddcf30d621 100644
--- a/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H
+++ b/src/thermophysicalModels/radiation/submodels/absorptionEmissionModel/wideBandAbsorptionEmission/wideBandAbsorptionEmission.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -40,11 +40,8 @@ Description
 
     The emission constant proportionality is specified per band (EhrrCoeff).
 
-    The coefficients for the species in the lookup table have to be specified
-    for use in moles x P [atm].i.e. (k[i] = species[i]*p*9.869231e-6).
-
-    The coefficients for CO and soot or any other added are multiplied by the
-    respective mass fraction being solved.
+    The coefficients for the species have to be specified for use in
+    moles x P [atm], i.e. (k[i] = species[i]*p*9.869231e-6).
 
     The look Up table file should be in the constant directory.
 
@@ -159,8 +156,8 @@ private:
         //- Proportion of the heat released rate emitted
         FixedList<scalar, maxBands_> iEhrrCoeffs_;
 
-        //- Lookup table of species related to ft
-        mutable interpolationLookUpTable<scalar> lookUpTable_;
+        //- Look-up table of species related to ft
+        mutable autoPtr<interpolationLookUpTable<scalar>> lookUpTablePtr_;
 
         //- Thermo package
         const fluidThermo& thermo_;
@@ -193,24 +190,14 @@ public:
 
     // Member Functions
 
-        // Access
+        //- Absorption coefficient for continuous phase
+        tmp<volScalarField> aCont(const label bandi = 0) const;
 
-            // Absorption coefficient
+        //- Emission coefficient for continuous phase
+        tmp<volScalarField> eCont(const label bandi = 0) const;
 
-                //- Absorption coefficient for continuous phase
-                tmp<volScalarField> aCont(const label bandI = 0) const;
-
-
-            // Emission coefficient
-
-                //- Emission coefficient for continuous phase
-                tmp<volScalarField> eCont(const label bandI = 0) const;
-
-
-            // Emission contribution
-
-                //- Emission contribution for continuous phase
-                tmp<volScalarField> ECont(const label bandI = 0) const;
+        //- Emission contribution for continuous phase
+        tmp<volScalarField> ECont(const label bandi = 0) const;
 
 
         inline bool isGrey() const
@@ -225,15 +212,15 @@ public:
         }
 
         //- Lower and upper limit of band i
-        inline const Vector2D<scalar>& bands(const label i) const
+        inline const Vector2D<scalar>& bands(const label bandi) const
         {
-            return iBands_[i];
+            return iBands_[bandi];
         }
 
         //- Correct rays
         void correct
         (
-            volScalarField& a_,
+            volScalarField& a,
             PtrList<volScalarField>& aLambda
         ) const;
 };

From a8261a5626bed657ba64dc97fc0b014e30b59df4 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://openfoam.org>
Date: Thu, 22 Mar 2018 11:09:42 +0000
Subject: [PATCH 3/8] STYLE: fvPatchFields: Corrected documentation

Resolves bug-report https://bugs.openfoam.org/view.php?id=2885
---
 .../outletMappedUniformInletFvPatchField.C    | 30 +++++++++----------
 .../outletMappedUniformInletFvPatchField.H    | 10 +++----
 .../timeVaryingMappedFixedValueFvPatchField.H |  2 +-
 3 files changed, 21 insertions(+), 21 deletions(-)

diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C
index 5490a8a930..ddc9fd7eea 100644
--- a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C
+++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.C
@@ -43,6 +43,21 @@ outletMappedUniformInletFvPatchField
 {}
 
 
+template<class Type>
+Foam::outletMappedUniformInletFvPatchField<Type>::
+outletMappedUniformInletFvPatchField
+(
+    const fvPatch& p,
+    const DimensionedField<Type, volMesh>& iF,
+    const dictionary& dict
+)
+:
+    fixedValueFvPatchField<Type>(p, iF, dict),
+    outletPatchName_(dict.lookup("outletPatch")),
+    phiName_(dict.lookupOrDefault<word>("phi", "phi"))
+{}
+
+
 template<class Type>
 Foam::outletMappedUniformInletFvPatchField<Type>::
 outletMappedUniformInletFvPatchField
@@ -59,21 +74,6 @@ outletMappedUniformInletFvPatchField
 {}
 
 
-template<class Type>
-Foam::outletMappedUniformInletFvPatchField<Type>::
-outletMappedUniformInletFvPatchField
-(
-    const fvPatch& p,
-    const DimensionedField<Type, volMesh>& iF,
-    const dictionary& dict
-)
-:
-    fixedValueFvPatchField<Type>(p, iF, dict),
-    outletPatchName_(dict.lookup("outletPatch")),
-    phiName_(dict.lookupOrDefault<word>("phi", "phi"))
-{}
-
-
 template<class Type>
 Foam::outletMappedUniformInletFvPatchField<Type>::
 outletMappedUniformInletFvPatchField
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H
index dd367ca3f8..7c4ea1ceaa 100644
--- a/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H
+++ b/src/finiteVolume/fields/fvPatchFields/derived/outletMappedUniformInlet/outletMappedUniformInletFvPatchField.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,14 +29,14 @@ Group
 
 Description
     This boundary conditon averages the field over the "outlet" patch specified
-    by name "outletPatch" and applies this as the uniform value of the
-    field over this patch.
+    by name "outletPatch" and applies this as the uniform value of the field
+    over this patch.
 
 Usage
     \table
         Property        | Description             | Required    | Default value
-        outletPatch     | name of outlet patch    | yes         |
-        phi             | flux field name         | no          | phi
+        outletPatch     | Name of outlet patch    | yes         |
+        phi             | Flux field name         | no          | phi
     \endtable
 
     Example of the boundary condition specification:
diff --git a/src/finiteVolume/fields/fvPatchFields/derived/timeVaryingMappedFixedValue/timeVaryingMappedFixedValueFvPatchField.H b/src/finiteVolume/fields/fvPatchFields/derived/timeVaryingMappedFixedValue/timeVaryingMappedFixedValueFvPatchField.H
index 2320ffe2ab..46f720a723 100644
--- a/src/finiteVolume/fields/fvPatchFields/derived/timeVaryingMappedFixedValue/timeVaryingMappedFixedValueFvPatchField.H
+++ b/src/finiteVolume/fields/fvPatchFields/derived/timeVaryingMappedFixedValue/timeVaryingMappedFixedValueFvPatchField.H
@@ -51,7 +51,7 @@ Usage
         setAverage   | Switch to activate setting of average value | no | false
         perturb      | Perturb points for regular geometries    | no | 1e-5
         points       | Name of points file        | no          | points
-        fieldTableName | Alternative field name to sample | no | this field name
+        fieldTable   | Alternative field name to sample | no | this field name
         mapMethod    | Type of mapping            | no | planarInterpolation
         offset       | Offset to mapped values    | no | Zero
     \endtable

From 54ddd14168a13d108d287e8f8afea8397dd90345 Mon Sep 17 00:00:00 2001
From: Henry Weller <http://openfoam.org>
Date: Sun, 18 Mar 2018 14:08:20 +0000
Subject: [PATCH 4/8] STYLE: src/atmosphericModels: New library containing
 models specific to atmospheric flows

Converted the atmBoundaryLayerInlet boundary conditions to inletOutlet to handle
changes in inflow orientation without the need for re-meshing.
---
 .../turbulenceModels/Make/files               |  7 -----
 src/atmosphericModels/Make/files              |  7 +++++
 src/atmosphericModels/Make/options            | 18 +++++++++++
 .../atmBoundaryLayer/atmBoundaryLayer.C       |  0
 .../atmBoundaryLayer/atmBoundaryLayer.H       |  0
 ...ndaryLayerInletEpsilonFvPatchScalarField.C | 30 ++++++++++++++-----
 ...ndaryLayerInletEpsilonFvPatchScalarField.H |  6 ++--
 ...atmBoundaryLayerInletKFvPatchScalarField.C | 30 ++++++++++++++-----
 ...atmBoundaryLayerInletKFvPatchScalarField.H |  6 ++--
 ...daryLayerInletVelocityFvPatchVectorField.C | 30 ++++++++++++++-----
 ...daryLayerInletVelocityFvPatchVectorField.H |  6 ++--
 ...tkAtmRoughWallFunctionFvPatchScalarField.C |  2 +-
 ...tkAtmRoughWallFunctionFvPatchScalarField.H |  0
 .../turbineSiting/system/controlDict          |  2 ++
 14 files changed, 103 insertions(+), 41 deletions(-)
 create mode 100644 src/atmosphericModels/Make/files
 create mode 100644 src/atmosphericModels/Make/options
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.C (100%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H (100%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C (82%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H (97%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C (82%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H (97%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C (83%)
 rename src/{TurbulenceModels/turbulenceModels/RAS => atmosphericModels}/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H (97%)
 rename src/{TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions => atmosphericModels/derivedFvPatchFields}/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C (98%)
 rename src/{TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions => atmosphericModels/derivedFvPatchFields}/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.H (100%)

diff --git a/src/TurbulenceModels/turbulenceModels/Make/files b/src/TurbulenceModels/turbulenceModels/Make/files
index 315ba24226..7a6ad66294 100644
--- a/src/TurbulenceModels/turbulenceModels/Make/files
+++ b/src/TurbulenceModels/turbulenceModels/Make/files
@@ -38,7 +38,6 @@ $(nutWallFunctions)/nutWallFunction/nutWallFunctionFvPatchScalarField.C
 
 $(nutWallFunctions)/nutkWallFunction/nutkWallFunctionFvPatchScalarField.C
 $(nutWallFunctions)/nutkRoughWallFunction/nutkRoughWallFunctionFvPatchScalarField.C
-$(nutWallFunctions)/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
 
 $(nutWallFunctions)/nutUBlendedWallFunction/nutUBlendedWallFunctionFvPatchScalarField.C
 $(nutWallFunctions)/nutUWallFunction/nutUWallFunctionFvPatchScalarField.C
@@ -70,11 +69,5 @@ RASBCs = RAS/derivedFvPatchFields
 $(RASBCs)/turbulentMixingLengthDissipationRateInlet/turbulentMixingLengthDissipationRateInletFvPatchScalarField.C
 $(RASBCs)/turbulentMixingLengthFrequencyInlet/turbulentMixingLengthFrequencyInletFvPatchScalarField.C
 
-/* Atmospheric boundary layer BCs */
-$(RASBCs)/atmBoundaryLayer/atmBoundaryLayer.C
-$(RASBCs)/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
-$(RASBCs)/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
-$(RASBCs)/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
-
 
 LIB = $(FOAM_LIBBIN)/libturbulenceModels
diff --git a/src/atmosphericModels/Make/files b/src/atmosphericModels/Make/files
new file mode 100644
index 0000000000..a5a1c27f22
--- /dev/null
+++ b/src/atmosphericModels/Make/files
@@ -0,0 +1,7 @@
+derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.C
+derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
+derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
+derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
+derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
+
+LIB = $(FOAM_LIBBIN)/libatmosphericModels
diff --git a/src/atmosphericModels/Make/options b/src/atmosphericModels/Make/options
new file mode 100644
index 0000000000..3ffea18bec
--- /dev/null
+++ b/src/atmosphericModels/Make/options
@@ -0,0 +1,18 @@
+EXE_INC = \
+    -I$(LIB_SRC)/TurbulenceModels/turbulenceModels/lnInclude \
+    -I$(LIB_SRC)/TurbulenceModels/incompressible/lnInclude \
+    -I$(LIB_SRC)/transportModels \
+    -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/triSurface/lnInclude \
+    -I$(LIB_SRC)/fvOptions/lnInclude
+
+LIB_LIBS = \
+    -lturbulenceModels \
+    -lincompressibleTurbulenceModels \
+    -lincompressibleTransportModels \
+    -lfiniteVolume \
+    -lmeshTools \
+    -ltriSurface \
+    -lfvOptions
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.C
similarity index 100%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.C
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.C
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
similarity index 100%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
similarity index 82%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
index 20a917ec33..65e748f3dc 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,7 +43,7 @@ atmBoundaryLayerInletEpsilonFvPatchScalarField
     const DimensionedField<scalar, volMesh>& iF
 )
 :
-    fixedValueFvPatchScalarField(p, iF),
+    inletOutletFvPatchScalarField(p, iF),
     atmBoundaryLayer()
 {}
 
@@ -56,10 +56,23 @@ atmBoundaryLayerInletEpsilonFvPatchScalarField
     const dictionary& dict
 )
 :
-    fixedValueFvPatchScalarField(p, iF, dict, false),
+    inletOutletFvPatchScalarField(p, iF),
     atmBoundaryLayer(patch().Cf(), dict)
 {
-    scalarField::operator=(epsilon(patch().Cf()));
+    phiName_ = dict.lookupOrDefault<word>("phi", "phi");
+
+    refValue() = epsilon(patch().Cf());
+    refGrad() = 0;
+    valueFraction() = 1;
+
+    if (dict.found("value"))
+    {
+        scalarField::operator=(scalarField("value", dict, p.size()));
+    }
+    else
+    {
+        scalarField::operator=(refValue());
+    }
 }
 
 
@@ -72,7 +85,7 @@ atmBoundaryLayerInletEpsilonFvPatchScalarField
     const fvPatchFieldMapper& mapper
 )
 :
-    fixedValueFvPatchScalarField(psf, p, iF, mapper),
+    inletOutletFvPatchScalarField(psf, p, iF, mapper),
     atmBoundaryLayer(psf, mapper)
 {}
 
@@ -84,7 +97,7 @@ atmBoundaryLayerInletEpsilonFvPatchScalarField
     const DimensionedField<scalar, volMesh>& iF
 )
 :
-    fixedValueFvPatchScalarField(psf, iF),
+    inletOutletFvPatchScalarField(psf, iF),
     atmBoundaryLayer(psf)
 {}
 
@@ -96,7 +109,7 @@ void atmBoundaryLayerInletEpsilonFvPatchScalarField::autoMap
     const fvPatchFieldMapper& m
 )
 {
-    fixedValueFvPatchScalarField::autoMap(m);
+    inletOutletFvPatchScalarField::autoMap(m);
     atmBoundaryLayer::autoMap(m);
 }
 
@@ -107,7 +120,7 @@ void atmBoundaryLayerInletEpsilonFvPatchScalarField::rmap
     const labelList& addr
 )
 {
-    fixedValueFvPatchScalarField::rmap(psf, addr);
+    inletOutletFvPatchScalarField::rmap(psf, addr);
 
     const atmBoundaryLayerInletEpsilonFvPatchScalarField& blpsf =
         refCast<const atmBoundaryLayerInletEpsilonFvPatchScalarField>(psf);
@@ -120,6 +133,7 @@ void atmBoundaryLayerInletEpsilonFvPatchScalarField::write(Ostream& os) const
 {
     fvPatchScalarField::write(os);
     atmBoundaryLayer::write(os);
+    writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
     writeEntry("value", os);
 }
 
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H
similarity index 97%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H
index a73e943975..ff347f966c 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -60,7 +60,7 @@ SourceFiles
 #define atmBoundaryLayerInletEpsilonFvPatchScalarField_H
 
 #include "fvPatchFields.H"
-#include "fixedValueFvPatchFields.H"
+#include "inletOutletFvPatchFields.H"
 #include "atmBoundaryLayer.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -74,7 +74,7 @@ namespace Foam
 
 class atmBoundaryLayerInletEpsilonFvPatchScalarField
 :
-    public fixedValueFvPatchScalarField,
+    public inletOutletFvPatchScalarField,
     public atmBoundaryLayer
 {
 
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
similarity index 82%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
index 65017e2520..4078d85482 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,7 +43,7 @@ atmBoundaryLayerInletKFvPatchScalarField
     const DimensionedField<scalar, volMesh>& iF
 )
 :
-    fixedValueFvPatchScalarField(p, iF),
+    inletOutletFvPatchScalarField(p, iF),
     atmBoundaryLayer()
 {}
 
@@ -56,10 +56,23 @@ atmBoundaryLayerInletKFvPatchScalarField
     const dictionary& dict
 )
 :
-    fixedValueFvPatchScalarField(p, iF, dict, false),
+    inletOutletFvPatchScalarField(p, iF),
     atmBoundaryLayer(patch().Cf(), dict)
 {
-    scalarField::operator=(k(patch().Cf()));
+    phiName_ = dict.lookupOrDefault<word>("phi", "phi");
+
+    refValue() = k(patch().Cf());
+    refGrad() = 0;
+    valueFraction() = 1;
+
+    if (dict.found("value"))
+    {
+        scalarField::operator=(scalarField("value", dict, p.size()));
+    }
+    else
+    {
+        scalarField::operator=(refValue());
+    }
 }
 
 
@@ -72,7 +85,7 @@ atmBoundaryLayerInletKFvPatchScalarField
     const fvPatchFieldMapper& mapper
 )
 :
-    fixedValueFvPatchScalarField(psf, p, iF, mapper),
+    inletOutletFvPatchScalarField(psf, p, iF, mapper),
     atmBoundaryLayer(psf, mapper)
 {}
 
@@ -84,7 +97,7 @@ atmBoundaryLayerInletKFvPatchScalarField
     const DimensionedField<scalar, volMesh>& iF
 )
 :
-    fixedValueFvPatchScalarField(psf, iF),
+    inletOutletFvPatchScalarField(psf, iF),
     atmBoundaryLayer(psf)
 {}
 
@@ -96,7 +109,7 @@ void atmBoundaryLayerInletKFvPatchScalarField::autoMap
     const fvPatchFieldMapper& m
 )
 {
-    fixedValueFvPatchScalarField::autoMap(m);
+    inletOutletFvPatchScalarField::autoMap(m);
     atmBoundaryLayer::autoMap(m);
 }
 
@@ -107,7 +120,7 @@ void atmBoundaryLayerInletKFvPatchScalarField::rmap
     const labelList& addr
 )
 {
-    fixedValueFvPatchScalarField::rmap(psf, addr);
+    inletOutletFvPatchScalarField::rmap(psf, addr);
 
     const atmBoundaryLayerInletKFvPatchScalarField& blpsf =
         refCast<const atmBoundaryLayerInletKFvPatchScalarField>(psf);
@@ -120,6 +133,7 @@ void atmBoundaryLayerInletKFvPatchScalarField::write(Ostream& os) const
 {
     fvPatchScalarField::write(os);
     atmBoundaryLayer::write(os);
+    writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
     writeEntry("value", os);
 }
 
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H
similarity index 97%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H
index b0cac8c4dc..b575ca3e33 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2014-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2014-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -60,7 +60,7 @@ SourceFiles
 #define atmBoundaryLayerInletKFvPatchScalarField_H
 
 #include "fvPatchFields.H"
-#include "fixedValueFvPatchFields.H"
+#include "inletOutletFvPatchFields.H"
 #include "atmBoundaryLayer.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -74,7 +74,7 @@ namespace Foam
 
 class atmBoundaryLayerInletKFvPatchScalarField
 :
-    public fixedValueFvPatchScalarField,
+    public inletOutletFvPatchScalarField,
     public atmBoundaryLayer
 {
 
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
similarity index 83%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
index 348d041cb0..ec9ee6e79c 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -43,7 +43,7 @@ atmBoundaryLayerInletVelocityFvPatchVectorField
     const DimensionedField<vector, volMesh>& iF
 )
 :
-    fixedValueFvPatchVectorField(p, iF),
+    inletOutletFvPatchVectorField(p, iF),
     atmBoundaryLayer()
 {}
 
@@ -56,10 +56,23 @@ atmBoundaryLayerInletVelocityFvPatchVectorField
     const dictionary& dict
 )
 :
-    fixedValueFvPatchVectorField(p, iF, dict, false),
+    inletOutletFvPatchVectorField(p, iF),
     atmBoundaryLayer(patch().Cf(), dict)
 {
-    vectorField::operator=(U(patch().Cf()));
+    phiName_ = dict.lookupOrDefault<word>("phi", "phi");
+
+    refValue() = U(patch().Cf());
+    refGrad() = Zero;
+    valueFraction() = 1;
+
+    if (dict.found("value"))
+    {
+        vectorField::operator=(vectorField("value", dict, p.size()));
+    }
+    else
+    {
+        vectorField::operator=(refValue());
+    }
 }
 
 
@@ -72,7 +85,7 @@ atmBoundaryLayerInletVelocityFvPatchVectorField
     const fvPatchFieldMapper& mapper
 )
 :
-    fixedValueFvPatchVectorField(pvf, p, iF, mapper),
+    inletOutletFvPatchVectorField(pvf, p, iF, mapper),
     atmBoundaryLayer(pvf, mapper)
 {}
 
@@ -84,7 +97,7 @@ atmBoundaryLayerInletVelocityFvPatchVectorField
     const DimensionedField<vector, volMesh>& iF
 )
 :
-    fixedValueFvPatchVectorField(pvf, iF),
+    inletOutletFvPatchVectorField(pvf, iF),
     atmBoundaryLayer(pvf)
 {}
 
@@ -96,7 +109,7 @@ void atmBoundaryLayerInletVelocityFvPatchVectorField::autoMap
     const fvPatchFieldMapper& m
 )
 {
-    fixedValueFvPatchVectorField::autoMap(m);
+    inletOutletFvPatchVectorField::autoMap(m);
     atmBoundaryLayer::autoMap(m);
 }
 
@@ -107,7 +120,7 @@ void atmBoundaryLayerInletVelocityFvPatchVectorField::rmap
     const labelList& addr
 )
 {
-    fixedValueFvPatchVectorField::rmap(pvf, addr);
+    inletOutletFvPatchVectorField::rmap(pvf, addr);
 
     const atmBoundaryLayerInletVelocityFvPatchVectorField& blpvf =
         refCast<const atmBoundaryLayerInletVelocityFvPatchVectorField>(pvf);
@@ -120,6 +133,7 @@ void atmBoundaryLayerInletVelocityFvPatchVectorField::write(Ostream& os) const
 {
     fvPatchVectorField::write(os);
     atmBoundaryLayer::write(os);
+    writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
     writeEntry("value", os);
 }
 
diff --git a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H
similarity index 97%
rename from src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H
rename to src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H
index 7a22c0ef32..22b421f4aa 100644
--- a/src/TurbulenceModels/turbulenceModels/RAS/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -61,7 +61,7 @@ SourceFiles
 #define atmBoundaryLayerInletVelocityFvPatchVectorField_H
 
 #include "fvPatchFields.H"
-#include "fixedValueFvPatchFields.H"
+#include "inletOutletFvPatchFields.H"
 #include "atmBoundaryLayer.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
@@ -75,7 +75,7 @@ namespace Foam
 
 class atmBoundaryLayerInletVelocityFvPatchVectorField
 :
-    public fixedValueFvPatchVectorField,
+    public inletOutletFvPatchVectorField,
     public atmBoundaryLayer
 {
 
diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C b/src/atmosphericModels/derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
similarity index 98%
rename from src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
rename to src/atmosphericModels/derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
index 668a07a2c0..7bf60ee50d 100644
--- a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
diff --git a/src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.H b/src/atmosphericModels/derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.H
similarity index 100%
rename from src/TurbulenceModels/turbulenceModels/derivedFvPatchFields/wallFunctions/nutWallFunctions/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.H
rename to src/atmosphericModels/derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.H
diff --git a/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict b/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict
index 2654bdd256..bedecb22e9 100644
--- a/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict
+++ b/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict
@@ -44,4 +44,6 @@ timePrecision   6;
 
 runTimeModifiable yes;
 
+libs ("libatmosphericModels.so");
+
 // ************************************************************************* //

From c812b4c6f763f2845b70ab26b0a36968e2ac3b0b Mon Sep 17 00:00:00 2001
From: Henry Weller <http://openfoam.org>
Date: Tue, 20 Mar 2018 22:26:07 +0000
Subject: [PATCH 5/8] ENH: atmosphericModels: Added Lopes da Costa porosity and
 turbulence models

Specialized variants of the power law porosity and k epsilon turbulence models
developed to simulate atmospheric flow over forested and non-forested complex
terrain.

Class
    Foam::powerLawLopesdaCosta

Description
    Variant of the power law porosity model with spatially varying
    drag coefficient

    given by:

        \f[
            S = -\rho C_d \Sigma |U|^{(C_1 - 1)} U
        \f]

    where
    \vartable
        \Sigma | Porosity surface area per unit volume
        C_d    | Model linear coefficient
        C_1    | Model exponent coefficient
    \endvartable

    Reference:
    \verbatim
        Costa, J. C. P. L. D. (2007).
        Atmospheric flow over forested and non-forested complex terrain.
    \endverbatim

Class
    Foam::RASModels::kEpsilonLopesdaCosta

Description
    Variant of the standard k-epsilon turbulence model with additional source
    terms to handle the changes in turbulence in porous regions represented by
    the powerLawLopesdaCosta porosity model.

    Reference:
    \verbatim
        Costa, J. C. P. L. D. (2007).
        Atmospheric flow over forested and non-forested complex terrain.
    \endverbatim

    The default model coefficients are
    \verbatim
        kEpsilonLopesdaCostaCoeffs
        {
            Cmu         0.09;
            C1          1.44;
            C2          1.92;
            sigmak      1.0;
            sigmaEps    1.3;
        }
    \endverbatim

Tutorial case to follow.
---
 src/Allwmake                                  |   1 +
 src/atmosphericModels/Make/files              |   3 +
 .../atmosphericTurbulentTransportModels.C     |  35 ++
 .../atmBoundaryLayer/atmBoundaryLayer.H       |   2 +-
 .../kEpsilonLopesdaCosta.C                    | 477 ++++++++++++++++++
 .../kEpsilonLopesdaCosta.H                    | 245 +++++++++
 .../powerLawLopesdaCosta.C                    | 419 +++++++++++++++
 .../powerLawLopesdaCosta.H                    | 229 +++++++++
 .../powerLawLopesdaCostaTemplates.C           |  87 ++++
 .../porosityModel/porosityModel.H             |   5 +-
 .../porosityModel/porosityModelI.H            |   8 +-
 .../explicitPorositySource.C                  |   9 +-
 .../explicitPorositySource.H                  |   7 +-
 .../triSurfaceTools/triSurfaceTools.C         |  91 +++-
 .../triSurfaceTools/triSurfaceTools.H         |  14 +-
 15 files changed, 1617 insertions(+), 15 deletions(-)
 create mode 100644 src/atmosphericModels/atmosphericTurbulentTransportModels.C
 create mode 100644 src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.C
 create mode 100644 src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H
 create mode 100644 src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C
 create mode 100644 src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H
 create mode 100644 src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C

diff --git a/src/Allwmake b/src/Allwmake
index 92a389d77e..300792371c 100755
--- a/src/Allwmake
+++ b/src/Allwmake
@@ -82,6 +82,7 @@ wmake $targetType sixDoFRigidBodyState
 wmake $targetType rigidBodyDynamics
 wmake $targetType rigidBodyMeshMotion
 wmake $targetType semiPermeableBaffle
+wmake $targetType atmosphericModels
 
 # Needs access to Turbulence
 
diff --git a/src/atmosphericModels/Make/files b/src/atmosphericModels/Make/files
index a5a1c27f22..644149f555 100644
--- a/src/atmosphericModels/Make/files
+++ b/src/atmosphericModels/Make/files
@@ -4,4 +4,7 @@ derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarF
 derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
 derivedFvPatchFields/nutkAtmRoughWallFunction/nutkAtmRoughWallFunctionFvPatchScalarField.C
 
+atmosphericTurbulentTransportModels.C
+porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C
+
 LIB = $(FOAM_LIBBIN)/libatmosphericModels
diff --git a/src/atmosphericModels/atmosphericTurbulentTransportModels.C b/src/atmosphericModels/atmosphericTurbulentTransportModels.C
new file mode 100644
index 0000000000..cf31e0c250
--- /dev/null
+++ b/src/atmosphericModels/atmosphericTurbulentTransportModels.C
@@ -0,0 +1,35 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 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 "turbulentTransportModels.H"
+
+// -------------------------------------------------------------------------- //
+// RAS models
+// -------------------------------------------------------------------------- //
+
+#include "kEpsilonLopesdaCosta.H"
+makeRASModel(kEpsilonLopesdaCosta);
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
index 9c37660b49..a7360f8d38 100644
--- a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2014-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2014-2018 OpenFOAM Foundation
      \\/     M anipulation  | Copyright (C) 2015 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
diff --git a/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.C b/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.C
new file mode 100644
index 0000000000..805d0ca16b
--- /dev/null
+++ b/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.C
@@ -0,0 +1,477 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 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 "kEpsilonLopesdaCosta.H"
+#include "fvOptions.H"
+#include "explicitPorositySource.H"
+#include "bound.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace RASModels
+{
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+template<class BasicTurbulenceModel>
+void kEpsilonLopesdaCosta<BasicTurbulenceModel>::setPorosityCoefficient
+(
+    volScalarField::Internal& C,
+    const porosityModels::powerLawLopesdaCosta& pm
+)
+{
+    if (pm.dict().found(C.name()))
+    {
+        const labelList& cellZoneIDs = pm.cellZoneIDs();
+
+        const scalar Cpm = readScalar(pm.dict().lookup(C.name()));
+
+        forAll(cellZoneIDs, zonei)
+        {
+            const labelList& cells =
+                this->mesh_.cellZones()[cellZoneIDs[zonei]];
+
+            forAll(cells, i)
+            {
+                const label celli = cells[i];
+                C[celli] = Cpm;
+            }
+        }
+    }
+}
+
+
+template<class BasicTurbulenceModel>
+void kEpsilonLopesdaCosta<BasicTurbulenceModel>::setCdSigma
+(
+    volScalarField::Internal& C,
+    const porosityModels::powerLawLopesdaCosta& pm
+)
+{
+    if (pm.dict().found(C.name()))
+    {
+        const labelList& cellZoneIDs = pm.cellZoneIDs();
+        const scalarField& Sigma = pm.Sigma();
+
+        const scalar Cpm = readScalar(pm.dict().lookup(C.name()));
+
+        forAll(cellZoneIDs, zonei)
+        {
+            const labelList& cells =
+                this->mesh_.cellZones()[cellZoneIDs[zonei]];
+
+            forAll(cells, i)
+            {
+                const label celli = cells[i];
+                C[celli] = Cpm*Sigma[celli];
+            }
+        }
+    }
+}
+
+
+template<class BasicTurbulenceModel>
+void kEpsilonLopesdaCosta<BasicTurbulenceModel>::setPorosityCoefficients()
+{
+    fv::options::optionList& fvOptions(fv::options::New(this->mesh_));
+
+    forAll(fvOptions, i)
+    {
+        if (isA<fv::explicitPorositySource>(fvOptions[i]))
+        {
+            const fv::explicitPorositySource& eps =
+                refCast<const fv::explicitPorositySource>(fvOptions[i]);
+
+            if (isA<porosityModels::powerLawLopesdaCosta>(eps.model()))
+            {
+                const porosityModels::powerLawLopesdaCosta& pm =
+                    refCast<const porosityModels::powerLawLopesdaCosta>
+                    (
+                        eps.model()
+                    );
+
+                setPorosityCoefficient(Cmu_, pm);
+                setPorosityCoefficient(C1_, pm);
+                setPorosityCoefficient(C2_, pm);
+                setPorosityCoefficient(sigmak_, pm);
+                setPorosityCoefficient(sigmaEps_, pm);
+
+                setCdSigma(CdSigma_, pm);
+                setPorosityCoefficient(betap_, pm);
+                setPorosityCoefficient(betad_, pm);
+                setPorosityCoefficient(C4_, pm);
+                setPorosityCoefficient(C5_, pm);
+            }
+        }
+    }
+}
+
+
+template<class BasicTurbulenceModel>
+void kEpsilonLopesdaCosta<BasicTurbulenceModel>::correctNut()
+{
+    this->nut_ = Cmu_*sqr(k_)/epsilon_;
+    this->nut_.correctBoundaryConditions();
+    fv::options::New(this->mesh_).correct(this->nut_);
+
+    BasicTurbulenceModel::correctNut();
+}
+
+
+template<class BasicTurbulenceModel>
+tmp<fvScalarMatrix> kEpsilonLopesdaCosta<BasicTurbulenceModel>::kSource
+(
+    const volScalarField::Internal& magU,
+    const volScalarField::Internal& magU3
+) const
+{
+    return fvm::Su(CdSigma_*(betap_*magU3 - betad_*magU*k_()), k_);
+}
+
+
+template<class BasicTurbulenceModel>
+tmp<fvScalarMatrix>
+kEpsilonLopesdaCosta<BasicTurbulenceModel>::epsilonSource
+(
+    const volScalarField::Internal& magU,
+    const volScalarField::Internal& magU3
+) const
+{
+    return fvm::Su
+    (
+        CdSigma_
+       *(C4_*betap_*epsilon_()/k_()*magU3 - C5_*betad_*magU*epsilon_()),
+        epsilon_
+    );
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class BasicTurbulenceModel>
+kEpsilonLopesdaCosta<BasicTurbulenceModel>::kEpsilonLopesdaCosta
+(
+    const alphaField& alpha,
+    const rhoField& rho,
+    const volVectorField& U,
+    const surfaceScalarField& alphaRhoPhi,
+    const surfaceScalarField& phi,
+    const transportModel& transport,
+    const word& propertiesName,
+    const word& type
+)
+:
+    eddyViscosity<RASModel<BasicTurbulenceModel>>
+    (
+        type,
+        alpha,
+        rho,
+        U,
+        alphaRhoPhi,
+        phi,
+        transport,
+        propertiesName
+    ),
+
+    Cmu_
+    (
+        IOobject
+        (
+            "Cmu",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensioned<scalar>::lookupOrAddToDict
+        (
+            "Cmu",
+            this->coeffDict_,
+            0.09
+        )
+    ),
+    C1_
+    (
+        IOobject
+        (
+            "C1",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensioned<scalar>::lookupOrAddToDict
+        (
+            "C1",
+            this->coeffDict_,
+            1.44
+        )
+    ),
+    C2_
+    (
+        IOobject
+        (
+            "C2",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensioned<scalar>::lookupOrAddToDict
+        (
+            "C2",
+            this->coeffDict_,
+            1.92
+        )
+    ),
+    sigmak_
+    (
+        IOobject
+        (
+            "sigmak",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensioned<scalar>::lookupOrAddToDict
+        (
+            "sigmak",
+            this->coeffDict_,
+            1.0
+        )
+    ),
+    sigmaEps_
+    (
+        IOobject
+        (
+            "sigmaEps",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensioned<scalar>::lookupOrAddToDict
+        (
+            "sigmaEps",
+            this->coeffDict_,
+            1.3
+        )
+    ),
+
+    CdSigma_
+    (
+        IOobject
+        (
+            "CdSigma",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensionedScalar("CdSigma", dimless/dimLength, 0)
+    ),
+    betap_
+    (
+        IOobject
+        (
+            "betap",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensionedScalar("betap", dimless, 0)
+    ),
+    betad_
+    (
+        IOobject
+        (
+            "betad",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensionedScalar("betad", dimless, 0)
+    ),
+    C4_
+    (
+        IOobject
+        (
+            "C4",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensionedScalar("C4", dimless, 0)
+    ),
+    C5_
+    (
+        IOobject
+        (
+            "C5",
+            this->runTime_.timeName(),
+            this->mesh_
+        ),
+        this->mesh_,
+        dimensionedScalar("C5", dimless, 0)
+    ),
+
+    k_
+    (
+        IOobject
+        (
+            IOobject::groupName("k", alphaRhoPhi.group()),
+            this->runTime_.timeName(),
+            this->mesh_,
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE
+        ),
+        this->mesh_
+    ),
+    epsilon_
+    (
+        IOobject
+        (
+            IOobject::groupName("epsilon", alphaRhoPhi.group()),
+            this->runTime_.timeName(),
+            this->mesh_,
+            IOobject::MUST_READ,
+            IOobject::AUTO_WRITE
+        ),
+        this->mesh_
+    )
+{
+    bound(k_, this->kMin_);
+    bound(epsilon_, this->epsilonMin_);
+
+    if (type == typeName)
+    {
+        this->printCoeffs(type);
+    }
+
+    setPorosityCoefficients();
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class BasicTurbulenceModel>
+bool kEpsilonLopesdaCosta<BasicTurbulenceModel>::read()
+{
+    if (eddyViscosity<RASModel<BasicTurbulenceModel>>::read())
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+template<class BasicTurbulenceModel>
+void kEpsilonLopesdaCosta<BasicTurbulenceModel>::correct()
+{
+    if (!this->turbulence_)
+    {
+        return;
+    }
+
+    // Local references
+    const alphaField& alpha = this->alpha_;
+    const rhoField& rho = this->rho_;
+    const surfaceScalarField& alphaRhoPhi = this->alphaRhoPhi_;
+    const volVectorField& U = this->U_;
+    volScalarField& nut = this->nut_;
+    fv::options& fvOptions(fv::options::New(this->mesh_));
+
+    eddyViscosity<RASModel<BasicTurbulenceModel>>::correct();
+
+    volScalarField::Internal divU
+    (
+        fvc::div(fvc::absolute(this->phi(), U))().v()
+    );
+
+    tmp<volTensorField> tgradU = fvc::grad(U);
+    volScalarField::Internal G
+    (
+        this->GName(),
+        nut.v()*(dev(twoSymm(tgradU().v())) && tgradU().v())
+    );
+    tgradU.clear();
+
+    // Update epsilon and G at the wall
+    epsilon_.boundaryFieldRef().updateCoeffs();
+
+    volScalarField::Internal magU(mag(U));
+    volScalarField::Internal magU3(pow3(magU));
+
+    // Dissipation equation
+    tmp<fvScalarMatrix> epsEqn
+    (
+        fvm::ddt(alpha, rho, epsilon_)
+      + fvm::div(alphaRhoPhi, epsilon_)
+      - fvm::laplacian(alpha*rho*DepsilonEff(), epsilon_)
+     ==
+        C1_*alpha()*rho()*G*epsilon_()/k_()
+      - fvm::SuSp(((2.0/3.0)*C1_)*alpha()*rho()*divU, epsilon_)
+      - fvm::Sp(C2_*alpha()*rho()*epsilon_()/k_(), epsilon_)
+      + epsilonSource(magU, magU3)
+      + fvOptions(alpha, rho, epsilon_)
+    );
+
+    epsEqn.ref().relax();
+    fvOptions.constrain(epsEqn.ref());
+    epsEqn.ref().boundaryManipulate(epsilon_.boundaryFieldRef());
+    solve(epsEqn);
+    fvOptions.correct(epsilon_);
+    bound(epsilon_, this->epsilonMin_);
+
+    // Turbulent kinetic energy equation
+    tmp<fvScalarMatrix> kEqn
+    (
+        fvm::ddt(alpha, rho, k_)
+      + fvm::div(alphaRhoPhi, k_)
+      - fvm::laplacian(alpha*rho*DkEff(), k_)
+     ==
+        alpha()*rho()*G
+      - fvm::SuSp((2.0/3.0)*alpha()*rho()*divU, k_)
+      - fvm::Sp(alpha()*rho()*epsilon_()/k_(), k_)
+      + kSource(magU, magU3)
+      + fvOptions(alpha, rho, k_)
+    );
+
+    kEqn.ref().relax();
+    fvOptions.constrain(kEqn.ref());
+    solve(kEqn);
+    fvOptions.correct(k_);
+    bound(k_, this->kMin_);
+
+    correctNut();
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace Foam
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H b/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H
new file mode 100644
index 0000000000..b72992d9a4
--- /dev/null
+++ b/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H
@@ -0,0 +1,245 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 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::RASModels::kEpsilonLopesdaCosta
+
+Group
+    grpRASTurbulence
+
+Description
+    Variant of the standard k-epsilon turbulence model with additional source
+    terms to handle the changes in turbulence in porous regions represented by
+    the powerLawLopesdaCosta porosity model.
+
+    Reference:
+    \verbatim
+        Costa, J. C. P. L. D. (2007).
+        Atmospheric flow over forested and non-forested complex terrain.
+    \endverbatim
+
+    The default model coefficients are
+    \verbatim
+        kEpsilonLopesdaCostaCoeffs
+        {
+            Cmu         0.09;
+            C1          1.44;
+            C2          1.92;
+            sigmak      1.0;
+            sigmaEps    1.3;
+        }
+    \endverbatim
+
+See also
+    Foam::RASModels::kEpsilon
+    Foam::porosityModels::powerLawLopesdaCosta
+
+SourceFiles
+    kEpsilonLopesdaCosta.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef kEpsilonLopesdaCosta_H
+#define kEpsilonLopesdaCosta_H
+
+#include "RASModel.H"
+#include "eddyViscosity.H"
+#include "powerLawLopesdaCosta.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace RASModels
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class kEpsilonLopesdaCosta Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class BasicTurbulenceModel>
+class kEpsilonLopesdaCosta
+:
+    public eddyViscosity<RASModel<BasicTurbulenceModel>>
+{
+    // Private Member Functions
+
+        // Disallow default bitwise copy construct and assignment
+        kEpsilonLopesdaCosta(const kEpsilonLopesdaCosta&);
+        void operator=(const kEpsilonLopesdaCosta&);
+
+
+protected:
+
+    // Protected data
+
+        // Standard model coefficients
+
+            volScalarField Cmu_;
+            volScalarField::Internal C1_;
+            volScalarField::Internal C2_;
+            volScalarField sigmak_;
+            volScalarField sigmaEps_;
+
+        // Lopes da Costa porosity coefficients
+
+            volScalarField::Internal CdSigma_;
+            volScalarField::Internal betap_;
+            volScalarField::Internal betad_;
+            volScalarField::Internal C4_;
+            volScalarField::Internal C5_;
+
+
+        // Fields
+
+            volScalarField k_;
+            volScalarField epsilon_;
+
+
+    // Protected Member Functions
+
+        void setPorosityCoefficient
+        (
+            volScalarField::Internal& C,
+            const porosityModels::powerLawLopesdaCosta& pm
+        );
+
+        void setCdSigma
+        (
+            volScalarField::Internal& C,
+            const porosityModels::powerLawLopesdaCosta& pm
+        );
+
+        void setPorosityCoefficients();
+
+        virtual void correctNut();
+
+        virtual tmp<fvScalarMatrix> kSource
+        (
+            const volScalarField::Internal& magU,
+            const volScalarField::Internal& magU3
+        ) const;
+
+        virtual tmp<fvScalarMatrix> epsilonSource
+        (
+            const volScalarField::Internal& magU,
+            const volScalarField::Internal& magU3
+        ) const;
+
+
+public:
+
+    typedef typename BasicTurbulenceModel::alphaField alphaField;
+    typedef typename BasicTurbulenceModel::rhoField rhoField;
+    typedef typename BasicTurbulenceModel::transportModel transportModel;
+
+
+    //- Runtime type information
+    TypeName("kEpsilonLopesdaCosta");
+
+
+    // Constructors
+
+        //- Construct from components
+        kEpsilonLopesdaCosta
+        (
+            const alphaField& alpha,
+            const rhoField& rho,
+            const volVectorField& U,
+            const surfaceScalarField& alphaRhoPhi,
+            const surfaceScalarField& phi,
+            const transportModel& transport,
+            const word& propertiesName = turbulenceModel::propertiesName,
+            const word& type = typeName
+        );
+
+
+    //- Destructor
+    virtual ~kEpsilonLopesdaCosta()
+    {}
+
+
+    // Member Functions
+
+        //- Re-read model coefficients if they have changed
+        virtual bool read();
+
+        //- Return the effective diffusivity for k
+        tmp<volScalarField> DkEff() const
+        {
+            return tmp<volScalarField>
+            (
+                new volScalarField
+                (
+                    "DkEff",
+                    (this->nut_/sigmak_ + this->nu())
+                )
+            );
+        }
+
+        //- Return the effective diffusivity for epsilon
+        tmp<volScalarField> DepsilonEff() const
+        {
+            return tmp<volScalarField>
+            (
+                new volScalarField
+                (
+                    "DepsilonEff",
+                    (this->nut_/sigmaEps_ + this->nu())
+                )
+            );
+        }
+
+        //- Return the turbulence kinetic energy
+        virtual tmp<volScalarField> k() const
+        {
+            return k_;
+        }
+
+        //- Return the turbulence kinetic energy dissipation rate
+        virtual tmp<volScalarField> epsilon() const
+        {
+            return epsilon_;
+        }
+
+        //- Solve the turbulence equations and correct the turbulence viscosity
+        virtual void correct();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace RASModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "kEpsilonLopesdaCosta.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C
new file mode 100644
index 0000000000..e5d13f50e9
--- /dev/null
+++ b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.C
@@ -0,0 +1,419 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 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 "addToRunTimeSelectionTable.H"
+#include "powerLawLopesdaCosta.H"
+#include "geometricOneField.H"
+#include "fvMatrices.H"
+#include "triSurfaceMesh.H"
+#include "triSurfaceTools.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    namespace porosityModels
+    {
+        defineTypeNameAndDebug(powerLawLopesdaCosta, 0);
+        addToRunTimeSelectionTable(porosityModel, powerLawLopesdaCosta, mesh);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::porosityModels::powerLawLopesdaCostaZone::powerLawLopesdaCostaZone
+(
+    const word& name,
+    const word& modelType,
+    const fvMesh& mesh,
+    const dictionary& dict
+)
+:
+    zoneName_(name + ":porousZone")
+{
+    dictionary coeffs(dict.optionalSubDict(modelType + "Coeffs"));
+
+    // Vertical direction within porous region
+    vector zDir(coeffs.lookup("zDir"));
+
+    // Span of the search for the cells in the porous region
+    scalar searchSpan(readScalar(coeffs.lookup("searchSpan")));
+
+    // Top surface file name defining the extent of the porous region
+    word topSurfaceFileName(coeffs.lookup("topSurface"));
+
+    // List of the ground patches defining the lower surface
+    // of the porous region
+    List<wordRe> groundPatches(coeffs.lookup("groundPatches"));
+
+    // Functional form of the porosity surface area per unit volume as a
+    // function of the normalized vertical position
+    autoPtr<Function1<scalar >> SigmaFunc
+    (
+        Function1<scalar>::New("Sigma", dict)
+    );
+
+    // Searchable triSurface for the top of the porous region
+    triSurfaceMesh searchSurf
+    (
+        IOobject
+        (
+            topSurfaceFileName,
+            mesh.time().constant(),
+            "triSurface",
+            mesh.time()
+        )
+    );
+
+    // Limit the porous cell search to those within the lateral and vertical
+    // extent of the top surface
+
+    boundBox surfBounds(searchSurf.points());
+    boundBox searchBounds
+    (
+        surfBounds.min() - searchSpan*zDir, surfBounds.max()
+    );
+
+    const pointField& C = mesh.cellCentres();
+
+    // List of cells within the porous region
+    labelList porousCells(C.size());
+    label porousCelli = 0;
+
+    forAll(C, celli)
+    {
+        if (searchBounds.contains(C[celli]))
+        {
+            porousCells[porousCelli++] = celli;
+        }
+    }
+
+    porousCells.setSize(porousCelli);
+
+    pointField start(porousCelli);
+    pointField end(porousCelli);
+
+    forAll(porousCells, porousCelli)
+    {
+        start[porousCelli] = C[porousCells[porousCelli]];
+        end[porousCelli] = start[porousCelli] + searchSpan*zDir;
+    }
+
+    // Field of distance between the cell centre and the corresponding point
+    // on the porous region top surface
+    scalarField zTop(porousCelli);
+
+    // Search the porous cells for those with a corresponding point on the
+    // porous region top surface
+    List<pointIndexHit> hitInfo;
+    searchSurf.findLine(start, end, hitInfo);
+
+    porousCelli = 0;
+
+    forAll(porousCells, celli)
+    {
+        const pointIndexHit& hit = hitInfo[celli];
+
+        if (hit.hit())
+        {
+            porousCells[porousCelli] = porousCells[celli];
+
+            zTop[porousCelli] =
+                (hit.hitPoint() - C[porousCells[porousCelli]]) & zDir;
+
+            porousCelli++;
+        }
+    }
+
+    // Resize the porous cells list and height field
+    porousCells.setSize(porousCelli);
+    zTop.setSize(porousCelli);
+
+    // Create a triSurface for the ground patch(s)
+    triSurface groundSurface
+    (
+        triSurfaceTools::triangulate
+        (
+            mesh.boundaryMesh(),
+            mesh.boundaryMesh().patchSet(groundPatches),
+            searchBounds
+        )
+    );
+
+    // Combine the ground triSurfaces across all processors
+    if (Pstream::parRun())
+    {
+        List<List<labelledTri>> groundSurfaceProcTris(Pstream::nProcs());
+        List<pointField> groundSurfaceProcPoints(Pstream::nProcs());
+
+        groundSurfaceProcTris[Pstream::myProcNo()] = groundSurface;
+        groundSurfaceProcPoints[Pstream::myProcNo()] = groundSurface.points();
+
+        Pstream::gatherList(groundSurfaceProcTris);
+        Pstream::scatterList(groundSurfaceProcTris);
+        Pstream::gatherList(groundSurfaceProcPoints);
+        Pstream::scatterList(groundSurfaceProcPoints);
+
+        label nTris = 0;
+        forAll(groundSurfaceProcTris, i)
+        {
+            nTris += groundSurfaceProcTris[i].size();
+        }
+
+        List<labelledTri> groundSurfaceTris(nTris);
+        label trii = 0;
+        label offset = 0;
+        forAll(groundSurfaceProcTris, i)
+        {
+            forAll(groundSurfaceProcTris[i], j)
+            {
+                groundSurfaceTris[trii] = groundSurfaceProcTris[i][j];
+                groundSurfaceTris[trii][0] += offset;
+                groundSurfaceTris[trii][1] += offset;
+                groundSurfaceTris[trii][2] += offset;
+                trii++;
+            }
+            offset += groundSurfaceProcPoints[i].size();
+        }
+
+        label nPoints = 0;
+        forAll(groundSurfaceProcPoints, i)
+        {
+            nPoints += groundSurfaceProcPoints[i].size();
+        }
+
+        pointField groundSurfacePoints(nPoints);
+
+        label pointi = 0;
+        forAll(groundSurfaceProcPoints, i)
+        {
+            forAll(groundSurfaceProcPoints[i], j)
+            {
+                groundSurfacePoints[pointi++] = groundSurfaceProcPoints[i][j];
+            }
+        }
+
+        groundSurface = triSurface(groundSurfaceTris, groundSurfacePoints);
+    }
+
+    // Create a searchable triSurface for the ground
+    triSurfaceSearch groundSearch(groundSurface);
+
+    // Search the porous cells for the corresponding point on the ground
+
+    start.setSize(porousCelli);
+    end.setSize(porousCelli);
+
+    forAll(porousCells, porousCelli)
+    {
+        start[porousCelli] = C[porousCells[porousCelli]];
+        end[porousCelli] = start[porousCelli] - searchSpan*zDir;
+    }
+
+    groundSearch.findLine(start, end, hitInfo);
+
+    scalarField zBottom(porousCelli);
+
+    forAll(porousCells, porousCelli)
+    {
+        const pointIndexHit& hit = hitInfo[porousCelli];
+
+        if (hit.hit())
+        {
+            zBottom[porousCelli] =
+                (C[porousCells[porousCelli]] - hit.hitPoint()) & zDir;
+        }
+    }
+
+    // Create the normalized height field
+    scalarField zNorm(zBottom/(zBottom + zTop));
+
+    // Create the porosity surface area per unit volume zone field
+    Sigma_ = SigmaFunc->value(zNorm);
+
+    // Create the porous region cellZone and add to the mesh cellZones
+
+    cellZoneMesh& cellZones = const_cast<cellZoneMesh&>(mesh.cellZones());
+
+    label zoneID = cellZones.findZoneID(zoneName_);
+
+    if (zoneID == -1)
+    {
+        zoneID = cellZones.size();
+        cellZones.setSize(zoneID + 1);
+
+        cellZones.set
+        (
+            zoneID,
+            new cellZone
+            (
+                zoneName_,
+                porousCells,
+                zoneID,
+                cellZones
+            )
+        );
+    }
+    else
+    {
+        FatalErrorInFunction
+            << "Unable to create porous cellZone " << zoneName_
+            << ": zone already exists"
+            << abort(FatalError);
+    }
+}
+
+
+Foam::porosityModels::powerLawLopesdaCosta::powerLawLopesdaCosta
+(
+    const word& name,
+    const word& modelType,
+    const fvMesh& mesh,
+    const dictionary& dict,
+    const word& dummyCellZoneName
+)
+:
+    powerLawLopesdaCostaZone(name, modelType, mesh, dict),
+    porosityModel
+    (
+        name,
+        modelType,
+        mesh,
+        dict,
+        powerLawLopesdaCostaZone::zoneName_
+    ),
+    Cd_(readScalar(coeffs_.lookup("Cd"))),
+    C1_(readScalar(coeffs_.lookup("C1"))),
+    rhoName_(coeffs_.lookupOrDefault<word>("rho", "rho"))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::porosityModels::powerLawLopesdaCosta::~powerLawLopesdaCosta()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+const Foam::scalarField&
+Foam::porosityModels::powerLawLopesdaCostaZone::Sigma() const
+{
+    return Sigma_;
+}
+
+
+void Foam::porosityModels::powerLawLopesdaCosta::calcTransformModelData()
+{}
+
+
+void Foam::porosityModels::powerLawLopesdaCosta::calcForce
+(
+    const volVectorField& U,
+    const volScalarField& rho,
+    const volScalarField& mu,
+    vectorField& force
+) const
+{
+    scalarField Udiag(U.size(), 0.0);
+    const scalarField& V = mesh_.V();
+
+    apply(Udiag, V, rho, U);
+
+    force = Udiag*U;
+}
+
+
+void Foam::porosityModels::powerLawLopesdaCosta::correct
+(
+    fvVectorMatrix& UEqn
+) const
+{
+    const vectorField& U = UEqn.psi();
+    const scalarField& V = mesh_.V();
+    scalarField& Udiag = UEqn.diag();
+
+    if (UEqn.dimensions() == dimForce)
+    {
+        const volScalarField& rho =
+            mesh_.lookupObject<volScalarField>(rhoName_);
+
+        apply(Udiag, V, rho, U);
+    }
+    else
+    {
+        apply(Udiag, V, geometricOneField(), U);
+    }
+}
+
+
+void Foam::porosityModels::powerLawLopesdaCosta::correct
+(
+    fvVectorMatrix& UEqn,
+    const volScalarField& rho,
+    const volScalarField& mu
+) const
+{
+    const vectorField& U = UEqn.psi();
+    const scalarField& V = mesh_.V();
+    scalarField& Udiag = UEqn.diag();
+
+    apply(Udiag, V, rho, U);
+}
+
+
+void Foam::porosityModels::powerLawLopesdaCosta::correct
+(
+    const fvVectorMatrix& UEqn,
+    volTensorField& AU
+) const
+{
+    const vectorField& U = UEqn.psi();
+
+    if (UEqn.dimensions() == dimForce)
+    {
+        const volScalarField& rho =
+            mesh_.lookupObject<volScalarField>(rhoName_);
+
+        apply(AU, rho, U);
+    }
+    else
+    {
+        apply(AU, geometricOneField(), U);
+    }
+}
+
+
+bool Foam::porosityModels::powerLawLopesdaCosta::writeData(Ostream& os) const
+{
+    os  << indent << name_ << endl;
+    dict_.write(os);
+
+    return true;
+}
+
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H
new file mode 100644
index 0000000000..f19677c808
--- /dev/null
+++ b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCosta.H
@@ -0,0 +1,229 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 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::powerLawLopesdaCosta
+
+Description
+    Variant of the power law porosity model with spatially varying
+    drag coefficient
+
+    given by:
+
+        \f[
+            S = -\rho C_d \Sigma |U|^{(C_1 - 1)} U
+        \f]
+
+    where
+    \vartable
+        \Sigma | Porosity surface area per unit volume
+        C_d    | Model linear coefficient
+        C_1    | Model exponent coefficient
+    \endvartable
+
+    Reference:
+    \verbatim
+        Costa, J. C. P. L. D. (2007).
+        Atmospheric flow over forested and non-forested complex terrain.
+    \endverbatim
+
+See also
+    Foam::RASModels::kEpsilonLopesdaCosta
+
+SourceFiles
+    powerLawLopesdaCosta.C
+    powerLawLopesdaCostaTemplates.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef powerLawLopesdaCosta_H
+#define powerLawLopesdaCosta_H
+
+#include "porosityModel.H"
+#include "Function1.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace porosityModels
+{
+
+/*---------------------------------------------------------------------------*\
+                  Class powerLawLopesdaCostaZone Declaration
+\*---------------------------------------------------------------------------*/
+
+class powerLawLopesdaCostaZone
+{
+protected:
+
+    // Protected data
+
+        //- Automatically generated zone name for this porous zone
+        const word zoneName_;
+
+        //- Porosity surface area per unit volume zone field
+        scalarField Sigma_;
+
+
+public:
+
+    //- Constructor
+    powerLawLopesdaCostaZone
+    (
+        const word& name,
+        const word& modelType,
+        const fvMesh& mesh,
+        const dictionary& dict
+    );
+
+    // Member Functions
+
+        //- Return the porosity surface area per unit volume zone field
+        const scalarField& Sigma() const;
+};
+
+
+/*---------------------------------------------------------------------------*\
+                      Class powerLawLopesdaCosta Declaration
+\*---------------------------------------------------------------------------*/
+
+class powerLawLopesdaCosta
+:
+    public powerLawLopesdaCostaZone,
+    public porosityModel
+{
+    // Private data
+
+        //- Cd coefficient
+        scalar Cd_;
+
+        //- C1 coefficient
+        scalar C1_;
+
+        //- Name of density field
+        word rhoName_;
+
+
+    // Private Member Functions
+
+        //- Apply resistance
+        template<class RhoFieldType>
+        void apply
+        (
+            scalarField& Udiag,
+            const scalarField& V,
+            const RhoFieldType& rho,
+            const vectorField& U
+        ) const;
+
+        //- Apply resistance
+        template<class RhoFieldType>
+        void apply
+        (
+            tensorField& AU,
+            const RhoFieldType& rho,
+            const vectorField& U
+        ) const;
+
+        //- Disallow default bitwise copy construct
+        powerLawLopesdaCosta(const powerLawLopesdaCosta&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const powerLawLopesdaCosta&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("powerLawLopesdaCosta");
+
+    //- Constructor
+    powerLawLopesdaCosta
+    (
+        const word& name,
+        const word& modelType,
+        const fvMesh& mesh,
+        const dictionary& dict,
+        const word& cellZoneName
+    );
+
+    //- Destructor
+    virtual ~powerLawLopesdaCosta();
+
+
+    // Member Functions
+
+        //- Transform the model data wrt mesh changes
+        virtual void calcTransformModelData();
+
+        //- Calculate the porosity force
+        virtual void calcForce
+        (
+            const volVectorField& U,
+            const volScalarField& rho,
+            const volScalarField& mu,
+            vectorField& force
+        ) const;
+
+        //- Add resistance
+        virtual void correct(fvVectorMatrix& UEqn) const;
+
+        //- Add resistance
+        virtual void correct
+        (
+            fvVectorMatrix& UEqn,
+            const volScalarField& rho,
+            const volScalarField& mu
+        ) const;
+
+        //- Add resistance
+        virtual void correct
+        (
+            const fvVectorMatrix& UEqn,
+            volTensorField& AU
+        ) const;
+
+
+    // I-O
+
+        //- Write
+        bool writeData(Ostream& os) const;
+};
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace porosityModels
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "powerLawLopesdaCostaTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C
new file mode 100644
index 0000000000..be9142d191
--- /dev/null
+++ b/src/atmosphericModels/porosityModels/powerLawLopesdaCosta/powerLawLopesdaCostaTemplates.C
@@ -0,0 +1,87 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 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 "volFields.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class RhoFieldType>
+void Foam::porosityModels::powerLawLopesdaCosta::apply
+(
+    scalarField& Udiag,
+    const scalarField& V,
+    const RhoFieldType& rho,
+    const vectorField& U
+) const
+{
+    const scalar C1m1b2 = (C1_ - 1.0)/2.0;
+
+    forAll(cellZoneIDs_, zonei)
+    {
+        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zonei]];
+
+        forAll(cells, i)
+        {
+            const label celli = cells[i];
+
+            Udiag[celli] +=
+                V[celli]*rho[celli]
+               *Cd_*Sigma_[i]*pow(magSqr(U[celli]), C1m1b2);
+        }
+    }
+}
+
+
+template<class RhoFieldType>
+void Foam::porosityModels::powerLawLopesdaCosta::apply
+(
+    tensorField& AU,
+    const RhoFieldType& rho,
+    const vectorField& U
+) const
+{
+    const scalar C1m1b2 = (C1_ - 1.0)/2.0;
+
+    forAll(cellZoneIDs_, zonei)
+    {
+        const labelList& cells = mesh_.cellZones()[cellZoneIDs_[zonei]];
+
+        forAll(cells, i)
+        {
+            const label celli = cells[i];
+
+            AU[celli] =
+                AU[celli]
+              + I
+               *(
+                   0.5*rho[celli]*Cd_*Sigma_[i]
+                  *pow(magSqr(U[celli]), C1m1b2)
+                );
+        }
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H
index bc7c496dd6..e3b00a535c 100644
--- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H
+++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModel.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2015 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -223,6 +223,9 @@ public:
         //- Return const access to the cell zone IDs
         inline const labelList& cellZoneIDs() const;
 
+        //- Return dictionary used for model construction
+        const dictionary& dict() const;
+
         //- Transform the model data wrt mesh changes
         virtual void transformModelData();
 
diff --git a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H
index 63f28091b3..ecf384a274 100644
--- a/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H
+++ b/src/finiteVolume/cfdTools/general/porosityModel/porosityModel/porosityModelI.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -35,6 +35,12 @@ inline bool Foam::porosityModel::active() const
 }
 
 
+inline const Foam::dictionary& Foam::porosityModel::dict() const
+{
+    return dict_;
+}
+
+
 inline const Foam::labelList& Foam::porosityModel::cellZoneIDs() const
 {
     return cellZoneIDs_;
diff --git a/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.C b/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.C
index 8b8c37762c..a90b1903b3 100644
--- a/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.C
+++ b/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.C
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -61,13 +61,6 @@ Foam::fv::explicitPorositySource::explicitPorositySource
 {
     read(dict);
 
-    if (selectionMode_ != smCellZone)
-    {
-        FatalErrorInFunction
-            << "selection mode is " << selectionModeTypeNames_[selectionMode_]
-            << exit(FatalError);
-    }
-
     porosityPtr_.reset
     (
         porosityModel::New
diff --git a/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.H b/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.H
index 576fac6843..c521774881 100644
--- a/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.H
+++ b/src/fvOptions/sources/derived/explicitPorositySource/explicitPorositySource.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2018 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -134,6 +134,11 @@ public:
 
     // Member Functions
 
+        const porosityModel& model() const
+        {
+            return porosityPtr_();
+        }
+
         // Add explicit and implicit contributions
 
             //- Add implicit contribution to momentum equation
diff --git a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C
index 1aaf5824f7..09bc91ed8d 100644
--- a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C
+++ b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.C
@@ -2276,7 +2276,6 @@ Foam::triSurfaceTools::sideType Foam::triSurfaceTools::surfaceSide
 }
 
 
-// triangulation of boundaryMesh
 Foam::triSurface Foam::triSurfaceTools::triangulate
 (
     const polyBoundaryMesh& bMesh,
@@ -2357,6 +2356,96 @@ Foam::triSurface Foam::triSurfaceTools::triangulate
 }
 
 
+Foam::triSurface Foam::triSurfaceTools::triangulate
+(
+    const polyBoundaryMesh& bMesh,
+    const labelHashSet& includePatches,
+    const boundBox& bBox,
+    const bool verbose
+)
+{
+    const polyMesh& mesh = bMesh.mesh();
+
+    // Storage for surfaceMesh. Size estimate.
+    DynamicList<labelledTri> triangles
+    (
+        mesh.nFaces() - mesh.nInternalFaces()
+    );
+
+    label newPatchi = 0;
+
+    forAllConstIter(labelHashSet, includePatches, iter)
+    {
+        const label patchi = iter.key();
+        const polyPatch& patch = bMesh[patchi];
+        const pointField& points = patch.points();
+
+        label nTriTotal = 0;
+
+        forAll(patch, patchFacei)
+        {
+            const face& f = patch[patchFacei];
+
+            if (bBox.containsAny(points, f))
+            {
+                faceList triFaces(f.nTriangles(points));
+
+                label nTri = 0;
+
+                f.triangles(points, nTri, triFaces);
+
+                forAll(triFaces, triFacei)
+                {
+                    const face& f = triFaces[triFacei];
+
+                    triangles.append(labelledTri(f[0], f[1], f[2], newPatchi));
+
+                    nTriTotal++;
+                }
+            }
+        }
+
+        if (verbose)
+        {
+            Pout<< patch.name() << " : generated " << nTriTotal
+                << " triangles from " << patch.size() << " faces with"
+                << " new patchid " << newPatchi << endl;
+        }
+
+        newPatchi++;
+    }
+    triangles.shrink();
+
+    // Create globally numbered tri surface
+    triSurface rawSurface(triangles, mesh.points());
+
+    // Create locally numbered tri surface
+    triSurface surface
+    (
+        rawSurface.localFaces(),
+        rawSurface.localPoints()
+    );
+
+    // Add patch names to surface
+    surface.patches().setSize(newPatchi);
+
+    newPatchi = 0;
+
+    forAllConstIter(labelHashSet, includePatches, iter)
+    {
+        const label patchi = iter.key();
+        const polyPatch& patch = bMesh[patchi];
+
+        surface.patches()[newPatchi].name() = patch.name();
+        surface.patches()[newPatchi].geometricType() = patch.type();
+
+        newPatchi++;
+    }
+
+    return surface;
+}
+
+
 // triangulation of boundaryMesh
 Foam::triSurface Foam::triSurfaceTools::triangulateFaceCentre
 (
diff --git a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.H b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.H
index 84067be52a..a761186b90 100644
--- a/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.H
+++ b/src/meshTools/triSurface/triSurfaceTools/triSurfaceTools.H
@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2018 OpenFOAM Foundation
      \\/     M anipulation  | Copyright (C) 2016 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
@@ -67,6 +67,7 @@ namespace Foam
 {
 
 // Forward declaration of classes
+class boundBox;
 class edge;
 class labelledTri;
 class polyBoundaryMesh;
@@ -76,7 +77,6 @@ class face;
 class Time;
 template<class Face> class MeshedSurface;
 
-
 /*---------------------------------------------------------------------------*\
                            Class triSurfaceTools Declaration
 \*---------------------------------------------------------------------------*/
@@ -484,6 +484,16 @@ public:
             const bool verbose = false
         );
 
+
+        static triSurface triangulate
+        (
+            const polyBoundaryMesh& bMesh,
+            const labelHashSet& includePatches,
+            const boundBox& bBox,
+            const bool verbose = false
+        );
+
+
         //- Face-centre triangulation of (selected patches of) boundaryMesh.
         //  Needs
         //  polyMesh (or polyBoundaryMesh) since only at this level are the

From 46b768628cec06ac9b851c12741acbc1a86b7a27 Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Wed, 16 May 2018 18:03:14 +0100
Subject: [PATCH 6/8] STYLE: use FieldBase for template invariant parts

- use refCount instead of the tmp<...>::refCount alias
---
 src/OpenFOAM/Make/files                       | 22 +++---
 .../DimensionedField/SubDimensionedFieldI.H   |  2 +-
 .../FieldFields/FieldField/FieldField.C       |  5 +-
 .../FieldFields/FieldField/FieldField.H       |  2 +-
 src/OpenFOAM/fields/Fields/Field/Field.C      |  6 --
 src/OpenFOAM/fields/Fields/Field/Field.H      |  8 +-
 src/OpenFOAM/fields/Fields/Field/FieldBase.C  | 33 ++++++++
 src/OpenFOAM/fields/Fields/Field/FieldBase.H  | 77 +++++++++++++++++++
 src/OpenFOAM/fields/Fields/Field/SubField.H   |  2 +-
 src/OpenFOAM/fields/Fields/Field/SubFieldI.H  |  2 +-
 .../functions/Function1/Function1/Function1.C | 17 ++--
 .../functions/Function1/Function1/Function1.H |  9 +--
 src/finiteArea/faMatrices/faMatrix/faMatrix.H |  2 +-
 .../convectionScheme/convectionScheme.C       |  2 +-
 .../convectionScheme/convectionScheme.H       |  2 +-
 .../d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.C    | 10 ---
 .../d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.H    | 14 ++--
 .../ddtSchemes/ddtScheme/ddtScheme.C          | 18 +----
 .../ddtSchemes/ddtScheme/ddtScheme.H          | 13 ++--
 .../divSchemes/divScheme/divScheme.C          | 10 ---
 .../divSchemes/divScheme/divScheme.H          | 12 +--
 .../gradSchemes/gradScheme/gradScheme.C       |  6 --
 .../gradSchemes/gradScheme/gradScheme.H       | 12 +--
 .../laplacianScheme/laplacianScheme.C         |  7 --
 .../laplacianScheme/laplacianScheme.H         | 12 +--
 .../snGradSchemes/snGradScheme/snGradScheme.C |  7 --
 .../snGradSchemes/snGradScheme/snGradScheme.H | 12 +--
 .../fvMatrices/fvMatrix/fvMatrix.C            |  2 +-
 .../fvMatrices/fvMatrix/fvMatrix.H            |  2 +-
 .../multivariateSurfaceInterpolationScheme.C  | 11 ---
 .../multivariateSurfaceInterpolationScheme.H  | 12 +--
 .../surfaceInterpolationScheme.C              | 10 +--
 .../surfaceInterpolationScheme.H              | 12 +--
 33 files changed, 196 insertions(+), 177 deletions(-)
 create mode 100644 src/OpenFOAM/fields/Fields/Field/FieldBase.C
 create mode 100644 src/OpenFOAM/fields/Fields/Field/FieldBase.H

diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index b0edac8799..96fc8ca422 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -627,36 +627,36 @@ fields/cloud/cloud.C
 
 Fields = fields/Fields
 
+$(Fields)/Field/FieldBase.C
 $(Fields)/labelField/labelField.C
-$(Fields)/scalarField/scalarField.C
-$(Fields)/vectorField/vectorField.C
-$(Fields)/vector2DField/vector2DField.C
-$(Fields)/sphericalTensorField/sphericalTensorField.C
-$(Fields)/diagTensorField/diagTensorField.C
-$(Fields)/symmTensorField/symmTensorField.C
-$(Fields)/tensorField/tensorField.C
-$(Fields)/quaternionField/quaternionField.C
-$(Fields)/triadField/triadField.C
-$(Fields)/complexFields/complexFields.C
-
 $(Fields)/labelField/labelIOField.C
 $(Fields)/labelField/labelFieldIOField.C
+$(Fields)/scalarField/scalarField.C
 $(Fields)/scalarField/scalarIOField.C
 $(Fields)/scalarField/scalarFieldIOField.C
+$(Fields)/vectorField/vectorField.C
 $(Fields)/vectorField/vectorIOField.C
 $(Fields)/vectorField/vectorFieldIOField.C
+$(Fields)/vector2DField/vector2DField.C
 $(Fields)/vector2DField/vector2DIOField.C
 $(Fields)/vector2DField/vector2DFieldIOField.C
+$(Fields)/sphericalTensorField/sphericalTensorField.C
 $(Fields)/sphericalTensorField/sphericalTensorIOField.C
 $(Fields)/sphericalTensorField/sphericalTensorFieldIOField.C
+$(Fields)/diagTensorField/diagTensorField.C
 $(Fields)/diagTensorField/diagTensorIOField.C
 $(Fields)/diagTensorField/diagTensorFieldIOField.C
+$(Fields)/symmTensorField/symmTensorField.C
 $(Fields)/symmTensorField/symmTensorIOField.C
 $(Fields)/symmTensorField/symmTensorFieldIOField.C
+$(Fields)/tensorField/tensorField.C
 $(Fields)/tensorField/tensorIOField.C
 $(Fields)/tensorField/tensorFieldIOField.C
+$(Fields)/quaternionField/quaternionField.C
 $(Fields)/quaternionField/quaternionIOField.C
+$(Fields)/triadField/triadField.C
 $(Fields)/triadField/triadIOField.C
+$(Fields)/complexFields/complexFields.C
 $(Fields)/transformField/transformField.C
 
 
diff --git a/src/OpenFOAM/fields/DimensionedFields/DimensionedField/SubDimensionedFieldI.H b/src/OpenFOAM/fields/DimensionedFields/DimensionedField/SubDimensionedFieldI.H
index 93147ece1d..215a525421 100644
--- a/src/OpenFOAM/fields/DimensionedFields/DimensionedField/SubDimensionedFieldI.H
+++ b/src/OpenFOAM/fields/DimensionedFields/DimensionedField/SubDimensionedFieldI.H
@@ -63,7 +63,7 @@ inline Foam::SubDimensionedField<Type, GeoMesh>::SubDimensionedField
     const SubDimensionedField<Type, GeoMesh>& sfield
 )
 :
-    tmp<SubDimensionedField<Type, GeoMesh>>::refCount(),
+    refCount(),
     SubField<Type>(sfield)
 {}
 
diff --git a/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.C b/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.C
index d9778371b4..ad98f6a02f 100644
--- a/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.C
+++ b/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.C
@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-    Generic fieldField type.
-
 \*---------------------------------------------------------------------------*/
 
 #include "FieldField.H"
@@ -137,7 +134,7 @@ FieldField<Field, Type>::FieldField
 template<template<class> class Field, class Type>
 FieldField<Field, Type>::FieldField(const FieldField<Field, Type>& f)
 :
-    tmp<FieldField<Field, Type>>::refCount(),
+    refCount(),
     PtrList<Field<Type>>(f)
 {}
 
diff --git a/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.H b/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.H
index 19d4ff5cdc..bae0375b68 100644
--- a/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.H
+++ b/src/OpenFOAM/fields/FieldFields/FieldField/FieldField.H
@@ -73,7 +73,7 @@ Ostream& operator<<
 template<template<class> class Field, class Type>
 class FieldField
 :
-    public tmp<FieldField<Field, Type>>::refCount,
+    public refCount,
     public PtrList<Field<Type>>
 {
 
diff --git a/src/OpenFOAM/fields/Fields/Field/Field.C b/src/OpenFOAM/fields/Fields/Field/Field.C
index 7cc7c55e6f..5422b19a13 100644
--- a/src/OpenFOAM/fields/Fields/Field/Field.C
+++ b/src/OpenFOAM/fields/Fields/Field/Field.C
@@ -29,12 +29,6 @@ License
 #include "contiguous.H"
 #include "mapDistributeBase.H"
 
-// * * * * * * * * * * * * * * * Static Members  * * * * * * * * * * * * * * //
-
-template<class Type>
-const char* const Foam::Field<Type>::typeName("Field");
-
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 template<class Type>
diff --git a/src/OpenFOAM/fields/Fields/Field/Field.H b/src/OpenFOAM/fields/Fields/Field/Field.H
index e8a91b0306..ae4f2bfcdb 100644
--- a/src/OpenFOAM/fields/Fields/Field/Field.H
+++ b/src/OpenFOAM/fields/Fields/Field/Field.H
@@ -47,6 +47,7 @@ SourceFiles
 #include "VectorSpace.H"
 #include "scalarList.H"
 #include "labelList.H"
+#include "FieldBase.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -73,7 +74,7 @@ Ostream& operator<<(Ostream&, const tmp<Field<Type>>&);
 template<class Type>
 class Field
 :
-    public tmp<Field<Type>>::refCount,
+    public FieldBase,
     public List<Type>
 {
 
@@ -86,11 +87,6 @@ public:
     typedef SubField<Type> subField;
 
 
-    // Static data members
-
-        static const char* const typeName;
-
-
     // Static Member Functions
 
         //- Return nullObject reference field
diff --git a/src/OpenFOAM/fields/Fields/Field/FieldBase.C b/src/OpenFOAM/fields/Fields/Field/FieldBase.C
new file mode 100644
index 0000000000..91157104f4
--- /dev/null
+++ b/src/OpenFOAM/fields/Fields/Field/FieldBase.C
@@ -0,0 +1,33 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     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 "FieldBase.H"
+
+// * * * * * * * * * * * * * * * Static Members  * * * * * * * * * * * * * * //
+
+const char* const Foam::FieldBase::typeName("Field");
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/fields/Fields/Field/FieldBase.H b/src/OpenFOAM/fields/Fields/Field/FieldBase.H
new file mode 100644
index 0000000000..5de4169ad3
--- /dev/null
+++ b/src/OpenFOAM/fields/Fields/Field/FieldBase.H
@@ -0,0 +1,77 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     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::FieldBase
+
+Description
+    Template invariant parts for Field
+
+SourceFiles
+    FieldBase.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef FieldBase_H
+#define FieldBase_H
+
+#include "refCount.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class FieldBase Declaration
+\*---------------------------------------------------------------------------*/
+
+class FieldBase
+:
+    public refCount
+{
+public:
+
+    // Static data members
+
+        //- Typename for Field
+        static const char* const typeName;
+
+
+    // Constructors
+
+        //- Construct null, with refCount zero
+        constexpr FieldBase() noexcept
+        :
+            refCount()
+        {}
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/fields/Fields/Field/SubField.H b/src/OpenFOAM/fields/Fields/Field/SubField.H
index 07ecd1443f..1a7fb9cb91 100644
--- a/src/OpenFOAM/fields/Fields/Field/SubField.H
+++ b/src/OpenFOAM/fields/Fields/Field/SubField.H
@@ -59,7 +59,7 @@ template<class Type> class SubField;
 template<class Type>
 class SubField
 :
-    public tmp<SubField<Type>>::refCount,
+    public refCount,
     public SubList<Type>
 {
 
diff --git a/src/OpenFOAM/fields/Fields/Field/SubFieldI.H b/src/OpenFOAM/fields/Fields/Field/SubFieldI.H
index 4919f78fe1..d0ed35b5dc 100644
--- a/src/OpenFOAM/fields/Fields/Field/SubFieldI.H
+++ b/src/OpenFOAM/fields/Fields/Field/SubFieldI.H
@@ -74,7 +74,7 @@ inline Foam::SubField<Type>::SubField
     const SubField<Type>& sfield
 )
 :
-    tmp<SubField<Type>>::refCount(),
+    refCount(),
     SubList<Type>(sfield)
 {}
 
diff --git a/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.C b/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.C
index d03e6e8425..0f8cb617c5 100644
--- a/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.C
+++ b/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.C
@@ -38,18 +38,11 @@ Foam::Function1<Type>::Function1(const word& entryName)
 template<class Type>
 Foam::Function1<Type>::Function1(const Function1<Type>& de)
 :
-    tmp<Function1<Type>>::refCount(),
+    refCount(),
     name_(de.name_)
 {}
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-Foam::Function1<Type>::~Function1()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type>
@@ -111,8 +104,8 @@ Foam::FieldFunction1<Function1Type>::value
     const scalarField& x
 ) const
 {
-    tmp<Field<Type>> tfld(new Field<Type>(x.size()));
-    Field<Type>& fld = tfld.ref();
+    auto tfld = tmp<Field<Type>>::New(x.size());
+    auto& fld = tfld.ref();
 
     forAll(x, i)
     {
@@ -152,8 +145,8 @@ Foam::FieldFunction1<Function1Type>::integrate
     const scalarField& x2
 ) const
 {
-    tmp<Field<Type>> tfld(new Field<Type>(x1.size()));
-    Field<Type>& fld = tfld.ref();
+    auto tfld = tmp<Field<Type>>::New(x1.size());
+    auto& fld = tfld.ref();
 
     forAll(x1, i)
     {
diff --git a/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.H b/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.H
index 327a045ec2..0784f12e17 100644
--- a/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.H
+++ b/src/OpenFOAM/primitives/functions/Function1/Function1/Function1.H
@@ -61,11 +61,11 @@ template<class Type> Ostream& operator<<(Ostream&, const Function1<Type>&);
 template<class Type>
 class Function1
 :
-    public tmp<Function1<Type>>::refCount
+    public refCount
 {
     // Private Member Functions
 
-        //- Disallow default bitwise assignment
+        //- No copy assignment
         void operator=(const Function1<Type>&) = delete;
 
 
@@ -119,7 +119,7 @@ public:
 
 
     //- Destructor
-    virtual ~Function1();
+    virtual ~Function1() = default;
 
 
     // Member Functions
@@ -194,8 +194,7 @@ public:
 
 
     //- Destructor
-    virtual ~FieldFunction1()
-    {}
+    virtual ~FieldFunction1() = default;
 
 
     // Member Functions
diff --git a/src/finiteArea/faMatrices/faMatrix/faMatrix.H b/src/finiteArea/faMatrices/faMatrix/faMatrix.H
index 82289005a6..702a366d74 100644
--- a/src/finiteArea/faMatrices/faMatrix/faMatrix.H
+++ b/src/finiteArea/faMatrices/faMatrix/faMatrix.H
@@ -71,7 +71,7 @@ Ostream& operator<<(Ostream&, const faMatrix<Type>&);
 template<class Type>
 class faMatrix
 :
-    public tmp<faMatrix<Type>>::refCount,
+    public refCount,
     public lduMatrix
 {
     // Private data
diff --git a/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.C b/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.C
index ca63afc788..c238481597 100644
--- a/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.C
+++ b/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.C
@@ -46,7 +46,7 @@ namespace fv
 template<class Type>
 convectionScheme<Type>::convectionScheme(const convectionScheme& cs)
 :
-    tmp<convectionScheme<Type>>::refCount(),
+    refCount(),
     mesh_(cs.mesh_)
 {}
 
diff --git a/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.H b/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.H
index aa409f3c26..543ec60ff0 100644
--- a/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.H
+++ b/src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.H
@@ -67,7 +67,7 @@ namespace fv
 template<class Type>
 class convectionScheme
 :
-    public tmp<convectionScheme<Type>>::refCount
+    public refCount
 {
     // Private data
 
diff --git a/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.C b/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.C
index ffdc4c25b8..70208e9790 100644
--- a/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.C
+++ b/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.C
@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-    Abstract base class for finite volume calculus d2dt2 schemes.
-
 \*---------------------------------------------------------------------------*/
 
 #include "fv.H"
@@ -84,13 +81,6 @@ tmp<d2dt2Scheme<Type>> d2dt2Scheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-d2dt2Scheme<Type>::~d2dt2Scheme()
-{}
-
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace fv
diff --git a/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.H b/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.H
index b00888c577..f52dab06ee 100644
--- a/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.H
+++ b/src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.H
@@ -25,7 +25,7 @@ Class
     Foam::fv::d2dt2Scheme
 
 Description
-    Abstract base class for d2dt2 schemes.
+    Abstract base class for finite volume d2dt2 schemes.
 
 SourceFiles
     d2dt2Scheme.C
@@ -64,7 +64,7 @@ namespace fv
 template<class Type>
 class d2dt2Scheme
 :
-    public tmp<d2dt2Scheme<Type>>::refCount
+    public refCount
 {
 
 protected:
@@ -76,11 +76,11 @@ protected:
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        d2dt2Scheme(const d2dt2Scheme&);
+        //- No copy construct
+        d2dt2Scheme(const d2dt2Scheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const d2dt2Scheme&);
+        //- No copy assignment
+        void operator=(const d2dt2Scheme&) = delete;
 
 
 public:
@@ -127,7 +127,7 @@ public:
 
 
     //- Destructor
-    virtual ~d2dt2Scheme();
+    virtual ~d2dt2Scheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
index fd4f502ed9..a7fa362d72 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.C
@@ -84,13 +84,6 @@ tmp<ddtScheme<Type>> ddtScheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-ddtScheme<Type>::~ddtScheme()
-{}
-
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 template<class Type>
@@ -120,14 +113,11 @@ tmp<fvMatrix<Type>> ddtScheme<Type>::fvmDdt
 {
     NotImplemented;
 
-    return tmp<fvMatrix<Type>>
+    return tmp<fvMatrix<Type>>::New
     (
-        new fvMatrix<Type>
-        (
-            vf,
-            alpha.dimensions()*rho.dimensions()
-            *vf.dimensions()*dimVol/dimTime
-        )
+        vf,
+        alpha.dimensions()*rho.dimensions()
+        *vf.dimensions()*dimVol/dimTime
     );
 }
 
diff --git a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H
index 8a1dfbaa8d..55fcc68ca1 100644
--- a/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H
+++ b/src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H
@@ -30,7 +30,6 @@ Group
 Description
     Abstract base class for ddt schemes.
 
-
 SourceFiles
     ddtScheme.C
 
@@ -68,7 +67,7 @@ namespace fv
 template<class Type>
 class ddtScheme
 :
-    public tmp<ddtScheme<Type>>::refCount
+    public refCount
 {
 
 protected:
@@ -83,11 +82,11 @@ protected:
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        ddtScheme(const ddtScheme&);
+        //- No copy construct
+        ddtScheme(const ddtScheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const ddtScheme&);
+        //- No copy assignment
+        void operator=(const ddtScheme&) = delete;
 
 
 public:
@@ -136,7 +135,7 @@ public:
 
 
     //- Destructor
-    virtual ~ddtScheme();
+    virtual ~ddtScheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.C b/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.C
index 8c71aafcb3..66b7f83df2 100644
--- a/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.C
+++ b/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.C
@@ -85,16 +85,6 @@ tmp<divScheme<Type>> divScheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-divScheme<Type>::~divScheme()
-{}
-
-
-// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
-
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace fv
diff --git a/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.H b/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.H
index 1f93d225c9..3b5e948a0e 100644
--- a/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.H
+++ b/src/finiteVolume/finiteVolume/divSchemes/divScheme/divScheme.H
@@ -67,7 +67,7 @@ namespace fv
 template<class Type>
 class divScheme
 :
-    public tmp<divScheme<Type>>::refCount
+    public refCount
 {
 
 protected:
@@ -80,11 +80,11 @@ protected:
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        divScheme(const divScheme&);
+        //- No copy construct
+        divScheme(const divScheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const divScheme&);
+        //- No copy assignment
+        void operator=(const divScheme&) = delete;
 
 
 public:
@@ -133,7 +133,7 @@ public:
 
 
     //- Destructor
-    virtual ~divScheme();
+    virtual ~divScheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.C b/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.C
index 1a055779f1..414fe7b8ce 100644
--- a/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.C
+++ b/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.C
@@ -72,12 +72,6 @@ Foam::tmp<Foam::fv::gradScheme<Type>> Foam::fv::gradScheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-Foam::fv::gradScheme<Type>::~gradScheme()
-{}
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 template<class Type>
diff --git a/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.H b/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.H
index 4e6ca98bf6..1ed984b213 100644
--- a/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.H
+++ b/src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.H
@@ -60,7 +60,7 @@ namespace fv
 template<class Type>
 class gradScheme
 :
-    public tmp<gradScheme<Type>>::refCount
+    public refCount
 {
     // Private data
 
@@ -69,11 +69,11 @@ class gradScheme
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        gradScheme(const gradScheme&);
+        //- No copy construct
+        gradScheme(const gradScheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const gradScheme&);
+        //- No copy assignment
+        void operator=(const gradScheme&) = delete;
 
 
 public:
@@ -114,7 +114,7 @@ public:
 
 
     //- Destructor
-    virtual ~gradScheme();
+    virtual ~gradScheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.C b/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.C
index 0a774bdf5d..fefd505876 100644
--- a/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.C
+++ b/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.C
@@ -83,13 +83,6 @@ tmp<laplacianScheme<Type, GType>> laplacianScheme<Type, GType>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type, class GType>
-laplacianScheme<Type, GType>::~laplacianScheme()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type, class GType>
diff --git a/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.H b/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.H
index cb9e02ab38..ea2d07bc95 100644
--- a/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.H
+++ b/src/finiteVolume/finiteVolume/laplacianSchemes/laplacianScheme/laplacianScheme.H
@@ -68,7 +68,7 @@ namespace fv
 template<class Type, class GType>
 class laplacianScheme
 :
-    public tmp<laplacianScheme<Type, GType>>::refCount
+    public refCount
 {
 
 protected:
@@ -84,11 +84,11 @@ private:
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        laplacianScheme(const laplacianScheme&);
+        //- No copy construct
+        laplacianScheme(const laplacianScheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const laplacianScheme&);
+        //- No copy assignment
+        void operator=(const laplacianScheme&) = delete;
 
 
 public:
@@ -162,7 +162,7 @@ public:
 
 
     //- Destructor
-    virtual ~laplacianScheme();
+    virtual ~laplacianScheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.C b/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.C
index d5422e729d..245a014065 100644
--- a/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.C
+++ b/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.C
@@ -85,13 +85,6 @@ tmp<snGradScheme<Type>> snGradScheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-snGradScheme<Type>::~snGradScheme()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type>
diff --git a/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.H b/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.H
index 14821c17c2..5c9ad02d4c 100644
--- a/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.H
+++ b/src/finiteVolume/finiteVolume/snGradSchemes/snGradScheme/snGradScheme.H
@@ -63,7 +63,7 @@ namespace fv
 template<class Type>
 class snGradScheme
 :
-    public tmp<snGradScheme<Type>>::refCount
+    public refCount
 {
     // Private data
 
@@ -73,11 +73,11 @@ class snGradScheme
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        snGradScheme(const snGradScheme&);
+        //- No copy construct
+        snGradScheme(const snGradScheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const snGradScheme&);
+        //- No copy assignment
+        void operator=(const snGradScheme&) = delete;
 
 
 public:
@@ -118,7 +118,7 @@ public:
 
 
     //- Destructor
-    virtual ~snGradScheme();
+    virtual ~snGradScheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
index a312cd97f9..8a18e7446a 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.C
@@ -321,7 +321,7 @@ Foam::fvMatrix<Type>::fvMatrix
 template<class Type>
 Foam::fvMatrix<Type>::fvMatrix(const fvMatrix<Type>& fvm)
 :
-    tmp<fvMatrix<Type>>::refCount(),
+    refCount(),
     lduMatrix(fvm),
     psi_(fvm.psi_),
     dimensions_(fvm.dimensions_),
diff --git a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
index 4cad056ab2..1acbb4300e 100644
--- a/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
+++ b/src/finiteVolume/fvMatrices/fvMatrix/fvMatrix.H
@@ -114,7 +114,7 @@ template<class T> class UIndirectList;
 template<class Type>
 class fvMatrix
 :
-    public tmp<fvMatrix<Type>>::refCount,
+    public refCount,
     public lduMatrix
 {
     // Private data
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.C b/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.C
index e9400fbfee..5a4a26285f 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.C
@@ -21,9 +21,6 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
 
-Description
-    Abstract base class for surface interpolation schemes.
-
 \*---------------------------------------------------------------------------*/
 
 #include "fv.H"
@@ -85,12 +82,4 @@ Foam::multivariateSurfaceInterpolationScheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-Foam::multivariateSurfaceInterpolationScheme<Type>::
-~multivariateSurfaceInterpolationScheme()
-{}
-
-
 // ************************************************************************* //
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.H b/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.H
index 4669b43d63..a3c20cc7dd 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.H
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/multivariateSchemes/multivariateSurfaceInterpolationScheme/multivariateSurfaceInterpolationScheme.H
@@ -50,7 +50,7 @@ namespace Foam
 template<class Type>
 class multivariateSurfaceInterpolationScheme
 :
-    public tmp<multivariateSurfaceInterpolationScheme<Type>>::refCount
+    public refCount
 {
 
 public:
@@ -85,14 +85,14 @@ private:
 
     // Private Member Functions
 
-        //- Disallow default bitwise copy construct
+        //- No copy construct
         multivariateSurfaceInterpolationScheme
         (
             const multivariateSurfaceInterpolationScheme&
-        );
+        ) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const multivariateSurfaceInterpolationScheme&);
+        //- No copy assignment
+        void operator=(const multivariateSurfaceInterpolationScheme&) = delete;
 
 
 public:
@@ -143,7 +143,7 @@ public:
 
 
     //- Destructor
-    virtual ~multivariateSurfaceInterpolationScheme();
+    virtual ~multivariateSurfaceInterpolationScheme() = default;
 
 
     // Member Functions
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
index b35c82f8a7..8e14179101 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.C
@@ -123,13 +123,6 @@ Foam::surfaceInterpolationScheme<Type>::New
 }
 
 
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-template<class Type>
-Foam::surfaceInterpolationScheme<Type>::~surfaceInterpolationScheme()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type>
@@ -147,8 +140,7 @@ Foam::surfaceInterpolationScheme<Type>::interpolate
             << "Interpolating "
             << vf.type() << " "
             << vf.name()
-            << " from cells to faces "
-               "without explicit correction"
+            << " from cells to faces without explicit correction"
             << endl;
     }
 
diff --git a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H
index bf233ef3f9..871578b626 100644
--- a/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H
+++ b/src/finiteVolume/interpolation/surfaceInterpolation/surfaceInterpolationScheme/surfaceInterpolationScheme.H
@@ -55,7 +55,7 @@ class fvMesh;
 template<class Type>
 class surfaceInterpolationScheme
 :
-    public tmp<surfaceInterpolationScheme<Type>>::refCount
+    public refCount
 {
     // Private data
 
@@ -65,11 +65,11 @@ class surfaceInterpolationScheme
 
     // Private Member Functions
 
-        //- Disallow copy construct
-        surfaceInterpolationScheme(const surfaceInterpolationScheme&);
+        //- No copy construct
+        surfaceInterpolationScheme(const surfaceInterpolationScheme&) = delete;
 
-        //- Disallow default bitwise assignment
-        void operator=(const surfaceInterpolationScheme&);
+        //- No copy assignment
+        void operator=(const surfaceInterpolationScheme&) = delete;
 
 
 public:
@@ -134,7 +134,7 @@ public:
 
 
     //- Destructor
-    virtual ~surfaceInterpolationScheme();
+    virtual ~surfaceInterpolationScheme() = default;
 
 
     // Member Functions

From 48d654cf19a1aa57009918c1865caa374a9c5e4f Mon Sep 17 00:00:00 2001
From: Mark Olesen <Mark.Olesen@esi-group.com>
Date: Thu, 17 May 2018 09:56:36 +0100
Subject: [PATCH 7/8] ENH: avoid memory leaks for HashPtrTable, PtrMap
 insertion (issue #749)

- disallow insert() of raw pointers, since a failed insertion
  (ie, entry already existed) results in an unmanaged pointer.

  Either insert using an autoPtr, or set() with raw pointers or autoPtr.

- IOobjectList::add() now takes an autoPtr instead of an object reference

- IOobjectList::remove() now returns an autoPtr instead of a raw pointer
---
 .../multiphaseSystem/multiphaseSystem.C       |   6 +-
 .../HeatAndMassTransferPhaseSystem.C          |  20 +--
 .../HeatTransferPhaseSystem.C                 |  14 +-
 ...terfaceCompositionPhaseChangePhaseSystem.C |  14 +-
 .../MomentumTransferPhaseSystem.C             |  33 ++---
 .../ThermalPhaseChangePhaseSystem.C           |  16 +--
 applications/test/HashPtrTable/Make/files     |   4 +-
 ...est-hashPtrTable.C => Test-HashPtrTable.C} |  34 ++++-
 applications/test/PtrMap/Make/files           |   3 +
 applications/test/PtrMap/Make/options         |   0
 applications/test/PtrMap/Test-PtrMap.C        | 130 ++++++++++++++++++
 .../utilities/mesh/advanced/PDRMesh/PDRMesh.C |  41 +++---
 modules/catalyst                              |   2 +-
 .../HashTables/HashPtrTable/HashPtrTable.C    |  61 ++++----
 .../HashTables/HashPtrTable/HashPtrTable.H    |  60 ++++++--
 .../HashTables/HashPtrTable/HashPtrTableIO.C  |   8 +-
 .../HashTables/HashTable/HashTable.C          |   2 +-
 .../HashTables/HashTable/HashTable.H          |  40 +++---
 .../HashTables/HashTable/HashTableI.H         |   4 +-
 .../containers/HashTables/PtrMap/PtrMap.H     |  25 +++-
 src/OpenFOAM/db/IOobjectList/IOobjectList.C   |  82 +++++++----
 src/OpenFOAM/db/IOobjectList/IOobjectList.H   |  47 ++++---
 .../masterUncollatedFileOperation.C           |   2 +-
 src/OpenFOAM/graph/graph.C                    |   8 +-
 .../graphics/runTimePostProcessing/scene.C    |   7 +-
 .../genericFaPatchField/genericFaPatchField.C |  53 ++++---
 .../genericFvPatchField/genericFvPatchField.C |  53 ++++---
 .../genericPointPatchField.C                  |  27 ++--
 src/sampling/probes/probes.C                  |  12 +-
 .../chemkinReader/chemkinLexer.L              |   2 +-
 30 files changed, 549 insertions(+), 261 deletions(-)
 rename applications/test/HashPtrTable/{Test-hashPtrTable.C => Test-HashPtrTable.C} (77%)
 create mode 100644 applications/test/PtrMap/Make/files
 create mode 100644 applications/test/PtrMap/Make/options
 create mode 100644 applications/test/PtrMap/Test-PtrMap.C

diff --git a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C
index 3d51ee57cd..2815e56690 100644
--- a/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C
+++ b/applications/solvers/multiphase/multiphaseEulerFoam/multiphaseSystem/multiphaseSystem.C
@@ -426,7 +426,7 @@ Foam::multiphaseSystem::multiphaseSystem
                 iter(),
                 *phases_.lookup(iter.key().first()),
                 *phases_.lookup(iter.key().second())
-            ).ptr()
+            )
         );
     }
 
@@ -664,7 +664,7 @@ Foam::tmp<Foam::volVectorField> Foam::multiphaseSystem::Svm
 Foam::autoPtr<Foam::multiphaseSystem::dragCoeffFields>
 Foam::multiphaseSystem::dragCoeffs() const
 {
-    autoPtr<dragCoeffFields> dragCoeffsPtr(new dragCoeffFields);
+    auto dragCoeffsPtr = autoPtr<dragCoeffFields>::New();
 
     forAllConstIter(dragModelTable, dragModels_, iter)
     {
@@ -706,7 +706,7 @@ Foam::multiphaseSystem::dragCoeffs() const
             }
         }
 
-        dragCoeffsPtr().insert(iter.key(), Kptr);
+        dragCoeffsPtr().set(iter.key(), Kptr);
     }
 
     return dragCoeffsPtr;
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C
index 16d01d3a71..b985ef2ad0 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatAndMassTransferPhaseSystem/HeatAndMassTransferPhaseSystem.C
@@ -70,7 +70,7 @@ HeatAndMassTransferPhaseSystem
 
         // Initially assume no mass transfer
 
-        dmdt_.insert
+        dmdt_.set
         (
             pair,
             new volScalarField
@@ -88,7 +88,7 @@ HeatAndMassTransferPhaseSystem
             )
         );
 
-        dmdtExplicit_.insert
+        dmdtExplicit_.set
         (
             pair,
             new volScalarField
@@ -107,7 +107,7 @@ HeatAndMassTransferPhaseSystem
         volScalarField H1(heatTransferModels_[pair][pair.first()]->K());
         volScalarField H2(heatTransferModels_[pair][pair.second()]->K());
 
-        Tf_.insert
+        Tf_.set
         (
             pair,
             new volScalarField
@@ -257,18 +257,12 @@ template<class BasePhaseSystem>
 Foam::autoPtr<Foam::phaseSystem::heatTransferTable>
 Foam::HeatAndMassTransferPhaseSystem<BasePhaseSystem>::heatTransfer() const
 {
-    autoPtr<phaseSystem::heatTransferTable> eqnsPtr
-    (
-        new phaseSystem::heatTransferTable()
-    );
+    auto eqnsPtr = autoPtr<phaseSystem::heatTransferTable>::New();
+    auto& eqns = *eqnsPtr;
 
-    phaseSystem::heatTransferTable& eqns = eqnsPtr();
-
-    forAll(this->phaseModels_, phasei)
+    for (const phaseModel& phase : this->phaseModels_)
     {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
-        eqns.insert
+        eqns.set
         (
             phase.name(),
             new fvScalarMatrix(phase.thermo().he(), dimEnergy/dimTime)
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C
index 2654fce6f9..7f8b768c94 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/HeatTransferPhaseSystem/HeatTransferPhaseSystem.C
@@ -107,18 +107,12 @@ template<class BasePhaseSystem>
 Foam::autoPtr<Foam::phaseSystem::heatTransferTable>
 Foam::HeatTransferPhaseSystem<BasePhaseSystem>::heatTransfer() const
 {
-    autoPtr<phaseSystem::heatTransferTable> eqnsPtr
-    (
-        new phaseSystem::heatTransferTable()
-    );
+    auto eqnsPtr = autoPtr<phaseSystem::heatTransferTable>::New();
+    auto& eqns = *eqnsPtr;
 
-    phaseSystem::heatTransferTable& eqns = eqnsPtr();
-
-    forAll(this->phaseModels_, phasei)
+    for (const phaseModel& phase : this->phaseModels_)
     {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
-        eqns.insert
+        eqns.set
         (
             phase.name(),
             new fvScalarMatrix(phase.thermo().he(), dimEnergy/dimTime)
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C
index 54fec2f8ec..a68bd91745 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/InterfaceCompositionPhaseChangePhaseSystem/InterfaceCompositionPhaseChangePhaseSystem.C
@@ -62,22 +62,16 @@ Foam::InterfaceCompositionPhaseChangePhaseSystem<BasePhaseSystem>::
 massTransfer() const
 {
     // Create a mass transfer matrix for each species of each phase
-    autoPtr<phaseSystem::massTransferTable> eqnsPtr
-    (
-        new phaseSystem::massTransferTable()
-    );
+    auto eqnsPtr = autoPtr<phaseSystem::massTransferTable>::New();
+    auto& eqns = *eqnsPtr;
 
-    phaseSystem::massTransferTable& eqns = eqnsPtr();
-
-    forAll(this->phaseModels_, phasei)
+    for (const phaseModel& phase : this->phaseModels_)
     {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
         const PtrList<volScalarField>& Yi = phase.Y();
 
         forAll(Yi, i)
         {
-            eqns.insert
+            eqns.set
             (
                 Yi[i].name(),
                 new fvScalarMatrix(Yi[i], dimMass/dimTime)
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C
index 812fbce2ac..619dbd1631 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/MomentumTransferPhaseSystem/MomentumTransferPhaseSystem.C
@@ -87,7 +87,7 @@ MomentumTransferPhaseSystem
         const phasePair& pair =
             *(this->phasePairs_[dragModelIter.key()]);
 
-        Kds_.insert
+        Kds_.set
         (
             pair,
             new volScalarField
@@ -103,7 +103,7 @@ MomentumTransferPhaseSystem
         const phasePair& pair =
             *(this->phasePairs_[virtualMassModelIter.key()]);
 
-        Vms_.insert
+        Vms_.set
         (
             pair,
             new volScalarField
@@ -373,18 +373,12 @@ Foam::autoPtr<Foam::phaseSystem::momentumTransferTable>
 Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::momentumTransfer() const
 {
     // Create a momentum transfer matrix for each phase
-    autoPtr<phaseSystem::momentumTransferTable> eqnsPtr
-    (
-        new phaseSystem::momentumTransferTable()
-    );
+    auto eqnsPtr = autoPtr<phaseSystem::momentumTransferTable>::New();
+    auto& eqns = *eqnsPtr;
 
-    phaseSystem::momentumTransferTable& eqns = eqnsPtr();
-
-    forAll(this->phaseModels_, phasei)
+    for (const phaseModel& phase : this->phaseModels_)
     {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
-        eqns.insert
+        eqns.set
         (
             phase.name(),
             new fvVectorMatrix(phase.U(), dimMass*dimVelocity/dimTime)
@@ -492,11 +486,8 @@ template<class BasePhaseSystem>
 Foam::autoPtr<Foam::PtrList<Foam::volVectorField>>
 Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::Fs() const
 {
-    autoPtr<PtrList<volVectorField>> tFs
-    (
-        new PtrList<volVectorField>(this->phases().size())
-    );
-    PtrList<volVectorField>& Fs = tFs();
+    auto tFs = autoPtr<PtrList<volVectorField>>::New(this->phases().size());
+    auto& Fs = *tFs;
 
     // Add the lift force
     forAllConstIters(liftModels_, modelIter)
@@ -570,11 +561,9 @@ Foam::MomentumTransferPhaseSystem<BasePhaseSystem>::phiDs
     const PtrList<volScalarField>& rAUs
 ) const
 {
-    autoPtr<PtrList<surfaceScalarField>> tphiDs
-    (
-        new PtrList<surfaceScalarField>(this->phases().size())
-    );
-    PtrList<surfaceScalarField>& phiDs = tphiDs();
+    auto tphiDs =
+        autoPtr<PtrList<surfaceScalarField>>::New(this->phases().size());
+    auto& phiDs = *tphiDs;
 
     // Add the turbulent dispersion force
     forAllConstIters(turbulentDispersionModels_, turbulentDispersionModelIter)
diff --git a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C
index 9ffa4a882f..d63ca7d058 100644
--- a/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C
+++ b/applications/solvers/multiphase/reactingEulerFoam/phaseSystems/PhaseSystems/ThermalPhaseChangePhaseSystem/ThermalPhaseChangePhaseSystem.C
@@ -52,7 +52,7 @@ ThermalPhaseChangePhaseSystem
         }
 
         // Initially assume no mass transfer
-        iDmdt_.insert
+        iDmdt_.set
         (
             pair,
             new volScalarField
@@ -186,22 +186,16 @@ Foam::autoPtr<Foam::phaseSystem::massTransferTable>
 Foam::ThermalPhaseChangePhaseSystem<BasePhaseSystem>::massTransfer() const
 {
     // Create a mass transfer matrix for each species of each phase
-    autoPtr<phaseSystem::massTransferTable> eqnsPtr
-    (
-        new phaseSystem::massTransferTable()
-    );
+    auto eqnsPtr = autoPtr<phaseSystem::massTransferTable>::New();
+    auto& eqns = *eqnsPtr;
 
-    phaseSystem::massTransferTable& eqns = eqnsPtr();
-
-    forAll(this->phaseModels_, phasei)
+    for (const phaseModel& phase : this->phaseModels_)
     {
-        const phaseModel& phase = this->phaseModels_[phasei];
-
         const PtrList<volScalarField>& Yi = phase.Y();
 
         forAll(Yi, i)
         {
-            eqns.insert
+            eqns.set
             (
                 Yi[i].name(),
                 new fvScalarMatrix(Yi[i], dimMass/dimTime)
diff --git a/applications/test/HashPtrTable/Make/files b/applications/test/HashPtrTable/Make/files
index 239bab4cd7..00d71a8230 100644
--- a/applications/test/HashPtrTable/Make/files
+++ b/applications/test/HashPtrTable/Make/files
@@ -1,3 +1,3 @@
-Test-hashPtrTable.C
+Test-HashPtrTable.C
 
-EXE = $(FOAM_USER_APPBIN)/Test-hashPtrTable
+EXE = $(FOAM_USER_APPBIN)/Test-HashPtrTable
diff --git a/applications/test/HashPtrTable/Test-hashPtrTable.C b/applications/test/HashPtrTable/Test-HashPtrTable.C
similarity index 77%
rename from applications/test/HashPtrTable/Test-hashPtrTable.C
rename to applications/test/HashPtrTable/Test-HashPtrTable.C
index e62d596869..79ef9f3d39 100644
--- a/applications/test/HashPtrTable/Test-hashPtrTable.C
+++ b/applications/test/HashPtrTable/Test-HashPtrTable.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -27,6 +27,7 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include <iostream>
+#include "autoPtr.H"
 #include "HashPtrTable.H"
 
 using namespace Foam;
@@ -42,7 +43,7 @@ void printTable(const HashPtrTable<T>& table)
         Info<< iter.key() << " = ";
         if (ptr)
         {
-            Info<< *ptr;
+            Info<< *ptr << " (" << long(ptr) << ")";
         }
         else
         {
@@ -77,11 +78,11 @@ void printTable(const HashPtrTable<T>& table)
 int main()
 {
     HashPtrTable<double> myTable;
-    myTable.insert("abc", new double(42.1));
-    myTable.insert("def", nullptr);
-    myTable.insert("pi", new double(3.14159));
-    myTable.insert("natlog", new double(2.718282));
-    myTable.insert("sqrt2", new double(1.414214));
+    myTable.set("abc", new double(42.1));
+    myTable.set("def", nullptr);
+    myTable.set("pi", new double(3.14159));
+    myTable.set("natlog", new double(2.718282));
+    myTable.insert("sqrt2", autoPtr<double>::New(1.414214));
 
     // Info<< myTable << endl;
     printTable(myTable);
@@ -105,6 +106,25 @@ int main()
 
     printTable(myTable);
 
+    HashPtrTable<double> moved(std::move(copy));
+
+    Info<< nl << "test movable" << nl;
+    Info<<"input:" << nl;
+    printTable(copy);
+
+    Info<<"output:" << nl;
+    printTable(moved);
+
+    HashPtrTable<double> other;
+
+    Info<<"move assign" << nl;
+
+    other = std::move(moved);
+    printTable(other);
+
+    Info<<"old" << nl;
+    printTable(moved);
+
     return 0;
 }
 
diff --git a/applications/test/PtrMap/Make/files b/applications/test/PtrMap/Make/files
new file mode 100644
index 0000000000..7b281cd0fd
--- /dev/null
+++ b/applications/test/PtrMap/Make/files
@@ -0,0 +1,3 @@
+Test-PtrMap.C
+
+EXE = $(FOAM_USER_APPBIN)/Test-PtrMap
diff --git a/applications/test/PtrMap/Make/options b/applications/test/PtrMap/Make/options
new file mode 100644
index 0000000000..e69de29bb2
diff --git a/applications/test/PtrMap/Test-PtrMap.C b/applications/test/PtrMap/Test-PtrMap.C
new file mode 100644
index 0000000000..6d5866bc6d
--- /dev/null
+++ b/applications/test/PtrMap/Test-PtrMap.C
@@ -0,0 +1,130 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2018 OpenCFD Ltd.
+     \\/     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/>.
+
+Description
+
+
+\*---------------------------------------------------------------------------*/
+
+#include <iostream>
+#include "PtrMap.H"
+
+using namespace Foam;
+
+template<class T>
+void printTable(const PtrMap<T>& table)
+{
+    Info<< table.size() << nl << "(" << nl;
+
+    forAllConstIters(table, iter)
+    {
+        const T* ptr = iter.object();
+        Info<< iter.key() << " = ";
+        if (ptr)
+        {
+            Info<< *ptr << " (" << long(ptr) << ")";
+        }
+        else
+        {
+            Info<< "nullptr";
+        }
+        Info<< nl;
+    }
+
+    Info<< ")" << endl;
+
+    // Values only, with for-range
+    Info<< "values (";
+    for (auto val : table)
+    {
+        Info<< ' ';
+        if (val)
+        {
+            Info<< *val;
+        }
+        else
+        {
+            Info<< "nullptr";
+        }
+    }
+    Info<< " )" << nl;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+//  Main program:
+
+int main()
+{
+    PtrMap<double> myTable;
+    myTable.set(1, new double(42.1));
+    myTable.set(2, nullptr);
+    myTable.set(3, new double(3.14159));
+    myTable.set(4, new double(2.718282));
+    myTable.set(4, new double(1.414214));
+
+    // Info<< myTable << endl;
+    printTable(myTable);
+
+    PtrMap<double> copy(myTable);
+
+    // Info<< copy << endl;
+    printTable(copy);
+    Info<< copy << endl;
+
+    Info<<"\nerase some existing and non-existing entries" << nl;
+
+    auto iter = myTable.find(3);
+    myTable.erase(iter);
+
+    iter = myTable.find(1000);  // unknown key
+    myTable.erase(iter);
+
+    myTable.erase(1);
+    iter = myTable.find(100000);  // unknown key
+
+    printTable(myTable);
+
+    PtrMap<double> moved(std::move(copy));
+
+    Info<< nl << "test movable" << nl;
+    Info<<"input:" << nl;
+    printTable(copy);
+
+    Info<<"output:" << nl;
+    printTable(moved);
+
+    PtrMap<double> other;
+
+    Info<<"move assign" << nl;
+
+    other = std::move(moved);
+    printTable(other);
+
+    Info<<"old" << nl;
+    printTable(moved);
+
+    return 0;
+}
+
+// ************************************************************************* //
diff --git a/applications/utilities/mesh/advanced/PDRMesh/PDRMesh.C b/applications/utilities/mesh/advanced/PDRMesh/PDRMesh.C
index 62d6a8cc4c..5fd3b34781 100644
--- a/applications/utilities/mesh/advanced/PDRMesh/PDRMesh.C
+++ b/applications/utilities/mesh/advanced/PDRMesh/PDRMesh.C
@@ -813,28 +813,33 @@ int main(int argc, char *argv[])
 
     // Read all fields in time and constant directories
     IOobjectList objects(mesh, runTime.timeName());
-    IOobjectList timeObjects(IOobjectList(mesh, mesh.facesInstance()));
-    forAllConstIter(IOobjectList, timeObjects, iter)
     {
-        if
-        (
-            iter()->headerClassName() == volScalarField::typeName
-         || iter()->headerClassName() == volVectorField::typeName
-         || iter()->headerClassName() == volSphericalTensorField::typeName
-         || iter()->headerClassName() == volTensorField::typeName
-         || iter()->headerClassName() == volSymmTensorField::typeName
-         || iter()->headerClassName() == surfaceScalarField::typeName
-         || iter()->headerClassName() == surfaceVectorField::typeName
-         || iter()->headerClassName()
-            == surfaceSphericalTensorField::typeName
-         || iter()->headerClassName() == surfaceSymmTensorField::typeName
-         || iter()->headerClassName() == surfaceTensorField::typeName
-        )
+        IOobjectList timeObjects(mesh, mesh.facesInstance());
+
+        // Transfer specific types
+        forAllIters(timeObjects, iter)
         {
-            objects.add(*iter());
+            autoPtr<IOobject> objPtr(timeObjects.remove(iter));
+            const auto& obj = *objPtr;
+
+            if
+            (
+                obj.headerClassName() == volScalarField::typeName
+             || obj.headerClassName() == volVectorField::typeName
+             || obj.headerClassName() == volSphericalTensorField::typeName
+             || obj.headerClassName() == volTensorField::typeName
+             || obj.headerClassName() == volSymmTensorField::typeName
+             || obj.headerClassName() == surfaceScalarField::typeName
+             || obj.headerClassName() == surfaceVectorField::typeName
+             || obj.headerClassName() == surfaceSphericalTensorField::typeName
+             || obj.headerClassName() == surfaceSymmTensorField::typeName
+             || obj.headerClassName() == surfaceTensorField::typeName
+            )
+            {
+                objects.add(objPtr);
+            }
         }
     }
-
     // Read vol fields and subset.
 
     wordList scalarNames(objects.names(volScalarField::typeName));
diff --git a/modules/catalyst b/modules/catalyst
index 5828d45108..d31e13e12b 160000
--- a/modules/catalyst
+++ b/modules/catalyst
@@ -1 +1 @@
-Subproject commit 5828d4510816948b034aa9afdf0b7b05659a9f59
+Subproject commit d31e13e12b4d9b727b28d39bc34086f2cb75acfa
diff --git a/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.C b/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.C
index e7ad07ec2f..57c4b9bd37 100644
--- a/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.C
+++ b/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -28,20 +28,6 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-template<class T, class Key, class Hash>
-Foam::HashPtrTable<T, Key, Hash>::HashPtrTable()
-:
-    parent_type()
-{}
-
-
-template<class T, class Key, class Hash>
-Foam::HashPtrTable<T, Key, Hash>::HashPtrTable(const label size)
-:
-    parent_type(size)
-{}
-
-
 template<class T, class Key, class Hash>
 Foam::HashPtrTable<T, Key, Hash>::HashPtrTable
 (
@@ -55,16 +41,26 @@ Foam::HashPtrTable<T, Key, Hash>::HashPtrTable
         const T* ptr = iter.object();
         if (ptr)
         {
-            this->insert(iter.key(), new T(*ptr));
+            this->set(iter.key(), new T(*ptr));
         }
         else
         {
-            this->insert(iter.key(), nullptr);
+            this->set(iter.key(), nullptr);
         }
     }
 }
 
 
+template<class T, class Key, class Hash>
+Foam::HashPtrTable<T, Key, Hash>::HashPtrTable
+(
+    HashPtrTable<T, Key, Hash>&& ht
+)
+:
+    parent_type(std::move(ht))
+{}
+
+
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 template<class T, class Key, class Hash>
@@ -77,13 +73,13 @@ Foam::HashPtrTable<T, Key, Hash>::~HashPtrTable()
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class T, class Key, class Hash>
-T* Foam::HashPtrTable<T, Key, Hash>::remove(iterator& iter)
+Foam::autoPtr<T> Foam::HashPtrTable<T, Key, Hash>::remove(iterator& iter)
 {
     if (iter.found())
     {
-        T* ptr = iter.object();
+        autoPtr<T> aptr(iter.object());
         this->parent_type::erase(iter);
-        return ptr;
+        return aptr;
     }
 
     return nullptr;
@@ -140,11 +136,10 @@ void Foam::HashPtrTable<T, Key, Hash>::operator=
     const HashPtrTable<T, Key, Hash>& rhs
 )
 {
-    // Check for assignment to self
     if (this == &rhs)
     {
         FatalErrorInFunction
-            << "attempted assignment to self"
+            << "attempted copy assignment to self"
             << abort(FatalError);
     }
 
@@ -155,16 +150,34 @@ void Foam::HashPtrTable<T, Key, Hash>::operator=
         const T* ptr = iter.object();
         if (ptr)
         {
-            this->insert(iter.key(), new T(*ptr));
+            this->set(iter.key(), new T(*ptr));
         }
         else
         {
-            this->insert(iter.key(), nullptr);
+            this->set(iter.key(), nullptr);
         }
     }
 }
 
 
+template<class T, class Key, class Hash>
+void Foam::HashPtrTable<T, Key, Hash>::operator=
+(
+    HashPtrTable<T, Key, Hash>&& rhs
+)
+{
+    if (this == &rhs)
+    {
+        FatalErrorInFunction
+            << "attempted move assignment to self"
+            << abort(FatalError);
+    }
+
+    this->clear();
+    this->transfer(rhs);
+}
+
+
 // * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
 
 #include "HashPtrTableIO.C"
diff --git a/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.H b/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.H
index b705c27c9f..9dedb2a0ff 100644
--- a/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.H
+++ b/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTable.H
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -43,11 +43,12 @@ SourceFiles
 namespace Foam
 {
 
+// Forward declarations
+
 class Istream;
 class Ostream;
 
-// Forward declaration of friend functions and operators
-
+template<class T> class autoPtr;
 template<class T, class Key, class Hash> class HashPtrTable;
 
 template<class T, class Key, class Hash>
@@ -92,10 +93,10 @@ public:
     // Constructors
 
         //- Construct null with default table capacity
-        HashPtrTable();
+        inline HashPtrTable();
 
         //- Construct given initial table capacity
-        explicit HashPtrTable(const label size);
+        explicit inline HashPtrTable(const label size);
 
         //- Construct from Istream using given Istream constructor class
         template<class INew>
@@ -107,9 +108,12 @@ public:
         //- Construct from dictionary with default dictionary constructor class
         explicit HashPtrTable(const dictionary& dict);
 
-        //- Construct as copy
+        //- Copy construct
         HashPtrTable(const this_type& ht);
 
+        //- Move construct
+        HashPtrTable(this_type&& ht);
+
 
     //- Destructor
     ~HashPtrTable();
@@ -117,20 +121,22 @@ public:
 
     // Member Functions
 
-      // Edit
+    // Edit
 
         //- Remove and return the pointer specified by given iterator.
         //  Includes a safeguard against the end-iterator.
-        T* remove(iterator& iter);
+        autoPtr<T> remove(iterator& iter);
 
-        //- Erase an entry specified by given iterator.
+        //- Erase an entry specified by given iterator and delete the
+        //- allocated pointer.
         //  Includes a safeguard against the end-iterator.
         bool erase(iterator& iter);
 
-        //- Erase an entry specified by the given key
+        //- Erase an entry specified by the given key and delete the
+        //- allocated pointer.
         bool erase(const Key& key);
 
-        //- Clear all entries from table and deleting any allocated pointers
+        //- Clear all entries from table and delete any allocated pointers
         void clear();
 
         //- Write
@@ -142,6 +148,9 @@ public:
         //- Copy assignment
         void operator=(const this_type& rhs);
 
+        //- Move assignment
+        void operator=(this_type&& rhs);
+
 
     // IOstream Operators
 
@@ -156,6 +165,31 @@ public:
             Ostream& os,
             const HashPtrTable<T, Key, Hash>& tbl
         );
+
+
+    // Housekeeping
+
+        //- No insert() with raw pointers (potential memory leaks)
+        //  Use insert() with autoPtr or set()
+        inline bool insert(const Key&, T*) = delete;
+
+        //- Insert a new entry, not overwriting existing entries.
+        //  \return True if the entry inserted, which means that it did
+        //  not previously exist in the table.
+        inline bool insert(const Key& key, autoPtr<T>& aptr);
+
+        //- Insert a new entry, not overwriting existing entries.
+        inline bool insert(const Key& key, autoPtr<T>&& aptr);
+
+        //- Use set(), not insert() to avoid potential memory leaks
+        //- Assign a new entry, overwriting existing entries.
+        inline bool set(const Key& key, T* ptr);
+
+        //- Assign a new entry, overwriting existing entries.
+        inline bool set(const Key& key, autoPtr<T>& aptr);
+
+        //- Assign a new entry, overwriting existing entries.
+        inline bool set(const Key& key, autoPtr<T>&& aptr);
 };
 
 
@@ -165,6 +199,10 @@ public:
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+#include "HashPtrTableI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
 #ifdef NoRepository
     #include "HashPtrTable.C"
 #endif
diff --git a/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTableIO.C b/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTableIO.C
index fdf3e554a0..fb123d50da 100644
--- a/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTableIO.C
+++ b/src/OpenFOAM/containers/HashTables/HashPtrTable/HashPtrTableIO.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2015 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -65,7 +65,7 @@ void Foam::HashPtrTable<T, Key, Hash>::read(Istream& is, const INew& inew)
                 {
                     Key key;
                     is >> key;
-                    this->insert(key, inew(key, is).ptr());
+                    this->set(key, inew(key, is).ptr());
 
                     is.fatalCheck
                     (
@@ -110,7 +110,7 @@ void Foam::HashPtrTable<T, Key, Hash>::read(Istream& is, const INew& inew)
             is.putBack(lastToken);
             Key key;
             is >> key;
-            this->insert(key, inew(key, is).ptr());
+            this->set(key, inew(key, is).ptr());
 
             is.fatalCheck
             (
@@ -147,7 +147,7 @@ void Foam::HashPtrTable<T, Key, Hash>::read
     {
         const word& k = iter().keyword();
 
-        this->insert(k, inew(dict.subDict(k)).ptr());
+        this->set(k, inew(dict.subDict(k)).ptr());
     }
 }
 
diff --git a/src/OpenFOAM/containers/HashTables/HashTable/HashTable.C b/src/OpenFOAM/containers/HashTables/HashTable/HashTable.C
index b57be56135..7405064ffd 100644
--- a/src/OpenFOAM/containers/HashTables/HashTable/HashTable.C
+++ b/src/OpenFOAM/containers/HashTables/HashTable/HashTable.C
@@ -298,7 +298,7 @@ Foam::label Foam::HashTable<T, Key, Hash>::countEntries
 
 
 template<class T, class Key, class Hash>
-bool Foam::HashTable<T, Key, Hash>::set
+bool Foam::HashTable<T, Key, Hash>::setEntry
 (
     const Key& key,
     const T& obj,
diff --git a/src/OpenFOAM/containers/HashTables/HashTable/HashTable.H b/src/OpenFOAM/containers/HashTables/HashTable/HashTable.H
index 6d67ba3264..cecd82083a 100644
--- a/src/OpenFOAM/containers/HashTables/HashTable/HashTable.H
+++ b/src/OpenFOAM/containers/HashTables/HashTable/HashTable.H
@@ -165,8 +165,10 @@ class HashTable
             }
 
         private:
-            //- Disallow default bitwise copy construct / assignment
+            //- No copy construct
             pair_entry(const pair_entry&) = delete;
+
+            //- No copy assignment
             void operator=(const pair_entry&) = delete;
         };
 
@@ -214,8 +216,10 @@ class HashTable
 
 
         private:
-            //- Disallow default bitwise copy construct / assignment
+            //- No copy construct
             unary_entry(const unary_entry&) = delete;
+
+            //- No copy assignment
             void operator=(const unary_entry&) = delete;
         };
 
@@ -249,7 +253,7 @@ class HashTable
 
         //- Assign a new hash-entry to a possibly already existing key.
         //  \return True if the new entry was set.
-        bool set(const Key& key, const T& obj, const bool overwrite);
+        bool setEntry(const Key& key, const T& obj, const bool overwrite);
 
 
 public:
@@ -606,10 +610,10 @@ public:
         //- Return existing entry or insert a new entry.
         inline T& operator()(const Key& key, const T& deflt);
 
-        //- Copy assignment
+        //- Copy assign
         void operator=(const HashTable<T, Key, Hash>& rhs);
 
-        //- Copy assignment from an initializer list
+        //- Copy assign from an initializer list
         void operator=(std::initializer_list<std::pair<Key, T>> rhs);
 
         //- Move assign
@@ -640,7 +644,7 @@ protected:
         class Iterator
         {
         public:
-          // Typedefs
+        // Typedefs
             using iterator_category = std::forward_iterator_tag;
             using difference_type = this_type::difference_type;
 
@@ -672,7 +676,7 @@ protected:
             >::type;
 
 
-          // Member Functions
+        // Member Functions
 
             //- True if iterator points to an entry
             //  This can be used directly instead of comparing to end()
@@ -681,7 +685,7 @@ protected:
             //- The key associated with the iterator
             inline const Key& key() const;
 
-          // Member Operators
+        // Member Operators
 
             //- Compare hash-entry element pointers.
             //  Independent of const/non-const access
@@ -695,7 +699,7 @@ protected:
         protected:
             friend class HashTable;  // For begin/find constructors
 
-          // Protected Data
+        // Protected Data
 
             //- The selected entry.
             //  MUST be the first member for easy comparison between iterators
@@ -712,7 +716,7 @@ protected:
             label index_;
 
 
-          // Protected Constructors
+        // Protected Constructors
 
             //- Construct null (end iterator)
             inline Iterator();
@@ -724,7 +728,7 @@ protected:
             Iterator(table_type* tbl, const Key& key);
 
 
-          // Protected Member Functions
+        // Protected Member Functions
 
             //- Increment to the next position
             inline void increment();
@@ -756,7 +760,7 @@ public:
             public Iterator<false>
         {
         public:
-          // Typedefs
+        // Typedefs
             using iterator_category = std::forward_iterator_tag;
             using difference_type   = this_type::difference_type;
 
@@ -766,7 +770,7 @@ public:
             using pointer     = this_type::pointer;
             using reference   = this_type::reference;
 
-          // Constructors
+        // Constructors
 
             //- Construct null (end iterator)
             inline iterator() = default;
@@ -778,7 +782,7 @@ public:
             {}
 
 
-          // Member functions/operators
+        // Member functions/operators
 
             //- Non-const access to referenced object
             using Iterator<false>::object;
@@ -800,7 +804,7 @@ public:
             public Iterator<true>
         {
         public:
-          // Typedefs
+        // Typedefs
             using iterator_category = std::forward_iterator_tag;
             using difference_type   = this_type::difference_type;
 
@@ -810,7 +814,7 @@ public:
             using pointer     = this_type::const_pointer;
             using reference   = this_type::const_reference;
 
-          // Constructors
+        // Constructors
 
             //- Construct null (end iterator)
             inline const_iterator() = default;
@@ -837,7 +841,7 @@ public:
             {}
 
 
-          // Member functions/operators
+        // Member functions/operators
 
             //- Const access to referenced object
             using Iterator<true>::object;
@@ -849,7 +853,7 @@ public:
             inline const_iterator& operator++();
             inline const_iterator operator++(int);
 
-          // Assignment
+        // Assignment
 
             const_iterator& operator=(const const_iterator&) = default;
 
diff --git a/src/OpenFOAM/containers/HashTables/HashTable/HashTableI.H b/src/OpenFOAM/containers/HashTables/HashTable/HashTableI.H
index 424d3db57d..47b1b0a91d 100644
--- a/src/OpenFOAM/containers/HashTables/HashTable/HashTableI.H
+++ b/src/OpenFOAM/containers/HashTables/HashTable/HashTableI.H
@@ -121,7 +121,7 @@ inline bool Foam::HashTable<T, Key, Hash>::insert
     const T& obj
 )
 {
-    return this->set(key, obj, false);  // No overwrite
+    return this->setEntry(key, obj, false);  // No overwrite
 }
 
 
@@ -132,7 +132,7 @@ inline bool Foam::HashTable<T, Key, Hash>::set
     const T& obj
 )
 {
-    return this->set(key, obj, true);   // Overwrite
+    return this->setEntry(key, obj, true);   // Overwrite
 }
 
 
diff --git a/src/OpenFOAM/containers/HashTables/PtrMap/PtrMap.H b/src/OpenFOAM/containers/HashTables/PtrMap/PtrMap.H
index 915edf49ad..ed639238d9 100644
--- a/src/OpenFOAM/containers/HashTables/PtrMap/PtrMap.H
+++ b/src/OpenFOAM/containers/HashTables/PtrMap/PtrMap.H
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2017-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -80,11 +80,32 @@ public:
             parent_type(is)
         {}
 
-        //- Construct as copy
+        //- Copy construct
         PtrMap(const this_type& map)
         :
             parent_type(map)
         {}
+
+        //- Move construct
+        PtrMap(this_type&& map)
+        :
+            parent_type(std::move(map))
+        {}
+
+
+    // Member Operators
+
+        //- Copy assignment
+        void operator=(const this_type& rhs)
+        {
+            parent_type::operator=(rhs);
+        }
+
+        //- Move assignment
+        void operator=(this_type&& rhs)
+        {
+            parent_type::operator=(std::move(rhs));
+        }
 };
 
 
diff --git a/src/OpenFOAM/db/IOobjectList/IOobjectList.C b/src/OpenFOAM/db/IOobjectList/IOobjectList.C
index a91134d91c..e8f77d9d0e 100644
--- a/src/OpenFOAM/db/IOobjectList/IOobjectList.C
+++ b/src/OpenFOAM/db/IOobjectList/IOobjectList.C
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2017 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2016-2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2016-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -52,7 +52,7 @@ namespace Foam
                     InfoInFunction << "Found " << iter.key() << endl;
                 }
 
-                results.insert
+                results.set
                 (
                     iter.key(),
                     new IOobject(*(iter.object()))
@@ -123,9 +123,27 @@ namespace Foam
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-Foam::IOobjectList::IOobjectList(const label nIoObjects)
+Foam::IOobjectList::IOobjectList()
 :
-    HashPtrTable<IOobject>(nIoObjects)
+    HashPtrTable<IOobject>()
+{}
+
+
+Foam::IOobjectList::IOobjectList(const label nObjects)
+:
+    HashPtrTable<IOobject>(nObjects)  // Could also use 2*nObjects instead
+{}
+
+
+Foam::IOobjectList::IOobjectList(const IOobjectList& list)
+:
+    HashPtrTable<IOobject>(list)
+{}
+
+
+Foam::IOobjectList::IOobjectList(IOobjectList&& list)
+:
+    HashPtrTable<IOobject>(std::move(list))
 {}
 
 
@@ -142,7 +160,7 @@ Foam::IOobjectList::IOobjectList
     HashPtrTable<IOobject>()
 {
     word newInstance;
-    fileNameList ObjectNames = fileHandler().readObjects
+    fileNameList objNames = fileHandler().readObjects
     (
         db,
         instance,
@@ -150,9 +168,9 @@ Foam::IOobjectList::IOobjectList
         newInstance
     );
 
-    for (const auto& objName : ObjectNames)
+    for (const auto& objName : objNames)
     {
-        IOobject* objectPtr = new IOobject
+        auto objectPtr = autoPtr<IOobject>::New
         (
             objName,
             newInstance,
@@ -183,35 +201,35 @@ Foam::IOobjectList::IOobjectList
         {
             insert(objectPtr->name(), objectPtr);
         }
-        else
-        {
-            delete objectPtr;
-        }
     }
 }
 
 
-Foam::IOobjectList::IOobjectList(const IOobjectList& iolist)
-:
-    HashPtrTable<IOobject>(iolist)
-{}
-
-
-// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
-
-Foam::IOobjectList::~IOobjectList()
-{}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-bool Foam::IOobjectList::add(IOobject& io)
+bool Foam::IOobjectList::add(autoPtr<IOobject>& objectPtr)
 {
-    return insert(io.name(), &io);
+    if (objectPtr.valid())
+    {
+        return insert(objectPtr->name(), objectPtr);
+    }
+
+    return false;
 }
 
 
-bool Foam::IOobjectList::remove(IOobject& io)
+bool Foam::IOobjectList::add(autoPtr<IOobject>&& objectPtr)
+{
+    if (objectPtr.valid())
+    {
+        return insert(objectPtr->name(), objectPtr);
+    }
+
+    return false;
+}
+
+
+bool Foam::IOobjectList::remove(const IOobject& io)
 {
     return erase(io.name());
 }
@@ -219,7 +237,7 @@ bool Foam::IOobjectList::remove(IOobject& io)
 
 Foam::IOobject* Foam::IOobjectList::lookup(const word& name) const
 {
-    const_iterator iter = find(name);
+    const_iterator iter = cfind(name);
 
     if (iter.found())
     {
@@ -267,7 +285,7 @@ Foam::IOobjectList Foam::IOobjectList::lookupClass(const word& clsName) const
                 InfoInFunction << "Found " << iter.key() << endl;
             }
 
-            results.insert
+            results.set
             (
                 iter.key(),
                 new IOobject(*(iter.object()))
@@ -369,6 +387,14 @@ Foam::wordList Foam::IOobjectList::sortedNames
 }
 
 
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+void Foam::IOobjectList::operator=(IOobjectList&& list)
+{
+    transfer(list);
+}
+
+
 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
 
 Foam::Ostream& Foam::operator<<(Ostream& os, const IOobjectList& list)
diff --git a/src/OpenFOAM/db/IOobjectList/IOobjectList.H b/src/OpenFOAM/db/IOobjectList/IOobjectList.H
index 3478473db1..247c6a0633 100644
--- a/src/OpenFOAM/db/IOobjectList/IOobjectList.H
+++ b/src/OpenFOAM/db/IOobjectList/IOobjectList.H
@@ -3,7 +3,7 @@
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
     \\  /    A nd           | Copyright (C) 2011-2016 OpenFOAM Foundation
-     \\/     M anipulation  | Copyright (C) 2016-2017 OpenCFD Ltd.
+     \\/     M anipulation  | Copyright (C) 2016-2018 OpenCFD Ltd.
 -------------------------------------------------------------------------------
 License
     This file is part of OpenFOAM.
@@ -45,7 +45,7 @@ SourceFiles
 namespace Foam
 {
 
-// Forward declaration of friend functions and operators
+// Forward declarations
 class IOobjectList;
 Ostream& operator<<(Ostream& os, const IOobjectList& list);
 
@@ -58,18 +58,21 @@ class IOobjectList
 :
     public HashPtrTable<IOobject>
 {
-    // Private Member Functions
-
-        //- No copy assignment
-        void operator=(const IOobjectList&) = delete;
-
-
 public:
 
     // Constructors
 
-        //- Construct given an initial estimate for the number of entries
-        explicit IOobjectList(const label nIoObjects = 128);
+        //- Construct null with default (128) table capacity
+        IOobjectList();
+
+        //- Construct given initial table capacity
+        explicit IOobjectList(const label nObjects);
+
+        //- Copy construct
+        IOobjectList(const IOobjectList& list);
+
+        //- Move construct
+        IOobjectList(IOobjectList&& list);
 
         //- Construct from objectRegistry and instance path
         IOobjectList
@@ -82,12 +85,9 @@ public:
             bool registerObject = true
         );
 
-        //- Construct as copy
-        IOobjectList(const IOobjectList& iolist);
-
 
     //- Destructor
-    ~IOobjectList();
+    ~IOobjectList() = default;
 
 
     // Member functions
@@ -95,10 +95,16 @@ public:
     // Basic methods
 
         //- Add an IOobject to the list
-        bool add(IOobject& io);
+        bool add(autoPtr<IOobject>& objectPtr);
+
+        //- Add an IOobject to the list
+        bool add(autoPtr<IOobject>&& objectPtr);
+
+        //- Remove an IOobject from the list, by iterator
+        using HashPtrTable<IOobject>::remove;
 
         //- Remove an IOobject from the list
-        bool remove(IOobject& io);
+        bool remove(const IOobject& io);
 
 
     // Lookup
@@ -230,6 +236,15 @@ public:
         wordList sortedNames(const word& clsName, const wordRes& matcher) const;
 
 
+    // Member Operators
+
+        //- No copy assignment
+        void operator=(const IOobjectList&) = delete;
+
+        //- Move assignment
+        void operator=(IOobjectList&& list);
+
+
     // Ostream Operator
 
         friend Ostream& operator<<(Ostream& os, const IOobjectList& list);
diff --git a/src/OpenFOAM/global/fileOperations/masterUncollatedFileOperation/masterUncollatedFileOperation.C b/src/OpenFOAM/global/fileOperations/masterUncollatedFileOperation/masterUncollatedFileOperation.C
index 4d97c0bdad..1a81109921 100644
--- a/src/OpenFOAM/global/fileOperations/masterUncollatedFileOperation/masterUncollatedFileOperation.C
+++ b/src/OpenFOAM/global/fileOperations/masterUncollatedFileOperation/masterUncollatedFileOperation.C
@@ -1460,7 +1460,7 @@ Foam::instantList Foam::fileOperations::masterUncollatedFileOperation::findTimes
 
         instantList* tPtr = new instantList(std::move(times));
 
-        times_.insert(directory, tPtr);
+        times_.set(directory, tPtr);
 
         if (debug)
         {
diff --git a/src/OpenFOAM/graph/graph.C b/src/OpenFOAM/graph/graph.C
index 26e0325129..5e0d11579f 100644
--- a/src/OpenFOAM/graph/graph.C
+++ b/src/OpenFOAM/graph/graph.C
@@ -66,7 +66,7 @@ void Foam::graph::readCurves(Istream& is)
         y[i] = xyData[i].y_;
     }
 
-    insert
+    set
     (
         wordify(yName_),
         new curve(wordify(yName_), curve::curveStyle::CONTINUOUS, y)
@@ -105,7 +105,11 @@ Foam::graph::graph
     yName_(yName),
     x_(x)
 {
-    insert(wordify(yName), new curve(yName, curve::curveStyle::CONTINUOUS, y));
+    set
+    (
+        wordify(yName),
+        new curve(yName, curve::curveStyle::CONTINUOUS, y)
+    );
 }
 
 
diff --git a/src/functionObjects/graphics/runTimePostProcessing/scene.C b/src/functionObjects/graphics/runTimePostProcessing/scene.C
index d984e9f03b..e39f3010e5 100644
--- a/src/functionObjects/graphics/runTimePostProcessing/scene.C
+++ b/src/functionObjects/graphics/runTimePostProcessing/scene.C
@@ -122,11 +122,10 @@ void Foam::functionObjects::runTimePostPro::scene::readColours
     const dictionary& dict
 )
 {
-    const wordList colours = dict.toc();
-    forAll(colours, i)
+    const wordList colours(dict.toc());
+    for (const word& c : colours)
     {
-        const word& c = colours[i];
-        colours_.insert(c, Function1<vector>::New(c, dict).ptr());
+        colours_.insert(c, Function1<vector>::New(c, dict));
     }
 }
 
diff --git a/src/genericPatchFields/genericFaPatchField/genericFaPatchField.C b/src/genericPatchFields/genericFaPatchField/genericFaPatchField.C
index 5d5e7ccc04..805e5c0e52 100644
--- a/src/genericPatchFields/genericFaPatchField/genericFaPatchField.C
+++ b/src/genericPatchFields/genericFaPatchField/genericFaPatchField.C
@@ -109,7 +109,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                             scalarFields_.insert
                             (
                                 iter().keyword(),
-                                new scalarField(0)
+                                autoPtr<scalarField>::New()
                             );
                         }
                         else
@@ -133,7 +133,8 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                      == token::Compound<List<scalar>>::typeName
                     )
                     {
-                        scalarField* fPtr = new scalarField;
+                        auto fPtr = autoPtr<scalarField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<scalar>>>
@@ -167,7 +168,8 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                      == token::Compound<List<vector>>::typeName
                     )
                     {
-                        vectorField* fPtr = new vectorField;
+                        auto fPtr = autoPtr<vectorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<vector>>>
@@ -201,7 +203,8 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                      == token::Compound<List<sphericalTensor>>::typeName
                     )
                     {
-                        sphericalTensorField* fPtr = new sphericalTensorField;
+                        auto fPtr = autoPtr<sphericalTensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast
@@ -238,7 +241,8 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                      == token::Compound<List<symmTensor>>::typeName
                     )
                     {
-                        symmTensorField* fPtr = new symmTensorField;
+                        auto fPtr = autoPtr<symmTensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast
@@ -275,7 +279,8 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                      == token::Compound<List<tensor>>::typeName
                     )
                     {
-                        tensorField* fPtr = new tensorField;
+                        auto fPtr = autoPtr<tensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<tensor>>>
@@ -331,7 +336,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                         scalarFields_.insert
                         (
                             iter().keyword(),
-                            new scalarField
+                            autoPtr<scalarField>::New
                             (
                                 this->size(),
                                 fieldToken.number()
@@ -352,7 +357,11 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                             vectorFields_.insert
                             (
                                 iter().keyword(),
-                                new vectorField(this->size(), vs)
+                                autoPtr<vectorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else if (l.size() == sphericalTensor::nComponents)
@@ -362,7 +371,11 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                             sphericalTensorFields_.insert
                             (
                                 iter().keyword(),
-                                new sphericalTensorField(this->size(), vs)
+                                autoPtr<sphericalTensorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else if (l.size() == symmTensor::nComponents)
@@ -372,7 +385,11 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                             symmTensorFields_.insert
                             (
                                 iter().keyword(),
-                                new symmTensorField(this->size(), vs)
+                                autoPtr<symmTensorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else if (l.size() == tensor::nComponents)
@@ -387,7 +404,11 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
                             tensorFields_.insert
                             (
                                 iter().keyword(),
-                                new tensorField(this->size(), vs)
+                                autoPtr<tensorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else
@@ -434,7 +455,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
         scalarFields_.insert
         (
             iter.key(),
-            new scalarField(*iter(), mapper)
+            autoPtr<scalarField>::New(*iter(), mapper)
         );
     }
 
@@ -448,7 +469,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
         vectorFields_.insert
         (
             iter.key(),
-            new vectorField(*iter(), mapper)
+            autoPtr<vectorField>::New(*iter(), mapper)
         );
     }
 
@@ -462,7 +483,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
         sphericalTensorFields_.insert
         (
             iter.key(),
-            new sphericalTensorField(*iter(), mapper)
+            autoPtr<sphericalTensorField>::New(*iter(), mapper)
         );
     }
 
@@ -476,7 +497,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
         symmTensorFields_.insert
         (
             iter.key(),
-            new symmTensorField(*iter(), mapper)
+            autoPtr<symmTensorField>::New(*iter(), mapper)
         );
     }
 
@@ -490,7 +511,7 @@ Foam::genericFaPatchField<Type>::genericFaPatchField
         tensorFields_.insert
         (
             iter.key(),
-            new tensorField(*iter(), mapper)
+            autoPtr<tensorField>::New(*iter(), mapper)
         );
     }
 }
diff --git a/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C b/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C
index d0a7f439d2..b47b219245 100644
--- a/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C
+++ b/src/genericPatchFields/genericFvPatchField/genericFvPatchField.C
@@ -109,7 +109,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                             scalarFields_.insert
                             (
                                 iter().keyword(),
-                                new scalarField(0)
+                                autoPtr<scalarField>::New()
                             );
                         }
                         else
@@ -133,7 +133,8 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                      == token::Compound<List<scalar>>::typeName
                     )
                     {
-                        scalarField* fPtr = new scalarField;
+                        auto fPtr = autoPtr<scalarField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<scalar>>>
@@ -167,7 +168,8 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                      == token::Compound<List<vector>>::typeName
                     )
                     {
-                        vectorField* fPtr = new vectorField;
+                        auto fPtr = autoPtr<vectorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<vector>>>
@@ -201,7 +203,8 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                      == token::Compound<List<sphericalTensor>>::typeName
                     )
                     {
-                        sphericalTensorField* fPtr = new sphericalTensorField;
+                        auto fPtr = autoPtr<sphericalTensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast
@@ -238,7 +241,8 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                      == token::Compound<List<symmTensor>>::typeName
                     )
                     {
-                        symmTensorField* fPtr = new symmTensorField;
+                        auto fPtr = autoPtr<symmTensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast
@@ -275,7 +279,8 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                      == token::Compound<List<tensor>>::typeName
                     )
                     {
-                        tensorField* fPtr = new tensorField;
+                        auto fPtr = autoPtr<tensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<tensor>>>
@@ -331,7 +336,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                         scalarFields_.insert
                         (
                             iter().keyword(),
-                            new scalarField
+                            autoPtr<scalarField>::New
                             (
                                 this->size(),
                                 fieldToken.number()
@@ -352,7 +357,11 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                             vectorFields_.insert
                             (
                                 iter().keyword(),
-                                new vectorField(this->size(), vs)
+                                autoPtr<vectorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else if (l.size() == sphericalTensor::nComponents)
@@ -362,7 +371,11 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                             sphericalTensorFields_.insert
                             (
                                 iter().keyword(),
-                                new sphericalTensorField(this->size(), vs)
+                                autoPtr<sphericalTensorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else if (l.size() == symmTensor::nComponents)
@@ -372,7 +385,11 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                             symmTensorFields_.insert
                             (
                                 iter().keyword(),
-                                new symmTensorField(this->size(), vs)
+                                autoPtr<symmTensorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else if (l.size() == tensor::nComponents)
@@ -387,7 +404,11 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
                             tensorFields_.insert
                             (
                                 iter().keyword(),
-                                new tensorField(this->size(), vs)
+                                autoPtr<tensorField>::New
+                                (
+                                    this->size(),
+                                    vs
+                                )
                             );
                         }
                         else
@@ -434,7 +455,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
         scalarFields_.insert
         (
             iter.key(),
-            new scalarField(*iter(), mapper)
+            autoPtr<scalarField>::New(*iter(), mapper)
         );
     }
 
@@ -448,7 +469,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
         vectorFields_.insert
         (
             iter.key(),
-            new vectorField(*iter(), mapper)
+            autoPtr<vectorField>::New(*iter(), mapper)
         );
     }
 
@@ -462,7 +483,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
         sphericalTensorFields_.insert
         (
             iter.key(),
-            new sphericalTensorField(*iter(), mapper)
+            autoPtr<sphericalTensorField>::New(*iter(), mapper)
         );
     }
 
@@ -476,7 +497,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
         symmTensorFields_.insert
         (
             iter.key(),
-            new symmTensorField(*iter(), mapper)
+            autoPtr<symmTensorField>::New(*iter(), mapper)
         );
     }
 
@@ -490,7 +511,7 @@ Foam::genericFvPatchField<Type>::genericFvPatchField
         tensorFields_.insert
         (
             iter.key(),
-            new tensorField(*iter(), mapper)
+            autoPtr<tensorField>::New(*iter(), mapper)
         );
     }
 }
diff --git a/src/genericPatchFields/genericPointPatchField/genericPointPatchField.C b/src/genericPatchFields/genericPointPatchField/genericPointPatchField.C
index 52744cdf3c..3ef8cf3b51 100644
--- a/src/genericPatchFields/genericPointPatchField/genericPointPatchField.C
+++ b/src/genericPatchFields/genericPointPatchField/genericPointPatchField.C
@@ -87,7 +87,7 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
                             scalarFields_.insert
                             (
                                 iter().keyword(),
-                                new scalarField(0)
+                                autoPtr<scalarField>::New()
                             );
                         }
                         else
@@ -111,7 +111,8 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
                      == token::Compound<List<scalar>>::typeName
                     )
                     {
-                        scalarField* fPtr = new scalarField;
+                        auto fPtr = autoPtr<scalarField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<scalar>>>
@@ -145,7 +146,8 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
                      == token::Compound<List<vector>>::typeName
                     )
                     {
-                        vectorField* fPtr = new vectorField;
+                        auto fPtr = autoPtr<vectorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<vector>>>
@@ -179,7 +181,8 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
                      == token::Compound<List<sphericalTensor>>::typeName
                     )
                     {
-                        sphericalTensorField* fPtr = new sphericalTensorField;
+                        auto fPtr = autoPtr<sphericalTensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast
@@ -216,7 +219,8 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
                      == token::Compound<List<symmTensor>>::typeName
                     )
                     {
-                        symmTensorField* fPtr = new symmTensorField;
+                        auto fPtr = autoPtr<symmTensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast
@@ -253,7 +257,8 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
                      == token::Compound<List<tensor>>::typeName
                     )
                     {
-                        tensorField* fPtr = new tensorField;
+                        auto fPtr = autoPtr<tensorField>::New();
+
                         fPtr->transfer
                         (
                             dynamicCast<token::Compound<List<tensor>>>
@@ -325,7 +330,7 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
         scalarFields_.insert
         (
             iter.key(),
-            new scalarField(*iter(), mapper)
+            autoPtr<scalarField>::New(*iter(), mapper)
         );
     }
 
@@ -339,7 +344,7 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
         vectorFields_.insert
         (
             iter.key(),
-            new vectorField(*iter(), mapper)
+            autoPtr<vectorField>::New(*iter(), mapper)
         );
     }
 
@@ -353,7 +358,7 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
         sphericalTensorFields_.insert
         (
             iter.key(),
-            new sphericalTensorField(*iter(), mapper)
+            autoPtr<sphericalTensorField>::New(*iter(), mapper)
         );
     }
 
@@ -367,7 +372,7 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
         symmTensorFields_.insert
         (
             iter.key(),
-            new symmTensorField(*iter(), mapper)
+            autoPtr<symmTensorField>::New(*iter(), mapper)
         );
     }
 
@@ -381,7 +386,7 @@ Foam::genericPointPatchField<Type>::genericPointPatchField
         tensorFields_.insert
         (
             iter.key(),
-            new tensorField(*iter(), mapper)
+            autoPtr<tensorField>::New(*iter(), mapper)
         );
     }
 }
diff --git a/src/sampling/probes/probes.C b/src/sampling/probes/probes.C
index 562702f8c9..427d038b9b 100644
--- a/src/sampling/probes/probes.C
+++ b/src/sampling/probes/probes.C
@@ -217,13 +217,13 @@ Foam::label Foam::probes::prepare()
         probeDir.clean();  // Remove unneeded ".."
 
         // ignore known fields, close streams for fields that no longer exist
-        forAllIter(HashPtrTable<OFstream>, probeFilePtrs_, iter)
+        forAllIters(probeFilePtrs_, iter)
         {
             if (!currentFields.erase(iter.key()))
             {
                 DebugInfo<< "close probe stream: " << iter()->name() << endl;
 
-                delete probeFilePtrs_.remove(iter);
+                probeFilePtrs_.remove(iter);
             }
         }
 
@@ -233,9 +233,8 @@ Foam::label Foam::probes::prepare()
             // Create directory if does not exist.
             mkDir(probeDir);
 
-            OFstream* fPtr = new OFstream(probeDir/fieldName);
-
-            OFstream& fout = *fPtr;
+            auto fPtr = autoPtr<OFstream>::New(probeDir/fieldName);
+            auto& fout = *fPtr;
 
             DebugInfo<< "open probe stream: " << fout.name() << endl;
 
@@ -426,9 +425,8 @@ void Foam::probes::updateMesh(const mapPolyMesh& mpm)
             DynamicList<label> elems(faceList_.size());
 
             const labelList& reverseMap = mpm.reverseFaceMap();
-            forAll(faceList_, i)
+            for (const label facei : faceList_)
             {
-                label facei = faceList_[i];
                 if (facei != -1)
                 {
                     label newFacei = reverseMap[facei];
diff --git a/src/thermophysicalModels/reactionThermo/chemistryReaders/chemkinReader/chemkinLexer.L b/src/thermophysicalModels/reactionThermo/chemistryReaders/chemkinReader/chemkinLexer.L
index 9db567ba47..784d38f20e 100644
--- a/src/thermophysicalModels/reactionThermo/chemistryReaders/chemkinReader/chemkinLexer.L
+++ b/src/thermophysicalModels/reactionThermo/chemistryReaders/chemkinReader/chemkinLexer.L
@@ -600,7 +600,7 @@ bool finishReaction = false;
         speciesThermo_.insert
         (
             currentSpecieName,
-            new gasHThermoPhysics
+            autoPtr<gasHThermoPhysics>::New
             (
                 janafThermo<perfectGas<specie>>
                 (

From 6180e0db8b7161f0cbb1b366c3cf618d29aaf751 Mon Sep 17 00:00:00 2001
From: Andrew Heather <a.heather@opencfd.co.uk>
Date: Fri, 18 May 2018 10:49:10 +0100
Subject: [PATCH 8/8] INT: Additional integration updates/clean-up

---
 applications/solvers/compressible/rhoPimpleFoam/Make/options | 3 ++-
 applications/solvers/compressible/rhoSimpleFoam/Make/options | 3 ++-
 .../heatTransfer/buoyantBoussinesqPimpleFoam/Make/options    | 3 ++-
 .../heatTransfer/buoyantBoussinesqSimpleFoam/Make/options    | 3 ++-
 .../solvers/heatTransfer/buoyantPimpleFoam/Make/options      | 3 ++-
 .../solvers/heatTransfer/buoyantSimpleFoam/Make/options      | 3 ++-
 applications/solvers/incompressible/pimpleFoam/Make/options  | 3 ++-
 applications/solvers/incompressible/simpleFoam/Make/options  | 3 ++-
 src/atmosphericModels/Make/options                           | 4 ++--
 .../derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H | 5 ++---
 .../atmBoundaryLayerInletEpsilonFvPatchScalarField.C         | 2 +-
 .../atmBoundaryLayerInletKFvPatchScalarField.C               | 2 +-
 .../atmBoundaryLayerInletVelocityFvPatchVectorField.C        | 2 +-
 .../kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H              | 2 +-
 .../simpleFoam/turbineSiting/system/controlDict              | 2 --
 15 files changed, 24 insertions(+), 19 deletions(-)

diff --git a/applications/solvers/compressible/rhoPimpleFoam/Make/options b/applications/solvers/compressible/rhoPimpleFoam/Make/options
index c3f820e3be..958a2d2fc8 100644
--- a/applications/solvers/compressible/rhoPimpleFoam/Make/options
+++ b/applications/solvers/compressible/rhoPimpleFoam/Make/options
@@ -22,4 +22,5 @@ EXE_LIBS = \
     -lfvOptions \
     -ldynamicFvMesh \
     -ltopoChangerFvMesh \
-    -ldynamicMesh
+    -ldynamicMesh \
+    -latmosphericModels
diff --git a/applications/solvers/compressible/rhoSimpleFoam/Make/options b/applications/solvers/compressible/rhoSimpleFoam/Make/options
index aecfa919e7..d3514cf84d 100644
--- a/applications/solvers/compressible/rhoSimpleFoam/Make/options
+++ b/applications/solvers/compressible/rhoSimpleFoam/Make/options
@@ -17,4 +17,5 @@ EXE_LIBS = \
     -lfiniteVolume \
     -lsampling \
     -lmeshTools \
-    -lfvOptions
+    -lfvOptions \
+    -latmosphericModels
diff --git a/applications/solvers/heatTransfer/buoyantBoussinesqPimpleFoam/Make/options b/applications/solvers/heatTransfer/buoyantBoussinesqPimpleFoam/Make/options
index 043be5d13c..767fc60e5b 100644
--- a/applications/solvers/heatTransfer/buoyantBoussinesqPimpleFoam/Make/options
+++ b/applications/solvers/heatTransfer/buoyantBoussinesqPimpleFoam/Make/options
@@ -18,4 +18,5 @@ EXE_LIBS = \
     -lfiniteVolume \
     -lfvOptions \
     -lsampling \
-    -lmeshTools
+    -lmeshTools \
+    -latmosphericModels
diff --git a/applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/Make/options b/applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/Make/options
index 9046994122..e805213895 100644
--- a/applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/Make/options
+++ b/applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/Make/options
@@ -16,4 +16,5 @@ EXE_LIBS = \
     -lfiniteVolume \
     -lsampling \
     -lmeshTools \
-    -lfvOptions
+    -lfvOptions \
+    -latmosphericModels
diff --git a/applications/solvers/heatTransfer/buoyantPimpleFoam/Make/options b/applications/solvers/heatTransfer/buoyantPimpleFoam/Make/options
index 934083f2bb..cc85f5e8cc 100644
--- a/applications/solvers/heatTransfer/buoyantPimpleFoam/Make/options
+++ b/applications/solvers/heatTransfer/buoyantPimpleFoam/Make/options
@@ -18,4 +18,5 @@ EXE_LIBS = \
     -lradiationModels \
     -lspecie \
     -lturbulenceModels \
-    -lcompressibleTurbulenceModels
+    -lcompressibleTurbulenceModels \
+    -latmosphericModels
diff --git a/applications/solvers/heatTransfer/buoyantSimpleFoam/Make/options b/applications/solvers/heatTransfer/buoyantSimpleFoam/Make/options
index 0d9433bac7..c5694b8e51 100644
--- a/applications/solvers/heatTransfer/buoyantSimpleFoam/Make/options
+++ b/applications/solvers/heatTransfer/buoyantSimpleFoam/Make/options
@@ -19,4 +19,5 @@ EXE_LIBS = \
     -lradiationModels \
     -lturbulenceModels \
     -lcompressibleTurbulenceModels \
-    -lmeshTools
+    -lmeshTools \
+    -latmosphericModels
diff --git a/applications/solvers/incompressible/pimpleFoam/Make/options b/applications/solvers/incompressible/pimpleFoam/Make/options
index 14a5a2d924..dc33faebf8 100644
--- a/applications/solvers/incompressible/pimpleFoam/Make/options
+++ b/applications/solvers/incompressible/pimpleFoam/Make/options
@@ -19,4 +19,5 @@ EXE_LIBS = \
     -ldynamicFvMesh \
     -ltopoChangerFvMesh \
     -ldynamicMesh \
-    -lmeshTools
+    -lmeshTools \
+    -latmosphericModels
diff --git a/applications/solvers/incompressible/simpleFoam/Make/options b/applications/solvers/incompressible/simpleFoam/Make/options
index 1d9ded80bf..13060603d1 100644
--- a/applications/solvers/incompressible/simpleFoam/Make/options
+++ b/applications/solvers/incompressible/simpleFoam/Make/options
@@ -15,4 +15,5 @@ EXE_LIBS = \
     -lfiniteVolume \
     -lmeshTools \
     -lfvOptions \
-    -lsampling
+    -lsampling \
+    -latmosphericModels
diff --git a/src/atmosphericModels/Make/options b/src/atmosphericModels/Make/options
index 3ffea18bec..93f257f6ee 100644
--- a/src/atmosphericModels/Make/options
+++ b/src/atmosphericModels/Make/options
@@ -5,7 +5,7 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/triSurface/lnInclude \
+    -I$(LIB_SRC)/surfMesh/lnInclude \
     -I$(LIB_SRC)/fvOptions/lnInclude
 
 LIB_LIBS = \
@@ -14,5 +14,5 @@ LIB_LIBS = \
     -lincompressibleTransportModels \
     -lfiniteVolume \
     -lmeshTools \
-    -ltriSurface \
+    -lsurfMesh \
     -lfvOptions
diff --git a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
index a7360f8d38..0f244ed59e 100644
--- a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayer/atmBoundaryLayer.H
@@ -123,7 +123,7 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-       Class atmBoundaryLayer Declaration
+                      Class atmBoundaryLayer Declaration
 \*---------------------------------------------------------------------------*/
 
 class atmBoundaryLayer
@@ -168,8 +168,7 @@ public:
         //- Construct from the coordinates field and dictionary
         atmBoundaryLayer(const vectorField& p, const dictionary&);
 
-        //- Construct by mapping given
-        // atmBoundaryLayer onto a new patch
+        //- Construct by mapping given atmBoundaryLayer onto a new patch
         atmBoundaryLayer
         (
             const atmBoundaryLayer&,
diff --git a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
index 65e748f3dc..8aeeed713d 100644
--- a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletEpsilon/atmBoundaryLayerInletEpsilonFvPatchScalarField.C
@@ -133,7 +133,7 @@ void atmBoundaryLayerInletEpsilonFvPatchScalarField::write(Ostream& os) const
 {
     fvPatchScalarField::write(os);
     atmBoundaryLayer::write(os);
-    writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
+    os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
     writeEntry("value", os);
 }
 
diff --git a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
index 4078d85482..06f366089d 100644
--- a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletK/atmBoundaryLayerInletKFvPatchScalarField.C
@@ -133,7 +133,7 @@ void atmBoundaryLayerInletKFvPatchScalarField::write(Ostream& os) const
 {
     fvPatchScalarField::write(os);
     atmBoundaryLayer::write(os);
-    writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
+    os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
     writeEntry("value", os);
 }
 
diff --git a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
index ec9ee6e79c..b1d423e545 100644
--- a/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
+++ b/src/atmosphericModels/derivedFvPatchFields/atmBoundaryLayerInletVelocity/atmBoundaryLayerInletVelocityFvPatchVectorField.C
@@ -133,7 +133,7 @@ void atmBoundaryLayerInletVelocityFvPatchVectorField::write(Ostream& os) const
 {
     fvPatchVectorField::write(os);
     atmBoundaryLayer::write(os);
-    writeEntryIfDifferent<word>(os, "phi", "phi", phiName_);
+    os.writeEntryIfDifferent<word>("phi", "phi", phiName_);
     writeEntry("value", os);
 }
 
diff --git a/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H b/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H
index b72992d9a4..a9f2affd6b 100644
--- a/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H
+++ b/src/atmosphericModels/kEpsilonLopesdaCosta/kEpsilonLopesdaCosta.H
@@ -74,7 +74,7 @@ namespace RASModels
 {
 
 /*---------------------------------------------------------------------------*\
-                           Class kEpsilonLopesdaCosta Declaration
+                    Class kEpsilonLopesdaCosta Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class BasicTurbulenceModel>
diff --git a/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict b/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict
index bedecb22e9..2654bdd256 100644
--- a/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict
+++ b/tutorials/incompressible/simpleFoam/turbineSiting/system/controlDict
@@ -44,6 +44,4 @@ timePrecision   6;
 
 runTimeModifiable yes;
 
-libs ("libatmosphericModels.so");
-
 // ************************************************************************* //