Merge branch 'master' of ssh://noisy/home/noisy3/OpenFOAM/OpenFOAM-dev

applications/solvers/combustion/PDRFoam/PDRModels/XiEqModels/basicXiSubXiEq/basicXiSubXiEq.H

@@ -29,6 +29,36 @@ Description
     Basic sub-grid obstacle flame-wrinking enhancement factor model.
     Details supplied by J Puttock 2/7/06.
 
+    <b> Sub-grid flame area generation </b>
+
+    \f$ n = N - \hat{\dwea{\vec{U}}}.n_{s}.\hat{\dwea{\vec{U}}} \f$
+    \f$ n_{r} = \sqrt{n} \f$
+
+    where:
+
+        \f$ \hat{\dwea{\vec{U}}} = \dwea{\vec{U}} / \vert \dwea{\vec{U}}
+        \vert \f$
+
+        \f$ b = \hat{\dwea{\vec{U}}}.B.\hat{\dwea{\vec{U}}} / n_{r} \f$
+
+    where:
+
+        \f$ B \f$ is the file "B".
+
+        \f$ N \f$ is the file "N".
+
+        \f$  n_{s} \f$ is the file "ns".
+
+    The flame area enhancement factor \f$ \Xi_{sub} \f$ is expected to
+    approach:
+
+    \f[
+        \Xi_{{sub}_{eq}} =
+            1 + max(2.2 \sqrt{b}, min(0.34 \frac{\vert \dwea{\vec{U}}
+            \vert}{{\vec{U}}^{'}}, 1.6)) \times min(\frac{n}{4}, 1)
+    \f]
+
+
 SourceFiles
     basicSubGrid.C
 

applications/solvers/combustion/PDRFoam/PDRModels/XiGModels/basicXiSubG/basicXiSubG.H

@@ -25,10 +25,28 @@ License
 Class
     basicSubGrid
 
+
 Description
+
     Basic sub-grid obstacle flame-wrinking generation rate coefficient model.
     Details supplied by J Puttock 2/7/06.
 
+    \f$ G_{sub} \f$ denotes the generation coefficient and it is given by
+
+    \f[
+        G_{sub} = k_{1} /frac{\vert \dwea{\vec{U}} \vert}{L_{obs}}
+                 \frac{/Xi_{{sub}_{eq}}-1}{/Xi_{sub}}
+    \f]
+
+    and the removal:
+
+    \f[ - k_{1} /frac{\vert \dwea{\vec{U}} \vert}{L_{sub}}
+    \frac{\Xi_{sub}-1}{\Xi_{sub}} \f]
+
+    Finally, \f$ G_{sub} \f$ is added to generation rate \f$ G_{in} \f$
+    due to the turbulence.
+
+
 SourceFiles
     basicSubGrid.C
 

applications/solvers/combustion/PDRFoam/PDRModels/dragModels/basic/basic.H

@@ -29,6 +29,50 @@ Description
     Basic sub-grid obstacle drag model.
     Details supplied by J Puttock 2/7/06.
 
+    <b> Sub-grid drag term </b>
+
+    The resistance term (force per unit of volume) is given by:
+
+    \f[
+        R = -\frac{1}{2} \rho \vert \dwea{\vec{U}} \vert \dwea{\vec{U}}.D
+    \f]
+
+    where:
+
+        \f$ D \f$ is the tensor field "CR" in \f$ m^{-1} \f$
+
+        This is term is treated implicitly in UEqn.H
+
+    <b> Sub-grid turbulence generation </b>
+
+    The turbulence source term \f$ G_{R} \f$ occurring in the
+    \f$ \kappa-\epsilon \f$ equations for the generation of turbulence due
+    to interaction with unresolved obstacles :
+
+    \f$ G_{R} = C_{s}\beta_{\nu}
+    \mu_{eff} A_{w}^{2}(\dwea{\vec{U}}-\dwea{\vec{U}_{s}})^2 + \frac{1}{2}
+    \rho \vert \dwea{\vec{U}} \vert \dwea{\vec{U}}.T.\dwea{\vec{U}} \f$
+
+    where:
+
+        \f$ C_{s} \f$ = 1
+
+        \f$ \beta_{\nu} \f$ is the volume porosity (file "betav").
+
+        \f$ \mu_{eff} \f$ is the effective viscosity.
+
+        \f$ A_{w}^{2}\f$ is the obstacle surface area per unit of volume
+        (file "Aw").
+
+        \f$ \dwea{\vec{U}_{s}} \f$ is the slip velocity and is considered
+        \f$ \frac{1}{2}. \dwea{\vec{U}} \f$.
+
+        \f$ T \f$ is a tensor in the file CT.
+
+    The term \f$ G_{R} \f$ is treated explicitly in the \f$ \kappa-\epsilon
+    \f$ Eqs in the PDRkEpsilon.C file.
+
+
 SourceFiles
     basic.C
 
@@ -40,7 +84,6 @@ SourceFiles
 #include "PDRDragModel.H"
 #include "XiEqModel.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {

applications/solvers/combustion/PDRFoam/PDRModels/turbulence/PDRkEpsilon/PDRkEpsilon.H

@@ -26,7 +26,17 @@ Class
     PDRkEpsilon
 
 Description
-    Standard k-epsilon turbulence model.
+    Standard k-epsilon turbulence model with additional source terms
+    corresponding to PDR basic drag model (basic.H)
+
+    The turbulence source term \f$ G_{R} \f$ appears in the
+    \f$ \kappa-\epsilon \f$ equation for the generation of turbulence due to
+    interaction with unresolved obstacles.
+
+    In the \f$ \epsilon  \f$ equation \f$ C_{1} G_{R} \f$ is added as a source
+    term.
+
+    In the \f$ \kappa \f$ equation \f$ G_{R} \f$ is added as a source term.
 
 SourceFiles
     PDRkEpsilon.C

applications/solvers/combustion/PDRFoam/XiModels/XiModel/XiModel.H

@@ -27,6 +27,57 @@ Class
 
 Description
     Base-class for all Xi models used by the b-Xi combustion model.
+    See Technical Report SH/RE/01R for details on the PDR modelling.
+
+    Xi is given through an algebraic expression (algebraic.H),
+    by solving a transport equation (transport.H) or a fixed value (fixed.H).
+    See report TR/HGW/10 for details on the Weller two equations model.
+
+    In the algebraic and transport methods \f$\Xi_{eq}\f$ is calculated in
+    similar way. In the algebraic approach, \f$\Xi_{eq}\f$ is the value used in
+    the \f$ b \f$ transport equation.
+
+    \f$\Xi_{eq}\f$ is calculated as follows:
+
+    \f$\Xi_{eq} = 1 + (1 + 2\Xi_{coeff}(0.5 - \dwea{b}))(\Xi^* - 1)\f$
+
+    where:
+
+        \f$ \dwea{b} \f$ is the regress variable.
+
+        \f$ \Xi^* \f$ is the total equilibrium wrinkling combining the effects
+        of the flame inestability and turbulence interaction and is given by
+
+        \f[
+            \Xi^* = \frac {R}{R - G_\eta - G_{in}}
+        \f]
+
+    where:
+
+        \f$ G_\eta \f$ is the generation rate of wrinkling due to turbulence
+        interaction.
+
+        \f$ G_{in} = \kappa \rho_{u}/\rho_{b} \f$ is the generation
+         rate due to the flame inestability.
+
+    By adding the removal rates of the two effects:
+
+        \f[
+            R = G_\eta \frac{\Xi_{\eta_{eq}}}{\Xi_{\eta_{eq}} - 1}
+              + G_{in} \frac{\Xi_{{in}_{eq}}}{\Xi_{{in}_{eq}} - 1}
+        \f]
+
+    where:
+
+        \f$ R \f$ is the total removal.
+
+        \f$ G_\eta \f$ is a model constant.
+
+        \f$ \Xi_{\eta_{eq}} \f$ is the flame wrinkling due to turbulence.
+
+        \f$ \Xi_{{in}_{eq}} \f$ is the equilibrium level of the flame wrinkling
+        generated by inestability. It is a constant (default 2.5).
+
 
 SourceFiles
     XiModel.C
@@ -51,6 +102,8 @@ namespace Foam
                           Class XiModel Declaration
 \*---------------------------------------------------------------------------*/
 
+
+
 class XiModel
 {
 

applications/solvers/combustion/PDRFoam/laminarFlameSpeed/SCOPE/SCOPELaminarFlameSpeed.H

@@ -28,6 +28,33 @@ Class
 Description
     Laminar flame speed obtained from the SCOPE correlation.
 
+    Seven parameters are specified in terms of polynomial functions of
+    stoichiometry. Two polynomials are fitted, covering different parts of the
+    flammable range. If the mixture is outside the fitted range, linear
+    interpolation is used between the extreme of the polynomio and the upper or
+    lower flammable limit with the Markstein number constant.
+
+    Variations of pressure and temperature from the reference values are taken
+    into account through \f$ pexp \f$ and \f$ texp \f$
+
+    The laminar burning velocity fitting polynomial is:
+
+    \f$ Su = a_{0}(1+a_{1}x+K+..a_{i}x^{i}..+a_{6}x^{6}) (p/p_{ref})^{pexp}
+    (T/T_{ref})^{texp} \f$
+
+    where:
+
+        \f$ a_{i} \f$ are the polinomial coefficients.
+
+        \f$ pexp \f$ and \f$ texp \f$ are the pressure and temperature factors
+        respectively.
+
+        \f$ x \f$ is the equivalence ratio.
+
+        \f$ T_{ref} \f$ and \f$ p_{ref} \f$ are the temperature and pressure
+        references for the laminar burning velocity.
+
+
 SourceFiles
     SCOPELaminarFlameSpeed.C
 
@@ -125,7 +152,7 @@ class SCOPE
         //  corrected for temperature and pressure dependence
         inline scalar Su0pTphi(scalar p, scalar Tu, scalar phi) const;
 
-        //- Laminar flame speed evaluated from the given uniform 
+        //- Laminar flame speed evaluated from the given uniform
         //  equivalence ratio corrected for temperature and pressure dependence
         tmp<volScalarField> Su0pTphi
         (
@@ -134,7 +161,7 @@ class SCOPE
             scalar phi
         ) const;
 
-        //- Laminar flame speed evaluated from the given equivalence ratio 
+        //- Laminar flame speed evaluated from the given equivalence ratio
         //  distribution corrected for temperature and pressure dependence
         tmp<volScalarField> Su0pTphi
         (
@@ -144,7 +171,7 @@ class SCOPE
         ) const;
 
         //- Return the Markstein number
-        //  evaluated from the given equivalence ratio 
+        //  evaluated from the given equivalence ratio
         tmp<volScalarField> Ma(const volScalarField& phi) const;
 
         //- Construct as copy (not implemented)

applications/test/DynamicList/DynamicListTest.C

@@ -43,16 +43,21 @@ int main(int argc, char *argv[])
     ldl[0](3) = 3;
     ldl[0](1) = 1;
 
-    ldl[0].setSize(5);     // increase allocated size
+    ldl[0].setCapacity(5);    // increase allocated size
-    ldl[1].setSize(10);    // increase allocated size
+    ldl[1].setCapacity(10);   // increase allocated size
-    ldl[1](2) = 2;
+    ldl[0].reserve(15);       // should increase allocated size
+    ldl[1].reserve(5);        // should not decrease allocated size
+    ldl[1](3) = 2;            // allocates space and sets value
+
+    // this works without a segfault, but doesn't change the list size
+    ldl[0][4] = 4;
 
     ldl[1] = 3;
 
     Info<< "<ldl>" << ldl << "</ldl>" << nl << "sizes: ";
     forAll(ldl, i)
     {
-        Info<< " " << ldl[i].size() << "/" << ldl[i].allocSize();
+        Info<< " " << ldl[i].size() << "/" << ldl[i].capacity();
     }
     Info<< endl;
 
@@ -63,7 +68,7 @@ int main(int argc, char *argv[])
     Info<< "<ldl>" << ldl << "</ldl>" << nl << "sizes: ";
     forAll(ldl, i)
     {
-        Info<< " " << ldl[i].size() << "/" << ldl[i].allocSize();
+        Info<< " " << ldl[i].size() << "/" << ldl[i].capacity();
     }
     Info<< endl;
 
@@ -78,10 +83,10 @@ int main(int argc, char *argv[])
     {
         dlA.append(i);
     }
-    dlA.setSize(10);
+    dlA.setCapacity(10);
 
     Info<< "<dlA>" << dlA << "</dlA>" << nl << "sizes: "
-        << " " << dlA.size() << "/" << dlA.allocSize() << endl;
+        << " " << dlA.size() << "/" << dlA.capacity() << endl;
 
     dlB.transfer(dlA);
 
@@ -91,10 +96,54 @@ int main(int argc, char *argv[])
 
     Info<< "Transferred to dlB" << endl;
     Info<< "<dlA>" << dlA << "</dlA>" << nl << "sizes: "
-        << " " << dlA.size() << "/" << dlA.allocSize() << endl;
+        << " " << dlA.size() << "/" << dlA.capacity() << endl;
     Info<< "<dlB>" << dlB << "</dlB>" << nl << "sizes: "
-        << " " << dlB.size() << "/" << dlB.allocSize() << endl;
+        << " " << dlB.size() << "/" << dlB.capacity() << endl;
 
+    // try with a normal list:
+    List<label> lstA;
+    lstA.transfer(dlB);
+    Info<< "Transferred to normal list" << endl;
+    Info<< "<lstA>" << lstA << "</lstA>" << nl << "sizes: "
+        << " " << lstA.size() << endl;
+    Info<< "<dlB>" << dlB << "</dlB>" << nl << "sizes: "
+        << " " << dlB.size() << "/" << dlB.capacity() << endl;
+
+    // Copy back and append a few time
+    for (label i=0; i < 3; i++)
+    {
+        dlB.append(lstA);
+    }
+
+    Info<< "appended list a few times" << endl;
+    Info<< "<dlB>" << dlB << "</dlB>" << nl << "sizes: "
+        << " " << dlB.size() << "/" << dlB.capacity() << endl;
+
+    // assign the list (should maintain allocated space)
+    dlB = lstA;
+    Info<< "assigned list" << endl;
+    Info<< "<dlB>" << dlB << "</dlB>" << nl << "sizes: "
+        << " " << dlB.size() << "/" << dlB.capacity() << endl;
+
+
+    // Copy back and append a few time
+    for (label i=0; i < 3; i++)
+    {
+        dlB.append(lstA);
+    }
+
+
+    // check allocation granularity
+    DynamicList<label, 6, 0> dlC;
+
+    Info<< "<dlC>" << dlC << "</dlC>" << nl << "sizes: "
+        << " " << dlC.size() << "/" << dlC.capacity() << endl;
+
+    dlC.reserve(dlB.size());
+    dlC = dlB;
+
+    Info<< "<dlC>" << dlC << "</dlC>" << nl << "sizes: "
+        << " " << dlC.size() << "/" << dlC.capacity() << endl;
 
     return 0;
 }

applications/test/dictionary/dictionaryTest.C

@@ -45,6 +45,27 @@ int main(int argc, char *argv[])
 
     Info<< testDict << endl;
 
+    {
+        dictionary someDict;
+        someDict.add(keyType("a.*", true), "subdictValue");
+
+        dictionary dict;
+        dict.add("someDict", someDict);
+        dict.add(keyType(".*", true), "parentValue");
+
+        Info<< "dict:" << dict << endl;
+
+        // Wildcard find.
+        Info<< "Wildcard find \"abc\" in top directory : "
+            << dict.lookup("abc") << endl;
+        Info<< "Wildcard find \"abc\" in sub directory : "
+            << dict.subDict("someDict").lookup("abc")
+            << endl;
+        Info<< "Recursive wildcard find \"def\" in sub directory : "
+            << dict.subDict("someDict").lookup("def", true)
+            << endl;
+    }
+
     return 0;
 }
 

applications/test/sort/sortListTest.C

@@ -27,6 +27,7 @@ Description
 \*---------------------------------------------------------------------------*/
 
 #include "SortableList.H"
+#include "ListOps.H"
 
 using namespace Foam;
 
@@ -35,42 +36,101 @@ using namespace Foam;
 
 int main(int argc, char *argv[])
 {
-    labelList orig(5);
+    labelList orig(8);
     orig[0] = 7;
-    orig[1] = 4;
+    orig[1] = 9;
     orig[2] = 1;
     orig[3] = 2;
-    orig[4] = 9;
+    orig[4] = 4;
+    orig[5] = 7;
+    orig[6] = 4;
+    orig[7] = 0;
+
+    labelList order;
 
     labelList a(orig);
-    Info << "before: " << a << endl;
+    sortedOrder(a, order);
+
+    Info<< "unsorted: " << a << endl;
     sort(a);
-    Info << "after:  " << a << endl;
+    Info<< "sorted:   " << a << endl;
+    Info<< "indices:  " << order << endl;
 
     SortableList<label> b(orig);
-    Info << "sorted:  " << b << endl;
+    Info<< "unsorted: " << orig << endl;
-    Info << "indices: " << b.indices() << endl;
+    Info<< "sorted:   " << b << endl;
+    Info<< "indices:  " << b.indices() << endl;
 
-    Info << "shrunk:  " << b.shrink() << endl;
+    Info<< "shrunk:   " << b.shrink() << endl;
-    Info << "indices: " << b.indices() << endl;
+    Info<< "indices:  " << b.indices() << endl;
 
     // repeat by assignment
     b = orig;
-    Info << "unsorted: " << b << endl;
+    Info<< "unsorted: " << b << endl;
     b.sort();
 
-    Info << "sorted:   " << b << endl;
+    Info<< "sorted:   " << b << endl;
-    Info << "indices:  " << b.indices() << endl;
+    Info<< "indices:  " << b.indices() << endl;
+
+    // find unique/duplicate values
+    b = orig;
+
+    Info<< "unsorted: " << b << endl;
+    uniqueOrder(b, order);
+    Info<< "unique:  " << order << endl;
+    duplicateOrder(b, order);
+    Info<< "duplicate:" << order << endl;
+
+    // sort reverse
+    Info<< "unsorted: " << b << endl;
+    b.reverseSort();
+    Info<< "rsort:    " << b << endl;
+    Info<< "indices:  " << b.indices() << endl;
 
     // transfer assignment
-    b.transfer(orig);
+    a = orig;
-    Info << "unsorted: " << b << endl;
+    b.transfer(a);
+    Info<< "unsorted: " << b << endl;
     b.sort();
 
-    Info << "sorted:   " << b << endl;
+    Info<< "sorted:   " << b << endl;
-    Info << "indices:  " << b.indices() << endl;
+    Info<< "indices:  " << b.indices() << endl;
 
-    Info << "End\n" << endl;
+    a.transfer(b);
+
+    Info<< "plain:    " << a << endl;
+    Info<< "sorted:   " << b << endl;
+    Info<< "indices:  " << b.indices() << endl;
+
+    // sort/duplicate/unique with identical values
+    b.setSize(8);
+    b = 5;
+
+    Info<< "unsorted: " << b << endl;
+
+    uniqueOrder(b, order);
+    Info<< "unique:  " << order << endl;
+    duplicateOrder(b, order);
+    Info<< "duplicate:" << order << endl;
+    b.sort();
+
+    Info<< "sorted:   " << b << endl;
+    Info<< "indices:  " << b.indices() << endl;
+
+    // with a single value
+    b.setSize(1);
+
+    Info<< "unsorted: " << b << endl;
+    uniqueOrder(b, order);
+    Info<< "unique:  " << order << endl;
+    duplicateOrder(b, order);
+    Info<< "duplicate:" << order << endl;
+    b.sort();
+
+    Info<< "sorted:   " << b << endl;
+    Info<< "indices:  " << b.indices() << endl;
+
+    Info<< "\nEnd\n" << endl;
 
     return 0;
 }

applications/test/xfer/xferListTest.C

@@ -108,13 +108,13 @@ int main(int argc, char *argv[])
     face f1(dl);
     face f2(xferCopy<labelList>(dl));
 
-    Info<< "dl[" << dl.size() << "/" << dl.allocSize() << "] " << dl << endl;
+    Info<< "dl[" << dl.size() << "/" << dl.capacity() << "] " << dl << endl;
     Info<< "f1: " << f1 << endl;
     Info<< "f2: " << f2 << endl;
 
     // note: using xferMoveTo to ensure the correct transfer() method is called
     face f3( xferMoveTo<labelList>(dl) );
-    Info<< "dl[" << dl.size() << "/" << dl.allocSize() << "] " << dl << endl;
+    Info<< "dl[" << dl.size() << "/" << dl.capacity() << "] " << dl << endl;
     Info<< "f3: " << f3 << endl;
 
     return 0;

applications/utilities/mesh/conversion/polyDualMesh/Make/files

@@ -1,3 +1,4 @@
-polyDualMeshApp.C
+meshDualiser.C
+makePolyDualMesh.C
 
 EXE = $(FOAM_APPBIN)/polyDualMesh

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

@@ -1,6 +1,9 @@
 EXE_INC = \
     -I$(LIB_SRC)/meshTools/lnInclude \
-    -I$(LIB_SRC)/conversion/lnInclude
+    -I$(LIB_SRC)/finiteVolume/lnInclude \
+    -I$(LIB_SRC)/dynamicMesh/lnInclude
 
 EXE_LIBS = \
-    -lmeshTools -lconversion
+    -lfiniteVolume \
+    -ldynamicMesh \
+    -lmeshTools

applications/utilities/mesh/conversion/polyDualMesh/makePolyDualMesh.C

@@ -0,0 +1,515 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2007 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Description
+    Calculate the dual of a polyMesh. Adheres to all the feature&patch edges.
+
+    Feature angle:
+    convex features : point becomes single boundary cell with multiple
+                      boundary faces.
+    concave features: point becomes multiple boundary cells.
+
+    -splitAllFaces:
+    Normally only constructs a single face between two cells. This single face
+    might be too distorted. splitAllFaces will create a single face for every
+    original cell the face passes through. The mesh will thus have
+    multiple faces inbetween two cells! (so is not strictly upper-triangular
+    anymore - checkMesh will complain)
+    -doNotPreserveFaceZones:
+    By default all faceZones are preserved by marking all faces, edges and
+    points on them as features. The -doNotPreserveFaceZones disables this
+    behaviour.
+
+\*---------------------------------------------------------------------------*/
+
+#include "argList.H"
+#include "Time.H"
+#include "fvMesh.H"
+#include "mathematicalConstants.H"
+#include "polyTopoChange.H"
+#include "mapPolyMesh.H"
+#include "PackedList.H"
+#include "meshTools.H"
+#include "OFstream.H"
+#include "meshDualiser.H"
+#include "ReadFields.H"
+#include "volFields.H"
+#include "surfaceFields.H"
+
+using namespace Foam;
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+// Naive feature detection. All boundary edges with angle > featureAngle become
+// feature edges. All points on feature edges become feature points. All
+// boundary faces become feature faces.
+void simpleMarkFeatures
+(
+    const polyMesh& mesh,
+    const PackedList<1>& isBoundaryEdge,
+    const scalar featureAngle,
+    const bool doNotPreserveFaceZones,
+
+    labelList& featureFaces,
+    labelList& featureEdges,
+    labelList& singleCellFeaturePoints,
+    labelList& multiCellFeaturePoints
+)
+{
+    scalar minCos = Foam::cos(featureAngle * mathematicalConstant::pi/180.0);
+
+    const polyBoundaryMesh& patches = mesh.boundaryMesh();
+
+    // Working sets
+    labelHashSet featureEdgeSet;
+    labelHashSet singleCellFeaturePointSet;
+    labelHashSet multiCellFeaturePointSet;
+
+
+    // 1. Mark all edges between patches
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    forAll(patches, patchI)
+    {
+        const polyPatch& pp = patches[patchI];
+        const labelList& meshEdges = pp.meshEdges();
+
+        // All patch corner edges. These need to be feature points & edges!
+        for (label edgeI = pp.nInternalEdges(); edgeI < pp.nEdges(); edgeI++)
+        {
+            label meshEdgeI = meshEdges[edgeI];
+            featureEdgeSet.insert(meshEdgeI);
+            singleCellFeaturePointSet.insert(mesh.edges()[meshEdgeI][0]);
+            singleCellFeaturePointSet.insert(mesh.edges()[meshEdgeI][1]);
+        }
+    }
+
+
+
+    // 2. Mark all geometric feature edges
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+    // Make distinction between convex features where the boundary point becomes
+    // a single cell and concave features where the boundary point becomes
+    // multiple 'half' cells.
+
+    // Addressing for all outside faces
+    primitivePatch allBoundary
+    (
+        SubList<face>
+        (
+            mesh.faces(),
+            mesh.nFaces()-mesh.nInternalFaces(),
+            mesh.nInternalFaces()
+        ),
+        mesh.points()
+    );
+
+    // Check for non-manifold points (surface pinched at point)
+    allBoundary.checkPointManifold(false, &singleCellFeaturePointSet);
+
+    // Check for non-manifold edges (surface pinched at edge)
+    const labelListList& edgeFaces = allBoundary.edgeFaces();
+    const labelList& meshPoints = allBoundary.meshPoints();
+
+    forAll(edgeFaces, edgeI)
+    {
+        const labelList& eFaces = edgeFaces[edgeI];
+
+        if (eFaces.size() > 2)
+        {
+            const edge& e = allBoundary.edges()[edgeI];
+
+            //Info<< "Detected non-manifold boundary edge:" << edgeI
+            //    << " coords:"
+            //    << allBoundary.points()[meshPoints[e[0]]]
+            //    << allBoundary.points()[meshPoints[e[1]]] << endl;
+
+            singleCellFeaturePointSet.insert(meshPoints[e[0]]);
+            singleCellFeaturePointSet.insert(meshPoints[e[1]]);
+        }
+    }
+
+    // Check for features.
+    forAll(edgeFaces, edgeI)
+    {
+        const labelList& eFaces = edgeFaces[edgeI];
+
+        if (eFaces.size() == 2)
+        {
+            label f0 = eFaces[0];
+            label f1 = eFaces[1];
+
+            // check angle
+            const vector& n0 = allBoundary.faceNormals()[f0];
+            const vector& n1 = allBoundary.faceNormals()[f1];
+
+            if ((n0 & n1) < minCos)
+            {
+                const edge& e = allBoundary.edges()[edgeI];
+                label v0 = meshPoints[e[0]];
+                label v1 = meshPoints[e[1]];
+
+                label meshEdgeI = meshTools::findEdge(mesh, v0, v1);
+                featureEdgeSet.insert(meshEdgeI);
+
+                // Check if convex or concave by looking at angle
+                // between face centres and normal
+                vector c1c0
+                (
+                    allBoundary[f1].centre(allBoundary.points())
+                  - allBoundary[f0].centre(allBoundary.points())
+                );
+
+                if ((c1c0 & n0) > SMALL)
+                {
+                    // Found concave edge. Make into multiCell features
+                    Info<< "Detected concave feature edge:" << edgeI
+                        << " cos:" << (c1c0 & n0)
+                        << " coords:"
+                        << allBoundary.points()[v0]
+                        << allBoundary.points()[v1]
+                        << endl;
+
+                    singleCellFeaturePointSet.erase(v0);
+                    multiCellFeaturePointSet.insert(v0);
+                    singleCellFeaturePointSet.erase(v1);
+                    multiCellFeaturePointSet.insert(v1);
+                }
+                else
+                {
+                    // Convex. singleCell feature.
+                    if (!multiCellFeaturePointSet.found(v0))
+                    {
+                        singleCellFeaturePointSet.insert(v0);
+                    }
+                    if (!multiCellFeaturePointSet.found(v1))
+                    {
+                        singleCellFeaturePointSet.insert(v1);
+                    }
+                }
+            }
+        }
+    }
+
+
+    // 3. Mark all feature faces
+    // ~~~~~~~~~~~~~~~~~~~~~~~~~
+
+    // Face centres that need inclusion in the dual mesh
+    labelHashSet featureFaceSet(mesh.nFaces()-mesh.nInternalFaces());
+    // A. boundary faces.
+    for (label faceI = mesh.nInternalFaces(); faceI < mesh.nFaces(); faceI++)
+    {
+        featureFaceSet.insert(faceI);
+    }
+
+    // B. face zones.
+    const faceZoneMesh& faceZones = mesh.faceZones();
+
+    if (doNotPreserveFaceZones)
+    {
+        if (faceZones.size() > 0)
+        {
+            WarningIn("simpleMarkFeatures(..)")
+                << "Detected " << faceZones.size()
+                << " faceZones. These will not be preserved."
+                << endl;
+        }
+    }
+    else
+    {
+        if (faceZones.size() > 0)
+        {
+            Info<< "Detected " << faceZones.size()
+                << " faceZones. Preserving these by marking their"
+                << " points, edges and faces as features." << endl;
+        }
+
+        forAll(faceZones, zoneI)
+        {
+            const faceZone& fz = faceZones[zoneI];
+
+            Info<< "Inserting all faces in faceZone " << fz.name()
+                << " as features." << endl;
+
+            forAll(fz, i)
+            {
+                label faceI = fz[i];
+                const face& f = mesh.faces()[faceI];
+                const labelList& fEdges = mesh.faceEdges()[faceI];
+
+                featureFaceSet.insert(faceI);
+                forAll(f, fp)
+                {
+                    // Mark point as multi cell point (since both sides of
+                    // face should have different cells)
+                    singleCellFeaturePointSet.erase(f[fp]);
+                    multiCellFeaturePointSet.insert(f[fp]);
+
+                    // Make sure there are points on the edges.
+                    featureEdgeSet.insert(fEdges[fp]);
+                }
+            }
+        }
+    }
+
+    // Transfer to arguments
+    featureFaces = featureFaceSet.toc();
+    featureEdges = featureEdgeSet.toc();
+    singleCellFeaturePoints = singleCellFeaturePointSet.toc();
+    multiCellFeaturePoints = multiCellFeaturePointSet.toc();
+}
+
+
+// Dump features to .obj files
+void dumpFeatures
+(
+    const polyMesh& mesh,
+    const labelList& featureFaces,
+    const labelList& featureEdges,
+    const labelList& singleCellFeaturePoints,
+    const labelList& multiCellFeaturePoints
+)
+{
+    {
+        OFstream str("featureFaces.obj");
+        Info<< "Dumping centres of featureFaces to obj file " << str.name()
+            << endl;
+        forAll(featureFaces, i)
+        {
+            meshTools::writeOBJ(str, mesh.faceCentres()[featureFaces[i]]);
+        }
+    }
+    {
+        OFstream str("featureEdges.obj");
+        Info<< "Dumping featureEdges to obj file " << str.name() << endl;
+        label vertI = 0;
+
+        forAll(featureEdges, i)
+        {
+            const edge& e = mesh.edges()[featureEdges[i]];
+            meshTools::writeOBJ(str, mesh.points()[e[0]]);
+            vertI++;
+            meshTools::writeOBJ(str, mesh.points()[e[1]]);
+            vertI++;
+            str<< "l " << vertI-1 << ' ' << vertI << nl;
+        }
+    }
+    {
+        OFstream str("singleCellFeaturePoints.obj");
+        Info<< "Dumping featurePoints that become a single cell to obj file "
+            << str.name() << endl;
+        forAll(singleCellFeaturePoints, i)
+        {
+            meshTools::writeOBJ(str, mesh.points()[singleCellFeaturePoints[i]]);
+        }
+    }
+    {
+        OFstream str("multiCellFeaturePoints.obj");
+        Info<< "Dumping featurePoints that become multiple cells to obj file "
+            << str.name() << endl;
+        forAll(multiCellFeaturePoints, i)
+        {
+            meshTools::writeOBJ(str, mesh.points()[multiCellFeaturePoints[i]]);
+        }
+    }
+}
+
+
+int main(int argc, char *argv[])
+{
+    argList::noParallel();
+#   include "addTimeOptions.H"
+    argList::validArgs.append("feature angle[0-180]");
+    argList::validOptions.insert("splitAllFaces", "");
+    argList::validOptions.insert("doNotPreserveFaceZones", "");
+    argList::validOptions.insert("overwrite", "");
+
+#   include "setRootCase.H"
+#   include "createTime.H"
+    // Get times list
+    instantList Times = runTime.times();
+#   include "checkTimeOptions.H"
+    runTime.setTime(Times[startTime], startTime);
+
+#   include "createMesh.H"
+
+    // Mark boundary edges and points.
+    // (Note: in 1.4.2 we can use the built-in mesh point ordering
+    //  facility instead)
+    PackedList<1> isBoundaryEdge(mesh.nEdges());
+    for (label faceI = mesh.nInternalFaces(); faceI < mesh.nFaces(); faceI++)
+    {
+        const labelList& fEdges = mesh.faceEdges()[faceI];
+
+        forAll(fEdges, i)
+        {
+            isBoundaryEdge.set(fEdges[i], 1);
+        }
+    }
+
+    scalar featureAngle(readScalar(IStringStream(args.additionalArgs()[0])()));
+
+    scalar minCos = Foam::cos(featureAngle * mathematicalConstant::pi/180.0);
+
+    Info<< "Feature:" << featureAngle << endl
+        << "minCos :" << minCos << endl
+        << endl;
+
+
+    const bool splitAllFaces = args.options().found("splitAllFaces");
+    const bool overwrite = args.options().found("overwrite");
+    const bool doNotPreserveFaceZones = args.options().found
+    (
+        "doNotPreserveFaceZones"
+    );
+
+    // Face(centre)s that need inclusion in the dual mesh
+    labelList featureFaces;
+    // Edge(centre)s  ,,
+    labelList featureEdges;
+    // Points (that become a single cell) that need inclusion in the dual mesh
+    labelList singleCellFeaturePoints;
+    // Points (that become a mulitple cells)        ,,
+    labelList multiCellFeaturePoints;
+
+    // Sample implementation of feature detection.
+    simpleMarkFeatures
+    (
+        mesh,
+        isBoundaryEdge,
+        featureAngle,
+        doNotPreserveFaceZones,
+
+        featureFaces,
+        featureEdges,
+        singleCellFeaturePoints,
+        multiCellFeaturePoints
+    );
+
+    // If we want to split all polyMesh faces into one dualface per cell
+    // we are passing through we also need a point
+    // at the polyMesh facecentre and edgemid of the faces we want to
+    // split.
+    if (splitAllFaces)
+    {
+        featureEdges = identity(mesh.nEdges());
+        featureFaces = identity(mesh.nFaces());
+    }
+
+    // Write obj files for debugging
+    dumpFeatures
+    (
+        mesh,
+        featureFaces,
+        featureEdges,
+        singleCellFeaturePoints,
+        multiCellFeaturePoints
+    );
+
+
+
+    // Read objects in time directory
+    IOobjectList objects(mesh, runTime.timeName());
+
+    // Read vol fields.
+    PtrList<volScalarField> vsFlds;
+    ReadFields(mesh, objects, vsFlds);
+
+    PtrList<volVectorField> vvFlds;
+    ReadFields(mesh, objects, vvFlds);
+
+    PtrList<volSphericalTensorField> vstFlds;
+    ReadFields(mesh, objects, vstFlds);
+
+    PtrList<volSymmTensorField> vsymtFlds;
+    ReadFields(mesh, objects, vsymtFlds);
+
+    PtrList<volTensorField> vtFlds;
+    ReadFields(mesh, objects, vtFlds);
+
+    // Read surface fields.
+    PtrList<surfaceScalarField> ssFlds;
+    ReadFields(mesh, objects, ssFlds);
+
+    PtrList<surfaceVectorField> svFlds;
+    ReadFields(mesh, objects, svFlds);
+
+    PtrList<surfaceSphericalTensorField> sstFlds;
+    ReadFields(mesh, objects, sstFlds);
+
+    PtrList<surfaceSymmTensorField> ssymtFlds;
+    ReadFields(mesh, objects, ssymtFlds);
+
+    PtrList<surfaceTensorField> stFlds;
+    ReadFields(mesh, objects, stFlds);
+
+
+    // Topo change container
+    polyTopoChange meshMod(mesh.boundaryMesh().size());
+
+    // Mesh dualiser engine
+    meshDualiser dualMaker(mesh);
+
+    // Insert all commands into polyTopoChange to create dual of mesh. This does
+    // all the hard work.
+    dualMaker.setRefinement
+    (
+        splitAllFaces,
+        featureFaces,
+        featureEdges,
+        singleCellFeaturePoints,
+        multiCellFeaturePoints,
+        meshMod
+    );
+
+    // Create mesh, return map from old to new mesh.
+    autoPtr<mapPolyMesh> map = meshMod.changeMesh(mesh, false);
+
+    // Update fields
+    mesh.updateMesh(map);
+
+    // Optionally inflate mesh
+    if (map().hasMotionPoints())
+    {
+        mesh.movePoints(map().preMotionPoints());
+    }
+
+    if (!overwrite)
+    {
+        runTime++;
+        mesh.setInstance(runTime.timeName());
+    }
+
+    Info<< "Writing dual mesh to " << runTime.timeName() << endl;
+
+    mesh.write();
+ 
+    Info<< "End\n" << endl;
+
+    return 0;
+}
+
+
+// ************************************************************************* //

applications/utilities/mesh/conversion/polyDualMesh/meshDualiser.H

@@ -0,0 +1,262 @@
+/*---------------------------------------------------------------------------*\
+  =========                 |
+  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+   \\    /   O peration     |
+    \\  /    A nd           | Copyright (C) 1991-2007 OpenCFD Ltd.
+     \\/     M anipulation  |
+-------------------------------------------------------------------------------
+License
+    This file is part of OpenFOAM.
+
+    OpenFOAM is free software; you can redistribute it and/or modify it
+    under the terms of the GNU General Public License as published by the
+    Free Software Foundation; either version 2 of the License, or (at your
+    option) any later version.
+
+    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+    for more details.
+
+    You should have received a copy of the GNU General Public License
+    along with OpenFOAM; if not, write to the Free Software Foundation,
+    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+
+Class
+    meshDualiser
+
+Description
+    Creates dual of polyMesh. Every point becomes a cell (or multiple cells
+    for feature points), a walk around every edge creates faces between them.
+
+    Put all points you want in the final mesh into featurePoints; all edge(mid)s
+    you want in the final mesh into featureEdges; all face(centre)s in
+    faceFaces.
+
+    Usually to preserve boundaries:
+        - all boundary faces are featureFaces
+        - all edges and points inbetween different patches are
+          featureEdges/points.
+
+    In same way you can also preserve internal faces (e.g. faceZones)
+
+SourceFiles
+    meshDualiser.C
+
+\*---------------------------------------------------------------------------*/
+
+#ifndef meshDualiser_H
+#define meshDualiser_H
+
+#include "DynamicList.H"
+#include "PackedList.H"
+#include "boolList.H"
+#include "typeInfo.H"
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+namespace Foam
+{
+
+class polyMesh;
+class polyTopoChange;
+
+/*---------------------------------------------------------------------------*\
+                           Class meshDualiser Declaration
+\*---------------------------------------------------------------------------*/
+
+class meshDualiser
+{
+    // Private data
+
+        const polyMesh& mesh_;
+
+        //- From point on cell to dual cell. Either single entry or
+        //  one entry per pointCells
+        labelListList pointToDualCells_;
+
+        //- From point to dual point (or -1 if not feature point).
+        labelList pointToDualPoint_;
+
+        //- From cell to dual point. All cells become point
+        labelList cellToDualPoint_;
+
+        //- From face to dual point (or -1 if not feature face)
+        labelList faceToDualPoint_;
+
+        //- From edge to dual point (or -1 if not feature edge)
+        labelList edgeToDualPoint_;
+
+
+    // Private Member Functions
+
+        static void checkPolyTopoChange(const polyTopoChange&);
+
+        static void dumpPolyTopoChange(const polyTopoChange&, const fileName&);
+
+        //- Find dual cell given point and cell
+        label findDualCell(const label cellI, const label pointI) const;
+
+        //- Helper function to generate dualpoints on all boundary edges
+        //  emanating from (boundary & feature) point
+        void generateDualBoundaryEdges
+        (
+            const PackedList<1>&,
+            const label pointI,
+            polyTopoChange&
+        );
+
+        //- Check that owner and neighbour of face have same dual cell
+        bool sameDualCell
+        (
+            const label faceI,
+            const label pointI
+        ) const;
+
+        //- Add internal face
+        label addInternalFace
+        (
+            const label masterPointI,
+            const label masterEdgeI,
+            const label masterFaceI,
+
+            const bool edgeOrder,
+            const label dualCell0,
+            const label dualCell1,
+            const DynamicList<label>& verts,
+            polyTopoChange& meshMod
+        ) const;
+
+        //- Add boundary face
+        label addBoundaryFace
+        (
+            const label masterPointI,
+            const label masterEdgeI,
+            const label masterFaceI,
+
+            const label dualCellI,
+            const label patchI,
+            const DynamicList<label>& verts,
+            polyTopoChange& meshMod
+        ) const;
+
+        //- Create internal faces walking around edge
+        void createFacesAroundEdge
+        (
+            const bool splitFace,
+            const PackedList<1>&,
+            const label edgeI,
+            const label startFaceI,
+            polyTopoChange&,
+            boolList& doneEFaces
+        ) const;
+
+        //- Create single internal face from internal face
+        void createFaceFromInternalFace
+        (
+            const label faceI,
+            label& fp,
+            polyTopoChange&
+        ) const;
+
+        //- Creates boundary faces walking around point on patchI.
+        void createFacesAroundBoundaryPoint
+        (
+            const label patchI,
+            const label patchPointI,
+            const label startFaceI,
+            polyTopoChange&,
+            boolList& donePFaces            // pFaces visited
+        ) const;
+
+        //- Disallow default bitwise copy construct
+        meshDualiser(const meshDualiser&);
+
+        //- Disallow default bitwise assignment
+        void operator=(const meshDualiser&);
+
+
+public:
+
+    //- Runtime type information
+    ClassName("meshDualiser");
+
+
+    // Constructors
+
+        //- Construct from mesh
+        meshDualiser(const polyMesh&);
+
+
+    // Member Functions
+
+        // Access
+
+            //- From point on cell to dual cell. Either single entry or
+            //  one entry per pointCells.
+            const labelListList& pointToDualCells() const
+            {
+                return pointToDualCells_;
+            }
+
+            //- From point to dual point (or -1 if not feature point).
+            const labelList& pointToDualPoint() const
+            {
+                return pointToDualPoint_;
+            }
+
+            //- From cell to dual point (at cell centre). All cells become
+            //  points.
+            const labelList& cellToDualPoint() const
+            {
+                return cellToDualPoint_;
+            }
+
+            //- From face to dual point (at face centre; or -1 if not
+            //  feature face).
+            const labelList& faceToDualPoint() const
+            {
+                return faceToDualPoint_;
+            }
+
+            //- From edge to dual point (at edge mid; or -1 if not feature
+            //  edge).
+            const labelList& edgeToDualPoint() const
+            {
+                return edgeToDualPoint_;
+            }
+
+
+        // Edit
+
+
+            //- Insert all changes into meshMod to convert the polyMesh into
+            //  its dual.
+            //  featureFaces  : faces where we want a point at the face centre
+            //  featureEdges  : edges           ,,                 edge mid
+            //  featurePoints : points          ,,    point. Two variants:
+            //      singleCellFeaturePoints : point becomes one dualcell.
+            //          Use this for e.g. convex boundary points.
+            //      multiCellFeaturePoints : one dualcell per original cell
+            //          around point. Use this for e.g. concave boundary points
+            //          since it prevents big concave boundary cells.
+            void setRefinement
+            (
+                const bool splitFace,
+                const labelList& featureFaces,
+                const labelList& featureEdges,
+                const labelList& singleCellFeaturePoints,
+                const labelList& multiCellFeaturePoints,
+                polyTopoChange& meshMod
+            );
+};
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+} // End namespace Foam
+
+// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
+
+#endif
+
+// ************************************************************************* //

applications/utilities/mesh/conversion/polyDualMesh/polyDualMeshApp.C

@@ -1,79 +0,0 @@
-/*---------------------------------------------------------------------------*\
-  =========                 |
-  \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-   \\    /   O peration     |
-    \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
-     \\/     M anipulation  |
--------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software; you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by the
-    Free Software Foundation; either version 2 of the License, or (at your
-    option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM; if not, write to the Free Software Foundation,
-    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-
-Application
-    polyDualMesh
-
-Description
-    Calculate the dual of a polyMesh. Adheres to all the feature&patch edges
-
-\*---------------------------------------------------------------------------*/
-
-#include "argList.H"
-#include "Time.H"
-#include "polyDualMesh.H"
-#include "mathematicalConstants.H"
-
-using namespace Foam;
-
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-int main(int argc, char *argv[])
-{
-    argList::noParallel();
-    argList::validArgs.append("feature angle[0-180]");
-    argList::validOptions.insert("overwrite", "");
-
-#   include "setRootCase.H"
-#   include "createTime.H"
-    runTime.functionObjects().off();
-#   include "createPolyMesh.H"
-
-    scalar featureAngle(readScalar(IStringStream(args.additionalArgs()[0])()));
-    bool overwrite = args.options().found("overwrite");
-
-    scalar minCos = Foam::cos(featureAngle * mathematicalConstant::pi/180.0);
-
-    Info<< "Feature:" << featureAngle << endl
-        << "minCos :" << minCos << endl
-        << endl;
-
-    polyDualMesh dualMesh(mesh, minCos);
-
-    if (!overwrite)
-    {
-        runTime++;
-    }
-
-    Pout<< "Writing dualMesh to " << runTime.timeName() << endl;
-
-    dualMesh.write();
-
-    Info << "End\n" << endl;
-
-    return 0;
-}
-
-
-// ************************************************************************* //

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

@@ -109,7 +109,7 @@ Foam::label Foam::checkTopology
         {
             Info<< "    Number of regions: " << rs.nRegions() << " (OK)."
                 << endl;
-        
+
         }
         else
         {
@@ -214,7 +214,7 @@ Foam::label Foam::checkTopology
                 const pointField& pts = pp.points();
                 const labelList& mp = pp.meshPoints();
 
-                boundBox bb(vector::zero, vector::zero);
+                boundBox bb;   // zero-sized
                 if (returnReduce(mp.size(), sumOp<label>()) > 0)
                 {
                     bb.min() = pts[mp[0]];

applications/utilities/mesh/manipulation/mergeMeshes/mergePolyMesh.C

@@ -144,7 +144,7 @@ Foam::mergePolyMesh::mergePolyMesh(const IOobject& io)
     wordList curPointZoneNames = pointZones().names();
     if (curPointZoneNames.size() > 0)
     {
-        pointZoneNames_.setSize(2*curPointZoneNames.size());
+        pointZoneNames_.setCapacity(2*curPointZoneNames.size());
     }
 
     forAll (curPointZoneNames, zoneI)
@@ -157,7 +157,7 @@ Foam::mergePolyMesh::mergePolyMesh(const IOobject& io)
 
     if (curFaceZoneNames.size() > 0)
     {
-        faceZoneNames_.setSize(2*curFaceZoneNames.size());
+        faceZoneNames_.setCapacity(2*curFaceZoneNames.size());
     }
     forAll (curFaceZoneNames, zoneI)
     {
@@ -169,7 +169,7 @@ Foam::mergePolyMesh::mergePolyMesh(const IOobject& io)
 
     if (curCellZoneNames.size() > 0)
     {
-        cellZoneNames_.setSize(2*curCellZoneNames.size());
+        cellZoneNames_.setCapacity(2*curCellZoneNames.size());
     }
     forAll (curCellZoneNames, zoneI)
     {

applications/utilities/parallelProcessing/decomposePar/domainDecomposition.C

@@ -354,7 +354,7 @@ bool domainDecomposition::writeDecomposition()
             // Estimate size
             forAll(zonePoints, zoneI)
             {
-                zonePoints[zoneI].setSize(pz[zoneI].size() / nProcs_);
+                zonePoints[zoneI].setCapacity(pz[zoneI].size() / nProcs_);
             }
 
             // Use the pointToZone map to find out the single zone (if any),
@@ -423,8 +423,8 @@ bool domainDecomposition::writeDecomposition()
             {
                 label procSize = fz[zoneI].size() / nProcs_;
 
-                zoneFaces[zoneI].setSize(procSize);
+                zoneFaces[zoneI].setCapacity(procSize);
-                zoneFaceFlips[zoneI].setSize(procSize);
+                zoneFaceFlips[zoneI].setCapacity(procSize);
             }
 
             // Go through all the zoned faces and find out if they
@@ -514,7 +514,7 @@ bool domainDecomposition::writeDecomposition()
             // Estimate size
             forAll(zoneCells, zoneI)
             {
-                zoneCells[zoneI].setSize(cz[zoneI].size() / nProcs_);
+                zoneCells[zoneI].setCapacity(cz[zoneI].size() / nProcs_);
             }
 
             forAll (curCellLabels, celli)

applications/utilities/parallelProcessing/reconstructParMesh/reconstructParMesh.C

@@ -273,7 +273,7 @@ autoPtr<mapPolyMesh> mergeSharedPoints
         }
     }
 
-    return map;    
+    return map;
 }
 
 
@@ -418,11 +418,7 @@ int main(int argc, char *argv[])
     // Read point on individual processors to determine merge tolerance
     // (otherwise single cell domains might give problems)
 
-    boundBox bb
+    boundBox bb = boundBox::invertedBox;
-    (
-        point(GREAT, GREAT, GREAT),
-        point(-GREAT, -GREAT, -GREAT)
-    );
 
     for (label procI = 0; procI < nProcs; procI++)
     {

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

@@ -617,8 +617,10 @@ void Foam::vtkPV3Foam::addPatchNames(vtkRenderer* renderer)
         {
             const labelList& eFaces = edgeFaces[edgeI];
 
-            if (eFaces.size() != 2)
+            if (eFaces.size() == 1)
             {
+                // Note: could also do ones with > 2 faces but this gives
+                // too many zones for baffles
                 featEdge[edgeI] = true;
             }
             else if (mag(n[eFaces[0]] & n[eFaces[1]]) < 0.5)

applications/utilities/postProcessing/wall/wallHeatFlux/wallHeatFlux.C

@@ -64,7 +64,7 @@ int main(int argc, char *argv[])
         Info<< "\nWall heat fluxes [W]" << endl;
         forAll(patchHeatFlux, patchi)
         {
-            if (typeid(mesh.boundary()) == typeid(wallFvPatch))
+            if (typeid(mesh.boundary()[patchi]) == typeid(wallFvPatch))
             {
                 Info<< mesh.boundary()[patchi].name()
                     << " "

applications/utilities/preProcessing/mapFields/MapVolFields.H

@@ -29,7 +29,6 @@ License
 
 #include "GeometricField.H"
 #include "meshToMesh.H"
-#include "cuttingPlane.H"
 #include "IOobjectList.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

applications/utilities/surface/surfaceCoordinateSystemTransform/surfaceCoordinateSystemTransform.C

@@ -199,7 +199,6 @@ int main(int argc, char *argv[])
         if (args.options().found("clean"))
         {
             surf.cleanup(true);
-            surf.checkOrientation(true);
         }
 
         if (fromCsys.valid())

applications/utilities/surface/surfaceMeshConvert/surfaceMeshConvert.C

@@ -104,13 +104,19 @@ int main(int argc, char *argv[])
     {
         triSurface surf(importName);
 
+        Info<< "Read surface:" << endl;
+        surf.writeStats(Info);
+        Info<< endl;
+
         if (args.options().found("clean"))
         {
+            Info<< "Cleaning up surface" << endl;
             surf.cleanup(true);
-            surf.checkOrientation(true);
+            surf.writeStats(Info);
+            Info<< endl;
         }
 
-        Info << "writing " << exportName;
+        Info<< "writing " << exportName;
         if (scaleFactor <= 0)
         {
             Info<< " without scaling" << endl;
@@ -119,6 +125,8 @@ int main(int argc, char *argv[])
         {
             Info<< " with scaling " << scaleFactor << endl;
             surf.scalePoints(scaleFactor);
+            surf.writeStats(Info);
+            Info<< endl;
         }
 
         // write sorted by region
@@ -128,34 +136,16 @@ int main(int argc, char *argv[])
     {
         UnsortedMeshedSurface<face> surf(importName);
 
-        if (args.options().found("clean"))
+        Info<< "Read surface:" << endl;
-        {
+        surf.writeStats(Info);
-            surf.cleanup(true);
+        Info<< endl;
-            surf.checkOrientation(true);
-        }
-
-        Info << "writing " << exportName;
-        if (scaleFactor <= 0)
-        {
-            Info<< " without scaling" << endl;
-        }
-        else
-        {
-            Info<< " with scaling " << scaleFactor << endl;
-            surf.scalePoints(scaleFactor);
-        }
-
-        surf.write(exportName);
-    }
-#if 1
-    else if (args.options().found("triFace"))
-    {
-        MeshedSurface<triFace> surf(importName);
 
         if (args.options().found("clean"))
         {
+            Info<< "Cleaning up surface" << endl;
             surf.cleanup(true);
-            surf.checkOrientation(true);
+            surf.writeStats(Info);
+            Info<< endl;
         }
 
         Info<< "writing " << exportName;
@@ -167,6 +157,39 @@ int main(int argc, char *argv[])
         {
             Info<< " with scaling " << scaleFactor << endl;
             surf.scalePoints(scaleFactor);
+            surf.writeStats(Info);
+            Info<< endl;
+        }
+        surf.write(exportName);
+    }
+#if 1
+    else if (args.options().found("triFace"))
+    {
+        MeshedSurface<triFace> surf(importName);
+
+        Info<< "Read surface:" << endl;
+        surf.writeStats(Info);
+        Info<< endl;
+
+        if (args.options().found("clean"))
+        {
+            Info<< "Cleaning up surface" << endl;
+            surf.cleanup(true);
+            surf.writeStats(Info);
+            Info<< endl;
+        }
+
+        Info<< "writing " << exportName;
+        if (scaleFactor <= 0)
+        {
+            Info<< " without scaling" << endl;
+        }
+        else
+        {
+            Info<< " with scaling " << scaleFactor << endl;
+            surf.scalePoints(scaleFactor);
+            surf.writeStats(Info);
+            Info<< endl;
         }
         surf.write(exportName);
     }
@@ -175,10 +198,16 @@ int main(int argc, char *argv[])
     {
         MeshedSurface<face> surf(importName);
 
+        Info<< "Read surface:" << endl;
+        surf.writeStats(Info);
+        Info<< endl;
+
         if (args.options().found("clean"))
         {
+            Info<< "Cleaning up surface" << endl;
             surf.cleanup(true);
-            surf.checkOrientation(true);
+            surf.writeStats(Info);
+            Info<< endl;
         }
 
         Info<< "writing " << exportName;
@@ -190,6 +219,8 @@ int main(int argc, char *argv[])
         {
             Info<< " with scaling " << scaleFactor << endl;
             surf.scalePoints(scaleFactor);
+            surf.writeStats(Info);
+            Info<< endl;
         }
         surf.write(exportName);
     }

bin/doxyFilt

@@ -50,9 +50,9 @@ then
    */applications/solvers/*.C | */applications/utilities/*.C )
       awkScript=$WM_PROJECT_DIR/bin/tools/doxyFilt-top.awk
       ;;
-   */applications/solvers/*.H | */applications/utilities/*.H )
+#   */applications/solvers/*.H | */applications/utilities/*.H )
-      awkScript=$WM_PROJECT_DIR/bin/tools/doxyFilt-ignore.awk
+#      awkScript=$WM_PROJECT_DIR/bin/tools/doxyFilt-ignore.awk
-      ;;
+#      ;;
    esac
 
    awk -f $awkScript $1 | \

bin/foamCheckSourceDeps

@@ -1,70 +0,0 @@
-#!/bin/sh
-#------------------------------------------------------------------------------
-# =========                 |
-# \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
-#  \\    /   O peration     |
-#   \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
-#    \\/     M anipulation  |
-#-------------------------------------------------------------------------------
-# License
-#     This file is part of OpenFOAM.
-#
-#     OpenFOAM is free software; you can redistribute it and/or modify it
-#     under the terms of the GNU General Public License as published by the
-#     Free Software Foundation; either version 2 of the License, or (at your
-#     option) any later version.
-#
-#     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-#     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-#     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-#     for more details.
-#
-#     You should have received a copy of the GNU General Public License
-#     along with OpenFOAM; if not, write to the Free Software Foundation,
-#     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-#
-# Script
-#     foamCheckSourceDeps
-#
-# Description
-#     Usage: foamCheckSourceDeps [dir1 .. dirN]
-#
-#     Search for *.dep files that are without a corresponding .C or .L file.
-#     These could indicate a directory that has been moved.
-#         - print questionable directory  and dep file
-#------------------------------------------------------------------------------
-if [ "$1" = "-h" -o "$1" = "-help" ]
-then
-   cat <<USAGE 1>&2
-Usage: ${0##*/} [dir1 .. dirN]
-
-    Search for .dep files that are without a corresponding .C or .L file.
-    This could indicate a directory that has been moved.
-        - print questionable directory and file
-USAGE
-   exit 1
-fi
-
-if [ "$#" -eq 0 ]
-then
-   set -- .
-fi
-
-for checkDir
-do
-   if [ -d $checkDir ]
-   then
-      find $checkDir -name '*.dep' -print | while read depFile;
-      do
-         Cfile=$(echo $depFile | sed -e 's/\.dep$/.C/')
-         # also check flex files
-         Lfile=$(echo $depFile | sed -e 's/\.C$/.L/')
-         if [ ! -f $Cfile -a ! -f $Lfile ]
-         then
-            echo "$(dirname $Cfile)  ${depFile##*/}"
-         fi
-      done
-   fi
-done
-
-# -----------------------------------------------------------------------------

bin/rmdepold

@@ -0,0 +1,109 @@
+#!/bin/sh
+#------------------------------------------------------------------------------
+# =========                 |
+# \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+#  \\    /   O peration     |
+#   \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+#    \\/     M anipulation  |
+#-------------------------------------------------------------------------------
+# License
+#     This file is part of OpenFOAM.
+#
+#     OpenFOAM is free software; you can redistribute it and/or modify it
+#     under the terms of the GNU General Public License as published by the
+#     Free Software Foundation; either version 2 of the License, or (at your
+#     option) any later version.
+#
+#     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+#     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+#     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+#     for more details.
+#
+#     You should have received a copy of the GNU General Public License
+#     along with OpenFOAM; if not, write to the Free Software Foundation,
+#     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+#
+# Script
+#     rmdepold
+#
+# Description
+#     Usage: rmdepold [dir1 .. dirN]
+#
+#     Remove *.dep files that are without a corresponding .C or .L file.
+#     This often occurs when a directory has been moved.
+#         - print questionable directory and the *.dep file
+#         - optionally remove empty directories
+#------------------------------------------------------------------------------
+usage() {
+    while [ "$#" -ge 1 ]; do echo "$1"; shift; done
+    cat<<USAGE 1>&2
+Usage: ${0##*/} [OPTION] [dir1 .. dirN]
+options:
+  -rmdir      find and remove empty directories (recursively)
+
+Remove *.dep files that are without a corresponding .C or .L file.
+This often occurs when a directory has been moved.
+  - print questionable directory and file
+  - optionally remove empty directories
+
+USAGE
+    exit 1
+}
+
+unset optRmdir
+
+# parse options
+while [ "$#" -gt 0 ]
+do
+    case "$1" in
+    -h | -help)
+        usage
+        ;;
+    -rmdir)
+        optRmdir=true
+        shift
+        ;;
+    -*)
+       usage "unknown option: '$*'"
+       ;;
+    *)
+       break
+       ;;
+    esac
+done
+
+# default is the current directory
+[ "$#" -gt 0 ] || set -- .
+
+for checkDir
+do
+    if [ -d $checkDir ]
+    then
+        echo "searching: $checkDir"
+    else
+        echo "skipping non-dir: $checkDir"
+        continue
+    fi
+
+    find $checkDir -name '*.dep' -print | while read depFile
+    do
+        # check C++ and Flex files
+        if [ ! -r "${depFile%dep}C" -a ! -r "${depFile%dep}L" ];
+        then
+            echo "rm $depFile"
+            rm -f $depFile 2>/dev/null
+        fi
+    done
+
+    # remove empty dirs
+    if [ "$optRmdir" ]
+    then
+        # get subdirs ourselves so we can avoid particular directories
+        for dir in $(find $checkDir -mindepth 1 -maxdepth 1 -type d \( -name .git -prune -o -print \) )
+        do
+            echo "check dir: $dir"
+            find $dir -depth -type d -empty -exec rmdir {} \; -print
+        done
+    fi
+done
+# -----------------------------------------------------------------------------

doc/Doxygen/Doxyfile

@@ -538,6 +538,7 @@ EXCLUDE_SYMBOLS        =
 
 EXAMPLE_PATH           =
 
+
 # If the value of the EXAMPLE_PATH tag contains directories, you can use the
 # EXAMPLE_PATTERNS tag to specify one or more wildcard pattern (like *.cpp
 # and *.h) to filter out the source-files in the directories. If left
@@ -787,6 +788,10 @@ TREEVIEW_WIDTH         = 250
 # configuration options related to the LaTeX output
 #---------------------------------------------------------------------------
 
+# Path for OpenCFD LaTeX macros
+
+@INCLUDE_PATH = $(WM_PROJECT_DIR)/doc/Doxygen/Macros/
+
 # If the GENERATE_LATEX tag is set to YES (the default) Doxygen will
 # generate Latex output.
 
@@ -824,7 +829,7 @@ PAPER_TYPE             = a4wide
 # The EXTRA_PACKAGES tag can be to specify one or more names of LaTeX
 # packages that should be included in the LaTeX output.
 
-EXTRA_PACKAGES         =
+EXTRA_PACKAGES         = $(WM_PROJECT_DIR)/doc/Doxygen/Macros/tensorOperator
 
 # The LATEX_HEADER tag can be used to specify a personal LaTeX header for
 # the generated latex document. The header should contain everything until

doc/Doxygen/Macros/tensorOperator.sty

@@ -0,0 +1,129 @@
+%-----------------------------------------------------------------------------
+% =========                 |
+% \\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox
+%  \\    /   O peration     |
+%   \\  /    A nd           | Copyright (C) 1991-2008 OpenCFD Ltd.
+%    \\/     M anipulation  |
+%------------------------------------------------------------------------------
+% License
+%     This file is part of OpenFOAM.
+%
+%     OpenFOAM is free software; you can redistribute it and/or modify it
+%     under the terms of the GNU General Public License as published by the
+%     Free Software Foundation; either version 2 of the License, or (at your
+%     option) any later version.
+%
+%     OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
+%     ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
+%     FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
+%     for more details.
+%
+%     You should have received a copy of the GNU General Public License
+%     along with OpenFOAM; if not, write to the Free Software Foundation,
+%     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
+%
+% LaTeX Style file
+%     tensorOperator.sty
+%
+% Description
+%     Standard OpenCFD LaTeX macros for typesetting tensor algebra.
+%
+%------------------------------------------------------------------------------
+
+% tensor style
+% ~~~~~~~~~~~~
+\renewcommand{\vec}[1] {\ensuremath{\mathbf #1}}
+\newcommand{\gvec}[1] {\ensuremath{\mbox{\boldmath$\bf#1$}}}
+
+% products
+% ~~~~~~~~
+\newcommand{\anyprod}{\star}
+\newcommand{\cprod} {\times}
+\newcommand{\dprod} {\,{\scriptscriptstyle \stackrel{\bullet}{{}}}\,}
+\newcommand{\ddprod} {\,{\scriptscriptstyle \stackrel{\bullet}{\bullet}}\,}
+\newcommand{\tdprod} {\,{\scriptscriptstyle \stackrel{3}{\bullet}}\,}
+\newcommand{\tprod} {\,{\scriptscriptstyle \stackrel{\otimes}{{}}}\,}
+
+% operations
+% ~~~~~~~~~~
+\newcommand{\adj} {\ensuremath{\operatorname{adj}}}
+\newcommand{\cof} {\ensuremath{\operatorname{cof}}}
+\newcommand{\diag} {\ensuremath{\operatorname{diag}}}
+\newcommand{\dev} {\ensuremath{\operatorname{dev}}}
+
+\newcommand{\Hodge} {\ensuremath{\operatorname{\stackrel{\displaystyle \ast}{}}}}
+\newcommand{\hyd} {\ensuremath{\operatorname{hyd}}}
+\renewcommand{\max} {\ensuremath{\operatorname{max}}}
+\renewcommand{\min} {\ensuremath{\operatorname{min}}}
+\newcommand{\inv} {\ensuremath{\operatorname{inv}}}
+\newcommand{\sym} {\ensuremath{\operatorname{symm}}}    % symm ?
+\newcommand{\skw} {\ensuremath{\operatorname{skew}}}    % skew already defined
+\newcommand{\tr} {\ensuremath{\operatorname{tr}}}
+\newcommand{\trans}[1] {\ensuremath{#1^{\operatorname{T}}}}
+
+% alternative tensor operators for hypersonics etc.
+% ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+\newcommand{\devs}[1] {\overset{\scriptscriptstyle\circ}{#1}}
+%\newcommand{\trans}[1] {\ensuremath{#1^{\operatorname{T}}}}
+\newcommand{\symms}[1] {\overleftrightarrow{#1}}
+\newlength{\skewslength}
+\newlength{\skewsheight}
+\newcommand{\skews}[1]{
+    \settowidth{\skewslength}{#1}%
+    \settoheight{\skewsheight}{#1}%
+    \addtolength{\skewsheight}{0.4mm}%
+    {\overleftrightarrow{#1}\hspace{-.5\skewslength}%
+    \rule[\skewsheight]{.4pt}{1.4mm}
+    \hspace{.5\skewslength}%
+}}
+%\newcommand{\skew}[1] {\ensuremath{#1^{\operatorname{A}}}}
+
+% spatial derivatives
+% ~~~~~~~~~~~~~~~~~~~
+\newcommand{\curl}{\ensuremath{\nabla\cprod}}
+\renewcommand{\div} {\ensuremath{\nabla\dprod}}
+\newcommand{\grad}{\ensuremath{\nabla}}
+\newcommand{\laplacian}{\ensuremath{\nabla^{2}}}
+
+% temporal derivatives
+% ~~~~~~~~~~~~~~~~~~~~
+\newcommand{\ddt}[1] {\ensuremath{\frac{\partial #1}{\partial t }}}
+\newcommand{\DDt}[1] {\ensuremath{\frac{D #1}{D t}}}
+\newcommand{\DpDt}[2] {\ensuremath{\frac{d_{#1} #2}{d t }}}
+\newcommand{\dsdts}[1] {\ensuremath{\frac{\partial ^2 #1}{\partial t^2}}}
+\newcommand{\rate}[1] {\ensuremath{\dot{#1}}}
+
+\newcommand{\genDer}{\mathcal{L}}
+
+% time average symbols
+% ~~~~~~~~~~~~~~~~~~~~
+\newcommand{\av}[1] {\ensuremath{\overline{#1}}}
+\newcommand{\corrtwo}[2] {{\dwea{\dprime{#1} \dprime{#2}}}}
+\newcommand{\curly}[1] {{\cal #1}}
+\newcommand{\dprime}[1] {\ensuremath{{#1}^{^{\prime \prime}}}}
+\newcommand{\dwea}[1] {\ensuremath{\widetilde{#1}}}
+\newcommand{\dweafluc}[1] {\ensuremath{\dprime{#1}}}
+\newcommand{\fluc}[1] {\ensuremath{#1^{\prime}}}
+
+% index style
+% ~~~~~~~~~~~
+\newcommand{\veci}[2][i] {\ensuremath{#2_{#1}}}
+\newcommand{\teni}[2][ij] {\ensuremath{#2_{#1}}}
+\newcommand{\tenTi}[2][ji] {\ensuremath{#2_{#1}}}
+
+% index operations
+% ~~~~~~~~~~~~~~~~
+\newcommand{\deltai}[1] {\ensuremath{\partial_{#1}}}
+
+% Sub-subscripts
+% ~~~~~~~~~~~~~~
+\newcommand{\eff} {{\scriptscriptstyle e\!f\!\!f\!}}
+
+% unknown use
+% ~~~~~~~~~~~
+%\font\bigtenrm=cmr12 scaled 1200
+%\newcommand{\eexp}[1]{{\hbox{$\textfont1=\bigtenrm e$}}^{\raise3pt
+%\hbox{$#1$}}}
+
+
+% ------------------------------------------------------------------------------

etc/settings.csh

@@ -187,7 +187,7 @@ case MPICH-GM:
     setenv FOAM_MPI_LIBBIN $FOAM_LIBBIN/mpich-gm
     breaksw
 
-case MPICH-GM:
+case HPMPI:
     setenv MPI_HOME /opt/hpmpi
     setenv MPI_ARCH_PATH $MPI_HOME
     setenv MPICH_ROOT=$MPI_ARCH_PATH

src/OSspecific/Unix/regularExpression.H

@@ -37,6 +37,8 @@ SourceFiles
 
 #include <sys/types.h>
 #include <regex.h>
+#include "string.H"
+#include "IOstreams.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -73,7 +75,7 @@ public:
         {
             preg_ = new regex_t;
 
-            if (regcomp(preg_, s.c_str(), REG_EXTENDED|REG_NOSUB) != 0)
+            if (regcomp(preg_, s.c_str(), REG_EXTENDED) != 0)
             {
                 FatalErrorIn
                 (
@@ -102,10 +104,12 @@ public:
         //- Matches?
         inline bool matches(const string& s) const
         {
-            size_t nmatch = 0;
+            size_t nmatch = 1;
-            regmatch_t *pmatch = NULL;
+            regmatch_t pmatch[1];
 
-            return regexec(preg_, s.c_str(), nmatch, pmatch, 0) == 0;
+            int errVal = regexec(preg_, s.c_str(), nmatch, pmatch, 0);
+
+            return (errVal == 0 && pmatch[0].rm_eo == label(s.size()));
         }
 };
 

src/OpenFOAM/containers/Lists/DynamicList/DynamicList.C

@@ -28,23 +28,23 @@ License
 
 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
 
-// Construct from Istream
+
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::DynamicList(Istream& is)
 :
     List<T>(is),
-    allocSize_(List<T>::size())
+    capacity_(List<T>::size())
 {}
 
 
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 Foam::Ostream& Foam::operator<<
 (
-    Foam::Ostream& os,
+    Ostream& os,
-    const Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>& DL
+    const DynamicList<T, SizeInc, SizeMult, SizeDiv>& lst
 )
 {
-    os << static_cast<const List<T>&>(DL);
+    os << static_cast<const List<T>&>(lst);
     return os;
 }
 
@@ -52,12 +52,12 @@ Foam::Ostream& Foam::operator<<
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 Foam::Istream& Foam::operator>>
 (
-    Foam::Istream& is,
+    Istream& is,
-    Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>& DL
+    DynamicList<T, SizeInc, SizeMult, SizeDiv>& lst
 )
 {
-    is >> static_cast<List<T>&>(DL);
+    is >> static_cast<List<T>&>(lst);
-    DL.allocSize_ = DL.List<T>::size();
+    lst.capacity_ = lst.List<T>::size();
 
     return is;
 }

src/OpenFOAM/containers/Lists/DynamicList/DynamicList.H

@@ -32,7 +32,7 @@ Description
     Internal storage is a compact array and the list can be shrunk to compact
     storage. The increase of list size is controlled by three template
     parameters, which allows the list storage to either increase by the given
-    increment or the given multiplier and divider (allowing non-integer
+    increment or by the given multiplier and divider (allowing non-integer
     multiples).
 
 SourceFiles
@@ -81,14 +81,11 @@ class DynamicList
 {
     // Private data
 
-        //- Allocated size for underlying List.
+        //- The capacity (allocated size) of the underlying list.
-        label allocSize_;
+        label capacity_;
 
     // Private Member Functions
 
-        // Disabled, since the usefulness and semantics are not quite clear
-        void setSize(const label, const T&);
-
 public:
 
     // Related types
@@ -116,20 +113,30 @@ public:
         // Access
 
         //- Size of the underlying storage.
-        inline label allocSize() const;
+        inline label capacity() const;
-
 
         // Edit
 
-        //- Alter the list size.
+        //- Alter the size of the underlying storage.
-        //  When the new size is greater than the addressed list size, the
+        //  The addressed size will be truncated if needed to fit, but will
-        //  allocated list sizes is adjusted and the
+        //  remain otherwise untouched.
-        //  addressed size does not change.
+        //  Use this or reserve() in combination with append().
-        //  Otherwise the addressed list size is just reduced and the
+        inline void setCapacity(const label);
-        //  allocated size does not change.
+
+        //- Alter the addressed list size.
+        //  New space will be allocated if required.
+        //  Use this to resize the list prior to using the operator[] for
+        //  setting values (as per List usage).
         inline void setSize(const label);
 
-        //- Clear the list, i.e. set the size to zero.
+        //- Alter the addressed list size and fill new space with a constant.
+        inline void setSize(const label, const T&);
+
+        //- Reserve allocation space for at least this size.
+        //  Never shrinks the allocated size, use setCapacity() for that.
+        inline void reserve(const label);
+
+        //- Clear the addressed list, i.e. set the size to zero.
         //  Allocated size does not change
         inline void clear();
 
@@ -152,6 +159,9 @@ public:
         //- Append an element at the end of the list
         inline void append(const T& e);
 
+        //- Append a List at the end of this list
+        inline void append(const UList<T>&);
+
         //- Remove and return the top element
         inline T remove();
 
@@ -162,7 +172,7 @@ public:
         inline void operator=(const T&);
 
         //- Assignment from List<T>. Also handles assignment from DynamicList.
-        inline void operator=(const List<T>&);
+        inline void operator=(const UList<T>&);
 
 
     // IOstream operators

src/OpenFOAM/containers/Lists/DynamicList/DynamicListI.H

@@ -30,7 +30,7 @@ template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::DynamicList()
 :
     List<T>(SizeInc),
-    allocSize_(SizeInc)
+    capacity_(SizeInc)
 {
     List<T>::size(0);
 }
@@ -39,11 +39,11 @@ inline Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::DynamicList()
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::DynamicList
 (
-    const label s
+    const label nElem
 )
 :
-    List<T>(s),
+    List<T>(nElem),
-    allocSize_(s)
+    capacity_(nElem)
 {
     List<T>::size(0);
 }
@@ -56,39 +56,101 @@ inline Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::DynamicList
 )
 :
     List<T>(lst),
-    allocSize_(lst.size())
+    capacity_(lst.size())
 {}
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
-inline Foam::label Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::allocSize()
+inline Foam::label Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::capacity()
 const
 {
-    return allocSize_;
+    return capacity_;
+}
+
+
+template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::setCapacity
+(
+    const label nElem
+)
+{
+    label nextFree = List<T>::size();
+    capacity_ = nElem;
+
+    if (nextFree > capacity_)
+    {
+        // truncate addressed sizes too
+        nextFree = capacity_;
+    }
+
+    List<T>::setSize(capacity_);
+    List<T>::size(nextFree);
+}
+
+
+template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::reserve
+(
+    const label nElem
+)
+{
+    // allocate more capacity?
+    if (nElem > capacity_)
+    {
+        capacity_ = max
+        (
+            nElem,
+            label(SizeInc + capacity_ * SizeMult / SizeDiv)
+        );
+
+        // adjust allocated size, leave addressed size untouched
+        label nextFree = List<T>::size();
+        List<T>::setSize(capacity_);
+        List<T>::size(nextFree);
+    }
 }
 
 
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::setSize
 (
-    const label s
+    const label nElem
+)
+{
+    // allocate more capacity?
+    if (nElem > capacity_)
+    {
+        capacity_ = max
+        (
+            nElem,
+            label(SizeInc + capacity_ * SizeMult / SizeDiv)
+        );
+
+        List<T>::setSize(capacity_);
+    }
+
+    // adjust addressed size
+    List<T>::size(nElem);
+}
+
+
+template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::setSize
+(
+    const label nElem,
+    const T& t
 )
 {
     label nextFree = List<T>::size();
-    if (s <= nextFree)
+    setSize(nElem);
+
+    // set new elements to constant value
+    while (nextFree < nElem)
     {
-        // adjust addressed size, leave allocated size untouched
+        this->operator[](nextFree++) = t;
-        nextFree = s;
     }
-    else
-    {
-        // adjust allocated size, leave addressed size untouched
-        allocSize_ = s;
-        List<T>::setSize(allocSize_);
-    }
-    List<T>::size(nextFree);
 }
 
 
@@ -103,7 +165,7 @@ template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::clearStorage()
 {
     List<T>::clear();
-    allocSize_ = 0;
+    capacity_ = 0;
 }
 
 
@@ -111,13 +173,15 @@ template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>&
 Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::shrink()
 {
-    if (allocSize_ > List<T>::size())
+    label nextFree = List<T>::size();
+    if (capacity_ > nextFree)
     {
-        allocSize_ = List<T>::size();
+        // use the full list when resizing
-        // force re-allocation/copying in List<T>::setSize() by temporarily
+        List<T>::size(capacity_);
-        // faking a larger list size that will be truncated
+        // the new size
-        List<T>::size(allocSize_+1);
+        capacity_ = nextFree;
-        List<T>::setSize(allocSize_);
+        List<T>::setSize(capacity_);
+        List<T>::size(nextFree);
     }
     return *this;
 }
@@ -127,7 +191,7 @@ template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline void
 Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::transfer(List<T>& lst)
 {
-    allocSize_ = lst.size();
+    capacity_ = lst.size();
     List<T>::transfer(lst);   // take over storage, clear addressing for lst.
 }
 
@@ -140,41 +204,58 @@ Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::transfer
 )
 {
     // take over storage as-is (without shrink), clear addressing for lst.
-    allocSize_ = lst.allocSize_;
+    capacity_ = lst.capacity_;
-    lst.allocSize_ = 0;
+    lst.capacity_ = 0;
 
     List<T>::transfer(static_cast<List<T>&>(lst));
 }
 
 
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
-inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::append(const T& e)
+inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::append
+(
+    const T& t
+)
+{
+    label elemI = List<T>::size();
+    setSize(elemI + 1);
+
+    this->operator[](elemI) = t;
+}
+
+
+template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::append
+(
+    const UList<T>& lst
+)
 {
-    // Work on copy free index since gets overwritten by setSize
     label nextFree = List<T>::size();
 
-    nextFree++;
+    if (this == &lst)
-
-    if (nextFree > allocSize_)
     {
-        allocSize_ = max
+        FatalErrorIn
         (
-            nextFree,
+            "DynamicList<T, SizeInc, SizeMult, SizeDiv>::append"
-            label(SizeMult*allocSize_/SizeDiv + SizeInc)
+            "(const UList<T>&)"
-        );
+        )   << "attempted appending to self" << abort(FatalError);
-        List<T>::setSize(allocSize_);
     }
 
-    List<T>::size(nextFree);
+    setSize(nextFree + lst.size());
 
-    this->operator[](nextFree - 1) = e;
+    forAll(lst, elemI)
+    {
+        this->operator[](nextFree++) = lst[elemI];
+    }
 }
 
 
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline T Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::remove()
 {
-    if (List<T>::size() == 0)
+    label elemI = List<T>::size() - 1;
+
+    if (elemI < 0)
     {
         FatalErrorIn
         (
@@ -182,11 +263,9 @@ inline T Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::remove()
         )   << "List is empty" << abort(FatalError);
     }
 
-    label nextFree = List<T>::size()-1;
+    const T& val = List<T>::operator[](elemI);
 
-    const T& val = List<T>::operator[](nextFree);
+    List<T>::size(elemI);
-
-    List<T>::size(nextFree);
 
     return val;
 }
@@ -197,26 +276,15 @@ inline T Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::remove()
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline T& Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::operator()
 (
-    const label i
+    const label elemI
 )
 {
-    label nextFree = List<T>::size();
+    if (elemI >= List<T>::size())
-
-    nextFree = max(nextFree, i + 1);
-
-    if (nextFree > allocSize_)
     {
-        allocSize_ = max
+        setSize(elemI + 1);
-        (
-            nextFree,
-            label(SizeMult*allocSize_/SizeDiv + SizeInc)
-        );
-        List<T>::setSize(allocSize_);
     }
 
-    List<T>::size(nextFree);
+    return this->operator[](elemI);
-
-    return this->operator[](i);
 }
 
 
@@ -226,21 +294,39 @@ inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::operator=
     const T& t
 )
 {
-    List<T>::operator=(t);
+    UList<T>::operator=(t);
 }
 
 
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
 inline void Foam::DynamicList<T, SizeInc, SizeMult, SizeDiv>::operator=
 (
-    const List<T>& lst
+    const UList<T>& lst
 )
 {
-    // make the entire storage available for the copy operation:
+    if (this == &lst)
-    List<T>::size(allocSize_);
+    {
+        FatalErrorIn
+        (
+            "DynamicList<T, SizeInc, SizeMult, SizeDiv>::operator="
+            "(const UList<T>&)"
+        )   << "attempted assignment to self" << abort(FatalError);
+    }
 
-    List<T>::operator=(lst);
+    if (capacity_ >= lst.size())
-    allocSize_ = List<T>::size();
+    {
+        // can copy w/o reallocating, match initial size to avoid reallocation
+        List<T>::size(lst.size());
+        List<T>::operator=(lst);
+    }
+    else
+    {
+        // make everything available for the copy operation
+        List<T>::size(capacity_);
+
+        List<T>::operator=(lst);
+        capacity_ = List<T>::size();
+    }
 }
 
 

src/OpenFOAM/containers/Lists/List/List.C

@@ -432,14 +432,19 @@ void Foam::List<T>::transfer(DynamicList<T, SizeInc, SizeMult, SizeDiv>& a)
 {
     // shrink the allocated space to the number of elements used
     a.shrink();
+    transfer(static_cast<List<T>&>(a));
+    a.clearStorage();
+}
 
-    if (this->v_) delete[] this->v_;
-    this->size_ = a.size_;
-    this->v_ = a.v_;
 
-    a.size_ = 0;
+// Transfer the contents of the argument SortableList into this List
-    a.v_ = 0;
+// and anull the argument list
-    a.allocSize_ = 0;
+template<class T>
+void Foam::List<T>::transfer(SortableList<T>& a)
+{
+    // shrink away the sort indices
+    a.shrink();
+    transfer(static_cast<List<T>&>(a));
 }
 
 

src/OpenFOAM/containers/Lists/List/List.H

@@ -63,7 +63,8 @@ template<class T, label Size> class FixedList;
 template<class T> class PtrList;
 template<class T> class SLList;
 template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
- class DynamicList;
+    class DynamicList;
+template<class T> class SortableList;
 template<class T> class IndirectList;
 template<class T> class BiIndirectList;
 
@@ -173,6 +174,10 @@ public:
             template<unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
             void transfer(DynamicList<T, SizeInc, SizeMult, SizeDiv>&);
 
+            //- Transfer the contents of the argument List into this List
+            //  and annull the argument list.
+            void transfer(SortableList<T>&);
+
             //- Return subscript-checked element of UList.
             inline T& newElmt(const label);
 

src/OpenFOAM/containers/Lists/ListOps/ListOps.H

@@ -82,6 +82,13 @@ void inplaceMapKey(const UList<label>& oldToNew, Container&);
 template<class T>
 void sortedOrder(const UList<T>&, labelList& order);
 
+//- Generate (sorted) indices corresponding to duplicate list values
+template<class T>
+void duplicateOrder(const UList<T>&, labelList& order);
+
+//- Generate (sorted) indices corresponding to unique list values
+template<class T>
+void uniqueOrder(const UList<T>&, labelList& order);
 
 //- Extract elements of List whose region is certain value.
 //  Use e.g. to extract all selected elements:

src/OpenFOAM/containers/Lists/ListOps/ListOpsTemplates.C

@@ -174,16 +174,65 @@ void Foam::sortedOrder
 )
 {
     order.setSize(lst.size());
-
     forAll(order, elemI)
     {
         order[elemI] = elemI;
     }
-
     Foam::stableSort(order, typename UList<T>::less(lst));
 }
 
 
+template<class T>
+void Foam::duplicateOrder
+(
+    const UList<T>& lst,
+    labelList& order
+)
+{
+    if (lst.size() < 2)
+    {
+        order.clear();
+        return;
+    }
+
+    sortedOrder(lst, order);
+
+    label n = 0;
+    for (label i = 0; i < order.size() - 1; ++i)
+    {
+        if (lst[order[i]] == lst[order[i+1]])
+        {
+            order[n++] = order[i];
+        }
+    }
+    order.setSize(n);
+}
+
+
+template<class T>
+void Foam::uniqueOrder
+(
+    const UList<T>& lst,
+    labelList& order
+)
+{
+    sortedOrder(lst, order);
+
+    if (order.size() > 1)
+    {
+        label n = 0;
+        for (label i = 0; i < order.size() - 1; ++i)
+        {
+            if (lst[order[i]] != lst[order[i+1]])
+            {
+                order[n++] = order[i];
+            }
+        }
+        order.setSize(n);
+    }
+}
+
+
 template<class T, class ListType>
 ListType Foam::subset
 (

src/OpenFOAM/containers/Lists/SortableList/SortableList.C

@@ -24,6 +24,23 @@ License
 
 \*---------------------------------------------------------------------------*/
 
+// * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
+
+template <class Type>
+void Foam::SortableList<Type>::sortIndices(List<label>& ind) const
+{
+    // list lengths must be identical
+    ind.setSize(this->size());
+
+    forAll(ind, i)
+    {
+        ind[i] = i;
+    }
+
+    Foam::stableSort(ind, typename UList<Type>::less(*this));
+}
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 template <class Type>
@@ -73,13 +90,6 @@ Foam::SortableList<Type>::SortableList(const SortableList<Type>& lst)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-template <class Type>
-void Foam::SortableList<Type>::setSize(const label newSize)
-{
-    List<Type>::setSize(newSize);
-    indices_.setSize(newSize, -1);
-}
-
 
 template <class Type>
 void Foam::SortableList<Type>::clear()
@@ -100,15 +110,7 @@ Foam::List<Type>& Foam::SortableList<Type>::shrink()
 template <class Type>
 void Foam::SortableList<Type>::sort()
 {
-    // list lengths must be identical
+    sortIndices(indices_);
-    indices_.setSize(this->size());
-
-    forAll(indices_, i)
-    {
-        indices_[i] = i;
-    }
-
-    Foam::stableSort(indices_, typename UList<Type>::less(*this));
 
     List<Type> lst(this->size());
     forAll(indices_, i)
@@ -120,10 +122,33 @@ void Foam::SortableList<Type>::sort()
 }
 
 
+template <class Type>
+void Foam::SortableList<Type>::reverseSort()
+{
+    sortIndices(indices_);
+
+    List<Type> lst(this->size());
+    label endI = indices_.size();
+    forAll(indices_, i)
+    {
+        lst[--endI] = this->operator[](indices_[i]);
+    }
+
+    List<Type>::transfer(lst);
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 template <class Type>
-void Foam::SortableList<Type>::operator=(const UList<Type>& rhs)
+inline void Foam::SortableList<Type>::operator=(const Type& t)
+{
+    UList<Type>::operator=(t);
+}
+
+
+template <class Type>
+inline void Foam::SortableList<Type>::operator=(const UList<Type>& rhs)
 {
     List<Type>::operator=(rhs);
     indices_.clear();
@@ -131,7 +156,7 @@ void Foam::SortableList<Type>::operator=(const UList<Type>& rhs)
 
 
 template <class Type>
-void Foam::SortableList<Type>::operator=(const SortableList<Type>& rhs)
+inline void Foam::SortableList<Type>::operator=(const SortableList<Type>& rhs)
 {
     List<Type>::operator=(rhs);
     indices_ = rhs.indices();

src/OpenFOAM/containers/Lists/SortableList/SortableList.H

@@ -60,6 +60,9 @@ class SortableList
         //- Original indices
         labelList indices_;
 
+        //- Resize and sort the parameter according to the list values
+        void sortIndices(List<label>&) const;
+
 public:
 
     // Constructors
@@ -99,9 +102,6 @@ public:
             return indices_;
         }
 
-        //- Size the list. Growing can cause undefined indices (until next sort)
-        void setSize(const label);
-
         //- Clear the list and the indices
         void clear();
 
@@ -112,13 +112,19 @@ public:
         //  also resizes the indices as required
         void sort();
 
+        //- Reverse (stable) sort the list
+        void reverseSort();
+
     // Member Operators
 
+        //- Assignment of all entries to the given value
+        inline void operator=(const Type&);
+
         //- Assignment from UList operator. Takes linear time.
-        void operator=(const UList<Type>&);
+        inline void operator=(const UList<Type>&);
 
         //- Assignment operator. Takes linear time.
-        void operator=(const SortableList<Type>&);
+        inline void operator=(const SortableList<Type>&);
 
 };
 

src/OpenFOAM/containers/Lists/SubList/SubList.H

@@ -87,6 +87,9 @@ public:
 
         //- Allow cast to a const List<T>&
         inline operator const Foam::List<T>&() const;
+
+        //- Assignment of all entries to the given value
+        inline void operator=(const T&);
 };
 
 

src/OpenFOAM/containers/Lists/SubList/SubListI.H

@@ -54,7 +54,7 @@ inline Foam::SubList<T>::SubList
 #   ifdef FULLDEBUG
 
     // Artificially allowing the start of a zero-sized subList to be
-    // one past the end of the original list.  
+    // one past the end of the original list.
     if (subSize > 0)
     {
         list.checkStart(startIndex);
@@ -73,7 +73,7 @@ inline Foam::SubList<T>::SubList
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class T>
-const Foam::SubList<T>& Foam::SubList<T>::null()
+inline const Foam::SubList<T>& Foam::SubList<T>::null()
 {
     SubList<T>* nullPtr = reinterpret_cast<SubList<T>*>(NULL);
     return *nullPtr;
@@ -89,4 +89,11 @@ inline Foam::SubList<T>::operator const Foam::List<T>&() const
 }
 
 
+template<class T>
+inline void Foam::SubList<T>::operator=(const T& t)
+{
+    UList<T>::operator=(t);
+}
+
+
 // ************************************************************************* //

src/OpenFOAM/containers/Lists/UList/UListI.H

@@ -50,14 +50,14 @@ inline Foam::UList<T>::UList(T* __restrict__ v, label size)
 template<class T>
 inline Foam::label Foam::UList<T>::fcIndex(const label i) const
 {
-    return (i == size()-1 ? 0 : i+1); 
+    return (i == size()-1 ? 0 : i+1);
 }
 
 
 template<class T>
 inline Foam::label Foam::UList<T>::rcIndex(const label i) const
 {
-    return (i == 0 ? size()-1 : i-1); 
+    return (i == 0 ? size()-1 : i-1);
 }
 
 

src/OpenFOAM/containers/Lists/UList/UListIO.C

@@ -45,7 +45,7 @@ void Foam::UList<T>::writeEntry(Ostream& os) const
     {
         os  << word("List<" + word(pTraits<T>::typeName) + '>') << " ";
     }
-    
+
     os << *this;
 }
 

src/OpenFOAM/db/Time/timeSelector.C

@@ -56,22 +56,12 @@ Foam::List<bool> Foam::timeSelector::selected(const List<instant>& Times) const
 {
     List<bool> lst(Times.size(), false);
 
-    // check ranges
+    // check ranges, avoid false positive on constant/
-    forAll(Times, i)
+    forAll(Times, timeI)
     {
-        if (selected(Times[i]))
+        if (Times[timeI].name() != "constant" && selected(Times[timeI]))
         {
-            lst[i] = true;
+            lst[timeI] = true;
-        }
-    }
-
-    // avoid false positive on "constant"
-    forAll(Times, i)
-    {
-        if (Times[i].name() == "constant")
-        {
-            lst[i] = false;
-            break;
         }
     }
 
@@ -85,15 +75,15 @@ Foam::List<bool> Foam::timeSelector::selected(const List<instant>& Times) const
             int nearestIndex = -1;
             scalar nearestDiff = Foam::GREAT;
 
-            forAll(Times, timeIndex)
+            forAll(Times, timeI)
             {
-                if (Times[timeIndex].name() == "constant") continue;
+                if (Times[timeI].name() == "constant") continue;
 
-                scalar diff = fabs(Times[timeIndex].value() - target);
+                scalar diff = fabs(Times[timeI].value() - target);
                 if (diff < nearestDiff)
                 {
                     nearestDiff = diff;
-                    nearestIndex = timeIndex;
+                    nearestIndex = timeI;
                 }
             }
 
@@ -156,6 +146,27 @@ Foam::List<Foam::instant> Foam::timeSelector::select
     {
         List<bool> selectTimes(timeDirs.size(), true);
 
+        // determine locations of constant/ and 0/ directories
+        label constantIdx = -1;
+        label zeroIdx = -1;
+
+        forAll(timeDirs, timeI)
+        {
+            if (timeDirs[timeI].name() == "constant")
+            {
+                constantIdx = timeI;
+            }
+            else if (timeDirs[timeI].value() == 0)
+            {
+                zeroIdx = timeI;
+            }
+
+            if (constantIdx >= 0 && zeroIdx >= 0)
+            {
+                break;
+            }
+        }
+
         if (args.options().found("time"))
         {
             selectTimes = timeSelector
@@ -166,25 +177,31 @@ Foam::List<Foam::instant> Foam::timeSelector::select
         else if (args.options().found("latestTime"))
         {
             selectTimes = false;
+            const label latestIdx = timeDirs.size() - 1;
 
             // avoid false match on constant/ or 0/
-            if (timeDirs.size() > 2)
+            if (latestIdx != constantIdx && latestIdx != zeroIdx)
             {
-                selectTimes[timeDirs.size() - 1] = true;
+                selectTimes[latestIdx] = true;
             }
         }
 
-        if (timeDirs.size() > 1)
+        // special treatment for constant/
+        if (constantIdx >= 0)
         {
-            selectTimes[0] = args.options().found("constant");
+            selectTimes[constantIdx] = args.options().found("constant");
+        }
 
+        // special treatment for 0/
+        if (zeroIdx >= 0)
+        {
             if (args.options().found("noZero"))
             {
-                selectTimes[1] = false;
+                selectTimes[zeroIdx] = false;
             }
             else if (argList::validOptions.found("zeroTime"))
             {
-                selectTimes[1] = args.options().found("zeroTime");
+                selectTimes[zeroIdx] = args.options().found("zeroTime");
             }
         }
 

src/OpenFOAM/db/dictionary/dictionary.C

@@ -27,6 +27,7 @@ License
 #include "dictionary.H"
 #include "primitiveEntry.H"
 #include "dictionaryEntry.H"
+#include "regularExpression.H"
 
 /* * * * * * * * * * * * * * * Static Member Data  * * * * * * * * * * * * * */
 

src/OpenFOAM/db/dictionary/dictionary.H

@@ -60,7 +60,6 @@ SourceFiles
 #include "HashTable.H"
 #include "wordList.H"
 #include "className.H"
-#include "regularExpression.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -68,7 +67,7 @@ namespace Foam
 {
 
 // Forward declaration of friend functions and operators
-
+class regularExpression;
 class dictionary;
 Istream& operator>>(Istream&, dictionary&);
 Ostream& operator<<(Ostream&, const dictionary&);

src/OpenFOAM/db/dictionary/dictionaryIO.C

@@ -27,6 +27,7 @@ License
 #include "dictionary.H"
 #include "IFstream.H"
 #include "inputModeEntry.H"
+#include "regularExpression.H"
 
 // * * * * * * * * * * * * * Private Member Functions * * * * * * * * * * * //
 

src/OpenFOAM/fields/Fields/transformField/transformField.C

@@ -26,6 +26,7 @@ License
 
 #include "transformField.H"
 #include "FieldM.H"
+#include "diagTensor.H"
 
 // * * * * * * * * * * * * * * * global functions  * * * * * * * * * * * * * //
 
@@ -75,7 +76,8 @@ void Foam::transform
 {
     vector T = tr.t();
 
-    if (mag(tr.r().w()) > SMALL)
+    // Check if any rotation
+    if (mag(tr.r().R() - I) > SMALL)
     {
         transform(rtf, tr.r(), tf);
 
@@ -90,6 +92,10 @@ void Foam::transform
         {
             TFOR_ALL_F_OP_S_OP_F(vector, rtf, =, vector, T, +, vector, tf);
         }
+        else
+        {
+            rtf = vector::zero;
+        }
     }
 }
 

src/OpenFOAM/meshes/boundBox/boundBox.C

@@ -27,18 +27,34 @@ License
 #include "boundBox.H"
 #include "PstreamReduceOps.H"
 
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+
+const Foam::boundBox Foam::boundBox::greatBox
+(
+    point(-VGREAT, -VGREAT, -VGREAT),
+    point(VGREAT, VGREAT, VGREAT)
+);
+
+
+const Foam::boundBox Foam::boundBox::invertedBox
+(
+    point(VGREAT, VGREAT, VGREAT),
+    point(-VGREAT, -VGREAT, -VGREAT)
+);
+
+
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 Foam::boundBox::boundBox(const pointField& points, const bool doReduce)
 :
-    min_(vector::zero),
+    min_(point::zero),
-    max_(vector::zero)
+    max_(point::zero)
 {
     if (points.size() == 0)
     {
         if (Pstream::parRun() && doReduce)
         {
-            // Use values which get overwritten by reduce minOp,maxOp below
+            // Use values that get overwritten by reduce minOp, maxOp below
             min_ = point(VGREAT, VGREAT, VGREAT);
             max_ = point(-VGREAT, -VGREAT, -VGREAT);
         }
@@ -89,7 +105,6 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const boundBox& bb)
 
     // Check state of Ostream
     os.check("Ostream& operator<<(Ostream&, const boundBox&)");
-
     return os;
 }
 
@@ -111,7 +126,6 @@ Foam::Istream& Foam::operator>>(Istream& is, boundBox& bb)
 
     // Check state of Istream
     is.check("Istream& operator>>(Istream&, boundBox&)");
-
     return is;
 }
 

src/OpenFOAM/meshes/boundBox/boundBox.H

@@ -48,7 +48,7 @@ Ostream& operator<<(Ostream& os, const boundBox& b);
 
 
 /*---------------------------------------------------------------------------*\
-                           Class boundBox Declaration
+                         Class boundBox Declaration
 \*---------------------------------------------------------------------------*/
 
 class boundBox
@@ -61,13 +61,22 @@ class boundBox
 
 public:
 
+    // Static data members
+
+        //- A very large boundBox: min/max == -/+ VGREAT
+        static const boundBox greatBox;
+
+        //- A very large inverted boundBox: min/max == +/- VGREAT
+        static const boundBox invertedBox;
+
+
     // Constructors
 
-        //- Construct null setting points to zero
+        //- Construct null, setting points to zero
         boundBox()
         :
-            min_(vector::zero),
+            min_(point::zero),
-            max_(vector::zero)
+            max_(point::zero)
         {}
 
         //- Construct from components
@@ -77,8 +86,8 @@ public:
             max_(max)
         {}
 
-        //- Construct as the bounding box of the given pointField. Does
+        //- Construct as the bounding box of the given pointField.
-        //  parallel communication (doReduce = true)
+        //  Does parallel communication (doReduce = true)
         boundBox(const pointField& points, const bool doReduce = true);
 
         //- Construct from Istream
@@ -117,39 +126,26 @@ public:
 
         // Query
 
-            //- Intersects other boundingbox?
+            //- Intersects other boundingBox?
             bool overlaps(const boundBox& bb) const
             {
-                if
+                return
                 (
-                    (min_.x() <= bb.max().x()) &&
+                    min_.x() <= bb.max().x() && max_.x() >= bb.min().x()
-                    (min_.y() <= bb.max().y()) &&
+                 && min_.y() <= bb.max().y() && max_.y() >= bb.min().y()
-                    (min_.z() <= bb.max().z()) &&
+                 && min_.z() <= bb.max().z() && max_.z() >= bb.min().z()
-
+                );
-                    (max_.x() >= bb.min().x()) &&
-                    (max_.y() >= bb.min().y()) &&
-                    (max_.z() >= bb.min().z())
-                )
-                {
-                    return true;
-                }
-                else
-                {
-                    return false;
-                }
             }
 
-
             //- Contains a point?
             bool contains(const point& pt) const
             {
                 return
-                    pt.x() >= min().x()
+                (
-                 && pt.y() >= min().y()
+                    pt.x() >= min().x() && pt.x() <= max().x()
-                 && pt.z() >= min().z()
+                 && pt.y() >= min().y() && pt.y() <= max().y()
-                 && pt.x() <= max().x()
+                 && pt.z() >= min().z() && pt.z() <= max().z()
-                 && pt.y() <= max().y()
+                );
-                 && pt.z() <= max().z();
             }
 
 
@@ -173,7 +169,7 @@ public:
 };
 
 
-//- Specify data associated with boundBox type is contiguous
+//- Specify data associated with boundBox type are contiguous
 template<>
 inline bool contiguous<boundBox>() {return contiguous<point>();}
 

src/OpenFOAM/meshes/meshShapes/face/face.C

@@ -119,7 +119,7 @@ Foam::label Foam::face::mostConcaveAngle
 }
 
 
-void Foam::face::split
+Foam::label Foam::face::split
 (
     const face::splitMode mode,
     const pointField& points,
@@ -129,6 +129,8 @@ void Foam::face::split
     faceList& quadFaces
 ) const
 {
+    label oldIndices = (triI + quadI);
+
     if (size() <= 2)
     {
         FatalErrorIn
@@ -143,7 +145,7 @@ void Foam::face::split
     if (size() == 3)
     {
         // Triangle. Just copy.
-        if ((mode == COUNTQUAD) || (mode == COUNTTRIANGLE))
+        if (mode == COUNTTRIANGLE || mode == COUNTQUAD)
         {
             triI++;
         }
@@ -250,7 +252,7 @@ void Foam::face::split
         }
         else
         {
-            // folded round
+            // folded around
             diff = minIndex + size() - startIndex;
         }
 
@@ -281,6 +283,8 @@ void Foam::face::split
         face1.split(mode, points, triI, quadI, triFaces, quadFaces);
         face2.split(mode, points, triI, quadI, triFaces, quadFaces);
     }
+
+    return (triI + quadI - oldIndices);
 }
 
 
@@ -749,30 +753,27 @@ int Foam::face::edgeDirection(const edge& e) const
 
 
 // Number of triangles directly known from number of vertices
-Foam::label Foam::face::nTriangles
+Foam::label Foam::face::nTriangles(const pointField&) const
-(
-    const pointField&
-) const
 {
-    return size() - 2;
+    return nTriangles();
 }
 
 
-void Foam::face::triangles
+Foam::label Foam::face::triangles
 (
     const pointField& points,
     label& triI,
     faceList& triFaces
 ) const
 {
-    faceList quadFaces;
     label quadI = 0;
+    faceList quadFaces;
 
-    split(SPLITTRIANGLE, points, triI, quadI, triFaces, quadFaces);
+    return split(SPLITTRIANGLE, points, triI, quadI, triFaces, quadFaces);
 }
 
 
-void Foam::face::nTrianglesQuads
+Foam::label Foam::face::nTrianglesQuads
 (
     const pointField& points,
     label& triI,
@@ -782,11 +783,11 @@ void Foam::face::nTrianglesQuads
     faceList triFaces;
     faceList quadFaces;
 
-    split(COUNTQUAD, points, triI, quadI, triFaces, quadFaces);
+    return split(COUNTQUAD, points, triI, quadI, triFaces, quadFaces);
 }
 
 
-void Foam::face::trianglesQuads
+Foam::label Foam::face::trianglesQuads
 (
     const pointField& points,
     label& triI,
@@ -795,7 +796,7 @@ void Foam::face::trianglesQuads
     faceList& quadFaces
 ) const
 {
-    split(SPLITQUAD, points, triI, quadI, triFaces, quadFaces);
+    return split(SPLITQUAD, points, triI, quadI, triFaces, quadFaces);
 }
 
 

src/OpenFOAM/meshes/meshShapes/face/face.H

@@ -57,11 +57,14 @@ namespace Foam
 
 class face;
 class triFace;
+
+template<class T, unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+class DynamicList;
+
 inline bool operator==(const face& a, const face& b);
 inline bool operator!=(const face& a, const face& b);
 inline Istream& operator>>(Istream&, face&);
 
-
 /*---------------------------------------------------------------------------*\
                            Class face Declaration
 \*---------------------------------------------------------------------------*/
@@ -108,9 +111,10 @@ class face
             SPLITQUAD       // split into triangles&quads
         };
 
-        //- Split face into triangles or triangles&quads. Stores results
+        //- Split face into triangles or triangles&quads.
-        //  quadFaces[quadI], triFaces[triI]
+        //  Stores results quadFaces[quadI], triFaces[triI]
-        void split
+        //  Returns number of new faces created
+        label split
         (
             const splitMode mode,
             const pointField& points,
@@ -270,30 +274,47 @@ public:
 
         // Face splitting utilities
 
+            //- Number of triangles after splitting
+            inline label nTriangles() const;
+
             //- Number of triangles after splitting
             label nTriangles(const pointField& points) const;
 
-            //- Split into triangles using existing points. Result in
+            //- Split into triangles using existing points.
-            //  triFaces[triI..triI+nTri]. Splits intelligently to maximize
+            //  Result in triFaces[triI..triI+nTri].
-            //  triangle quality.
+            //  Splits intelligently to maximize triangle quality.
-            void triangles
+            //  Returns number of faces created.
+            label triangles
             (
                 const pointField& points,
                 label& triI,
                 faceList& triFaces
             ) const;
 
+            //- Split into triangles using existing points.
+            //  Append to DynamicList.
+            //  Returns number of faces created.
+            template<unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+            label triangles
+            (
+                const pointField& points,
+                DynamicList<face, SizeInc, SizeMult, SizeDiv>& triFaces
+            ) const;
+
             //- Number of triangles and quads after splitting
-            void nTrianglesQuads
+            //  Returns the sum of both
+            label nTrianglesQuads
             (
                 const pointField& points,
                 label& nTris,
                 label& nQuads
             ) const;
 
-            //- Split into triangles and quads. Result in triFaces (starting at
+            //- Split into triangles and quads.
-            //  triI and quadFaces (starting at quadI)
+            //  Results in triFaces (starting at triI) and quadFaces
-            void trianglesQuads
+            //  (starting at quadI).
+            //  Returns number of new faces created.
+            label trianglesQuads
             (
                 const pointField& points,
                 label& triI,

src/OpenFOAM/meshes/meshShapes/face/faceI.H

@@ -127,6 +127,11 @@ inline Foam::label Foam::face::prevLabel(const label i) const
     return operator[](rcIndex(i));
 }
 
+// Number of triangles directly known from number of vertices
+inline Foam::label Foam::face::nTriangles() const
+{
+    return size() - 2;
+}
 
 // * * * * * * * * * * * * * * * Friend Operators  * * * * * * * * * * * * * //
 

src/OpenFOAM/meshes/meshShapes/face/faceTemplates.C

@@ -25,9 +25,28 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "face.H"
+#include "DynamicList.H"
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<unsigned SizeInc, unsigned SizeMult, unsigned SizeDiv>
+Foam::label Foam::face::triangles
+(
+    const pointField& points,
+    DynamicList<face, SizeInc, SizeMult, SizeDiv>& triFaces
+) const
+{
+    label triI = triFaces.size();
+    label quadI = 0;
+    faceList quadFaces;
+
+    // adjust the addressable size (and allocate space if needed)
+    triFaces.setSize(triI + nTriangles());
+
+    return split(SPLITTRIANGLE, points, triI, quadI, triFaces, quadFaces);
+}
+
+
 template<class Type>
 Type Foam::face::average
 (

src/OpenFOAM/meshes/meshShapes/triFace/triFace.H

@@ -124,6 +124,12 @@ public:
             //- Return vector normal
             inline vector normal(const pointField&) const;
 
+            //- Number of triangles after splitting
+            inline label nTriangles() const;
+
+            //- Return face with reverse direction
+            inline triFace reverseFace() const;
+
             //- Return swept-volume
             inline scalar sweptVol
             (

src/OpenFOAM/meshes/meshShapes/triFace/triFaceI.H

@@ -78,9 +78,9 @@ inline Foam::triFace::triFace
 }
 
 
-inline Foam::triFace::triFace(const UList<label>& l)
+inline Foam::triFace::triFace(const UList<label>& lst)
 :
-    FixedList<label, 3>(l)
+    FixedList<label, 3>(lst)
 {}
 
 
@@ -96,7 +96,7 @@ inline Foam::label Foam::triFace::collapse()
 {
     // we cannot resize a FixedList, so mark duplicates with '-1'
     // (the lower vertex is retained)
-    // catch any '-1' - ie, if called twice
+    // catch any '-1' (eg, if called twice)
 
     label n = 3;
     if (operator[](0) == operator[](1) || operator[](1) == -1)
@@ -154,13 +154,13 @@ inline Foam::edgeList Foam::triFace::edges() const
     edgeList e(3);
 
     e[0].start() = operator[](0);
-    e[0].end() = operator[](1);
+    e[0].end()   = operator[](1);
 
     e[1].start() = operator[](1);
-    e[1].end() = operator[](2);
+    e[1].end()   = operator[](2);
 
     e[2].start() = operator[](2);
-    e[2].end() = operator[](0);
+    e[2].end()   = operator[](0);
 
     return e;
 }
@@ -213,7 +213,7 @@ inline Foam::scalar Foam::triFace::mag(const pointField& points) const
     return ::Foam::mag(normal(points));
 }
 
-
+// could also delegate to triPointRef(...).normal()
 inline Foam::vector Foam::triFace::normal(const pointField& points) const
 {
     return 0.5*
@@ -224,6 +224,19 @@ inline Foam::vector Foam::triFace::normal(const pointField& points) const
 }
 
 
+inline Foam::label Foam::triFace::nTriangles() const
+{
+    return 1;
+}
+
+
+inline Foam::triFace Foam::triFace::reverseFace() const
+{
+    // The starting points of the original and reverse face are identical.
+    return triFace(operator[](0), operator[](2), operator[](1));
+}
+
+
 inline Foam::scalar Foam::triFace::sweptVol
 (
     const pointField& opts,

src/OpenFOAM/meshes/primitiveMesh/primitiveMeshCellPoints.C

@@ -27,14 +27,9 @@ License
 #include "primitiveMesh.H"
 #include "ListOps.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-const labelListList& primitiveMesh::cellPoints() const
+const Foam::labelListList& Foam::primitiveMesh::cellPoints() const
 {
     if (!cpPtr_)
     {
@@ -56,12 +51,12 @@ const labelListList& primitiveMesh::cellPoints() const
         cpPtr_ = new labelListList(nCells());
         invertManyToMany(nCells(), pointCells(), *cpPtr_);
     }
-    
+
     return *cpPtr_;
 }
 
 
-const labelList& primitiveMesh::cellPoints
+const Foam::labelList& Foam::primitiveMesh::cellPoints
 (
     const label cellI,
     DynamicList<label>& storage
@@ -84,14 +79,14 @@ const labelList& primitiveMesh::cellPoints
 
             forAll(f, fp)
             {
-                labelSet_.insert(f[fp]);    
+                labelSet_.insert(f[fp]);
             }
         }
 
         storage.clear();
-        if (labelSet_.size() > storage.allocSize())
+        if (labelSet_.size() > storage.capacity())
         {
-            storage.setSize(labelSet_.size());
+            storage.setCapacity(labelSet_.size());
         }
 
         forAllConstIter(labelHashSet, labelSet_, iter)
@@ -104,7 +99,7 @@ const labelList& primitiveMesh::cellPoints
 }
 
 
-const labelList& primitiveMesh::cellPoints(const label cellI) const
+const Foam::labelList& Foam::primitiveMesh::cellPoints(const label cellI) const
 {
     return cellPoints(cellI, labels_);
 }
@@ -112,6 +107,4 @@ const labelList& primitiveMesh::cellPoints(const label cellI) const
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace Foam
-
 // ************************************************************************* //

src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C

@@ -30,15 +30,10 @@ License
 #include "SortableList.H"
 #include "ListOps.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 // Returns edgeI between two points.
-Foam::label primitiveMesh::getEdge
+Foam::label Foam::primitiveMesh::getEdge
 (
     List<DynamicList<label> >& pe,
     DynamicList<edge>& es,
@@ -76,7 +71,7 @@ Foam::label primitiveMesh::getEdge
 }
 
 
-void primitiveMesh::calcEdges(const bool doFaceEdges) const
+void Foam::primitiveMesh::calcEdges(const bool doFaceEdges) const
 {
     if (debug)
     {
@@ -113,7 +108,7 @@ void primitiveMesh::calcEdges(const bool doFaceEdges) const
         List<DynamicList<label> > pe(nPoints());
         forAll(pe, pointI)
         {
-            pe[pointI].setSize(primitiveMesh::edgesPerPoint_);
+            pe[pointI].setCapacity(primitiveMesh::edgesPerPoint_);
         }
 
         // Estimate edges storage
@@ -337,7 +332,7 @@ void primitiveMesh::calcEdges(const bool doFaceEdges) const
 
                         oldToNew[edgeI] = internal0EdgeI++;
                     }
-                }        
+                }
             }
             else
             {
@@ -460,8 +455,8 @@ void primitiveMesh::calcEdges(const bool doFaceEdges) const
 }
 
 
-label primitiveMesh::findFirstCommonElementFromSortedLists
+Foam::label Foam::primitiveMesh::findFirstCommonElementFromSortedLists
-( 
+(
     const labelList& list1,
     const labelList& list2
 )
@@ -494,7 +489,7 @@ label primitiveMesh::findFirstCommonElementFromSortedLists
             "primitiveMesh::findFirstCommonElementFromSortedLists"
             "(const labelList&, const labelList&)"
         )   << "No common elements in lists " << list1 << " and " << list2
-            << abort(FatalError);        
+            << abort(FatalError);
     }
     return result;
 }
@@ -502,7 +497,7 @@ label primitiveMesh::findFirstCommonElementFromSortedLists
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-const edgeList& primitiveMesh::edges() const
+const Foam::edgeList& Foam::primitiveMesh::edges() const
 {
     if (!edgesPtr_)
     {
@@ -513,7 +508,7 @@ const edgeList& primitiveMesh::edges() const
     return *edgesPtr_;
 }
 
-const labelListList& primitiveMesh::pointEdges() const
+const Foam::labelListList& Foam::primitiveMesh::pointEdges() const
 {
     if (!pePtr_)
     {
@@ -525,7 +520,7 @@ const labelListList& primitiveMesh::pointEdges() const
 }
 
 
-const labelListList& primitiveMesh::faceEdges() const
+const Foam::labelListList& Foam::primitiveMesh::faceEdges() const
 {
     if (!fePtr_)
     {
@@ -576,7 +571,7 @@ const labelListList& primitiveMesh::faceEdges() const
 }
 
 
-void primitiveMesh::clearOutEdges()
+void Foam::primitiveMesh::clearOutEdges()
 {
     deleteDemandDrivenData(edgesPtr_);
     deleteDemandDrivenData(pePtr_);
@@ -586,7 +581,7 @@ void primitiveMesh::clearOutEdges()
 }
 
 
-const labelList& primitiveMesh::faceEdges
+const Foam::labelList& Foam::primitiveMesh::faceEdges
 (
     const label faceI,
     DynamicList<label>& storage
@@ -602,9 +597,9 @@ const labelList& primitiveMesh::faceEdges
         const face& f = faces()[faceI];
 
         storage.clear();
-        if (f.size() > storage.allocSize())
+        if (f.size() > storage.capacity())
         {
-            storage.setSize(f.size());
+            storage.setCapacity(f.size());
         }
 
         forAll(f, fp)
@@ -624,13 +619,13 @@ const labelList& primitiveMesh::faceEdges
 }
 
 
-const labelList& primitiveMesh::faceEdges(const label faceI) const
+const Foam::labelList& Foam::primitiveMesh::faceEdges(const label faceI) const
 {
     return faceEdges(faceI, labels_);
 }
 
 
-const labelList& primitiveMesh::cellEdges
+const Foam::labelList& Foam::primitiveMesh::cellEdges
 (
     const label cellI,
     DynamicList<label>& storage
@@ -652,15 +647,15 @@ const labelList& primitiveMesh::cellEdges
 
             forAll(fe, feI)
             {
-                labelSet_.insert(fe[feI]);    
+                labelSet_.insert(fe[feI]);
             }
         }
 
         storage.clear();
 
-        if (labelSet_.size() > storage.allocSize())
+        if (labelSet_.size() > storage.capacity())
         {
-            storage.setSize(labelSet_.size());
+            storage.setCapacity(labelSet_.size());
         }
 
         forAllConstIter(labelHashSet, labelSet_, iter)
@@ -673,7 +668,7 @@ const labelList& primitiveMesh::cellEdges
 }
 
 
-const labelList& primitiveMesh::cellEdges(const label cellI) const
+const Foam::labelList& Foam::primitiveMesh::cellEdges(const label cellI) const
 {
     return cellEdges(cellI, labels_);;
 }
@@ -681,6 +676,4 @@ const labelList& primitiveMesh::cellEdges(const label cellI) const
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace Foam
-
 // ************************************************************************* //

src/OpenFOAM/meshes/primitiveMesh/primitiveMeshPointCells.C

@@ -27,14 +27,9 @@ License
 #include "primitiveMesh.H"
 #include "cell.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-void primitiveMesh::calcPointCells() const
+void Foam::primitiveMesh::calcPointCells() const
 {
     // Loop through cells and mark up points
 
@@ -111,7 +106,7 @@ void primitiveMesh::calcPointCells() const
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-const labelListList& primitiveMesh::pointCells() const
+const Foam::labelListList& Foam::primitiveMesh::pointCells() const
 {
     if (!pcPtr_)
     {
@@ -122,7 +117,7 @@ const labelListList& primitiveMesh::pointCells() const
 }
 
 
-const labelList& primitiveMesh::pointCells
+const Foam::labelList& Foam::primitiveMesh::pointCells
 (
     const label pointI,
     DynamicList<label>& storage
@@ -155,26 +150,29 @@ const labelList& primitiveMesh::pointCells
         }
 
         // Filter duplicates
-        sort(storage);
+        if (storage.size() > 1)
-
-        label n = 1;
-
-        for (label i = 1; i < storage.size(); i++)
         {
-            if (storage[i] != storage[i-1])
+            sort(storage);
-            {
-                storage[n++] = storage[i];
-            }
-        }
 
-        storage.setSize(n);
+            label n = 1;
+            for (label i = 1; i < storage.size(); i++)
+            {
+                if (storage[i-1] != storage[i])
+                {
+                    storage[n++] = storage[i];
+                }
+            }
+
+            // truncate addressed list
+            storage.setSize(n);
+        }
 
         return storage;
     }
 }
 
 
-const labelList& primitiveMesh::pointCells(const label pointI) const
+const Foam::labelList& Foam::primitiveMesh::pointCells(const label pointI) const
 {
     return pointCells(pointI, labels_);
 }
@@ -182,6 +180,4 @@ const labelList& primitiveMesh::pointCells(const label pointI) const
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace Foam
-
 // ************************************************************************* //

src/OpenFOAM/meshes/primitiveMesh/primitiveMeshPointPoints.C

@@ -26,14 +26,9 @@ License
 
 #include "primitiveMesh.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
-void primitiveMesh::calcPointPoints() const
+void Foam::primitiveMesh::calcPointPoints() const
 {
     if (debug)
     {
@@ -94,7 +89,7 @@ void primitiveMesh::calcPointPoints() const
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-const labelListList& primitiveMesh::pointPoints() const
+const Foam::labelListList& Foam::primitiveMesh::pointPoints() const
 {
     if (!ppPtr_)
     {
@@ -105,7 +100,7 @@ const labelListList& primitiveMesh::pointPoints() const
 }
 
 
-const labelList& primitiveMesh::pointPoints
+const Foam::labelList& Foam::primitiveMesh::pointPoints
 (
     const label pointI,
     DynamicList<label>& storage
@@ -122,9 +117,9 @@ const labelList& primitiveMesh::pointPoints
 
         storage.clear();
 
-        if (pEdges.size() > storage.allocSize())
+        if (pEdges.size() > storage.capacity())
         {
-            storage.setSize(pEdges.size());
+            storage.setCapacity(pEdges.size());
         }
 
         forAll(pEdges, i)
@@ -137,7 +132,10 @@ const labelList& primitiveMesh::pointPoints
 }
 
 
-const labelList& primitiveMesh::pointPoints(const label pointI) const
+const Foam::labelList& Foam::primitiveMesh::pointPoints
+(
+    const label pointI
+) const
 {
     return pointPoints(pointI, labels_);
 }
@@ -145,6 +143,4 @@ const labelList& primitiveMesh::pointPoints(const label pointI) const
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace Foam
-
 // ************************************************************************* //

src/autoMesh/autoHexMesh/autoHexMeshDriver/autoHexMeshDriver.C

@@ -96,11 +96,7 @@ Foam::scalar Foam::autoHexMeshDriver::getMergeDistance(const scalar mergeTol)
 //)
 //{
 //    // Determine outside point.
-//    boundBox overallBb
+//    boundBox overallBb = boundBox::invertedBox;
-//    (
-//        point(GREAT, GREAT, GREAT),
-//        point(-GREAT, -GREAT, -GREAT)
-//    );
 //
 //    bool hasSurface = false;
 //

src/autoMesh/autoHexMesh/shellSurfaces/shellSurfaces.C

@@ -158,11 +158,7 @@ void Foam::shellSurfaces::setAndCheckLevels
 void Foam::shellSurfaces::orient()
 {
     // Determine outside point.
-    boundBox overallBb
+    boundBox overallBb = boundBox::invertedBox;
-    (
-        point(GREAT, GREAT, GREAT),
-        point(-GREAT, -GREAT, -GREAT)
-    );
 
     bool hasSurface = false;
 

src/conversion/polyDualMesh/polyDualMesh.C

@@ -236,13 +236,13 @@ Foam::labelList Foam::polyDualMesh::collectPatchSideFace
     if (pointToDualPoint[meshPointI] >= 0)
     {
         // Number of pFaces + 2 boundary edge + feature point
-        dualFace.setSize(pFaces.size()+2+1);
+        dualFace.setCapacity(pFaces.size()+2+1);
         // Store dualVertex for feature edge
         dualFace.append(pointToDualPoint[meshPointI]);
     }
     else
     {
-        dualFace.setSize(pFaces.size()+2);
+        dualFace.setCapacity(pFaces.size()+2);
     }
 
     // Store dual vertex for starting edge.
@@ -958,13 +958,13 @@ void Foam::polyDualMesh::calcDual
         if (edgeToDualPoint[edgeI] >= 0)
         {
             // Number of cells + 2 boundary faces + feature edge point
-            dualFace.setSize(mesh.edgeCells()[edgeI].size()+2+1);
+            dualFace.setCapacity(mesh.edgeCells()[edgeI].size()+2+1);
             // Store dualVertex for feature edge
             dualFace.append(edgeToDualPoint[edgeI]);
         }
         else
         {
-            dualFace.setSize(mesh.edgeCells()[edgeI].size()+2);
+            dualFace.setCapacity(mesh.edgeCells()[edgeI].size()+2);
         }
 
         // Store dual vertex for starting face.

src/dynamicMesh/boundaryMesh/boundaryMesh.C

@@ -1250,7 +1250,7 @@ void Foam::boundaryMesh::patchify
 
         forAll(patchFaces, newPatchI)
         {
-            patchFaces[newPatchI].setSize(nAvgFaces);
+            patchFaces[newPatchI].setCapacity(nAvgFaces);
         }
 
         //

src/dynamicMesh/meshCut/cellLooper/topoCellLooper.C

@@ -37,16 +37,11 @@ License
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
+// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 namespace Foam
 {
-
+   defineTypeNameAndDebug(topoCellLooper, 0);
-// * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
+   addToRunTimeSelectionTable(cellLooper, topoCellLooper, word);
-
-defineTypeNameAndDebug(topoCellLooper, 0);
-
-addToRunTimeSelectionTable(cellLooper, topoCellLooper, word);
-
-
 }
 
 // Angle for polys to be considered splitHexes.
@@ -67,7 +62,7 @@ void Foam::topoCellLooper::subsetList
 {
     if (startI == 0)
     {
-        // Truncate (setsize decides itself not to do anything if nothing
+        // Truncate (setSize decides itself not to do anything if nothing
         // changed)
         if (freeI < 0)
         {
@@ -228,7 +223,7 @@ Foam::labelList Foam::topoCellLooper::getSuperEdge
 
     do
     {
-        vertI = mesh().edges()[edgeI].otherVertex(vertI);    
+        vertI = mesh().edges()[edgeI].otherVertex(vertI);
 
         superVerts[superVertI++] = vertI;
 
@@ -237,7 +232,7 @@ Foam::labelList Foam::topoCellLooper::getSuperEdge
         edgeI = meshTools::otherEdge(mesh(), fEdges, edgeI, vertI);
     }
     while (!features.isFeaturePoint(prevEdgeI, edgeI));
-  
+
     superVerts.setSize(superVertI);
 
     return superVerts;
@@ -500,7 +495,7 @@ void Foam::topoCellLooper::walkSplitHex
                 nextEdgeI,
                 nextVertI
             );
-            
+
             edgeI = nextEdgeI;
             vertI = nextVertI;
         }
@@ -615,7 +610,7 @@ void Foam::topoCellLooper::walkSplitHex
                 labelList nextFaces =
                     getVertFacesNonEdge
                     (
-                        cellI, 
+                        cellI,
                         edgeI,
                         vertI
                     );
@@ -720,7 +715,7 @@ bool Foam::topoCellLooper::cut
     if (mesh().cellShapes()[cellI].model() == hex_)
     {
         // Let parent handle hex case.
-        return 
+        return
             hexCellLooper::cut
             (
                 refDir,
@@ -752,7 +747,7 @@ bool Foam::topoCellLooper::cut
 
             if (edgeI != -1)
             {
-                // Found non-feature edge. Start walking from vertex on edge. 
+                // Found non-feature edge. Start walking from vertex on edge.
                 vertI = mesh().edges()[edgeI].start();
             }
             else

src/dynamicMesh/polyMeshAdder/polyMeshAdder.C

@@ -47,7 +47,7 @@ void Foam::polyMeshAdder::append
     DynamicList<label>& dynLst
 )
 {
-    dynLst.setSize(dynLst.size() + lst.size());
+    dynLst.setCapacity(dynLst.size() + lst.size());
 
     forAll(lst, i)
     {
@@ -902,7 +902,7 @@ void Foam::polyMeshAdder::mergePointZones
     List<DynamicList<label> >& pzPoints
 )
 {
-    zoneNames.setSize(pz0.size() + pz1.size());
+    zoneNames.setCapacity(pz0.size() + pz1.size());
 
     // Names
     append(pz0.names(), zoneNames);
@@ -922,7 +922,7 @@ void Foam::polyMeshAdder::mergePointZones
     {
         DynamicList<label>& newZone = pzPoints[zoneI];
 
-        newZone.setSize(pz0[zoneI].size());
+        newZone.setCapacity(pz0[zoneI].size());
 
         append(from0ToAllPoints, pz0[zoneI], newZone);
     }
@@ -933,7 +933,7 @@ void Foam::polyMeshAdder::mergePointZones
         // Relabel all points of zone and add to correct pzPoints.
         DynamicList<label>& newZone = pzPoints[from1ToAll[zoneI]];
 
-        newZone.setSize(newZone.size() + pz1[zoneI].size());
+        newZone.setCapacity(newZone.size() + pz1[zoneI].size());
 
         append(from1ToAllPoints, pz1[zoneI], newZone);
     }
@@ -958,7 +958,7 @@ void Foam::polyMeshAdder::mergeFaceZones
     List<DynamicList<bool> >& fzFlips
 )
 {
-    zoneNames.setSize(fz0.size() + fz1.size());
+    zoneNames.setCapacity(fz0.size() + fz1.size());
 
     append(fz0.names(), zoneNames);
 
@@ -979,8 +979,8 @@ void Foam::polyMeshAdder::mergeFaceZones
         DynamicList<label>& newZone = fzFaces[zoneI];
         DynamicList<bool>& newFlip = fzFlips[zoneI];
 
-        newZone.setSize(fz0[zoneI].size());
+        newZone.setCapacity(fz0[zoneI].size());
-        newFlip.setSize(newZone.size());
+        newFlip.setCapacity(newZone.size());
 
         const labelList& addressing = fz0[zoneI];
         const boolList& flipMap = fz0[zoneI].flipMap();
@@ -1003,8 +1003,8 @@ void Foam::polyMeshAdder::mergeFaceZones
         DynamicList<label>& newZone = fzFaces[from1ToAll[zoneI]];
         DynamicList<bool>& newFlip = fzFlips[from1ToAll[zoneI]];
 
-        newZone.setSize(newZone.size() + fz1[zoneI].size());
+        newZone.setCapacity(newZone.size() + fz1[zoneI].size());
-        newFlip.setSize(newZone.size());
+        newFlip.setCapacity(newZone.size());
 
         const labelList& addressing = fz1[zoneI];
         const boolList& flipMap = fz1[zoneI].flipMap();
@@ -1040,7 +1040,7 @@ void Foam::polyMeshAdder::mergeCellZones
     List<DynamicList<label> >& czCells
 )
 {
-    zoneNames.setSize(cz0.size() + cz1.size());
+    zoneNames.setCapacity(cz0.size() + cz1.size());
 
     append(cz0.names(), zoneNames);
 
@@ -1056,7 +1056,7 @@ void Foam::polyMeshAdder::mergeCellZones
     czCells.setSize(zoneNames.size());
     forAll(cz0, zoneI)
     {
-        czCells[zoneI].setSize(cz0[zoneI].size());
+        czCells[zoneI].setCapacity(cz0[zoneI].size());
         // Insert mesh0 cells
         append(cz0[zoneI], czCells[zoneI]);
     }
@@ -1067,7 +1067,7 @@ void Foam::polyMeshAdder::mergeCellZones
     {
         DynamicList<label>& newZone = czCells[from1ToAll[zoneI]];
 
-        newZone.setSize(newZone.size() + cz1[zoneI].size());
+        newZone.setCapacity(newZone.size() + cz1[zoneI].size());
 
         append(from1ToAllCells, cz1[zoneI], newZone);
     }

src/dynamicMesh/polyMeshAdder/polyMeshAdderTemplates.C

@@ -31,7 +31,7 @@ License
 template<class T>
 void Foam::polyMeshAdder::append(const List<T>& lst, DynamicList<T>& dynLst)
 {
-    dynLst.setSize(dynLst.size() + lst.size());
+    dynLst.setCapacity(dynLst.size() + lst.size());
 
     forAll(lst, i)
     {

src/dynamicMesh/polyTopoChange/polyTopoChange/hexRef8.C

@@ -2432,8 +2432,8 @@ Foam::labelList Foam::hexRef8::consistentSlowRefinement
 
 
         // Transfer into seedFaces, seedFacesInfo
-        seedFaces.setSize(changedFacesInfo.size());
+        seedFaces.setCapacity(changedFacesInfo.size());
-        seedFacesInfo.setSize(changedFacesInfo.size());
+        seedFacesInfo.setCapacity(changedFacesInfo.size());
 
         forAllConstIter(Map<refinementData>, changedFacesInfo, iter)
         {

src/dynamicMesh/polyTopoChange/polyTopoChange/polyTopoChange.C

@@ -726,7 +726,7 @@ void Foam::polyTopoChange::getFaceOrder
 }
 
 
-// Compact and reorder faces according to map.
+// Reorder and compact faces according to map.
 void Foam::polyTopoChange::reorderCompactFaces
 (
     const label newSize,
@@ -2999,9 +2999,6 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::polyTopoChange::changeMesh
             patchStarts,
             syncParallel
         );
-        // Invalidate new points to go into map.
-        newPoints.clear();
-
         mesh.changing(true);
     }
 
@@ -3010,14 +3007,8 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::polyTopoChange::changeMesh
         retiredPoints_.clear();
         retiredPoints_.resize(0);
 
-        faces_.clear();
-        faces_.setSize(0);
         region_.clear();
         region_.setSize(0);
-        faceOwner_.clear();
-        faceOwner_.setSize(0);
-        faceNeighbour_.clear();
-        faceNeighbour_.setSize(0);
     }
 
 
@@ -3237,27 +3228,20 @@ Foam::autoPtr<Foam::mapPolyMesh> Foam::polyTopoChange::makeMesh
         new fvMesh
         (
             io,
-            xferMove<pointField>(newPoints),
+            xferMove(newPoints),
-            xferMove<faceList>(faces_),
+            xferMoveTo<faceList>(faces_),
-            xferMove<labelList>(faceOwner_),
+            xferMoveTo<labelList>(faceOwner_),
-            xferMove<labelList>(faceNeighbour_)
+            xferMoveTo<labelList>(faceNeighbour_)
         )
     );
     fvMesh& newMesh = newMeshPtr();
 
     // Clear out primitives
     {
-        newPoints.clear();
         retiredPoints_.clear();
         retiredPoints_.resize(0);
-        faces_.clear();
-        faces_.setSize(0);
         region_.clear();
         region_.setSize(0);
-        faceOwner_.clear();
-        faceOwner_.setSize(0);
-        faceNeighbour_.clear();
-        faceNeighbour_.setSize(0);
     }
 
 

src/dynamicMesh/polyTopoChange/polyTopoChange/refinementHistory.C

@@ -411,7 +411,7 @@ Foam::refinementHistory::refinementHistory
     else
     {
         visibleCells_.setSize(nCells);
-        splitCells_.setSize(nCells);
+        splitCells_.setCapacity(nCells);
 
         for (label cellI = 0; cellI < nCells; cellI++)
         {

src/engine/ignition/ignitionSite.C

@@ -37,14 +37,17 @@ namespace Foam
 
 void ignitionSite::findIgnitionCells(const fvMesh& mesh)
 {
+    // Bit tricky: generate C and V before shortcutting if cannot find
+    // cell locally. mesh.C generation uses parallel communication.
+    const volVectorField& centres = mesh.C();
+    const scalarField& vols = mesh.V();
+
     label ignCell = mesh.findCell(location_);
     if (ignCell == -1)
     {
         return;
     }
 
-    const volVectorField& centres = mesh.C();
-    const scalarField& vols = mesh.V();
     scalar radius = diameter_/2.0;
 
     cells_.setSize(1);

src/finiteVolume/cfdTools/general/findRefCell/findRefCell.C

@@ -49,6 +49,22 @@ void Foam::setRefCell
             if (Pstream::master())
             {
                 refCelli = readLabel(dict.lookup(refCellName));
+
+                if (refCelli < 0 || refCelli >= field.mesh().nCells())
+                {
+                    FatalErrorIn
+                    (
+                        "void Foam::setRefCell"
+                         "("
+                         "    const volScalarField&,"
+                         "    const dictionary&,"
+                         "    label& scalar&,"
+                         "    bool"
+                         ")"
+                    )   << "Illegal master cellID " << refCelli
+                        << ". Should be 0.." << field.mesh().nCells()
+                        << exit(FatalError);
+                }
             }
             else
             {
@@ -75,7 +91,7 @@ void Foam::setRefCell
                 )
                   << "Unable to set reference cell for field " << field.name()
                   << nl << "    Reference point " << refPointName
-                  << " found on multiple domains" << nl << abort(FatalError);
+                  << " found on multiple domains" << nl << exit(FatalError);
             }
         }
         else
@@ -92,7 +108,7 @@ void Foam::setRefCell
             )
               << "Unable to set reference cell for field" << field.name() << nl
               << "    Please supply either " << refCellName
-              << " or " << refPointName << nl << abort(FatalError);
+              << " or " << refPointName << nl << exit(FatalError);
         }
 
         refValue = readScalar(dict.lookup(refValueName));

src/finiteVolume/fields/fvPatchFields/derived/timeVaryingMappedFixedValue/timeVaryingMappedFixedValueFvPatchField.C

@@ -207,8 +207,11 @@ void timeVaryingMappedFixedValueFvPatchField<Type>::autoMap
 )
 {
     fixedValueFvPatchField<Type>::autoMap(m);
-    startSampledValues_.autoMap(m);
+    if (startSampledValues_.size() > 0)
-    endSampledValues_.autoMap(m);
+    {
+        startSampledValues_.autoMap(m);
+        endSampledValues_.autoMap(m);
+    }
 }
 
 

src/lagrangian/dieselSpray/Make/options

@@ -1,6 +1,7 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude \
+    -I$(LIB_SRC)/turbulenceModels \
     -I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
     -I$(LIB_SRC)/turbulenceModels/compressible/RAS/lnInclude \
     -I$(LIB_SRC)/turbulenceModels/LES/LESdeltas/lnInclude \

src/lagrangian/intermediate/Make/options

@@ -10,6 +10,7 @@ EXE_INC = \
     -I$(LIB_SRC)/thermophysicalModels/solidMixture/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/combustion/lnInclude \
     -I$(LIB_SRC)/thermophysicalModels/radiation/lnInclude \
+    -I$(LIB_SRC)/turbulenceModels \
     -I$(LIB_SRC)/turbulenceModels/compressible/turbulenceModel \
     -I$(LIB_SRC)/turbulenceModels/compressible/RAS/lnInclude \
     -I$(LIB_SRC)/turbulenceModels/LES/LESdeltas/lnInclude \

src/meshTools/indexedOctree/indexedOctree.C

@@ -157,7 +157,7 @@ void indexedOctree<Type>::divide
     List<DynamicList<label> > subIndices(8);
     for (direction octant = 0; octant < subIndices.size(); octant++)
     {
-        subIndices[octant].setSize(indices.size()/8);
+        subIndices[octant].setCapacity(indices.size()/8);
     }
 
     // Precalculate bounding boxes.

src/meshTools/octree/treeBoundBox.C

@@ -110,8 +110,10 @@ Foam::treeBoundBox::treeBoundBox(const UList<point>& points)
 {
     if (points.size() == 0)
     {
-        WarningIn("treeBoundBox::treeBoundBox(const UList<point>&)")
+        WarningIn
-            << "cannot find bounding box for zero sized pointField"
+        (
+            "treeBoundBox::treeBoundBox(const UList<point>&)"
+        )   << "cannot find bounding box for zero-sized pointField"
             << "returning zero" << endl;
 
         return;
@@ -132,17 +134,18 @@ Foam::treeBoundBox::treeBoundBox(const UList<point>& points)
 Foam::treeBoundBox::treeBoundBox
 (
     const UList<point>& points,
-    const labelList& meshPoints
+    const UList<label>& meshPoints
 )
 :
     boundBox()
 {
-    if (meshPoints.size() == 0)
+    if (points.size() == 0 || meshPoints.size() == 0)
     {
         WarningIn
         (
-            "treeBoundBox::treeBoundBox(const UList<point>&, const labelList)"
+            "treeBoundBox::treeBoundBox"
-        )   << "cannot find bounding box for zero sized pointField"
+            "(const UList<point>&, const UList<label>&)"
+        )   << "cannot find bounding box for zero-sized pointField"
             << "returning zero" << endl;
 
         return;

src/meshTools/octree/treeBoundBox.H

@@ -87,6 +87,7 @@ public:
 
     // Static data members
 
+        //- As per boundBox::greatBox, but with GREAT instead of VGREAT
         static const treeBoundBox greatBox;
 
         //- Bits used for octant/point coding. Every octant/corner point
@@ -171,7 +172,7 @@ public:
         treeBoundBox(const UList<point>& points);
 
         //- Construct as subset of points
-        treeBoundBox(const UList<point>&, const labelList& meshPoints);
+        treeBoundBox(const UList<point>&, const UList<label>& meshPoints);
 
         //- Construct from Istream
         treeBoundBox(Istream&);

src/meshTools/octree/treeBoundBoxI.H

@@ -27,29 +27,21 @@ License
 #include "treeBoundBox.H"
 #include "Random.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-// Construct null setting points to zero
+inline Foam::treeBoundBox::treeBoundBox()
-inline treeBoundBox::treeBoundBox()
 :
     boundBox()
 {}
 
 
-// Construct from components
+inline Foam::treeBoundBox::treeBoundBox(const point& min, const point& max)
-inline treeBoundBox::treeBoundBox(const point& min, const point& max)
 :
     boundBox(min, max)
 {}
 
 
-// Construct from components
+inline Foam::treeBoundBox::treeBoundBox(const boundBox& bb)
-inline treeBoundBox::treeBoundBox(const boundBox& bb)
 :
     boundBox(bb)
 {}
@@ -57,7 +49,7 @@ inline treeBoundBox::treeBoundBox(const boundBox& bb)
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline scalar treeBoundBox::minDim() const
+inline Foam::scalar Foam::treeBoundBox::minDim() const
 {
     return ::Foam::min
     (
@@ -71,7 +63,7 @@ inline scalar treeBoundBox::minDim() const
 }
 
 
-inline scalar treeBoundBox::maxDim() const
+inline Foam::scalar Foam::treeBoundBox::maxDim() const
 {
     return ::Foam::max
     (
@@ -85,7 +77,7 @@ inline scalar treeBoundBox::maxDim() const
 }
 
 
-inline scalar treeBoundBox::avgDim() const
+inline Foam::scalar Foam::treeBoundBox::avgDim() const
 {
     return
     (
@@ -96,19 +88,19 @@ inline scalar treeBoundBox::avgDim() const
 }
 
 
-inline scalar treeBoundBox::typDim() const
+inline Foam::scalar Foam::treeBoundBox::typDim() const
 {
     return avgDim();
 }
 
 
-inline point treeBoundBox::mid() const
+inline Foam::point Foam::treeBoundBox::mid() const
 {
     return 0.5*(min() + max());
 }
 
 
-inline point treeBoundBox::corner(const direction octant) const
+inline Foam::point Foam::treeBoundBox::corner(const direction octant) const
 {
     return point
     (
@@ -119,7 +111,7 @@ inline point treeBoundBox::corner(const direction octant) const
 }
 
 // Returns octant in which sample resides. Reverse of subBbox.
-inline direction treeBoundBox::subOctant(const point& sample) const
+inline Foam::direction Foam::treeBoundBox::subOctant(const point& sample) const
 {
     point mid = 0.5*(max() + min());
 
@@ -146,7 +138,7 @@ inline direction treeBoundBox::subOctant(const point& sample) const
 
 // Returns octant in which sample resides. Reverse of subBbox. Precalculated
 // midpoint
-inline direction treeBoundBox::subOctant
+inline Foam::direction Foam::treeBoundBox::subOctant
 (
     const point& mid,
     const point& sample
@@ -175,8 +167,11 @@ inline direction treeBoundBox::subOctant
 
 // Returns octant in which sample resides. Reverse of subBbox. Flags sample
 // exactly on edge.
-inline direction treeBoundBox::subOctant(const point& sample, bool& onEdge)
+inline Foam::direction Foam::treeBoundBox::subOctant
- const
+(
+    const point& sample,
+    bool& onEdge
+) const
 {
     point mid = 0.5*(max() + min());
 
@@ -216,7 +211,7 @@ inline direction treeBoundBox::subOctant(const point& sample, bool& onEdge)
 
 // Returns octant in which sample resides. Reverse of subBbox. Precalculated
 // midpoint
-inline direction treeBoundBox::subOctant
+inline Foam::direction Foam::treeBoundBox::subOctant
 (
     const point& mid,
     const point& sample,
@@ -261,7 +256,7 @@ inline direction treeBoundBox::subOctant
 // Precalculated midpoint. If the sample is on the dividing line between
 // the octants the direction vector determines which octant to use
 // (i.e. in which octant the sample would be if it were moved along dir)
-inline direction treeBoundBox::subOctant
+inline Foam::direction Foam::treeBoundBox::subOctant
 (
     const point& mid,
     const vector& dir,
@@ -319,9 +314,9 @@ inline direction treeBoundBox::subOctant
 }
 
 
-// Returns reference to octantOrder which defines the 
+// Returns reference to octantOrder which defines the
 // order to do the search.
-inline void treeBoundBox::searchOrder
+inline void Foam::treeBoundBox::searchOrder
 (
     const point& sample,
     FixedList<direction,8>& octantOrder
@@ -380,7 +375,7 @@ inline void treeBoundBox::searchOrder
 	{
 	    min = treeBoundBox::FRONTHALF;
 	    mid = treeBoundBox::TOPHALF;
-	    max = treeBoundBox::RIGHTHALF;	    
+	    max = treeBoundBox::RIGHTHALF;
 	}
 	else if( dist.x() < dist.z())
 	{
@@ -412,13 +407,13 @@ inline void treeBoundBox::searchOrder
 
 // true if bb's intersect or overlap.
 // Note: <= to make sure we catch all.
-inline bool treeBoundBox::overlaps(const treeBoundBox& bb) const
+inline bool Foam::treeBoundBox::overlaps(const treeBoundBox& bb) const
 {
     return boundBox::overlaps(bb);
 }
 
 
-inline bool treeBoundBox::contains(const point& sample) const
+inline bool Foam::treeBoundBox::contains(const point& sample) const
 {
     return
     (
@@ -433,7 +428,11 @@ inline bool treeBoundBox::contains(const point& sample) const
 
 
 //- Return slightly wider bounding box
-inline treeBoundBox treeBoundBox::extend(Random& rndGen, const scalar s) const
+inline Foam::treeBoundBox Foam::treeBoundBox::extend
+(
+    Random& rndGen,
+    const scalar s
+) const
 {
     treeBoundBox bb(*this);
 
@@ -456,6 +455,4 @@ inline treeBoundBox treeBoundBox::extend(Random& rndGen, const scalar s) const
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
-} // End namespace Foam
-
 // ************************************************************************* //

src/meshTools/searchableSurface/distributedTriSurfaceMesh.C

@@ -97,7 +97,11 @@ void Foam::distributedTriSurfaceMesh::splitSegment
     point clipPt0, clipPt1;
 
 
-    // 1. Fully local already handled outside
+    // 1. Fully local already handled outside. Note: retest is cheap.
+    if (isLocal(procBb_[Pstream::myProcNo()], start, end))
+    {
+        return;
+    }
 
 
     // 2. Check if fully inside other processor. Rare occurrence
@@ -109,17 +113,14 @@ void Foam::distributedTriSurfaceMesh::splitSegment
         {
             const List<treeBoundBox>& bbs = procBb_[procI];
 
-            forAll(bbs, bbI)
+            if (isLocal(bbs, start, end))
             {
-                if (bbs[bbI].contains(start) && bbs[bbI].contains(end))
+                //Pout<< "    Completely remote segment:"
-                {
+                //    << start << end << " on proc:" << procI << endl;
-                    //Pout<< "    Completely remote segment:"
+                sendMap[procI].append(allSegments.size());
-                    //    << start << end << " on proc:" << procI << endl;
+                allSegmentMap.append(segmentI);
-                    sendMap[procI].append(allSegments.size());
+                allSegments.append(segment(start, end));
-                    allSegmentMap.append(segmentI);
+                return;
-                    allSegments.append(segment(start, end));
-                    return;
-                }
             }
         }
     }
@@ -776,8 +777,7 @@ void Foam::distributedTriSurfaceMesh::calcBounds
     pointIsUsed = 0U;
 
     nPoints = 0;
-    bb.min() = point(VGREAT, VGREAT, VGREAT);
+    bb = boundBox::invertedBox;
-    bb.max() = point(-VGREAT, -VGREAT, -VGREAT);
 
     const triSurface& s = static_cast<const triSurface&>(*this);
 
@@ -2123,7 +2123,7 @@ void Foam::distributedTriSurfaceMesh::distribute
            }
         }
     }
-    
+
 
     faceMap.reset
     (

src/meshTools/searchableSurface/searchableSphere.C

@@ -292,7 +292,8 @@ void Foam::searchableSphere::getNormal
         if (info[i].hit())
         {
             normal[i] = info[i].hitPoint() - centre_;
-            normal[i] /= mag(normal[i]);
+
+            normal[i] /= mag(normal[i])+VSMALL;
         }
         else
         {

src/meshTools/triSurface/surfaceFeatures/surfaceFeatures.C

@@ -247,7 +247,7 @@ Foam::label Foam::surfaceFeatures::nextFeatEdge
 // prevPointI. Marks feature edges visited in featVisited by assigning them
 // the current feature line number. Returns cumulative length of edges walked.
 // Works in one of two modes:
-// - mark : step to edges with featVisited = -1. 
+// - mark : step to edges with featVisited = -1.
 //          Mark edges visited with currentFeatI.
 // - clear : step to edges with featVisited = currentFeatI
 //           Mark edges visited with -2 and erase from feature edges.
@@ -257,7 +257,7 @@ Foam::surfaceFeatures::labelScalar Foam::surfaceFeatures::walkSegment
     const List<edgeStatus>& edgeStat,
     const label startEdgeI,
     const label startPointI,
-    const label currentFeatI,     
+    const label currentFeatI,
     labelList& featVisited
 )
 {
@@ -360,7 +360,7 @@ Foam::surfaceFeatures::surfaceFeatures
     const labelList& featureEdges,
     const label externalStart,
     const label internalStart
-)    
+)
 :
     surf_(surf),
     featurePoints_(featurePoints),
@@ -457,7 +457,7 @@ Foam::labelList Foam::surfaceFeatures::selectFeatureEdges
 
     if (regionEdges)
     {
-        selectedEdges.setSize(selectedEdges.size() + nRegionEdges());
+        selectedEdges.setCapacity(selectedEdges.size() + nRegionEdges());
 
         for (label i = 0; i < externalStart_; i++)
         {
@@ -467,7 +467,7 @@ Foam::labelList Foam::surfaceFeatures::selectFeatureEdges
 
     if (externalEdges)
     {
-        selectedEdges.setSize(selectedEdges.size() + nExternalEdges());
+        selectedEdges.setCapacity(selectedEdges.size() + nExternalEdges());
 
         for (label i = externalStart_; i < internalStart_; i++)
         {
@@ -477,7 +477,7 @@ Foam::labelList Foam::surfaceFeatures::selectFeatureEdges
 
     if (internalEdges)
     {
-        selectedEdges.setSize(selectedEdges.size() + nInternalEdges());
+        selectedEdges.setCapacity(selectedEdges.size() + nInternalEdges());
 
         for (label i = internalStart_; i < featureEdges_.size(); i++)
         {
@@ -530,8 +530,8 @@ void Foam::surfaceFeatures::findFeatures(const scalar includedAngle)
 
                     // Check if convex or concave by looking at angle
                     // between face centres and normal
-                    vector f0Tof1 = 
+                    vector f0Tof1 =
-                        surf_[face1].centre(points) 
+                        surf_[face1].centre(points)
                         - surf_[face0].centre(points);
 
                     if ((f0Tof1 & faceNormals[face0]) > 0.0)
@@ -683,11 +683,11 @@ void Foam::surfaceFeatures::writeDict(Ostream& writeFile) const
 {
 
     dictionary featInfoDict;
-    featInfoDict.add("externalStart", externalStart_); 
+    featInfoDict.add("externalStart", externalStart_);
     featInfoDict.add("internalStart", internalStart_);
     featInfoDict.add("featureEdges", featureEdges_);
     featInfoDict.add("featurePoints", featurePoints_);
-    
+
     featInfoDict.write(writeFile);
 }
 
@@ -1152,7 +1152,7 @@ void Foam::surfaceFeatures::nearestSurfEdge
 (
     const labelList& selectedEdges,
     const pointField& samples,
-    const vector& searchSpan,   // Search span 
+    const vector& searchSpan,   // Search span
     labelList& edgeLabel,
     labelList& edgeEndPoint,
     pointField& edgePoint
@@ -1163,7 +1163,7 @@ void Foam::surfaceFeatures::nearestSurfEdge
     edgePoint.setSize(samples.size());
 
     const pointField& localPoints = surf_.localPoints();
-    
+
     octree<octreeDataEdges> ppTree
     (
         treeBoundBox(localPoints),  // overall search domain
@@ -1232,7 +1232,7 @@ void Foam::surfaceFeatures::nearestSurfEdge
     const labelList& selectedSampleEdges,
     const pointField& samplePoints,
 
-    const vector& searchSpan,   // Search span 
+    const vector& searchSpan,   // Search span
     labelList& edgeLabel,       // label of surface edge or -1
     pointField& pointOnEdge,    // point on above edge
     pointField& pointOnFeature  // point on sample edge
@@ -1242,7 +1242,6 @@ void Foam::surfaceFeatures::nearestSurfEdge
     pointOnEdge.setSize(selectedSampleEdges.size());
     pointOnFeature.setSize(selectedSampleEdges.size());
 
-    
 
     octree<octreeDataEdges> ppTree
     (

src/sampling/Make/files

@@ -12,11 +12,13 @@ sampledSet/midPointAndFace/midPointAndFaceSet.C
 sampledSet/sampledSets/sampledSets.C
 sampledSet/sampledSetsFunctionObject/sampledSetsFunctionObject.C
 
-sampledSet/writers/writer/writers.C
+setWriters = sampledSet/writers
-sampledSet/writers/raw/rawWriters.C
+
-sampledSet/writers/xmgr/xmgrWriters.C
+$(setWriters)/writers.C
-sampledSet/writers/gnuplot/gnuplotWriters.C
+$(setWriters)/gnuplot/gnuplotSetWriterRunTime.C
-sampledSet/writers/jplot/jplotWriters.C
+$(setWriters)/jplot/jplotSetWriterRunTime.C
+$(setWriters)/raw/rawSetWriterRunTime.C
+$(setWriters)/xmgrace/xmgraceSetWriterRunTime.C
 
 cuttingPlane/cuttingPlane.C
 
@@ -31,13 +33,15 @@ sampledSurface/sampledSurface/sampledSurface.C
 sampledSurface/sampledSurfaces/sampledSurfaces.C
 sampledSurface/sampledSurfacesFunctionObject/sampledSurfacesFunctionObject.C
 
-sampledSurface/writers/surfaceWriters.C
+surfWriters = sampledSurface/writers
-sampledSurface/writers/foamFile/foamFileWriters.C
+
-sampledSurface/writers/dx/dxWriters.C
+$(surfWriters)/surfaceWriters.C
-sampledSurface/writers/raw/rawSurfaceWriters.C
+$(surfWriters)/dx/dxSurfaceWriterRunTime.C
-sampledSurface/writers/vtk/vtkWriters.C
+$(surfWriters)/foamFile/foamFileSurfaceWriterRunTime.C
-sampledSurface/writers/stl/stlWriters.C
+$(surfWriters)/null/nullSurfaceWriterRunTime.C
-sampledSurface/writers/null/nullWriters.C
+$(surfWriters)/raw/rawSurfaceWriterRunTime.C
+$(surfWriters)/stl/stlSurfaceWriterRunTime.C
+$(surfWriters)/vtk/vtkSurfaceWriterRunTime.C
 
 graphField/writePatchGraph.C
 graphField/writeCellGraph.C

src/sampling/Make/options

@@ -1,10 +1,12 @@
 EXE_INC = \
     -I$(LIB_SRC)/finiteVolume/lnInclude \
     -I$(LIB_SRC)/meshTools/lnInclude \
+    -I$(LIB_SRC)/surfMesh/lnInclude \
     -I$(LIB_SRC)/triSurface/lnInclude \
     -I$(LIB_SRC)/lagrangian/basic/lnInclude
 
 LIB_LIBS = \
     -lfiniteVolume \
     -lmeshTools \
+    -lsurfMesh \
     -ltriSurface

src/sampling/cuttingPlane/cuttingPlane.C

@@ -45,7 +45,7 @@ void Foam::cuttingPlane::calcCutCells
 (
     const primitiveMesh& mesh,
     const scalarField& dotProducts,
-    const labelList& cellIdLabels
+    const UList<label>& cellIdLabels
 )
 {
     const labelListList& cellEdges = mesh.cellEdges();
@@ -103,11 +103,12 @@ void Foam::cuttingPlane::calcCutCells
 
 // Determine for each edge the intersection point. Calculates
 // - cutPoints_ : coordinates of all intersection points
-// - edgePoint  : per edge -1 or the index into cutPoints_
+// - edgePoint  : per edge -1 or the index into cutPoints
-Foam::labelList Foam::cuttingPlane::intersectEdges
+void Foam::cuttingPlane::intersectEdges
 (
     const primitiveMesh& mesh,
-    const scalarField& dotProducts
+    const scalarField& dotProducts,
+    List<label>& edgePoint
 )
 {
     // Use the dotProducts to find out the cut edges.
@@ -115,7 +116,7 @@ Foam::labelList Foam::cuttingPlane::intersectEdges
     const pointField& points = mesh.points();
 
     // Per edge -1 or the label of the intersection point
-    labelList edgePoint(edges.size(), -1);
+    edgePoint.setSize(edges.size());
 
     DynamicList<point> dynCuttingPoints(4*cutCells_.size());
 
@@ -129,7 +130,9 @@ Foam::labelList Foam::cuttingPlane::intersectEdges
          || (dotProducts[e[1]] < zeroish && dotProducts[e[0]] > positive)
         )
         {
-            // Edge is cut.
+            // Edge is cut
+            edgePoint[edgeI] = dynCuttingPoints.size();
+
             const point& p0 = points[e[0]];
             const point& p1 = points[e[1]];
 
@@ -139,7 +142,7 @@ Foam::labelList Foam::cuttingPlane::intersectEdges
             {
                 dynCuttingPoints.append(p0);
             }
-            else if (alpha > 1.0)
+            else if (alpha >= 1.0)
             {
                 dynCuttingPoints.append(p1);
             }
@@ -147,15 +150,14 @@ Foam::labelList Foam::cuttingPlane::intersectEdges
             {
                 dynCuttingPoints.append((1-alpha)*p0 + alpha*p1);
             }
-
+        }
-            edgePoint[edgeI] = dynCuttingPoints.size() - 1;
+        else
+        {
+            edgePoint[edgeI] = -1;
         }
     }
 
-    dynCuttingPoints.shrink();
+    this->storedPoints().transfer(dynCuttingPoints);
-    cutPoints_.transfer(dynCuttingPoints);
-
-    return edgePoint;
 }
 
 
@@ -164,7 +166,7 @@ Foam::labelList Foam::cuttingPlane::intersectEdges
 bool Foam::cuttingPlane::walkCell
 (
     const primitiveMesh& mesh,
-    const labelList& edgePoint,
+    const UList<label>& edgePoint,
     const label cellI,
     const label startEdgeI,
     DynamicList<label>& faceVerts
@@ -175,6 +177,7 @@ bool Foam::cuttingPlane::walkCell
 
     label nIter = 0;
 
+    faceVerts.clear();
     do
     {
         faceVerts.append(edgePoint[edgeI]);
@@ -255,11 +258,17 @@ bool Foam::cuttingPlane::walkCell
 void Foam::cuttingPlane::walkCellCuts
 (
     const primitiveMesh& mesh,
-    const labelList& edgePoint
+    const UList<label>& edgePoint
 )
 {
-    cutFaces_.setSize(cutCells_.size());
+    const pointField& cutPoints = this->points();
-    label cutFaceI = 0;
+
+    // use dynamic lists to handle triangulation and/or missed cuts
+    DynamicList<face>  dynCutFaces(cutCells_.size());
+    DynamicList<label> dynCutCells(cutCells_.size());
+
+    // scratch space for calculating the face vertices
+    DynamicList<label> faceVerts(10);
 
     forAll(cutCells_, i)
     {
@@ -290,7 +299,6 @@ void Foam::cuttingPlane::walkCellCuts
         }
 
         // Walk from starting edge around the circumference of the cell.
-        DynamicList<label> faceVerts(2*mesh.faces()[cellI].size());
         bool okCut = walkCell
         (
             mesh,
@@ -302,54 +310,46 @@ void Foam::cuttingPlane::walkCellCuts
 
         if (okCut)
         {
-            faceVerts.shrink();
+            face f(faceVerts);
 
-            face cutFace(faceVerts);
+            // Orient face to point in the same direction as the plane normal
-
+            if ((f.normal(cutPoints) && normal()) < 0)
-            // Orient face.
-            if ((cutFace.normal(cutPoints_) && normal()) < 0)
             {
-                cutFace = cutFace.reverseFace();
+                f = f.reverseFace();
             }
 
-            cutFaces_[cutFaceI++] = cutFace;
+            // the cut faces are usually quite ugly, so always triangulate
+            label nTri = f.triangles(cutPoints, dynCutFaces);
+            while (nTri--)
+            {
+                dynCutCells.append(cellI);
+            }
         }
     }
 
-    cutFaces_.setSize(cutFaceI);
+    this->storedFaces().transfer(dynCutFaces);
+    cutCells_.transfer(dynCutCells);
 }
 
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
 // Construct without cutting
-Foam::cuttingPlane::cuttingPlane(const plane& newPlane)
+Foam::cuttingPlane::cuttingPlane(const plane& pln)
 :
-    plane(newPlane)
+    plane(pln)
 {}
 
 
-// Construct from components
-Foam::cuttingPlane::cuttingPlane
-(
-    const primitiveMesh& mesh,
-    const plane& newPlane
-)
-:
-    plane(newPlane)
-{
-    reCut(mesh);
-}
-
 // Construct from mesh reference and plane, restricted to a list of cells
 Foam::cuttingPlane::cuttingPlane
 (
     const primitiveMesh& mesh,
-    const plane& newPlane,
+    const plane& pln,
-    const labelList& cellIdLabels
+    const UList<label>& cellIdLabels
 )
 :
-    plane(newPlane)
+    plane(pln)
 {
     reCut(mesh, cellIdLabels);
 }
@@ -358,101 +358,51 @@ Foam::cuttingPlane::cuttingPlane
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-// reCut mesh with existing planeDesc
-void Foam::cuttingPlane::reCut
-(
-    const primitiveMesh& mesh
-)
-{
-    cutCells_.clear();
-    cutPoints_.clear();
-    cutFaces_.clear();
-
-    scalarField dotProducts = (mesh.points() - refPoint()) & normal();
-
-    //// Perturb points cuts so edges are cut properly.
-    //const tmp<scalarField> tdotProducts = stabilise(rawDotProducts, SMALL);
-    //const scalarField& dotProducts = tdotProducts();
-
-    // Determine cells that are (probably) cut.
-    calcCutCells(mesh, dotProducts);
-
-    // Determine cutPoints and return list of edge cuts. (per edge -1 or
-    // the label of the intersection point (in cutPoints_)
-    labelList edgePoint(intersectEdges(mesh, dotProducts));
-
-    // Do topological walk around cell to find closed loop.
-    walkCellCuts(mesh, edgePoint);
-}
-
-
 // recut mesh with existing planeDesc
 void Foam::cuttingPlane::reCut
 (
     const primitiveMesh& mesh,
-    const labelList& cellIdLabels
+    const UList<label>& cellIdLabels
 )
 {
+    MeshStorage::clear();
     cutCells_.clear();
-    cutPoints_.clear();
-    cutFaces_.clear();
 
     scalarField dotProducts = (mesh.points() - refPoint()) & normal();
 
     // Determine cells that are (probably) cut.
     calcCutCells(mesh, dotProducts, cellIdLabels);
 
-    // Determine cutPoints and return list of edge cuts. (per edge -1 or
+    // Determine cutPoints and return list of edge cuts.
-    // the label of the intersection point (in cutPoints_)
+    // per edge -1 or the label of the intersection point
-    labelList edgePoint(intersectEdges(mesh, dotProducts));
+    labelList edgePoint;
+    intersectEdges(mesh, dotProducts, edgePoint);
 
     // Do topological walk around cell to find closed loop.
     walkCellCuts(mesh, edgePoint);
 }
 
 
-
+// remap action on triangulation
-// Return plane used
+void Foam::cuttingPlane::remapFaces
-const Foam::plane& Foam::cuttingPlane::planeDesc() const
+(
+    const UList<label>& faceMap
+)
 {
-    return static_cast<const plane&>(*this);
+    // recalculate the cells cut
-}
+    if (&faceMap && faceMap.size())
-
-
-// Return vectorField of cutting points
-const Foam::pointField& Foam::cuttingPlane::points() const
-{
-    return cutPoints_;
-}
-
-
-// Return unallocFaceList of points in cells
-const Foam::faceList& Foam::cuttingPlane::faces() const
-{
-    return cutFaces_;
-}
-
-
-// Return labelList of cut cells
-const Foam::labelList& Foam::cuttingPlane::cells() const
-{
-    return cutCells_;
-}
-
-
-bool Foam::cuttingPlane::cut()
-{
-    if (cutCells_.size() > 0)
     {
-        return true;
+        MeshStorage::remapFaces(faceMap);
-    }
+
-    else
+        List<label> newCutCells(faceMap.size());
-    {
+        forAll(faceMap, faceI)
-        return false;
+        {
+            newCutCells[faceI] = cutCells_[faceMap[faceI]];
+        }
+        cutCells_.transfer(newCutCells);
     }
 }
 
-
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
 void Foam::cuttingPlane::operator=(const cuttingPlane& rhs)
@@ -465,11 +415,9 @@ void Foam::cuttingPlane::operator=(const cuttingPlane& rhs)
             << abort(FatalError);
     }
 
+    static_cast<MeshStorage&>(*this) = rhs;
     static_cast<plane&>(*this) = rhs;
-
+    cutCells_ = rhs.cutCells();
-    cutCells_ = rhs.cells();
-    cutPoints_ = rhs.points();
-    cutFaces_ = rhs.faces();
 }
 
 

src/sampling/cuttingPlane/cuttingPlane.H

@@ -28,7 +28,8 @@ Class
 Description
     Constructs plane through mesh.
 
-    No attempt at resolving degenerate cases.
+    No attempt at resolving degenerate cases. Since the cut faces are
+    usually quite ugly, they will always be triangulated.
 
 Note
     When the cutting plane coincides with a mesh face, the cell edge on the
@@ -45,6 +46,7 @@ SourceFiles
 #include "plane.H"
 #include "pointField.H"
 #include "faceList.H"
+#include "MeshedSurface.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -54,25 +56,22 @@ namespace Foam
 class primitiveMesh;
 
 /*---------------------------------------------------------------------------*\
-                           Class cuttingPlane Declaration
+                       Class cuttingPlane Declaration
 \*---------------------------------------------------------------------------*/
 
 class cuttingPlane
 :
+    public PrimitiveMeshedSurface<face>,
     public plane
 {
+    //- Private typedefs for convenience
+        typedef PrimitiveMeshedSurface<face> MeshStorage;
+
     // Private data
 
         //- List of cells cut by the plane
         labelList cutCells_;
 
-        //- Intersection points
-        pointField cutPoints_;
-
-        //- Cut faces in terms of intersection points
-        faceList cutFaces_;
-
-
     // Private Member Functions
 
         //- Determine cut cells, possibly restricted to a list of cells
@@ -80,22 +79,23 @@ class cuttingPlane
         (
             const primitiveMesh&,
             const scalarField& dotProducts,
-            const labelList& cellIdLabels = labelList::null()
+            const UList<label>& cellIdLabels = UList<label>::null()
         );
 
-        //- Determine intersection points (cutPoints_).
+        //- Determine intersection points (cutPoints).
-        labelList intersectEdges
+        void intersectEdges
         (
             const primitiveMesh&,
-            const scalarField& dotProducts
+            const scalarField& dotProducts,
+            List<label>& edgePoint
         );
 
-        //- Walk around circumference of cell starting from startEdgeI crossing
+        //- Walk circumference of cell, starting from startEdgeI crossing
         //  only cut edges. Record cutPoint labels in faceVerts.
         static bool walkCell
         (
             const primitiveMesh&,
-            const labelList& edgePoint,
+            const UList<label>& edgePoint,
             const label cellI,
             const label startEdgeI,
             DynamicList<label>& faceVerts
@@ -105,7 +105,7 @@ class cuttingPlane
         void walkCellCuts
         (
             const primitiveMesh& mesh,
-            const labelList& edgePoint
+            const UList<label>& edgePoint
         );
 
 
@@ -116,53 +116,51 @@ protected:
         //- Construct plane description without cutting
         cuttingPlane(const plane&);
 
-
     // Protected Member Functions
 
-        //- recut mesh with existing planeDesc
-        void reCut(const primitiveMesh&);
-
         //- recut mesh with existing planeDesc, restricted to a list of cells
         void reCut
         (
             const primitiveMesh&,
-            const labelList& cellIdLabels
+            const UList<label>& cellIdLabels = UList<label>::null()
         );
 
+        //- remap action on triangulation or cleanup
+        virtual void remapFaces(const UList<label>& faceMap);
 
 public:
 
     // Constructors
 
-        //- Construct from components: Mesh reference and plane
-        cuttingPlane(const primitiveMesh&, const plane&);
-
         //- Construct from mesh reference and plane,
-        //  restricted to a list of cells
+        //  possibly restricted to a list of cells
         cuttingPlane
         (
             const primitiveMesh&,
             const plane&,
-            const labelList& cellIdLabels
+            const UList<label>& cellIdLabels = UList<label>::null()
         );
 
 
     // Member Functions
 
         //- Return plane used
-        const plane& planeDesc() const;
+        const plane& planeDesc() const
+        {
+            return static_cast<const plane&>(*this);
+        }
 
-        //- Return pointField of cutting points
+        //- Return List of cells cut by the plane
-        const pointField& points() const;
+        const labelList& cutCells() const
-
+        {
-        //- Return faceList of points in cells
+            return cutCells_;
-        const faceList& faces() const;
+        }
-
-        //- Return labelList of cut cells
-        const labelList& cells() const;
 
         //- Return true or false to question: have any cells been cut?
-        bool cut();
+        bool cut() const
+        {
+            return (cutCells_.size() != 0);
+        }
 
         //- Sample the cell field
         template<class Type>
@@ -185,7 +183,7 @@ public:
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #ifdef NoRepository
-#   include "cuttingPlaneSample.C"
+#   include "cuttingPlaneTemplates.C"
 #endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/sampling/cuttingPlane/cuttingPlaneTemplates.C

@@ -39,7 +39,7 @@ Foam::tmp<Foam::Field<Type> > Foam::cuttingPlane::sample
     const Field<Type>& sf
 ) const
 {
-    return tmp<Field<Type> >(new Field<Type>(sf, cells()));
+    return tmp<Field<Type> >(new Field<Type>(sf, cutCells()));
 }
 
 

src/sampling/sampledSet/writers/gnuplot/gnuplotSetWriter.C

@@ -24,7 +24,7 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "gnuplot.H"
+#include "gnuplotSetWriter.H"
 #include "clock.H"
 #include "coordSet.H"
 #include "fileName.H"
@@ -34,9 +34,8 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-// Construct from components
 template<class Type>
-Foam::gnuplot<Type>::gnuplot()
+Foam::gnuplotSetWriter<Type>::gnuplotSetWriter()
 :
     writer<Type>()
 {}
@@ -44,14 +43,14 @@ Foam::gnuplot<Type>::gnuplot()
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 template<class Type>
-Foam::gnuplot<Type>::~gnuplot()
+Foam::gnuplotSetWriter<Type>::~gnuplotSetWriter()
 {}
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type>
-Foam::fileName Foam::gnuplot<Type>::getFileName
+Foam::fileName Foam::gnuplotSetWriter<Type>::getFileName
 (
     const coordSet& points,
     const wordList& valueSetNames
@@ -62,7 +61,7 @@ Foam::fileName Foam::gnuplot<Type>::getFileName
 
 
 template<class Type>
-void Foam::gnuplot<Type>::write
+void Foam::gnuplotSetWriter<Type>::write
 (
     const coordSet& points,
     const wordList& valueSetNames,

src/sampling/sampledSet/writers/gnuplot/gnuplotSetWriter.H

@@ -23,17 +23,17 @@ License
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
 Class
-    Foam::gnuplot
+    Foam::gnuplotSetWriter
 
 Description
 
 SourceFiles
-    gnuplot.C
+    gnuplotSetWriter.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef gnuplot_H
+#ifndef gnuplotSetWriter_H
-#define gnuplot_H
+#define gnuplotSetWriter_H
 
 #include "writer.H"
 
@@ -43,11 +43,11 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-                           Class gnuplot Declaration
+                     Class gnuplotSetWriter Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class Type>
-class gnuplot
+class gnuplotSetWriter
 :
     public writer<Type>
 {
@@ -61,11 +61,11 @@ public:
     // Constructors
 
         //- Construct null
-        gnuplot();
+        gnuplotSetWriter();
 
 
     //- Destructor
-    virtual ~gnuplot();
+    virtual ~gnuplotSetWriter();
 
 
     // Member Functions
@@ -93,7 +93,7 @@ public:
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #ifdef NoRepository
-#   include "gnuplot.C"
+#   include "gnuplotSetWriter.C"
 #endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/sampling/sampledSet/writers/gnuplot/gnuplotSetWriterRunTime.C

@@ -24,14 +24,15 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "xmgrWriters.H"
+#include "gnuplotSetWriter.H"
+#include "writers.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
-    makeWriters(xmgr)
+    makeSetWriters(gnuplotSetWriter)
 }
 
 // ************************************************************************* //

src/sampling/sampledSet/writers/gnuplot/gnuplotWriters.H

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

src/sampling/sampledSet/writers/jplot/jplotSetWriter.C

@@ -24,7 +24,7 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "jplot.H"
+#include "jplotSetWriter.H"
 #include "clock.H"
 #include "coordSet.H"
 #include "fileName.H"
@@ -34,21 +34,19 @@ License
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 template<class Type>
-Foam::Ostream& Foam::jplot<Type>::writeHeader(Ostream& os) const
+Foam::Ostream& Foam::jplotSetWriter<Type>::writeHeader(Ostream& os) const
 {
     return os
-        << "# JPlot input file" << endl
+        << "# JPlot input file" << nl
-        << "#" << endl
+        << "#" << nl << nl
-        << endl
         << "# Generated by sample on " << clock::date().c_str() << endl;
 }
 
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-// Construct from components
 template<class Type>
-Foam::jplot<Type>::jplot()
+Foam::jplotSetWriter<Type>::jplotSetWriter()
 :
     writer<Type>()
 {}
@@ -57,14 +55,14 @@ Foam::jplot<Type>::jplot()
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 template<class Type>
-Foam::jplot<Type>::~jplot()
+Foam::jplotSetWriter<Type>::~jplotSetWriter()
 {}
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type>
-Foam::fileName Foam::jplot<Type>::getFileName
+Foam::fileName Foam::jplotSetWriter<Type>::getFileName
 (
     const coordSet& points,
     const wordList& valueSetNames
@@ -75,7 +73,7 @@ Foam::fileName Foam::jplot<Type>::getFileName
 
 
 template<class Type>
-void Foam::jplot<Type>::write
+void Foam::jplotSetWriter<Type>::write
 (
     const coordSet& points,
     const wordList& valueSetNames,
@@ -83,7 +81,7 @@ void Foam::jplot<Type>::write
     Ostream& os
 ) const
 {
-    os  << "# JPlot file" << endl
+    os  << "# JPlot file" << nl
         << "# column 1: " << points.name() << endl;
 
     forAll(valueSets, i)

src/sampling/sampledSet/writers/jplot/jplotSetWriter.H

@@ -23,15 +23,15 @@ License
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
 Class
-    Foam::jplot
+    Foam::jplotSetWriter
 
 SourceFiles
-    jplot.C
+    jplotSetWriter.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef jplot_H
+#ifndef jplotSetWriter_H
-#define jplot_H
+#define jplotSetWriter_H
 
 #include "writer.H"
 #include "vector.H"
@@ -42,11 +42,11 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-                           Class jplot Declaration
+                      Class jplotSetWriter Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class Type>
-class jplot
+class jplotSetWriter
 :
     public writer<Type>
 {
@@ -65,11 +65,11 @@ public:
     // Constructors
 
         //- Construct null
-        jplot();
+        jplotSetWriter();
 
 
     //- Destructor
-    virtual ~jplot();
+    virtual ~jplotSetWriter();
 
 
     // Member Functions
@@ -97,7 +97,7 @@ public:
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #ifdef NoRepository
-#   include "jplot.C"
+#   include "jplotSetWriter.C"
 #endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/sampling/sampledSet/writers/jplot/jplotSetWriterRunTime.C

@@ -24,14 +24,15 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "rawWriters.H"
+#include "jplotSetWriter.H"
+#include "writers.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
-    makeWriters(raw)
+    makeSetWriters(jplotSetWriter)
 }
 
 // ************************************************************************* //

src/sampling/sampledSet/writers/jplot/jplotWriters.H

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

src/sampling/sampledSet/writers/raw/rawSetWriter.C

@@ -24,16 +24,15 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "raw.H"
+#include "rawSetWriter.H"
 #include "coordSet.H"
 #include "fileName.H"
 #include "OFstream.H"
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-// Construct from components
 template<class Type>
-Foam::raw<Type>::raw()
+Foam::rawSetWriter<Type>::rawSetWriter()
 :
     writer<Type>()
 {}
@@ -42,14 +41,14 @@ Foam::raw<Type>::raw()
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 template<class Type>
-Foam::raw<Type>::~raw()
+Foam::rawSetWriter<Type>::~rawSetWriter()
 {}
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template<class Type>
-Foam::fileName Foam::raw<Type>::getFileName
+Foam::fileName Foam::rawSetWriter<Type>::getFileName
 (
     const coordSet& points,
     const wordList& valueSetNames
@@ -60,7 +59,7 @@ Foam::fileName Foam::raw<Type>::getFileName
 
 
 template<class Type>
-void Foam::raw<Type>::write
+void Foam::rawSetWriter<Type>::write
 (
     const coordSet& points,
     const wordList& valueSetNames,

src/sampling/sampledSet/writers/raw/rawSetWriter.H

@@ -23,17 +23,17 @@ License
     Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
 
 Class
-    Foam::raw
+    Foam::rawSetWriter
 
 Description
 
 SourceFiles
-    raw.C
+    rawSetWriter.C
 
 \*---------------------------------------------------------------------------*/
 
-#ifndef raw_H
+#ifndef rawSetWriter_H
-#define raw_H
+#define rawSetWriter_H
 
 #include "writer.H"
 
@@ -43,11 +43,11 @@ namespace Foam
 {
 
 /*---------------------------------------------------------------------------*\
-                           Class raw Declaration
+                       Class rawSetWriter Declaration
 \*---------------------------------------------------------------------------*/
 
 template<class Type>
-class raw
+class rawSetWriter
 :
     public writer<Type>
 {
@@ -61,11 +61,11 @@ public:
     // Constructors
 
         //- Construct null
-        raw();
+        rawSetWriter();
 
 
     //- Destructor
-    virtual ~raw();
+    virtual ~rawSetWriter();
 
 
     // Member Functions
@@ -93,7 +93,7 @@ public:
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 #ifdef NoRepository
-#   include "raw.C"
+#   include "rawSetWriter.C"
 #endif
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //

src/sampling/sampledSet/writers/raw/rawSetWriterRunTime.C

@@ -24,14 +24,15 @@ License
 
 \*---------------------------------------------------------------------------*/
 
-#include "jplotWriters.H"
+#include "rawSetWriter.H"
+#include "writers.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 namespace Foam
 {
-    makeWriters(jplot)
+    makeSetWriters(rawSetWriter)
 }
 
 // ************************************************************************* //