From 6b1fdc5132c0e39c699fcfd1de8d32dd5fc3494d Mon Sep 17 00:00:00 2001 From: Jurgen Van Gael Date: Sat, 15 May 2010 17:45:31 +0800 Subject: [PATCH] Added Cholesky factorization. --- .../Double/Decomposition/Cholesky.cs | 78 ++++++++++++ .../Double/Decomposition/DenseCholesky.cs | 117 ++++++++++++++++++ .../Double/Decomposition/ExtensionMethods.cs | 51 ++++++++ src/Numerics/Numerics.csproj | 4 + .../Double/Decomposition/CholeskyTests.cs | 98 +++++++++++++++ src/UnitTests/UnitTests.csproj | 1 + 6 files changed, 349 insertions(+) create mode 100644 src/Numerics/LinearAlgebra/Double/Decomposition/Cholesky.cs create mode 100644 src/Numerics/LinearAlgebra/Double/Decomposition/DenseCholesky.cs create mode 100644 src/Numerics/LinearAlgebra/Double/Decomposition/ExtensionMethods.cs create mode 100644 src/UnitTests/LinearAlgebraTests/Double/Decomposition/CholeskyTests.cs diff --git a/src/Numerics/LinearAlgebra/Double/Decomposition/Cholesky.cs b/src/Numerics/LinearAlgebra/Double/Decomposition/Cholesky.cs new file mode 100644 index 00000000..3246596c --- /dev/null +++ b/src/Numerics/LinearAlgebra/Double/Decomposition/Cholesky.cs @@ -0,0 +1,78 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// +// Copyright (c) 2009-2010 Math.NET +// +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double.Decomposition +{ + using System; + using Properties; + + /// + /// A class which encapsulates the functionality of a Cholesky decomposition. + /// For a symmetric, positive definite matrix A, the Cholesky decomposition + /// is an lower triangular matrix L so that A = L*L'. + /// + /// + /// The computation of the Cholesky decomposition is done at construction time. If the matrix is not symmetric + /// or positive definite, the constructor will throw an exception. + /// + public abstract class Cholesky + { + /// + /// Internal method which routes the call to perform the Cholesky factorization to the appropriate class. + /// + /// The matrix to factor. + /// A cholesky factorization object. + internal static Cholesky Create(Matrix matrix) + { + var dense = matrix as DenseMatrix; + if (dense != null) + { + return new DenseCholesky(dense); + } + + throw new NotImplementedException(); + } + + /// + /// A reference to the lower triangular form of the Cholesky matrix. + /// + public abstract Matrix Factor { get; } + + /// + /// The determinant of the matrix for which the Cholesky matrix was computed. + /// + public abstract double Determinant { get; } + + /// + /// The log determinant of the matrix for which the Cholesky matrix was computed. + /// + public abstract double DeterminantLn { get; } + } +} diff --git a/src/Numerics/LinearAlgebra/Double/Decomposition/DenseCholesky.cs b/src/Numerics/LinearAlgebra/Double/Decomposition/DenseCholesky.cs new file mode 100644 index 00000000..1f03a90c --- /dev/null +++ b/src/Numerics/LinearAlgebra/Double/Decomposition/DenseCholesky.cs @@ -0,0 +1,117 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// +// Copyright (c) 2009-2010 Math.NET +// +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double.Decomposition +{ + using System; + using Properties; + + /// + /// A class which encapsulates the functionality of a Cholesky decomposition for dense matrices. + /// For a symmetric, positive definite matrix A, the Cholesky decomposition + /// is an lower triangular matrix L so that A = L*L'. + /// + /// + /// The computation of the Cholesky decomposition is done at construction time. If the matrix is not symmetric + /// or positive definite, the constructor will throw an exception. + /// + public class DenseCholesky : Cholesky + { + /// + /// Stores the Cholesky factor for the decomposition. + /// + private DenseMatrix _factor; + + /// + /// Initializes a new instance of the class. This object will compute the + /// Cholesky decomposition when the constructor is called and cache it's factorization. + /// + /// The matrix to factor. + /// If is null. + /// If is not a square matrix. + /// If is not positive definite. + public DenseCholesky(DenseMatrix matrix) + { + if (matrix == null) + { + throw new ArgumentNullException("matrix"); + } + + if (matrix.RowCount != matrix.ColumnCount) + { + throw new ArgumentException(Resources.ArgumentMatrixSquare); + } + + // Create a new matrix for the Cholesky factor, then perform factorization (while overwriting). + _factor = (DenseMatrix)matrix.Clone(); + Control.LinearAlgebraProvider.CholeskyFactor(_factor.Data, _factor.RowCount); + } + + /// + /// Returns the lower triangular form of the Cholesky matrix. + /// + public override Matrix Factor + { + get { return _factor; } + } + + /// + /// The determinant of the matrix for which the Cholesky matrix was computed. + /// + public override double Determinant + { + get + { + double det = 1.0; + for (int j = 0; j < _factor.RowCount; j++) + { + det *= (_factor[j, j] * _factor[j, j]); + } + return det; + } + } + + /// + /// The log determinant of the matrix for which the Cholesky matrix was computed. + /// + public override double DeterminantLn + { + get + { + double det = 0.0; + for (int j = 0; j < _factor.RowCount; j++) + { + det += (2.0 * Math.Log(_factor[j, j])); + } + return det; + } + } + } +} diff --git a/src/Numerics/LinearAlgebra/Double/Decomposition/ExtensionMethods.cs b/src/Numerics/LinearAlgebra/Double/Decomposition/ExtensionMethods.cs new file mode 100644 index 00000000..d2653edb --- /dev/null +++ b/src/Numerics/LinearAlgebra/Double/Decomposition/ExtensionMethods.cs @@ -0,0 +1,51 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// +// Copyright (c) 2009-2010 Math.NET +// +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double.Decomposition +{ + using System; + using Properties; + + /// + /// Extension methods which return decompositions for the various matrix classes. + /// + public static class ExtensionMethods + { + /// + /// Computes the cholesky decomposition for a matrix. + /// + /// The matrix to factor. + /// The Cholesky decomposition object. + public static Cholesky Cholesky(this Matrix matrix) + { + return Decomposition.Cholesky.Create(matrix); + } + } +} diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj index 17b808ed..79f6d96f 100644 --- a/src/Numerics/Numerics.csproj +++ b/src/Numerics/Numerics.csproj @@ -139,6 +139,9 @@ + + + @@ -246,6 +249,7 @@ true +