diff --git a/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.C b/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.C
index e16018024036947c704cdfa723be431a93ff4265..2e8f7421d951405d3dc7826cfd494adf99ca6837 100644
--- a/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.C
+++ b/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.C
@@ -27,6 +27,11 @@ License
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
+template<class Type>
+Foam::LLTMatrix<Type>::LLTMatrix()
+{}
+
+
 template<class Type>
 Foam::LLTMatrix<Type>::LLTMatrix(const SquareMatrix<Type>& M)
 {
diff --git a/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.H b/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.H
index 9a48c38cd708dc30b265fec72a3390ba7b7276ed..622efc5a62f5a263819039e8fdc14fa33dce20c9 100644
--- a/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.H
+++ b/src/OpenFOAM/matrices/LLTMatrix/LLTMatrix.H
@@ -61,6 +61,9 @@ public:
 
     // Constructors
 
+        //- Construct null
+        LLTMatrix();
+
         //- Construct from a square matrix and perform the decomposition
         LLTMatrix(const SquareMatrix<Type>& M);
 
diff --git a/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.C b/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.C
index 5a4076633d57313f9193c4f0d8f2012b2ff63d94..ecff7ea54c66b7ebee95bc0e9b54bce106fbf3d1 100644
--- a/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.C
+++ b/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.C
@@ -39,6 +39,12 @@ namespace Foam
 
 // * * * * * * * * * * * * * * * * Constructors  * * * * * * * * * * * * * * //
 
+Foam::LUscalarMatrix::LUscalarMatrix()
+:
+    comm_(Pstream::worldComm)
+{}
+
+
 Foam::LUscalarMatrix::LUscalarMatrix(const scalarSquareMatrix& matrix)
 :
     scalarSquareMatrix(matrix),
@@ -249,8 +255,6 @@ void Foam::LUscalarMatrix::convert
             }
         }
     }
-
-    //printDiagonalDominance();
 }
 
 
@@ -380,8 +384,6 @@ void Foam::LUscalarMatrix::convert
             }
         }
     }
-
-    //printDiagonalDominance();
 }
 
 
@@ -402,4 +404,12 @@ void Foam::LUscalarMatrix::printDiagonalDominance() const
 }
 
 
+void Foam::LUscalarMatrix::decompose(const scalarSquareMatrix& M)
+{
+    scalarSquareMatrix::operator=(M);
+    pivotIndices_.setSize(m());
+    LUDecompose(*this, pivotIndices_);
+}
+
+
 // ************************************************************************* //
diff --git a/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.H b/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.H
index fc289fb3a5cf0a1c87fba5411fb9d2ae944b5414..c8a83e5e5a8bc27d7000e86b3726f26e8da581b5 100644
--- a/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.H
+++ b/src/OpenFOAM/matrices/LUscalarMatrix/LUscalarMatrix.H
@@ -25,7 +25,7 @@ Class
     Foam::LUscalarMatrix
 
 Description
-    Foam::LUscalarMatrix
+    Class to perform the LU decomposition on a symmetric matrix.
 
 SourceFiles
     LUscalarMatrix.C
@@ -93,10 +93,14 @@ public:
     // Declare name of the class and its debug switch
     ClassName("LUscalarMatrix");
 
+
     // Constructors
 
-        //- Construct from scalarSquareMatrix and perform LU decomposition
-        LUscalarMatrix(const scalarSquareMatrix&);
+        //- Construct null
+        LUscalarMatrix();
+
+        //- Construct from and perform LU decomposition of the matrix M
+        LUscalarMatrix(const scalarSquareMatrix& M);
 
         //- Construct from lduMatrix and perform LU decomposition
         LUscalarMatrix
@@ -109,6 +113,9 @@ public:
 
     // Member Functions
 
+        //- Perform the LU decomposition of the matrix M
+        void decompose(const scalarSquareMatrix& M);
+
         //- Solve the matrix using the LU decomposition with pivoting
         //  returning the solution in the source
         template<class T>