diff --git a/src/OpenFOAM/Make/files b/src/OpenFOAM/Make/files
index b67459678..9c2ce22ba 100644
--- a/src/OpenFOAM/Make/files
+++ b/src/OpenFOAM/Make/files
@@ -414,6 +414,7 @@ $(primitiveMesh)/primitiveMeshFaceCentresAndAreas.C
 $(primitiveMesh)/primitiveMeshFindCell.C
 $(primitiveMesh)/primitiveMeshPointCells.C
 $(primitiveMesh)/primitiveMeshPointFaces.C
+$(primitiveMesh)/primitiveMeshPointEdges.C
 $(primitiveMesh)/primitiveMeshPointPoints.C
 $(primitiveMesh)/primitiveMeshCellPoints.C
 $(primitiveMesh)/primitiveMeshCalcCellShapes.C
diff --git a/src/OpenFOAM/coordinateSystems/coordinateRotation/EulerCoordinateRotation.C b/src/OpenFOAM/coordinateSystems/coordinateRotation/EulerCoordinateRotation.C
index 76ca321ae..99037d65b 100644
--- a/src/OpenFOAM/coordinateSystems/coordinateRotation/EulerCoordinateRotation.C
+++ b/src/OpenFOAM/coordinateSystems/coordinateRotation/EulerCoordinateRotation.C
@@ -28,6 +28,7 @@ License
 
 #include "Switch.H"
 #include "mathematicalConstants.H"
+#include "dictionary.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/coordinateSystems/coordinateRotation/STARCDCoordinateRotation.C b/src/OpenFOAM/coordinateSystems/coordinateRotation/STARCDCoordinateRotation.C
index 8a6d42914..524eb18c5 100644
--- a/src/OpenFOAM/coordinateSystems/coordinateRotation/STARCDCoordinateRotation.C
+++ b/src/OpenFOAM/coordinateSystems/coordinateRotation/STARCDCoordinateRotation.C
@@ -28,6 +28,7 @@ License
 
 #include "Switch.H"
 #include "mathematicalConstants.H"
+#include "dictionary.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/coordinateSystems/coordinateRotation/axisCoordinateRotation.C b/src/OpenFOAM/coordinateSystems/coordinateRotation/axisCoordinateRotation.C
index 9932acaa6..01ed4bdcb 100644
--- a/src/OpenFOAM/coordinateSystems/coordinateRotation/axisCoordinateRotation.C
+++ b/src/OpenFOAM/coordinateSystems/coordinateRotation/axisCoordinateRotation.C
@@ -28,6 +28,7 @@ License
 
 #include "Switch.H"
 #include "mathematicalConstants.H"
+#include "dictionary.H"
 #include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/coordinateSystems/cylindricalCS.C b/src/OpenFOAM/coordinateSystems/cylindricalCS.C
index 82c4989ff..468bca0a9 100644
--- a/src/OpenFOAM/coordinateSystems/cylindricalCS.C
+++ b/src/OpenFOAM/coordinateSystems/cylindricalCS.C
@@ -26,7 +26,6 @@ License
 
 #include "cylindricalCS.H"
 
-#include "one.H"
 #include "Switch.H"
 #include "mathematicalConstants.H"
 #include "addToRunTimeSelectionTable.H"
@@ -133,7 +132,7 @@ Foam::vector Foam::cylindricalCS::localToGlobal
 {
     scalar theta
     (
-        local.y() * ( inDegrees_ ? mathematicalConstant::pi/180.0 : 1.0 )  
+        local.y() * ( inDegrees_ ? mathematicalConstant::pi/180.0 : 1.0 )
     );
 
     return coordinateSystem::localToGlobal
diff --git a/src/OpenFOAM/coordinateSystems/cylindricalCS.H b/src/OpenFOAM/coordinateSystems/cylindricalCS.H
index b84acd093..d8e2eda9f 100644
--- a/src/OpenFOAM/coordinateSystems/cylindricalCS.H
+++ b/src/OpenFOAM/coordinateSystems/cylindricalCS.H
@@ -95,13 +95,13 @@ public:
     // Constructors
 
         //- Construct null
-        cylindricalCS(const bool inDegrees=true);
+        cylindricalCS(const bool inDegrees = true);
 
         //- Construct copy
         cylindricalCS
         (
             const coordinateSystem&,
-            const bool inDegrees=true
+            const bool inDegrees = true
         );
 
         //- Construct copy with a different name
@@ -109,7 +109,7 @@ public:
         (
             const word& name,
             const coordinateSystem&,
-            const bool inDegrees=true
+            const bool inDegrees = true
         );
 
         //- Construct from origin and rotation
@@ -118,7 +118,7 @@ public:
             const word& name,
             const point& origin,
             const coordinateRotation&,
-            const bool inDegrees=true
+            const bool inDegrees = true
         );
 
         //- Construct from origin and 2 axes
@@ -128,7 +128,7 @@ public:
             const point& origin,
             const vector& axis,
             const vector& dirn,
-            const bool inDegrees=true
+            const bool inDegrees = true
         );
 
         //- Construct from dictionary
diff --git a/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.C b/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.C
index 16d531c8c..d5a668184 100644
--- a/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.C
+++ b/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.C
@@ -25,8 +25,10 @@ License
 \*---------------------------------------------------------------------------*/
 
 #include "ellipticCylindricalCS.H"
-#include "addToRunTimeSelectionTable.H"
+
+#include "Switch.H"
 #include "mathematicalConstants.H"
+#include "addToRunTimeSelectionTable.H"
 
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
@@ -59,6 +61,7 @@ Foam::ellipticCylindricalCS::ellipticCylindricalCS
 )
 :
     coordinateSystem(cs),
+    a_(0),
     inDegrees_(inDegrees)
 {}
 
@@ -71,6 +74,22 @@ Foam::ellipticCylindricalCS::ellipticCylindricalCS
 )
 :
     coordinateSystem(name, cs),
+    a_(0),
+    inDegrees_(inDegrees)
+{}
+
+
+Foam::ellipticCylindricalCS::ellipticCylindricalCS
+(
+    const word& name,
+    const point& origin,
+    const coordinateRotation& cr,
+    const scalar a,
+    const bool inDegrees
+)
+:
+    coordinateSystem(name, origin, cr),
+    a_(a),
     inDegrees_(inDegrees)
 {}
 
@@ -91,21 +110,6 @@ Foam::ellipticCylindricalCS::ellipticCylindricalCS
 {}
 
 
-Foam::ellipticCylindricalCS::ellipticCylindricalCS
-(
-    const word& name,
-    const point& origin,
-    const coordinateRotation& cr,
-    const scalar a,
-    const bool inDegrees
-)
-:
-    coordinateSystem(name, origin, cr),
-    a_(a),
-    inDegrees_(inDegrees)
-{}
-
-
 Foam::ellipticCylindricalCS::ellipticCylindricalCS
 (
     const word& name,
@@ -128,8 +132,7 @@ Foam::vector Foam::ellipticCylindricalCS::localToGlobal
 {
     // Notation: u = local.x() v = local.y() z = local.z();
     scalar theta =
-        local.y()
-       *( inDegrees_ ? mathematicalConstant::pi/180.0 : 1.0 )
+        local.y()*( inDegrees_ ? mathematicalConstant::pi/180.0 : 1.0 );
 
     return coordinateSystem::localToGlobal
     (
@@ -150,8 +153,8 @@ Foam::tmp<Foam::vectorField> Foam::ellipticCylindricalCS::localToGlobal
 ) const
 {
     scalarField theta =
-        local.component(vector::Y)
-       *( inDegrees_ ? mathematicalConstant::pi/180.0 : 1.0 )
+        local.component(vector::Y)*
+        ( inDegrees_ ? mathematicalConstant::pi/180.0 : 1.0 );
 
     vectorField lc(local.size());
     lc.replace
diff --git a/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.H b/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.H
index 8b4c673fa..4fc4fea2d 100644
--- a/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.H
+++ b/src/OpenFOAM/coordinateSystems/ellipticCylindricalCS.H
@@ -108,14 +108,14 @@ public:
         ellipticCylindricalCS(const bool inDegrees = true);
 
         //- Construct copy
-        sphericalCS
+        ellipticCylindricalCS
         (
             const coordinateSystem&,
             const bool inDegrees = true
         );
 
         //- Construct copy with a different name
-        sphericalCS
+        ellipticCylindricalCS
         (
             const word& name,
             const coordinateSystem&,
@@ -132,6 +132,17 @@ public:
             const bool inDegrees = true
         );
 
+        //- Construct from origin and 2 axes
+        ellipticCylindricalCS
+        (
+            const word& name,
+            const point& origin,
+            const vector& axis,
+            const vector& dirn,
+            const scalar a,
+            const bool inDegrees = true
+        );
+
         //- Construct from dictionary
         ellipticCylindricalCS
         (
diff --git a/src/OpenFOAM/db/IOobject/IOobjectI.H b/src/OpenFOAM/db/IOobject/IOobjectI.H
index f1de065c7..a8e3cb822 100644
--- a/src/OpenFOAM/db/IOobject/IOobjectI.H
+++ b/src/OpenFOAM/db/IOobject/IOobjectI.H
@@ -33,16 +33,20 @@ inline Stream& Foam::IOobject::writeBanner(Stream& os, bool noHint)
 {
     static bool spacesSet = false;
     static char spaces[40];
-    static char spacesRev[40];
 
     if (!spacesSet)
     {
         memset(spaces, ' ', 40);
-        spaces[38 - strlen(Foam::FOAMversion)] = '\0';
-
-        memset(spacesRev, ' ', 40);
-        spacesRev[38 - strlen(Foam::FOAMDevRevisionNumber)] = '\0';
 
+        size_t len = strlen(Foam::FOAMversion);
+        if (len < 38)
+        {
+            spaces[38 - len] = '\0';
+        }
+        else
+        {
+            spaces[0] = '\0';
+        }
         spacesSet = true;
     }
 
@@ -59,10 +63,10 @@ inline Stream& Foam::IOobject::writeBanner(Stream& os, bool noHint)
 
     os  <<
         "| =========                 |                                                 |\n"
-        "| \\\\      /  F ield         | OpenFOAM: The Open Source CFD Toolbox           |\n"
-        "|  \\\\    /   O peration     | Version:  " << FOAMversion           << spaces << "|\n"
-        "|   \\\\  /    A nd           | Revision: " << FOAMDevRevisionNumber << spacesRev << "|\n"
-        "|    \\\\/     M anipulation  | Web:             www.OpenFOAM.org               |\n"
+        "| \\\\      /  F ield         | OpenFOAM Extend Project: Open source CFD        |\n"
+        "|  \\\\    /   O peration     | Version:  " << FOAMversion << spaces << "|\n"
+        "|   \\\\  /    A nd           | Web:      www.extend-project.de                 |\n"
+        "|    \\\\/     M anipulation  |                                                 |\n"
         "\\*---------------------------------------------------------------------------*/\n";
 
     return os;
diff --git a/src/OpenFOAM/db/dlLibraryTable/dlLibraryTable.H b/src/OpenFOAM/db/dlLibraryTable/dlLibraryTable.H
index 4ddc223af..47f33c866 100644
--- a/src/OpenFOAM/db/dlLibraryTable/dlLibraryTable.H
+++ b/src/OpenFOAM/db/dlLibraryTable/dlLibraryTable.H
@@ -49,28 +49,6 @@ namespace Foam
                        Class dlLibraryTable Declaration
 \*---------------------------------------------------------------------------*/
 
-//- A means of hashing pointer addresses
-template<>
-class Hash<void*>
-{
-
-public:
-
-    Hash()
-    {}
-
-    long operator()(const void* const& p) const
-    {
-        return long(p);
-    }
-
-    label operator()(const void* const& p, const label tableSize) const
-    {
-        return abs(operator()(p)) % tableSize;
-    }
-};
-
-
 class dlLibraryTable
 :
     public HashTable<fileName, void*, Hash<void*> >
diff --git a/src/OpenFOAM/db/regIOobject/regIOobjectWrite.C b/src/OpenFOAM/db/regIOobject/regIOobjectWrite.C
index 5c6f1e16d..e803481f4 100644
--- a/src/OpenFOAM/db/regIOobject/regIOobjectWrite.C
+++ b/src/OpenFOAM/db/regIOobject/regIOobjectWrite.C
@@ -84,7 +84,15 @@ bool Foam::regIOobject::writeObject
 
     {
         // Try opening an OFstream for object
-        OFstream os(objectPath(), fmt, ver, cmp);
+        // Stream open for over-write.  HJ, 17/Aug/2010
+        OFstream os
+        (
+            objectPath(),
+            ios_base::out|ios_base::trunc,
+            fmt,
+            ver,
+            cmp
+        );
 
         // If any of these fail, return (leave error handling to Ostream class)
         if (!os.good())
diff --git a/src/OpenFOAM/fields/Fields/Field/SubFieldI.H b/src/OpenFOAM/fields/Fields/Field/SubFieldI.H
index 9992ded59..07d57a15b 100644
--- a/src/OpenFOAM/fields/Fields/Field/SubFieldI.H
+++ b/src/OpenFOAM/fields/Fields/Field/SubFieldI.H
@@ -42,7 +42,7 @@ inline Foam::SubField<Type>::SubField
     const UList<Type>& list
 )
 :
-q    SubList<Type>(list, list.size())
+    SubList<Type>(list, list.size())
 {}
 
 
diff --git a/src/OpenFOAM/matrices/Matrix/Matrix.C b/src/OpenFOAM/matrices/Matrix/Matrix.C
index 231862557..c45dd2cbb 100644
--- a/src/OpenFOAM/matrices/Matrix/Matrix.C
+++ b/src/OpenFOAM/matrices/Matrix/Matrix.C
@@ -28,13 +28,13 @@ License
 
 // * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * * //
 
-template<class T>
-void Foam::Matrix<T>::allocate()
+template<class Form, class Type>
+void Foam::Matrix<Form, Type>::allocate()
 {
     if (n_ && m_)
     {
-        v_ = new T*[n_];
-        v_[0] = new T[n_*m_];
+        v_ = new Type*[n_];
+        v_[0] = new Type[n_*m_];
 
         for (register label i=1; i<n_; i++)
         {
@@ -46,8 +46,8 @@ void Foam::Matrix<T>::allocate()
 
 // * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * * //
 
-template<class T>
-Foam::Matrix<T>::~Matrix()
+template<class Form, class Type>
+Foam::Matrix<Form, Type>::~Matrix()
 {
     if (v_)
     {
@@ -59,34 +59,8 @@ Foam::Matrix<T>::~Matrix()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-template<class T>
-const Foam::Matrix<T>& Foam::Matrix<T>::null()
-{
-    Matrix<T>* nullPtr = reinterpret_cast<Matrix<T>*>(NULL);
-    return *nullPtr;
-}
-
-
-template<class T>
-Foam::Matrix<T>::Matrix(const label n, const label m)
-:
-    v_(NULL),
-    n_(n),
-    m_(m)
-{
-    if (n_ < 0 || m_ < 0)
-    {
-        FatalErrorIn("Matrix<T>::Matrix(const label n, const label m)")
-            << "bad n, m " << n_ << ", " << m_
-            << abort(FatalError);
-    }
-
-    allocate();
-}
-
-
-template<class T>
-Foam::Matrix<T>::Matrix(const label n, const label m, const T& a)
+template<class Form, class Type>
+Foam::Matrix<Form, Type>::Matrix(const label n, const label m)
 :
     v_(NULL),
     n_(n),
@@ -96,7 +70,28 @@ Foam::Matrix<T>::Matrix(const label n, const label m, const T& a)
     {
         FatalErrorIn
         (
-            "Matrix<T>::Matrix(const label n, const label m, const T&)"
+            "Matrix<Form, Type>::Matrix(const label n, const label m)"
+        )   << "bad n, m " << n_ << ", " << m_
+            << abort(FatalError);
+    }
+
+    allocate();
+}
+
+
+template<class Form, class Type>
+Foam::Matrix<Form, Type>::Matrix(const label n, const label m, const Type& a)
+:
+    v_(NULL),
+    n_(n),
+    m_(m)
+{
+    if (n_ < 0 || m_ < 0)
+    {
+        FatalErrorIn
+        (
+            "Matrix<Form, Type>::Matrix"
+            "(const label n, const label m, const T&)"
         )   << "bad n, m " << n_ << ", " << m_
             << abort(FatalError);
     }
@@ -105,7 +100,7 @@ Foam::Matrix<T>::Matrix(const label n, const label m, const T& a)
 
     if (v_)
     {
-        T* v = v_[0];
+        Type* v = v_[0];
 
         label nm = n_*m_;
 
@@ -117,8 +112,8 @@ Foam::Matrix<T>::Matrix(const label n, const label m, const T& a)
 }
 
 
-template<class T>
-Foam::Matrix<T>::Matrix(const Matrix<T>& a)
+template<class Form, class Type>
+Foam::Matrix<Form, Type>::Matrix(const Matrix<Form, Type>& a)
 :
     v_(NULL),
     n_(a.n_),
@@ -127,8 +122,8 @@ Foam::Matrix<T>::Matrix(const Matrix<T>& a)
     if (a.v_)
     {
         allocate();
-        T* v = v_[0];
-        const T* av = a.v_[0];
+        Type* v = v_[0];
+        const Type* av = a.v_[0];
 
         label nm = n_*m_;
         for (register label i=0; i<nm; i++)
@@ -138,9 +133,9 @@ Foam::Matrix<T>::Matrix(const Matrix<T>& a)
     }
 }
 
-        
-template<class T>
-void Foam::Matrix<T>::clear()
+
+template<class Form, class Type>
+void Foam::Matrix<Form, Type>::clear()
 {
     if (v_)
     {
@@ -153,8 +148,8 @@ void Foam::Matrix<T>::clear()
 }
 
 
-template<class T>
-void Foam::Matrix<T>::transfer(Matrix<T>& a)
+template<class Form, class Type>
+void Foam::Matrix<Form, Type>::transfer(Matrix<Form, Type>& a)
 {
     clear();
 
@@ -169,14 +164,32 @@ void Foam::Matrix<T>::transfer(Matrix<T>& a)
 }
 
 
+template<class Form, class Type>
+Form Foam::Matrix<Form, Type>::T() const
+{
+    const Matrix<Form, Type>& A = *this;
+    Form At(m(), n());
+
+    for (register label i=0; i<n(); i++)
+    {
+        for (register label j=0; j<m(); j++)
+        {
+            At[j][i] = A[i][j];
+        }
+    }
+
+    return At;
+}
+
+
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
-template<class T>
-void Foam::Matrix<T>::operator=(const T& t)
+template<class Form, class Type>
+void Foam::Matrix<Form, Type>::operator=(const Type& t)
 {
     if (v_)
     {
-        T* v = v_[0];
+        Type* v = v_[0];
 
         label nm = n_*m_;
         for (register label i=0; i<nm; i++)
@@ -188,12 +201,12 @@ void Foam::Matrix<T>::operator=(const T& t)
 
 
 // Assignment operator. Takes linear time.
-template<class T>
-void Foam::Matrix<T>::operator=(const Matrix<T>& a)
+template<class Form, class Type>
+void Foam::Matrix<Form, Type>::operator=(const Matrix<Form, Type>& a)
 {
     if (this == &a)
     {
-        FatalErrorIn("Matrix<T>::operator=(const Matrix<T>&)")
+        FatalErrorIn("Matrix<Form, Type>::operator=(const Matrix<Form, Type>&)")
             << "attempted assignment to self"
             << abort(FatalError);
     }
@@ -204,12 +217,12 @@ void Foam::Matrix<T>::operator=(const Matrix<T>& a)
         n_ = a.n_;
         m_ = a.m_;
         allocate();
-    } 
+    }
 
     if (v_)
     {
-        T* v = v_[0];
-        const T* av = a.v_[0];
+        Type* v = v_[0];
+        const Type* av = a.v_[0];
 
         label nm = n_*m_;
         for (register label i=0; i<nm; i++)
@@ -222,15 +235,15 @@ void Foam::Matrix<T>::operator=(const Matrix<T>& a)
 
 // * * * * * * * * * * * * * * * Global Functions  * * * * * * * * * * * * * //
 
-template<class T>
-const T& Foam::max(const Matrix<T>& a)
+template<class Form, class Type>
+const Type& Foam::max(const Matrix<Form, Type>& a)
 {
-    label nm = a.n_*a.m_;
+    label nm = a.n()*a.m();
 
     if (nm)
     {
         label curMaxI = 0;
-        const T* v = a.v_[0];
+        const Type* v = a[0];
 
         for (register label i=1; i<nm; i++)
         {
@@ -244,7 +257,7 @@ const T& Foam::max(const Matrix<T>& a)
     }
     else
     {
-        FatalErrorIn("max(const Matrix<T>&)")
+        FatalErrorIn("max(const Matrix<Form, Type>&)")
             << "matrix is empty"
             << abort(FatalError);
 
@@ -254,15 +267,15 @@ const T& Foam::max(const Matrix<T>& a)
 }
 
 
-template<class T>
-const T& Foam::min(const Matrix<T>& a)
+template<class Form, class Type>
+const Type& Foam::min(const Matrix<Form, Type>& a)
 {
-    label nm = a.n_*a.m_;
+    label nm = a.n()*a.m();
 
     if (nm)
     {
         label curMinI = 0;
-        const T* v = a.v_[0];
+        const Type* v = a[0];
 
         for (register label i=1; i<nm; i++)
         {
@@ -276,7 +289,7 @@ const T& Foam::min(const Matrix<T>& a)
     }
     else
     {
-        FatalErrorIn("min(const Matrix<T>&)")
+        FatalErrorIn("min(const Matrix<Form, Type>&)")
             << "matrix is empty"
             << abort(FatalError);
 
@@ -288,17 +301,17 @@ const T& Foam::min(const Matrix<T>& a)
 
 // * * * * * * * * * * * * * * * Global Operators  * * * * * * * * * * * * * //
 
-template<class T>
-Foam::Matrix<T> Foam::operator-(const Matrix<T>& a)
+template<class Form, class Type>
+Form Foam::operator-(const Matrix<Form, Type>& a)
 {
-    Matrix<T> na(a.n_, a.m_);
+    Form na(a.n(), a.m());
 
-    if (a.n_ && a.m_)
+    if (a.n() && a.m())
     {
-        T* nav = na.v_[0];
-        const T* av = a.v_[0];
+        Type* nav = na[0];
+        const Type* av = a[0];
 
-        label nm = a.n_*a.m_;
+        label nm = a.n()*a.m();
         for (register label i=0; i<nm; i++)
         {
             nav[i] = -av[i];
@@ -309,32 +322,36 @@ Foam::Matrix<T> Foam::operator-(const Matrix<T>& a)
 }
 
 
-template<class T>
-Foam::Matrix<T> Foam::operator+(const Matrix<T>& a, const Matrix<T>& b)
+template<class Form, class Type>
+Form Foam::operator+(const Matrix<Form, Type>& a, const Matrix<Form, Type>& b)
 {
-    if (a.n_ != b.n_)
+    if (a.n() != b.n())
     {
-        FatalErrorIn("Matrix<T>::operator+(const Matrix<T>&, const Matrix<T>&)")
-            << "attempted add matrices with different number of rows: "
-            << a.n_ << ", " << b.n_
+        FatalErrorIn
+        (
+            "Matrix<Form, Type>::operator+(const Matrix<Form, Type>&, const Matrix<Form, Type>&)"
+        )   << "attempted add matrices with different number of rows: "
+            << a.n() << ", " << b.n()
             << abort(FatalError);
     }
 
-    if (a.m_ != b.m_)
+    if (a.m() != b.m())
     {
-        FatalErrorIn("Matrix<T>::operator+(const Matrix<T>&, const Matrix<T>&)")
-            << "attempted add matrices with different number of columns: "
-            << a.m_ << ", " << b.m_
+        FatalErrorIn
+        (
+            "Matrix<Form, Type>::operator+(const Matrix<Form, Type>&, const Matrix<Form, Type>&)"
+        )   << "attempted add matrices with different number of columns: "
+            << a.m() << ", " << b.m()
             << abort(FatalError);
     }
 
-    Matrix<T> ab(a.n_, a.m_);
+    Form ab(a.n(), a.m());
 
-    T* abv = ab.v_[0];
-    const T* av = a.v_[0];
-    const T* bv = b.v_[0];
+    Type* abv = ab[0];
+    const Type* av = a[0];
+    const Type* bv = b[0];
 
-    label nm = a.n_*a.m_;;
+    label nm = a.n()*a.m();
     for (register label i=0; i<nm; i++)
     {
         abv[i] = av[i] + bv[i];
@@ -344,32 +361,36 @@ Foam::Matrix<T> Foam::operator+(const Matrix<T>& a, const Matrix<T>& b)
 }
 
 
-template<class T>
-Foam::Matrix<T> Foam::operator-(const Matrix<T>& a, const Matrix<T>& b)
+template<class Form, class Type>
+Form Foam::operator-(const Matrix<Form, Type>& a, const Matrix<Form, Type>& b)
 {
-    if (a.n_ != b.n_)
+    if (a.n() != b.n())
     {
-        FatalErrorIn("Matrix<T>::operator-(const Matrix<T>&, const Matrix<T>&)")
-            << "attempted add matrices with different number of rows: "
-            << a.n_ << ", " << b.n_
+        FatalErrorIn
+        (
+            "Matrix<Form, Type>::operator-(const Matrix<Form, Type>&, const Matrix<Form, Type>&)"
+        )   << "attempted add matrices with different number of rows: "
+            << a.n() << ", " << b.n()
             << abort(FatalError);
     }
 
-    if (a.m_ != b.m_)
+    if (a.m() != b.m())
     {
-        FatalErrorIn("Matrix<T>::operator-(const Matrix<T>&, const Matrix<T>&)")
-            << "attempted add matrices with different number of columns: "
-            << a.m_ << ", " << b.m_
+        FatalErrorIn
+        (
+            "Matrix<Form, Type>::operator-(const Matrix<Form, Type>&, const Matrix<Form, Type>&)"
+        )   << "attempted add matrices with different number of columns: "
+            << a.m() << ", " << b.m()
             << abort(FatalError);
     }
 
-    Matrix<T> ab(a.n_, a.m_);
+    Form ab(a.n(), a.m());
 
-    T* abv = ab.v_[0];
-    const T* av = a.v_[0];
-    const T* bv = b.v_[0];
+    Type* abv = ab[0];
+    const Type* av = a[0];
+    const Type* bv = b[0];
 
-    label nm = a.n_*a.m_;;
+    label nm = a.n()*a.m();
     for (register label i=0; i<nm; i++)
     {
         abv[i] = av[i] - bv[i];
@@ -379,17 +400,17 @@ Foam::Matrix<T> Foam::operator-(const Matrix<T>& a, const Matrix<T>& b)
 }
 
 
-template<class T>
-Foam::Matrix<T> Foam::operator*(const scalar s, const Matrix<T>& a)
+template<class Form, class Type>
+Form Foam::operator*(const scalar s, const Matrix<Form, Type>& a)
 {
-    Matrix<T> sa(a.n_, a.m_);
+    Form sa(a.n(), a.m());
 
-    if (a.n_ && a.m_)
+    if (a.n() && a.m())
     {
-        T* sav = sa.v_[0];
-        const T* av = a.v_[0];
+        Type* sav = sa[0];
+        const Type* av = a[0];
 
-        label nm = a.n_*a.m_;;
+        label nm = a.n()*a.m();
         for (register label i=0; i<nm; i++)
         {
             sav[i] = s*av[i];
diff --git a/src/OpenFOAM/matrices/Matrix/Matrix.H b/src/OpenFOAM/matrices/Matrix/Matrix.H
index 5e1f64b48..694aa0b74 100644
--- a/src/OpenFOAM/matrices/Matrix/Matrix.H
+++ b/src/OpenFOAM/matrices/Matrix/Matrix.H
@@ -51,31 +51,32 @@ namespace Foam
 
 // Forward declaration of friend functions and operators
 
-template<class T> class Matrix;
+template<class Form, class Type> class Matrix;
 
-template<class T> const T& max(const Matrix<T>&);
-template<class T> const T& min(const Matrix<T>&);
+template<class Form, class Type> Istream& operator>>
+(
+    Istream&,
+    Matrix<Form, Type>&
+);
 
-template<class T> Matrix<T> operator-(const Matrix<T>&);
-template<class T> Matrix<T> operator+(const Matrix<T>&, const Matrix<T>&);
-template<class T> Matrix<T> operator-(const Matrix<T>&, const Matrix<T>&);
-template<class T> Matrix<T> operator*(const scalar, const Matrix<T>&);
-
-template<class T> Istream& operator>>(Istream&, Matrix<T>&);
-template<class T> Ostream& operator<<(Ostream&, const Matrix<T>&);
+template<class Form, class Type> Ostream& operator<<
+(
+    Ostream&,
+    const Matrix<Form, Type>&
+);
 
 
 /*---------------------------------------------------------------------------*\
                            Class Matrix Declaration
 \*---------------------------------------------------------------------------*/
 
-template<class T>
+template<class Form, class Type>
 class Matrix
 {
     // Private data
 
         //- Row pointers
-        T** __restrict__ v_;
+        Type** __restrict__ v_;
 
         //- Number of rows and columns in Matrix.
         label n_, m_;
@@ -87,6 +88,12 @@ class Matrix
 
 public:
 
+    // Static Member Functions
+
+        //- Return a null Matrix
+        inline static const Matrix<Form, Type>& null();
+
+
     // Constructors
 
         //- Null constructor.
@@ -97,16 +104,16 @@ public:
 
         //- Construct with given number of rows and columns
         //  and value for all elements.
-        Matrix(const label n, const label m, const T&);
+        Matrix(const label n, const label m, const Type&);
 
         //- Copy constructor.
-        Matrix(const Matrix<T>&);
+    Matrix(const Matrix<Form, Type>&);
 
         //- Construct from Istream.
         Matrix(Istream&);
 
         //- Clone
-        inline autoPtr<Matrix<T> > clone() const;
+        inline autoPtr<Matrix<Form, Type> > clone() const;
 
 
     // Destructor
@@ -116,10 +123,6 @@ public:
 
     // Member functions
 
-        //- Return a null Matrix
-        static const Matrix<T>& null();
-
-
         // Access
 
             //- Return the number of rows
@@ -148,48 +151,72 @@ public:
 
             //- Transfer the contents of the argument Matrix into this Matrix
             //  and annull the argument Matrix.
-            void transfer(Matrix<T>&);
+            void transfer(Matrix<Form, Type>&);
+
+
+        //- Return the transpose of the matrix
+        Form T() const;
 
 
     // Member operators
 
         //- Return subscript-checked element of Matrix.
-        inline T* operator[](const label);
+        inline Type* operator[](const label);
 
         //- Return subscript-checked element of constant Matrix.
-        inline const T* operator[](const label) const;
+        inline const Type* operator[](const label) const;
 
         //- Assignment operator. Takes linear time.
-        void operator=(const Matrix<T>&);
+        void operator=(const Matrix<Form, Type>&);
 
         //- Assignment of all entries to the given value
-        void operator=(const T&);
-
-
-    // Friend functions
-
-        friend const T& max<T>(const Matrix<T>&);
-        friend const T& min<T>(const Matrix<T>&);
-
-
-    // Friend operators
-
-        friend Matrix<T> operator-<T>(const Matrix<T>&);
-        friend Matrix<T> operator+<T>(const Matrix<T>&, const Matrix<T>&);
-        friend Matrix<T> operator-<T>(const Matrix<T>&, const Matrix<T>&);
-        friend Matrix<T> operator*<T>(const scalar, const Matrix<T>&);
+        void operator=(const Type&);
 
 
     // IOstream operators
 
         //- Read Matrix from Istream, discarding contents of existing Matrix.
-        friend Istream& operator>> <T>(Istream&, Matrix<T>&);
+        friend Istream& operator>> <Form, Type>
+        (
+            Istream&,
+            Matrix<Form, Type>&
+        );
 
         // Write Matrix to Ostream.
-        friend Ostream& operator<< <T>(Ostream&, const Matrix<T>&);
+        friend Ostream& operator<< <Form, Type>
+        (
+            Ostream&,
+            const Matrix<Form, Type>&
+        );
 };
 
 
+// Global functions and operators
+
+template<class Form, class Type> const Type& max(const Matrix<Form, Type>&);
+template<class Form, class Type> const Type& min(const Matrix<Form, Type>&);
+
+template<class Form, class Type> Form operator-(const Matrix<Form, Type>&);
+
+template<class Form, class Type> Form operator+
+(
+    const Matrix<Form, Type>&,
+    const Matrix<Form, Type>&
+);
+
+template<class Form, class Type> Form operator-
+(
+    const Matrix<Form, Type>&,
+    const Matrix<Form, Type>&
+);
+
+template<class Form, class Type> Form operator*
+(
+    const scalar,
+    const Matrix<Form, Type>&
+);
+
+
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
 } // End namespace Foam
diff --git a/src/OpenFOAM/matrices/Matrix/MatrixI.H b/src/OpenFOAM/matrices/Matrix/MatrixI.H
index acb2d8928..4273a8d24 100644
--- a/src/OpenFOAM/matrices/Matrix/MatrixI.H
+++ b/src/OpenFOAM/matrices/Matrix/MatrixI.H
@@ -26,8 +26,8 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
-template<class T>
-inline Foam::Matrix<T>::Matrix()
+template<class Form, class Type>
+inline Foam::Matrix<Form, Type>::Matrix()
 :
     v_(NULL),
     n_(0),
@@ -35,71 +35,75 @@ inline Foam::Matrix<T>::Matrix()
 {}
 
 
-template<class T>
-inline Foam::autoPtr<Foam::Matrix<T> > Foam::Matrix<T>::clone() const
+template<class Form, class Type>
+inline Foam::autoPtr<Foam::Matrix<Form, Type> >
+Foam::Matrix<Form, Type>::clone() const
 {
-    return autoPtr<Matrix<T> >(new Matrix<T>(*this));
+    return autoPtr<Matrix<Form, Type> >(new Matrix<Form, Type>(*this));
 }
 
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
+template<class Form, class Type>
+inline const Foam::Matrix<Form, Type>& Foam::Matrix<Form, Type>::null()
+{
+    return *reinterpret_cast< Matrix<Form, Type>* >(0);
+}
+
+
 //- Return the number of rows
-template<class T>
-inline Foam::label Foam::Matrix<T>::n() const
+template<class Form, class Type>
+inline Foam::label Foam::Matrix<Form, Type>::n() const
 {
     return n_;
 }
 
 
-//- Return the number of columns
-template<class T>
-inline Foam::label Foam::Matrix<T>::m() const
+template<class Form, class Type>
+inline Foam::label Foam::Matrix<Form, Type>::m() const
 {
     return m_;
 }
 
 
-//- Return the number of columns
-template<class T>
-inline Foam::label Foam::Matrix<T>::size() const
+template<class Form, class Type>
+inline Foam::label Foam::Matrix<Form, Type>::size() const
 {
     return n_*m_;
 }
 
 
-// Check index i is within valid range (0 ... n-1).
-template<class T>
-inline void Foam::Matrix<T>::checki(const label i) const
+template<class Form, class Type>
+inline void Foam::Matrix<Form, Type>::checki(const label i) const
 {
     if (!n_)
     {
-        FatalErrorIn("Matrix<T>::checki(const label)")
+        FatalErrorIn("Matrix<Form, Type>::checki(const label)")
             << "attempt to access element from zero sized row"
             << abort(FatalError);
     }
     else if (i<0 || i>=n_)
     {
-        FatalErrorIn("Matrix<T>::checki(const label)")
+        FatalErrorIn("Matrix<Form, Type>::checki(const label)")
             << "index " << i << " out of range 0 ... " << n_-1
             << abort(FatalError);
     }
 }
 
 
-// Check index j is within valid range (0 ... n-1).
-template<class T>
-inline void Foam::Matrix<T>::checkj(const label j) const
+template<class Form, class Type>
+inline void Foam::Matrix<Form, Type>::checkj(const label j) const
 {
     if (!m_)
     {
-        FatalErrorIn("Matrix<T>::checkj(const label)")
+        FatalErrorIn("Matrix<Form, Type>::checkj(const label)")
             << "attempt to access element from zero sized column"
             << abort(FatalError);
     }
     else if (j<0 || j>=m_)
     {
-        FatalErrorIn("Matrix<T>::checkj(const label)")
+        FatalErrorIn("Matrix<Form, Type>::checkj(const label)")
             << "index " << j << " out of range 0 ... " << m_-1
             << abort(FatalError);
     }
@@ -108,9 +112,8 @@ inline void Foam::Matrix<T>::checkj(const label j) const
 
 // * * * * * * * * * * * * * * * Member Operators  * * * * * * * * * * * * * //
 
-// Return subscript-checked element access
-template<class T>
-inline T* Foam::Matrix<T>::operator[](const label i)
+template<class Form, class Type>
+inline Type* Foam::Matrix<Form, Type>::operator[](const label i)
 {
 #   ifdef FULLDEBUG
     checki(i);
@@ -119,9 +122,8 @@ inline T* Foam::Matrix<T>::operator[](const label i)
 }
 
 
-// Return subscript-checked const element access
-template<class T>
-inline const T* Foam::Matrix<T>::operator[](const label i) const
+template<class Form, class Type>
+inline const Type* Foam::Matrix<Form, Type>::operator[](const label i) const
 {
 #   ifdef FULLDEBUG
     checki(i);
diff --git a/src/OpenFOAM/matrices/Matrix/MatrixIO.C b/src/OpenFOAM/matrices/Matrix/MatrixIO.C
index 99a3df9da..d97ae8803 100644
--- a/src/OpenFOAM/matrices/Matrix/MatrixIO.C
+++ b/src/OpenFOAM/matrices/Matrix/MatrixIO.C
@@ -32,8 +32,8 @@ License
 
 // * * * * * * * * * * * * * * * Ostream Operator *  * * * * * * * * * * * * //
 
-template<class T>
-Foam::Matrix<T>::Matrix(Istream& is)
+template<class Form, class Type>
+Foam::Matrix<Form, Type>::Matrix(Istream& is)
 :
     v_(NULL),
     n_(0),
@@ -43,17 +43,17 @@ Foam::Matrix<T>::Matrix(Istream& is)
 }
 
 
-template<class T>
-Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
+template<class Form, class Type>
+Foam::Istream& Foam::operator>>(Istream& is, Matrix<Form, Type>& M)
 {
     // Anull matrix
     M.clear();
 
-    is.fatalCheck("operator>>(Istream&, Matrix<T>&)");
+    is.fatalCheck("operator>>(Istream&, Matrix<Form, Type>&)");
 
     token firstToken(is);
 
-    is.fatalCheck("operator>>(Istream&, Matrix<T>&) : reading first token");
+    is.fatalCheck("operator>>(Istream&, Matrix<Form, Type>&) : reading first token");
 
     if (firstToken.isLabel())
     {
@@ -63,7 +63,7 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
         label nm = M.n_*M.m_;
 
         // Read list contents depending on data format
-        if (is.format() == IOstream::ASCII || !contiguous<T>())
+        if (is.format() == IOstream::ASCII || !contiguous<Type>())
         {
             // Read beginning of contents
             char listDelimiter = is.readBeginList("Matrix");
@@ -71,7 +71,7 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
             if (nm)
             {
                 M.allocate();
-                T* v = M.v_[0];
+                Type* v = M.v_[0];
 
                 if (listDelimiter == token::BEGIN_LIST)
                 {
@@ -88,7 +88,7 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
 
                             is.fatalCheck
                             (
-                                "operator>>(Istream&, Matrix<T>&) : "
+                                "operator>>(Istream&, Matrix<Form, Type>&) : "
                                 "reading entry"
                             );
                         }
@@ -98,12 +98,12 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
                 }
                 else
                 {
-                    T element;
+                    Type element;
                     is >> element;
 
                     is.fatalCheck
                     (
-                        "operator>>(Istream&, Matrix<T>&) : "
+                        "operator>>(Istream&, Matrix<Form, Type>&) : "
                         "reading the single entry"
                     );
 
@@ -122,13 +122,13 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
             if (nm)
             {
                 M.allocate();
-                T* v = M.v_[0];
+                Type* v = M.v_[0];
 
-                is.read(reinterpret_cast<char*>(v), nm*sizeof(T));
+                is.read(reinterpret_cast<char*>(v), nm*sizeof(Type));
 
                 is.fatalCheck
                 (
-                    "operator>>(Istream&, Matrix<T>&) : "
+                    "operator>>(Istream&, Matrix<Form, Type>&) : "
                     "reading the binary block"
                 );
             }
@@ -136,7 +136,7 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
     }
     else
     {
-        FatalIOErrorIn("operator>>(Istream&, Matrix<T>&)", is)
+        FatalIOErrorIn("operator>>(Istream&, Matrix<Form, Type>&)", is)
             << "incorrect first token, expected <int>, found "
             << firstToken.info()
             << exit(FatalIOError);
@@ -146,23 +146,23 @@ Foam::Istream& Foam::operator>>(Istream& is, Matrix<T>& M)
 }
 
 
-template<class T>
-Foam::Ostream& Foam::operator<<(Ostream& os, const Matrix<T>& M)
+template<class Form, class Type>
+Foam::Ostream& Foam::operator<<(Ostream& os, const Matrix<Form, Type>& M)
 {
     label nm = M.n_*M.m_;
 
     os  << M.n() << token::SPACE << M.m();
 
     // Write list contents depending on data format
-    if (os.format() == IOstream::ASCII || !contiguous<T>())
+    if (os.format() == IOstream::ASCII || !contiguous<Type>())
     {
         if (nm)
         {
             bool uniform = false;
 
-            const T* v = M.v_[0];
+            const Type* v = M.v_[0];
 
-            if (nm > 1 && contiguous<T>())
+            if (nm > 1 && contiguous<Type>())
             {
                 uniform = true;
 
@@ -189,7 +189,7 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const Matrix<T>& M)
             }
             // Fix: matrices smaller than 20x20 will be written square.
             // HJ, 22/Jan/2009
-            else if (nm < 400 && contiguous<T>())
+            else if (nm < 400 && contiguous<Type>())
             {
                 // Write size of list and start contents delimiter
                 os  << token::BEGIN_LIST;
@@ -248,7 +248,7 @@ Foam::Ostream& Foam::operator<<(Ostream& os, const Matrix<T>& M)
     {
         if (nm)
         {
-            os.write(reinterpret_cast<const char*>(M.v_[0]), nm*sizeof(T));
+            os.write(reinterpret_cast<const char*>(M.v_[0]), nm*sizeof(Type));
         }
     }
 
diff --git a/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.C b/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.C
index e8797a1ab..28333e6e6 100644
--- a/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.C
+++ b/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.C
@@ -25,7 +25,7 @@ License
 Class
     DenseMatrixTools
 
-\*----------------------------------------------------------------------------*/
+\*---------------------------------------------------------------------------*/
 
 #include "DenseMatrixTools.H"
 #include "Swap.H"
@@ -36,7 +36,7 @@ Class
 template<class T>
 static void Foam::DenseMatrixTools::solve
 (
-    Matrix<T>& A,
+    SquareMatrix<T>& A,
     List<T>& x,
     List<T>& b
 )
@@ -107,12 +107,12 @@ static void Foam::DenseMatrixTools::solve
 }
 
 
-template<class T>
+template<class Form, class Type>
 static void Foam::DenseMatrixTools::qrDecompose
 (
     const label nCols,
-    FieldField<Field, T>& A,
-    Matrix<T>& R
+    FieldField<Field, Type>& A,
+    Matrix<Form, T>& R
 )
 {
     // Note: consider Arnoldi algorithm for speed-up.  HJ, 14/Sep/2006
@@ -128,9 +128,9 @@ static void Foam::DenseMatrixTools::qrDecompose
 
         for (label k = j + 1; k < nCols; k++)
         {
-            const Field<T>& Aj = A[j];
-            Field<T>& Ak = A[k];
-            T& Rjk = R[j][k];
+            const Field<Type>& Aj = A[j];
+            Field<Type>& Ak = A[k];
+            Type& Rjk = R[j][k];
 
             Rjk = gSumProd(Aj, Ak);
 
diff --git a/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.H b/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.H
index 9545a6b99..178f330d5 100644
--- a/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.H
+++ b/src/OpenFOAM/matrices/Matrix/tools/DenseMatrixTools.H
@@ -33,7 +33,7 @@ SourceFiles
 #ifndef DenseMatrixTools_H
 #define DenseMatrixTools_H
 
-#include "Matrix.H"
+#include "SquareMatrix.H"
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
@@ -47,16 +47,16 @@ namespace Foam
 namespace DenseMatrixTools
 {
     //- Solve using Gauss Elimination with pivoting.  Destroys matrix and b
-    template<class T>
-    static void solve(Matrix<T>& A, List<T>& x, List<T>& b);
+    template<class Type>
+    static void solve(SquareMatrix<Type>& A, List<Type>& x, List<Type>& b);
 
     //- Q-R decomposition
-    template<class T>
+    template<class Form, class T>
     static void qrDecompose
     (
         const label nCols,
         FieldField<Field, T>& A,
-        Matrix<T>& R
+        Matrix<class Form, class Type>& R
     );
 };
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
index f6dd24e9b..f41cbc56e 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrix.H
@@ -250,6 +250,11 @@ public:
             //- Maximum number of iterations
             label maxIter_;
 
+
+    protected:
+
+        // Protected data
+
             //- Matrix reference
             const lduMatrix& matrix_;
 
@@ -451,7 +456,9 @@ public:
     //- Abstract base-class for lduMatrix smoothers
     class smoother
     {
-        // Private data
+    protected:
+
+        // Protected data
 
             //- Matrix reference
             const lduMatrix& matrix_;
@@ -585,7 +592,9 @@ public:
     //- Abstract base-class for lduMatrix preconditioners
     class preconditioner
     {
-        // Private data
+    protected:
+
+        // Protected data
 
             //- Matrix reference
             const lduMatrix& matrix_;
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixPreconditioner.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixPreconditioner.C
index 627c4e1d3..6e006291b 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixPreconditioner.C
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixPreconditioner.C
@@ -84,7 +84,7 @@ Foam::lduPreconditioner::New
 
     const dictionary& controls = e.isDict() ? e.dict() : dictionary::null;
 
-    if (sol.matrix().symmetric())
+    if (matrix.symmetric())
     {
         symMatrixConstructorTable::iterator constructorIter =
             symMatrixConstructorTablePtr_->find(preconName);
@@ -121,7 +121,7 @@ Foam::lduPreconditioner::New
             )
         );
     }
-    else if (sol.matrix().asymmetric())
+    else if (matrix.asymmetric())
     {
         asymMatrixConstructorTable::iterator constructorIter =
             asymMatrixConstructorTablePtr_->find(preconName);
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSmoother.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSmoother.C
index c5a919aea..b2ecabe12 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSmoother.C
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSmoother.C
@@ -172,6 +172,8 @@ Foam::autoPtr<Foam::lduMatrix::smoother> Foam::lduMatrix::smoother::New
         )   << "cannot solve incomplete matrix, "
                "no diagonal or off-diagonal coefficient"
             << exit(FatalIOError);
+
+        return autoPtr<lduSmoother>(NULL);
     }
 }
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C
index 37b0537fe..7e48e1fc3 100644
--- a/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/lduMatrix/lduMatrixSolver.C
@@ -48,7 +48,7 @@ Foam::autoPtr<Foam::lduMatrix::solver> Foam::lduMatrix::solver::New
     const dictionary& dict
 )
 {
-    word name(dict.lookup("solver"));
+    word solverName(dict.lookup("solver"));
 
     if (matrix.diagonal())
     {
@@ -243,10 +243,11 @@ Foam::lduMatrix::solver::solver
 
 void Foam::lduMatrix::solver::readControls()
 {
-    minIter_   = dict_.lookupOrDefault<label>("minIter", 0);
-    maxIter_   = dict_.lookupOrDefault<label>("maxIter", 1000);
     tolerance_ = dict_.lookupOrDefault<scalar>("tolerance", 1e-6);
-    relTol_    = dict_.lookupOrDefault<scalar>("relTol", 0);
+    relTolerance_ = dict_.lookupOrDefault<scalar>("relTol", 0);
+
+    minIter_ = dict_.lookupOrDefault<label>("minIter", 0);
+    maxIter_ = dict_.lookupOrDefault<label>("maxIter", 1000);
 }
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/preconditioners/DICPreconditioner/DICPreconditioner.C b/src/OpenFOAM/matrices/lduMatrix/preconditioners/DICPreconditioner/DICPreconditioner.C
index e2b85efb0..527613bfd 100644
--- a/src/OpenFOAM/matrices/lduMatrix/preconditioners/DICPreconditioner/DICPreconditioner.C
+++ b/src/OpenFOAM/matrices/lduMatrix/preconditioners/DICPreconditioner/DICPreconditioner.C
@@ -38,21 +38,6 @@ namespace Foam
 }
 
 
-// * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
-
-Foam::DICPreconditioner::DICPreconditioner
-(
-    const lduMatrix::solver& sol,
-    const dictionary&
-)
-:
-    lduMatrix::preconditioner(sol),
-    rD_(sol.matrix().diag())
-{
-    calcReciprocalD(rD_, sol.matrix());
-}
-
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 void Foam::DICPreconditioner::calcReciprocalD
@@ -93,7 +78,7 @@ Foam::DICPreconditioner::DICPreconditioner
     const FieldField<Field, scalar>& coupleBouCoeffs,
     const FieldField<Field, scalar>& coupleIntCoeffs,
     const lduInterfaceFieldPtrsList& interfaces,
-    const dictionary& 
+    const dictionary&
 )
 :
     lduPreconditioner
diff --git a/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.C b/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.C
index b16e3f045..a70c39b5e 100644
--- a/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.C
+++ b/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.C
@@ -84,7 +84,7 @@ Foam::GAMGPreconditioner::~GAMGPreconditioner()
 void Foam::GAMGPreconditioner::readControls()
 {
     GAMG_.readControls();
-    nVcycles_ = controlDict_.lookupOrDefault<label>("nVcycles", 2);
+    nVcycles_ = GAMG_.dict().lookupOrDefault<label>("nVcycles", 2);
 }
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.H b/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.H
index 39007d83a..5dc497d31 100644
--- a/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.H
+++ b/src/OpenFOAM/matrices/lduMatrix/preconditioners/GAMGPreconditioner/GAMGPreconditioner.H
@@ -62,7 +62,7 @@ class GAMGPreconditioner
         //- Number of V-cycles to perform
         label nVcycles_;
 
-        //- Read the control parameters from the controlDict_
+        //- Read the control parameters from the dictionary
         virtual void readControls();
 
 public:
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/BICCG/BICCG.H b/src/OpenFOAM/matrices/lduMatrix/solvers/BICCG/BICCG.H
index 6556b363c..f7d17c4de 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/BICCG/BICCG.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/BICCG/BICCG.H
@@ -94,7 +94,7 @@ public:
             const FieldField<Field, scalar>& coupleBouCoeffs,
             const FieldField<Field, scalar>& coupleIntCoeffs,
             const lduInterfaceFieldPtrsList& interfaces,
-            const dictionary& solverControls
+            const dictionary& dict
         );
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
index 2562d8f31..a43eb3780 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolver.C
@@ -108,7 +108,7 @@ Foam::GAMGSolver::GAMGSolver
     nFinestSweeps_(2),
     scaleCorrection_(matrix.symmetric()),
     directSolveCoarsest_(false),
-    agglomeration_(GAMGAgglomeration::New(matrix_, dict())),
+    agglomeration_(GAMGAgglomeration::New(matrix_, dict)),
 
     matrixLevels_(agglomeration_.size()),
     interfaceLevels_(agglomeration_.size()),
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
index a6dd97da3..09bc4bced 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/GAMG/GAMGSolverSolve.C
@@ -75,7 +75,7 @@ Foam::lduMatrix::solverPerformance Foam::GAMGSolver::solve
         PtrList<scalarField> coarseB;
 
         // Create the smoothers for all levels
-        PtrList<lduMatrix::smoother> smoothers;
+        PtrList<lduSmoother> smoothers;
 
         // Initialise the above data structures
         initVcycle(coarseCorrX, coarseB, smoothers);
@@ -115,7 +115,7 @@ Foam::lduMatrix::solverPerformance Foam::GAMGSolver::solve
 
 void Foam::GAMGSolver::Vcycle
 (
-    const PtrList<lduMatrix::smoother>& smoothers,
+    const PtrList<lduSmoother>& smoothers,
     scalarField& x,
     const scalarField& b,
     scalarField& Ax,
@@ -352,7 +352,7 @@ void Foam::GAMGSolver::initVcycle
 (
     PtrList<scalarField>& coarseCorrX,
     PtrList<scalarField>& coarseB,
-    PtrList<lduMatrix::smoother>& smoothers
+    PtrList<lduSmoother>& smoothers
 ) const
 {
     coarseCorrX.setSize(matrixLevels_.size());
@@ -363,13 +363,13 @@ void Foam::GAMGSolver::initVcycle
     smoothers.set
     (
         0,
-        lduMatrix::smoother::New
+        lduSmoother::New
         (
-            dict().lookup("smoother"),
             matrix_,
             coupleBouCoeffs_,
             coupleIntCoeffs_,
-            interfaces_
+            interfaces_,
+            dict()
         )
     );
 
@@ -396,13 +396,13 @@ void Foam::GAMGSolver::initVcycle
         smoothers.set
         (
             leveli + 1,
-            lduMatrix::smoother::New
+            lduSmoother::New
             (
-                dict().lookup("smoother"),
                 matrixLevels_[leveli],
                 coupleLevelsBouCoeffs_[leveli],
                 coupleLevelsIntCoeffs_[leveli],
-                interfaceLevels_[leveli]
+                interfaceLevels_[leveli],
+                dict()
             )
         );
     }
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.C
index 810f83864..5775ce8b7 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.C
@@ -72,7 +72,7 @@ Foam::ICCG::ICCG
     const FieldField<Field, scalar>& coupleBouCoeffs,
     const FieldField<Field, scalar>& coupleIntCoeffs,
     const lduInterfaceFieldPtrsList& interfaces,
-    const dictionary& solverControls
+    const dictionary& dict
 )
 :
     PCG
@@ -82,7 +82,7 @@ Foam::ICCG::ICCG
         coupleBouCoeffs,
         coupleIntCoeffs,
         interfaces,
-        solverControls
+        dict
     )
 {}
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.H b/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.H
index 585548a76..99f343d85 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/ICCG/ICCG.H
@@ -94,7 +94,7 @@ public:
             const FieldField<Field, scalar>& coupleBouCoeffs,
             const FieldField<Field, scalar>& coupleIntCoeffs,
             const lduInterfaceFieldPtrsList& interfaces,
-            const dictionary& solverControls
+            const dictionary& dict
         );
 
         //- Construct from matrix components and tolerances
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
index 0127fd30c..edc0caa73 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PBiCG/PBiCG.C
@@ -46,7 +46,7 @@ Foam::PBiCG::PBiCG
     const FieldField<Field, scalar>& coupleBouCoeffs,
     const FieldField<Field, scalar>& coupleIntCoeffs,
     const lduInterfaceFieldPtrsList& interfaces,
-    const dictionary& solverControls
+    const dictionary& dict
 )
 :
     lduSolver
@@ -56,7 +56,7 @@ Foam::PBiCG::PBiCG
         coupleBouCoeffs,
         coupleIntCoeffs,
         interfaces,
-        solverControls
+        dict
     )
 {}
 
@@ -71,10 +71,10 @@ Foam::lduSolverPerformance Foam::PBiCG::solve
 ) const
 {
     // --- Setup class containing solver performance data
-    lduMatrix::solverPerformance solverPerf
+    lduSolverPerformance solverPerf
     (
-        lduMatrix::preconditioner::getName(controlDict_) + typeName,
-        fieldName_
+        lduMatrix::preconditioner::getName(dict()) + typeName,
+        fieldName()
     );
 
 
@@ -206,7 +206,7 @@ Foam::lduSolverPerformance Foam::PBiCG::solve
 
             for (register label cell=0; cell<nCells; cell++)
             {
-                psiPtr[cell] += alpha*pAPtr[cell];
+                xPtr[cell] += alpha*pAPtr[cell];
                 rAPtr[cell] -= alpha*wAPtr[cell];
                 rTPtr[cell] -= alpha*wTPtr[cell];
             }
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
index 43741a6e1..2e8d23b22 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.C
@@ -46,7 +46,7 @@ Foam::PCG::PCG
     const FieldField<Field, scalar>& coupleBouCoeffs,
     const FieldField<Field, scalar>& coupleIntCoeffs,
     const lduInterfaceFieldPtrsList& interfaces,
-    const dictionary& solverControls
+    const dictionary& dict
 )
 :
     lduSolver
@@ -56,7 +56,7 @@ Foam::PCG::PCG
         coupleBouCoeffs,
         coupleIntCoeffs,
         interfaces,
-        solverControls
+        dict
     )
 {}
 
@@ -71,10 +71,10 @@ Foam::lduSolverPerformance Foam::PCG::solve
 ) const
 {
     // --- Setup class containing solver performance data
-    lduMatrix::solverPerformance solverPerf
+    lduSolverPerformance solverPerf
     (
-        lduMatrix::preconditioner::getName(controlDict_) + typeName,
-        fieldName_
+        lduMatrix::preconditioner::getName(dict()) + typeName,
+        fieldName()
     );
 
     register label nCells = x.size();
@@ -194,7 +194,7 @@ Foam::lduSolverPerformance Foam::PCG::solve
 
             for (register label cell=0; cell<nCells; cell++)
             {
-                psiPtr[cell] += alpha*pAPtr[cell];
+                xPtr[cell] += alpha*pAPtr[cell];
                 rAPtr[cell] -= alpha*wAPtr[cell];
             }
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.H b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.H
index 7e9764e00..7c5537b78 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/PCG/PCG.H
@@ -77,7 +77,7 @@ public:
             const FieldField<Field, scalar>& coupleBouCoeffs,
             const FieldField<Field, scalar>& coupleIntCoeffs,
             const lduInterfaceFieldPtrsList& interfaces,
-            const dictionary& solverControls
+            const dictionary& dict
         );
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.C
index ebccfc168..f6bb0b55a 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.C
@@ -43,7 +43,7 @@ Foam::diagonalSolver::diagonalSolver
     const FieldField<Field, scalar>& coupleBouCoeffs,
     const FieldField<Field, scalar>& coupleIntCoeffs,
     const lduInterfaceFieldPtrsList& interfaces,
-    const dictionary& solverControls
+    const dictionary& dict
 )
 :
     lduMatrix::solver
@@ -53,7 +53,7 @@ Foam::diagonalSolver::diagonalSolver
         coupleBouCoeffs,
         coupleIntCoeffs,
         interfaces,
-        solverControls
+        dict
     )
 {}
 
@@ -72,7 +72,7 @@ Foam::lduMatrix::solverPerformance Foam::diagonalSolver::solve
     return lduSolverPerformance
     (
         typeName,
-        fieldName_,
+        fieldName(),
         0,
         0,
         0,
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.H b/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.H
index dbf2679f8..92328116e 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.H
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/diagonalSolver/diagonalSolver.H
@@ -76,7 +76,7 @@ public:
             const FieldField<Field, scalar>& coupleBouCoeffs,
             const FieldField<Field, scalar>& coupleIntCoeffs,
             const lduInterfaceFieldPtrsList& interfaces,
-            const dictionary& solverControls
+            const dictionary& dict
         );
 
 
diff --git a/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C b/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
index 8a8c33b76..90ad1dcd0 100644
--- a/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
+++ b/src/OpenFOAM/matrices/lduMatrix/solvers/smoothSolver/smoothSolver.C
@@ -49,7 +49,7 @@ Foam::smoothSolver::smoothSolver
     const FieldField<Field, scalar>& coupleBouCoeffs,
     const FieldField<Field, scalar>& coupleIntCoeffs,
     const lduInterfaceFieldPtrsList& interfaces,
-    const dictionary& solverControls
+    const dictionary& dict
 )
 :
     lduSolver
@@ -59,7 +59,7 @@ Foam::smoothSolver::smoothSolver
         coupleBouCoeffs,
         coupleIntCoeffs,
         interfaces,
-        solverControls
+        dict
     )
 {
     readControls();
@@ -71,7 +71,7 @@ Foam::smoothSolver::smoothSolver
 void Foam::smoothSolver::readControls()
 {
     lduSolver::readControls();
-    nSweeps_ = controlDict_.lookupOrDefault<label>("nSweeps", 1);
+    nSweeps_ = dict().lookupOrDefault<label>("nSweeps", 1);
 }
 
 
@@ -93,6 +93,7 @@ Foam::lduSolverPerformance Foam::smoothSolver::solve
             matrix_,
             coupleBouCoeffs_,
             coupleIntCoeffs_,
+            interfaces_,
             dict()
         );
 
@@ -140,6 +141,7 @@ Foam::lduSolverPerformance Foam::smoothSolver::solve
                 matrix_,
                 coupleBouCoeffs_,
                 coupleIntCoeffs_,
+                interfaces_,
                 dict()
             );
 
diff --git a/src/OpenFOAM/meshes/pointMesh/pointMesh.H b/src/OpenFOAM/meshes/pointMesh/pointMesh.H
index b74aad95c..ed82b8608 100644
--- a/src/OpenFOAM/meshes/pointMesh/pointMesh.H
+++ b/src/OpenFOAM/meshes/pointMesh/pointMesh.H
@@ -106,7 +106,7 @@ public:
         //  HJ, 16/Oct/2008
         label nPoints() const
         {
-            return mesh_.nPoints();
+            return mesh.GeoMesh<polyMesh>::mesh_.nPoints();
         }
 
         //- Return number of points
@@ -115,7 +115,7 @@ public:
         //  HJ, 16/Oct/2008
         label nCells() const
         {
-            return mesh_.nCells();
+            return mesh.GeoMesh<polyMesh>::mesh_.nCells();
         }
 
         //- Return reference to boundary mesh
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.C
index 6ba9864e1..9b2d5bd18 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.C
@@ -26,11 +26,6 @@ License
 
 #include "Map.H"
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-namespace Foam
-{
-
 // * * * * * * * * * * * * * * Static Data Members * * * * * * * * * * * * * //
 
 template
@@ -41,7 +36,8 @@ template
     class PointType
 >
 const bool
-PrimitivePatch<Face, FaceList, PointField, PointType>::nSquaredProjection_
+Foam::PrimitivePatch<Face, FaceList, PointField, PointType>::
+nSquaredProjection_
 (
     debug::optimisationSwitch("nSquaredProjection", 0) > 0
 );
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.H b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.H
index 357c1d142..8d6c75915 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.H
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatch.H
@@ -69,14 +69,14 @@ class objectHit;
 template<class T> class Map;
 
 /*---------------------------------------------------------------------------*\
-                        Class PrimitivePatchName Declaration
+                      Class PrimitivePatchName Declaration
 \*---------------------------------------------------------------------------*/
 
 TemplateName(PrimitivePatch);
 
 
 /*---------------------------------------------------------------------------*\
-                           Class PrimitivePatch Declaration
+                        Class PrimitivePatch Declaration
 \*---------------------------------------------------------------------------*/
 
 template
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchBdryPoints.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchBdryPoints.C
index 1989b59ad..43ce8c2eb 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchBdryPoints.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchBdryPoints.C
@@ -27,7 +27,6 @@ License
 #include "PrimitivePatch.H"
 #include "HashSet.H"
 
-
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
 template
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchCheck.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchCheck.C
index c9355e623..4d35a47f1 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchCheck.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchCheck.C
@@ -30,7 +30,7 @@ Description
 #include "PrimitivePatch.H"
 #include "Map.H"
 #include "ListOps.H"
-
+#include "OFstream.H"
 
 // * * * * * * * * * * * * * Private Member Functions  * * * * * * * * * * * //
 
@@ -497,8 +497,4 @@ writeVTKNormals
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
-} // End namespace Foam
-
 // ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchClear.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchClear.C
index c92fd5bd7..2603b065e 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchClear.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchClear.C
@@ -27,7 +27,6 @@ License
 #include "PrimitivePatch.H"
 #include "demandDrivenData.H"
 
-
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 template
diff --git a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchMeshData.C b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchMeshData.C
index c0b0612e7..9ef9b7aa1 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchMeshData.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/PrimitivePatch/PrimitivePatchMeshData.C
@@ -91,31 +91,29 @@ calcMeshData() const
     forAll(*this, facei)
     {
         const Face& curPoints = this->operator[](facei);
-    
+
         forAll(curPoints, pointi)
         {
             markedPoints.insert(curPoints[pointi], -1);
         }
     }
- 
+
     // Create the storage and store the meshPoints.  Mesh points are
     // the ones marked by the usage loop above
     meshPointsPtr_ = new labelList(markedPoints.toc());
     labelList& pointPatch = *meshPointsPtr_;
-    
+
     // Sort the list to preserve compatibility with the old ordering
     sort(pointPatch);
-    
+
     // For every point in map give it its label in mesh points
     forAll(pointPatch, pointi)
     {
         markedPoints.find(pointPatch[pointi])() = pointi;
     }
 
-    //- Unsorted version:
-    //DynamicList<label> meshPoints(2*this->size());
-    //forAll(*this, facei)
-    //{
+    forAll(*this, faceI)
+    {
         const Face& curPoints = this->operator[](faceI);
 
         forAll (curPoints, pointI)
@@ -124,20 +122,6 @@ calcMeshData() const
         }
     }
 
-    // Create the storage and store the meshPoints.  Mesh points are
-    // the ones marked by the usage loop above
-    meshPointsPtr_ = new labelList(markedPoints.toc());
-    labelList& pointPatch = *meshPointsPtr_;
-
-    // Sort the list to preserve compatibility with the old ordering
-    sort(pointPatch);
-
-    // For every point in map give it its label in mesh points
-    forAll (pointPatch, pointI)
-    {
-        markedPoints.find(pointPatch[pointI])() = pointI;
-    }
-
     // Create local faces. Note that we start off from copy of original face
     // list (even though vertices are overwritten below). This is done so
     // additional data gets copied (e.g. region number of labelledTri)
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C
index d70e0e96c..5a5dcf312 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.C
@@ -125,85 +125,6 @@ Foam::primitiveMesh::~primitiveMesh()
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-bool Foam::primitiveMesh::calcPointOrder
-(
-    label& nInternalPoints,
-    labelList& oldToNew,
-    const faceList& faces,
-    const label nInternalFaces,
-    const label nPoints
-)
-{
-    // Internal points are points that are not used by a boundary face.
-
-    // Map from old to new position
-    oldToNew.setSize(nPoints);
-    oldToNew = -1;
-
-
-    // 1. Create compact addressing for boundary points. Start off by indexing
-    // from 0 inside oldToNew. (shifted up later on)
-
-    label nBoundaryPoints = 0;
-    for (label faceI = nInternalFaces; faceI < faces.size(); faceI++)
-    {
-        const face& f = faces[faceI];
-
-        forAll(f, fp)
-        {
-            label pointI = f[fp];
-
-            if (oldToNew[pointI] == -1)
-            {
-                oldToNew[pointI] = nBoundaryPoints++;
-            }
-        }
-    }
-
-    // Now we know the number of boundary and internal points
-
-    nInternalPoints = nPoints - nBoundaryPoints;
-
-    // Move the boundary addressing up
-    forAll(oldToNew, pointI)
-    {
-        if (oldToNew[pointI] != -1)
-        {
-            oldToNew[pointI] += nInternalPoints;
-        }
-    }
-
-
-    // 2. Compact the internal points. Detect whether internal and boundary
-    // points are mixed.
-
-    label internalPointI = 0;
-
-    bool ordered = true;
-
-    for (label faceI = 0; faceI < nInternalFaces; faceI++)
-    {
-        const face& f = faces[faceI];
-
-        forAll(f, fp)
-        {
-            label pointI = f[fp];
-
-            if (oldToNew[pointI] == -1)
-            {
-                if (pointI >= nInternalPoints)
-                {
-                    ordered = false;
-                }
-                oldToNew[pointI] = internalPointI++;
-            }
-        }
-    }
-
-    return ordered;
-}
-
-
 void Foam::primitiveMesh::reset
 (
     const label nPoints,
@@ -216,40 +137,14 @@ void Foam::primitiveMesh::reset
 
     nPoints_ = nPoints;
     nEdges_ = -1;
-    nInternal0Edges_ = -1;
-    nInternal1Edges_ = -1;
-    nInternalEdges_ = -1;
 
     nInternalFaces_ = nInternalFaces;
     nFaces_ = nFaces;
     nCells_ = nCells;
 
-    // Check if points are ordered
-    label nInternalPoints;
-    labelList pointMap;
-
-    bool isOrdered = calcPointOrder
-    (
-        nInternalPoints,
-        pointMap,
-        faces(),
-        nInternalFaces_,
-        nPoints_
-    );
-
-    if (isOrdered)
-    {
-        nInternalPoints_ = nInternalPoints;
-    }
-    else
-    {
-        nInternalPoints_ = -1;
-    }
-
     if (debug)
     {
         Pout<< "primitiveMesh::reset : mesh reset to"
-            << " nInternalPoints:" << nInternalPoints_
             << " nPoints:" << nPoints_
             << " nEdges:" << nEdges_
             << " nInternalFaces:" << nInternalFaces_
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H
index 9f570967b..cf3a48bd6 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMesh.H
@@ -200,12 +200,15 @@ class primitiveMesh
             //- Calculate point-point addressing
             void calcPointPoints() const;
 
-            //- Calculate edges, pointEdges and faceEdges (if doFaceEdges=true)
-            //  During edge calculation, a larger set of data is assembled.
-            //  Create and destroy as a set, using clearOutEdges()
-            void calcEdges(const bool doFaceEdges) const;
+            //- Calculate edges, pointEdges and faceEdges
+            void calcEdges(const bool) const;
+
+            //- During edge calculation, a larger set of data is assembled.
+            // Create and destroy as a set, using clearOutEdges()
             void clearOutEdges();
-            //- Helper: return (after optional creation) edge between two points
+
+            //- Helper:
+            //  return (after optional creation) edge between two points
             static label getEdge
             (
                 List<DynamicList<label> >&,
@@ -432,18 +435,8 @@ public:
                     const label nCells = -1
                 );
 
-                //- Helper function to calculate point ordering. Returns true
-                //  if points already ordered, false and fills pointMap (old to
-                //  new). Map splits points into those not used by any boundary
-                //  face and those that are.
-                static bool calcPointOrder
-                (
-                    label& nInternalPoints,
-                    labelList& pointMap,
-                    const faceList&,
-                    const label nInternalFaces,
-                    const label nPoints
-                );
+            // Removed calcPointOrder - garbage.  HJ, 27/Aug/2010
+
 
             // Return mesh connectivity
 
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshClear.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshClear.C
index be6971e0e..0e9a5f978 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshClear.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshClear.C
@@ -94,7 +94,7 @@ void Foam::primitiveMesh::printAllocated() const
     {
         Pout<< "    Point-point" << endl;
     }
-    
+
     if (cpPtr_)
     {
         Pout<< "    Cell-point" << endl;
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C
index f063693c3..8ebe95fb7 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshEdges.C
@@ -71,386 +71,232 @@ Foam::label Foam::primitiveMesh::getEdge
 }
 
 
-void Foam::primitiveMesh::calcEdges(const bool doFaceEdges) const
+void Foam::primitiveMesh::calcEdges() const
 {
     if (debug)
     {
-        Pout<< "primitiveMesh::calcEdges(const bool) : "
-            << "calculating edges, pointEdges and optionally faceEdges"
+        Pout<< "primitiveMesh::calcEdges() : "
+            << "calculating edges, pointEdges and faceEdges"
             << endl;
     }
 
     // It is an error to attempt to recalculate edges
     // if the pointer is already set
-    if ((edgesPtr_ || pePtr_) || (doFaceEdges && fePtr_))
+    if (edgesPtr_ || fePtr_)
     {
         FatalErrorIn("primitiveMesh::calcEdges(const bool) const")
-            << "edges or pointEdges or faceEdges already calculated"
+            << "edges or faceEdges already calculated"
             << abort(FatalError);
     }
     else
     {
+        // Note: replaced with the original and correct algorithm
+        // which preserved correct ordering for edges list
+        // Version 1.5.x, 1.5-dev were wrong, with failures in parallel
+        // point-based data exchange.  Fixed in 1.6-ext and
+        // consecutive versions
+        // HJ, 27/Aug/2010
+
         // ALGORITHM:
-        // Go through the faces list. Search pointEdges for existing edge.
-        // If not found create edge and add to pointEdges.
-        // In second loop sort edges in order of increasing neighbouring
-        // point.
-        // This algorithm replaces the one using pointFaces which used more
-        // allocations but less memory and was on practical cases
-        // quite a bit slower.
+        // Go through the pointFace list.  Go through the list of faces for that
+        // point and ask for edges.  If the edge has got the point in question
+        // AND the second point in the edge is larger than the first, add the
+        // edge to the list.  At the same time, add the edge label to the list
+        // of edges for the current face (faceEdges) and log the other face as
+        // the neighbour of this face.
 
-        const faceList& fcs = faces();
+        const faceList& f = faces();
 
-        // Size up lists
-        // ~~~~~~~~~~~~~
+        const labelListList& pf = pointFaces();
 
-        // Estimate pointEdges storage
-        List<DynamicList<label> > pe(nPoints());
-        forAll(pe, pointI)
+        fePtr_ = new labelListList(nFaces());
+        labelListList& fe = *fePtr_;
+
+        // count the maximum number of edges
+        label maxEdges = 0;
+
+        // create a list of edges for each face and store for efficiency
+        edgeListList edgesOfFace(nFaces());
+
+        forAll (f, faceI)
         {
-            pe[pointI].setCapacity(primitiveMesh::edgesPerPoint_);
-        }
+            edgesOfFace[faceI] = f[faceI].edges();
 
-        // Estimate edges storage
-        DynamicList<edge> es(pe.size()*primitiveMesh::edgesPerPoint_/2);
+            maxEdges += f[faceI].nEdges();
 
-        // Estimate faceEdges storage
-        if (doFaceEdges)
-        {
-            fePtr_ = new labelListList(fcs.size());
-            labelListList& faceEdges = *fePtr_;
-            forAll(fcs, faceI)
+            labelList& curFE = fe[faceI];
+
+            curFE.setSize(f[faceI].nEdges());
+
+            forAll (curFE, curFEI)
             {
-                faceEdges[faceI].setSize(fcs[faceI].size());
+                curFE[curFEI] = -1;
             }
         }
 
+        // EDGE CALCULATION
 
-        // Find consecutive face points in edge list
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
+        edgesPtr_ = new edgeList(maxEdges);
+        edgeList& e = *edgesPtr_;
+        label nEdges = 0;
 
-        // Edges on boundary faces
-        label nExtEdges = 0;
-        // Edges using no boundary point
-        nInternal0Edges_ = 0;
-        // Edges using 1 boundary point
-        label nInt1Edges = 0;
-        // Edges using two boundary points but not on boundary face:
-        // edges.size()-nExtEdges-nInternal0Edges_-nInt1Edges
-
-        // Ordering is different if points are ordered.
-        if (nInternalPoints_ == -1)
+        forAll (pf, pointI)
         {
-            // No ordering. No distinction between types.
-            forAll(fcs, faceI)
+            const labelList& curFaces = pf[pointI];
+
+            // create a list of labels to keep the neighbours that
+            // have already been added
+            DynamicList<label, edgesPerPoint_> addedNeighbours;
+            DynamicList<DynamicList<label, edgesPerPoint_> >
+                faceGivingNeighbour;
+            DynamicList<DynamicList<label, edgesPerPoint_> >
+                edgeOfFaceGivingNeighbour;
+
+            forAll (curFaces, faceI)
             {
-                const face& f = fcs[faceI];
+                // get the edges
+                const edgeList& fEdges = edgesOfFace[curFaces[faceI]];
 
-                forAll(f, fp)
+                // for every edge
+                forAll(fEdges, edgeI)
                 {
-                    label pointI = f[fp];
-                    label nextPointI = f[f.fcIndex(fp)];
+                    const edge& ends = fEdges[edgeI];
 
-                    label edgeI = getEdge(pe, es, pointI, nextPointI);
+                    // does the edge has got the point in question
+                    bool found = false;
+                    label secondPoint = -1;
 
-                    if (doFaceEdges)
+                    if (ends.start() == pointI)
                     {
-                        (*fePtr_)[faceI][fp] = edgeI;
+                        found = true;
+                        secondPoint = ends.end();
                     }
-                }
-            }
-            // Assume all edges internal
-            nExtEdges = 0;
-            nInternal0Edges_ = es.size();
-            nInt1Edges = 0;
-        }
-        else
-        {
-            // 1. Do external faces first. This creates external edges.
-            for (label faceI = nInternalFaces_; faceI < fcs.size(); faceI++)
-            {
-                const face& f = fcs[faceI];
 
-                forAll(f, fp)
-                {
-                    label pointI = f[fp];
-                    label nextPointI = f[f.fcIndex(fp)];
-
-                    label oldNEdges = es.size();
-                    label edgeI = getEdge(pe, es, pointI, nextPointI);
-
-                    if (es.size() > oldNEdges)
+                    if (ends.end() == pointI)
                     {
-                        nExtEdges++;
+                        found = true;
+                        secondPoint = ends.start();
                     }
-                    if (doFaceEdges)
+
+                    // if the edge has got the point and second label is larger
+                    // than first, it is a candidate for adding
+                    if (found && (secondPoint > pointI))
                     {
-                        (*fePtr_)[faceI][fp] = edgeI;
-                    }
-                }
-            }
+                        // check if the edge has already been added
+                        bool added = false;
 
-            // 2. Do internal faces. This creates internal edges.
-            for (label faceI = 0; faceI < nInternalFaces_; faceI++)
-            {
-                const face& f = fcs[faceI];
-
-                forAll(f, fp)
-                {
-                    label pointI = f[fp];
-                    label nextPointI = f[f.fcIndex(fp)];
-
-                    label oldNEdges = es.size();
-                    label edgeI = getEdge(pe, es, pointI, nextPointI);
-
-                    if (es.size() > oldNEdges)
-                    {
-                        if (pointI < nInternalPoints_)
+                        forAll (addedNeighbours, eopI)
                         {
-                            if (nextPointI < nInternalPoints_)
+                            if (secondPoint == addedNeighbours[eopI])
                             {
-                                nInternal0Edges_++;
-                            }
-                            else
-                            {
-                                nInt1Edges++;
+                                // Edge is already added. New face sharing it
+                                added = true;
+
+                                // Remember the face and edge giving neighbour
+                                faceGivingNeighbour[eopI].append
+                                    (curFaces[faceI]);
+
+                                edgeOfFaceGivingNeighbour[eopI].append(edgeI);
+
+                                break;
                             }
                         }
-                        else
+
+                        // If not added, add the edge to the list
+                        if (!added)
                         {
-                            if (nextPointI < nInternalPoints_)
-                            {
-                                nInt1Edges++;
-                            }
-                            else
-                            {
-                                // Internal edge with two points on boundary
-                            }
-                        }
-                    }
-                    if (doFaceEdges)
-                    {
-                        (*fePtr_)[faceI][fp] = edgeI;
-                    }
-                }
-            }
-        }
+                            addedNeighbours.append(secondPoint);
 
-
-        // For unsorted meshes the edges will be in order of occurrence inside
-        // the faces. For sorted meshes the first nExtEdges will be the external
-        // edges.
-
-        if (nInternalPoints_ != -1)
-        {
-            nInternalEdges_ = es.size()-nExtEdges;
-            nInternal1Edges_ = nInternal0Edges_+nInt1Edges;
-
-            //Pout<< "Edge overview:" << nl
-            //    << "    total number of edges           : " << es.size()
-            //    << nl
-            //    << "    boundary edges                  : " << nExtEdges
-            //    << nl
-            //    << "    edges using no boundary point   : "
-            //    << nInternal0Edges_
-            //    << nl
-            //    << "    edges using one boundary points : " << nInt1Edges
-            //   << nl
-            //    << "    edges using two boundary points : "
-            //    << es.size()-nExtEdges-nInternal0Edges_-nInt1Edges << endl;
-        }
-
-
-        // Like faces sort edges in order of increasing neigbouring point.
-        // ~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~~
-        // Automatically if points are sorted into internal and external points
-        // the edges will be sorted into
-        // - internal point to internal point
-        // - internal point to external point
-        // - external point to external point
-        // The problem is that the last one mixes external edges with internal
-        // edges with two points on the boundary.
-
-
-        // Map to sort into new upper-triangular order
-        labelList oldToNew(es.size());
-        if (debug)
-        {
-            oldToNew = -1;
-        }
-
-        // start of edges with 0 boundary points
-        label internal0EdgeI = 0;
-
-        // start of edges with 1 boundary points
-        label internal1EdgeI = nInternal0Edges_;
-
-        // start of edges with 2 boundary points
-        label internal2EdgeI = nInternal1Edges_;
-
-        // start of external edges
-        label externalEdgeI = nInternalEdges_;
-
-
-        // To sort neighbouring points in increasing order. Defined outside
-        // for optimisation reasons: if all connectivity size same will need
-        // no reallocations
-        SortableList<label> nbrPoints(primitiveMesh::edgesPerPoint_);
-
-        forAll(pe, pointI)
-        {
-            const DynamicList<label>& pEdges = pe[pointI];
-
-            nbrPoints.setSize(pEdges.size());
-
-            forAll(pEdges, i)
-            {
-                const edge& e = es[pEdges[i]];
-
-                label nbrPointI = e.otherVertex(pointI);
-
-                if (nbrPointI < pointI)
-                {
-                    nbrPoints[i] = -1;
-                }
-                else
-                {
-                    nbrPoints[i] = nbrPointI;
-                }
-            }
-            nbrPoints.sort();
-
-
-            if (nInternalPoints_ == -1)
-            {
-                // Sort all upper-triangular
-                forAll(nbrPoints, i)
-                {
-                    if (nbrPoints[i] != -1)
-                    {
-                        label edgeI = pEdges[nbrPoints.indices()[i]];
-
-                        oldToNew[edgeI] = internal0EdgeI++;
-                    }
-                }
-            }
-            else
-            {
-                if (pointI < nInternalPoints_)
-                {
-                    forAll(nbrPoints, i)
-                    {
-                        label nbrPointI = nbrPoints[i];
-
-                        label edgeI = pEdges[nbrPoints.indices()[i]];
-
-                        if (nbrPointI != -1)
-                        {
-                            if (edgeI < nExtEdges)
-                            {
-                                // External edge
-                                oldToNew[edgeI] = externalEdgeI++;
-                            }
-                            else if (nbrPointI < nInternalPoints_)
-                            {
-                                // Both points inside
-                                oldToNew[edgeI] = internal0EdgeI++;
-                            }
-                            else
-                            {
-                                // One points inside, one outside
-                                oldToNew[edgeI] = internal1EdgeI++;
-                            }
-                        }
-                    }
-                }
-                else
-                {
-                    forAll(nbrPoints, i)
-                    {
-                        label nbrPointI = nbrPoints[i];
-
-                        label edgeI = pEdges[nbrPoints.indices()[i]];
-
-                        if (nbrPointI != -1)
-                        {
-                            if (edgeI < nExtEdges)
-                            {
-                                // External edge
-                                oldToNew[edgeI] = externalEdgeI++;
-                            }
-                            else if (nbrPointI < nInternalPoints_)
-                            {
-                                // Not possible!
-                                FatalErrorIn("primitiveMesh::calcEdges(..)")
-                                    << abort(FatalError);
-                            }
-                            else
-                            {
-                                // Both points outside
-                                oldToNew[edgeI] = internal2EdgeI++;
-                            }
+                            // Remember the face and subShape giving neighbour
+                            faceGivingNeighbour(addedNeighbours.size() - 1)
+                                .append(curFaces[faceI]);
+                            edgeOfFaceGivingNeighbour
+                                (addedNeighbours.size() - 1).append(edgeI);
                         }
                     }
                 }
             }
-        }
 
+            // All edges for the current point found. Before adding them to the
+            // list, it is necessary to sort them in the increasing order of the
+            // neighbouring point.
 
-        if (debug)
-        {
-            label edgeI = findIndex(oldToNew, -1);
+            // Make real list out of SLList to simplify the manipulation.
+            // Also, make another list to "remember" how the original list was
+            // reshuffled.
+            labelList shuffleList(addedNeighbours.size());
 
-            if (edgeI != -1)
+            forAll (shuffleList, i)
             {
-                const edge& e = es[edgeI];
+                shuffleList[i] = i;
+            }
 
-                FatalErrorIn("primitiveMesh::calcEdges(..)")
-                    << "Did not sort edge " << edgeI << " points:" << e
-                    << " coords:" << points()[e[0]] << points()[e[1]]
-                    << endl
-                    << "Current buckets:" << endl
-                    << "    internal0EdgeI:" << internal0EdgeI << endl
-                    << "    internal1EdgeI:" << internal1EdgeI << endl
-                    << "    internal2EdgeI:" << internal2EdgeI << endl
-                    << "    externalEdgeI:" << externalEdgeI << endl
-                    << abort(FatalError);
+            // Use a simple sort to sort the addedNeighbours list.
+            //  Other two lists mimic the same behaviour
+            label i, j, a, b;
+
+            for (j = 1; j <= addedNeighbours.size() - 1; j++)
+            {
+                a = addedNeighbours[j];
+                b = shuffleList[j];
+
+                i = j - 1;
+
+                while (i >= 0 && addedNeighbours[i] > a)
+                {
+                    addedNeighbours[i + 1] = addedNeighbours[i];
+                    shuffleList[i + 1] = shuffleList[i];
+                    i--;
+                }
+
+                addedNeighbours[i + 1] = a;
+                shuffleList[i + 1] = b;
+            }
+
+            labelList reshuffleList(shuffleList.size());
+
+            forAll(shuffleList, i)
+            {
+                reshuffleList[shuffleList[i]] = i;
+            }
+
+            // Reshuffle other lists
+
+            labelListList fgn(faceGivingNeighbour.size());
+
+            forAll (faceGivingNeighbour, i)
+            {
+                fgn[reshuffleList[i]].transfer(faceGivingNeighbour[i].shrink());
+            }
+
+            labelListList eofgn(edgeOfFaceGivingNeighbour.size());
+
+            forAll (edgeOfFaceGivingNeighbour, i)
+            {
+                eofgn[reshuffleList[i]].transfer
+                (
+                    edgeOfFaceGivingNeighbour[i].shrink()
+                );
+            }
+
+            // adding the edges
+            forAll(addedNeighbours, edgeI)
+            {
+                const labelList& curFgn = fgn[edgeI];
+                const labelList& curEofgn = eofgn[edgeI];
+
+                forAll (curFgn, fgnI)
+                {
+                    fe[curFgn[fgnI]][curEofgn[fgnI]] = nEdges;
+                }
+
+                e[nEdges] = edge(pointI, addedNeighbours[edgeI]);
+                nEdges++;
             }
         }
 
-
-
-        // Renumber and transfer.
-
-        // Edges
-        edgesPtr_ = new edgeList(es.size());
-        edgeList& edges = *edgesPtr_;
-        forAll(es, edgeI)
-        {
-            edges[oldToNew[edgeI]] = es[edgeI];
-        }
-
-        // pointEdges
-        pePtr_ = new labelListList(nPoints());
-        labelListList& pointEdges = *pePtr_;
-        forAll(pe, pointI)
-        {
-            DynamicList<label>& pEdges = pe[pointI];
-            pEdges.shrink();
-            inplaceRenumber(oldToNew, pEdges);
-            pointEdges[pointI].transfer(pEdges);
-            Foam::sort(pointEdges[pointI]);
-        }
-
-        // faceEdges
-        if (doFaceEdges)
-        {
-            labelListList& faceEdges = *fePtr_;
-            forAll(faceEdges, faceI)
-            {
-                inplaceRenumber(oldToNew, faceEdges[faceI]);
-            }
-        }
+        // reset the size
+        e.setSize(nEdges);
     }
 }
 
@@ -505,70 +351,22 @@ const Foam::edgeList& Foam::primitiveMesh::edges() const
 {
     if (!edgesPtr_)
     {
-        //calcEdges(true);
-        calcEdges(false);
+        // 1.6.x merge: reverting to correct code
+        // HJ, 17/Aug/2010
+        calcEdges();
     }
 
     return *edgesPtr_;
 }
 
-const Foam::labelListList& Foam::primitiveMesh::pointEdges() const
-{
-    if (!pePtr_)
-    {
-        //calcEdges(true);
-        calcEdges(false);
-    }
-
-    return *pePtr_;
-}
-
 
 const Foam::labelListList& Foam::primitiveMesh::faceEdges() const
 {
     if (!fePtr_)
     {
-        if (debug)
-        {
-            Pout<< "primitiveMesh::faceEdges() : "
-                << "calculating faceEdges" << endl;
-        }
-
-        //calcEdges(true);
-        const faceList& fcs = faces();
-        const labelListList& pe = pointEdges();
-        const edgeList& es = edges();
-
-        fePtr_ = new labelListList(fcs.size());
-        labelListList& faceEdges = *fePtr_;
-
-        forAll(fcs, faceI)
-        {
-            const face& f = fcs[faceI];
-
-            labelList& fEdges = faceEdges[faceI];
-            fEdges.setSize(f.size());
-
-            forAll(f, fp)
-            {
-                label pointI = f[fp];
-                label nextPointI = f[f.fcIndex(fp)];
-
-                // Find edge between pointI, nextPontI
-                const labelList& pEdges = pe[pointI];
-
-                forAll(pEdges, i)
-                {
-                    label edgeI = pEdges[i];
-
-                    if (es[edgeI].otherVertex(pointI) == nextPointI)
-                    {
-                        fEdges[fp] = edgeI;
-                        break;
-                    }
-                }
-            }
-        }
+        // 1.6.x merge: reverting to correct code
+        // HJ, 17/Aug/2010
+        calcEdges();
     }
 
     return *fePtr_;
@@ -578,7 +376,6 @@ const Foam::labelListList& Foam::primitiveMesh::faceEdges() const
 void Foam::primitiveMesh::clearOutEdges()
 {
     deleteDemandDrivenData(edgesPtr_);
-    deleteDemandDrivenData(pePtr_);
     deleteDemandDrivenData(fePtr_);
     labels_.clear();
     labelSet_.clear();
@@ -678,6 +475,4 @@ const Foam::labelList& Foam::primitiveMesh::cellEdges(const label cellI) const
 }
 
 
-// * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
-
 // ************************************************************************* //
diff --git a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshI.H b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshI.H
index d98ca9a02..bd554f32d 100644
--- a/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshI.H
+++ b/src/OpenFOAM/meshes/primitiveMesh/primitiveMeshI.H
@@ -33,42 +33,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
-inline label primitiveMesh::nInternalPoints() const
-{
-    return nInternalPoints_;
-}
-
-
 inline label primitiveMesh::nPoints() const
 {
     return nPoints_;
 }
 
 
-inline label primitiveMesh::nInternal0Edges() const
-{
-    // Force edge calculation
-    nEdges();
-    return nInternal0Edges_;
-}
-
-
-inline label primitiveMesh::nInternal1Edges() const
-{
-    // Force edge calculation
-    nEdges();
-    return nInternal1Edges_;
-}
-
-
-inline label primitiveMesh::nInternalEdges() const
-{
-    // Force edge calculation
-    nEdges();
-    return nInternalEdges_;
-}
-
-
 inline label primitiveMesh::nEdges() const
 {
     if (nEdges_ < 0)
diff --git a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
index 98df63f21..035fb0052 100644
--- a/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
+++ b/src/OpenFOAM/meshes/primitiveShapes/triangle/triangleI.H
@@ -624,46 +624,45 @@ inline pointHit triangle<Point, PointRef>::fastIntersection
             // Ray parallel to triangle. Return miss
             return intersection;
         }
-    }
 
-    const scalar inv_det = 1/det;
+        const scalar inv_det = 1/det;
 
         /* calculate distance from a_ to ray origin */
         const vector tVec = orig - a_;
 
-    /* calculate U parameter and test bounds */
-    const scalar u = (tVec & pVec)*inv_det;
+        /* calculate U parameter and test bounds */
+        const scalar u = (tVec & pVec)*inv_det;
 
-    if (u < -tol || u > 1.0 + tol)
-    {
-        // return miss
-        return intersection;
-    }
+        if (u < -tol || u > 1.0 + tol)
+        {
+            // return miss
+            return intersection;
+        }
 
-    /* prepare to test V parameter */
-    const vector qVec = tVec ^ edge1;
+        /* prepare to test V parameter */
+        const vector qVec = tVec ^ edge1;
 
-    /* calculate V parameter and test bounds */
-    const scalar v = (dir & qVec) * inv_det;
+        /* calculate V parameter and test bounds */
+        const scalar v = (dir & qVec) * inv_det;
 
-    if (v < -tol || u + v > 1.0 + tol)
-    {
-        // return miss
-        return intersection;
-    }
+        if (v < -tol || u + v > 1.0 + tol)
+        {
+            // return miss
+            return intersection;
+        }
 
-    /* calculate t, scale parameters, ray intersects triangle */
-    const scalar t = (edge2 & qVec)*inv_det;
+        /* calculate t, scale parameters, ray intersects triangle */
+        const scalar t = (edge2 & qVec)*inv_det;
 
-    if (alg == intersection::HALF_RAY && t < -tol)
-    {
-        // Wrong side of orig. Return miss
-        return intersection;
-    }
+        if (alg == intersection::HALF_RAY && t < -tol)
+        {
+            // Wrong side of orig. Return miss
+            return intersection;
+        }
 
-    intersection.setHit();
-    intersection.setPoint(a_ + u*edge1 + v*edge2);
-    intersection.setDistance(t);
+        intersection.setHit();
+        intersection.setPoint(a_ + u*edge1 + v*edge2);
+        intersection.setDistance(t);
     }
     else
     {
diff --git a/src/OpenFOAM/primitives/Scalar/Scalar.H b/src/OpenFOAM/primitives/Scalar/Scalar.H
index 4a91584ac..18bd14b32 100644
--- a/src/OpenFOAM/primitives/Scalar/Scalar.H
+++ b/src/OpenFOAM/primitives/Scalar/Scalar.H
@@ -208,6 +208,11 @@ inline Scalar cmptDivide(const Scalar s1, const Scalar s2)
     return s1/s2;
 }
 
+inline Scalar cmptSumMultiply(const Scalar s1, const Scalar s2)
+{
+    return s1*s2;
+}
+
 inline Scalar cmptMax(const Scalar s)
 {
     return s;
@@ -275,6 +280,24 @@ inline Scalar stabilise(const Scalar s, const Scalar small)
 }
 
 
+inline Scalar cmptStabilise
+(
+    const Scalar s,
+    const Scalar small,
+    const Scalar value
+)
+{
+    if (mag(s) < small)
+    {
+        return sign(s)*value;
+    }
+    else
+    {
+        return s;
+    }
+}
+
+
 // * * * * * * * * * * * * * * * IOstream Operators  * * * * * * * * * * * * //
 
 Scalar readScalar(Istream&);
diff --git a/wmake/src/lex.yy.c b/wmake/src/lex.yy.c
deleted file mode 100644
index 5f3ae14ac..000000000
--- a/wmake/src/lex.yy.c
+++ /dev/null
@@ -1,2394 +0,0 @@
-
-#line 3 "lex.yy.c"
-
-#define  YY_INT_ALIGNED short int
-
-/* A lexical scanner generated by flex */
-
-#define FLEX_SCANNER
-#define YY_FLEX_MAJOR_VERSION 2
-#define YY_FLEX_MINOR_VERSION 5
-#define YY_FLEX_SUBMINOR_VERSION 33
-#if YY_FLEX_SUBMINOR_VERSION > 0
-#define FLEX_BETA
-#endif
-
-/* First, we deal with  platform-specific or compiler-specific issues. */
-
-/* begin standard C headers. */
-#include <stdio.h>
-#include <string.h>
-#include <errno.h>
-#include <stdlib.h>
-
-/* end standard C headers. */
-
-/* flex integer type definitions */
-
-#ifndef FLEXINT_H
-#define FLEXINT_H
-
-/* C99 systems have <inttypes.h>. Non-C99 systems may or may not. */
-
-#if __STDC_VERSION__ >= 199901L
-
-/* C99 says to define __STDC_LIMIT_MACROS before including stdint.h,
- * if you want the limit (max/min) macros for int types. 
- */
-#ifndef __STDC_LIMIT_MACROS
-#define __STDC_LIMIT_MACROS 1
-#endif
-
-#include <inttypes.h>
-typedef int8_t flex_int8_t;
-typedef uint8_t flex_uint8_t;
-typedef int16_t flex_int16_t;
-typedef uint16_t flex_uint16_t;
-typedef int32_t flex_int32_t;
-typedef uint32_t flex_uint32_t;
-#else
-typedef signed char flex_int8_t;
-typedef short int flex_int16_t;
-typedef int flex_int32_t;
-typedef unsigned char flex_uint8_t; 
-typedef unsigned short int flex_uint16_t;
-typedef unsigned int flex_uint32_t;
-#endif /* ! C99 */
-
-/* Limits of integral types. */
-#ifndef INT8_MIN
-#define INT8_MIN               (-128)
-#endif
-#ifndef INT16_MIN
-#define INT16_MIN              (-32767-1)
-#endif
-#ifndef INT32_MIN
-#define INT32_MIN              (-2147483647-1)
-#endif
-#ifndef INT8_MAX
-#define INT8_MAX               (127)
-#endif
-#ifndef INT16_MAX
-#define INT16_MAX              (32767)
-#endif
-#ifndef INT32_MAX
-#define INT32_MAX              (2147483647)
-#endif
-#ifndef UINT8_MAX
-#define UINT8_MAX              (255U)
-#endif
-#ifndef UINT16_MAX
-#define UINT16_MAX             (65535U)
-#endif
-#ifndef UINT32_MAX
-#define UINT32_MAX             (4294967295U)
-#endif
-
-#endif /* ! FLEXINT_H */
-
-#ifdef __cplusplus
-
-/* The "const" storage-class-modifier is valid. */
-#define YY_USE_CONST
-
-#else	/* ! __cplusplus */
-
-#if __STDC__
-
-#define YY_USE_CONST
-
-#endif	/* __STDC__ */
-#endif	/* ! __cplusplus */
-
-#ifdef YY_USE_CONST
-#define yyconst const
-#else
-#define yyconst
-#endif
-
-/* Returned upon end-of-file. */
-#define YY_NULL 0
-
-/* Promotes a possibly negative, possibly signed char to an unsigned
- * integer for use as an array index.  If the signed char is negative,
- * we want to instead treat it as an 8-bit unsigned char, hence the
- * double cast.
- */
-#define YY_SC_TO_UI(c) ((unsigned int) (unsigned char) c)
-
-/* Enter a start condition.  This macro really ought to take a parameter,
- * but we do it the disgusting crufty way forced on us by the ()-less
- * definition of BEGIN.
- */
-#define BEGIN (yy_start) = 1 + 2 *
-
-/* Translate the current start state into a value that can be later handed
- * to BEGIN to return to the state.  The YYSTATE alias is for lex
- * compatibility.
- */
-#define YY_START (((yy_start) - 1) / 2)
-#define YYSTATE YY_START
-
-/* Action number for EOF rule of a given start state. */
-#define YY_STATE_EOF(state) (YY_END_OF_BUFFER + state + 1)
-
-/* Special action meaning "start processing a new file". */
-#define YY_NEW_FILE yyrestart(yyin  )
-
-#define YY_END_OF_BUFFER_CHAR 0
-
-/* Size of default input buffer. */
-#ifndef YY_BUF_SIZE
-#define YY_BUF_SIZE 16384
-#endif
-
-/* The state buf must be large enough to hold one state per character in the main buffer.
- */
-#define YY_STATE_BUF_SIZE   ((YY_BUF_SIZE + 2) * sizeof(yy_state_type))
-
-#ifndef YY_TYPEDEF_YY_BUFFER_STATE
-#define YY_TYPEDEF_YY_BUFFER_STATE
-typedef struct yy_buffer_state *YY_BUFFER_STATE;
-#endif
-
-extern int yyleng;
-
-extern FILE *yyin, *yyout;
-
-#define EOB_ACT_CONTINUE_SCAN 0
-#define EOB_ACT_END_OF_FILE 1
-#define EOB_ACT_LAST_MATCH 2
-
-    #define YY_LESS_LINENO(n)
-    
-/* Return all but the first "n" matched characters back to the input stream. */
-#define yyless(n) \
-	do \
-		{ \
-		/* Undo effects of setting up yytext. */ \
-        int yyless_macro_arg = (n); \
-        YY_LESS_LINENO(yyless_macro_arg);\
-		*yy_cp = (yy_hold_char); \
-		YY_RESTORE_YY_MORE_OFFSET \
-		(yy_c_buf_p) = yy_cp = yy_bp + yyless_macro_arg - YY_MORE_ADJ; \
-		YY_DO_BEFORE_ACTION; /* set up yytext again */ \
-		} \
-	while ( 0 )
-
-#define unput(c) yyunput( c, (yytext_ptr)  )
-
-/* The following is because we cannot portably get our hands on size_t
- * (without autoconf's help, which isn't available because we want
- * flex-generated scanners to compile on their own).
- */
-
-#ifndef YY_TYPEDEF_YY_SIZE_T
-#define YY_TYPEDEF_YY_SIZE_T
-typedef unsigned int yy_size_t;
-#endif
-
-#ifndef YY_STRUCT_YY_BUFFER_STATE
-#define YY_STRUCT_YY_BUFFER_STATE
-struct yy_buffer_state
-	{
-	FILE *yy_input_file;
-
-	char *yy_ch_buf;		/* input buffer */
-	char *yy_buf_pos;		/* current position in input buffer */
-
-	/* Size of input buffer in bytes, not including room for EOB
-	 * characters.
-	 */
-	yy_size_t yy_buf_size;
-
-	/* Number of characters read into yy_ch_buf, not including EOB
-	 * characters.
-	 */
-	int yy_n_chars;
-
-	/* Whether we "own" the buffer - i.e., we know we created it,
-	 * and can realloc() it to grow it, and should free() it to
-	 * delete it.
-	 */
-	int yy_is_our_buffer;
-
-	/* Whether this is an "interactive" input source; if so, and
-	 * if we're using stdio for input, then we want to use getc()
-	 * instead of fread(), to make sure we stop fetching input after
-	 * each newline.
-	 */
-	int yy_is_interactive;
-
-	/* Whether we're considered to be at the beginning of a line.
-	 * If so, '^' rules will be active on the next match, otherwise
-	 * not.
-	 */
-	int yy_at_bol;
-
-    int yy_bs_lineno; /**< The line count. */
-    int yy_bs_column; /**< The column count. */
-    
-	/* Whether to try to fill the input buffer when we reach the
-	 * end of it.
-	 */
-	int yy_fill_buffer;
-
-	int yy_buffer_status;
-
-#define YY_BUFFER_NEW 0
-#define YY_BUFFER_NORMAL 1
-	/* When an EOF's been seen but there's still some text to process
-	 * then we mark the buffer as YY_EOF_PENDING, to indicate that we
-	 * shouldn't try reading from the input source any more.  We might
-	 * still have a bunch of tokens to match, though, because of
-	 * possible backing-up.
-	 *
-	 * When we actually see the EOF, we change the status to "new"
-	 * (via yyrestart()), so that the user can continue scanning by
-	 * just pointing yyin at a new input file.
-	 */
-#define YY_BUFFER_EOF_PENDING 2
-
-	};
-#endif /* !YY_STRUCT_YY_BUFFER_STATE */
-
-/* Stack of input buffers. */
-static size_t yy_buffer_stack_top = 0; /**< index of top of stack. */
-static size_t yy_buffer_stack_max = 0; /**< capacity of stack. */
-static YY_BUFFER_STATE * yy_buffer_stack = 0; /**< Stack as an array. */
-
-/* We provide macros for accessing buffer states in case in the
- * future we want to put the buffer states in a more general
- * "scanner state".
- *
- * Returns the top of the stack, or NULL.
- */
-#define YY_CURRENT_BUFFER ( (yy_buffer_stack) \
-                          ? (yy_buffer_stack)[(yy_buffer_stack_top)] \
-                          : NULL)
-
-/* Same as previous macro, but useful when we know that the buffer stack is not
- * NULL or when we need an lvalue. For internal use only.
- */
-#define YY_CURRENT_BUFFER_LVALUE (yy_buffer_stack)[(yy_buffer_stack_top)]
-
-/* yy_hold_char holds the character lost when yytext is formed. */
-static char yy_hold_char;
-static int yy_n_chars;		/* number of characters read into yy_ch_buf */
-int yyleng;
-
-/* Points to current character in buffer. */
-static char *yy_c_buf_p = (char *) 0;
-static int yy_init = 0;		/* whether we need to initialize */
-static int yy_start = 0;	/* start state number */
-
-/* Flag which is used to allow yywrap()'s to do buffer switches
- * instead of setting up a fresh yyin.  A bit of a hack ...
- */
-static int yy_did_buffer_switch_on_eof;
-
-void yyrestart (FILE *input_file  );
-void yy_switch_to_buffer (YY_BUFFER_STATE new_buffer  );
-YY_BUFFER_STATE yy_create_buffer (FILE *file,int size  );
-void yy_delete_buffer (YY_BUFFER_STATE b  );
-void yy_flush_buffer (YY_BUFFER_STATE b  );
-void yypush_buffer_state (YY_BUFFER_STATE new_buffer  );
-void yypop_buffer_state (void );
-
-static void yyensure_buffer_stack (void );
-static void yy_load_buffer_state (void );
-static void yy_init_buffer (YY_BUFFER_STATE b,FILE *file  );
-
-#define YY_FLUSH_BUFFER yy_flush_buffer(YY_CURRENT_BUFFER )
-
-YY_BUFFER_STATE yy_scan_buffer (char *base,yy_size_t size  );
-YY_BUFFER_STATE yy_scan_string (yyconst char *yy_str  );
-YY_BUFFER_STATE yy_scan_bytes (yyconst char *bytes,int len  );
-
-void *yyalloc (yy_size_t  );
-void *yyrealloc (void *,yy_size_t  );
-void yyfree (void *  );
-
-#define yy_new_buffer yy_create_buffer
-
-#define yy_set_interactive(is_interactive) \
-	{ \
-	if ( ! YY_CURRENT_BUFFER ){ \
-        yyensure_buffer_stack (); \
-		YY_CURRENT_BUFFER_LVALUE =    \
-            yy_create_buffer(yyin,YY_BUF_SIZE ); \
-	} \
-	YY_CURRENT_BUFFER_LVALUE->yy_is_interactive = is_interactive; \
-	}
-
-#define yy_set_bol(at_bol) \
-	{ \
-	if ( ! YY_CURRENT_BUFFER ){\
-        yyensure_buffer_stack (); \
-		YY_CURRENT_BUFFER_LVALUE =    \
-            yy_create_buffer(yyin,YY_BUF_SIZE ); \
-	} \
-	YY_CURRENT_BUFFER_LVALUE->yy_at_bol = at_bol; \
-	}
-
-#define YY_AT_BOL() (YY_CURRENT_BUFFER_LVALUE->yy_at_bol)
-
-/* Begin user sect3 */
-
-typedef unsigned char YY_CHAR;
-
-FILE *yyin = (FILE *) 0, *yyout = (FILE *) 0;
-
-typedef int yy_state_type;
-
-extern int yylineno;
-
-int yylineno = 1;
-
-extern char *yytext;
-#define yytext_ptr yytext
-
-static yy_state_type yy_get_previous_state (void );
-static yy_state_type yy_try_NUL_trans (yy_state_type current_state  );
-static int yy_get_next_buffer (void );
-static void yy_fatal_error (yyconst char msg[]  );
-
-/* Done after the current pattern has been matched and before the
- * corresponding action - sets up yytext.
- */
-#define YY_DO_BEFORE_ACTION \
-	(yytext_ptr) = yy_bp; \
-	yyleng = (size_t) (yy_cp - yy_bp); \
-	(yy_hold_char) = *yy_cp; \
-	*yy_cp = '\0'; \
-	(yy_c_buf_p) = yy_cp;
-
-#define YY_NUM_RULES 19
-#define YY_END_OF_BUFFER 20
-/* This struct is not used in this scanner,
-   but its presence is necessary. */
-struct yy_trans_info
-	{
-	flex_int32_t yy_verify;
-	flex_int32_t yy_nxt;
-	};
-static yyconst flex_int16_t yy_accept[152] =
-    {   0,
-        0,    0,    0,    0,    6,    6,    0,    0,    0,    0,
-       17,   17,   20,   18,   18,   18,   18,   18,   18,   18,
-       18,    3,    3,    6,   19,   19,   19,   14,   19,   19,
-       19,   19,   19,   17,    0,    2,    0,    0,    0,    0,
-        0,    0,    0,    0,    4,    6,    0,    0,   14,    0,
-       15,    0,    0,    0,    0,    0,   17,    0,    0,    1,
-        0,    0,    0,    0,    0,    0,    0,    0,    0,    0,
-        0,    0,    0,    0,    0,   10,    0,    0,    0,    0,
-        0,    9,    0,    0,    0,   12,    0,    0,    0,    0,
-        0,   10,    0,   10,    0,    0,    0,    8,    0,    9,
-
-        0,    9,   12,    0,   12,    0,    0,    0,   11,    0,
-        0,    0,    0,   10,    8,    0,    8,    0,    9,   12,
-       11,    0,   11,    0,    0,    0,    7,    8,    0,   11,
-        0,    0,    0,    0,    0,    0,    0,    0,   13,    0,
-        5,   13,    0,   13,    0,   13,    0,    0,    0,   16,
-        0
-    } ;
-
-static yyconst flex_int32_t yy_ec[256] =
-    {   0,
-        1,    1,    1,    1,    1,    1,    1,    1,    2,    3,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    4,    1,    5,    6,    1,    1,    1,    7,    1,
-        1,    8,    1,    1,    1,    1,    9,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,   10,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,   11,    1,   12,   13,
-
-       14,    1,   15,   16,   17,   18,    1,   19,   20,   21,
-       22,   23,    1,   24,   25,   26,   27,   28,   29,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1
-    } ;
-
-static yyconst flex_int32_t yy_meta[30] =
-    {   0,
-        1,    1,    2,    3,    4,    1,    5,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1,    1,
-        1,    1,    1,    1,    1,    1,    1,    1,    1
-    } ;
-
-static yyconst flex_int16_t yy_base[181] =
-    {   0,
-        0,    8,  353,  351,    2,   15,    0,    0,   23,   46,
-      349,  348,  354,  467,  467,  317,    7,   30,   53,   20,
-      293,  467,  300,    0,  467,    3,   13,   19,   30,   52,
-       64,   59,   68,    0,  301,  467,  298,   75,   80,  268,
-       83,   84,  261,  253,  467,    0,   36,   83,   86,  262,
-      245,   94,   95,   98,  101,  105,    0,  245,  242,  467,
-      115,  229,  213,  102,  119,  110,  128,  129,  224,  138,
-      202,  191,  138,  204,  142,  146,  152,  151,  159,  198,
-      160,  166,  167,  189,  173,  174,  180,  190,  185,  159,
-      154,  163,  187,  188,  195,  156,  196,  202,  201,  143,
-
-      208,  214,  139,  215,  221,  222,  135,  228,  229,  124,
-      238,  122,  120,  238,  118,  242,  246,  249,  256,  257,
-       80,  263,  264,   59,   78,   59,  256,  270,  269,  276,
-       54,  281,  282,   42,  292,  290,   43,  294,  298,   23,
-      467,   33,  302,  306,   26,  310,   16,  287,  315,  467,
-      467,  322,  327,  332,  337,  342,  347,  352,  357,  362,
-      367,  371,  376,  381,  386,  391,  396,  401,  406,  411,
-      416,  421,  426,  431,  436,  441,  446,  451,  456,  461
-    } ;
-
-static yyconst flex_int16_t yy_def[181] =
-    {   0,
-      152,  152,  153,  153,  154,  154,  155,  155,  156,  156,
-      157,  157,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  158,  151,  159,  160,  159,  159,  159,
-      159,  159,  159,  161,  151,  151,  162,  151,  151,  151,
-      151,  151,  151,  151,  151,  158,  159,  160,  159,  160,
-      160,  159,  159,  159,  159,  159,  161,  151,  162,  151,
-      151,  151,  151,  163,  159,  159,  164,  165,  151,  151,
-      151,  151,  163,  166,  167,  163,  168,  159,  164,  169,
-      170,  164,  165,  171,  172,  165,  173,  151,  151,  151,
-      151,  166,  167,  167,  168,  174,  175,  168,  159,  169,
-
-      170,  170,  171,  172,  172,  173,  176,  177,  173,  151,
-      151,  151,  151,  167,  174,  175,  175,  159,  170,  172,
-      176,  177,  177,  151,  151,  151,  151,  175,  159,  177,
-      151,  151,  178,  151,  151,  178,  179,  180,  178,  151,
-      151,  179,  180,  180,  151,  180,  151,  151,  151,  151,
-        0,  151,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151,  151,  151,  151,  151
-    } ;
-
-static yyconst flex_int16_t yy_nxt[497] =
-    {   0,
-      151,   15,  151,   16,   25,   25,   25,  151,   17,   18,
-       48,   19,   49,   20,   36,   37,   17,   25,   25,   25,
-       48,   42,   51,   42,   21,   25,   48,   25,   49,  148,
-       27,   38,   28,   38,   29,   39,   43,   48,  147,   49,
-       30,   31,  142,   48,   32,   49,   40,   33,   25,  145,
-       25,   52,  142,   27,   38,   28,   41,   29,   39,   48,
-      140,   49,   53,   30,   31,  134,   48,   32,   49,   40,
-       33,   48,  132,   49,   54,   48,   38,   49,   38,  131,
-       39,   42,   55,   42,   38,   42,   61,   42,   39,  121,
-       48,   40,   51,   48,   56,   49,   43,   44,  131,   40,
-
-       43,   48,   48,   49,   49,   48,   74,   49,   48,   75,
-       49,   76,   48,   64,   49,   67,   38,   48,   70,   49,
-       39,  127,   65,  127,   66,   68,   48,  115,   49,   77,
-       78,   40,   80,   84,  126,   81,   85,   82,   86,   38,
-      124,   89,   74,   39,  121,   75,   74,   76,  103,   75,
-       74,   94,  100,   75,   40,   76,   96,   87,   48,   97,
-       49,   98,   99,   80,   80,  115,   81,   81,   82,  102,
-       80,   84,   92,   81,   85,   82,   86,   84,   84,  113,
-       85,   85,  105,   86,  107,  112,   38,  108,  111,  109,
-       39,   74,   74,  110,   75,   75,  114,  114,  103,   96,
-
-       96,   40,   97,   97,   98,  117,   96,  100,   48,   97,
-       49,   98,   80,   92,   91,   81,  118,  119,   80,   84,
-       90,   81,   85,  119,  120,   84,  107,   88,   85,  108,
-      120,  109,  107,  107,   72,  108,  108,  123,  109,   38,
-       71,   38,   74,   39,   60,   75,   96,  114,   69,   97,
-       96,  128,   48,   97,  125,  128,   48,  127,   49,  127,
-       80,   84,  129,   81,   85,  119,  120,  107,  107,   48,
-      108,  108,  130,  130,   96,   63,   48,   97,   49,  128,
-      107,   62,  135,  108,  135,  130,  137,   44,  149,  138,
-      149,  139,  133,  135,  137,  135,  141,  138,  137,  139,
-
-       60,  138,  137,  144,   58,  138,  137,  139,   45,  138,
-      137,  146,   44,  138,  137,  146,  149,  138,  149,  146,
-       35,  150,   14,   14,   14,   14,   14,   22,   22,   22,
-       22,   22,   24,   24,   24,   24,   24,   25,   25,   25,
-       25,   25,   26,   26,   26,   26,   26,   34,   34,   34,
-       34,   34,   46,  151,   25,   25,   46,   47,   23,   47,
-       23,   47,   50,  151,   50,  151,   50,   57,   57,   57,
-       57,   59,   59,   59,   59,   59,   73,  151,   73,   73,
-       73,   79,  151,   79,   79,   79,   83,  151,   83,   83,
-       83,   74,  151,   74,   74,   74,   93,  151,   93,   93,
-
-       93,   95,  151,   95,   95,   95,   80,  151,   80,   80,
-       80,  101,  151,  101,  101,  101,   84,  151,   84,   84,
-       84,  104,  151,  104,  104,  104,  106,  151,  106,  106,
-      106,   96,  151,   96,   96,   96,  116,  151,  116,  116,
-      116,  107,  151,  107,  107,  107,  122,  151,  122,  122,
-      122,  136,  151,  136,  136,  136,  137,  151,  137,  137,
-      137,  143,  151,  143,  143,  143,   13,  151,  151,  151,
-      151,  151,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151
-
-    } ;
-
-static yyconst flex_int16_t yy_chk[497] =
-    {   0,
-        0,    1,    0,    1,    5,    5,    5,    0,    1,    2,
-       26,    2,   26,    2,   17,   17,    2,    6,    6,    6,
-       27,   20,   27,   20,    2,    9,   28,    9,   28,  147,
-        9,   18,    9,   18,    9,   18,   20,   29,  145,   29,
-        9,    9,  142,   47,    9,   47,   18,    9,   10,  140,
-       10,   29,  137,   10,   19,   10,   19,   10,   19,   30,
-      134,   30,   30,   10,   10,  131,   32,   10,   32,   19,
-       10,   31,  126,   31,   31,   33,   38,   33,   38,  124,
-       38,   39,   32,   39,   41,   42,   41,   42,   41,  121,
-       48,   38,   48,   49,   33,   49,   39,  125,  125,   41,
-
-       42,   52,   53,   52,   53,   54,   64,   54,   55,   64,
-       55,   64,   56,   52,   56,   55,   61,   66,   61,   66,
-       61,  113,   53,  113,   54,   56,   65,  115,   65,   65,
-       66,   61,   67,   68,  112,   67,   68,   67,   68,   70,
-      110,   70,   73,   70,  107,   73,   75,   73,  103,   75,
-       76,   75,  100,   76,   70,   76,   77,   68,   78,   77,
-       78,   77,   78,   79,   81,   96,   79,   81,   79,   81,
-       82,   83,   92,   82,   83,   82,   83,   85,   86,   91,
-       85,   86,   85,   86,   87,   90,   89,   87,   89,   87,
-       89,   93,   94,   88,   93,   94,   93,   94,   84,   95,
-
-       97,   89,   95,   97,   95,   97,   98,   80,   99,   98,
-       99,   98,  101,   74,   72,  101,   99,  101,  102,  104,
-       71,  102,  104,  102,  104,  105,  106,   69,  105,  106,
-      105,  106,  108,  109,   63,  108,  109,  108,  109,  111,
-       62,  111,  114,  111,   59,  114,  116,  114,   58,  116,
-      117,  116,   51,  117,  111,  117,  118,  127,  118,  127,
-      119,  120,  118,  119,  120,  119,  120,  122,  123,   50,
-      122,  123,  122,  123,  128,   44,  129,  128,  129,  128,
-      130,   43,  132,  130,  132,  130,  133,   40,  148,  133,
-      148,  133,  129,  135,  136,  135,  135,  136,  138,  136,
-
-       37,  138,  139,  138,   35,  139,  143,  139,   23,  143,
-      144,  143,   21,  144,  146,  144,  149,  146,  149,  146,
-       16,  149,  152,  152,  152,  152,  152,  153,  153,  153,
-      153,  153,  154,  154,  154,  154,  154,  155,  155,  155,
-      155,  155,  156,  156,  156,  156,  156,  157,  157,  157,
-      157,  157,  158,   13,   12,   11,  158,  159,    4,  159,
-        3,  159,  160,    0,  160,    0,  160,  161,  161,  161,
-      161,  162,  162,  162,  162,  162,  163,    0,  163,  163,
-      163,  164,    0,  164,  164,  164,  165,    0,  165,  165,
-      165,  166,    0,  166,  166,  166,  167,    0,  167,  167,
-
-      167,  168,    0,  168,  168,  168,  169,    0,  169,  169,
-      169,  170,    0,  170,  170,  170,  171,    0,  171,  171,
-      171,  172,    0,  172,  172,  172,  173,    0,  173,  173,
-      173,  174,    0,  174,  174,  174,  175,    0,  175,  175,
-      175,  176,    0,  176,  176,  176,  177,    0,  177,  177,
-      177,  178,    0,  178,  178,  178,  179,    0,  179,  179,
-      179,  180,    0,  180,  180,  180,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151,  151,  151,  151,  151,
-      151,  151,  151,  151,  151,  151
-
-    } ;
-
-static yy_state_type yy_last_accepting_state;
-static char *yy_last_accepting_cpos;
-
-extern int yy_flex_debug;
-int yy_flex_debug = 0;
-
-/* The intent behind this definition is that it'll catch
- * any uses of REJECT which flex missed.
- */
-#define REJECT reject_used_but_not_detected
-#define yymore() yymore_used_but_not_detected
-#define YY_MORE_ADJ 0
-#define YY_RESTORE_YY_MORE_OFFSET
-char *yytext;
-#line 1 "wmkdep.l"
-#line 2 "wmkdep.l"
-/*---------------------------------------------------------------------------*\
- =========                   |
- \\      /   F ield          | OpenFOAM: The Open Source CFD Toolbox
-  \\    /    O peration      |
-   \\  /     A nd            | Copyright held by original author
-    \\/      M anipulation   |
-------------------------------------------------------------------------------
-License
-    This file is part of OpenFOAM.
-
-    OpenFOAM is free software; you can redistribute it and/or modify it
-    under the terms of the GNU General Public License as published by the
-    Free Software Foundation; either version 2 of the License, or (at your
-    option) any later version.
-
-    OpenFOAM is distributed in the hope that it will be useful, but WITHOUT
-    ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or
-    FITNESS FOR A PARTICULAR PURPOSE.  See the GNU General Public License
-    for more details.
-
-    You should have received a copy of the GNU General Public License
-    along with OpenFOAM; if not, write to the Free Software Foundation,
-    Inc., 51 Franklin St, Fifth Floor, Boston, MA 02110-1301 USA
-
-Application
-    wmkdep
-
-Description
-     A fast dependency list generator which emulates the behaviour and
-     output of cpp -M. However, the output contains no duplications and 
-     is ~40% faster than cpp.
-
-     The algorithm uses flex to scan for includes and searches the files
-     found.  Each file is entererd into a hash table so that files are
-     scaned only once.  This is why this program is faster than cpp.
-
-Usage
-    wmkdep [ -Idirectory ... -Idirectory] filename
-
-\*---------------------------------------------------------------------------*/
-
-#define FILE_STACK_SIZE 300
-#define HASH_TABLE_SIZE 500
-
-#include <stdlib.h>
-#include <string.h>
-#include <sys/types.h>
-#include <dirent.h>
-
-void nextFile(const char* fileName);
-void importFile(const char* fileName);
-void importDirectory(const char* dirName);
-
-#   undef yywrap	/* sometimes a macro by default */
-
-
-#line 672 "lex.yy.c"
-
-#define INITIAL 0
-#define CMNT 1
-#define CFNAME 2
-#define SCFNAME 3
-#define JFNAME 4
-#define FFNAME 5
-
-#ifndef YY_NO_UNISTD_H
-/* Special case for "unistd.h", since it is non-ANSI. We include it way
- * down here because we want the user's section 1 to have been scanned first.
- * The user has a chance to override it with an option.
- */
-#include <unistd.h>
-#endif
-
-#ifndef YY_EXTRA_TYPE
-#define YY_EXTRA_TYPE void *
-#endif
-
-static int yy_init_globals (void );
-
-/* Macros after this point can all be overridden by user definitions in
- * section 1.
- */
-
-#ifndef YY_SKIP_YYWRAP
-#ifdef __cplusplus
-extern "C" int yywrap (void );
-#else
-extern int yywrap (void );
-#endif
-#endif
-
-    static void yyunput (int c,char *buf_ptr  );
-    
-#ifndef yytext_ptr
-static void yy_flex_strncpy (char *,yyconst char *,int );
-#endif
-
-#ifdef YY_NEED_STRLEN
-static int yy_flex_strlen (yyconst char * );
-#endif
-
-#ifndef YY_NO_INPUT
-
-#ifdef __cplusplus
-static int yyinput (void );
-#else
-static int input (void );
-#endif
-
-#endif
-
-/* Amount of stuff to slurp up with each read. */
-#ifndef YY_READ_BUF_SIZE
-#define YY_READ_BUF_SIZE 8192
-#endif
-
-/* Copy whatever the last rule matched to the standard output. */
-#ifndef ECHO
-/* This used to be an fputs(), but since the string might contain NUL's,
- * we now use fwrite().
- */
-#define ECHO (void) fwrite( yytext, yyleng, 1, yyout )
-#endif
-
-/* Gets input and stuffs it into "buf".  number of characters read, or YY_NULL,
- * is returned in "result".
- */
-#ifndef YY_INPUT
-#define YY_INPUT(buf,result,max_size) \
-	if ( YY_CURRENT_BUFFER_LVALUE->yy_is_interactive ) \
-		{ \
-		int c = '*'; \
-		size_t n; \
-		for ( n = 0; n < max_size && \
-			     (c = getc( yyin )) != EOF && c != '\n'; ++n ) \
-			buf[n] = (char) c; \
-		if ( c == '\n' ) \
-			buf[n++] = (char) c; \
-		if ( c == EOF && ferror( yyin ) ) \
-			YY_FATAL_ERROR( "input in flex scanner failed" ); \
-		result = n; \
-		} \
-	else \
-		{ \
-		errno=0; \
-		while ( (result = fread(buf, 1, max_size, yyin))==0 && ferror(yyin)) \
-			{ \
-			if( errno != EINTR) \
-				{ \
-				YY_FATAL_ERROR( "input in flex scanner failed" ); \
-				break; \
-				} \
-			errno=0; \
-			clearerr(yyin); \
-			} \
-		}\
-\
-
-#endif
-
-/* No semi-colon after return; correct usage is to write "yyterminate();" -
- * we don't want an extra ';' after the "return" because that will cause
- * some compilers to complain about unreachable statements.
- */
-#ifndef yyterminate
-#define yyterminate() return YY_NULL
-#endif
-
-/* Number of entries by which start-condition stack grows. */
-#ifndef YY_START_STACK_INCR
-#define YY_START_STACK_INCR 25
-#endif
-
-/* Report a fatal error. */
-#ifndef YY_FATAL_ERROR
-#define YY_FATAL_ERROR(msg) yy_fatal_error( msg )
-#endif
-
-/* end tables serialization structures and prototypes */
-
-/* Default declaration of generated scanner - a define so the user can
- * easily add parameters.
- */
-#ifndef YY_DECL
-#define YY_DECL_IS_OURS 1
-
-extern int yylex (void);
-
-#define YY_DECL int yylex (void)
-#endif /* !YY_DECL */
-
-/* Code executed at the beginning of each rule, after yytext and yyleng
- * have been set up.
- */
-#ifndef YY_USER_ACTION
-#define YY_USER_ACTION
-#endif
-
-/* Code executed at the end of each rule. */
-#ifndef YY_BREAK
-#define YY_BREAK break;
-#endif
-
-#define YY_RULE_SETUP \
-	if ( yyleng > 0 ) \
-		YY_CURRENT_BUFFER_LVALUE->yy_at_bol = \
-				(yytext[yyleng - 1] == '\n'); \
-	YY_USER_ACTION
-
-/** The main scanner function which does all the work.
- */
-YY_DECL
-{
-	register yy_state_type yy_current_state;
-	register char *yy_cp, *yy_bp;
-	register int yy_act;
-    
-#line 60 "wmkdep.l"
-
-
-#line 836 "lex.yy.c"
-
-	if ( !(yy_init) )
-		{
-		(yy_init) = 1;
-
-#ifdef YY_USER_INIT
-		YY_USER_INIT;
-#endif
-
-		if ( ! (yy_start) )
-			(yy_start) = 1;	/* first start state */
-
-		if ( ! yyin )
-			yyin = stdin;
-
-		if ( ! yyout )
-			yyout = stdout;
-
-		if ( ! YY_CURRENT_BUFFER ) {
-			yyensure_buffer_stack ();
-			YY_CURRENT_BUFFER_LVALUE =
-				yy_create_buffer(yyin,YY_BUF_SIZE );
-		}
-
-		yy_load_buffer_state( );
-		}
-
-	while ( 1 )		/* loops until end-of-file is reached */
-		{
-		yy_cp = (yy_c_buf_p);
-
-		/* Support of yytext. */
-		*yy_cp = (yy_hold_char);
-
-		/* yy_bp points to the position in yy_ch_buf of the start of
-		 * the current run.
-		 */
-		yy_bp = yy_cp;
-
-		yy_current_state = (yy_start);
-		yy_current_state += YY_AT_BOL();
-yy_match:
-		do
-			{
-			register YY_CHAR yy_c = yy_ec[YY_SC_TO_UI(*yy_cp)];
-			if ( yy_accept[yy_current_state] )
-				{
-				(yy_last_accepting_state) = yy_current_state;
-				(yy_last_accepting_cpos) = yy_cp;
-				}
-			while ( yy_chk[yy_base[yy_current_state] + yy_c] != yy_current_state )
-				{
-				yy_current_state = (int) yy_def[yy_current_state];
-				if ( yy_current_state >= 152 )
-					yy_c = yy_meta[(unsigned int) yy_c];
-				}
-			yy_current_state = yy_nxt[yy_base[yy_current_state] + (unsigned int) yy_c];
-			++yy_cp;
-			}
-		while ( yy_base[yy_current_state] != 467 );
-
-yy_find_action:
-		yy_act = yy_accept[yy_current_state];
-		if ( yy_act == 0 )
-			{ /* have to back up */
-			yy_cp = (yy_last_accepting_cpos);
-			yy_current_state = (yy_last_accepting_state);
-			yy_act = yy_accept[yy_current_state];
-			}
-
-		YY_DO_BEFORE_ACTION;
-
-do_action:	/* This label is used only to access EOF actions. */
-
-		switch ( yy_act )
-	{ /* beginning of action switch */
-			case 0: /* must back up */
-			/* undo the effects of YY_DO_BEFORE_ACTION */
-			*yy_cp = (yy_hold_char);
-			yy_cp = (yy_last_accepting_cpos);
-			yy_current_state = (yy_last_accepting_state);
-			goto yy_find_action;
-
-case 1:
-/* rule 1 can match eol */
-YY_RULE_SETUP
-#line 62 "wmkdep.l"
-;            /* remove c++ style one line comments */
-	YY_BREAK
-case 2:
-YY_RULE_SETUP
-#line 64 "wmkdep.l"
-BEGIN(CMNT);    /* start removing c style comment */
-	YY_BREAK
-case 3:
-/* rule 3 can match eol */
-YY_RULE_SETUP
-#line 65 "wmkdep.l"
-;
-	YY_BREAK
-case 4:
-YY_RULE_SETUP
-#line 66 "wmkdep.l"
-BEGIN(INITIAL); /* end removing c style comment */
-	YY_BREAK
-case 5:
-YY_RULE_SETUP
-#line 68 "wmkdep.l"
-BEGIN(CFNAME);  /* c-file name */
-	YY_BREAK
-case 6:
-YY_RULE_SETUP
-#line 69 "wmkdep.l"
-{ BEGIN(INITIAL); nextFile(yytext); } /*"*/
-	YY_BREAK
-case 7:
-YY_RULE_SETUP
-#line 71 "wmkdep.l"
-BEGIN(JFNAME);  /* java-file name */
-	YY_BREAK
-case 8:
-YY_RULE_SETUP
-#line 72 "wmkdep.l"
-BEGIN(INITIAL);
-	YY_BREAK
-case 9:
-YY_RULE_SETUP
-#line 73 "wmkdep.l"
-BEGIN(INITIAL);
-	YY_BREAK
-case 10:
-YY_RULE_SETUP
-#line 74 "wmkdep.l"
-BEGIN(INITIAL);
-	YY_BREAK
-case 11:
-YY_RULE_SETUP
-#line 75 "wmkdep.l"
-BEGIN(INITIAL);
-	YY_BREAK
-case 12:
-YY_RULE_SETUP
-#line 76 "wmkdep.l"
-BEGIN(INITIAL);
-	YY_BREAK
-case 13:
-YY_RULE_SETUP
-#line 77 "wmkdep.l"
-BEGIN(INITIAL);
-	YY_BREAK
-case 14:
-YY_RULE_SETUP
-#line 78 "wmkdep.l"
-{ BEGIN(INITIAL); importFile(yytext); } /*"*/
-	YY_BREAK
-case 15:
-YY_RULE_SETUP
-#line 79 "wmkdep.l"
-{ BEGIN(INITIAL); importDirectory(yytext); } /*"*/
-	YY_BREAK
-case 16:
-YY_RULE_SETUP
-#line 81 "wmkdep.l"
-BEGIN(FFNAME);  /* FORTRAN-file name */
-	YY_BREAK
-case 17:
-/* rule 17 can match eol */
-YY_RULE_SETUP
-#line 82 "wmkdep.l"
-{ BEGIN(INITIAL); nextFile(yytext); } /*'*/
-	YY_BREAK
-case 18:
-/* rule 18 can match eol */
-YY_RULE_SETUP
-#line 84 "wmkdep.l"
-;
-	YY_BREAK
-case 19:
-YY_RULE_SETUP
-#line 86 "wmkdep.l"
-ECHO;
-	YY_BREAK
-#line 1019 "lex.yy.c"
-case YY_STATE_EOF(INITIAL):
-case YY_STATE_EOF(CMNT):
-case YY_STATE_EOF(CFNAME):
-case YY_STATE_EOF(SCFNAME):
-case YY_STATE_EOF(JFNAME):
-case YY_STATE_EOF(FFNAME):
-	yyterminate();
-
-	case YY_END_OF_BUFFER:
-		{
-		/* Amount of text matched not including the EOB char. */
-		int yy_amount_of_matched_text = (int) (yy_cp - (yytext_ptr)) - 1;
-
-		/* Undo the effects of YY_DO_BEFORE_ACTION. */
-		*yy_cp = (yy_hold_char);
-		YY_RESTORE_YY_MORE_OFFSET
-
-		if ( YY_CURRENT_BUFFER_LVALUE->yy_buffer_status == YY_BUFFER_NEW )
-			{
-			/* We're scanning a new file or input source.  It's
-			 * possible that this happened because the user
-			 * just pointed yyin at a new source and called
-			 * yylex().  If so, then we have to assure
-			 * consistency between YY_CURRENT_BUFFER and our
-			 * globals.  Here is the right place to do so, because
-			 * this is the first action (other than possibly a
-			 * back-up) that will match for the new input source.
-			 */
-			(yy_n_chars) = YY_CURRENT_BUFFER_LVALUE->yy_n_chars;
-			YY_CURRENT_BUFFER_LVALUE->yy_input_file = yyin;
-			YY_CURRENT_BUFFER_LVALUE->yy_buffer_status = YY_BUFFER_NORMAL;
-			}
-
-		/* Note that here we test for yy_c_buf_p "<=" to the position
-		 * of the first EOB in the buffer, since yy_c_buf_p will
-		 * already have been incremented past the NUL character
-		 * (since all states make transitions on EOB to the
-		 * end-of-buffer state).  Contrast this with the test
-		 * in input().
-		 */
-		if ( (yy_c_buf_p) <= &YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[(yy_n_chars)] )
-			{ /* This was really a NUL. */
-			yy_state_type yy_next_state;
-
-			(yy_c_buf_p) = (yytext_ptr) + yy_amount_of_matched_text;
-
-			yy_current_state = yy_get_previous_state(  );
-
-			/* Okay, we're now positioned to make the NUL
-			 * transition.  We couldn't have
-			 * yy_get_previous_state() go ahead and do it
-			 * for us because it doesn't know how to deal
-			 * with the possibility of jamming (and we don't
-			 * want to build jamming into it because then it
-			 * will run more slowly).
-			 */
-
-			yy_next_state = yy_try_NUL_trans( yy_current_state );
-
-			yy_bp = (yytext_ptr) + YY_MORE_ADJ;
-
-			if ( yy_next_state )
-				{
-				/* Consume the NUL. */
-				yy_cp = ++(yy_c_buf_p);
-				yy_current_state = yy_next_state;
-				goto yy_match;
-				}
-
-			else
-				{
-				yy_cp = (yy_c_buf_p);
-				goto yy_find_action;
-				}
-			}
-
-		else switch ( yy_get_next_buffer(  ) )
-			{
-			case EOB_ACT_END_OF_FILE:
-				{
-				(yy_did_buffer_switch_on_eof) = 0;
-
-				if ( yywrap( ) )
-					{
-					/* Note: because we've taken care in
-					 * yy_get_next_buffer() to have set up
-					 * yytext, we can now set up
-					 * yy_c_buf_p so that if some total
-					 * hoser (like flex itself) wants to
-					 * call the scanner after we return the
-					 * YY_NULL, it'll still work - another
-					 * YY_NULL will get returned.
-					 */
-					(yy_c_buf_p) = (yytext_ptr) + YY_MORE_ADJ;
-
-					yy_act = YY_STATE_EOF(YY_START);
-					goto do_action;
-					}
-
-				else
-					{
-					if ( ! (yy_did_buffer_switch_on_eof) )
-						YY_NEW_FILE;
-					}
-				break;
-				}
-
-			case EOB_ACT_CONTINUE_SCAN:
-				(yy_c_buf_p) =
-					(yytext_ptr) + yy_amount_of_matched_text;
-
-				yy_current_state = yy_get_previous_state(  );
-
-				yy_cp = (yy_c_buf_p);
-				yy_bp = (yytext_ptr) + YY_MORE_ADJ;
-				goto yy_match;
-
-			case EOB_ACT_LAST_MATCH:
-				(yy_c_buf_p) =
-				&YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[(yy_n_chars)];
-
-				yy_current_state = yy_get_previous_state(  );
-
-				yy_cp = (yy_c_buf_p);
-				yy_bp = (yytext_ptr) + YY_MORE_ADJ;
-				goto yy_find_action;
-			}
-		break;
-		}
-
-	default:
-		YY_FATAL_ERROR(
-			"fatal flex scanner internal error--no action found" );
-	} /* end of action switch */
-		} /* end of scanning one token */
-} /* end of yylex */
-
-/* yy_get_next_buffer - try to read in a new buffer
- *
- * Returns a code representing an action:
- *	EOB_ACT_LAST_MATCH -
- *	EOB_ACT_CONTINUE_SCAN - continue scanning from current position
- *	EOB_ACT_END_OF_FILE - end of file
- */
-static int yy_get_next_buffer (void)
-{
-    	register char *dest = YY_CURRENT_BUFFER_LVALUE->yy_ch_buf;
-	register char *source = (yytext_ptr);
-	register int number_to_move, i;
-	int ret_val;
-
-	if ( (yy_c_buf_p) > &YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[(yy_n_chars) + 1] )
-		YY_FATAL_ERROR(
-		"fatal flex scanner internal error--end of buffer missed" );
-
-	if ( YY_CURRENT_BUFFER_LVALUE->yy_fill_buffer == 0 )
-		{ /* Don't try to fill the buffer, so this is an EOF. */
-		if ( (yy_c_buf_p) - (yytext_ptr) - YY_MORE_ADJ == 1 )
-			{
-			/* We matched a single character, the EOB, so
-			 * treat this as a final EOF.
-			 */
-			return EOB_ACT_END_OF_FILE;
-			}
-
-		else
-			{
-			/* We matched some text prior to the EOB, first
-			 * process it.
-			 */
-			return EOB_ACT_LAST_MATCH;
-			}
-		}
-
-	/* Try to read more data. */
-
-	/* First move last chars to start of buffer. */
-	number_to_move = (int) ((yy_c_buf_p) - (yytext_ptr)) - 1;
-
-	for ( i = 0; i < number_to_move; ++i )
-		*(dest++) = *(source++);
-
-	if ( YY_CURRENT_BUFFER_LVALUE->yy_buffer_status == YY_BUFFER_EOF_PENDING )
-		/* don't do the read, it's not guaranteed to return an EOF,
-		 * just force an EOF
-		 */
-		YY_CURRENT_BUFFER_LVALUE->yy_n_chars = (yy_n_chars) = 0;
-
-	else
-		{
-			int num_to_read =
-			YY_CURRENT_BUFFER_LVALUE->yy_buf_size - number_to_move - 1;
-
-		while ( num_to_read <= 0 )
-			{ /* Not enough room in the buffer - grow it. */
-
-			/* just a shorter name for the current buffer */
-			YY_BUFFER_STATE b = YY_CURRENT_BUFFER;
-
-			int yy_c_buf_p_offset =
-				(int) ((yy_c_buf_p) - b->yy_ch_buf);
-
-			if ( b->yy_is_our_buffer )
-				{
-				int new_size = b->yy_buf_size * 2;
-
-				if ( new_size <= 0 )
-					b->yy_buf_size += b->yy_buf_size / 8;
-				else
-					b->yy_buf_size *= 2;
-
-				b->yy_ch_buf = (char *)
-					/* Include room in for 2 EOB chars. */
-					yyrealloc((void *) b->yy_ch_buf,b->yy_buf_size + 2  );
-				}
-			else
-				/* Can't grow it, we don't own it. */
-				b->yy_ch_buf = 0;
-
-			if ( ! b->yy_ch_buf )
-				YY_FATAL_ERROR(
-				"fatal error - scanner input buffer overflow" );
-
-			(yy_c_buf_p) = &b->yy_ch_buf[yy_c_buf_p_offset];
-
-			num_to_read = YY_CURRENT_BUFFER_LVALUE->yy_buf_size -
-						number_to_move - 1;
-
-			}
-
-		if ( num_to_read > YY_READ_BUF_SIZE )
-			num_to_read = YY_READ_BUF_SIZE;
-
-		/* Read in more data. */
-		YY_INPUT( (&YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[number_to_move]),
-			(yy_n_chars), num_to_read );
-
-		YY_CURRENT_BUFFER_LVALUE->yy_n_chars = (yy_n_chars);
-		}
-
-	if ( (yy_n_chars) == 0 )
-		{
-		if ( number_to_move == YY_MORE_ADJ )
-			{
-			ret_val = EOB_ACT_END_OF_FILE;
-			yyrestart(yyin  );
-			}
-
-		else
-			{
-			ret_val = EOB_ACT_LAST_MATCH;
-			YY_CURRENT_BUFFER_LVALUE->yy_buffer_status =
-				YY_BUFFER_EOF_PENDING;
-			}
-		}
-
-	else
-		ret_val = EOB_ACT_CONTINUE_SCAN;
-
-	(yy_n_chars) += number_to_move;
-	YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[(yy_n_chars)] = YY_END_OF_BUFFER_CHAR;
-	YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[(yy_n_chars) + 1] = YY_END_OF_BUFFER_CHAR;
-
-	(yytext_ptr) = &YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[0];
-
-	return ret_val;
-}
-
-/* yy_get_previous_state - get the state just before the EOB char was reached */
-
-    static yy_state_type yy_get_previous_state (void)
-{
-	register yy_state_type yy_current_state;
-	register char *yy_cp;
-    
-	yy_current_state = (yy_start);
-	yy_current_state += YY_AT_BOL();
-
-	for ( yy_cp = (yytext_ptr) + YY_MORE_ADJ; yy_cp < (yy_c_buf_p); ++yy_cp )
-		{
-		register YY_CHAR yy_c = (*yy_cp ? yy_ec[YY_SC_TO_UI(*yy_cp)] : 1);
-		if ( yy_accept[yy_current_state] )
-			{
-			(yy_last_accepting_state) = yy_current_state;
-			(yy_last_accepting_cpos) = yy_cp;
-			}
-		while ( yy_chk[yy_base[yy_current_state] + yy_c] != yy_current_state )
-			{
-			yy_current_state = (int) yy_def[yy_current_state];
-			if ( yy_current_state >= 152 )
-				yy_c = yy_meta[(unsigned int) yy_c];
-			}
-		yy_current_state = yy_nxt[yy_base[yy_current_state] + (unsigned int) yy_c];
-		}
-
-	return yy_current_state;
-}
-
-/* yy_try_NUL_trans - try to make a transition on the NUL character
- *
- * synopsis
- *	next_state = yy_try_NUL_trans( current_state );
- */
-    static yy_state_type yy_try_NUL_trans  (yy_state_type yy_current_state )
-{
-	register int yy_is_jam;
-    	register char *yy_cp = (yy_c_buf_p);
-
-	register YY_CHAR yy_c = 1;
-	if ( yy_accept[yy_current_state] )
-		{
-		(yy_last_accepting_state) = yy_current_state;
-		(yy_last_accepting_cpos) = yy_cp;
-		}
-	while ( yy_chk[yy_base[yy_current_state] + yy_c] != yy_current_state )
-		{
-		yy_current_state = (int) yy_def[yy_current_state];
-		if ( yy_current_state >= 152 )
-			yy_c = yy_meta[(unsigned int) yy_c];
-		}
-	yy_current_state = yy_nxt[yy_base[yy_current_state] + (unsigned int) yy_c];
-	yy_is_jam = (yy_current_state == 151);
-
-	return yy_is_jam ? 0 : yy_current_state;
-}
-
-    static void yyunput (int c, register char * yy_bp )
-{
-	register char *yy_cp;
-    
-    yy_cp = (yy_c_buf_p);
-
-	/* undo effects of setting up yytext */
-	*yy_cp = (yy_hold_char);
-
-	if ( yy_cp < YY_CURRENT_BUFFER_LVALUE->yy_ch_buf + 2 )
-		{ /* need to shift things up to make room */
-		/* +2 for EOB chars. */
-		register int number_to_move = (yy_n_chars) + 2;
-		register char *dest = &YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[
-					YY_CURRENT_BUFFER_LVALUE->yy_buf_size + 2];
-		register char *source =
-				&YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[number_to_move];
-
-		while ( source > YY_CURRENT_BUFFER_LVALUE->yy_ch_buf )
-			*--dest = *--source;
-
-		yy_cp += (int) (dest - source);
-		yy_bp += (int) (dest - source);
-		YY_CURRENT_BUFFER_LVALUE->yy_n_chars =
-			(yy_n_chars) = YY_CURRENT_BUFFER_LVALUE->yy_buf_size;
-
-		if ( yy_cp < YY_CURRENT_BUFFER_LVALUE->yy_ch_buf + 2 )
-			YY_FATAL_ERROR( "flex scanner push-back overflow" );
-		}
-
-	*--yy_cp = (char) c;
-
-	(yytext_ptr) = yy_bp;
-	(yy_hold_char) = *yy_cp;
-	(yy_c_buf_p) = yy_cp;
-}
-
-#ifndef YY_NO_INPUT
-#ifdef __cplusplus
-    static int yyinput (void)
-#else
-    static int input  (void)
-#endif
-
-{
-	int c;
-    
-	*(yy_c_buf_p) = (yy_hold_char);
-
-	if ( *(yy_c_buf_p) == YY_END_OF_BUFFER_CHAR )
-		{
-		/* yy_c_buf_p now points to the character we want to return.
-		 * If this occurs *before* the EOB characters, then it's a
-		 * valid NUL; if not, then we've hit the end of the buffer.
-		 */
-		if ( (yy_c_buf_p) < &YY_CURRENT_BUFFER_LVALUE->yy_ch_buf[(yy_n_chars)] )
-			/* This was really a NUL. */
-			*(yy_c_buf_p) = '\0';
-
-		else
-			{ /* need more input */
-			int offset = (yy_c_buf_p) - (yytext_ptr);
-			++(yy_c_buf_p);
-
-			switch ( yy_get_next_buffer(  ) )
-				{
-				case EOB_ACT_LAST_MATCH:
-					/* This happens because yy_g_n_b()
-					 * sees that we've accumulated a
-					 * token and flags that we need to
-					 * try matching the token before
-					 * proceeding.  But for input(),
-					 * there's no matching to consider.
-					 * So convert the EOB_ACT_LAST_MATCH
-					 * to EOB_ACT_END_OF_FILE.
-					 */
-
-					/* Reset buffer status. */
-					yyrestart(yyin );
-
-					/*FALLTHROUGH*/
-
-				case EOB_ACT_END_OF_FILE:
-					{
-					if ( yywrap( ) )
-						return EOF;
-
-					if ( ! (yy_did_buffer_switch_on_eof) )
-						YY_NEW_FILE;
-#ifdef __cplusplus
-					return yyinput();
-#else
-					return input();
-#endif
-					}
-
-				case EOB_ACT_CONTINUE_SCAN:
-					(yy_c_buf_p) = (yytext_ptr) + offset;
-					break;
-				}
-			}
-		}
-
-	c = *(unsigned char *) (yy_c_buf_p);	/* cast for 8-bit char's */
-	*(yy_c_buf_p) = '\0';	/* preserve yytext */
-	(yy_hold_char) = *++(yy_c_buf_p);
-
-	YY_CURRENT_BUFFER_LVALUE->yy_at_bol = (c == '\n');
-
-	return c;
-}
-#endif	/* ifndef YY_NO_INPUT */
-
-/** Immediately switch to a different input stream.
- * @param input_file A readable stream.
- * 
- * @note This function does not reset the start condition to @c INITIAL .
- */
-    void yyrestart  (FILE * input_file )
-{
-    
-	if ( ! YY_CURRENT_BUFFER ){
-        yyensure_buffer_stack ();
-		YY_CURRENT_BUFFER_LVALUE =
-            yy_create_buffer(yyin,YY_BUF_SIZE );
-	}
-
-	yy_init_buffer(YY_CURRENT_BUFFER,input_file );
-	yy_load_buffer_state( );
-}
-
-/** Switch to a different input buffer.
- * @param new_buffer The new input buffer.
- * 
- */
-    void yy_switch_to_buffer  (YY_BUFFER_STATE  new_buffer )
-{
-    
-	/* TODO. We should be able to replace this entire function body
-	 * with
-	 *		yypop_buffer_state();
-	 *		yypush_buffer_state(new_buffer);
-     */
-	yyensure_buffer_stack ();
-	if ( YY_CURRENT_BUFFER == new_buffer )
-		return;
-
-	if ( YY_CURRENT_BUFFER )
-		{
-		/* Flush out information for old buffer. */
-		*(yy_c_buf_p) = (yy_hold_char);
-		YY_CURRENT_BUFFER_LVALUE->yy_buf_pos = (yy_c_buf_p);
-		YY_CURRENT_BUFFER_LVALUE->yy_n_chars = (yy_n_chars);
-		}
-
-	YY_CURRENT_BUFFER_LVALUE = new_buffer;
-	yy_load_buffer_state( );
-
-	/* We don't actually know whether we did this switch during
-	 * EOF (yywrap()) processing, but the only time this flag
-	 * is looked at is after yywrap() is called, so it's safe
-	 * to go ahead and always set it.
-	 */
-	(yy_did_buffer_switch_on_eof) = 1;
-}
-
-static void yy_load_buffer_state  (void)
-{
-    	(yy_n_chars) = YY_CURRENT_BUFFER_LVALUE->yy_n_chars;
-	(yytext_ptr) = (yy_c_buf_p) = YY_CURRENT_BUFFER_LVALUE->yy_buf_pos;
-	yyin = YY_CURRENT_BUFFER_LVALUE->yy_input_file;
-	(yy_hold_char) = *(yy_c_buf_p);
-}
-
-/** Allocate and initialize an input buffer state.
- * @param file A readable stream.
- * @param size The character buffer size in bytes. When in doubt, use @c YY_BUF_SIZE.
- * 
- * @return the allocated buffer state.
- */
-    YY_BUFFER_STATE yy_create_buffer  (FILE * file, int  size )
-{
-	YY_BUFFER_STATE b;
-    
-	b = (YY_BUFFER_STATE) yyalloc(sizeof( struct yy_buffer_state )  );
-	if ( ! b )
-		YY_FATAL_ERROR( "out of dynamic memory in yy_create_buffer()" );
-
-	b->yy_buf_size = size;
-
-	/* yy_ch_buf has to be 2 characters longer than the size given because
-	 * we need to put in 2 end-of-buffer characters.
-	 */
-	b->yy_ch_buf = (char *) yyalloc(b->yy_buf_size + 2  );
-	if ( ! b->yy_ch_buf )
-		YY_FATAL_ERROR( "out of dynamic memory in yy_create_buffer()" );
-
-	b->yy_is_our_buffer = 1;
-
-	yy_init_buffer(b,file );
-
-	return b;
-}
-
-/** Destroy the buffer.
- * @param b a buffer created with yy_create_buffer()
- * 
- */
-    void yy_delete_buffer (YY_BUFFER_STATE  b )
-{
-    
-	if ( ! b )
-		return;
-
-	if ( b == YY_CURRENT_BUFFER ) /* Not sure if we should pop here. */
-		YY_CURRENT_BUFFER_LVALUE = (YY_BUFFER_STATE) 0;
-
-	if ( b->yy_is_our_buffer )
-		yyfree((void *) b->yy_ch_buf  );
-
-	yyfree((void *) b  );
-}
-
-#ifndef __cplusplus
-extern int isatty (int );
-#endif /* __cplusplus */
-    
-/* Initializes or reinitializes a buffer.
- * This function is sometimes called more than once on the same buffer,
- * such as during a yyrestart() or at EOF.
- */
-    static void yy_init_buffer  (YY_BUFFER_STATE  b, FILE * file )
-
-{
-	int oerrno = errno;
-    
-	yy_flush_buffer(b );
-
-	b->yy_input_file = file;
-	b->yy_fill_buffer = 1;
-
-    /* If b is the current buffer, then yy_init_buffer was _probably_
-     * called from yyrestart() or through yy_get_next_buffer.
-     * In that case, we don't want to reset the lineno or column.
-     */
-    if (b != YY_CURRENT_BUFFER){
-        b->yy_bs_lineno = 1;
-        b->yy_bs_column = 0;
-    }
-
-        b->yy_is_interactive = file ? (isatty( fileno(file) ) > 0) : 0;
-    
-	errno = oerrno;
-}
-
-/** Discard all buffered characters. On the next scan, YY_INPUT will be called.
- * @param b the buffer state to be flushed, usually @c YY_CURRENT_BUFFER.
- * 
- */
-    void yy_flush_buffer (YY_BUFFER_STATE  b )
-{
-    	if ( ! b )
-		return;
-
-	b->yy_n_chars = 0;
-
-	/* We always need two end-of-buffer characters.  The first causes
-	 * a transition to the end-of-buffer state.  The second causes
-	 * a jam in that state.
-	 */
-	b->yy_ch_buf[0] = YY_END_OF_BUFFER_CHAR;
-	b->yy_ch_buf[1] = YY_END_OF_BUFFER_CHAR;
-
-	b->yy_buf_pos = &b->yy_ch_buf[0];
-
-	b->yy_at_bol = 1;
-	b->yy_buffer_status = YY_BUFFER_NEW;
-
-	if ( b == YY_CURRENT_BUFFER )
-		yy_load_buffer_state( );
-}
-
-/** Pushes the new state onto the stack. The new state becomes
- *  the current state. This function will allocate the stack
- *  if necessary.
- *  @param new_buffer The new state.
- *  
- */
-void yypush_buffer_state (YY_BUFFER_STATE new_buffer )
-{
-    	if (new_buffer == NULL)
-		return;
-
-	yyensure_buffer_stack();
-
-	/* This block is copied from yy_switch_to_buffer. */
-	if ( YY_CURRENT_BUFFER )
-		{
-		/* Flush out information for old buffer. */
-		*(yy_c_buf_p) = (yy_hold_char);
-		YY_CURRENT_BUFFER_LVALUE->yy_buf_pos = (yy_c_buf_p);
-		YY_CURRENT_BUFFER_LVALUE->yy_n_chars = (yy_n_chars);
-		}
-
-	/* Only push if top exists. Otherwise, replace top. */
-	if (YY_CURRENT_BUFFER)
-		(yy_buffer_stack_top)++;
-	YY_CURRENT_BUFFER_LVALUE = new_buffer;
-
-	/* copied from yy_switch_to_buffer. */
-	yy_load_buffer_state( );
-	(yy_did_buffer_switch_on_eof) = 1;
-}
-
-/** Removes and deletes the top of the stack, if present.
- *  The next element becomes the new top.
- *  
- */
-void yypop_buffer_state (void)
-{
-    	if (!YY_CURRENT_BUFFER)
-		return;
-
-	yy_delete_buffer(YY_CURRENT_BUFFER );
-	YY_CURRENT_BUFFER_LVALUE = NULL;
-	if ((yy_buffer_stack_top) > 0)
-		--(yy_buffer_stack_top);
-
-	if (YY_CURRENT_BUFFER) {
-		yy_load_buffer_state( );
-		(yy_did_buffer_switch_on_eof) = 1;
-	}
-}
-
-/* Allocates the stack if it does not exist.
- *  Guarantees space for at least one push.
- */
-static void yyensure_buffer_stack (void)
-{
-	int num_to_alloc;
-    
-	if (!(yy_buffer_stack)) {
-
-		/* First allocation is just for 2 elements, since we don't know if this
-		 * scanner will even need a stack. We use 2 instead of 1 to avoid an
-		 * immediate realloc on the next call.
-         */
-		num_to_alloc = 1;
-		(yy_buffer_stack) = (struct yy_buffer_state**)yyalloc
-								(num_to_alloc * sizeof(struct yy_buffer_state*)
-								);
-		
-		memset((yy_buffer_stack), 0, num_to_alloc * sizeof(struct yy_buffer_state*));
-				
-		(yy_buffer_stack_max) = num_to_alloc;
-		(yy_buffer_stack_top) = 0;
-		return;
-	}
-
-	if ((yy_buffer_stack_top) >= ((yy_buffer_stack_max)) - 1){
-
-		/* Increase the buffer to prepare for a possible push. */
-		int grow_size = 8 /* arbitrary grow size */;
-
-		num_to_alloc = (yy_buffer_stack_max) + grow_size;
-		(yy_buffer_stack) = (struct yy_buffer_state**)yyrealloc
-								((yy_buffer_stack),
-								num_to_alloc * sizeof(struct yy_buffer_state*)
-								);
-
-		/* zero only the new slots.*/
-		memset((yy_buffer_stack) + (yy_buffer_stack_max), 0, grow_size * sizeof(struct yy_buffer_state*));
-		(yy_buffer_stack_max) = num_to_alloc;
-	}
-}
-
-/** Setup the input buffer state to scan directly from a user-specified character buffer.
- * @param base the character buffer
- * @param size the size in bytes of the character buffer
- * 
- * @return the newly allocated buffer state object. 
- */
-YY_BUFFER_STATE yy_scan_buffer  (char * base, yy_size_t  size )
-{
-	YY_BUFFER_STATE b;
-    
-	if ( size < 2 ||
-	     base[size-2] != YY_END_OF_BUFFER_CHAR ||
-	     base[size-1] != YY_END_OF_BUFFER_CHAR )
-		/* They forgot to leave room for the EOB's. */
-		return 0;
-
-	b = (YY_BUFFER_STATE) yyalloc(sizeof( struct yy_buffer_state )  );
-	if ( ! b )
-		YY_FATAL_ERROR( "out of dynamic memory in yy_scan_buffer()" );
-
-	b->yy_buf_size = size - 2;	/* "- 2" to take care of EOB's */
-	b->yy_buf_pos = b->yy_ch_buf = base;
-	b->yy_is_our_buffer = 0;
-	b->yy_input_file = 0;
-	b->yy_n_chars = b->yy_buf_size;
-	b->yy_is_interactive = 0;
-	b->yy_at_bol = 1;
-	b->yy_fill_buffer = 0;
-	b->yy_buffer_status = YY_BUFFER_NEW;
-
-	yy_switch_to_buffer(b  );
-
-	return b;
-}
-
-/** Setup the input buffer state to scan a string. The next call to yylex() will
- * scan from a @e copy of @a str.
- * @param str a NUL-terminated string to scan
- * 
- * @return the newly allocated buffer state object.
- * @note If you want to scan bytes that may contain NUL values, then use
- *       yy_scan_bytes() instead.
- */
-YY_BUFFER_STATE yy_scan_string (yyconst char * yystr )
-{
-    
-	return yy_scan_bytes(yystr,strlen(yystr) );
-}
-
-/** Setup the input buffer state to scan the given bytes. The next call to yylex() will
- * scan from a @e copy of @a bytes.
- * @param bytes the byte buffer to scan
- * @param len the number of bytes in the buffer pointed to by @a bytes.
- * 
- * @return the newly allocated buffer state object.
- */
-YY_BUFFER_STATE yy_scan_bytes  (yyconst char * yybytes, int  _yybytes_len )
-{
-	YY_BUFFER_STATE b;
-	char *buf;
-	yy_size_t n;
-	int i;
-    
-	/* Get memory for full buffer, including space for trailing EOB's. */
-	n = _yybytes_len + 2;
-	buf = (char *) yyalloc(n  );
-	if ( ! buf )
-		YY_FATAL_ERROR( "out of dynamic memory in yy_scan_bytes()" );
-
-	for ( i = 0; i < _yybytes_len; ++i )
-		buf[i] = yybytes[i];
-
-	buf[_yybytes_len] = buf[_yybytes_len+1] = YY_END_OF_BUFFER_CHAR;
-
-	b = yy_scan_buffer(buf,n );
-	if ( ! b )
-		YY_FATAL_ERROR( "bad buffer in yy_scan_bytes()" );
-
-	/* It's okay to grow etc. this buffer, and we should throw it
-	 * away when we're done.
-	 */
-	b->yy_is_our_buffer = 1;
-
-	return b;
-}
-
-#ifndef YY_EXIT_FAILURE
-#define YY_EXIT_FAILURE 2
-#endif
-
-static void yy_fatal_error (yyconst char* msg )
-{
-    	(void) fprintf( stderr, "%s\n", msg );
-	exit( YY_EXIT_FAILURE );
-}
-
-/* Redefine yyless() so it works in section 3 code. */
-
-#undef yyless
-#define yyless(n) \
-	do \
-		{ \
-		/* Undo effects of setting up yytext. */ \
-        int yyless_macro_arg = (n); \
-        YY_LESS_LINENO(yyless_macro_arg);\
-		yytext[yyleng] = (yy_hold_char); \
-		(yy_c_buf_p) = yytext + yyless_macro_arg; \
-		(yy_hold_char) = *(yy_c_buf_p); \
-		*(yy_c_buf_p) = '\0'; \
-		yyleng = yyless_macro_arg; \
-		} \
-	while ( 0 )
-
-/* Accessor  methods (get/set functions) to struct members. */
-
-/** Get the current line number.
- * 
- */
-int yyget_lineno  (void)
-{
-        
-    return yylineno;
-}
-
-/** Get the input stream.
- * 
- */
-FILE *yyget_in  (void)
-{
-        return yyin;
-}
-
-/** Get the output stream.
- * 
- */
-FILE *yyget_out  (void)
-{
-        return yyout;
-}
-
-/** Get the length of the current token.
- * 
- */
-int yyget_leng  (void)
-{
-        return yyleng;
-}
-
-/** Get the current token.
- * 
- */
-
-char *yyget_text  (void)
-{
-        return yytext;
-}
-
-/** Set the current line number.
- * @param line_number
- * 
- */
-void yyset_lineno (int  line_number )
-{
-    
-    yylineno = line_number;
-}
-
-/** Set the input stream. This does not discard the current
- * input buffer.
- * @param in_str A readable stream.
- * 
- * @see yy_switch_to_buffer
- */
-void yyset_in (FILE *  in_str )
-{
-        yyin = in_str ;
-}
-
-void yyset_out (FILE *  out_str )
-{
-        yyout = out_str ;
-}
-
-int yyget_debug  (void)
-{
-        return yy_flex_debug;
-}
-
-void yyset_debug (int  bdebug )
-{
-        yy_flex_debug = bdebug ;
-}
-
-static int yy_init_globals (void)
-{
-        /* Initialization is the same as for the non-reentrant scanner.
-     * This function is called from yylex_destroy(), so don't allocate here.
-     */
-
-    (yy_buffer_stack) = 0;
-    (yy_buffer_stack_top) = 0;
-    (yy_buffer_stack_max) = 0;
-    (yy_c_buf_p) = (char *) 0;
-    (yy_init) = 0;
-    (yy_start) = 0;
-
-/* Defined in main.c */
-#ifdef YY_STDINIT
-    yyin = stdin;
-    yyout = stdout;
-#else
-    yyin = (FILE *) 0;
-    yyout = (FILE *) 0;
-#endif
-
-    /* For future reference: Set errno on error, since we are called by
-     * yylex_init()
-     */
-    return 0;
-}
-
-/* yylex_destroy is for both reentrant and non-reentrant scanners. */
-int yylex_destroy  (void)
-{
-    
-    /* Pop the buffer stack, destroying each element. */
-	while(YY_CURRENT_BUFFER){
-		yy_delete_buffer(YY_CURRENT_BUFFER  );
-		YY_CURRENT_BUFFER_LVALUE = NULL;
-		yypop_buffer_state();
-	}
-
-	/* Destroy the stack itself. */
-	yyfree((yy_buffer_stack) );
-	(yy_buffer_stack) = NULL;
-
-    /* Reset the globals. This is important in a non-reentrant scanner so the next time
-     * yylex() is called, initialization will occur. */
-    yy_init_globals( );
-
-    return 0;
-}
-
-/*
- * Internal utility routines.
- */
-
-#ifndef yytext_ptr
-static void yy_flex_strncpy (char* s1, yyconst char * s2, int n )
-{
-	register int i;
-	for ( i = 0; i < n; ++i )
-		s1[i] = s2[i];
-}
-#endif
-
-#ifdef YY_NEED_STRLEN
-static int yy_flex_strlen (yyconst char * s )
-{
-	register int n;
-	for ( n = 0; s[n]; ++n )
-		;
-
-	return n;
-}
-#endif
-
-void *yyalloc (yy_size_t  size )
-{
-	return (void *) malloc( size );
-}
-
-void *yyrealloc  (void * ptr, yy_size_t  size )
-{
-	/* The cast to (char *) in the following accommodates both
-	 * implementations that use char* generic pointers, and those
-	 * that use void* generic pointers.  It works with the latter
-	 * because both ANSI C and C++ allow castless assignment from
-	 * any pointer type to void*, and deal with argument conversions
-	 * as though doing an assignment.
-	 */
-	return (void *) realloc( (char *) ptr, size );
-}
-
-void yyfree (void * ptr )
-{
-	free( (char *) ptr );	/* see yyrealloc() for (char *) cast */
-}
-
-#define YYTABLES_NAME "yytables"
-
-#line 86 "wmkdep.l"
-
-
-
-
-int nDirectories;
-char** directories;
-char* sourceFile = NULL;
-char* sourceExt = NULL;
-char* objectFile = NULL;
-char* classFile = NULL;
-char* depFile = NULL;
-
-int main(int argc, char* argv[])
-{
-    char *dotPtr, *slashPtr;
-    int i;
-
-    if (argc == 1)
-    {
-        fprintf(stderr,"input file not supplied\n");
-        exit(1);
-    }
-
-    sourceFile = (char*)malloc(strlen(argv[argc-1]) + 1);
-    strcpy(sourceFile, argv[argc-1]);
-    fprintf(stderr,"Making dependency list for source file %s\n", sourceFile);
-    
-    
-    /* Get list of -I directories. */
-    nDirectories = 0;
-
-    for (i = 1; i < argc; i++)
-    {
-        if (strncmp(argv[i], "-I", 2) == 0)
-        {
-            nDirectories++;
-        }
-    }
-
-    directories = (char**)malloc(sizeof(char*)*nDirectories);
-
-    nDirectories = 0;
-
-    for (i = 1; i < argc; i++)
-    {
-        if (strncmp(argv[i], "-I", 2) == 0)
-        {
-            directories[nDirectories++] = strdup(argv[i] + 2);
-        }
-    }
-
-
-
-    if ((dotPtr = strrchr(sourceFile, '.')) == NULL)
-    {
-        fprintf
-        (
-            stderr,
-            "Cannot find extension in source file name %s\n",
-            sourceFile
-        );
-        exit(1);
-    }
-
-    if ((slashPtr = strrchr(sourceFile, '/')) == NULL)
-    {
-        slashPtr = sourceFile;
-    }
-    else
-    {
-        slashPtr++;
-    }
-
-    sourceExt = (char*)malloc(strlen(sourceFile));
-    sourceExt[0] = 0;
-    strncat
-    (
-        sourceExt,
-        dotPtr+1,
-        (&sourceFile[strlen(sourceFile) - 1] - dotPtr)/sizeof(char)
-    );
-
-    /* initialise depFile to zero and use strncat rather than strncpy
-       because there is a bug in the SGI strncat that if 0 preceeds the .
-       it inserts a space */
-    depFile = (char*)malloc(strlen(sourceFile) + 3);
-    depFile[0] = 0;
-    strncat(depFile, sourceFile, (dotPtr - sourceFile)/sizeof(char));
-    strcat(depFile, ".dep");
-
-    if (strcmp(sourceExt, "java") == 0)
-    {
-        classFile = (char*)malloc(strlen(sourceFile) + 17);
-        strcpy(classFile, "$(CLASSES_DIR)/");
-        strncat(classFile, sourceFile, (dotPtr - sourceFile)/sizeof(char));
-        strcat(classFile, ".class");
-
-        printf("%s: %s\n", classFile, depFile);
-    }
-    else
-    {
-        objectFile = (char*)malloc(strlen(sourceFile) + 16);
-        strcpy(objectFile, "$(OBJECTS_DIR)/");
-        strncat(objectFile, slashPtr, (dotPtr - slashPtr)/sizeof(char));
-        strcat(objectFile, ".o");
-
-        printf("%s: %s\n", objectFile, depFile);
-    }
-
-    nextFile(sourceFile);
-    yylex();
-
-
-    for (i = 0; i < nDirectories; i++)
-    {
-        free(directories[i]);
-    }
-    free(directories);
-
-    free(sourceFile);
-    free(sourceExt);
-    free(objectFile);
-    free(classFile);
-    free(depFile);
-
-    return 0;
-}
-
-
-int currentBuffer = 0;                     /* Buffer pointer stack counter */
-YY_BUFFER_STATE buffers[FILE_STACK_SIZE];  /* Buffer pointer stack */
-
-/* file name entry in hash table */
-
-struct FileName
-{
-    char* name;
-    struct FileName* next;
-};
-
-struct FileName* fileHashTable[HASH_TABLE_SIZE];   /* File hash table */
-struct FileName* dirHashTable[HASH_TABLE_SIZE];    /* Directory hash table */
-
-
-/* lookup name in hash table, if not found insert in table */
-
-int lookUp(struct FileName** hashTable, const char* p)
-{
-    int ii = 0;
-    struct FileName* n;
-    struct FileName* nn;
-
-    /* hash */
-    const char* pp = p;
-    while (*pp) ii = ii<<1 ^ *pp++;
-    if (ii < 0) ii = -ii;
-    ii %= HASH_TABLE_SIZE;
-
-    /* search */
-    for (n=hashTable[ii]; n; n=n->next)
-    {
-        if (strcmp(p, n->name) == 0)
-        {
-            /* entry found so return true */
-            return 1;
-        }
-    }
-
-    /* insert */
-    nn = (struct FileName*)malloc(sizeof(struct FileName));
-    nn->name = (char*)malloc(strlen(p)+1);
-    strcpy(nn->name, p);
-    nn->next = hashTable[ii];
-    hashTable[ii] = nn;
-
-    /* entry not found, and therefore added. return false */
-    return 0;
-}
-
-
-
-/* Add a directory name to the file name */
-
-char* addDirectoryName(const char* directoryName, const char* fileName)
-{
-    char* pathName;
-    pathName = (char*)malloc(strlen(directoryName) + strlen(fileName) + 2);
-
-    strcpy(pathName, directoryName);
-
-    if (directoryName[strlen(directoryName)-1] != '/')
-    {
-        strcat(pathName, "/");
-    }
-
-    strcat(pathName, fileName);
-
-    return pathName;
-}
-
-
-/* open a file and create buffer and put on stack stack */
-
-void nextFile(const char* fileName)
-{
-    int d;
-    char* pathName;
-
-    if (lookUp(fileHashTable, fileName)) return;
-
-    if (currentBuffer >= FILE_STACK_SIZE)
-    {
-        fprintf
-        (
-            stderr,
-            "depth of file search exceeds stack size %d "
-            "while opening %s for file %s\n",
-            FILE_STACK_SIZE, fileName, sourceFile
-        );
-        exit(1);
-    }
-
-    /* Pointer to new file which is set if the file is successfully opened */
-    FILE* newyyin = NULL;
-
-    if (!(newyyin = fopen(fileName, "r")))
-    {
-        for (d=0; d<nDirectories; d++)
-        {
-            pathName = addDirectoryName(directories[d], fileName);
-
-            if ((newyyin = fopen(pathName, "r")))
-            { 
-                printf("%s: %s\n", depFile, pathName);
-                buffers[currentBuffer++] = YY_CURRENT_BUFFER;
-                yy_switch_to_buffer(yy_create_buffer(newyyin,YY_BUF_SIZE));
-
-                free(pathName);
-
-                return;
-            }
-
-            free(pathName);
-        }
-
-        fprintf
-        (
-            stderr, 
-            "could not open file %s for source file %s\n", 
-            fileName, sourceFile
-        );
-
-        fflush(stdout);
-    }
-    else
-    {
-        printf("%s: %s\n", depFile, fileName);
-        fflush(stdout);
-
-        buffers[currentBuffer++] = YY_CURRENT_BUFFER;
-        yy_switch_to_buffer(yy_create_buffer(newyyin,YY_BUF_SIZE));
-    }
-}
-
-
-void dotToSlash(char* fileName)
-{
-    int i, len;
-
-    len = strlen(fileName);
-
-    for (i=0; i<len; i++)
-    {
-        if (fileName[i] == '.') fileName[i] = '/';
-    }
-}
-
-
-void importFile(const char* fileName)
-{
-    char* javaFileName;
-    int fileNameLen;
-
-    fileNameLen = strlen(fileName);
-    javaFileName = (char*)malloc(fileNameLen + 6);
-    javaFileName[0] = 0;
-    strncat(javaFileName, fileName, fileNameLen - 1);
-    dotToSlash(javaFileName);
-    strcat(javaFileName, ".java");
-
-    nextFile(javaFileName);
-
-    free(javaFileName);
-}
-
-
-void importDirectory(const char* dirName)
-{
-    int dirNameLen;
-    char *uDirName, *path;
-    DIR *source;
-    struct dirent *list;
-
-    if (lookUp(dirHashTable, dirName)) return;
-
-    dirNameLen = strlen(dirName);
-    uDirName = strdup(dirName);
-    dotToSlash(uDirName);
-    uDirName[dirNameLen-2] = 0;
-
-    if ((source = opendir(uDirName)) == NULL)
-    {
-        fprintf
-        (
-            stderr, 
-            "could not open directory %s\n", 
-            uDirName
-        );
-        free(uDirName);
-
-        return;
-    }
-    else
-    {
-        /* Read and parse all the entries in the directory */
-        while ((list = readdir(source)) != NULL)
-        {
-            if
-            (
-                strstr(list->d_name, ".java")
-             && !strstr(list->d_name, ".java~")
-            )
-            {
-                path = addDirectoryName(uDirName, list->d_name);
-                nextFile(path);
-                free(path);
-            }
-        }
-
-        closedir(source);
-        free(uDirName);
-    }
-}
-
-
-/* The lexer calls yywrap to handle EOF conditions */
-
-int yywrap()
-{
-    /* Close the file for the buffer which has just reached EOF */
-    /* This causes strange problems
-    fclose(yyin);
-    yyin = 0;
-    */
-
-    /* Delete the buffer */
-    yy_delete_buffer(YY_CURRENT_BUFFER);
-
-    /* Set buffer counter to previous buffer */
-    currentBuffer--;
-
-    if (currentBuffer >= 0) /* if buffer counter refers to a valid file */
-    {
-        /* reset input buffer to the previous buffer on the stack */
-        yy_switch_to_buffer(buffers[currentBuffer]);
-
-        /* Return to the normal state for the previous buffer on the stack */
-        BEGIN(INITIAL);
-
-        /* return 0 to inform lex to continue reading */
-        return 0;
-    }
-    else                 /* else there are no more buffers on the stack */
-    {
-        /* return 1 to inform lex finish now that all buffers have been read */
-        return 1;
-    }
-}
-
-
-/*****************************************************************************/
-