diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduInterface/BlockLduInterface.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduInterface/BlockLduInterface.H
index bbeaf215a..e0f34d61a 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduInterface/BlockLduInterface.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduInterface/BlockLduInterface.H
@@ -86,7 +86,7 @@ public:
         // Coupled interface matrix update
 
             //- Initialise matrix update
-            virtual void initMatrixUpdate
+            virtual void initInterfaceMatrixUpdate
             (
                 const BlockLduMatrix<Type>& matrix,
                 Field<Type>& Ax,
@@ -95,7 +95,7 @@ public:
             {}
 
             //- Update result based on interface functionality
-            virtual void updateMatrix
+            virtual void updateInterfaceMatrix
             (
                 const BlockLduMatrix<Type>& matrix,
                 const FieldField<CoeffField, Type>& coeffs,
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.C
index 56e550609..bf1e1e4e6 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.C
@@ -54,8 +54,19 @@ Foam::BlockLduMatrix<Type>::BlockLduMatrix(const lduMesh& ldu)
     diagPtr_(NULL),
     upperPtr_(NULL),
     lowerPtr_(NULL),
+    interfaces_(),
+    coupleUpper_(ldu.lduAddr().nPatches()),
+    coupleLower_(ldu.lduAddr().nPatches()),
     fixedEqns_(ldu.lduAddr().size()/fixFillIn)
-{}
+{
+    const lduAddressing& addr = ldu.lduAddr();
+
+    forAll (coupleUpper_, i)
+    {
+        coupleUpper_.set(i, new CoeffField<Type>(addr.patchAddr(i).size()));
+        coupleLower_.set(i, new CoeffField<Type>(addr.patchAddr(i).size()));
+    }
+}
 
 
 template<class Type>
@@ -66,6 +77,9 @@ Foam::BlockLduMatrix<Type>::BlockLduMatrix(const BlockLduMatrix<Type>& A)
     diagPtr_(NULL),
     upperPtr_(NULL),
     lowerPtr_(NULL),
+    interfaces_(),
+    coupleUpper_(ldu.lduAddr().nPatches()),
+    coupleLower_(ldu.lduAddr().nPatches()),
     fixedEqns_(A.fixedEqns_)
 {
     if (A.diagPtr_)
@@ -82,8 +96,19 @@ Foam::BlockLduMatrix<Type>::BlockLduMatrix(const BlockLduMatrix<Type>& A)
     {
         lowerPtr_ = new TypeCoeffField(*(A.lowerPtr_));
     }
+
+    const lduAddressing& addr = ldu.lduAddr();
+
+    forAll (coupleUpper_, i)
+    {
+        coupleUpper_.set(i, A.coupleUpper_.clone());
+        coupleLower_.set(i, , A.coupleLower_.clone());
+    }
 }
 
+
+//HJ, problematic: memmory management.
+// Reconsider.  HJ, 7/Nov/2010
 template<class Type>
 Foam::BlockLduMatrix<Type>::BlockLduMatrix(BlockLduMatrix<Type>& A, bool reUse)
 :
@@ -92,6 +117,9 @@ Foam::BlockLduMatrix<Type>::BlockLduMatrix(BlockLduMatrix<Type>& A, bool reUse)
     diagPtr_(NULL),
     upperPtr_(NULL),
     lowerPtr_(NULL),
+    interfaces_(),
+    coupleUpper_(),
+    coupleLower_(),
     fixedEqns_(A.fixedEqns_)
 {
     if (reUse)
@@ -131,8 +159,13 @@ Foam::BlockLduMatrix<Type>::BlockLduMatrix(BlockLduMatrix<Type>& A, bool reUse)
             lowerPtr_ = new TypeCoeffField(*(A.lowerPtr_));
         }
     }
+
+    // Interface data
+    coupleUpper_ = A.coupleUpper_;
+    coupleLower_ = A.coupleUpper_;
 }
 
+
 // * * * * * * * * * * * * * * * * Destructor  * * * * * * * * * * * * * * * //
 
 template<class Type>
@@ -143,6 +176,7 @@ Foam::BlockLduMatrix<Type>::~BlockLduMatrix()
     deleteDemandDrivenData(lowerPtr_);
 }
 
+
 // * * * * * * * * * * * * * * * Member Functions  * * * * * * * * * * * * * //
 
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.H
index 1cf200cba..79f324d58 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrix.H
@@ -51,6 +51,7 @@ SourceFiles
 #include "coeffFields.H"
 #include "FieldField.H"
 #include "lduMesh.H"
+#include "BlockLduInterface.H"
 #include "BlockLduInterfaceFieldPtrsList.H"
 #include "Map.H"
 
@@ -97,7 +98,7 @@ private:
         const lduMesh& lduMesh_;
 
 
-        // Block matrix elements (excluding interfaces)
+        // Block matrix elements
 
             //- Diagonal coefficients
             CoeffField<Type>* diagPtr_;
@@ -109,6 +110,18 @@ private:
             CoeffField<Type> *lowerPtr_;
 
 
+        // Coupling
+
+            //- List of coupled interfaces
+            PtrList<BlockLduInterface<Type> > interfaces_;
+
+            //- Coupled interface coefficients, upper
+            FieldField<CoeffField, Type> coupleUpper_;
+
+            //- Coupled interface coefficients, lower
+            FieldField<CoeffField, Type> coupleLower_;
+
+
         // Constraints
 
             //- Equation triangle map
@@ -209,6 +222,7 @@ public:
         //- Construct as copy or re-use as specified.
         BlockLduMatrix(BlockLduMatrix<Type>&, bool reUse);
 
+
     // Destructor
 
         virtual ~BlockLduMatrix();
@@ -259,6 +273,30 @@ public:
             //- Return lower coefficients
             const TypeCoeffField& lower() const;
 
+            //- Return access to coupled interface coefficients, upper
+            FieldField<CoeffField, Type>& coupleUpper()
+            {
+                return coupleUpper_;
+            }
+
+            //- Return coupled interface coefficients, upper
+            const FieldField<CoeffField, Type>& coupleUpper() const
+            {
+                return coupleUpper_;
+            }
+
+            //- Return access to coupled interface coefficients, lower
+            FieldField<CoeffField, Type>& coupleLower()
+            {
+                return coupleLower_;
+            }
+
+            //- Return coupled interface coefficients, lower
+            const FieldField<CoeffField, Type>& coupleLower() const
+            {
+                return coupleLower_;
+            }
+
 
         // Matrix structure
 
@@ -319,10 +357,7 @@ public:
             void Amul
             (
                 TypeField& Ax,
-                const TypeField& x,
-                const FieldField<CoeffField, Type>& boundaryCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                    interfaces
+                const TypeField& x
             ) const;
 
             //- Matrix multiplication without coupled interface update
@@ -336,10 +371,7 @@ public:
             void Tmul
             (
                 TypeField& Ax,
-                const TypeField& x,
-                const FieldField<CoeffField, Type>& boundaryCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                    interfaces
+                const TypeField& x
             ) const;
 
             //- Matrix transpose multiplication without coupled
@@ -366,7 +398,6 @@ public:
             void initInterfaces
             (
                 const FieldField<CoeffField, Type>& interfaceCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
                 TypeField& result,
                 const TypeField& psi
             ) const;
@@ -375,7 +406,6 @@ public:
             void updateInterfaces
             (
                 const FieldField<CoeffField, Type>& interfaceCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
                 TypeField& result,
                 const TypeField& psi
             ) const;
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixATmul.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixATmul.C
index a3e18b5cb..2d8c79231 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixATmul.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixATmul.C
@@ -35,20 +35,18 @@ template<class Type>
 void Foam::BlockLduMatrix<Type>::Amul
 (
     TypeField& Ax,
-    const TypeField& x,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces
+    const TypeField& x
 ) const
 {
     Ax = pTraits<Type>::zero;
 
     // Initialise the update of coupled interfaces
-    initInterfaces(boundaryCoeffs, interfaces, Ax, x);
+    initInterfaces(coupleUpper_, Ax, x);
 
     AmulCore(Ax, x);
 
     // Update coupled interfaces
-    updateInterfaces(boundaryCoeffs, interfaces, Ax, x);
+    updateInterfaces(coupleUpper_, Ax, x);
 }
 
 
@@ -179,20 +177,18 @@ template<class Type>
 void Foam::BlockLduMatrix<Type>::Tmul
 (
     TypeField& Ax,
-    const TypeField& x,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces
+    const TypeField& x
 ) const
 {
     Ax = pTraits<Type>::zero;
 
     // Initialise the update of coupled interfaces
-    initInterfaces(boundaryCoeffs, interfaces, Ax, x);
+    initInterfaces(coupleLower_, Ax, x);
 
     TmulCore(Ax, x);
 
     // Update coupled interfaces
-    updateInterfaces(boundaryCoeffs, interfaces, Ax, x);
+    updateInterfaces(coupleLower_, Ax, x);
 }
 
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixUpdateInterfaces.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixUpdateInterfaces.C
index a5b494996..780e17001 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixUpdateInterfaces.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduMatrix/BlockLduMatrixUpdateInterfaces.C
@@ -35,7 +35,6 @@ template<class Type>
 void Foam::BlockLduMatrix<Type>::initInterfaces
 (
     const FieldField<CoeffField, Type>& interfaceCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     TypeField& result,
     const TypeField& psi
 ) const
@@ -46,18 +45,18 @@ void Foam::BlockLduMatrix<Type>::initInterfaces
      || Pstream::defaultCommsType == Pstream::nonBlocking
     )
     {
-        forAll (interfaces, interfaceI)
+        forAll (interfaces_, interfaceI)
         {
-            if (interfaces.set(interfaceI))
+            if (interfaces_.set(interfaceI))
             {
-                interfaces[interfaceI].initInterfaceMatrixUpdate
-                (
-                    psi,
-                    result,
-                    *this,
-                    interfaceCoeffs[interfaceI],
-                    Pstream::defaultCommsType
-                );
+//                 interfaces_[interfaceI].initInterfaceMatrixUpdate
+//                 (
+//                     psi,
+//                     result,
+//                     *this,
+//                     interfaceCoeffs[interfaceI],
+//                     Pstream::defaultCommsType
+//                 );
             }
         }
     }
@@ -69,24 +68,24 @@ void Foam::BlockLduMatrix<Type>::initInterfaces
         // beyond the end of the schedule which only handles "normal" patches
         for
         (
-            label interfaceI=patchSchedule.size()/2;
-            interfaceI<interfaces.size();
+            label interfaceI = patchSchedule.size()/2;
+            interfaceI < interfaces_.size();
             interfaceI++
         )
         {
-            if (interfaces.set(interfaceI))
+            if (interfaces_.set(interfaceI))
             {
-                interfaces[interfaceI].initInterfaceMatrixUpdate
-                (
-                    psi,
-                    result,
-                    *this,
-                    interfaceCoeffs[interfaceI],
-                    Pstream::blocking
-                );
+//                 interfaces_[interfaceI].initInterfaceMatrixUpdate
+//                 (
+//                     psi,
+//                     result,
+//                     *this,
+//                     interfaceCoeffs[interfaceI],
+//                     Pstream::blocking
+//                 );
             }
         }
-    }    
+    }
     else
     {
         FatalErrorIn("BlockLduMatrix<Type>::initMatrixInterfaces")
@@ -101,7 +100,6 @@ template<class Type>
 void Foam::BlockLduMatrix<Type>::updateInterfaces
 (
     const FieldField<CoeffField, Type>& interfaceCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     TypeField& result,
     const TypeField& psi
 ) const
@@ -119,18 +117,18 @@ void Foam::BlockLduMatrix<Type>::updateInterfaces
             OPstream::waitRequests();
         }
 
-        forAll (interfaces, interfaceI)
+        forAll (interfaces_, interfaceI)
         {
-            if (interfaces.set(interfaceI))
+            if (interfaces_.set(interfaceI))
             {
-                interfaces[interfaceI].updateInterfaceMatrix
-                (
-                    psi,
-                    result,
-                    *this,
-                    interfaceCoeffs[interfaceI],
-                    Pstream::defaultCommsType
-                );
+//                 interfaces_[interfaceI].updateInterfaceMatrix
+//                 (
+//                     psi,
+//                     result,
+//                     *this,
+//                     interfaceCoeffs[interfaceI],
+//                     Pstream::defaultCommsType
+//                 );
             }
         }
     }
@@ -143,29 +141,29 @@ void Foam::BlockLduMatrix<Type>::updateInterfaces
         {
             label interfaceI = patchSchedule[i].patch;
 
-            if (interfaces.set(interfaceI))
+            if (interfaces_.set(interfaceI))
             {
                 if (patchSchedule[i].init)
                 {
-                    interfaces[interfaceI].initInterfaceMatrixUpdate
-                    (
-                        psi,
-                        result,
-                        *this,
-                        interfaceCoeffs[interfaceI],
-                        Pstream::scheduled
-                    );
+//                     interfaces_[interfaceI].initInterfaceMatrixUpdate
+//                     (
+//                         psi,
+//                         result,
+//                         *this,
+//                         interfaceCoeffs[interfaceI],
+//                         Pstream::scheduled
+//                     );
                 }
                 else
                 {
-                    interfaces[interfaceI].updateInterfaceMatrix
-                    (
-                        psi,
-                        result,
-                        *this,
-                        interfaceCoeffs[interfaceI],
-                        Pstream::scheduled
-                    );
+//                     interfaces_[interfaceI].updateInterfaceMatrix
+//                     (
+//                         psi,
+//                         result,
+//                         *this,
+//                         interfaceCoeffs[interfaceI],
+//                         Pstream::scheduled
+//                     );
                 }
             }
         }
@@ -174,21 +172,21 @@ void Foam::BlockLduMatrix<Type>::updateInterfaces
         // beyond the end of the schedule which only handles "normal" patches
         for
         (
-            label interfaceI=patchSchedule.size()/2;
-            interfaceI<interfaces.size();
+            label interfaceI = patchSchedule.size()/2;
+            interfaceI < interfaces_.size();
             interfaceI++
         )
         {
-            if (interfaces.set(interfaceI))
+            if (interfaces_.set(interfaceI))
             {
-                interfaces[interfaceI].updateInterfaceMatrix
-                (
-                    psi,
-                    result,
-                    *this,
-                    interfaceCoeffs[interfaceI],
-                    Pstream::blocking
-                );
+//                 interfaces_[interfaceI].updateInterfaceMatrix
+//                 (
+//                     psi,
+//                     result,
+//                     *this,
+//                     interfaceCoeffs[interfaceI],
+//                     Pstream::blocking
+//                 );
             }
         }
     }
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockGaussSeidelSmoother/BlockGaussSeidelSmoother.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockGaussSeidelSmoother/BlockGaussSeidelSmoother.H
index 250daf470..c6172f2d4 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockGaussSeidelSmoother/BlockGaussSeidelSmoother.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockGaussSeidelSmoother/BlockGaussSeidelSmoother.H
@@ -83,12 +83,10 @@ public:
         BlockGaussSeidelSmoother
         (
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
             const dictionary& dict
         )
         :
-            BlockLduSmoother<Type>(matrix, boundaryCoeffs, interfaces),
+            BlockLduSmoother<Type>(matrix),
             gs_(matrix)
         {}
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockILUSmoother/BlockILUSmoother.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockILUSmoother/BlockILUSmoother.H
index e88f7dd07..7bda49efb 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockILUSmoother/BlockILUSmoother.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockILUSmoother/BlockILUSmoother.H
@@ -86,12 +86,10 @@ public:
         BlockILUSmoother
         (
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
             const dictionary& dict
         )
         :
-            BlockLduSmoother<Type>(matrix, boundaryCoeffs, interfaces),
+            BlockLduSmoother<Type>(matrix),
             precon_(matrix),
             xCorr_(matrix.lduAddr().size()),
             residual_(matrix.lduAddr().size())
@@ -117,13 +115,7 @@ public:
             for (label sweep = 0; sweep < nSweeps; sweep++)
             {
                 // Calculate residual
-                this-> matrix_.Amul
-                (
-                    residual_,
-                    x,
-                    BlockLduSmoother<Type>::boundaryCoeffs_,
-                    BlockLduSmoother<Type>::interfaces_
-                );
+                this-> matrix_.Amul(residual_, x);
 
                 // residual = b - Ax
                 forAll (b, i)
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.C
index bc1aae095..26a54736d 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.C
@@ -41,8 +41,6 @@ template<class Type>
 Foam::autoPtr<Foam::BlockLduSmoother<Type> > Foam::BlockLduSmoother<Type>::New
 (
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 {
@@ -72,8 +70,6 @@ Foam::autoPtr<Foam::BlockLduSmoother<Type> > Foam::BlockLduSmoother<Type>::New
             "autoPtr<BlockLduSmoother> BlockLduSmoother::New\n"
             "(\n"
             "    const BlockLduMatrix<Type>& matrix,\n"
-            "    const FieldField<CoeffField, Type>& boundaryCoeffs,\n"
-            "    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,\n"
             "    const dictionary& dict\n"
             ")",
             dict
@@ -89,8 +85,6 @@ Foam::autoPtr<Foam::BlockLduSmoother<Type> > Foam::BlockLduSmoother<Type>::New
         constructorIter()
         (
             matrix,
-            boundaryCoeffs,
-            interfaces,
             dict
         )
     );
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.H
index 5da53b7a9..35f9ef0c7 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/BlockLduSmoother.H
@@ -66,12 +66,6 @@ protected:
         //- Matrix reference
         const BlockLduMatrix<Type>& matrix_;
 
-        //- Reference to boundary coefficients
-        const FieldField<CoeffField, Type>& boundaryCoeffs_;
-
-        //- Reference to interfaces
-        const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces_;
-
 
     // Protected member functions
 
@@ -94,14 +88,10 @@ public:
             dictionary,
             (
                 const BlockLduMatrix<Type>& matrix,
-                const FieldField<CoeffField, Type>& boundaryCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
                 const dictionary& dict
             ),
             (
                 matrix,
-                boundaryCoeffs,
-                interfaces,
                 dict
             )
         );
@@ -110,16 +100,9 @@ public:
     // Constructors
 
         //- Construct from matrix
-        BlockLduSmoother
-        (
-            const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces
-        )
+        BlockLduSmoother(const BlockLduMatrix<Type>& matrix)
         :
-            matrix_(matrix),
-            boundaryCoeffs_(boundaryCoeffs),
-            interfaces_(interfaces)
+            matrix_(matrix)
         {}
 
 
@@ -129,8 +112,6 @@ public:
         static autoPtr<BlockLduSmoother<Type> > New
         (
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
             const dictionary& dict
         );
 
@@ -172,8 +153,8 @@ addToRunTimeSelectionTable(SmootherType, typeSmootherType, dictionary);
 #define makeBlockSmoothers(smootherType)                                       \
                                                                              \
 makeBlockSmoother(blockScalarSmoother, smootherType##Scalar);                  \
-makeBlockSmoother(blockVectorSmoother, smootherType##Vector);
-// makeBlockSmoother(blockTensorSmoother, smootherType##Tensor);
+makeBlockSmoother(blockVectorSmoother, smootherType##Vector);                  \
+makeBlockSmoother(blockTensorSmoother, smootherType##Tensor);
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.C
index d9d22d8af..8e7cfec9a 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.C
@@ -38,11 +38,11 @@ namespace Foam
 
 defineNamedTemplateTypeNameAndDebug(blockScalarSmoother, 0);
 defineNamedTemplateTypeNameAndDebug(blockVectorSmoother, 0);
-// defineNamedTemplateTypeNameAndDebug(blockTensorSmoother, 0);
+defineNamedTemplateTypeNameAndDebug(blockTensorSmoother, 0);
 
 defineTemplateRunTimeSelectionTable(blockScalarSmoother, dictionary);
 defineTemplateRunTimeSelectionTable(blockVectorSmoother, dictionary);
-// defineTemplateRunTimeSelectionTable(blockTensorSmoother, dictionary);
+defineTemplateRunTimeSelectionTable(blockTensorSmoother, dictionary);
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.H
index 15a3a4f7a..377184e50 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSmoothers/BlockLduSmoother/blockLduSmoothers.H
@@ -50,7 +50,7 @@ namespace Foam
 
 typedef BlockLduSmoother<scalar> blockScalarSmoother;
 typedef BlockLduSmoother<vector> blockVectorSmoother;
-// typedef BlockLduSmoother<tensor> blockTensorSmoother;
+typedef BlockLduSmoother<tensor> blockTensorSmoother;
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.C
index 28b532c01..c3acf2ac3 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.C
@@ -37,8 +37,6 @@ Foam::BlockBiCGStabSolver<Type>::BlockBiCGStabSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
@@ -46,8 +44,6 @@ Foam::BlockBiCGStabSolver<Type>::BlockBiCGStabSolver
     (
         fieldName,
         matrix,
-        boundaryCoeffs,
-        interfaces,
         dict
     ),
     preconPtr_
@@ -89,13 +85,7 @@ Foam::BlockBiCGStabSolver<Type>::solve
     Field<Type> p(x.size());
 
     // Calculate initial residual
-    matrix.Amul
-    (
-        p,
-        x,
-        BlockLduSolver<Type>::boundaryCoeffs_,
-        BlockLduSolver<Type>::interfaces_
-    );
+    matrix.Amul(p, x);
     Field<Type> r(b - p);
 
     solverPerf.initialResidual() = gSum(cmptMag(r))/norm;
@@ -148,13 +138,7 @@ Foam::BlockBiCGStabSolver<Type>::solve
             }
 
             preconPtr_->precondition(ph, p);
-            matrix.Amul
-            (
-                v,
-                ph,
-                BlockLduSolver<Type>::boundaryCoeffs_,
-                BlockLduSolver<Type>::interfaces_
-            );
+            matrix.Amul(v, ph);
             alpha = rho/gSumProd(rw, v);
 
             forAll (s, i)
@@ -163,13 +147,7 @@ Foam::BlockBiCGStabSolver<Type>::solve
             }
 
             preconPtr_->preconditionT(sh, s);
-            matrix.Amul
-            (
-                t,
-                sh,
-                BlockLduSolver<Type>::boundaryCoeffs_,
-                BlockLduSolver<Type>::interfaces_
-            );
+            matrix.Amul(t, sh);
             omega = gSumProd(t, s)/gSumProd(t, t);
 
             forAll (x, i)
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.H
index 04a992c5b..8d676de8c 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockBiCGStab/BlockBiCGStabSolver.H
@@ -85,10 +85,7 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
-            const dictionary& solver
+            const dictionary& dict
         );
 
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.C
index 610e0a9fd..90f3015be 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.C
@@ -37,8 +37,6 @@ Foam::BlockCGSolver<Type>::BlockCGSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
@@ -46,8 +44,6 @@ Foam::BlockCGSolver<Type>::BlockCGSolver
     (
         fieldName,
         matrix,
-        boundaryCoeffs,
-        interfaces,
         dict
     ),
     preconPtr_
@@ -88,13 +84,7 @@ typename Foam::BlockSolverPerformance<Type> Foam::BlockCGSolver<Type>::solve
     Field<Type> wA(x.size());
 
     // Calculate initial residual
-    matrix.Amul
-    (
-        wA,
-        x,
-        BlockLduSolver<Type>::boundaryCoeffs_,
-        BlockLduSolver<Type>::interfaces_
-    );
+    matrix.Amul(wA, x);
     Field<Type> rA(b - wA);
 
     solverPerf.initialResidual() = gSum(cmptMag(rA))/norm;
@@ -129,13 +119,7 @@ typename Foam::BlockSolverPerformance<Type> Foam::BlockCGSolver<Type>::solve
             }
 
             // Update preconditioner residual
-            matrix.Amul
-            (
-                wA,
-                pA,
-                BlockLduSolver<Type>::boundaryCoeffs_,
-                BlockLduSolver<Type>::interfaces_
-            );
+            matrix.Amul(wA, pA);
 
             wApA = gSumProd(wA, pA);
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.H
index b7d459fa6..a18d0a012 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockCG/BlockCGSolver.H
@@ -85,9 +85,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.C
index abb22057f..663f62175 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.C
@@ -38,19 +38,10 @@ Foam::BlockDiagonalSolver<Type>::BlockDiagonalSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
-    BlockLduSolver<Type>
-    (
-        fieldName,
-        matrix,
-        boundaryCoeffs,
-        interfaces,
-        dict
-    )
+    BlockLduSolver<Type>(fieldName, matrix, dict)
 {}
 
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.H
index e5a7cd70e..64d2951c9 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockDiagonal/BlockDiagonalSolver.H
@@ -87,9 +87,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.C
index 255704687..b91e1e56a 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.C
@@ -71,8 +71,6 @@ Foam::BlockGMRESSolver<Type>::BlockGMRESSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
@@ -80,8 +78,6 @@ Foam::BlockGMRESSolver<Type>::BlockGMRESSolver
     (
         fieldName,
         matrix,
-        boundaryCoeffs,
-        interfaces,
         dict
     ),
     preconPtr_
@@ -124,13 +120,7 @@ Foam::BlockGMRESSolver<Type>::solve
     Field<Type> wA(x.size());
 
     // Calculate initial residual
-    matrix.Amul
-    (
-        wA,
-        x,
-        BlockLduSolver<Type>::boundaryCoeffs_,
-        BlockLduSolver<Type>::interfaces_
-    );
+    matrix.Amul(wA, x);
     Field<Type> rA(b - wA);
 
     solverPerf.initialResidual() = gSum(cmptMag(rA))/norm;
@@ -178,13 +168,7 @@ Foam::BlockGMRESSolver<Type>::solve
                 V[i] /= beta;
 
                 // Arnoldi's method
-                matrix.Amul
-                (
-                    rA,
-                    V[i],
-                    BlockLduSolver<Type>::boundaryCoeffs_,
-                    BlockLduSolver<Type>::interfaces_
-                );
+                matrix.Amul(rA, V[i]);
 
                 // Execute preconditioning
                 preconPtr_->precondition(wA, rA);
@@ -247,13 +231,7 @@ Foam::BlockGMRESSolver<Type>::solve
             }
 
             // Re-calculate the residual
-            matrix.Amul
-            (
-                wA,
-                x,
-                BlockLduSolver<Type>::boundaryCoeffs_,
-                BlockLduSolver<Type>::interfaces_
-            );
+            matrix.Amul(wA, x);
 
             forAll (rA, raI)
             {
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.H
index 94615a6bd..2844cbffa 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGMRES/BlockGMRESSolver.H
@@ -98,9 +98,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.C
index ee887a798..0819dc9d3 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.C
@@ -39,8 +39,6 @@ Foam::BlockGaussSeidelSolver<Type>::BlockGaussSeidelSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
@@ -48,8 +46,6 @@ Foam::BlockGaussSeidelSolver<Type>::BlockGaussSeidelSolver
     (
         fieldName,
         matrix,
-        boundaryCoeffs,
-        interfaces,
         dict
     ),
     gs_(matrix),
@@ -82,13 +78,7 @@ Foam::BlockGaussSeidelSolver<Type>::solve
     Field<Type> wA(x.size());
 
     // Calculate residual.  Note: sign of residual swapped for efficiency
-    matrix.Amul
-    (
-        wA,
-        x,
-        BlockLduSolver<Type>::boundaryCoeffs_,
-        BlockLduSolver<Type>::interfaces_
-    );
+    matrix.Amul(wA, x);
     wA -= b;
 
     solverPerf.initialResidual() = gSum(cmptMag(wA))/norm;
@@ -111,13 +101,7 @@ Foam::BlockGaussSeidelSolver<Type>::solve
 
             // Re-calculate residual.  Note: sign of residual swapped
             // for efficiency
-            matrix.Amul
-            (
-                wA,
-                x,
-                BlockLduSolver<Type>::boundaryCoeffs_,
-                BlockLduSolver<Type>::interfaces_
-            );
+            matrix.Amul(wA, x);
             wA -= b;
 
             solverPerf.finalResidual() = gSum(cmptMag(wA))/norm;
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.H
index b7a4bf7f5..a29fbcc0c 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockGaussSeidel/BlockGaussSeidelSolver.H
@@ -90,9 +90,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.C
index b33a4910b..d4e8dd0e4 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.C
@@ -37,12 +37,10 @@ Foam::BlockIterativeSolver<Type>::BlockIterativeSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
-    BlockLduSolver<Type>(fieldName, matrix, boundaryCoeffs, interfaces, dict),
+    BlockLduSolver<Type>(fieldName, matrix, dict),
     tolerance_(readScalar(this->dict().lookup("tolerance"))),
     relTolerance_(readScalar(this->dict().lookup("relTol"))),
     minIter_(readLabel(this->dict().lookup("minIter"))),
@@ -71,21 +69,13 @@ Foam::scalar Foam::BlockIterativeSolver<Type>::normFactor
     Type xRef = gAverage(x);
 
     // Calculate A.x
-    matrix.Amul
-    (
-        wA,
-        x,
-        BlockLduSolver<Type>::boundaryCoeffs_,
-        BlockLduSolver<Type>::interfaces_
-    );
+    matrix.Amul(wA, x);
 
     // Calculate A.xRef, temporarily using pA for storage
     matrix.Amul
     (
         pA,
-        Field<Type>(nRows, xRef),
-        BlockLduSolver<Type>::boundaryCoeffs_,
-        BlockLduSolver<Type>::interfaces_
+        Field<Type>(nRows, xRef)
     );
 
     scalar normFactor = gSum(mag(wA - pA) + mag(b - pA)) + this->small_;
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.H
index 6a1f46ad9..c34033c48 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockIterativeSolver/BlockIterativeSolver.H
@@ -104,9 +104,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.C
index 9bf52ada6..bd3ec1be7 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.C
@@ -76,16 +76,12 @@ template<class Type>
 Foam::BlockLduSolver<Type>::BlockLduSolver
 (
     const word& fieldName,
-    const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces
+    const BlockLduMatrix<Type>& matrix
 )
 :
     fieldName_(fieldName),
     dict_(),
-    matrix_(matrix),
-    boundaryCoeffs_(boundaryCoeffs),
-    interfaces_(interfaces)
+    matrix_(matrix)
 {}
 
 
@@ -94,16 +90,12 @@ Foam::BlockLduSolver<Type>::BlockLduSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
     fieldName_(fieldName),
     dict_(dict),
-    matrix_(matrix),
-    boundaryCoeffs_(boundaryCoeffs),
-    interfaces_(interfaces)
+    matrix_(matrix)
 {}
 
 
@@ -115,8 +107,6 @@ Foam::BlockLduSolver<Type>::New
 (
     const word& fieldName,
     BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 {
@@ -126,14 +116,7 @@ Foam::BlockLduSolver<Type>::New
     {
         return autoPtr<BlockLduSolver<Type> >
         (
-            new BlockDiagonalSolver<Type>
-            (
-                fieldName,
-                matrix,
-                boundaryCoeffs,
-                interfaces,
-                dict
-            )
+            new BlockDiagonalSolver<Type>(fieldName, matrix, dict)
         );
     }
     else if (matrix.symmetric())
@@ -146,8 +129,6 @@ Foam::BlockLduSolver<Type>::New
                 "(\n"
                 "    const word& fieldName,\n"
                 "    BlockLduMatrix<Type>& matrix,\n"
-                "    const FieldField<CoeffField, Type>& boundaryCoeffs,\n"
-                "    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,\n"
                 "    const dictionary& dict\n"
                 ")"
             )   << "Initialization problem." << endl;
@@ -164,8 +145,6 @@ Foam::BlockLduSolver<Type>::New
                 "(\n"
                 "    const word& fieldName,\n"
                 "    BlockLduMatrix<Type>& matrix,\n"
-                "    const FieldField<CoeffField, Type>& boundaryCoeffs,\n"
-                "    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,\n"
                 "    const dictionary& dict\n"
                 ")",
                 dict
@@ -182,8 +161,6 @@ Foam::BlockLduSolver<Type>::New
             (
                 fieldName,
                 matrix,
-                boundaryCoeffs,
-                interfaces,
                 dict
             )
         );
@@ -198,8 +175,6 @@ Foam::BlockLduSolver<Type>::New
                 "(\n"
                 "    const word& fieldName,\n"
                 "    BlockLduMatrix<Type>& matrix,\n"
-                "    const FieldField<CoeffField, Type>& boundaryCoeffs,\n"
-                "    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,\n"
                 "    const dictionary& dict\n"
                 ")"
             )   << "Initialization problem." << endl;
@@ -216,8 +191,6 @@ Foam::BlockLduSolver<Type>::New
                 "(\n"
                 "    const word& fieldName,\n"
                 "    BlockLduMatrix<Type>& matrix,\n"
-                "    const FieldField<CoeffField, Type>& boundaryCoeffs,\n"
-                "    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,\n"
                 "    const dictionary& dict\n"
                 ")",
                 dict
@@ -234,8 +207,6 @@ Foam::BlockLduSolver<Type>::New
             (
                 fieldName,
                 matrix,
-                boundaryCoeffs,
-                interfaces,
                 dict
             )
         );
@@ -248,9 +219,6 @@ Foam::BlockLduSolver<Type>::New
             "(\n"
             "    const word& fieldName,\n"
             "    BlockLduMatrix<Type>& matrix,\n"
-            "    const FieldField<CoeffField, Type>& boundaryCoeffs,\n"
-            "    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& "
-            "interfaces,\n"
             "    const dictionary& dict\n"
             ")"
         )   << "cannot solve incomplete matrix, "
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.H
index ef51e784f..62aa78379 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/BlockLduSolver.H
@@ -67,6 +67,7 @@ class BlockLduSolver
         //- Control data dictionary
         dictionary dict_;
 
+
 protected:
 
     // Protected data types
@@ -79,11 +80,6 @@ protected:
         //- Matrix
         const BlockLduMatrix<Type>& matrix_;
 
-        // Reference to boundary coefficients
-        const FieldField<CoeffField, Type>& boundaryCoeffs_;
-
-        //- Reference to interfaces
-        const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces_;
 
     // Protected Member Functions
 
@@ -128,16 +124,11 @@ public:
             (
                 const word& fieldName,
                 BlockLduMatrix<Type>& matrix,
-                const FieldField<CoeffField, Type>& boundaryCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                    interfaces,
                 const dictionary& dict
             ),
             (
                 fieldName,
                 matrix,
-                boundaryCoeffs,
-                interfaces,
                 dict
             )
         );
@@ -150,16 +141,11 @@ public:
             (
                 const word& fieldName,
                 BlockLduMatrix<Type>& matrix,
-                const FieldField<CoeffField, Type>& boundaryCoeffs,
-                const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                    interfaces,
                 const dictionary& dict
             ),
             (
                 fieldName,
                 matrix,
-                boundaryCoeffs,
-                interfaces,
                 dict
             )
         );
@@ -172,10 +158,7 @@ public:
         BlockLduSolver
         (
             const word& fieldName,
-            const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces
+            const BlockLduMatrix<Type>& matrix
         );
 
         //- Construct from dictionary
@@ -183,9 +166,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
@@ -197,9 +177,6 @@ public:
         (
             const word& fieldName,
             BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
@@ -245,8 +222,8 @@ defineNamedTemplateTypeNameAndDebug(typeSolverType, 0);
 #define makeBlockSolvers(solverType)                                         \
                                                                              \
 makeBlockSolverTypeName(solverType##Scalar);                                 \
-makeBlockSolverTypeName(solverType##Vector);
-// makeBlockSolverTypeName(solverType##Tensor);
+makeBlockSolverTypeName(solverType##Vector);                                 \
+makeBlockSolverTypeName(solverType##Tensor);
 
 
 #define addSolverToBlockMatrix(Type, typeSolverType, typeMatrix)             \
@@ -256,8 +233,8 @@ addToRunTimeSelectionTable(block##Type##Solver, typeSolverType, typeMatrix);
 #define addSolversToBlockMatrix(solverType, typeMatrix)                      \
                                                                              \
 addSolverToBlockMatrix(Scalar, solverType##Scalar, typeMatrix);              \
-addSolverToBlockMatrix(Vector, solverType##Vector, typeMatrix);
-// addSolverToBlockMatrix(Tensor, solverType##Tensor, typeMatrix);
+addSolverToBlockMatrix(Vector, solverType##Vector, typeMatrix);              \
+addSolverToBlockMatrix(Tensor, solverType##Tensor, typeMatrix);
 
 #define addSymSolverToBlockMatrix(solverType)                                \
                                                                              \
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.C
index 6fba557ce..b40e7907f 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.C
@@ -38,7 +38,7 @@ namespace Foam
 
 defineNamedTemplateTypeNameAndDebug(blockScalarSolver, 0);
 defineNamedTemplateTypeNameAndDebug(blockVectorSolver, 0);
-// defineNamedTemplateTypeNameAndDebug(blockTensorSolver, 0);
+defineNamedTemplateTypeNameAndDebug(blockTensorSolver, 0);
 
 
 // Define the constructor function hash tables for symmetric solvers
@@ -55,11 +55,11 @@ defineTemplateRunTimeSelectionTable
     symMatrix
 );
 
-// defineTemplateRunTimeSelectionTable
-// (
-//     blockTensorSolver,
-//     symMatrix
-// );
+defineTemplateRunTimeSelectionTable
+(
+    blockTensorSolver,
+    symMatrix
+);
 
 
 // Define the constructor function hash tables for asymmetric solvers
@@ -76,11 +76,11 @@ defineTemplateRunTimeSelectionTable
     asymMatrix
 );
 
-// defineTemplateRunTimeSelectionTable
-// (
-//     blockTensorSolver,
-//     asymMatrix
-// );
+defineTemplateRunTimeSelectionTable
+(
+    blockTensorSolver,
+    asymMatrix
+);
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.H
index 964a10fea..4d0f526cf 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/BlockLduSolver/blockLduSolvers.H
@@ -44,7 +44,7 @@ namespace Foam
 
 typedef BlockLduSolver<scalar> blockScalarSolver;
 typedef BlockLduSolver<vector> blockVectorSolver;
-// typedef BlockLduSolver<tensor> blockTensorSolver;
+typedef BlockLduSolver<tensor> blockTensorSolver;
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.C
index 0ab5397d1..46aa30edf 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.C
@@ -38,19 +38,10 @@ Foam::SegregatedSolver<Type>::SegregatedSolver
 (
     const word& fieldName,
     const BlockLduMatrix<Type>& matrix,
-    const FieldField<CoeffField, Type>& boundaryCoeffs,
-    const typename BlockLduInterfaceFieldPtrsList<Type>::Type& interfaces,
     const dictionary& dict
 )
 :
-    BlockLduSolver<Type>
-    (
-        fieldName,
-        matrix,
-        boundaryCoeffs,
-        interfaces,
-        dict
-    ),
+    BlockLduSolver<Type>(fieldName, matrix, dict),
     scalarX_(matrix.lduAddr().size()),
     scalarMatrix_(matrix.mesh()),
     scalarB_(matrix.lduAddr().size())
@@ -181,8 +172,6 @@ Foam::BlockSolverPerformance<Type> Foam::SegregatedSolver<Type>::solve
             (
                 this->fieldName(),
                 scalarMatrix_,
-                boundaryCoeffs,
-                interfaces_,
                 this->dict()
             )->solve(scalarX_, scalarB_);
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.H b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.H
index f88d16e1b..4ea0d51be 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.H
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/SegregatedSolver.H
@@ -92,9 +92,6 @@ public:
         (
             const word& fieldName,
             const BlockLduMatrix<Type>& matrix,
-            const FieldField<CoeffField, Type>& boundaryCoeffs,
-            const typename BlockLduInterfaceFieldPtrsList<Type>::Type&
-                interfaces,
             const dictionary& dict
         );
 
diff --git a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/segregatedSolvers.C b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/segregatedSolvers.C
index add717747..57b39747d 100644
--- a/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/segregatedSolvers.C
+++ b/src/OpenFOAM/matrices/blockLduMatrix/BlockLduSolvers/Segregated/segregatedSolvers.C
@@ -40,10 +40,10 @@ makeBlockSolverTypeName(segregatedSolverVector);
 makeBlockSolverTypeName(segregatedSolverTensor);
 
 addSolverToBlockMatrix(Vector, segregatedSolverVector, symMatrix);
-// addSolverToBlockMatrix(Tensor, segregatedSolverTensor, symMatrix);
+addSolverToBlockMatrix(Tensor, segregatedSolverTensor, symMatrix);
 
 addSolverToBlockMatrix(Vector, segregatedSolverVector, asymMatrix);
-// addSolverToBlockMatrix(Tensor, segregatedSolverTensor, asymMatrix);
+addSolverToBlockMatrix(Tensor, segregatedSolverTensor, asymMatrix);
 
 
 // * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * * //