Merge branch 'feature/cvMesh' of ssh://opencfd-laurence.no-ip.org:8004/home/dm4/OpenFOAM/OpenFOAM-dev into feature/cvMesh

applications/solvers/heatTransfer/buoyantBoussinesqPimpleFoam/TEqn.H

@@ -1,14 +1,14 @@
 {
-    kappat = turbulence->nut()/Prt;
+    alphat = turbulence->nut()/Prt;
-    kappat.correctBoundaryConditions();
+    alphat.correctBoundaryConditions();
 
-    volScalarField kappaEff("kappaEff", turbulence->nu()/Pr + kappat);
+    volScalarField alphaEff("alphaEff", turbulence->nu()/Pr + alphat);
 
     fvScalarMatrix TEqn
     (
         fvm::ddt(T)
       + fvm::div(phi, T)
-      - fvm::laplacian(kappaEff, T)
+      - fvm::laplacian(alphaEff, T)
      ==
         radiation->ST(rhoCpRef, T)
     );

applications/solvers/heatTransfer/buoyantBoussinesqPimpleFoam/createFields.H

@@ -65,12 +65,12 @@
     );
 
     // kinematic turbulent thermal thermal conductivity m2/s
-    Info<< "Reading field kappat\n" << endl;
+    Info<< "Reading field alphat\n" << endl;
-    volScalarField kappat
+    volScalarField alphat
     (
         IOobject
         (
-            "kappat",
+            "alphat",
             runTime.timeName(),
             mesh,
             IOobject::MUST_READ,
@@ -116,4 +116,3 @@
             pRefValue - getRefCellValue(p, pRefCell)
         );
     }
-

applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/TEqn.H

@@ -1,13 +1,13 @@
 {
-    kappat = turbulence->nut()/Prt;
+    alphat = turbulence->nut()/Prt;
-    kappat.correctBoundaryConditions();
+    alphat.correctBoundaryConditions();
 
-    volScalarField kappaEff("kappaEff", turbulence->nu()/Pr + kappat);
+    volScalarField alphaEff("alphaEff", turbulence->nu()/Pr + alphat);
 
     fvScalarMatrix TEqn
     (
         fvm::div(phi, T)
-      - fvm::laplacian(kappaEff, T)
+      - fvm::laplacian(alphaEff, T)
     );
 
     TEqn.relax();

applications/solvers/heatTransfer/buoyantBoussinesqSimpleFoam/createFields.H

@@ -65,12 +65,12 @@
     );
 
     // kinematic turbulent thermal thermal conductivity m2/s
-    Info<< "Reading field kappat\n" << endl;
+    Info<< "Reading field alphat\n" << endl;
-    volScalarField kappat
+    volScalarField alphat
     (
         IOobject
         (
-            "kappat",
+            "alphat",
             runTime.timeName(),
             mesh,
             IOobject::MUST_READ,

applications/solvers/heatTransfer/chtMultiRegionFoam/solid/setRegionSolidFields.H

@@ -12,7 +12,7 @@
 
     if (!thermo.isotropic())
     {
-        tmp<volVectorField> tkappaCp = thermo.Kappa()/cp;
+        tmp<volVectorField> tkappaByCp = thermo.Kappa()/cp;
 
         const coordinateSystem& coodSys = coordinates[i];
         tAnialpha =
@@ -29,13 +29,13 @@
                         IOobject::NO_WRITE
                     ),
                     mesh,
-                    tkappaCp().dimensions(),
+                    tkappaByCp().dimensions(),
                     zeroGradientFvPatchVectorField::typeName
                 )
             );
 
         volSymmTensorField& Anialpha = tAnialpha();
-        Anialpha.internalField() = coodSys.R().transformVector(tkappaCp());
+        Anialpha.internalField() = coodSys.R().transformVector(tkappaByCp());
         Anialpha.correctBoundaryConditions();
     }
 

applications/solvers/incompressible/pimpleFoam/pimpleDyMFoam/UEqn.H

@@ -17,7 +17,7 @@ rAU = 1.0/UEqn().A();
 
 if (pimple.momentumPredictor())
 {
-    solve(UEqn() == -fvc::grad(p) + fvOptions(U));
+    solve(UEqn() == -fvc::grad(p));
 
     fvOptions.correct(U);
 }

applications/utilities/mesh/generation/cvMesh/conformalVoronoiMesh/cellSizeControlSurfaces/cellSizeFunction/linearDistance/linearDistance.C

@@ -129,7 +129,7 @@ bool linearDistance::sizeLocations
         shapePts.resize(1);
         shapeSizes.resize(1);
 
-        shapePts[0] = pt - n*distance_;
+        shapePts[0] = pt + n*distance_;
         shapeSizes[0] = distanceCellSize_;
     }
 

applications/utilities/mesh/generation/snappyHexMesh/snappyHexMesh.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -47,6 +47,8 @@ Description
 #include "snapParameters.H"
 #include "layerParameters.H"
 #include "vtkSetWriter.H"
+#include "faceSet.H"
+#include "motionSmoother.H"
 
 using namespace Foam;
 
@@ -675,6 +677,31 @@ int main(int argc, char *argv[])
     }
 
 
+    {
+        // Check final mesh
+        Info<< "Checking final mesh ..." << endl;
+        faceSet wrongFaces(mesh, "wrongFaces", mesh.nFaces()/100);
+        motionSmoother::checkMesh(false, mesh, motionDict, wrongFaces);
+        const label nErrors = returnReduce
+        (
+            wrongFaces.size(),
+            sumOp<label>()
+        );
+
+        if (nErrors > 0)
+        {
+            Info<< "Finished meshing with " << nErrors << " illegal faces"
+                << " (concave, zero area or negative cell pyramid volume)"
+                << endl;
+            wrongFaces.write();
+        }
+        else
+        {
+            Info<< "Finished meshing without any errors" << endl;
+        }
+    }
+
+
     Info<< "Finished meshing in = "
         << runTime.elapsedCpuTime() << " s." << endl;
 

applications/utilities/mesh/manipulation/checkMesh/checkGeometry.C

@@ -834,5 +834,39 @@ Foam::label Foam::checkGeometry(const polyMesh& mesh, const bool allGeometry)
         }
     }
 
+    if (allGeometry)
+    {
+        faceSet faces(mesh, "lowWeightFaces", mesh.nFaces()/100);
+        if (mesh.checkFaceWeight(true, 0.05, &faces))
+        {
+            noFailedChecks++;
+
+            label nFaces = returnReduce(faces.size(), sumOp<label>());
+
+            Info<< "  <<Writing " << nFaces
+                << " faces with low interpolation weights to set "
+                << faces.name() << endl;
+            faces.instance() = mesh.pointsInstance();
+            faces.write();
+        }
+    }
+
+    if (allGeometry)
+    {
+        faceSet faces(mesh, "lowVolRatioFaces", mesh.nFaces()/100);
+        if (mesh.checkVolRatio(true, 0.05, &faces))
+        {
+            noFailedChecks++;
+
+            label nFaces = returnReduce(faces.size(), sumOp<label>());
+
+            Info<< "  <<Writing " << nFaces
+                << " faces with low volume ratio cells to set "
+                << faces.name() << endl;
+            faces.instance() = mesh.pointsInstance();
+            faces.write();
+        }
+    }
+
     return noFailedChecks;
 }

applications/utilities/mesh/manipulation/createBaffles/createBafflesDict

@@ -138,6 +138,8 @@ baffles
         //- Select faces and orientation through a searchableSurface
         type        searchableSurface;
         surface     searchablePlate;
+        //name sphere.stl;  // name if surface=triSurfaceMesh
+
         origin      (0.099 -0.006 0.004);
         span        (0 0.012 0.012);
 

applications/utilities/mesh/manipulation/createBaffles/faceSelection/searchableSurfaceSelection.C

@@ -28,6 +28,7 @@ License
 #include "syncTools.H"
 #include "searchableSurface.H"
 #include "fvMesh.H"
+#include "Time.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -61,7 +62,15 @@ Foam::faceSelections::searchableSurfaceSelection::searchableSurfaceSelection
         searchableSurface::New
         (
             word(dict.lookup("surface")),
+            IOobject
+            (
+                dict.lookupOrDefault("name", mesh.objectRegistry::db().name()),
+                mesh.time().constant(),
+                "triSurface",
                 mesh.objectRegistry::db(),
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE
+            ),
             dict
         )
     )

applications/utilities/mesh/manipulation/topoSet/topoSetDict

@@ -378,6 +378,7 @@ FoamFile
 //            surface searchableSphere;
 //            centre  (0.05 0.05 0.005);
 //            radius  0.025;
+//            //name  sphere.stl; // Optional name if surface triSurfaceMesh
 //        }
 //    }
 //

applications/utilities/parallelProcessing/decomposePar/decomposePar.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -539,7 +539,7 @@ int main(int argc, char *argv[])
                     cloud::prefix/cloudDirs[i]
                 );
 
-                IOobject* positionsPtr = sprayObjs.lookup("positions");
+                IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
 
                 if (positionsPtr)
                 {

applications/utilities/parallelProcessing/reconstructPar/reconstructPar.C

@@ -529,10 +529,8 @@ int main(int argc, char *argv[])
                                 cloud::prefix/cloudDirs[i]
                             );
 
-                            IOobject* positionsPtr = sprayObjs.lookup
+                            IOobject* positionsPtr =
-                            (
+                                sprayObjs.lookup(word("positions"));
-                                "positions"
-                            );
 
                             if (positionsPtr)
                             {

applications/utilities/postProcessing/dataConversion/foamToEnsight/ensightMesh.C

@@ -1026,7 +1026,7 @@ void Foam::ensightMesh::write
 
     if (timeIndex == 0)
     {
-        timeFile += "000.";
+        timeFile += "0000.";
     }
     else if (meshMoving)
     {

applications/utilities/postProcessing/dataConversion/foamToEnsight/foamToEnsight.C

@@ -275,7 +275,7 @@ int main(int argc, char *argv[])
                 cloud::prefix/cloudDirs[cloudI]
             );
 
-            IOobject* positionsPtr = cloudObjs.lookup("positions");
+            IOobject* positionsPtr = cloudObjs.lookup(word("positions"));
 
             if (positionsPtr)
             {

applications/utilities/postProcessing/dataConversion/foamToEnsightParts/foamToEnsightParts.C

@@ -380,7 +380,8 @@ int main(int argc, char *argv[])
             );
 
             // check that the positions field is present for this time
-            if (cloudObjs.lookup("positions"))
+            IOobject* positionPtr = cloudObjs.lookup(word("positions"));
+            if (positionPtr != NULL)
             {
                 ensightParticlePositions
                 (

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

@@ -448,7 +448,7 @@ int main(int argc, char *argv[])
                 cloud::prefix/cloudDirs[i]
             );
 
-            IOobject* positionsPtr = sprayObjs.lookup("positions");
+            IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
 
             if (positionsPtr)
             {
@@ -1094,7 +1094,7 @@ int main(int argc, char *argv[])
                 cloud::prefix/cloudName
             );
 
-            IOobject* positionsPtr = sprayObjs.lookup("positions");
+            IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
 
             if (positionsPtr)
             {

applications/utilities/postProcessing/graphics/PV398Readers/Allwmake

@@ -2,6 +2,8 @@
 cd ${0%/*} || exit 1    # run from this directory
 #set -x
 
+if [ "$ParaView_VERSION" == "3.98.1" ]
+    then
     if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
     then
         [ -n "$PV_PLUGIN_PATH" ] || {
@@ -18,5 +20,8 @@ then
     else
         echo "ERROR: ParaView not found in $ParaView_DIR"
     fi
+else
+    echo "WARN: PV398 readers not building: ParaView_VERSION=$ParaView_VERSION"
+fi
 
 # ----------------------------------------------------------------- end-of-file

applications/utilities/postProcessing/graphics/PV398Readers/PV398FoamReader/vtkPV398Foam/vtkPV398FoamMeshLagrangian.C

@@ -63,7 +63,7 @@ vtkPolyData* Foam::vtkPV398Foam::lagrangianVTKMesh
         cloud::prefix/cloudName
     );
 
-    IOobject* positionsPtr = sprayObjs.lookup("positions");
+    IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
     if (positionsPtr)
     {
         Cloud<passiveParticle> parcels(mesh, cloudName, false);

applications/utilities/postProcessing/graphics/PV3Readers/Allwmake

@@ -2,6 +2,8 @@
 cd ${0%/*} || exit 1    # run from this directory
 #set -x
 
+if [ "$ParaView_VERSION" != "3.98.1" ]
+    then
     if [ -d "$ParaView_DIR" -a -r "$ParaView_DIR" ]
     then
         [ -n "$PV_PLUGIN_PATH" ] || {
@@ -18,5 +20,8 @@ then
     else
         echo "ERROR: ParaView not found in $ParaView_DIR"
     fi
+else
+    echo "WARN: PV3 readers not building: ParaView_VERSION=$ParaView_VERSION"
+fi
 
 # ----------------------------------------------------------------- end-of-file

applications/utilities/postProcessing/graphics/PV3Readers/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMesh.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -179,7 +179,10 @@ void Foam::vtkPV3Foam::convertMeshPatches
 
         const word patchName = getPartName(partId);
 
-        labelHashSet patchIds(patches.patchSet(List<wordRe>(1, patchName)));
+        labelHashSet patchIds
+        (
+            patches.patchSet(List<wordRe>(1, wordRe(patchName)))
+        );
 
         if (debug)
         {

applications/utilities/postProcessing/graphics/PV3Readers/PV3FoamReader/vtkPV3Foam/vtkPV3FoamMeshLagrangian.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -63,7 +63,7 @@ vtkPolyData* Foam::vtkPV3Foam::lagrangianVTKMesh
         cloud::prefix/cloudName
     );
 
-    IOobject* positionsPtr = sprayObjs.lookup("positions");
+    IOobject* positionsPtr = sprayObjs.lookup(word("positions"));
     if (positionsPtr)
     {
         Cloud<passiveParticle> parcels(mesh, cloudName, false);

applications/utilities/postProcessing/wall/wallGradU/wallGradU.C

@@ -30,6 +30,7 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "fvCFD.H"
+#include "wallFvPatch.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -83,11 +84,18 @@ int main(int argc, char *argv[])
                 )
             );
 
+            const fvPatchList& patches = mesh.boundary();
+
             forAll(wallGradU.boundaryField(), patchi)
+            {
+                const fvPatch& currPatch = patches[patchi];
+
+                if (isA<wallFvPatch>(currPatch))
                 {
                     wallGradU.boundaryField()[patchi] =
                         -U.boundaryField()[patchi].snGrad();
                 }
+            }
 
             wallGradU.write();
         }

applications/utilities/preProcessing/mapFields/mapLagrangian.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -120,7 +120,7 @@ void mapLagrangian(const meshToMesh& meshToMeshInterp)
             cloud::prefix/cloudDirs[cloudI]
         );
 
-        IOobject* positionsPtr = objects.lookup("positions");
+        IOobject* positionsPtr = objects.lookup(word("positions"));
 
         if (positionsPtr)
         {

src/OpenFOAM/db/IOobjectList/IOobjectList.C

@@ -148,7 +148,7 @@ Foam::IOobject* Foam::IOobjectList::lookup(const word& name) const
 }
 
 
-Foam::IOobjectList Foam::IOobjectList::lookupRe(const wordRe& name) const
+Foam::IOobjectList Foam::IOobjectList::lookup(const wordRe& name) const
 {
     IOobjectList objectsOfName(size());
 
@@ -169,6 +169,29 @@ Foam::IOobjectList Foam::IOobjectList::lookupRe(const wordRe& name) const
 }
 
 
+Foam::IOobjectList Foam::IOobjectList::lookup(const wordReList& patterns) const
+{
+    wordReListMatcher names(patterns);
+
+    IOobjectList objectsOfName(size());
+
+    forAllConstIter(HashPtrTable<IOobject>, *this, iter)
+    {
+        if (names.match(iter()->name()))
+        {
+            if (IOobject::debug)
+            {
+                Info<< "IOobjectList::lookupRe : found " << iter.key() << endl;
+            }
+
+            objectsOfName.insert(iter.key(), new IOobject(*iter()));
+        }
+    }
+
+    return objectsOfName;
+}
+
+
 Foam::IOobjectList Foam::IOobjectList::lookupClass(const word& ClassName) const
 {
     IOobjectList objectsOfClass(size());

src/OpenFOAM/db/IOobjectList/IOobjectList.H

@@ -37,6 +37,7 @@ SourceFiles
 
 #include "HashPtrTable.H"
 #include "IOobject.H"
+#include "wordReList.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -92,7 +93,10 @@ public:
         IOobject* lookup(const word& name) const;
 
         //- Return the list for all IOobects whose name matches name
-        IOobjectList lookupRe(const wordRe& name) const;
+        IOobjectList lookup(const wordRe& name) const;
+
+        //- Return the list for all IOobects whose name matches name
+        IOobjectList lookup(const wordReList& patterns) const;
 
         //- Return the list for all IOobjects of a given class
         IOobjectList lookupClass(const word& className) const;

src/OpenFOAM/db/objectRegistry/objectRegistry.H

@@ -37,6 +37,7 @@ SourceFiles
 
 #include "HashTable.H"
 #include "regIOobject.H"
+#include "wordReList.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -147,6 +148,14 @@ public:
             template<class Type>
             wordList names() const;
 
+            //- Return the list of objects whose name matches the input regExp
+            template<class Type>
+            wordList names(const wordRe& name) const;
+
+            //- Return the list of objects whose name matches the input regExp
+            template<class Type>
+            wordList names(const wordReList& name) const;
+
             //- Lookup and return a const sub-objectRegistry. Optionally create
             //  it if it does not exist.
             const objectRegistry& subRegistry
@@ -163,10 +172,6 @@ public:
             template<class Type>
             HashTable<Type*> lookupClass(const bool strict = false);
 
-            //- Return the list of objects whose name matches the input regExp
-            template<class Type>
-            wordList foundObjectRe(const wordRe& name) const;
-
             //- Is the named Type found?
             template<class Type>
             bool foundObject(const word& name) const;

src/OpenFOAM/db/objectRegistry/objectRegistryTemplates.C

@@ -24,7 +24,7 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "objectRegistry.H"
-
+#include "stringListOps.H"
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
@@ -48,6 +48,40 @@ Foam::wordList Foam::objectRegistry::names() const
 }
 
 
+template<class Type>
+Foam::wordList Foam::objectRegistry::names(const wordRe& name) const
+{
+    wordList objectNames(size());
+
+    label count = 0;
+    forAllConstIter(HashTable<regIOobject*>, *this, iter)
+    {
+        if (isA<Type>(*iter()))
+        {
+            const word& objectName = iter()->name();
+
+            if (name.match(objectName))
+            {
+                objectNames[count++] = objectName;
+            }
+        }
+    }
+
+    objectNames.setSize(count);
+
+    return objectNames;
+}
+
+
+template<class Type>
+Foam::wordList Foam::objectRegistry::names(const wordReList& patterns) const
+{
+    wordList names(this->names<Type>());
+
+    return wordList(names, findStrings(patterns, names));
+}
+
+
 template<class Type>
 Foam::HashTable<const Type*> Foam::objectRegistry::lookupClass
 (
@@ -127,31 +161,6 @@ bool Foam::objectRegistry::foundObject(const word& name) const
 }
 
 
-template<class Type>
-Foam::wordList Foam::objectRegistry::foundObjectRe(const wordRe& name) const
-{
-    wordList objectNames(size());
-
-    label count = 0;
-    forAllConstIter(HashTable<regIOobject*>, *this, iter)
-    {
-        if (isA<Type>(*iter()))
-        {
-            const word& objectName = iter()->name();
-
-            if (name.match(objectName))
-            {
-                objectNames[count++] = objectName;
-            }
-        }
-    }
-
-    objectNames.setSize(count);
-
-    return objectNames;
-}
-
-
 template<class Type>
 const Type& Foam::objectRegistry::lookupObject(const word& name) const
 {

src/OpenFOAM/meshes/polyMesh/polyBoundaryMesh/polyBoundaryMesh.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -737,7 +737,7 @@ Foam::labelHashSet Foam::polyBoundaryMesh::patchSet
 
     forAll(patchNames, i)
     {
-         const word& patchName = patchNames[i];
+         const wordRe& patchName = patchNames[i];
 
         // Treat the given patch names as wild-cards and search the set
         // of all patch names for matches

src/OpenFOAM/meshes/polyMesh/polyMesh.H

@@ -265,6 +265,24 @@ private:
                 const Vector<label>& meshD
             ) const;
 
+            bool checkFaceWeight
+            (
+                const vectorField& fCtrs,
+                const vectorField& fAreas,
+                const vectorField& cellCtrs,
+                const bool report,
+                const scalar minWeight,
+                labelHashSet* setPtr
+            ) const;
+
+            bool checkVolRatio
+            (
+                const scalarField& cellVols,
+                const bool report,
+                const scalar minRatio,
+                labelHashSet* setPtr
+            ) const;
+
 public:
 
     // Public typedefs
@@ -612,6 +630,22 @@ public:
                 const bool detailedReport = false
             ) const;
 
+            //- Check for face weights
+            virtual bool checkFaceWeight
+            (
+                const bool report,
+                const scalar minWeight = 0.05,
+                labelHashSet* setPtr = NULL
+            ) const;
+
+            //- Check for neighbouring cell volumes
+            virtual bool checkVolRatio
+            (
+                const bool report,
+                const scalar minRatio = 0.01,
+                labelHashSet* setPtr = NULL
+            ) const;
+
 
         // Helper functions
 

src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshCheck.C

@@ -500,6 +500,190 @@ bool Foam::polyMesh::checkCellDeterminant
 }
 
 
+bool Foam::polyMesh::checkFaceWeight
+(
+    const vectorField& fCtrs,
+    const vectorField& fAreas,
+    const vectorField& cellCtrs,
+    const bool report,
+    const scalar minWeight,
+    labelHashSet* setPtr
+) const
+{
+    if (debug)
+    {
+        Info<< "bool polyMesh::checkFaceWeight(const bool"
+            << ", labelHashSet*) const: "
+            << "checking for low face interpolation weights" << endl;
+    }
+
+    tmp<scalarField> tfaceWght = polyMeshTools::faceWeights
+    (
+        *this,
+        fCtrs,
+        fAreas,
+        cellCtrs
+    );
+    scalarField& faceWght = tfaceWght();
+
+
+    label nErrorFaces = 0;
+    scalar minDet = GREAT;
+    scalar sumDet = 0.0;
+    label nSummed = 0;
+
+    // Statistics only for internal and masters of coupled faces
+    PackedBoolList isMasterFace(syncTools::getInternalOrMasterFaces(*this));
+
+    forAll (faceWght, faceI)
+    {
+        if (faceWght[faceI] < minWeight)
+        {
+            // Note: insert both sides of coupled faces
+            if (setPtr)
+            {
+                setPtr->insert(faceI);
+            }
+
+            nErrorFaces++;
+        }
+
+        // Note: statistics only on master of coupled faces
+        if (isMasterFace[faceI])
+        {
+            minDet = min(minDet, faceWght[faceI]);
+            sumDet += faceWght[faceI];
+            nSummed++;
+        }
+    }
+
+    reduce(nErrorFaces, sumOp<label>());
+    reduce(minDet, minOp<scalar>());
+    reduce(sumDet, sumOp<scalar>());
+    reduce(nSummed, sumOp<label>());
+
+    if (debug || report)
+    {
+        if (nSummed > 0)
+        {
+            Info<< "    Face interpolation weight : minimum: " << minDet
+                << " average: " << sumDet/nSummed
+                << endl;
+        }
+    }
+
+    if (nErrorFaces > 0)
+    {
+        if (debug || report)
+        {
+            Info<< " ***Faces with small interpolation found, number of faces: "
+                << nErrorFaces << endl;
+        }
+
+        return true;
+    }
+    else
+    {
+        if (debug || report)
+        {
+            Info<< "    Face interpolation weight check OK." << endl;
+        }
+
+        return false;
+    }
+
+    return false;
+}
+
+
+bool Foam::polyMesh::checkVolRatio
+(
+    const scalarField& cellVols,
+    const bool report,
+    const scalar minRatio,
+    labelHashSet* setPtr
+) const
+{
+    if (debug)
+    {
+        Info<< "bool polyMesh::checkVolRatio(const bool"
+            << ", labelHashSet*) const: "
+            << "checking for low volume ratio" << endl;
+    }
+
+    tmp<scalarField> tvolRatio = polyMeshTools::volRatio(*this, cellVols);
+    scalarField& volRatio = tvolRatio();
+
+
+    label nErrorFaces = 0;
+    scalar minDet = GREAT;
+    scalar sumDet = 0.0;
+    label nSummed = 0;
+
+    // Statistics only for internal and masters of coupled faces
+    PackedBoolList isMasterFace(syncTools::getInternalOrMasterFaces(*this));
+
+    forAll (volRatio, faceI)
+    {
+        if (volRatio[faceI] < minRatio)
+        {
+            // Note: insert both sides of coupled faces
+            if (setPtr)
+            {
+                setPtr->insert(faceI);
+            }
+
+            nErrorFaces++;
+        }
+
+        // Note: statistics only on master of coupled faces
+        if (isMasterFace[faceI])
+        {
+            minDet = min(minDet, volRatio[faceI]);
+            sumDet += volRatio[faceI];
+            nSummed++;
+        }
+    }
+
+    reduce(nErrorFaces, sumOp<label>());
+    reduce(minDet, minOp<scalar>());
+    reduce(sumDet, sumOp<scalar>());
+    reduce(nSummed, sumOp<label>());
+
+    if (debug || report)
+    {
+        if (nSummed > 0)
+        {
+            Info<< "    Face volume ratio : minimum: " << minDet
+                << " average: " << sumDet/nSummed
+                << endl;
+        }
+    }
+
+    if (nErrorFaces > 0)
+    {
+        if (debug || report)
+        {
+            Info<< " ***Faces with small volume ratio found, number of faces: "
+                << nErrorFaces << endl;
+        }
+
+        return true;
+    }
+    else
+    {
+        if (debug || report)
+        {
+            Info<< "    Face volume ratio check OK." << endl;
+        }
+
+        return false;
+    }
+
+    return false;
+}
+
+
 //- Could override checkClosedBoundary to not look at (collocated!) coupled
 //  faces
 //bool Foam::polyMesh::checkClosedBoundary(const bool report) const
@@ -582,6 +766,36 @@ bool Foam::polyMesh::checkCellDeterminant
 }
 
 
+bool Foam::polyMesh::checkFaceWeight
+(
+    const bool report,
+    const scalar minWeight,
+    labelHashSet* setPtr
+) const
+{
+    return checkFaceWeight
+    (
+        faceCentres(),
+        faceAreas(),
+        cellCentres(),
+        report,
+        minWeight,
+        setPtr
+    );
+}
+
+
+bool Foam::polyMesh::checkVolRatio
+(
+    const bool report,
+    const scalar minRatio,
+    labelHashSet* setPtr
+) const
+{
+    return checkVolRatio(cellVolumes(), report, minRatio, setPtr);
+}
+
+
 bool Foam::polyMesh::checkMeshMotion
 (
     const pointField& newPoints,

src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.C

@@ -187,4 +187,112 @@ Foam::tmp<Foam::scalarField> Foam::polyMeshTools::faceSkewness
 }
 
 
+Foam::tmp<Foam::scalarField> Foam::polyMeshTools::faceWeights
+(
+    const polyMesh& mesh,
+    const vectorField& fCtrs,
+    const vectorField& fAreas,
+    const vectorField& cellCtrs
+)
+{
+    const labelList& own = mesh.faceOwner();
+    const labelList& nei = mesh.faceNeighbour();
+    const polyBoundaryMesh& pbm = mesh.boundaryMesh();
+
+    tmp<scalarField> tweight(new scalarField(mesh.nFaces(), 1.0));
+    scalarField& weight = tweight();
+
+    // Internal faces
+    forAll(nei, faceI)
+    {
+        const point& fc = fCtrs[faceI];
+        const vector& fa = fAreas[faceI];
+
+        scalar dOwn = mag(fa & (fc-cellCtrs[own[faceI]]));
+        scalar dNei = mag(fa & (cellCtrs[nei[faceI]]-fc));
+
+        weight[faceI] = min(dNei,dOwn)/(dNei+dOwn+VSMALL);
+    }
+
+
+    // Coupled faces
+
+    pointField neiCc;
+    syncTools::swapBoundaryCellPositions(mesh, cellCtrs, neiCc);
+
+    forAll(pbm, patchI)
+    {
+        const polyPatch& pp = pbm[patchI];
+        if (pp.coupled())
+        {
+            forAll(pp, i)
+            {
+                label faceI = pp.start() + i;
+                label bFaceI = faceI - mesh.nInternalFaces();
+
+                const point& fc = fCtrs[faceI];
+                const vector& fa = fAreas[faceI];
+
+                scalar dOwn = mag(fa & (fc-cellCtrs[own[faceI]]));
+                scalar dNei = mag(fa & (neiCc[bFaceI]-fc));
+
+                weight[faceI] = min(dNei,dOwn)/(dNei+dOwn+VSMALL);
+            }
+        }
+    }
+
+    return tweight;
+}
+
+
+Foam::tmp<Foam::scalarField> Foam::polyMeshTools::volRatio
+(
+    const polyMesh& mesh,
+    const scalarField& vol
+)
+{
+    const labelList& own = mesh.faceOwner();
+    const labelList& nei = mesh.faceNeighbour();
+    const polyBoundaryMesh& pbm = mesh.boundaryMesh();
+
+    tmp<scalarField> tratio(new scalarField(mesh.nFaces(), 1.0));
+    scalarField& ratio = tratio();
+
+    // Internal faces
+    forAll(nei, faceI)
+    {
+        scalar volOwn = vol[own[faceI]];
+        scalar volNei = vol[nei[faceI]];
+
+        ratio[faceI] = min(volOwn,volNei)/(volOwn+volNei+VSMALL);
+    }
+
+
+    // Coupled faces
+
+    scalarField neiVol;
+    syncTools::swapBoundaryCellList(mesh, vol, neiVol);
+
+    forAll(pbm, patchI)
+    {
+        const polyPatch& pp = pbm[patchI];
+        if (pp.coupled())
+        {
+            forAll(pp, i)
+            {
+                label faceI = pp.start() + i;
+                label bFaceI = faceI - mesh.nInternalFaces();
+
+                scalar volOwn = vol[own[faceI]];
+                scalar volNei = neiVol[bFaceI];
+
+                ratio[faceI] = min(volOwn,volNei)/(volOwn+volNei+VSMALL);
+            }
+        }
+    }
+
+    return tratio;
+}
+
+
 // ************************************************************************* //

src/OpenFOAM/meshes/polyMesh/polyMeshCheck/polyMeshTools.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -95,6 +95,22 @@ public:
 //        );
 //    }
 
+    //- Generate interpolation factors field
+    static tmp<scalarField> faceWeights
+    (
+        const polyMesh& mesh,
+        const vectorField& fCtrs,
+        const vectorField& fAreas,
+        const vectorField& cellCtrs
+    );
+
+    //- Generate volume ratio field
+    static tmp<scalarField> volRatio
+    (
+        const polyMesh& mesh,
+        const scalarField& vol
+    );
+
 };
 
 

src/OpenFOAM/meshes/polyMesh/polyPatches/constraint/processorCyclic/processorCyclicPolyPatch.H

@@ -63,8 +63,6 @@ class processorCyclicPolyPatch
         //- Index of originating patch
         mutable label referPatchID_;
 
-    // Private member functions
-
 
 protected:
 
@@ -223,7 +221,6 @@ public:
 
 
     // Destructor
-
     virtual ~processorCyclicPolyPatch();
 
 

src/OpenFOAM/primitives/strings/wordRe/wordRe.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -120,22 +120,25 @@ public:
         inline wordRe(const wordRe&);
 
         //- Construct from keyType
-        inline wordRe(const keyType&, const compOption=LITERAL);
+        inline explicit wordRe(const keyType&);
+
+        //- Construct from keyType
+        inline wordRe(const keyType&, const compOption);
 
         //- Construct as copy of word
-        inline wordRe(const word&);
+        inline explicit wordRe(const word&);
 
         //- Construct as copy of character array
         //  Optionally specify how it should be treated.
-        inline wordRe(const char*, const compOption = LITERAL);
+        inline explicit wordRe(const char*, const compOption = LITERAL);
 
         //- Construct as copy of string.
         //  Optionally specify how it should be treated.
-        inline wordRe(const string&, const compOption = LITERAL);
+        inline explicit wordRe(const string&, const compOption = LITERAL);
 
         //- Construct as copy of std::string
         //  Optionally specify how it should be treated.
-        inline wordRe(const std::string&, const compOption = LITERAL);
+        inline explicit wordRe(const std::string&, const compOption = LITERAL);
 
         //- Construct from Istream
         //  Words are treated as literals, strings with an auto-test

src/OpenFOAM/primitives/strings/wordRe/wordReI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -65,6 +65,18 @@ inline Foam::wordRe::wordRe(const word& str)
 {}
 
 
+inline Foam::wordRe::wordRe(const keyType& str)
+:
+    word(str, false),
+    re_()
+{
+    if (str.isPattern())
+    {
+        compile();
+    }
+}
+
+
 inline Foam::wordRe::wordRe(const keyType& str, const compOption opt)
 :
     word(str, false),

src/dynamicMesh/motionSmoother/motionSmoother.C

@@ -678,34 +678,8 @@ void Foam::motionSmoother::correct()
 }
 
 
-void Foam::motionSmoother::setDisplacement(pointField& patchDisp)
+void Foam::motionSmoother::setDisplacementPatchFields()
 {
-    // See comment in .H file about shared points.
-    const polyBoundaryMesh& patches = mesh_.boundaryMesh();
-
-    forAll(patches, patchI)
-    {
-        const polyPatch& pp = patches[patchI];
-
-        if (pp.coupled())
-        {
-            const labelList& meshPoints = pp.meshPoints();
-
-            forAll(meshPoints, i)
-            {
-                displacement_[meshPoints[i]] = vector::zero;
-            }
-        }
-    }
-
-    const labelList& ppMeshPoints = pp_.meshPoints();
-
-    // Set internal point data from displacement on combined patch points.
-    forAll(ppMeshPoints, patchPointI)
-    {
-        displacement_[ppMeshPoints[patchPointI]] = patchDisp[patchPointI];
-    }
-
     // Adapt the fixedValue bc's (i.e. copy internal point data to
     // boundaryField for all affected patches)
     forAll(adaptPatchIDs_, i)
@@ -765,6 +739,42 @@ void Foam::motionSmoother::setDisplacement(pointField& patchDisp)
         displacement_.boundaryField()[patchI] ==
             displacement_.boundaryField()[patchI].patchInternalField();
     }
+}
+
+
+void Foam::motionSmoother::setDisplacement(pointField& patchDisp)
+{
+    // See comment in .H file about shared points.
+    const polyBoundaryMesh& patches = mesh_.boundaryMesh();
+
+    forAll(patches, patchI)
+    {
+        const polyPatch& pp = patches[patchI];
+
+        if (pp.coupled())
+        {
+            const labelList& meshPoints = pp.meshPoints();
+
+            forAll(meshPoints, i)
+            {
+                displacement_[meshPoints[i]] = vector::zero;
+            }
+        }
+    }
+
+    const labelList& ppMeshPoints = pp_.meshPoints();
+
+    // Set internal point data from displacement on combined patch points.
+    forAll(ppMeshPoints, patchPointI)
+    {
+        displacement_[ppMeshPoints[patchPointI]] = patchDisp[patchPointI];
+    }
+
+
+    // Adapt the fixedValue bc's (i.e. copy internal point data to
+    // boundaryField for all affected patches)
+    setDisplacementPatchFields();
+
 
     if (debug)
     {

src/dynamicMesh/motionSmoother/motionSmoother.H

@@ -385,6 +385,10 @@ public:
             //- Take over existing mesh position.
             void correct();
 
+            //- Set patch fields on displacement to be consistent with
+            //  internal values.
+            void setDisplacementPatchFields();
+
             //- Set displacement field from displacement on patch points.
             //  Modify provided displacement to be consistent with actual
             //  boundary conditions on displacement. Note: resets the

src/edgeMesh/extendedFeatureEdgeMesh/extendedFeatureEdgeMeshTemplates.C

@@ -176,18 +176,18 @@ void Foam::extendedFeatureEdgeMesh::sortPointsAndEdges
         forAll(fpfe, eI)
         {
             const label oldEdgeIndex = fpfe[eI];
-            const label newEdgeIndex = edgeMap[oldEdgeIndex];
 
-            if (newEdgeIndex != -1)
+            const label newFeatureEdgeIndex = edgeMap[oldEdgeIndex];
+
+            if (newFeatureEdgeIndex != -1)
             {
-                newFeatureEdges.append(newEdgeIndex);
+                newFeatureEdges.append(newFeatureEdgeIndex);
             }
         }
 
         featurePointFeatureEdges[pI].transfer(newFeatureEdges);
     }
 
-
     // Reorder the edges by classification
     List<DynamicList<label> > allEds(nEdgeTypes);
 
@@ -272,7 +272,6 @@ void Foam::extendedFeatureEdgeMesh::sortPointsAndEdges
 
     // Initialise sorted edge related data
     edgeDirections_ = edgeDirections/(mag(edgeDirections) + VSMALL);
-
     edgeNormals_ = edgeNormals;
     regionEdges_ = regionEdges;
 

src/finiteVolume/Make/files

@@ -356,6 +356,8 @@ $(snGradSchemes)/faceCorrectedSnGrad/faceCorrectedSnGrads.C
 $(snGradSchemes)/limitedSnGrad/limitedSnGrads.C
 $(snGradSchemes)/uncorrectedSnGrad/uncorrectedSnGrads.C
 $(snGradSchemes)/orthogonalSnGrad/orthogonalSnGrads.C
+$(snGradSchemes)/quadraticFitSnGrad/quadraticFitSnGrads.C
+$(snGradSchemes)/linearFitSnGrad/linearFitSnGrads.C
 
 convectionSchemes = finiteVolume/convectionSchemes
 $(convectionSchemes)/convectionScheme/convectionSchemes.C

src/finiteVolume/fields/fvPatchFields/derived/mappedField/mappedPatchFieldBase.H

@@ -25,7 +25,18 @@ Class
     Foam::mappedPatchFieldBase
 
 Description
-    Functionality for sampling fields using mappedPatchBase.
+    Functionality for sampling fields using mappedPatchBase. Every call to
+    mappedField() returns a sampled field, optionally scaled to maintain an
+    area-weighted average.
+
+    Example usage:
+
+    {
+        fieldName           T;          // default is same as fvPatchField
+        setAverage          false;
+        average             1.0;        // only if setAverage=true
+        interpolationScheme cellPoint;  // default is cell
+    }
 
 SourceFiles
     mappedPatchFieldBase.C

src/finiteVolume/fields/fvPatchFields/derived/pressureInletOutletVelocity/pressureInletOutletVelocityFvPatchVectorField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -29,7 +29,7 @@ Group
 
 Description
     This velocity inlet/outlet boundary condition is applied to pressure
-    boundaries where the pressure is specified.  A zero-gradient condtion is
+    boundaries where the pressure is specified.  A zero-gradient condition is
     applied for outflow (as defined by the flux); for inflow, the velocity is
     obtained from the patch-face normal component of the internal-cell value.
 

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

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -62,7 +62,11 @@ Foam::fv::fourthGrad<Type>::calcGrad
     // Assemble the second-order least-square gradient
     // Calculate the second-order least-square gradient
     tmp<GeometricField<GradType, fvPatchField, volMesh> > tsecondfGrad
-        = leastSquaresGrad<Type>(mesh).grad(vsf);
+      = leastSquaresGrad<Type>(mesh).grad
+        (
+            vsf,
+            "leastSquaresGrad(" + vsf.name() + ")"
+        );
     const GeometricField<GradType, fvPatchField, volMesh>& secondfGrad =
         tsecondfGrad();
 

src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.C

@@ -0,0 +1,270 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "CentredFitSnGradData.H"
+#include "surfaceFields.H"
+#include "SVD.H"
+#include "extendedCentredCellToFaceStencil.H"
+
+// * * * * * * * * * * * * * * * * Constructors * * * * * * * * * * * * * * //
+
+template<class Polynomial>
+Foam::CentredFitSnGradData<Polynomial>::CentredFitSnGradData
+(
+    const fvMesh& mesh,
+    const extendedCentredCellToFaceStencil& stencil,
+    const scalar linearLimitFactor,
+    const scalar centralWeight
+)
+:
+    FitData
+    <
+        CentredFitSnGradData<Polynomial>,
+        extendedCentredCellToFaceStencil,
+        Polynomial
+    >
+    (
+        mesh, stencil, true, linearLimitFactor, centralWeight
+    ),
+    coeffs_(mesh.nFaces())
+{
+    if (debug)
+    {
+        Info<< "Contructing CentredFitSnGradData<Polynomial>" << endl;
+    }
+
+    calcFit();
+
+    if (debug)
+    {
+        Info<< "CentredFitSnGradData<Polynomial>::CentredFitSnGradData() :"
+            << "Finished constructing polynomialFit data"
+            << endl;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class Polynomial>
+void Foam::CentredFitSnGradData<Polynomial>::calcFit
+(
+    scalarList& coeffsi,
+    const List<point>& C,
+    const scalar wLin,
+    const scalar deltaCoeff,
+    const label facei
+)
+{
+    vector idir(1,0,0);
+    vector jdir(0,1,0);
+    vector kdir(0,0,1);
+    this->findFaceDirs(idir, jdir, kdir, facei);
+
+    // Setup the point weights
+    scalarList wts(C.size(), scalar(1));
+    wts[0] = this->centralWeight();
+    wts[1] = this->centralWeight();
+
+    // Reference point
+    point p0 = this->mesh().faceCentres()[facei];
+
+    // p0 -> p vector in the face-local coordinate system
+    vector d;
+
+    // Local coordinate scaling
+    scalar scale = 1;
+
+    // Matrix of the polynomial components
+    scalarRectangularMatrix B(C.size(), this->minSize(), scalar(0));
+
+    forAll(C, ip)
+    {
+        const point& p = C[ip];
+        const vector p0p = p - p0;
+
+        d.x() = p0p & idir;
+        d.y() = p0p & jdir;
+        d.z() = p0p & kdir;
+
+        if (ip == 0)
+        {
+            scale = cmptMax(cmptMag((d)));
+        }
+
+        // Scale the radius vector
+        d /= scale;
+
+        Polynomial::addCoeffs(B[ip], d, wts[ip], this->dim());
+    }
+
+    // Additional weighting for constant and linear terms
+    for (label i = 0; i < B.n(); i++)
+    {
+        B[i][0] *= wts[0];
+        B[i][1] *= wts[0];
+    }
+
+    // Set the fit
+    label stencilSize = C.size();
+    coeffsi.setSize(stencilSize);
+
+    bool goodFit = false;
+    for (int iIt = 0; iIt < 8 && !goodFit; iIt++)
+    {
+        SVD svd(B, SMALL);
+
+        for (label i=0; i<stencilSize; i++)
+        {
+            coeffsi[i] = wts[1]*wts[i]*svd.VSinvUt()[1][i]/scale;
+        }
+
+        goodFit =
+        (
+            mag(wts[0]*wts[0]*svd.VSinvUt()[0][0] - wLin)
+          < this->linearLimitFactor()*wLin)
+         && (mag(wts[0]*wts[1]*svd.VSinvUt()[0][1] - (1 - wLin)
+        ) < this->linearLimitFactor()*(1 - wLin))
+         && coeffsi[0] < 0 && coeffsi[1] > 0
+         && mag(coeffsi[0] + deltaCoeff) < 0.5*deltaCoeff
+         && mag(coeffsi[1] - deltaCoeff) < 0.5*deltaCoeff;
+
+        if (!goodFit)
+        {
+            // (not good fit so increase weight in the centre and weight
+            //  for constant and linear terms)
+
+            WarningIn
+            (
+                "CentredFitSnGradData<Polynomial>::calcFit"
+                "(const List<point>& C, const label facei"
+            )   << "Cannot fit face " << facei << " iteration " << iIt
+                << " with sum of weights " << sum(coeffsi) << nl
+                << "    Weights " << coeffsi << nl
+                << "    Linear weights " << wLin << " " << 1 - wLin << nl
+                << "    deltaCoeff " << deltaCoeff << nl
+                << "    sing vals " << svd.S() << nl
+                << "Components of goodFit:\n"
+                << "    wts[0]*wts[0]*svd.VSinvUt()[0][0] = "
+                << wts[0]*wts[0]*svd.VSinvUt()[0][0] << nl
+                << "    wts[0]*wts[1]*svd.VSinvUt()[0][1] = "
+                << wts[0]*wts[1]*svd.VSinvUt()[0][1]
+                << " dim = " << this->dim() << endl;
+
+            wts[0] *= 10;
+            wts[1] *= 10;
+
+            for (label j = 0; j < B.m(); j++)
+            {
+                B[0][j] *= 10;
+                B[1][j] *= 10;
+            }
+
+            for (label i = 0; i < B.n(); i++)
+            {
+                B[i][0] *= 10;
+                B[i][1] *= 10;
+            }
+        }
+    }
+
+    if (goodFit)
+    {
+        // Remove the uncorrected coefficients
+        coeffsi[0] += deltaCoeff;
+        coeffsi[1] -= deltaCoeff;
+    }
+    else
+    {
+        WarningIn
+        (
+            "CentredFitSnGradData<Polynomial>::calcFit(..)"
+        )   << "Could not fit face " << facei
+            << "    Coefficients = " << coeffsi
+            << ", reverting to uncorrected." << endl;
+
+        coeffsi = 0;
+    }
+}
+
+
+template<class Polynomial>
+void Foam::CentredFitSnGradData<Polynomial>::calcFit()
+{
+    const fvMesh& mesh = this->mesh();
+
+    // Get the cell/face centres in stencil order.
+    // Centred face stencils no good for triangles or tets.
+    // Need bigger stencils
+    List<List<point> > stencilPoints(mesh.nFaces());
+    this->stencil().collectData(mesh.C(), stencilPoints);
+
+    // find the fit coefficients for every face in the mesh
+
+    const surfaceScalarField& w = mesh.surfaceInterpolation::weights();
+    const surfaceScalarField& dC = mesh.nonOrthDeltaCoeffs();
+
+    for (label facei = 0; facei < mesh.nInternalFaces(); facei++)
+    {
+        calcFit
+        (
+            coeffs_[facei],
+            stencilPoints[facei],
+            w[facei],
+            dC[facei],
+            facei
+        );
+    }
+
+    const surfaceScalarField::GeometricBoundaryField& bw = w.boundaryField();
+    const surfaceScalarField::GeometricBoundaryField& bdC = dC.boundaryField();
+
+    forAll(bw, patchi)
+    {
+        const fvsPatchScalarField& pw = bw[patchi];
+        const fvsPatchScalarField& pdC = bdC[patchi];
+
+        if (pw.coupled())
+        {
+            label facei = pw.patch().start();
+
+            forAll(pw, i)
+            {
+                calcFit
+                (
+                    coeffs_[facei],
+                    stencilPoints[facei],
+                    pw[i],
+                    pdC[i],
+                    facei
+                );
+                facei++;
+            }
+        }
+    }
+}
+
+
+// ************************************************************************* //

src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradData.H

@@ -0,0 +1,127 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::CentredFitSnGradData
+
+Description
+    Data for centred fit snGrad schemes
+
+SourceFiles
+    CentredFitSnGradData.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef CentredFitSnGradData_H
+#define CentredFitSnGradData_H
+
+#include "FitData.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class extendedCentredCellToFaceStencil;
+
+/*---------------------------------------------------------------------------*\
+                    Class CentredFitSnGradData Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Polynomial>
+class CentredFitSnGradData
+:
+    public FitData
+    <
+        CentredFitSnGradData<Polynomial>,
+        extendedCentredCellToFaceStencil,
+        Polynomial
+    >
+{
+    // Private data
+
+        //- For each cell in the mesh store the values which multiply the
+        //  values of the stencil to obtain the gradient for each direction
+        List<scalarList> coeffs_;
+
+
+public:
+
+    TypeName("CentredFitSnGradData");
+
+
+    // Constructors
+
+        //- Construct from components
+        CentredFitSnGradData
+        (
+            const fvMesh& mesh,
+            const extendedCentredCellToFaceStencil& stencil,
+            const scalar linearLimitFactor,
+            const scalar centralWeight
+        );
+
+
+    //- Destructor
+    virtual ~CentredFitSnGradData()
+    {}
+
+
+    // Member functions
+
+        //- Return reference to fit coefficients
+        const List<scalarList>& coeffs() const
+        {
+            return coeffs_;
+        }
+
+        //- Calculate the fit for the specified face and set the coefficients
+        void calcFit
+        (
+            scalarList& coeffsi, // coefficients to be set
+            const List<point>&,  // Stencil points
+            const scalar wLin,   // Weight for linear approximation (weights
+                                 // nearest neighbours)
+            const scalar deltaCoeff, // uncorrected delta coefficient
+            const label faci     // Current face index
+        );
+
+        void calcFit();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "CentredFitSnGradData.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/finiteVolume/finiteVolume/snGradSchemes/CentredFitSnGrad/CentredFitSnGradScheme.H

@@ -0,0 +1,183 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::CentredFitSnGradScheme
+
+Description
+    Centred fit snGrad scheme which applies an explicit correction to snGrad
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef CentredFitSnGradScheme_H
+#define CentredFitSnGradScheme_H
+
+#include "CentredFitSnGradData.H"
+#include "snGradScheme.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class fvMesh;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace fv
+{
+/*---------------------------------------------------------------------------*\
+                           Class CentredFitSnGradSnGradScheme Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class Type, class Polynomial, class Stencil>
+class CentredFitSnGradScheme
+:
+    public snGradScheme<Type>
+{
+    // Private Data
+
+        //- Factor the fit is allowed to deviate from linear.
+        //  This limits the amount of high-order correction and increases
+        //  stability on bad meshes
+        const scalar linearLimitFactor_;
+
+        //- Weights for central stencil
+        const scalar centralWeight_;
+
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        CentredFitSnGradScheme(const CentredFitSnGradScheme&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const CentredFitSnGradScheme&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName("CentredFitSnGradScheme");
+
+
+    // Constructors
+
+        //- Construct from mesh and Istream
+        CentredFitSnGradScheme(const fvMesh& mesh, Istream& is)
+        :
+            snGradScheme<Type>(mesh),
+            linearLimitFactor_(readScalar(is)),
+            centralWeight_(1000)
+        {}
+
+
+    // Member Functions
+
+        //- Return the interpolation weighting factors for the given field
+        virtual tmp<surfaceScalarField> deltaCoeffs
+        (
+            const GeometricField<Type, fvPatchField, volMesh>&
+        ) const
+        {
+            return this->mesh().nonOrthDeltaCoeffs();
+        }
+
+        //- Return true if this scheme uses an explicit correction
+        virtual bool corrected() const
+        {
+            return true;
+        }
+
+        //- Return the explicit correction to the face-interpolate
+        virtual tmp<GeometricField<Type, fvsPatchField, surfaceMesh> >
+        correction
+        (
+            const GeometricField<Type, fvPatchField, volMesh>& vf
+        ) const
+        {
+            const fvMesh& mesh = this->mesh();
+
+            const extendedCentredCellToFaceStencil& stencil = Stencil::New
+            (
+                mesh
+            );
+
+            const CentredFitSnGradData<Polynomial>& cfd =
+                CentredFitSnGradData<Polynomial>::New
+                (
+                    mesh,
+                    stencil,
+                    linearLimitFactor_,
+                    centralWeight_
+                );
+
+            tmp<GeometricField<Type, fvsPatchField, surfaceMesh> > sft
+            (
+                stencil.weightedSum(vf, cfd.coeffs())
+            );
+
+            sft().dimensions() /= dimLength;
+
+            return sft;
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace fv
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// Add the patch constructor functions to the hash tables
+
+#define makeCentredFitSnGradTypeScheme(SS, POLYNOMIAL, STENCIL, TYPE) \
+                                                                      \
+typedef CentredFitSnGradScheme<TYPE, POLYNOMIAL, STENCIL>             \
+    CentredFitSnGradScheme##TYPE##POLYNOMIAL##STENCIL##_;             \
+defineTemplateTypeNameAndDebugWithName                                \
+    (CentredFitSnGradScheme##TYPE##POLYNOMIAL##STENCIL##_, #SS, 0);   \
+                                                                      \
+snGradScheme<TYPE>::addMeshConstructorToTable                         \
+<CentredFitSnGradScheme<TYPE, POLYNOMIAL, STENCIL> >                  \
+    add##SS##STENCIL##TYPE##MeshConstructorToTable_;
+
+#define makeCentredFitSnGradScheme(SS, POLYNOMIAL, STENCIL)           \
+                                                                      \
+makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,scalar)          \
+makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,vector)          \
+makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,sphericalTensor) \
+makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,symmTensor)      \
+makeCentredFitSnGradTypeScheme(SS,POLYNOMIAL,STENCIL,tensor)
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/finiteVolume/finiteVolume/snGradSchemes/linearFitSnGrad/linearFitSnGrads.C

@@ -0,0 +1,51 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "CentredFitSnGradScheme.H"
+#include "linearFitPolynomial.H"
+#include "centredFECCellToFaceStencilObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTemplateTypeNameAndDebug
+    (
+        CentredFitSnGradData<linearFitPolynomial>,
+        0
+    );
+
+    namespace fv
+    {
+        makeCentredFitSnGradScheme
+        (
+            linearFit,
+            linearFitPolynomial,
+            centredFECCellToFaceStencilObject
+        );
+    }
+}
+
+// ************************************************************************* //

src/finiteVolume/finiteVolume/snGradSchemes/quadraticFitSnGrad/quadraticFitSnGrads.C

@@ -0,0 +1,51 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "CentredFitSnGradScheme.H"
+#include "quadraticFitPolynomial.H"
+#include "centredCFCCellToFaceStencilObject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTemplateTypeNameAndDebug
+    (
+        CentredFitSnGradData<quadraticFitPolynomial>,
+        0
+    );
+
+    namespace fv
+    {
+        makeCentredFitSnGradScheme
+        (
+            quadraticFit,
+            quadraticFitPolynomial,
+            centredCFCCellToFaceStencilObject
+        );
+    }
+}
+
+// ************************************************************************* //

src/finiteVolume/interpolation/surfaceInterpolation/schemes/FitData/FitData.H

@@ -80,7 +80,7 @@ class FitData
         const label minSize_;
 
 
-    // Private Member Functions
+protected:
 
         //- Find the normal direction (i) and j and k directions for face faci
         void findFaceDirs
@@ -93,9 +93,6 @@ class FitData
 
 public:
 
-    //TypeName("FitData");
-
-
     // Constructors
 
         //- Construct from components
@@ -122,6 +119,30 @@ public:
             return stencil_;
         }
 
+        //- Factor the fit is allowed to deviate from the base scheme
+        scalar linearLimitFactor() const
+        {
+            return linearLimitFactor_;
+        }
+
+        //- Return weight for central stencil
+        scalar centralWeight() const
+        {
+            return centralWeight_;
+        }
+
+        //- Dimensionality of the geometry
+        label dim() const
+        {
+            return dim_;
+        }
+
+        //- Minimum stencil size
+        label minSize() const
+        {
+            return minSize_;
+        }
+
         bool linearCorrection() const
         {
             return linearCorrection_;
@@ -140,7 +161,6 @@ public:
         //- Calculate the fit for all the faces
         virtual void calcFit() = 0;
 
-
         //- Recalculate weights (but not stencil) when the mesh moves
         bool movePoints();
 };

src/lagrangian/intermediate/parcels/Templates/ReactingMultiphaseParcel/ReactingMultiphaseParcel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -477,7 +477,7 @@ void Foam::ReactingMultiphaseParcel<ParcelType>::calc
             const scalar T4 = pow4(this->T_);
             td.cloud().radAreaP()[cellI] += dt*np0*ap;
             td.cloud().radT4()[cellI] += dt*np0*T4;
-            td.cloud().radAreaP()[cellI] += dt*np0*ap*T4;
+            td.cloud().radAreaPT4()[cellI] += dt*np0*ap*T4;
         }
     }
 }

src/lagrangian/intermediate/parcels/Templates/ReactingParcel/ReactingParcel.C

@@ -466,7 +466,7 @@ void Foam::ReactingParcel<ParcelType>::calc
             const scalar T4 = pow4(this->T_);
             td.cloud().radAreaP()[cellI] += dt*np0*ap;
             td.cloud().radT4()[cellI] += dt*np0*T4;
-            td.cloud().radAreaP()[cellI] += dt*np0*ap*T4;
+            td.cloud().radAreaPT4()[cellI] += dt*np0*ap*T4;
         }
     }
 }

src/lagrangian/intermediate/parcels/Templates/ThermoParcel/ThermoParcel.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -260,7 +260,7 @@ void Foam::ThermoParcel<ParcelType>::calc
             const scalar T4 = pow4(this->T_);
             td.cloud().radAreaP()[cellI] += dt*np0*ap;
             td.cloud().radT4()[cellI] += dt*np0*T4;
-            td.cloud().radAreaP()[cellI] += dt*np0*ap*T4;
+            td.cloud().radAreaPT4()[cellI] += dt*np0*ap*T4;
         }
     }
 }

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/autoSnapDriverFeature.C

@@ -767,7 +767,7 @@ void Foam::autoSnapDriver::binFeatureFace
 ) const
 {
     // What to do with very far attraction? For now just ignore the face
-    if (magSqr(faceDisp) < sqr(snapDist))
+    if (magSqr(faceDisp) < sqr(snapDist) && mag(faceSurfaceNormal) > VSMALL)
     {
         const point pt = fc + faceDisp;
 

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -138,6 +138,7 @@ Foam::layerParameters::layerParameters
         readLabel(dict.lookup("nSmoothSurfaceNormals"))
     ),
     nSmoothNormals_(readLabel(dict.lookup("nSmoothNormals"))),
+    nSmoothDisplacement_(dict.lookupOrDefault("nSmoothDisplacement", 0)),
     nSmoothThickness_(readLabel(dict.lookup("nSmoothThickness"))),
     maxFaceThicknessRatio_
     (
@@ -278,7 +279,7 @@ Foam::layerParameters::layerParameters
             const keyType& key = iter().keyword();
             const labelHashSet patchIDs
             (
-                boundaryMesh.patchSet(List<wordRe>(1, key))
+                boundaryMesh.patchSet(List<wordRe>(1, wordRe(key)))
             );
 
             if (patchIDs.size() == 0)

src/mesh/autoMesh/autoHexMesh/autoHexMeshDriver/layerParameters/layerParameters.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -119,6 +119,8 @@ private:
 
         label nSmoothNormals_;
 
+        label nSmoothDisplacement_;
+
         label nSmoothThickness_;
 
         scalar maxFaceThicknessRatio_;
@@ -275,6 +277,12 @@ public:
             return layerTerminationCos_;
         }
 
+        //- Smooth internal displacement
+        label nSmoothDisplacement() const
+        {
+            return nSmoothDisplacement_;
+        }
+
         //- Smooth layer thickness over surface patches
         label nSmoothThickness() const
         {

src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolation.H

@@ -50,14 +50,11 @@ SourceFiles
 #define AMIInterpolation_H
 
 #include "className.H"
-#include "DynamicList.H"
 #include "searchableSurface.H"
+#include "treeBoundBoxList.H"
 #include "boolList.H"
 #include "primitivePatch.H"
 #include "faceAreaIntersect.H"
-#include "indexedOctree.H"
-#include "treeDataPrimitivePatch.H"
-#include "treeBoundBoxList.H"
 #include "globalIndex.H"
 #include "ops.H"
 
@@ -82,12 +79,38 @@ class AMIInterpolation
 :
     public AMIInterpolationName
 {
-    //- local typedef to octree tree-type
+public:
-    typedef treeDataPrimitivePatch<TargetPatch> treeType;
 
+    // Public data types
+
+        //- Enumeration specifying interpolation method
+        enum interpolationMethod
+        {
+            imDirect,
+            imMapNearest,
+            imFaceAreaWeight
+        };
+
+        //- Convert interpolationMethod to word representation
+        static word interpolationMethodToWord
+        (
+            const interpolationMethod& method
+        );
+
+        //- Convert word to interpolationMethod
+        static interpolationMethod wordTointerpolationMethod
+        (
+            const word& method
+        );
+
+
+private:
 
     // Private data
 
+        //- Interpolation method
+        interpolationMethod method_;
+
         //- Flag to indicate that the two patches are co-directional and
         //  that the orientation of the target patch should be reversed
         const bool reverseTarget_;
@@ -96,6 +119,7 @@ class AMIInterpolation
         //  cases
         label singlePatchProc_;
 
+
         // Source patch
 
             //- Source face areas
@@ -110,10 +134,6 @@ class AMIInterpolation
             //- Sum of weights of target faces per source face
             scalarField srcWeightsSum_;
 
-            //- Labels of faces that are not overlapped by any target faces
-            //  (should be empty for correct functioning)
-            labelList srcNonOverlap_;
-
 
         // Target patch
 
@@ -130,9 +150,6 @@ class AMIInterpolation
             scalarField tgtWeightsSum_;
 
 
-        //- Octree used to find face seeds
-        autoPtr<indexedOctree<treeType> > treePtr_;
-
         //- Face triangulation mode
         const faceAreaIntersect::triangulationMode triMode_;
 
@@ -152,30 +169,6 @@ class AMIInterpolation
         void operator=(const AMIInterpolation&);
 
 
-        // Helper functions
-
-            //- Write triangle intersection to OBJ file
-            void writeIntersectionOBJ
-            (
-                const scalar area,
-                const face& f1,
-                const face& f2,
-                const pointField& f1Points,
-                const pointField& f2Points
-            ) const;
-
-            //- Check that patches are valid
-            void checkPatches
-            (
-                const SourcePatch& srcPatch,
-                const TargetPatch& tgtPatch
-            ) const;
-
-            //- Reset the octree for the traget patch face search
-            void resetTree(const TargetPatch& tgtPatch);
-
-
-
         // Parallel functionality
 
             //- Calculate if patches are on multiple processors
@@ -229,83 +222,8 @@ class AMIInterpolation
             ) const;
 
 
-        // Marching front
-
-            //- Find face on target patch that overlaps source face
-            label findTargetFace
-            (
-                const label srcFaceI,
-                const SourcePatch& srcPatch
-            ) const;
-
-            //- Add faces neighbouring faceI to the ID list
-            void appendNbrFaces
-            (
-                const label faceI,
-                const TargetPatch& patch,
-                const DynamicList<label>& visitedFaces,
-                DynamicList<label>& faceIDs
-            ) const;
-
-            bool processSourceFace
-            (
-                const SourcePatch& srcPatch,
-                const TargetPatch& tgtPatch,
-                const label srcFaceI,
-                const label tgtStartFaceI,
-
-                DynamicList<label>& nbrFaces,
-                DynamicList<label>& visitedFaces,
-                List<DynamicList<label> >& srcAddr,
-                List<DynamicList<scalar> >& srcWght,
-                List<DynamicList<label> >& tgtAddr,
-                List<DynamicList<scalar> >& tgtWght
-            );
-
-            void restartUncoveredSourceFace
-            (
-                const SourcePatch& srcPatch,
-                const TargetPatch& tgtPatch,
-                List<DynamicList<label> >& srcAddr,
-                List<DynamicList<scalar> >& srcWght,
-                List<DynamicList<label> >& tgtAddr,
-                List<DynamicList<scalar> >& tgtWght
-            );
-
-            //- Set the source and target seed faces
-            void setNextFaces
-            (
-                label& startSeedI,
-                label& srcFaceI,
-                label& tgtFaceI,
-                const SourcePatch& srcPatch0,
-                const TargetPatch& tgtPatch0,
-                const boolList& mapFlag,
-                labelList& seedFaces,
-                const DynamicList<label>& visitedFaces
-            ) const;
-
-
         // Evaluation
 
-            //- Area of intersection between source and target faces
-            scalar interArea
-            (
-                const label srcFaceI,
-                const label tgtFaceI,
-                const SourcePatch& srcPatch,
-                const TargetPatch& tgtPatch
-            ) const;
-
-            //- Calculate addressing
-            void calcAddressing
-            (
-                const SourcePatch& srcPatch,
-                const TargetPatch& tgtPatch,
-                label srcFaceI = -1,
-                label tgtFaceI = -1
-            );
-
             //- Normalise the (area) weights - suppresses numerical error in
             //  weights calculation
             //  NOTE: if area weights are incorrect by 'a significant amount'
@@ -352,6 +270,7 @@ public:
             const SourcePatch& srcPatch,
             const TargetPatch& tgtPatch,
             const faceAreaIntersect::triangulationMode& triMode,
+            const interpolationMethod& method = imFaceAreaWeight,
             const bool reverseTarget = false
         );
 
@@ -362,6 +281,7 @@ public:
             const TargetPatch& tgtPatch,
             const autoPtr<searchableSurface>& surf,
             const faceAreaIntersect::triangulationMode& triMode,
+            const interpolationMethod& method = imFaceAreaWeight,
             const bool reverseTarget = false
         );
 
@@ -378,6 +298,12 @@ public:
     //- Destructor
     ~AMIInterpolation();
 
+    // Typedef to SourcePatch type this AMIInterplation is instantiated on
+    typedef SourcePatch sourcePatchType;
+
+    // Typedef to TargetPatch type this AMIInterplation is instantiated on
+    typedef TargetPatch targetPatchType;
+
 
     // Member Functions
 
@@ -403,10 +329,6 @@ public:
                 //  patch weights (i.e. the sum before normalisation)
                 inline const scalarField& srcWeightsSum() const;
 
-                //- Labels of faces that are not overlapped by any target faces
-                //  (should be empty for correct functioning)
-                inline const labelList& srcNonOverlap() const;
-
                 //- Source map pointer - valid only if singlePatchProc = -1
                 //  This gets source data into a form to be consumed by
                 //  tgtAddress, tgtWeights

src/meshTools/AMIInterpolation/AMIInterpolation/AMIInterpolationI.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -55,14 +55,6 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeights() const
 }
 
 
-template<class SourcePatch, class TargetPatch>
-inline const Foam::labelList&
-Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcNonOverlap() const
-{
-    return srcNonOverlap_;
-}
-
-
 template<class SourcePatch, class TargetPatch>
 inline const Foam::scalarField&
 Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeightsSum() const

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.C

@@ -0,0 +1,339 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "AMIMethod.H"
+#include "meshTools.H"
+#include "mapDistribute.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::checkPatches() const
+{
+    if (debug && (!srcPatch_.size() || !tgtPatch_.size()))
+    {
+        Pout<< "AMI: Patches not on processor: Source faces = "
+            << srcPatch_.size() << ", target faces = " << tgtPatch_.size()
+            << endl;
+    }
+
+
+    const scalar maxBoundsError = 0.05;
+
+    // check bounds of source and target
+    boundBox bbSrc(srcPatch_.points(), srcPatch_.meshPoints(), true);
+    boundBox bbTgt(tgtPatch_.points(), tgtPatch_.meshPoints(), true);
+
+    boundBox bbTgtInf(bbTgt);
+    bbTgtInf.inflate(maxBoundsError);
+
+    if (!bbTgtInf.contains(bbSrc))
+    {
+        WarningIn("AMIMethod<SourcePatch, TargetPatch>::checkPatches()")
+            << "Source and target patch bounding boxes are not similar" << nl
+            << "    source box span     : " << bbSrc.span() << nl
+            << "    target box span     : " << bbTgt.span() << nl
+            << "    source box          : " << bbSrc << nl
+            << "    target box          : " << bbTgt << nl
+            << "    inflated target box : " << bbTgtInf << endl;
+    }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+bool Foam::AMIMethod<SourcePatch, TargetPatch>::initialise
+(
+    labelListList& srcAddress,
+    scalarListList& srcWeights,
+    labelListList& tgtAddress,
+    scalarListList& tgtWeights,
+    label& srcFaceI,
+    label& tgtFaceI
+)
+{
+    // set initial sizes for weights and addressing - must be done even if
+    // returns false below
+    srcAddress.setSize(srcPatch_.size());
+    srcWeights.setSize(srcPatch_.size());
+    tgtAddress.setSize(tgtPatch_.size());
+    tgtWeights.setSize(tgtPatch_.size());
+
+    // check that patch sizes are valid
+    if (!srcPatch_.size())
+    {
+        return false;
+    }
+    else if (!tgtPatch_.size())
+    {
+        WarningIn
+        (
+            "void Foam::AMIMethod<SourcePatch, TargetPatch>::initialise"
+            "("
+                "label&, "
+                "label&"
+            ")"
+        )
+            << srcPatch_.size() << " source faces but no target faces" << endl;
+
+        return false;
+    }
+
+    // reset the octree
+    resetTree();
+
+    // find initial face match using brute force/octree search
+    if ((srcFaceI == -1) || (tgtFaceI == -1))
+    {
+        srcFaceI = 0;
+        tgtFaceI = 0;
+        bool foundFace = false;
+        forAll(srcPatch_, faceI)
+        {
+            tgtFaceI = findTargetFace(faceI);
+            if (tgtFaceI >= 0)
+            {
+                srcFaceI = faceI;
+                foundFace = true;
+                break;
+            }
+        }
+
+        if (!foundFace)
+        {
+            FatalErrorIn
+            (
+                "void Foam::AMIMethod<SourcePatch, TargetPatch>::initialise"
+                "("
+                    "label&, "
+                    "label&"
+                ")"
+            )   << "Unable to find initial target face" << abort(FatalError);
+        }
+    }
+
+    if (debug)
+    {
+        Pout<< "AMI: initial target face = " << tgtFaceI << endl;
+    }
+
+    return true;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::writeIntersectionOBJ
+(
+    const scalar area,
+    const face& f1,
+    const face& f2,
+    const pointField& f1Points,
+    const pointField& f2Points
+) const
+{
+    static label count = 1;
+
+    const pointField f1pts = f1.points(f1Points);
+    const pointField f2pts = f2.points(f2Points);
+
+    Pout<< "Face intersection area (" << count <<  "):" << nl
+        << "    f1 face = " << f1 << nl
+        << "    f1 pts  = " << f1pts << nl
+        << "    f2 face = " << f2 << nl
+        << "    f2 pts  = " << f2pts << nl
+        << "    area    = " << area
+        << endl;
+
+    OFstream os("areas" + name(count) + ".obj");
+
+    forAll(f1pts, i)
+    {
+        meshTools::writeOBJ(os, f1pts[i]);
+    }
+    os<< "l";
+    forAll(f1pts, i)
+    {
+        os<< " " << i + 1;
+    }
+    os<< " 1" << endl;
+
+
+    forAll(f2pts, i)
+    {
+        meshTools::writeOBJ(os, f2pts[i]);
+    }
+    os<< "l";
+    forAll(f2pts, i)
+    {
+        os<< " " << f1pts.size() + i + 1;
+    }
+    os<< " " << f1pts.size() + 1 << endl;
+
+    count++;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::resetTree()
+{
+    // Clear the old octree
+    treePtr_.clear();
+
+    treeBoundBox bb(tgtPatch_.points());
+    bb.inflate(0.01);
+
+    if (!treePtr_.valid())
+    {
+        treePtr_.reset
+        (
+            new indexedOctree<treeType>
+            (
+                treeType(false, tgtPatch_),
+                bb,                         // overall search domain
+                8,                          // maxLevel
+                10,                         // leaf size
+                3.0                         // duplicity
+            )
+        );
+    }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+Foam::label Foam::AMIMethod<SourcePatch, TargetPatch>::findTargetFace
+(
+    const label srcFaceI
+) const
+{
+    label targetFaceI = -1;
+
+    const pointField& srcPts = srcPatch_.points();
+    const face& srcFace = srcPatch_[srcFaceI];
+    const point srcPt = srcFace.centre(srcPts);
+    const scalar srcFaceArea = srcMagSf_[srcFaceI];
+
+    pointIndexHit sample = treePtr_->findNearest(srcPt, 10.0*srcFaceArea);
+
+
+    if (debug)
+    {
+        Pout<< "Source point = " << srcPt << ", Sample point = "
+            << sample.hitPoint() << ", Sample index = " << sample.index()
+            << endl;
+    }
+
+    if (sample.hit())
+    {
+        targetFaceI = sample.index();
+    }
+
+    return targetFaceI;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIMethod<SourcePatch, TargetPatch>::appendNbrFaces
+(
+    const label faceI,
+    const TargetPatch& patch,
+    const DynamicList<label>& visitedFaces,
+    DynamicList<label>& faceIDs
+) const
+{
+    const labelList& nbrFaces = patch.faceFaces()[faceI];
+
+    // filter out faces already visited from src face neighbours
+    forAll(nbrFaces, i)
+    {
+        label nbrFaceI = nbrFaces[i];
+        bool valid = true;
+        forAll(visitedFaces, j)
+        {
+            if (nbrFaceI == visitedFaces[j])
+            {
+                valid = false;
+                break;
+            }
+        }
+
+        if (valid)
+        {
+            forAll(faceIDs, j)
+            {
+                if (nbrFaceI == faceIDs[j])
+                {
+                    valid = false;
+                    break;
+                }
+            }
+        }
+
+        if (valid)
+        {
+            faceIDs.append(nbrFaceI);
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::AMIMethod<SourcePatch, TargetPatch>::AMIMethod
+(
+    const SourcePatch& srcPatch,
+    const TargetPatch& tgtPatch,
+    const scalarField& srcMagSf,
+    const scalarField& tgtMagSf,
+    const faceAreaIntersect::triangulationMode& triMode,
+    const bool reverseTarget
+)
+:
+    srcPatch_(srcPatch),
+    tgtPatch_(tgtPatch),
+    reverseTarget_(reverseTarget),
+    srcMagSf_(srcMagSf),
+    tgtMagSf_(tgtMagSf),
+    srcNonOverlap_(),
+    triMode_(triMode)
+{
+    checkPatches();
+
+    label srcSize = returnReduce(srcPatch_.size(), sumOp<label>());
+    label tgtSize = returnReduce(tgtPatch_.size(), sumOp<label>());
+
+    IInfo<< "AMI: Creating addressing and weights between "
+        << srcSize << " source faces and " << tgtSize << " target faces"
+        << endl;
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::AMIMethod<SourcePatch, TargetPatch>::~AMIMethod()
+{}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethod.H

@@ -0,0 +1,268 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::AMIMethod
+
+Description
+    Base class for Arbitrary Mesh Interface (AMI) methods
+
+SourceFiles
+    AMIMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef AMIMethod_H
+#define AMIMethod_H
+
+#include "className.H"
+#include "DynamicList.H"
+#include "faceAreaIntersect.H"
+#include "indexedOctree.H"
+#include "treeDataPrimitivePatch.H"
+#include "treeBoundBoxList.H"
+#include "runTimeSelectionTables.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class AMIMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class AMIMethod
+{
+
+private:
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        AMIMethod(const AMIMethod&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const AMIMethod&);
+
+
+protected:
+
+    //- local typedef to octree tree-type
+    typedef treeDataPrimitivePatch<TargetPatch> treeType;
+
+
+    // Protected data
+
+        //- Reference to source patch
+        const SourcePatch& srcPatch_;
+
+        //- Reference to target patch
+        const SourcePatch& tgtPatch_;
+
+        //- Flag to indicate that the two patches are co-directional and
+        //  that the orientation of the target patch should be reversed
+        const bool reverseTarget_;
+
+        //- Source face areas
+        const scalarField& srcMagSf_;
+
+        //- Target face areas
+        const scalarField& tgtMagSf_;
+
+        //- Labels of faces that are not overlapped by any target faces
+        //  (should be empty for correct functioning)
+        labelList srcNonOverlap_;
+
+        //- Octree used to find face seeds
+        autoPtr<indexedOctree<treeType> > treePtr_;
+
+        //- Face triangulation mode
+        const faceAreaIntersect::triangulationMode triMode_;
+
+
+    // Protected Member Functions
+
+        // Helper functions
+
+            //- Check AMI patch coupling
+            void checkPatches() const;
+
+            //- Initialise and return true if all ok
+            bool initialise
+            (
+                labelListList& srcAddress,
+                scalarListList& srcWeights,
+                labelListList& tgtAddress,
+                scalarListList& tgtWeights,
+                label& srcFaceI,
+                label& tgtFaceI
+            );
+
+            //- Write triangle intersection to OBJ file
+            void writeIntersectionOBJ
+            (
+                const scalar area,
+                const face& f1,
+                const face& f2,
+                const pointField& f1Points,
+                const pointField& f2Points
+            ) const;
+
+
+        // Common AMI method functions
+
+            //- Reset the octree for the target patch face search
+            void resetTree();
+
+            //- Find face on target patch that overlaps source face
+            label findTargetFace(const label srcFaceI) const;
+
+            //- Add faces neighbouring faceI to the ID list
+            void appendNbrFaces
+            (
+                const label faceI,
+                const TargetPatch& patch,
+                const DynamicList<label>& visitedFaces,
+                DynamicList<label>& faceIDs
+            ) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("AMIMethod");
+
+    //- Declare runtime constructor selection table
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        AMIMethod,
+        components,
+        (
+            const SourcePatch& srcPatch,
+            const TargetPatch& tgtPatch,
+            const scalarField& srcMagSf,
+            const scalarField& tgtMagSf,
+            const faceAreaIntersect::triangulationMode& triMode,
+            const bool reverseTarget
+        ),
+        (srcPatch, tgtPatch, srcMagSf, tgtMagSf, triMode, reverseTarget)
+    );
+
+    //- Construct from components
+    AMIMethod
+    (
+        const SourcePatch& srcPatch,
+        const TargetPatch& tgtPatch,
+        const scalarField& srcMagSf,
+        const scalarField& tgtMagSf,
+        const faceAreaIntersect::triangulationMode& triMode,
+        const bool reverseTarget
+    );
+
+    //- Selector
+    static autoPtr<AMIMethod> New
+    (
+        const word& methodName,
+        const SourcePatch& srcPatch,
+        const TargetPatch& tgtPatch,
+        const scalarField& srcMagSf,
+        const scalarField& tgtMagSf,
+        const faceAreaIntersect::triangulationMode& triMode,
+        const bool reverseTarget
+    );
+
+
+    //- Destructor
+    virtual ~AMIMethod();
+
+
+    // Member Functions
+
+        // Access
+
+            //- Labels of faces that are not overlapped by any target faces
+            //  Note: this should be empty for correct functioning
+            inline const labelList& srcNonOverlap() const;
+
+
+        // Manipulation
+
+            //- Update addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcAddress,
+                scalarListList& srcWeights,
+                labelListList& tgtAddress,
+                scalarListList& tgtWeights,
+                label srcFaceI = -1,
+                label tgtFaceI = -1
+            ) = 0;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#define makeAMIMethod(AMIType)                                                \
+                                                                              \
+    typedef AMIMethod<AMIType::sourcePatchType,AMIType::targetPatchType>      \
+        AMIMethod##AMIType;                                                   \
+                                                                              \
+    defineNamedTemplateTypeNameAndDebug(AMIMethod##AMIType, 0);               \
+    defineTemplateRunTimeSelectionTable(AMIMethod##AMIType, components);
+
+
+#define makeAMIMethodType(AMIType, Method)                                    \
+                                                                              \
+    typedef Method<AMIType::sourcePatchType,AMIType::targetPatchType>         \
+        Method##AMIType;                                                      \
+                                                                              \
+    defineNamedTemplateTypeNameAndDebug(Method##AMIType, 0);                  \
+                                                                              \
+    AMIMethod<AMIType::sourcePatchType,AMIType::targetPatchType>::            \
+        addcomponentsConstructorToTable<Method##AMIType>                      \
+        add##Method##AMIType##ConstructorToTable_;
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "AMIMethodI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "AMIMethod.C"
+#   include "AMIMethodNew.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodI.H

@@ -0,0 +1,34 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+inline const Foam::labelList&
+Foam::AMIMethod<SourcePatch, TargetPatch>::srcNonOverlap() const
+{
+    return srcNonOverlap_;
+}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/AMIMethod/AMIMethodNew.C

@@ -0,0 +1,84 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::autoPtr<Foam::AMIMethod<SourcePatch, TargetPatch> >
+Foam::AMIMethod<SourcePatch, TargetPatch>::New
+(
+    const word& methodName,
+    const SourcePatch& srcPatch,
+    const TargetPatch& tgtPatch,
+    const scalarField& srcMagSf,
+    const scalarField& tgtMagSf,
+    const faceAreaIntersect::triangulationMode& triMode,
+    const bool reverseTarget
+)
+{
+    if (debug)
+    {
+        Info<< "Selecting AMIMethod " << methodName << endl;
+    }
+
+    typename componentsConstructorTable::iterator cstrIter =
+        componentsConstructorTablePtr_->find(methodName);
+
+    if (cstrIter == componentsConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "AMIMethod<SourcePatch, TargetPatch>::New"
+            "("
+                "const word&, "
+                "const SourcePatch&, "
+                "const TargetPatch&, "
+                "const scalarField&, "
+                "const scalarField&, "
+                "const faceAreaIntersect::triangulationMode&, "
+                "const bool"
+            ")"
+        )   << "Unknown AMIMethod type "
+            << methodName << nl << nl
+            << "Valid AMIMethod types are:" << nl
+            << componentsConstructorTablePtr_->sortedToc() << exit(FatalError);
+    }
+
+    return autoPtr<AMIMethod<SourcePatch, TargetPatch> >
+    (
+        cstrIter()
+        (
+            srcPatch,
+            tgtPatch,
+            srcMagSf,
+            tgtMagSf,
+            triMode,
+            reverseTarget
+        )
+    );
+}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.C

@@ -0,0 +1,224 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "directAMI.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::directAMI<SourcePatch, TargetPatch>::appendToDirectSeeds
+(
+    boolList& mapFlag,
+    labelList& srcTgtSeed,
+    DynamicList<label>& srcSeeds,
+    DynamicList<label>& nonOverlapFaces,
+    label& srcFaceI,
+    label& tgtFaceI
+) const
+{
+    const labelList& srcNbr = this->srcPatch_.faceFaces()[srcFaceI];
+    const labelList& tgtNbr = this->tgtPatch_.faceFaces()[tgtFaceI];
+
+    const pointField& srcPoints = this->srcPatch_.points();
+    const pointField& tgtPoints = this->tgtPatch_.points();
+
+    const vectorField& srcCf = this->srcPatch_.faceCentres();
+
+    forAll(srcNbr, i)
+    {
+        label srcI = srcNbr[i];
+
+        if (mapFlag[srcI] && srcTgtSeed[srcI] == -1)
+        {
+            const face& srcF = this->srcPatch_[srcI];
+            const vector srcN = srcF.normal(srcPoints);
+
+            bool found = false;
+            forAll(tgtNbr, j)
+            {
+                label tgtI = tgtNbr[j];
+                const face& tgtF = this->tgtPatch_[tgtI];
+                pointHit ray = tgtF.ray(srcCf[srcI], srcN, tgtPoints);
+
+                if (ray.hit())
+                {
+                    // new match - append to lists
+                    found = true;
+
+                    srcTgtSeed[srcI] = tgtI;
+                    srcSeeds.append(srcI);
+
+                    break;
+                }
+            }
+
+            if (!found)
+            {
+                // no match available for source face srcI
+                mapFlag[srcI] = false;
+                nonOverlapFaces.append(srcI);
+            }
+        }
+    }
+
+    if (srcSeeds.size())
+    {
+        srcFaceI = srcSeeds.remove();
+        tgtFaceI = srcTgtSeed[srcFaceI];
+    }
+    else
+    {
+        srcFaceI = -1;
+        tgtFaceI = -1;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::directAMI<SourcePatch, TargetPatch>::directAMI
+(
+    const SourcePatch& srcPatch,
+    const TargetPatch& tgtPatch,
+    const scalarField& srcMagSf,
+    const scalarField& tgtMagSf,
+    const faceAreaIntersect::triangulationMode& triMode,
+    const bool reverseTarget
+)
+:
+    AMIMethod<SourcePatch, TargetPatch>
+    (
+        srcPatch,
+        tgtPatch,
+        srcMagSf,
+        tgtMagSf,
+        triMode,
+        reverseTarget
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::directAMI<SourcePatch, TargetPatch>::~directAMI()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::directAMI<SourcePatch, TargetPatch>::calculate
+(
+    labelListList& srcAddress,
+    scalarListList& srcWeights,
+    labelListList& tgtAddress,
+    scalarListList& tgtWeights,
+    label srcFaceI,
+    label tgtFaceI
+)
+{
+    bool ok =
+        this->initialise
+        (
+            srcAddress,
+            srcWeights,
+            tgtAddress,
+            tgtWeights,
+            srcFaceI,
+            tgtFaceI
+        );
+
+    if (!ok)
+    {
+        return;
+    }
+
+
+    // temporary storage for addressing and weights
+    List<DynamicList<label> > srcAddr(this->srcPatch_.size());
+    List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
+
+
+    // construct weights and addressing
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // list of faces currently visited for srcFaceI to avoid multiple hits
+    DynamicList<label> srcSeeds(10);
+
+    // list to keep track of tgt faces used to seed src faces
+    labelList srcTgtSeed(srcAddr.size(), -1);
+    srcTgtSeed[srcFaceI] = tgtFaceI;
+
+    // list to keep track of whether src face can be mapped
+    boolList mapFlag(srcAddr.size(), true);
+
+    DynamicList<label> nonOverlapFaces;
+    do
+    {
+        srcAddr[srcFaceI].append(tgtFaceI);
+        tgtAddr[tgtFaceI].append(srcFaceI);
+
+        mapFlag[srcFaceI] = false;
+
+        // Do advancing front starting from srcFaceI, tgtFaceI
+        appendToDirectSeeds
+        (
+            mapFlag,
+            srcTgtSeed,
+            srcSeeds,
+            nonOverlapFaces,
+            srcFaceI,
+            tgtFaceI
+        );
+    } while (srcFaceI >= 0);
+
+    if (nonOverlapFaces.size() != 0)
+    {
+        Pout<< "    AMI: " << nonOverlapFaces.size()
+            << " non-overlap faces identified"
+            << endl;
+
+        this->srcNonOverlap_.transfer(nonOverlapFaces);
+    }
+
+    // transfer data to persistent storage
+    forAll(srcAddr, i)
+    {
+        scalar magSf = this->srcMagSf_[i];
+        srcAddress[i].transfer(srcAddr[i]);
+        srcWeights[i] = scalarList(1, magSf);
+    }
+    forAll(tgtAddr, i)
+    {
+        scalar magSf = this->tgtMagSf_[i];
+        tgtAddress[i].transfer(tgtAddr[i]);
+        tgtWeights[i] = scalarList(1, magSf);
+    }
+}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/directAMI/directAMI.H

@@ -0,0 +1,144 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::directAMI
+
+Description
+    Direct mapped Arbitrary Mesh Interface (AMI) method
+
+SourceFiles
+    directAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef directAMI_H
+#define directAMI_H
+
+#include "AMIMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                          Class directAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class directAMI
+:
+    public AMIMethod<SourcePatch, TargetPatch>
+{
+
+private:
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        directAMI(const directAMI&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const directAMI&);
+
+        // Marching front
+
+            //- Append to list of src face seed indices
+            void appendToDirectSeeds
+            (
+                boolList& mapFlag,
+                labelList& srcTgtSeed,
+                DynamicList<label>& srcSeeds,
+                DynamicList<label>& nonOverlapFaces,
+                label& srcFaceI,
+                label& tgtFaceI
+            ) const;
+
+
+        // Evaluation
+
+            //- Area of intersection between source and target faces
+            scalar interArea
+            (
+                const label srcFaceI,
+                const label tgtFaceI
+            ) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("directAMI");
+
+
+    // Constructors
+
+        //- Construct from components
+        directAMI
+        (
+            const SourcePatch& srcPatch,
+            const TargetPatch& tgtPatch,
+            const scalarField& srcMagSf,
+            const scalarField& tgtMagSf,
+            const faceAreaIntersect::triangulationMode& triMode,
+            const bool reverseTarget = false
+        );
+
+
+    //- Destructor
+    virtual ~directAMI();
+
+
+    // Member Functions
+
+        // Manipulation
+
+            //- Update addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcAddress,
+                scalarListList& srcWeights,
+                labelListList& tgtAddress,
+                scalarListList& tgtWeights,
+                label srcFaceI = -1,
+                label tgtFaceI = -1
+            );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "directAMI.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.C

@@ -0,0 +1,553 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "faceAreaWeightAMI.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+bool Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::processSourceFace
+(
+    const label srcFaceI,
+    const label tgtStartFaceI,
+
+    // list of tgt face neighbour faces
+    DynamicList<label>& nbrFaces,
+    // list of faces currently visited for srcFaceI to avoid multiple hits
+    DynamicList<label>& visitedFaces,
+
+    // temporary storage for addressing and weights
+    List<DynamicList<label> >& srcAddr,
+    List<DynamicList<scalar> >& srcWght,
+    List<DynamicList<label> >& tgtAddr,
+    List<DynamicList<scalar> >& tgtWght
+)
+{
+    nbrFaces.clear();
+    visitedFaces.clear();
+
+    // append initial target face and neighbours
+    nbrFaces.append(tgtStartFaceI);
+    this->appendNbrFaces
+    (
+        tgtStartFaceI,
+        this->tgtPatch_,
+        visitedFaces,
+        nbrFaces
+    );
+
+    bool faceProcessed = false;
+
+    do
+    {
+        // process new target face
+        label tgtFaceI = nbrFaces.remove();
+        visitedFaces.append(tgtFaceI);
+        scalar area = interArea(srcFaceI, tgtFaceI);
+
+        // store when intersection area > 0
+        if (area > 0)
+        {
+            srcAddr[srcFaceI].append(tgtFaceI);
+            srcWght[srcFaceI].append(area);
+
+            tgtAddr[tgtFaceI].append(srcFaceI);
+            tgtWght[tgtFaceI].append(area);
+
+            this->appendNbrFaces
+            (
+                tgtFaceI,
+                this->tgtPatch_,
+                visitedFaces,
+                nbrFaces
+            );
+
+            faceProcessed = true;
+        }
+
+    } while (nbrFaces.size() > 0);
+
+    return faceProcessed;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::setNextFaces
+(
+    label& startSeedI,
+    label& srcFaceI,
+    label& tgtFaceI,
+    const boolList& mapFlag,
+    labelList& seedFaces,
+    const DynamicList<label>& visitedFaces
+) const
+{
+    const labelList& srcNbrFaces = this->srcPatch_.faceFaces()[srcFaceI];
+
+    // set possible seeds for later use
+    bool valuesSet = false;
+    forAll(srcNbrFaces, i)
+    {
+        label faceS = srcNbrFaces[i];
+
+        if (mapFlag[faceS] && seedFaces[faceS] == -1)
+        {
+            forAll(visitedFaces, j)
+            {
+                label faceT = visitedFaces[j];
+                scalar area = interArea(faceS, faceT);
+                scalar areaTotal = this->srcMagSf_[srcFaceI];
+
+                // Check that faces have enough overlap for robust walking
+                if (area/areaTotal > faceAreaIntersect::tolerance())
+                {
+                    // TODO - throwing area away - re-use in next iteration?
+
+                    seedFaces[faceS] = faceT;
+
+                    if (!valuesSet)
+                    {
+                        srcFaceI = faceS;
+                        tgtFaceI = faceT;
+                        valuesSet = true;
+                    }
+                }
+            }
+        }
+    }
+
+    // set next src and tgt faces if not set above
+    if (valuesSet)
+    {
+        return;
+    }
+    else
+    {
+        // try to use existing seed
+        bool foundNextSeed = false;
+        for (label faceI = startSeedI; faceI < mapFlag.size(); faceI++)
+        {
+            if (mapFlag[faceI])
+            {
+                if (!foundNextSeed)
+                {
+                    startSeedI = faceI;
+                    foundNextSeed = true;
+                }
+
+                if (seedFaces[faceI] != -1)
+                {
+                    srcFaceI = faceI;
+                    tgtFaceI = seedFaces[faceI];
+
+                    return;
+                }
+            }
+        }
+
+        // perform new search to find match
+        if (debug)
+        {
+            Pout<< "Advancing front stalled: searching for new "
+                << "target face" << endl;
+        }
+
+        foundNextSeed = false;
+        for (label faceI = startSeedI; faceI < mapFlag.size(); faceI++)
+        {
+            if (mapFlag[faceI])
+            {
+                if (!foundNextSeed)
+                {
+                    startSeedI = faceI + 1;
+                    foundNextSeed = true;
+                }
+
+                srcFaceI = faceI;
+                tgtFaceI = this->findTargetFace(srcFaceI);
+
+                if (tgtFaceI >= 0)
+                {
+                    return;
+                }
+            }
+        }
+
+        FatalErrorIn
+        (
+            "void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::"
+            "setNextFaces"
+            "("
+                "label&, "
+                "label&, "
+                "label&, "
+                "const boolList&, "
+                "labelList&, "
+                "const DynamicList<label>&"
+            ") const"
+        )  << "Unable to set source and target faces" << abort(FatalError);
+    }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+Foam::scalar Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::interArea
+(
+    const label srcFaceI,
+    const label tgtFaceI
+) const
+{
+    const pointField& srcPoints = this->srcPatch_.points();
+    const pointField& tgtPoints = this->tgtPatch_.points();
+
+    // references to candidate faces
+    const face& src = this->srcPatch_[srcFaceI];
+    const face& tgt = this->tgtPatch_[tgtFaceI];
+
+    // quick reject if either face has zero area
+    // Note: do not used stored face areas for target patch
+    if
+    (
+        (this->srcMagSf_[srcFaceI] < ROOTVSMALL)
+     || (tgt.mag(tgtPoints) < ROOTVSMALL)
+    )
+    {
+        return 0.0;
+    }
+
+    // create intersection object
+    faceAreaIntersect inter(srcPoints, tgtPoints, this->reverseTarget_);
+
+    // crude resultant norm
+    vector n(-src.normal(srcPoints));
+    if (this->reverseTarget_)
+    {
+        n -= tgt.normal(tgtPoints);
+    }
+    else
+    {
+        n += tgt.normal(tgtPoints);
+    }
+    n *= 0.5;
+
+    scalar area = 0;
+    if (mag(n) > ROOTVSMALL)
+    {
+        area = inter.calc(src, tgt, n, this->triMode_);
+    }
+    else
+    {
+        WarningIn
+        (
+            "void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::"
+            "interArea"
+            "("
+                "const label, "
+                "const label"
+            ") const"
+        )   << "Invalid normal for source face " << srcFaceI
+            << " points " << UIndirectList<point>(srcPoints, src)
+            << " target face " << tgtFaceI
+            << " points " << UIndirectList<point>(tgtPoints, tgt)
+            << endl;
+    }
+
+
+    if ((debug > 1) && (area > 0))
+    {
+        this->writeIntersectionOBJ(area, src, tgt, srcPoints, tgtPoints);
+    }
+
+    return area;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::
+restartUncoveredSourceFace
+(
+    List<DynamicList<label> >& srcAddr,
+    List<DynamicList<scalar> >& srcWght,
+    List<DynamicList<label> >& tgtAddr,
+    List<DynamicList<scalar> >& tgtWght
+)
+{
+    // Collect all src faces with a low weight
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    labelHashSet lowWeightFaces(100);
+    forAll(srcWght, srcFaceI)
+    {
+        scalar s = sum(srcWght[srcFaceI]);
+        scalar t = s/this->srcMagSf_[srcFaceI];
+
+        if (t < 0.5)
+        {
+            lowWeightFaces.insert(srcFaceI);
+        }
+    }
+
+    if (debug)
+    {
+        Pout<< "faceAreaWeightAMI: restarting search on "
+            << lowWeightFaces.size() << " faces since sum of weights < 0.5"
+            << endl;
+    }
+
+    if (lowWeightFaces.size() > 0)
+    {
+        // Erase all the lowWeight source faces from the target
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+        DynamicList<label> okSrcFaces(10);
+        DynamicList<scalar> okSrcWeights(10);
+        forAll(tgtAddr, tgtFaceI)
+        {
+            okSrcFaces.clear();
+            okSrcWeights.clear();
+            DynamicList<label>& srcFaces = tgtAddr[tgtFaceI];
+            DynamicList<scalar>& srcWeights = tgtWght[tgtFaceI];
+            forAll(srcFaces, i)
+            {
+                if (!lowWeightFaces.found(srcFaces[i]))
+                {
+                    okSrcFaces.append(srcFaces[i]);
+                    okSrcWeights.append(srcWeights[i]);
+                }
+            }
+            if (okSrcFaces.size() < srcFaces.size())
+            {
+                srcFaces.transfer(okSrcFaces);
+                srcWeights.transfer(okSrcWeights);
+            }
+        }
+
+
+
+        // Restart search from best hit
+        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+        // list of tgt face neighbour faces
+        DynamicList<label> nbrFaces(10);
+
+        // list of faces currently visited for srcFaceI to avoid multiple hits
+        DynamicList<label> visitedFaces(10);
+
+        forAllConstIter(labelHashSet, lowWeightFaces, iter)
+        {
+            label srcFaceI = iter.key();
+            label tgtFaceI = this->findTargetFace(srcFaceI);
+            if (tgtFaceI != -1)
+            {
+                //bool faceProcessed =
+                processSourceFace
+                (
+                    srcFaceI,
+                    tgtFaceI,
+
+                    nbrFaces,
+                    visitedFaces,
+
+                    srcAddr,
+                    srcWght,
+                    tgtAddr,
+                    tgtWght
+                );
+                // ? Check faceProcessed to see if restarting has worked.
+            }
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::faceAreaWeightAMI
+(
+    const SourcePatch& srcPatch,
+    const TargetPatch& tgtPatch,
+    const scalarField& srcMagSf,
+    const scalarField& tgtMagSf,
+    const faceAreaIntersect::triangulationMode& triMode,
+    const bool reverseTarget
+)
+:
+    AMIMethod<SourcePatch, TargetPatch>
+    (
+        srcPatch,
+        tgtPatch,
+        srcMagSf,
+        tgtMagSf,
+        triMode,
+        reverseTarget
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::~faceAreaWeightAMI()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::faceAreaWeightAMI<SourcePatch, TargetPatch>::calculate
+(
+    labelListList& srcAddress,
+    scalarListList& srcWeights,
+    labelListList& tgtAddress,
+    scalarListList& tgtWeights,
+    label srcFaceI,
+    label tgtFaceI
+)
+{
+    bool ok =
+        this->initialise
+        (
+            srcAddress,
+            srcWeights,
+            tgtAddress,
+            tgtWeights,
+            srcFaceI,
+            tgtFaceI
+        );
+
+    if (!ok)
+    {
+        return;
+    }
+
+    // temporary storage for addressing and weights
+    List<DynamicList<label> > srcAddr(this->srcPatch_.size());
+    List<DynamicList<scalar> > srcWght(srcAddr.size());
+    List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
+    List<DynamicList<scalar> > tgtWght(tgtAddr.size());
+
+
+    // construct weights and addressing
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    label nFacesRemaining = srcAddr.size();
+
+    // list of tgt face neighbour faces
+    DynamicList<label> nbrFaces(10);
+
+    // list of faces currently visited for srcFaceI to avoid multiple hits
+    DynamicList<label> visitedFaces(10);
+
+    // list to keep track of tgt faces used to seed src faces
+    labelList seedFaces(nFacesRemaining, -1);
+    seedFaces[srcFaceI] = tgtFaceI;
+
+    // list to keep track of whether src face can be mapped
+    boolList mapFlag(nFacesRemaining, true);
+
+    // reset starting seed
+    label startSeedI = 0;
+
+    DynamicList<label> nonOverlapFaces;
+    do
+    {
+        // Do advancing front starting from srcFaceI,tgtFaceI
+        bool faceProcessed = processSourceFace
+        (
+            srcFaceI,
+            tgtFaceI,
+
+            nbrFaces,
+            visitedFaces,
+
+            srcAddr,
+            srcWght,
+            tgtAddr,
+            tgtWght
+        );
+
+        mapFlag[srcFaceI] = false;
+
+        nFacesRemaining--;
+
+        if (!faceProcessed)
+        {
+            nonOverlapFaces.append(srcFaceI);
+        }
+
+        // choose new src face from current src face neighbour
+        if (nFacesRemaining > 0)
+        {
+            setNextFaces
+            (
+                startSeedI,
+                srcFaceI,
+                tgtFaceI,
+                mapFlag,
+                seedFaces,
+                visitedFaces
+            );
+        }
+    } while (nFacesRemaining > 0);
+
+    if (nonOverlapFaces.size() != 0)
+    {
+        Pout<< "    AMI: " << nonOverlapFaces.size()
+            << " non-overlap faces identified"
+            << endl;
+
+        this->srcNonOverlap_.transfer(nonOverlapFaces);
+    }
+
+
+    // Check for badly covered faces
+    if (debug)
+    {
+        restartUncoveredSourceFace
+        (
+            srcAddr,
+            srcWght,
+            tgtAddr,
+            tgtWght
+        );
+    }
+
+
+    // transfer data to persistent storage
+    forAll(srcAddr, i)
+    {
+        srcAddress[i].transfer(srcAddr[i]);
+        srcWeights[i].transfer(srcWght[i]);
+    }
+    forAll(tgtAddr, i)
+    {
+        tgtAddress[i].transfer(tgtAddr[i]);
+        tgtWeights[i].transfer(tgtWght[i]);
+    }
+}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/faceAreaWeightAMI/faceAreaWeightAMI.H

@@ -0,0 +1,166 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::faceAreaWeightAMI
+
+Description
+    Face area weighted Arbitrary Mesh Interface (AMI) method
+
+SourceFiles
+    faceAreaWeightAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef faceAreaWeightAMI_H
+#define faceAreaWeightAMI_H
+
+#include "AMIMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class faceAreaWeightAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class faceAreaWeightAMI
+:
+    public AMIMethod<SourcePatch, TargetPatch>
+{
+
+private:
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        faceAreaWeightAMI(const faceAreaWeightAMI&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const faceAreaWeightAMI&);
+
+        // Marching front
+
+            //- Determine overlap contributions for source face srcFaceI
+            bool processSourceFace
+            (
+                const label srcFaceI,
+                const label tgtStartFaceI,
+                DynamicList<label>& nbrFaces,
+                DynamicList<label>& visitedFaces,
+                List<DynamicList<label> >& srcAddr,
+                List<DynamicList<scalar> >& srcWght,
+                List<DynamicList<label> >& tgtAddr,
+                List<DynamicList<scalar> >& tgtWght
+            );
+
+            //- Attempt to re-evaluate source faces that have not been included
+            void restartUncoveredSourceFace
+            (
+                List<DynamicList<label> >& srcAddr,
+                List<DynamicList<scalar> >& srcWght,
+                List<DynamicList<label> >& tgtAddr,
+                List<DynamicList<scalar> >& tgtWght
+            );
+
+            //- Set the source and target seed faces
+            void setNextFaces
+            (
+                label& startSeedI,
+                label& srcFaceI,
+                label& tgtFaceI,
+                const boolList& mapFlag,
+                labelList& seedFaces,
+                const DynamicList<label>& visitedFaces
+            ) const;
+
+
+        // Evaluation
+
+            //- Area of intersection between source and target faces
+            scalar interArea
+            (
+                const label srcFaceI,
+                const label tgtFaceI
+            ) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("faceAreaWeightAMI");
+
+
+    // Constructors
+
+        //- Construct from components
+        faceAreaWeightAMI
+        (
+            const SourcePatch& srcPatch,
+            const TargetPatch& tgtPatch,
+            const scalarField& srcMagSf,
+            const scalarField& tgtMagSf,
+            const faceAreaIntersect::triangulationMode& triMode,
+            const bool reverseTarget = false
+        );
+
+
+    //- Destructor
+    virtual ~faceAreaWeightAMI();
+
+
+    // Member Functions
+
+        // Manipulation
+
+            //- Update addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcAddress,
+                scalarListList& srcWeights,
+                labelListList& tgtAddress,
+                scalarListList& tgtWeights,
+                label srcFaceI = -1,
+                label tgtFaceI = -1
+            );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "faceAreaWeightAMI.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.C

@@ -0,0 +1,343 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "mapNearestAMI.H"
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::mapNearestAMI<SourcePatch, TargetPatch>::findNearestFace
+(
+    const SourcePatch& srcPatch,
+    const TargetPatch& tgtPatch,
+    const label& srcFaceI,
+    label& tgtFaceI
+) const
+{
+    const vectorField& srcCf = srcPatch.faceCentres();
+    const vectorField& tgtCf = tgtPatch.faceCentres();
+
+    const vector srcP = srcCf[srcFaceI];
+
+    DynamicList<label> tgtFaces(10);
+    tgtFaces.append(tgtFaceI);
+
+    DynamicList<label> visitedFaces(10);
+
+    scalar d = GREAT;
+
+    do
+    {
+        label tgtI = tgtFaces.remove();
+        visitedFaces.append(tgtI);
+
+        scalar dTest = magSqr(tgtCf[tgtI] - srcP);
+        if (dTest < d)
+        {
+            tgtFaceI = tgtI;
+            d = dTest;
+
+            this->appendNbrFaces
+            (
+                tgtFaceI,
+                tgtPatch,
+                visitedFaces,
+                tgtFaces
+            );
+        }
+
+    } while (tgtFaces.size() > 0);
+}
+
+
+template<class SourcePatch, class TargetPatch>
+void Foam::mapNearestAMI<SourcePatch, TargetPatch>::setNextNearestFaces
+(
+    boolList& mapFlag,
+    label& startSeedI,
+    label& srcFaceI,
+    label& tgtFaceI
+) const
+{
+    const labelList& srcNbr = this->srcPatch_.faceFaces()[srcFaceI];
+
+    srcFaceI = -1;
+
+    forAll(srcNbr, i)
+    {
+        label faceI = srcNbr[i];
+        if (mapFlag[faceI])
+        {
+            srcFaceI = faceI;
+            startSeedI = faceI + 1;
+
+            return;
+        }
+    }
+
+    forAll(mapFlag, srcFaceI)
+    {
+        if (!mapFlag[srcFaceI])
+        {
+            tgtFaceI = this->findTargetFace(srcFaceI);
+            if (tgtFaceI == -1)
+            {
+                const vectorField& srcCf = this->srcPatch_.faceCentres();
+
+                FatalErrorIn
+                (
+                    "void Foam::mapNearestAMI<SourcePatch, TargetPatch>::"
+                    "setNextNearestFaces"
+                    "("
+                        "boolList&, "
+                        "label&, "
+                        "label&, "
+                        "label&"
+                    ") const"
+                )
+                    << "Unable to find target face for source face "
+                    << srcFaceI << " with face centre " << srcCf[srcFaceI]
+                    << abort(FatalError);
+            }
+
+            break;
+        }
+    }
+}
+
+
+template<class SourcePatch, class TargetPatch>
+Foam::label Foam::mapNearestAMI<SourcePatch, TargetPatch>::findMappedSrcFace
+(
+    const label tgtFaceI,
+    const List<DynamicList<label> >& tgtToSrc
+) const
+{
+    DynamicList<label> testFaces(10);
+    DynamicList<label> visitedFaces(10);
+
+    testFaces.append(tgtFaceI);
+
+    do
+    {
+        // search target tgtFaceI neighbours for match with source face
+        label tgtI = testFaces.remove();
+
+        if (findIndex(visitedFaces, tgtI) == -1)
+        {
+            visitedFaces.append(tgtI);
+
+            if (tgtToSrc[tgtI].size())
+            {
+                return tgtToSrc[tgtI][0];
+            }
+            else
+            {
+                const labelList& nbrFaces = this->tgtPatch_.faceFaces()[tgtI];
+
+                forAll(nbrFaces, i)
+                {
+                    if (findIndex(visitedFaces, nbrFaces[i]) == -1)
+                    {
+                        testFaces.append(nbrFaces[i]);
+                    }
+                }
+            }
+        }
+    } while (testFaces.size());
+
+    // did not find any match - should not be possible to get here!
+    return -1;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::mapNearestAMI<SourcePatch, TargetPatch>::mapNearestAMI
+(
+    const SourcePatch& srcPatch,
+    const TargetPatch& tgtPatch,
+    const scalarField& srcMagSf,
+    const scalarField& tgtMagSf,
+    const faceAreaIntersect::triangulationMode& triMode,
+    const bool reverseTarget
+)
+:
+    AMIMethod<SourcePatch, TargetPatch>
+    (
+        srcPatch,
+        tgtPatch,
+        srcMagSf,
+        tgtMagSf,
+        triMode,
+        reverseTarget
+    )
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+Foam::mapNearestAMI<SourcePatch, TargetPatch>::~mapNearestAMI()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+template<class SourcePatch, class TargetPatch>
+void Foam::mapNearestAMI<SourcePatch, TargetPatch>::calculate
+(
+    labelListList& srcAddress,
+    scalarListList& srcWeights,
+    labelListList& tgtAddress,
+    scalarListList& tgtWeights,
+    label srcFaceI,
+    label tgtFaceI
+)
+{
+    bool ok =
+        this->initialise
+        (
+            srcAddress,
+            srcWeights,
+            tgtAddress,
+            tgtWeights,
+            srcFaceI,
+            tgtFaceI
+        );
+
+    if (!ok)
+    {
+        return;
+    }
+
+
+    // temporary storage for addressing and weights
+    List<DynamicList<label> > srcAddr(this->srcPatch_.size());
+    List<DynamicList<label> > tgtAddr(this->tgtPatch_.size());
+
+
+    // construct weights and addressing
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // list to keep track of whether src face can be mapped
+    boolList mapFlag(srcAddr.size(), true);
+
+    // reset starting seed
+    label startSeedI = 0;
+
+    DynamicList<label> nonOverlapFaces;
+    do
+    {
+        findNearestFace(this->srcPatch_, this->tgtPatch_, srcFaceI, tgtFaceI);
+
+        srcAddr[srcFaceI].append(tgtFaceI);
+        tgtAddr[tgtFaceI].append(srcFaceI);
+
+        mapFlag[srcFaceI] = false;
+
+        // Do advancing front starting from srcFaceI, tgtFaceI
+        setNextNearestFaces
+        (
+            mapFlag,
+            startSeedI,
+            srcFaceI,
+            tgtFaceI
+        );
+    } while (srcFaceI >= 0);
+
+
+    // for the case of multiple source faces per target face, select the
+    // nearest source face only and discard the others
+    const vectorField& srcCf = this->srcPatch_.faceCentres();
+    const vectorField& tgtCf = this->tgtPatch_.faceCentres();
+
+    forAll(tgtAddr, targetFaceI)
+    {
+        if (tgtAddr[targetFaceI].size() > 1)
+        {
+            const vector& tgtC = tgtCf[tgtFaceI];
+
+            DynamicList<label>& srcFaces = tgtAddr[targetFaceI];
+
+            label srcFaceI = srcFaces[0];
+            scalar d = magSqr(tgtC - srcCf[srcFaceI]);
+
+            for (label i = 1; i < srcFaces.size(); i++)
+            {
+                label srcI = srcFaces[i];
+                scalar dNew = magSqr(tgtC - srcCf[srcI]);
+                if (dNew < d)
+                {
+                    d = dNew;
+                    srcFaceI = srcI;
+                }
+            }
+
+            srcFaces.clear();
+            srcFaces.append(srcFaceI);
+        }
+    }
+
+    // If there are more target faces than source faces, some target faces
+    // might not yet be mapped
+    forAll(tgtAddr, tgtFaceI)
+    {
+        if (tgtAddr[tgtFaceI].empty())
+        {
+            label srcFaceI = findMappedSrcFace(tgtFaceI, tgtAddr);
+
+            // note - reversed search from src->tgt to tgt->src
+            findNearestFace
+            (
+                this->tgtPatch_,
+                this->srcPatch_,
+                tgtFaceI,
+                srcFaceI
+            );
+
+            tgtAddr[tgtFaceI].append(srcFaceI);
+        }
+    }
+
+
+    // transfer data to persistent storage
+    forAll(srcAddr, i)
+    {
+        scalar magSf = this->srcMagSf_[i];
+        srcAddress[i].transfer(srcAddr[i]);
+        srcWeights[i] = scalarList(1, magSf);
+    }
+    forAll(tgtAddr, i)
+    {
+        scalar magSf = this->tgtMagSf_[i];
+        tgtAddress[i].transfer(tgtAddr[i]);
+        tgtWeights[i] = scalarList(1, magSf);
+    }
+}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIMethod/mapNearestAMI/mapNearestAMI.H

@@ -0,0 +1,158 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::mapNearestAMI
+
+Description
+    Nearest-mapping Arbitrary Mesh Interface (AMI) method
+
+SourceFiles
+    mapNearestAMI.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef mapNearestAMI_H
+#define mapNearestAMI_H
+
+#include "AMIMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class mapNearestAMI Declaration
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+class mapNearestAMI
+:
+    public AMIMethod<SourcePatch, TargetPatch>
+{
+
+private:
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        mapNearestAMI(const mapNearestAMI&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const mapNearestAMI&);
+
+        // Marching front
+
+            //- Find nearest target face for source face srcFaceI
+            void findNearestFace
+            (
+                const SourcePatch& srcPatch,
+                const TargetPatch& tgtPatch,
+                const label& srcFaceI,
+                label& tgtFaceI
+            ) const;
+
+            //- Determine next source-target face pair
+            void setNextNearestFaces
+            (
+                boolList& mapFlag,
+                label& startSeedI,
+                label& srcFaceI,
+                label& tgtFaceI
+            ) const;
+
+            //- Find mapped source face
+            label findMappedSrcFace
+            (
+                const label tgtFaceI,
+                const List<DynamicList<label> >& tgtToSrc
+            ) const;
+
+
+        // Evaluation
+
+            //- Area of intersection between source and target faces
+            scalar interArea
+            (
+                const label srcFaceI,
+                const label tgtFaceI
+            ) const;
+
+
+public:
+
+    //- Runtime type information
+    TypeName("mapNearestAMI");
+
+
+    // Constructors
+
+        //- Construct from components
+        mapNearestAMI
+        (
+            const SourcePatch& srcPatch,
+            const TargetPatch& tgtPatch,
+            const scalarField& srcMagSf,
+            const scalarField& tgtMagSf,
+            const faceAreaIntersect::triangulationMode& triMode,
+            const bool reverseTarget = false
+        );
+
+
+    //- Destructor
+    virtual ~mapNearestAMI();
+
+
+    // Member Functions
+
+        // Manipulation
+
+            //- Update addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcAddress,
+                scalarListList& srcWeights,
+                labelListList& tgtAddress,
+                scalarListList& tgtWeights,
+                label srcFaceI = -1,
+                label tgtFaceI = -1
+            );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+#   include "mapNearestAMI.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/AMIInterpolation/AMIPatchToPatchInterpolation.C

@@ -0,0 +1,44 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "AMIPatchToPatchInterpolation.H"
+#include "AMIMethod.H"
+#include "directAMI.H"
+#include "mapNearestAMI.H"
+#include "faceAreaWeightAMI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    makeAMIMethod(AMIPatchToPatchInterpolation);
+
+    makeAMIMethodType(AMIPatchToPatchInterpolation, directAMI);
+    makeAMIMethodType(AMIPatchToPatchInterpolation, mapNearestAMI);
+    makeAMIMethodType(AMIPatchToPatchInterpolation, faceAreaWeightAMI);
+}
+
+
+// ************************************************************************* //

src/meshTools/AMIInterpolation/patches/cyclic/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C

@@ -283,6 +283,7 @@ void Foam::cyclicAMIPolyPatch::resetAMI() const
                 nbrPatch0,
                 surfPtr(),
                 faceAreaIntersect::tmMesh,
+                AMIPatchToPatchInterpolation::imFaceAreaWeight,
                 AMIReverse_
             )
         );

src/meshTools/Make/files

@@ -171,6 +171,7 @@ twoDPointCorrector/twoDPointCorrector.C
 
 AMI=AMIInterpolation
 $(AMI)/AMIInterpolation/AMIInterpolationName.C
+$(AMI)/AMIInterpolation/AMIPatchToPatchInterpolation.C
 $(AMI)/faceAreaIntersect/faceAreaIntersect.C
 $(AMI)/GAMG/interfaces/cyclicAMIGAMGInterface/cyclicAMIGAMGInterface.C
 $(AMI)/GAMG/interfaceFields/cyclicAMIGAMGInterfaceField/cyclicAMIGAMGInterfaceField.C

src/meshTools/mappedPatches/mappedPolyPatch/mappedPatchBase.C

@@ -849,6 +849,7 @@ void Foam::mappedPatchBase::calcAMI() const
             samplePolyPatch(), // nbrPatch0,
             surfPtr(),
             faceAreaIntersect::tmMesh,
+            AMIPatchToPatchInterpolation::imFaceAreaWeight,
             AMIReverse_
         )
     );

src/meshTools/regionCoupled/patches/regionCoupledPolyPatch/regionCoupledBase.C

@@ -91,6 +91,7 @@ void Foam::regionCoupledBase::resetAMI() const
                 nbrPatch0,
                 surfPtr(),
                 faceAreaIntersect::tmMesh,
+                AMIPatchToPatchInterpolation::imFaceAreaWeight,
                 AMIReverse_
             )
         );

src/meshTools/sets/faceZoneSources/searchableSurfaceToFaceZone/searchableSurfaceToFaceZone.C

@@ -28,6 +28,7 @@ License
 #include "faceZoneSet.H"
 #include "searchableSurface.H"
 #include "syncTools.H"
+#include "Time.H"
 
 #include "addToRunTimeSelectionTable.H"
 
@@ -69,7 +70,15 @@ Foam::searchableSurfaceToFaceZone::searchableSurfaceToFaceZone
         searchableSurface::New
         (
             word(dict.lookup("surface")),
+            IOobject
+            (
+                dict.lookupOrDefault("name", mesh.objectRegistry::db().name()),
+                mesh.time().constant(),
+                "triSurface",
                 mesh.objectRegistry::db(),
+                IOobject::MUST_READ,
+                IOobject::NO_WRITE
+            ),
             dict
         )
     )

src/regionModels/regionModel/regionModel/regionModel.C

@@ -221,6 +221,7 @@ Foam::regionModels::regionModel::interRegionAMI
                     p,
                     nbrP,
                     faceAreaIntersect::tmMesh,
+                    AMIPatchToPatchInterpolation::imFaceAreaWeight,
                     flip
                 )
             );
@@ -261,6 +262,7 @@ Foam::regionModels::regionModel::interRegionAMI
                 p,
                 nbrP,
                 faceAreaIntersect::tmMesh,
+                AMIPatchToPatchInterpolation::imFaceAreaWeight,
                 flip
             )
         );

src/sampling/Make/files

@@ -60,10 +60,14 @@ $(meshToMesh)/calculateMeshToMeshAddressing.C
 $(meshToMesh)/calculateMeshToMeshWeights.C
 
 meshToMeshNew = meshToMeshInterpolation/meshToMeshNew
-$(meshToMeshNew)/calcDirect.C
-$(meshToMeshNew)/calcMapNearest.C
-$(meshToMeshNew)/calcCellVolumeWeight.C
 $(meshToMeshNew)/meshToMeshNew.C
 $(meshToMeshNew)/meshToMeshNewParallelOps.C
+meshToMeshNewMethods = meshToMeshInterpolation/meshToMeshNew/calcMethod
+$(meshToMeshNewMethods)/meshToMeshMethod/meshToMeshMethod.C
+$(meshToMeshNewMethods)/meshToMeshMethod/meshToMeshMethodNew.C
+$(meshToMeshNewMethods)/cellVolumeWeight/cellVolumeWeightMethod.C
+$(meshToMeshNewMethods)/direct/directMethod.C
+$(meshToMeshNewMethods)/mapNearest/mapNearestMethod.C
+
 
 LIB = $(FOAM_LIBBIN)/libsampling

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcDirect.C

@@ -1,159 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "meshToMeshNew.H"
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::meshToMeshNew::calcDirect
-(
-    const polyMesh& src,
-    const polyMesh& tgt,
-    const label srcSeedI,
-    const label tgtSeedI
-)
-{
-    // store a list of src cells already mapped
-    boolList srcSeedFlag(src.nCells(), true);
-    labelList srcTgtSeed(src.nCells(), -1);
-
-    List<DynamicList<label> > srcToTgt(src.nCells());
-    List<DynamicList<label> > tgtToSrc(tgt.nCells());
-
-    DynamicList<label> srcSeeds;
-
-    const scalarField& srcVc = src.cellVolumes();
-    const scalarField& tgtVc = tgt.cellVolumes();
-
-    label srcCellI = srcSeedI;
-    label tgtCellI = tgtSeedI;
-
-    do
-    {
-        // store src/tgt cell pair
-        srcToTgt[srcCellI].append(tgtCellI);
-        tgtToSrc[tgtCellI].append(srcCellI);
-
-        // mark source cell srcSeedI as matched
-        srcSeedFlag[srcCellI] = false;
-
-        // accumulate intersection volume
-        V_ += srcVc[srcCellI];
-
-        // find new source seed cell
-        appendToDirectSeeds
-        (
-            src,
-            tgt,
-            srcSeedFlag,
-            srcTgtSeed,
-            srcSeeds,
-            srcCellI,
-            tgtCellI
-        );
-    }
-    while (srcCellI >= 0);
-
-    // transfer addressing into persistent storage
-    forAll(srcToTgtCellAddr_, i)
-    {
-        scalar v = srcVc[i];
-        srcToTgtCellAddr_[i].transfer(srcToTgt[i]);
-        srcToTgtCellWght_[i] = scalarList(srcToTgtCellAddr_[i].size(), v);
-    }
-
-    forAll(tgtToSrcCellAddr_, i)
-    {
-        scalar v = tgtVc[i];
-        tgtToSrcCellAddr_[i].transfer(tgtToSrc[i]);
-        tgtToSrcCellWght_[i] = scalarList(tgtToSrcCellAddr_[i].size(), v);
-    }
-}
-
-
-void Foam::meshToMeshNew::appendToDirectSeeds
-(
-    const polyMesh& src,
-    const polyMesh& tgt,
-    boolList& mapFlag,
-    labelList& srcTgtSeed,
-    DynamicList<label>& srcSeeds,
-    label& srcSeedI,
-    label& tgtSeedI
-) const
-{
-    const labelList& srcNbr = src.cellCells()[srcSeedI];
-    const labelList& tgtNbr = tgt.cellCells()[tgtSeedI];
-
-    const vectorField& srcCentre = src.cellCentres();
-
-    forAll(srcNbr, i)
-    {
-        label srcI = srcNbr[i];
-
-        if (mapFlag[srcI] && (srcTgtSeed[srcI] == -1))
-        {
-            // source cell srcI not yet mapped
-
-            // identfy if target cell exists for source cell srcI
-            bool found = false;
-            forAll(tgtNbr, j)
-            {
-                label tgtI = tgtNbr[j];
-
-                if (tgt.pointInCell(srcCentre[srcI], tgtI))
-                {
-                    // new match - append to lists
-                    found = true;
-
-                    srcTgtSeed[srcI] = tgtI;
-                    srcSeeds.append(srcI);
-
-                    break;
-                }
-            }
-
-            if (!found)
-            {
-                // no match available for source cell srcI
-                mapFlag[srcI] = false;
-            }
-        }
-    }
-
-    if (srcSeeds.size())
-    {
-        srcSeedI = srcSeeds.remove();
-        tgtSeedI = srcTgtSeed[srcSeedI];
-    }
-    else
-    {
-        srcSeedI = -1;
-        tgtSeedI = -1;
-    }
-}
-
-
-// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMapNearest.C

@@ -1,265 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software: you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by
-    the Free Software Foundation, either version 3 of the License, or
-    (at your option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
-
-\*---------------------------------------------------------------------------*/
-
-#include "meshToMeshNew.H"
-#include "ListOps.H"
-
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
-
-void Foam::meshToMeshNew::calcMapNearest
-(
-    const polyMesh& src,
-    const polyMesh& tgt,
-    const label srcSeedI,
-    const label tgtSeedI,
-    const labelList& srcCellIDs,
-    boolList& mapFlag,
-    label& startSeedI
-)
-{
-    List<DynamicList<label> > srcToTgt(src.nCells());
-    List<DynamicList<label> > tgtToSrc(tgt.nCells());
-
-    const scalarField& srcVc = src.cellVolumes();
-    const scalarField& tgtVc = tgt.cellVolumes();
-
-    label srcCellI = srcSeedI;
-    label tgtCellI = tgtSeedI;
-
-    do
-    {
-        // find nearest tgt cell
-        findNearestCell(src, tgt, srcCellI, tgtCellI);
-
-        // store src/tgt cell pair
-        srcToTgt[srcCellI].append(tgtCellI);
-        tgtToSrc[tgtCellI].append(srcCellI);
-
-        // mark source cell srcCellI and tgtCellI as matched
-        mapFlag[srcCellI] = false;
-
-        // accumulate intersection volume
-        V_ += srcVc[srcCellI];
-
-        // find new source cell
-        setNextNearestCells
-        (
-            startSeedI,
-            srcCellI,
-            tgtCellI,
-            mapFlag,
-            src,
-            tgt,
-            srcCellIDs
-        );
-    }
-    while (srcCellI >= 0);
-
-
-    // for the case of multiple source cells per target cell, select the
-    // nearest source cell only and discard the others
-    const vectorField& srcCc = src.cellCentres();
-    const vectorField& tgtCc = tgt.cellCentres();
-
-    forAll(tgtToSrc, targetCellI)
-    {
-        if (tgtToSrc[targetCellI].size() > 1)
-        {
-            const vector& tgtC = tgtCc[tgtCellI];
-
-            DynamicList<label>& srcCells = tgtToSrc[targetCellI];
-
-            label srcCellI = srcCells[0];
-            scalar d = magSqr(tgtC - srcCc[srcCellI]);
-
-            for (label i = 1; i < srcCells.size(); i++)
-            {
-                label srcI = srcCells[i];
-                scalar dNew = magSqr(tgtC - srcCc[srcI]);
-                if (dNew < d)
-                {
-                    d = dNew;
-                    srcCellI = srcI;
-                }
-            }
-
-            srcCells.clear();
-            srcCells.append(srcCellI);
-        }
-    }
-
-    // If there are more target cells than source cells, some target cells
-    // might not yet be mapped
-    forAll(tgtToSrc, tgtCellI)
-    {
-        if (tgtToSrc[tgtCellI].empty())
-        {
-            label srcCellI = findMappedSrcCell(tgt, tgtCellI, tgtToSrc);
-
-            findNearestCell(tgt, src, tgtCellI, srcCellI);
-
-            tgtToSrc[tgtCellI].append(srcCellI);
-        }
-    }
-
-    // transfer addressing into persistent storage
-    forAll(srcToTgtCellAddr_, i)
-    {
-        scalar v = srcVc[i];
-        srcToTgtCellWght_[i] = scalarList(srcToTgt[i].size(), v);
-        srcToTgtCellAddr_[i].transfer(srcToTgt[i]);
-    }
-
-    forAll(tgtToSrcCellAddr_, i)
-    {
-        scalar v = tgtVc[i];
-        tgtToSrcCellWght_[i] = scalarList(tgtToSrc[i].size(), v);
-        tgtToSrcCellAddr_[i].transfer(tgtToSrc[i]);
-    }
-}
-
-
-void Foam::meshToMeshNew::findNearestCell
-(
-    const polyMesh& src,
-    const polyMesh& tgt,
-    const label srcCellI,
-    label& tgtCellI
-)
-{
-    const vectorField& srcC = src.cellCentres();
-    const vectorField& tgtC = tgt.cellCentres();
-
-    const vector& srcP = srcC[srcCellI];
-
-    DynamicList<label> tgtCells(10);
-    tgtCells.append(tgtCellI);
-
-    DynamicList<label> visitedCells(10);
-
-    scalar d = GREAT;
-
-    do
-    {
-        label tgtI = tgtCells.remove();
-        visitedCells.append(tgtI);
-
-        scalar dTest = magSqr(tgtC[tgtI] - srcP);
-        if (dTest < d)
-        {
-            tgtCellI = tgtI;
-            d = dTest;
-            appendNbrCells(tgtCellI, tgt, visitedCells, tgtCells);
-        }
-
-    } while (tgtCells.size() > 0);
-}
-
-
-void Foam::meshToMeshNew::setNextNearestCells
-(
-    label& startSeedI,
-    label& srcCellI,
-    label& tgtCellI,
-    boolList& mapFlag,
-    const polyMesh& src,
-    const polyMesh& tgt,
-    const labelList& srcCellIDs
-)
-{
-    const labelList& srcNbr = src.cellCells()[srcCellI];
-
-    srcCellI = -1;
-    forAll(srcNbr, i)
-    {
-        label cellI = srcNbr[i];
-        if (mapFlag[cellI])
-        {
-            srcCellI = cellI;
-            startSeedI = cellI + 1;
-
-            return;
-        }
-    }
-
-    (void)findInitialSeeds
-    (
-        src,
-        tgt,
-        srcCellIDs,
-        mapFlag,
-        startSeedI,
-        srcCellI,
-        tgtCellI
-    );
-}
-
-
-Foam::label Foam::meshToMeshNew::findMappedSrcCell
-(
-    const polyMesh& tgt,
-    const label tgtCellI,
-    const List<DynamicList<label> >& tgtToSrc
-) const
-{
-    DynamicList<label> testCells(10);
-    DynamicList<label> visitedCells(10);
-
-    testCells.append(tgtCellI);
-
-    do
-    {
-        // search target tgtCellI neighbours for match with source cell
-        label tgtI = testCells.remove();
-
-        if (findIndex(visitedCells, tgtI) == -1)
-        {
-            visitedCells.append(tgtI);
-
-            if (tgtToSrc[tgtI].size())
-            {
-                return tgtToSrc[tgtI][0];
-            }
-            else
-            {
-                const labelList& nbrCells = tgt.cellCells()[tgtI];
-
-                forAll(nbrCells, i)
-                {
-                    if (findIndex(visitedCells, nbrCells[i]) == -1)
-                    {
-                        testCells.append(nbrCells[i]);
-                    }
-                }
-            }
-        }
-    } while (testCells.size());
-
-    // did not find any match - should not be possible to get here!
-    return -1;
-}
-
-
-// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.C

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2012-2013 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -23,15 +23,79 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "meshToMeshNew.H"
+#include "cellVolumeWeightMethod.H"
-#include "tetOverlapVolume.H"
+#include "indexedOctree.H"
+#include "treeDataCell.H"
+#include "addToRunTimeSelectionTable.H"
 
-// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
-void Foam::meshToMeshNew::calcCellVolumeWeight
+namespace Foam
+{
+    defineTypeNameAndDebug(cellVolumeWeightMethod, 0);
+    addToRunTimeSelectionTable
     (
-    const polyMesh& src,
+        meshToMeshMethod,
-    const polyMesh& tgt,
+        cellVolumeWeightMethod,
+        components
+    );
+}
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+bool Foam::cellVolumeWeightMethod::findInitialSeeds
+(
+    const labelList& srcCellIDs,
+    const boolList& mapFlag,
+    const label startSeedI,
+    label& srcSeedI,
+    label& tgtSeedI
+) const
+{
+    const cellList& srcCells = src_.cells();
+    const faceList& srcFaces = src_.faces();
+    const pointField& srcPts = src_.points();
+
+    for (label i = startSeedI; i < srcCellIDs.size(); i++)
+    {
+        label srcI = srcCellIDs[i];
+
+        if (mapFlag[srcI])
+        {
+            const pointField
+                pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
+
+            forAll(pts, ptI)
+            {
+                const point& pt = pts[ptI];
+                label tgtI = tgt_.cellTree().findInside(pt);
+
+                if (tgtI != -1 && intersect(srcI, tgtI))
+                {
+                    srcSeedI = srcI;
+                    tgtSeedI = tgtI;
+
+                    return true;
+                }
+            }
+        }
+    }
+
+    if (debug)
+    {
+        Pout<< "could not find starting seed" << endl;
+    }
+
+    return false;
+}
+
+
+void Foam::cellVolumeWeightMethod::calculateAddressing
+(
+    labelListList& srcToTgtCellAddr,
+    scalarListList& srcToTgtCellWght,
+    labelListList& tgtToSrcCellAddr,
+    scalarListList& tgtToSrcCellWght,
     const label srcSeedI,
     const label tgtSeedI,
     const labelList& srcCellIDs,
@@ -42,11 +106,11 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
     label srcCellI = srcSeedI;
     label tgtCellI = tgtSeedI;
 
-    List<DynamicList<label> > srcToTgtAddr(src.nCells());
+    List<DynamicList<label> > srcToTgtAddr(src_.nCells());
-    List<DynamicList<scalar> > srcToTgtWght(src.nCells());
+    List<DynamicList<scalar> > srcToTgtWght(src_.nCells());
 
-    List<DynamicList<label> > tgtToSrcAddr(tgt.nCells());
+    List<DynamicList<label> > tgtToSrcAddr(tgt_.nCells());
-    List<DynamicList<scalar> > tgtToSrcWght(tgt.nCells());
+    List<DynamicList<scalar> > tgtToSrcWght(tgt_.nCells());
 
     // list of tgt cell neighbour cells
     DynamicList<label> nbrTgtCells(10);
@@ -55,10 +119,10 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
     DynamicList<label> visitedTgtCells(10);
 
     // list to keep track of tgt cells used to seed src cells
-    labelList seedCells(src.nCells(), -1);
+    labelList seedCells(src_.nCells(), -1);
     seedCells[srcCellI] = tgtCellI;
 
-    const scalarField& srcVol = src.cellVolumes();
+    const scalarField& srcVol = src_.cellVolumes();
 
     do
     {
@@ -67,14 +131,14 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
 
         // append initial target cell and neighbours
         nbrTgtCells.append(tgtCellI);
-        appendNbrCells(tgtCellI, tgt, visitedTgtCells, nbrTgtCells);
+        appendNbrCells(tgtCellI, tgt_, visitedTgtCells, nbrTgtCells);
 
         do
         {
             tgtCellI = nbrTgtCells.remove();
             visitedTgtCells.append(tgtCellI);
 
-            scalar vol = interVol(src, tgt, srcCellI, tgtCellI);
+            scalar vol = interVol(srcCellI, tgtCellI);
 
             // accumulate addressing and weights for valid intersection
             if (vol/srcVol[srcCellI] > tolerance_)
@@ -86,7 +150,7 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
                 tgtToSrcAddr[tgtCellI].append(srcCellI);
                 tgtToSrcWght[tgtCellI].append(vol);
 
-                appendNbrCells(tgtCellI, tgt, visitedTgtCells, nbrTgtCells);
+                appendNbrCells(tgtCellI, tgt_, visitedTgtCells, nbrTgtCells);
 
                 // accumulate intersection volume
                 V_ += vol;
@@ -102,8 +166,6 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
             startSeedI,
             srcCellI,
             tgtCellI,
-            src,
-            tgt,
             srcCellIDs,
             mapFlag,
             visitedTgtCells,
@@ -113,34 +175,32 @@ void Foam::meshToMeshNew::calcCellVolumeWeight
     while (srcCellI != -1);
 
     // transfer addressing into persistent storage
-    forAll(srcToTgtCellAddr_, i)
+    forAll(srcToTgtCellAddr, i)
     {
-        srcToTgtCellAddr_[i].transfer(srcToTgtAddr[i]);
+        srcToTgtCellAddr[i].transfer(srcToTgtAddr[i]);
-        srcToTgtCellWght_[i].transfer(srcToTgtWght[i]);
+        srcToTgtCellWght[i].transfer(srcToTgtWght[i]);
     }
 
-    forAll(tgtToSrcCellAddr_, i)
+    forAll(tgtToSrcCellAddr, i)
     {
-        tgtToSrcCellAddr_[i].transfer(tgtToSrcAddr[i]);
+        tgtToSrcCellAddr[i].transfer(tgtToSrcAddr[i]);
-        tgtToSrcCellWght_[i].transfer(tgtToSrcWght[i]);
+        tgtToSrcCellWght[i].transfer(tgtToSrcWght[i]);
     }
 }
 
 
-void Foam::meshToMeshNew::setNextCells
+void Foam::cellVolumeWeightMethod::setNextCells
 (
     label& startSeedI,
     label& srcCellI,
     label& tgtCellI,
-    const polyMesh& src,
-    const polyMesh& tgt,
     const labelList& srcCellIDs,
     const boolList& mapFlag,
     const DynamicList<label>& visitedCells,
     labelList& seedCells
 ) const
 {
-    const labelList& srcNbrCells = src.cellCells()[srcCellI];
+    const labelList& srcNbrCells = src_.cellCells()[srcCellI];
 
     // set possible seeds for later use by querying all src cell neighbours
     // with all visited target cells
@@ -155,7 +215,7 @@ void Foam::meshToMeshNew::setNextCells
             {
                 label cellT = visitedCells[j];
 
-                if (intersect(src, tgt, cellS, cellT))
+                if (intersect(cellS, cellT))
                 {
                     seedCells[cellS] = cellT;
 
@@ -211,8 +271,6 @@ void Foam::meshToMeshNew::setNextCells
         bool restart =
             findInitialSeeds
             (
-                src,
-                tgt,
                 srcCellIDs,
                 mapFlag,
                 startSeedI,
@@ -233,68 +291,91 @@ void Foam::meshToMeshNew::setNextCells
 }
 
 
-bool Foam::meshToMeshNew::intersect
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::cellVolumeWeightMethod::cellVolumeWeightMethod
 (
     const polyMesh& src,
-    const polyMesh& tgt,
+    const polyMesh& tgt
-    const label srcCellI,
-    const label tgtCellI
-) const
-{
-    scalar threshold = tolerance_*src.cellVolumes()[srcCellI];
-
-    tetOverlapVolume overlapEngine;
-
-    treeBoundBox bbTgtCell
-    (
-        pointField
-        (
-            tgt.points(),
-            tgt.cellPoints()[tgtCellI]
 )
+:
+    meshToMeshMethod(src, tgt)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::cellVolumeWeightMethod::~cellVolumeWeightMethod()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::cellVolumeWeightMethod::calculate
+(
+    labelListList& srcToTgtAddr,
+    scalarListList& srcToTgtWght,
+    labelListList& tgtToSrcAddr,
+    scalarListList& tgtToSrcWght
+)
+{
+    bool ok = initialise
+    (
+        srcToTgtAddr,
+        srcToTgtWght,
+        tgtToSrcAddr,
+        tgtToSrcWght
     );
 
-    return overlapEngine.cellCellOverlapMinDecomp
+    if (!ok)
-    (
+    {
-        src,
+        return;
-        srcCellI,
-        tgt,
-        tgtCellI,
-        bbTgtCell,
-        threshold
-    );
     }
 
+    // (potentially) participating source mesh cells
+    const labelList srcCellIDs(maskCells());
 
-Foam::scalar Foam::meshToMeshNew::interVol
+    // list to keep track of whether src cell can be mapped
+    boolList mapFlag(src_.nCells(), false);
+    UIndirectList<bool>(mapFlag, srcCellIDs) = true;
+
+    // find initial point in tgt mesh
+    label srcSeedI = -1;
+    label tgtSeedI = -1;
+    label startSeedI = 0;
+
+    bool startWalk =
+        findInitialSeeds
         (
-    const polyMesh& src,
+            srcCellIDs,
-    const polyMesh& tgt,
+            mapFlag,
-    const label srcCellI,
+            startSeedI,
-    const label tgtCellI
+            srcSeedI,
-) const
+            tgtSeedI
+        );
+
+    if (startWalk)
     {
-    tetOverlapVolume overlapEngine;
+        calculateAddressing
-
-    treeBoundBox bbTgtCell
         (
-        pointField
+            srcToTgtAddr,
-        (
+            srcToTgtWght,
-            tgt.points(),
+            tgtToSrcAddr,
-            tgt.cellPoints()[tgtCellI]
+            tgtToSrcWght,
-        )
+            srcSeedI,
+            tgtSeedI,
+            srcCellIDs,
+            mapFlag,
+            startSeedI
         );
-
+    }
-    scalar vol = overlapEngine.cellCellOverlapVolumeMinDecomp
+    else
-    (
+    {
-        src,
+        // if meshes are collocated, after inflating the source mesh bounding
-        srcCellI,
+        // box tgt mesh cells may be transferred, but may still not overlap
-        tgt,
+        // with the source mesh
-        tgtCellI,
+        return;
-        bbTgtCell
+    }
-    );
-
-    return vol;
 }
 
 

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/cellVolumeWeight/cellVolumeWeightMethod.H

@@ -0,0 +1,138 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::cellVolumeWeightMethod
+
+Description
+    Cell-volume-weighted mesh-to-mesh interpolation class
+
+    Volume conservative.
+
+SourceFiles
+    cellVolumeWeightMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cellVolumeWeightMethod_H
+#define cellVolumeWeightMethod_H
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                   Class cellVolumeWeightMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class cellVolumeWeightMethod
+:
+    public meshToMeshMethod
+{
+protected:
+
+    // Protected Member Functions
+
+        //- Find indices of overlapping cells in src and tgt meshes - returns
+        //  true if found a matching pair
+        bool findInitialSeeds
+        (
+            const labelList& srcCellIDs,
+            const boolList& mapFlag,
+            const label startSeedI,
+            label& srcSeedI,
+            label& tgtSeedI
+        ) const;
+
+        //- Calculate the mesh-to-mesh addressing and weights
+        void calculateAddressing
+        (
+            labelListList& srcToTgtCellAddr,
+            scalarListList& srcToTgtCellWght,
+            labelListList& tgtToSrcCellAddr,
+            scalarListList& tgtToSrcCellWght,
+            const label srcSeedI,
+            const label tgtSeedI,
+            const labelList& srcCellIDs,
+            boolList& mapFlag,
+            label& startSeedI
+        );
+
+        //- Set the next cells in the advancing front algorithm
+        void setNextCells
+        (
+            label& startSeedI,
+            label& srcCellI,
+            label& tgtCellI,
+            const labelList& srcCellIDs,
+            const boolList& mapFlag,
+            const DynamicList<label>& visitedCells,
+            labelList& seedCells
+        ) const;
+
+        //- Disallow default bitwise copy construct
+        cellVolumeWeightMethod(const cellVolumeWeightMethod&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const cellVolumeWeightMethod&);
+
+
+public:
+
+    //- Run-time type information
+    TypeName("cellVolumeWeight");
+
+    //- Construct from source and target meshes
+    cellVolumeWeightMethod(const polyMesh& src, const polyMesh& tgt);
+
+    //- Destructor
+    virtual ~cellVolumeWeightMethod();
+
+
+    // Member Functions
+
+        // Evaluate
+
+            //- Calculate addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcToTgtAddr,
+                scalarListList& srcToTgtWght,
+                labelListList& tgtToTgtAddr,
+                scalarListList& tgtToTgtWght
+            );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.C

@@ -0,0 +1,303 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "directMethod.H"
+#include "indexedOctree.H"
+#include "treeDataCell.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(directMethod, 0);
+    addToRunTimeSelectionTable(meshToMeshMethod, directMethod, components);
+}
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+bool Foam::directMethod::findInitialSeeds
+(
+    const labelList& srcCellIDs,
+    const boolList& mapFlag,
+    const label startSeedI,
+    label& srcSeedI,
+    label& tgtSeedI
+) const
+{
+    const cellList& srcCells = src_.cells();
+    const faceList& srcFaces = src_.faces();
+    const pointField& srcPts = src_.points();
+
+    for (label i = startSeedI; i < srcCellIDs.size(); i++)
+    {
+        label srcI = srcCellIDs[i];
+
+        if (mapFlag[srcI])
+        {
+            const pointField
+                pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
+
+            forAll(pts, ptI)
+            {
+                const point& pt = pts[ptI];
+                label tgtI = tgt_.cellTree().findInside(pt);
+
+                if (tgtI != -1 && intersect(srcI, tgtI))
+                {
+                    srcSeedI = srcI;
+                    tgtSeedI = tgtI;
+
+                    return true;
+                }
+            }
+        }
+    }
+
+    if (debug)
+    {
+        Pout<< "could not find starting seed" << endl;
+    }
+
+    return false;
+}
+
+
+void Foam::directMethod::calculateAddressing
+(
+    labelListList& srcToTgtCellAddr,
+    scalarListList& srcToTgtCellWght,
+    labelListList& tgtToSrcCellAddr,
+    scalarListList& tgtToSrcCellWght,
+    const label srcSeedI,
+    const label tgtSeedI,
+    const labelList& srcCellIDs, // not used
+    boolList& mapFlag,
+    label& startSeedI
+)
+{
+    // store a list of src cells already mapped
+    labelList srcTgtSeed(src_.nCells(), -1);
+
+    List<DynamicList<label> > srcToTgt(src_.nCells());
+    List<DynamicList<label> > tgtToSrc(tgt_.nCells());
+
+    DynamicList<label> srcSeeds(10);
+
+    const scalarField& srcVc = src_.cellVolumes();
+    const scalarField& tgtVc = tgt_.cellVolumes();
+
+    label srcCellI = srcSeedI;
+    label tgtCellI = tgtSeedI;
+
+    do
+    {
+        // store src/tgt cell pair
+        srcToTgt[srcCellI].append(tgtCellI);
+        tgtToSrc[tgtCellI].append(srcCellI);
+
+        // mark source cell srcSeedI as matched
+        mapFlag[srcCellI] = false;
+
+        // accumulate intersection volume
+        V_ += srcVc[srcCellI];
+
+        // find new source seed cell
+        appendToDirectSeeds
+        (
+            mapFlag,
+            srcTgtSeed,
+            srcSeeds,
+            srcCellI,
+            tgtCellI
+        );
+    }
+    while (srcCellI >= 0);
+
+    // transfer addressing into persistent storage
+    forAll(srcToTgtCellAddr, i)
+    {
+        scalar v = srcVc[i];
+        srcToTgtCellAddr[i].transfer(srcToTgt[i]);
+        srcToTgtCellWght[i] = scalarList(1, v);
+    }
+
+    forAll(tgtToSrcCellAddr, i)
+    {
+        scalar v = tgtVc[i];
+        tgtToSrcCellAddr[i].transfer(tgtToSrc[i]);
+        tgtToSrcCellWght[i] = scalarList(1, v);
+    }
+}
+
+
+void Foam::directMethod::appendToDirectSeeds
+(
+    boolList& mapFlag,
+    labelList& srcTgtSeed,
+    DynamicList<label>& srcSeeds,
+    label& srcSeedI,
+    label& tgtSeedI
+) const
+{
+    const labelList& srcNbr = src_.cellCells()[srcSeedI];
+    const labelList& tgtNbr = tgt_.cellCells()[tgtSeedI];
+
+    const vectorField& srcCentre = src_.cellCentres();
+
+    forAll(srcNbr, i)
+    {
+        label srcI = srcNbr[i];
+
+        if (mapFlag[srcI] && (srcTgtSeed[srcI] == -1))
+        {
+            // source cell srcI not yet mapped
+
+            // identfy if target cell exists for source cell srcI
+            bool found = false;
+            forAll(tgtNbr, j)
+            {
+                label tgtI = tgtNbr[j];
+
+                if (tgt_.pointInCell(srcCentre[srcI], tgtI))
+                {
+                    // new match - append to lists
+                    found = true;
+
+                    srcTgtSeed[srcI] = tgtI;
+                    srcSeeds.append(srcI);
+
+                    break;
+                }
+            }
+
+            if (!found)
+            {
+                // no match available for source cell srcI
+                mapFlag[srcI] = false;
+            }
+        }
+    }
+
+    if (srcSeeds.size())
+    {
+        srcSeedI = srcSeeds.remove();
+        tgtSeedI = srcTgtSeed[srcSeedI];
+    }
+    else
+    {
+        srcSeedI = -1;
+        tgtSeedI = -1;
+    }
+}
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::directMethod::directMethod
+(
+    const polyMesh& src,
+    const polyMesh& tgt
+)
+:
+    meshToMeshMethod(src, tgt)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::directMethod::~directMethod()
+{}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::directMethod::calculate
+(
+    labelListList& srcToTgtAddr,
+    scalarListList& srcToTgtWght,
+    labelListList& tgtToSrcAddr,
+    scalarListList& tgtToSrcWght
+)
+{
+    bool ok = initialise
+    (
+        srcToTgtAddr,
+        srcToTgtWght,
+        tgtToSrcAddr,
+        tgtToSrcWght
+    );
+
+    if (!ok)
+    {
+        return;
+    }
+
+    // (potentially) participating source mesh cells
+    const labelList srcCellIDs(maskCells());
+
+    // list to keep track of whether src cell can be mapped
+    boolList mapFlag(src_.nCells(), false);
+    UIndirectList<bool>(mapFlag, srcCellIDs) = true;
+
+    // find initial point in tgt mesh
+    label srcSeedI = -1;
+    label tgtSeedI = -1;
+    label startSeedI = 0;
+
+    bool startWalk =
+        findInitialSeeds
+        (
+            srcCellIDs,
+            mapFlag,
+            startSeedI,
+            srcSeedI,
+            tgtSeedI
+        );
+
+    if (startWalk)
+    {
+        calculateAddressing
+        (
+            srcToTgtAddr,
+            srcToTgtWght,
+            tgtToSrcAddr,
+            tgtToSrcWght,
+            srcSeedI,
+            tgtSeedI,
+            srcCellIDs,
+            mapFlag,
+            startSeedI
+        );
+    }
+    else
+    {
+        // if meshes are collocated, after inflating the source mesh bounding
+        // box tgt mesh cells may be transferred, but may still not overlap
+        // with the source mesh
+        return;
+    }
+}
+
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/direct/directMethod.H

@@ -0,0 +1,136 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::directMethod
+
+Description
+    Direct (one-to-one cell correspondence) mesh-to-mesh interpolation class
+
+    Volume conservative.
+
+SourceFiles
+    directMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef directMethod_H
+#define directMethod_H
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                        Class directMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class directMethod
+:
+    public meshToMeshMethod
+{
+protected:
+
+    // Protected Member Functions
+
+        //- Find indices of overlapping cells in src and tgt meshes - returns
+        //  true if found a matching pair
+        bool findInitialSeeds
+        (
+            const labelList& srcCellIDs,
+            const boolList& mapFlag,
+            const label startSeedI,
+            label& srcSeedI,
+            label& tgtSeedI
+        ) const;
+
+        //- Calculate the mesh-to-mesh addressing and weights
+        void calculateAddressing
+        (
+            labelListList& srcToTgtCellAddr,
+            scalarListList& srcToTgtCellWght,
+            labelListList& tgtToSrcCellAddr,
+            scalarListList& tgtToSrcCellWght,
+            const label srcSeedI,
+            const label tgtSeedI,
+            const labelList& srcCellIDs,
+            boolList& mapFlag,
+            label& startSeedI
+        );
+
+        //- Append to list of src mesh seed indices
+        void appendToDirectSeeds
+        (
+            boolList& mapFlag,
+            labelList& srcTgtSeed,
+            DynamicList<label>& srcSeeds,
+            label& srcSeedI,
+            label& tgtSeedI
+        ) const;
+
+        //- Disallow default bitwise copy construct
+        directMethod(const directMethod&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const directMethod&);
+
+
+public:
+
+    //- Run-time type information
+    TypeName("direct");
+
+    //- Construct from source and target meshes
+    directMethod(const polyMesh& src, const polyMesh& tgt);
+
+    //- Destructor
+    virtual ~directMethod();
+
+
+    // Member Functions
+
+        // Evaluate
+
+            //- Calculate addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcToTgtAddr,
+                scalarListList& srcToTgtWght,
+                labelListList& tgtToTgtAddr,
+                scalarListList& tgtToTgtWght
+            );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.C

@@ -0,0 +1,424 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "mapNearestMethod.H"
+#include "pointIndexHit.H"
+#include "indexedOctree.H"
+#include "treeDataCell.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(mapNearestMethod, 0);
+    addToRunTimeSelectionTable(meshToMeshMethod, mapNearestMethod, components);
+}
+
+// * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * * //
+
+bool Foam::mapNearestMethod::findInitialSeeds
+(
+    const labelList& srcCellIDs,
+    const boolList& mapFlag,
+    const label startSeedI,
+    label& srcSeedI,
+    label& tgtSeedI
+) const
+{
+    const cellList& srcCells = src_.cells();
+    const faceList& srcFaces = src_.faces();
+    const pointField& srcPts = src_.points();
+
+    for (label i = startSeedI; i < srcCellIDs.size(); i++)
+    {
+        label srcI = srcCellIDs[i];
+
+        if (mapFlag[srcI])
+        {
+            const pointField
+                pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
+
+            const point& pt = pts[0];
+            pointIndexHit hit = tgt_.cellTree().findNearest(pt, GREAT);
+
+            if (hit.hit())
+            {
+                srcSeedI = srcI;
+                tgtSeedI = hit.index();
+
+                return true;
+            }
+            else
+            {
+                FatalErrorIn
+                (
+                    "bool Foam::mapNearestMethod::findInitialSeeds"
+                    "("
+                        "const labelList&, "
+                        "const boolList&, "
+                        "const label, "
+                        "label&, "
+                        "label&"
+                    ") const"
+                )
+                    << "Unable to find nearest target cell"
+                    << " for source cell " << srcI
+                    << " with centre "
+                    << srcCells[srcI].centre(srcPts, srcFaces)
+                    << abort(FatalError);
+            }
+
+            break;
+        }
+    }
+
+    if (debug)
+    {
+        Pout<< "could not find starting seed" << endl;
+    }
+
+    return false;
+}
+
+
+void Foam::mapNearestMethod::calculateAddressing
+(
+    labelListList& srcToTgtCellAddr,
+    scalarListList& srcToTgtCellWght,
+    labelListList& tgtToSrcCellAddr,
+    scalarListList& tgtToSrcCellWght,
+    const label srcSeedI,
+    const label tgtSeedI,
+    const labelList& srcCellIDs,
+    boolList& mapFlag,
+    label& startSeedI
+)
+{
+    List<DynamicList<label> > srcToTgt(src_.nCells());
+    List<DynamicList<label> > tgtToSrc(tgt_.nCells());
+
+    const scalarField& srcVc = src_.cellVolumes();
+    const scalarField& tgtVc = tgt_.cellVolumes();
+
+    label srcCellI = srcSeedI;
+    label tgtCellI = tgtSeedI;
+
+    do
+    {
+        // find nearest tgt cell
+        findNearestCell(src_, tgt_, srcCellI, tgtCellI);
+
+        // store src/tgt cell pair
+        srcToTgt[srcCellI].append(tgtCellI);
+        tgtToSrc[tgtCellI].append(srcCellI);
+
+        // mark source cell srcCellI and tgtCellI as matched
+        mapFlag[srcCellI] = false;
+
+        // accumulate intersection volume
+        V_ += srcVc[srcCellI];
+
+        // find new source cell
+        setNextNearestCells
+        (
+            startSeedI,
+            srcCellI,
+            tgtCellI,
+            mapFlag,
+            srcCellIDs
+        );
+    }
+    while (srcCellI >= 0);
+
+
+    // for the case of multiple source cells per target cell, select the
+    // nearest source cell only and discard the others
+    const vectorField& srcCc = src_.cellCentres();
+    const vectorField& tgtCc = tgt_.cellCentres();
+
+    forAll(tgtToSrc, targetCellI)
+    {
+        if (tgtToSrc[targetCellI].size() > 1)
+        {
+            const vector& tgtC = tgtCc[tgtCellI];
+
+            DynamicList<label>& srcCells = tgtToSrc[targetCellI];
+
+            label srcCellI = srcCells[0];
+            scalar d = magSqr(tgtC - srcCc[srcCellI]);
+
+            for (label i = 1; i < srcCells.size(); i++)
+            {
+                label srcI = srcCells[i];
+                scalar dNew = magSqr(tgtC - srcCc[srcI]);
+                if (dNew < d)
+                {
+                    d = dNew;
+                    srcCellI = srcI;
+                }
+            }
+
+            srcCells.clear();
+            srcCells.append(srcCellI);
+        }
+    }
+
+    // If there are more target cells than source cells, some target cells
+    // might not yet be mapped
+    forAll(tgtToSrc, tgtCellI)
+    {
+        if (tgtToSrc[tgtCellI].empty())
+        {
+            label srcCellI = findMappedSrcCell(tgtCellI, tgtToSrc);
+
+            findNearestCell(tgt_, src_, tgtCellI, srcCellI);
+
+            tgtToSrc[tgtCellI].append(srcCellI);
+        }
+    }
+
+    // transfer addressing into persistent storage
+    forAll(srcToTgtCellAddr, i)
+    {
+        scalar v = srcVc[i];
+        srcToTgtCellAddr[i].transfer(srcToTgt[i]);
+        srcToTgtCellWght[i] = scalarList(1, v);
+    }
+
+    forAll(tgtToSrcCellAddr, i)
+    {
+        scalar v = tgtVc[i];
+        tgtToSrcCellAddr[i].transfer(tgtToSrc[i]);
+        tgtToSrcCellWght[i] = scalarList(1, v);
+    }
+}
+
+
+void Foam::mapNearestMethod::findNearestCell
+(
+    const polyMesh& mesh1,
+    const polyMesh& mesh2,
+    const label cell1,
+    label& cell2
+) const
+{
+    const vectorField& Cc1 = mesh1.cellCentres();
+    const vectorField& Cc2 = mesh2.cellCentres();
+
+    const vector& p1 = Cc1[cell1];
+
+    DynamicList<label> cells2(10);
+    cells2.append(cell2);
+
+    DynamicList<label> visitedCells(10);
+
+    scalar d = GREAT;
+
+    do
+    {
+        label c2 = cells2.remove();
+        visitedCells.append(c2);
+
+        scalar dTest = magSqr(Cc2[c2] - p1);
+        if (dTest < d)
+        {
+            cell2 = c2;
+            d = dTest;
+            appendNbrCells(cell2, mesh2, visitedCells, cells2);
+        }
+
+    } while (cells2.size() > 0);
+}
+
+
+void Foam::mapNearestMethod::setNextNearestCells
+(
+    label& startSeedI,
+    label& srcCellI,
+    label& tgtCellI,
+    boolList& mapFlag,
+    const labelList& srcCellIDs
+) const
+{
+    const labelList& srcNbr = src_.cellCells()[srcCellI];
+
+    srcCellI = -1;
+    forAll(srcNbr, i)
+    {
+        label cellI = srcNbr[i];
+        if (mapFlag[cellI])
+        {
+            srcCellI = cellI;
+            startSeedI = cellI + 1;
+
+            return;
+        }
+    }
+
+    (void)findInitialSeeds
+    (
+        srcCellIDs,
+        mapFlag,
+        startSeedI,
+        srcCellI,
+        tgtCellI
+    );
+}
+
+
+Foam::label Foam::mapNearestMethod::findMappedSrcCell
+(
+    const label tgtCellI,
+    const List<DynamicList<label> >& tgtToSrc
+) const
+{
+    DynamicList<label> testCells(10);
+    DynamicList<label> visitedCells(10);
+
+    testCells.append(tgtCellI);
+
+    do
+    {
+        // search target tgtCellI neighbours for match with source cell
+        label tgtI = testCells.remove();
+
+        if (findIndex(visitedCells, tgtI) == -1)
+        {
+            visitedCells.append(tgtI);
+
+            if (tgtToSrc[tgtI].size())
+            {
+                return tgtToSrc[tgtI][0];
+            }
+            else
+            {
+                const labelList& nbrCells = tgt_.cellCells()[tgtI];
+
+                forAll(nbrCells, i)
+                {
+                    if (findIndex(visitedCells, nbrCells[i]) == -1)
+                    {
+                        testCells.append(nbrCells[i]);
+                    }
+                }
+            }
+        }
+    } while (testCells.size());
+
+    // did not find any match - should not be possible to get here!
+    return -1;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::mapNearestMethod::mapNearestMethod
+(
+    const polyMesh& src,
+    const polyMesh& tgt
+)
+:
+    meshToMeshMethod(src, tgt)
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * * //
+
+Foam::mapNearestMethod::~mapNearestMethod()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::mapNearestMethod::calculate
+(
+    labelListList& srcToTgtAddr,
+    scalarListList& srcToTgtWght,
+    labelListList& tgtToSrcAddr,
+    scalarListList& tgtToSrcWght
+)
+{
+    bool ok = initialise
+    (
+        srcToTgtAddr,
+        srcToTgtWght,
+        tgtToSrcAddr,
+        tgtToSrcWght
+    );
+
+    if (!ok)
+    {
+        return;
+    }
+
+    // (potentially) participating source mesh cells
+    const labelList srcCellIDs(maskCells());
+
+    // list to keep track of whether src cell can be mapped
+    boolList mapFlag(src_.nCells(), false);
+    UIndirectList<bool>(mapFlag, srcCellIDs) = true;
+
+    // find initial point in tgt mesh
+    label srcSeedI = -1;
+    label tgtSeedI = -1;
+    label startSeedI = 0;
+
+    bool startWalk =
+        findInitialSeeds
+        (
+            srcCellIDs,
+            mapFlag,
+            startSeedI,
+            srcSeedI,
+            tgtSeedI
+        );
+
+    if (startWalk)
+    {
+        calculateAddressing
+        (
+            srcToTgtAddr,
+            srcToTgtWght,
+            tgtToSrcAddr,
+            tgtToSrcWght,
+            srcSeedI,
+            tgtSeedI,
+            srcCellIDs,
+            mapFlag,
+            startSeedI
+        );
+    }
+    else
+    {
+        // if meshes are collocated, after inflating the source mesh bounding
+        // box tgt mesh cells may be transferred, but may still not overlap
+        // with the source mesh
+        return;
+    }
+}
+
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/mapNearest/mapNearestMethod.H

@@ -0,0 +1,152 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::mapNearestMethod
+
+Description
+    Map nearest mesh-to-mesh interpolation class
+
+    Not volume conservative.
+
+SourceFiles
+    mapNearestMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef mapNearestMethod_H
+#define mapNearestMethod_H
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class mapNearestMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class mapNearestMethod
+:
+    public meshToMeshMethod
+{
+protected:
+
+    // Protected Member Functions
+
+        //- Find indices of overlapping cells in src and tgt meshes - returns
+        //  true if found a matching pair
+        bool findInitialSeeds
+        (
+            const labelList& srcCellIDs,
+            const boolList& mapFlag,
+            const label startSeedI,
+            label& srcSeedI,
+            label& tgtSeedI
+        ) const;
+
+        //- Calculate the mesh-to-mesh addressing and weights
+        void calculateAddressing
+        (
+            labelListList& srcToTgtCellAddr,
+            scalarListList& srcToTgtCellWght,
+            labelListList& tgtToSrcCellAddr,
+            scalarListList& tgtToSrcCellWght,
+            const label srcSeedI,
+            const label tgtSeedI,
+            const labelList& srcCellIDs,
+            boolList& mapFlag,
+            label& startSeedI
+        );
+
+        //- Find the nearest cell on mesh2 for cell1 on mesh1
+        void findNearestCell
+        (
+            const polyMesh& mesh1,
+            const polyMesh& mesh2,
+            const label cell1,
+            label& cell2
+        ) const;
+
+        //- Set the next cells for the marching front algorithm
+        void setNextNearestCells
+        (
+            label& startSeedI,
+            label& srcCellI,
+            label& tgtCellI,
+            boolList& mapFlag,
+            const labelList& srcCellIDs
+        ) const;
+
+        //- Find a source cell mapped to target cell tgtCellI
+        label findMappedSrcCell
+        (
+            const label tgtCellI,
+            const List<DynamicList<label> >& tgtToSrc
+        ) const;
+
+        //- Disallow default bitwise copy construct
+        mapNearestMethod(const mapNearestMethod&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const mapNearestMethod&);
+
+
+public:
+
+    //- Run-time type information
+    TypeName("mapNearest");
+
+    //- Construct from source and target meshes
+    mapNearestMethod(const polyMesh& src, const polyMesh& tgt);
+
+    //- Destructor
+    virtual ~mapNearestMethod();
+
+
+    // Member Functions
+
+        // Evaluate
+
+            //- Calculate addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcToTgtAddr,
+                scalarListList& srcToTgtWght,
+                labelListList& tgtToTgtAddr,
+                scalarListList& tgtToTgtWght
+            );
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.C

@@ -0,0 +1,274 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "meshToMeshMethod.H"
+#include "tetOverlapVolume.H"
+#include "OFstream.H"
+#include "Time.H"
+#include "treeBoundBox.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(meshToMeshMethod, 0);
+    defineRunTimeSelectionTable(meshToMeshMethod, components);
+}
+
+Foam::scalar Foam::meshToMeshMethod::tolerance_ = 1e-6;
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+Foam::labelList Foam::meshToMeshMethod::maskCells() const
+{
+    boundBox intersectBb
+    (
+        max(src_.bounds().min(), tgt_.bounds().min()),
+        min(src_.bounds().max(), tgt_.bounds().max())
+    );
+
+    intersectBb.inflate(0.01);
+
+    const cellList& srcCells = src_.cells();
+    const faceList& srcFaces = src_.faces();
+    const pointField& srcPts = src_.points();
+
+    DynamicList<label> cells(src_.nCells());
+    forAll(srcCells, srcI)
+    {
+        boundBox cellBb(srcCells[srcI].points(srcFaces, srcPts), false);
+        if (intersectBb.overlaps(cellBb))
+        {
+            cells.append(srcI);
+        }
+    }
+
+    if (debug)
+    {
+        Pout<< "participating source mesh cells: " << cells.size() << endl;
+    }
+
+    return cells;
+}
+
+
+bool Foam::meshToMeshMethod::intersect
+(
+    const label srcCellI,
+    const label tgtCellI
+) const
+{
+    scalar threshold = tolerance_*src_.cellVolumes()[srcCellI];
+
+    tetOverlapVolume overlapEngine;
+
+    treeBoundBox bbTgtCell
+    (
+        pointField
+        (
+            tgt_.points(),
+            tgt_.cellPoints()[tgtCellI]
+        )
+    );
+
+    return overlapEngine.cellCellOverlapMinDecomp
+    (
+        src_,
+        srcCellI,
+        tgt_,
+        tgtCellI,
+        bbTgtCell,
+        threshold
+    );
+}
+
+
+Foam::scalar Foam::meshToMeshMethod::interVol
+(
+    const label srcCellI,
+    const label tgtCellI
+) const
+{
+    tetOverlapVolume overlapEngine;
+
+    treeBoundBox bbTgtCell
+    (
+        pointField
+        (
+            tgt_.points(),
+            tgt_.cellPoints()[tgtCellI]
+        )
+    );
+
+    scalar vol = overlapEngine.cellCellOverlapVolumeMinDecomp
+    (
+        src_,
+        srcCellI,
+        tgt_,
+        tgtCellI,
+        bbTgtCell
+    );
+
+    return vol;
+}
+
+
+void Foam::meshToMeshMethod::appendNbrCells
+(
+    const label cellI,
+    const polyMesh& mesh,
+    const DynamicList<label>& visitedCells,
+    DynamicList<label>& nbrCellIDs
+) const
+{
+    const labelList& nbrCells = mesh.cellCells()[cellI];
+
+    // filter out cells already visited from cell neighbours
+    forAll(nbrCells, i)
+    {
+        label nbrCellI = nbrCells[i];
+
+        if
+        (
+            (findIndex(visitedCells, nbrCellI) == -1)
+         && (findIndex(nbrCellIDs, nbrCellI) == -1)
+        )
+        {
+            nbrCellIDs.append(nbrCellI);
+        }
+    }
+}
+
+
+bool Foam::meshToMeshMethod::initialise
+(
+    labelListList& srcToTgtAddr,
+    scalarListList& srcToTgtWght,
+    labelListList& tgtToSrcAddr,
+    scalarListList& tgtToSrcWght
+) const
+{
+    srcToTgtAddr.setSize(src_.nCells());
+    srcToTgtWght.setSize(src_.nCells());
+    tgtToSrcAddr.setSize(tgt_.nCells());
+    tgtToSrcWght.setSize(tgt_.nCells());
+
+    if (!src_.nCells())
+    {
+        return false;
+    }
+    else if (!tgt_.nCells())
+    {
+        if (debug)
+        {
+            Pout<< "mesh interpolation: hhave " << src_.nCells() << " source "
+                << " cells but no target cells" << endl;
+        }
+
+        return false;
+    }
+
+    return true;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::meshToMeshMethod::meshToMeshMethod
+(
+    const polyMesh& src,
+    const polyMesh& tgt
+)
+:
+    src_(src),
+    tgt_(tgt),
+    V_(0.0)
+{
+    if (!src_.nCells() || !tgt_.nCells())
+    {
+        if (debug)
+        {
+            Pout<< "mesh interpolation: cells not on processor: Source cells = "
+                << src_.nCells() << ", target cells = " << tgt_.nCells()
+                << endl;
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::meshToMeshMethod::~meshToMeshMethod()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+void Foam::meshToMeshMethod::writeConnectivity
+(
+    const polyMesh& mesh1,
+    const polyMesh& mesh2,
+    const labelListList& mesh1ToMesh2Addr
+) const
+{
+    Pout<< "Source size = " << mesh1.nCells() << endl;
+    Pout<< "Target size = " << mesh2.nCells() << endl;
+
+    word fName("addressing_" + mesh1.name() + "_to_" + mesh2.name());
+
+    if (Pstream::parRun())
+    {
+        fName = fName +  "_proc" + Foam::name(Pstream::myProcNo());
+    }
+
+    OFstream os(src_.time().path()/fName + ".obj");
+
+    label vertI = 0;
+    forAll(mesh1ToMesh2Addr, i)
+    {
+        const labelList& addr = mesh1ToMesh2Addr[i];
+        forAll(addr, j)
+        {
+            label cellI = addr[j];
+            const vector& c0 = mesh1.cellCentres()[i];
+
+            const cell& c = mesh2.cells()[cellI];
+            const pointField pts(c.points(mesh2.faces(), mesh2.points()));
+            forAll(pts, j)
+            {
+                const point& p = pts[j];
+                os  << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
+                vertI++;
+                os  << "v " << c0.x() << ' ' << c0.y() << ' ' << c0.z()
+                    << nl;
+                vertI++;
+                os  << "l " << vertI - 1 << ' ' << vertI << nl;
+            }
+        }
+    }
+}
+
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethod.H

@@ -0,0 +1,188 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::meshToMeshMethod
+
+Description
+    Base class for mesh-to-mesh calculation methods
+
+SourceFiles
+    meshToMeshMethod.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef meshToMeshMethod_H
+#define meshToMeshMethod_H
+
+#include "polyMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class meshToMeshMethod Declaration
+\*---------------------------------------------------------------------------*/
+
+class meshToMeshMethod
+{
+
+protected:
+
+    // Protected data
+
+        //- Reference to the source mesh
+        const polyMesh& src_;
+
+        //- Reference to the target mesh
+        const polyMesh& tgt_;
+
+        //- Cell total volume in overlap region [m3]
+        scalar V_;
+
+        //- Tolerance used in volume overlap calculations
+        static scalar tolerance_;
+
+
+    // Protected Member Functions
+
+        //- Return src cell IDs for the overlap region
+        labelList maskCells() const;
+
+        //- Return the true if cells intersect
+        bool intersect(const label srcCellI, const label tgtCellI) const;
+
+        //- Return the intersection volume between two cells
+        scalar interVol(const label srcCellI, const label tgtCellI) const;
+
+        //- Append target cell neihgbour cells to cellIDs list
+        void appendNbrCells
+        (
+            const label tgtCellI,
+            const polyMesh& mesh,
+            const DynamicList<label>& visitedTgtCells,
+            DynamicList<label>& nbrTgtCellIDs
+        ) const;
+
+        bool initialise
+        (
+            labelListList& srcToTgtAddr,
+            scalarListList& srcToTgtWght,
+            labelListList& tgtToTgtAddr,
+            scalarListList& tgtToTgtWght
+        ) const;
+
+        //- Disallow default bitwise copy construct
+        meshToMeshMethod(const meshToMeshMethod&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const meshToMeshMethod&);
+
+
+public:
+
+    //- Run-time type information
+    TypeName("meshToMeshMethod");
+
+    //- Declare runtime constructor selection table
+    declareRunTimeSelectionTable
+    (
+        autoPtr,
+        meshToMeshMethod,
+        components,
+        (
+            const polyMesh& src,
+            const polyMesh& tgt
+        ),
+        (src, tgt)
+    );
+
+    //- Construct from source and target meshes
+    meshToMeshMethod(const polyMesh& src, const polyMesh& tgt);
+
+    //- Selector
+    static autoPtr<meshToMeshMethod> New
+    (
+        const word& methodName,
+        const polyMesh& src,
+        const polyMesh& tgt
+    );
+
+
+    //- Destructor
+    virtual ~meshToMeshMethod();
+
+
+    // Member Functions
+
+        // Evaluate
+
+            //- Calculate addressing and weights
+            virtual void calculate
+            (
+                labelListList& srcToTgtAddr,
+                scalarListList& srcToTgtWght,
+                labelListList& tgtToTgtAddr,
+                scalarListList& tgtToTgtWght
+            ) = 0;
+
+
+        // Access
+
+            //- Return const access to the source mesh
+            inline const polyMesh& src() const;
+
+            //- Return const access to the target mesh
+            inline const polyMesh& tgt() const;
+
+            //- Return const access to the overlap volume
+            inline scalar V() const;
+
+
+        // Check
+
+            //- Write the connectivity (debugging)
+            void writeConnectivity
+            (
+                const polyMesh& mesh1,
+                const polyMesh& mesh2,
+                const labelListList& mesh1ToMesh2Addr
+            ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "meshToMeshMethodI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodI.H

@@ -0,0 +1,44 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+const Foam::polyMesh& Foam::meshToMeshMethod::src() const
+{
+    return src_;
+}
+
+
+const Foam::polyMesh& Foam::meshToMeshMethod::tgt() const
+{
+    return tgt_;
+}
+
+
+Foam::scalar Foam::meshToMeshMethod::V() const
+{
+    return V_;
+}
+
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/calcMethod/meshToMeshMethod/meshToMeshMethodNew.C

@@ -0,0 +1,65 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2013 OpenFOAM Foundation
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "meshToMeshMethod.H"
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+Foam::autoPtr<Foam::meshToMeshMethod> Foam::meshToMeshMethod::New
+(
+    const word& methodName,
+    const polyMesh& src,
+    const polyMesh& tgt
+)
+{
+    if (debug)
+    {
+        Info<< "Selecting AMIMethod " << methodName << endl;
+    }
+
+    componentsConstructorTable::iterator cstrIter =
+        componentsConstructorTablePtr_->find(methodName);
+
+    if (cstrIter == componentsConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "Foam::autoPtr<Foam::meshToMeshMethod> Foam::meshToMeshMethod::New"
+            "("
+                "const word&, "
+                "const polyMesh&, "
+                "const polyMesh&"
+            ")"
+        )   << "Unknown meshToMesh type "
+            << methodName << nl << nl
+            << "Valid meshToMesh types are:" << nl
+            << componentsConstructorTablePtr_->sortedToc() << exit(FatalError);
+    }
+
+    return autoPtr<meshToMeshMethod>(cstrIter()(src, tgt));
+}
+
+
+// ************************************************************************* //

src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.C

@@ -24,14 +24,9 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "meshToMeshNew.H"
-#include "OFstream.H"
 #include "Time.H"
 #include "globalIndex.H"
-#include "mergePoints.H"
+#include "meshToMeshMethod.H"
-#include "treeBoundBox.H"
-#include "indexedOctree.H"
-#include "treeDataCell.H"
-#include "ListOps.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -55,55 +50,9 @@ namespace Foam
         meshToMeshNew::interpolationMethodNames_;
 }
 
-Foam::scalar Foam::meshToMeshNew::tolerance_ = 1e-6;
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-void Foam::meshToMeshNew::writeConnectivity
-(
-    const polyMesh& src,
-    const polyMesh& tgt,
-    const labelListList& srcToTargetAddr
-) const
-{
-    Pout<< "Source size = " << src.nCells() << endl;
-    Pout<< "Target size = " << tgt.nCells() << endl;
-
-    word fName("addressing_" + src.name() + "_to_" + tgt.name());
-
-    if (Pstream::parRun())
-    {
-        fName = fName +  "_proc" + Foam::name(Pstream::myProcNo());
-    }
-
-    OFstream os(src.time().path()/fName + ".obj");
-
-    label vertI = 0;
-    forAll(srcToTargetAddr, i)
-    {
-        const labelList& tgtAddress = srcToTargetAddr[i];
-        forAll(tgtAddress, j)
-        {
-            label tgtI = tgtAddress[j];
-            const vector& c0 = src.cellCentres()[i];
-
-            const cell& c = tgt.cells()[tgtI];
-            const pointField pts(c.points(tgt.faces(), tgt.points()));
-            forAll(pts, j)
-            {
-                const point& p = pts[j];
-                os  << "v " << p.x() << ' ' << p.y() << ' ' << p.z() << nl;
-                vertI++;
-                os  << "v " << c0.x() << ' ' << c0.y() << ' ' << c0.z()
-                    << nl;
-                vertI++;
-                os  << "l " << vertI - 1 << ' ' << vertI << nl;
-            }
-        }
-    }
-}
-
-
 Foam::labelList Foam::meshToMeshNew::maskCells
 (
     const polyMesh& src,
@@ -141,82 +90,6 @@ Foam::labelList Foam::meshToMeshNew::maskCells
 }
 
 
-bool Foam::meshToMeshNew::findInitialSeeds
-(
-    const polyMesh& src,
-    const polyMesh& tgt,
-    const labelList& srcCellIDs,
-    const boolList& mapFlag,
-    const label startSeedI,
-    label& srcSeedI,
-    label& tgtSeedI
-) const
-{
-    const cellList& srcCells = src.cells();
-    const faceList& srcFaces = src.faces();
-    const pointField& srcPts = src.points();
-
-    for (label i = startSeedI; i < srcCellIDs.size(); i++)
-    {
-        label srcI = srcCellIDs[i];
-
-        if (mapFlag[srcI])
-        {
-            const pointField
-                pts(srcCells[srcI].points(srcFaces, srcPts).xfer());
-
-            forAll(pts, ptI)
-            {
-                const point& pt = pts[ptI];
-                label tgtI = tgt.cellTree().findInside(pt);
-
-                if (tgtI != -1 && intersect(src, tgt, srcI, tgtI))
-                {
-                    srcSeedI = srcI;
-                    tgtSeedI = tgtI;
-
-                    return true;
-                }
-            }
-        }
-    }
-
-    if (debug)
-    {
-        Pout<< "could not find starting seed" << endl;
-    }
-
-    return false;
-}
-
-
-void Foam::meshToMeshNew::appendNbrCells
-(
-    const label cellI,
-    const polyMesh& mesh,
-    const DynamicList<label>& visitedCells,
-    DynamicList<label>& nbrCellIDs
-) const
-{
-    const labelList& nbrCells = mesh.cellCells()[cellI];
-
-    // filter out cells already visited from cell neighbours
-    forAll(nbrCells, i)
-    {
-        label nbrCellI = nbrCells[i];
-
-        if
-        (
-            (findIndex(visitedCells, nbrCellI) == -1)
-         && (findIndex(nbrCellIDs, nbrCellI) == -1)
-        )
-        {
-            nbrCellIDs.append(nbrCellI);
-        }
-    }
-}
-
-
 void Foam::meshToMeshNew::normaliseWeights
 (
     const word& descriptor,
@@ -260,124 +133,29 @@ void Foam::meshToMeshNew::calcAddressing
     const polyMesh& tgt
 )
 {
-    srcToTgtCellAddr_.setSize(src.nCells());
+    autoPtr<meshToMeshMethod> methodPtr
-    srcToTgtCellWght_.setSize(src.nCells());
-
-    tgtToSrcCellAddr_.setSize(tgt.nCells());
-    tgtToSrcCellWght_.setSize(tgt.nCells());
-
-    if (!src.nCells() || !tgt.nCells())
-    {
-        if (debug)
-        {
-            Pout<< "mesh interpolation: cells not on processor: Source cells = "
-                << src.nCells() << ", target cells = " << tgt.nCells()
-                << endl;
-        }
-    }
-
-    if (!src.nCells())
-    {
-        return;
-    }
-    else if (!tgt.nCells())
-    {
-        if (debug)
-        {
-            Pout<< "mesh interpolation: hhave " << src.nCells() << " source "
-                << " cells but no target cells" << endl;
-        }
-
-        return;
-    }
-
-    // (potentially) participating source mesh cells
-    const labelList srcCellIDs = maskCells(src, tgt);
-
-    // list to keep track of whether src cell can be mapped
-    boolList mapFlag(src.nCells(), false);
-    UIndirectList<bool>(mapFlag, srcCellIDs) = true;
-
-    // find initial point in tgt mesh
-    label srcSeedI = -1;
-    label tgtSeedI = -1;
-    label startSeedI = 0;
-
-    bool startWalk =
-        findInitialSeeds
     (
+        meshToMeshMethod::New
+        (
+            interpolationMethodNames_[method_],
             src,
-            tgt,
+            tgt
-            srcCellIDs,
-            mapFlag,
-            startSeedI,
-            srcSeedI,
-            tgtSeedI
-        );
-
-    if (!startWalk)
-    {
-        // if meshes are collocated, after inflating the source mesh bounding
-        // box tgt mesh cells may be transferred, but may still not overlap
-        // with the source mesh
-        return;
-    }
-
-
-    switch (method_)
-    {
-        case imDirect:
-        {
-            calcDirect(src, tgt, srcSeedI, tgtSeedI);
-            break;
-        }
-        case imMapNearest:
-        {
-            calcMapNearest
-            (
-                src,
-                tgt,
-                srcSeedI,
-                tgtSeedI,
-                srcCellIDs,
-                mapFlag,
-                startSeedI
-            );
-            break;
-        }
-        case imCellVolumeWeight:
-        {
-            calcCellVolumeWeight
-            (
-                src,
-                tgt,
-                srcSeedI,
-                tgtSeedI,
-                srcCellIDs,
-                mapFlag,
-                startSeedI
-            );
-            break;
-        }
-        default:
-        {
-            FatalErrorIn
-            (
-                "void Foam::meshToMeshNew::calcAddressing"
-                "("
-                    "const polyMesh&, "
-                    "const polyMesh&"
-                ")"
         )
-                << "Unknown interpolation method"
+    );
-                << abort(FatalError);
-        }
-    }
 
+    methodPtr->calculate
+    (
+        srcToTgtCellAddr_,
+        srcToTgtCellWght_,
+        tgtToSrcCellAddr_,
+        tgtToSrcCellWght_
+    );
+
+    V_ = methodPtr->V();
 
     if (debug)
     {
-        writeConnectivity(src, tgt, srcToTgtCellAddr_);
+        methodPtr->writeConnectivity(src, tgt, srcToTgtCellAddr_);
     }
 }
 
@@ -577,11 +355,59 @@ void Foam::meshToMeshNew::calculate()
 }
 
 
+Foam::AMIPatchToPatchInterpolation::interpolationMethod
+Foam::meshToMeshNew::interpolationMethodAMI
+(
+    const interpolationMethod method
+) const
+{
+    switch (method_)
+    {
+        case imDirect:
+        {
+            return AMIPatchToPatchInterpolation::imDirect;
+            break;
+        }
+        case imMapNearest:
+        {
+            return AMIPatchToPatchInterpolation::imMapNearest;
+            break;
+        }
+        case imCellVolumeWeight:
+        {
+            return AMIPatchToPatchInterpolation::imFaceAreaWeight;
+            break;
+        }
+        default:
+        {
+            FatalErrorIn
+            (
+                "Foam::AMIPatchToPatchInterpolation::interpolationMethod"
+                "Foam::meshToMeshNew::interpolationMethodAMI"
+                "("
+                    "const interpolationMethod method"
+                ") const"
+            )
+                << "Unhandled enumeration " << method_
+                << abort(FatalError);
+        }
+    }
+
+    return AMIPatchToPatchInterpolation::imDirect;
+}
+
+
 const Foam::PtrList<Foam::AMIPatchToPatchInterpolation>&
 Foam::meshToMeshNew::patchAMIs() const
 {
     if (patchAMIs_.empty())
     {
+        const word amiMethod =
+            AMIPatchToPatchInterpolation::interpolationMethodToWord
+            (
+                interpolationMethodAMI(method_)
+            );
+
         patchAMIs_.setSize(srcPatchID_.size());
 
         forAll(srcPatchID_, i)
@@ -593,7 +419,9 @@ Foam::meshToMeshNew::patchAMIs() const
             const polyPatch& tgtPP = tgtRegion_.boundaryMesh()[tgtPatchI];
 
             Info<< "Creating AMI between source patch " << srcPP.name()
-                << " and target patch " << tgtPP.name() << endl;
+                << " and target patch " << tgtPP.name()
+                << " using " << amiMethod
+                << endl;
 
             Info<< incrIndent;
 
@@ -605,6 +433,7 @@ Foam::meshToMeshNew::patchAMIs() const
                     srcPP,
                     tgtPP,
                     faceAreaIntersect::tmMesh,
+                    interpolationMethodAMI(method_),
                     true // flip target patch since patch normals are aligned
                 )
             );

src/sampling/meshToMeshInterpolation/meshToMeshNew/meshToMeshNew.H

@@ -27,20 +27,14 @@ Class
 Description
     Class to calculate the cell-addressing between two overlapping meshes
 
-    Three methods are currently available:
+    Mapping is performed using a run-time selectable interpolation mothod
-    - direct           : 1-to-1 mapping between meshes
-    - mapNearest       : assign nearest cell values without interpolation
-    - cellVolumeWeight : volume consistent mapping
-
-    The \c direct and \c cellVolumeWeight options are volume conservative,
-    whereas mapNearest is non-conservative.
 
+SeeAlso
+    meshToMeshMethod
 
 SourceFiles
-    calcDirect.C
-    calcMapNearest.C
-    calcCellVolumeWeight.C
     meshToMeshNew.C
+    meshToMeshNewParallelOps.C
     meshToMeshNewTemplates.C
 
 \*---------------------------------------------------------------------------*/
@@ -70,7 +64,7 @@ public:
 
     // Public data types
 
-        //- Enumeration specifying required accuracy
+        //- Enumeration specifying interpolation method
         enum interpolationMethod
         {
             imDirect,
@@ -128,9 +122,6 @@ private:
         //- Target map pointer - parallel running only
         autoPtr<mapDistribute> tgtMapPtr_;
 
-        //- Tolerance used in volume overlap calculations
-        static scalar tolerance_;
-
 
     // Private Member Functions
 
@@ -138,39 +129,9 @@ private:
         template<class Type>
         void add(UList<Type>& fld, const label offset) const;
 
-        //- Write the connectivity (debugging)
-        void writeConnectivity
-        (
-            const polyMesh& src,
-            const polyMesh& tgt,
-            const labelListList& srcToTargetAddr
-        ) const;
-
         //- Return src cell IDs for the overlap region
         labelList maskCells(const polyMesh& src, const polyMesh& tgt) const;
 
-        //- Find indices of overlapping cells in src and tgt meshes - returns
-        //  true if found a matching pair
-        bool findInitialSeeds
-        (
-            const polyMesh& src,
-            const polyMesh& tgt,
-            const labelList& srcCellIDs,
-            const boolList& mapFlag,
-            const label startSeedI,
-            label& srcSeedI,
-            label& tgtSeedI
-        ) const;
-
-        //- Append target cell neihgbour cells to cellIDs list
-        void appendNbrCells
-        (
-            const label tgtCellI,
-            const polyMesh& tgt,
-            const DynamicList<label>& visitedTgtCells,
-            DynamicList<label>& nbrTgtCellIDs
-        ) const;
-
         //- Normalise the interpolation weights
         void normaliseWeights
         (
@@ -180,121 +141,6 @@ private:
             scalarListList& wght
         ) const;
 
-
-        // Direct (one-to-one) mapping
-
-            //- Main driver routine for direct mapping
-            void calcDirect
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const label srcSeedI,
-                const label tgtSeedI
-            );
-
-            //- Append to list of src mesh seed indices
-            void appendToDirectSeeds
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                boolList& mapFlag,
-                labelList& srcTgtSeed,
-                DynamicList<label>& srcSeeds,
-                label& srcSeedI,
-                label& tgtSeedI
-            ) const;
-
-        // Nearest (non-conformal) mapping
-
-            //- Main driver routine for nearest-mapping routine
-            void calcMapNearest
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const label srcSeedI,
-                const label tgtSeedI,
-                const labelList& srcCellIDs,
-                boolList& mapFlag,
-                label& startSeedI
-            );
-
-            //- Find target cell index of cell closest to source cell
-            void findNearestCell
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const label srcCellI,
-                label& tgtCellI
-            );
-
-            //- Set the next pair of cells
-            void setNextNearestCells
-            (
-                label& startSeedI,
-                label& srcCellI,
-                label& tgtCellI,
-                boolList& mapFlag,
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const labelList& srcCellIDs
-            );
-
-            //- Find source cell for target cell
-            label findMappedSrcCell
-            (
-                const polyMesh& tgt,
-                const label tgtCellI,
-                const List<DynamicList<label> >& tgtToSrc
-            ) const;
-
-
-        // Cell volume weighted (non-conformal) interpolation
-
-            //- Main driver routine for cell volume weighted interpolation
-            void calcCellVolumeWeight
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const label srcSeedI,
-                const label tgtSeedI,
-                const labelList& srcCellIDs,
-                boolList& mapFlag,
-                label& startSeedI
-            );
-
-            //- Set the next cells in the advancing front algorithm
-            void setNextCells
-            (
-                label& startSeedI,
-                label& srcCellI,
-                label& tgtCellI,
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const labelList& srcCellIDs,
-                const boolList& mapFlag,
-                const DynamicList<label>& visitedCells,
-                labelList& seedCells
-            ) const;
-
-            //- Return the true if cells intersect
-            bool intersect
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const label srcCellI,
-                const label tgtCellI
-            ) const;
-
-            //- Return the intersection volume between two cells
-            scalar interVol
-            (
-                const polyMesh& src,
-                const polyMesh& tgt,
-                const label srcCellI,
-                const label tgtCellI
-            ) const;
-
-
         //- Calculate the addressing between overalping regions of src and tgt
         //  meshes
         void calcAddressing(const polyMesh& src, const polyMesh& tgt);
@@ -302,6 +148,13 @@ private:
         //- Calculate - main driver function
         void calculate();
 
+        //- Conversion between mesh and patch interpolation methods
+        AMIPatchToPatchInterpolation::interpolationMethod
+        interpolationMethodAMI
+        (
+            const interpolationMethod method
+        ) const;
+
         //- Return the list of AMIs between source and target patches
         const PtrList<AMIPatchToPatchInterpolation>& patchAMIs() const;
 

src/sampling/sampledSurface/sampledSurfaces/sampledSurfacesTemplates.C

@@ -27,6 +27,7 @@ License
 #include "volFields.H"
 #include "surfaceFields.H"
 #include "ListListOps.H"
+#include "stringListOps.H"
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
@@ -146,7 +147,6 @@ void Foam::sampledSurfaces::sampleAndWrite
 }
 
 
-
 template<class Type>
 void Foam::sampledSurfaces::sampleAndWrite
 (
@@ -168,19 +168,32 @@ void Foam::sampledSurfaces::sampleAndWrite
 template<class GeoField>
 void Foam::sampledSurfaces::sampleAndWrite(const IOobjectList& objects)
 {
+    wordList names;
     if (loadFromFiles_)
     {
         IOobjectList fieldObjects(objects.lookupClass(GeoField::typeName));
-
+        names = fieldObjects.names();
-        forAll(fieldSelection_, fieldI)
+    }
+    else
     {
-            const wordRe fieldNameRe = fieldSelection_[fieldI];
+        names = mesh_.thisDb().names<GeoField>();
-            IOobjectList fieldIO = fieldObjects.lookupRe(fieldNameRe);
+    }
 
-            forAllConstIter(IOobjectList, fieldIO, iter)
+    labelList nameIDs(findStrings(fieldSelection_, names));
+
+    wordHashSet fieldNames(wordList(names, nameIDs));
+
+    forAllConstIter(wordHashSet, fieldNames, iter)
     {
-                const word& fieldName = iter()->name();
+        const word& fieldName = iter.key();
 
+        if ((Pstream::master()) && verbose_)
+        {
+            Pout<< "sampleAndWrite: " << fieldName << endl;
+        }
+
+        if (loadFromFiles_)
+        {
             const GeoField fld
             (
                 IOobject
@@ -193,35 +206,10 @@ void Foam::sampledSurfaces::sampleAndWrite(const IOobjectList& objects)
                 mesh_
             );
 
-                if ((Pstream::master()) && verbose_)
-                {
-                    Pout<< "sampleAndWrite: " << fieldName << endl;
-                }
-
             sampleAndWrite(fld);
         }
-        }
-    }
         else
         {
-        forAll(fieldSelection_, fieldI)
-        {
-            const wordRe& fieldNameRe = fieldSelection_[fieldI];
-
-            const wordList dbFields
-            (
-                mesh_.thisDb().foundObjectRe<GeoField>(fieldNameRe)
-            );
-
-            forAll(dbFields, i)
-            {
-                const word& fieldName = dbFields[i];
-
-                if ((Pstream::master()) && verbose_)
-                {
-                    Pout<< "sampleAndWrite: " << fieldName << endl;
-                }
-
             sampleAndWrite
             (
                 mesh_.thisDb().lookupObject<GeoField>(fieldName)
@@ -229,7 +217,6 @@ void Foam::sampledSurfaces::sampleAndWrite(const IOobjectList& objects)
         }
     }
 }
-}
 
 
 // ************************************************************************* //

src/thermophysicalModels/radiationModels/radiationModel/P1/P1.C

@@ -241,10 +241,14 @@ void Foam::radiation::P1::calculate()
     );
 
     // Calculate radiative heat flux on boundaries.
-    forAll(mesh_.boundaryMesh(), patchI)
+    forAll(mesh_.boundaryMesh(), patchi)
     {
-        Qr_.boundaryField()[patchI] =
+        if (!G_.boundaryField()[patchi].coupled())
-            -gamma.boundaryField()[patchI]*G_.boundaryField()[patchI].snGrad();
+        {
+            Qr_.boundaryField()[patchi] =
+                -gamma.boundaryField()[patchi]
+                *G_.boundaryField()[patchi].snGrad();
+        }
     }
 }
 

src/thermophysicalModels/solidChemistryModel/basicSolidChemistryModel/basicSolidChemistryModelNew.C

@@ -53,11 +53,12 @@ New
 
     Info<< "Selecting chemistry type " << chemistryTypeDict << endl;
 
-    const int nCmpt = 12;
+    const int nCmpt = 13;
     const char* cmptNames[nCmpt] =
     {
         "chemistrySolver",
         "chemistryThermo",
+        "baseChemistry",
         "transport",
         "thermo",
         "equationOfState",
@@ -107,8 +108,9 @@ New
     word chemistryTypeName
     (
         word(chemistryTypeDict.lookup("chemistrySolver")) + '<'
-      + word(chemistryTypeDict.lookup("chemistryThermo")) + ','
+      + word(chemistryTypeDict.lookup("chemistryThermo")) + '<'
-      + solidThermoTypeName + ',' + gasThermoTypeName + ">"
+      + typeName + ','
+      + solidThermoTypeName + ',' + gasThermoTypeName + ">>"
     );
 
     Info<< "chemistryTypeName " << chemistryTypeName << endl;

src/thermophysicalModels/solidChemistryModel/makeSolidChemistryModel.H

@@ -48,7 +48,8 @@ namespace Foam
     defineTemplateTypeNameAndDebugWithName                                    \
     (                                                                         \
         sChemistryl##Comp##SThermo,                                           \
-        (#sChemistry"<"#Comp"," + SThermo::typeName() + ">").c_str(),         \
+        (word(sChemistry::typeName_()) + "<"#Comp"," + SThermo::typeName() +  \
+        + ">").c_str(),                                                       \
         0                                                                     \
     );                                                                        \
                                                                               \

src/thermophysicalModels/solidChemistryModel/solidChemistrySolver/makeSolidChemistrySolverType.H

@@ -47,8 +47,8 @@ namespace Foam
     defineTemplateTypeNameAndDebugWithName                                    \
     (                                                                         \
         SS##Schem##Comp##SThermo##GThermo,                                    \
-        (#SS"<"#Schem"<"#Comp"," + SThermo::typeName() + ","                  \
+        (#SS"<" + word(Schem::typeName_()) +"<"#Comp"," + SThermo::typeName() \
-      + GThermo::typeName() + ">>").c_str(),                                  \
+        + ","  + GThermo::typeName() + ">>").c_str(),                         \
         0                                                                     \
     );                                                                        \
                                                                               \

src/transportModels/incompressible/incompressibleTwoPhaseMixture/incompressibleTwoPhaseMixture.C

@@ -175,8 +175,14 @@ bool Foam::incompressibleTwoPhaseMixture::read()
     {
         if
         (
-            nuModel1_().read(subDict(phase1Name_))
+            nuModel1_().read
-         && nuModel2_().read(subDict(phase2Name_))
+            (
+                subDict(phase1Name_ == "1" ? "phase1": phase1Name_)
+            )
+         && nuModel2_().read
+            (
+                subDict(phase2Name_ == "2" ? "phase2": phase2Name_)
+            )
         )
         {
             nuModel1_->viscosityProperties().lookup("rho") >> rho1_;

src/triSurface/triSurface/triSurface.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -175,6 +175,7 @@ class triSurface
         void writeDXGeometry(const bool, Ostream&) const;
         void writeDXTrailer(Ostream&) const;
 
+
     // Static private functions
 
         //- Convert faces to labelledTri. All get same region.
@@ -203,6 +204,7 @@ class triSurface
         //- read non-comment line
         static string getLineNoComment(IFstream&);
 
+
 protected:
 
     // Protected Member Functions
@@ -219,6 +221,7 @@ protected:
             return static_cast<List<Face>&>(*this);
         }
 
+
 public:
 
     // Public typedefs
@@ -287,7 +290,6 @@ public:
         triSurface(const triSurface&);
 
 
-
     //- Destructor
     ~triSurface();
 
@@ -324,6 +326,7 @@ public:
             //  If >2 neighbours: undetermined.
             const labelList& edgeOwner() const;
 
+
         // Edit
 
             //- Move points
@@ -381,6 +384,7 @@ public:
                 labelList& faceMap
             ) const;
 
+
         // Write
 
             //- Write to Ostream in simple FOAM format

src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/alphatWallFunctions/alphatJayatillekeWallFunction/alphatJayatillekeWallFunctionFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -122,7 +122,7 @@ class alphatJayatillekeWallFunctionFvPatchScalarField
 public:
 
     //- Runtime type information
-    TypeName("alphatJayatillekeWallFunction");
+    TypeName("compressible::alphatJayatillekeWallFunction");
 
 
     // Constructors

src/turbulenceModels/compressible/RAS/derivedFvPatchFields/wallFunctions/alphatWallFunctions/alphatWallFunction/alphatWallFunctionFvPatchScalarField.H

@@ -2,7 +2,7 @@
   =========                 |
   \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
    \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 2011-2012 OpenFOAM Foundation
+    \\  /    A nd           | Copyright (C) 2011-2013 OpenFOAM Foundation
      \\/     M anipulation  |
 -------------------------------------------------------------------------------
 License
@@ -105,7 +105,7 @@ class alphatWallFunctionFvPatchScalarField
 public:
 
     //- Runtime type information
-    TypeName("alphatWallFunction");
+    TypeName("compressible::alphatWallFunction");
 
 
     // Constructors