diff --git a/src/Numerics/LinearAlgebra/Complex/SquareMatrix.cs b/src/Numerics/LinearAlgebra/Complex/SquareMatrix.cs new file mode 100644 index 00000000..2e95ea27 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Complex/SquareMatrix.cs @@ -0,0 +1,62 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Complex +{ + using System; + using System.Numerics; + using Properties; + using Storage; + + /// + /// Abstract class for square matrices. + /// + [Serializable] + public abstract class SquareMatrix : Matrix + { + /// + /// Number of rows or columns. + /// + protected readonly int Order; + + /// + /// Initializes a new instance of the class. + /// + /// + /// If the matrix is not square. + /// + protected SquareMatrix(MatrixStorage storage) + : base(storage) + { + if (storage.RowCount != storage.ColumnCount) + { + throw new ArgumentException(Resources.ArgumentMatrixSquare); + } + + Order = storage.RowCount; + } + } +} diff --git a/src/Numerics/LinearAlgebra/Complex/SymmetricDenseMatrix.cs b/src/Numerics/LinearAlgebra/Complex/SymmetricDenseMatrix.cs new file mode 100644 index 00000000..df3305a0 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Complex/SymmetricDenseMatrix.cs @@ -0,0 +1,767 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Complex +{ + using System; + using System.Numerics; + using Generic; + using MathNet.Numerics.Distributions; + using MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static; + using Properties; + using Storage; + + /// + /// A Symmetric Matrix class with dense storage. + /// + /// The underlying storage is a one dimensional array in column-major order. + /// The Upper Triangle is stored(it is equal to the Lower Triangle) + [Serializable] + public class SymmetricDenseMatrix : SymmetricMatrix + { + readonly DenseColumnMajorSymmetricMatrixStorage _storage; + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + readonly Complex[] _data; + + internal SymmetricDenseMatrix(DenseColumnMajorSymmetricMatrixStorage storage) + : base(storage) + { + _storage = storage; + _data = _storage.Data; + } + + /// + /// Initializes a new instance of the class. This matrix is square with a given size. + /// + /// The order of the matrix. + /// + /// If is less than one. + /// + public SymmetricDenseMatrix(int order) + : this(new DenseColumnMajorSymmetricMatrixStorage(order)) + { + } + + /// + /// Initializes a new instance of the class with all entries set to a particular value. + /// + /// + /// The order of the matrix. + /// + /// The value which we assign to each element of the matrix. + public SymmetricDenseMatrix(int order, Complex value) + : this(order) + { + for (var i = 0; i < _data.Length; i++) + { + _data[i] = value; + } + } + + /// + /// Initializes a new instance of the class from a one dimensional array. This constructor + /// will reference the one dimensional array and not copy it. + /// + /// The size of the square matrix. + /// + /// The one dimensional array to create this matrix from. Column-major and row-major order is identical on a symmetric matrix: http://en.wikipedia.org/wiki/Row-major_order + /// + /// + /// If does not represent a packed array. + /// + public SymmetricDenseMatrix(int order, Complex[] array) + : this(new DenseColumnMajorSymmetricMatrixStorage(order, array)) + { + } + + /// + /// Initializes a new instance of the class from a 2D array. This constructor + /// will allocate a completely new memory block for storing the symmetric dense matrix. + /// + /// The 2D array to create this matrix from. + /// + /// If is not a square array. + /// + /// + /// If is not a symmetric array. + /// + public SymmetricDenseMatrix(Complex[,] array) + : this(array.GetLength(0)) + { + if (!CheckIfSymmetric(array)) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = array[row, column]; + } + } + } + + /// + /// Initializes a new instance of the class, copying + /// the values from the given matrix. Matrix must be Symmetric. + /// + /// The matrix to copy. + /// + /// If is not a square matrix. + /// + /// + /// If is not a symmetric matrix. + /// + public SymmetricDenseMatrix(Matrix matrix) + : this(matrix.RowCount) + { + var symmetricMatrix = matrix as SymmetricDenseMatrix; + + if (!matrix.IsSymmetric) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + if (symmetricMatrix == null) + { + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = matrix[row, column]; + } + } + } + else + { + matrix.CopyTo(this); + } + } + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + public Complex[] Data + { + get { return _data; } + } + + /// + /// Creates a SymmetricDenseMatrix for the given number of rows and columns. + /// If rows and columns are not equal, returns a DenseMatrix instead. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// True if all fields must be mutable (e.g. not a diagonal matrix). + /// + /// A DenseMatrix or SymmetricDenseMatrix with the given dimensions. + /// + /// /// + /// If is not equal to . + /// Symmetric arrays are always square + /// + public override Matrix CreateMatrix(int numberOfRows, int numberOfColumns, bool fullyMutable = false) + { + if (numberOfRows != numberOfColumns || fullyMutable) + { + return new DenseMatrix(numberOfRows, numberOfColumns); + } + + return new SymmetricDenseMatrix(numberOfRows); + } + + /// + /// Creates a with a the given dimension. + /// + /// The size of the vector. + /// True if all fields must be mutable. + /// + /// A with the given dimension. + /// + public override Vector CreateVector(int size, bool fullyMutable = false) + { + return new DenseVector(size); + } + + #region Static constructors for special matrices. + + /// + /// Initializes a square with all zero's except for ones on the diagonal. + /// + /// the size of the square matrix. + /// A symmetric dense identity matrix. + /// + /// If is less than one. + /// + public static SymmetricDenseMatrix Identity(int order) + { + var m = new SymmetricDenseMatrix(order); + for (var i = 0; i < order; i++) + { + m.At(i, i, Complex.One); + } + + return m; + } + + #endregion + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The matrix to store the result of add + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + protected override void DoAdd(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoAdd(other, result); + } + else + { + Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The matrix to store the result of the subtraction. + protected override void DoSubtract(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoSubtract(other, result); + } + else + { + Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// The scalar to multiply the matrix with. + /// The matrix to store the result of the multiplication. + protected override void DoMultiply(Complex scalar, Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + if (denseResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + } + } + + /// + /// Multiplies this matrix with a vector and places the results into the result vector. + /// + /// The vector to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Vector rightSide, Vector result) + { + var denseRight = rightSide as DenseVector; + var denseResult = result as DenseVector; + + if (denseRight == null || denseResult == null) + { + base.DoMultiply(rightSide, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(rightSide, result); + } + } + + /// + /// Multiplies this matrix with another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(other, result); + } + } + + /// + /// Multiplies this matrix with transpose of another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoTransposeAndMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoTransposeAndMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoTransposeAndMultiply(other, result); + } + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// The result of the negation. + protected override void DoNegate(Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + + if (denseResult == null) + { + base.DoNegate(result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + } + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise multiply with this one. + /// The matrix to store the result of the pointwise multiplication. + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise divide this one by. + /// The matrix to store the result of the pointwise division. + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix LowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column <= row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix StrictlyLowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column < row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix UpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix StrictlyUpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row + 1; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Computes the trace of this matrix. + /// + /// The trace of this matrix + public override Complex Trace() + { + // Matrix is always square. + var sum = Complex.Zero; + for (var i = 0; i < RowCount; i++) + { + sum += At(i, i); + } + + return sum; + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = distribution.Sample(); + } + } + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = distribution.Sample(); + } + } + } + + /// + /// Adds two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to add. + /// The right matrix to add. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Add(rightSide); + } + + /// + /// Returns a Matrix containing the same values of . + /// + /// The matrix to get the values from. + /// A matrix containing a the same values as . + /// If is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Clone(); + } + + /// + /// Subtracts two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to subtract. + /// The right matrix to subtract. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Subtract(rightSide); + } + + /// + /// Negates each element of the matrix. + /// + /// The matrix to negate. + /// A matrix containing the negated values. + /// If is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Negate(); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, Complex rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(Complex leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Multiply(leftSide); + } + + /// + /// Multiplies two matrices. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to multiply. + /// The right matrix to multiply. + /// The result of multiplication. + /// If or is . + /// If the dimensions of or don't conform. + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide.ColumnCount != rightSide.RowCount) + { + throw new ArgumentException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix and a Vector. + /// + /// The matrix to multiply. + /// The vector to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(SymmetricDenseMatrix leftSide, DenseVector rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (DenseVector)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Vector and a Matrix. + /// + /// The vector to multiply. + /// The matrix to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(DenseVector leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (DenseVector)rightSide.LeftMultiply(leftSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator %(SymmetricDenseMatrix leftSide, Complex rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Modulus(rightSide); + } + } +} \ No newline at end of file diff --git a/src/Numerics/LinearAlgebra/Complex/SymmetricMatrix.cs b/src/Numerics/LinearAlgebra/Complex/SymmetricMatrix.cs new file mode 100644 index 00000000..6f1b5975 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Complex/SymmetricMatrix.cs @@ -0,0 +1,657 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Complex +{ + using System; + using System.Numerics; + using Generic; + using Distributions; + using Properties; + using Storage; + + /// + /// Abstract class for symmetric matrices. + /// + [Serializable] + public abstract class SymmetricMatrix : SquareMatrix + { + /// + /// Initializes a new instance of the class. + /// + protected SymmetricMatrix(MatrixStorage storage) + : base(storage) + { + } + + /// + /// Returns a value indicating whether the array is symmetric. + /// + /// + /// The array to check for symmetry. + /// + /// + /// True is array is symmetric, false if not symmetric. + /// + public static bool CheckIfSymmetric(Complex[,] array) + { + var rows = array.GetLength(0); + var columns = array.GetLength(1); + + if (rows != columns) + { + return false; + } + + for (var row = 0; row < rows; row++) + { + for (var column = 0; column < columns; column++) + { + if (column >= row) + { + continue; + } + + if (!array[row, column].Equals(array[column, row])) + { + return false; + } + } + } + + return true; + } + + /// + /// Gets a value indicating whether this matrix is symmetric. + /// + public override sealed bool IsSymmetric + { + get + { + return true; + } + } + + /// + /// Returns the transpose of this matrix. The transpose is equal and this method returns a reference to this matrix. + /// + /// + /// The transpose of this matrix. + /// + public override sealed Matrix Transpose() + { + return this.Clone(); + } + + /// + /// Returns the conjugate transpose of this matrix. + /// + /// The conjugate transpose of this matrix. + public override Matrix ConjugateTranspose() + { + var ret = CreateMatrix(ColumnCount, RowCount); + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + ret.At(row, column, At(column, row).Conjugate()); + } + } + + return ret; + } + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The result of the addition. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Add(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoAdd(other, result); + return result; + } + + /// + /// Adds another matrix to this matrix. + /// + /// + /// The matrix to add to this matrix. + /// + /// + /// The matrix to store the result of the addition. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoAdd(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric + non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoAdd(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) + symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The result of the subtraction. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Subtract(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoSubtract(other, result); + return result; + } + + + /// + /// Subtracts another matrix from this matrix. + /// + /// + /// The matrix to subtract to this matrix. + /// + /// + /// The matrix to store the result of subtraction. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoSubtract(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric - non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoSubtract(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) - symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// + /// The scalar to multiply the matrix with. + /// + /// + /// The matrix to store the result of the multiplication. + /// + protected override void DoMultiply(Complex scalar, Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * scalar); + } + } + } + } + + /// + /// Multiplies the transpose of this matrix with another matrix and places the results into the result matrix. + /// + /// + /// The matrix to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Matrix other, Matrix result) + { + DoMultiply(other, result); + } + + /// + /// Multiplies the transpose of this matrix with a vector and places the results into the result vector. + /// + /// + /// The vector to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Vector rightSide, Vector result) + { + DoMultiply(rightSide, result); + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// + /// The result of the negation. + /// + protected override void DoNegate(Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoNegate(result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column != ColumnCount; column++) + { + symmetricResult[row, column] = -At(row, column); + } + } + } + } + + /// + /// Pointwise multiplies this matrix with another matrix. + /// + /// The matrix to pointwise multiply with this one. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise multiplication of this matrix and . + public override Matrix PointwiseMultiply(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other, "other"); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseMultiply(other, result); + return result; + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise multiply with this one. + /// + /// + /// The matrix to store the result of the pointwise multiplication. + /// + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise* non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Pointwise divide this matrix by another matrix. + /// + /// The matrix to pointwise subtract this one by. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise division of this matrix and . + public override Matrix PointwiseDivide(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseDivide(other, result); + return result; + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise divide this one by. + /// + /// + /// The matrix to store the result of the pointwise division. + /// + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise/ non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) / symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Computes the modulus for each element of the matrix. + /// + /// + /// The divisor to use. + /// + /// + /// Matrix to store the results in. + /// + protected override void DoModulus(Complex divisor, Matrix result) + { + throw new NotImplementedException(); + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, distribution.Sample()); + } + } + } + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, distribution.Sample()); + } + } + } + } + + /// + /// Creates a new matrix and inserts the given column at the given index. + /// + /// The index of where to insert the column. + /// The column to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of columns. + /// If the size of != the number of rows. + public override Matrix InsertColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Inserting a column is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given array to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, Complex[] column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given Vector to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Creates a new matrix and inserts the given row at the given index. + /// + /// The index of where to insert the row. + /// The row to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of rows. + /// If the size of != the number of columns. + public override Matrix InsertRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Inserting a row is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given Vector to the specified row. + /// + /// The row to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given array to the specified row. + /// + /// The row to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, Complex[] row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + } +} diff --git a/src/Numerics/LinearAlgebra/Complex32/SquareMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/SquareMatrix.cs new file mode 100644 index 00000000..337a5118 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Complex32/SquareMatrix.cs @@ -0,0 +1,62 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Complex32 +{ + using System; + using Numerics; + using Properties; + using Storage; + + /// + /// Abstract class for square matrices. + /// + [Serializable] + public abstract class SquareMatrix : Matrix + { + /// + /// Number of rows or columns. + /// + protected readonly int Order; + + /// + /// Initializes a new instance of the class. + /// + /// + /// If the matrix is not square. + /// + protected SquareMatrix(MatrixStorage storage) + : base(storage) + { + if (storage.RowCount != storage.ColumnCount) + { + throw new ArgumentException(Resources.ArgumentMatrixSquare); + } + + Order = storage.RowCount; + } + } +} diff --git a/src/Numerics/LinearAlgebra/Complex32/SymmetricDenseMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/SymmetricDenseMatrix.cs new file mode 100644 index 00000000..38cbc136 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Complex32/SymmetricDenseMatrix.cs @@ -0,0 +1,767 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Complex32 +{ + using System; + using Numerics; + using Generic; + using MathNet.Numerics.Distributions; + using MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static; + using Properties; + using Storage; + + /// + /// A Symmetric Matrix class with dense storage. + /// + /// The underlying storage is a one dimensional array in column-major order. + /// The Upper Triangle is stored(it is equal to the Lower Triangle) + [Serializable] + public class SymmetricDenseMatrix : SymmetricMatrix + { + readonly DenseColumnMajorSymmetricMatrixStorage _storage; + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + readonly Complex32[] _data; + + internal SymmetricDenseMatrix(DenseColumnMajorSymmetricMatrixStorage storage) + : base(storage) + { + _storage = storage; + _data = _storage.Data; + } + + /// + /// Initializes a new instance of the class. This matrix is square with a given size. + /// + /// The order of the matrix. + /// + /// If is less than one. + /// + public SymmetricDenseMatrix(int order) + : this(new DenseColumnMajorSymmetricMatrixStorage(order)) + { + } + + /// + /// Initializes a new instance of the class with all entries set to a particular value. + /// + /// + /// The order of the matrix. + /// + /// The value which we assign to each element of the matrix. + public SymmetricDenseMatrix(int order, Complex32 value) + : this(order) + { + for (var i = 0; i < _data.Length; i++) + { + _data[i] = value; + } + } + + /// + /// Initializes a new instance of the class from a one dimensional array. This constructor + /// will reference the one dimensional array and not copy it. + /// + /// The size of the square matrix. + /// + /// The one dimensional array to create this matrix from. Column-major and row-major order is identical on a symmetric matrix: http://en.wikipedia.org/wiki/Row-major_order + /// + /// + /// If does not represent a packed array. + /// + public SymmetricDenseMatrix(int order, Complex32[] array) + : this(new DenseColumnMajorSymmetricMatrixStorage(order, array)) + { + } + + /// + /// Initializes a new instance of the class from a 2D array. This constructor + /// will allocate a completely new memory block for storing the symmetric dense matrix. + /// + /// The 2D array to create this matrix from. + /// + /// If is not a square array. + /// + /// + /// If is not a symmetric array. + /// + public SymmetricDenseMatrix(Complex32[,] array) + : this(array.GetLength(0)) + { + if (!CheckIfSymmetric(array)) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = array[row, column]; + } + } + } + + /// + /// Initializes a new instance of the class, copying + /// the values from the given matrix. Matrix must be Symmetric. + /// + /// The matrix to copy. + /// + /// If is not a square matrix. + /// + /// + /// If is not a symmetric matrix. + /// + public SymmetricDenseMatrix(Matrix matrix) + : this(matrix.RowCount) + { + var symmetricMatrix = matrix as SymmetricDenseMatrix; + + if (!matrix.IsSymmetric) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + if (symmetricMatrix == null) + { + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = matrix[row, column]; + } + } + } + else + { + matrix.CopyTo(this); + } + } + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + public Complex32[] Data + { + get { return _data; } + } + + /// + /// Creates a SymmetricDenseMatrix for the given number of rows and columns. + /// If rows and columns are not equal, returns a DenseMatrix instead. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// True if all fields must be mutable (e.g. not a diagonal matrix). + /// + /// A DenseMatrix or SymmetricDenseMatrix with the given dimensions. + /// + /// /// + /// If is not equal to . + /// Symmetric arrays are always square + /// + public override Matrix CreateMatrix(int numberOfRows, int numberOfColumns, bool fullyMutable = false) + { + if (numberOfRows != numberOfColumns || fullyMutable) + { + return new DenseMatrix(numberOfRows, numberOfColumns); + } + + return new SymmetricDenseMatrix(numberOfRows); + } + + /// + /// Creates a with a the given dimension. + /// + /// The size of the vector. + /// True if all fields must be mutable. + /// + /// A with the given dimension. + /// + public override Vector CreateVector(int size, bool fullyMutable = false) + { + return new DenseVector(size); + } + + #region Static constructors for special matrices. + + /// + /// Initializes a square with all zero's except for ones on the diagonal. + /// + /// the size of the square matrix. + /// A symmetric dense identity matrix. + /// + /// If is less than one. + /// + public static SymmetricDenseMatrix Identity(int order) + { + var m = new SymmetricDenseMatrix(order); + for (var i = 0; i < order; i++) + { + m.At(i, i, Complex32.One); + } + + return m; + } + + #endregion + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The matrix to store the result of add + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + protected override void DoAdd(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoAdd(other, result); + } + else + { + Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The matrix to store the result of the subtraction. + protected override void DoSubtract(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoSubtract(other, result); + } + else + { + Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// The scalar to multiply the matrix with. + /// The matrix to store the result of the multiplication. + protected override void DoMultiply(Complex32 scalar, Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + if (denseResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + } + } + + /// + /// Multiplies this matrix with a vector and places the results into the result vector. + /// + /// The vector to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Vector rightSide, Vector result) + { + var denseRight = rightSide as DenseVector; + var denseResult = result as DenseVector; + + if (denseRight == null || denseResult == null) + { + base.DoMultiply(rightSide, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(rightSide, result); + } + } + + /// + /// Multiplies this matrix with another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(other, result); + } + } + + /// + /// Multiplies this matrix with transpose of another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoTransposeAndMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoTransposeAndMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoTransposeAndMultiply(other, result); + } + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// The result of the negation. + protected override void DoNegate(Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + + if (denseResult == null) + { + base.DoNegate(result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + } + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise multiply with this one. + /// The matrix to store the result of the pointwise multiplication. + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise divide this one by. + /// The matrix to store the result of the pointwise division. + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix LowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column <= row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix StrictlyLowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column < row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix UpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix StrictlyUpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row + 1; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Computes the trace of this matrix. + /// + /// The trace of this matrix + public override Complex32 Trace() + { + // Matrix is always square. + var sum = Complex32.Zero; + for (var i = 0; i < RowCount; i++) + { + sum += At(i, i); + } + + return sum; + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = Convert.ToSingle(distribution.Sample()); + } + } + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = distribution.Sample(); + } + } + } + + /// + /// Adds two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to add. + /// The right matrix to add. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Add(rightSide); + } + + /// + /// Returns a Matrix containing the same values of . + /// + /// The matrix to get the values from. + /// A matrix containing a the same values as . + /// If is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Clone(); + } + + /// + /// Subtracts two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to subtract. + /// The right matrix to subtract. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Subtract(rightSide); + } + + /// + /// Negates each element of the matrix. + /// + /// The matrix to negate. + /// A matrix containing the negated values. + /// If is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Negate(); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, Complex32 rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(Complex32 leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Multiply(leftSide); + } + + /// + /// Multiplies two matrices. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to multiply. + /// The right matrix to multiply. + /// The result of multiplication. + /// If or is . + /// If the dimensions of or don't conform. + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide.ColumnCount != rightSide.RowCount) + { + throw new ArgumentException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix and a Vector. + /// + /// The matrix to multiply. + /// The vector to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(SymmetricDenseMatrix leftSide, DenseVector rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (DenseVector)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Vector and a Matrix. + /// + /// The vector to multiply. + /// The matrix to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(DenseVector leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (DenseVector)rightSide.LeftMultiply(leftSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator %(SymmetricDenseMatrix leftSide, Complex32 rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Modulus(rightSide); + } + } +} \ No newline at end of file diff --git a/src/Numerics/LinearAlgebra/Complex32/SymmetricMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/SymmetricMatrix.cs new file mode 100644 index 00000000..e283fd63 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Complex32/SymmetricMatrix.cs @@ -0,0 +1,657 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Complex32 +{ + using System; + using Numerics; + using Generic; + using Distributions; + using Properties; + using Storage; + + /// + /// Abstract class for symmetric matrices. + /// + [Serializable] + public abstract class SymmetricMatrix : SquareMatrix + { + /// + /// Initializes a new instance of the class. + /// + protected SymmetricMatrix(MatrixStorage storage) + : base(storage) + { + } + + /// + /// Returns a value indicating whether the array is symmetric. + /// + /// + /// The array to check for symmetry. + /// + /// + /// True is array is symmetric, false if not symmetric. + /// + public static bool CheckIfSymmetric(Complex32[,] array) + { + var rows = array.GetLength(0); + var columns = array.GetLength(1); + + if (rows != columns) + { + return false; + } + + for (var row = 0; row < rows; row++) + { + for (var column = 0; column < columns; column++) + { + if (column >= row) + { + continue; + } + + if (!array[row, column].Equals(array[column, row])) + { + return false; + } + } + } + + return true; + } + + /// + /// Gets a value indicating whether this matrix is symmetric. + /// + public override sealed bool IsSymmetric + { + get + { + return true; + } + } + + /// + /// Returns the transpose of this matrix. The transpose is equal and this method returns a reference to this matrix. + /// + /// + /// The transpose of this matrix. + /// + public override sealed Matrix Transpose() + { + return this.Clone(); + } + + /// + /// Returns the conjugate transpose of this matrix. + /// + /// The conjugate transpose of this matrix. + public override Matrix ConjugateTranspose() + { + var ret = CreateMatrix(ColumnCount, RowCount); + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + ret.At(row, column, At(column, row).Conjugate()); + } + } + + return ret; + } + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The result of the addition. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Add(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoAdd(other, result); + return result; + } + + /// + /// Adds another matrix to this matrix. + /// + /// + /// The matrix to add to this matrix. + /// + /// + /// The matrix to store the result of the addition. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoAdd(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric + non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoAdd(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) + symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The result of the subtraction. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Subtract(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoSubtract(other, result); + return result; + } + + + /// + /// Subtracts another matrix from this matrix. + /// + /// + /// The matrix to subtract to this matrix. + /// + /// + /// The matrix to store the result of subtraction. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoSubtract(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric - non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoSubtract(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) - symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// + /// The scalar to multiply the matrix with. + /// + /// + /// The matrix to store the result of the multiplication. + /// + protected override void DoMultiply(Complex32 scalar, Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * scalar); + } + } + } + } + + /// + /// Multiplies the transpose of this matrix with another matrix and places the results into the result matrix. + /// + /// + /// The matrix to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Matrix other, Matrix result) + { + DoMultiply(other, result); + } + + /// + /// Multiplies the transpose of this matrix with a vector and places the results into the result vector. + /// + /// + /// The vector to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Vector rightSide, Vector result) + { + DoMultiply(rightSide, result); + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// + /// The result of the negation. + /// + protected override void DoNegate(Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoNegate(result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column != ColumnCount; column++) + { + symmetricResult[row, column] = -At(row, column); + } + } + } + } + + /// + /// Pointwise multiplies this matrix with another matrix. + /// + /// The matrix to pointwise multiply with this one. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise multiplication of this matrix and . + public override Matrix PointwiseMultiply(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other, "other"); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseMultiply(other, result); + return result; + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise multiply with this one. + /// + /// + /// The matrix to store the result of the pointwise multiplication. + /// + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise* non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Pointwise divide this matrix by another matrix. + /// + /// The matrix to pointwise subtract this one by. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise division of this matrix and . + public override Matrix PointwiseDivide(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseDivide(other, result); + return result; + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise divide this one by. + /// + /// + /// The matrix to store the result of the pointwise division. + /// + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise/ non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) / symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Computes the modulus for each element of the matrix. + /// + /// + /// The divisor to use. + /// + /// + /// Matrix to store the results in. + /// + protected override void DoModulus(Complex32 divisor, Matrix result) + { + throw new NotImplementedException(); + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, Convert.ToSingle(distribution.Sample())); + } + } + } + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, distribution.Sample()); + } + } + } + } + + /// + /// Creates a new matrix and inserts the given column at the given index. + /// + /// The index of where to insert the column. + /// The column to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of columns. + /// If the size of != the number of rows. + public override Matrix InsertColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Inserting a column is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given array to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, Complex32[] column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given Vector to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Creates a new matrix and inserts the given row at the given index. + /// + /// The index of where to insert the row. + /// The row to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of rows. + /// If the size of != the number of columns. + public override Matrix InsertRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Inserting a row is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given Vector to the specified row. + /// + /// The row to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given array to the specified row. + /// + /// The row to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, Complex32[] row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + } +} diff --git a/src/Numerics/LinearAlgebra/Double/SquareMatrix.cs b/src/Numerics/LinearAlgebra/Double/SquareMatrix.cs new file mode 100644 index 00000000..f7a72e8c --- /dev/null +++ b/src/Numerics/LinearAlgebra/Double/SquareMatrix.cs @@ -0,0 +1,61 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double +{ + using System; + using Properties; + using Storage; + + /// + /// Abstract class for square matrices. + /// + [Serializable] + public abstract class SquareMatrix : Matrix + { + /// + /// Number of rows or columns. + /// + protected readonly int Order; + + /// + /// Initializes a new instance of the class. + /// + /// + /// If the matrix is not square. + /// + protected SquareMatrix(MatrixStorage storage) + : base(storage) + { + if (storage.RowCount != storage.ColumnCount) + { + throw new ArgumentException(Resources.ArgumentMatrixSquare); + } + + Order = storage.RowCount; + } + } +} diff --git a/src/Numerics/LinearAlgebra/Double/SymmetricDenseMatrix.cs b/src/Numerics/LinearAlgebra/Double/SymmetricDenseMatrix.cs new file mode 100644 index 00000000..1a952572 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Double/SymmetricDenseMatrix.cs @@ -0,0 +1,794 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double +{ + using System; + using Generic; + using MathNet.Numerics.Distributions; + using MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static; + using MathNet.Numerics.Threading; + using Properties; + using Storage; + + /// + /// A Symmetric Matrix class with dense storage. + /// + /// The underlying storage is a one dimensional array in column-major order. + /// The Upper Triangle is stored(it is equal to the Lower Triangle) + [Serializable] + public class SymmetricDenseMatrix : SymmetricMatrix + { + readonly DenseColumnMajorSymmetricMatrixStorage _storage; + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + readonly double[] _data; + + internal SymmetricDenseMatrix(DenseColumnMajorSymmetricMatrixStorage storage) + : base(storage) + { + _storage = storage; + _data = _storage.Data; + } + + /// + /// Initializes a new instance of the class. This matrix is square with a given size. + /// + /// The order of the matrix. + /// + /// If is less than one. + /// + public SymmetricDenseMatrix(int order) + : this(new DenseColumnMajorSymmetricMatrixStorage(order)) + { + } + + /// + /// Initializes a new instance of the class with all entries set to a particular value. + /// + /// + /// The order of the matrix. + /// + /// The value which we assign to each element of the matrix. + public SymmetricDenseMatrix(int order, double value) + : this(order) + { + for (var i = 0; i < _data.Length; i++) + { + _data[i] = value; + } + } + + /// + /// Initializes a new instance of the class from a one dimensional array. This constructor + /// will reference the one dimensional array and not copy it. + /// + /// The size of the square matrix. + /// + /// The one dimensional array to create this matrix from. Column-major and row-major order is identical on a symmetric matrix: http://en.wikipedia.org/wiki/Row-major_order + /// + /// + /// If does not represent a packed array. + /// + public SymmetricDenseMatrix(int order, double[] array) + : this(new DenseColumnMajorSymmetricMatrixStorage(order, array)) + { + } + + /// + /// Initializes a new instance of the class from a 2D array. This constructor + /// will allocate a completely new memory block for storing the symmetric dense matrix. + /// + /// The 2D array to create this matrix from. + /// + /// If is not a square array. + /// + /// + /// If is not a symmetric array. + /// + public SymmetricDenseMatrix(double[,] array) + : this(array.GetLength(0)) + { + if (!CheckIfSymmetric(array)) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = array[row, column]; + } + } + } + + /// + /// Initializes a new instance of the class, copying + /// the values from the given matrix. Matrix must be Symmetric. + /// + /// The matrix to copy. + /// + /// If is not a square matrix. + /// + /// + /// If is not a symmetric matrix. + /// + public SymmetricDenseMatrix(Matrix matrix) + : this(matrix.RowCount) + { + var symmetricMatrix = matrix as SymmetricDenseMatrix; + + if (!matrix.IsSymmetric) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + if (symmetricMatrix == null) + { + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = matrix[row, column]; + } + } + } + else + { + matrix.CopyTo(this); + } + } + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + public double[] Data + { + get { return _data; } + } + + /// + /// Creates a SymmetricDenseMatrix for the given number of rows and columns. + /// If rows and columns are not equal, returns a DenseMatrix instead. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// True if all fields must be mutable (e.g. not a diagonal matrix). + /// + /// A DenseMatrix or SymmetricDenseMatrix with the given dimensions. + /// + /// /// + /// If is not equal to . + /// Symmetric arrays are always square + /// + public override Matrix CreateMatrix(int numberOfRows, int numberOfColumns, bool fullyMutable = false) + { + if (numberOfRows != numberOfColumns || fullyMutable) + { + return new DenseMatrix(numberOfRows, numberOfColumns); + } + + return new SymmetricDenseMatrix(numberOfRows); + } + + /// + /// Creates a with a the given dimension. + /// + /// The size of the vector. + /// True if all fields must be mutable. + /// + /// A with the given dimension. + /// + public override Vector CreateVector(int size, bool fullyMutable = false) + { + return new DenseVector(size); + } + + #region Static constructors for special matrices. + + /// + /// Initializes a square with all zero's except for ones on the diagonal. + /// + /// the size of the square matrix. + /// A symmetric dense identity matrix. + /// + /// If is less than one. + /// + public static SymmetricDenseMatrix Identity(int order) + { + var m = new SymmetricDenseMatrix(order); + for (var i = 0; i < order; i++) + { + m.At(i, i, 1.0); + } + + return m; + } + + #endregion + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The matrix to store the result of add + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + protected override void DoAdd(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoAdd(other, result); + } + else + { + Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The matrix to store the result of the subtraction. + protected override void DoSubtract(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoSubtract(other, result); + } + else + { + Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// The scalar to multiply the matrix with. + /// The matrix to store the result of the multiplication. + protected override void DoMultiply(double scalar, Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + if (denseResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + } + } + + /// + /// Multiplies this matrix with a vector and places the results into the result vector. + /// + /// The vector to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Vector rightSide, Vector result) + { + var denseRight = rightSide as DenseVector; + var denseResult = result as DenseVector; + + if (denseRight == null || denseResult == null) + { + base.DoMultiply(rightSide, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(rightSide, result); + } + } + + /// + /// Multiplies this matrix with another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(other, result); + } + } + + /// + /// Multiplies this matrix with transpose of another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoTransposeAndMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoTransposeAndMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoTransposeAndMultiply(other, result); + } + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// The result of the negation. + protected override void DoNegate(Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + + if (denseResult == null) + { + base.DoNegate(result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + } + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise multiply with this one. + /// The matrix to store the result of the pointwise multiplication. + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise divide this one by. + /// The matrix to store the result of the pointwise division. + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix LowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column <= row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix StrictlyLowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column < row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix UpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix StrictlyUpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row + 1; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Computes the modulus for each element of the matrix. + /// + /// The divisor to use. + /// Matrix to store the results in. + protected override void DoModulus(double divisor, Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + + if (denseResult == null) + { + base.DoModulus(divisor, result); + } + else + { + if (!ReferenceEquals(this, result)) + { + CopyTo(result); + } + + CommonParallel.For( + 0, + _data.Length, + index => denseResult._data[index] %= divisor); + } + } + + /// + /// Computes the trace of this matrix. + /// + /// The trace of this matrix + public override double Trace() + { + // Matrix is always square. + var sum = 0.0; + for (var i = 0; i < RowCount; i++) + { + sum += At(i, i); + } + + return sum; + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = distribution.Sample(); + } + } + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = distribution.Sample(); + } + } + } + + /// + /// Adds two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to add. + /// The right matrix to add. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Add(rightSide); + } + + /// + /// Returns a Matrix containing the same values of . + /// + /// The matrix to get the values from. + /// A matrix containing a the same values as . + /// If is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Clone(); + } + + /// + /// Subtracts two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to subtract. + /// The right matrix to subtract. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Subtract(rightSide); + } + + /// + /// Negates each element of the matrix. + /// + /// The matrix to negate. + /// A matrix containing the negated values. + /// If is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Negate(); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, double rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(double leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Multiply(leftSide); + } + + /// + /// Multiplies two matrices. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to multiply. + /// The right matrix to multiply. + /// The result of multiplication. + /// If or is . + /// If the dimensions of or don't conform. + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide.ColumnCount != rightSide.RowCount) + { + throw new ArgumentException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix and a Vector. + /// + /// The matrix to multiply. + /// The vector to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(SymmetricDenseMatrix leftSide, DenseVector rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (DenseVector)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Vector and a Matrix. + /// + /// The vector to multiply. + /// The matrix to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(DenseVector leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (DenseVector)rightSide.LeftMultiply(leftSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator %(SymmetricDenseMatrix leftSide, double rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Modulus(rightSide); + } + } +} \ No newline at end of file diff --git a/src/Numerics/LinearAlgebra/Double/SymmetricMatrix.cs b/src/Numerics/LinearAlgebra/Double/SymmetricMatrix.cs new file mode 100644 index 00000000..a00a0206 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Double/SymmetricMatrix.cs @@ -0,0 +1,653 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double +{ + using System; + using Generic; + using Distributions; + using Properties; + using Storage; + + /// + /// Abstract class for symmetric matrices. + /// + [Serializable] + public abstract class SymmetricMatrix : SquareMatrix + { + /// + /// Initializes a new instance of the class. + /// + protected SymmetricMatrix(MatrixStorage storage) + : base(storage) + { + } + + /// + /// Returns a value indicating whether the array is symmetric. + /// + /// + /// The array to check for symmetry. + /// + /// + /// True is array is symmetric, false if not symmetric. + /// + public static bool CheckIfSymmetric(double[,] array) + { + var rows = array.GetLength(0); + var columns = array.GetLength(1); + + if (rows != columns) + { + return false; + } + + for (var row = 0; row < rows; row++) + { + for (var column = 0; column < columns; column++) + { + if (column >= row) + { + continue; + } + + if (!array[row, column].Equals(array[column, row])) + { + return false; + } + } + } + + return true; + } + + /// + /// Gets a value indicating whether this matrix is symmetric. + /// + public override sealed bool IsSymmetric + { + get + { + return true; + } + } + + /// + /// Returns the transpose of this matrix. The transpose is equal and this method returns a reference to this matrix. + /// + /// + /// The transpose of this matrix. + /// + public override sealed Matrix Transpose() + { + return this.Clone(); + } + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The result of the addition. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Add(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoAdd(other, result); + return result; + } + + /// + /// Adds another matrix to this matrix. + /// + /// + /// The matrix to add to this matrix. + /// + /// + /// The matrix to store the result of the addition. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoAdd(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric + non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoAdd(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) + symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The result of the subtraction. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Subtract(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoSubtract(other, result); + return result; + } + + + /// + /// Subtracts another matrix from this matrix. + /// + /// + /// The matrix to subtract to this matrix. + /// + /// + /// The matrix to store the result of subtraction. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoSubtract(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric - non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoSubtract(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) - symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// + /// The scalar to multiply the matrix with. + /// + /// + /// The matrix to store the result of the multiplication. + /// + protected override void DoMultiply(double scalar, Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * scalar); + } + } + } + } + + /// + /// Multiplies the transpose of this matrix with another matrix and places the results into the result matrix. + /// + /// + /// The matrix to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Matrix other, Matrix result) + { + DoMultiply(other, result); + } + + /// + /// Multiplies the transpose of this matrix with a vector and places the results into the result vector. + /// + /// + /// The vector to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Vector rightSide, Vector result) + { + DoMultiply(rightSide, result); + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// + /// The result of the negation. + /// + protected override void DoNegate(Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoNegate(result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column != ColumnCount; column++) + { + symmetricResult[row, column] = -At(row, column); + } + } + } + } + + /// + /// Pointwise multiplies this matrix with another matrix. + /// + /// The matrix to pointwise multiply with this one. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise multiplication of this matrix and . + public override Matrix PointwiseMultiply(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other, "other"); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseMultiply(other, result); + return result; + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise multiply with this one. + /// + /// + /// The matrix to store the result of the pointwise multiplication. + /// + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise* non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Pointwise divide this matrix by another matrix. + /// + /// The matrix to pointwise subtract this one by. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise division of this matrix and . + public override Matrix PointwiseDivide(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseDivide(other, result); + return result; + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise divide this one by. + /// + /// + /// The matrix to store the result of the pointwise division. + /// + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise/ non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) / symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Computes the modulus for each element of the matrix. + /// + /// + /// The divisor to use. + /// + /// + /// Matrix to store the results in. + /// + protected override void DoModulus(double divisor, Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoModulus(divisor, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) % divisor); + } + } + } + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, distribution.Sample()); + } + } + } + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, distribution.Sample()); + } + } + } + } + + /// + /// Creates a new matrix and inserts the given column at the given index. + /// + /// The index of where to insert the column. + /// The column to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of columns. + /// If the size of != the number of rows. + public override Matrix InsertColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Inserting a column is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given array to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, double[] column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given Vector to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Creates a new matrix and inserts the given row at the given index. + /// + /// The index of where to insert the row. + /// The row to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of rows. + /// If the size of != the number of columns. + public override Matrix InsertRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Inserting a row is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given Vector to the specified row. + /// + /// The row to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given array to the specified row. + /// + /// The row to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, double[] row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + } +} diff --git a/src/Numerics/LinearAlgebra/Generic/Matrix.Arithmetic.cs b/src/Numerics/LinearAlgebra/Generic/Matrix.Arithmetic.cs index bd8d7db8..e23d2b51 100644 --- a/src/Numerics/LinearAlgebra/Generic/Matrix.Arithmetic.cs +++ b/src/Numerics/LinearAlgebra/Generic/Matrix.Arithmetic.cs @@ -481,7 +481,7 @@ namespace MathNet.Numerics.LinearAlgebra.Generic throw DimensionsDontMatch(this, other); } - var result = CreateMatrix(RowCount, other.ColumnCount); + var result = CreateMatrix(RowCount, other.ColumnCount, true); Multiply(other, result); return result; } @@ -550,7 +550,7 @@ namespace MathNet.Numerics.LinearAlgebra.Generic throw DimensionsDontMatch(this, other); } - var result = CreateMatrix(RowCount, other.RowCount); + var result = CreateMatrix(RowCount, other.RowCount, true); TransposeAndMultiply(other, result); return result; } @@ -683,7 +683,7 @@ namespace MathNet.Numerics.LinearAlgebra.Generic throw DimensionsDontMatch(this, other); } - var result = CreateMatrix(ColumnCount, other.ColumnCount); + var result = CreateMatrix(ColumnCount, other.ColumnCount, true); TransposeThisAndMultiply(other, result); return result; } diff --git a/src/Numerics/LinearAlgebra/Single/SquareMatrix.cs b/src/Numerics/LinearAlgebra/Single/SquareMatrix.cs new file mode 100644 index 00000000..c51efdfb --- /dev/null +++ b/src/Numerics/LinearAlgebra/Single/SquareMatrix.cs @@ -0,0 +1,61 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Single +{ + using System; + using Properties; + using Storage; + + /// + /// Abstract class for square matrices. + /// + [Serializable] + public abstract class SquareMatrix : Matrix + { + /// + /// Number of rows or columns. + /// + protected readonly int Order; + + /// + /// Initializes a new instance of the class. + /// + /// + /// If the matrix is not square. + /// + protected SquareMatrix(MatrixStorage storage) + : base(storage) + { + if (storage.RowCount != storage.ColumnCount) + { + throw new ArgumentException(Resources.ArgumentMatrixSquare); + } + + Order = storage.RowCount; + } + } +} diff --git a/src/Numerics/LinearAlgebra/Single/SymmetricDenseMatrix.cs b/src/Numerics/LinearAlgebra/Single/SymmetricDenseMatrix.cs new file mode 100644 index 00000000..cef55387 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Single/SymmetricDenseMatrix.cs @@ -0,0 +1,794 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Single +{ + using System; + using Generic; + using MathNet.Numerics.Distributions; + using MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static; + using MathNet.Numerics.Threading; + using Properties; + using Storage; + + /// + /// A Symmetric Matrix class with dense storage. + /// + /// The underlying storage is a one dimensional array in column-major order. + /// The Upper Triangle is stored(it is equal to the Lower Triangle) + [Serializable] + public class SymmetricDenseMatrix : SymmetricMatrix + { + readonly DenseColumnMajorSymmetricMatrixStorage _storage; + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + readonly float[] _data; + + internal SymmetricDenseMatrix(DenseColumnMajorSymmetricMatrixStorage storage) + : base(storage) + { + _storage = storage; + _data = _storage.Data; + } + + /// + /// Initializes a new instance of the class. This matrix is square with a given size. + /// + /// The order of the matrix. + /// + /// If is less than one. + /// + public SymmetricDenseMatrix(int order) + : this(new DenseColumnMajorSymmetricMatrixStorage(order)) + { + } + + /// + /// Initializes a new instance of the class with all entries set to a particular value. + /// + /// + /// The order of the matrix. + /// + /// The value which we assign to each element of the matrix. + public SymmetricDenseMatrix(int order, float value) + : this(order) + { + for (var i = 0; i < _data.Length; i++) + { + _data[i] = value; + } + } + + /// + /// Initializes a new instance of the class from a one dimensional array. This constructor + /// will reference the one dimensional array and not copy it. + /// + /// The size of the square matrix. + /// + /// The one dimensional array to create this matrix from. Column-major and row-major order is identical on a symmetric matrix: http://en.wikipedia.org/wiki/Row-major_order + /// + /// + /// If does not represent a packed array. + /// + public SymmetricDenseMatrix(int order, float[] array) + : this(new DenseColumnMajorSymmetricMatrixStorage(order, array)) + { + } + + /// + /// Initializes a new instance of the class from a 2D array. This constructor + /// will allocate a completely new memory block for storing the symmetric dense matrix. + /// + /// The 2D array to create this matrix from. + /// + /// If is not a square array. + /// + /// + /// If is not a symmetric array. + /// + public SymmetricDenseMatrix(float[,] array) + : this(array.GetLength(0)) + { + if (!CheckIfSymmetric(array)) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = array[row, column]; + } + } + } + + /// + /// Initializes a new instance of the class, copying + /// the values from the given matrix. Matrix must be Symmetric. + /// + /// The matrix to copy. + /// + /// If is not a square matrix. + /// + /// + /// If is not a symmetric matrix. + /// + public SymmetricDenseMatrix(Matrix matrix) + : this(matrix.RowCount) + { + var symmetricMatrix = matrix as SymmetricDenseMatrix; + + if (!matrix.IsSymmetric) + { + throw new ArgumentException(Resources.ArgumentMatrixSymmetric); + } + + if (symmetricMatrix == null) + { + var indexer = new PackedStorageIndexerUpper(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + _data[indexer.Of(row, column)] = matrix[row, column]; + } + } + } + else + { + matrix.CopyTo(this); + } + } + + /// + /// Gets the matrix's data. + /// + /// The matrix's data. + public float[] Data + { + get { return _data; } + } + + /// + /// Creates a SymmetricDenseMatrix for the given number of rows and columns. + /// If rows and columns are not equal, returns a DenseMatrix instead. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// True if all fields must be mutable (e.g. not a diagonal matrix). + /// + /// A DenseMatrix or SymmetricDenseMatrix with the given dimensions. + /// + /// /// + /// If is not equal to . + /// Symmetric arrays are always square + /// + public override Matrix CreateMatrix(int numberOfRows, int numberOfColumns, bool fullyMutable = false) + { + if (numberOfRows != numberOfColumns || fullyMutable) + { + return new DenseMatrix(numberOfRows, numberOfColumns); + } + + return new SymmetricDenseMatrix(numberOfRows); + } + + /// + /// Creates a with a the given dimension. + /// + /// The size of the vector. + /// True if all fields must be mutable. + /// + /// A with the given dimension. + /// + public override Vector CreateVector(int size, bool fullyMutable = false) + { + return new DenseVector(size); + } + + #region Static constructors for special matrices. + + /// + /// Initializes a square with all zero's except for ones on the diagonal. + /// + /// the size of the square matrix. + /// A symmetric dense identity matrix. + /// + /// If is less than one. + /// + public static SymmetricDenseMatrix Identity(int order) + { + var m = new SymmetricDenseMatrix(order); + for (var i = 0; i < order; i++) + { + m.At(i, i, 1.0f); + } + + return m; + } + + #endregion + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The matrix to store the result of add + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + protected override void DoAdd(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoAdd(other, result); + } + else + { + Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The matrix to store the result of the subtraction. + protected override void DoSubtract(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + if (denseOther == null || denseResult == null) + { + base.DoSubtract(other, result); + } + else + { + Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// The scalar to multiply the matrix with. + /// The matrix to store the result of the multiplication. + protected override void DoMultiply(float scalar, Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + if (denseResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + } + } + + /// + /// Multiplies this matrix with a vector and places the results into the result vector. + /// + /// The vector to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Vector rightSide, Vector result) + { + var denseRight = rightSide as DenseVector; + var denseResult = result as DenseVector; + + if (denseRight == null || denseResult == null) + { + base.DoMultiply(rightSide, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(rightSide, result); + } + } + + /// + /// Multiplies this matrix with another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoMultiply(other, result); + } + } + + /// + /// Multiplies this matrix with transpose of another matrix and places the results into the result matrix. + /// + /// The matrix to multiply with. + /// The result of the multiplication. + protected override void DoTransposeAndMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoTransposeAndMultiply(other, result); + } + else + { + // TODO: Change this when symmetric methods are implemented in the Linear Algebra Providers. + base.DoTransposeAndMultiply(other, result); + } + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// The result of the negation. + protected override void DoNegate(Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + + if (denseResult == null) + { + base.DoNegate(result); + } + else + { + Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + } + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise multiply with this one. + /// The matrix to store the result of the pointwise multiplication. + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// The matrix to pointwise divide this one by. + /// The matrix to store the result of the pointwise division. + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var denseOther = other as SymmetricDenseMatrix; + var denseResult = result as SymmetricDenseMatrix; + + if (denseOther == null || denseResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + } + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix LowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column <= row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the lower triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The lower triangle of this matrix. + public override Matrix StrictlyLowerTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = 0; column < row; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix UpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Returns a new matrix containing the upper triangle of this matrix. The new matrix + /// does not contain the diagonal elements of this matrix. + /// + /// The upper triangle of this matrix. + public override Matrix StrictlyUpperTriangle() + { + var ret = new DenseMatrix(Order); + for (var row = 0; row < Order; row++) + { + for (var column = row + 1; column < Order; column++) + { + ret[row, column] = At(row, column); + } + } + + return ret; + } + + /// + /// Computes the modulus for each element of the matrix. + /// + /// The divisor to use. + /// Matrix to store the results in. + protected override void DoModulus(float divisor, Matrix result) + { + var denseResult = result as SymmetricDenseMatrix; + + if (denseResult == null) + { + base.DoModulus(divisor, result); + } + else + { + if (!ReferenceEquals(this, result)) + { + CopyTo(result); + } + + CommonParallel.For( + 0, + _data.Length, + index => denseResult._data[index] %= divisor); + } + } + + /// + /// Computes the trace of this matrix. + /// + /// The trace of this matrix + public override float Trace() + { + // Matrix is always square. + var sum = 0.0f; + for (var i = 0; i < RowCount; i++) + { + sum += At(i, i); + } + + return sum; + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = Convert.ToSingle(distribution.Sample()); + } + } + } + + /// + /// Populates a symmetric matrix with random elements. + /// + /// The symmetric matrix to populate. + /// Continuous Random Distribution to generate elements from. + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var denseMatrix = matrix as SymmetricDenseMatrix; + + if (denseMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var i = 0; i < denseMatrix._data.Length; i++) + { + denseMatrix._data[i] = distribution.Sample(); + } + } + } + + /// + /// Adds two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to add. + /// The right matrix to add. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Add(rightSide); + } + + /// + /// Returns a Matrix containing the same values of . + /// + /// The matrix to get the values from. + /// A matrix containing a the same values as . + /// If is . + public static SymmetricDenseMatrix operator +(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Clone(); + } + + /// + /// Subtracts two matrices together and returns the results. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to subtract. + /// The right matrix to subtract. + /// The result of the addition. + /// If and don't have the same dimensions. + /// If or is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (leftSide.RowCount != rightSide.RowCount) + { + throw new ArgumentOutOfRangeException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Subtract(rightSide); + } + + /// + /// Negates each element of the matrix. + /// + /// The matrix to negate. + /// A matrix containing the negated values. + /// If is . + public static SymmetricDenseMatrix operator -(SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Negate(); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, float rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator *(float leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (SymmetricDenseMatrix)rightSide.Multiply(leftSide); + } + + /// + /// Multiplies two matrices. + /// + /// This operator will allocate new memory for the result. It will + /// choose the representation of either or depending on which + /// is denser. + /// The left matrix to multiply. + /// The right matrix to multiply. + /// The result of multiplication. + /// If or is . + /// If the dimensions of or don't conform. + public static SymmetricDenseMatrix operator *(SymmetricDenseMatrix leftSide, SymmetricDenseMatrix rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + if (leftSide.ColumnCount != rightSide.RowCount) + { + throw new ArgumentException(Resources.ArgumentMatrixDimensions); + } + + return (SymmetricDenseMatrix)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Matrix and a Vector. + /// + /// The matrix to multiply. + /// The vector to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(SymmetricDenseMatrix leftSide, DenseVector rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (DenseVector)leftSide.Multiply(rightSide); + } + + /// + /// Multiplies a Vector and a Matrix. + /// + /// The vector to multiply. + /// The matrix to multiply. + /// The result of multiplication. + /// If or is . + public static DenseVector operator *(DenseVector leftSide, SymmetricDenseMatrix rightSide) + { + if (rightSide == null) + { + throw new ArgumentNullException("rightSide"); + } + + return (DenseVector)rightSide.LeftMultiply(leftSide); + } + + /// + /// Multiplies a Matrix by a constant and returns the result. + /// + /// The matrix to multiply. + /// The constant to multiply the matrix by. + /// The result of the multiplication. + /// If is . + public static SymmetricDenseMatrix operator %(SymmetricDenseMatrix leftSide, float rightSide) + { + if (leftSide == null) + { + throw new ArgumentNullException("leftSide"); + } + + return (SymmetricDenseMatrix)leftSide.Modulus(rightSide); + } + } +} \ No newline at end of file diff --git a/src/Numerics/LinearAlgebra/Single/SymmetricMatrix.cs b/src/Numerics/LinearAlgebra/Single/SymmetricMatrix.cs new file mode 100644 index 00000000..1ec1d112 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Single/SymmetricMatrix.cs @@ -0,0 +1,653 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Single +{ + using System; + using Generic; + using Distributions; + using Properties; + using Storage; + + /// + /// Abstract class for symmetric matrices. + /// + [Serializable] + public abstract class SymmetricMatrix : SquareMatrix + { + /// + /// Initializes a new instance of the class. + /// + protected SymmetricMatrix(MatrixStorage storage) + : base(storage) + { + } + + /// + /// Returns a value indicating whether the array is symmetric. + /// + /// + /// The array to check for symmetry. + /// + /// + /// True is array is symmetric, false if not symmetric. + /// + public static bool CheckIfSymmetric(float[,] array) + { + var rows = array.GetLength(0); + var columns = array.GetLength(1); + + if (rows != columns) + { + return false; + } + + for (var row = 0; row < rows; row++) + { + for (var column = 0; column < columns; column++) + { + if (column >= row) + { + continue; + } + + if (!array[row, column].Equals(array[column, row])) + { + return false; + } + } + } + + return true; + } + + /// + /// Gets a value indicating whether this matrix is symmetric. + /// + public override sealed bool IsSymmetric + { + get + { + return true; + } + } + + /// + /// Returns the transpose of this matrix. The transpose is equal and this method returns a reference to this matrix. + /// + /// + /// The transpose of this matrix. + /// + public override sealed Matrix Transpose() + { + return this.Clone(); + } + + /// + /// Adds another matrix to this matrix. + /// + /// The matrix to add to this matrix. + /// The result of the addition. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Add(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoAdd(other, result); + return result; + } + + /// + /// Adds another matrix to this matrix. + /// + /// + /// The matrix to add to this matrix. + /// + /// + /// The matrix to store the result of the addition. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoAdd(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric + non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoAdd(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) + symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Subtracts another matrix from this matrix. + /// + /// The matrix to subtract. + /// The result of the subtraction. + /// If the other matrix is . + /// If the two matrices don't have the same dimensions. + public override Matrix Subtract(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (other.RowCount != RowCount || other.ColumnCount != ColumnCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + DoSubtract(other, result); + return result; + } + + + /// + /// Subtracts another matrix from this matrix. + /// + /// + /// The matrix to subtract to this matrix. + /// + /// + /// The matrix to store the result of subtraction. + /// + /// + /// If the other matrix is . + /// + /// + /// If the two matrices don't have the same dimensions. + /// + protected override void DoSubtract(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric - non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoSubtract(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) - symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Multiplies each element of the matrix by a scalar and places results into the result matrix. + /// + /// + /// The scalar to multiply the matrix with. + /// + /// + /// The matrix to store the result of the multiplication. + /// + protected override void DoMultiply(float scalar, Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoMultiply(scalar, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * scalar); + } + } + } + } + + /// + /// Multiplies the transpose of this matrix with another matrix and places the results into the result matrix. + /// + /// + /// The matrix to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Matrix other, Matrix result) + { + DoMultiply(other, result); + } + + /// + /// Multiplies the transpose of this matrix with a vector and places the results into the result vector. + /// + /// + /// The vector to multiply with. + /// + /// + /// The result of the multiplication. + /// + protected override sealed void DoTransposeThisAndMultiply(Vector rightSide, Vector result) + { + DoMultiply(rightSide, result); + } + + /// + /// Negate each element of this matrix and place the results into the result matrix. + /// + /// + /// The result of the negation. + /// + protected override void DoNegate(Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoNegate(result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column != ColumnCount; column++) + { + symmetricResult[row, column] = -At(row, column); + } + } + } + } + + /// + /// Pointwise multiplies this matrix with another matrix. + /// + /// The matrix to pointwise multiply with this one. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise multiplication of this matrix and . + public override Matrix PointwiseMultiply(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other, "other"); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseMultiply(other, result); + return result; + } + + /// + /// Pointwise multiplies this matrix with another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise multiply with this one. + /// + /// + /// The matrix to store the result of the pointwise multiplication. + /// + protected override void DoPointwiseMultiply(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise* non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseMultiply(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) * symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Pointwise divide this matrix by another matrix. + /// + /// The matrix to pointwise subtract this one by. + /// If the other matrix is . + /// If this matrix and are not the same size. + /// A new matrix that is the pointwise division of this matrix and . + public override Matrix PointwiseDivide(Matrix other) + { + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (ColumnCount != other.ColumnCount || RowCount != other.RowCount) + { + throw DimensionsDontMatch(this, other); + } + + Matrix result; + if (other is SymmetricMatrix) + { + result = CreateMatrix(RowCount, ColumnCount); + } + else + { + result = CreateMatrix(RowCount, ColumnCount, true); + } + + PointwiseDivide(other, result); + return result; + } + + /// + /// Pointwise divide this matrix by another matrix and stores the result into the result matrix. + /// + /// + /// The matrix to pointwise divide this one by. + /// + /// + /// The matrix to store the result of the pointwise division. + /// + protected override void DoPointwiseDivide(Matrix other, Matrix result) + { + var symmetricOther = other as SymmetricMatrix; + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult != null && !other.IsSymmetric) + { + throw new InvalidOperationException("Symmetric pointwise/ non-symmetric matrix cannot be a symmetric matrix"); + } + + if (symmetricOther == null || symmetricResult == null) + { + base.DoPointwiseDivide(other, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) / symmetricOther.At(row, column)); + } + } + } + } + + /// + /// Computes the modulus for each element of the matrix. + /// + /// + /// The divisor to use. + /// + /// + /// Matrix to store the results in. + /// + protected override void DoModulus(float divisor, Matrix result) + { + var symmetricResult = result as SymmetricMatrix; + + if (symmetricResult == null) + { + base.DoModulus(divisor, result); + } + else + { + for (var row = 0; row < RowCount; row++) + { + for (var column = row; column < ColumnCount; column++) + { + symmetricResult.At(row, column, At(row, column) % divisor); + } + } + } + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IContinuousDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, Convert.ToSingle(distribution.Sample())); + } + } + } + } + + /// + /// Populates a matrix with random elements. + /// + /// + /// The matrix to populate. + /// + /// + /// Continuous Random Distribution to generate elements from. + /// + protected override void DoRandom(Matrix matrix, IDiscreteDistribution distribution) + { + var symmetricMatrix = matrix as SymmetricMatrix; + if (symmetricMatrix == null) + { + base.DoRandom(matrix, distribution); + } + else + { + for (var row = 0; row < matrix.RowCount; row++) + { + for (var column = row; column < matrix.ColumnCount; column++) + { + symmetricMatrix.At(row, column, distribution.Sample()); + } + } + } + } + + /// + /// Creates a new matrix and inserts the given column at the given index. + /// + /// The index of where to insert the column. + /// The column to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of columns. + /// If the size of != the number of rows. + public override Matrix InsertColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Inserting a column is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given array to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, float[] column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given Vector to the specified column. The changes retain the symmetry of the matrix. + /// + /// The column to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of columns. + /// If the size of does not + /// equal the number of rows of this Matrix. + public override void SetColumn(int columnIndex, Vector column) + { + throw new InvalidOperationException("Setting a column is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Creates a new matrix and inserts the given row at the given index. + /// + /// The index of where to insert the row. + /// The row to insert. + /// A new matrix with the inserted column. + /// If is . + /// If is < zero or > the number of rows. + /// If the size of != the number of columns. + public override Matrix InsertRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Inserting a row is not supported on a symmetric matrix. Symmetric matrices are square"); + } + + /// + /// Copies the values of the given Vector to the specified row. + /// + /// The row to copy the values to. + /// The vector to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, Vector row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + + /// + /// Copies the values of the given array to the specified row. + /// + /// The row to copy the values to. + /// The array to copy the values from. + /// If is . + /// If is less than zero, + /// or greater than or equal to the number of rows. + /// If the size of does not + /// equal the number of columns of this Matrix. + public override void SetRow(int rowIndex, float[] row) + { + throw new InvalidOperationException("Setting a row is not supported on a symmetric matrix. It will violate symmetry"); + } + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorSymmetricMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorSymmetricMatrixStorage.cs new file mode 100644 index 00000000..3d0ee7f7 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorSymmetricMatrixStorage.cs @@ -0,0 +1,84 @@ +using System; + +namespace MathNet.Numerics.LinearAlgebra.Storage +{ + using MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static; + using MathNet.Numerics.Properties; + + public class DenseColumnMajorSymmetricMatrixStorage : SymmetricMatrixStorage + where T : struct, IEquatable, IFormattable + { + // [ruegg] public fields are OK here + + public readonly T[] Data; + + public readonly PackedStorageIndexerUpper Indexer; + + internal DenseColumnMajorSymmetricMatrixStorage(int order) + : base(order) + { + Indexer = new PackedStorageIndexerUpper(order); + Data = new T[Indexer.DataLength]; + } + + internal DenseColumnMajorSymmetricMatrixStorage(int order, T[] data) + : base(order) + { + if (data == null) + { + throw new ArgumentNullException("data"); + } + + Indexer = new PackedStorageIndexerUpper(order); + + if (data.Length != Indexer.DataLength) + { + throw new ArgumentOutOfRangeException("data", string.Format(Resources.ArgumentArrayWrongLength, Indexer.DataLength)); + } + + Data = data; + } + + /// + /// Retrieves the requested element without range checking. + /// + /// + /// The row of the element. + /// + /// + /// The column of the element. + /// + /// + /// The requested element. + /// + /// Not range-checked. + public override T At(int row, int column) + { + var r = Math.Min(row, column); + var c = Math.Max(row, column); + return Data[Indexer.Of(r, c)]; + } + + /// + /// Sets the element without range checking. + /// + /// The row of the element. + /// The column of the element. + /// The value to set the element to. + /// WARNING: This method is not thread safe. Use "lock" with it and be sure to avoid deadlocks. + public override void At(int row, int column, T value) + { + if (row > column) + { + throw new IndexOutOfRangeException("Setting an element in the strictly lower triangle of a symmetric matrix is disabled to avoid errors"); + } + + Data[Indexer.Of(row, column)] = value; + } + + public override void Clear() + { + Array.Clear(Data, 0, Data.Length); + } + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/Indexers/IStorageIndexer.cs b/src/Numerics/LinearAlgebra/Storage/Indexers/IStorageIndexer.cs new file mode 100644 index 00000000..7a7e5975 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/Indexers/IStorageIndexer.cs @@ -0,0 +1,28 @@ +namespace MathNet.Numerics.LinearAlgebra.Storage.Indexers +{ + /// + /// Abstract class that defines common features for all storage schemes. + /// + public interface IStorageIndexer + { + /// + /// Retrieves the index of the requested element without parameter checking. + /// + /// The row of the element. + /// + /// The column of the element. + /// + /// The requested index. + /// + int Of(int row, int column); + + /// + /// Retrieves the index of the requested diagonal element without parameter checking. + /// + /// The row=column of the diagonal element. + /// + /// The requested index. + /// + int OfDiagonal(int row); + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/Indexers/Static/PackedStorageIndexer.cs b/src/Numerics/LinearAlgebra/Storage/Indexers/Static/PackedStorageIndexer.cs new file mode 100644 index 00000000..56a22c4a --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/Indexers/Static/PackedStorageIndexer.cs @@ -0,0 +1,50 @@ +namespace MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static +{ + using System; + using Properties; + + /// + /// A class for managing indexing when using Packed Storage, which is a column-major packing scheme for dense Symmetric, Hermitian or Triangular square matrices. + /// + public abstract class PackedStorageIndexer : StaticStorageIndexer + { + /// + /// Number of rows or columns. + /// + protected readonly int Order; + + /// + /// Length of the stored data. + /// + private readonly int _dataLength; + + /// + /// Initializes a new instance of the class. + /// + /// + /// The order of the matrix. + /// + /// is out of range. + protected PackedStorageIndexer(int order) + { + if (order <= 0) + { + throw new ArgumentOutOfRangeException(Resources.MatrixRowsOrColumnsMustBePositive); + } + + Order = order; + _dataLength = order * (order + 1) / 2; + } + + /// + /// Gets the length of the stored data. + /// + public override int DataLength + { + get + { + return _dataLength; + } + } + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/Indexers/Static/PackedStorageIndexerUpper.cs b/src/Numerics/LinearAlgebra/Storage/Indexers/Static/PackedStorageIndexerUpper.cs new file mode 100644 index 00000000..9989b864 --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/Indexers/Static/PackedStorageIndexerUpper.cs @@ -0,0 +1,88 @@ +namespace MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static +{ + using System; + + /// + /// A class for managing indexes when using Packed Storage, which is a column-major packing scheme for Symmetric, Hermitian or Triangular square matrices. + /// This variation provides indexes for storing the upper triangle of a matrix (row less than or equal to column). + /// + public class PackedStorageIndexerUpper : PackedStorageIndexer + { + /// + /// Initializes a new instance of the class. + /// + /// + /// The order of the matrix. + /// + internal PackedStorageIndexerUpper(int order) + : base(order) + { + } + + /// + /// Gets the index of the given element. + /// + /// + /// The row of the element. + /// + /// + /// The column of the element. + /// + /// + /// This method is parameter checked. and to get values without parameter checking. + /// + public override int this[int row, int column] + { + get + { + if (row < 0 || row >= Order) + { + throw new ArgumentOutOfRangeException("row"); + } + + if (column < 0 || column >= Order) + { + throw new ArgumentOutOfRangeException("column"); + } + + if (row > column) + { + throw new ArgumentException("Row must be less than or equal to column"); + } + + return this.Of(row, column); + } + } + + /// + /// Retrieves the index of the requested element without parameter checking. Row must be less than or equal to column. + /// + /// + /// The row of the element. + /// + /// + /// The column of the element. + /// + /// + /// The requested index. + /// + public override int Of(int row, int column) + { + return row + ((column * (column + 1)) / 2); + } + + /// + /// Retrieves the index of the requested diagonal element without parameter checking. + /// + /// + /// The row=column of the diagonal element. + /// + /// + /// The requested index. + /// + public override int OfDiagonal(int row) + { + return (row * (row + 3)) / 2; + } + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/Indexers/Static/StaticStorageIndexer.cs b/src/Numerics/LinearAlgebra/Storage/Indexers/Static/StaticStorageIndexer.cs new file mode 100644 index 00000000..d447fbee --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/Indexers/Static/StaticStorageIndexer.cs @@ -0,0 +1,56 @@ +namespace MathNet.Numerics.LinearAlgebra.Storage.Indexers.Static +{ + /// + /// Classes that contain indexing information of a static storage scheme. + /// + /// + /// A static storage scheme is always the same and only depends on the size of the matrix. + /// + public abstract class StaticStorageIndexer : IStorageIndexer + { + /// + /// Gets the index of the given element. + /// + /// + /// The row of the element. + /// + /// + /// The column of the element. + /// + /// + /// This method is parameter checked. and to get values without parameter checking. + /// + public abstract int this[int row, int column] + { + get; + } + + /// + /// Gets the length of the stored data. + /// + public abstract int DataLength + { + get; + } + + /// + /// Retrieves the index of the requested element without parameter checking. + /// + /// The row of the element. + /// + /// The column of the element. + /// + /// The requested index. + /// + public abstract int Of(int row, int column); + + /// + /// Retrieves the index of the requested diagonal element without parameter checking. + /// + /// The row=column of the diagonal element. + /// + /// The requested index. + /// + public abstract int OfDiagonal(int row); + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/SymmetricMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/SymmetricMatrixStorage.cs new file mode 100644 index 00000000..744c286b --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/SymmetricMatrixStorage.cs @@ -0,0 +1,64 @@ +using System; +using System.Collections.Generic; +using System.Linq; +using System.Text; + +namespace MathNet.Numerics.LinearAlgebra.Storage +{ + using MathNet.Numerics.Properties; + + public abstract class SymmetricMatrixStorage : MatrixStorage + where T : struct, IEquatable, IFormattable + { + // [ruegg] public fields are OK here + + protected SymmetricMatrixStorage(int order) + : base(order, order) + { + } + + public override bool IsFullyMutable + { + get { return false; } + } + + public override bool IsMutable(int row, int column) + { + return row <= column; + } + + public override void Clear() + { + for (var i = 0; i < RowCount; i++) + { + for (var j = i; j < ColumnCount; j++) + { + At(i, j, default(T)); + } + } + } + + public override void Clear(int rowIndex, int rowCount, int columnIndex, int columnCount) + { + for (var i = rowIndex; i < rowIndex + rowCount; i++) + { + for (var j = Math.Max(columnIndex, i); j < columnIndex + columnCount; j++) + { + At(i, j, default(T)); + } + } + } + + /// Parameters assumed to be validated already. + public override void CopyTo(MatrixStorage target, bool skipClearing = false) + { + for (int j = 0; j < ColumnCount; j++) + { + for (int i = 0; i <= j; i++) + { + target.At(i, j, At(i, j)); + } + } + } + } +} diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj index 8a6a9de7..5b5429a6 100644 --- a/src/Numerics/Numerics.csproj +++ b/src/Numerics/Numerics.csproj @@ -125,7 +125,24 @@ + + + + + + + + + + + + + + + + + @@ -339,6 +356,7 @@ + @@ -465,6 +483,7 @@ + diff --git a/src/UnitTests/LinearAlgebraTests/Complex/SymmetricDenseMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Complex/SymmetricDenseMatrixTests.cs new file mode 100644 index 00000000..08d8b3f1 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Complex/SymmetricDenseMatrixTests.cs @@ -0,0 +1,217 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex +{ + using System; + using System.Collections.Generic; + using System.Numerics; + + using MathNet.Numerics.LinearAlgebra.Complex; + + using NUnit.Framework; + + /// + /// Symmetric Dense matrix tests. + /// + public class SymmetricDenseMatrixTests : SymmetricMatrixTests + { + /// + /// Creates a matrix for the given number of rows and columns. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// + /// A matrix with the given dimensions. + /// + protected override Matrix CreateMatrix(int rows, int columns) + { + return new DenseMatrix(rows, columns); + } + + /// + /// Creates a matrix from a 2D array. + /// + /// + /// The 2D array to create this matrix from. + /// + /// + /// A matrix with the given values. + /// + protected override Matrix CreateMatrix(Complex[,] data) + { + if (SymmetricMatrix.CheckIfSymmetric(data)) + { + return new SymmetricDenseMatrix(data); + } + + return new DenseMatrix(data); + } + + /// + /// Creates a vector of the given size. + /// + /// + /// The size of the vector to create. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(int size) + { + return new DenseVector(size); + } + + /// + /// Creates a vector from an array. + /// + /// + /// The array to create this vector from. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(Complex[] data) + { + return new DenseVector(data); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom1DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(3, new[] { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(0.0, 1), new Complex(3.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) }) }, + { "Square3x3", new SymmetricDenseMatrix(3, new[] { new Complex(-1.1, 1), new Complex(2.0, 1), new Complex(1.1, 1), new Complex(3.0, 1), new Complex(0.0, 1), new Complex(6.6, 1) }) }, + { "Square4x4", new SymmetricDenseMatrix(4, new[] { new Complex(1.1, 1), new Complex(2.0, 1), new Complex(5.0, 1), new Complex(-3.0, 1), new Complex(-6.0, 1), new Complex(8.0, 1), new Complex(4.4, 1), new Complex(7.0, 1), new Complex(9.0, 1), new Complex(10.0, 1) }) }, + { "Singular4x4", new SymmetricDenseMatrix(4, new[] { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(5.0, 1), new Complex(0.0, 1), new Complex(0.0, 1), new Complex(0.0, 1), new Complex(4.0, 1), new Complex(7.0, 1), new Complex(0.0, 1), new Complex(10.0, 1) }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(3, new[] { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(2.0, 1), new Complex(3.0, 1), new Complex(0.0, 1), new Complex(3.0, 1) }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(4, new [] { new Complex(0, 1), new Complex(1, 1), new Complex(2, 1), new Complex(3, 1), new Complex(4, 1), new Complex(5, 1), new Complex(6, 1), new Complex(7, 1), new Complex(8, 1), new Complex(9, 1) }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from array is a reference. + /// + [Test] + public void MatrixFrom1DArrayIsReference() + { + var data = new Complex[] { new Complex(1, 1), new Complex(1, 1), new Complex(1, 1), new Complex(1, 1), new Complex(1, 1), new Complex(1, 1) }; + var matrix = new SymmetricDenseMatrix(3, data); + matrix[0, 0] = new Complex(10.0, 2); + Assert.AreEqual(new Complex(10.0, 2), data[0]); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom2DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(new[,] { { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(3.0, 1) }, { new Complex(2.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) }, { new Complex(3.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) } }) }, + { "Square3x3", new SymmetricDenseMatrix(new[,] { { new Complex(-1.1, 1), new Complex(2.0, 1), new Complex(3.0, 1) }, { new Complex(2.0, 1), new Complex(1.1, 1), new Complex(0.0, 1) }, { new Complex(3.0, 1), new Complex(0.0, 1), new Complex(6.6, 1) } }) }, + { "Square4x4", new SymmetricDenseMatrix(new[,] { { new Complex(1.1, 1), new Complex(2.0, 1), new Complex(-3.0, 1), new Complex(4.4, 1) }, { new Complex(2.0, 1), new Complex(5.0, 1), new Complex(-6.0, 1), new Complex(7.0, 1) }, { new Complex(-3.0, 1), new Complex(-6.0, 1), new Complex(8.0, 1), new Complex(9.0, 1) }, { new Complex(4.4, 1), new Complex(7.0, 1), new Complex(9.0, 1), new Complex(10.0, 1) } }) }, + { "Singular4x4", new SymmetricDenseMatrix(new[,] { { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(0.0, 1), new Complex(4.0, 1) }, { new Complex(2.0, 1), new Complex(5.0, 1), new Complex(0.0, 1), new Complex(7.0, 1) }, { new Complex(0.0, 1), new Complex(0.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) }, { new Complex(4.0, 1), new Complex(7.0, 1), new Complex(0.0, 1), new Complex(10.0, 1) } }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(new[,] { { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(3.0, 1) }, { new Complex(2.0, 1), new Complex(2.0, 1), new Complex(0.0, 1) }, { new Complex(3.0, 1), new Complex(0.0, 1), new Complex(3.0, 1) } }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(new [,] { { new Complex(0, 1), new Complex(1, 1), new Complex(3, 1), new Complex(6, 1) }, { new Complex(1, 1), new Complex(2, 1), new Complex(4, 1), new Complex(7, 1) }, { new Complex(3, 1), new Complex(4, 1), new Complex(5, 1), new Complex(8, 1) }, { new Complex(6, 1), new Complex(7, 1), new Complex(8, 1), new Complex(9, 1) } }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from two-dimensional array is a copy. + /// + [Test] + public void MatrixFrom2DArrayIsCopy() + { + var matrix = new DenseMatrix(TestData2D["Singular3x3"]); + matrix[0, 0] = new Complex(10.0, 2); + Assert.AreEqual(new Complex(1.0, 1), TestData2D["Singular3x3"][0, 0]); + } + + /// + /// Can create a matrix with uniform values. + /// + [Test] + public void CanCreateMatrixWithUniformValues() + { + var matrix = new SymmetricDenseMatrix(10, new Complex(10.0, 2)); + var value = new Complex(10.0, 2); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], value); + } + } + } + + /// + /// Can create an identity matrix. + /// + [Test] + public void CanCreateIdentity() + { + var matrix = SymmetricDenseMatrix.Identity(5); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(i == j ? Complex.One : Complex.Zero, matrix[i, j]); + } + } + } + + /// + /// Identity with wrong order throws ArgumentOutOfRangeException. + /// + /// The size of the square matrix + [TestCase(0)] + [TestCase(-1)] + public void IdentityWithWrongOrderThrowsArgumentOutOfRangeException(int order) + { + Assert.Throws(() => SymmetricDenseMatrix.Identity(order)); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Complex/SymmetricMatrixTests.Arithmetic.cs b/src/UnitTests/LinearAlgebraTests/Complex/SymmetricMatrixTests.Arithmetic.cs new file mode 100644 index 00000000..d3c59ac3 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Complex/SymmetricMatrixTests.Arithmetic.cs @@ -0,0 +1,192 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex +{ + using System.Collections.Generic; + using LinearAlgebra.Complex; + using NUnit.Framework; + using System.Numerics; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices + /// + public abstract partial class SymmetricMatrixTests + { + /// + /// Setup test matrices. + /// Singular and Square matrices are overridden here with symmetric ones so that calls to base methods work as intended. + /// Additional NonSymmetric matrices are defined for some tests. + /// + [SetUp] + public override void SetupMatrices() + { + TestData2D = new Dictionary + { + { "Singular3x3", new[,] { { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(3.0, 1) }, { new Complex(2.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) }, { new Complex(3.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) } } }, + { "Square3x3", new[,] { { new Complex(-1.1, 1), new Complex(2.0, 1), new Complex(3.0, 1) }, { new Complex(2.0, 1), new Complex(1.1, 1), new Complex(0.0, 1) }, { new Complex(3.0, 1), new Complex(0.0, 1), new Complex(6.6, 1) } } }, + { "Square4x4", new[,] { { new Complex(1.1, 1), new Complex(2.0, 1), new Complex(-3.0, 1), new Complex(4.4, 1) }, { new Complex(2.0, 1), new Complex(5.0, 1), new Complex(-6.0, 1), new Complex(7.0, 1) }, { new Complex(-3.0, 1), new Complex(-6.0, 1), new Complex(8.0, 1), new Complex(9.0, 1) }, { new Complex(4.4, 1), new Complex(7.0, 1), new Complex(9.0, 1), new Complex(10.0, 1) } } }, + { "Singular4x4", new[,] { { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(0.0, 1), new Complex(4.0, 1) }, { new Complex(2.0, 1), new Complex(5.0, 1), new Complex(0.0, 1), new Complex(7.0, 1) }, { new Complex(0.0, 1), new Complex(0.0, 1), new Complex(0.0, 1), new Complex(0.0, 1) }, { new Complex(4.0, 1), new Complex(7.0, 1), new Complex(0.0, 1), new Complex(10.0, 1) } } }, + { "Tall3x2", new[,] { { new Complex(-1.1, 1), new Complex(-2.2, 1) }, { new Complex(0.0, 1), new Complex(1.1, 1) }, { new Complex(-4.4, 1), new Complex(5.5, 1) } } }, + { "Wide2x3", new[,] { { new Complex(-1.1, 1), new Complex(-2.2, 1), new Complex(-3.3, 1) }, { new Complex(0.0, 1), new Complex(1.1, 1), new Complex(2.2, 1) } } }, + { "Symmetric3x3", new[,] { { new Complex(1.0, 1), new Complex(2.0, 1), new Complex(3.0, 1) }, { new Complex(2.0, 1), new Complex(2.0, 1), new Complex(0.0, 1) }, { new Complex(3.0, 1), new Complex(0.0, 1), new Complex(3.0, 1) } } }, + { "NonSymmetric3x3", new[,] { { new Complex(-1.1, 1), new Complex(-2.2, 1), new Complex(-3.3, 1) }, { new Complex(0.0, 1), new Complex(1.1, 1), new Complex(2.2, 1) }, { new Complex(-4.4, 1), new Complex(5.5, 1), new Complex(6.6, 1) } } }, + { "NonSymmetric4x4", new[,] { { new Complex(-1.1, 1), new Complex(-2.2, 1), new Complex(-3.3, 1), new Complex(-4.4, 1) }, { new Complex(0.0, 1), new Complex(1.1, 1), new Complex(2.2, 1), new Complex(3.3, 1) }, { new Complex(1.0, 1), new Complex(2.1, 1), new Complex(6.2, 1), new Complex(4.3, 1) }, { new Complex(-4.4, 1), new Complex(5.5, 1), new Complex(6.6, 1), new Complex(-7.7, 1) } } }, + { "IndexTester4x4", new [,] { { new Complex(0, 1), new Complex(1, 1), new Complex(3, 1), new Complex(6, 1) }, { new Complex(1, 1), new Complex(2, 1), new Complex(4, 1), new Complex(7, 1) }, { new Complex(3, 1), new Complex(4, 1), new Complex(5, 1), new Complex(8, 1) }, { new Complex(6, 1), new Complex(7, 1), new Complex(8, 1), new Complex(9, 1) } } } + }; + + TestMatrices = new Dictionary(); + + foreach (var name in TestData2D.Keys) + { + TestMatrices.Add(name, CreateMatrix(TestData2D[name])); + } + } + + /// + /// Can add a non-symmetric matrix to this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanAddNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Add(matrixB); + + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] + matrixB[i, j]); + } + } + } + + /// + /// Can subtract a non-symmetric matrix from this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanSubtractNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Subtract(matrixB); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] - matrixB[i, j]); + } + } + } + + /// + /// Can compute Frobenius norm. + /// + public override void CanComputeFrobeniusNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + } + + + + /// + /// Can compute Infinity norm. + /// + public override void CanComputeInfinityNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + } + + + + /// + /// Can compute L1 norm. + /// + public override void CanComputeL1Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + } + + + + /// + /// Can compute L2 norm. + /// + public override void CanComputeL2Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Complex/SymmetricMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Complex/SymmetricMatrixTests.cs new file mode 100644 index 00000000..3121e9b8 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Complex/SymmetricMatrixTests.cs @@ -0,0 +1,124 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex +{ + using System.Numerics; + using MathNet.Numerics.LinearAlgebra.Complex; + using NUnit.Framework; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices. + /// + public abstract partial class SymmetricMatrixTests : MatrixTests + { + /// + /// Can check if a matrix is symmetric. + /// + [Test] + public override void CanCheckIfMatrixIsSymmetric() + { + var matrix = TestMatrices["Square3x3"]; + Assert.IsTrue(matrix.IsSymmetric); + + matrix = TestMatrices["NonSymmetric3x3"]; + Assert.IsFalse(matrix.IsSymmetric); + } + + /// + /// Can check if a [,] array is symmetric. + /// + [Test] + public void CanCheckIfArrayIsSymmetric() + { + Assert.IsTrue(SymmetricMatrix.CheckIfSymmetric(TestData2D["Square3x3"])); + Assert.IsFalse(SymmetricMatrix.CheckIfSymmetric(TestData2D["NonSymmetric3x3"])); + } + + /// + /// Test whether the index enumerator returns the correct values. + /// + [Test] + public void CanUseIndexedEnumerator() + { + var matrix = TestMatrices["Singular3x3"]; + var enumerator = matrix.IndexedEnumerator().GetEnumerator(); + enumerator.MoveNext(); + var item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(new Complex(1.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(new Complex(2.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(new Complex(3.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(new Complex(2.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(new Complex(0.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(new Complex(0.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(new Complex(3.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(new Complex(0.0, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(new Complex(0.0, 1), item.Item3); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricDenseMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricDenseMatrixTests.cs new file mode 100644 index 00000000..75747b33 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricDenseMatrixTests.cs @@ -0,0 +1,217 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 +{ + using System; + using System.Collections.Generic; + using Numerics; + + using MathNet.Numerics.LinearAlgebra.Complex32; + + using NUnit.Framework; + + /// + /// Symmetric Dense matrix tests. + /// + public class SymmetricDenseMatrixTests : SymmetricMatrixTests + { + /// + /// Creates a matrix for the given number of rows and columns. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// + /// A matrix with the given dimensions. + /// + protected override Matrix CreateMatrix(int rows, int columns) + { + return new DenseMatrix(rows, columns); + } + + /// + /// Creates a matrix from a 2D array. + /// + /// + /// The 2D array to create this matrix from. + /// + /// + /// A matrix with the given values. + /// + protected override Matrix CreateMatrix(Complex32[,] data) + { + if (SymmetricMatrix.CheckIfSymmetric(data)) + { + return new SymmetricDenseMatrix(data); + } + + return new DenseMatrix(data); + } + + /// + /// Creates a vector of the given size. + /// + /// + /// The size of the vector to create. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(int size) + { + return new DenseVector(size); + } + + /// + /// Creates a vector from an array. + /// + /// + /// The array to create this vector from. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(Complex32[] data) + { + return new DenseVector(data); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom1DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(3, new[] { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(0.0f, 1), new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) }) }, + { "Square3x3", new SymmetricDenseMatrix(3, new[] { new Complex32(-1.1f, 1), new Complex32(2.0f, 1), new Complex32(1.1f, 1), new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(6.6f, 1) }) }, + { "Square4x4", new SymmetricDenseMatrix(4, new[] { new Complex32(1.1f, 1), new Complex32(2.0f, 1), new Complex32(5.0f, 1), new Complex32(-3.0f, 1), new Complex32(-6.0f, 1), new Complex32(8.0f, 1), new Complex32(4.4f, 1), new Complex32(7.0f, 1), new Complex32(9.0f, 1), new Complex32(10.0f, 1) }) }, + { "Singular4x4", new SymmetricDenseMatrix(4, new[] { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(5.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1), new Complex32(4.0f, 1), new Complex32(7.0f, 1), new Complex32(0.0f, 1), new Complex32(10.0f, 1) }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(3, new[] { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(3.0f, 1) }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(4, new [] { new Complex32(0, 1), new Complex32(1, 1), new Complex32(2, 1), new Complex32(3, 1), new Complex32(4, 1), new Complex32(5, 1), new Complex32(6, 1), new Complex32(7, 1), new Complex32(8, 1), new Complex32(9, 1) }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from array is a reference. + /// + [Test] + public void MatrixFrom1DArrayIsReference() + { + var data = new Complex32[] { new Complex32(1, 1), new Complex32(1, 1), new Complex32(1, 1), new Complex32(1, 1), new Complex32(1, 1), new Complex32(1, 1) }; + var matrix = new SymmetricDenseMatrix(3, data); + matrix[0, 0] = new Complex32(10.0f, 2); + Assert.AreEqual(new Complex32(10.0f, 2), data[0]); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom2DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(new[,] { { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) }, { new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) } }) }, + { "Square3x3", new SymmetricDenseMatrix(new[,] { { new Complex32(-1.1f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(1.1f, 1), new Complex32(0.0f, 1) }, { new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(6.6f, 1) } }) }, + { "Square4x4", new SymmetricDenseMatrix(new[,] { { new Complex32(1.1f, 1), new Complex32(2.0f, 1), new Complex32(-3.0f, 1), new Complex32(4.4f, 1) }, { new Complex32(2.0f, 1), new Complex32(5.0f, 1), new Complex32(-6.0f, 1), new Complex32(7.0f, 1) }, { new Complex32(-3.0f, 1), new Complex32(-6.0f, 1), new Complex32(8.0f, 1), new Complex32(9.0f, 1) }, { new Complex32(4.4f, 1), new Complex32(7.0f, 1), new Complex32(9.0f, 1), new Complex32(10.0f, 1) } }) }, + { "Singular4x4", new SymmetricDenseMatrix(new[,] { { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(0.0f, 1), new Complex32(4.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(5.0f, 1), new Complex32(0.0f, 1), new Complex32(7.0f, 1) }, { new Complex32(0.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) }, { new Complex32(4.0f, 1), new Complex32(7.0f, 1), new Complex32(0.0f, 1), new Complex32(10.0f, 1) } }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(new[,] { { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(2.0f, 1), new Complex32(0.0f, 1) }, { new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(3.0f, 1) } }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(new [,] { { new Complex32(0, 1), new Complex32(1, 1), new Complex32(3, 1), new Complex32(6, 1) }, { new Complex32(1, 1), new Complex32(2, 1), new Complex32(4, 1), new Complex32(7, 1) }, { new Complex32(3, 1), new Complex32(4, 1), new Complex32(5, 1), new Complex32(8, 1) }, { new Complex32(6, 1), new Complex32(7, 1), new Complex32(8, 1), new Complex32(9, 1) } }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from two-dimensional array is a copy. + /// + [Test] + public void MatrixFrom2DArrayIsCopy() + { + var matrix = new DenseMatrix(TestData2D["Singular3x3"]); + matrix[0, 0] = new Complex32(10.0f, 2); + Assert.AreEqual(new Complex32(1.0f, 1), TestData2D["Singular3x3"][0, 0]); + } + + /// + /// Can create a matrix with uniform values. + /// + [Test] + public void CanCreateMatrixWithUniformValues() + { + var matrix = new SymmetricDenseMatrix(10, new Complex32(10.0f, 2)); + var value = new Complex32(10.0f, 2); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], value); + } + } + } + + /// + /// Can create an identity matrix. + /// + [Test] + public void CanCreateIdentity() + { + var matrix = SymmetricDenseMatrix.Identity(5); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(i == j ? Complex32.One : Complex32.Zero, matrix[i, j]); + } + } + } + + /// + /// Identity with wrong order throws ArgumentOutOfRangeException. + /// + /// The size of the square matrix + [TestCase(0)] + [TestCase(-1)] + public void IdentityWithWrongOrderThrowsArgumentOutOfRangeException(int order) + { + Assert.Throws(() => SymmetricDenseMatrix.Identity(order)); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricMatrixTests.Arithmetic.cs b/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricMatrixTests.Arithmetic.cs new file mode 100644 index 00000000..3d924a53 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricMatrixTests.Arithmetic.cs @@ -0,0 +1,192 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 +{ + using System.Collections.Generic; + using LinearAlgebra.Complex32; + using NUnit.Framework; + using Numerics; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices + /// + public abstract partial class SymmetricMatrixTests + { + /// + /// Setup test matrices. + /// Singular and Square matrices are overridden here with symmetric ones so that calls to base methods work as intended. + /// Additional NonSymmetric matrices are defined for some tests. + /// + [SetUp] + public override void SetupMatrices() + { + TestData2D = new Dictionary + { + { "Singular3x3", new[,] { { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) }, { new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) } } }, + { "Square3x3", new[,] { { new Complex32(-1.1f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(1.1f, 1), new Complex32(0.0f, 1) }, { new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(6.6f, 1) } } }, + { "Square4x4", new[,] { { new Complex32(1.1f, 1), new Complex32(2.0f, 1), new Complex32(-3.0f, 1), new Complex32(4.4f, 1) }, { new Complex32(2.0f, 1), new Complex32(5.0f, 1), new Complex32(-6.0f, 1), new Complex32(7.0f, 1) }, { new Complex32(-3.0f, 1), new Complex32(-6.0f, 1), new Complex32(8.0f, 1), new Complex32(9.0f, 1) }, { new Complex32(4.4f, 1), new Complex32(7.0f, 1), new Complex32(9.0f, 1), new Complex32(10.0f, 1) } } }, + { "Singular4x4", new[,] { { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(0.0f, 1), new Complex32(4.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(5.0f, 1), new Complex32(0.0f, 1), new Complex32(7.0f, 1) }, { new Complex32(0.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1), new Complex32(0.0f, 1) }, { new Complex32(4.0f, 1), new Complex32(7.0f, 1), new Complex32(0.0f, 1), new Complex32(10.0f, 1) } } }, + { "Tall3x2", new[,] { { new Complex32(-1.1f, 1), new Complex32(-2.2f, 1) }, { new Complex32(0.0f, 1), new Complex32(1.1f, 1) }, { new Complex32(-4.4f, 1), new Complex32(5.5f, 1) } } }, + { "Wide2x3", new[,] { { new Complex32(-1.1f, 1), new Complex32(-2.2f, 1), new Complex32(-3.3f, 1) }, { new Complex32(0.0f, 1), new Complex32(1.1f, 1), new Complex32(2.2f, 1) } } }, + { "Symmetric3x3", new[,] { { new Complex32(1.0f, 1), new Complex32(2.0f, 1), new Complex32(3.0f, 1) }, { new Complex32(2.0f, 1), new Complex32(2.0f, 1), new Complex32(0.0f, 1) }, { new Complex32(3.0f, 1), new Complex32(0.0f, 1), new Complex32(3.0f, 1) } } }, + { "NonSymmetric3x3", new[,] { { new Complex32(-1.1f, 1), new Complex32(-2.2f, 1), new Complex32(-3.3f, 1) }, { new Complex32(0.0f, 1), new Complex32(1.1f, 1), new Complex32(2.2f, 1) }, { new Complex32(-4.4f, 1), new Complex32(5.5f, 1), new Complex32(6.6f, 1) } } }, + { "NonSymmetric4x4", new[,] { { new Complex32(-1.1f, 1), new Complex32(-2.2f, 1), new Complex32(-3.3f, 1), new Complex32(-4.4f, 1) }, { new Complex32(0.0f, 1), new Complex32(1.1f, 1), new Complex32(2.2f, 1), new Complex32(3.3f, 1) }, { new Complex32(1.0f, 1), new Complex32(2.1f, 1), new Complex32(6.2f, 1), new Complex32(4.3f, 1) }, { new Complex32(-4.4f, 1), new Complex32(5.5f, 1), new Complex32(6.6f, 1), new Complex32(-7.7f, 1) } } }, + { "IndexTester4x4", new [,] { { new Complex32(0, 1), new Complex32(1, 1), new Complex32(3, 1), new Complex32(6, 1) }, { new Complex32(1, 1), new Complex32(2, 1), new Complex32(4, 1), new Complex32(7, 1) }, { new Complex32(3, 1), new Complex32(4, 1), new Complex32(5, 1), new Complex32(8, 1) }, { new Complex32(6, 1), new Complex32(7, 1), new Complex32(8, 1), new Complex32(9, 1) } } } + }; + + TestMatrices = new Dictionary(); + + foreach (var name in TestData2D.Keys) + { + TestMatrices.Add(name, CreateMatrix(TestData2D[name])); + } + } + + /// + /// Can add a non-symmetric matrix to this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanAddNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Add(matrixB); + + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] + matrixB[i, j]); + } + } + } + + /// + /// Can subtract a non-symmetric matrix from this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanSubtractNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Subtract(matrixB); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] - matrixB[i, j]); + } + } + } + + /// + /// Can compute Frobenius norm. + /// + public override void CanComputeFrobeniusNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + } + + + + /// + /// Can compute Infinity norm. + /// + public override void CanComputeInfinityNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + } + + + + /// + /// Can compute L1 norm. + /// + public override void CanComputeL1Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + } + + + + /// + /// Can compute L2 norm. + /// + public override void CanComputeL2Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricMatrixTests.cs new file mode 100644 index 00000000..03b4877d --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Complex32/SymmetricMatrixTests.cs @@ -0,0 +1,125 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 +{ + using MathNet.Numerics.LinearAlgebra.Complex32; + + using Numerics; + using NUnit.Framework; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices. + /// + public abstract partial class SymmetricMatrixTests : MatrixTests + { + /// + /// Can check if a matrix is symmetric. + /// + [Test] + public override void CanCheckIfMatrixIsSymmetric() + { + var matrix = TestMatrices["Square3x3"]; + Assert.IsTrue(matrix.IsSymmetric); + + matrix = TestMatrices["NonSymmetric3x3"]; + Assert.IsFalse(matrix.IsSymmetric); + } + + /// + /// Can check if a [,] array is symmetric. + /// + [Test] + public void CanCheckIfArrayIsSymmetric() + { + Assert.IsTrue(SymmetricMatrix.CheckIfSymmetric(TestData2D["Square3x3"])); + Assert.IsFalse(SymmetricMatrix.CheckIfSymmetric(TestData2D["NonSymmetric3x3"])); + } + + /// + /// Test whether the index enumerator returns the correct values. + /// + [Test] + public void CanUseIndexedEnumerator() + { + var matrix = TestMatrices["Singular3x3"]; + var enumerator = matrix.IndexedEnumerator().GetEnumerator(); + enumerator.MoveNext(); + var item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(new Complex32(1.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(new Complex32(2.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(new Complex32(3.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(new Complex32(2.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(new Complex32(0.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(new Complex32(0.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(new Complex32(3.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(new Complex32(0.0f, 1), item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(new Complex32(0.0f, 1), item.Item3); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Double/MatrixTests.Arithmetic.cs b/src/UnitTests/LinearAlgebraTests/Double/MatrixTests.Arithmetic.cs index b4354ebe..3287e4ae 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/MatrixTests.Arithmetic.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/MatrixTests.Arithmetic.cs @@ -948,10 +948,10 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double } /// - /// Can calculate Kronecker product. + /// Can calculate Kronecker product into a result matrix. /// [Test] - public void CanKroneckerProduct() + public void CanKroneckerProductIntoResult() { var matrixA = TestMatrices["Wide2x3"]; var matrixB = TestMatrices["Square3x3"]; @@ -972,11 +972,12 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double } } + /// - /// Can calculate Kronecker product into a result matrix. + /// Can calculate Kronecker product. /// [Test] - public void CanKroneckerProductIntoResult() + public void CanKroneckerProduct() { var matrixA = TestMatrices["Wide2x3"]; var matrixB = TestMatrices["Square3x3"]; diff --git a/src/UnitTests/LinearAlgebraTests/Double/SymmetricDenseMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Double/SymmetricDenseMatrixTests.cs new file mode 100644 index 00000000..b0078994 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Double/SymmetricDenseMatrixTests.cs @@ -0,0 +1,215 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double +{ + using System; + using System.Collections.Generic; + + using MathNet.Numerics.LinearAlgebra.Double; + + using NUnit.Framework; + + /// + /// Symmetric Dense matrix tests. + /// + public class SymmetricDenseMatrixTests : SymmetricMatrixTests + { + /// + /// Creates a matrix for the given number of rows and columns. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// + /// A matrix with the given dimensions. + /// + protected override Matrix CreateMatrix(int rows, int columns) + { + return new DenseMatrix(rows, columns); + } + + /// + /// Creates a matrix from a 2D array. + /// + /// + /// The 2D array to create this matrix from. + /// + /// + /// A matrix with the given values. + /// + protected override Matrix CreateMatrix(double[,] data) + { + if (SymmetricMatrix.CheckIfSymmetric(data)) + { + return new SymmetricDenseMatrix(data); + } + + return new DenseMatrix(data); + } + + /// + /// Creates a vector of the given size. + /// + /// + /// The size of the vector to create. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(int size) + { + return new DenseVector(size); + } + + /// + /// Creates a vector from an array. + /// + /// + /// The array to create this vector from. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(double[] data) + { + return new DenseVector(data); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom1DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(3, new[] { 1.0, 2.0, 0.0, 3.0, 0.0, 0.0 }) }, + { "Square3x3", new SymmetricDenseMatrix(3, new[] { -1.1, 2.0, 1.1, 3.0, 0.0, 6.6 }) }, + { "Square4x4", new SymmetricDenseMatrix(4, new[] { 1.1, 2.0, 5.0, -3.0, -6.0, 8.0, 4.4, 7.0, 9.0, 10.0 }) }, + { "Singular4x4", new SymmetricDenseMatrix(4, new[] { 1.0, 2.0, 5.0, 0.0, 0.0, 0.0, 4.0, 7.0, 0.0, 10.0 }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(3, new[] { 1.0, 2.0, 2.0, 3.0, 0.0, 3.0 }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(4, new double[] { 0, 1, 2, 3, 4, 5, 6, 7, 8, 9 }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from array is a reference. + /// + [Test] + public void MatrixFrom1DArrayIsReference() + { + var data = new double[] { 1, 1, 1, 1, 1, 1 }; + var matrix = new SymmetricDenseMatrix(3, data); + matrix[0, 0] = 10.0; + Assert.AreEqual(10.0, data[0]); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom2DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(new[,] { { 1.0, 2.0, 3.0 }, { 2.0, 0.0, 0.0 }, { 3.0, 0.0, 0.0 } }) }, + { "Square3x3", new SymmetricDenseMatrix(new[,] { { -1.1, 2.0, 3.0 }, { 2.0, 1.1, 0.0 }, { 3.0, 0.0, 6.6 } }) }, + { "Square4x4", new SymmetricDenseMatrix(new[,] { { 1.1, 2.0, -3.0, 4.4 }, { 2.0, 5.0, -6.0, 7.0 }, { -3.0, -6.0, 8.0, 9.0 }, { 4.4, 7.0, 9.0, 10.0 } }) }, + { "Singular4x4", new SymmetricDenseMatrix(new[,] { { 1.0, 2.0, 0.0, 4.0 }, { 2.0, 5.0, 0.0, 7.0 }, { 0.0, 0.0, 0.0, 0.0 }, { 4.0, 7.0, 0.0, 10.0 } }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(new[,] { { 1.0, 2.0, 3.0 }, { 2.0, 2.0, 0.0 }, { 3.0, 0.0, 3.0 } }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(new double[,] { { 0, 1, 3, 6 }, { 1, 2, 4, 7 }, { 3, 4, 5, 8 }, { 6, 7, 8, 9 } }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from two-dimensional array is a copy. + /// + [Test] + public void MatrixFrom2DArrayIsCopy() + { + var matrix = new DenseMatrix(TestData2D["Singular3x3"]); + matrix[0, 0] = 10.0; + Assert.AreEqual(1.0, TestData2D["Singular3x3"][0, 0]); + } + + /// + /// Can create a matrix with uniform values. + /// + [Test] + public void CanCreateMatrixWithUniformValues() + { + var matrix = new SymmetricDenseMatrix(10, 10.0); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], 10.0); + } + } + } + + /// + /// Can create an identity matrix. + /// + [Test] + public void CanCreateIdentity() + { + var matrix = SymmetricDenseMatrix.Identity(5); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(i == j ? 1.0 : 0.0, matrix[i, j]); + } + } + } + + /// + /// Identity with wrong order throws ArgumentOutOfRangeException. + /// + /// The size of the square matrix + [TestCase(0)] + [TestCase(-1)] + public void IdentityWithWrongOrderThrowsArgumentOutOfRangeException(int order) + { + Assert.Throws(() => SymmetricDenseMatrix.Identity(order)); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Double/SymmetricMatrixTests.Arithmetic.cs b/src/UnitTests/LinearAlgebraTests/Double/SymmetricMatrixTests.Arithmetic.cs new file mode 100644 index 00000000..cb7accb7 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Double/SymmetricMatrixTests.Arithmetic.cs @@ -0,0 +1,191 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double +{ + using System.Collections.Generic; + using LinearAlgebra.Double; + using NUnit.Framework; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices + /// + public abstract partial class SymmetricMatrixTests + { + /// + /// Setup test matrices. + /// Singular and Square matrices are overridden here with symmetric ones so that calls to base methods work as intended. + /// Additional NonSymmetric matrices are defined for some tests. + /// + [SetUp] + public override void SetupMatrices() + { + TestData2D = new Dictionary + { + { "Singular3x3", new[,] { { 1.0, 2.0, 3.0 }, { 2.0, 0.0, 0.0 }, { 3.0, 0.0, 0.0 } } }, + { "Square3x3", new[,] { { -1.1, 2.0, 3.0 }, { 2.0, 1.1, 0.0 }, { 3.0, 0.0, 6.6 } } }, + { "Square4x4", new[,] { { 1.1, 2.0, -3.0, 4.4 }, { 2.0, 5.0, -6.0, 7.0 }, { -3.0, -6.0, 8.0, 9.0 }, { 4.4, 7.0, 9.0, 10.0 } } }, + { "Singular4x4", new[,] { { 1.0, 2.0, 0.0, 4.0 }, { 2.0, 5.0, 0.0, 7.0 }, { 0.0, 0.0, 0.0, 0.0 }, { 4.0, 7.0, 0.0, 10.0 } } }, + { "Tall3x2", new[,] { { -1.1, -2.2 }, { 0.0, 1.1 }, { -4.4, 5.5 } } }, + { "Wide2x3", new[,] { { -1.1, -2.2, -3.3 }, { 0.0, 1.1, 2.2 } } }, + { "Symmetric3x3", new[,] { { 1.0, 2.0, 3.0 }, { 2.0, 2.0, 0.0 }, { 3.0, 0.0, 3.0 } } }, + { "NonSymmetric3x3", new[,] { { -1.1, -2.2, -3.3 }, { 0.0, 1.1, 2.2 }, { -4.4, 5.5, 6.6 } } }, + { "NonSymmetric4x4", new[,] { { -1.1, -2.2, -3.3, -4.4 }, { 0.0, 1.1, 2.2, 3.3 }, { 1.0, 2.1, 6.2, 4.3 }, { -4.4, 5.5, 6.6, -7.7 } } }, + { "IndexTester4x4", new double[,] { { 0, 1, 3, 6 }, { 1, 2, 4, 7 }, { 3, 4, 5, 8 }, { 6, 7, 8, 9 } } } + }; + + TestMatrices = new Dictionary(); + + foreach (var name in TestData2D.Keys) + { + TestMatrices.Add(name, CreateMatrix(TestData2D[name])); + } + } + + /// + /// Can add a non-symmetric matrix to this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanAddNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Add(matrixB); + + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] + matrixB[i, j]); + } + } + } + + /// + /// Can subtract a non-symmetric matrix from this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanSubtractNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Subtract(matrixB); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] - matrixB[i, j]); + } + } + } + + /// + /// Can compute Frobenius norm. + /// + public override void CanComputeFrobeniusNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + } + + + + /// + /// Can compute Infinity norm. + /// + public override void CanComputeInfinityNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + } + + + + /// + /// Can compute L1 norm. + /// + public override void CanComputeL1Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + } + + + + /// + /// Can compute L2 norm. + /// + public override void CanComputeL2Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Double/SymmetricMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Double/SymmetricMatrixTests.cs new file mode 100644 index 00000000..d95359f0 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Double/SymmetricMatrixTests.cs @@ -0,0 +1,123 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double +{ + using MathNet.Numerics.LinearAlgebra.Double; + using NUnit.Framework; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices. + /// + public abstract partial class SymmetricMatrixTests : MatrixTests + { + /// + /// Can check if a matrix is symmetric. + /// + [Test] + public override void CanCheckIfMatrixIsSymmetric() + { + var matrix = TestMatrices["Square3x3"]; + Assert.IsTrue(matrix.IsSymmetric); + + matrix = TestMatrices["NonSymmetric3x3"]; + Assert.IsFalse(matrix.IsSymmetric); + } + + /// + /// Can check if a [,] array is symmetric. + /// + [Test] + public void CanCheckIfArrayIsSymmetric() + { + Assert.IsTrue(SymmetricMatrix.CheckIfSymmetric(TestData2D["Square3x3"])); + Assert.IsFalse(SymmetricMatrix.CheckIfSymmetric(TestData2D["NonSymmetric3x3"])); + } + + /// + /// Test whether the index enumerator returns the correct values. + /// + [Test] + public void CanUseIndexedEnumerator() + { + var matrix = TestMatrices["Singular3x3"]; + var enumerator = matrix.IndexedEnumerator().GetEnumerator(); + enumerator.MoveNext(); + var item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(1.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(2.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(3.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(2.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(0.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(0.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(3.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(0.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(0.0, item.Item3); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Single/SymmetricDenseMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Single/SymmetricDenseMatrixTests.cs new file mode 100644 index 00000000..5ec190c2 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Single/SymmetricDenseMatrixTests.cs @@ -0,0 +1,215 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single +{ + using System; + using System.Collections.Generic; + + using MathNet.Numerics.LinearAlgebra.Single; + + using NUnit.Framework; + + /// + /// Symmetric Dense matrix tests. + /// + public class SymmetricDenseMatrixTests : SymmetricMatrixTests + { + /// + /// Creates a matrix for the given number of rows and columns. + /// + /// + /// The number of rows. + /// + /// + /// The number of columns. + /// + /// + /// A matrix with the given dimensions. + /// + protected override Matrix CreateMatrix(int rows, int columns) + { + return new DenseMatrix(rows, columns); + } + + /// + /// Creates a matrix from a 2D array. + /// + /// + /// The 2D array to create this matrix from. + /// + /// + /// A matrix with the given values. + /// + protected override Matrix CreateMatrix(float[,] data) + { + if (SymmetricMatrix.CheckIfSymmetric(data)) + { + return new SymmetricDenseMatrix(data); + } + + return new DenseMatrix(data); + } + + /// + /// Creates a vector of the given size. + /// + /// + /// The size of the vector to create. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(int size) + { + return new DenseVector(size); + } + + /// + /// Creates a vector from an array. + /// + /// + /// The array to create this vector from. + /// + /// + /// The new vector. + /// + protected override Vector CreateVector(float[] data) + { + return new DenseVector(data); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom1DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(3, new[] { 1.0f, 2.0f, 0.0f, 3.0f, 0.0f, 0.0f }) }, + { "Square3x3", new SymmetricDenseMatrix(3, new[] { -1.1f, 2.0f, 1.1f, 3.0f, 0.0f, 6.6f }) }, + { "Square4x4", new SymmetricDenseMatrix(4, new[] { 1.1f, 2.0f, 5.0f, -3.0f, -6.0f, 8.0f, 4.4f, 7.0f, 9.0f, 10.0f }) }, + { "Singular4x4", new SymmetricDenseMatrix(4, new[] { 1.0f, 2.0f, 5.0f, 0.0f, 0.0f, 0.0f, 4.0f, 7.0f, 0.0f, 10.0f }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(3, new[] { 1.0f, 2.0f, 2.0f, 3.0f, 0.0f, 3.0f }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(4, new [] { 0f, 1f, 2f, 3f, 4f, 5f, 6f, 7f, 8f, 9f }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from array is a reference. + /// + [Test] + public void MatrixFrom1DArrayIsReference() + { + var data = new float[] { 1, 1, 1, 1, 1, 1 }; + var matrix = new SymmetricDenseMatrix(3, data); + matrix[0, 0] = 10.0f; + Assert.AreEqual(10.0f, data[0]); + } + + /// + /// Can create a matrix form array. + /// + [Test] + public void CanCreateMatrixFrom2DArray() + { + var testData = new Dictionary + { + { "Singular3x3", new SymmetricDenseMatrix(new[,] { { 1.0f, 2.0f, 3.0f }, { 2.0f, 0.0f, 0.0f }, { 3.0f, 0.0f, 0.0f } }) }, + { "Square3x3", new SymmetricDenseMatrix(new[,] { { -1.1f, 2.0f, 3.0f }, { 2.0f, 1.1f, 0.0f }, { 3.0f, 0.0f, 6.6f } }) }, + { "Square4x4", new SymmetricDenseMatrix(new[,] { { 1.1f, 2.0f, -3.0f, 4.4f }, { 2.0f, 5.0f, -6.0f, 7.0f }, { -3.0f, -6.0f, 8.0f, 9.0f }, { 4.4f, 7.0f, 9.0f, 10.0f } }) }, + { "Singular4x4", new SymmetricDenseMatrix(new[,] { { 1.0f, 2.0f, 0.0f, 4.0f }, { 2.0f, 5.0f, 0.0f, 7.0f }, { 0.0f, 0.0f, 0.0f, 0.0f }, { 4.0f, 7.0f, 0.0f, 10.0f } }) }, + { "Symmetric3x3", new SymmetricDenseMatrix(new[,] { { 1.0f, 2.0f, 3.0f }, { 2.0f, 2.0f, 0.0f }, { 3.0f, 0.0f, 3.0f } }) }, + { "IndexTester4x4", new SymmetricDenseMatrix(new [,] { { 0f, 1f, 3f, 6f }, { 1f, 2f, 4f, 7f }, { 3f, 4f, 5f, 8f }, { 6f, 7f, 8f, 9f } }) } + }; + + foreach (var name in testData.Keys) + { + Assert.AreEqual(TestMatrices[name], testData[name]); + } + } + + /// + /// Matrix from two-dimensional array is a copy. + /// + [Test] + public void MatrixFrom2DArrayIsCopy() + { + var matrix = new DenseMatrix(TestData2D["Singular3x3"]); + matrix[0, 0] = 10.0f; + Assert.AreEqual(1.0f, TestData2D["Singular3x3"][0, 0]); + } + + /// + /// Can create a matrix with uniform values. + /// + [Test] + public void CanCreateMatrixWithUniformValues() + { + var matrix = new SymmetricDenseMatrix(10, 10.0f); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], 10.0f); + } + } + } + + /// + /// Can create an identity matrix. + /// + [Test] + public void CanCreateIdentity() + { + var matrix = SymmetricDenseMatrix.Identity(5); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(i == j ? 1.0f : 0.0f, matrix[i, j]); + } + } + } + + /// + /// Identity with wrong order throws ArgumentOutOfRangeException. + /// + /// The size of the square matrix + [TestCase(0)] + [TestCase(-1)] + public void IdentityWithWrongOrderThrowsArgumentOutOfRangeException(int order) + { + Assert.Throws(() => SymmetricDenseMatrix.Identity(order)); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Single/SymmetricMatrixTests.Arithmetic.cs b/src/UnitTests/LinearAlgebraTests/Single/SymmetricMatrixTests.Arithmetic.cs new file mode 100644 index 00000000..331ac876 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Single/SymmetricMatrixTests.Arithmetic.cs @@ -0,0 +1,191 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single +{ + using System.Collections.Generic; + using LinearAlgebra.Single; + using NUnit.Framework; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices + /// + public abstract partial class SymmetricMatrixTests + { + /// + /// Setup test matrices. + /// Singular and Square matrices are overridden here with symmetric ones so that calls to base methods work as intended. + /// Additional NonSymmetric matrices are defined for some tests. + /// + [SetUp] + public override void SetupMatrices() + { + TestData2D = new Dictionary + { + { "Singular3x3", new[,] { { 1.0f, 2.0f, 3.0f }, { 2.0f, 0.0f, 0.0f }, { 3.0f, 0.0f, 0.0f } } }, + { "Square3x3", new[,] { { -1.1f, 2.0f, 3.0f }, { 2.0f, 1.1f, 0.0f }, { 3.0f, 0.0f, 6.6f } } }, + { "Square4x4", new[,] { { 1.1f, 2.0f, -3.0f, 4.4f }, { 2.0f, 5.0f, -6.0f, 7.0f }, { -3.0f, -6.0f, 8.0f, 9.0f }, { 4.4f, 7.0f, 9.0f, 10.0f } } }, + { "Singular4x4", new[,] { { 1.0f, 2.0f, 0.0f, 4.0f }, { 2.0f, 5.0f, 0.0f, 7.0f }, { 0.0f, 0.0f, 0.0f, 0.0f }, { 4.0f, 7.0f, 0.0f, 10.0f } } }, + { "Tall3x2", new[,] { { -1.1f, -2.2f }, { 0.0f, 1.1f }, { -4.4f, 5.5f } } }, + { "Wide2x3", new[,] { { -1.1f, -2.2f, -3.3f }, { 0.0f, 1.1f, 2.2f } } }, + { "Symmetric3x3", new[,] { { 1.0f, 2.0f, 3.0f }, { 2.0f, 2.0f, 0.0f }, { 3.0f, 0.0f, 3.0f } } }, + { "NonSymmetric3x3", new[,] { { -1.1f, -2.2f, -3.3f }, { 0.0f, 1.1f, 2.2f }, { -4.4f, 5.5f, 6.6f } } }, + { "NonSymmetric4x4", new[,] { { -1.1f, -2.2f, -3.3f, -4.4f }, { 0.0f, 1.1f, 2.2f, 3.3f }, { 1.0f, 2.1f, 6.2f, 4.3f }, { -4.4f, 5.5f, 6.6f, -7.7f } } }, + { "IndexTester4x4", new [,] { { 0f, 1f, 3f, 6f }, { 1f, 2f, 4f, 7f }, { 3f, 4f, 5f, 8f }, { 6f, 7f, 8f, 9f } } } + }; + + TestMatrices = new Dictionary(); + + foreach (var name in TestData2D.Keys) + { + TestMatrices.Add(name, CreateMatrix(TestData2D[name])); + } + } + + /// + /// Can add a non-symmetric matrix to this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanAddNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Add(matrixB); + + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] + matrixB[i, j]); + } + } + } + + /// + /// Can subtract a non-symmetric matrix from this symmetric matrix. + /// + /// Matrix A name. + /// Matrix B name. + [Test, Sequential] + public void CanSubtractNonSymmetricMatrix([Values("Square3x3", "Square4x4")] string mtxA, [Values("NonSymmetric3x3", "NonSymmetric4x4")] string mtxB) + { + var matrixA = TestMatrices[mtxA]; + var matrixB = TestMatrices[mtxB]; + + var matrix = matrixA.Clone(); + matrix = matrix.Subtract(matrixB); + for (var i = 0; i < matrix.RowCount; i++) + { + for (var j = 0; j < matrix.ColumnCount; j++) + { + Assert.AreEqual(matrix[i, j], matrixA[i, j] - matrixB[i, j]); + } + } + } + + /// + /// Can compute Frobenius norm. + /// + public override void CanComputeFrobeniusNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.FrobeniusNorm(), matrix.FrobeniusNorm(), 14); + } + + + + /// + /// Can compute Infinity norm. + /// + public override void CanComputeInfinityNorm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.InfinityNorm(), matrix.InfinityNorm(), 14); + } + + + + /// + /// Can compute L1 norm. + /// + public override void CanComputeL1Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L1Norm(), matrix.L1Norm(), 14); + } + + + + /// + /// Can compute L2 norm. + /// + public override void CanComputeL2Norm() + { + var matrix = TestMatrices["Square3x3"]; + var denseMatrix = new DenseMatrix(TestData2D["Square3x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Wide2x3"]; + denseMatrix = new DenseMatrix(TestData2D["Wide2x3"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + + matrix = TestMatrices["Tall3x2"]; + denseMatrix = new DenseMatrix(TestData2D["Tall3x2"]); + AssertHelpers.AlmostEqual(denseMatrix.L2Norm(), matrix.L2Norm(), 14); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/LinearAlgebraTests/Single/SymmetricMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Single/SymmetricMatrixTests.cs new file mode 100644 index 00000000..a0b9d4e5 --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/Single/SymmetricMatrixTests.cs @@ -0,0 +1,123 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// Copyright (c) 2009-2010 Math.NET +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single +{ + using MathNet.Numerics.LinearAlgebra.Single; + using NUnit.Framework; + + /// + /// Abstract class with the common set of matrix tests for symmetric matrices. + /// + public abstract partial class SymmetricMatrixTests : MatrixTests + { + /// + /// Can check if a matrix is symmetric. + /// + [Test] + public override void CanCheckIfMatrixIsSymmetric() + { + var matrix = TestMatrices["Square3x3"]; + Assert.IsTrue(matrix.IsSymmetric); + + matrix = TestMatrices["NonSymmetric3x3"]; + Assert.IsFalse(matrix.IsSymmetric); + } + + /// + /// Can check if a [,] array is symmetric. + /// + [Test] + public void CanCheckIfArrayIsSymmetric() + { + Assert.IsTrue(SymmetricMatrix.CheckIfSymmetric(TestData2D["Square3x3"])); + Assert.IsFalse(SymmetricMatrix.CheckIfSymmetric(TestData2D["NonSymmetric3x3"])); + } + + /// + /// Test whether the index enumerator returns the correct values. + /// + [Test] + public void CanUseIndexedEnumerator() + { + var matrix = TestMatrices["Singular3x3"]; + var enumerator = matrix.IndexedEnumerator().GetEnumerator(); + enumerator.MoveNext(); + var item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(1.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(2.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(0, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(3.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(2.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(0.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(1, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(0.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(0, item.Item2); + Assert.AreEqual(3.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(1, item.Item2); + Assert.AreEqual(0.0, item.Item3); + + enumerator.MoveNext(); + item = enumerator.Current; + Assert.AreEqual(2, item.Item1); + Assert.AreEqual(2, item.Item2); + Assert.AreEqual(0.0, item.Item3); + } + } +} \ No newline at end of file diff --git a/src/UnitTests/UnitTests.csproj b/src/UnitTests/UnitTests.csproj index 06f03500..7bb0cd41 100644 --- a/src/UnitTests/UnitTests.csproj +++ b/src/UnitTests/UnitTests.csproj @@ -272,6 +272,9 @@ Code + + + Code @@ -420,6 +423,9 @@ Code + + + Code @@ -568,6 +574,9 @@ Code + + + Code @@ -718,6 +727,9 @@ Code + + + Code