Merge branch 'master' of /home/hunt2/mattijs/OpenFOAM/OpenFOAM-dev/.

applications/test/BSpline/BSplineTest.C

@@ -0,0 +1,67 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+#include "argList.H"
+
+#include "vector.H"
+#include "IFstream.H"
+#include "BSpline.H"
+
+using namespace Foam;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+// Main program:
+
+int main(int argc, char *argv[])
+{
+    argList::noParallel();
+    argList::validArgs.insert("file .. fileN");
+
+    argList args(argc, argv, false, true);
+
+    forAll(args.additionalArgs(), argI)
+    {
+        const string& srcFile = args.additionalArgs()[argI];
+        Info<< nl << "reading " << srcFile << nl;
+        IFstream ifs(srcFile);
+
+        List<pointField> splinePointFields(ifs);
+
+        forAll(splinePointFields, splineI)
+        {
+            Info<<"convert " << splinePointFields[splineI] << " to bspline" << endl;
+
+            BSpline spl(splinePointFields[splineI], vector::zero, vector::zero);
+
+            Info<< "1/2 = " << spl.position(0.5) << endl;
+        }
+    }
+
+
+    return 0;
+}
+
+
+// ************************************************************************* //

applications/test/BSpline/Make/files

@@ -0,0 +1,3 @@
+BSplineTest.C
+
+EXE = $(FOAM_USER_APPBIN)/BSplineTest

applications/test/BSpline/Make/options

@@ -0,0 +1,5 @@
+EXE_INC = \
+    -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/mesh/blockMesh/lnInclude
+
+EXE_LIBS = -lblockMesh

applications/test/BSpline/test-splines

@@ -0,0 +1,69 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.6                                   |
+|   \\  /    A nd           | Web:      http://www.OpenFOAM.org               |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+(
+    // Upper body longitudinal splines.
+    (
+        (-0.22685 -0.01125166 0)       // 7
+        (-0.21685 -0.01340204 0)
+        (-0.20685 -0.01529684 0)
+        (-0.19685 -0.01694748 0)
+        (-0.18685 -0.01836538 0)
+        (-0.17685 -0.01956197 0)
+        (-0.16685 -0.02054868 0)
+        (-0.15685 -0.02133693 0)
+        (-0.14685 -0.02193816 0)
+        (-0.13685 -0.02236377 0)
+        (-0.12685 -0.02262521 0)
+        (-0.11685 -0.02273389 0)       // 2
+    )
+
+    (
+        (-0.22685 0 0.01125166)        // 8
+        (-0.21685 0 0.01340204)
+        (-0.20685 0 0.01529684)
+        (-0.19685 0 0.01694748)
+        (-0.18685 0 0.01836538)
+        (-0.17685 0 0.01956197)
+        (-0.16685 0 0.02054868)
+        (-0.15685 0 0.02133693)
+        (-0.14685 0 0.02193816)
+        (-0.13685 0 0.02236377)
+        (-0.12685 0 0.02262521)
+        (-0.11685 0 0.02273389)        // 3
+    )
+
+    (
+        (-0.22685 0.01125166 0)        // 9
+        (-0.21685 0.01340204 0)
+        (-0.20685 0.01529684 0)
+        (-0.19685 0.01694748 0)
+        (-0.18685 0.01836538 0)
+        (-0.17685 0.01956197 0)
+        (-0.16685 0.02054868 0)
+        (-0.15685 0.02133693 0)
+        (-0.14685 0.02193816 0)
+        (-0.13685 0.02236377 0)
+        (-0.12685 0.02262521 0)
+        (-0.11685 0.02273389 0)        // 4
+    )
+
+    (
+        (-0.22685 0 -0.01125166)       // 6
+        (-0.21685 0 -0.01340204)
+        (-0.20685 0 -0.01529684)
+        (-0.19685 0 -0.01694748)
+        (-0.18685 0 -0.01836538)
+        (-0.17685 0 -0.01956197)
+        (-0.16685 0 -0.02054868)
+        (-0.15685 0 -0.02133693)
+        (-0.14685 0 -0.02193816)
+        (-0.13685 0 -0.02236377)
+        (-0.12685 0 -0.02262521)
+        (-0.11685 0 -0.02273389)       // 1
+    )
+);

applications/utilities/postProcessing/dataConversion/foamToVTK/writeFuns.C

@@ -220,7 +220,7 @@ void Foam::writeFuns::writePointDataHeader
 }
 
 
-void Foam::writeFuns::insert(const scalar& pt, DynamicList<floatScalar>& dest)
+void Foam::writeFuns::insert(const scalar pt, DynamicList<floatScalar>& dest)
 {
     dest.append(float(pt));
 }

applications/utilities/postProcessing/dataConversion/foamToVTK/writeFuns.H

@@ -86,7 +86,7 @@ public:
 
     // Convert to VTK and store
 
-        static void insert(const scalar&, DynamicList<floatScalar>&);
+        static void insert(const scalar, DynamicList<floatScalar>&);
         static void insert(const point&, DynamicList<floatScalar>&);
         static void insert(const sphericalTensor&, DynamicList<floatScalar>&);
         static void insert(const symmTensor&, DynamicList<floatScalar>&);

applications/utilities/preProcessing/wallFunctionTable/Allwclean

@@ -0,0 +1,5 @@
+#!/bin/sh
+
+wclean libso tabulatedWallFunction
+wclean
+

applications/utilities/preProcessing/wallFunctionTable/Allwmake

@@ -0,0 +1,5 @@
+#!/bin/sh
+
+wmake libso tabulatedWallFunction
+wmake
+

applications/utilities/preProcessing/wallFunctionTable/Make/files

@@ -0,0 +1,3 @@
+wallFunctionTableApp.C
+
+EXE = $(FOAM_APPBIN)/wallFunctionTable

applications/utilities/preProcessing/wallFunctionTable/Make/options

@@ -0,0 +1,8 @@
+EXE_INC = \
+    -I./tabulatedWallFunction/lnInclude \
+    -I$(LIB_SRC)/turbulenceModels/incompressible/RAS/lnInclude \
+    -I$(LIB_SRC)/finiteVolume/lnInclude
+
+EXE_LIBS = \
+    -ltabulatedWallFunctions \
+    -lfiniteVolume

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/Make/files

@@ -0,0 +1,7 @@
+tabulatedWallFunction/tabulatedWallFunction.C
+tabulatedWallFunction/newTabulatedWallFunction.C
+
+SpaldingsLaw/SpaldingsLaw.C
+general/general.C
+
+LIB = $(FOAM_LIBBIN)/libtabulatedWallFunctions

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/Make/options

@@ -0,0 +1,9 @@
+
+EXE_INC = \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/turbulenceModels \
+    -I$(LIB_SRC)/turbulenceModels/incompressible/RAS/lnInclude \
+    -I$(LIB_SRC)/transportModels
+
+LIB_LIBS = \
+    -lfiniteVolume

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/SpaldingsLaw/SpaldingsLaw.C

@@ -0,0 +1,197 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "SpaldingsLaw.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    namespace tabulatedWallFunctions
+    {
+        defineTypeNameAndDebug(SpaldingsLaw, 0);
+        addToRunTimeSelectionTable
+        (
+            tabulatedWallFunction,
+            SpaldingsLaw,
+            dictionary
+        );
+    }
+}
+
+const Foam::label Foam::tabulatedWallFunctions::SpaldingsLaw::maxIters_ = 1000;
+
+const Foam::scalar
+    Foam::tabulatedWallFunctions::SpaldingsLaw::tolerance_ = 1e-4;
+
+
+// * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
+
+void Foam::tabulatedWallFunctions::SpaldingsLaw::invertFunction()
+{
+    // Initialise u+ and R
+    scalar Re = 0.0;
+    scalar uPlus = 1;
+
+    // Populate the table
+    for (label i=0; i<invertedTable_.size(); i++)
+    {
+        if (invertedTable_.log10())
+        {
+            Re = pow(10, (i*invertedTable_.dx() + invertedTable_.x0()));
+        }
+        else
+        {
+            Re = i*invertedTable_.dx() + invertedTable_.x0();
+        }
+
+        // Use latest available u+ estimate
+        if (i > 0)
+        {
+            uPlus = invertedTable_[i-1];
+        }
+
+        // Newton iterations to determine u+
+        label iter = 0;
+        scalar error = GREAT;
+        do
+        {
+            scalar kUPlus = min(kappa_*uPlus, 50);
+            scalar A =
+                E_*sqr(uPlus)
+              + uPlus
+               *(exp(kUPlus) - pow3(kUPlus)/6 - 0.5*sqr(kUPlus) - kUPlus - 1);
+
+            scalar f = - Re + A/E_;
+
+            scalar df =
+                (
+                    2*E_*uPlus
+                  + exp(kUPlus)*(kUPlus + 1)
+                  - 2/3*pow3(kUPlus)
+                  - 1.5*sqr(kUPlus)
+                  - 2*kUPlus
+                  - 1
+                )/E_;
+
+            scalar uPlusNew = uPlus - f/(df + ROOTVSMALL);
+            error = mag((uPlus - uPlusNew)/uPlusNew);
+            uPlus = uPlusNew;
+        } while (error > tolerance_ && ++iter < maxIters_);
+
+        if (iter == maxIters_)
+        {
+            WarningIn("SpaldingsLaw::invertFunction()")
+                << "Newton iterations not converged:" << nl
+                << "    iters = " << iter << ", error = " << error << endl;
+        }
+
+        // Set new values - constrain u+ >= 0
+        invertedTable_[i] = max(0, uPlus);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::tabulatedWallFunctions::SpaldingsLaw::SpaldingsLaw
+(
+    const dictionary& dict,
+    const polyMesh& mesh
+)
+:
+    tabulatedWallFunction(dict, mesh, typeName),
+    kappa_(readScalar(coeffDict_.lookup("kappa"))),
+    E_(readScalar(coeffDict_.lookup("E")))
+{
+    invertFunction();
+
+    if (debug)
+    {
+        writeData(Info);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::tabulatedWallFunctions::SpaldingsLaw::~SpaldingsLaw()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::scalar Foam::tabulatedWallFunctions::SpaldingsLaw::yPlus
+(
+    const scalar uPlus
+) const
+{
+    scalar kUPlus = min(kappa_*uPlus, 50);
+
+    return
+        uPlus
+      + 1/E_*(exp(kUPlus) - pow3(kUPlus)/6 - 0.5*sqr(kUPlus) - kUPlus - 1);
+}
+
+
+Foam::scalar Foam::tabulatedWallFunctions::SpaldingsLaw::Re
+(
+    const scalar uPlus
+) const
+{
+    return uPlus*yPlus(uPlus);
+}
+
+
+void Foam::tabulatedWallFunctions::SpaldingsLaw::writeData(Ostream& os) const
+{
+    if (invertedTable_.log10())
+    {
+        os  << "log10(Re), y+, u+:" << endl;
+        forAll(invertedTable_, i)
+        {
+            scalar uPlus = invertedTable_[i];
+            scalar Re = ::log10(this->Re(uPlus));
+            scalar yPlus = this->yPlus(uPlus);
+            os  << Re << ", " << yPlus << ", " << uPlus << endl;
+        }
+    }
+    else
+    {
+        os  << "Re, y+, u+:" << endl;
+        forAll(invertedTable_, i)
+        {
+            scalar uPlus = invertedTable_[i];
+            scalar Re = this->Re(uPlus);
+            scalar yPlus = this->yPlus(uPlus);
+            os  << Re << ", " << yPlus << ", " << uPlus << endl;
+        }
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/SpaldingsLaw/SpaldingsLaw.H

@@ -0,0 +1,140 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Foam::tabulatedWallFunctions::SpaldingsLaw
+
+Description
+    Computes U+ as a function of Reynolds number by inverting Spaldings law.
+
+    Example dictionary specification:
+
+        tabulatedWallFunction SpaldingsLaw;
+
+        // Output table info
+        tableName       uPlusWallFunctionData; // Output table name
+        log10           yes;                // Rey interpreted as log10(Rey)
+        dx              0.2;                // Interval log10(Rey)
+        x0              -3;                 // Minimum log10(Rey)
+        xMax            7;                  // Maximum log10(Rey)
+
+        SpaldingsLawCoeffs
+        {
+            kappa           0.41;           // Von Karman constant
+            E               9.8;            // Law-of-the-wall E coefficient
+        }
+
+
+SourceFiles
+    SpaldingsLaw.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef SpaldingsLaw_H
+#define SpaldingsLaw_H
+
+#include "tabulatedWallFunction.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace tabulatedWallFunctions
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class SpaldingsLaw Declaration
+\*---------------------------------------------------------------------------*/
+
+class SpaldingsLaw
+:
+    public tabulatedWallFunction
+{
+protected:
+
+    // Protected data
+
+        //- Von Karman constant
+        scalar kappa_;
+
+        //- Law-of-the-wall E coefficient
+        scalar E_;
+
+
+        // Newton iteration solution properties
+
+            //- Maximum number of iterations
+            static const label maxIters_;
+
+            //- Tolerance
+            static const scalar tolerance_;
+
+
+    // Protected member functions
+
+        //- Invert the function
+        virtual void invertFunction();
+
+
+public:
+
+    //- Run-time type information
+    TypeName("SpaldingsLaw");
+
+
+    // Constructors
+    SpaldingsLaw(const dictionary& dict, const polyMesh& mesh);
+
+    //- Destructor
+    virtual ~SpaldingsLaw();
+
+
+    // Member Functions
+
+        // Access
+
+            //- Return y+ as a function of u+
+            virtual scalar yPlus(const scalar uPlus) const;
+
+            //- Return Reynolds number as a function of u+
+            virtual scalar Re(const scalar uPlus) const;
+
+
+        // I-O
+
+            //- Write to Ostream
+            virtual void writeData(Ostream& os) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace tabulatedWallFunctions
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/general/general.C

@@ -0,0 +1,254 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "general.H"
+#include "addToRunTimeSelectionTable.H"
+#include "Tuple2.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    namespace tabulatedWallFunctions
+    {
+        defineTypeNameAndDebug(general, 0);
+        addToRunTimeSelectionTable
+        (
+            tabulatedWallFunction,
+            general,
+            dictionary
+        );
+    }
+}
+
+template<>
+const char*
+Foam::NamedEnum<Foam::tabulatedWallFunctions::general::interpolationType, 1>::
+    names[] =
+    {
+        "linear"
+    };
+
+const
+Foam::NamedEnum<Foam::tabulatedWallFunctions::general::interpolationType, 1>
+    Foam::tabulatedWallFunctions::general::interpolationTypeNames_;
+
+
+// * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
+
+void Foam::tabulatedWallFunctions::general::invertTable()
+{
+    scalarList Rey(uPlus_.size(), 0.0);
+
+    // Calculate Reynolds number
+    forAll(uPlus_, i)
+    {
+        Rey[i] = yPlus_[i]*uPlus_[i];
+        if (invertedTable_.log10())
+        {
+            Rey[i] = ::log10(max(ROOTVSMALL, Rey[i]));
+        }
+    }
+
+    // Populate the U+ vs Re table
+    forAll(invertedTable_, i)
+    {
+        scalar Re = i*invertedTable_.dx() + invertedTable_.x0();
+        invertedTable_[i] = max(0, interpolate(Re, Rey, uPlus_));
+    }
+}
+
+
+Foam::scalar Foam::tabulatedWallFunctions::general::interpolate
+(
+    const scalar xi,
+    const scalarList& x,
+    const scalarList& fx
+) const
+{
+    switch (interpType_)
+    {
+        case itLinear:
+        {
+            if (xi < x[0])
+            {
+                return fx[0];
+            }
+            else if (xi > x[x.size()-1])
+            {
+                return fx[x.size()-1];
+            }
+            else
+            {
+                label i2 = 0;
+                while (x[i2] < xi)
+                {
+                    i2++;
+                }
+                label i1 = i2 - 1;
+
+                return (xi - x[i1])/(x[i2] - x[i1])*(fx[i2] - fx[i1]) + fx[i1];
+            }
+
+            break;
+        }
+        default:
+        {
+            FatalErrorIn
+            (
+                "tabulatedWallFunctions::general::interpolate"
+                "("
+                    "const scalar, "
+                    "const scalarList&, "
+                    "const scalarList&"
+                ")"
+            )   << "Unknown interpolation method" << nl
+                << abort(FatalError);
+        }
+    }
+
+    return 0.0;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::tabulatedWallFunctions::general::general
+(
+    const dictionary& dict,
+    const polyMesh& mesh
+)
+:
+    tabulatedWallFunction(dict, mesh, typeName),
+    interpType_(interpolationTypeNames_[coeffDict_.lookup("interpType")]),
+    yPlus_(),
+    uPlus_(),
+    log10YPlus_(coeffDict_.lookup("log10YPlus")),
+    log10UPlus_(coeffDict_.lookup("log10UPlus"))
+{
+    List<Tuple2<scalar, scalar> > inputTable = coeffDict_.lookup("inputTable");
+    if (inputTable.size() < 2)
+    {
+        FatalErrorIn
+        (
+            "tabulatedWallFunctions::general::general"
+            "("
+                "const dictionary&, "
+                "const polyMesh&"
+            ")"
+        )   << "Input table must have at least 2 values" << nl
+            << exit(FatalError);
+    }
+
+    yPlus_.setSize(inputTable.size());
+    uPlus_.setSize(inputTable.size());
+
+    forAll(inputTable, i)
+    {
+        if (log10YPlus_)
+        {
+            yPlus_[i] = pow(10, inputTable[i].first());
+        }
+        else
+        {
+            yPlus_[i] = inputTable[i].first();
+        }
+
+        if (log10UPlus_)
+        {
+            uPlus_[i] = pow(10, inputTable[i].second());
+        }
+        else
+        {
+            uPlus_[i] = inputTable[i].second();
+        }
+    }
+
+    invertTable();
+
+    if (debug)
+    {
+        writeData(Info);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::tabulatedWallFunctions::general::~general()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+Foam::scalar Foam::tabulatedWallFunctions::general::yPlus
+(
+    const scalar uPlus
+) const
+{
+    return interpolate(uPlus, uPlus_, yPlus_);
+}
+
+
+Foam::scalar Foam::tabulatedWallFunctions::general::Re
+(
+    const scalar uPlus
+) const
+{
+    return uPlus*yPlus(uPlus);
+}
+
+
+// * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * * //
+
+void Foam::tabulatedWallFunctions::general::writeData(Ostream& os) const
+{
+    if (invertedTable_.log10())
+    {
+        os  << "log10(Re), y+, u+:" << endl;
+        forAll(invertedTable_, i)
+        {
+            scalar uPlus = invertedTable_[i];
+            scalar Re = ::log10(this->Re(uPlus));
+            scalar yPlus = this->yPlus(uPlus);
+            os  << Re << ", " << yPlus << ", " << uPlus << endl;
+        }
+    }
+    else
+    {
+        os  << "Re, y+, u+:" << endl;
+        forAll(invertedTable_, i)
+        {
+            scalar uPlus = invertedTable_[i];
+            scalar Re = this->Re(uPlus);
+            scalar yPlus = this->yPlus(uPlus);
+            os  << Re << ", " << yPlus << ", " << uPlus << endl;
+        }
+    }
+}
+
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/general/general.H

@@ -0,0 +1,172 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Foam::tabulatedWallFunctions::general
+
+Description
+    Computes U+ as a function of Reynolds number by inverting table of
+    y+ vs U+
+
+    Example dictionary specification:
+
+        tabulatedWallFunction general;
+
+        // Output table info
+        tableName       uPlusWallFunctionData; // Output table name
+        log10           yes;                // Re interpreted as log10(Rey)
+        dx              0.2;                // Interval log10(Rey)
+        x0              -3;                 // Minimum log10(Rey)
+        xMax            7;                  // Maximum log10(Rey)
+
+        generalCoeffs
+        {
+            interpType      linear;         // Interpolation method
+            log10YPlus      true;           // y+ values defined as log10(y+)
+            log10UPlus      true;           // U+ values defined as log10(y+)
+            inputTable
+            (
+                (yPlusValue0    uPlusValue0)
+                ...
+                (yPlusValueN    uPlusValueN)
+            );
+
+        }
+
+
+SourceFiles
+    general.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef general_H
+#define general_H
+
+#include "tabulatedWallFunction.H"
+#include "NamedEnum.H"
+#include "Switch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace tabulatedWallFunctions
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class general Declaration
+\*---------------------------------------------------------------------------*/
+
+class general
+:
+    public tabulatedWallFunction
+{
+public:
+
+    // Public data types
+
+        //- Enumeration listing available interpolation types
+        enum interpolationType
+        {
+            itLinear
+        };
+
+        static const NamedEnum<interpolationType, 1> interpolationTypeNames_;
+
+
+protected:
+
+    // Protected data
+
+        //- Type of interpolation to apply when inverting the data set
+        interpolationType interpType_;
+
+        //- Input y+ values
+        List<scalar> yPlus_;
+
+        //- Input U+ values
+        List<scalar> uPlus_;
+
+        //- Are y+ values entered as log10(y+)?
+        Switch log10YPlus_;
+
+        //- Are U+ values entered as log10(U+)?
+        Switch log10UPlus_;
+
+
+    // Protected member functions
+
+        //- Invert the table
+        virtual void invertTable();
+
+        //- Interpolate
+        virtual scalar interpolate
+        (
+            const scalar xi,
+            const scalarList& x,
+            const scalarList& fx
+        ) const;
+
+
+public:
+
+    //- Run-time type information
+    TypeName("general");
+
+
+    // Constructors
+    general(const dictionary& dict, const polyMesh& mesh);
+
+    //- Destructor
+    virtual ~general();
+
+
+    // Member Functions
+
+        // Access
+
+            //- Return y+ as a function of u+
+            virtual scalar yPlus(const scalar uPlus) const;
+
+            //- Return Reynolds number as a function of u+
+            virtual scalar Re(const scalar uPlus) const;
+
+
+        // I-O
+
+            //- Write to Ostream
+            virtual void writeData(Ostream& os) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace tabulatedWallFunctions
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/tabulatedWallFunction/newTabulatedWallFunction.C

@@ -0,0 +1,71 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "tabulatedWallFunction.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace tabulatedWallFunctions
+{
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+autoPtr<tabulatedWallFunction> tabulatedWallFunction::New
+(
+    const dictionary& dict,
+    const polyMesh& mesh
+)
+{
+    word twfTypeName = dict.lookup("tabulatedWallFunction");
+
+    Info<< "Selecting tabulatedWallFunction " << twfTypeName << endl;
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(twfTypeName);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "tabulatedWallFunction::New(const dictionary&, const polyMesh&)"
+        )   << "Unknown tabulatedWallFunction type " << twfTypeName
+            << nl << nl << "Valid tabulatedWallFunction types are:" << nl
+            << dictionaryConstructorTablePtr_->toc()
+            << exit(FatalError);
+    }
+
+    return autoPtr<tabulatedWallFunction>(cstrIter()(dict, mesh));
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace tabulatedWallFunctions
+} // End namespace Foam
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/tabulatedWallFunction/tabulatedWallFunction.C

@@ -0,0 +1,88 @@
+
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "tabulatedWallFunction.H"
+#include "Time.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    namespace tabulatedWallFunctions
+    {
+        defineTypeNameAndDebug(tabulatedWallFunction, 0);
+        defineRunTimeSelectionTable(tabulatedWallFunction, dictionary);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::tabulatedWallFunctions::tabulatedWallFunction::tabulatedWallFunction
+(
+    const dictionary& dict,
+    const polyMesh& mesh,
+    const word& name
+)
+:
+    dict_(dict),
+    mesh_(mesh),
+    coeffDict_(dict.subDict(name + "Coeffs")),
+    invertedTableName_(dict.lookup("invertedTableName")),
+    invertedTable_(invertedTableName_, mesh_, dict, true)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::tabulatedWallFunctions::tabulatedWallFunction::~tabulatedWallFunction()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+void Foam::tabulatedWallFunctions::tabulatedWallFunction::write()
+{
+    if (invertedTable_.log10())
+    {
+        invertedTable_.note() =
+            "U+ as a function of log10(Re) computed using " + type();
+    }
+    else
+    {
+        invertedTable_.note() =
+            "U+ as a function of Re computed using " + type();
+    }
+
+    Info<< "Writing inverted table to\n    " << invertedTable_.objectPath()
+        << endl;
+
+    invertedTable_.write();
+}
+
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/tabulatedWallFunction/tabulatedWallFunction.H

@@ -0,0 +1,148 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Foam::tabulatedWallFunctions::tabulatedWallFunction
+
+Description
+    Base class for models that generate tabulated wall function data.
+
+SourceFiles
+    tabulatedWallFunction.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef tabulatedWallFunction_H
+#define tabulatedWallFunction_H
+
+#include "dictionary.H"
+#include "polyMesh.H"
+#include "runTimeSelectionTables.H"
+#include "Switch.H"
+#include "uniformInterpolationTable.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+namespace tabulatedWallFunctions
+{
+
+/*---------------------------------------------------------------------------*\
+                   Class tabulatedWallFunction Declaration
+\*---------------------------------------------------------------------------*/
+
+class tabulatedWallFunction
+{
+protected:
+
+    // Proteced data
+
+        //- Main dictionary
+        const dictionary dict_;
+
+        //- Reference to the mesh database
+        const polyMesh& mesh_;
+
+        //- Coefficients dictionary
+        const dictionary coeffDict_;
+
+        //- Name of inverted table
+        word invertedTableName_;
+
+        //- Inverted table
+        uniformInterpolationTable<scalar> invertedTable_;
+
+
+public:
+
+    //- Run-time type information
+    TypeName("tabulatedWallFunction");
+
+    //- Declare runtime constructor selection table
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        tabulatedWallFunction,
+        dictionary,
+        (
+            const dictionary& dict,
+            const polyMesh& mesh
+        ),
+        (dict, mesh)
+    );
+
+    //- Constructor
+    tabulatedWallFunction
+    (
+        const dictionary& dict,
+        const polyMesh& mesh,
+        const word& name
+    );
+
+    //- Destructor
+    virtual ~tabulatedWallFunction();
+
+
+    //- Selector
+    static autoPtr<tabulatedWallFunction> New
+    (
+        const dictionary& dict,
+        const polyMesh& mesh
+    );
+
+
+    // Member Functions
+
+        // Access
+
+            //- Return the inverted table name
+            inline const word& invertedTableName() const;
+
+            //- Return the inverted table
+            inline const uniformInterpolationTable<scalar>&
+                invertedTable() const;
+
+
+        // I-O
+
+            //- Write
+            virtual void write();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace tabulatedWallFunctions
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "tabulatedWallFunctionI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/tabulatedWallFunction/tabulatedWallFunction/tabulatedWallFunctionI.H

@@ -0,0 +1,45 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "tabulatedWallFunction.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+inline const Foam::word&
+Foam::tabulatedWallFunctions::tabulatedWallFunction::invertedTableName() const
+{
+    return invertedTableName_;
+}
+
+
+inline const Foam::uniformInterpolationTable<Foam::scalar>&
+Foam::tabulatedWallFunctions::tabulatedWallFunction::invertedTable() const
+{
+    return invertedTable_;
+}
+
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/wallFunctionDict

@@ -0,0 +1,39 @@
+/*--------------------------------*- C++ -*----------------------------------*\
+| =========                 |                                                 |
+| \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |
+|  \\    /   O peration     | Version:  1.6                                   |
+|   \\  /    A nd           | Web:      www.OpenFOAM.org                      |
+|    \\/     M anipulation  |                                                 |
+\*---------------------------------------------------------------------------*/
+FoamFile
+{
+    version     2.0;
+    format      ascii;
+    class       dictionary;
+    location    "constant";
+    object      wallFunctionDict;
+}
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+tabulatedWallFunction SpaldingsLaw;
+
+invertedTableName uPlusWallFunctionData;
+
+dx              0.2;
+
+x0              -3;
+
+xMax            7;
+
+log10           yes;
+
+bound           yes;
+
+SpaldingsLawCoeffs
+{
+    kappa           0.41;
+    E               9.8;
+}
+
+
+// ************************************************************************* //

applications/utilities/preProcessing/wallFunctionTable/wallFunctionTableApp.C

@@ -0,0 +1,68 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Application
+    wallFunctionTable
+
+Description
+    Generates a table suitable for use by tabulated wall functions
+
+\*---------------------------------------------------------------------------*/
+
+#include "fvCFD.H"
+#include "tabulatedWallFunction.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+int main(int argc, char *argv[])
+{
+    #include "setRootCase.H"
+    #include "createTime.H"
+    #include "createMesh.H"
+
+    IOdictionary dict
+    (
+        IOobject
+        (
+            "wallFunctionDict",
+            runTime.constant(),
+            mesh,
+            IOobject::MUST_READ
+        )
+    );
+
+    autoPtr<tabulatedWallFunctions::tabulatedWallFunction>
+        twf(tabulatedWallFunctions::tabulatedWallFunction::New(dict, mesh));
+
+//    twf->writeData(Info);
+
+    twf->write();
+
+    Info << "End\n" << endl;
+
+    return 0;
+}
+
+
+// ************************************************************************* //

src/OpenFOAM/db/scalarRange/scalarRange.C

@@ -42,7 +42,7 @@ Foam::scalarRange::scalarRange()
 {}
 
 
-Foam::scalarRange::scalarRange(const scalar& lower, const scalar& upper)
+Foam::scalarRange::scalarRange(const scalar lower, const scalar upper)
 :
     type_(RANGE),
     value_(lower),
@@ -123,7 +123,7 @@ Foam::scalar Foam::scalarRange::upper() const
 }
 
 
-bool Foam::scalarRange::selected(const scalar& value) const
+bool Foam::scalarRange::selected(const scalar value) const
 {
     switch (type_)
     {

src/OpenFOAM/db/scalarRange/scalarRange.H

@@ -92,7 +92,7 @@ public:
         scalarRange();
 
         //- Construct a Range
-        scalarRange(const scalar& lower, const scalar& upper);
+        scalarRange(const scalar lower, const scalar upper);
 
         //- Construct from Istream.
         //  Since commas can be used as list delimiters,
@@ -119,7 +119,7 @@ public:
         scalar upper() const;
 
         //- Return true if the value is within the range
-        bool selected(const scalar&) const;
+        bool selected(const scalar) const;
 
 
     // Member Operators

src/OpenFOAM/db/scalarRange/scalarRanges.C

@@ -57,7 +57,7 @@ Foam::scalarRanges::scalarRanges(Istream& is)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-bool Foam::scalarRanges::selected(const scalar& value) const
+bool Foam::scalarRanges::selected(const scalar value) const
 {
     forAll(*this, i)
     {

src/OpenFOAM/db/scalarRange/scalarRanges.H

@@ -66,7 +66,7 @@ public:
     // Member Functions
 
         //- Return true if the given value is within the ranges
-        bool selected(const scalar&) const;
+        bool selected(const scalar) const;
 
         //- Return the set of selected entries in the given list
         //  that are within the ranges

src/OpenFOAM/global/unitConversion/unitConversion.H

@@ -22,7 +22,7 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Namespace
+InNamespace
     Foam
 
 Description
@@ -35,8 +35,6 @@ Description
 
 #include "mathematicalConstants.H"
 
-using namespace Foam::constant::mathematical;
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
@@ -45,15 +43,15 @@ namespace Foam
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 //- Conversion from degrees to radians
-inline scalar degToRad(const scalar& deg)
+inline scalar degToRad(const scalar deg)
 {
-    return (deg*pi/180.0);
+    return (deg * constant::mathematical::pi/180.0);
 }
 
 //- Conversion from radians to degrees
-inline scalar radToDeg(const scalar& rad)
+inline scalar radToDeg(const scalar rad)
 {
-    return (rad*180.0/pi);
+    return (rad * 180.0/constant::mathematical::pi);
 }
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/OpenFOAM/interpolations/uniformInterpolationTable/uniformInterpolationTable.C

@@ -0,0 +1,236 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+#include "uniformInterpolationTable.H"
+#include "Time.H"
+
+// * * * * * * * * * * * *  Private Member Functions * * * * * * * * * * * * //
+
+template <class Type>
+void Foam::uniformInterpolationTable<Type>::checkTable() const
+{
+    if (size() < 2)
+    {
+        FatalErrorIn("uniformInterpolationTable<Type>::checkTable()")
+            << "Table " << name() << ": must have at least 2 values." << nl
+            << "Table size = " << size() << nl
+            << "    min, interval width = " << x0_ << ", " << dx_ << nl
+            << exit(FatalError);
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template <class Type>
+Foam::uniformInterpolationTable<Type>::uniformInterpolationTable
+(
+    const IOobject& io,
+    bool readFields
+)
+:
+    IOobject(io),
+    List<scalar>(2, 0.0),
+    x0_(0.0),
+    dx_(1.0),
+    log10_(false),
+    bound_(false)
+{
+    if (readFields)
+    {
+        IOdictionary dict(io);
+
+        dict.lookup("data") >> *this;
+        dict.lookup("x0") >> x0_;
+        dict.lookup("dx") >> dx_;
+        dict.lookup("log10") >> log10_;
+        dict.lookup("bound") >> bound_;
+    }
+
+    checkTable();
+}
+
+
+template <class Type>
+Foam::uniformInterpolationTable<Type>::uniformInterpolationTable
+(
+    const word& tableName,
+    const objectRegistry& db,
+    const dictionary& dict,
+    const bool initialiseOnly
+)
+:
+    IOobject
+    (
+        tableName,
+        db.time().constant(),
+        db,
+        IOobject::NO_READ,
+        IOobject::NO_WRITE,
+        false // if used in BCs, could be used by multiple patches
+    ),
+    List<scalar>(2, 0.0),
+    x0_(readScalar(dict.lookup("x0"))),
+    dx_(readScalar(dict.lookup("dx"))),
+    log10_(dict.lookup("log10")),
+    bound_(dict.lookup("bound"))
+{
+    if (initialiseOnly)
+    {
+        scalar xMax = readScalar(dict.lookup("xMax"));
+        label nIntervals = static_cast<label>(xMax - x0_)/dx_ + 1;
+        this->setSize(nIntervals);
+    }
+    else
+    {
+        dict.lookup("data") >> *this;
+    }
+
+    checkTable();
+}
+
+
+template <class Type>
+Foam::uniformInterpolationTable<Type>::uniformInterpolationTable
+(
+    const uniformInterpolationTable& uit
+)
+:
+    IOobject(uit),
+    List<scalar>(uit),
+    x0_(uit.x0_),
+    dx_(uit.dx_),
+    log10_(uit.log10_),
+    bound_(uit.bound_)
+{
+    checkTable();
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+template <class Type>
+Foam::uniformInterpolationTable<Type>::~uniformInterpolationTable()
+{}
+
+
+// * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * * //
+
+template <class Type>
+Type Foam::uniformInterpolationTable<Type>::interpolate(scalar x) const
+{
+    if (bound_)
+    {
+        x = max(min(xMax() - SMALL*dx_, x), x0_);
+    }
+    else
+    {
+        if (x < x0_)
+        {
+            FatalErrorIn
+            (
+                "uniformInterpolationTable<Type>::interpolate(scalar x)"
+            )   << "Supplied value is less than minimum table value:" << nl
+                << "xMin=" << x0_ << ", xMax=" << xMax() << ", x=" << x << nl
+                << exit(FatalError);
+        }
+
+        if (x > xMax())
+        {
+            FatalErrorIn
+            (
+                "uniformInterpolationTable<Type>::interpolate(scalar x)"
+            )   << "Supplied value is greater than maximum table value:" << nl
+                << "xMin=" << x0_ << ", xMax=" << xMax() << ", x=" << x << nl
+                << exit(FatalError);
+        }
+    }
+
+    label i = static_cast<label>((x - x0_)/dx_);
+
+    scalar xLo = x0_ + i*dx_;
+
+    Type fx = (x - xLo)/dx_*(operator[](i+1) - operator[](i)) + operator[](i);
+
+    if (debug)
+    {
+        Info<< "Table: " << name() << ", x=" << x
+            << ", x_lo=" << xLo << ", x_hi=" << xLo + dx_
+            << ", f(x_lo)=" << operator[](i) << ", f(x_hi)=" << operator[](i+1)
+            << ", f(x)=" << fx << endl;
+    }
+
+    return fx;
+}
+
+
+template <class Type>
+Type Foam::uniformInterpolationTable<Type>::interpolateLog10
+(
+    scalar x
+) const
+{
+    if (log10_)
+    {
+        if (x > 0)
+        {
+            x = ::log10(x);
+        }
+        else if (bound_ && (x <= 0))
+        {
+            x = x0_;
+        }
+        else
+        {
+            FatalErrorIn
+            (
+                "uniformInterpolationTable<Type>::interpolateLog10(scalar x)"
+            )   << "Table " << name() << nl
+                << "Supplied value must be greater than 0 when in log10 mode"
+                << nl << "x=" << x << nl << exit(FatalError);
+        }
+    }
+
+    return interpolate(x);
+}
+
+
+template <class Type>
+void Foam::uniformInterpolationTable<Type>::write() const
+{
+    IOdictionary dict(*this);
+
+    dict.add("data", static_cast<const List<scalar>&>(*this));
+    dict.add("x0", x0_);
+    dict.add("dx", dx_);
+    dict.add("log10", log10_);
+    dict.add("bound", bound_);
+
+    dict.regIOobject::write();
+}
+
+
+// ************************************************************************* //

src/OpenFOAM/interpolations/uniformInterpolationTable/uniformInterpolationTable.H

@@ -0,0 +1,221 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Foam::uniformInterpolationTable
+
+Description
+    Table with uniform interval in independant variable, with linear
+    interpolation
+
+    Example usage (scalar): values specified in a dictionary:
+
+    {
+        x0              0;      // lower limit
+        dx              0.2;    // fixed interval
+        log10           true;   // take log(10) when interpolating?
+        data                    // list of dependent data values
+        (
+            7870                // value at x0
+            7870                // value at x0 + dx
+            ...
+            7870                // value at x0 + n*dx
+        );
+    }
+
+SourceFiles
+    uniformInterpolationTable.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef uniformInterpolationTable_H
+#define uniformInterpolationTable_H
+
+#include "List.H"
+#include "Switch.H"
+#include "IOobject.H"
+#include "objectRegistry.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                 Class uniformInterpolationTable Declaration
+\*---------------------------------------------------------------------------*/
+
+template <class Type>
+class uniformInterpolationTable
+:
+    public IOobject,
+    public List<Type>
+{
+    // Private data
+
+        // Control parameetrs
+
+            //- Lower limit
+            scalar x0_;
+
+            //- Fixed interval
+            scalar dx_;
+
+            //- Flag to indicate that x data is given in log10(x) form
+            Switch log10_;
+
+            //- Bound x values
+            Switch bound_;
+
+
+    // Private Member Functions
+
+        //- Check that the table is valid
+        void checkTable() const;
+
+        //- Disallow default bitwise assignment
+        void operator=(const uniformInterpolationTable&);
+
+
+public:
+
+    // Constructors
+
+        //- Construct from IOobject and readFields flag. Creates a null object
+        //  if readFields = false
+        uniformInterpolationTable(const IOobject& io, const bool readFields);
+
+        //- Construct from name, objectRegistry and dictionary.
+        //  If initialiseOnly flag is set, control parameters are read from
+        //  the dictionary, but not the data table
+        uniformInterpolationTable
+        (
+            const word& tableName,
+            const objectRegistry& db,
+            const dictionary& dict,
+            const bool initialiseOnly = false
+        );
+
+        //- Construct as copy
+        uniformInterpolationTable(const uniformInterpolationTable& uit);
+
+
+    //- Destructor
+    ~uniformInterpolationTable();
+
+
+    // Member Functions
+
+        // Access
+
+            //- Return the lower limit
+            inline scalar x0() const;
+
+            //- Return the fixed interval
+            inline scalar dx() const;
+
+            //- Return the log10(x) flag
+            inline const Switch& log10() const;
+
+            //- Return the bound flag
+            inline const Switch& bound() const;
+
+
+        // Edit
+
+            //- Return the lower limit
+            inline scalar& x0();
+
+            //- Return the fixed interval
+            inline scalar& dx();
+
+            //- Return the log10(x) flag
+            inline Switch& log10();
+
+            //- Return the bound flag
+            inline Switch& bound();
+
+
+        // Evaluation
+
+            //- Return the minimum x value
+            inline scalar xMin() const;
+
+            //- Return the maximum x value
+            inline scalar xMax() const;
+
+            //- Interpolate
+            Type interpolate(scalar x) const;
+
+            //- Interpolate - takes log10 flag into account
+            Type interpolateLog10(scalar x) const;
+
+
+        // Override ancestor size() function and [] operator
+
+            //- Return the size of the table
+            label size() const
+            {
+                return List<Type>::size();
+            }
+
+            //- Use List[] operator for read access
+            Type operator[](label x) const
+            {
+                return List<Type>::operator[](x);
+            }
+
+            //- Use List[] operator for write access
+            Type& operator[](label x)
+            {
+                return List<Type>::operator[](x);
+            }
+
+
+        // I-O
+
+            //- Write
+            void write() const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "uniformInterpolationTableI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "uniformInterpolationTable.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/OpenFOAM/interpolations/uniformInterpolationTable/uniformInterpolationTableI.H

@@ -0,0 +1,97 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2009-2009 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 2 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, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+template <class Type>
+Foam::scalar Foam::uniformInterpolationTable<Type>::x0() const
+{
+    return x0_;
+}
+
+
+template <class Type>
+Foam::scalar Foam::uniformInterpolationTable<Type>::dx() const
+{
+    return dx_;
+}
+
+
+template <class Type>
+const Foam::Switch& Foam::uniformInterpolationTable<Type>::log10() const
+{
+    return log10_;
+}
+
+
+template <class Type>
+const Foam::Switch& Foam::uniformInterpolationTable<Type>::bound() const
+{
+    return bound_;
+}
+
+
+template <class Type>
+Foam::scalar& Foam::uniformInterpolationTable<Type>::x0()
+{
+    return x0_;
+}
+
+
+template <class Type>
+Foam::scalar& Foam::uniformInterpolationTable<Type>::dx()
+{
+    return dx_;
+}
+
+
+template <class Type>
+Foam::Switch& Foam::uniformInterpolationTable<Type>::log10()
+{
+    return log10_;
+}
+
+
+template <class Type>
+Foam::Switch& Foam::uniformInterpolationTable<Type>::bound()
+{
+    return bound_;
+}
+
+
+template <class Type>
+Foam::scalar Foam::uniformInterpolationTable<Type>::xMin() const
+{
+    return x0_;
+}
+
+
+template <class Type>
+Foam::scalar Foam::uniformInterpolationTable<Type>::xMax() const
+{
+    return x0_ + dx_*(size() - 1);
+}
+
+
+// ************************************************************************* //

src/OpenFOAM/matrices/simpleMatrix/simpleMatrix.C

@@ -32,7 +32,20 @@ template<class Type>
 Foam::simpleMatrix<Type>::simpleMatrix(const label mSize)
 :
     scalarSquareMatrix(mSize),
-    source_(mSize, pTraits<Type>::zero)
+    source_(mSize)
+{}
+
+
+template<class Type>
+Foam::simpleMatrix<Type>::simpleMatrix
+(
+    const label mSize,
+    const scalar coeffVal,
+    const Type& sourceVal
+)
+:
+    scalarSquareMatrix(mSize, mSize, coeffVal),
+    source_(mSize, sourceVal)
 {}
 
 

src/OpenFOAM/matrices/simpleMatrix/simpleMatrix.H

@@ -26,7 +26,7 @@ Class
     Foam::simpleMatrix
 
 Description
-    Foam::simpleMatrix
+    A simple square matrix solver with scalar coefficients.
 
 SourceFiles
     simpleMatrix.C
@@ -75,8 +75,12 @@ public:
     // Constructors
 
         //- Construct given size
+        //  Note: this does not initialise the coefficients or the source.
         simpleMatrix(const label);
 
+        //- Construct given size and initial values for coefficients and source
+        simpleMatrix(const label, const scalar, const Type&);
+
         //- Construct from components
         simpleMatrix(const scalarSquareMatrix&, const Field<Type>&);
 
@@ -91,11 +95,13 @@ public:
 
         // Access
 
+            //- Return access to the source
             Field<Type>& source()
             {
                 return source_;
             }
 
+            //- Return const-access to the source
             const Field<Type>& source() const
             {
                 return source_;

src/OpenFOAM/matrices/solution/solution.C

@@ -60,6 +60,7 @@ Foam::solution::solution(const objectRegistry& obr, const fileName& dictName)
         )
     ),
     cache_(ITstream("cache", tokenList())()),
+    caching_(false),
     relaxationFactors_(ITstream("relaxationFactors", tokenList())()),
     defaultRelaxationFactor_(0),
     solvers_(ITstream("solvers", tokenList())())
@@ -150,12 +151,19 @@ Foam::label Foam::solution::upgradeSolverDict
 
 bool Foam::solution::cache(const word& name) const
 {
-    if (debug)
+    if (caching_)
     {
-        Info<< "Find cache entry for " << name << endl;
-    }
+        if (debug)
+        {
+            Info<< "Cache: find entry for " << name << endl;
+        }
 
-    return cache_.found(name);
+        return cache_.found(name);
+    }
+    else
+    {
+        return false;
+    }
 }
 
 
@@ -248,6 +256,7 @@ bool Foam::solution::read()
         if (dict.found("cache"))
         {
             cache_ = dict.subDict("cache");
+            caching_ = cache_.lookupOrDefault<Switch>("active", true);
         }
 
         if (dict.found("relaxationFactors"))

src/OpenFOAM/matrices/solution/solution.H

@@ -37,6 +37,7 @@ SourceFiles
 #define solution_H
 
 #include "IOdictionary.H"
+#include "Switch.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -56,6 +57,9 @@ class solution
         //- Dictionary of temporary fields to cache
         dictionary cache_;
 
+        //- Switch for the caching mechanism
+        Switch caching_;
+
         //- Dictionary of relaxation factors for all the fields
         dictionary relaxationFactors_;
 

src/OpenFOAM/meshes/primitiveShapes/objectHit/objectHit.H

@@ -76,7 +76,7 @@ public:
         {}
 
         //- Construct from components
-        objectHit(const bool success, const label& obj)
+        objectHit(const bool success, const label obj)
         :
             hit_(success),
             hitObject_(obj)
@@ -111,7 +111,7 @@ public:
         {
             return ((a.hit_ == b.hit_) && (a.hitObject_ == b.hitObject_));
         }
-            
+
         friend bool operator!=(const objectHit& a, const objectHit& b)
         {
             return (!(a == b));

src/OpenFOAM/primitives/hashes/Hash/Hash.H

@@ -107,7 +107,7 @@ public:
     //- Incrementally hash a label.
     //  This will necessarily return a different value than the
     //  non-incremental version.
-    unsigned operator()(const label& p, unsigned seed) const
+    unsigned operator()(const label p, unsigned seed) const
     {
         return Hasher(&p, sizeof(label), seed);
     }
@@ -115,11 +115,10 @@ public:
     //- Return the unsigned representation of a label.
     //  This helps if people have relied on the hash value corresponding to
     //  the natural order.
-    unsigned operator()(const label& p) const
+    unsigned operator()(const label p) const
     {
         return p;
     }
-
 };
 
 

src/OpenFOAM/primitives/random/Random.C

@@ -47,8 +47,7 @@ namespace Foam
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-// construct given seed
-Random::Random(const label& seed)
+Random::Random(const label seed)
 {
     if (seed > 1)
     {

src/OpenFOAM/primitives/random/Random.H

@@ -61,7 +61,7 @@ public:
     // Constructors
 
         //- Construct given seed
-        Random(const label&);
+        Random(const label);
 
 
     // Member functions
@@ -70,14 +70,19 @@ public:
 
         //- scalar [0..1] (so including 0,1)
         scalar scalar01();
+
         //- vector with every component scalar01
         vector vector01();
+
         //- sphericalTensor with every component scalar01
         sphericalTensor sphericalTensor01();
+
         //- symmTensor with every component scalar01
         symmTensor symmTensor01();
+
         //- tensor with every component scalar01
         tensor tensor01();
+
         //- label [lower..upper]
         label integer(const label lower, const label upper);
 

src/conversion/ensight/file/ensightFile.C

@@ -74,7 +74,7 @@ bool Foam::ensightFile::allowUndef(bool value)
 }
 
 
-Foam::scalar Foam::ensightFile::undefValue(const scalar& value)
+Foam::scalar Foam::ensightFile::undefValue(const scalar value)
 {
     // enable its use too
     allowUndef_ = true;
@@ -133,7 +133,7 @@ Foam::Ostream& Foam::ensightFile::write(const string& value)
 }
 
 
-Foam::Ostream& Foam::ensightFile::write(const label& value)
+Foam::Ostream& Foam::ensightFile::write(const label value)
 {
     if (format() == IOstream::BINARY)
     {
@@ -157,7 +157,7 @@ Foam::Ostream& Foam::ensightFile::write(const label& value)
 
 Foam::Ostream& Foam::ensightFile::write
 (
-    const label& value,
+    const label value,
     const label fieldWidth
 )
 {
@@ -181,7 +181,7 @@ Foam::Ostream& Foam::ensightFile::write
 }
 
 
-Foam::Ostream& Foam::ensightFile::write(const scalar& value)
+Foam::Ostream& Foam::ensightFile::write(const scalar value)
 {
     if (format() == IOstream::BINARY)
     {

src/conversion/ensight/file/ensightFile.H

@@ -103,7 +103,7 @@ public:
         //- Assign the value to represent undef in the results
         //  Returns the previous value
         //  NB: do not use values larger than floatScalarVGREAT
-        static scalar undefValue(const scalar&);
+        static scalar undefValue(const scalar);
 
 
     // Output
@@ -124,13 +124,13 @@ public:
         Ostream& write(const string& value);
 
         //- write integer as "%10d" or as binary
-        Ostream& write(const label& value);
+        Ostream& write(const label value);
 
         //- write integer with specified width or as binary
-        Ostream& write(const label& value, const label fieldWidth);
+        Ostream& write(const label value, const label fieldWidth);
 
         //- write float as "%12.5e" or as binary
-        Ostream& write(const scalar& value);
+        Ostream& write(const scalar value);
 
         //- Add carriage return to ascii stream
         void newline();

src/conversion/meshTables/cellTable.C

@@ -22,8 +22,6 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Description
-
 \*---------------------------------------------------------------------------*/
 
 #include "cellTable.H"
@@ -81,7 +79,7 @@ void Foam::cellTable::addDefaults()
 
 void Foam::cellTable::setEntry
 (
-    const label& id,
+    const label id,
     const word& keyWord,
     const word& value
 )
@@ -192,7 +190,7 @@ Foam::Map<Foam::word> Foam::cellTable::names
 }
 
 
-Foam::word Foam::cellTable::name(const label& id) const
+Foam::word Foam::cellTable::name(const label id) const
 {
     word theName("cellTable_" + Foam::name(id));
 
@@ -289,19 +287,19 @@ Foam::Map<Foam::word> Foam::cellTable::shells() const
 
 
 
-void Foam::cellTable::setMaterial(const label& id, const word& matlType)
+void Foam::cellTable::setMaterial(const label id, const word& matlType)
 {
     setEntry(id, "MaterialType", matlType);
 }
 
 
-void Foam::cellTable::setName(const label& id, const word& name)
+void Foam::cellTable::setName(const label id, const word& name)
 {
     setEntry(id, "Label", name);
 }
 
 
-void Foam::cellTable::setName(const label& id)
+void Foam::cellTable::setName(const label id)
 {
     iterator iter = find(id);
 

src/conversion/meshTables/cellTable.H

@@ -96,7 +96,7 @@ class cellTable
         //- Add required entries - MaterialType
         void addDefaults();
 
-        void setEntry(const label& id, const word& keyWord, const word& value);
+        void setEntry(const label id, const word& keyWord, const word& value);
 
         //- Disallow default bitwise copy construct
         cellTable(const cellTable&);
@@ -133,7 +133,7 @@ public:
 
         //- Return the name corresponding to id
         //  returns cellTable_ID if not otherwise defined
-        word name(const label& id) const;
+        word name(const label id) const;
 
         //- Return a Map of (id => name)
         Map<word> names() const;
@@ -157,13 +157,13 @@ public:
         Map<word> materialTypes() const;
 
         //- Assign material Type
-        void setMaterial(const label&, const word&);
+        void setMaterial(const label, const word&);
 
         //- Assign name
-        void setName(const label&, const word&);
+        void setName(const label, const word&);
 
         //- Assign default name if not already set
-        void setName(const label&);
+        void setName(const label);
 
         //- Read constant/cellTable
         void readDict

src/conversion/meshWriter/meshWriter.H

@@ -182,7 +182,7 @@ public:
             virtual bool writeSurface
             (
                 const fileName& timeName = fileName::null,
-                const bool& triangulate = false
+                const bool triangulate = false
             ) const
             {
                 return false;

src/conversion/meshWriter/starcd/STARCDMeshWriter.C

@@ -530,7 +530,7 @@ bool Foam::meshWriters::STARCD::write(const fileName& meshName) const
 bool Foam::meshWriters::STARCD::writeSurface
 (
     const fileName& meshName,
-    const bool& triangulate
+    const bool triangulate
 ) const
 {
     fileName baseName(meshName);

src/conversion/meshWriter/starcd/STARCDMeshWriter.H

@@ -135,7 +135,7 @@ public:
             virtual bool writeSurface
             (
                 const fileName& meshName = fileName::null,
-                const bool& triangulate = false
+                const bool triangulate = false
             ) const;
 };
 

src/dynamicMesh/polyMeshAdder/faceCoupleInfo.C

@@ -1536,7 +1536,7 @@ void Foam::faceCoupleInfo::perfectPointMatch
         FatalErrorIn
         (
             "faceCoupleInfo::perfectPointMatch"
-            "(const scalar&, const bool)"
+            "(const scalar, const bool)"
         )   << "Did not match all of the master faces to the slave faces"
             << endl
             << "This usually means that the slave patch and master patch"

src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C

@@ -55,9 +55,9 @@ namespace Foam
     class ifEqEqOp
     {
         public:
-        void operator()(label& x, const label& y) const
+        void operator()(label& x, const label y) const
         {
-            x = (x==y) ? x : value;
+            x = (x == y) ? x : value;
         }
     };
 }

src/finiteVolume/Make/files

@@ -285,6 +285,7 @@ $(divSchemes)/gaussDivScheme/gaussDivSchemes.C
 gradSchemes = finiteVolume/gradSchemes
 $(gradSchemes)/gradScheme/gradSchemes.C
 $(gradSchemes)/gaussGrad/gaussGrads.C
+
 $(gradSchemes)/leastSquaresGrad/leastSquaresVectors.C
 $(gradSchemes)/leastSquaresGrad/leastSquaresGrads.C
 $(gradSchemes)/extendedLeastSquaresGrad/extendedLeastSquaresVectors.C

src/finiteVolume/cfdTools/general/porousMedia/porousZone.H

@@ -273,7 +273,7 @@ public:
             }
 
             //- Return porosity
-            const scalar& porosity() const
+            scalar porosity() const
             {
                 return porosity_;
             }

src/finiteVolume/fields/fvPatchFields/derived/rotatingWallVelocity/rotatingWallVelocityFvPatchVectorField.H

@@ -149,7 +149,7 @@ public:
             }
 
             //- Return the rotational speed
-            const scalar& omega() const
+            scalar omega() const
             {
                 return omega_;
             }

src/finiteVolume/finiteVolume/convectionSchemes/convectionScheme/convectionScheme.H

@@ -177,7 +177,7 @@ public:
         virtual tmp<fvMatrix<Type> > fvmDiv
         (
             const surfaceScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) const = 0;
 
         virtual tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDiv

src/finiteVolume/finiteVolume/convectionSchemes/gaussConvectionScheme/gaussConvectionScheme.C

@@ -69,7 +69,7 @@ tmp<fvMatrix<Type> >
 gaussConvectionScheme<Type>::fvmDiv
 (
     const surfaceScalarField& faceFlux,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 ) const
 {
     tmp<surfaceScalarField> tweights = tinterpScheme_().weights(vf);
@@ -89,9 +89,9 @@ gaussConvectionScheme<Type>::fvmDiv
     fvm.upper() = fvm.lower() + faceFlux.internalField();
     fvm.negSumDiag();
 
-    forAll(fvm.psi().boundaryField(), patchI)
+    forAll(vf.boundaryField(), patchI)
     {
-        const fvPatchField<Type>& psf = fvm.psi().boundaryField()[patchI];
+        const fvPatchField<Type>& psf = vf.boundaryField()[patchI];
         const fvsPatchScalarField& patchFlux = faceFlux.boundaryField()[patchI];
         const fvsPatchScalarField& pw = weights.boundaryField()[patchI];
 

src/finiteVolume/finiteVolume/convectionSchemes/gaussConvectionScheme/gaussConvectionScheme.H

@@ -124,7 +124,7 @@ public:
         tmp<fvMatrix<Type> > fvmDiv
         (
             const surfaceScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) const;
 
         tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDiv

src/finiteVolume/finiteVolume/convectionSchemes/multivariateGaussConvectionScheme/multivariateGaussConvectionScheme.C

@@ -22,8 +22,6 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Description
-    
 \*---------------------------------------------------------------------------*/
 
 #include "multivariateGaussConvectionScheme.H"
@@ -81,7 +79,7 @@ tmp<fvMatrix<Type> >
 multivariateGaussConvectionScheme<Type>::fvmDiv
 (
     const surfaceScalarField& faceFlux,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 ) const
 {
     return gaussConvectionScheme<Type>

src/finiteVolume/finiteVolume/convectionSchemes/multivariateGaussConvectionScheme/multivariateGaussConvectionScheme.H

@@ -114,7 +114,7 @@ public:
         tmp<fvMatrix<Type> > fvmDiv
         (
             const surfaceScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) const;
 
         tmp<GeometricField<Type, fvPatchField, volMesh> > fvcDiv

src/finiteVolume/finiteVolume/d2dt2Schemes/EulerD2dt2Scheme/EulerD2dt2Scheme.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "EulerD2dt2Scheme.H"
@@ -181,7 +181,7 @@ EulerD2dt2Scheme<Type>::fvcD2dt2
                     coefft
                    *(rho.boundaryField() + rho.oldTime().boundaryField())
                    *vf.boundaryField()
-                        
+
                   - (
                         coefft
                        *(
@@ -232,7 +232,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 EulerD2dt2Scheme<Type>::fvmD2dt2
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -292,7 +292,7 @@ tmp<fvMatrix<Type> >
 EulerD2dt2Scheme<Type>::fvmD2dt2
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -353,7 +353,7 @@ tmp<fvMatrix<Type> >
 EulerD2dt2Scheme<Type>::fvmD2dt2
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm

src/finiteVolume/finiteVolume/d2dt2Schemes/EulerD2dt2Scheme/EulerD2dt2Scheme.H

@@ -109,19 +109,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmD2dt2
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmD2dt2
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmD2dt2
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 };
 

src/finiteVolume/finiteVolume/d2dt2Schemes/d2dt2Scheme/d2dt2Scheme.H

@@ -153,19 +153,19 @@ public:
 
         virtual tmp<fvMatrix<Type> > fvmD2dt2
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
         virtual tmp<fvMatrix<Type> > fvmD2dt2
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
         virtual tmp<fvMatrix<Type> > fvmD2dt2
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 };
 

src/finiteVolume/finiteVolume/d2dt2Schemes/steadyStateD2dt2Scheme/steadyStateD2dt2Scheme.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "steadyStateD2dt2Scheme.H"
@@ -107,7 +107,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 steadyStateD2dt2Scheme<Type>::fvmD2dt2
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -128,7 +128,7 @@ tmp<fvMatrix<Type> >
 steadyStateD2dt2Scheme<Type>::fvmD2dt2
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -149,7 +149,7 @@ tmp<fvMatrix<Type> >
 steadyStateD2dt2Scheme<Type>::fvmD2dt2
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm

src/finiteVolume/finiteVolume/d2dt2Schemes/steadyStateD2dt2Scheme/steadyStateD2dt2Scheme.H

@@ -108,19 +108,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmD2dt2
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmD2dt2
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmD2dt2
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 };
 

src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "CoEulerDdtScheme.H"
@@ -69,10 +69,10 @@ tmp<volScalarField> CoEulerDdtScheme<Type>::CorDeltaT() const
 
     forAll(owner, faceI)
     {
-        corDeltaT[owner[faceI]] = 
+        corDeltaT[owner[faceI]] =
             max(corDeltaT[owner[faceI]], cofrDeltaT[faceI]);
 
-        corDeltaT[neighbour[faceI]] = 
+        corDeltaT[neighbour[faceI]] =
             max(corDeltaT[neighbour[faceI]], cofrDeltaT[faceI]);
     }
 
@@ -127,7 +127,7 @@ tmp<surfaceScalarField> CoEulerDdtScheme<Type>::CofrDeltaT() const
         const volScalarField& rho =
             static_cast<const objectRegistry&>(mesh())
            .lookupObject<volScalarField>(rhoName_).oldTime();
- 
+
         surfaceScalarField Co
         (
             mesh().surfaceInterpolation::deltaCoeffs()
@@ -369,7 +369,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 CoEulerDdtScheme<Type>::fvmDdt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -386,7 +386,7 @@ CoEulerDdtScheme<Type>::fvmDdt
     scalarField rDeltaT = CorDeltaT()().internalField();
 
     fvm.diag() = rDeltaT*mesh().V();
-    
+
     if (mesh().moving())
     {
         fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V0();
@@ -405,7 +405,7 @@ tmp<fvMatrix<Type> >
 CoEulerDdtScheme<Type>::fvmDdt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -421,7 +421,7 @@ CoEulerDdtScheme<Type>::fvmDdt
     scalarField rDeltaT = CorDeltaT()().internalField();
 
     fvm.diag() = rDeltaT*rho.value()*mesh().V();
-    
+
     if (mesh().moving())
     {
         fvm.source() = rDeltaT
@@ -442,7 +442,7 @@ tmp<fvMatrix<Type> >
 CoEulerDdtScheme<Type>::fvmDdt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -588,7 +588,7 @@ CoEulerDdtScheme<Type>::fvcDdtPhiCorr
                 )
             );
         }
-        else if 
+        else if
         (
             U.dimensions() == dimVelocity
          && phi.dimensions() == dimDensity*dimVelocity*dimArea
@@ -617,7 +617,7 @@ CoEulerDdtScheme<Type>::fvcDdtPhiCorr
                 )
             );
         }
-        else if 
+        else if
         (
             U.dimensions() == dimDensity*dimVelocity
          && phi.dimensions() == dimDensity*dimVelocity*dimArea

src/finiteVolume/finiteVolume/ddtSchemes/CoEulerDdtScheme/CoEulerDdtScheme.H

@@ -142,19 +142,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;

src/finiteVolume/finiteVolume/ddtSchemes/CrankNicholsonDdtScheme/CrankNicholsonDdtScheme.C

@@ -625,7 +625,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 CrankNicholsonDdtScheme<Type>::fvmDdt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     DDt0Field<GeometricField<Type, fvPatchField, volMesh> >& ddt0 =
@@ -709,7 +709,7 @@ tmp<fvMatrix<Type> >
 CrankNicholsonDdtScheme<Type>::fvmDdt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     DDt0Field<GeometricField<Type, fvPatchField, volMesh> >& ddt0 =
@@ -791,7 +791,7 @@ tmp<fvMatrix<Type> >
 CrankNicholsonDdtScheme<Type>::fvmDdt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     DDt0Field<GeometricField<Type, fvPatchField, volMesh> >& ddt0 =

src/finiteVolume/finiteVolume/ddtSchemes/CrankNicholsonDdtScheme/CrankNicholsonDdtScheme.H

@@ -90,7 +90,7 @@ class CrankNicholsonDdtScheme
 
             //- Return the start-time index
             label startTimeIndex() const;
-            
+
             //- Cast to the underlying GeoField
             GeoField& operator()();
 
@@ -213,19 +213,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;

src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.C

@@ -262,7 +262,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 EulerDdtScheme<Type>::fvmDdt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -298,7 +298,7 @@ tmp<fvMatrix<Type> >
 EulerDdtScheme<Type>::fvmDdt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -335,7 +335,7 @@ tmp<fvMatrix<Type> >
 EulerDdtScheme<Type>::fvmDdt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm

src/finiteVolume/finiteVolume/ddtSchemes/EulerDdtScheme/EulerDdtScheme.H

@@ -120,19 +120,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;

src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "SLTSDdtScheme.H"
@@ -369,7 +369,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 SLTSDdtScheme<Type>::fvmDdt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -388,7 +388,7 @@ SLTSDdtScheme<Type>::fvmDdt
     Info<< "max/min rDeltaT " << max(rDeltaT) << " " << min(rDeltaT) << endl;
 
     fvm.diag() = rDeltaT*mesh().V();
-    
+
     if (mesh().moving())
     {
         fvm.source() = rDeltaT*vf.oldTime().internalField()*mesh().V0();
@@ -407,7 +407,7 @@ tmp<fvMatrix<Type> >
 SLTSDdtScheme<Type>::fvmDdt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -423,7 +423,7 @@ SLTSDdtScheme<Type>::fvmDdt
     scalarField rDeltaT = SLrDeltaT()().internalField();
 
     fvm.diag() = rDeltaT*rho.value()*mesh().V();
-    
+
     if (mesh().moving())
     {
         fvm.source() = rDeltaT
@@ -444,7 +444,7 @@ tmp<fvMatrix<Type> >
 SLTSDdtScheme<Type>::fvmDdt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -590,7 +590,7 @@ SLTSDdtScheme<Type>::fvcDdtPhiCorr
                 )
             );
         }
-        else if 
+        else if
         (
             U.dimensions() == dimVelocity
          && phi.dimensions() == dimDensity*dimVelocity*dimArea
@@ -619,7 +619,7 @@ SLTSDdtScheme<Type>::fvcDdtPhiCorr
                 )
             );
         }
-        else if 
+        else if
         (
             U.dimensions() == dimDensity*dimVelocity
          && phi.dimensions() == dimDensity*dimVelocity*dimArea

src/finiteVolume/finiteVolume/ddtSchemes/SLTSDdtScheme/SLTSDdtScheme.H

@@ -143,19 +143,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;

src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.C

@@ -361,7 +361,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 backwardDdtScheme<Type>::fvmDdt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -413,7 +413,7 @@ tmp<fvMatrix<Type> >
 backwardDdtScheme<Type>::fvmDdt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -464,7 +464,7 @@ tmp<fvMatrix<Type> >
 backwardDdtScheme<Type>::fvmDdt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm

src/finiteVolume/finiteVolume/ddtSchemes/backwardDdtScheme/backwardDdtScheme.H

@@ -131,19 +131,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;

src/finiteVolume/finiteVolume/ddtSchemes/boundedBackwardDdtScheme/boundedBackwardDdtScheme.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "boundedBackwardDdtScheme.H"
@@ -413,7 +413,7 @@ boundedBackwardDdtScheme::fvcDdt
 tmp<fvScalarMatrix>
 boundedBackwardDdtScheme::fvmDdt
 (
-    volScalarField& vf
+    const volScalarField& vf
 )
 {
     tmp<fvScalarMatrix> tfvm
@@ -484,7 +484,7 @@ tmp<fvScalarMatrix>
 boundedBackwardDdtScheme::fvmDdt
 (
     const dimensionedScalar& rho,
-    volScalarField& vf
+    const volScalarField& vf
 )
 {
     tmp<fvScalarMatrix> tfvm
@@ -554,7 +554,7 @@ tmp<fvScalarMatrix>
 boundedBackwardDdtScheme::fvmDdt
 (
     const volScalarField& rho,
-    volScalarField& vf
+    const volScalarField& vf
 )
 {
     tmp<fvScalarMatrix> tfvm

src/finiteVolume/finiteVolume/ddtSchemes/boundedBackwardDdtScheme/boundedBackwardDdtScheme.H

@@ -142,19 +142,19 @@ public:
 
         tmp<fvScalarMatrix> fvmDdt
         (
-            volScalarField&
+            const volScalarField&
         );
 
         tmp<fvScalarMatrix> fvmDdt
         (
             const dimensionedScalar&,
-            volScalarField&
+            const volScalarField&
         );
 
         tmp<fvScalarMatrix> fvmDdt
         (
             const volScalarField&,
-            volScalarField&
+            const volScalarField&
         );
 
         tmp<surfaceScalarField> fvcDdtPhiCorr

src/finiteVolume/finiteVolume/ddtSchemes/boundedBackwardDdtScheme/boundedBackwardDdtSchemes.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "boundedBackwardDdtScheme.H"

src/finiteVolume/finiteVolume/ddtSchemes/ddtScheme/ddtScheme.H

@@ -164,19 +164,19 @@ public:
 
         virtual tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
         virtual tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
         virtual tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         ) = 0;
 
 

src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "steadyStateDdtScheme.H"
@@ -162,7 +162,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 steadyStateDdtScheme<Type>::fvmDdt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -183,7 +183,7 @@ tmp<fvMatrix<Type> >
 steadyStateDdtScheme<Type>::fvmDdt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm
@@ -204,7 +204,7 @@ tmp<fvMatrix<Type> >
 steadyStateDdtScheme<Type>::fvmDdt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm

src/finiteVolume/finiteVolume/ddtSchemes/steadyStateDdtScheme/steadyStateDdtScheme.H

@@ -119,19 +119,19 @@ public:
 
         tmp<fvMatrix<Type> > fvmDdt
         (
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         tmp<fvMatrix<Type> > fvmDdt
         (
             const volScalarField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         typedef typename ddtScheme<Type>::fluxFieldType fluxFieldType;

src/finiteVolume/finiteVolume/fvc/fvcGrad.C

@@ -54,7 +54,7 @@ grad
     const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf
 )
 {
-    return fv::gaussGrad<Type>::grad(ssf);
+    return fv::gaussGrad<Type>::gradf(ssf, "grad(" + ssf.name() + ')');
 }
 
 
@@ -99,7 +99,7 @@ grad
     (
         vf.mesh(),
         vf.mesh().gradScheme(name)
-    )().grad(vf);
+    )().grad(vf, name);
 }
 
 

src/finiteVolume/finiteVolume/fvm/fvmD2dt2.C

@@ -22,9 +22,6 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Description
-    
-
 \*---------------------------------------------------------------------------*/
 
 #include "volFields.H"
@@ -48,7 +45,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 d2dt2
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fv::d2dt2Scheme<Type>::New
@@ -64,7 +61,7 @@ tmp<fvMatrix<Type> >
 d2dt2
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fv::d2dt2Scheme<Type>::New
@@ -80,7 +77,7 @@ tmp<fvMatrix<Type> >
 d2dt2
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fv::d2dt2Scheme<Type>::New

src/finiteVolume/finiteVolume/fvm/fvmD2dt2.H

@@ -54,20 +54,20 @@ namespace fvm
     tmp<fvMatrix<Type> > d2dt2
     (
         const dimensionedScalar&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > d2dt2
     (
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > d2dt2
     (
         const volScalarField&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 }
 

src/finiteVolume/finiteVolume/fvm/fvmDdt.C

@@ -45,7 +45,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 ddt
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fv::ddtScheme<Type>::New
@@ -61,7 +61,7 @@ tmp<fvMatrix<Type> >
 ddt
 (
     const oneField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return ddt(vf);
@@ -73,7 +73,7 @@ tmp<fvMatrix<Type> >
 ddt
 (
     const dimensionedScalar& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fv::ddtScheme<Type>::New
@@ -89,7 +89,7 @@ tmp<fvMatrix<Type> >
 ddt
 (
     const volScalarField& rho,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fv::ddtScheme<Type>::New

src/finiteVolume/finiteVolume/fvm/fvmDdt.H

@@ -54,28 +54,28 @@ namespace fvm
     template<class Type>
     tmp<fvMatrix<Type> > ddt
     (
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > ddt
     (
         const oneField&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > ddt
     (
         const dimensionedScalar&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > ddt
     (
         const volScalarField&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 }
 

src/finiteVolume/finiteVolume/fvm/fvmDiv.C

@@ -22,9 +22,6 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Description
-    
-
 \*---------------------------------------------------------------------------*/
 
 #include "fvmDiv.H"
@@ -49,7 +46,7 @@ tmp<fvMatrix<Type> >
 div
 (
     const surfaceScalarField& flux,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -66,7 +63,7 @@ tmp<fvMatrix<Type> >
 div
 (
     const tmp<surfaceScalarField>& tflux,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -81,7 +78,7 @@ tmp<fvMatrix<Type> >
 div
 (
     const surfaceScalarField& flux,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fvm::div(flux, vf, "div("+flux.name()+','+vf.name()+')');
@@ -92,7 +89,7 @@ tmp<fvMatrix<Type> >
 div
 (
     const tmp<surfaceScalarField>& tflux,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > Div(fvm::div(tflux(), vf));

src/finiteVolume/finiteVolume/fvm/fvmDiv.H

@@ -56,7 +56,7 @@ namespace fvm
     tmp<fvMatrix<Type> > div
     (
         const surfaceScalarField&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word& name
     );
 
@@ -64,7 +64,7 @@ namespace fvm
     tmp<fvMatrix<Type> > div
     (
         const tmp<surfaceScalarField>&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word& name
     );
 
@@ -73,14 +73,14 @@ namespace fvm
     tmp<fvMatrix<Type> > div
     (
         const surfaceScalarField&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > div
     (
         const tmp<surfaceScalarField>&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 }
 

src/finiteVolume/finiteVolume/fvm/fvmLaplacian.C

@@ -45,7 +45,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 laplacian
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -70,7 +70,7 @@ template<class Type>
 tmp<fvMatrix<Type> >
 laplacian
 (
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     surfaceScalarField Gamma
@@ -100,7 +100,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const zeroField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -116,7 +116,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const zeroField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return tmp<fvMatrix<Type> >
@@ -131,7 +131,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const oneField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -144,7 +144,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const oneField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fvm::laplacian(vf);
@@ -156,11 +156,11 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const dimensioned<GType>& gamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
-    GeometricField<GType, fvsPatchField, surfaceMesh> Gamma
+    const GeometricField<GType, fvsPatchField, surfaceMesh> Gamma
     (
         IOobject
         (
@@ -182,10 +182,10 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const dimensioned<GType>& gamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
-    GeometricField<GType, fvsPatchField, surfaceMesh> Gamma
+    const GeometricField<GType, fvsPatchField, surfaceMesh> Gamma
     (
         IOobject
         (
@@ -209,7 +209,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const GeometricField<GType, fvPatchField, volMesh>& gamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -226,7 +226,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const tmp<GeometricField<GType, fvPatchField, volMesh> >& tgamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -241,7 +241,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const GeometricField<GType, fvPatchField, volMesh>& gamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fvm::laplacian
@@ -258,7 +258,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const tmp<GeometricField<GType, fvPatchField, volMesh> >& tgamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > Laplacian(fvm::laplacian(tgamma(), vf));
@@ -274,7 +274,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const GeometricField<GType, fvsPatchField, surfaceMesh>& gamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -291,7 +291,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const tmp<GeometricField<GType, fvsPatchField, surfaceMesh> >& tgamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf,
+    const GeometricField<Type, fvPatchField, volMesh>& vf,
     const word& name
 )
 {
@@ -306,7 +306,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const GeometricField<GType, fvsPatchField, surfaceMesh>& gamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return fvm::laplacian
@@ -323,7 +323,7 @@ tmp<fvMatrix<Type> >
 laplacian
 (
     const tmp<GeometricField<GType, fvsPatchField, surfaceMesh> >& tGamma,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm(fvm::laplacian(tGamma(), vf));

src/finiteVolume/finiteVolume/fvm/fvmLaplacian.H

@@ -55,14 +55,14 @@ namespace fvm
     template<class Type>
     tmp<fvMatrix<Type> > laplacian
     (
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
     template<class Type>
     tmp<fvMatrix<Type> > laplacian
     (
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
 
@@ -70,7 +70,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const zeroField&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -78,7 +78,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const zeroField&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
 
@@ -86,7 +86,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const oneField&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -94,7 +94,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const oneField&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
 
@@ -102,7 +102,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const dimensioned<GType>&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -110,7 +110,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const dimensioned<GType>&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
 
@@ -118,7 +118,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const GeometricField<GType, fvPatchField, volMesh>&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -126,7 +126,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const GeometricField<GType, fvPatchField, volMesh>&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
 
@@ -134,7 +134,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const tmp<GeometricField<GType, fvPatchField, volMesh> >&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -142,7 +142,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const tmp<GeometricField<GType, fvPatchField, volMesh> >&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
 
@@ -150,7 +150,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const GeometricField<GType, fvsPatchField, surfaceMesh>&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -158,7 +158,7 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const tmp<GeometricField<GType, fvsPatchField, surfaceMesh> >&,
-        GeometricField<Type, fvPatchField, volMesh>&,
+        const GeometricField<Type, fvPatchField, volMesh>&,
         const word&
     );
 
@@ -166,14 +166,14 @@ namespace fvm
     tmp<fvMatrix<Type> > laplacian
     (
         const GeometricField<GType, fvsPatchField, surfaceMesh>&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 
     template<class Type, class GType>
     tmp<fvMatrix<Type> > laplacian
     (
         const tmp<GeometricField<GType, fvsPatchField, surfaceMesh> >&,
-        GeometricField<Type, fvPatchField, volMesh>&
+        const GeometricField<Type, fvPatchField, volMesh>&
     );
 }
 

src/finiteVolume/finiteVolume/fvm/fvmSup.C

@@ -35,7 +35,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Su
 (
     const DimensionedField<Type, volMesh>& su,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     const fvMesh& mesh = vf.mesh();
@@ -60,7 +60,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Su
 (
     const tmp<DimensionedField<Type, volMesh> >& tsu,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm = fvm::Su(tsu(), vf);
@@ -73,7 +73,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Su
 (
     const tmp<GeometricField<Type, fvPatchField, volMesh> >& tsu,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm = fvm::Su(tsu(), vf);
@@ -86,7 +86,7 @@ Foam::zeroField
 Foam::fvm::Su
 (
     const zeroField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return zeroField();
@@ -98,7 +98,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Sp
 (
     const DimensionedField<scalar, volMesh>& sp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     const fvMesh& mesh = vf.mesh();
@@ -123,7 +123,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Sp
 (
     const tmp<DimensionedField<scalar, volMesh> >& tsp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm = fvm::Sp(tsp(), vf);
@@ -136,7 +136,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Sp
 (
     const tmp<volScalarField>& tsp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm = fvm::Sp(tsp(), vf);
@@ -150,7 +150,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::Sp
 (
     const dimensionedScalar& sp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     const fvMesh& mesh = vf.mesh();
@@ -175,7 +175,7 @@ Foam::zeroField
 Foam::fvm::Sp
 (
     const zeroField&,
-    GeometricField<Type, fvPatchField, volMesh>&
+    const GeometricField<Type, fvPatchField, volMesh>&
 )
 {
     return zeroField();
@@ -187,7 +187,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::SuSp
 (
     const DimensionedField<scalar, volMesh>& susp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     const fvMesh& mesh = vf.mesh();
@@ -215,7 +215,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::SuSp
 (
     const tmp<DimensionedField<scalar, volMesh> >& tsusp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm = fvm::SuSp(tsusp(), vf);
@@ -228,7 +228,7 @@ Foam::tmp<Foam::fvMatrix<Type> >
 Foam::fvm::SuSp
 (
     const tmp<volScalarField>& tsusp,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     tmp<fvMatrix<Type> > tfvm = fvm::SuSp(tsusp(), vf);
@@ -241,7 +241,7 @@ Foam::zeroField
 Foam::fvm::SuSp
 (
     const zeroField&,
-    GeometricField<Type, fvPatchField, volMesh>& vf
+    const GeometricField<Type, fvPatchField, volMesh>& vf
 )
 {
     return zeroField();

src/finiteVolume/finiteVolume/fvm/fvmSup.H

@@ -56,28 +56,28 @@ namespace fvm
         tmp<fvMatrix<Type> > Su
         (
             const DimensionedField<Type, volMesh>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<fvMatrix<Type> > Su
         (
             const tmp<DimensionedField<Type, volMesh> >&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<fvMatrix<Type> > Su
         (
             const tmp<GeometricField<Type, fvPatchField, volMesh> >&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         zeroField Su
         (
             const zeroField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -87,21 +87,21 @@ namespace fvm
         tmp<fvMatrix<Type> > Sp
         (
             const DimensionedField<scalar, volMesh>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<fvMatrix<Type> > Sp
         (
             const tmp<DimensionedField<scalar, volMesh> >&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<fvMatrix<Type> > Sp
         (
             const tmp<volScalarField>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -109,7 +109,7 @@ namespace fvm
         tmp<fvMatrix<Type> > Sp
         (
             const dimensionedScalar&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -117,7 +117,7 @@ namespace fvm
         zeroField Sp
         (
             const zeroField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
 
@@ -127,28 +127,28 @@ namespace fvm
         tmp<fvMatrix<Type> > SuSp
         (
             const DimensionedField<scalar, volMesh>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<fvMatrix<Type> > SuSp
         (
             const tmp<DimensionedField<scalar, volMesh> >&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         tmp<fvMatrix<Type> > SuSp
         (
             const tmp<volScalarField>&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 
         template<class Type>
         zeroField SuSp
         (
             const zeroField&,
-            GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>&
         );
 }
 

src/finiteVolume/finiteVolume/gradSchemes/extendedLeastSquaresGrad/extendedLeastSquaresGrad.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "extendedLeastSquaresGrad.H"
@@ -35,27 +35,20 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-namespace Foam
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace fv
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 template<class Type>
-tmp
+Foam::tmp
 <
-    GeometricField
+    Foam::GeometricField
     <
-        typename outerProduct<vector, Type>::type, fvPatchField, volMesh
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+    Foam::volMesh
     >
 >
-extendedLeastSquaresGrad<Type>::grad
+Foam::fv::extendedLeastSquaresGrad<Type>::calcGrad
 (
-    const GeometricField<Type, fvPatchField, volMesh>& vsf
+    const GeometricField<Type, fvPatchField, volMesh>& vsf,
+    const word& name
 ) const
 {
     typedef typename outerProduct<vector, Type>::type GradType;
@@ -68,7 +61,7 @@ extendedLeastSquaresGrad<Type>::grad
         (
             IOobject
             (
-                "grad("+vsf.name()+')',
+                name,
                 vsf.instance(),
                 mesh,
                 IOobject::NO_READ,
@@ -120,7 +113,7 @@ extendedLeastSquaresGrad<Type>::grad
 
         if (vsf.boundaryField()[patchi].coupled())
         {
-            Field<Type> neiVsf = 
+            Field<Type> neiVsf =
                 vsf.boundaryField()[patchi].patchNeighbourField();
 
             forAll(neiVsf, patchFaceI)
@@ -162,12 +155,4 @@ extendedLeastSquaresGrad<Type>::grad
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace fv
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
 // ************************************************************************* //

src/finiteVolume/finiteVolume/gradSchemes/extendedLeastSquaresGrad/extendedLeastSquaresGrad.H

@@ -102,13 +102,16 @@ public:
 
     // Member Functions
 
-        tmp
+        //- Return the gradient of the given field to the gradScheme::grad
+        //  for optional caching
+        virtual tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
-        > grad
+        > calcGrad
         (
-            const GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>& vsf,
+            const word& name
         ) const;
 };
 

src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "fourthGrad.H"
@@ -35,27 +35,20 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-namespace Foam
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace fv
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 template<class Type>
-tmp
+Foam::tmp
 <
-    GeometricField
+    Foam::GeometricField
     <
-        typename outerProduct<vector, Type>::type, fvPatchField, volMesh
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
     >
 >
-fourthGrad<Type>::grad
+Foam::fv::fourthGrad<Type>::calcGrad
 (
-    const GeometricField<Type, fvPatchField, volMesh>& vsf
+    const GeometricField<Type, fvPatchField, volMesh>& vsf,
+    const word& name
 ) const
 {
     // The fourth-order gradient is calculated in two passes.  First,
@@ -80,7 +73,7 @@ fourthGrad<Type>::grad
         (
             IOobject
             (
-                "grad("+vsf.name()+')',
+                name,
                 vsf.instance(),
                 mesh,
                 IOobject::NO_READ,
@@ -130,7 +123,7 @@ fourthGrad<Type>::grad
             const scalarField& lambdap = lambda.boundaryField()[patchi];
 
             // Build the d-vectors
-            vectorField pd = 
+            vectorField pd =
                 mesh.Sf().boundaryField()[patchi]
                /(
                    mesh.magSf().boundaryField()[patchi]
@@ -171,12 +164,4 @@ fourthGrad<Type>::grad
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace fv
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
 // ************************************************************************* //

src/finiteVolume/finiteVolume/gradSchemes/fourthGrad/fourthGrad.H

@@ -83,13 +83,16 @@ public:
 
     // Member Functions
 
-        tmp
+        //- Return the gradient of the given field to the gradScheme::grad
+        //  for optional caching
+        virtual tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
-        > grad
+        > calcGrad
         (
-            const GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>& vsf,
+            const word& name
         ) const;
 };
 

src/finiteVolume/finiteVolume/gradSchemes/gaussGrad/gaussGrad.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "gaussGrad.H"
@@ -29,27 +29,20 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-namespace Foam
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace fv
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 template<class Type>
-tmp
+Foam::tmp
 <
-    GeometricField
+    Foam::GeometricField
     <
-        typename outerProduct<vector, Type>::type, fvPatchField, volMesh
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
     >
 >
-gaussGrad<Type>::grad
+Foam::fv::gaussGrad<Type>::gradf
 (
-    const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf
+    const GeometricField<Type, fvsPatchField, surfaceMesh>& ssf,
+    const word& name
 )
 {
     typedef typename outerProduct<vector, Type>::type GradType;
@@ -62,7 +55,7 @@ gaussGrad<Type>::grad
         (
             IOobject
             (
-                "grad("+ssf.name()+')',
+                name,
                 ssf.instance(),
                 mesh,
                 IOobject::NO_READ,
@@ -119,27 +112,29 @@ gaussGrad<Type>::grad
 
 
 template<class Type>
-tmp
+Foam::tmp
 <
-    GeometricField
+    Foam::GeometricField
     <
-        typename outerProduct<vector, Type>::type, fvPatchField, volMesh
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
     >
 >
-gaussGrad<Type>::grad
+Foam::fv::gaussGrad<Type>::calcGrad
 (
-    const GeometricField<Type, fvPatchField, volMesh>& vsf
+    const GeometricField<Type, fvPatchField, volMesh>& vsf,
+    const word& name
 ) const
 {
     typedef typename outerProduct<vector, Type>::type GradType;
 
     tmp<GeometricField<GradType, fvPatchField, volMesh> > tgGrad
     (
-        grad(tinterpScheme_().interpolate(vsf))
+        gradf(tinterpScheme_().interpolate(vsf), name)
     );
     GeometricField<GradType, fvPatchField, volMesh>& gGrad = tgGrad();
 
-    gGrad.rename("grad(" + vsf.name() + ')');
     correctBoundaryConditions(vsf, gGrad);
 
     return tgGrad;
@@ -147,7 +142,7 @@ gaussGrad<Type>::grad
 
 
 template<class Type>
-void gaussGrad<Type>::correctBoundaryConditions
+void Foam::fv::gaussGrad<Type>::correctBoundaryConditions
 (
     const GeometricField<Type, fvPatchField, volMesh>& vsf,
     GeometricField
@@ -174,12 +169,4 @@ void gaussGrad<Type>::correctBoundaryConditions
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace fv
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
 // ************************************************************************* //

src/finiteVolume/finiteVolume/gradSchemes/gaussGrad/gaussGrad.H

@@ -89,7 +89,7 @@ public:
             tinterpScheme_(new linear<Type>(mesh))
         {}
 
-        //- Construct from Istream
+        //- Construct from mesh and Istream
         gaussGrad(const fvMesh& mesh, Istream& is)
         :
             gradScheme<Type>(mesh),
@@ -116,31 +116,31 @@ public:
 
     // Member Functions
 
-        //- Return the gradient of the given field 
+        //- Return the gradient of the given field
         //  calculated using Gauss' theorem on the given surface field
         static
         tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
-        > grad
+        > gradf
         (
-            const GeometricField<Type, fvsPatchField, surfaceMesh>&
+            const GeometricField<Type, fvsPatchField, surfaceMesh>&,
+            const word& name
         );
 
-
-        //- Return the gradient of the given field calculated
-        //  using Gauss' theorem on the interpolated field
-        tmp
+        //- Return the gradient of the given field to the gradScheme::grad
+        //  for optional caching
+        virtual tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
-        > grad
+        > calcGrad
         (
-            const GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>& vsf,
+            const word& name
         ) const;
 
-
         //- Correct the boundary values of the gradient using the patchField
         // snGrad functions
         static void correctBoundaryConditions

src/finiteVolume/finiteVolume/gradSchemes/gaussGrad/gaussGrads.C

@@ -22,8 +22,6 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Description
-    
 \*---------------------------------------------------------------------------*/
 
 #include "fvMesh.H"

src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.C

@@ -22,28 +22,16 @@ License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
-Description
-    Abstract base class for finite volume calculus gradient schemes.
-
 \*---------------------------------------------------------------------------*/
 
 #include "fv.H"
-#include "HashTable.H"
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace fv
-{
+#include "objectRegistry.H"
+#include "solution.H"
 
 // * * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * //
 
 template<class Type>
-tmp<gradScheme<Type> > gradScheme<Type>::New
+Foam::tmp<Foam::fv::gradScheme<Type> > Foam::fv::gradScheme<Type>::New
 (
     const fvMesh& mesh,
     Istream& schemeData
@@ -51,7 +39,8 @@ tmp<gradScheme<Type> > gradScheme<Type>::New
 {
     if (fv::debug)
     {
-        Info<< "gradScheme<Type>::New(Istream& schemeData) : "
+        Info<< "gradScheme<Type>::New"
+               "(const fvMesh& mesh, Istream& schemeData) : "
                "constructing gradScheme<Type>"
             << endl;
     }
@@ -60,7 +49,8 @@ tmp<gradScheme<Type> > gradScheme<Type>::New
     {
         FatalIOErrorIn
         (
-            "gradScheme<Type>::New(Istream& schemeData)",
+            "gradScheme<Type>::New"
+            "(const fvMesh& mesh, Istream& schemeData)",
             schemeData
         )   << "Grad scheme not specified" << endl << endl
             << "Valid grad schemes are :" << endl
@@ -77,7 +67,8 @@ tmp<gradScheme<Type> > gradScheme<Type>::New
     {
         FatalIOErrorIn
         (
-            "gradScheme<Type>::New(Istream& schemeData)",
+            "gradScheme<Type>::New"
+            "(const fvMesh& mesh, Istream& schemeData)",
             schemeData
         )   << "unknown grad scheme " << schemeName << endl << endl
             << "Valid grad schemes are :" << endl
@@ -92,16 +83,151 @@ tmp<gradScheme<Type> > gradScheme<Type>::New
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 template<class Type>
-gradScheme<Type>::~gradScheme()
+Foam::fv::gradScheme<Type>::~gradScheme()
 {}
 
-
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace fv
+namespace Foam
+{
+    template<class Type>
+    inline void cachePrintMessage
+    (
+        const char* message,
+        const word& name,
+        const GeometricField<Type, fvPatchField, volMesh>& vf
+    )
+    {
+        if (solution::debug)
+        {
+            Info<< "Cache: " << message << token::SPACE << name
+                << ", " << vf.name() << " event No. " << vf.eventNo()
+                << endl;
+        }
+    }
+}
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+template<class Type>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
+    >
+>
+Foam::fv::gradScheme<Type>::grad
+(
+    const GeometricField<Type, fvPatchField, volMesh>& vsf,
+    const word& name
+) const
+{
+    typedef typename outerProduct<vector, Type>::type GradType;
+    typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
+
+    if (!this->mesh().changing() && this->mesh().cache(name))
+    {
+        if (!mesh().objectRegistry::foundObject<GradFieldType>(name))
+        {
+            cachePrintMessage("Calculating and caching", name, vsf);
+            tmp<GradFieldType> tgGrad = calcGrad(vsf, name);
+            regIOobject::store(tgGrad.ptr());
+        }
+
+        cachePrintMessage("Retreiving", name, vsf);
+        GradFieldType& gGrad = const_cast<GradFieldType&>
+        (
+            mesh().objectRegistry::lookupObject<GradFieldType>(name)
+        );
+
+        if (gGrad.upToDate(vsf))
+        {
+            return gGrad;
+        }
+        else
+        {
+            cachePrintMessage("Deleting", name, vsf);
+            gGrad.release();
+            delete &gGrad;
+
+            cachePrintMessage("Recalculating", name, vsf);
+            tmp<GradFieldType> tgGrad = calcGrad(vsf, name);
+
+            cachePrintMessage("Storing", name, vsf);
+            regIOobject::store(tgGrad.ptr());
+            GradFieldType& gGrad = const_cast<GradFieldType&>
+            (
+                mesh().objectRegistry::lookupObject<GradFieldType>(name)
+            );
+
+            return gGrad;
+        }
+    }
+    else
+    {
+        if (mesh().objectRegistry::foundObject<GradFieldType>(name))
+        {
+            GradFieldType& gGrad = const_cast<GradFieldType&>
+            (
+                mesh().objectRegistry::lookupObject<GradFieldType>(name)
+            );
+
+            if (gGrad.ownedByRegistry())
+            {
+                cachePrintMessage("Deleting", name, vsf);
+                gGrad.release();
+                delete &gGrad;
+            }
+        }
+
+        cachePrintMessage("Calculating", name, vsf);
+        return calcGrad(vsf, name);
+    }
+}
+
+
+template<class Type>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
+    >
+>
+Foam::fv::gradScheme<Type>::grad
+(
+    const GeometricField<Type, fvPatchField, volMesh>& vsf
+) const
+{
+    return grad(vsf, "grad(" + vsf.name() + ')');
+}
+
+
+template<class Type>
+Foam::tmp
+<
+    Foam::GeometricField
+    <
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
+    >
+>
+Foam::fv::gradScheme<Type>::grad
+(
+    const tmp<GeometricField<Type, fvPatchField, volMesh> >& tvsf
+) const
+{
+    typedef typename outerProduct<vector, Type>::type GradType;
+    typedef GeometricField<GradType, fvPatchField, volMesh> GradFieldType;
+
+    tmp<GradFieldType> tgrad = grad(tvsf());
+    tvsf.clear();
+    return tgrad;
+}
 
-} // End namespace Foam
 
 // ************************************************************************* //

src/finiteVolume/finiteVolume/gradSchemes/gradScheme/gradScheme.H

@@ -114,9 +114,8 @@ public:
         );
 
 
-    // Destructor
-
-        virtual ~gradScheme();
+    //- Destructor
+    virtual ~gradScheme();
 
 
     // Member Functions
@@ -127,15 +126,54 @@ public:
             return mesh_;
         }
 
-        //- Calculate and return the grad of the given field
+        //- Calculate and return the grad of the given field.
+        //  Used by grad either to recalculate the cached gradient when it is
+        //  out of date with respect to the field or when it is not cached.
         virtual tmp
+        <
+            GeometricField
+            <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
+        > calcGrad
+        (
+            const GeometricField<Type, fvPatchField, volMesh>&,
+            const word& name
+        ) const = 0;
+
+        //- Calculate and return the grad of the given field
+        //  which may have been cached
+        tmp
+        <
+            GeometricField
+            <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
+        > grad
+        (
+            const GeometricField<Type, fvPatchField, volMesh>&,
+            const word& name
+        ) const;
+
+        //- Calculate and return the grad of the given field
+        //  with the default name
+        //  which may have been cached
+        tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
         > grad
         (
             const GeometricField<Type, fvPatchField, volMesh>&
-        ) const = 0;
+        ) const;
+
+        //- Calculate and return the grad of the given field
+        //  with the default name
+        //  which may have been cached
+        tmp
+        <
+            GeometricField
+            <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
+        > grad
+        (
+            const tmp<GeometricField<Type, fvPatchField, volMesh> >&
+        ) const;
 };
 
 

src/finiteVolume/finiteVolume/gradSchemes/leastSquaresGrad/leastSquaresGrad.C

@@ -21,7 +21,7 @@ License
     You should have received a copy of the GNU General Public License
     along with OpenFOAM; if not, write to the Free Software Foundation,
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-    
+
 \*---------------------------------------------------------------------------*/
 
 #include "leastSquaresGrad.H"
@@ -35,27 +35,20 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-namespace Foam
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace fv
-{
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 template<class Type>
-tmp
+Foam::tmp
 <
-    GeometricField
+    Foam::GeometricField
     <
-        typename outerProduct<vector, Type>::type, fvPatchField, volMesh
+        typename Foam::outerProduct<Foam::vector, Type>::type,
+        Foam::fvPatchField,
+        Foam::volMesh
     >
 >
-leastSquaresGrad<Type>::grad
+Foam::fv::leastSquaresGrad<Type>::calcGrad
 (
-    const GeometricField<Type, fvPatchField, volMesh>& vsf
+    const GeometricField<Type, fvPatchField, volMesh>& vsf,
+    const word& name
 ) const
 {
     typedef typename outerProduct<vector, Type>::type GradType;
@@ -68,7 +61,7 @@ leastSquaresGrad<Type>::grad
         (
             IOobject
             (
-                "grad("+vsf.name()+')',
+                name,
                 vsf.instance(),
                 mesh,
                 IOobject::NO_READ,
@@ -116,7 +109,7 @@ leastSquaresGrad<Type>::grad
 
         if (vsf.boundaryField()[patchi].coupled())
         {
-            Field<Type> neiVsf = 
+            Field<Type> neiVsf =
                 vsf.boundaryField()[patchi].patchNeighbourField();
 
             forAll(neiVsf, patchFaceI)
@@ -147,12 +140,4 @@ leastSquaresGrad<Type>::grad
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace fv
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
 // ************************************************************************* //

src/finiteVolume/finiteVolume/gradSchemes/leastSquaresGrad/leastSquaresGrad.H

@@ -89,13 +89,16 @@ public:
 
     // Member Functions
 
-        tmp
+        //- Return the gradient of the given field to the gradScheme::grad
+        //  for optional caching
+        virtual tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
-        > grad
+        > calcGrad
         (
-            const GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>& vsf,
+            const word& name
         ) const;
 };
 

src/finiteVolume/finiteVolume/gradSchemes/limitedGradSchemes/cellLimitedGrad/cellLimitedGrad.H

@@ -117,18 +117,82 @@ public:
             const Type& extrapolate
         );
 
-
-        tmp
+        //- Return the gradient of the given field to the gradScheme::grad
+        //  for optional caching
+        virtual tmp
         <
             GeometricField
             <typename outerProduct<vector, Type>::type, fvPatchField, volMesh>
-        > grad
+        > calcGrad
         (
-            const GeometricField<Type, fvPatchField, volMesh>&
+            const GeometricField<Type, fvPatchField, volMesh>& vsf,
+            const word& name
         ) const;
 };
 
 
+// * * * * * * * * * * * * Inline Member Function  * * * * * * * * * * * * * //
+
+template<>
+inline void cellLimitedGrad<scalar>::limitFace
+(
+    scalar& limiter,
+    const scalar& maxDelta,
+    const scalar& minDelta,
+    const scalar& extrapolate
+)
+{
+    if (extrapolate > maxDelta + VSMALL)
+    {
+        limiter = min(limiter, maxDelta/extrapolate);
+    }
+    else if (extrapolate < minDelta - VSMALL)
+    {
+        limiter = min(limiter, minDelta/extrapolate);
+    }
+}
+
+
+template<class Type>
+inline void cellLimitedGrad<Type>::limitFace
+(
+    Type& limiter,
+    const Type& maxDelta,
+    const Type& minDelta,
+    const Type& extrapolate
+)
+{
+    for(direction cmpt=0; cmpt<Type::nComponents; cmpt++)
+    {
+        cellLimitedGrad<scalar>::limitFace
+        (
+            limiter.component(cmpt),
+            maxDelta.component(cmpt),
+            minDelta.component(cmpt),
+            extrapolate.component(cmpt)
+        );
+    }
+}
+
+
+// * * * * * * * * Template Member Function Specialisations  * * * * * * * * //
+
+template<>
+tmp<volVectorField> cellLimitedGrad<scalar>::calcGrad
+(
+    const volScalarField& vsf,
+    const word& name
+) const;
+
+
+template<>
+tmp<volTensorField> cellLimitedGrad<vector>::calcGrad
+(
+    const volVectorField& vsf,
+    const word& name
+) const;
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace fv

src/finiteVolume/finiteVolume/gradSchemes/limitedGradSchemes/cellLimitedGrad/cellLimitedGrads.C

@@ -36,70 +36,25 @@ License
 
 namespace Foam
 {
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 namespace fv
 {
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-makeFvGradScheme(cellLimitedGrad)
+    makeFvGradScheme(cellLimitedGrad)
+}
+}
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 template<>
-inline void cellLimitedGrad<scalar>::limitFace
+Foam::tmp<Foam::volVectorField>
+Foam::fv::cellLimitedGrad<Foam::scalar>::calcGrad
 (
-    scalar& limiter,
-    const scalar& maxDelta,
-    const scalar& minDelta,
-    const scalar& extrapolate
-)
-{
-    if (extrapolate > maxDelta + VSMALL)
-    {
-        limiter = min(limiter, maxDelta/extrapolate);
-    }
-    else if (extrapolate < minDelta - VSMALL)
-    {
-        limiter = min(limiter, minDelta/extrapolate);
-    }
-}
-
-template<class Type>
-inline void cellLimitedGrad<Type>::limitFace
-(
-    Type& limiter,
-    const Type& maxDelta,
-    const Type& minDelta,
-    const Type& extrapolate
-)
-{
-    for(direction cmpt=0; cmpt<Type::nComponents; cmpt++)
-    {
-        cellLimitedGrad<scalar>::limitFace
-        (
-            limiter.component(cmpt),
-            maxDelta.component(cmpt),
-            minDelta.component(cmpt),
-            extrapolate.component(cmpt)
-        );
-    }
-}
-
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-template<>
-tmp<volVectorField> cellLimitedGrad<scalar>::grad
-(
-    const volScalarField& vsf
+    const volScalarField& vsf,
+    const word& name
 ) const
 {
     const fvMesh& mesh = vsf.mesh();
 
-    tmp<volVectorField> tGrad = basicGradScheme_().grad(vsf);
+    tmp<volVectorField> tGrad = basicGradScheme_().calcGrad(vsf, name);
 
     if (k_ < SMALL)
     {
@@ -244,14 +199,16 @@ tmp<volVectorField> cellLimitedGrad<scalar>::grad
 
 
 template<>
-tmp<volTensorField> cellLimitedGrad<vector>::grad
+Foam::tmp<Foam::volTensorField>
+Foam::fv::cellLimitedGrad<Foam::vector>::calcGrad
 (
-    const volVectorField& vsf
+    const volVectorField& vsf,
+    const word& name
 ) const
 {
     const fvMesh& mesh = vsf.mesh();
 
-    tmp<volTensorField> tGrad = basicGradScheme_().grad(vsf);
+    tmp<volTensorField> tGrad = basicGradScheme_().calcGrad(vsf, name);
 
     if (k_ < SMALL)
     {
@@ -402,12 +359,4 @@ tmp<volTensorField> cellLimitedGrad<vector>::grad
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace fv
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
 // ************************************************************************* //