Merge remote branch 'remotes/origin/testLoop/PhilipposeRajanAdaptiveMeshingMerge'

applications/solvers/compressible/sonicFoamAutoMotion/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \

applications/solvers/engine/icoDyMEngineFoam/Make/options

@@ -1,8 +1,8 @@
 EXE_INC = \
     -I$(LIB_SRC)/engine/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/topoChangerFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/topoChangerFvMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 

applications/solvers/engine/sonicTurbDyMEngineFoam/Make/options

@@ -1,8 +1,8 @@
 EXE_INC = \
     -I$(LIB_SRC)/engine/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/basic/lnInclude \
     -I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel
 

applications/solvers/engine/turbDyMEngineFoam/Make/options

@@ -1,8 +1,8 @@
 EXE_INC = \
     -I$(LIB_SRC)/engine/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/topoChangerFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/topoChangerFvMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/transportModels \
     -I$(LIB_SRC)/transportModels/incompressible/singlePhaseTransportModel \

applications/solvers/incompressible/icoDyMFoam/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 

applications/solvers/incompressible/icoDyMSimpleFoam/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 

applications/solvers/incompressible/pimpleDyMFoam/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
     -I$(LIB_SRC)/transportModels \

applications/solvers/multiphase/compressibleInterDyMFoam/Make/options

@@ -4,9 +4,9 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
     -I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude
 
 EXE_LIBS = \
     -linterfaceProperties \

applications/solvers/multiphase/interDyMFoam/Make/options

@@ -5,9 +5,9 @@ EXE_INC = \
     -I$(LIB_SRC)/transportModels/interfaceProperties/lnInclude \
     -I$(LIB_SRC)/turbulenceModels/incompressible/turbulenceModel \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude
 
 EXE_LIBS = \
     -linterfaceProperties \

applications/solvers/stressAnalysis/icoFsiFoam/Make/options

@@ -1,10 +1,10 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     $(WM_DECOMP_INC) \
     -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
 
 
 EXE_LIBS = \

applications/solvers/surfaceTracking/bubbleInterTrackFoam/Make/options

@@ -4,9 +4,9 @@ EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/finiteArea/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
     $(WM_DECOMP_INC)
 
 EXE_LIBS = \

applications/solvers/surfaceTracking/freeSurface/Make/options

@@ -2,8 +2,8 @@ EXE_INC = \
     -I$(LIB_SRC)/finiteArea/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
     -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
     $(WM_DECOMP_INC)
 

applications/solvers/surfaceTracking/interTrackFoam/Make/options

@@ -4,8 +4,8 @@ EXE_INC = \
     -I$(LIB_SRC)/finiteArea/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     $(WM_DECOMP_INC)
 
 EXE_LIBS = \

applications/solvers/surfaceTracking/sloshingInterTrackFoam/Make/options

@@ -3,9 +3,9 @@ EXE_INC = \
     -I$(LIB_SRC)/finiteArea/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
     $(WM_DECOMP_INC)
 
 EXE_LIBS = \

applications/utilities/mesh/advanced/autoRefineMesh/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude

applications/utilities/mesh/advanced/collapseEdges/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/mesh/advanced/combinePatchFaces/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
 
 
 EXE_LIBS = \

applications/utilities/mesh/advanced/modifyMesh/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ltriSurface \

applications/utilities/mesh/advanced/refineHexMesh/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 

applications/utilities/mesh/advanced/refineWallLayer/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/mesh/advanced/removeFaces/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 
 EXE_LIBS = \

applications/utilities/mesh/advanced/selectCells/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude

applications/utilities/mesh/advanced/splitCells/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/mesh/conversion/fluent3DMeshToFoam/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lmeshTools \

applications/utilities/mesh/conversion/fluentMeshToFoam/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lmeshTools \

applications/utilities/mesh/conversion/gmshToFoam/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \
-    -lmeshTools 
+    -lmeshTools

applications/utilities/mesh/conversion/polyDualMesh/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \

applications/utilities/mesh/generation/blockMesh/Make/options

@@ -2,7 +2,7 @@ EXE_INC = \
     -I$(curvedEdges) \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/ODE/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lmeshTools \

applications/utilities/mesh/generation/extrude2DMesh/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     /* -DFULLDEBUG -g -O0 */ \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 
 EXE_LIBS = \

applications/utilities/mesh/generation/extrudeMesh/Make/options

@@ -2,7 +2,7 @@ EXE_INC = \
     -IextrudedMesh \
     -IextrudeModel/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lmeshTools \

applications/utilities/mesh/generation/snappyHexMesh/Make/options

@@ -3,7 +3,7 @@ EXE_INC = \
     -I$(LIB_SRC)/autoMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/edgeMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 

applications/utilities/mesh/manipulation/attachMesh/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/mesh/manipulation/autoPatch/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 
 EXE_LIBS = \

applications/utilities/mesh/manipulation/createBaffles/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \

applications/utilities/mesh/manipulation/createPatch/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \

applications/utilities/mesh/manipulation/mergeMeshes/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/mesh/manipulation/mergeOrSplitBaffles/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \

applications/utilities/mesh/manipulation/moveDyMEngineMesh/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/engine/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude

applications/utilities/mesh/manipulation/moveDynamicMesh/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicFvMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicFvMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude
 

applications/utilities/mesh/manipulation/moveMesh/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \

applications/utilities/mesh/manipulation/refineMesh/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 

applications/utilities/mesh/manipulation/renumberMesh/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/decompositionMethods/decompositionMethods/lnInclude
 

applications/utilities/mesh/manipulation/splitMesh/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/meshTools/lnInclude 
+    -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/mesh/manipulation/splitMeshRegions/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 
 EXE_LIBS = \

applications/utilities/mesh/manipulation/stitchMesh/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/finiteVolume/lnInclude 
+    -I$(LIB_SRC)/finiteVolume/lnInclude
 
 EXE_LIBS = \
     -ldynamicMesh \

applications/utilities/parallelProcessing/reconstructParMesh/Make/options

@@ -1,5 +1,5 @@
 EXE_INC = \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude
 

applications/utilities/parallelProcessing/redistributeMeshPar/Make/options

@@ -2,7 +2,7 @@ EXE_INC = \
     -I$(LIB_SRC)/decompositionMethods/decompositionMethods/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 EXE_LIBS = \
     -lfiniteVolume \

applications/utilities/surface/surfaceToPatch/Make/options

@@ -1,7 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude
 
 
 EXE_LIBS = \

src/Allwmake

@@ -38,20 +38,15 @@ wmake libso finiteArea
 wmake libso lduSolvers
 wmake libso blockMatrix
 
-wmake libso dynamicMesh
-
 (cd tetDecompositionFiniteElement ; ./Allwmake)
 
-wmake libso fvMotionSolver
+# Make dynamicMesh tools
-(cd tetDecompositionMotionSolver ; ./Allwmake)
+dynamicMesh/AllwmakeLnInclude
+dynamicMesh/Allwmake
 
 wmake libso coupledMatrix
 wmake libso sampling
 
-wmake libso dynamicFvMesh
-wmake libso topoChangerFvMesh
-wmake libso engine
-
 wmake libso ODE
 wmake libso POD
 wmake libso randomProcesses
@@ -68,7 +63,6 @@ wmake libso errorEstimation
 
 fvAgglomerationMethods/Allwmake
 
-wmake libso fvMotionSolver
 wmake libso engine
 
 # ----------------------------------------------------------------- end-of-file

src/OpenFOAM/containers/Lists/Stack/Stack.H

@@ -0,0 +1,115 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Stack
+
+Description
+    Thread-safe stack implementation using an internal DynamicList
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    StackI.H
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef Stack_H
+#define Stack_H
+
+#include "label.H"
+#include "DynamicList.H"
+#include "multiThreader.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class Stack Declaration
+\*---------------------------------------------------------------------------*/
+
+class Stack
+{
+    // Private data
+
+        //- Internal list
+        DynamicList<label> stack_;
+
+        //- Mutex for multi-threading
+        Mutex stackMutex_;
+
+public:
+
+    // Constructor
+    Stack(){}
+
+    // Member functions for access to the stack
+
+        //- Push items on to the stack
+        inline void push(const label index);
+
+        //- Insert item onto stack
+        inline void insert(const label index);
+
+        //- Pop an item off the stack
+        inline label pop();
+
+        //- Remove a specific index off the stack
+        inline void remove(const label index);
+
+        //- Return if a stack is empty or not
+        inline bool empty();
+
+        //- Return the size of the stack
+        inline label size();
+
+        //- Clear out the stack
+        inline void clear();
+
+        //- Print out the stack
+        inline void print();
+
+        //- Access the stack as a list
+        inline label operator[](const label index);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "StackI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/OpenFOAM/containers/Lists/Stack/StackI.H

@@ -0,0 +1,158 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    Stack
+
+Description
+    Member functions of the Stack class
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*---------------------------------------------------------------------------*/
+
+#include "ListOps.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+// Push items on to the stack
+inline void Stack::push(const label index)
+{
+    stackMutex_.lock();
+
+    if (findIndex(stack_, index) == -1)
+    {
+        stack_.append(index);
+    }
+
+    stackMutex_.unlock();
+}
+
+
+//- Insert item onto stack (no checking)
+inline void Stack::insert(const label index)
+{
+    stack_.append(index);
+}
+
+
+// Pop an item off the stack
+inline label Stack::pop()
+{
+    stackMutex_.lock();
+
+    const label index = stack_.remove();
+
+    stackMutex_.unlock();
+
+    return index;
+}
+
+
+// Remove a specific item off the stack
+inline void Stack::remove(const label index)
+{
+    stackMutex_.lock();
+
+    label loc = findIndex(stack_, index);
+
+    if (loc != -1)
+    {
+        // Create a new list
+        labelList newList(stack_.size() - 1);
+
+        label n = 0;
+
+        // Copy items upto location
+        for(label i = 0; i < loc; i++)
+        {
+            newList[n++] = stack_[i];
+        }
+
+        // Copy items from location
+        for(label i = (loc + 1); i < stack_.size(); i++)
+        {
+            newList[n++] = stack_[i];
+        }
+
+        // Overwrite
+        stack_ = newList;
+    }
+
+    stackMutex_.unlock();
+}
+
+
+// Return if the stack is empty or not
+inline bool Stack::empty()
+{
+    return (stack_.size() == 0);
+}
+
+
+//- Return the size of the stack
+inline label Stack::size()
+{
+    return stack_.size();
+}
+
+
+//- Clear out the stack
+inline void Stack::clear()
+{
+    stackMutex_.lock();
+
+    stack_.clear();
+
+    stackMutex_.unlock();
+}
+
+
+//- Print out the stack
+inline void Stack::print()
+{
+    Info << stack_ << endl;
+}
+
+
+//- Access the stack as a list
+inline label Stack::operator[]
+(
+    const label index
+)
+{
+    return stack_[index];
+}
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/autoMesh/Make/options

@@ -1,6 +1,6 @@
 EXE_INC = \
     -I$(LIB_SRC)/decompositionMethods/decompositionMethods/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \

src/dynamicFvMesh/Make/options

@@ -1,17 +0,0 @@
-EXE_INC = \
-    -I$(LIB_SRC)/triSurface/lnInclude \
-    -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/dynamicMesh/lnInclude \
-    -I$(LIB_SRC)/finiteVolume/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionMotionSolver/lnInclude \
-    -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
-    $(WM_DECOMP_INC) \
-    -I$(LIB_SRC)/fvMotionSolver/lnInclude
-
-LIB_LIBS = \
-    -ltriSurface \
-    -lmeshTools \
-    -ldynamicMesh \
-    -lfiniteVolume \
-    $(WM_DECOMP_LIBS) \
-    -lfvMotionSolver

src/dynamicMesh/Allwmake

@@ -0,0 +1,15 @@
+#!/bin/sh
+cd ${0%/*} || exit 1    # run from this directory
+set -x
+
+wmakeLnInclude dynamicMesh
+
+wmake libso dynamicMesh
+
+# Make meshMotion solvers
+meshMotion/Allwmake
+
+wmake libso dynamicFvMesh
+wmake libso topoChangerFvMesh
+
+# ----------------------------------------------------------------- end-of-file

src/dynamicMesh/AllwmakeLnInclude

@@ -0,0 +1,10 @@
+#!/bin/sh
+cd ${0%/*} || exit 1    # run from this directory
+set -x
+
+wmakeLnInclude dynamicMesh
+wmakeLnInclude dynamicFvMesh
+wmakeLnInclude topoChangerFvMesh
+
+# Make lnInclude directories for motionSolvers
+meshMotion/AllwmakeLnInclude

src/dynamicMesh/dynamicFvMesh/Make/files

@@ -20,7 +20,27 @@ $(solidBodyMotionFunctions)/translation/translation.C
 mixerGgiFvMesh/mixerGgiFvMesh.C
 turboFvMesh/turboFvMesh.C
 
-tetMetrics/tetMetric.C
+tetMetrics = dynamicTopoFvMesh/tetMetrics
-tetMetrics/tetMetrics.C
+$(tetMetrics)/tetMetric.C
+$(tetMetrics)/tetMetrics.C
+
+lengthScaleEstimator = dynamicTopoFvMesh/lengthScaleEstimator
+$(lengthScaleEstimator)/lengthScaleEstimator.C
+
+eMesh = dynamicTopoFvMesh/eMesh
+$(eMesh)/eMesh.C
+$(eMesh)/eMeshDemandDrivenData.C
+$(eMesh)/eBoundaryMesh/eBoundaryMesh.C
+ePatches = $(eMesh)/ePatches
+$(ePatches)/ePatch/ePatch.C
+$(ePatches)/ePatch/newEPatch.C
+dynamicTopoFvMesh/meshOps.C
+dynamicTopoFvMesh/dynamicTopoFvMesh.C
+dynamicTopoFvMesh/dynamicTopoFvMeshCheck.C
+dynamicTopoFvMesh/dynamicTopoFvMeshReOrder.C
+dynamicTopoFvMesh/dynamicTopoFvMeshMapping.C
+dynamicTopoFvMesh/edgeSwap.C
+dynamicTopoFvMesh/edgeBisect.C
+dynamicTopoFvMesh/edgeCollapse.C
 
 LIB = $(FOAM_LIBBIN)/libdynamicFvMesh

src/dynamicMesh/dynamicFvMesh/Make/options

@@ -0,0 +1,21 @@
+EXE_INC = \
+    -I$(LIB_SRC)/triSurface/lnInclude \
+    -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/dynamicMesh/lnInclude \
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/tetDecompositionMotionSolver/lnInclude \
+    -I$(LIB_SRC)/tetDecompositionFiniteElement/lnInclude \
+    $(WM_DECOMP_INC) \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/fvMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/RBFMotionSolver/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/meshMotion/mesquiteMotionSolver/lnInclude
+
+LIB_LIBS = \
+    -ltriSurface \
+    -lmeshTools \
+    -ldynamicMesh \
+    -lfiniteVolume \
+    $(WM_DECOMP_LIBS) \
+    -lfvMotionSolver \
+    -lRBFMotionSolver \
+    -lmesquiteMotionSolver

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/changeMap.H

@@ -0,0 +1,162 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    changeMap
+
+Description
+    Accumulate topology change statistics
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    changeMapI.H
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef changeMap_H
+#define changeMap_H
+
+#include "labelList.H"
+#include "FixedList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                           Class changeMap Declaration
+\*---------------------------------------------------------------------------*/
+
+class changeMap
+{
+    // Sliver type index.
+    //  Type 1: Sliver
+    //  Type 2: Cap
+    //  Type 3: Spade
+    //  Type 4: Wedge
+    label type_;
+
+    // Data specific to sliver-type cells
+    label firstEdge_;
+    label secondEdge_;
+
+    // Data specific to cap-type cells
+    label apexPoint_;
+    label opposingFace_;
+
+    // Entities that were added during the operation.
+    //  - Add a master mapping for slaves.
+    //  - Masters get a map of -1.
+    List<FixedList<label,2> > addedPoints_;
+    List<FixedList<label,2> > addedEdges_;
+    List<FixedList<label,2> > addedFaces_;
+    List<FixedList<label,2> > addedCells_;
+
+public:
+
+    // Constructor
+    changeMap()
+    :
+        type_(-1),
+        firstEdge_(-1),
+        secondEdge_(-1),
+        apexPoint_(-1),
+        opposingFace_(-1),
+        addedPoints_(0),
+        addedEdges_(0),
+        addedFaces_(0),
+        addedCells_(0)
+    {}
+
+    //- Access
+
+        // Type
+        inline label& type();
+        inline label type() const;
+
+        // For sliver-type cells, opposite edges
+        // are identified for removal.
+        inline label& firstEdge();
+        inline label& secondEdge();
+
+        // For cap-type cells, the face requiring splitting
+        // is identified for removal.
+        inline label& apexPoint();
+        inline label& opposingFace();
+
+        // Added entities
+        inline void addPoint
+        (
+            const label pIndex,
+            const label master = -1
+        );
+
+        inline void addEdge
+        (
+            const label eIndex,
+            const label master = -1
+        );
+
+        inline void addFace
+        (
+            const label fIndex,
+            const label master = -1
+        );
+
+        inline void addCell
+        (
+            const label cIndex,
+            const label master = -1
+        );
+
+        // Return the list of added entities
+        inline const List<FixedList<label,2> >& addedPointList() const;
+        inline const List<FixedList<label,2> >& addedEdgeList() const;
+        inline const List<FixedList<label,2> >& addedFaceList() const;
+        inline const List<FixedList<label,2> >& addedCellList() const;
+
+    //- Operators
+
+        inline void operator=(const changeMap& rhs);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "changeMapI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/changeMapI.H

@@ -0,0 +1,181 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+\*---------------------------------------------------------------------------*/
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+// Type
+inline label& changeMap::type()
+{
+    return type_;
+}
+
+
+inline label changeMap::type() const
+{
+    return type_;
+}
+
+
+// For sliver-type cells, opposite edges
+// are identified for removal.
+inline label& changeMap::firstEdge()
+{
+    return firstEdge_;
+}
+
+
+inline label& changeMap::secondEdge()
+{
+    return secondEdge_;
+}
+
+
+// For cap-type cells, the face requiring splitting
+// is identified for removal.
+inline label& changeMap::apexPoint()
+{
+    return apexPoint_;
+}
+
+
+inline label& changeMap::opposingFace()
+{
+    return opposingFace_;
+}
+
+
+// Added entities
+inline void changeMap::addPoint
+(
+    const label pIndex,
+    const label master
+)
+{
+    label curSize = addedPoints_.size();
+
+    addedPoints_.setSize(curSize + 1);
+
+    addedPoints_[curSize][0] = pIndex;
+    addedPoints_[curSize][1] = master;
+}
+
+
+inline void changeMap::addEdge
+(
+    const label eIndex,
+    const label master
+)
+{
+    label curSize = addedEdges_.size();
+
+    addedEdges_.setSize(curSize + 1);
+
+    addedEdges_[curSize][0] = eIndex;
+    addedEdges_[curSize][1] = master;
+}
+
+
+inline void changeMap::addFace
+(
+    const label fIndex,
+    const label master
+)
+{
+    label curSize = addedFaces_.size();
+
+    addedFaces_.setSize(curSize + 1);
+
+    addedFaces_[curSize][0] = fIndex;
+    addedFaces_[curSize][1] = master;
+}
+
+
+inline void changeMap::addCell
+(
+    const label cIndex,
+    const label master
+)
+{
+    label curSize = addedCells_.size();
+
+    addedCells_.setSize(curSize + 1);
+
+    addedCells_[curSize][0] = cIndex;
+    addedCells_[curSize][1] = master;
+}
+
+
+// Return an added point
+inline const List<FixedList<label,2> >&
+changeMap::addedPointList() const
+{
+    return addedPoints_;
+}
+
+
+// Return the list of added entities
+inline const List<FixedList<label,2> >&
+changeMap::addedEdgeList() const
+{
+    return addedEdges_;
+}
+
+
+inline const List<FixedList<label,2> >&
+changeMap::addedFaceList() const
+{
+    return addedFaces_;
+}
+
+
+inline const List<FixedList<label,2> >&
+changeMap::addedCellList() const
+{
+    return addedCells_;
+}
+
+
+inline void changeMap::operator=(const changeMap& rhs)
+{
+    type_ = rhs.type_;
+
+    firstEdge_    = rhs.firstEdge_;
+    secondEdge_   = rhs.secondEdge_;
+    apexPoint_    = rhs.apexPoint_;
+    opposingFace_ = rhs.opposingFace_;
+
+    addedPoints_  = rhs.addedPoints_;
+    addedEdges_   = rhs.addedEdges_;
+    addedFaces_   = rhs.addedFaces_;
+    addedCells_   = rhs.addedCells_;
+}
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/dynamicTopoFvMesh.H

@@ -0,0 +1,864 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    dynamicTopoFvMesh
+
+Description
+    An implementation of dynamic changes to mesh-topology
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    dynamicTopoFvMesh.C
+    dynamicTopoFvMeshI.H
+    dynamicTopoFvMeshCheck.C
+    dynamicTopoFvMeshReOrder.C
+    dynamicTopoFvMeshMapping.C
+    edgeBisect.C
+    edgeCollapse.C
+    edgeSwap.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef dynamicTopoFvMesh_H
+#define dynamicTopoFvMesh_H
+
+#include "Switch.H"
+#include "tetMetric.H"
+#include "DynamicField.H"
+#include "threadHandler.H"
+#include "dynamicFvMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Class forward declarations
+class eMesh;
+class Stack;
+class changeMap;
+class objectMap;
+class topoMapper;
+class motionSolver;
+class lengthScaleEstimator;
+
+/*---------------------------------------------------------------------------*\
+                      Class dynamicTopoFvMesh Declaration
+\*---------------------------------------------------------------------------*/
+
+class dynamicTopoFvMesh
+:
+    public dynamicFvMesh
+{
+    // Private data
+
+        //- Topology change flag
+        bool topoChangeFlag_;
+
+        //- Dynamic mesh dictionary
+        IOdictionary dict_;
+
+        //- Should all options be mandatorily specified?
+        //  Handy for first-time use.
+        Switch mandatory_;
+
+        //- Mesh characteristics [2D/3D]
+        Switch twoDMesh_;
+
+        //- Edge refinement switch
+        Switch edgeRefinement_;
+
+        //- Switch for cell-bandwidth reduction
+        Switch bandWidthReduction_;
+
+        //- Specify the re-meshing interval
+        label interval_;
+
+        //- Edge-mesh
+        autoPtr<eMesh> eMeshPtr_;
+
+        //- Field mapper
+        autoPtr<topoMapper> mapper_;
+
+        //- Mesh motion solver
+        autoPtr<motionSolver> motionSolver_;
+
+        //- Length scale estimator
+        autoPtr<lengthScaleEstimator> lengthEstimator_;
+
+        //- Lists that dynamically resize during topo-changes
+        //   - Since resizes happen infrequently,
+        //     scale up by 10% to save memory.
+        //   - Current C++ standard doesn't
+        //     support template typedefs,
+        //     so this is a work-around
+        template <class T>
+        class resizableList
+        {
+        public:
+            typedef DynamicList<T, 0, 11, 10> Type;
+        };
+
+        template <class T>
+        class resizableField
+        {
+        public:
+            typedef DynamicField<T> Type;
+        };
+
+        resizableField<point>::Type oldPoints_, points_;
+        resizableList<face>::Type faces_;
+        resizableList<label>::Type owner_, neighbour_;
+        resizableList<cell>::Type cells_;
+        resizableList<edge>::Type edges_;
+        resizableList<labelList>::Type pointEdges_, edgePoints_;
+        resizableList<labelList>::Type edgeFaces_, faceEdges_;
+        resizableList<scalar>::Type lengthScale_;
+
+        //- Size information
+        labelList oldPatchSizes_, patchSizes_;
+        labelList oldPatchStarts_, patchStarts_;
+        labelList oldEdgePatchSizes_, edgePatchSizes_;
+        labelList oldEdgePatchStarts_, edgePatchStarts_;
+        labelList oldPatchNMeshPoints_, patchNMeshPoints_;
+        label nOldPoints_, nPoints_;
+        label nOldEdges_, nEdges_;
+        label nOldFaces_, nFaces_;
+        label nOldCells_, nCells_;
+        label nOldInternalFaces_, nInternalFaces_;
+        label nOldInternalEdges_, nInternalEdges_;
+
+        //- Maps [Old-to-New]
+        labelList reversePointMap_;
+        labelList reverseEdgeMap_;
+        labelList reverseFaceMap_;
+        labelList reverseCellMap_;
+
+        //- Maps [New-to-Old]
+        labelList pointMap_;
+        labelList edgeMap_;
+        labelList faceMap_;
+        labelList cellMap_;
+
+        //- Maps for the renumbering of added entities
+        Map<label> addedPointRenumbering_;
+        Map<label> addedEdgeRenumbering_;
+        Map<label> addedFaceRenumbering_;
+        Map<label> addedCellRenumbering_;
+        Map<label> addedFacePatches_;
+        Map<label> addedEdgePatches_;
+        Map<label> addedPointZones_;
+        Map<label> addedFaceZones_;
+        Map<label> addedCellZones_;
+
+        // Information for field-mapping
+        Map<labelList> faceParents_;
+        List<scalarField> faceWeights_;
+        List<vectorField> faceCentres_;
+        Map<labelList> cellParents_;
+        List<scalarField> cellWeights_;
+        List<vectorField> cellCentres_;
+
+        // Information for mapPolyMesh
+        List<objectMap> pointsFromPoints_;
+        List<objectMap> facesFromPoints_;
+        List<objectMap> facesFromEdges_;
+        List<objectMap> facesFromFaces_;
+        List<objectMap> cellsFromPoints_;
+        List<objectMap> cellsFromEdges_;
+        List<objectMap> cellsFromFaces_;
+        List<objectMap> cellsFromCells_;
+
+        //- Maps to keep track of entities deleted after addition
+        labelHashSet deletedPoints_;
+        labelHashSet deletedEdges_;
+        labelHashSet deletedFaces_;
+        labelHashSet deletedCells_;
+
+        //- List of flipped faces / modified old-points
+        labelHashSet flipFaces_;
+        Map<labelList> modPoints_;
+
+        //- Run-time statistics
+        label maxModifications_;
+        FixedList<label, 8> statistics_;
+
+        //- Sliver exudation
+        scalar sliverThreshold_;
+        Map<scalar> thresholdSlivers_;
+
+        //- Specific to proximity-based refinement
+        List<labelPair> slicePairs_;
+
+        //- 3D Mesh Flipping data
+        label maxTetsPerEdge_;
+        scalar swapDeviation_;
+        Switch allowTableResize_;
+        labelList noSwapPatchIDs_;
+
+        //- Stack-list of entities to be checked for topo-changes.
+        List<Stack> entityStack_;
+
+        //- Support for multithreading
+        autoPtr<IOmultiThreader> threader_;
+
+        // Pointer-list of thread-handlers
+        typedef threadHandler<dynamicTopoFvMesh> meshHandler;
+
+        PtrList<meshHandler> handlerPtr_;
+
+        // Entity mutexes used for synchronization
+        // in multi-threaded reOrdering
+        FixedList<Mutex, 4> entityMutex_;
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        dynamicTopoFvMesh(const dynamicTopoFvMesh&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const dynamicTopoFvMesh&);
+
+        // Perform a Delaunay test on an internal face
+        bool testDelaunay(const label fIndex) const;
+
+        //- Quality metric for tetrahedra in 3D
+        tetMetric::tetMetricReturnType tetMetric_;
+
+        // Compute mapping weights for modified entities
+        void computeMapping
+        (
+            const scalar matchTol,
+            const bool skipMapping,
+            const label faceStart,
+            const label faceSize,
+            const label cellStart,
+            const label cellSize
+        );
+
+        // Static equivalent for multiThreading
+        static void computeMappingThread(void *argument);
+
+        // Compute parents for inverse-distance weighting
+        void computeParents
+        (
+            const label index,
+            const labelList& mapCandidates,
+            const labelListList& oldNeighbourList,
+            const vectorField& oldCentres,
+            const label dimension,
+            labelList& parents
+        ) const;
+
+        // Routine to invoke threaded mapping
+        void threadedMapping(scalar matchTol, bool skipMapping);
+
+        // Initialize mesh edges and related connectivity lists
+        void initEdges();
+
+        // Load the mesh-quality metric from the library
+        void loadMetric();
+
+        // Load the mesh-motion solver
+        void loadMotionSolver();
+
+        // Load the field mapper
+        void loadFieldMapper();
+
+        // Load the length scale estimator
+        void loadLengthScaleEstimator();
+
+        // Initialize the threading environment
+        void initializeThreadingEnvironment(const label specThreads = -1);
+
+        // Return a non-const reference to the lengthScaleEstimator
+        inline lengthScaleEstimator& lengthEstimator();
+
+        // Return a const reference to the lengthScaleEstimator
+        inline const lengthScaleEstimator& lengthEstimator() const;
+
+        // Return a const reference to the multiThreader
+        inline const multiThreader& threader() const;
+
+        // Return a reference to the entity mutexes
+        inline const Mutex& entityMutex(const label entity) const;
+
+        // Return the edge index for a provided edge
+        inline label getEdgeIndex(const edge& edgeToCheck) const;
+
+        // Given a boundary face, pick out a boundary edge that
+        // contains a triangular face.
+        inline label getTriBoundaryEdge(const label fIndex) const;
+
+        // Return length-scale at an face-location in the mesh [2D]
+        inline scalar faceLengthScale(const label fIndex) const;
+
+        // Return length-scale at an edge-location in the mesh [3D]
+        inline scalar edgeLengthScale(const label eIndex) const;
+
+        // Check for bisection
+        inline bool checkBisection(const label index) const;
+
+        // Check for collapse
+        inline bool checkCollapse(const label index) const;
+
+        // MultiThreaded topology modifier
+        void threadedTopoModifier();
+
+        // 2D Edge-swapping engine
+        static void swap2DEdges(void *argument);
+
+        // 3D Edge-swapping engine
+        static void swap3DEdges(void *argument);
+
+        // Edge refinement engine
+        static void edgeRefinementEngine(void *argument);
+
+        // Return the entity stack for a particular thread
+        inline Stack& stack(const label threadID);
+
+        // Return the integer ID for a given thread
+        inline label self() const;
+
+        // Initialize stacks
+        inline void initStacks(const labelHashSet& entities);
+
+        // Method to determine the boundary patch index for a given face
+        inline label whichPatch(const label index) const;
+
+        // Method to determine the boundary patch index for a given edge
+        inline label whichEdgePatch(const label index) const;
+
+        // Report topo-modification status
+        inline label status(const label type) const;
+
+        // Method for the swapping of a quad-face in 2D
+        const changeMap
+        swapQuadFace
+        (
+            const label fIndex
+        );
+
+        // Method for the bisection of a quad-face in 2D
+        const changeMap
+        bisectQuadFace
+        (
+            const label fIndex,
+            const changeMap& masterMap,
+            bool checkOnly = false,
+            bool forceOp = false
+        );
+
+        // Method for the collapse of a quad-face in 2D
+        const changeMap
+        collapseQuadFace
+        (
+            const label fIndex,
+            label overRideCase = -1,
+            bool checkOnly = false,
+            bool forceOp = false
+        );
+
+        // Method for the swapping of an edge in 3D
+        void swapEdge
+        (
+            const label eIndex,
+            bool forceOp = false
+        );
+
+        // Method for the bisection of an edge in 3D
+        const changeMap
+        bisectEdge
+        (
+            const label eIndex,
+            bool checkOnly = false,
+            bool forceOp = false
+        );
+
+        // Method for the trisection of a face in 3D
+        const changeMap
+        trisectFace
+        (
+            const label fIndex,
+            bool checkOnly = false,
+            bool forceOp = false
+        );
+
+        // Method for the collapse of an edge in 3D
+        const changeMap
+        collapseEdge
+        (
+            const label eIndex,
+            label overRideCase = -1,
+            bool checkOnly = false,
+            bool forceOp = false
+        );
+
+        // Utility method to check for invalid face-bisection.
+        bool checkBisection
+        (
+            const label fIndex,
+            const label bFaceIndex,
+            bool forceOp = false
+        ) const;
+
+        // Utility method to check for invalid face-collapse.
+        bool checkCollapse
+        (
+            const labelList& triFaces,
+            const FixedList<label,2>& c0BdyIndex,
+            const FixedList<label,2>& c1BdyIndex,
+            const FixedList<label,2>& pointIndex,
+            const FixedList<point,2>& newPoint,
+            const FixedList<point,2>& oldPoint,
+            scalar& collapseQuality,
+            const bool checkNeighbour,
+            bool forceOp = false
+        ) const;
+
+        // Utility method to check for invalid edge-collapse.
+        bool checkCollapse
+        (
+            const point& newPoint,
+            const point& oldPoint,
+            const label pointIndex,
+            const label cellIndex,
+            labelHashSet& cellsChecked,
+            scalar& collapseQuality,
+            bool forceOp = false
+        ) const;
+
+        // Remove 2D sliver cells from the mesh
+        void remove2DSlivers();
+
+        // Identify the sliver type in 3D
+        const changeMap identifySliverType(const label cIndex) const;
+
+        // Remove sliver cells
+        void removeSlivers();
+
+        // Insert the specified cell to the mesh
+        label insertCell
+        (
+            const cell& newCell,
+            const scalar lengthScale,
+            const label zoneID = -1
+        );
+
+        // Set fill-in mapping information for a particular cell
+        void setCellMapping
+        (
+            const label cIndex,
+            const labelList& mapCells,
+            bool addEntry = true
+        );
+
+        // Set fill-in mapping information for a particular face
+        void setFaceMapping
+        (
+            const label fIndex,
+            const labelList& mapFaces = labelList(0)
+        );
+
+        // Remove the specified cells from the mesh
+        const changeMap
+        removeCells
+        (
+            const labelList& cList,
+            const label patch
+        );
+
+        // Remove the specified cell from the mesh
+        void removeCell
+        (
+            const label cIndex
+        );
+
+        // Insert the specified face to the mesh
+        label insertFace
+        (
+            const label patch,
+            const face& newFace,
+            const label newOwner,
+            const label newNeighbour,
+            const label zoneID = -1
+        );
+
+        // Remove the specified face from the mesh
+        void removeFace
+        (
+            const label fIndex
+        );
+
+        // Split a set of internal faces into boundary faces,
+        // and add them to the specified patch.
+        void splitInternalFaces
+        (
+            const label patchIndex,
+            const labelList& internalFaces,
+            const Map<bool>& cellColors
+        );
+
+        // Insert the specified edge to the mesh
+        label insertEdge
+        (
+            const label patch,
+            const edge& newEdge,
+            const labelList& edgeFaces,
+            const labelList& edgePoints = labelList::null()
+        );
+
+        // Remove the specified edge from the mesh
+        void removeEdge
+        (
+            const label eIndex
+        );
+
+        // Insert the specified point to the mesh
+        label insertPoint
+        (
+            const point& newPoint,
+            const point& oldPoint,
+            const labelList& mapPoints,
+            const label zoneID = -1
+        );
+
+        // Remove the specified point from the mesh
+        void removePoint
+        (
+            const label pIndex
+        );
+
+        // Utility method to build edgePoints for an edge [3D]
+        void buildEdgePoints
+        (
+            const label eIndex,
+            const label checkIndex = 0
+        );
+
+        // Utility to check whether points of an edge lie on a boundary
+        const FixedList<bool,2> checkEdgeBoundary(const label eIndex) const;
+
+        // Set a particular face index as flipped.
+        inline void setFlip(const label fIndex);
+
+        // Utility method to compute the minimum quality of a vertex hull
+        scalar computeMinQuality(const label eIndex) const;
+
+        // Utility method to compute the quality of a vertex hull
+        // around an edge after bisection.
+        scalar computeBisectionQuality(const label eIndex) const;
+
+        // Utility method to compute the quality of cells
+        // around a face after trisection.
+        scalar computeTrisectionQuality(const label fIndex) const;
+
+        // Check whether the given edge is on a bounding curve
+        bool checkBoundingCurve(const label eIndex) const;
+
+        // Allocate dynamic programming tables
+        void initTables
+        (
+            labelList& m,
+            PtrList<scalarListList>& Q,
+            PtrList<labelListList>& K,
+            PtrList<labelListList>& triangulations,
+            const label checkIndex = -1
+        ) const;
+
+        // Check triangulation quality for an edge index
+        bool checkQuality
+        (
+            const label eIndex,
+            const labelList& m,
+            const PtrList<scalarListList>& Q,
+            const scalar minQuality,
+            const label checkIndex = 0
+        ) const;
+
+        // Utility method to fill the dynamic programming tables
+        bool fillTables
+        (
+            const label eIndex,
+            const scalar minQuality,
+            labelList& m,
+            PtrList<scalarListList>& Q,
+            PtrList<labelListList>& K,
+            PtrList<labelListList>& triangulations,
+            const label checkIndex = 0
+        ) const;
+
+        // Print out tables for debugging
+        void printTables
+        (
+            const labelList& m,
+            const PtrList<scalarListList>& Q,
+            const PtrList<labelListList>& K,
+            const label checkIndex = 0
+        ) const;
+
+        // Remove the edge according to the swap sequence
+        const changeMap
+        removeEdgeFlips
+        (
+            const label eIndex,
+            const scalar minQuality,
+            const PtrList<labelListList>& K,
+            PtrList<labelListList>& triangulations,
+            const label checkIndex = 0
+        );
+
+        // Extract triangulations from the programming table
+        void extractTriangulation
+        (
+            const label i,
+            const label j,
+            const labelListList& K,
+            label& numTriangulations,
+            labelListList& triangulations
+        ) const;
+
+        // Identify the 3-2 swap from the triangulation sequence
+        label identify32Swap
+        (
+            const label eIndex,
+            const labelList& hullVertices,
+            const labelListList& triangulations
+        ) const;
+
+        // Check old-volumes for the configuration.
+        bool checkTriangulationVolumes
+        (
+            const label eIndex,
+            const labelList& hullVertices,
+            const labelListList& triangulations
+        ) const;
+
+        // Routine to check whether the triangulation at the
+        // index lies on the boundary of the vertex ring.
+        bool boundaryTriangulation
+        (
+            const label index,
+            label& isolatedVertex,
+            labelListList& triangulations
+        ) const;
+
+        // Routine to perform 2-3 swaps
+        const changeMap
+        swap23
+        (
+            const label isolatedVertex,
+            const label eIndex,
+            const label triangulationIndex,
+            const label numTriangulations,
+            const labelListList& triangulations,
+            const labelList& hullVertices,
+            const labelList& hullFaces,
+            const labelList& hullCells
+        );
+
+        // Routine to perform 3-2 or 2-2 swaps
+        const changeMap
+        swap32
+        (
+            const label eIndex,
+            const label triangulationIndex,
+            const label numTriangulations,
+            const labelListList& triangulations,
+            const labelList& hullVertices,
+            const labelList& hullFaces,
+            const labelList& hullCells
+        );
+
+        // Reorder points after a topology change
+        void reOrderPoints
+        (
+            pointField& points,
+            pointField& preMotionPoints,
+            labelListList& pointZoneMap,
+            bool threaded = false
+        );
+
+        // Static equivalent for multi-threading
+        static void reOrderPointsThread(void *argument);
+
+        // Reorder edges after a topology change
+        void reOrderEdges
+        (
+            edgeList& edges,
+            labelListList& edgeFaces,
+            labelListList& faceEdges,
+            bool threaded = false
+        );
+
+        // Static equivalent for multi-threading
+        static void reOrderEdgesThread(void *argument);
+
+        // Reorder faces in upper-triangular order after a topology change
+        void reOrderFaces
+        (
+            faceList& faces,
+            labelList& owner,
+            labelList& neighbour,
+            labelListList& faceEdges,
+            labelListList& faceZoneFaceMap,
+            bool threaded = false
+        );
+
+        // Static equivalent for multi-threading
+        static void reOrderFacesThread(void *argument);
+
+        // Reorder & renumber cells with bandwidth
+        // reduction after a topology change
+        void reOrderCells
+        (
+            labelListList& cellZoneMap,
+            bool threaded = false
+        );
+
+        // Static equivalent for multi-threading
+        static void reOrderCellsThread(void *argument);
+
+        // Reorder the mesh in upper-triangular order,
+        // and generate mapping information
+        void reOrderMesh
+        (
+            pointField& points,
+            pointField& preMotionPoints,
+            edgeList& edges,
+            faceList& faces,
+            labelList& owner,
+            labelList& neighbour,
+            labelListList& faceEdges,
+            labelListList& edgeFaces,
+            labelListList& pointZoneMap,
+            labelListList& faceZoneFaceMap,
+            labelListList& cellZoneMap
+        );
+
+        // Invoke reOrdering with multiple threads
+        void threadedMeshReOrdering
+        (
+            pointField& points,
+            pointField& preMotionPoints,
+            edgeList& edges,
+            faceList& faces,
+            labelList& owner,
+            labelList& neighbour,
+            labelListList& faceEdges,
+            labelListList& edgeFaces,
+            labelListList& pointZoneMap,
+            labelListList& faceZoneFaceMap,
+            labelListList& cellZoneMap
+        );
+
+        // Reset the mesh and generate mapping information
+        bool resetMesh();
+
+        // Output an entity as a VTK file
+        void writeVTK
+        (
+            const word& name,
+            const label entity,
+            const label primitiveType = 3,
+            const bool useOldConnectivity = false,
+            const bool useOldPoints = false
+        ) const;
+
+        // Output a list of entities as a VTK file
+        void writeVTK
+        (
+            const word& name,
+            const labelList& cList,
+            const label primitiveType = 3,
+            const bool useOldConnectivity = false,
+            const bool useOldPoints = false
+        ) const;
+
+        // Return the status report interval
+        scalar reportInterval() const;
+
+        // Check the state of local connectivity lists
+        void checkConnectivity(const label maxErrors = 0) const;
+
+        // Test an edge / face for proximity with other non-neighbouring faces.
+        // Return the scalar distance to the nearest face.
+        scalar testProximity
+        (
+            const label index,
+            label& proximityFace
+        ) const;
+
+        // Calculate the edge length-scale for the mesh
+        void calculateLengthScale(bool dump = false);
+
+        // Read optional dictionary parameters
+        void readOptionalParameters(bool reRead = false);
+
+        // Dump cell-quality statistics
+        bool meshQuality(bool outputOption);
+
+public:
+
+    // Declare the name of the class and its debug switch
+        TypeName("dynamicTopoFvMesh");
+
+    // Constructors
+
+        //- Construct from existing IOobject
+        explicit dynamicTopoFvMesh(const IOobject& io);
+
+    // Destructor
+
+        virtual ~dynamicTopoFvMesh();
+
+    // Member Functions
+
+        // Map all fields in time using a customized mapper
+        virtual void mapFields(const mapPolyMesh& meshMap) const;
+
+        // Update the mesh for motion / topology changes
+        // Return true if topology changes have occurred
+        virtual bool update();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#include "dynamicTopoFvMeshI.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/dynamicTopoFvMeshI.H

@@ -0,0 +1,818 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    dynamicTopoFvMesh
+
+Description
+    An implementation of dynamic changes to mesh-topology
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*---------------------------------------------------------------------------*/
+
+#include "Stack.H"
+#include "meshOps.H"
+#include "tetPointRef.H"
+#include "linePointRef.H"
+#include "lengthScaleEstimator.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+// Return a non-const reference to the lengthScaleEstimator
+inline lengthScaleEstimator&
+dynamicTopoFvMesh::lengthEstimator()
+{
+    if (!lengthEstimator_.valid())
+    {
+        FatalErrorIn
+        (
+            "inline lengthScaleEstimator& "
+            "dynamicTopoFvMesh::lengthEstimator()"
+        ) << nl
+          << " Invalid request. Length scale estimator was not loaded. "
+          << abort(FatalError);
+    }
+
+    return lengthEstimator_();
+}
+
+
+// Return a const reference to the lengthScaleEstimator
+inline const lengthScaleEstimator&
+dynamicTopoFvMesh::lengthEstimator() const
+{
+    if (!lengthEstimator_.valid())
+    {
+        FatalErrorIn
+        (
+            "inline const lengthScaleEstimator& "
+            "dynamicTopoFvMesh::lengthEstimator() const"
+        ) << nl
+          << " Invalid request. Length scale estimator was not loaded. "
+          << abort(FatalError);
+    }
+
+    return lengthEstimator_();
+}
+
+
+// Return a const reference to the multiThreader
+inline const multiThreader&
+dynamicTopoFvMesh::threader() const
+{
+    if (!threader_.valid())
+    {
+        FatalErrorIn
+        (
+            "inline const multiThreader& "
+            "dynamicTopoFvMesh::threader() const"
+        ) << nl
+          << " Invalid request. multiThreader was not loaded. "
+          << abort(FatalError);
+    }
+
+    return threader_();
+}
+
+
+// Return a reference to the entity mutexes.
+// The index 'entity' ranges from 0 to 3 for point/edge/face/cell.
+inline const Mutex& dynamicTopoFvMesh::entityMutex
+(
+    const label entity
+) const
+{
+    return entityMutex_[entity];
+}
+
+
+// Return the edge index for a provided edge
+inline label dynamicTopoFvMesh::getEdgeIndex
+(
+    const edge& edgeToCheck
+) const
+{
+    if (twoDMesh_)
+    {
+        // No efficient search method. Loop through all edges.
+        forAll(edges_, edgeI)
+        {
+            if (edges_[edgeI] == edgeToCheck)
+            {
+                return edgeI;
+            }
+        }
+    }
+    else
+    {
+        // Look througg pointEdges list
+        const labelList& pEdges = pointEdges_[edgeToCheck.start()];
+
+        forAll(pEdges, edgeI)
+        {
+            if (edges_[pEdges[edgeI]] == edgeToCheck)
+            {
+                return pEdges[edgeI];
+            }
+        }
+    }
+
+    // Probably not a valid edge.
+    FatalErrorIn
+    (
+        "inline label dynamicTopoFvMesh::getEdgeIndex"
+        "(const edge& edgeToCheck) const"
+    )
+        << "Could not find an appropriate edge index for edge:"
+        << edgeToCheck
+        << abort(FatalError);
+
+    return -1;
+}
+
+
+// Given any quad face, pick out a boundary edge that
+// contains a triangular face. For 2D simplical meshes only.
+inline label dynamicTopoFvMesh::getTriBoundaryEdge
+(
+    const label fIndex
+) const
+{
+    const labelList& fEdges = faceEdges_[fIndex];
+
+    forAll(fEdges, edgeI)
+    {
+        // Obtain edgeFaces for this edge
+        const labelList& eFaces = edgeFaces_[fEdges[edgeI]];
+
+        forAll(eFaces, faceI)
+        {
+            if (faces_[eFaces[faceI]].size() == 3)
+            {
+                // Found a triangular face. Return this edge.
+                return fEdges[edgeI];
+            }
+        }
+    }
+
+    // This bit should never happen.
+    FatalErrorIn
+    (
+        "inline label dynamicTopoFvMesh::getTriBoundaryEdge"
+        "(const label fIndex) const"
+    )
+        << "Cannot find a triangular face bordering face: "
+        << fIndex << " :: " << faces_[fIndex]
+        << abort(FatalError);
+
+    return -1;
+}
+
+
+// Return length-scale at an face-location in the mesh [2D]
+inline scalar dynamicTopoFvMesh::faceLengthScale
+(
+    const label fIndex
+) const
+{
+    // Reset the scale first
+    scalar scale = 0.0;
+
+    label facePatch = whichPatch(fIndex);
+
+    // Determine whether the face is internal
+    if (facePatch < 0)
+    {
+#       ifdef FULLDEBUG
+        // Check whether neighbour is valid
+        if (neighbour_[fIndex] == -1)
+        {
+            FatalErrorIn
+            (
+                "inline scalar dynamicTopoFvMesh::faceLengthScale"
+                "(const label fIndex) const"
+            )
+                << nl << "Face: " << fIndex
+                << ": " << faces_[fIndex]
+                << " is not internal."
+                << abort(FatalError);
+        }
+#       endif
+
+        scale =
+        (
+            0.5 *
+            (
+                lengthScale_[owner_[fIndex]]
+              + lengthScale_[neighbour_[fIndex]]
+            )
+        );
+    }
+    else
+    {
+        // Fetch the fixed-length scale
+        scale = lengthEstimator().fixedLengthScale(fIndex, facePatch);
+
+        // If this is a floating face, pick the owner length-scale
+        if (lengthEstimator().isFreePatch(facePatch))
+        {
+            scale = lengthScale_[owner_[fIndex]];
+        }
+
+        // If proximity-based refinement is requested,
+        // test the proximity to the nearest non-neighbour.
+        if (lengthEstimator().isProximityPatch(facePatch))
+        {
+            label proximityFace = -1;
+
+            // Perform a proximity-check.
+            scalar distance = testProximity(fIndex, proximityFace);
+
+            if (debug > 3 && self() == 0)
+            {
+                if
+                (
+                    (proximityFace > -1) &&
+                    ((distance / 5.0) < scale)
+                )
+                {
+                    Info << " Closest opposing face detected for face: " << nl
+                         << '\t' << fIndex
+                         << " :: " << faces_[fIndex]
+                         << " was face:\n"
+                         << '\t' << proximityFace
+                         << " :: " << polyMesh::faces()[proximityFace] << nl
+                         << " with distance: " << distance
+                         << endl;
+                }
+            }
+
+            scale =
+            (
+                Foam::min
+                (
+                    scale,
+                    ((distance / 3.0) - SMALL)/lengthEstimator().ratioMax()
+                )
+            );
+        }
+
+        // Limit scales if necessary
+        lengthEstimator().limitScale(scale);
+    }
+
+    return scale;
+}
+
+
+// Compute length-scale at an edge-location in the mesh [3D]
+inline scalar dynamicTopoFvMesh::edgeLengthScale
+(
+    const label eIndex
+) const
+{
+    // Reset the scale first
+    scalar scale = 0.0;
+
+    const labelList& eFaces = edgeFaces_[eIndex];
+
+    label edgePatch = whichEdgePatch(eIndex);
+
+    // Determine whether the edge is internal
+    if (edgePatch < 0)
+    {
+        forAll(eFaces, faceI)
+        {
+#           ifdef FULLDEBUG
+            // Check whether neighbour is valid
+            if (neighbour_[eFaces[faceI]] == -1)
+            {
+                FatalErrorIn
+                (
+                    "inline scalar dynamicTopoFvMesh::edgeLengthScale"
+                    "(const label eIndex) const"
+                )
+                    << nl << "Face: " << eFaces[faceI]
+                    << ": " << faces_[eFaces[faceI]]
+                    << " is not internal, while edge: "
+                    << eIndex << ": " << edges_[eIndex] << " is."
+                    << abort(FatalError);
+            }
+#           endif
+
+            scale += lengthScale_[owner_[eFaces[faceI]]];
+            scale += lengthScale_[neighbour_[eFaces[faceI]]];
+        }
+
+        scale /= (2.0*eFaces.size());
+    }
+    else
+    {
+        // Search for boundary faces, and average their scale
+        forAll(eFaces, faceI)
+        {
+            if (neighbour_[eFaces[faceI]] == -1)
+            {
+                scale +=
+                (
+                    lengthEstimator().fixedLengthScale
+                    (
+                        eFaces[faceI],
+                        edgePatch
+                    )
+                );
+            }
+        }
+
+        scale *= 0.5;
+
+        // If proximity-based refinement is requested,
+        // test the proximity to the nearest non-neighbour.
+        if (lengthEstimator().isProximityPatch(edgePatch))
+        {
+            label proximityFace = -1;
+
+            // Perform a proximity-check.
+            scalar distance = testProximity(eIndex, proximityFace);
+
+            if (debug > 3 && self() == 0)
+            {
+                if
+                (
+                    (proximityFace > -1) &&
+                    ((distance / 5.0) < scale)
+                )
+                {
+                    Info << " Closest opposing face detected for edge: " << nl
+                         << '\t' << eIndex
+                         << " :: " << edges_[eIndex]
+                         << " was face:\n"
+                         << '\t' << proximityFace
+                         << " :: " << polyMesh::faces()[proximityFace] << nl
+                         << " with distance: " << distance
+                         << endl;
+                }
+            }
+
+            scale =
+            (
+                Foam::min
+                (
+                    scale,
+                    ((distance / 3.0) - SMALL)/lengthEstimator().ratioMax()
+                )
+            );
+        }
+
+        // If curvature-based refinement is requested,
+        // test the variation in face-normal directions.
+        if (lengthEstimator().isCurvaturePatch(edgePatch))
+        {
+            // Obtain face-normals for both faces.
+            label count = 0;
+            FixedList<vector, 2> fNorm;
+
+            forAll(eFaces, faceI)
+            {
+                if (neighbour_[eFaces[faceI]] == -1)
+                {
+                    // Obtain the normal.
+                    fNorm[count] = faces_[eFaces[faceI]].normal(points_);
+
+                    // Normalize it.
+                    fNorm[count] /= mag(fNorm[count]);
+
+                    count++;
+                }
+            }
+
+            scalar deviation = (fNorm[0] & fNorm[1]);
+            scalar refDeviation = lengthEstimator().curvatureDeviation();
+
+            if (mag(deviation) < refDeviation)
+            {
+                // Fetch the edge
+                const edge& edgeToCheck = edges_[eIndex];
+
+                // Get the edge-length.
+                scalar length =
+                (
+                    linePointRef
+                    (
+                        points_[edgeToCheck.start()],
+                        points_[edgeToCheck.end()]
+                    ).mag()
+                );
+
+                if (debug > 3 && self() == 0)
+                {
+                    Info << "Deviation: " << deviation << nl
+                         << "curvatureDeviation: " << refDeviation
+                         << ", Edge: " << eIndex << ", Length: " << length
+                         << ", Scale: " << scale << nl
+                         << " Half-length: " << (0.5*length) << nl
+                         << " MinRatio: "
+                         << (lengthEstimator().ratioMin()*scale)
+                         << endl;
+                }
+
+                scale =
+                (
+                    Foam::min
+                    (
+                        scale,
+                        ((length - SMALL)/lengthEstimator().ratioMax())
+                    )
+                );
+            }
+        }
+
+        // Limit scales if necessary
+        lengthEstimator().limitScale(scale);
+    }
+
+    return scale;
+}
+
+
+// Check for edge bisection
+inline bool dynamicTopoFvMesh::checkBisection
+(
+    const label index
+) const
+{
+    scalar length = 0.0, scale = 0.0;
+
+    if (twoDMesh_)
+    {
+        // Fetch the edge
+        const edge& edgeToCheck = edges_[getTriBoundaryEdge(index)];
+
+        // Measure the boundary edge-length of the face in question
+        length =
+        (
+            linePointRef
+            (
+                points_[edgeToCheck.start()],
+                points_[edgeToCheck.end()]
+            ).mag()
+        );
+
+        // Determine the length-scale at this face
+        scale = faceLengthScale(index);
+    }
+    else
+    {
+        // Fetch the edge
+        const edge& edgeToCheck = edges_[index];
+
+        // Measure the edge-length
+        length =
+        (
+            linePointRef
+            (
+                points_[edgeToCheck.start()],
+                points_[edgeToCheck.end()]
+            ).mag()
+        );
+
+        // Determine the length-scale at this point in the mesh
+        scale = edgeLengthScale(index);
+    }
+
+    if (length > lengthEstimator().ratioMax() * scale)
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// Check for edge collapse
+inline bool dynamicTopoFvMesh::checkCollapse
+(
+    const label index
+) const
+{
+    scalar length = 0.0, scale = 0.0;
+
+    if (twoDMesh_)
+    {
+        // Fetch the edge
+        const edge& edgeToCheck = edges_[getTriBoundaryEdge(index)];
+
+        // Measure the boundary edge-length of the face in question
+        length =
+        (
+            linePointRef
+            (
+                points_[edgeToCheck.start()],
+                points_[edgeToCheck.end()]
+            ).mag()
+        );
+
+        // Determine the length-scale at this face
+        scale = faceLengthScale(index);
+    }
+    else
+    {
+        // Fetch the edge
+        const edge& edgeToCheck = edges_[index];
+
+        // Measure the edge-length
+        length =
+        (
+            linePointRef
+            (
+                points_[edgeToCheck.start()],
+                points_[edgeToCheck.end()]
+            ).mag()
+        );
+
+        // Determine the length-scale at this point in the mesh
+        scale = edgeLengthScale(index);
+    }
+
+    if (length < lengthEstimator().ratioMin() * scale)
+    {
+        return true;
+    }
+    else
+    {
+        return false;
+    }
+}
+
+
+// Return the entity stack
+inline Stack& dynamicTopoFvMesh::stack
+(
+    const label threadID
+)
+{
+    return entityStack_[threadID];
+}
+
+
+// Return the integer ID for a given thread
+// Return zero for single-threaded operation
+inline label dynamicTopoFvMesh::self() const
+{
+    if (threader_->multiThreaded())
+    {
+        for (label i = 1; i <= threader_->getNumThreads(); i++)
+        {
+            if (handlerPtr_[i].self())
+            {
+                return i;
+            }
+        }
+    }
+
+    return 0;
+}
+
+
+// Initialize edge-stacks
+inline void dynamicTopoFvMesh::initStacks
+(
+    const labelHashSet& entities
+)
+{
+    forAll(entityStack_, stackI)
+    {
+        entityStack_[stackI].clear();
+    }
+
+    // Prepare a filling sequence based on threading operation
+    label tIndex = 0;
+    labelList tID(threader().getNumThreads());
+
+    if (threader_->multiThreaded())
+    {
+        forAll(tID, tI)
+        {
+            tID[tI] = (tI + 1);
+        }
+    }
+    else
+    {
+        tID = 0;
+    }
+
+    if (twoDMesh_)
+    {
+        forAll(faces_, faceI)
+        {
+            if (faces_[faceI].size() == 4)
+            {
+                stack(tID[tIndex]).insert(faceI);
+
+                tIndex = tID.fcIndex(tIndex);
+            }
+        }
+    }
+    else
+    {
+        forAll(edges_, edgeI)
+        {
+            if (edgeFaces_[edgeI].size())
+            {
+                stack(tID[tIndex]).insert(edgeI);
+
+                tIndex = tID.fcIndex(tIndex);
+            }
+        }
+    }
+}
+
+
+// Method to determine the old boundary patch index for a given face
+// Similar to the polyBoundaryMesh routine, but works on local information
+inline label dynamicTopoFvMesh::whichPatch
+(
+    const label index
+) const
+{
+    if (index < nOldInternalFaces_)
+    {
+        return -1;
+    }
+
+    for (label i = 0; i < boundaryMesh().size(); i++)
+    {
+        if
+        (
+            index >= oldPatchStarts_[i]
+         && index < oldPatchStarts_[i] + oldPatchSizes_[i]
+        )
+        {
+            return i;
+        }
+    }
+
+    // If not in any of the above, it's possible that the face was added
+    // at the end of the list. Check addedFacePatches_ for the patch info
+    if (addedFacePatches_.found(index))
+    {
+        return addedFacePatches_[index];
+    }
+    else
+    {
+        FatalErrorIn
+        (
+            "label dynamicTopoFvMesh::whichPatch"
+            "(const label index) const"
+        )
+            << "Cannot find patch information for face index: "
+            << index << nl
+            << " It appears that face ordering is"
+            << " inconsistent with patch information." << nl
+            << " Mesh info: " << nl
+            << "    nOldInternalFaces: " << nOldInternalFaces_ << nl
+            << "    oldPatchStarts: " << oldPatchStarts_ << nl
+            << "    oldPatchSizes: " << oldPatchSizes_ << nl
+            << abort(FatalError);
+    }
+
+    return -2;
+}
+
+
+// Method to determine the old boundary patch index for a given edge
+inline label dynamicTopoFvMesh::whichEdgePatch
+(
+    const label index
+) const
+{
+    if (index < nOldInternalEdges_)
+    {
+        return -1;
+    }
+
+    for (label i = 0; i < boundaryMesh().size(); i++)
+    {
+        if
+        (
+            index >= oldEdgePatchStarts_[i]
+         && index < oldEdgePatchStarts_[i] + oldEdgePatchSizes_[i]
+        )
+        {
+            return i;
+        }
+    }
+
+    // If not in any of the above, it's possible that the edge was added
+    // at the end of the list. Check addedEdgePatches_ for the patch info
+    if (addedEdgePatches_.found(index))
+    {
+        return addedEdgePatches_[index];
+    }
+    else
+    {
+        FatalErrorIn
+        (
+            "label dynamicTopoFvMesh::whichEdgePatch"
+            "(const label index) const"
+        )
+            << "Cannot find patch information for edge index: "
+            << index << nl
+            << " It appears that edge ordering is"
+            << " inconsistent with patch information." << nl
+            << " Mesh info: " << nl
+            << "    nOldInternalEdges: " << nOldInternalEdges_ << nl
+            << "    oldEdgePatchStarts: " << oldEdgePatchStarts_ << nl
+            << "    oldEdgePatchSizes: " << oldEdgePatchSizes_ << nl
+            << abort(FatalError);
+    }
+
+    return -2;
+}
+
+
+// Report topo-modification status
+//  - Enumerants are as follows
+//    [0] nModifications
+//    [1] Total swaps (internal and surface)
+//    [2] Surface swaps
+//    [3] Total bisections
+//    [4] Total collapses
+//    [5] Surface bisections
+//    [6] Surface collapses
+//    [7] Slivers removed
+inline label dynamicTopoFvMesh::status(const label type) const
+{
+    if (type < statistics_.size())
+    {
+        return statistics_[type];
+    }
+    else
+    {
+        FatalErrorIn
+        (
+            "inline label dynamicTopoFvMesh::status"
+            "(const label type) const"
+        )
+            << " Unknown type: " << type << nl
+            << abort(FatalError);
+    }
+
+    return 0;
+}
+
+
+// Set a particular face index as flipped.
+inline void dynamicTopoFvMesh::setFlip(const label fIndex)
+{
+    if (fIndex < nOldFaces_)
+    {
+        if (flipFaces_.found(fIndex))
+        {
+            flipFaces_.erase(fIndex);
+        }
+        else
+        {
+            flipFaces_.insert(fIndex);
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/dynamicTopoFvMeshMapping.C

@@ -0,0 +1,597 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    dynamicTopoFvMesh
+
+Description
+    Functions specific to conservative mapping
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*----------------------------------------------------------------------------*/
+
+#include "dynamicTopoFvMesh.H"
+
+#include "meshOps.H"
+#include "IOmanip.H"
+#include "triFace.H"
+#include "objectMap.H"
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+// Compute mapping weights for modified entities
+void dynamicTopoFvMesh::computeMapping
+(
+    const scalar matchTol,
+    const bool skipMapping,
+    const label faceStart,
+    const label faceSize,
+    const label cellStart,
+    const label cellSize
+)
+{
+    // Compute cell mapping
+    for (label cellI = cellStart; cellI < (cellStart + cellSize); cellI++)
+    {
+        label cIndex = cellsFromCells_[cellI].index();
+
+        if (skipMapping)
+        {
+            // Set empty mapping parameters
+            const labelList& mo = cellParents_[cIndex];
+
+            cellsFromCells_[cellI].masterObjects() = mo;
+            cellWeights_[cellI].setSize(mo.size(), (1.0/(mo.size() + VSMALL)));
+            cellCentres_[cellI].setSize(mo.size(), vector::zero);
+        }
+        else
+        {
+            // Compute master objects for inverse-distance weighting
+            computeParents
+            (
+                cIndex,
+                cellParents_[cIndex],
+                polyMesh::cellCells(),
+                polyMesh::cellCentres(),
+                3,
+                cellsFromCells_[cellI].masterObjects()
+            );
+        }
+    }
+
+    // Compute face mapping
+    for (label faceI = faceStart; faceI < (faceStart + faceSize); faceI++)
+    {
+        label fIndex = facesFromFaces_[faceI].index();
+        label patchIndex = whichPatch(fIndex);
+
+        // Skip mapping for internal faces.
+        if (patchIndex == -1)
+        {
+            // Set dummy masters, so that the conventional
+            // faceMapper doesn't incur a seg-fault.
+            facesFromFaces_[faceI].masterObjects() = labelList(1, 0);
+            continue;
+        }
+
+        if (skipMapping)
+        {
+            // Set empty mapping parameters
+            const labelList& mo = faceParents_[fIndex];
+
+            facesFromFaces_[faceI].masterObjects() = mo;
+            faceWeights_[faceI].setSize(mo.size(), (1.0/(mo.size() + VSMALL)));
+            faceCentres_[faceI].setSize(mo.size(), vector::zero);
+        }
+        else
+        {
+            // Compute master objects for inverse-distance weighting
+            computeParents
+            (
+                fIndex,
+                faceParents_[fIndex],
+                boundaryMesh()[patchIndex].faceFaces(),
+                boundaryMesh()[patchIndex].faceCentres(),
+                2,
+                facesFromFaces_[faceI].masterObjects()
+            );
+        }
+    }
+}
+
+
+// Static equivalent for multiThreading
+void dynamicTopoFvMesh::computeMappingThread(void *argument)
+{
+    // Recast the argument
+    meshHandler *thread = static_cast<meshHandler*>(argument);
+
+    if (thread->slave())
+    {
+        thread->sendSignal(meshHandler::START);
+    }
+
+    dynamicTopoFvMesh& mesh = thread->reference();
+
+    // Recast the pointers for the argument
+    scalar& matchTol  = *(static_cast<scalar*>(thread->operator()(0)));
+    bool& skipMapping = *(static_cast<bool*>(thread->operator()(1)));
+    label& faceStart  = *(static_cast<label*>(thread->operator()(2)));
+    label& faceSize   = *(static_cast<label*>(thread->operator()(3)));
+    label& cellStart  = *(static_cast<label*>(thread->operator()(4)));
+    label& cellSize   = *(static_cast<label*>(thread->operator()(5)));
+
+    // Now calculate addressing
+    mesh.computeMapping
+    (
+        matchTol,
+        skipMapping,
+        faceStart, faceSize,
+        cellStart, cellSize
+    );
+
+    if (thread->slave())
+    {
+        thread->sendSignal(meshHandler::STOP);
+    }
+}
+
+
+// Routine to invoke threaded mapping
+void dynamicTopoFvMesh::threadedMapping
+(
+    scalar matchTol,
+    bool skipMapping
+)
+{
+    label nThreads = threader_->getNumThreads();
+
+    // If mapping is being skipped, issue a warning.
+    if (skipMapping)
+    {
+        Info << " *** Mapping is being skipped *** " << endl;
+    }
+
+    // Check if single-threaded
+    if (nThreads == 1)
+    {
+        computeMapping
+        (
+            matchTol,
+            skipMapping,
+            0, facesFromFaces_.size(),
+            0, cellsFromCells_.size()
+        );
+
+        return;
+    }
+
+    // Set one handler per thread
+    PtrList<meshHandler> hdl(nThreads);
+
+    forAll(hdl, i)
+    {
+        hdl.set(i, new meshHandler(*this, threader()));
+    }
+
+    // Simple load-balancing scheme
+    FixedList<label, 2> index(-1);
+    FixedList<labelList, 2> tStarts(labelList(nThreads, 0));
+    FixedList<labelList, 2> tSizes(labelList(nThreads, 0));
+
+    index[0] = facesFromFaces_.size();
+    index[1] = cellsFromCells_.size();
+
+    if (debug > 2)
+    {
+        Info << " Mapping Faces: " << index[0] << endl;
+        Info << " Mapping Cells: " << index[1] << endl;
+    }
+
+    forAll(index, indexI)
+    {
+        label j = 0, total = 0;
+
+        while (index[indexI]--)
+        {
+            tSizes[indexI][(j = tSizes[indexI].fcIndex(j))]++;
+        }
+
+        for (label i = 1; i < tStarts[indexI].size(); i++)
+        {
+            tStarts[indexI][i] = tSizes[indexI][i-1] + total;
+
+            total += tSizes[indexI][i-1];
+        }
+
+        if (debug > 2)
+        {
+            Info << " Load starts: " << tStarts[indexI] << endl;
+            Info << " Load sizes: " << tSizes[indexI] << endl;
+        }
+    }
+
+    // Set the argument list for each thread
+    forAll(hdl, i)
+    {
+        // Size up the argument list
+        hdl[i].setSize(6);
+
+        // Set match tolerance
+        hdl[i].set(0, &matchTol);
+
+        // Set the skipMapping flag
+        hdl[i].set(1, &skipMapping);
+
+        // Set the start/size indices
+        hdl[i].set(2, &(tStarts[0][i]));
+        hdl[i].set(3, &(tSizes[0][i]));
+        hdl[i].set(4, &(tStarts[1][i]));
+        hdl[i].set(5, &(tSizes[1][i]));
+    }
+
+    // Prior to multi-threaded operation,
+    // force calculation of demand-driven data.
+    polyMesh::cells();
+    primitiveMesh::cellCells();
+    primitiveMesh::cellCentres();
+
+    const polyBoundaryMesh& boundary = boundaryMesh();
+
+    forAll(boundary, patchI)
+    {
+        boundary[patchI].faceFaces();
+        boundary[patchI].faceCentres();
+    }
+
+    // Execute threads in linear sequence
+    executeThreads(identity(nThreads), hdl, &computeMappingThread);
+}
+
+
+// Compute parents for inverse-distance weighting
+void dynamicTopoFvMesh::computeParents
+(
+    const label index,
+    const labelList& mapCandidates,
+    const labelListList& oldNeighbourList,
+    const vectorField& oldCentres,
+    const label dimension,
+    labelList& parents
+) const
+{
+    if (parents.size())
+    {
+        FatalErrorIn
+        (
+            "\n\n"
+            "void dynamicTopoFvMesh::computeParents\n"
+            "(\n"
+            "    const label index,\n"
+            "    const labelList& mapCandidates,\n"
+            "    const labelListList& oldNeighbourList,\n"
+            "    const vectorField& oldCentres,\n"
+            "    const label dimension,\n"
+            "    labelList& parents\n"
+            ") const\n"
+        )
+            << " Addressing has already been calculated." << nl
+            << " Index: " << index << nl
+            << " Type: " << (dimension == 2 ? "Face" : "Cell") << nl
+            << " mapCandidates: " << mapCandidates << nl
+            << " Parents: " << parents << nl
+            << abort(FatalError);
+    }
+
+    // Figure out the patch offset and centre
+    label offset = -1;
+    vector centre = vector::zero;
+
+    if (dimension == 2)
+    {
+        offset = boundaryMesh()[whichPatch(index)].start();
+
+        centre = faces_[index].centre(oldPoints_);
+    }
+    else
+    if (dimension == 3)
+    {
+        offset = 0;
+        scalar dummyVol = 0.0;
+
+        meshOps::cellCentreAndVolume
+        (
+            index,
+            oldPoints_,
+            faces_,
+            cells_,
+            owner_,
+            centre,
+            dummyVol
+        );
+    }
+
+    // Maintain a check-list
+    labelHashSet checked;
+
+    // Insert candidates first
+    forAll(mapCandidates, indexI)
+    {
+        checked.insert(mapCandidates[indexI] - offset);
+    }
+
+    // Loop for three outward levels from candidates
+    for (label i = 0; i < 3; i++)
+    {
+        // Fetch the set of candidates
+        const labelList checkList = checked.toc();
+
+        forAll(checkList, indexI)
+        {
+            const labelList& oldNei = oldNeighbourList[checkList[indexI]];
+
+            forAll(oldNei, entityI)
+            {
+                if (!checked.found(oldNei[entityI]))
+                {
+                    checked.insert(oldNei[entityI]);
+                }
+            }
+        }
+    }
+
+    // Loop through accumulated candidates and fetch the nearest one.
+    scalar minDist = GREAT;
+    label minLocation = -1;
+
+    const labelList checkList = checked.toc();
+
+    forAll(checkList, indexI)
+    {
+        scalar dist = magSqr(oldCentres[checkList[indexI]] - centre);
+
+        if (dist < minDist)
+        {
+            minDist = dist;
+            minLocation = checkList[indexI];
+        }
+    }
+
+    // Set the final list
+    const labelList& neiList = oldNeighbourList[minLocation];
+
+    parents.setSize(neiList.size() + 1, -1);
+
+    // Fill indices
+    forAll(neiList, indexI)
+    {
+        parents[indexI] = neiList[indexI] + offset;
+    }
+
+    // Set final index
+    parents[neiList.size()] = minLocation + offset;
+}
+
+
+// Set fill-in mapping information for a particular cell
+void dynamicTopoFvMesh::setCellMapping
+(
+    const label cIndex,
+    const labelList& mapCells,
+    bool addEntry
+)
+{
+    if (addEntry)
+    {
+        if (debug > 3)
+        {
+            Info << "Inserting mapping cell: " << cIndex << nl
+                 << " mapCells: " << mapCells
+                 << endl;
+        }
+
+        // Insert index into the list, and overwrite if necessary
+        label index = -1;
+
+        forAll(cellsFromCells_, indexI)
+        {
+            if (cellsFromCells_[indexI].index() == cIndex)
+            {
+                index = indexI;
+                break;
+            }
+        }
+
+        if (index == -1)
+        {
+            meshOps::sizeUpList
+            (
+                objectMap(cIndex, labelList(0)),
+                cellsFromCells_
+            );
+        }
+        else
+        {
+            cellsFromCells_[index].masterObjects() = labelList(0);
+        }
+    }
+
+    // Update cell-parents information
+    labelHashSet masterCells;
+
+    forAll(mapCells, cellI)
+    {
+        if (mapCells[cellI] < 0)
+        {
+            continue;
+        }
+
+        if (mapCells[cellI] < nOldCells_)
+        {
+            masterCells.insert(mapCells[cellI]);
+        }
+        else
+        if (cellParents_.found(mapCells[cellI]))
+        {
+            const labelList& nParents = cellParents_[mapCells[cellI]];
+
+            forAll(nParents, cI)
+            {
+                masterCells.insert(nParents[cI]);
+            }
+        }
+    }
+
+    cellParents_.set(cIndex, masterCells.toc());
+}
+
+
+// Set fill-in mapping information for a particular face
+void dynamicTopoFvMesh::setFaceMapping
+(
+    const label fIndex,
+    const labelList& mapFaces
+)
+{
+    label patch = whichPatch(fIndex);
+
+    if (debug > 3)
+    {
+        Info << "Inserting mapping face: " << fIndex << nl
+             << " patch: " << patch << nl
+             << " mapFaces: " << mapFaces
+             << endl;
+    }
+
+    bool foundError = false;
+
+    // Check to ensure that internal faces are not mapped
+    // from any faces in the mesh
+    if (patch == -1 && mapFaces.size())
+    {
+        foundError = true;
+    }
+
+    // Check to ensure that boundary faces map
+    // only from other faces on the same patch
+    if (patch > -1 && mapFaces.empty())
+    {
+        foundError = true;
+    }
+
+    if (foundError)
+    {
+        FatalErrorIn
+        (
+            "\n"
+            "void dynamicTopoFvMesh::setFaceMapping\n"
+            "(\n"
+            "    const label fIndex,\n"
+            "    const labelList& mapFaces\n"
+            ")"
+        )
+            << nl << " Incompatible mapping. " << nl
+            << "  Possible reasons: " << nl
+            << "    1. No mapping specified for a boundary face; " << nl
+            << "    2. Mapping specified for an internal face, " << nl
+            << "       when none was expected." << nl << nl
+            << " Face: " << fIndex << nl
+            << " Patch: " << patch << nl
+            << " Owner: " << owner_[fIndex] << nl
+            << " Neighbour: " << neighbour_[fIndex] << nl
+            << " mapFaces: " << mapFaces << nl
+            << abort(FatalError);
+    }
+
+    // Insert addressing into the list, and overwrite if necessary
+    label index = -1;
+
+    forAll(facesFromFaces_, indexI)
+    {
+        if (facesFromFaces_[indexI].index() == fIndex)
+        {
+            index = indexI;
+            break;
+        }
+    }
+
+    if (index == -1)
+    {
+        meshOps::sizeUpList
+        (
+            objectMap(fIndex, labelList(0)),
+            facesFromFaces_
+        );
+    }
+    else
+    {
+        facesFromFaces_[index].masterObjects() = labelList(0);
+    }
+
+    // For internal faces, set dummy maps / weights, and bail out
+    if (patch == -1)
+    {
+        return;
+    }
+
+    // Update face-parents information
+    labelHashSet masterFaces;
+
+    forAll(mapFaces, faceI)
+    {
+        if (mapFaces[faceI] < 0)
+        {
+            continue;
+        }
+
+        if (mapFaces[faceI] < nOldFaces_)
+        {
+            masterFaces.insert(mapFaces[faceI]);
+        }
+        else
+        if (faceParents_.found(mapFaces[faceI]))
+        {
+            const labelList& nParents = faceParents_[mapFaces[faceI]];
+
+            forAll(nParents, fI)
+            {
+                masterFaces.insert(nParents[fI]);
+            }
+        }
+    }
+
+    faceParents_.set(fIndex, masterFaces.toc());
+}
+
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/eBoundaryMesh/eBoundaryMesh.C

@@ -0,0 +1,295 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+
+\*---------------------------------------------------------------------------*/
+
+#include "eBoundaryMesh.H"
+#include "eMesh.H"
+#include "primitiveMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(eBoundaryMesh, 0);
+
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+// Construct from dictionary
+eBoundaryMesh::eBoundaryMesh
+(
+    const IOobject& io,
+    const eMesh& mesh
+)
+:
+    ePatchList(),
+    regIOobject(io),
+    mesh_(mesh)
+{
+    if
+    (
+        (readOpt() == IOobject::MUST_READ) ||
+        (readOpt() == IOobject::READ_IF_PRESENT && headerOk())
+    )
+    {
+        readFromInputStream();
+    }
+}
+
+
+// Construct given size. Patches will be set later
+eBoundaryMesh::eBoundaryMesh
+(
+    const IOobject& io,
+    const eMesh& em,
+    const label size
+)
+:
+    ePatchList(size),
+    regIOobject(io),
+    mesh_(em)
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+//- Read from input stream
+void eBoundaryMesh::readFromInputStream()
+{
+    ePatchList& patches = *this;
+
+    // Read polyPatchList
+    Istream& is = readStream(typeName);
+
+    PtrList<entry> patchEntries(is);
+    patches.setSize(patchEntries.size());
+
+    forAll(patches, patchI)
+    {
+        patches.set
+        (
+            patchI,
+            ePatch::New
+            (
+                patchEntries[patchI].keyword(),
+                patchEntries[patchI].dict(),
+                patchI,
+                *this
+            )
+        );
+    }
+
+    // Check state of IOstream
+    is.check
+    (
+        "eBoundaryMesh::eBoundaryMesh"
+        "(const IOobject&, const faMesh&)"
+    );
+
+    close();
+}
+
+// Return the mesh reference
+const eMesh& eBoundaryMesh::mesh() const
+{
+    return mesh_;
+}
+
+
+// Return a list of patch types
+wordList eBoundaryMesh::types() const
+{
+    const ePatchList& patches = *this;
+
+    wordList t(patches.size());
+
+    forAll (patches, patchI)
+    {
+        t[patchI] = patches[patchI].type();
+    }
+
+    return t;
+}
+
+
+// Return a list of patch names
+wordList eBoundaryMesh::names() const
+{
+    const ePatchList& patches = *this;
+
+    wordList t(patches.size());
+
+    forAll (patches, patchI)
+    {
+        t[patchI] = patches[patchI].name();
+    }
+
+    return t;
+}
+
+
+// Return a list of patch sizes
+labelList eBoundaryMesh::patchSizes() const
+{
+    const ePatchList& patches = *this;
+
+    labelList t(patches.size());
+
+    forAll (patches, patchI)
+    {
+        t[patchI] = patches[patchI].size();
+    }
+
+    return t;
+}
+
+
+// Return a list of patch starts
+labelList eBoundaryMesh::patchStarts() const
+{
+    const ePatchList& patches = *this;
+
+    labelList t(patches.size());
+
+    forAll (patches, patchI)
+    {
+        t[patchI] = patches[patchI].start();
+    }
+
+    return t;
+}
+
+
+//- Return patch index for a given edge label
+label eBoundaryMesh::whichPatch(const label edgeIndex) const
+{
+    // Find out which patch the current edge belongs to by comparing label
+    // with patch start labels.
+    // If the face is internal, return -1;
+    // if it is off the end of the list, abort
+    if (edgeIndex >= mesh().nEdges())
+    {
+        FatalErrorIn
+        (
+            "eBoundaryMesh::whichPatch(const label edgeIndex) const"
+        )   << "given label greater than the number of geometric edges"
+            << abort(FatalError);
+    }
+
+    if (edgeIndex < mesh().nInternalEdges())
+    {
+        return -1;
+    }
+
+    forAll (*this, patchI)
+    {
+        const ePatch& bp = operator[](patchI);
+
+        if
+        (
+            edgeIndex >= bp.start()
+         && edgeIndex < bp.start() + bp.size()
+        )
+        {
+            return patchI;
+        }
+    }
+
+    // If not in any of above, it is trouble!
+    FatalErrorIn
+    (
+        "label eBoundaryMesh::whichPatch(const label edgeIndex) const"
+    )   << "Cannot find edge " << edgeIndex << " in any of the patches "
+        << names() << nl
+        << "It seems your patches are not consistent with the mesh :"
+        << " internalEdges:" << mesh().nInternalEdges()
+        << "  total number of edges:" << mesh().nEdges()
+        << abort(FatalError);
+
+    return -1;
+}
+
+
+label eBoundaryMesh::findPatchID(const word& patchName) const
+{
+    const ePatchList& patches = *this;
+
+    forAll (patches, patchI)
+    {
+        if (patches[patchI].name() == patchName)
+        {
+            return patchI;
+        }
+    }
+
+    // Patch not found
+    return -1;
+}
+
+
+// writeData member function required by regIOobject
+bool eBoundaryMesh::writeData(Ostream& os) const
+{
+    const ePatchList& patches = *this;
+
+    os  << patches.size() << nl << token::BEGIN_LIST << incrIndent << nl;
+
+    forAll(patches, patchi)
+    {
+        os  << indent << patches[patchi].name() << nl
+            << indent << token::BEGIN_BLOCK << nl
+            << incrIndent << patches[patchi] << decrIndent
+            << indent << token::END_BLOCK << endl;
+    }
+
+    os  << decrIndent << token::END_LIST;
+
+    // Check state of IOstream
+    os.check("eBoundaryMesh::writeData(Ostream& os) const");
+
+    return os.good();
+}
+
+
+Ostream& operator<<(Ostream& os, const eBoundaryMesh& bm)
+{
+    bm.writeData(os);
+    return os;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/eBoundaryMesh/eBoundaryMesh.H

@@ -0,0 +1,144 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    eBoundaryMesh
+
+SourceFiles
+    eBoundaryMesh.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef eBoundaryMesh_H
+#define eBoundaryMesh_H
+
+#include "ePatchList.H"
+#include "wordList.H"
+#include "pointField.H"
+#include "regIOobject.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+class eMesh;
+
+/*---------------------------------------------------------------------------*\
+                      Class eBoundaryMesh Declaration
+\*---------------------------------------------------------------------------*/
+
+class eBoundaryMesh
+:
+    public ePatchList,
+    public regIOobject
+{
+    // private data
+
+        //- Reference to mesh
+        const eMesh& mesh_;
+
+        //- Disallow construct as copy
+        eBoundaryMesh(const eBoundaryMesh&);
+
+        //- Disallow assignment
+        void operator=(const eBoundaryMesh&);
+
+public:
+
+    //- Runtime type information
+    TypeName("eBoundaryMesh");
+
+
+    // Constructors
+
+        //- Construct from dictionary
+        eBoundaryMesh
+        (
+            const IOobject& io,
+            const eMesh& mesh
+        );
+
+        //- Construct given size
+        eBoundaryMesh
+        (
+            const IOobject& io,
+            const eMesh& mesh,
+            const label size
+        );
+
+
+    // Destructor - default
+
+
+    // Member functions
+
+        // Access
+
+            //- Return the mesh reference
+            const eMesh& mesh() const;
+
+            //- Return a list of patch types
+            wordList types() const;
+
+            //- Return a list of patch names
+            wordList names() const;
+
+            //- Return a list of patch sizes
+            labelList patchSizes() const;
+
+            //- Return a list of patch starts
+            labelList patchStarts() const;
+
+            //- Return patch index for a given edge label
+            label whichPatch(const label edgeIndex) const;
+
+            //- Find patch index given a name
+            label findPatchID(const word& patchName) const;
+
+            //- writeData member function required by regIOobject
+            bool writeData(Ostream&) const;
+
+        // Input / Output
+
+            //- Read from input stream
+            void readFromInputStream();
+
+    // Ostream operator
+
+        friend Ostream& operator<<(Ostream&, const eBoundaryMesh&);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/eMesh.C

@@ -0,0 +1,482 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+
+    Mesh needed to do edge-based addressing.
+
+\*---------------------------------------------------------------------------*/
+
+#include "eMesh.H"
+#include "demandDrivenData.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(eMesh, 0);
+
+word eMesh::meshSubDir = "eMesh";
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+void eMesh::clearGeom() const
+{
+    if (debug)
+    {
+        Info<< "void eMesh::clearGeom() const : "
+            << "Clearing geometry" << endl;
+    }
+}
+
+
+void eMesh::clearAddressing() const
+{
+    if (debug)
+    {
+        Info<< "void eMesh::clearAddressing() const : "
+            << "Clearing addressing" << endl;
+    }
+
+    edges_.clear();
+
+    deleteDemandDrivenData(pePtr_);
+    deleteDemandDrivenData(epPtr_);
+    deleteDemandDrivenData(fePtr_);
+    deleteDemandDrivenData(efPtr_);
+}
+
+
+// Helper function to isolate points on triangular faces
+label eMesh::findIsolatedPoint(const face& f, const edge& e) const
+{
+    // Check the first point
+    if ( f[0] != e.start() && f[0] != e.end() )
+    {
+        return f[0];
+    }
+
+    // Check the second point
+    if ( f[1] != e.start() && f[1] != e.end() )
+    {
+        return f[1];
+    }
+
+    // Check the third point
+    if ( f[2] != e.start() && f[2] != e.end() )
+    {
+        return f[2];
+    }
+
+    return -1;
+}
+
+
+//- Helper function to determine the orientation of a triangular face
+label eMesh::edgeDirection(const face& f, const edge& e) const
+{
+    if
+    (
+        (f[0] == e.start() && f[1] == e.end())
+     || (f[1] == e.start() && f[2] == e.end())
+     || (f[2] == e.start() && f[0] == e.end())
+    )
+    {
+        // Return counter-clockwise
+        return 1;
+    }
+    else
+    {
+        // Return clockwise
+        return -1;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+eMesh::eMesh(const polyMesh& pMesh, const word& subDir)
+:
+    objectRegistry(pMesh.time()),
+    mesh_(pMesh),
+    edges_
+    (
+        IOobject
+        (
+            "edges",
+            mesh_.facesInstance(),
+            subDir,
+            mesh_,
+            IOobject::READ_IF_PRESENT,
+            IOobject::NO_WRITE
+        )
+    ),
+    boundary_
+    (
+        IOobject
+        (
+            "edgeBoundary",
+            mesh_.facesInstance(),
+            subDir,
+            mesh_,
+            IOobject::READ_IF_PRESENT,
+            IOobject::NO_WRITE
+        ),
+        *this
+    ),
+    pePtr_(NULL),
+    epPtr_(NULL),
+    fePtr_(NULL),
+    efPtr_(NULL)
+{
+    if (debug)
+    {
+        Info << "eMesh::eMesh(const polyMesh&, const word&) : "
+             << "Creating eMesh from polyMesh"
+             << endl;
+    }
+
+    // Re-initialize / override meshSubDir
+    meshSubDir = subDir;
+
+    // Try to read from disk.
+    if (edges_.headerOk() && boundary_.headerOk())
+    {
+        // Set sizes
+        nEdges_ = edges_.size();
+        nInternalEdges_ = boundary_[0].start();
+    }
+    else
+    {
+        // Could not read ordered edges, so calculate it instead.
+        calcOrderedEdgeList();
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+eMesh::~eMesh()
+{
+    clearOut();
+}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+fileName eMesh::meshDir() const
+{
+    return mesh_.dbDir();
+}
+
+
+fileName eMesh::meshSubDirectory() const
+{
+    return mesh_.dbDir()/meshSubDir;
+}
+
+
+const Time& eMesh::time() const
+{
+    return mesh_.time();
+}
+
+
+label eMesh::nEdges() const
+{
+    return nEdges_;
+}
+
+
+label eMesh::nInternalEdges() const
+{
+    return nInternalEdges_;
+}
+
+
+const edgeList& eMesh::edges() const
+{
+    return edges_;
+}
+
+
+void eMesh::addEdgePatches(const List<ePatch*>& p)
+{
+    if (debug)
+    {
+        Info << "void eMesh::addEdgePatches(const List<ePatch*>& p) : "
+             << "Adding patches to eMesh" << endl;
+    }
+
+    if (boundary().size() > 0)
+    {
+        FatalErrorIn
+        (
+            "void eMesh::addEdgePatches(const List<ePatch*>& p)"
+        )
+            << "Boundary already exists"
+            << abort(FatalError);
+    }
+
+    boundary_.setSize(p.size());
+
+    forAll(p, patchI)
+    {
+        boundary_.set(patchI, p[patchI]);
+    }
+}
+
+
+const objectRegistry& eMesh::db() const
+{
+    return mesh_.db();
+}
+
+
+const eBoundaryMesh& eMesh::boundary() const
+{
+    return boundary_;
+}
+
+
+void eMesh::resetPrimitives
+(
+    edgeList& edges,
+    labelListList& faceEdges,
+    labelListList& edgeFaces,
+    const labelList& patchSizes,
+    const labelList& patchStarts,
+    const bool reUse,
+    const bool storePrimitives
+)
+{
+    // Clear out geometry and addressing
+    clearOut();
+
+    // Initialize pointers for storage
+    if (storePrimitives)
+    {
+        fePtr_ =
+        (
+            new labelListIOList
+            (
+                IOobject
+                (
+                    "faceEdges",
+                    mesh_.facesInstance(),
+                    meshSubDir,
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                )
+            )
+        );
+
+        efPtr_ =
+        (
+            new labelListIOList
+            (
+                IOobject
+                (
+                    "edgeFaces",
+                    mesh_.facesInstance(),
+                    meshSubDir,
+                    mesh_,
+                    IOobject::NO_READ,
+                    IOobject::NO_WRITE
+                )
+            )
+        );
+
+        if (reUse)
+        {
+            edges_.transfer(edges);
+            fePtr_->transfer(faceEdges);
+            efPtr_->transfer(edgeFaces);
+        }
+        else
+        {
+            edges_ = edges;
+            fePtr_->operator=(faceEdges);
+            efPtr_->operator=(edgeFaces);
+        }
+    }
+
+    // Reset patch sizes and starts
+    forAll(boundary_, patchI)
+    {
+        boundary_[patchI] = ePatch
+        (
+            boundary_[patchI].name(),
+            patchSizes[patchI],
+            patchStarts[patchI],
+            patchI,
+            boundary_
+        );
+    }
+
+    // Reset size information
+    nEdges_ = edges.size();
+    nInternalEdges_ = boundary_[0].start();
+
+    // Set mesh files as changed
+    setInstance(time().timeName());
+}
+
+
+//- Clear demand-driven data
+void eMesh::clearOut() const
+{
+    clearGeom();
+    clearAddressing();
+}
+
+
+//- Set the instance for mesh files
+void eMesh::setInstance(const fileName& inst)
+{
+    if (debug)
+    {
+        Info<< "void eMesh::setInstance(const fileName& inst) : "
+            << "Resetting file instance to " << inst << endl;
+    }
+
+    if (edges_.size())
+    {
+        edges_.writeOpt() = IOobject::AUTO_WRITE;
+        edges_.instance() = inst;
+    }
+
+    if (efPtr_)
+    {
+        efPtr_->writeOpt() = IOobject::AUTO_WRITE;
+        efPtr_->instance() = inst;
+    }
+
+    if (fePtr_)
+    {
+        fePtr_->writeOpt() = IOobject::AUTO_WRITE;
+        fePtr_->instance() = inst;
+    }
+
+    // Write out boundary information only
+    // if all others are being written out
+    if (edges_.size() && efPtr_ && fePtr_)
+    {
+        boundary_.writeOpt() = IOobject::AUTO_WRITE;
+        boundary_.instance() = inst;
+    }
+}
+
+
+const labelListList& eMesh::pointEdges() const
+{
+    if (!pePtr_)
+    {
+        calcPointEdges();
+    }
+
+    return *pePtr_;
+}
+
+
+const labelListList& eMesh::edgePoints() const
+{
+    if (!epPtr_)
+    {
+        calcEdgePoints();
+    }
+
+    return *epPtr_;
+}
+
+
+const labelListList& eMesh::faceEdges() const
+{
+    if (!fePtr_)
+    {
+        calcFaceEdges();
+    }
+
+    return *fePtr_;
+}
+
+
+const labelListList& eMesh::edgeFaces() const
+{
+    if (!efPtr_)
+    {
+        calcEdgeFaces();
+    }
+
+    return *efPtr_;
+}
+
+
+bool eMesh::write() const
+{
+    if (edges_.size())
+    {
+        edges_.write();
+    }
+
+    if (efPtr_)
+    {
+        efPtr_->write();
+    }
+
+    if (fePtr_)
+    {
+        fePtr_->write();
+    }
+
+    if (edges_.size() && efPtr_ && fePtr_)
+    {
+        boundary_.write();
+    }
+
+    return false;
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+bool eMesh::operator!=(const eMesh& m) const
+{
+    return &m != this;
+}
+
+bool eMesh::operator==(const eMesh& m) const
+{
+    return &m == this;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/eMesh.H

@@ -0,0 +1,251 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    eMesh
+
+Description
+    Mesh class to store edge-based connectivity structures.
+
+SourceFiles
+    eMesh.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef eMesh_H
+#define eMesh_H
+
+#include "Time.H"
+#include "polyMesh.H"
+#include "edgeIOList.H"
+#include "labelIOList.H"
+#include "eBoundaryMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Class forward declarations
+
+/*---------------------------------------------------------------------------*\
+                           Class eMesh Declaration
+\*---------------------------------------------------------------------------*/
+
+class eMesh
+:
+    public objectRegistry
+{
+    // Private data
+
+        const polyMesh& mesh_;
+
+        //- Ordered edge-list
+        mutable edgeIOList edges_;
+
+        //- Boundary mesh
+        mutable eBoundaryMesh boundary_;
+
+        // Primitive size data
+
+            //- Number of edges
+            mutable label nEdges_;
+
+            //- Number of internal edges
+            mutable label nInternalEdges_;
+
+    // Demand-driven data
+
+            //- Mapping between ordered edge-data and regular edges.
+            labelList reverseEdgeMap_;
+
+            //- Point-edges
+            mutable labelListList* pePtr_;
+
+            //- Edge-points
+            mutable labelListList* epPtr_;
+
+            //- Face-edges
+            mutable labelListIOList* fePtr_;
+
+            //- Edge-faces
+            mutable labelListIOList* efPtr_;
+
+    // Private Member Functions
+
+        //- Helper function to isolate points on triangular faces
+        label findIsolatedPoint(const face& f, const edge& e) const;
+
+        //- Helper function to determine the orientation of a triangular face
+        label edgeDirection(const face& f, const edge& e) const;
+
+        //- Disallow default bitwise copy construct
+        eMesh(const eMesh&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const eMesh&);
+
+        // Private member functions to calculate demand driven data
+
+            //- Calculate ordered edges (and edgeFaces)
+            void calcOrderedEdgeList();
+
+            //- Calculate ordered edgePoints
+            void calcEdgePoints() const;
+
+            //- Calculate point-edges
+            void calcPointEdges() const;
+
+            //- Calculate face-edges
+            void calcFaceEdges() const;
+
+            //- Calculate edge-faces
+            void calcEdgeFaces() const;
+
+            //- Clear geometry
+            void clearGeom() const;
+
+            //- Clear addressing
+            void clearAddressing() const;
+
+public:
+
+    // Public typedefs
+
+        typedef eMesh Mesh;
+        typedef eBoundaryMesh BoundaryMesh;
+
+
+    //- Runtime type information
+    TypeName("eMesh");
+
+    //- Return the mesh sub-directory name (usually "eMesh")
+    static word meshSubDir;
+
+
+    // Constructors
+
+        //- Construct from IOobject and polyMesh reference
+        eMesh(const polyMesh& m, const word& subDir = eMesh::meshSubDir);
+
+
+    // Destructor
+
+        virtual ~eMesh();
+
+
+    // Member Functions
+
+        // Helpers
+
+            //- Add boundary patches. Constructor helper
+            void addEdgePatches(const List<ePatch*> &);
+
+
+        // Database
+
+            //- Return reference to the mesh database
+            virtual const objectRegistry& db() const;
+
+            //- Return the base mesh directory (dbDir())
+            fileName meshDir() const;
+
+            //- Return the local mesh directory (dbDir()/meshSubDir)
+            fileName meshSubDirectory() const;
+
+            //- Return reference to time
+            const Time& time() const;
+
+
+            //- Mesh size parameters
+
+                label nEdges() const;
+
+                label nInternalEdges() const;
+
+            // Primitive mesh data
+
+                //- Return constant reference to the ordered edge-list
+                const edgeList& edges() const;
+
+        // Access
+
+            //- Return constant reference to boundary mesh
+            const eBoundaryMesh& boundary() const;
+
+
+        // Demand-driven data
+
+            //- Return constant reference to the pointEdges list
+            const labelListList& pointEdges() const;
+
+            //- Return constant reference to the edgePoints list
+            //- Valid for 3D tetrahedral meshes only.
+            const labelListList& edgePoints() const;
+
+            //- Return constant reference to the faceEdges list
+            const labelListList& faceEdges() const;
+
+            //- Return constant reference to the edgeFaces list
+            const labelListList& edgeFaces() const;
+
+        // Reset primitive data
+        void resetPrimitives
+        (
+            edgeList& edges,
+            labelListList& faceEdges,
+            labelListList& edgeFaces,
+            const labelList& patchSizes,
+            const labelList& patchStarts,
+            const bool reUse,
+            const bool storePrimitives
+        );
+
+        //- Clear demand-driven data
+        void clearOut() const;
+
+        //- Set the instance for mesh files
+        void setInstance(const fileName& inst);
+
+        //- Write mesh
+        virtual bool write() const;
+
+
+    // Member Operators
+
+        bool operator!=(const eMesh& m) const;
+
+        bool operator==(const eMesh& m) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/eMeshDemandDrivenData.C

@@ -0,0 +1,467 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+
+\*---------------------------------------------------------------------------*/
+
+#include "eMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+//- Calculate ordered edges
+void eMesh::calcOrderedEdgeList()
+{
+    if (debug)
+    {
+        Info<< "void eMesh::calcOrderedEdges() const : "
+            << "Calculating ordered edges" << endl;
+    }
+
+    if (edges_.size() || boundary_.size())
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcOrderedEdges() const"
+        )   << "Ordered edges already allocated."
+            << abort(FatalError);
+    }
+
+    // Set size first.
+    nEdges_ = mesh_.nEdges();
+    nInternalEdges_ = 0;
+    boundary_.setSize(mesh_.boundaryMesh().size());
+
+    // Allocate lists for re-ordering
+    labelList edgePatch(nEdges_, -1);
+    labelList edgePatchStarts(mesh_.boundaryMesh().size(), -1);
+    labelList edgePatchSizes(mesh_.boundaryMesh().size(), 0);
+
+    reverseEdgeMap_.setSize(nEdges_);
+
+    // Obtain connectivity from primitive mesh
+    const edgeList& edges = mesh_.edges();
+    const labelListList& fEdges = mesh_.faceEdges();
+
+    // Edge-patches are the same as faces
+    for (label i = mesh_.nInternalFaces(); i < mesh_.nFaces(); i++)
+    {
+        const labelList& fEdge = fEdges[i];
+
+        forAll(fEdge, edgeI)
+        {
+            edgePatch[fEdge[edgeI]] = mesh_.boundaryMesh().whichPatch(i);
+        }
+    }
+
+    // Loop through edgePatch and renumber internal edges
+    forAll(edgePatch, edgeI)
+    {
+        if (edgePatch[edgeI] == -1)
+        {
+            reverseEdgeMap_[edgeI] = nInternalEdges_++;
+        }
+        else
+        {
+            edgePatchSizes[edgePatch[edgeI]]++;
+        }
+    }
+
+    // Calculate patch-starts
+    label startCount = nInternalEdges_;
+
+    forAll(edgePatchStarts, patchI)
+    {
+        edgePatchStarts[patchI] = startCount;
+        startCount += edgePatchSizes[patchI];
+    }
+
+    // Now renumber boundary edges
+    labelList patchCount(edgePatchStarts);
+
+    forAll(edgePatch, edgeI)
+    {
+        if (edgePatch[edgeI] > -1)
+        {
+            reverseEdgeMap_[edgeI] = patchCount[edgePatch[edgeI]]++;
+        }
+    }
+
+    // Renumber and fill in edges
+    edges_.setSize(nEdges_, edge(-1,-1));
+
+    forAll(edges, edgeI)
+    {
+        edges_[reverseEdgeMap_[edgeI]] = edges[edgeI];
+    }
+
+    // Now set the boundary, copy name. (type is default)
+    forAll(boundary_, patchI)
+    {
+        boundary_.set
+        (
+            patchI,
+            ePatch::New
+            (
+                ePatch::typeName_(),
+                mesh_.boundaryMesh()[patchI].name(),
+                edgePatchSizes[patchI],
+                edgePatchStarts[patchI],
+                patchI,
+                boundary_
+            )
+        );
+    }
+
+    // Force calculation of other data...
+    calcEdgeFaces();
+    calcFaceEdges();
+}
+
+void eMesh::calcEdgePoints() const
+{
+    if (debug)
+    {
+        Info<< "void eMesh::calcEdgePoints() const : "
+            << "Calculating ordered edgePoints" << endl;
+    }
+
+    if (epPtr_)
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcEdgePoints() const"
+        )   << "edgePoints already allocated."
+            << abort(FatalError);
+    }
+
+    if (!efPtr_->size())
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcEdgePoints() const"
+        )   << "edgeFaces doesn't exist."
+            << abort(FatalError);
+    }
+
+    // Size the list
+    epPtr_ = new labelListList(nEdges_);
+    labelListList& edgePoints = *epPtr_;
+
+    // EdgePoints are ordered such that points appear in
+    // counter-clockwise direction around point[0] of the edge.
+    // If a boundary point is present, that must begin the list.
+    // NOTE: Will work only on tetrahedral meshes!
+
+    bool found;
+    label faceIndex = -1, cellIndex = -1;
+    const labelList& owner = mesh_.faceOwner();
+    const labelList& neighbour = mesh_.faceNeighbour();
+    const cellList& cells = mesh_.cells();
+    const faceList& faces = mesh_.faces();
+    const labelListList& eFaces = this->edgeFaces();
+
+    for (label eIndex=0; eIndex < nEdges_; eIndex++)
+    {
+        const edge& e = edges_[eIndex];
+        const labelList& eFace = eFaces[eIndex];
+
+        // Size the list
+        edgePoints[eIndex].setSize(eFace.size(), -1);
+
+        if (eIndex < nInternalEdges_)
+        {
+            // Pick the first face and start with that
+            faceIndex = eFace[0];
+        }
+        else
+        {
+            // Need to find a properly oriented start-face
+            forAll(eFace, faceI)
+            {
+                if
+                (
+                    (!mesh_.isInternalFace(eFace[faceI]))
+                 && (edgeDirection(faces[eFace[faceI]], e) == 1)
+                )
+                {
+                    faceIndex = eFace[faceI];
+                    break;
+                }
+            }
+        }
+
+        // Shuffle vertices to appear in CCW order
+        for (label j=0; j < eFace.size(); j++)
+        {
+            const face& f = faces[faceIndex];
+
+            // Add the isolated point
+            edgePoints[eIndex][j] = findIsolatedPoint(f, e);
+
+            // Figure out how this edge is oriented.
+            if (edgeDirection(f, e) == 1)
+            {
+                // Counter-clockwise. Pick the owner.
+                cellIndex = owner[faceIndex];
+            }
+            else
+            if (mesh_.isInternalFace(faceIndex))
+            {
+                // Clockwise. Pick the neighbour.
+                cellIndex = neighbour[faceIndex];
+            }
+            else
+            {
+                // Looks like we've hit a boundary face. Break out.
+                break;
+            }
+
+            const cell& cellToCheck = cells[cellIndex];
+
+            found = false;
+
+            // Assuming tet-cells,
+            // Loop through edgeFaces and get the next face
+            forAll(eFace, faceI)
+            {
+                if
+                (
+                    eFace[faceI] != faceIndex
+                 && eFace[faceI] == cellToCheck[0]
+                )
+                {
+                    faceIndex = cellToCheck[0];
+                    found = true; break;
+                }
+
+                if
+                (
+                    eFace[faceI] != faceIndex
+                 && eFace[faceI] == cellToCheck[1]
+                )
+                {
+                    faceIndex = cellToCheck[1];
+                    found = true; break;
+                }
+
+                if
+                (
+                    eFace[faceI] != faceIndex
+                 && eFace[faceI] == cellToCheck[2]
+                )
+                {
+                    faceIndex = cellToCheck[2];
+                    found = true; break;
+                }
+
+                if
+                (
+                    eFace[faceI] != faceIndex
+                 && eFace[faceI] == cellToCheck[3]
+                )
+                {
+                    faceIndex = cellToCheck[3];
+                    found = true; break;
+                }
+            }
+
+#           ifdef FULLDEBUG
+            if (!found)
+            {
+                // Something's terribly wrong
+                FatalErrorIn
+                (
+                    "void eMesh::calcEdgePoints() const"
+                )
+                << " Failed to determine a vertex ring. " << nl
+                << " edgeFaces connectivity is inconsistent. " << nl
+                << "Edge: " << eIndex << ":: " << e << nl
+                << "edgeFaces: " << eFace
+                << abort(FatalError);
+            }
+#           endif
+        }
+    }
+}
+
+
+void eMesh::calcPointEdges() const
+{
+    if (debug)
+    {
+        Info<< "void eMesh::calcPointEdges() const : "
+            << "Calculating PointEdges" << endl;
+    }
+
+    if (pePtr_)
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcPointEdges() const"
+        )   << "pePtr_ already allocated"
+            << abort(FatalError);
+    }
+
+    pePtr_ = new labelListList(mesh_.nPoints());
+    invertManyToMany(mesh_.nPoints(), edges(), *pePtr_);
+}
+
+
+void eMesh::calcFaceEdges() const
+{
+    if (debug)
+    {
+        Info<< "void eMesh::calcFaceEdges() const : "
+            << "Calculating FaceEdges" << endl;
+    }
+
+    if (fePtr_)
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcFaceEdges() const"
+        )   << "fePtr_ already allocated"
+            << abort(FatalError);
+    }
+
+    fePtr_ =
+    (
+        new labelListIOList
+        (
+            IOobject
+            (
+                "faceEdges",
+                mesh_.facesInstance(),
+                meshSubDir,
+                mesh_,
+                IOobject::READ_IF_PRESENT,
+                IOobject::NO_WRITE
+            ),
+            mesh_.nFaces()
+        )
+    );
+
+    labelListIOList& faceEdges = *fePtr_;
+
+    // If the read was successful, return.
+    if (faceEdges.headerOk())
+    {
+        return;
+    }
+
+    // If edgeFaces exists, use that.
+    if (efPtr_)
+    {
+        invertManyToMany(mesh_.nFaces(), edgeFaces(), faceEdges);
+    }
+    else
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcFaceEdges() const"
+        )   << "Cannot calculate faceEdges."
+            << abort(FatalError);
+    }
+}
+
+void eMesh::calcEdgeFaces() const
+{
+    if (debug)
+    {
+        Info<< "void eMesh::calcEdgeFaces() const : "
+            << "Calculating EdgeFaces" << endl;
+    }
+
+    if (efPtr_)
+    {
+        FatalErrorIn
+        (
+            "void eMesh::calcEdgeFaces() const"
+        )   << "efPtr_ already allocated."
+            << abort(FatalError);
+    }
+
+    efPtr_ =
+    (
+        new labelListIOList
+        (
+            IOobject
+            (
+                "edgeFaces",
+                mesh_.facesInstance(),
+                meshSubDir,
+                mesh_,
+                IOobject::READ_IF_PRESENT,
+                IOobject::NO_WRITE
+            ),
+            nEdges_
+        )
+    );
+
+    labelListIOList& edgeFaces = *efPtr_;
+
+    // If the read was successful, return.
+    if (edgeFaces.headerOk())
+    {
+        return;
+    }
+
+    if (fePtr_)
+    {
+        // If faceEdges exists, use that.
+        invertManyToMany(nEdges_, faceEdges(), edgeFaces);
+    }
+    else
+    {
+        // Obtain connectivity from primitive mesh.
+        if (!reverseEdgeMap_.size())
+        {
+            FatalErrorIn
+            (
+                "void eMesh::calcEdgeFaces() const"
+            )   << "reverseEdgeMap has not been allocated."
+                << abort(FatalError);
+        }
+
+        const labelListList& eFaces = mesh_.edgeFaces();
+
+        forAll(eFaces, edgeI)
+        {
+            edgeFaces[reverseEdgeMap_[edgeI]] = eFaces[edgeI];
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/ePatches/ePatch/ePatch.C

@@ -0,0 +1,138 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+
+\*---------------------------------------------------------------------------*/
+
+#include "ePatch.H"
+#include "addToRunTimeSelectionTable.H"
+#include "eBoundaryMesh.H"
+#include "eMesh.H"
+#include "polyMesh.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(ePatch, 0);
+
+defineRunTimeSelectionTable(ePatch, word);
+defineRunTimeSelectionTable(ePatch, dictionary);
+
+addToRunTimeSelectionTable(ePatch, ePatch, word);
+addToRunTimeSelectionTable(ePatch, ePatch, dictionary);
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+Foam::ePatch::ePatch
+(
+    const word& name,
+    const label size,
+    const label start,
+    const label index,
+    const eBoundaryMesh& bm
+)
+:
+    patchIdentifier(name, index),
+    boundaryMesh_(bm),
+    start_(start),
+    size_(size)
+{}
+
+
+// Construct from dictionary
+ePatch::ePatch
+(
+    const word& name,
+    const dictionary& dict,
+    const label index,
+    const eBoundaryMesh& bm
+)
+:
+    patchIdentifier(name, dict, index),
+    boundaryMesh_(bm),
+    start_(readLabel(dict.lookup("start"))),
+    size_(readLabel(dict.lookup("size")))
+{}
+
+ePatch::ePatch(const ePatch& p, const eBoundaryMesh& bm)
+:
+    patchIdentifier(p, p.index()),
+    boundaryMesh_(bm),
+    start_(p.start()),
+    size_(p.size())
+{}
+
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+ePatch::~ePatch()
+{}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+const eBoundaryMesh& ePatch::boundaryMesh() const
+{
+    return boundaryMesh_;
+}
+
+
+void ePatch::write(Ostream& os) const
+{
+    os.writeKeyword("type") << type() << token::END_STATEMENT << nl;
+    os.writeKeyword("start") << start() << token::END_STATEMENT << nl;
+    os.writeKeyword("size") << size() << token::END_STATEMENT << nl;
+    patchIdentifier::write(os);
+}
+
+
+//- Assignment operator
+void ePatch::operator=(const ePatch& p)
+{
+    patchIdentifier::operator=(p);
+    start_ = p.start_;
+    size_ = p.size_;
+}
+
+
+// * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
+
+Ostream& operator<<(Ostream& os, const ePatch& p)
+{
+    p.write(os);
+    os.check("Ostream& operator<<(Ostream& f, const faPatch& p)");
+    return os;
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/ePatches/ePatch/ePatch.H

@@ -0,0 +1,265 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    ePatch
+
+SourceFiles
+    ePatch.C
+    newEPatch.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef ePatch_H
+#define ePatch_H
+
+#include "patchIdentifier.H"
+#include "labelList.H"
+#include "pointField.H"
+#include "typeInfo.H"
+#include "autoPtr.H"
+#include "runTimeSelectionTables.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class eBoundaryMesh;
+
+/*---------------------------------------------------------------------------*\
+                           Class ePatch Declaration
+\*---------------------------------------------------------------------------*/
+
+class ePatch
+:
+    public patchIdentifier
+{
+private:
+
+    // Private data
+
+        //- Reference to boundary mesh
+        const eBoundaryMesh& boundaryMesh_;
+
+        //- Start index of the patch
+        label start_;
+
+        //- Size of the patch
+        label size_;
+
+        // Demand-driven private data
+
+
+    // Private Member Functions
+
+        //- Disallow construct as copy
+        ePatch(const ePatch&);
+
+protected:
+
+        // The ePatch geometry initialisation is called by eBoundaryMesh
+        friend class eBoundaryMesh;
+
+        //- Initialise the calculation of the patch geometry
+        virtual void initGeometry()
+        {}
+
+        //- Calculate the patch geometry
+        virtual void calcGeometry()
+        {}
+
+        //- Initialise the patches for moving points
+        virtual void initMovePoints(const pointField&)
+        {}
+
+        //- Correct patch after moving points
+        virtual void movePoints(const pointField&)
+        {}
+
+        //- Initialise the update of the patch topology
+        virtual void initUpdateMesh()
+        {}
+
+        //- Update of the patch topology
+        virtual void updateMesh()
+        {}
+
+
+public:
+
+    typedef eBoundaryMesh BoundaryMesh;
+
+
+    //- Runtime type information
+    TypeName("patch");
+
+    // Declare run-time constructor selection tables
+
+        declareRunTimeSelectionTable
+        (
+            autoPtr,
+            ePatch,
+            word,
+            (
+                const word& name,
+                const label size,
+                const label start,
+                const label index,
+                const eBoundaryMesh& bm
+            ),
+            (name, size, start, index, bm)
+        );
+
+
+        declareRunTimeSelectionTable
+        (
+            autoPtr,
+            ePatch,
+            dictionary,
+            (
+                const word& name,
+                const dictionary& dict,
+                const label index,
+                const eBoundaryMesh& bm
+            ),
+            (name, dict, index, bm)
+        );
+
+
+    // Constructors
+
+        //- Construct from components
+        ePatch
+        (
+            const word& name,
+            const label size,
+            const label start,
+            const label index,
+            const eBoundaryMesh& bm
+        );
+
+
+        //- Construct from dictionary
+        ePatch
+        (
+            const word& name,
+            const dictionary& dict,
+            const label index,
+            const eBoundaryMesh& bm
+        );
+
+        //- Construct as copy, resetting the boundary mesh
+        ePatch(const ePatch&, const eBoundaryMesh&);
+
+
+    // Selectors
+
+        //- Return a pointer to a new patch created on freestore from
+        //  components
+        static autoPtr<ePatch> New
+        (
+            const word& patchType,
+            const word& name,
+            const label size,
+            const label start,
+            const label index,
+            const eBoundaryMesh& bm
+        );
+
+
+        //- Return a pointer to a new patch created 
+        // on freestore from dictionary
+        static autoPtr<ePatch> New
+        (
+            const word& name,
+            const dictionary& dict,
+            const label index,
+            const eBoundaryMesh& bm
+        );
+
+
+    // Destructor
+
+        virtual ~ePatch();
+
+
+    // Member Functions
+
+        //- Return the index of this patch in the boundaryMesh
+        label index() const
+        {
+            return patchIdentifier::index();
+        }
+
+        //- Return boundaryMesh reference
+        const eBoundaryMesh& boundaryMesh() const;
+
+        //- Return true if this patch is coupled
+        virtual bool coupled() const
+        {
+            return false;
+        }
+
+        //- Patch start in edge list
+        label start() const
+        {
+            return start_;
+        }
+
+        //- Patch size
+        virtual label size() const
+        {
+            return size_;
+        }
+
+        //- Slice list to patch
+        template<class T>
+        typename List<T>::subList patchSlice(const List<T>& l) const
+        {
+            return typename List<T>::subList(l, size(), start());
+        }
+
+        //- Write
+        virtual void write(Ostream&) const;
+
+        //- Assignment operator
+        void operator=(const ePatch&);
+
+
+    // Ostream Operator
+
+        friend Ostream& operator<<(Ostream&, const ePatch&);
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/ePatches/ePatch/ePatchList.H

@@ -0,0 +1,56 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Type
+    ePatchList
+
+Description
+    Container classes for ePatch
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef ePatchList_H
+#define ePatchList_H
+
+#include "ePatch.H"
+#include "PtrList.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+typedef PtrList<ePatch> ePatchList;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/eMesh/ePatches/ePatch/newEPatch.C

@@ -0,0 +1,118 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+
+\*---------------------------------------------------------------------------*/
+
+#include "ePatch.H"
+#include "dictionary.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+autoPtr<ePatch> ePatch::New
+(
+    const word& patchType,
+    const word& name,
+    const label size,
+    const label start,
+    const label index,
+    const eBoundaryMesh& bm
+)
+{
+    if (debug)
+    {
+        Info<< "ePatch::New(const word&, const word&, const label, "
+               "const label, const label, const eBoundaryMesh&) : "
+               "constructing ePatch"
+            << endl;
+    }
+
+    wordConstructorTable::iterator cstrIter =
+        wordConstructorTablePtr_->find(patchType);
+
+    if (cstrIter == wordConstructorTablePtr_->end())
+    {
+        FatalErrorIn
+        (
+            "ePatch::New(const word&, const word&, const label, "
+            "const label, const label, const eBoundaryMesh&) "
+        )   << "Unknown ePatch type " << patchType << " for patch " << name
+            << endl << endl
+            << "Valid ePatch types are :" << endl
+            << wordConstructorTablePtr_->toc()
+            << exit(FatalError);
+    }
+
+    return autoPtr<ePatch>(cstrIter()(name, size, start, index, bm));
+}
+
+
+autoPtr<ePatch> ePatch::New
+(
+    const word& name,
+    const dictionary& dict,
+    const label index,
+    const eBoundaryMesh& bm
+)
+{
+    if (debug)
+    {
+        Info<< "ePatch::New(const word&, const dictionary&, const label, "
+               "const eBoundaryMesh&) : constructing ePatch"
+            << endl;
+    }
+
+    word patchType(dict.lookup("type"));
+
+    dictionaryConstructorTable::iterator cstrIter =
+        dictionaryConstructorTablePtr_->find(patchType);
+
+    if (cstrIter == dictionaryConstructorTablePtr_->end())
+    {
+        FatalIOErrorIn
+        (
+            "ePatch::New(const word&, const dictionary&, "
+            "const label, const eBoundaryMesh&)",
+            dict
+        )   << "Unknown ePatch type " << patchType << endl << endl
+            << "Valid ePatch types are :" << endl
+            << dictionaryConstructorTablePtr_->toc()
+            << exit(FatalIOError);
+    }
+
+    return autoPtr<ePatch>(cstrIter()(name, dict, index, bm));
+}
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/fluxCorrector.C

@@ -0,0 +1,173 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    fluxCorrector
+
+Description
+    Implementation of the fluxCorrector base class
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*----------------------------------------------------------------------------*/
+
+#include "fluxCorrector.H"
+#include "dlLibraryTable.H"
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+defineTypeNameAndDebug(fluxCorrector, 0);
+defineRunTimeSelectionTable(fluxCorrector, mesh);
+
+// * * * * * * * * * * * * * * * * Selectors * * * * * * * * * * * * * * * * //
+
+autoPtr<fluxCorrector> fluxCorrector::New
+(
+    const fvMesh& mesh,
+    const dictionary& dict
+)
+{
+    // Check if an optional entry was specified
+    if (dict.found("fluxCorrector"))
+    {
+        word correctorTypeName(dict.lookup("fluxCorrector"));
+
+        // Open any supplied libraries in dictionary
+        dlLibraryTable::open
+        (
+            dict,
+            "fluxCorrectorLibs",
+            meshConstructorTablePtr_
+        );
+
+        if (!meshConstructorTablePtr_)
+        {
+            FatalErrorIn
+            (
+                "autoPtr<fluxCorrector> fluxCorrector::New"
+                "(const fvMesh& mesh, const dictionary& dict)"
+            )   << "fluxCorrector table is empty"
+                << exit(FatalError);
+        }
+
+        correctorTypeName = word(dict.lookup("fluxCorrector"));
+
+        meshConstructorTable::iterator cstrIter =
+            meshConstructorTablePtr_->find(correctorTypeName);
+
+        if (cstrIter == meshConstructorTablePtr_->end())
+        {
+            FatalErrorIn
+            (
+                "autoPtr<fluxCorrector> fluxCorrector::New"
+                "(const fvMesh& mesh, const dictionary& dict)"
+            )   << "Unknown fluxCorrector type " << correctorTypeName
+                << endl << endl
+                << "Valid fluxCorrector types are: " << endl
+                << meshConstructorTablePtr_->toc()
+                << exit(FatalError);
+        }
+
+        return autoPtr<fluxCorrector>(cstrIter()(mesh, dict));
+    }
+
+    // Return the default fluxCorrector
+    return autoPtr<fluxCorrector>(new fluxCorrector(mesh, dict));
+}
+
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+//- Return reference to mesh
+const fvMesh& fluxCorrector::mesh() const
+{
+    return mesh_;
+}
+
+
+//- Return reference to dictionary
+const dictionary& fluxCorrector::dict() const
+{
+    return dict_;
+}
+
+
+//- Is flux-correction required?
+bool fluxCorrector::required() const
+{
+    // Support for cases which do not involve a flow-solver.
+    // Notify the user, just in case.
+    Info << " ~~~ No flux correction ~~~ " << endl;
+
+    return false;
+}
+
+
+//- Interpolate fluxes to a specified list of faces
+void fluxCorrector::interpolateFluxes(const labelList& faces) const
+{
+    if (required())
+    {
+        // Throw an error stating that the scheme hasn't been implemented
+        notImplemented
+        (
+            "void fluxCorrector::interpolateFluxes"
+            "(const labelList& faces) const"
+        );
+    }
+}
+
+
+//- Update fluxes in the registry, if required
+void fluxCorrector::updateFluxes() const
+{
+    if (required())
+    {
+        // Throw an error stating that the scheme hasn't been implemented
+        notImplemented("void fluxCorrector::updateFluxes() const");
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+void fluxCorrector::operator=(const fluxCorrector& rhs)
+{
+    // Check for assignment to self
+    if (this == &rhs)
+    {
+        FatalErrorIn("fluxCorrector::operator=(const fluxCorrector&)")
+            << "Attempted assignment to self"
+            << abort(FatalError);
+    }
+}
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/fluxCorrector.H

@@ -0,0 +1,153 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    fluxCorrector
+
+Description
+    Virtual base class that deals with flux-correction after topo-changes.
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    fluxCorrector.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef fluxCorrector_H
+#define fluxCorrector_H
+
+#include "fvMesh.H"
+#include "runTimeSelectionTables.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Forward class declarations
+
+/*---------------------------------------------------------------------------*\
+                         Class fluxCorrector Declaration
+\*---------------------------------------------------------------------------*/
+
+class fluxCorrector
+{
+    // Private data
+
+        //- Reference to fvMesh
+        const fvMesh& mesh_;
+
+        //- Reference to dictionary
+        const dictionary& dict_;
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        fluxCorrector(const fluxCorrector&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const fluxCorrector&);
+
+public:
+
+    //- Runtime type information
+    TypeName("fluxCorrector");
+
+
+    // Declare run-time constructor selection tables
+
+        declareRunTimeSelectionTable
+        (
+            autoPtr,
+            fluxCorrector,
+            mesh,
+            (
+                const fvMesh& mesh,
+                const dictionary& dict
+            ),
+            (mesh, dict)
+        );
+
+
+    // Selectors
+
+        //- Select constructed from fvMesh
+        static autoPtr<fluxCorrector> New
+        (
+            const fvMesh& mesh,
+            const dictionary& dict
+        );
+
+
+    // Constructors
+
+        //- Construct from fvMesh and dictionary
+        fluxCorrector(const fvMesh& mesh, const dictionary& dict)
+        :
+            mesh_(mesh),
+            dict_(dict)
+        {}
+
+
+    // Destructor
+
+        virtual ~fluxCorrector()
+        {}
+
+
+    // Member Functions
+
+        //- Return reference to mesh
+        const fvMesh& mesh() const;
+
+        //- Return reference to dictionary
+        const dictionary& dict() const;
+
+        //- Is flux-correction required?
+        virtual bool required() const;
+
+        //- Interpolate fluxes to a specified list of faces
+        virtual void interpolateFluxes(const labelList& faces) const;
+
+        //- Update fluxes in the registry
+        virtual void updateFluxes() const;
+};
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+#ifdef NoRepository
+#include "fluxCorrector.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoBoundaryMeshMapper.H

@@ -0,0 +1,109 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoBoundaryMeshMapper
+
+Description
+    This object provides mapping and fill-in information for boundary data
+    between the two meshes after the topological change.  It is
+    constructed from mapPolyMesh.
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef topoBoundaryMeshMapper_H
+#define topoBoundaryMeshMapper_H
+
+#include "PtrList.H"
+#include "topoPatchMapper.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+/*---------------------------------------------------------------------------*\
+                    Class topoBoundaryMeshMapper Declaration
+\*---------------------------------------------------------------------------*/
+
+class topoBoundaryMeshMapper
+:
+    public PtrList<topoPatchMapper>
+{
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        topoBoundaryMeshMapper(const topoBoundaryMeshMapper&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const topoBoundaryMeshMapper&);
+
+
+public:
+
+    // Constructors
+
+        //- Construct from components
+        topoBoundaryMeshMapper
+        (
+            const fvMesh& mesh,
+            const mapPolyMesh& mpm,
+            const topoMapper& mapper
+        )
+        :
+            PtrList<topoPatchMapper>(mesh.boundary().size())
+        {
+            const fvBoundaryMesh& patches = mesh.boundary();
+
+            forAll (patches, patchI)
+            {
+                set
+                (
+                    patchI,
+                    new topoPatchMapper
+                    (
+                        patches[patchI],
+                        mpm,
+                        mapper
+                    )
+                );
+            }
+        }
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoCellMapper.C

@@ -0,0 +1,570 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoCellMapper
+
+Description
+    Implementation of the topoCellMapper class
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*----------------------------------------------------------------------------*/
+
+#include "IOmanip.H"
+#include "topoMapper.H"
+#include "mapPolyMesh.H"
+#include "topoCellMapper.H"
+#include "demandDrivenData.H"
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+//- Clear out local storage
+void topoCellMapper::clearOut()
+{
+    deleteDemandDrivenData(directAddrPtr_);
+    deleteDemandDrivenData(interpolationAddrPtr_);
+    deleteDemandDrivenData(weightsPtr_);
+    deleteDemandDrivenData(insertedCellLabelsPtr_);
+    deleteDemandDrivenData(volumesPtr_);
+    deleteDemandDrivenData(centresPtr_);
+}
+
+
+//- Calculate addressing for interpolative mapping
+void topoCellMapper::calcAddressing() const
+{
+    if (directAddrPtr_ || interpolationAddrPtr_)
+    {
+        FatalErrorIn("void topoCellMapper::calcAddressing() const")
+            << "Addressing already calculated."
+            << abort(FatalError);
+    }
+
+    // Allocate for inserted cell labels
+    label nInsertedCells = 0;
+
+    insertedCellLabelsPtr_ = new labelList(mesh_.nCells(), -1);
+    labelList& insertedCells = *insertedCellLabelsPtr_;
+
+    if (direct())
+    {
+        // Direct addressing, no weights
+        directAddrPtr_ = new labelList(mpm_.cellMap());
+    }
+    else
+    {
+        // Interpolative addressing
+        interpolationAddrPtr_ = new labelListList(mesh_.nCells());
+        labelListList& addr = *interpolationAddrPtr_;
+
+        const List<objectMap>& cfc = mpm_.cellsFromCellsMap();
+
+        forAll (cfc, cfcI)
+        {
+            // Get addressing
+            const labelList& mo = cfc[cfcI].masterObjects();
+
+            label cellI = cfc[cfcI].index();
+
+            if (addr[cellI].size() > 0)
+            {
+                FatalErrorIn("void topoCellMapper::calcAddressing() const")
+                    << "Master cell " << cellI
+                    << " mapped from cells " << mo
+                    << " is already destination for mapping."
+                    << abort(FatalError);
+            }
+
+            // Set master objects
+            addr[cellI] = mo;
+        }
+
+        // Do mapped cells. Note that this can already be set by cellsFromCells
+        // so check if addressing size still zero.
+        const labelList& cm = mpm_.cellMap();
+
+        forAll (cm, cellI)
+        {
+            // Mapped from a single cell
+            if (cm[cellI] > -1 && addr[cellI].empty())
+            {
+                addr[cellI] = labelList(1, cm[cellI]);
+            }
+
+            // Check for inserted cells without any addressing
+            if (cm[cellI] < 0 && addr[cellI].empty())
+            {
+                insertedCells[nInsertedCells++] = cellI;
+            }
+        }
+    }
+
+    // Shorten inserted cells to actual size
+    insertedCells.setSize(nInsertedCells);
+
+    if (nInsertedCells)
+    {
+        FatalErrorIn("void topoCellMapper::calcAddressing() const")
+            << " Found " << nInsertedCells << " which are"
+            << " not mapped from any parent cells." << nl
+            << " List: " << nl
+            << insertedCells
+            << abort(FatalError);
+    }
+}
+
+
+//- Calculate inverse-distance weights for interpolative mapping
+void topoCellMapper::calcInverseDistanceWeights() const
+{
+    if (weightsPtr_)
+    {
+        FatalErrorIn
+        (
+            "void topoCellMapper::calcInverseDistanceWeights() const"
+        )
+            << "Weights already calculated."
+            << abort(FatalError);
+    }
+
+    // Fetch interpolative addressing
+    const labelListList& addr = addressing();
+
+    // Allocate memory
+    weightsPtr_ = new scalarListList(size());
+    scalarListList& w = *weightsPtr_;
+
+    // Obtain cell-centre information from old/new meshes
+    const vectorField& oldCentres = tMapper_.internalCentres();
+    const vectorField& newCentres = mesh_.cellCentres();
+
+    forAll(addr, cellI)
+    {
+        const labelList& mo = addr[cellI];
+
+        // Do mapped cells
+        if (mo.size() == 1)
+        {
+            w[cellI] = scalarList(1, 1.0);
+        }
+        else
+        {
+            // Map from masters, inverse-distance weights
+            scalar totalWeight = 0.0;
+            w[cellI] = scalarList(mo.size(), 0.0);
+
+            forAll (mo, oldCellI)
+            {
+                w[cellI][oldCellI] =
+                (
+                    1.0/stabilise
+                    (
+                        magSqr
+                        (
+                            newCentres[cellI]
+                          - oldCentres[mo[oldCellI]]
+                        ),
+                        VSMALL
+                    )
+                );
+
+                totalWeight += w[cellI][oldCellI];
+            }
+
+            // Normalize weights
+            scalar normFactor = (1.0/totalWeight);
+
+            forAll (mo, oldCellI)
+            {
+                w[cellI][oldCellI] *= normFactor;
+            }
+        }
+    }
+}
+
+
+//- Calculate intersection weights for conservative mapping
+void topoCellMapper::calcIntersectionWeightsAndCentres() const
+{
+    if (volumesPtr_ || centresPtr_)
+    {
+        FatalErrorIn
+        (
+            "void topoCellMapper::"
+            "calcIntersectionWeightsAndCentres() const"
+        )
+            << "Weights already calculated."
+            << abort(FatalError);
+    }
+
+    // Fetch interpolative addressing
+    const labelListList& addr = addressing();
+
+    // Allocate memory
+    volumesPtr_ = new List<scalarField>(size(), scalarField(0));
+    List<scalarField>& v = *volumesPtr_;
+
+    centresPtr_ = new List<vectorField>(size(), vectorField(0));
+    List<vectorField>& x = *centresPtr_;
+
+    // Obtain stored cell-centres
+    const vectorField& cellCentres = tMapper_.internalCentres();
+
+    // Fetch maps
+    const List<objectMap>& cfc = mpm_.cellsFromCellsMap();
+    const List<vectorField>& mapCellCentres = tMapper_.cellCentres();
+    const List<scalarField>& mapCellWeights = tMapper_.cellWeights();
+
+    // Fill in maps first
+    forAll(cfc, indexI)
+    {
+        x[cfc[indexI].index()] = mapCellCentres[indexI];
+        v[cfc[indexI].index()] = mapCellWeights[indexI];
+    }
+
+    // Now do mapped cells
+    forAll(addr, cellI)
+    {
+        const labelList& mo = addr[cellI];
+
+        // Check if this is indeed a mapped cell
+        if (mo.size() == 1 && x[cellI].empty() && v[cellI].empty())
+        {
+            x[cellI] = vectorField(1, cellCentres[mo[0]]);
+            v[cellI] = scalarField(1, 1.0);
+        }
+    }
+}
+
+
+const List<scalarField>& topoCellMapper::intersectionWeights() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const List<scalarField>& "
+            "topoCellMapper::intersectionWeights() const"
+        )   << "Requested interpolative weights for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!volumesPtr_)
+    {
+        calcIntersectionWeightsAndCentres();
+    }
+
+    return *volumesPtr_;
+}
+
+
+const List<vectorField>& topoCellMapper::intersectionCentres() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const List<vectorField>& "
+            "topoCellMapper::intersectionCentres() const"
+        )   << "Requested interpolative weights for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!centresPtr_)
+    {
+        calcIntersectionWeightsAndCentres();
+    }
+
+    return *centresPtr_;
+}
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+// Construct from components
+topoCellMapper::topoCellMapper
+(
+    const mapPolyMesh& mpm,
+    const topoMapper& mapper
+)
+:
+    mesh_(mpm.mesh()),
+    mpm_(mpm),
+    tMapper_(mapper),
+    direct_(false),
+    directAddrPtr_(NULL),
+    interpolationAddrPtr_(NULL),
+    weightsPtr_(NULL),
+    insertedCellLabelsPtr_(NULL),
+    volumesPtr_(NULL),
+    centresPtr_(NULL)
+{
+    // Check for possibility of direct mapping
+    if
+    (
+        (min(mpm_.cellMap()) > -1)
+     && mpm_.cellsFromPointsMap().empty()
+     && mpm_.cellsFromEdgesMap().empty()
+     && mpm_.cellsFromFacesMap().empty()
+     && mpm_.cellsFromCellsMap().empty()
+    )
+    {
+        direct_ = true;
+    }
+    else
+    {
+        direct_ = false;
+    }
+}
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+topoCellMapper::~topoCellMapper()
+{
+    clearOut();
+}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+//- Return size
+label topoCellMapper::size() const
+{
+    return mesh_.nCells();
+}
+
+
+//- Return size before mapping
+label topoCellMapper::sizeBeforeMapping() const
+{
+    return mpm_.nOldCells();
+}
+
+
+//- Is the mapping direct
+bool topoCellMapper::direct() const
+{
+    return direct_;
+}
+
+
+//- Return direct addressing
+const unallocLabelList& topoCellMapper::directAddressing() const
+{
+    if (!direct())
+    {
+        FatalErrorIn
+        (
+            "const unallocLabelList& "
+            "topoCellMapper::directAddressing() const"
+        )   << "Requested direct addressing for an interpolative mapper."
+            << abort(FatalError);
+    }
+
+    if (!directAddrPtr_)
+    {
+        calcAddressing();
+    }
+
+    return *directAddrPtr_;
+}
+
+
+//- Return interpolation addressing
+const labelListList& topoCellMapper::addressing() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const labelListList& "
+            "topoCellMapper::addressing() const"
+        )   << "Requested interpolative addressing for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!interpolationAddrPtr_)
+    {
+        calcAddressing();
+    }
+
+    return *interpolationAddrPtr_;
+}
+
+
+//- Return weights
+const scalarListList& topoCellMapper::weights() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const scalarListList& "
+            "topoCellMapper::weights() const"
+        )   << "Requested interpolative weights for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!weightsPtr_)
+    {
+        calcInverseDistanceWeights();
+    }
+
+    return *weightsPtr_;
+}
+
+
+//- Are there any inserted cells
+bool topoCellMapper::insertedObjects() const
+{
+    return insertedObjectLabels().size();
+}
+
+
+//- Return list of inserted cells
+const labelList& topoCellMapper::insertedObjectLabels() const
+{
+    if (!insertedCellLabelsPtr_)
+    {
+        calcAddressing();
+    }
+
+    return *insertedCellLabelsPtr_;
+}
+
+
+//- Conservatively map the internal field
+template <class Type, class gradType>
+void topoCellMapper::mapInternalField
+(
+    const word& fieldName,
+    const Field<gradType>& gF,
+    Field<Type>& iF
+) const
+{
+    if (iF.size() != sizeBeforeMapping() || gF.size() != sizeBeforeMapping())
+    {
+        FatalErrorIn
+        (
+            "\n\n"
+            "void topoCellMapper::mapInternalField<Type>\n"
+            "(\n"
+            "    const word& fieldName,\n"
+            "    const Field<gradType>& gF,\n"
+            "    Field<Type>& iF\n"
+            ") const\n"
+        )  << "Incompatible size before mapping." << nl
+           << " Field: " << fieldName << nl
+           << " Field size: " << iF.size() << nl
+           << " Gradient Field size: " << gF.size() << nl
+           << " map size: " << sizeBeforeMapping() << nl
+           << abort(FatalError);
+    }
+
+    // Fetch addressing
+    const labelListList& cAddressing = addressing();
+    const List<scalarField>& wC = intersectionWeights();
+    const List<vectorField>& xC = intersectionCentres();
+
+    // Fetch geometry
+    const vectorField& centres = tMapper_.internalCentres();
+
+    // Compute the integral of the source field
+    Type intSource = sum(iF * tMapper_.cellVolumes());
+
+    // Copy the original field
+    Field<Type> fieldCpy(iF);
+
+    // Resize to current dimensions
+    iF.setSize(size());
+
+    // Map the internal field
+    forAll(iF, cellI)
+    {
+        const labelList& addr = cAddressing[cellI];
+
+        iF[cellI] = pTraits<Type>::zero;
+
+        forAll(addr, cellJ)
+        {
+            // Accumulate volume-weighted Taylor-series interpolate
+            iF[cellI] +=
+            (
+                wC[cellI][cellJ] *
+                (
+                    fieldCpy[addr[cellJ]]
+                  + (
+                        gF[addr[cellJ]] &
+                        (
+                            xC[cellI][cellJ]
+                          - centres[addr[cellJ]]
+                        )
+                    )
+                )
+            );
+        }
+    }
+
+    // Compute the integral of the target field
+    Type intTarget = sum(iF * mesh_.cellVolumes());
+
+    if (polyMesh::debug)
+    {
+        int oldP = Info().precision();
+
+        // Compare the global integral
+        Info << " Field : " << fieldName
+             << " integral errors : "
+             << setprecision(15)
+             << " source : " << mag(intSource)
+             << " target : " << mag(intTarget)
+             << " norm : "
+             << (mag(intTarget - intSource) / (mag(intSource) + VSMALL))
+             << setprecision(oldP)
+             << endl;
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+void topoCellMapper::operator=(const topoCellMapper& rhs)
+{
+    // Check for assignment to self
+    if (this == &rhs)
+    {
+        FatalErrorIn("void topoCellMapper::operator=")
+            << "Attempted assignment to self"
+            << abort(FatalError);
+    }
+}
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoCellMapper.H

@@ -0,0 +1,188 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoCellMapper
+
+Description
+    This object provides mapping and fill-in information for internal
+    cell data between the two meshes after topological changes.
+    It is constructed from mapPolyMesh.
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    topoCellMapper.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef topoCellMapper_H
+#define topoCellMapper_H
+
+#include "morphFieldMapper.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Class forward declarations
+class topoMapper;
+class mapPolyMesh;
+
+/*---------------------------------------------------------------------------*\
+                        Class topoCellMapper Declaration
+\*---------------------------------------------------------------------------*/
+
+class topoCellMapper
+:
+    public morphFieldMapper
+{
+    // Private data
+
+        //- Reference to polyMesh
+        const polyMesh& mesh_;
+
+        //- Reference to mapPolyMesh
+        const mapPolyMesh& mpm_;
+
+        //- Reference to the topoMapper
+        const topoMapper& tMapper_;
+
+        //- Is the mapping direct
+        bool direct_;
+
+    // Demand-driven private data
+
+        //- Direct addressing
+        mutable labelList* directAddrPtr_;
+
+        //- Interpolated addressing
+        mutable labelListList* interpolationAddrPtr_;
+
+        //- Inverse-distance weights
+        mutable scalarListList* weightsPtr_;
+
+        //- Inserted cells
+        mutable labelList* insertedCellLabelsPtr_;
+
+        //- Interpolation volumes
+        mutable List<scalarField>* volumesPtr_;
+
+        //- Interpolation centres
+        mutable List<vectorField>* centresPtr_;
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        topoCellMapper(const topoCellMapper&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const topoCellMapper&);
+
+        //- Calculate addressing for interpolative mapping
+        void calcAddressing() const;
+
+        //- Calculate inverse-distance weights for interpolative mapping
+        void calcInverseDistanceWeights() const;
+
+        //- Calculate intersection weights for conservative mapping
+        void calcIntersectionWeightsAndCentres() const;
+
+        //- Return intersection volume weights
+        const List<scalarField>& intersectionWeights() const;
+
+        //- Return intersection volume centres
+        const List<vectorField>& intersectionCentres() const;
+
+        //- Clear out local storage
+        void clearOut();
+
+public:
+
+    // Constructors
+
+        //- Construct from mapPolyMesh
+        topoCellMapper
+        (
+            const mapPolyMesh& mpm,
+            const topoMapper& mapper
+        );
+
+    // Destructor
+
+        virtual ~topoCellMapper();
+
+    // Member Functions
+
+        //- Return size
+        virtual label size() const;
+
+        //- Return size before mapping
+        virtual label sizeBeforeMapping() const;
+
+        //- Is the mapping direct
+        virtual bool direct() const;
+
+        //- Return direct addressing
+        virtual const unallocLabelList& directAddressing() const;
+
+        //- Return interpolation addressing
+        virtual const labelListList& addressing() const;
+
+        //- Return weights
+        virtual const scalarListList& weights() const;
+
+        //- Are there any inserted cells
+        virtual bool insertedObjects() const;
+
+        //- Return list of inserted cells
+        virtual const labelList& insertedObjectLabels() const;
+
+        //- Conservatively map the internal field
+        template <class Type, class gradType>
+        void mapInternalField
+        (
+            const word& fieldName,
+            const Field<gradType>& gF,
+            Field<Type>& iF
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+#ifdef NoRepository
+#    include "topoCellMapper.C"
+#endif
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoMapper.C

@@ -0,0 +1,637 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoMapper
+
+Description
+    Implementation of topoMapper
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*----------------------------------------------------------------------------*/
+
+#include "fvc.H"
+#include "topoMapper.H"
+#include "fluxCorrector.H"
+#include "topoCellMapper.H"
+#include "leastSquaresGrad.H"
+#include "topoSurfaceMapper.H"
+#include "topoBoundaryMeshMapper.H"
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+// Store gradients of fields on the mesh prior to topology changes
+template <class Type, class gradType>
+void topoMapper::storeGradients
+(
+    HashTable<autoPtr<gradType> >& gradTable
+) const
+{
+    // Define a few typedefs for convenience
+    typedef typename outerProduct<vector, Type>::type gCmptType;
+    typedef GeometricField<Type, fvPatchField, volMesh> volType;
+    typedef GeometricField<gCmptType, fvPatchField, volMesh> gVolType;
+
+    // Fetch all fields from registry
+    HashTable<const volType*> fields
+    (
+        mesh_.objectRegistry::lookupClass<volType>()
+    );
+
+    forAllConstIter(typename HashTable<const volType*>, fields, fIter)
+    {
+        const volType& field = *fIter();
+
+        // Compute the gradient.
+        tmp<gVolType> tGrad;
+
+        // If the fvSolution dictionary contains an entry,
+        // use that, otherwise, default to leastSquares
+        word gradName("grad(" + field.name() + ')');
+
+        if (mesh_.schemesDict().subDict("gradSchemes").found(gradName))
+        {
+            tGrad = fvc::grad(field);
+        }
+        else
+        {
+            tGrad = fv::leastSquaresGrad<Type>(mesh_).grad(field);
+        }
+
+        // Make a new entry, but don't register the field.
+        gradTable.insert
+        (
+            field.name(),
+            autoPtr<gradType>
+            (
+                new gradType
+                (
+                    IOobject
+                    (
+                        tGrad().name(),
+                        mesh_.time().timeName(),
+                        mesh_,
+                        IOobject::NO_READ,
+                        IOobject::NO_WRITE,
+                        false
+                    ),
+                    tGrad()
+                )
+            )
+        );
+    }
+}
+
+
+//- Fetch the gradient field (template specialisation)
+template <>
+const volVectorField& topoMapper::gradient(const word& name) const
+{
+    if (!sGrads_.found(name))
+    {
+        FatalErrorIn
+        (
+            "const volVectorField& "
+            "topoMapper::gradient(const word& name) const"
+        ) << nl << " Gradient for: " << name
+          << " has not been stored."
+          << abort(FatalError);
+    }
+
+    return sGrads_[name]();
+}
+
+
+//- Fetch the gradient field (template specialisation)
+template <>
+const volTensorField& topoMapper::gradient(const word& name) const
+{
+    if (!vGrads_.found(name))
+    {
+        FatalErrorIn
+        (
+            "const volTensorField& "
+            "topoMapper::gradient(const word& name) const"
+        ) << nl << " Gradient for: " << name
+          << " has not been stored."
+          << abort(FatalError);
+    }
+
+    return vGrads_[name]();
+}
+
+
+//- Store gradients prior to mesh reset
+void topoMapper::storeGradients() const
+{
+    storeGradients<scalar>(sGrads_);
+    storeGradients<vector>(vGrads_);
+
+    if (fvMesh::debug)
+    {
+        Info << "Registered volScalarFields: " << endl;
+        Info << sGrads_.toc() << endl;
+
+        Info << "Registered volVectorFields: " << endl;
+        Info << vGrads_.toc() << endl;
+    }
+}
+
+
+//- Store geometric information
+void topoMapper::storeGeometry() const
+{
+    if (cellCentresPtr_)
+    {
+        deleteDemandDrivenData(cellVolumesPtr_);
+        deleteDemandDrivenData(cellCentresPtr_);
+
+        patchAreasPtr_.clear();
+        patchCentresPtr_.clear();
+    }
+
+    // Set the cell-volumes pointer.
+    cellVolumesPtr_ = new scalarField(mesh_.cellVolumes());
+
+    // Set the cell-centres pointer.
+    cellCentresPtr_ = new vectorField(mesh_.cellCentres());
+
+    // Set patch-areas
+    patchAreasPtr_.setSize(mesh_.boundaryMesh().size());
+
+    forAll(mesh_.boundaryMesh(), patchI)
+    {
+        patchAreasPtr_.set
+        (
+            patchI,
+            new scalarField
+            (
+                mag(mesh_.boundaryMesh()[patchI].faceAreas())
+            )
+        );
+    }
+
+    // Set patch-centres.
+    patchCentresPtr_.setSize(mesh_.boundaryMesh().size());
+
+    forAll(mesh_.boundaryMesh(), patchI)
+    {
+        patchCentresPtr_.set
+        (
+            patchI,
+            new vectorField
+            (
+                mesh_.boundaryMesh()[patchI].faceCentres()
+            )
+        );
+    }
+}
+
+
+// * * * * * * * * * * * * * * * * Destructor * * * * * * * * * * * * * * *  //
+
+topoMapper::~topoMapper()
+{
+    clear();
+}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+//- Return reference to the mesh
+const fvMesh& topoMapper::mesh() const
+{
+    return mesh_;
+}
+
+
+//- Return reference to objectRegistry storing fields.
+const objectRegistry& topoMapper::thisDb() const
+{
+    return mesh_;
+}
+
+
+//- Set mapping information
+void topoMapper::setMapper(const mapPolyMesh& mpm) const
+{
+    if
+    (
+        cellMap_.valid() ||
+        surfaceMap_.valid() ||
+        boundaryMap_.valid()
+    )
+    {
+        FatalErrorIn
+        (
+            "void topoMapper::setMapper() const"
+        ) << nl << " Mapper has already been set. "
+          << abort(FatalError);
+    }
+
+    // Set pointers
+    cellMap_.set(new topoCellMapper(mpm, *this));
+    surfaceMap_.set(new topoSurfaceMapper(mpm, *this));
+    boundaryMap_.set(new topoBoundaryMeshMapper(mesh(), mpm, *this));
+}
+
+
+//- Set face weighting information
+void topoMapper::setFaceWeights
+(
+    const Xfer<List<scalarField> >& weights,
+    const Xfer<List<vectorField> >& centres
+) const
+{
+    faceWeights_.transfer(weights());
+    faceCentres_.transfer(centres());
+}
+
+
+//- Set cell weighting information
+void topoMapper::setCellWeights
+(
+    const Xfer<List<scalarField> >& weights,
+    const Xfer<List<vectorField> >& centres
+) const
+{
+    cellWeights_.transfer(weights());
+    cellCentres_.transfer(centres());
+}
+
+
+//- Fetch face weights
+const List<scalarField>& topoMapper::faceWeights() const
+{
+    return faceWeights_;
+}
+
+
+//- Fetch cell weights
+const List<scalarField>& topoMapper::cellWeights() const
+{
+    return cellWeights_;
+}
+
+
+//- Fetch face centres
+const List<vectorField>& topoMapper::faceCentres() const
+{
+    return faceCentres_;
+}
+
+
+//- Fetch cell centres
+const List<vectorField>& topoMapper::cellCentres() const
+{
+    return cellCentres_;
+}
+
+
+//- Store mesh information for the mapping stage
+void topoMapper::storeMeshInformation() const
+{
+    // Store field-gradients
+    storeGradients();
+
+    // Store geometry
+    storeGeometry();
+}
+
+
+//- Return stored cell volume information
+const scalarField& topoMapper::cellVolumes() const
+{
+    if (!cellVolumesPtr_)
+    {
+        FatalErrorIn
+        (
+            "const scalarField& topoMapper::cellVolumes()"
+        ) << nl << " Pointer has not been set. "
+          << abort(FatalError);
+    }
+
+    return *cellVolumesPtr_;
+}
+
+
+//- Return stored cell centre information
+const vectorField& topoMapper::internalCentres() const
+{
+    if (!cellCentresPtr_)
+    {
+        FatalErrorIn
+        (
+            "const vectorField& topoMapper::internalCentres()"
+        ) << nl << " Pointer has not been set. "
+          << abort(FatalError);
+    }
+
+    return *cellCentresPtr_;
+}
+
+
+//- Return stored patch areas information
+const scalarField& topoMapper::patchAreas(const label i) const
+{
+    if (!patchAreasPtr_.set(i))
+    {
+        FatalErrorIn
+        (
+            "const scalarField& topoMapper::patchAreas"
+            "(const label i) const"
+        ) << nl << " Pointer has not been set at index: " << i
+          << abort(FatalError);
+    }
+
+    return patchAreasPtr_[i];
+}
+
+
+//- Return stored patch centre information
+const vectorField& topoMapper::patchCentres(const label i) const
+{
+    if (!patchCentresPtr_.set(i))
+    {
+        FatalErrorIn
+        (
+            "const vectorField& topoMapper::patchCentres"
+            "(const label i) const"
+        ) << nl << " Pointer has not been set at index: " << i
+          << abort(FatalError);
+    }
+
+    return patchCentresPtr_[i];
+}
+
+
+// Conservatively map all volFields in the registry
+template <class Type>
+void topoMapper::conservativeMapVolFields() const
+{
+    // Define a few typedefs for convenience
+    typedef typename outerProduct<vector, Type>::type gCmptType;
+    typedef GeometricField<Type, fvPatchField, volMesh> volType;
+    typedef GeometricField<gCmptType, fvPatchField, volMesh> gradVolType;
+
+    HashTable<const volType*> fields(mesh_.lookupClass<volType>());
+
+    // Store old-times before mapping
+    forAllIter(typename HashTable<const volType*>, fields, fIter)
+    {
+        volType& field = const_cast<volType&>(*fIter());
+
+        field.storeOldTimes();
+    }
+
+    // Fetch internal/boundary mappers
+    const topoCellMapper& fMap = volMap();
+    const topoBoundaryMeshMapper& bMap = boundaryMap();
+
+    // Now map all fields
+    forAllIter(typename HashTable<const volType*>, fields, fIter)
+    {
+        volType& field = const_cast<volType&>(*fIter());
+
+        if (fvMesh::debug)
+        {
+            Info << "Conservatively mapping "
+                 << field.typeName
+                 << ' ' << field.name()
+                 << endl;
+        }
+
+        // Map the internal field
+        fMap.mapInternalField
+        (
+            field.name(),
+            gradient<gradVolType>(field.name()).internalField(),
+            field.internalField()
+        );
+
+        // Map patch fields
+        forAll(bMap, patchI)
+        {
+            bMap[patchI].mapPatchField
+            (
+                field.name(),
+                field.boundaryField()[patchI]
+            );
+        }
+
+        // Set the field instance
+        field.instance() = field.mesh().thisDb().time().timeName();
+    }
+}
+
+
+// Conservatively map all surfaceFields in the registry
+template <class Type>
+void topoMapper::conservativeMapSurfaceFields() const
+{
+    // Define a few typedefs for convenience
+    typedef GeometricField<Type, fvsPatchField, surfaceMesh> surfType;
+
+    HashTable<const surfType*> fields(mesh_.lookupClass<surfType>());
+
+    // Store old-times before mapping
+    forAllIter(typename HashTable<const surfType*>, fields, fIter)
+    {
+        surfType& field = const_cast<surfType&>(*fIter());
+
+        field.storeOldTimes();
+    }
+
+    // Fetch internal/boundary mappers
+    const topoSurfaceMapper& fMap = surfaceMap();
+    const topoBoundaryMeshMapper& bMap = boundaryMap();
+
+    // Now map all fields
+    forAllIter(typename HashTable<const surfType*>, fields, fIter)
+    {
+        surfType& field = const_cast<surfType&>(*fIter());
+
+        if (fvMesh::debug)
+        {
+            Info << "Conservatively mapping "
+                 << field.typeName
+                 << ' ' << field.name()
+                 << endl;
+        }
+
+        // Map the internal field
+        fMap.mapInternalField
+        (
+            field.name(),
+            field.internalField()
+        );
+
+        // Map patch fields
+        forAll(bMap, patchI)
+        {
+            bMap[patchI].mapPatchField
+            (
+                field.name(),
+                field.boundaryField()[patchI]
+            );
+        }
+
+        // Set the field instance
+        field.instance() = field.mesh().thisDb().time().timeName();
+    }
+}
+
+
+//- Correct fluxes after topology changes, if required
+void topoMapper::correctFluxes() const
+{
+    if (surfaceFluxCorrector().required())
+    {
+        // Supply a list of inserted faces for interpolation
+        surfaceFluxCorrector().interpolateFluxes
+        (
+            surfaceMap().insertedObjectLabels()
+        );
+
+        // Update fluxes
+        surfaceFluxCorrector().updateFluxes();
+    }
+}
+
+
+//- Return volume mapper
+const topoCellMapper& topoMapper::volMap() const
+{
+    if (!cellMap_.valid())
+    {
+        FatalErrorIn
+        (
+            "const topoCellMapper& topoMapper::volMap()"
+        ) << nl << " Volume mapper has not been set. "
+          << abort(FatalError);
+    }
+
+    return cellMap_();
+}
+
+
+//- Return surface mapper
+const topoSurfaceMapper& topoMapper::surfaceMap() const
+{
+    if (!surfaceMap_.valid())
+    {
+        FatalErrorIn
+        (
+            "const topoSurfaceMapper& topoMapper::surfaceMap()"
+        ) << nl << " Surface mapper has not been set. "
+          << abort(FatalError);
+    }
+
+    return surfaceMap_();
+}
+
+
+//- Return boundary mapper
+const topoBoundaryMeshMapper& topoMapper::boundaryMap() const
+{
+    if (!boundaryMap_.valid())
+    {
+        FatalErrorIn
+        (
+            "const topoBoundaryMeshMapper& topoMapper::boundaryMap()"
+        ) << nl << " Boundary mapper has not been set. "
+          << abort(FatalError);
+    }
+
+    return boundaryMap_();
+}
+
+
+//- Return flux-corrector
+const fluxCorrector& topoMapper::surfaceFluxCorrector() const
+{
+    if (!fluxCorrector_.valid())
+    {
+        FatalErrorIn
+        (
+            "const fluxCorrector& topoMapper::surfaceFluxCorrector()"
+        ) << nl << " fluxCorrector has not been set. "
+          << abort(FatalError);
+    }
+
+    return fluxCorrector_();
+}
+
+
+//- Clear out member data
+void topoMapper::clear()
+{
+    // Clear out mappers
+    cellMap_.clear();
+    surfaceMap_.clear();
+    boundaryMap_.clear();
+
+    // Clear stored gradients
+    sGrads_.clear();
+    vGrads_.clear();
+
+    // Wipe out geomtry information
+    deleteDemandDrivenData(cellVolumesPtr_);
+    deleteDemandDrivenData(cellCentresPtr_);
+    patchAreasPtr_.clear();
+    patchCentresPtr_.clear();
+
+    // Clear maps
+    faceWeights_.clear();
+    cellWeights_.clear();
+
+    faceCentres_.clear();
+    cellCentres_.clear();
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+void topoMapper::operator=(const topoMapper& rhs)
+{
+    // Check for assignment to self
+    if (this == &rhs)
+    {
+        FatalErrorIn
+        (
+            "topoMapper::operator=(const topoMapper&)"
+        )
+            << "Attempted assignment to self"
+            << abort(FatalError);
+    }
+}
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoMapper.H

@@ -0,0 +1,245 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoMapper
+
+Description
+    Class holds all necessary information for mapping fields associated with
+    dynamicTopoFvMesh and fvMesh.
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    topoMapper.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef topoMapper_H
+#define topoMapper_H
+
+#include "fluxCorrector.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Class forward declarations
+class topoCellMapper;
+class topoSurfaceMapper;
+class topoBoundaryMeshMapper;
+class mapPolyMesh;
+
+/*---------------------------------------------------------------------------*\
+                         Class topoMapper Declaration
+\*---------------------------------------------------------------------------*/
+
+class topoMapper
+{
+    // Private data
+
+        //- Reference to fvMesh
+        const fvMesh& mesh_;
+
+        //- Reference to the options dictionary
+        const dictionary& dict_;
+
+    // Demand-driven private data
+
+        //- Cell mapper
+        mutable autoPtr<topoCellMapper> cellMap_;
+
+        //- Surface mapper
+        mutable autoPtr<topoSurfaceMapper> surfaceMap_;
+
+        //- Boundary mapper
+        mutable autoPtr<topoBoundaryMeshMapper> boundaryMap_;
+
+        //- Flux corrector
+        mutable autoPtr<fluxCorrector> fluxCorrector_;
+
+        // Stored gradients for mapping
+        mutable HashTable<autoPtr<volVectorField> > sGrads_;
+        mutable HashTable<autoPtr<volTensorField> > vGrads_;
+
+        //- Geometric information on the old mesh
+        mutable scalarField* cellVolumesPtr_;
+        mutable vectorField* cellCentresPtr_;
+        mutable PtrList<scalarField> patchAreasPtr_;
+        mutable PtrList<vectorField> patchCentresPtr_;
+
+        //- Intersection volume weights
+        mutable List<scalarField> faceWeights_;
+        mutable List<scalarField> cellWeights_;
+
+        //- Intersection centre weights
+        mutable List<vectorField> faceCentres_;
+        mutable List<vectorField> cellCentres_;
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        topoMapper(const topoMapper&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const topoMapper&);
+
+        // Store gradients of volFields on the mesh
+        // prior to topology changes
+        template <class Type, class gradType>
+        void storeGradients
+        (
+            HashTable<autoPtr<gradType> >& gradTable
+        ) const;
+
+        //- Fetch the gradient field
+        template <class Type>
+        const Type& gradient(const word& name) const;
+
+        //- Store gradients prior to mesh reset
+        void storeGradients() const;
+
+        //- Set geometric information
+        void storeGeometry() const;
+
+public:
+
+    // Constructors
+
+        //- Construct from mesh and dictionary
+        topoMapper(const fvMesh& mesh, const dictionary& dict)
+        :
+            mesh_(mesh),
+            dict_(dict),
+            cellMap_(NULL),
+            surfaceMap_(NULL),
+            boundaryMap_(NULL),
+            fluxCorrector_(fluxCorrector::New(mesh, dict)),
+            cellVolumesPtr_(NULL),
+            cellCentresPtr_(NULL),
+            patchAreasPtr_(mesh.boundary().size()),
+            patchCentresPtr_(mesh.boundary().size())
+        {}
+
+    // Destructor
+
+        ~topoMapper();
+
+    // Member Functions
+
+        //- Return reference to the mesh
+        const fvMesh& mesh() const;
+
+        //- Return reference to objectRegistry storing fields.
+        const objectRegistry& thisDb() const;
+
+        //- Set mapping information
+        void setMapper(const mapPolyMesh& mpm) const;
+
+        //- Set face weighting information
+        void setFaceWeights
+        (
+            const Xfer<List<scalarField> >& weights,
+            const Xfer<List<vectorField> >& centres
+        ) const;
+
+        //- Set cell weighting information
+        void setCellWeights
+        (
+            const Xfer<List<scalarField> >& weights,
+            const Xfer<List<vectorField> >& centres
+        ) const;
+
+        //- Fetch face weights
+        const List<scalarField>& faceWeights() const;
+
+        //- Fetch cell weights
+        const List<scalarField>& cellWeights() const;
+
+        //- Fetch face centres
+        const List<vectorField>& faceCentres() const;
+
+        //- Fetch cell centres
+        const List<vectorField>& cellCentres() const;
+
+        //- Store mesh information for the mapping stage
+        void storeMeshInformation() const;
+
+        //- Return stored cell volume information
+        const scalarField& cellVolumes() const;
+
+        //- Return stored cell centre information
+        const vectorField& internalCentres() const;
+
+        //- Return stored patch area information
+        const scalarField& patchAreas(const label i) const;
+
+        //- Return stored patch centre information
+        const vectorField& patchCentres(const label i) const;
+
+        // Conservatively map all volFields in the registry
+        template <class Type>
+        void conservativeMapVolFields() const;
+
+        // Conservatively map all surfaceFields in the registry
+        template <class Type>
+        void conservativeMapSurfaceFields() const;
+
+        //- Correct fluxes after topology change
+        void correctFluxes() const;
+
+        //- Return volume mapper
+        const topoCellMapper& volMap() const;
+
+        //- Return surface mapper
+        const topoSurfaceMapper& surfaceMap() const;
+
+        //- Return boundary mapper
+        const topoBoundaryMeshMapper& boundaryMap() const;
+
+        //- Return flux-corrector
+        const fluxCorrector& surfaceFluxCorrector() const;
+
+        //- Clear out member data
+        void clear();
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+#ifdef NoRepository
+#include "topoMapper.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoPatchMapper.C

@@ -0,0 +1,776 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoPatchMapper
+
+Description
+    Implementation of the topoPatchMapper class
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+\*----------------------------------------------------------------------------*/
+
+#include "IOmanip.H"
+#include "topoMapper.H"
+#include "mapPolyMesh.H"
+#include "topoPatchMapper.H"
+#include "demandDrivenData.H"
+
+namespace Foam
+{
+
+// * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
+
+//- Clear out local storage
+void topoPatchMapper::clearOut()
+{
+    deleteDemandDrivenData(directAddrPtr_);
+    deleteDemandDrivenData(interpolationAddrPtr_);
+    deleteDemandDrivenData(weightsPtr_);
+    deleteDemandDrivenData(insertedFaceLabelsPtr_);
+    deleteDemandDrivenData(insertedFaceIndexMapPtr_);
+    deleteDemandDrivenData(insertedFaceAddressingPtr_);
+    deleteDemandDrivenData(areasPtr_);
+    deleteDemandDrivenData(centresPtr_);
+}
+
+
+//- Calculate the insertedFaceLabels list
+void topoPatchMapper::calcInsertedFaceAddressing() const
+{
+    if (insertedFaceLabelsPtr_ || insertedFaceAddressingPtr_)
+    {
+        FatalErrorIn
+        (
+            "void topoPatchMapper::calcInsertedFaceAddressing() const"
+        )   << " Inserted labels has already been calculated."
+            << abort(FatalError);
+    }
+
+    // Information from the old patch
+    const label oldPatchSize = mpm_.oldPatchSizes()[patch_.index()];
+    const label oldPatchStart = mpm_.oldPatchStarts()[patch_.index()];
+    const label oldPatchEnd = oldPatchStart + oldPatchSize;
+
+    // Allocate for inserted face labels and addressing
+    label nInsertedFaces = 0;
+
+    insertedFaceLabelsPtr_ = new labelList(patchSize(), -1);
+    labelList& insertedFaces = *insertedFaceLabelsPtr_;
+
+    insertedFaceIndexMapPtr_ = new labelList(patchSize(), -1);
+    labelList& insertedFacesMap = *insertedFaceIndexMapPtr_;
+
+    insertedFaceAddressingPtr_ = new labelListList(patchSize(), labelList(0));
+    labelListList& insertedAddressing = *insertedFaceAddressingPtr_;
+
+    // Fetch current patch start
+    label pStart = patch_.patch().start();
+
+    // Fetch the current boundary
+    const polyBoundaryMesh& boundary = mpm_.mesh().boundaryMesh();
+
+    // Loop through the facesFromFaces map, and ensure that
+    // inserted faces are only mapped from faces on the same patch.
+    const List<objectMap>& fff = mpm_.facesFromFacesMap();
+
+    forAll(fff, objectI)
+    {
+        const objectMap& fffI = fff[objectI];
+
+        // Only pick boundary faces in this patch
+        if (boundary.whichPatch(fffI.index()) == patch_.index())
+        {
+            if (fffI.masterObjects().empty())
+            {
+                FatalErrorIn
+                (
+                    "void topoPatchMapper::"
+                    "calcInsertedFaceAddressing() const"
+                )   << " Mapping for inserted boundary face is incorrect."
+                    << " Found an empty masterObjects list."
+                    << nl << " Face: " << fffI.index()
+                    << nl << " Patch: " << patch_.name()
+                    << abort(FatalError);
+            }
+            else
+            {
+                // Make an entry for the inserted label,
+                // and renumber addressing to patch.
+                insertedFaces[nInsertedFaces] = fffI.index() - pStart;
+
+                // Add a mapping entry for facesFromFaces indices
+                insertedFacesMap[nInsertedFaces] = objectI;
+
+                // Make an entry for addressing
+                labelList& addr = insertedAddressing[nInsertedFaces];
+
+                // Renumber addressing to patch.
+                // Also, check mapping for hits into
+                // other patches / internal faces.
+                addr = fffI.masterObjects();
+
+                forAll(addr, faceI)
+                {
+                    if
+                    (
+                        addr[faceI] >= oldPatchStart
+                     && addr[faceI] < oldPatchEnd
+                    )
+                    {
+                        addr[faceI] -= oldPatchStart;
+                    }
+                    else
+                    {
+                        FatalErrorIn
+                        (
+                            "void topoPatchMapper::"
+                            "calcInsertedFaceAddressing() const"
+                        )
+                            << "Addressing into another patch is not allowed."
+                            << nl << " Patch face index: " << faceI
+                            << nl << " addr[faceI]: " << addr[faceI]
+                            << nl << " oldPatchStart: " << oldPatchStart
+                            << nl << " oldPatchSize: " << oldPatchSize
+                            << nl << " oldPatchEnd: " << oldPatchEnd
+                            << abort(FatalError);
+                    }
+                }
+
+                nInsertedFaces++;
+            }
+        }
+    }
+
+    // Shorten inserted faces to actual size
+    insertedFaces.setSize(nInsertedFaces);
+    insertedFacesMap.setSize(nInsertedFaces);
+    insertedAddressing.setSize(nInsertedFaces);
+}
+
+
+//- Calculate addressing for mapping
+void topoPatchMapper::calcAddressing() const
+{
+    if (directAddrPtr_ || interpolationAddrPtr_)
+    {
+        FatalErrorIn("void topoPatchMapper::calcAddressing() const")
+            << "Addressing already calculated."
+            << abort(FatalError);
+    }
+
+    // Information from the old patch
+    const label oldPatchSize = mpm_.oldPatchSizes()[patch_.index()];
+    const label oldPatchStart = mpm_.oldPatchStarts()[patch_.index()];
+    const label oldPatchEnd = oldPatchStart + oldPatchSize;
+
+    // Assemble the maps: slice to patch
+    if (direct())
+    {
+        // Direct mapping - slice to size
+        directAddrPtr_ = new labelList(patch_.patch().patchSlice(mpm_.faceMap()));
+
+        labelList& addr = *directAddrPtr_;
+
+        // Shift to local patch indices.
+        // Also, check mapping for hits into other patches / internal faces.
+        forAll (addr, faceI)
+        {
+            if
+            (
+                addr[faceI] >= oldPatchStart
+             && addr[faceI] < oldPatchEnd
+            )
+            {
+                addr[faceI] -= oldPatchStart;
+            }
+            else
+            {
+                FatalErrorIn
+                (
+                    "void topoPatchMapper::calcAddressing() const"
+                )
+                    << "Addressing into another patch is not allowed."
+                    << nl << " Patch face index: " << faceI
+                    << nl << " addr[faceI]: " << addr[faceI]
+                    << nl << " oldPatchStart: " << oldPatchStart
+                    << nl << " oldPatchSize: " << oldPatchSize
+                    << nl << " oldPatchEnd: " << oldPatchEnd
+                    << abort(FatalError);
+            }
+        }
+    }
+    else
+    {
+        // Interpolative addressing
+        interpolationAddrPtr_ = new labelListList(patchSize(), labelList(0));
+        labelListList& addr = *interpolationAddrPtr_;
+
+        // Fetch the list of inserted faces / addressing
+        const labelList& insertedFaces = insertedObjectLabels();
+        const labelListList& insertedAddressing = insertedFaceAddressing();
+
+        // Make entries
+        forAll(insertedFaces, faceI)
+        {
+            addr[insertedFaces[faceI]] = insertedAddressing[faceI];
+        }
+
+        // Do mapped faces. Note that this can already be set by insertedFaces
+        // so check if addressing size still zero.
+        const labelList& fm = patch_.patch().patchSlice(mpm_.faceMap());
+
+        forAll(fm, faceI)
+        {
+            if (fm[faceI] > -1 && addr[faceI].size() == 0)
+            {
+                // Mapped from a single face
+                label oldFace = fm[faceI];
+
+                if
+                (
+                    oldFace >= oldPatchStart
+                 && oldFace < oldPatchEnd
+                )
+                {
+                    oldFace -= oldPatchStart;
+                }
+                else
+                {
+                    FatalErrorIn
+                    (
+                        "void topoPatchMapper::calcAddressing() const"
+                    )
+                        << "Addressing into another patch is not allowed."
+                        << nl << " Patch face index: " << faceI
+                        << nl << " faceMap[faceI]: " << oldFace
+                        << nl << " oldPatchStart: " << oldPatchStart
+                        << nl << " oldPatchSize: " << oldPatchSize
+                        << nl << " oldPatchEnd: " << oldPatchEnd
+                        << abort(FatalError);
+                }
+
+                addr[faceI] = labelList(1, oldFace);
+            }
+        }
+
+        // Check if we missed anything
+        forAll(addr, faceI)
+        {
+            if (addr[faceI].empty())
+            {
+                FatalErrorIn
+                (
+                    "void topoPatchMapper::calcAddressing() const"
+                )
+                    << "Addressing is missing." << nl
+                    << " Patch face index: " << faceI << nl
+                    << " nInsertedFaces: " << insertedFaces.size() << nl
+                    << " faceMap: " << fm[faceI] << nl
+                    << " Patch: " << patch_.name() << nl
+                    << abort(FatalError);
+            }
+        }
+    }
+}
+
+
+//- Calculate inverse-distance weights for interpolative mapping
+void topoPatchMapper::calcInverseDistanceWeights() const
+{
+    if (weightsPtr_)
+    {
+        FatalErrorIn
+        (
+            "void topoPatchMapper::calcInverseDistanceWeights() const"
+        )
+            << "Weights already calculated."
+            << abort(FatalError);
+    }
+
+    // Fetch interpolative addressing
+    const labelListList& addr = addressing();
+
+    // Allocate memory
+    weightsPtr_ = new scalarListList(patchSize());
+    scalarListList& w = *weightsPtr_;
+
+    // Obtain face-centre information from old/new meshes
+    const vectorField& oldCentres = tMapper_.patchCentres(patch_.index());
+    const vectorField& newCentres = patch_.patch().faceCentres();
+
+    forAll(addr, faceI)
+    {
+        const labelList& mo = addr[faceI];
+
+        // Do mapped faces
+        if (mo.size() == 1)
+        {
+            w[faceI] = scalarList(1, 1.0);
+        }
+        else
+        {
+            // Map from masters, inverse-distance weights
+            scalar totalWeight = 0.0;
+            w[faceI] = scalarList(mo.size(), 0.0);
+
+            forAll (mo, oldFaceI)
+            {
+                w[faceI][oldFaceI] =
+                (
+                    1.0/stabilise
+                    (
+                        magSqr
+                        (
+                            newCentres[faceI]
+                          - oldCentres[mo[oldFaceI]]
+                        ),
+                        VSMALL
+                    )
+                );
+
+                totalWeight += w[faceI][oldFaceI];
+            }
+
+            // Normalize weights
+            scalar normFactor = (1.0/totalWeight);
+
+            forAll (mo, oldFaceI)
+            {
+                w[faceI][oldFaceI] *= normFactor;
+            }
+        }
+    }
+}
+
+
+//- Calculate intersection weights for conservative mapping
+void topoPatchMapper::calcIntersectionWeightsAndCentres() const
+{
+    if (areasPtr_ || centresPtr_)
+    {
+        FatalErrorIn
+        (
+            "void topoPatchMapper::"
+            "calcIntersectionWeightsAndCentres() const"
+        )
+            << "Weights already calculated."
+            << abort(FatalError);
+    }
+
+    // Fetch interpolative addressing
+    const labelListList& addr = addressing();
+
+    // Allocate memory
+    areasPtr_ = new List<scalarField>(patchSize(), scalarField(0));
+    List<scalarField>& a = *areasPtr_;
+
+    centresPtr_ = new List<vectorField>(patchSize(), vectorField(0));
+    List<vectorField>& x = *centresPtr_;
+
+    // Obtain stored face-centres
+    const vectorField& faceCentres = tMapper_.patchCentres(patch_.index());
+
+    // Fetch maps
+    const List<vectorField>& mapFaceCentres = tMapper_.faceCentres();
+    const List<scalarField>& mapFaceWeights = tMapper_.faceWeights();
+
+    // Fill in maps first
+    const labelList& insertedFaces = insertedObjectLabels();
+    const labelList& insertedFacesMap = insertedObjectMap();
+
+    forAll(insertedFaces, faceI)
+    {
+        a[insertedFaces[faceI]] = mapFaceWeights[insertedFacesMap[faceI]];
+        x[insertedFaces[faceI]] = mapFaceCentres[insertedFacesMap[faceI]];
+    }
+
+    // Now do mapped faces
+    forAll(addr, faceI)
+    {
+        const labelList& mo = addr[faceI];
+
+        // Check if this is indeed a mapped face
+        if (mo.size() == 1 && x[faceI].empty() && a[faceI].empty())
+        {
+            x[faceI] = vectorField(1, faceCentres[mo[0]]);
+            a[faceI] = scalarField(1, 1.0);
+        }
+    }
+
+    // Final check to ensure everything went okay
+    forAll(a, faceI)
+    {
+        if (a[faceI].empty())
+        {
+            FatalErrorIn
+            (
+                "void topoPatchMapper::"
+                "calcIntersectionWeightsAndCentres() const"
+            )
+                << "Weights / centres is missing."
+                << nl << " Patch face index: " << faceI
+                << abort(FatalError);
+        }
+    }
+}
+
+
+const List<scalarField>& topoPatchMapper::intersectionWeights() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const List<scalarField>& "
+            "topoPatchMapper::intersectionWeights() const"
+        )   << "Requested interpolative weights for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!areasPtr_)
+    {
+        calcIntersectionWeightsAndCentres();
+    }
+
+    return *areasPtr_;
+}
+
+
+const List<vectorField>& topoPatchMapper::intersectionCentres() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const List<vectorField>& "
+            "topoPatchMapper::intersectionCentres() const"
+        )   << "Requested interpolative weights for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!centresPtr_)
+    {
+        calcIntersectionWeightsAndCentres();
+    }
+
+    return *centresPtr_;
+}
+
+
+// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
+
+// Construct from components
+topoPatchMapper::topoPatchMapper
+(
+    const fvPatch& patch,
+    const mapPolyMesh& mpm,
+    const topoMapper& mapper
+)
+:
+    patch_(patch),
+    mpm_(mpm),
+    tMapper_(mapper),
+    direct_(false),
+    directAddrPtr_(NULL),
+    interpolationAddrPtr_(NULL),
+    weightsPtr_(NULL),
+    insertedFaceLabelsPtr_(NULL),
+    insertedFaceIndexMapPtr_(NULL),
+    insertedFaceAddressingPtr_(NULL),
+    areasPtr_(NULL),
+    centresPtr_(NULL)
+{
+    // Check for the possibility of direct mapping
+    if (insertedObjects())
+    {
+        direct_ = false;
+    }
+    else
+    {
+        direct_ = true;
+    }
+}
+
+// * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
+
+topoPatchMapper::~topoPatchMapper()
+{
+    clearOut();
+}
+
+// * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
+
+//- Return the polyPatch size
+label topoPatchMapper::patchSize() const
+{
+    return patch_.patch().size();
+}
+
+
+//- Return size
+label topoPatchMapper::size() const
+{
+    return patch_.size();
+}
+
+
+//- Return size before mapping
+label topoPatchMapper::sizeBeforeMapping() const
+{
+    if (patch_.type() == "empty")
+    {
+        return 0;
+    }
+
+    return mpm_.oldPatchSizes()[patch_.index()];
+}
+
+
+//- Is the mapping direct
+bool topoPatchMapper::direct() const
+{
+    return direct_;
+}
+
+
+//- Return direct addressing
+const unallocLabelList& topoPatchMapper::directAddressing() const
+{
+    if (!direct())
+    {
+        FatalErrorIn
+        (
+            "const unallocLabelList& "
+            "topoPatchMapper::directAddressing() const"
+        )   << "Requested direct addressing for an interpolative mapper."
+            << abort(FatalError);
+    }
+
+    if (!directAddrPtr_)
+    {
+        calcAddressing();
+    }
+
+    return *directAddrPtr_;
+}
+
+
+//- Return interpolation addressing
+const labelListList& topoPatchMapper::addressing() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const labelListList& "
+            "topoPatchMapper::addressing() const"
+        )   << "Requested interpolative addressing for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!interpolationAddrPtr_)
+    {
+        calcAddressing();
+    }
+
+    return *interpolationAddrPtr_;
+}
+
+
+//- Return weights
+const scalarListList& topoPatchMapper::weights() const
+{
+    if (direct())
+    {
+        FatalErrorIn
+        (
+            "const scalarListList& "
+            "topoPatchMapper::weights() const"
+        )   << "Requested interpolative weights for a direct mapper."
+            << abort(FatalError);
+    }
+
+    if (!weightsPtr_)
+    {
+        calcInverseDistanceWeights();
+    }
+
+    return *weightsPtr_;
+}
+
+
+//- Are there any inserted faces
+bool topoPatchMapper::insertedObjects() const
+{
+    return insertedObjectLabels().size();
+}
+
+
+//- Return list of inserted faces
+const labelList& topoPatchMapper::insertedObjectLabels() const
+{
+    if (!insertedFaceLabelsPtr_)
+    {
+        calcInsertedFaceAddressing();
+    }
+
+    return *insertedFaceLabelsPtr_;
+}
+
+
+//- Return addressing map for inserted faces
+const labelList& topoPatchMapper::insertedObjectMap() const
+{
+    if (!insertedFaceIndexMapPtr_)
+    {
+        calcInsertedFaceAddressing();
+    }
+
+    return *insertedFaceIndexMapPtr_;
+}
+
+
+//- Return addressing for inserted faces
+const labelListList& topoPatchMapper::insertedFaceAddressing() const
+{
+    if (!insertedFaceAddressingPtr_)
+    {
+        calcInsertedFaceAddressing();
+    }
+
+    return *insertedFaceAddressingPtr_;
+}
+
+
+//- Map the patch field
+template <class Type>
+void topoPatchMapper::mapPatchField
+(
+    const word& fieldName,
+    Field<Type>& pF
+) const
+{
+    // To invoke inverse-distance weighting, use this:
+    // pF.autoMap(*this);
+
+    // Check for possibility of direct mapping
+    if (direct())
+    {
+        pF.autoMap(*this);
+    }
+    else
+    {
+        if (pF.size() != sizeBeforeMapping())
+        {
+            FatalErrorIn
+            (
+                "\n\n"
+                "void topoCellMapper::mapPatchField<Type>\n"
+                "(\n"
+                "    const word& fieldName,\n"
+                "    Field<Type>& iF\n"
+                ") const\n"
+            )  << "Incompatible size before mapping." << nl
+               << " Field: " << fieldName << nl
+               << " Field size: " << pF.size() << nl
+               << " map size: " << sizeBeforeMapping() << nl
+               << abort(FatalError);
+        }
+
+        // Fetch addressing
+        const labelListList& pAddressing = addressing();
+        const List<scalarField>& wF = intersectionWeights();
+
+        // Compute the integral of the source field
+        Type intSource = sum(pF * tMapper_.patchAreas(patch_.index()));
+
+        // Copy the original field
+        Field<Type> fieldCpy(pF);
+
+        // Resize to current dimensions
+        pF.setSize(size());
+
+        // Map the patch field
+        forAll(pF, faceI)
+        {
+            const labelList& addr = pAddressing[faceI];
+
+            pF[faceI] = pTraits<Type>::zero;
+
+            // Accumulate area-weighted interpolate
+            forAll(addr, faceJ)
+            {
+                pF[faceI] +=
+                (
+                    wF[faceI][faceJ] * fieldCpy[addr[faceJ]]
+                );
+            }
+        }
+
+        // Compute the integral of the target field
+        const polyPatch& ppI = mpm_.mesh().boundaryMesh()[patch_.index()];
+
+        Type intTarget = sum(pF * mag(ppI.faceAreas()));
+
+        if (polyMesh::debug)
+        {
+            int oldP = Info().precision();
+
+            // Compare the global integral
+            Info << " Field : " << fieldName
+                 << " Patch : " << ppI.name()
+                 << " integral errors : "
+                 << setprecision(15)
+                 << " source : " << mag(intSource)
+                 << " target : " << mag(intTarget)
+                 << " norm : "
+                 << (mag(intTarget - intSource) / (mag(intSource) + VSMALL))
+                 << setprecision(oldP)
+                 << endl;
+        }
+    }
+}
+
+
+// * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
+
+void topoPatchMapper::operator=(const topoPatchMapper& rhs)
+{
+    // Check for assignment to self
+    if (this == &rhs)
+    {
+        FatalErrorIn("void topoPatchMapper::operator=")
+            << "Attempted assignment to self"
+            << abort(FatalError);
+    }
+}
+
+} // End namespace Foam
+
+// ************************************************************************* //

src/dynamicMesh/dynamicFvMesh/dynamicTopoFvMesh/fieldMapping/topoPatchMapper.H

@@ -0,0 +1,207 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright held by original author
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    topoPatchMapper
+
+Description
+    This object provides mapping and fill-in information for patch data
+    between the two meshes after the topological change.  It is
+    constructed from mapPolyMesh.
+
+Author
+    Sandeep Menon
+    University of Massachusetts Amherst
+    All rights reserved
+
+SourceFiles
+    topoPatchMapper.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef topoPatchMapper_H
+#define topoPatchMapper_H
+
+#include "morphFieldMapper.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+// Class forward declarations
+class topoMapper;
+class mapPolyMesh;
+
+/*---------------------------------------------------------------------------*\
+                    Class topoPatchMapper Declaration
+\*---------------------------------------------------------------------------*/
+
+class topoPatchMapper
+:
+    public fvPatchFieldMapper
+{
+    // Private data
+
+        //- Reference to patch
+        const fvPatch& patch_;
+
+        //- Reference to mapPolyMesh
+        const mapPolyMesh& mpm_;
+
+        //- Reference to the topoMapper
+        const topoMapper& tMapper_;
+
+        //- Is the mapping direct
+        bool direct_;
+
+    // Demand-driven private data
+
+        //- Direct addressing
+        mutable labelList* directAddrPtr_;
+
+        //- Interpolated addressing
+        mutable labelListList* interpolationAddrPtr_;
+
+        //- Inverse-distance weights
+        mutable scalarListList* weightsPtr_;
+
+        //- Inserted faces
+        mutable labelList* insertedFaceLabelsPtr_;
+
+        //- Inserted face index map
+        mutable labelList* insertedFaceIndexMapPtr_;
+
+        //- Inserted face addressing
+        mutable labelListList* insertedFaceAddressingPtr_;
+
+        //- Interpolation areas
+        mutable List<scalarField>* areasPtr_;
+
+        //- Interpolation centres
+        mutable List<vectorField>* centresPtr_;
+
+    // Private Member Functions
+
+        //- Disallow default bitwise copy construct
+        topoPatchMapper(const topoPatchMapper&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const topoPatchMapper&);
+
+        //- Calculate the insertedFace addressing list
+        void calcInsertedFaceAddressing() const;
+
+        //- Calculate addressing for mapping
+        void calcAddressing() const;
+
+        //- Calculate inverse-distance weights for interpolative mapping
+        void calcInverseDistanceWeights() const;
+
+        //- Calculate intersection weights for conservative mapping
+        void calcIntersectionWeightsAndCentres() const;
+
+        //- Return intersection area weights
+        const List<scalarField>& intersectionWeights() const;
+
+        //- Return intersection area centres
+        const List<vectorField>& intersectionCentres() const;
+
+        //- Clear out local storage
+        void clearOut();
+
+public:
+
+    // Constructors
+
+        //- Construct from mapPolyMesh
+        topoPatchMapper
+        (
+            const fvPatch& patch,
+            const mapPolyMesh& mpm,
+            const topoMapper& mapper
+        );
+
+    // Destructor
+
+        virtual ~topoPatchMapper();
+
+    // Member Functions
+
+        //- Return the polyPatch size
+        label patchSize() const;
+
+        //- Return size
+        virtual label size() const;
+
+        //- Return size of field before mapping
+        virtual label sizeBeforeMapping() const;
+
+        //- Is the mapping direct
+        virtual bool direct() const;
+
+        //- Return direct addressing
+        virtual const unallocLabelList& directAddressing() const;
+
+        //- Return interpolated addressing
+        virtual const labelListList& addressing() const;
+
+        //- Return interpolaion weights
+        virtual const scalarListList& weights() const;
+
+        //- Are there any inserted faces
+        virtual bool insertedObjects() const;
+
+        //- Return list of inserted faces
+        virtual const labelList& insertedObjectLabels() const;
+
+        //- Return addressing map for inserted faces
+        const labelList& insertedObjectMap() const;
+
+        //- Return addressing for inserted faces
+        const labelListList& insertedFaceAddressing() const;
+
+        //- Map the patch field
+        template <class Type>
+        void mapPatchField
+        (
+            const word& fieldName,
+            Field<Type>& pF
+        ) const;
+};
+
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+#ifdef NoRepository
+#include "topoPatchMapper.C"
+#endif
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //