diff --git a/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs index 8a687a38..c1e591b4 100644 --- a/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs @@ -1021,7 +1021,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override QR QR(QRMethod method = QRMethod.Thin) { - return new DenseQR(this, method); + return DenseQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs index a05ceda7..c40393a8 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs @@ -164,7 +164,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); } /// @@ -199,7 +199,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); } } } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs index bacb264a..267832e8 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,9 +28,9 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; -using System; namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { @@ -39,6 +39,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization using Numerics; #else using System.Numerics; + #endif /// @@ -50,16 +51,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class DenseQR : QR + public sealed class DenseQR : QR { /// /// Gets or sets Tau vector. Contains additional information on Q - used for native solver. /// - public Complex[] Tau - { - get; - set; - } + Complex[] Tau { get; set; } /// /// Initializes a new instance of the class. This object will compute the @@ -69,35 +66,37 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The type of QR factorization to perform. /// If is null. /// If row count is less then column count - public DenseQR(DenseMatrix matrix, QRMethod method = QRMethod.Full) + public static DenseQR Create(DenseMatrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; - Tau = new Complex[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + var tau = new Complex[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + Matrix q; + Matrix r; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = new DenseMatrix(matrix.RowCount); - Control.LinearAlgebraProvider.QRFactor(((DenseMatrix)MatrixR).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)Q).Values, Tau); + r = matrix.Clone(); + q = new DenseMatrix(matrix.RowCount); + Control.LinearAlgebraProvider.QRFactor(((DenseMatrix) r).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) q).Values, tau); } else { - Q = matrix.Clone(); - MatrixR = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix)Q).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)MatrixR).Values, Tau); + q = matrix.Clone(); + r = new DenseMatrix(matrix.ColumnCount); + Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix) q).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) r).Values, tau); } + + return new DenseQR(q, r, method, tau); + } + + DenseQR(Matrix q, Matrix rFull, QRMethod method, Complex[] tau) + : base(q, rFull, method) + { + Tau = tau; } /// @@ -107,17 +106,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -131,7 +119,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -148,7 +136,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do QR factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, Method); } /// @@ -158,16 +146,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries if (Q.RowCount != input.Count) @@ -176,9 +154,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var dinput = input as DenseVector; @@ -193,7 +171,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, Method); } } } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/GramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/GramSchmidt.cs index f5bfeda2..0b7afb8a 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/GramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/GramSchmidt.cs @@ -62,16 +62,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = Complex.One; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0)) + det *= FullR.At(i, i); + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0)) { return 0; } @@ -89,9 +89,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0)) + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/QR.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/QR.cs index 129083c7..344320d6 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/QR.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/QR.cs @@ -3,7 +3,9 @@ // http://numerics.mathdotnet.com // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com -// Copyright (c) 2009-2010 Math.NET +// +// Copyright (c) 2009-2013 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 @@ -12,8 +14,10 @@ // 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 @@ -24,17 +28,18 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; +using MathNet.Numerics.Properties; namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { - using System; - using Properties; #if NOSYSNUMERICS using Complex = Numerics.Complex; #else using Complex = System.Numerics.Complex; + #endif /// @@ -51,6 +56,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// public abstract class QR : QR { + protected QR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { + } + /// /// Gets the absolute determinant value of the matrix for which the QR matrix was computed. /// @@ -58,16 +68,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = Complex.One; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0)) + det *= FullR.At(i, i); + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0)) { return 0; } @@ -85,9 +95,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0)) + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/UserGramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/UserGramSchmidt.cs index 6e7e873a..a34723d1 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/UserGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/UserGramSchmidt.cs @@ -159,19 +159,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j) / FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < input.ColumnCount; j++) { @@ -223,14 +223,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization // Solve R*X = Y; for (var k = Q.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k] * FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/UserQR.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/UserQR.cs index 050dd37e..eb0ce390 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/UserQR.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/UserQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,11 +28,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; +using System.Linq; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; using MathNet.Numerics.Threading; -using System; -using System.Linq; namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { @@ -41,6 +41,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization using Numerics; #else using System.Numerics; + #endif /// @@ -52,7 +53,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class UserQR : QR + public sealed class UserQR : QR { /// /// Initializes a new instance of the class. This object will compute the @@ -61,71 +62,70 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The matrix to factor. /// The QR factorization method to use. /// If is null. - public UserQR(Matrix matrix, QRMethod method = QRMethod.Full) + public static UserQR Create(Matrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; + Matrix q; + Matrix r; + var minmn = Math.Min(matrix.RowCount, matrix.ColumnCount); var u = new Complex[minmn][]; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); + r = matrix.Clone(); + q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); for (var i = 0; i < matrix.RowCount; i++) { - Q.At(i, i, 1.0f); + q.At(i, i, 1.0f); } for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(MatrixR, i, i); - ComputeQR(u[i], MatrixR, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(r, i, i); + ComputeQR(u[i], r, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.RowCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.RowCount, Control.NumberOfParallelWorkerThreads); } } else { - MatrixR = matrix.CreateMatrix(matrix.ColumnCount, matrix.ColumnCount); - Q = matrix.Clone(); + q = matrix.Clone(); for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(Q, i, i); - ComputeQR(u[i], Q, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(q, i, i); + ComputeQR(u[i], q, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } - MatrixR = Q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); - Q.Clear(); + r = q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); + q.Clear(); for (var i = 0; i < matrix.ColumnCount; i++) { - Q.At(i, i, 1.0f); + q.At(i, i, 1.0f); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } } + + return new UserQR(q, r, method); + } + + UserQR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { } /// @@ -135,7 +135,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The first row /// Column index /// Generated vector - private static Complex[] GenerateColumn(Matrix a, int row, int column) + static Complex[] GenerateColumn(Matrix a, int row, int column) { var ru = a.RowCount - row; var u = new Complex[ru]; @@ -146,7 +146,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization a.At(i, column, 0.0); } - var norm = u.Aggregate(Complex.Zero, (current, t) => current + (t.Magnitude * t.Magnitude)); + var norm = u.Aggregate(Complex.Zero, (current, t) => current + (t.Magnitude*t.Magnitude)); norm = norm.SquareRoot(); if (row == a.RowCount - 1 || norm.Magnitude == 0) @@ -158,7 +158,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization if (u[0].Magnitude != 0.0) { - norm = norm.Magnitude * (u[0] / u[0].Magnitude); + norm = norm.Magnitude*(u[0]/u[0].Magnitude); } a.At(row, column, -norm); @@ -170,10 +170,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization u[0] += 1.0; - var s = (1.0 / u[0]).SquareRoot(); + var s = (1.0/u[0]).SquareRoot(); for (var i = 0; i < ru; i++) { - u[i] = u[i].Conjugate() * s; + u[i] = u[i].Conjugate()*s; } return u; @@ -189,7 +189,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The first column /// The last column /// Number of available CPUs - private static void ComputeQR(Complex[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) + static void ComputeQR(Complex[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) { if (rowDim < rowStart || columnDim < columnStart) { @@ -200,8 +200,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization if ((availableCores > 1) && (tmpColCount > 200)) { - var tmpSplit = columnStart + (tmpColCount / 2); - var tmpCores = availableCores / 2; + var tmpSplit = columnStart + (tmpColCount/2); + var tmpCores = availableCores/2; CommonParallel.Invoke( () => ComputeQR(u, a, rowStart, rowDim, columnStart, tmpSplit, tmpCores), @@ -214,12 +214,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization var scale = Complex.Zero; for (var i = rowStart; i < rowDim; i++) { - scale += u[i - rowStart] * a.At(i, j); + scale += u[i - rowStart]*a.At(i, j); } for (var i = rowStart; i < rowDim; i++) { - a.At(i, j, a.At(i, j) - (u[i - rowStart].Conjugate() * scale)); + a.At(i, j, a.At(i, j) - (u[i - rowStart].Conjugate()*scale)); } } } @@ -232,17 +232,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -250,13 +239,13 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization } // The dimension compatibility conditions for X = A\B require the two matrices A and B to have the same number of rows - if (MatrixR.RowCount != input.RowCount) + if (FullR.RowCount != input.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameRowDimension); } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -264,20 +253,20 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new Complex[MatrixR.RowCount]; + var column = new Complex[FullR.RowCount]; for (var j = 0; j < input.ColumnCount; j++) { - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy.At(k, j); } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { var s = Complex.Zero; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i).Conjugate() * column[k]; + s += Q.At(k, i).Conjugate()*column[k]; } inputCopy.At(i, j, s); @@ -285,23 +274,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j)/FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j)*FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < inputCopy.ColumnCount; j++) { @@ -317,60 +306,50 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries - if (MatrixR.RowCount != input.Count) + if (FullR.RowCount != input.Count) { throw new ArgumentException(Resources.ArgumentVectorsSameLength); } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new Complex[MatrixR.RowCount]; - for (var k = 0; k < MatrixR.RowCount; k++) + var column = new Complex[FullR.RowCount]; + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy[k]; } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { var s = Complex.Zero; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i).Conjugate() * column[k]; + s += Q.At(k, i).Conjugate()*column[k]; } inputCopy[i] = s; } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k]*FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Complex/Matrix.cs b/src/Numerics/LinearAlgebra/Complex/Matrix.cs index a95f76bc..01ac95ad 100644 --- a/src/Numerics/LinearAlgebra/Complex/Matrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/Matrix.cs @@ -470,7 +470,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override QR QR(QRMethod method = QRMethod.Thin) { - return new UserQR(this, method); + return UserQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs index 4f247fee..1ec044b3 100644 --- a/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs @@ -1016,7 +1016,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override QR QR(QRMethod method = QRMethod.Thin) { - return new DenseQR(this, method); + return DenseQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs index 8c638910..3167094e 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs @@ -159,7 +159,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); } /// @@ -194,7 +194,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs index 7cd40ed1..32eb9c1d 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,9 +28,9 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; -using System; namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { @@ -45,16 +45,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class DenseQR : QR + public sealed class DenseQR : QR { /// /// Gets or sets Tau vector. Contains additional information on Q - used for native solver. /// - public Complex32[] Tau - { - get; - set; - } + Complex32[] Tau { get; set; } /// /// Initializes a new instance of the class. This object will compute the @@ -64,35 +60,37 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The QR factorization method to use. /// If is null. /// If row count is less then column count - public DenseQR(DenseMatrix matrix, QRMethod method = QRMethod.Full) + public static DenseQR Create(DenseMatrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; - Tau = new Complex32[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + var tau = new Complex32[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + Matrix q; + Matrix r; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = new DenseMatrix(matrix.RowCount); - Control.LinearAlgebraProvider.QRFactor(((DenseMatrix)MatrixR).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)Q).Values, Tau); + r = matrix.Clone(); + q = new DenseMatrix(matrix.RowCount); + Control.LinearAlgebraProvider.QRFactor(((DenseMatrix) r).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) q).Values, tau); } else { - Q = matrix.Clone(); - MatrixR = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix)Q).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)MatrixR).Values, Tau); + q = matrix.Clone(); + r = new DenseMatrix(matrix.ColumnCount); + Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix) q).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) r).Values, tau); } + + return new DenseQR(q, r, method, tau); + } + + DenseQR(Matrix q, Matrix rFull, QRMethod method, Complex32[] tau) + : base(q, rFull, method) + { + Tau = tau; } /// @@ -102,17 +100,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -126,7 +113,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -143,7 +130,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do QR factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, Method); } /// @@ -153,16 +140,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries if (Q.RowCount != input.Count) @@ -171,9 +148,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var dinput = input as DenseVector; @@ -188,7 +165,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, Method); } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/GramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/GramSchmidt.cs index 33bcd896..9a9cea79 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/GramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/GramSchmidt.cs @@ -57,16 +57,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = Complex32.One; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0f)) + det *= FullR.At(i, i); + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0f)) { return 0; } @@ -84,9 +84,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0f)) + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0f)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/QR.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/QR.cs index fd99b31e..6690843e 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/QR.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/QR.cs @@ -3,7 +3,9 @@ // http://numerics.mathdotnet.com // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com -// Copyright (c) 2009-2010 Math.NET +// +// Copyright (c) 2009-2013 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 @@ -12,8 +14,10 @@ // 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 @@ -24,13 +28,13 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; +using MathNet.Numerics.Properties; namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { - using System; using Numerics; - using Properties; /// /// A class which encapsulates the functionality of the QR decomposition. @@ -46,6 +50,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// public abstract class QR : QR { + protected QR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { + } + /// /// Gets the absolute determinant value of the matrix for which the QR matrix was computed. /// @@ -53,16 +62,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = Complex32.One; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0f)) + det *= FullR.At(i, i); + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0f)) { return 0; } @@ -80,9 +89,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (MatrixR.At(i, i).Magnitude.AlmostEqual(0.0f)) + if (FullR.At(i, i).Magnitude.AlmostEqual(0.0f)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/UserGramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/UserGramSchmidt.cs index e9614d9b..601391a7 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/UserGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/UserGramSchmidt.cs @@ -154,19 +154,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j) / FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < input.ColumnCount; j++) { @@ -218,14 +218,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization // Solve R*X = Y; for (var k = Q.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k] * FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/UserQR.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/UserQR.cs index b0108220..b595989b 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/UserQR.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/UserQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,11 +28,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; +using System.Linq; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; using MathNet.Numerics.Threading; -using System; -using System.Linq; namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { @@ -47,7 +47,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class UserQR : QR + public sealed class UserQR : QR { /// /// Initializes a new instance of the class. This object will compute the @@ -56,71 +56,70 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The matrix to factor. /// The QR factorization method to use. /// If is null. - public UserQR(Matrix matrix, QRMethod method = QRMethod.Full) + public static UserQR Create(Matrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; + Matrix q; + Matrix r; + var minmn = Math.Min(matrix.RowCount, matrix.ColumnCount); var u = new Complex32[minmn][]; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); + r = matrix.Clone(); + q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); for (var i = 0; i < matrix.RowCount; i++) { - Q.At(i, i, 1.0f); + q.At(i, i, 1.0f); } for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(MatrixR, i, i); - ComputeQR(u[i], MatrixR, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(r, i, i); + ComputeQR(u[i], r, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.RowCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.RowCount, Control.NumberOfParallelWorkerThreads); } } else { - MatrixR = matrix.CreateMatrix(matrix.ColumnCount, matrix.ColumnCount); - Q = matrix.Clone(); + q = matrix.Clone(); for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(Q, i, i); - ComputeQR(u[i], Q, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(q, i, i); + ComputeQR(u[i], q, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } - MatrixR = Q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); - Q.Clear(); + r = q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); + q.Clear(); for (var i = 0; i < matrix.ColumnCount; i++) { - Q.At(i, i, 1.0f); + q.At(i, i, 1.0f); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } } + + return new UserQR(q, r, method); + } + + UserQR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { } /// @@ -130,7 +129,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The first row /// Column index /// Generated vector - private static Complex32[] GenerateColumn(Matrix a, int row, int column) + static Complex32[] GenerateColumn(Matrix a, int row, int column) { var ru = a.RowCount - row; var u = new Complex32[ru]; @@ -141,19 +140,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization a.At(i, column, 0.0f); } - var norm = u.Aggregate(Complex32.Zero, (current, t) => current + (t.Magnitude * t.Magnitude)); + var norm = u.Aggregate(Complex32.Zero, (current, t) => current + (t.Magnitude*t.Magnitude)); norm = norm.SquareRoot(); if (row == a.RowCount - 1 || norm.Magnitude == 0) { a.At(row, column, -u[0]); - u[0] = (float)Constants.Sqrt2; + u[0] = (float) Constants.Sqrt2; return u; } if (u[0].Magnitude != 0.0f) { - norm = norm.Magnitude * (u[0] / u[0].Magnitude); + norm = norm.Magnitude*(u[0]/u[0].Magnitude); } a.At(row, column, -norm); @@ -165,10 +164,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization u[0] += 1.0f; - var s = (1.0f / u[0]).SquareRoot(); + var s = (1.0f/u[0]).SquareRoot(); for (var i = 0; i < ru; i++) { - u[i] = u[i].Conjugate() * s; + u[i] = u[i].Conjugate()*s; } return u; @@ -184,7 +183,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The first column /// The last column /// Number of available CPUs - private static void ComputeQR(Complex32[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) + static void ComputeQR(Complex32[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) { if (rowDim < rowStart || columnDim < columnStart) { @@ -195,8 +194,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization if ((availableCores > 1) && (tmpColCount > 200)) { - var tmpSplit = columnStart + (tmpColCount / 2); - var tmpCores = availableCores / 2; + var tmpSplit = columnStart + (tmpColCount/2); + var tmpCores = availableCores/2; CommonParallel.Invoke( () => ComputeQR(u, a, rowStart, rowDim, columnStart, tmpSplit, tmpCores), @@ -209,12 +208,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization var scale = Complex32.Zero; for (var i = rowStart; i < rowDim; i++) { - scale += u[i - rowStart] * a.At(i, j); + scale += u[i - rowStart]*a.At(i, j); } for (var i = rowStart; i < rowDim; i++) { - a.At(i, j, a.At(i, j) - (u[i - rowStart].Conjugate() * scale)); + a.At(i, j, a.At(i, j) - (u[i - rowStart].Conjugate()*scale)); } } } @@ -227,17 +226,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -245,13 +233,13 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization } // The dimension compatibility conditions for X = A\B require the two matrices A and B to have the same number of rows - if (MatrixR.RowCount != input.RowCount) + if (FullR.RowCount != input.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameRowDimension); } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -259,20 +247,20 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new Complex32[MatrixR.RowCount]; + var column = new Complex32[FullR.RowCount]; for (var j = 0; j < input.ColumnCount; j++) { - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy.At(k, j); } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { var s = Complex32.Zero; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i).Conjugate() * column[k]; + s += Q.At(k, i).Conjugate()*column[k]; } inputCopy.At(i, j, s); @@ -280,23 +268,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j)/FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j)*FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < inputCopy.ColumnCount; j++) { @@ -312,60 +300,50 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries - if (MatrixR.RowCount != input.Count) + if (FullR.RowCount != input.Count) { throw new ArgumentException(Resources.ArgumentVectorsSameLength); } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new Complex32[MatrixR.RowCount]; - for (var k = 0; k < MatrixR.RowCount; k++) + var column = new Complex32[FullR.RowCount]; + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy[k]; } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { var s = Complex32.Zero; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i).Conjugate() * column[k]; + s += Q.At(k, i).Conjugate()*column[k]; } inputCopy[i] = s; } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k]*FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Complex32/Matrix.cs b/src/Numerics/LinearAlgebra/Complex32/Matrix.cs index 189f9aa3..2f39afcb 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Matrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Matrix.cs @@ -465,7 +465,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override QR QR(QRMethod method = QRMethod.Thin) { - return new UserQR(this, method); + return UserQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs index de183c98..85bfa7e5 100644 --- a/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs @@ -1047,7 +1047,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override QR QR(QRMethod method = QRMethod.Thin) { - return new DenseQR(this, method); + return DenseQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs index 70bf4cb8..46e3a430 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs @@ -157,7 +157,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); } /// @@ -192,7 +192,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); } } } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs index dcf3391b..a4ac0c35 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,9 +28,9 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; -using System; namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { @@ -43,17 +43,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class DenseQR : QR + public sealed class DenseQR : QR { /// /// Gets or sets Tau vector. Contains additional information on Q - used for native solver. /// - public double[] Tau - { - get; - set; - } - + double[] Tau { get; set; } + /// /// Initializes a new instance of the class. This object will compute the /// QR factorization when the constructor is called and cache it's factorization. @@ -62,36 +58,37 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The type of QR factorization to perform. /// If is null. /// If row count is less then column count - public DenseQR(DenseMatrix matrix, QRMethod method = QRMethod.Full) + public static DenseQR Create(DenseMatrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; - Tau = new double[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + var tau = new double[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + Matrix q; + Matrix r; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = new DenseMatrix(matrix.RowCount); - Control.LinearAlgebraProvider.QRFactor(((DenseMatrix)MatrixR).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)Q).Values, Tau); + r = matrix.Clone(); + q = new DenseMatrix(matrix.RowCount); + Control.LinearAlgebraProvider.QRFactor(((DenseMatrix) r).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) q).Values, tau); } else { - Q = matrix.Clone(); - MatrixR = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix) Q).Values, matrix.RowCount, - matrix.ColumnCount, - ((DenseMatrix) MatrixR).Values, Tau); + q = matrix.Clone(); + r = new DenseMatrix(matrix.ColumnCount); + Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix) q).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) r).Values, tau); } + + return new DenseQR(q, r, method, tau); + } + + DenseQR(Matrix q, Matrix rFull, QRMethod method, double[] tau) + : base(q, rFull, method) + { + Tau = tau; } /// @@ -101,17 +98,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -125,7 +111,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -142,7 +128,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do QR factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, Method); } /// @@ -152,16 +138,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries if (Q.RowCount != input.Count) @@ -170,9 +146,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var dinput = input as DenseVector; @@ -187,7 +163,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, Method); } } } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/GramSchmidt.cs b/src/Numerics/LinearAlgebra/Double/Factorization/GramSchmidt.cs index 1e1b34a2..f0572a36 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/GramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/GramSchmidt.cs @@ -55,16 +55,16 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = 1.0; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0)) + det *= FullR.At(i, i); + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0)) { return 0; } @@ -82,9 +82,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0)) + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/QR.cs b/src/Numerics/LinearAlgebra/Double/Factorization/QR.cs index 0db05ed0..0eda8b05 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/QR.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/QR.cs @@ -3,7 +3,9 @@ // http://numerics.mathdotnet.com // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com -// Copyright (c) 2009-2010 Math.NET +// +// Copyright (c) 2009-2013 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 @@ -12,8 +14,10 @@ // 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 @@ -24,13 +28,12 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; +using MathNet.Numerics.Properties; namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { - using System; - using Properties; - /// /// A class which encapsulates the functionality of the QR decomposition. /// Any real square matrix A (m x n) may be decomposed as A = QR where Q is an orthogonal matrix @@ -45,6 +48,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// public abstract class QR : QR { + protected QR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { + } + /// /// Gets the absolute determinant value of the matrix for which the QR matrix was computed. /// @@ -52,16 +60,16 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = 1.0; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0)) + det *= FullR.At(i, i); + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0)) { return 0; } @@ -79,9 +87,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0)) + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/UserGramSchmidt.cs b/src/Numerics/LinearAlgebra/Double/Factorization/UserGramSchmidt.cs index ba902841..dfba5f0f 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/UserGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/UserGramSchmidt.cs @@ -147,19 +147,19 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j) / FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < input.ColumnCount; j++) { @@ -211,14 +211,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization // Solve R*X = Y; for (var k = Q.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k] * FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/UserQR.cs b/src/Numerics/LinearAlgebra/Double/Factorization/UserQR.cs index 93433589..fe355c09 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/UserQR.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/UserQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,11 +28,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; +using System.Linq; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; using MathNet.Numerics.Threading; -using System; -using System.Linq; namespace MathNet.Numerics.LinearAlgebra.Double.Factorization { @@ -45,7 +45,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class UserQR : QR + public sealed class UserQR : QR { /// /// Initializes a new instance of the class. This object will compute the @@ -54,74 +54,72 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The matrix to factor. /// The QR factorization method to use. /// If is null. - public UserQR(Matrix matrix, QRMethod method = QRMethod.Full) + public static UserQR Create(Matrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; + Matrix q; + Matrix r; + var minmn = Math.Min(matrix.RowCount, matrix.ColumnCount); var u = new double[minmn][]; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); + r = matrix.Clone(); + q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); for (var i = 0; i < matrix.RowCount; i++) { - Q.At(i, i, 1.0); + q.At(i, i, 1.0); } for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(MatrixR, i, i); - ComputeQR(u[i], MatrixR, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(r, i, i); + ComputeQR(u[i], r, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.RowCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.RowCount, Control.NumberOfParallelWorkerThreads); } } else { - MatrixR = matrix.CreateMatrix(matrix.ColumnCount, matrix.ColumnCount); - Q = matrix.Clone(); + q = matrix.Clone(); for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(Q, i, i); - ComputeQR(u[i], Q, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(q, i, i); + ComputeQR(u[i], q, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } - MatrixR = Q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); - Q.Clear(); - + r = q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); + q.Clear(); + for (var i = 0; i < matrix.ColumnCount; i++) { - Q.At(i, i, 1.0); + q.At(i, i, 1.0); } - + for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } - } + + return new UserQR(q, r, method); + } + + UserQR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { } - + /// /// Generate column from initial matrix to work array /// @@ -129,7 +127,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The first row /// Column index /// Generated vector - private static double[] GenerateColumn(Matrix a, int row, int column) + static double[] GenerateColumn(Matrix a, int row, int column) { var ru = a.RowCount - row; var u = new double[ru]; @@ -140,7 +138,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization a.At(i, row, 0.0); } - var norm = u.Sum(t => t * t); + var norm = u.Sum(t => t*t); norm = Math.Sqrt(norm); if (row == a.RowCount - 1 || norm == 0) @@ -150,13 +148,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization return u; } - var scale = 1.0 / norm; + var scale = 1.0/norm; if (u[0] < 0.0) { scale *= -1.0; } - a.At(row, column, -1.0 / scale); + a.At(row, column, -1.0/scale); for (var i = 0; i < ru; i++) { @@ -164,7 +162,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } u[0] += 1.0; - var s = Math.Sqrt(1.0 / u[0]); + var s = Math.Sqrt(1.0/u[0]); for (var i = 0; i < ru; i++) { @@ -184,7 +182,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The first column /// The last column /// Number of available CPUs - private static void ComputeQR(double[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) + static void ComputeQR(double[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) { if (rowDim < rowStart || columnDim < columnStart) { @@ -195,8 +193,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization if ((availableCores > 1) && (tmpColCount > 200)) { - var tmpSplit = columnStart + (tmpColCount / 2); - var tmpCores = availableCores / 2; + var tmpSplit = columnStart + (tmpColCount/2); + var tmpCores = availableCores/2; CommonParallel.Invoke( () => ComputeQR(u, a, rowStart, rowDim, columnStart, tmpSplit, tmpCores), @@ -209,12 +207,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization var scale = 0.0; for (var i = rowStart; i < rowDim; i++) { - scale += u[i - rowStart] * a.At(i, j); + scale += u[i - rowStart]*a.At(i, j); } for (var i = rowStart; i < rowDim; i++) { - a.At(i, j, a.At(i, j) - (u[i - rowStart] * scale)); + a.At(i, j, a.At(i, j) - (u[i - rowStart]*scale)); } } } @@ -227,17 +225,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -245,13 +232,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } // The dimension compatibility conditions for X = A\B require the two matrices A and B to have the same number of rows - if (MatrixR.RowCount != input.RowCount) + if (FullR.RowCount != input.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameRowDimension); } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -259,20 +246,20 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new double[MatrixR.RowCount]; + var column = new double[FullR.RowCount]; for (var j = 0; j < input.ColumnCount; j++) { - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy.At(k, j); } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { double s = 0; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i) * column[k]; + s += Q.At(k, i)*column[k]; } inputCopy.At(i, j, s); @@ -280,23 +267,23 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j)/FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j)*FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < inputCopy.ColumnCount; j++) { @@ -312,60 +299,50 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries - if (MatrixR.RowCount != input.Count) + if (FullR.RowCount != input.Count) { throw new ArgumentException(Resources.ArgumentVectorsSameLength); } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new double[MatrixR.RowCount]; - for (var k = 0; k < MatrixR.RowCount; k++) + var column = new double[FullR.RowCount]; + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy[k]; } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { double s = 0; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i) * column[k]; + s += Q.At(k, i)*column[k]; } inputCopy[i] = s; } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k]*FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Double/Matrix.cs b/src/Numerics/LinearAlgebra/Double/Matrix.cs index 336850d2..bb0b4fdc 100644 --- a/src/Numerics/LinearAlgebra/Double/Matrix.cs +++ b/src/Numerics/LinearAlgebra/Double/Matrix.cs @@ -471,7 +471,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override QR QR(QRMethod method = QRMethod.Thin) { - return new UserQR(this, method); + return UserQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Factorization/GramSchmidt.cs b/src/Numerics/LinearAlgebra/Factorization/GramSchmidt.cs index 1cb59e1b..04e4fd92 100644 --- a/src/Numerics/LinearAlgebra/Factorization/GramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Factorization/GramSchmidt.cs @@ -44,9 +44,8 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization where T : struct, IEquatable, IFormattable { protected GramSchmidt(Matrix q, Matrix rFull) + : base(q, rFull, QRMethod.Full) { - Q = q; - MatrixR = rFull; } } } diff --git a/src/Numerics/LinearAlgebra/Factorization/QR.cs b/src/Numerics/LinearAlgebra/Factorization/QR.cs index d75fe570..6535b0c6 100644 --- a/src/Numerics/LinearAlgebra/Factorization/QR.cs +++ b/src/Numerics/LinearAlgebra/Factorization/QR.cs @@ -3,7 +3,9 @@ // http://numerics.mathdotnet.com // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com -// Copyright (c) 2009-2010 Math.NET +// +// Copyright (c) 2009-2013 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 @@ -12,8 +14,10 @@ // 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 @@ -58,29 +62,26 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization /// /// Supported data types are double, single, , and . public abstract class QR : ISolver - where T : struct, IEquatable, IFormattable + where T : struct, IEquatable, IFormattable { readonly Lazy> _lazyR; - protected QR() + protected readonly Matrix FullR; + protected readonly QRMethod Method; + + protected QR(Matrix q, Matrix rFull, QRMethod method) { - _lazyR = new Lazy>(ComputeR); + Q = q; + FullR = rFull; + Method = method; + + _lazyR = new Lazy>(FullR.UpperTriangle); } /// /// Gets or sets orthogonal Q matrix /// - public Matrix Q { get; protected set; } - - /// - /// Gets or sets upper triangular factor R - /// - protected Matrix MatrixR { get; set; } - - /// - /// The QR factorization method. - /// - protected QRMethod QrMethod { get; set; } + public Matrix Q { get; private set; } /// /// Gets the upper triangular factor R. @@ -101,11 +102,6 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization /// true if the matrix is full rank; otherwise false. public abstract bool IsFullRank { get; } - private Matrix ComputeR() - { - return MatrixR.UpperTriangle(); - } - /// /// Solves a system of linear equations, AX = B, with A QR factorized. /// @@ -113,13 +109,7 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization /// The left hand side , X. public virtual Matrix Solve(Matrix input) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - var matrixX = input.CreateMatrix(MatrixR.ColumnCount, input.ColumnCount); + var matrixX = input.CreateMatrix(FullR.ColumnCount, input.ColumnCount); Solve(input, matrixX); return matrixX; } @@ -138,13 +128,7 @@ namespace MathNet.Numerics.LinearAlgebra.Factorization /// The left hand side , x. public virtual Vector Solve(Vector input) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - var x = input.CreateVector(MatrixR.ColumnCount); + var x = input.CreateVector(FullR.ColumnCount); Solve(input, x); return x; } diff --git a/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs index 42865d30..6d9ae7dd 100644 --- a/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs @@ -1047,7 +1047,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override QR QR(QRMethod method = QRMethod.Thin) { - return new DenseQR(this, method); + return DenseQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs index c5b803ff..b44c43db 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs @@ -157,7 +157,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, input.ColumnCount, dresult.Values, QRMethod.Thin); } /// @@ -192,7 +192,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)FullR).Values, Q.RowCount, FullR.ColumnCount, null, dinput.Values, 1, dresult.Values, QRMethod.Thin); } } } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs index 06800375..fc58e8c8 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,9 +28,9 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; -using System; namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { @@ -43,16 +43,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class DenseQR : QR + public sealed class DenseQR : QR { /// /// Gets or sets Tau vector. Contains additional information on Q - used for native solver. /// - internal float[] Tau - { - get; - set; - } + float[] Tau { get; set; } /// /// Initializes a new instance of the class. This object will compute the @@ -62,35 +58,37 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The QR factorization method to use. /// If is null. /// If row count is less then column count - public DenseQR(DenseMatrix matrix, QRMethod method = QRMethod.Full) + public static DenseQR Create(DenseMatrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; - Tau = new float[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + var Tau = new float[Math.Min(matrix.RowCount, matrix.ColumnCount)]; + Matrix Q; + Matrix MatrixR; if (method == QRMethod.Full) { MatrixR = matrix.Clone(); Q = new DenseMatrix(matrix.RowCount); - Control.LinearAlgebraProvider.QRFactor(((DenseMatrix)MatrixR).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)Q).Values, Tau); + Control.LinearAlgebraProvider.QRFactor(((DenseMatrix) MatrixR).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) Q).Values, Tau); } else { Q = matrix.Clone(); MatrixR = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix)Q).Values, matrix.RowCount, matrix.ColumnCount, - ((DenseMatrix)MatrixR).Values, Tau); + Control.LinearAlgebraProvider.ThinQRFactor(((DenseMatrix) Q).Values, matrix.RowCount, matrix.ColumnCount, ((DenseMatrix) MatrixR).Values, Tau); } + + return new DenseQR(Q, MatrixR, method, Tau); + } + + DenseQR(Matrix q, Matrix rFull, QRMethod method, float[] tau) + : base(q, rFull, method) + { + Tau = tau; } /// @@ -100,17 +98,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -124,7 +111,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -141,7 +128,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do QR factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, input.ColumnCount, dresult.Values, Method); } /// @@ -151,16 +138,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries if (Q.RowCount != input.Count) @@ -169,9 +146,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var dinput = input as DenseVector; @@ -186,7 +163,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)Q).Values, ((DenseMatrix)MatrixR).Values, Q.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, QrMethod); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix) Q).Values, ((DenseMatrix) FullR).Values, Q.RowCount, FullR.ColumnCount, Tau, dinput.Values, 1, dresult.Values, Method); } } } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/GramSchmidt.cs b/src/Numerics/LinearAlgebra/Single/Factorization/GramSchmidt.cs index f58ae735..2978389a 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/GramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/GramSchmidt.cs @@ -55,16 +55,16 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = 1.0; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0f)) + det *= FullR.At(i, i); + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0f)) { return 0; } @@ -82,9 +82,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0f)) + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0f)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/QR.cs b/src/Numerics/LinearAlgebra/Single/Factorization/QR.cs index dfdf6b5e..fff91b25 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/QR.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/QR.cs @@ -3,7 +3,9 @@ // http://numerics.mathdotnet.com // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com -// Copyright (c) 2009-2010 Math.NET +// +// Copyright (c) 2009-2013 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 @@ -12,8 +14,10 @@ // 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 @@ -24,13 +28,12 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; +using MathNet.Numerics.Properties; namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { - using System; - using Properties; - /// /// A class which encapsulates the functionality of the QR decomposition. /// Any real square matrix A (m x n) may be decomposed as A = QR where Q is an orthogonal matrix @@ -45,6 +48,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// public abstract class QR : QR { + protected QR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { + } + /// /// Gets the absolute determinant value of the matrix for which the QR matrix was computed. /// @@ -52,16 +60,16 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { get { - if (MatrixR.RowCount != MatrixR.ColumnCount) + if (FullR.RowCount != FullR.ColumnCount) { throw new ArgumentException(Resources.ArgumentMatrixSquare); } var det = 1.0; - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - det *= MatrixR.At(i, i); - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0f)) + det *= FullR.At(i, i); + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0f)) { return 0; } @@ -79,9 +87,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { get { - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { - if (Math.Abs(MatrixR.At(i, i)).AlmostEqual(0.0f)) + if (Math.Abs(FullR.At(i, i)).AlmostEqual(0.0f)) { return false; } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/UserGramSchmidt.cs b/src/Numerics/LinearAlgebra/Single/Factorization/UserGramSchmidt.cs index bf983bdd..d3c62a1d 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/UserGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/UserGramSchmidt.cs @@ -147,19 +147,19 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j) / FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < input.ColumnCount; j++) { @@ -211,14 +211,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization // Solve R*X = Y; for (var k = Q.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k] * FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/UserQR.cs b/src/Numerics/LinearAlgebra/Single/Factorization/UserQR.cs index 65e81d95..74e7b64f 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/UserQR.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/UserQR.cs @@ -4,7 +4,7 @@ // http://github.com/mathnet/mathnet-numerics // http://mathnetnumerics.codeplex.com // -// Copyright (c) 2009-2010 Math.NET +// Copyright (c) 2009-2013 Math.NET // // Permission is hereby granted, free of charge, to any person // obtaining a copy of this software and associated documentation @@ -28,11 +28,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; +using System.Linq; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.Properties; using MathNet.Numerics.Threading; -using System; -using System.Linq; namespace MathNet.Numerics.LinearAlgebra.Single.Factorization { @@ -45,7 +45,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// /// The computation of the QR decomposition is done at construction time by Householder transformation. /// - public class UserQR : QR + public sealed class UserQR : QR { /// /// Initializes a new instance of the class. This object will compute the @@ -54,71 +54,70 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The matrix to factor. /// The QR factorization method to use. /// If is null. - public UserQR(Matrix matrix, QRMethod method = QRMethod.Full) + public static UserQR Create(Matrix matrix, QRMethod method = QRMethod.Full) { - if (matrix == null) - { - throw new ArgumentNullException("matrix"); - } - if (matrix.RowCount < matrix.ColumnCount) { throw Matrix.DimensionsDontMatch(matrix); } - QrMethod = method; + Matrix q; + Matrix r; + var minmn = Math.Min(matrix.RowCount, matrix.ColumnCount); var u = new float[minmn][]; if (method == QRMethod.Full) { - MatrixR = matrix.Clone(); - Q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); + r = matrix.Clone(); + q = matrix.CreateMatrix(matrix.RowCount, matrix.RowCount); for (var i = 0; i < matrix.RowCount; i++) { - Q.At(i, i, 1.0f); + q.At(i, i, 1.0f); } for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(MatrixR, i, i); - ComputeQR(u[i], MatrixR, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(r, i, i); + ComputeQR(u[i], r, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.RowCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.RowCount, Control.NumberOfParallelWorkerThreads); } } else { - MatrixR = matrix.CreateMatrix(matrix.ColumnCount, matrix.ColumnCount); - Q = matrix.Clone(); + q = matrix.Clone(); for (var i = 0; i < minmn; i++) { - u[i] = GenerateColumn(Q, i, i); - ComputeQR(u[i], Q, i, matrix.RowCount, i + 1, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + u[i] = GenerateColumn(q, i, i); + ComputeQR(u[i], q, i, matrix.RowCount, i + 1, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } - MatrixR = Q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); - Q.Clear(); + r = q.SubMatrix(0, matrix.ColumnCount, 0, matrix.ColumnCount); + q.Clear(); for (var i = 0; i < matrix.ColumnCount; i++) { - Q.At(i, i, 1.0f); + q.At(i, i, 1.0f); } for (var i = minmn - 1; i >= 0; i--) { - ComputeQR(u[i], Q, i, matrix.RowCount, i, matrix.ColumnCount, - Control.NumberOfParallelWorkerThreads); + ComputeQR(u[i], q, i, matrix.RowCount, i, matrix.ColumnCount, Control.NumberOfParallelWorkerThreads); } } + + return new UserQR(q, r, method); + } + + UserQR(Matrix q, Matrix rFull, QRMethod method) + : base(q, rFull, method) + { } /// @@ -128,7 +127,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The first row /// Column index /// Generated vector - private static float[] GenerateColumn(Matrix a, int row, int column) + static float[] GenerateColumn(Matrix a, int row, int column) { var ru = a.RowCount - row; var u = new float[ru]; @@ -139,23 +138,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization a.At(i, row, 0.0f); } - var norm = u.Sum(t => t * t); - norm = (float)Math.Sqrt(norm); + var norm = u.Sum(t => t*t); + norm = (float) Math.Sqrt(norm); if (row == a.RowCount - 1 || norm == 0) { a.At(row, column, -u[0]); - u[0] = (float)Constants.Sqrt2; + u[0] = (float) Constants.Sqrt2; return u; } - var scale = 1.0f / norm; + var scale = 1.0f/norm; if (u[0] < 0.0) { scale *= -1.0f; } - a.At(row, column, -1.0f / scale); + a.At(row, column, -1.0f/scale); for (var i = 0; i < ru; i++) { @@ -163,7 +162,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } u[0] += 1.0f; - var s = (float)Math.Sqrt(1.0 / u[0]); + var s = (float) Math.Sqrt(1.0/u[0]); for (var i = 0; i < ru; i++) { @@ -183,7 +182,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The first column /// The last column /// Number of available CPUs - private static void ComputeQR(float[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) + static void ComputeQR(float[] u, Matrix a, int rowStart, int rowDim, int columnStart, int columnDim, int availableCores) { if (rowDim < rowStart || columnDim < columnStart) { @@ -194,8 +193,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization if ((availableCores > 1) && (tmpColCount > 200)) { - var tmpSplit = columnStart + (tmpColCount / 2); - var tmpCores = availableCores / 2; + var tmpSplit = columnStart + (tmpColCount/2); + var tmpCores = availableCores/2; CommonParallel.Invoke( () => ComputeQR(u, a, rowStart, rowDim, columnStart, tmpSplit, tmpCores), @@ -208,12 +207,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization var scale = 0.0f; for (var i = rowStart; i < rowDim; i++) { - scale += u[i - rowStart] * a.At(i, j); + scale += u[i - rowStart]*a.At(i, j); } for (var i = rowStart; i < rowDim; i++) { - a.At(i, j, a.At(i, j) - (u[i - rowStart] * scale)); + a.At(i, j, a.At(i, j) - (u[i - rowStart]*scale)); } } } @@ -226,17 +225,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The left hand side , X. public override void Solve(Matrix input, Matrix result) { - // Check for proper arguments. - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // The solution X should have the same number of columns as B if (input.ColumnCount != result.ColumnCount) { @@ -244,13 +232,13 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } // The dimension compatibility conditions for X = A\B require the two matrices A and B to have the same number of rows - if (MatrixR.RowCount != input.RowCount) + if (FullR.RowCount != input.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameRowDimension); } // The solution X row dimension is equal to the column dimension of A - if (MatrixR.ColumnCount != result.RowCount) + if (FullR.ColumnCount != result.RowCount) { throw new ArgumentException(Resources.ArgumentMatrixSameColumnDimension); } @@ -258,20 +246,20 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new float[MatrixR.RowCount]; + var column = new float[FullR.RowCount]; for (var j = 0; j < input.ColumnCount; j++) { - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy.At(k, j); } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { float s = 0; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i) * column[k]; + s += Q.At(k, i)*column[k]; } inputCopy.At(i, j, s); @@ -279,23 +267,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(k, j, inputCopy.At(k, j) / MatrixR.At(k, k)); + inputCopy.At(k, j, inputCopy.At(k, j)/FullR.At(k, k)); } for (var i = 0; i < k; i++) { for (var j = 0; j < input.ColumnCount; j++) { - inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j) * MatrixR.At(i, k))); + inputCopy.At(i, j, inputCopy.At(i, j) - (inputCopy.At(k, j)*FullR.At(i, k))); } } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { for (var j = 0; j < inputCopy.ColumnCount; j++) { @@ -311,60 +299,50 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization /// The left hand side , x. public override void Solve(Vector input, Vector result) { - if (input == null) - { - throw new ArgumentNullException("input"); - } - - if (result == null) - { - throw new ArgumentNullException("result"); - } - // Ax=b where A is an m x n matrix // Check that b is a column vector with m entries - if (MatrixR.RowCount != input.Count) + if (FullR.RowCount != input.Count) { throw new ArgumentException(Resources.ArgumentVectorsSameLength); } // Check that x is a column vector with n entries - if (MatrixR.ColumnCount != result.Count) + if (FullR.ColumnCount != result.Count) { - throw Matrix.DimensionsDontMatch(MatrixR, result); + throw Matrix.DimensionsDontMatch(FullR, result); } var inputCopy = input.Clone(); // Compute Y = transpose(Q)*B - var column = new float[MatrixR.RowCount]; - for (var k = 0; k < MatrixR.RowCount; k++) + var column = new float[FullR.RowCount]; + for (var k = 0; k < FullR.RowCount; k++) { column[k] = inputCopy[k]; } - for (var i = 0; i < MatrixR.RowCount; i++) + for (var i = 0; i < FullR.RowCount; i++) { float s = 0; - for (var k = 0; k < MatrixR.RowCount; k++) + for (var k = 0; k < FullR.RowCount; k++) { - s += Q.At(k, i) * column[k]; + s += Q.At(k, i)*column[k]; } inputCopy[i] = s; } // Solve R*X = Y; - for (var k = MatrixR.ColumnCount - 1; k >= 0; k--) + for (var k = FullR.ColumnCount - 1; k >= 0; k--) { - inputCopy[k] /= MatrixR.At(k, k); + inputCopy[k] /= FullR.At(k, k); for (var i = 0; i < k; i++) { - inputCopy[i] -= inputCopy[k] * MatrixR.At(i, k); + inputCopy[i] -= inputCopy[k]*FullR.At(i, k); } } - for (var i = 0; i < MatrixR.ColumnCount; i++) + for (var i = 0; i < FullR.ColumnCount; i++) { result[i] = inputCopy[i]; } diff --git a/src/Numerics/LinearAlgebra/Single/Matrix.cs b/src/Numerics/LinearAlgebra/Single/Matrix.cs index 759605f2..cce9f85d 100644 --- a/src/Numerics/LinearAlgebra/Single/Matrix.cs +++ b/src/Numerics/LinearAlgebra/Single/Matrix.cs @@ -471,7 +471,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override QR QR(QRMethod method = QRMethod.Thin) { - return new UserQR(this, method); + return UserQR.Create(this, method); } public override GramSchmidt GramSchmidt() diff --git a/src/UnitTests/LinearAlgebraTests/Complex/Factorization/QRTests.cs b/src/UnitTests/LinearAlgebraTests/Complex/Factorization/QRTests.cs index 1019b5fc..fde3ea30 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex/Factorization/QRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex/Factorization/QRTests.cs @@ -24,11 +24,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Complex; using MathNet.Numerics.LinearAlgebra.Complex.Factorization; using MathNet.Numerics.LinearAlgebra.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex.Factorization { @@ -39,22 +39,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex.Factorization /// public class QRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new DenseQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new DenseQR(new DenseMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new DenseMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex/Factorization/UserQRTests.cs b/src/UnitTests/LinearAlgebraTests/Complex/Factorization/UserQRTests.cs index 691f1ba4..281ac205 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex/Factorization/UserQRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex/Factorization/UserQRTests.cs @@ -24,10 +24,10 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Complex.Factorization; using MathNet.Numerics.LinearAlgebra.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex.Factorization { @@ -38,22 +38,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex.Factorization /// public class UserQRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new UserQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new UserQR(new UserDefinedMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new UserDefinedMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/QRTests.cs b/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/QRTests.cs index 991aae05..26623a4b 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/QRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/QRTests.cs @@ -24,11 +24,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Complex32; using MathNet.Numerics.LinearAlgebra.Complex32.Factorization; using MathNet.Numerics.LinearAlgebra.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Factorization { @@ -39,22 +39,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Factorization /// public class QRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new DenseQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new DenseQR(new DenseMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new DenseMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/UserQRTests.cs b/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/UserQRTests.cs index c6459845..3efce8a2 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/UserQRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex32/Factorization/UserQRTests.cs @@ -24,10 +24,10 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Complex32.Factorization; using MathNet.Numerics.LinearAlgebra.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Factorization { @@ -38,22 +38,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32.Factorization /// public class UserQRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new UserQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new UserQR(new UserDefinedMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new UserDefinedMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Double/Factorization/QRTests.cs b/src/UnitTests/LinearAlgebraTests/Double/Factorization/QRTests.cs index 80afd8ea..1d46f336 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/Factorization/QRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/Factorization/QRTests.cs @@ -24,11 +24,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Double; using MathNet.Numerics.LinearAlgebra.Double.Factorization; using MathNet.Numerics.LinearAlgebra.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double.Factorization { @@ -37,22 +37,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double.Factorization /// public class QRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new DenseQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new DenseQR(new DenseMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new DenseMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Double/Factorization/UserQRTests.cs b/src/UnitTests/LinearAlgebraTests/Double/Factorization/UserQRTests.cs index 8cf820bb..ac3210b2 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/Factorization/UserQRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/Factorization/UserQRTests.cs @@ -24,10 +24,10 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Double.Factorization; using MathNet.Numerics.LinearAlgebra.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double.Factorization { @@ -36,22 +36,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double.Factorization /// public class UserQRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new UserQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new UserQR(new UserDefinedMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new UserDefinedMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Single/Factorization/QRTests.cs b/src/UnitTests/LinearAlgebraTests/Single/Factorization/QRTests.cs index e5c823ba..e876373f 100644 --- a/src/UnitTests/LinearAlgebraTests/Single/Factorization/QRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Single/Factorization/QRTests.cs @@ -24,11 +24,11 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.LinearAlgebra.Single; using MathNet.Numerics.LinearAlgebra.Single.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Factorization { @@ -37,22 +37,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Factorization /// public class QRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new DenseQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new DenseQR(new DenseMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new DenseMatrix(3, 4))); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Single/Factorization/UserQRTests.cs b/src/UnitTests/LinearAlgebraTests/Single/Factorization/UserQRTests.cs index f4d034ff..ccdede6e 100644 --- a/src/UnitTests/LinearAlgebraTests/Single/Factorization/UserQRTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Single/Factorization/UserQRTests.cs @@ -24,10 +24,10 @@ // OTHER DEALINGS IN THE SOFTWARE. // +using System; using MathNet.Numerics.LinearAlgebra.Factorization; using MathNet.Numerics.LinearAlgebra.Single.Factorization; using NUnit.Framework; -using System; namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Factorization { @@ -36,22 +36,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single.Factorization /// public class UserQRTests { - /// - /// Constructor with null throws ArgumentNullException. - /// - [Test] - public void ConstructorNull() - { - Assert.Throws(() => new UserQR(null)); - } - /// /// Constructor with wide matrix throws ArgumentException. /// [Test] public void ConstructorWideMatrixThrowsInvalidMatrixOperationException() { - Assert.Throws(() => new UserQR(new UserDefinedMatrix(3, 4))); + Assert.Throws(() => UserQR.Create(new UserDefinedMatrix(3, 4))); } ///