From a372dbfb5391b7c620240d026de7e2ed4ffb356c Mon Sep 17 00:00:00 2001
From: andy <andy>
Date: Thu, 21 Jul 2011 17:48:20 +0100
Subject: [PATCH] ENH: Further updates to AMI functionality

---
 src/AMIInterpolation/AMIInterpolation.C       | 214 +++++++++++---
 src/AMIInterpolation/AMIInterpolation.H       | 172 ++++++++----
 src/AMIInterpolation/AMIInterpolationI.H      |  58 ++++
 src/AMIInterpolation/Make/files               |   1 +
 src/AMIInterpolation/Make/options             |   7 +-
 .../faceAreaIntersect/faceAreaIntersect.C     | 120 ++++----
 .../faceAreaIntersect/faceAreaIntersect.H     |  17 +-
 .../cyclicAMIFvPatch/cyclicAMIFvPatch.C       |   7 +-
 .../cyclicAMIFvPatch/cyclicAMIFvPatch.H       |   5 -
 .../cyclicAMIPointPatch/cyclicAMIPointPatch.C |  99 +++++++
 .../cyclicAMIPointPatch/cyclicAMIPointPatch.H | 121 ++++++++
 .../cyclicAMIPolyPatch/cyclicAMIPolyPatch.C   | 262 +++++++++++-------
 .../cyclicAMIPolyPatch/cyclicAMIPolyPatch.H   | 100 ++++---
 .../cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H  |  67 +++++
 .../cyclicAMIPolyPatchTemplates.C             |  64 +++++
 src/OpenFOAM/meshes/boundBox/boundBox.C       |  12 +-
 16 files changed, 1000 insertions(+), 326 deletions(-)
 create mode 100644 src/AMIInterpolation/AMIInterpolationI.H
 create mode 100644 src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.C
 create mode 100644 src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.H
 create mode 100644 src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H
 create mode 100644 src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchTemplates.C

diff --git a/src/AMIInterpolation/AMIInterpolation.C b/src/AMIInterpolation/AMIInterpolation.C
index 0e5c6348e5..99fae0b857 100644
--- a/src/AMIInterpolation/AMIInterpolation.C
+++ b/src/AMIInterpolation/AMIInterpolation.C
@@ -24,13 +24,14 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "AMIInterpolation.H"
-#include "faceAreaIntersect.H"
 #include "Random.H"
 #include "treeDataPrimitivePatch.H"
 #include "indexedOctree.H"
 #include "primitivePatch.H"
 #include "meshTools.H"
 
+#include "vtkSurfaceWriter.H"
+
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 template<class SourcePatch, class TargetPatch>
@@ -95,9 +96,9 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkPatches
     const scalar maxBoundsError = 0.05;
 
     // sanity checks
-    boundBox bbSrc(srcPatch.points());
-    boundBox bbTgt(tgtPatch.points());
-    boundBox bbSurf(srcPatch.points());
+    boundBox bbSrc(srcPatch.points(), srcPatch.meshPoints());
+    boundBox bbTgt(tgtPatch.points(), tgtPatch.meshPoints());
+    boundBox bbSurf(bbTgt); // USE POINTS OF SURFACE!!!!!!!!!!!!!!!!!!
 
 
     // projection surface bounds - check against bounds of source and target
@@ -132,7 +133,6 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkPatches
 
 
     // check bounds of source and target
-
     bbTgt.inflate(maxBoundsError);
 
     if (!bbTgt.contains(bbSrc))
@@ -144,8 +144,11 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::checkPatches
                 "const primitivePatch&, "
                 "const primitivePatch&"
             ")"
-        )   << "Source and target patch bounding boxes are not similar"
-            << endl;
+        )   << "Source and target patch bounding boxes are not similar" << nl
+            << "    src span : " << bbSrc.span() << nl
+            << "    tgt span : " << bbTgt.span() << nl
+            << "    source: " << bbSrc << nl
+            << "    target: " << bbTgt << endl;
     }
 }
 
@@ -157,6 +160,16 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::projectPointsToSurface
     pointField& pts
 ) const
 {
+    if (!projectPoints_)
+    {
+        return;
+    }
+
+    if (debug)
+    {
+        Info<< "AMI: projecting points to surface" << endl;
+    }
+
     List<pointIndexHit> nearInfo;
 
     surf.findNearest(pts, scalarField(pts.size(), GREAT), nearInfo);
@@ -313,7 +326,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::setNextFaces
     const boolList& mapFlag,
     labelList& seedFaces,
     const DynamicList<label>& visitedFaces
-) const
+)
 {
 //    const labelList& srcNbrFaces = srcPatch0.pointFaces()[srcFaceI];
     const labelList& srcNbrFaces = srcPatch0.faceFaces()[srcFaceI];
@@ -356,14 +369,24 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::setNextFaces
     else
     {
         // try to use existing seed
-        forAll(mapFlag, faceI)
+        bool foundNextSeed = false;
+        for (label faceI = startSeedI_; faceI < mapFlag.size(); faceI++)
         {
-            if (mapFlag[faceI] && seedFaces[faceI] != -1)
+            if (mapFlag[faceI])
             {
-                srcFaceI = faceI;
-                tgtFaceI = seedFaces[faceI];
+                if (!foundNextSeed)
+                {
+                    startSeedI_ = faceI + 1;
+                    foundNextSeed = true;
+                }
 
-                return;
+                if (seedFaces[faceI] != -1)
+                {
+                    srcFaceI = faceI;
+                    tgtFaceI = seedFaces[faceI];
+
+                    return;
+                }
             }
         }
 
@@ -374,10 +397,12 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::setNextFaces
                 << "target face" << endl;
         }
 
-        forAll(mapFlag, faceI)
+        for (label faceI = startSeedI_; faceI < mapFlag.size(); faceI++)
         {
             if (mapFlag[faceI])
             {
+                startSeedI_ = faceI + 1;
+
                 srcFaceI = faceI;
                 tgtFaceI = findTargetFace(srcFaceI, srcPatch0, tgtPatch0);
 
@@ -429,7 +454,7 @@ Foam::scalar Foam::AMIInterpolation<SourcePatch, TargetPatch>::interArea
     // crude resultant norm
     const vector n = 0.5*(tgt.normal(tgtPoints) - src.normal(srcPoints));
 
-    scalar area = inter.calc(src, tgt, n);
+    scalar area = inter.calc(src, tgt, n, triMode_);
 
     if ((debug > 1) && (area > 0))
     {
@@ -447,10 +472,11 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights
     const primitivePatch& patch,
     const pointField& points,
     const List<DynamicList<label> >& addr,
-    List<DynamicList<scalar> >& wght
+    List<DynamicList<scalar> >& wght,
+    const bool output
 )
 {
-    scalarList wghtSum(patch.size(), 0.0);
+    scalarField wghtSum(patch.size(), 0.0);
 
     // normalise weights by face areas
     forAll(wght, faceI)
@@ -464,8 +490,37 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::normaliseWeights
         }
     }
 
-    Info<< "Cumulative " << patchType <<  "source weights min/max = "
-        << min(wghtSum) << ", " << max(wghtSum) << endl;
+//    if (debug)
+    {
+        vtkSurfaceWriter writer;
+
+        writer.write
+        (
+            "VTK",
+            patchType,
+            patch.localPoints(),
+            patch.localFaces(),
+            "weights",
+            wghtSum,
+            false
+        );
+    }
+
+    if (output)
+    {
+        Info<< "AMI: Patch " << patchType << " weights min/max = "
+            << min(wghtSum) << ", " << max(wghtSum) << endl;
+    }
+
+
+    forAll(wght, faceI)
+    {
+        const DynamicList<label>& addressing = addr[faceI];
+        forAll(addressing, addrI)
+        {
+            wght[faceI][addrI] /= wghtSum[faceI];
+        }
+    }
 }
 
 
@@ -476,17 +531,29 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::AMIInterpolation
 (
     const SourcePatch& srcPatch,
     const TargetPatch& tgtPatch,
-    const searchableSurface& surf
+    const searchableSurface& surf,
+    const faceAreaIntersect::triangulationMode& triMode,
+    const bool projectPoints
 )
 :
     srcAddress_(),
     srcWeights_(),
     tgtAddress_(),
-    tgtWeights_()
+    tgtWeights_(),
+    startSeedI_(0),
+    triMode_(triMode),
+    projectPoints_(projectPoints)
 {
-    checkPatches(srcPatch, tgtPatch);
+    if (srcPatch.size() && tgtPatch.size())
+    {
+        checkPatches(srcPatch, tgtPatch);
 
-    calcAddressing(srcPatch, tgtPatch, surf);
+        Info<< "AMI: Creating interpolation addressing and weights" << nl
+            << "    Source faces = " << srcPatch.size() << nl
+            << "    Target faces = " << tgtPatch.size() << endl;
+
+        calcAddressing(srcPatch, tgtPatch, surf);
+    }
 }
 
 
@@ -557,7 +624,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
 
     if (debug)
     {
-        OFstream os("osTgtPoints.obj");
+        OFstream os("amiTgtPoints.obj");
         forAll(tgtPoints, i)
         {
             meshTools::writeOBJ(os, tgtPoints[i]);
@@ -565,11 +632,6 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
     }
 
 
-    if (debug)
-    {
-        Info<< "AMI: projecting points to surface" << endl;
-    }
-
     // map source and target patches onto projection surface
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
     projectPointsToSurface(surf, srcPoints);
@@ -589,7 +651,7 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
 
     // construct weights and addressing
     // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-    label facesRemaining = srcPatch0.size();
+    label nFacesRemaining = srcPatch0.size();
 
     // list of tgt face neighbour faces
     DynamicList<label> nbrFaces(10);
@@ -598,11 +660,13 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
     DynamicList<label> visitedFaces(10);
 
     // list to keep track of tgt faces used to seed src faces
-    labelList seedFaces(facesRemaining, -1);
+    labelList seedFaces(nFacesRemaining, -1);
     seedFaces[srcFaceI] = tgtFaceI;
 
     // list to keep track of whether src face can be mapped
-    boolList mapFlag(facesRemaining, true);
+    boolList mapFlag(nFacesRemaining, true);
+
+    label nFail = 0;
 
     do
     {
@@ -613,6 +677,8 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
         nbrFaces.append(tgtFaceI);
         appendNbrFaces(tgtFaceI, tgtPatch0, visitedFaces, nbrFaces);
 
+        bool faceProcessed = false;
+
         do
         {
             // process new target face
@@ -630,16 +696,30 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
                 tgtWght[tgtFaceI].append(area);
 
                 appendNbrFaces(tgtFaceI, tgtPatch0, visitedFaces, nbrFaces);
+
+                faceProcessed = true;
             }
 
         } while (nbrFaces.size() > 0);
 
-        mapFlag[srcFaceI] = false;
+        if (faceProcessed)
+        {
+            mapFlag[srcFaceI] = false;
 
-        facesRemaining--;
+            nFacesRemaining--;
+        }
+        else
+        {
+            nFail++;
+
+            if (nFail > 1000)
+            {
+                Info<< "nFail = " << nFail << endl;
+            }
+        }
 
         // choose new src face from current src face neighbour
-        if (facesRemaining > 0)
+        if (nFacesRemaining > 0)
         {
             setNextFaces
             (
@@ -652,12 +732,12 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
                 visitedFaces
             );
         }
-    } while (facesRemaining > 0);
+    } while (nFacesRemaining > 0);
 
 
     // weights normalisation
-    normaliseWeights("source", srcPatch0, srcPoints, srcAddr, srcWght);
-    normaliseWeights("target", tgtPatch0, tgtPoints, tgtAddr, tgtWght);
+    normaliseWeights("source", srcPatch0, srcPoints, srcAddr, srcWght, true);
+    normaliseWeights("target", tgtPatch0, tgtPoints, tgtAddr, tgtWght, true);
 
     // transfer data to persistent storage
     srcAddress_.setSize(srcAddr.size());
@@ -675,6 +755,32 @@ void Foam::AMIInterpolation<SourcePatch, TargetPatch>::calcAddressing
         tgtAddress_[i].transfer(tgtAddr[i]);
         tgtWeights_[i].transfer(tgtWght[i]);
     }
+
+// ADD PARALLEL BITS
+/*
+    if (Pstream::parRun())
+    {
+        // Source fields - addesses and weights
+        const label nSrc = srcAddess.size();
+        const globalIndex globalSrcFaces(nSrc);
+
+        List<Map<label> > compactMapSrcAddress;
+        labelListList transSrcAddress;
+
+        mapDistribute srcMapDist
+        (
+            globalSrcFaces,
+            srcAddress_,
+            globalIndexAndTransform(mesh), // mesh not available!!!!!!!!
+            transSrcAddress,
+            compactMapSrcAddress
+        );
+
+
+    }
+*/
+    // reset starting seed
+    startSeedI_ = 0;
 }
 
 
@@ -790,4 +896,36 @@ Foam::AMIInterpolation<SourcePatch, TargetPatch>::interpolateToTarget
 }
 
 
+template<class SourcePatch, class TargetPatch>
+void Foam::AMIInterpolation<SourcePatch, TargetPatch>::writeFaceConnectivity
+(
+    const primitivePatch& srcPatch,
+    const primitivePatch& tgtPatch
+)
+const
+{
+    OFstream os("faceConnectivity.obj");
+
+    label ptI = 1;
+
+    forAll(srcAddress_, i)
+    {
+        const labelList& addr = srcAddress_[i];
+        const point& srcPt = srcPatch.faceCentres()[i];
+        forAll(addr, j)
+        {
+            label tgtPtI = addr[j];
+            const point& tgtPt = tgtPatch.faceCentres()[tgtPtI];
+
+            meshTools::writeOBJ(os, srcPt);
+            meshTools::writeOBJ(os, tgtPt);
+
+            os  << "l " << ptI << " " << ptI + 1 << endl;
+
+            ptI += 2;
+        }
+    }
+}
+
+
 // ************************************************************************* //
diff --git a/src/AMIInterpolation/AMIInterpolation.H b/src/AMIInterpolation/AMIInterpolation.H
index 246c1660a6..8dd8b72c1b 100644
--- a/src/AMIInterpolation/AMIInterpolation.H
+++ b/src/AMIInterpolation/AMIInterpolation.H
@@ -26,7 +26,14 @@ Class
 
 Description
     Interpolation class dealing with transfer of data between two
-    primitive patches with an arbitrary mesh interface (AMI)
+    primitive patches with an arbitrary mesh interface (AMI).
+
+    Based on the algorithm given in:
+
+        Conservative interpolation between volume meshes by local Galerkin
+        projection, Farrell PE and Maddison JR, 2011, Comput. Methods Appl.
+        Mech Engrg, Volume 200, Issues 1-4, pp 89-100
+
 
 SourceFiles
     AMIInterpolation.C
@@ -42,6 +49,7 @@ SourceFiles
 #include "searchableSurface.H"
 #include "boolList.H"
 #include "primitivePatch.H"
+#include "faceAreaIntersect.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -84,6 +92,16 @@ class AMIInterpolation
             scalarListList tgtWeights_;
 
 
+        //- Starting face seed index
+        label startSeedI_;
+
+        //- Face triangulation mode
+        const faceAreaIntersect::triangulationMode triMode_;
+
+        //- Flag to indicate whether or not to project points
+        const bool projectPoints_;
+
+
     // Private Member Functions
 
         //- Disallow default bitwise copy construct
@@ -92,70 +110,94 @@ class AMIInterpolation
         //- Disallow default bitwise assignment
         void operator=(const AMIInterpolation&);
 
-        //- Write triangle intersection to OBJ file
-        void writeIntersectionOBJ
-        (
-            const scalar area,
-            const face& f1,
-            const face& f2,
-            const pointField& f1Points,
-            const pointField& f2Points
-        ) const;
 
-        void checkPatches
-        (
-            const primitivePatch& srcPatch,
-            const primitivePatch& tgtPatch
-        );
+        // Helper functions
 
-        void projectPointsToSurface
-        (
-            const searchableSurface& surf,
-            pointField& pts
-        ) const;
+            //- Write triangle intersection to OBJ file
+            void writeIntersectionOBJ
+            (
+                const scalar area,
+                const face& f1,
+                const face& f2,
+                const pointField& f1Points,
+                const pointField& f2Points
+            ) const;
 
-        label findTargetFace
-        (
-            const label srcFaceI,
-            const primitivePatch& srcPatch,
-            const primitivePatch& tgtPatch
-        ) const;
+            //- Check that patches are valid
+            void checkPatches
+            (
+                const primitivePatch& srcPatch,
+                const primitivePatch& tgtPatch
+            );
 
-        void appendNbrFaces
-        (
-            const label faceI,
-            const primitivePatch& patch,
-            const DynamicList<label>& visitedFaces,
-            DynamicList<label>& faceIDs
-        ) const;
 
-        void setNextFaces
-        (
-            label& srcFaceI,
-            label& tgtFaceI,
-            const primitivePatch& srcPatch0,
-            const primitivePatch& tgtPatch0,
-            const boolList& mapFlag,
-            labelList& seedFaces,
-            const DynamicList<label>& visitedFaces
-        ) const;
+        // Initialisation
 
-        scalar interArea
-        (
-            const label srcFaceI,
-            const label tgtFaceI,
-            const primitivePatch& srcPatch,
-            const primitivePatch& tgtPatch
-        ) const;
+            //- Project points to surface
+            void projectPointsToSurface
+            (
+                const searchableSurface& surf,
+                pointField& pts
+            ) const;
 
-        void normaliseWeights
-        (
-            const word& patchType,
-            const primitivePatch& patch,
-            const pointField& points,
-            const List<DynamicList<label> >& addr,
-            List<DynamicList<scalar> >& wght
-        );
+
+        // Marching front
+
+            //- Find face on target patch that overlaps source face
+            label findTargetFace
+            (
+                const label srcFaceI,
+                const primitivePatch& srcPatch,
+                const primitivePatch& tgtPatch
+            ) const;
+
+            //- Add faces neighbouring faceI to the ID list
+            void appendNbrFaces
+            (
+                const label faceI,
+                const primitivePatch& patch,
+                const DynamicList<label>& visitedFaces,
+                DynamicList<label>& faceIDs
+            ) const;
+
+            //- Set the source and target seed faces
+            void setNextFaces
+            (
+                label& srcFaceI,
+                label& tgtFaceI,
+                const primitivePatch& srcPatch0,
+                const primitivePatch& tgtPatch0,
+                const boolList& mapFlag,
+                labelList& seedFaces,
+                const DynamicList<label>& visitedFaces
+            );
+
+
+        // Evaluation
+
+            //- Area of intersection between source and target faces
+            scalar interArea
+            (
+                const label srcFaceI,
+                const label tgtFaceI,
+                const primitivePatch& srcPatch,
+                const primitivePatch& tgtPatch
+            ) const;
+
+            //- Normalise the (area) weights - suppresses numerical error in
+            //  weights calculation
+            //  NOTE: if area weights are incorrect by 'a significant amount'
+            //     normalisation may stabilise the solution, but will introduce
+            //     numerical error!
+            void normaliseWeights
+            (
+                const word& patchType,
+                const primitivePatch& patch,
+                const pointField& points,
+                const List<DynamicList<label> >& addr,
+                List<DynamicList<scalar> >& wght,
+                const bool output
+            );
 
 
 public:
@@ -167,7 +209,9 @@ public:
         (
             const SourcePatch& srcPatch,
             const TargetPatch& tgtPatch,
-            const searchableSurface& surf
+            const searchableSurface& surf,
+            const faceAreaIntersect::triangulationMode& triMode,
+            const bool projectPoints
         );
 
 
@@ -234,6 +278,16 @@ public:
             (
                 const tmp<Field<Type> >& tFld
             ) const;
+
+
+        // Checks
+
+            //- Write face connectivity as OBJ file
+            void writeFaceConnectivity
+            (
+                const primitivePatch& srcPatch,
+                const primitivePatch& tgtPatch
+            ) const;
 };
 
 
diff --git a/src/AMIInterpolation/AMIInterpolationI.H b/src/AMIInterpolation/AMIInterpolationI.H
new file mode 100644
index 0000000000..a38c83ebaf
--- /dev/null
+++ b/src/AMIInterpolation/AMIInterpolationI.H
@@ -0,0 +1,58 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+template<class SourcePatch, class TargetPatch>
+inline const Foam::labelListList&
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcAddress()
+{
+    return srcAddress_;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+inline const Foam::scalarListList&
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::srcWeights()
+{
+    return srcWeights_;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+inline const Foam::labelListList&
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtAddress()
+{
+    return tgtAddress_;
+}
+
+
+template<class SourcePatch, class TargetPatch>
+inline const Foam::scalarListList&
+Foam::AMIInterpolation<SourcePatch, TargetPatch>::tgtWeights()
+{
+    return tgtWeights_;
+}
+
+
+// ************************************************************************* //
diff --git a/src/AMIInterpolation/Make/files b/src/AMIInterpolation/Make/files
index 7a47cad93c..1c4f231be8 100644
--- a/src/AMIInterpolation/Make/files
+++ b/src/AMIInterpolation/Make/files
@@ -8,5 +8,6 @@ patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
 patches/cyclicAMIFvPatch/cyclicAMIFvPatch.C
 patches/cyclicAMIFvPatchField/cyclicAMIFvPatchFields.C
 patches/cyclicAMIFvsPatchField/cyclicAMIFvsPatchFields.C
+patches/cyclicAMIPointPatch/cyclicAMIPointPatch.C
 
 LIB = $(FOAM_USER_LIBBIN)/libAMIInterpolation
diff --git a/src/AMIInterpolation/Make/options b/src/AMIInterpolation/Make/options
index 431de9fbf4..b1bf11173c 100644
--- a/src/AMIInterpolation/Make/options
+++ b/src/AMIInterpolation/Make/options
@@ -1,8 +1,11 @@
 EXE_INC = \
     -DFULLDEBUG -g -O0 \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/meshTools/lnInclude
+    -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/sampling/lnInclude
+
 
 LIB_LIBS = \
     -lfiniteVolume \
-    -lmeshTools
+    -lmeshTools \
+    -lsampling
diff --git a/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C b/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C
index 174fcdd461..e6294ba63e 100644
--- a/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C
+++ b/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.C
@@ -32,9 +32,12 @@ void Foam::faceAreaIntersect::triSliceWithPlane
     const triPoints& tri,
     const plane& p,
     FixedList<triPoints, 10>& tris,
-    label& nTris
+    label& nTris,
+    const scalar len
 )
 {
+    const scalar matchTol = 1e-4;
+
     // distance to cutting plane
     FixedList<scalar, 3> d;
 
@@ -43,22 +46,16 @@ void Foam::faceAreaIntersect::triSliceWithPlane
     label nPos = 0;
     label posI = -1;
     label negI = -1;
-    bool coPlanar;
     forAll(tri, i)
     {
         d[i] = ((tri[i] - p.refPoint()) & p.normal());
 
-        if (mag(d[i]) < 1e-10)
+        if (mag(d[i]) < matchTol*len)
         {
-            coPlanar = true;
+            d[i] = 0.0;
             nCoPlanar++;
         }
         else
-        {
-            coPlanar = false;
-        }
-
-        if (!coPlanar)
         {
             if (d[i] > 0)
             {
@@ -133,44 +130,48 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
 )
 {
     // Work storage
-    FixedList<triPoints, 10> cutTris;
-    label nCutTris = 0;
+    FixedList<triPoints, 10> workTris1;
+    label nWorkTris1 = 0;
 
-    FixedList<triPoints, 10> tris;
-    label nTris = 0;
+    FixedList<triPoints, 10> workTris2;
+    label nWorkTris2 = 0;
 
     // cut source triangle with all inwards pointing faces of target triangle
-    // - triangles in cutTris are inside target triangle
+    // - triangles in workTris1 are inside target triangle
+
+    scalar t = sqrt(triArea(src));
 
     // edge 0
     {
         // cut triangle src with plane and put resulting sub-triangles in
-        // cutTris list
+        // workTris1 list
 
-        plane pl0(tgt[0], tgt[1], tgt[1] + n);
-        triSliceWithPlane(src, pl0, cutTris, nCutTris);
+        scalar s = mag(tgt[1] - tgt[0]);
+        plane pl0(tgt[0], tgt[1], tgt[1] + s*n);
+        triSliceWithPlane(src, pl0, workTris1, nWorkTris1, t);
     }
 
-    if (nCutTris == 0)
+    if (nWorkTris1 == 0)
     {
         return 0.0;
     }
 
     // edge1
     {
-        // cut cutTris with plane and put resulting sub-triangles in
-        // tris list (re-use tris storage)
+        // cut workTris1 with plane and put resulting sub-triangles in
+        // workTris2 list (re-use tris storage)
 
-        plane pl1(tgt[1], tgt[2], tgt[2] + n);
+        scalar s = mag(tgt[2] - tgt[1]);
+        plane pl1(tgt[1], tgt[2], tgt[2] + s*n);
 
-        nTris = 0;
+        nWorkTris2 = 0;
 
-        for (label i = 0; i < nCutTris; i++)
+        for (label i = 0; i < nWorkTris1; i++)
         {
-            triSliceWithPlane(cutTris[i], pl1, tris, nTris);
+            triSliceWithPlane(workTris1[i], pl1, workTris2, nWorkTris2, t);
         }
 
-        if (nTris == 0)
+        if (nWorkTris2 == 0)
         {
             return 0.0;
         }
@@ -178,19 +179,20 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
 
     // edge2
     {
-        // cut tris with plane and put resulting sub-triangles in
-        // cutTris list (re-use cutTris storage)
+        // cut workTris2 with plane and put resulting sub-triangles in
+        // workTris1 list (re-use workTris1 storage)
 
-        plane pl2(tgt[2], tgt[0], tgt[0] + n);
+        scalar s = mag(tgt[2] - tgt[0]);
+        plane pl2(tgt[2], tgt[0], tgt[0] + s*n);
 
-        nCutTris = 0;
+        nWorkTris1 = 0;
 
-        for (label i = 0; i < nTris; i++)
+        for (label i = 0; i < nWorkTris2; i++)
         {
-            triSliceWithPlane(tris[i], pl2, cutTris, nCutTris);
+            triSliceWithPlane(workTris2[i], pl2, workTris1, nWorkTris1, t);
         }
 
-        if (nCutTris == 0)
+        if (nWorkTris1 == 0)
         {
             return 0.0;
         }
@@ -198,10 +200,9 @@ Foam::scalar Foam::faceAreaIntersect::triangleIntersect
         {
             // calculate area of sub-triangles
             scalar area = 0.0;
-            for (label i = 0; i < nCutTris; i++)
+            for (label i = 0; i < nWorkTris1; i++)
             {
-                const triPoints& t = cutTris[i];
-                area += mag(0.5*((t[1] - t[0])^(t[2] - t[0])));
+                area += triArea(workTris1[i]);
             }
 
             return area;
@@ -229,26 +230,45 @@ Foam::scalar Foam::faceAreaIntersect::calc
 (
     const face& faceA,
     const face& faceB,
-    const vector& n
+    const vector& n,
+    const triangulationMode& triMode
 )
 {
     // split faces into triangles
     DynamicList<face> trisA;
     DynamicList<face> trisB;
 
-//    if (useTriangleFan)
-//    {
-//        triangleFan(faceA, trisA);
-//        triangleFan(faceB, trisB);
-//    }
-//    else
-//    {
-//        faceA.triangles(pointsA_, trisA);
-//        faceB.triangles(pointsB_, trisB);
-//    }
+    switch (triMode)
+    {
+        case tmFan:
+        {
+            triangleFan(faceA, trisA);
+            triangleFan(faceB, trisB);
 
-    faceA.triangles(pointsA_, trisA);
-    faceB.triangles(pointsB_, trisB);
+            break;
+        }
+        case tmMesh:
+        {
+            faceA.triangles(pointsA_, trisA);
+            faceB.triangles(pointsB_, trisB);
+
+            break;
+        }
+        default:
+        {
+            FatalErrorIn
+            (
+                "Foam::scalar Foam::faceAreaIntersect::calc"
+                "("
+                    "const face&, "
+                    "const face&, "
+                    "const vector&, "
+                    "const triangulationMode&"
+                ")"
+            )   << "Unknown triangulation mode enumeration"
+                << abort(FatalError);
+        }
+    }
 
     // intersect triangles
     scalar totalArea = 0.0;
@@ -256,13 +276,13 @@ Foam::scalar Foam::faceAreaIntersect::calc
     {
         triPoints tpA = getTriPoints(pointsA_, trisA[tA], false);
 
-        if (triArea(tpA) > ROOTVSMALL)
+//        if (triArea(tpA) > ROOTVSMALL)
         {
             forAll(trisB, tB)
             {
                 triPoints tpB = getTriPoints(pointsB_, trisB[tB], true);
 
-                if (triArea(tpB) > ROOTVSMALL)
+//                if (triArea(tpB) > ROOTVSMALL)
                 {
                     totalArea += triangleIntersect(tpA, tpB, n);
                 }
diff --git a/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H b/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H
index 77638e51dd..fd464b8458 100644
--- a/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H
+++ b/src/AMIInterpolation/faceAreaIntersect/faceAreaIntersect.H
@@ -53,8 +53,19 @@ namespace Foam
 
 class faceAreaIntersect
 {
+public:
+
     typedef FixedList<point, 3> triPoints;
 
+    enum triangulationMode
+    {
+        tmFan,
+        tmMesh
+    };
+
+
+private:
+
     // Private data
 
         //- Reference to the points of sideA
@@ -110,7 +121,8 @@ class faceAreaIntersect
             const triPoints& tri,
             const plane& p,
             FixedList<triPoints, 10>& tris,
-            label& nTris
+            label& nTris,
+            const scalar len
         );
 
         //- Return area of intersection of triangles src and tgt
@@ -141,7 +153,8 @@ public:
         (
             const face& faceA,
             const face& faceB,
-            const vector& n
+            const vector& n,
+            const triangulationMode& triMode
         );
 };
 
diff --git a/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.C b/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.C
index d000c52c31..7c11307f6d 100644
--- a/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.C
+++ b/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.C
@@ -47,7 +47,7 @@ void Foam::cyclicAMIFvPatch::makeWeights(scalarField& w) const
 
     const scalarField nbrDeltas
     (
-        interpolateToSource(nbrPatch.nf() & nbrPatchDelta())
+        interpolateToSource(nbrPatch.nf() & nbrPatch.fvPatch::delta())
     );
 
     forAll(deltas, faceI)
@@ -68,7 +68,7 @@ void Foam::cyclicAMIFvPatch::makeDeltaCoeffs(scalarField& dc) const
 
     const scalarField nbrDeltas
     (
-        interpolateToSource(nbrPatch.nf() & nbrPatchDelta())
+        interpolateToSource(nbrPatch.nf() & nbrPatch.fvPatch::delta())
     );
 
     forAll(deltas, faceI)
@@ -85,7 +85,8 @@ Foam::tmp<Foam::vectorField> Foam::cyclicAMIFvPatch::delta() const
 {
     const vectorField patchD(fvPatch::delta());
 
-    const vectorField& nbrPatchD(interpolateToSource(nbrPatchDelta()));
+    const cyclicAMIFvPatch& nbrPatch = neighbFvPatch();
+    const vectorField nbrPatchD(interpolateToSource(nbrPatch.fvPatch::delta()));
 
     tmp<vectorField> tpdv(new vectorField(patchD.size()));
     vectorField& pdv = tpdv();
diff --git a/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.H b/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.H
index 5ccdd01d89..c0d32162e8 100644
--- a/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.H
+++ b/src/AMIInterpolation/patches/cyclicAMIFvPatch/cyclicAMIFvPatch.H
@@ -158,11 +158,6 @@ public:
                 return cyclicAMIPolyPatch_.interpolateToSource(tFld);
             }
 
-            const vectorField& nbrPatchDelta() const
-            {
-                return cyclicAMIPolyPatch_.nbrPatchDelta();
-            }
-
 
         // Interface transfer functions
 
diff --git a/src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.C b/src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.C
new file mode 100644
index 0000000000..c41779078a
--- /dev/null
+++ b/src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.C
@@ -0,0 +1,99 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+#include "cyclicAMIPointPatch.H"
+#include "pointBoundaryMesh.H"
+#include "addToRunTimeSelectionTable.H"
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+namespace Foam
+{
+    defineTypeNameAndDebug(cyclicAMIPointPatch, 0);
+    addToRunTimeSelectionTable
+    (
+        facePointPatch,
+        cyclicAMIPointPatch,
+        polyPatch
+    );
+}
+
+
+// * * * * * * * * * * * * * Protected Member Functions  * * * * * * * * * * //
+
+void Foam::cyclicAMIPointPatch::initGeometry(PstreamBuffers&)
+{}
+
+
+void Foam::cyclicAMIPointPatch::calcGeometry(PstreamBuffers&)
+{}
+
+
+void Foam::cyclicAMIPointPatch::initMovePoints
+(
+    PstreamBuffers&,
+    const pointField&
+)
+{}
+
+
+void Foam::cyclicAMIPointPatch::movePoints(PstreamBuffers&, const pointField&)
+{}
+
+
+void Foam::cyclicAMIPointPatch::initUpdateMesh(PstreamBuffers& pBufs)
+{
+    facePointPatch::initUpdateMesh(pBufs);
+//    cyclicAMIPointPatch::initGeometry(pBufs);
+}
+
+
+void Foam::cyclicAMIPointPatch::updateMesh(PstreamBuffers& pBufs)
+{
+    facePointPatch::updateMesh(pBufs);
+//    cyclicAMIPointPatch::calcGeometry(pBufs);
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::cyclicAMIPointPatch::cyclicAMIPointPatch
+(
+    const polyPatch& patch,
+    const pointBoundaryMesh& bm
+)
+:
+    coupledFacePointPatch(patch, bm),
+    cyclicAMIPolyPatch_(refCast<const cyclicAMIPolyPatch>(patch))
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+Foam::cyclicAMIPointPatch::~cyclicAMIPointPatch()
+{}
+
+
+// ************************************************************************* //
diff --git a/src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.H b/src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.H
new file mode 100644
index 0000000000..ddfa71d7cd
--- /dev/null
+++ b/src/AMIInterpolation/patches/cyclicAMIPointPatch/cyclicAMIPointPatch.H
@@ -0,0 +1,121 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+Class
+    Foam::cyclicAMIPointPatch
+
+Description
+    Cyclic AMI point patch - place holder only
+
+SourceFiles
+    cyclicAMIPointPatch.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef cyclicAMIPointPatch_H
+#define cyclicAMIPointPatch_H
+
+#include "coupledFacePointPatch.H"
+#include "cyclicAMIPolyPatch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                      Class cyclicAMIPointPatch Declaration
+\*---------------------------------------------------------------------------*/
+
+class cyclicAMIPointPatch
+:
+    public coupledFacePointPatch
+{
+    // Private data
+
+        //- Local reference cast into the cyclic AMI patch
+        const cyclicAMIPolyPatch& cyclicAMIPolyPatch_;
+
+
+    // Private Member Functions
+
+        //- Disallow default construct as copy
+        cyclicAMIPointPatch(const cyclicAMIPointPatch&);
+
+        //- Disallow default assignment
+        void operator=(const cyclicAMIPointPatch&);
+
+
+protected:
+
+    // Protected Member Functions
+
+        //- Initialise the calculation of the patch geometry
+        virtual void initGeometry(PstreamBuffers&);
+
+        //- Calculate the patch geometry
+        virtual void calcGeometry(PstreamBuffers&);
+
+        //- Initialise the patches for moving points
+        virtual void initMovePoints(PstreamBuffers&, const pointField&);
+
+        //- Correct patches after moving points
+        virtual void movePoints(PstreamBuffers&, const pointField&);
+
+        //- Initialise the update of the patch topology
+        virtual void initUpdateMesh(PstreamBuffers&);
+
+        //- Update of the patch topology
+        virtual void updateMesh(PstreamBuffers&);
+
+
+public:
+
+    //- Runtime type information
+    TypeName(cyclicAMIPolyPatch::typeName_());
+
+
+    // Constructors
+
+        //- Construct from components
+        cyclicAMIPointPatch
+        (
+            const polyPatch& patch,
+            const pointBoundaryMesh& bm
+        );
+
+
+    //- Destructor
+    virtual ~cyclicAMIPointPatch();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //
diff --git a/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
index 7ae69325b7..ba2bf63630 100644
--- a/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
+++ b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.C
@@ -29,6 +29,7 @@ License
 #include "polyMesh.H"
 #include "Time.H"
 #include "addToRunTimeSelectionTable.H"
+#include "faceAreaIntersect.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -72,7 +73,6 @@ Foam::label Foam::cyclicAMIPolyPatch::findFaceMaxRadius
 
 // * * * * * * * * * * * Protecetd Member Functions  * * * * * * * * * * * * //
 
-
 void Foam::cyclicAMIPolyPatch::calcTransforms()
 {
     if (size())
@@ -121,7 +121,7 @@ void Foam::cyclicAMIPolyPatch::calcTransforms
             << " has transform type " << transformTypeNames[transform_]
             << ", neighbour patch " << nbrPatchName_
             << " has transform type "
-            << nbrPatch().transformTypeNames[transform_]
+            << nbrPatch().transformTypeNames[nbrPatch().transform_]
             << exit(FatalError);
     }
 
@@ -176,8 +176,7 @@ void Foam::cyclicAMIPolyPatch::calcTransforms
         }
         case TRANSLATIONAL:
         {
-            // Transform 0 points.
-
+            // Transform 0 points
             if (debug)
             {
                 Pout<< "cyclicAMIPolyPatch::calcTransforms :"
@@ -185,64 +184,33 @@ void Foam::cyclicAMIPolyPatch::calcTransforms
                     << endl;
             }
 
-                const_cast<tensorField&>(forwardT()).clear();
-                const_cast<tensorField&>(reverseT()).clear();
-                const_cast<vectorField&>(separation()) = vectorField
-                (
-                    1,
-                    separationVector_
-                );
-                const_cast<boolList&>(collocated()) = boolList(1, false);
+            const_cast<tensorField&>(forwardT()).clear();
+            const_cast<tensorField&>(reverseT()).clear();
+            const_cast<vectorField&>(separation()) = vectorField
+            (
+                1,
+                separationVector_
+            );
+            const_cast<boolList&>(collocated()) = boolList(1, false);
+
             break;
         }
         default:
         {
-/*
-            // Assumes that cyclic is planar. This is also the initial
-            // condition for patches without faces.
-
-            // Determine the face with max area on both halves. These
-            // two faces are used to determine the transformation tensors
-            label max0I = findMaxArea(pp0.points(), pp0);
-            vector n0 = pp0[max0I].normal(pp0.points());
-            n0 /= mag(n0) + VSMALL;
-
-            label max1I = findMaxArea(pp1.points(), pp1);
-            vector n1 = pp1[max1I].normal(pp1.points());
-            n1 /= mag(n1) + VSMALL;
-
-            if (mag(n0 & n1) < 1 - matchTolerance())
+            if (debug)
             {
-                if (debug)
-                {
-                    Pout<< "cyclicAMIPolyPatch::calcTransforms :"
-                        << " Detected rotation :"
-                        << " n0:" << n0 << " n1:" << n1 << endl;
-                }
-
-                // Rotation (around origin)
-                const tensor revT(rotationTensor(n0, -n1));
-                const_cast<tensorField&>(forwardT()) = tensorField(1, revT.T());
-                const_cast<tensorField&>(reverseT()) = tensorField(1, revT);
-                const_cast<vectorField&>(separation()).setSize(0);
-                const_cast<boolList&>(collocated()) = boolList(1, false);
+                Pout<< "Assuming cyclic AMI pairs are colocated" << endl;
             }
-            else
-            {
-                // Parallel translation
-                const_cast<tensorField&>(forwardT()).clear();
-                const_cast<tensorField&>(reverseT()).clear();
-                const_cast<vectorField&>(separation()) = vectorField
-                (
-                    1,
-                    separationVector_
-                );
-                const_cast<boolList&>(collocated()) = boolList(1, false);
 
-            }
-*/
+            const_cast<tensorField&>(forwardT()).clear();
+            const_cast<tensorField&>(reverseT()).clear();
+            const_cast<vectorField&>(separation()) = vectorField
+            (
+                1,
+                vector::zero
+            );
+            const_cast<boolList&>(collocated()) = boolList(1, true);
 
-            notImplemented("calcTransforms - unknown transform clause");
             break;
         }
     }
@@ -253,19 +221,49 @@ void Foam::cyclicAMIPolyPatch::calcTransforms
 }
 
 
+const Foam::searchableSurface& Foam::cyclicAMIPolyPatch::surf()
+{
+    if (projectionSurfaceType_ == "patch")
+    {
+        notImplemented("projectionSurfaceType_ == patch")
+    }
+    else
+    {
+        surfPtr_ =
+            searchableSurface::New
+            (
+                projectionSurfaceType_,
+                IOobject
+                (
+                    projectionName_,
+                    boundaryMesh().mesh().time().constant(),
+                    "triSurface",
+                    boundaryMesh().mesh(),
+                    IOobject::MUST_READ,
+                    IOobject::NO_WRITE
+                ),
+                dict_
+            );
+    }
+
+    return surfPtr_();
+}
+
+
 const Foam::AMIPatchToPatchInterpolation& Foam::cyclicAMIPolyPatch::AMI()
 {
-    if (!AMIPtr_)
+    if (!AMIPtr_.valid() && owner())
     {
         const polyPatch& nbr = nbrPatch();
         pointField nbrPoints = nbrPatch().localPoints();
 
+        if (debug)
+        {
+            OFstream os(name() + "_neighbourPatch-org.obj");
+            meshTools::writeOBJ(os, nbrPatch().localFaces(), nbrPoints);
+        }
 
-const word myName = nbrPatch().nbrPatch().nbrPatchName();
-
-OFstream osNorg(myName + "_neighbourPatch-org.obj");
-meshTools::writeOBJ(osNorg, nbrPatch().localFaces(), nbrPoints);
-
+        // transform neighbour patch to local system
         transformPosition(nbrPoints);
         primitivePatch nbrPatch0
         (
@@ -277,31 +275,30 @@ meshTools::writeOBJ(osNorg, nbrPatch().localFaces(), nbrPoints);
             nbrPoints
         );
 
-OFstream osN(myName + "_neighbourPatch-trans.obj");
-meshTools::writeOBJ(osN, nbrPatch0.localFaces(), nbrPoints);
+        if (debug)
+        {
+            OFstream osN(name() + "_neighbourPatch-trans.obj");
+            meshTools::writeOBJ(osN, nbrPatch0.localFaces(), nbrPoints);
 
-OFstream osO(myName + "_ownerPatch.obj");
-meshTools::writeOBJ(osO, this->localFaces(), localPoints());
+            OFstream osO(name() + "_ownerPatch.obj");
+            meshTools::writeOBJ(osO, this->localFaces(), localPoints());
+        }
 
         // Construct/apply AMI interpolation to determine addressing and weights
-        AMIPtr_ = new AMIPatchToPatchInterpolation
+        AMIPtr_.reset
         (
-            *this,
-            nbrPatch0,
-            surfPtr_()
+            new AMIPatchToPatchInterpolation
+            (
+                *this,
+                nbrPatch0,
+                surf(),
+                faceAreaIntersect::tmMesh,
+                projectPoints_
+            )
         );
-
-
-        vectorField nbrPatchNf
-        (
-            nbrPatch().faceAreas()/mag(nbrPatch().faceAreas())
-        );
-        nbrPatchDelta_ =
-            nbrPatch0.faceCentres() - nbrPatch().faceCellCentres();
-        nbrPatchDelta_ = nbrPatchNf*(nbrPatchNf & nbrPatchDelta_);
     }
 
-    return *AMIPtr_;
+    return AMIPtr_();
 }
 
 
@@ -348,9 +345,11 @@ void Foam::cyclicAMIPolyPatch::movePoints
 
     calcTransforms();
 
-    deleteDemandDrivenData(AMIPtr_);
-
-    AMI();
+    if (owner())
+    {
+        AMIPtr_.clear();
+        AMI();
+    }
 }
 
 
@@ -363,8 +362,6 @@ void Foam::cyclicAMIPolyPatch::initUpdateMesh(PstreamBuffers& pBufs)
 void Foam::cyclicAMIPolyPatch::updateMesh(PstreamBuffers& pBufs)
 {
     polyPatch::updateMesh(pBufs);
-    deleteDemandDrivenData(AMIPtr_);
-//    deleteDemandDrivenData(surfPtr_);
 }
 
 
@@ -380,6 +377,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
 )
 :
     coupledPolyPatch(name, size, start, index, bm),
+    dict_(dictionary::null),
     nbrPatchName_(word::null),
     nbrPatchID_(-1),
     transform_(UNKNOWN),
@@ -387,7 +385,10 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
     rotationCentre_(point::zero),
     separationVector_(vector::zero),
     AMIPtr_(NULL),
-    surfPtr_(NULL)
+    surfPtr_(NULL),
+    projectionSurfaceType_(word::null),
+    projectPoints_(false),
+    projectionName_(word::null)
 {
     // Neighbour patch might not be valid yet so no transformation
     // calculation possible
@@ -403,6 +404,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
 )
 :
     coupledPolyPatch(name, dict, index, bm),
+    dict_(dict),
     nbrPatchName_(dict.lookup("neighbourPatch")),
     nbrPatchID_(-1),
     transform_(UNKNOWN),
@@ -410,23 +412,10 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
     rotationCentre_(point::zero),
     separationVector_(vector::zero),
     AMIPtr_(NULL),
-    surfPtr_
-    (
-        searchableSurface::New
-        (
-            dict.lookup("projectionSurfaceType"),
-            IOobject
-            (
-                dict.lookup("projectionSurfaceName"),
-                bm.mesh().time().constant(),
-                "triSurface",
-                bm.mesh(),
-                IOobject::MUST_READ,
-                IOobject::NO_WRITE
-            ),
-            dict
-        )
-    )
+    surfPtr_(NULL),
+    projectionSurfaceType_(dict.lookup("projectionSurfaceType")),
+    projectPoints_(readBool(dict.lookup("projectPoints"))),
+    projectionName_(dict.lookupOrDefault("surfaceName", name))
 {
     if (nbrPatchName_ == name)
     {
@@ -500,6 +489,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
 )
 :
     coupledPolyPatch(pp, bm),
+    dict_(dictionary::null),
     nbrPatchName_(pp.nbrPatchName_),
     nbrPatchID_(-1),
     transform_(UNKNOWN),
@@ -507,7 +497,10 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
     rotationCentre_(point::zero),
     separationVector_(vector::zero),
     AMIPtr_(NULL),
-    surfPtr_(NULL)
+    surfPtr_(NULL),
+    projectionSurfaceType_(word::null),
+    projectPoints_(false),
+    projectionName_(word::null)
 {
     // Neighbour patch might not be valid yet so no transformation
     // calculation possible
@@ -525,6 +518,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
 )
 :
     coupledPolyPatch(pp, bm, index, newSize, newStart),
+    dict_(dictionary::null),
     nbrPatchName_(nbrPatchName),
     nbrPatchID_(-1),
     transform_(UNKNOWN),
@@ -532,7 +526,10 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
     rotationCentre_(point::zero),
     separationVector_(vector::zero),
     AMIPtr_(NULL),
-    surfPtr_(NULL)
+    surfPtr_(NULL),
+    projectionSurfaceType_(word::null),
+    projectPoints_(false),
+    projectionName_(word::null)
 {
     if (nbrPatchName_ == name())
     {
@@ -564,6 +561,7 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
 )
 :
     coupledPolyPatch(pp, bm, index, mapAddressing, newStart),
+    dict_(pp.dict_),
     nbrPatchName_(pp.nbrPatchName_),
     nbrPatchID_(-1),
     transform_(pp.transform_),
@@ -571,7 +569,10 @@ Foam::cyclicAMIPolyPatch::cyclicAMIPolyPatch
     rotationCentre_(pp.rotationCentre_),
     separationVector_(pp.separationVector_),
     AMIPtr_(NULL),
-    surfPtr_(NULL)
+    surfPtr_(NULL),
+    projectionSurfaceType_(pp.projectionSurfaceType_),
+    projectPoints_(pp.projectPoints_),
+    projectionName_(pp.projectionName_)
 {}
 
 
@@ -619,6 +620,12 @@ Foam::label Foam::cyclicAMIPolyPatch::nbrPatchID() const
 }
 
 
+bool Foam::cyclicAMIPolyPatch::owner() const
+{
+    return index() < nbrPatchID();
+}
+
+
 void Foam::cyclicAMIPolyPatch::transformPosition(pointField& l) const
 {
     if (!parallel())
@@ -712,9 +719,11 @@ void Foam::cyclicAMIPolyPatch::calcGeometry
         nbrAreas
     );
 
-    deleteDemandDrivenData(AMIPtr_);
-
-    AMI();
+    if (owner())
+    {
+        AMIPtr_.clear();
+        AMI();
+    }
 }
 
 
@@ -748,6 +757,45 @@ bool Foam::cyclicAMIPolyPatch::order
 void Foam::cyclicAMIPolyPatch::write(Ostream& os) const
 {
     coupledPolyPatch::write(os);
+    os.writeKeyword("neighbourPatch") << nbrPatchName_
+        << token::END_STATEMENT << nl;
+    os.writeKeyword("projectPoints") << projectPoints_
+        << token::END_STATEMENT << nl;
+    os.writeKeyword("projectionSurfaceType") << projectionSurfaceType_
+        << token::END_STATEMENT << nl;
+    os.writeKeyword("surfaceName") << projectionName_
+        << token::END_STATEMENT << nl;
+    switch (transform_)
+    {
+        case ROTATIONAL:
+        {
+            os.writeKeyword("transform") << transformTypeNames[transform_]
+                << token::END_STATEMENT << nl;
+            os.writeKeyword("rotationAxis") << rotationAxis_
+                << token::END_STATEMENT << nl;
+            os.writeKeyword("rotationCentre") << rotationCentre_
+                << token::END_STATEMENT << nl;
+            break;
+        }
+        case TRANSLATIONAL:
+        {
+            os.writeKeyword("transform") << transformTypeNames[transform_]
+                << token::END_STATEMENT << nl;
+            os.writeKeyword("separationVector") << separationVector_
+                << token::END_STATEMENT << nl;
+            break;
+        }
+        case NOORDERING:
+        {
+            os.writeKeyword("transform") << transformTypeNames[transform_]
+                << token::END_STATEMENT << nl;
+            break;
+        }
+        default:
+        {
+            // no additional info to write
+        }
+    }
 }
 
 
diff --git a/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
index c54ae0f853..cfd20c040f 100644
--- a/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
+++ b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatch.H
@@ -57,15 +57,15 @@ private:
 
     // Private data
 
+        //- Dictionary used during construction
+        const dictionary dict_;
+
         //- Name of other half
         const word nbrPatchName_;
 
         //- Index of other half
         mutable label nbrPatchID_;
 
-        //- Neighbour delta field to *this patch
-        vectorField nbrPatchDelta_;
-
 
         // Transformations
 
@@ -89,11 +89,20 @@ private:
 
 
         //- AMI interpolation class
-        AMIPatchToPatchInterpolation* AMIPtr_;
+        autoPtr<AMIPatchToPatchInterpolation> AMIPtr_;
 
         //- Projection surface
         autoPtr<searchableSurface> surfPtr_;
 
+        //- Projection surface type
+        word projectionSurfaceType_;
+
+        //- Flag to activate projection
+        bool projectPoints_;
+
+        //- Name of projection object (if read)
+        const word projectionName_;
+
 
     // Private Member Functions
 
@@ -138,9 +147,6 @@ protected:
 
 public:
 
-//    //- Declare friendship with processorCyclicPolyPatch
-//    friend class processorCyclicPolyPatch;
-
     //- Runtime type information
     TypeName("cyclicAMI");
 
@@ -255,60 +261,37 @@ public:
         // Access
 
             //- Neighbour patch name
-            const word& nbrPatchName() const
-            {
-                return nbrPatchName_;
-            }
+            inline const word& nbrPatchName() const;
 
             //- Neighbour patch ID
             virtual label nbrPatchID() const;
 
-            //- Does this side own the patch ?
-            virtual bool owner() const
-            {
-                return index() < nbrPatchID();
-            }
+            //- Does this side own the patch?
+            virtual bool owner() const;
 
-            const cyclicAMIPolyPatch& nbrPatch() const
-            {
-                const polyPatch& pp = this->boundaryMesh()[nbrPatchID()];
-                return refCast<const cyclicAMIPolyPatch>(pp);
-            }
+            //- Return a reference to the neihgjbour patch
+            inline const cyclicAMIPolyPatch& nbrPatch() const;
 
+            //- Return a reference to the projection surface
+            const searchableSurface& surf();
 
-            const vectorField& nbrPatchDelta() const
-            {
-                return nbrPatchDelta_;
-            }
-
+            //- Return a reference to the AMI interpolator
             const AMIPatchToPatchInterpolation& AMI();
 
 
             // Transformations
 
                 //- Type of transform
-                transformType transform() const
-                {
-                    return transform_;
-                }
+                inline transformType transform() const;
 
                 //- Axis of rotation for rotational cyclic AMI
-                const vector& rotationAxis() const
-                {
-                    return rotationAxis_;
-                }
+                inline const vector& rotationAxis() const;
 
                 //- point on axis of rotation for rotational cyclic AMI
-                const point& rotationCentre() const
-                {
-                    return rotationCentre_;
-                }
+                inline const point& rotationCentre() const;
 
                 //- Translation vector for translational cyclic AMI
-                const vector& separationVector() const
-                {
-                    return separationVector_;
-                }
+                inline const vector& separationVector() const;
 
                 //- Transform patch-based positions from nbr side to this side
                 virtual void transformPosition(pointField&) const;
@@ -320,17 +303,22 @@ public:
                     const label faceI
                 ) const;
 
-                template<class Type>
-                tmp<Field<Type> > interpolateToSource(const Field<Type>& fld) const
-                {
-                    return AMIPtr_->interpolateToSource(fld);
-                }
 
+            // Interpolations
+
+                //- Interpolate from target field to source
                 template<class Type>
-                tmp<Field<Type> > interpolateToSource(const tmp<Field<Type> >& tFld) const
-                {
-                    return AMIPtr_->interpolateToSource(tFld);
-                }
+                tmp<Field<Type> > interpolateToSource
+                (
+                    const Field<Type>& fld
+                ) const;
+
+                //- Interpolate from target tmp field to source
+                template<class Type>
+                tmp<Field<Type> > interpolateToSource
+                (
+                    const tmp<Field<Type> >& tFld
+                ) const;
 
 
         //- Calculate the patch geometry
@@ -377,6 +365,16 @@ public:
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+#include "cyclicAMIPolyPatchI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#ifdef NoRepository
+    #include "cyclicAMIPolyPatchTemplates.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
 #endif
 
 // ************************************************************************* //
diff --git a/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H
new file mode 100644
index 0000000000..8bb8badc77
--- /dev/null
+++ b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchI.H
@@ -0,0 +1,67 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+inline const Foam::word& Foam::cyclicAMIPolyPatch::nbrPatchName() const
+{
+    return nbrPatchName_;
+}
+
+
+inline const Foam::cyclicAMIPolyPatch&
+Foam::cyclicAMIPolyPatch::nbrPatch() const
+{
+    const polyPatch& pp = this->boundaryMesh()[nbrPatchID()];
+    return refCast<const cyclicAMIPolyPatch>(pp);
+}
+
+
+inline Foam::coupledPolyPatch::transformType
+Foam::cyclicAMIPolyPatch::transform() const
+{
+    return transform_;
+}
+
+
+inline const Foam::vector& Foam::cyclicAMIPolyPatch::rotationAxis() const
+{
+    return rotationAxis_;
+}
+
+
+inline const Foam::point& Foam::cyclicAMIPolyPatch::rotationCentre() const
+{
+    return rotationCentre_;
+}
+
+
+inline const Foam::vector& Foam::cyclicAMIPolyPatch::separationVector() const
+{
+    return separationVector_;
+}
+
+
+// ************************************************************************* //
diff --git a/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchTemplates.C b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchTemplates.C
new file mode 100644
index 0000000000..0e0243c5ef
--- /dev/null
+++ b/src/AMIInterpolation/patches/cyclicAMIPolyPatch/cyclicAMIPolyPatchTemplates.C
@@ -0,0 +1,64 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 2011-2011 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software: you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by
+    the Free Software Foundation, either version 3 of the License, or
+    (at your option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM.  If not, see <http://www.gnu.org/licenses/>.
+
+\*---------------------------------------------------------------------------*/
+
+//#include "cyclicAMIPolyPatch.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+template<class Type>
+Foam::tmp<Foam::Field<Type> > Foam::cyclicAMIPolyPatch::interpolateToSource
+(
+    const Field<Type>& fld
+) const
+{
+    if (owner())
+    {
+        return AMIPtr_->interpolateToSource(fld);
+    }
+    else
+    {
+        return nbrPatch().AMIPtr_->interpolateToTarget(fld);
+    }
+}
+
+
+template<class Type>
+Foam::tmp<Foam::Field<Type> > Foam::cyclicAMIPolyPatch::interpolateToSource
+(
+    const tmp<Field<Type> >& tFld
+) const
+{
+    if (owner())
+    {
+        return AMIPtr_->interpolateToSource(tFld);
+    }
+    else
+    {
+        return nbrPatch().AMIPtr_->interpolateToTarget(tFld);
+    }
+}
+
+
+// ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/boundBox/boundBox.C b/src/OpenFOAM/meshes/boundBox/boundBox.C
index 823b82b56a..63795f9145 100644
--- a/src/OpenFOAM/meshes/boundBox/boundBox.C
+++ b/src/OpenFOAM/meshes/boundBox/boundBox.C
@@ -165,16 +165,10 @@ Foam::tmp<Foam::pointField> Foam::boundBox::points() const
 
 void Foam::boundBox::inflate(const scalar s)
 {
-    boundBox bb(*this);
+    vector ext = vector::one*s*mag();
 
-    vector dir = bb.span();
-    scalar magDir = Foam::mag(dir);
-    dir /= magDir;
-
-    scalar ext = s*magDir;
-
-    min_ -= dir*ext;
-    max_ += dir*ext;
+    min_ -= ext;
+    max_ += ext;
 }