diff --git a/src/Examples/LinearAlgebra/VectorDataAccessor.cs b/src/Examples/LinearAlgebra/VectorDataAccessor.cs index 76fa8b40..9f7030a6 100644 --- a/src/Examples/LinearAlgebra/VectorDataAccessor.cs +++ b/src/Examples/LinearAlgebra/VectorDataAccessor.cs @@ -30,6 +30,8 @@ using MathNet.Numerics.LinearAlgebra.Double; namespace Examples.LinearAlgebraExamples { + using MathNet.Numerics.LinearAlgebra.Generic; + /// /// Vector data access, copying and conversion examples /// @@ -124,7 +126,7 @@ namespace Examples.LinearAlgebraExamples Console.WriteLine(); // 8. Copy part of vector into another vector. If you need to copy all data then use CopoTy(vector) method. - vectorA.CopyTo(cloneA, 3, 3, 4); + vectorA.CopySubVectorTo(cloneA, 3, 3, 4); Console.WriteLine(@"8. Copy part of vector into another vector"); Console.WriteLine(cloneA.ToString("#0.00\t", formatProvider)); Console.WriteLine(); diff --git a/src/Numerics.IO/LinearAlgebra/IO/MatlabWriter.cs b/src/Numerics.IO/LinearAlgebra/IO/MatlabWriter.cs index 815c0669..8120a79e 100644 --- a/src/Numerics.IO/LinearAlgebra/IO/MatlabWriter.cs +++ b/src/Numerics.IO/LinearAlgebra/IO/MatlabWriter.cs @@ -37,6 +37,7 @@ namespace MathNet.Numerics.LinearAlgebra.IO using Generic; using Matlab; using Properties; + using Storage; /// /// Writes matrices to a Matlab file. @@ -500,7 +501,7 @@ namespace MathNet.Numerics.LinearAlgebra.IO int count = 0; foreach (var column in matrix.ColumnEnumerator()) { - count += ((Double.SparseVector)column.Item2).NonZerosCount; + count += ((SparseVectorStorage)column.Item2.Storage).ValueCount; dataWriter.Write(count); } @@ -568,7 +569,7 @@ namespace MathNet.Numerics.LinearAlgebra.IO int count = 0; foreach (var column in matrix.ColumnEnumerator()) { - count += ((Single.SparseVector)column.Item2).NonZerosCount; + count += ((SparseVectorStorage)column.Item2.Storage).ValueCount; dataWriter.Write(count); } @@ -639,7 +640,7 @@ namespace MathNet.Numerics.LinearAlgebra.IO int count = 0; foreach (var column in matrix.ColumnEnumerator()) { - count += ((Complex.SparseVector)column.Item2).NonZerosCount; + count += ((SparseVectorStorage)column.Item2.Storage).ValueCount; dataWriter.Write(count); } @@ -718,7 +719,7 @@ namespace MathNet.Numerics.LinearAlgebra.IO int count = 0; foreach (var column in matrix.ColumnEnumerator()) { - count += ((Complex32.SparseVector)column.Item2).NonZerosCount; + count += ((SparseVectorStorage)column.Item2.Storage).ValueCount; dataWriter.Write(count); } diff --git a/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs index df43af87..89a51e71 100644 --- a/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/DenseMatrix.cs @@ -39,8 +39,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex [Serializable] public class DenseMatrix : Matrix { - readonly DenseColumnMajorMatrixStorage _storage; - /// /// Number of rows. /// @@ -59,15 +57,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// Gets the matrix's data. /// /// The matrix's data. - readonly Complex[] _data; + readonly Complex[] _values; internal DenseMatrix(DenseColumnMajorMatrixStorage storage) : base(storage) { - _storage = storage; - _rowCount = _storage.RowCount; - _columnCount = _storage.ColumnCount; - _data = _storage.Data; + _rowCount = storage.RowCount; + _columnCount = storage.ColumnCount; + _values = storage.Data; } /// @@ -109,9 +106,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public DenseMatrix(int rows, int columns, Complex value) : this(rows, columns) { - for (var i = 0; i < _data.Length; i++) + for (var i = 0; i < _values.Length; i++) { - _data[i] = value; + _values[i] = value; } } @@ -139,7 +136,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { for (var j = 0; j < _columnCount; j++) { - _data[(j * _rowCount) + i] = array[i, j]; + _values[(j * _rowCount) + i] = array[i, j]; } } } @@ -161,7 +158,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The matrix's data. public Complex[] Data { - get { return _data; } + get { return _values; } } /// @@ -203,7 +200,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var index = j * _rowCount; for (var i = 0; i < _rowCount; i++) { - ret._data[(i * _columnCount) + j] = _data[index + i]; + ret._values[(i * _columnCount) + j] = _values[index + i]; } } @@ -214,21 +211,21 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The L1 norm of the matrix. public override Complex L1Norm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _values); } /// Calculates the Frobenius norm of this matrix. /// The Frobenius norm of this matrix. public override Complex FrobeniusNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _values); } /// Calculates the infinity norm of this matrix. /// The infinity norm of this matrix. public override Complex InfinityNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _values); } #region Static constructors for special matrices. @@ -246,7 +243,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var m = new DenseMatrix(order); for (var i = 0; i < order; i++) { - m._data[(i * order) + i] = 1.0; + m._values[(i * order) + i] = 1.0; } return m; @@ -268,7 +265,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult._values); } } @@ -292,14 +289,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0, - denseResult.Data); + denseResult.Values); } } @@ -323,14 +320,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0, - denseResult._data); + denseResult._values); } } @@ -354,14 +351,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.Transpose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0, - denseResult._data); + denseResult._values); } } @@ -385,14 +382,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0, - denseResult.Data); + denseResult.Values); } } @@ -416,14 +413,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0, - denseResult._data); + denseResult._values); } } @@ -441,7 +438,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(-1, _values, denseResult._values); } } @@ -461,7 +458,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_values, denseOther._values, denseResult._values); } } @@ -481,7 +478,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseDivideArrays(_values, denseOther._values, denseResult._values); } } @@ -502,7 +499,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.AddArrays(_values, denseOther._values, denseResult._values); } } @@ -521,7 +518,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.SubtractArrays(_values, denseOther._values, denseResult._values); } } @@ -537,7 +534,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var index = j * _rowCount; for (var i = 0; i < _rowCount; i++) { - ret._data[(i * _columnCount) + j] = _data[index + i].Conjugate(); + ret._values[(i * _columnCount) + j] = _values[index + i].Conjugate(); } } @@ -559,7 +556,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var sum = Complex.Zero; for (var i = 0; i < _rowCount; i++) { - sum += _data[(i * _rowCount) + i]; + sum += _values[(i * _rowCount) + i]; } return sum; diff --git a/src/Numerics/LinearAlgebra/Complex/DenseVector.cs b/src/Numerics/LinearAlgebra/Complex/DenseVector.cs index eab389ba..9562ced0 100644 --- a/src/Numerics/LinearAlgebra/Complex/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex/DenseVector.cs @@ -33,6 +33,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex using Generic; using NumberTheory; using Properties; + using Storage; using Threading; /// @@ -41,6 +42,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex [Serializable] public class DenseVector : Vector { + /// + /// Number of elements + /// + readonly int _length; + + /// + /// Gets the vector's data. + /// + readonly Complex[] _values; + + internal DenseVector(DenseVectorStorage storage) + : base(storage) + { + _length = storage.Length; + _values = storage.Data; + } + /// /// Initializes a new instance of the class with a given size. /// @@ -50,9 +68,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// If is less than one. /// - public DenseVector(int size) : base(size) + public DenseVector(int size) + : this(new DenseVectorStorage(size)) { - Data = new Complex[size]; } /// @@ -68,11 +86,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// If is less than one. /// - public DenseVector(int size, Complex value) : this(size) + public DenseVector(int size, Complex value) + : this(size) { - for (var index = 0; index < Data.Length; index++) + for (var index = 0; index < _values.Length; index++) { - Data[index] = value; + _values[index] = value; } } @@ -83,27 +102,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// The vector to create the new vector from. /// - public DenseVector(Vector other) : this(other.Count) - { - CommonParallel.For( - 0, - Data.Length, - index => this[index] = other[index]); - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public DenseVector(DenseVector other) : this(other.Count) + public DenseVector(Vector other) + : this(other.Count) { - CommonParallel.For( - 0, - Data.Length, - index => Data[index] = other.Data[index]); + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -112,20 +114,18 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The array to create this vector from. /// The vector does not copy the array, but keeps a reference to it. Any /// changes to the vector will also change the array. - public DenseVector(Complex[] array) : base(array.Length) + public DenseVector(Complex[] array) + : this(new DenseVectorStorage(array.Length, array)) { - Data = array; } /// - /// Gets the vector's internal data. + /// Gets the vector's data. /// - /// The vector's internal data. - /// Changing values in the array also changes the corresponding value in vector. Use with care. - internal Complex[] Data + /// The vector's data. + public Complex[] Values { - get; - private set; + get { return _values; } } /// @@ -136,14 +136,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// A reference to the internal date of the given vector. /// - public static implicit operator Complex[](DenseVector vector) + public static explicit operator Complex[](DenseVector vector) { if (vector == null) { throw new ArgumentNullException(); } - return vector.Data; + return vector.Values; } /// @@ -169,10 +169,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// This vector as a column matrix. public override Matrix ToColumnMatrix() { - var matrix = new DenseMatrix(Count, 1); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(_length, 1); + for (var i = 0; i < _values.Length; i++) { - matrix.At(i, 0, Data[i]); + matrix.At(i, 0, _values[i]); } return matrix; @@ -184,33 +184,15 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// This vector as a row matrix. public override Matrix ToRowMatrix() { - var matrix = new DenseMatrix(1, Count); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(1, _length); + for (var i = 0; i < _values.Length; i++) { - matrix.At(0, i, Data[i]); + matrix.At(0, i, _values[i]); } return matrix; } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - /// If is negative or - /// greater than the size of the vector. - public override Complex this[int index] - { - get - { - return Data[index]; - } - - set - { - Data[index] = value; - } - } - /// /// Creates a matrix with the given dimensions using the same storage type /// as this vector. @@ -244,49 +226,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex return new DenseVector(size); } - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as DenseVector; - if (otherVector == null) - { - CommonParallel.For( - 0, - Data.Length, - index => target[index] = Data[index]); - } - else - { - Array.Copy(Data, 0, otherVector.Data, 0, Data.Length); - } - } - /// /// Adds a scalar to each element of the vector and stores the result in the result vector. /// @@ -303,8 +242,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] + scalar); + _values.Length, + index => dense._values[index] = _values[index] + scalar); } } @@ -319,7 +258,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, Complex.One, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, Complex.One, odense.Values, rdense.Values); } else { @@ -388,8 +327,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] - scalar); + _values.Length, + index => dense._values[index] = _values[index] - scalar); } } @@ -404,7 +343,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, -1.0, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, -1.0, odense.Values, rdense.Values); } else { @@ -463,11 +402,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new DenseVector(Count); + var result = new DenseVector(_length); CommonParallel.For( 0, - Data.Length, - index => result[index] = -Data[index]); + _values.Length, + index => result[index] = -_values[index]); return result; } @@ -487,7 +426,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, Data, denseResult.Data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult.Values); } } @@ -501,7 +440,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { var denseVector = other as DenseVector; - return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(Data, denseVector.Data); + return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(_values, denseVector.Values); } /// @@ -563,7 +502,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide"); } - return Control.LinearAlgebraProvider.DotProduct(leftSide.Data, rightSide.Data); + return Control.LinearAlgebraProvider.DotProduct(leftSide.Values, rightSide.Values); } /// @@ -590,10 +529,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override int AbsoluteMinimumIndex() { var index = 0; - var min = Data[index].Magnitude; - for (var i = 1; i < Count; i++) + var min = _values[index].Magnitude; + for (var i = 1; i < _length; i++) { - var test = Data[i].Magnitude; + var test = _values[i].Magnitude; if (test < min) { index = i; @@ -610,7 +549,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The value of the absolute minimum element. public override Complex AbsoluteMinimum() { - return Data[AbsoluteMinimumIndex()].Magnitude; + return _values[AbsoluteMinimumIndex()].Magnitude; } /// @@ -619,7 +558,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The value of the absolute maximum element. public override Complex AbsoluteMaximum() { - return Data[AbsoluteMaximumIndex()].Magnitude; + return _values[AbsoluteMaximumIndex()].Magnitude; } /// @@ -629,10 +568,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override int AbsoluteMaximumIndex() { var index = 0; - var max = Data[index].Magnitude; - for (var i = 1; i < Count; i++) + var max = _values[index].Magnitude; + for (var i = 1; i < _length; i++) { - var test = Data[i].Magnitude; + var test = _values[i].Magnitude; if (test > max) { index = i; @@ -656,7 +595,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// If is not positive. public override Vector SubVector(int index, int length) { - if (index < 0 || index >= Count) + if (index < 0 || index >= _length) { throw new ArgumentOutOfRangeException("index"); } @@ -666,7 +605,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex throw new ArgumentOutOfRangeException("length"); } - if (index + length > Count) + if (index + length > _length) { throw new ArgumentOutOfRangeException("length"); } @@ -675,8 +614,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex CommonParallel.For( index, - index + length, - i => result.Data[i - index] = Data[i]); + index + length, + i => result._values[i - index] = _values[i]); return result; } @@ -693,7 +632,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex throw new ArgumentNullException("values"); } - if (values.Length != Count) + if (values.Length != _length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "values"); } @@ -701,7 +640,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex CommonParallel.For( 0, values.Length, - i => Data[i] = values[i]); + i => _values[i] = values[i]); } /// @@ -712,9 +651,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { var sum = Complex.Zero; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Data[i]; + sum += _values[i]; } return sum; @@ -728,9 +667,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { var sum = Complex.Zero; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Data[i].Magnitude; + sum += _values[i].Magnitude; } return sum; @@ -752,8 +691,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] * other[index]); + _values.Length, + index => dense._values[index] = _values[index] * other[index]); } } @@ -774,8 +713,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] / other[index]); + _values.Length, + index => dense._values[index] = _values[index] / other[index]); } } @@ -807,7 +746,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { for (var j = 0; j < v.Count; j++) { - matrix.At(i, j, u.Data[i] * v.Data[j]); + matrix.At(i, j, u._values[i] * v._values[j]); } }); return matrix; @@ -845,19 +784,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex if (2.0 == p) { - return Data.Aggregate(Complex.Zero, SpecialFunctions.Hypotenuse).Magnitude; + return _values.Aggregate(Complex.Zero, SpecialFunctions.Hypotenuse).Magnitude; } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(Data, (i, v) => v.Magnitude, Math.Max, 0d); + return CommonParallel.Aggregate(_values, (i, v) => v.Magnitude, Math.Max, 0d); } var sum = 0.0; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Math.Pow(Data[i].Magnitude, p); + sum += Math.Pow(_values[i].Magnitude, p); } return Math.Pow(sum, 1.0 / p); @@ -1022,14 +961,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex #endregion - /// - /// Resets all values to zero. - /// - public override void Clear() - { - Array.Clear(Data, 0, Data.Length); - } - /// /// Conjugates vector and save result to /// @@ -1046,24 +977,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { CommonParallel.For( 0, - Count, - index => denseTarget.Data[index] = Data[index].Conjugate()); + _length, + index => denseTarget._values[index] = _values[index].Conjugate()); } } - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - internal protected override Complex At(int index) - { - return Data[index]; - } - - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - internal protected override void At(int index, Complex value) - { - Data[index] = value; - } } } diff --git a/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs index 6991fd30..ffb81b0e 100644 --- a/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/DiagonalMatrix.cs @@ -378,12 +378,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex return; } - if (_data.Length != denseSource.Data.Length) + if (_data.Length != denseSource.Values.Length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "source"); } - Array.Copy(denseSource.Data, _data, denseSource.Data.Length); + Array.Copy(denseSource.Values, _data, denseSource.Values.Length); } /// diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseCholesky.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseCholesky.cs index fdf78d58..ad8faee8 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseCholesky.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseCholesky.cs @@ -169,11 +169,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization } // Copy the contents of input to result. - CommonParallel.For(0, dinput.Data.Length, index => dresult.Data[index] = dinput.Data[index]); + CommonParallel.For(0, dinput.Values.Length, index => dresult.Values[index] = dinput.Values[index]); // Cholesky solve by overwriting result. var dfactor = (DenseMatrix)CholeskyFactor; - Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Data, 1); + Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Values, 1); } } } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseEvd.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseEvd.cs index 01bd6350..1bf067dd 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseEvd.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseEvd.cs @@ -107,7 +107,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization { var matrixH = matrix.ToArray(); NonsymmetricReduceToHessenberg(((DenseMatrix)MatrixEv).Data, matrixH, order); - NonsymmetricReduceHessenberToRealSchur(((DenseVector)VectorEv).Data, ((DenseMatrix)MatrixEv).Data, matrixH, order); + NonsymmetricReduceHessenberToRealSchur(((DenseVector)VectorEv).Values, ((DenseMatrix)MatrixEv).Data, matrixH, order); } MatrixD.SetDiagonal(VectorEv); diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs index eebf5f39..387ceba3 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseGramSchmidt.cs @@ -217,7 +217,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Data, 1, dresult.Data); + _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseLU.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseLU.cs index 49e3c1d2..bc35e77a 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseLU.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseLU.cs @@ -170,11 +170,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization } // Copy the contents of input to result. - CommonParallel.For(0, dinput.Data.Length, index => dresult.Data[index] = dinput.Data[index]); + CommonParallel.For(0, dinput.Values.Length, index => dresult.Values[index] = dinput.Values[index]); // LU solve by overwriting result. var dfactors = (DenseMatrix)Factors; - Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Data); + Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Values); } /// diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs index 341b9e9e..66ce7ce7 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseQR.cs @@ -188,7 +188,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseSvd.cs b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseSvd.cs index 6571eb6c..57c1c32e 100644 --- a/src/Numerics/LinearAlgebra/Complex/Factorization/DenseSvd.cs +++ b/src/Numerics/LinearAlgebra/Complex/Factorization/DenseSvd.cs @@ -70,7 +70,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization VectorS = new DenseVector(nm); MatrixU = new DenseMatrix(matrix.RowCount); MatrixVT = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); + Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); } /// @@ -126,7 +126,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do SVD factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); } /// @@ -176,7 +176,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex.Factorization throw new NotSupportedException("Can only do SVD factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs index 49bf71f5..c5764779 100644 --- a/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex/SparseMatrix.cs @@ -98,6 +98,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The number of columns. /// /// The value which we assign to each element of the matrix. + [Obsolete("Use a dense matrix instead.")] public SparseMatrix(int rows, int columns, Complex value) : this(rows, columns) { diff --git a/src/Numerics/LinearAlgebra/Complex/SparseVector.cs b/src/Numerics/LinearAlgebra/Complex/SparseVector.cs index f568e086..3afe6132 100644 --- a/src/Numerics/LinearAlgebra/Complex/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex/SparseVector.cs @@ -37,6 +37,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex using Generic; using NumberTheory; using Properties; + using Storage; using Threading; /// @@ -46,24 +47,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex [Serializable] public class SparseVector : Vector { - /// - /// Gets the vector's internal data. The array containing the actual values; only the non-zero values are stored. - /// - private Complex[] _nonZeroValues = new Complex[0]; + readonly SparseVectorStorage _storage; - /// - /// The indices of the non-zero entries. - /// - private int[] _nonZeroIndices = new int[0]; - - /// - /// Gets the number of non zero elements in the vector. - /// - /// The number of non zero elements. - public int NonZerosCount + internal SparseVector(SparseVectorStorage storage) + : base(storage) { - get; - private set; + _storage = storage; } #region Constructors @@ -77,7 +66,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// If is less than one. /// - public SparseVector(int size) : base(size) + public SparseVector(int size) + : this(new SparseVectorStorage(size)) { } @@ -94,58 +84,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// If is less than one. /// - public SparseVector(int size, Complex value) : this(size) + [Obsolete("Use a dense vector instead.")] + public SparseVector(int size, Complex value) + : this(new SparseVectorStorage(size)) { if (value == Complex.Zero) { - // Skip adding values return; } - // We already know that this vector is "full", let's allocate all needed memory - _nonZeroValues = new Complex[size]; - _nonZeroIndices = new int[size]; - NonZerosCount = size; + var valueCount = _storage.ValueCount = size; + var indices = _storage.Indices = new int[valueCount]; + var values = _storage.Values = new Complex[valueCount]; - CommonParallel.For( - 0, - Count, - index => - { - _nonZeroValues[index] = value; - _nonZeroIndices[index] = index; - }); - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public SparseVector(Vector other) : this(other.Count) - { - var vector = other as SparseVector; - if (vector == null) + for (int i = 0; i < values.Length; i++) { - for (var i = 0; i < other.Count; i++) - { - this[i] = other.At(i); - } - } - else - { - _nonZeroValues = new Complex[vector.NonZerosCount]; - _nonZeroIndices = new int[vector.NonZerosCount]; - NonZerosCount = vector.NonZerosCount; - - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - if (vector.NonZerosCount != 0) - { - CommonParallel.For(0, vector.NonZerosCount, index => _nonZeroValues[index] = vector._nonZeroValues[index]); - Buffer.BlockCopy(vector._nonZeroIndices, 0, _nonZeroIndices, 0, vector.NonZerosCount * Constants.SizeOfInt); - } + values[i] = value; + indices[i] = i; } } @@ -156,18 +111,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// The vector to create the new vector from. /// - public SparseVector(SparseVector other) : this(other.Count) + public SparseVector(Vector other) + : this(new SparseVectorStorage(other.Count)) { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - _nonZeroValues = new Complex[other.NonZerosCount]; - _nonZeroIndices = new int[other.NonZerosCount]; - NonZerosCount = other.NonZerosCount; - - if (other.NonZerosCount != 0) - { - CommonParallel.For(0, other.NonZerosCount, index => _nonZeroValues[index] = other._nonZeroValues[index]); - Buffer.BlockCopy(other._nonZeroIndices, 0, _nonZeroIndices, 0, other.NonZerosCount * Constants.SizeOfInt); - } + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -175,11 +122,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// /// The array to create this vector from. /// The vector copy the array. Any changes to the vector will NOT change the array. - public SparseVector(IList array) : this(array.Count) + public SparseVector(IList array) + : this(new SparseVectorStorage(array.Count)) { for (var i = 0; i < array.Count; i++) { - this[i] = array[i]; + Storage.At(i, array[i]); } } @@ -192,9 +140,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override Matrix ToColumnMatrix() { var matrix = new SparseMatrix(Count, 1); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(_nonZeroIndices[i], 0, _nonZeroValues[i]); + matrix.At(_storage.Indices[i], 0, _storage.Values[i]); } return matrix; @@ -207,9 +155,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override Matrix ToRowMatrix() { var matrix = new SparseMatrix(1, Count); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(0, _nonZeroIndices[i], _nonZeroValues[i]); + matrix.At(0, _storage.Indices[i], _storage.Values[i]); } return matrix; @@ -248,68 +196,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex return new SparseVector(size); } - /// - /// Resets all values to zero. - /// - public override void Clear() - { - NonZerosCount = 0; - } - - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as SparseVector; - if (otherVector == null) - { - target.Clear(); - - for (var index = 0; index < NonZerosCount; index++) - { - target.At(_nonZeroIndices[index], _nonZeroValues[index]); - } - } - else - { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - otherVector._nonZeroValues = new Complex[NonZerosCount]; - otherVector._nonZeroIndices = new int[NonZerosCount]; - otherVector.NonZerosCount = NonZerosCount; - - if (NonZerosCount != 0) - { - Array.Copy(_nonZeroValues, 0, otherVector._nonZeroValues, 0, NonZerosCount); - Buffer.BlockCopy(_nonZeroIndices, 0, otherVector._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); - } - } - } - /// /// Conjugates vector and save result to /// @@ -341,14 +227,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex else { // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - otherVector._nonZeroValues = new Complex[NonZerosCount]; - otherVector._nonZeroIndices = new int[NonZerosCount]; - otherVector.NonZerosCount = NonZerosCount; + otherVector._storage.Values = new Complex[_storage.ValueCount]; + otherVector._storage.Indices = new int[_storage.ValueCount]; + otherVector._storage.ValueCount = _storage.ValueCount; - if (NonZerosCount != 0) + if (_storage.ValueCount != 0) { - CommonParallel.For(0, NonZerosCount, index => otherVector._nonZeroValues[index] = _nonZeroValues[index].Conjugate()); - Buffer.BlockCopy(_nonZeroIndices, 0, otherVector._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); + CommonParallel.For(0, _storage.ValueCount, index => otherVector._storage.Values[index] = _storage.Values[index].Conjugate()); + Buffer.BlockCopy(_storage.Indices, 0, otherVector._storage.Indices, 0, _storage.ValueCount * Constants.SizeOfInt); } } } @@ -390,15 +276,15 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } //populate the non zero values from this - for (int j = 0; j < NonZerosCount; j++) + for (int j = 0; j < _storage.ValueCount; j++) { - vnonZeroValues[_nonZeroIndices[j]] = _nonZeroValues[j] + scalar; + vnonZeroValues[_storage.Indices[j]] = _storage.Values[j] + scalar; } //assign this vectors arrary to the new arrays. - _nonZeroValues = vnonZeroValues; - _nonZeroIndices = vnonZeroIndices; - NonZerosCount = Count; + _storage.Values = vnonZeroValues; + _storage.Indices = vnonZeroIndices; + _storage.ValueCount = Count; } else { @@ -439,22 +325,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] += otherSparse._nonZeroValues[j++]; + _storage.Values[i++] += otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != Complex.Zero) { - InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], otherValue); + InsertAtUnchecked(i++, otherSparse._storage.Indices[j], otherValue); } j++; } @@ -464,25 +350,25 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] + otherSparse.At(next)); + result.At(next, _storage.Values[i] + otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) + otherSparse._nonZeroValues[j]); + result.At(next, At(next) + otherSparse._storage.Values[j]); } j++; } @@ -585,22 +471,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] -= otherSparse._nonZeroValues[j++]; + _storage.Values[i++] -= otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != Complex.Zero) { - InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue); + InsertAtUnchecked(i++, otherSparse._storage.Indices[j], -otherValue); } j++; } @@ -610,25 +496,25 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] - otherSparse.At(next)); + result.At(next, _storage.Values[i] - otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) - otherSparse._nonZeroValues[j]); + result.At(next, At(next) - otherSparse._storage.Values[j]); } j++; } @@ -687,23 +573,18 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new SparseVector(Count) - { - _nonZeroValues = new Complex[NonZerosCount], - _nonZeroIndices = new int[NonZerosCount], - NonZerosCount = NonZerosCount - }; + var result = new SparseVectorStorage(Count); + var valueCount = result.ValueCount = _storage.ValueCount; + var indices = result.Indices = new int[valueCount]; + var values = result.Values = new Complex[valueCount]; - if (NonZerosCount != 0) + if (valueCount != 0) { - CommonParallel.For( - 0, - NonZerosCount, - index => result._nonZeroValues[index] = -_nonZeroValues[index]); - Buffer.BlockCopy(_nonZeroIndices, 0, result._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); + CommonParallel.For(0, valueCount, index => values[index] = -_storage.Values[index]); + Buffer.BlockCopy(_storage.Indices, 0, indices, 0, valueCount * Constants.SizeOfInt); } - return result; + return new SparseVector(result); } /// @@ -724,7 +605,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } var copy = new SparseVector(this); - Control.LinearAlgebraProvider.ScaleArray(scalar, copy._nonZeroValues, copy._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, copy._storage.Values, copy._storage.Values); return copy; } @@ -759,23 +640,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex if (sparseResult == null) { result.Clear(); - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - result.At(_nonZeroIndices[index], scalar * _nonZeroValues[index]); + result.At(_storage.Indices[index], scalar * _storage.Values[index]); } } else { if (!ReferenceEquals(this, result)) { - sparseResult.NonZerosCount = NonZerosCount; - sparseResult._nonZeroIndices = new int[NonZerosCount]; - Buffer.BlockCopy(_nonZeroIndices, 0, sparseResult._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); - sparseResult._nonZeroValues = new Complex[NonZerosCount]; - Array.Copy(_nonZeroValues, sparseResult._nonZeroValues, NonZerosCount); + sparseResult._storage.ValueCount = _storage.ValueCount; + sparseResult._storage.Indices = new int[_storage.ValueCount]; + Buffer.BlockCopy(_storage.Indices, 0, sparseResult._storage.Indices, 0, _storage.ValueCount * Constants.SizeOfInt); + sparseResult._storage.Values = new Complex[_storage.ValueCount]; + Array.Copy(_storage.Values, sparseResult._storage.Values, _storage.ValueCount); } - Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._nonZeroValues, sparseResult._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._storage.Values, sparseResult._storage.Values); } } @@ -794,16 +675,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * _nonZeroValues[i]; + result += _storage.Values[i] * _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * other.At(_nonZeroIndices[i]); + result += _storage.Values[i] * other.At(_storage.Indices[i]); } } @@ -895,17 +776,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The index of absolute minimum element. public override int AbsoluteMinimumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { // No non-zero elements. Return 0 return 0; } var index = 0; - var min = _nonZeroValues[index].Magnitude; - for (var i = 1; i < NonZerosCount; i++) + var min = _storage.Values[index].Magnitude; + for (var i = 1; i < _storage.ValueCount; i++) { - var test = _nonZeroValues[i].Magnitude; + var test = _storage.Values[i].Magnitude; if (test < min) { index = i; @@ -913,7 +794,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -984,9 +865,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override Complex Sum() { var result = Complex.Zero; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i]; + result += _storage.Values[i]; } return result; @@ -999,9 +880,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override Complex SumMagnitudes() { double result = 0; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i].Magnitude; + result += _storage.Values[i].Magnitude; } return result; @@ -1016,17 +897,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] *= _nonZeroValues[i]; + _storage.Values[i] *= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, other.At(index) * _nonZeroValues[i]); + var index = _storage.Indices[i]; + result.At(index, other.At(index) * _storage.Values[i]); } } } @@ -1040,17 +921,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] /= _nonZeroValues[i]; + _storage.Values[i] /= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, _nonZeroValues[i] / other.At(index)); + var index = _storage.Indices[i]; + result.At(index, _storage.Values[i] / other.At(index)); } } } @@ -1076,13 +957,13 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } var matrix = new SparseMatrix(u.Count, v.Count); - for (var i = 0; i < u.NonZerosCount; i++) + for (var i = 0; i < u._storage.ValueCount; i++) { - for (var j = 0; j < v.NonZerosCount; j++) + for (var j = 0; j < v._storage.ValueCount; j++) { - if (u._nonZeroIndices[i] == v._nonZeroIndices[j]) + if (u._storage.Indices[i] == v._storage.Indices[j]) { - matrix.At(i, j, u._nonZeroValues[i] * v._nonZeroValues[j]); + matrix.At(i, j, u._storage.Values[i] * v._storage.Values[j]); } } } @@ -1118,25 +999,25 @@ namespace MathNet.Numerics.LinearAlgebra.Complex throw new ArgumentOutOfRangeException("p"); } - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0.0; } if (2.0 == p) { - return _nonZeroValues.Aggregate(Complex.Zero, SpecialFunctions.Hypotenuse).Magnitude; + return _storage.Values.Aggregate(Complex.Zero, SpecialFunctions.Hypotenuse).Magnitude; } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(0, NonZerosCount, i => _nonZeroValues[i].Magnitude, Math.Max, 0d); + return CommonParallel.Aggregate(0, _storage.ValueCount, i => _storage.Values[i].Magnitude, Math.Max, 0d); } var sum = 0.0; - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - sum += Math.Pow(_nonZeroValues[index].Magnitude, p); + sum += Math.Pow(_storage.Values[index].Magnitude, p); } return Math.Pow(sum, 1.0 / p); @@ -1302,97 +1183,33 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } #endregion - /// - /// Gets the value at the given index. - /// - /// Value real index in array - /// The value at the given index. - internal protected override Complex At(int index) - { - // Search if item idex exists in NonZeroIndices array in range "0 - real nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - return itemIndex >= 0 ? _nonZeroValues[itemIndex] : Complex.Zero; - } - - /// - /// Delete, Add or Update the value in NonZeroValues and NonZeroIndices - /// - /// Value real index in array - /// The value to set. - /// This method assume that index is between 0 and Array Size - internal protected override void At(int index, Complex value) - { - // Search if "index" already exists in range "0 - complex nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - - if (itemIndex >= 0) - { - // Item already exist at itemIndex - if (value == Complex.Zero) - { - RemoveAtUnchecked(itemIndex); - } - else - { - _nonZeroValues[itemIndex] = value; - } - } - else - { - if (value != Complex.Zero) - { - InsertAtUnchecked(~itemIndex, index, value); - } - } - } - private void InsertAtUnchecked(int itemIndex, int index, Complex value) { // Check if the storage needs to be increased - if ((NonZerosCount == _nonZeroValues.Length) && (NonZerosCount < Count)) + if ((_storage.ValueCount == _storage.Values.Length) && (_storage.ValueCount < Count)) { // Value and Indices arrays are completely full so we increase the size - var size = Math.Min(_nonZeroValues.Length + GrowthSize(), Count); - Array.Resize(ref _nonZeroValues, size); - Array.Resize(ref _nonZeroIndices, size); + var size = Math.Min(_storage.Values.Length + GrowthSize(), Count); + Array.Resize(ref _storage.Values, size); + Array.Resize(ref _storage.Indices, size); } // Move all values (with a position larger than index) in the value array // to the next position // Move all values (with a position larger than index) in the columIndices // array to the next position - for (var i = NonZerosCount - 1; i > itemIndex - 1; i--) + for (var i = _storage.ValueCount - 1; i > itemIndex - 1; i--) { - _nonZeroValues[i + 1] = _nonZeroValues[i]; - _nonZeroIndices[i + 1] = _nonZeroIndices[i]; + _storage.Values[i + 1] = _storage.Values[i]; + _storage.Indices[i + 1] = _storage.Indices[i]; } // Add the value and the column index - _nonZeroValues[itemIndex] = value; - _nonZeroIndices[itemIndex] = index; + _storage.Values[itemIndex] = value; + _storage.Indices[itemIndex] = index; // increase the number of non-zero numbers by one - NonZerosCount += 1; - } - - private void RemoveAtUnchecked(int itemIndex) - { - // Value is zero. Let's delete it from Values and Indices array - for (var i = itemIndex + 1; i < NonZerosCount; i++) - { - _nonZeroValues[i - 1] = _nonZeroValues[i]; - _nonZeroIndices[i - 1] = _nonZeroIndices[i]; - } - - NonZerosCount -= 1; - - // Check whether we need to shrink the arrays. This is reasonable to do if - // there are a lot of non-zero elements and storage is two times bigger - if ((NonZerosCount > 1024) && (NonZerosCount < _nonZeroIndices.Length / 2)) - { - Array.Resize(ref _nonZeroValues, NonZerosCount); - Array.Resize(ref _nonZeroIndices, NonZerosCount); - } + _storage.ValueCount += 1; } /// @@ -1403,19 +1220,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex private int GrowthSize() { int delta; - if (_nonZeroValues.Length > 1024) + if (_storage.Values.Length > 1024) { - delta = _nonZeroValues.Length / 4; + delta = _storage.Values.Length / 4; } else { - if (_nonZeroValues.Length > 256) + if (_storage.Values.Length > 256) { delta = 512; } else { - delta = _nonZeroValues.Length > 64 ? 128 : 32; + delta = _storage.Values.Length > 64 ? 128 : 32; } } @@ -1428,7 +1245,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { if (Count > 20) { - return String.Format("SparseVectorOfComplex({0},{1},{2})", Count, NonZerosCount, GetHashCode()); + return String.Format("SparseVectorOfComplex({0},{1},{2})", Count, _storage.ValueCount, GetHashCode()); } return base.ToString(format, formatProvider); @@ -1442,14 +1259,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// public override int GetHashCode() { - var hashNum = Math.Min(NonZerosCount, 20); + var hashNum = Math.Min(_storage.ValueCount, 20); long hash = 0; for (var i = 0; i < hashNum; i++) { #if PORTABLE - hash ^= Precision.DoubleToInt64Bits(this._nonZeroValues[i].GetHashCode()); + hash ^= Precision.DoubleToInt64Bits(this._storage.Values[i].GetHashCode()); #else - hash ^= BitConverter.DoubleToInt64Bits(_nonZeroValues[i].GetHashCode()); + hash ^= BitConverter.DoubleToInt64Bits(_storage.Values[i].GetHashCode()); #endif } @@ -1493,27 +1310,27 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { - if (_nonZeroValues[i++] != Complex.Zero) + if (_storage.Values[i++] != Complex.Zero) { return false; } continue; } - if (i >= NonZerosCount || j < otherSparse.NonZerosCount && otherSparse._nonZeroIndices[j] < _nonZeroIndices[i]) + if (i >= _storage.ValueCount || j < otherSparse._storage.ValueCount && otherSparse._storage.Indices[j] < _storage.Indices[i]) { - if (otherSparse._nonZeroValues[j++] != Complex.Zero) + if (otherSparse._storage.Values[j++] != Complex.Zero) { return false; } continue; } - if (!_nonZeroValues[i].AlmostEqual(otherSparse._nonZeroValues[j])) + if (!_storage.Values[i].AlmostEqual(otherSparse._storage.Values[j])) { return false; } @@ -1541,9 +1358,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// public override IEnumerable> GetIndexedEnumerator() { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - yield return new Tuple(_nonZeroIndices[i], _nonZeroValues[i]); + yield return new Tuple(_storage.Indices[i], _storage.Values[i]); } } @@ -1556,9 +1373,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex public override Complex[] ToArray() { var ret = new Complex[Count]; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - ret[_nonZeroIndices[i]] = _nonZeroValues[i]; + ret[_storage.Indices[i]] = _storage.Values[i]; } return ret; diff --git a/src/Numerics/LinearAlgebra/Complex/Vector.cs b/src/Numerics/LinearAlgebra/Complex/Vector.cs index 21eb5ea5..6459df39 100644 --- a/src/Numerics/LinearAlgebra/Complex/Vector.cs +++ b/src/Numerics/LinearAlgebra/Complex/Vector.cs @@ -31,6 +31,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex using Distributions; using Generic; using Properties; + using Storage; using Threading; /// @@ -41,15 +42,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex { /// /// Initializes a new instance of the Vector class. - /// Constructs a Vector with the given size. /// - /// - /// The size of the Vector to construct. - /// - /// - /// If is less than one. - /// - protected Vector(int size) : base(size) + protected Vector(VectorStorage storage) + : base(storage) { } diff --git a/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs index ba6fcaf8..c5c64796 100644 --- a/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/DenseMatrix.cs @@ -39,8 +39,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 [Serializable] public class DenseMatrix : Matrix { - readonly DenseColumnMajorMatrixStorage _storage; - /// /// Number of rows. /// @@ -59,15 +57,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// Gets the matrix's data. /// /// The matrix's data. - readonly Complex32[] _data; + readonly Complex32[] _values; internal DenseMatrix(DenseColumnMajorMatrixStorage storage) : base(storage) { - _storage = storage; - _rowCount = _storage.RowCount; - _columnCount = _storage.ColumnCount; - _data = _storage.Data; + _rowCount = storage.RowCount; + _columnCount = storage.ColumnCount; + _values = storage.Data; } /// @@ -109,9 +106,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public DenseMatrix(int rows, int columns, Complex32 value) : this(rows, columns) { - for (var i = 0; i < _data.Length; i++) + for (var i = 0; i < _values.Length; i++) { - _data[i] = value; + _values[i] = value; } } @@ -139,7 +136,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { for (var j = 0; j < _columnCount; j++) { - _data[(j * _rowCount) + i] = array[i, j]; + _values[(j * _rowCount) + i] = array[i, j]; } } } @@ -161,7 +158,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The matrix's data. public Complex32[] Data { - get { return _data; } + get { return _values; } } /// @@ -203,7 +200,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var index = j * _rowCount; for (var i = 0; i < _rowCount; i++) { - ret._data[(i * _columnCount) + j] = _data[index + i]; + ret._values[(i * _columnCount) + j] = _values[index + i]; } } @@ -214,21 +211,21 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The L1 norm of the matrix. public override Complex32 L1Norm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _values); } /// Calculates the Frobenius norm of this matrix. /// The Frobenius norm of this matrix. public override Complex32 FrobeniusNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _values); } /// Calculates the infinity norm of this matrix. /// The infinity norm of this matrix. public override Complex32 InfinityNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _values); } #region Static constructors for special matrices. @@ -246,7 +243,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var m = new DenseMatrix(order); for (var i = 0; i < order; i++) { - m._data[(i * order) + i] = 1.0f; + m._values[(i * order) + i] = 1.0f; } return m; @@ -268,7 +265,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult._values); } } @@ -292,14 +289,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0f, - denseResult.Data); + denseResult.Values); } } @@ -323,14 +320,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0f, - denseResult._data); + denseResult._values); } } @@ -354,14 +351,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.Transpose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0f, - denseResult._data); + denseResult._values); } } @@ -385,14 +382,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0f, - denseResult.Data); + denseResult.Values); } } @@ -416,14 +413,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0f, - denseResult._data); + denseResult._values); } } @@ -441,7 +438,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(-1, _values, denseResult._values); } } @@ -461,7 +458,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_values, denseOther._values, denseResult._values); } } @@ -481,7 +478,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseDivideArrays(_values, denseOther._values, denseResult._values); } } @@ -502,7 +499,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.AddArrays(_values, denseOther._values, denseResult._values); } } @@ -521,7 +518,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.SubtractArrays(_values, denseOther._values, denseResult._values); } } @@ -537,7 +534,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var index = j * _rowCount; for (var i = 0; i < _rowCount; i++) { - ret._data[(i * _columnCount) + j] = _data[index + i].Conjugate(); + ret._values[(i * _columnCount) + j] = _values[index + i].Conjugate(); } } @@ -559,7 +556,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var sum = Complex32.Zero; for (var i = 0; i < _rowCount; i++) { - sum += _data[(i * _rowCount) + i]; + sum += _values[(i * _rowCount) + i]; } return sum; diff --git a/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs b/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs index 0273a55b..0886cc90 100644 --- a/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex32/DenseVector.cs @@ -34,6 +34,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 using NumberTheory; using Numerics; using Properties; + using Storage; using Threading; /// @@ -42,6 +43,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 [Serializable] public class DenseVector : Vector { + /// + /// Number of elements + /// + readonly int _length; + + /// + /// Gets the vector's data. + /// + readonly Complex32[] _values; + + internal DenseVector(DenseVectorStorage storage) + : base(storage) + { + _length = storage.Length; + _values = storage.Data; + } + /// /// Initializes a new instance of the class with a given size. /// @@ -51,9 +69,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// If is less than one. /// - public DenseVector(int size) : base(size) + public DenseVector(int size) + : this(new DenseVectorStorage(size)) { - Data = new Complex32[size]; } /// @@ -69,11 +87,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// If is less than one. /// - public DenseVector(int size, Complex32 value) : this(size) + public DenseVector(int size, Complex32 value) + : this(size) { - for (var index = 0; index < Data.Length; index++) + for (var index = 0; index < _values.Length; index++) { - Data[index] = value; + _values[index] = value; } } @@ -84,27 +103,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// The vector to create the new vector from. /// - public DenseVector(Vector other) : this(other.Count) + public DenseVector(Vector other) + : this(other.Count) { - CommonParallel.For( - 0, - Data.Length, - index => this[index] = other[index]); - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public DenseVector(DenseVector other) : this(other.Count) - { - CommonParallel.For( - 0, - Data.Length, - index => Data[index] = other.Data[index]); + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -113,20 +115,18 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The array to create this vector from. /// The vector does not copy the array, but keeps a reference to it. Any /// changes to the vector will also change the array. - public DenseVector(Complex32[] array) : base(array.Length) + public DenseVector(Complex32[] array) + : this(new DenseVectorStorage(array.Length, array)) { - Data = array; } /// - /// Gets the vector's internal data. + /// Gets the vector's data. /// - /// The vector's internal data. - /// Changing values in the array also changes the corresponding value in vector. Use with care. - internal Complex32[] Data + /// The vector's data. + public Complex32[] Values { - get; - private set; + get { return _values; } } /// @@ -137,14 +137,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// A reference to the internal date of the given vector. /// - public static implicit operator Complex32[](DenseVector vector) + public static explicit operator Complex32[](DenseVector vector) { if (vector == null) { throw new ArgumentNullException(); } - return vector.Data; + return vector.Values; } /// @@ -170,10 +170,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// This vector as a column matrix. public override Matrix ToColumnMatrix() { - var matrix = new DenseMatrix(Count, 1); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(_length, 1); + for (var i = 0; i < _values.Length; i++) { - matrix.At(i, 0, Data[i]); + matrix.At(i, 0, _values[i]); } return matrix; @@ -185,33 +185,15 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// This vector as a row matrix. public override Matrix ToRowMatrix() { - var matrix = new DenseMatrix(1, Count); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(1, _length); + for (var i = 0; i < _values.Length; i++) { - matrix.At(0, i, Data[i]); + matrix.At(0, i, _values[i]); } return matrix; } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - /// If is negative or - /// greater than the size of the vector. - public override Complex32 this[int index] - { - get - { - return Data[index]; - } - - set - { - Data[index] = value; - } - } - /// /// Creates a matrix with the given dimensions using the same storage type /// as this vector. @@ -245,49 +227,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 return new DenseVector(size); } - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as DenseVector; - if (otherVector == null) - { - CommonParallel.For( - 0, - Data.Length, - index => target[index] = Data[index]); - } - else - { - Array.Copy(Data, 0, otherVector.Data, 0, Data.Length); - } - } - /// /// Adds a scalar to each element of the vector and stores the result in the result vector. /// @@ -304,8 +243,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] + scalar); + _values.Length, + index => dense._values[index] = _values[index] + scalar); } } @@ -320,7 +259,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, Complex32.One, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, Complex32.One, odense.Values, rdense.Values); } else { @@ -389,8 +328,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] - scalar); + _values.Length, + index => dense._values[index] = _values[index] - scalar); } } @@ -405,7 +344,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, -1.0f, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, -1.0f, odense.Values, rdense.Values); } else { @@ -464,11 +403,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new DenseVector(Count); + var result = new DenseVector(_length); CommonParallel.For( 0, - Data.Length, - index => result[index] = -Data[index]); + _values.Length, + index => result[index] = -_values[index]); return result; } @@ -488,7 +427,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, Data, denseResult.Data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult.Values); } } @@ -502,7 +441,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { var denseVector = other as DenseVector; - return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(Data, denseVector.Data); + return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(_values, denseVector.Values); } /// @@ -564,7 +503,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide"); } - return Control.LinearAlgebraProvider.DotProduct(leftSide.Data, rightSide.Data); + return Control.LinearAlgebraProvider.DotProduct(leftSide.Values, rightSide.Values); } /// @@ -591,10 +530,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override int AbsoluteMinimumIndex() { var index = 0; - var min = Data[index].Magnitude; - for (var i = 1; i < Count; i++) + var min = _values[index].Magnitude; + for (var i = 1; i < _length; i++) { - var test = Data[i].Magnitude; + var test = _values[i].Magnitude; if (test < min) { index = i; @@ -611,7 +550,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The value of the absolute minimum element. public override Complex32 AbsoluteMinimum() { - return Data[AbsoluteMinimumIndex()].Magnitude; + return _values[AbsoluteMinimumIndex()].Magnitude; } /// @@ -620,7 +559,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The value of the absolute maximum element. public override Complex32 AbsoluteMaximum() { - return Data[AbsoluteMaximumIndex()].Magnitude; + return _values[AbsoluteMaximumIndex()].Magnitude; } /// @@ -630,10 +569,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override int AbsoluteMaximumIndex() { var index = 0; - var max = Data[index].Magnitude; - for (var i = 1; i < Count; i++) + var max = _values[index].Magnitude; + for (var i = 1; i < _length; i++) { - var test = Data[i].Magnitude; + var test = _values[i].Magnitude; if (test > max) { index = i; @@ -657,7 +596,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// If is not positive. public override Vector SubVector(int index, int length) { - if (index < 0 || index >= Count) + if (index < 0 || index >= _length) { throw new ArgumentOutOfRangeException("index"); } @@ -667,7 +606,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 throw new ArgumentOutOfRangeException("length"); } - if (index + length > Count) + if (index + length > _length) { throw new ArgumentOutOfRangeException("length"); } @@ -676,8 +615,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 CommonParallel.For( index, - index + length, - i => result.Data[i - index] = Data[i]); + index + length, + i => result._values[i - index] = _values[i]); return result; } @@ -694,7 +633,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 throw new ArgumentNullException("values"); } - if (values.Length != Count) + if (values.Length != _length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "values"); } @@ -702,7 +641,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 CommonParallel.For( 0, values.Length, - i => Data[i] = values[i]); + i => _values[i] = values[i]); } /// @@ -713,9 +652,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { var sum = Complex32.Zero; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Data[i]; + sum += _values[i]; } return sum; @@ -729,9 +668,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { var sum = Complex32.Zero; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Data[i].Magnitude; + sum += _values[i].Magnitude; } return sum; @@ -753,8 +692,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] * other[index]); + _values.Length, + index => dense._values[index] = _values[index] * other[index]); } } @@ -775,8 +714,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] / other[index]); + _values.Length, + index => dense._values[index] = _values[index] / other[index]); } } @@ -808,7 +747,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { for (var j = 0; j < v.Count; j++) { - matrix.At(i, j, u.Data[i] * v.Data[j]); + matrix.At(i, j, u._values[i] * v._values[j]); } }); return matrix; @@ -832,9 +771,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } var v = (DenseVector)CreateVector(length); - for (var index = 0; index < v.Data.Length; index++) + for (var index = 0; index < v._values.Length; index++) { - v.Data[index] = new Complex32((float)randomDistribution.Sample(), (float)randomDistribution.Sample()); + v._values[index] = new Complex32((float)randomDistribution.Sample(), (float)randomDistribution.Sample()); } return v; @@ -858,9 +797,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } var v = (DenseVector)CreateVector(length); - for (var index = 0; index < v.Data.Length; index++) + for (var index = 0; index < v._values.Length; index++) { - v.Data[index] = new Complex32(randomDistribution.Sample(), randomDistribution.Sample()); + v._values[index] = new Complex32(randomDistribution.Sample(), randomDistribution.Sample()); } return v; @@ -898,19 +837,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 if (2.0 == p) { - return Data.Aggregate(Complex32.Zero, SpecialFunctions.Hypotenuse).Magnitude; + return _values.Aggregate(Complex32.Zero, SpecialFunctions.Hypotenuse).Magnitude; } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(Data, (i, v) => v.Magnitude, Math.Max, 0f); + return CommonParallel.Aggregate(_values, (i, v) => v.Magnitude, Math.Max, 0f); } var sum = 0.0; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Math.Pow(Data[i].Magnitude, p); + sum += Math.Pow(_values[i].Magnitude, p); } return (float)Math.Pow(sum, 1.0 / p); @@ -1075,14 +1014,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 #endregion - /// - /// Resets all values to zero. - /// - public override void Clear() - { - Array.Clear(Data, 0, Data.Length); - } - /// /// Conjugates vector and save result to /// @@ -1099,25 +1030,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { CommonParallel.For( 0, - Count, - index => denseTarget.Data[index] = Data[index].Conjugate()); + _length, + index => denseTarget._values[index] = _values[index].Conjugate()); } } - - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - internal protected override Complex32 At(int index) - { - return Data[index]; - } - - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - internal protected override void At(int index, Complex32 value) - { - Data[index] = value; - } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs index fbeff2e1..28eb1a36 100644 --- a/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/DiagonalMatrix.cs @@ -378,12 +378,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 return; } - if (_data.Length != denseSource.Data.Length) + if (_data.Length != denseSource.Values.Length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "source"); } - Array.Copy(denseSource.Data, _data, denseSource.Data.Length); + Array.Copy(denseSource.Values, _data, denseSource.Values.Length); } /// diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseCholesky.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseCholesky.cs index 4365c98e..ee96174d 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseCholesky.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseCholesky.cs @@ -169,11 +169,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization } // Copy the contents of input to result. - CommonParallel.For(0, dinput.Data.Length, index => dresult.Data[index] = dinput.Data[index]); + CommonParallel.For(0, dinput.Values.Length, index => dresult.Values[index] = dinput.Values[index]); // Cholesky solve by overwriting result. var dfactor = (DenseMatrix)CholeskyFactor; - Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Data, 1); + Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Values, 1); } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseEvd.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseEvd.cs index d62047f1..9dca246b 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseEvd.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseEvd.cs @@ -108,7 +108,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization { var matrixH = matrix.ToArray(); NonsymmetricReduceToHessenberg(((DenseMatrix)MatrixEv).Data, matrixH, order); - NonsymmetricReduceHessenberToRealSchur(((LinearAlgebra.Complex.DenseVector)VectorEv).Data, ((DenseMatrix)MatrixEv).Data, matrixH, order); + NonsymmetricReduceHessenberToRealSchur(((LinearAlgebra.Complex.DenseVector)VectorEv).Values, ((DenseMatrix)MatrixEv).Data, matrixH, order); } for (var i = 0; i < VectorEv.Count; i++) diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs index 36919fc8..ccc50c73 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseGramSchmidt.cs @@ -217,7 +217,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Data, 1, dresult.Data); + _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseLU.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseLU.cs index 30572c09..b37f97bf 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseLU.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseLU.cs @@ -170,11 +170,11 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization } // Copy the contents of input to result. - CommonParallel.For(0, dinput.Data.Length, index => dresult.Data[index] = dinput.Data[index]); + CommonParallel.For(0, dinput.Values.Length, index => dresult.Values[index] = dinput.Values[index]); // LU solve by overwriting result. var dfactors = (DenseMatrix)Factors; - Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Data); + Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Values); } /// diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs index 04f9d8cb..e54af3dc 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseQR.cs @@ -188,7 +188,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseSvd.cs b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseSvd.cs index 515a5f0a..d2b5896e 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseSvd.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Factorization/DenseSvd.cs @@ -70,7 +70,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization VectorS = new DenseVector(nm); MatrixU = new DenseMatrix(matrix.RowCount); MatrixVT = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); + Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); } /// @@ -126,7 +126,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do SVD factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); } /// @@ -176,7 +176,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32.Factorization throw new NotSupportedException("Can only do SVD factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs index 1c46f3ac..5ec92abe 100644 --- a/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Complex32/SparseMatrix.cs @@ -98,6 +98,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The number of columns. /// /// The value which we assign to each element of the matrix. + [Obsolete("Use a dense matrix instead.")] public SparseMatrix(int rows, int columns, Complex32 value) : this(rows, columns) { diff --git a/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs b/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs index b03a3be8..a5c47215 100644 --- a/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs @@ -37,6 +37,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 using NumberTheory; using Numerics; using Properties; + using Storage; using Threading; /// @@ -46,24 +47,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 [Serializable] public class SparseVector : Vector { - /// - /// Gets the vector's internal data. The array containing the actual values; only the non-zero values are stored. - /// - private Complex32[] _nonZeroValues = new Complex32[0]; - - /// - /// The indices of the non-zero entries. - /// - private int[] _nonZeroIndices = new int[0]; + readonly SparseVectorStorage _storage; - /// - /// Gets the number of non zero elements in the vector. - /// - /// The number of non zero elements. - public int NonZerosCount + internal SparseVector(SparseVectorStorage storage) + : base(storage) { - get; - private set; + _storage = storage; } #region Constructors @@ -77,7 +66,8 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// If is less than one. /// - public SparseVector(int size) : base(size) + public SparseVector(int size) + : this(new SparseVectorStorage(size)) { } @@ -94,58 +84,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// If is less than one. /// - public SparseVector(int size, Complex32 value) : this(size) + [Obsolete("Use a dense vector instead.")] + public SparseVector(int size, Complex32 value) + : this(new SparseVectorStorage(size)) { if (value == Complex32.Zero) { - // Skip adding values return; } - // We already know that this vector is "full", let's allocate all needed memory - _nonZeroValues = new Complex32[size]; - _nonZeroIndices = new int[size]; - NonZerosCount = size; + var valueCount = _storage.ValueCount = size; + var indices = _storage.Indices = new int[valueCount]; + var values = _storage.Values = new Complex32[valueCount]; - CommonParallel.For( - 0, - Count, - index => - { - _nonZeroValues[index] = value; - _nonZeroIndices[index] = index; - }); - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public SparseVector(Vector other) : this(other.Count) - { - var vector = other as SparseVector; - if (vector == null) + for (int i = 0; i < values.Length; i++) { - for (var i = 0; i < other.Count; i++) - { - this[i] = other.At(i); - } - } - else - { - _nonZeroValues = new Complex32[vector.NonZerosCount]; - _nonZeroIndices = new int[vector.NonZerosCount]; - NonZerosCount = vector.NonZerosCount; - - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - if (vector.NonZerosCount != 0) - { - CommonParallel.For(0, vector.NonZerosCount, index => _nonZeroValues[index] = vector._nonZeroValues[index]); - Buffer.BlockCopy(vector._nonZeroIndices, 0, _nonZeroIndices, 0, vector.NonZerosCount * Constants.SizeOfInt); - } + values[i] = value; + indices[i] = i; } } @@ -156,18 +111,10 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// The vector to create the new vector from. /// - public SparseVector(SparseVector other) : this(other.Count) + public SparseVector(Vector other) + : this(new SparseVectorStorage(other.Count)) { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - _nonZeroValues = new Complex32[other.NonZerosCount]; - _nonZeroIndices = new int[other.NonZerosCount]; - NonZerosCount = other.NonZerosCount; - - if (other.NonZerosCount != 0) - { - CommonParallel.For(0, other.NonZerosCount, index => _nonZeroValues[index] = other._nonZeroValues[index]); - Buffer.BlockCopy(other._nonZeroIndices, 0, _nonZeroIndices, 0, other.NonZerosCount * Constants.SizeOfInt); - } + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -175,11 +122,12 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// /// The array to create this vector from. /// The vector copy the array. Any changes to the vector will NOT change the array. - public SparseVector(IList array) : this(array.Count) + public SparseVector(IList array) + : this(new SparseVectorStorage(array.Count)) { for (var i = 0; i < array.Count; i++) { - this[i] = array[i]; + Storage.At(i, array[i]); } } @@ -192,9 +140,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override Matrix ToColumnMatrix() { var matrix = new SparseMatrix(Count, 1); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(_nonZeroIndices[i], 0, _nonZeroValues[i]); + matrix.At(_storage.Indices[i], 0, _storage.Values[i]); } return matrix; @@ -207,44 +155,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override Matrix ToRowMatrix() { var matrix = new SparseMatrix(1, Count); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(0, _nonZeroIndices[i], _nonZeroValues[i]); + matrix.At(0, _storage.Indices[i], _storage.Values[i]); } return matrix; } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - /// If is negative or - /// greater than the size of the vector. - public override Complex32 this[int index] - { - get - { - // If index is out of bounds - if ((index < 0) || (index >= Count)) - { - throw new IndexOutOfRangeException(); - } - - return At(index); - } - - set - { - // If index is out of bounds - if ((index < 0) || (index >= Count)) - { - throw new IndexOutOfRangeException(); - } - - At(index, value); - } - } - /// /// Creates a matrix with the given dimensions using the same storage type /// as this vector. @@ -278,68 +196,6 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 return new SparseVector(size); } - /// - /// Resets all values to zero. - /// - public override void Clear() - { - NonZerosCount = 0; - } - - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as SparseVector; - if (otherVector == null) - { - target.Clear(); - - for (var index = 0; index < NonZerosCount; index++) - { - target.At(_nonZeroIndices[index], _nonZeroValues[index]); - } - } - else - { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - otherVector._nonZeroValues = new Complex32[NonZerosCount]; - otherVector._nonZeroIndices = new int[NonZerosCount]; - otherVector.NonZerosCount = NonZerosCount; - - if (NonZerosCount != 0) - { - Array.Copy(_nonZeroValues, 0, otherVector._nonZeroValues, 0, NonZerosCount); - Buffer.BlockCopy(_nonZeroIndices, 0, otherVector._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); - } - } - } - /// /// Conjugates vector and save result to /// @@ -371,14 +227,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 else { // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - otherVector._nonZeroValues = new Complex32[NonZerosCount]; - otherVector._nonZeroIndices = new int[NonZerosCount]; - otherVector.NonZerosCount = NonZerosCount; + otherVector._storage.Values = new Complex32[_storage.ValueCount]; + otherVector._storage.Indices = new int[_storage.ValueCount]; + otherVector._storage.ValueCount = _storage.ValueCount; - if (NonZerosCount != 0) + if (_storage.ValueCount != 0) { - CommonParallel.For(0, NonZerosCount, index => otherVector._nonZeroValues[index] = _nonZeroValues[index].Conjugate()); - Buffer.BlockCopy(_nonZeroIndices, 0, otherVector._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); + CommonParallel.For(0, _storage.ValueCount, index => otherVector._storage.Values[index] = _storage.Values[index].Conjugate()); + Buffer.BlockCopy(_storage.Indices, 0, otherVector._storage.Indices, 0, _storage.ValueCount * Constants.SizeOfInt); } } } @@ -420,15 +276,15 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } //populate the non zero values from this - for (int j = 0; j < NonZerosCount; j++) + for (int j = 0; j < _storage.ValueCount; j++) { - vnonZeroValues[_nonZeroIndices[j]] = _nonZeroValues[j] + scalar; + vnonZeroValues[_storage.Indices[j]] = _storage.Values[j] + scalar; } //assign this vectors arrary to the new arrays. - _nonZeroValues = vnonZeroValues; - _nonZeroIndices = vnonZeroIndices; - NonZerosCount = Count; + _storage.Values = vnonZeroValues; + _storage.Indices = vnonZeroIndices; + _storage.ValueCount = Count; } else { @@ -469,22 +325,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] += otherSparse._nonZeroValues[j++]; + _storage.Values[i++] += otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != Complex32.Zero) { - InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], otherValue); + InsertAtUnchecked(i++, otherSparse._storage.Indices[j], otherValue); } j++; } @@ -494,25 +350,25 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] + otherSparse.At(next)); + result.At(next, _storage.Values[i] + otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) + otherSparse._nonZeroValues[j]); + result.At(next, At(next) + otherSparse._storage.Values[j]); } j++; } @@ -615,22 +471,22 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] -= otherSparse._nonZeroValues[j++]; + _storage.Values[i++] -= otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != Complex32.Zero) { - InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue); + InsertAtUnchecked(i++, otherSparse._storage.Indices[j], -otherValue); } j++; } @@ -640,25 +496,25 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] - otherSparse.At(next)); + result.At(next, _storage.Values[i] - otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) - otherSparse._nonZeroValues[j]); + result.At(next, At(next) - otherSparse._storage.Values[j]); } j++; } @@ -717,23 +573,18 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new SparseVector(Count) - { - _nonZeroValues = new Complex32[NonZerosCount], - _nonZeroIndices = new int[NonZerosCount], - NonZerosCount = NonZerosCount - }; + var result = new SparseVectorStorage(Count); + var valueCount = result.ValueCount = _storage.ValueCount; + var indices = result.Indices = new int[valueCount]; + var values = result.Values = new Complex32[valueCount]; - if (NonZerosCount != 0) + if (valueCount != 0) { - CommonParallel.For( - 0, - NonZerosCount, - index => result._nonZeroValues[index] = -_nonZeroValues[index]); - Buffer.BlockCopy(_nonZeroIndices, 0, result._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); + CommonParallel.For(0, valueCount, index => values[index] = -_storage.Values[index]); + Buffer.BlockCopy(_storage.Indices, 0, indices, 0, valueCount * Constants.SizeOfInt); } - return result; + return new SparseVector(result); } /// @@ -754,7 +605,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } var copy = new SparseVector(this); - Control.LinearAlgebraProvider.ScaleArray(scalar, copy._nonZeroValues, copy._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, copy._storage.Values, copy._storage.Values); return copy; } @@ -789,23 +640,23 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 if (sparseResult == null) { result.Clear(); - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - result.At(_nonZeroIndices[index], scalar * _nonZeroValues[index]); + result.At(_storage.Indices[index], scalar * _storage.Values[index]); } } else { if (!ReferenceEquals(this, result)) { - sparseResult.NonZerosCount = NonZerosCount; - sparseResult._nonZeroIndices = new int[NonZerosCount]; - Buffer.BlockCopy(_nonZeroIndices, 0, sparseResult._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); - sparseResult._nonZeroValues = new Complex32[NonZerosCount]; - Array.Copy(_nonZeroValues, sparseResult._nonZeroValues, NonZerosCount); + sparseResult._storage.ValueCount = _storage.ValueCount; + sparseResult._storage.Indices = new int[_storage.ValueCount]; + Buffer.BlockCopy(_storage.Indices, 0, sparseResult._storage.Indices, 0, _storage.ValueCount * Constants.SizeOfInt); + sparseResult._storage.Values = new Complex32[_storage.ValueCount]; + Array.Copy(_storage.Values, sparseResult._storage.Values, _storage.ValueCount); } - Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._nonZeroValues, sparseResult._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._storage.Values, sparseResult._storage.Values); } } @@ -824,16 +675,16 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * _nonZeroValues[i]; + result += _storage.Values[i] * _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * other.At(_nonZeroIndices[i]); + result += _storage.Values[i] * other.At(_storage.Indices[i]); } } @@ -925,17 +776,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The index of absolute minimum element. public override int AbsoluteMinimumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { // No non-zero elements. Return 0 return 0; } var index = 0; - var min = _nonZeroValues[index].Magnitude; - for (var i = 1; i < NonZerosCount; i++) + var min = _storage.Values[index].Magnitude; + for (var i = 1; i < _storage.ValueCount; i++) { - var test = _nonZeroValues[i].Magnitude; + var test = _storage.Values[i].Magnitude; if (test < min) { index = i; @@ -943,7 +794,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -1014,9 +865,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override Complex32 Sum() { var result = Complex32.Zero; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i]; + result += _storage.Values[i]; } return result; @@ -1029,9 +880,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override Complex32 SumMagnitudes() { var result = 0.0f; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i].Magnitude; + result += _storage.Values[i].Magnitude; } return result; @@ -1046,17 +897,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] *= _nonZeroValues[i]; + _storage.Values[i] *= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, other.At(index) * _nonZeroValues[i]); + var index = _storage.Indices[i]; + result.At(index, other.At(index) * _storage.Values[i]); } } } @@ -1070,17 +921,17 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] /= _nonZeroValues[i]; + _storage.Values[i] /= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, _nonZeroValues[i] / other.At(index)); + var index = _storage.Indices[i]; + result.At(index, _storage.Values[i] / other.At(index)); } } } @@ -1106,13 +957,13 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } var matrix = new SparseMatrix(u.Count, v.Count); - for (var i = 0; i < u.NonZerosCount; i++) + for (var i = 0; i < u._storage.ValueCount; i++) { - for (var j = 0; j < v.NonZerosCount; j++) + for (var j = 0; j < v._storage.ValueCount; j++) { - if (u._nonZeroIndices[i] == v._nonZeroIndices[j]) + if (u._storage.Indices[i] == v._storage.Indices[j]) { - matrix.At(i, j, u._nonZeroValues[i] * v._nonZeroValues[j]); + matrix.At(i, j, u._storage.Values[i] * v._storage.Values[j]); } } } @@ -1148,25 +999,25 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 throw new ArgumentOutOfRangeException("p"); } - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return Complex32.Zero; } if (2.0 == p) { - return _nonZeroValues.Aggregate(Complex32.Zero, SpecialFunctions.Hypotenuse).Magnitude; + return _storage.Values.Aggregate(Complex32.Zero, SpecialFunctions.Hypotenuse).Magnitude; } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(0, NonZerosCount, i => _nonZeroValues[i].Magnitude, Math.Max, 0f); + return CommonParallel.Aggregate(0, _storage.ValueCount, i => _storage.Values[i].Magnitude, Math.Max, 0f); } var sum = 0.0; - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - sum += Math.Pow(_nonZeroValues[index].Magnitude, p); + sum += Math.Pow(_storage.Values[index].Magnitude, p); } return (float)Math.Pow(sum, 1.0 / p); @@ -1332,97 +1183,33 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } #endregion - /// - /// Gets the value at the given index. - /// - /// Value real index in array - /// The value at the given index. - internal protected override Complex32 At(int index) - { - // Search if item idex exists in NonZeroIndices array in range "0 - real nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - return itemIndex >= 0 ? _nonZeroValues[itemIndex] : Complex32.Zero; - } - - /// - /// Delete, Add or Update the value in NonZeroValues and NonZeroIndices - /// - /// Value real index in array - /// The value to set. - /// This method assume that index is between 0 and Array Size - internal protected override void At(int index, Complex32 value) - { - // Search if "index" already exists in range "0 - complex nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - - if (itemIndex >= 0) - { - // Item already exist at itemIndex - if (value == Complex32.Zero) - { - RemoveAtUnchecked(itemIndex); - } - else - { - _nonZeroValues[itemIndex] = value; - } - } - else - { - if (value != Complex32.Zero) - { - InsertAtUnchecked(~itemIndex, index, value); - } - } - } - private void InsertAtUnchecked(int itemIndex, int index, Complex32 value) { // Check if the storage needs to be increased - if ((NonZerosCount == _nonZeroValues.Length) && (NonZerosCount < Count)) + if ((_storage.ValueCount == _storage.Values.Length) && (_storage.ValueCount < Count)) { // Value and Indices arrays are completely full so we increase the size - var size = Math.Min(_nonZeroValues.Length + GrowthSize(), Count); - Array.Resize(ref _nonZeroValues, size); - Array.Resize(ref _nonZeroIndices, size); + var size = Math.Min(_storage.Values.Length + GrowthSize(), Count); + Array.Resize(ref _storage.Values, size); + Array.Resize(ref _storage.Indices, size); } // Move all values (with a position larger than index) in the value array // to the next position // Move all values (with a position larger than index) in the columIndices // array to the next position - for (var i = NonZerosCount - 1; i > itemIndex - 1; i--) + for (var i = _storage.ValueCount - 1; i > itemIndex - 1; i--) { - _nonZeroValues[i + 1] = _nonZeroValues[i]; - _nonZeroIndices[i + 1] = _nonZeroIndices[i]; + _storage.Values[i + 1] = _storage.Values[i]; + _storage.Indices[i + 1] = _storage.Indices[i]; } // Add the value and the column index - _nonZeroValues[itemIndex] = value; - _nonZeroIndices[itemIndex] = index; + _storage.Values[itemIndex] = value; + _storage.Indices[itemIndex] = index; // increase the number of non-zero numbers by one - NonZerosCount += 1; - } - - private void RemoveAtUnchecked(int itemIndex) - { - // Value is zero. Let's delete it from Values and Indices array - for (var i = itemIndex + 1; i < NonZerosCount; i++) - { - _nonZeroValues[i - 1] = _nonZeroValues[i]; - _nonZeroIndices[i - 1] = _nonZeroIndices[i]; - } - - NonZerosCount -= 1; - - // Check whether we need to shrink the arrays. This is reasonable to do if - // there are a lot of non-zero elements and storage is two times bigger - if ((NonZerosCount > 1024) && (NonZerosCount < _nonZeroIndices.Length / 2)) - { - Array.Resize(ref _nonZeroValues, NonZerosCount); - Array.Resize(ref _nonZeroIndices, NonZerosCount); - } + _storage.ValueCount += 1; } /// @@ -1433,19 +1220,19 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 private int GrowthSize() { int delta; - if (_nonZeroValues.Length > 1024) + if (_storage.Values.Length > 1024) { - delta = _nonZeroValues.Length / 4; + delta = _storage.Values.Length / 4; } else { - if (_nonZeroValues.Length > 256) + if (_storage.Values.Length > 256) { delta = 512; } else { - delta = _nonZeroValues.Length > 64 ? 128 : 32; + delta = _storage.Values.Length > 64 ? 128 : 32; } } @@ -1458,7 +1245,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { if (Count > 20) { - return String.Format("SparseVectorOfComplex32({0},{1},{2})", Count, NonZerosCount, GetHashCode()); + return String.Format("SparseVectorOfComplex32({0},{1},{2})", Count, _storage.ValueCount, GetHashCode()); } return base.ToString(format, formatProvider); @@ -1472,14 +1259,14 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// public override int GetHashCode() { - var hashNum = Math.Min(NonZerosCount, 20); + var hashNum = Math.Min(_storage.ValueCount, 20); long hash = 0; for (var i = 0; i < hashNum; i++) { #if PORTABLE - hash ^= Precision.DoubleToInt64Bits(this._nonZeroValues[i].GetHashCode()); + hash ^= Precision.DoubleToInt64Bits(this._storage.Values[i].GetHashCode()); #else - hash ^= BitConverter.DoubleToInt64Bits(_nonZeroValues[i].GetHashCode()); + hash ^= BitConverter.DoubleToInt64Bits(_storage.Values[i].GetHashCode()); #endif } @@ -1523,27 +1310,27 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { - if (_nonZeroValues[i++] != Complex32.Zero) + if (_storage.Values[i++] != Complex32.Zero) { return false; } continue; } - if (i >= NonZerosCount || j < otherSparse.NonZerosCount && otherSparse._nonZeroIndices[j] < _nonZeroIndices[i]) + if (i >= _storage.ValueCount || j < otherSparse._storage.ValueCount && otherSparse._storage.Indices[j] < _storage.Indices[i]) { - if (otherSparse._nonZeroValues[j++] != Complex32.Zero) + if (otherSparse._storage.Values[j++] != Complex32.Zero) { return false; } continue; } - if (!_nonZeroValues[i].AlmostEqual(otherSparse._nonZeroValues[j])) + if (!_storage.Values[i].AlmostEqual(otherSparse._storage.Values[j])) { return false; } @@ -1571,9 +1358,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// public override IEnumerable> GetIndexedEnumerator() { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - yield return new Tuple(_nonZeroIndices[i], _nonZeroValues[i]); + yield return new Tuple(_storage.Indices[i], _storage.Values[i]); } } @@ -1586,9 +1373,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 public override Complex32[] ToArray() { var ret = new Complex32[Count]; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - ret[_nonZeroIndices[i]] = _nonZeroValues[i]; + ret[_storage.Indices[i]] = _storage.Values[i]; } return ret; diff --git a/src/Numerics/LinearAlgebra/Complex32/Vector.cs b/src/Numerics/LinearAlgebra/Complex32/Vector.cs index c76eecbd..3beb431e 100644 --- a/src/Numerics/LinearAlgebra/Complex32/Vector.cs +++ b/src/Numerics/LinearAlgebra/Complex32/Vector.cs @@ -30,6 +30,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 using Distributions; using Generic; using Properties; + using Storage; using Threading; using Complex32 = Numerics.Complex32; @@ -41,15 +42,9 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 { /// /// Initializes a new instance of the Vector class. - /// Constructs a Vector with the given size. /// - /// - /// The size of the Vector to construct. - /// - /// - /// If is less than one. - /// - protected Vector(int size) : base(size) + protected Vector(VectorStorage storage) + : base(storage) { } diff --git a/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs index f59fdc41..78caef7d 100644 --- a/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/DenseMatrix.cs @@ -39,8 +39,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double [Serializable] public class DenseMatrix : Matrix { - readonly DenseColumnMajorMatrixStorage _storage; - /// /// Number of rows. /// @@ -59,15 +57,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// Gets the matrix's data. /// /// The matrix's data. - readonly double[] _data; + readonly double[] _values; internal DenseMatrix(DenseColumnMajorMatrixStorage storage) : base(storage) { - _storage = storage; - _rowCount = _storage.RowCount; - _columnCount = _storage.ColumnCount; - _data = _storage.Data; + _rowCount = storage.RowCount; + _columnCount = storage.ColumnCount; + _values = storage.Data; } /// @@ -109,9 +106,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double public DenseMatrix(int rows, int columns, double value) : this(rows, columns) { - for (var i = 0; i < _data.Length; i++) + for (var i = 0; i < _values.Length; i++) { - _data[i] = value; + _values[i] = value; } } @@ -139,7 +136,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double { for (var j = 0; j < _columnCount; j++) { - _data[(j * _rowCount) + i] = array[i, j]; + _values[(j * _rowCount) + i] = array[i, j]; } } } @@ -161,7 +158,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The matrix's data. public double[] Data { - get { return _data; } + get { return _values; } } /// @@ -203,7 +200,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double var index = j * _rowCount; for (var i = 0; i < _rowCount; i++) { - ret._data[(i * _columnCount) + j] = _data[index + i]; + ret._values[(i * _columnCount) + j] = _values[index + i]; } } @@ -214,21 +211,21 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The L1 norm of the matrix. public override double L1Norm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _values); } /// Calculates the Frobenius norm of this matrix. /// The Frobenius norm of this matrix. public override double FrobeniusNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _values); } /// Calculates the infinity norm of this matrix. /// The infinity norm of this matrix. public override double InfinityNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _values); } #region Static constructors for special matrices. @@ -246,7 +243,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double var m = new DenseMatrix(order); for (var i = 0; i < order; i++) { - m._data[(i * order) + i] = 1.0; + m._values[(i * order) + i] = 1.0; } return m; @@ -271,7 +268,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.AddArrays(_values, denseOther._values, denseResult._values); } } @@ -290,7 +287,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.SubtractArrays(_values, denseOther._values, denseResult._values); } } @@ -308,7 +305,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult._values); } } @@ -332,14 +329,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0, - denseResult.Data); + denseResult.Values); } } @@ -363,14 +360,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0, - denseResult._data); + denseResult._values); } } @@ -394,14 +391,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.Transpose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0, - denseResult._data); + denseResult._values); } } @@ -425,14 +422,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0, - denseResult.Data); + denseResult.Values); } } @@ -456,14 +453,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0, - denseResult._data); + denseResult._values); } } @@ -481,7 +478,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(-1, _values, denseResult._values); } } @@ -501,7 +498,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_values, denseOther._values, denseResult._values); } } @@ -521,7 +518,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseDivideArrays(_values, denseOther._values, denseResult._values); } } @@ -547,8 +544,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double CommonParallel.For( 0, - _data.Length, - index => denseResult._data[index] %= divisor); + _values.Length, + index => denseResult._values[index] %= divisor); } } @@ -567,7 +564,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double var sum = 0.0; for (var i = 0; i < _rowCount; i++) { - sum += _data[(i * _rowCount) + i]; + sum += _values[(i * _rowCount) + i]; } return sum; diff --git a/src/Numerics/LinearAlgebra/Double/DenseVector.cs b/src/Numerics/LinearAlgebra/Double/DenseVector.cs index ecfbe523..5aded844 100644 --- a/src/Numerics/LinearAlgebra/Double/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Double/DenseVector.cs @@ -33,6 +33,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double using Generic; using NumberTheory; using Properties; + using Storage; using Threading; /// @@ -41,6 +42,23 @@ namespace MathNet.Numerics.LinearAlgebra.Double [Serializable] public class DenseVector : Vector { + /// + /// Number of elements + /// + readonly int _length; + + /// + /// Gets the vector's data. + /// + readonly double[] _values; + + internal DenseVector(DenseVectorStorage storage) + : base(storage) + { + _length = storage.Length; + _values = storage.Data; + } + /// /// Initializes a new instance of the class with a given size. /// @@ -51,9 +69,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// If is less than one. /// public DenseVector(int size) - : base(size) + : this(new DenseVectorStorage(size)) { - Data = new double[size]; } /// @@ -72,9 +89,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double public DenseVector(int size, double value) : this(size) { - for (var index = 0; index < Data.Length; index++) + for (var index = 0; index < _values.Length; index++) { - Data[index] = value; + _values[index] = value; } } @@ -88,31 +105,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double public DenseVector(Vector other) : this(other.Count) { - var vector = other as DenseVector; - if (vector == null) - { - CommonParallel.For( - 0, - Data.Length, - index => this[index] = other[index]); - } - else - { - Buffer.BlockCopy(vector.Data, 0, Data, 0, Data.Length * Constants.SizeOfDouble); - } - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public DenseVector(DenseVector other) - : this(other.Count) - { - Buffer.BlockCopy(other.Data, 0, Data, 0, Data.Length * Constants.SizeOfDouble); + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -121,20 +114,18 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The array to create this vector from. /// The vector does not copy the array, but keeps a reference to it. Any /// changes to the vector will also change the array. - public DenseVector(double[] array) : base(array.Length) + public DenseVector(double[] array) + : this(new DenseVectorStorage(array.Length, array)) { - Data = array; } /// - /// Gets the vector's internal data. + /// Gets the vector's data. /// - /// The vector's internal data. - /// Changing values in the array also changes the corresponding value in vector. Use with care. - internal double[] Data + /// The vector's data. + public double[] Values { - get; - private set; + get { return _values; } } /// @@ -145,14 +136,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// /// A reference to the internal date of the given vector. /// - public static implicit operator double[](DenseVector vector) + public static explicit operator double[](DenseVector vector) { if (vector == null) { throw new ArgumentNullException(); } - return vector.Data; + return vector.Values; } /// @@ -178,10 +169,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// This vector as a column matrix. public override Matrix ToColumnMatrix() { - var matrix = new DenseMatrix(Count, 1); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(_length, 1); + for (var i = 0; i < _values.Length; i++) { - matrix.At(i, 0, Data[i]); + matrix.At(i, 0, _values[i]); } return matrix; @@ -193,33 +184,15 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// This vector as a row matrix. public override Matrix ToRowMatrix() { - var matrix = new DenseMatrix(1, Count); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(1, _length); + for (var i = 0; i < _values.Length; i++) { - matrix.At(0, i, Data[i]); + matrix.At(0, i, _values[i]); } return matrix; } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - /// If is negative or - /// greater than the size of the vector. - public override double this[int index] - { - get - { - return Data[index]; - } - - set - { - Data[index] = value; - } - } - /// /// Creates a matrix with the given dimensions using the same storage type /// as this vector. @@ -253,49 +226,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double return new DenseVector(size); } - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as DenseVector; - if (otherVector == null) - { - CommonParallel.For( - 0, - Data.Length, - index => target[index] = Data[index]); - } - else - { - Buffer.BlockCopy(Data, 0, otherVector.Data, 0, Data.Length * Constants.SizeOfDouble); - } - } - /// /// Adds a scalar to each element of the vector and stores the result in the result vector. /// @@ -312,8 +242,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] + scalar); + _values.Length, + index => dense._values[index] = _values[index] + scalar); } } @@ -328,7 +258,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, 1.0, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, 1.0, odense.Values, rdense.Values); } else { @@ -397,8 +327,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] - scalar); + _values.Length, + index => dense._values[index] = _values[index] - scalar); } } @@ -413,7 +343,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, -1.0, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, -1.0, odense.Values, rdense.Values); } else { @@ -472,11 +402,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new DenseVector(Count); + var result = new DenseVector(_length); CommonParallel.For( 0, - Data.Length, - index => result[index] = -Data[index]); + _values.Length, + index => result[index] = -_values[index]); return result; } @@ -496,7 +426,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, Data, denseResult.Data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult.Values); } } @@ -510,7 +440,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double { var denseVector = other as DenseVector; - return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(Data, denseVector.Data); + return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(_values, denseVector.Values); } /// @@ -572,7 +502,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide"); } - return Control.LinearAlgebraProvider.DotProduct(leftSide.Data, rightSide.Data); + return Control.LinearAlgebraProvider.DotProduct(leftSide.Values, rightSide.Values); } /// @@ -602,16 +532,16 @@ namespace MathNet.Numerics.LinearAlgebra.Double var denseResult = result as DenseVector; if (denseResult == null) { - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - result.At(index, Data[index] % divisor); + result.At(index, _values[index] % divisor); } } else { - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - denseResult.Data[index] = Data[index] % divisor; + denseResult._values[index] = _values[index] % divisor; } } } @@ -640,10 +570,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override int AbsoluteMinimumIndex() { var index = 0; - var min = Math.Abs(Data[index]); - for (var i = 1; i < Count; i++) + var min = Math.Abs(_values[index]); + for (var i = 1; i < _length; i++) { - var test = Math.Abs(Data[i]); + var test = Math.Abs(_values[i]); if (test < min) { index = i; @@ -660,7 +590,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The value of the absolute minimum element. public override double AbsoluteMinimum() { - return Math.Abs(Data[AbsoluteMinimumIndex()]); + return Math.Abs(_values[AbsoluteMinimumIndex()]); } /// @@ -669,7 +599,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The value of the absolute maximum element. public override double AbsoluteMaximum() { - return Math.Abs(Data[AbsoluteMaximumIndex()]); + return Math.Abs(_values[AbsoluteMaximumIndex()]); } /// @@ -679,10 +609,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override int AbsoluteMaximumIndex() { var index = 0; - var max = Math.Abs(Data[index]); - for (var i = 1; i < Count; i++) + var max = Math.Abs(_values[index]); + for (var i = 1; i < _length; i++) { - var test = Math.Abs(Data[i]); + var test = Math.Abs(_values[i]); if (test > max) { index = i; @@ -706,7 +636,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// If is not positive. public override Vector SubVector(int index, int length) { - if (index < 0 || index >= Count) + if (index < 0 || index >= _length) { throw new ArgumentOutOfRangeException("index"); } @@ -716,7 +646,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentOutOfRangeException("length"); } - if (index + length > Count) + if (index + length > _length) { throw new ArgumentOutOfRangeException("length"); } @@ -725,8 +655,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double CommonParallel.For( index, - index + length, - i => result.Data[i - index] = Data[i]); + index + length, + i => result._values[i - index] = _values[i]); return result; } @@ -743,7 +673,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentNullException("values"); } - if (values.Length != Count) + if (values.Length != _length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "values"); } @@ -751,7 +681,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double CommonParallel.For( 0, values.Length, - i => Data[i] = values[i]); + i => _values[i] = values[i]); } /// @@ -761,13 +691,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override int MaximumIndex() { var index = 0; - var max = Data[0]; - for (var i = 1; i < Count; i++) + var max = _values[0]; + for (var i = 1; i < _length; i++) { - if (max < Data[i]) + if (max < _values[i]) { index = i; - max = Data[i]; + max = _values[i]; } } @@ -781,13 +711,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override int MinimumIndex() { var index = 0; - var min = Data[0]; - for (var i = 1; i < Count; i++) + var min = _values[0]; + for (var i = 1; i < _length; i++) { - if (min > Data[i]) + if (min > _values[i]) { index = i; - min = Data[i]; + min = _values[i]; } } @@ -802,9 +732,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double { var sum = 0.0; - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - sum += Data[index]; + sum += _values[index]; } return sum; @@ -818,9 +748,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double { var sum = 0.0; - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - sum += Math.Abs(Data[index]); + sum += Math.Abs(_values[index]); } return sum; @@ -842,8 +772,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] * other[index]); + _values.Length, + index => dense._values[index] = _values[index] * other[index]); } } @@ -864,8 +794,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] / other[index]); + _values.Length, + index => dense._values[index] = _values[index] / other[index]); } } @@ -897,7 +827,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double { for (var j = 0; j < v.Count; j++) { - matrix.At(i, j, u.Data[i] * v.Data[j]); + matrix.At(i, j, u._values[i] * v._values[j]); } }); return matrix; @@ -937,18 +867,18 @@ namespace MathNet.Numerics.LinearAlgebra.Double if (2.0 == p) { - return Data.Aggregate(0.0, SpecialFunctions.Hypotenuse); + return _values.Aggregate(0.0, SpecialFunctions.Hypotenuse); } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(Data, (i, v) => Math.Abs(v), Math.Max, 0d); + return CommonParallel.Aggregate(_values, (i, v) => Math.Abs(v), Math.Max, 0d); } var sum = 0.0; - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - sum += Math.Pow(Math.Abs(Data[index]), p); + sum += Math.Pow(Math.Abs(_values[index]), p); } return Math.Pow(sum, 1.0 / p); @@ -1114,29 +1044,5 @@ namespace MathNet.Numerics.LinearAlgebra.Double } #endregion - - /// - /// Resets all values to zero. - /// - public override void Clear() - { - Array.Clear(Data, 0, Data.Length); - } - - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - internal protected override double At(int index) - { - return Data[index]; - } - - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - internal protected override void At(int index, double value) - { - Data[index] = value; - } } } diff --git a/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs index de0461b3..32737edb 100644 --- a/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/DiagonalMatrix.cs @@ -377,12 +377,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double return; } - if (_data.Length != denseSource.Data.Length) + if (_data.Length != denseSource.Values.Length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "source"); } - Buffer.BlockCopy(denseSource.Data, 0, _data, 0, denseSource.Data.Length * Constants.SizeOfDouble); + Buffer.BlockCopy(denseSource.Values, 0, _data, 0, denseSource.Values.Length * Constants.SizeOfDouble); } /// diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseCholesky.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseCholesky.cs index 60fcedd5..a75a3dc5 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseCholesky.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseCholesky.cs @@ -167,11 +167,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } // Copy the contents of input to result. - Buffer.BlockCopy(dinput.Data, 0, dresult.Data, 0, dinput.Data.Length * Constants.SizeOfDouble); + Buffer.BlockCopy(dinput.Values, 0, dresult.Values, 0, dinput.Values.Length * Constants.SizeOfDouble); // Cholesky solve by overwriting result. var dfactor = (DenseMatrix)CholeskyFactor; - Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Data, 1); + Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Values, 1); } } } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs index 63d6da6c..377001fa 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseGramSchmidt.cs @@ -216,7 +216,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Data, 1, dresult.Data); + _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseLU.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseLU.cs index 12cdf536..e71fc58e 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseLU.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseLU.cs @@ -168,11 +168,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization } // Copy the contents of input to result. - Buffer.BlockCopy(dinput.Data, 0, dresult.Data, 0, dinput.Data.Length * Constants.SizeOfDouble); + Buffer.BlockCopy(dinput.Values, 0, dresult.Values, 0, dinput.Values.Length * Constants.SizeOfDouble); // LU solve by overwriting result. var dfactors = (DenseMatrix)Factors; - Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Data); + Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Values); } /// diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs index 3cb58563..fca7385e 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseQR.cs @@ -188,7 +188,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Double/Factorization/DenseSvd.cs b/src/Numerics/LinearAlgebra/Double/Factorization/DenseSvd.cs index a0440f5b..7efd47fc 100644 --- a/src/Numerics/LinearAlgebra/Double/Factorization/DenseSvd.cs +++ b/src/Numerics/LinearAlgebra/Double/Factorization/DenseSvd.cs @@ -69,7 +69,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization VectorS = new DenseVector(nm); MatrixU = new DenseMatrix(matrix.RowCount); MatrixVT = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); + Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); } /// @@ -125,7 +125,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do SVD factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); } /// @@ -175,7 +175,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double.Factorization throw new NotSupportedException("Can only do SVD factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs index 24b83055..335847ec 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseMatrix.cs @@ -97,6 +97,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The number of columns. /// /// The value which we assign to each element of the matrix. + [Obsolete("Use a dense matrix instead.")] public SparseMatrix(int rows, int columns, double value) : this(rows, columns) { diff --git a/src/Numerics/LinearAlgebra/Double/SparseVector.cs b/src/Numerics/LinearAlgebra/Double/SparseVector.cs index 81afccc6..4d6d7c72 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseVector.cs @@ -37,6 +37,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double using Generic; using NumberTheory; using Properties; + using Storage; using Threading; /// @@ -46,24 +47,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double [Serializable] public class SparseVector : Vector { - /// - /// Gets the vector's internal data. The array containing the actual values; only the non-zero values are stored. - /// - private double[] _nonZeroValues = new double[0]; - - /// - /// The indices of the non-zero entries. - /// - private int[] _nonZeroIndices = new int[0]; + readonly SparseVectorStorage _storage; - /// - /// Gets the number of non zero elements in the vector. - /// - /// The number of non zero elements. - public int NonZerosCount + internal SparseVector(SparseVectorStorage storage) + : base(storage) { - get; - private set; + _storage = storage; } #region Constructors @@ -77,7 +66,8 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// /// If is less than one. /// - public SparseVector(int size) : base(size) + public SparseVector(int size) + : this(new SparseVectorStorage(size)) { } @@ -94,55 +84,23 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// /// If is less than one. /// - public SparseVector(int size, double value) : this(size) + [Obsolete("Use a dense vector instead.")] + public SparseVector(int size, double value) + : this(new SparseVectorStorage(size)) { if (value == 0.0) { - // Skip adding values return; } - // We already know that this vector is "full", let's allocate all needed memory - _nonZeroValues = new double[size]; - _nonZeroIndices = new int[size]; - NonZerosCount = size; - - CommonParallel.For( - 0, - Count, - index => - { - _nonZeroValues[index] = value; - _nonZeroIndices[index] = index; - }); - } + var valueCount = _storage.ValueCount = size; + var indices = _storage.Indices = new int[valueCount]; + var values = _storage.Values = new double[valueCount]; - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public SparseVector(Vector other) : this(other.Count) - { - var vector = other as SparseVector; - if (vector == null) + for (int i = 0; i < values.Length; i++) { - for (var i = 0; i < other.Count; i++) - { - this[i] = other.At(i); - } - } - else - { - _nonZeroValues = new double[vector.NonZerosCount]; - _nonZeroIndices = new int[vector.NonZerosCount]; - NonZerosCount = vector.NonZerosCount; - - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - Buffer.BlockCopy(vector._nonZeroValues, 0, _nonZeroValues, 0, vector.NonZerosCount * Constants.SizeOfDouble); - Buffer.BlockCopy(vector._nonZeroIndices, 0, _nonZeroIndices, 0, vector.NonZerosCount * Constants.SizeOfInt); + values[i] = value; + indices[i] = i; } } @@ -153,15 +111,10 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// /// The vector to create the new vector from. /// - public SparseVector(SparseVector other) : this(other.Count) + public SparseVector(Vector other) + : this(new SparseVectorStorage(other.Count)) { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - _nonZeroValues = new double[other.NonZerosCount]; - _nonZeroIndices = new int[other.NonZerosCount]; - NonZerosCount = other.NonZerosCount; - - Buffer.BlockCopy(other._nonZeroValues, 0, _nonZeroValues, 0, other.NonZerosCount * Constants.SizeOfDouble); - Buffer.BlockCopy(other._nonZeroIndices, 0, _nonZeroIndices, 0, other.NonZerosCount * Constants.SizeOfInt); + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -169,11 +122,12 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// /// The array to create this vector from. /// The vector copy the array. Any changes to the vector will NOT change the array. - public SparseVector(IList array) : this(array.Count) + public SparseVector(IList array) + : this(new SparseVectorStorage(array.Count)) { for (var i = 0; i < array.Count; i++) { - this[i] = array[i]; + Storage.At(i, array[i]); } } @@ -185,10 +139,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// This vector as a column matrix. public override Matrix ToColumnMatrix() { + var indices = _storage.Indices; + var values = _storage.Values; + var matrix = new SparseMatrix(Count, 1); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(_nonZeroIndices[i], 0, _nonZeroValues[i]); + matrix.At(indices[i], 0, values[i]); } return matrix; @@ -200,10 +157,13 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// This vector as a row matrix. public override Matrix ToRowMatrix() { + var indices = _storage.Indices; + var values = _storage.Values; + var matrix = new SparseMatrix(1, Count); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(0, _nonZeroIndices[i], _nonZeroValues[i]); + matrix.At(0, indices[i], values[i]); } return matrix; @@ -242,68 +202,6 @@ namespace MathNet.Numerics.LinearAlgebra.Double return new SparseVector(size); } - /// - /// Resets all values to zero. - /// - public override void Clear() - { - NonZerosCount = 0; - } - - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as SparseVector; - if (otherVector == null) - { - target.Clear(); - - for (var index = 0; index < NonZerosCount; index++) - { - target.At(_nonZeroIndices[index], _nonZeroValues[index]); - } - } - else - { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - otherVector._nonZeroValues = new double[NonZerosCount]; - otherVector._nonZeroIndices = new int[NonZerosCount]; - otherVector.NonZerosCount = NonZerosCount; - - if (NonZerosCount != 0) - { - Buffer.BlockCopy(_nonZeroValues, 0, otherVector._nonZeroValues, 0, NonZerosCount*Constants.SizeOfDouble); - Buffer.BlockCopy(_nonZeroIndices, 0, otherVector._nonZeroIndices, 0, NonZerosCount*Constants.SizeOfInt); - } - } - } - #region Operators and supplementary functions /// @@ -340,16 +238,19 @@ namespace MathNet.Numerics.LinearAlgebra.Double vnonZeroValues[index] = scalar; } + var indices = _storage.Indices; + var values = _storage.Values; + //populate the non zero values from this - for (int j = 0; j < NonZerosCount; j++) + for (int j = 0; j < _storage.ValueCount; j++) { - vnonZeroValues[_nonZeroIndices[j]] = _nonZeroValues[j] + scalar; + vnonZeroValues[indices[j]] = values[j] + scalar; } - //assign this vectors arrary to the new arrays. - _nonZeroValues = vnonZeroValues; - _nonZeroIndices = vnonZeroIndices; - NonZerosCount = Count; + //assign this vectors arrary to the new arrays. + _storage.Values = vnonZeroValues; + _storage.Indices = vnonZeroIndices; + _storage.ValueCount = Count; } else { @@ -390,22 +291,22 @@ namespace MathNet.Numerics.LinearAlgebra.Double if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] += otherSparse._nonZeroValues[j++]; + _storage.Values[i++] += otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != 0.0) { - InsertAtIndexUnchecked(i++, otherSparse._nonZeroIndices[j], otherValue); + InsertAtIndexUnchecked(i++, otherSparse._storage.Indices[j], otherValue); } j++; } @@ -415,25 +316,25 @@ namespace MathNet.Numerics.LinearAlgebra.Double { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] + otherSparse.At(next)); + result.At(next, _storage.Values[i] + otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) + otherSparse._nonZeroValues[j]); + result.At(next, At(next) + otherSparse._storage.Values[j]); } j++; } @@ -536,22 +437,22 @@ namespace MathNet.Numerics.LinearAlgebra.Double if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] -= otherSparse._nonZeroValues[j++]; + _storage.Values[i++] -= otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != 0.0) { - InsertAtIndexUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue); + InsertAtIndexUnchecked(i++, otherSparse._storage.Indices[j], -otherValue); } j++; } @@ -561,25 +462,25 @@ namespace MathNet.Numerics.LinearAlgebra.Double { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] - otherSparse.At(next)); + result.At(next, _storage.Values[i] - otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) - otherSparse._nonZeroValues[j]); + result.At(next, At(next) - otherSparse._storage.Values[j]); } j++; } @@ -638,24 +539,18 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new SparseVector(Count) - { - _nonZeroValues = new double[NonZerosCount], - _nonZeroIndices = new int[NonZerosCount], - NonZerosCount = NonZerosCount - }; + var result = new SparseVectorStorage(Count); + var valueCount = result.ValueCount = _storage.ValueCount; + var indices = result.Indices = new int[valueCount]; + var values = result.Values = new double[valueCount]; - if (NonZerosCount != 0) + if (valueCount != 0) { - CommonParallel.For( - 0, - NonZerosCount, - index => result._nonZeroValues[index] = -_nonZeroValues[index]); - - Buffer.BlockCopy(_nonZeroIndices, 0, result._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); + CommonParallel.For(0, valueCount, index => values[index] = -_storage.Values[index]); + Buffer.BlockCopy(_storage.Indices, 0, indices, 0, valueCount * Constants.SizeOfInt); } - return result; + return new SparseVector(result); } /// @@ -676,7 +571,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } var copy = new SparseVector(this); - Control.LinearAlgebraProvider.ScaleArray(scalar, copy._nonZeroValues, copy._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, copy._storage.Values, copy._storage.Values); return copy; } @@ -711,23 +606,23 @@ namespace MathNet.Numerics.LinearAlgebra.Double if (sparseResult == null) { result.Clear(); - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - result.At(_nonZeroIndices[index], scalar * _nonZeroValues[index]); + result.At(_storage.Indices[index], scalar * _storage.Values[index]); } } else { if (!ReferenceEquals(this, result)) { - sparseResult.NonZerosCount = NonZerosCount; - sparseResult._nonZeroIndices = new int[NonZerosCount]; - Buffer.BlockCopy(_nonZeroIndices, 0, sparseResult._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); - sparseResult._nonZeroValues = new double[NonZerosCount]; - Buffer.BlockCopy(_nonZeroValues, 0, sparseResult._nonZeroValues, 0, NonZerosCount * Constants.SizeOfDouble); + sparseResult._storage.ValueCount = _storage.ValueCount; + sparseResult._storage.Indices = new int[_storage.ValueCount]; + Buffer.BlockCopy(_storage.Indices, 0, sparseResult._storage.Indices, 0, _storage.ValueCount * Constants.SizeOfInt); + sparseResult._storage.Values = new double[_storage.ValueCount]; + Buffer.BlockCopy(_storage.Values, 0, sparseResult._storage.Values, 0, _storage.ValueCount * Constants.SizeOfDouble); } - Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._nonZeroValues, sparseResult._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._storage.Values, sparseResult._storage.Values); } } @@ -746,16 +641,16 @@ namespace MathNet.Numerics.LinearAlgebra.Double if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * _nonZeroValues[i]; + result += _storage.Values[i] * _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * other.At(_nonZeroIndices[i]); + result += _storage.Values[i] * other.At(_storage.Indices[i]); } } @@ -850,17 +745,17 @@ namespace MathNet.Numerics.LinearAlgebra.Double { if (ReferenceEquals(this, result)) { - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - _nonZeroValues[index] %= divisor; + _storage.Values[index] %= divisor; } } else { result.Clear(); - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - result.At(_nonZeroIndices[index], _nonZeroValues[index] % divisor); + result.At(_storage.Indices[index], _storage.Values[index] % divisor); } } } @@ -888,17 +783,17 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The index of absolute minimum element. public override int AbsoluteMinimumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { // No non-zero elements. Return 0 return 0; } var index = 0; - var min = Math.Abs(_nonZeroValues[index]); - for (var i = 1; i < NonZerosCount; i++) + var min = Math.Abs(_storage.Values[index]); + for (var i = 1; i < _storage.ValueCount; i++) { - var test = Math.Abs(_nonZeroValues[i]); + var test = Math.Abs(_storage.Values[i]); if (test < min) { index = i; @@ -906,7 +801,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -976,23 +871,23 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The index of absolute maximum element. public override int MaximumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0; } var index = 0; - var max = _nonZeroValues[0]; - for (var i = 1; i < NonZerosCount; i++) + var max = _storage.Values[0]; + for (var i = 1; i < _storage.ValueCount; i++) { - if (max < _nonZeroValues[i]) + if (max < _storage.Values[i]) { index = i; - max = _nonZeroValues[i]; + max = _storage.Values[i]; } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -1001,23 +896,23 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The index of minimum element. public override int MinimumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0; } var index = 0; - var min = _nonZeroValues[0]; - for (var i = 1; i < NonZerosCount; i++) + var min = _storage.Values[0]; + for (var i = 1; i < _storage.ValueCount; i++) { - if (min > _nonZeroValues[i]) + if (min > _storage.Values[i]) { index = i; - min = _nonZeroValues[i]; + min = _storage.Values[i]; } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -1027,9 +922,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override double Sum() { double result = 0; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i]; + result += _storage.Values[i]; } return result; @@ -1042,9 +937,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override double SumMagnitudes() { double result = 0; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += Math.Abs(_nonZeroValues[i]); + result += Math.Abs(_storage.Values[i]); } return result; @@ -1059,17 +954,17 @@ namespace MathNet.Numerics.LinearAlgebra.Double { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] *= _nonZeroValues[i]; + _storage.Values[i] *= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, other.At(index) * _nonZeroValues[i]); + var index = _storage.Indices[i]; + result.At(index, other.At(index) * _storage.Values[i]); } } } @@ -1083,17 +978,17 @@ namespace MathNet.Numerics.LinearAlgebra.Double { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] /= _nonZeroValues[i]; + _storage.Values[i] /= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, _nonZeroValues[i] / other.At(index)); + var index = _storage.Indices[i]; + result.At(index, _storage.Values[i] / other.At(index)); } } } @@ -1119,11 +1014,11 @@ namespace MathNet.Numerics.LinearAlgebra.Double } var matrix = new SparseMatrix(u.Count, v.Count); - for (var i = 0; i < u.NonZerosCount; i++) + for (var i = 0; i < u._storage.ValueCount; i++) { - for (var j = 0; j < v.NonZerosCount; j++) + for (var j = 0; j < v._storage.ValueCount; j++) { - matrix.At(i, j, u._nonZeroValues[i] * v._nonZeroValues[j]); + matrix.At(i, j, u._storage.Values[i] * v._storage.Values[j]); } } @@ -1158,25 +1053,25 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentOutOfRangeException("p"); } - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0.0; } if (2.0 == p) { - return _nonZeroValues.Aggregate(0.0, SpecialFunctions.Hypotenuse); + return _storage.Values.Aggregate(0.0, SpecialFunctions.Hypotenuse); } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(0, NonZerosCount, i => Math.Abs(_nonZeroValues[i]), Math.Max, 0d); + return CommonParallel.Aggregate(0, _storage.ValueCount, i => Math.Abs(_storage.Values[i]), Math.Max, 0d); } var sum = 0.0; - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - sum += Math.Pow(Math.Abs(_nonZeroValues[index]), p); + sum += Math.Pow(Math.Abs(_storage.Values[index]), p); } return Math.Pow(sum, 1.0 / p); @@ -1343,97 +1238,33 @@ namespace MathNet.Numerics.LinearAlgebra.Double #endregion - /// - /// Gets the value at the given index. - /// - /// Value real index in array - /// The value at the given index. - internal protected override double At(int index) - { - // Search if item idex exists in NonZeroIndices array in range "0 - real nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - return itemIndex >= 0 ? _nonZeroValues[itemIndex] : 0.0; - } - - /// - /// Delete, Add or Update the value in NonZeroValues and NonZeroIndices - /// - /// Value real index in array - /// The value to set. - /// This method assume that index is between 0 and Array Size - internal protected override void At(int index, double value) - { - // Search if "index" already exists in range "0 - real nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - - if (itemIndex >= 0) - { - // Item already exist at itemIndex - if (value == 0.0) - { - RemoveAtIndexUnchecked(itemIndex); - } - else - { - _nonZeroValues[itemIndex] = value; - } - } - else - { - if (value != 0.0) - { - InsertAtIndexUnchecked(~itemIndex, index, value); - } - } - } - private void InsertAtIndexUnchecked(int itemIndex, int index, double value) { // Check if the storage needs to be increased - if ((NonZerosCount == _nonZeroValues.Length) && (NonZerosCount < Count)) + if ((_storage.ValueCount == _storage.Values.Length) && (_storage.ValueCount < Count)) { // Value and Indices arrays are completely full so we increase the size - var size = Math.Min(_nonZeroValues.Length + GrowthSize(), Count); - Array.Resize(ref _nonZeroValues, size); - Array.Resize(ref _nonZeroIndices, size); + var size = Math.Min(_storage.Values.Length + GrowthSize(), Count); + Array.Resize(ref _storage.Values, size); + Array.Resize(ref _storage.Indices, size); } // Move all values (with a position larger than index) in the value array // to the next position // Move all values (with a position larger than index) in the columIndices // array to the next position - for (var i = NonZerosCount - 1; i > itemIndex - 1; i--) + for (var i = _storage.ValueCount - 1; i > itemIndex - 1; i--) { - _nonZeroValues[i + 1] = _nonZeroValues[i]; - _nonZeroIndices[i + 1] = _nonZeroIndices[i]; + _storage.Values[i + 1] = _storage.Values[i]; + _storage.Indices[i + 1] = _storage.Indices[i]; } // Add the value and the column index - _nonZeroValues[itemIndex] = value; - _nonZeroIndices[itemIndex] = index; + _storage.Values[itemIndex] = value; + _storage.Indices[itemIndex] = index; // increase the number of non-zero numbers by one - NonZerosCount += 1; - } - - private void RemoveAtIndexUnchecked(int itemIndex) - { - // Value is zero. Let's delete it from Values and Indices array - for (var i = itemIndex + 1; i < NonZerosCount; i++) - { - _nonZeroValues[i - 1] = _nonZeroValues[i]; - _nonZeroIndices[i - 1] = _nonZeroIndices[i]; - } - - NonZerosCount -= 1; - - // Check whether we need to shrink the arrays. This is reasonable to do if - // there are a lot of non-zero elements and storage is two times bigger - if ((NonZerosCount > 1024) && (NonZerosCount < _nonZeroIndices.Length / 2)) - { - Array.Resize(ref _nonZeroValues, NonZerosCount); - Array.Resize(ref _nonZeroIndices, NonZerosCount); - } + _storage.ValueCount += 1; } /// @@ -1444,19 +1275,19 @@ namespace MathNet.Numerics.LinearAlgebra.Double private int GrowthSize() { int delta; - if (_nonZeroValues.Length > 1024) + if (_storage.Values.Length > 1024) { - delta = _nonZeroValues.Length / 4; + delta = _storage.Values.Length / 4; } else { - if (_nonZeroValues.Length > 256) + if (_storage.Values.Length > 256) { delta = 512; } else { - delta = _nonZeroValues.Length > 64 ? 128 : 32; + delta = _storage.Values.Length > 64 ? 128 : 32; } } @@ -1469,7 +1300,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double { if (Count > 20) { - return String.Format("SparseVectorOfDouble({0},{1},{2})", Count, NonZerosCount, GetHashCode()); + return String.Format("SparseVectorOfDouble({0},{1},{2})", Count, _storage.ValueCount, GetHashCode()); } return base.ToString(format, formatProvider); @@ -1483,14 +1314,14 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// public override int GetHashCode() { - var hashNum = Math.Min(NonZerosCount, 20); + var hashNum = Math.Min(_storage.ValueCount, 20); long hash = 0; for (var i = 0; i < hashNum; i++) { #if PORTABLE - hash ^= Precision.DoubleToInt64Bits(this._nonZeroValues[i].GetHashCode()); + hash ^= Precision.DoubleToInt64Bits(this._storage.Values[i].GetHashCode()); #else - hash ^= BitConverter.DoubleToInt64Bits(_nonZeroValues[i].GetHashCode()); + hash ^= BitConverter.DoubleToInt64Bits(_storage.Values[i].GetHashCode()); #endif } @@ -1534,27 +1365,27 @@ namespace MathNet.Numerics.LinearAlgebra.Double } int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { - if (_nonZeroValues[i++] != 0d) + if (_storage.Values[i++] != 0d) { return false; } continue; } - if (i >= NonZerosCount || j < otherSparse.NonZerosCount && otherSparse._nonZeroIndices[j] < _nonZeroIndices[i]) + if (i >= _storage.ValueCount || j < otherSparse._storage.ValueCount && otherSparse._storage.Indices[j] < _storage.Indices[i]) { - if (otherSparse._nonZeroValues[j++] != 0d) + if (otherSparse._storage.Values[j++] != 0d) { return false; } continue; } - if (!_nonZeroValues[i].AlmostEqual(otherSparse._nonZeroValues[j])) + if (!_storage.Values[i].AlmostEqual(otherSparse._storage.Values[j])) { return false; } @@ -1582,9 +1413,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// public override IEnumerable> GetIndexedEnumerator() { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - yield return new Tuple(_nonZeroIndices[i], _nonZeroValues[i]); + yield return new Tuple(_storage.Indices[i], _storage.Values[i]); } } @@ -1597,9 +1428,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double public override double[] ToArray() { var ret = new double[Count]; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - ret[_nonZeroIndices[i]] = _nonZeroValues[i]; + ret[_storage.Indices[i]] = _storage.Values[i]; } return ret; diff --git a/src/Numerics/LinearAlgebra/Double/Vector.cs b/src/Numerics/LinearAlgebra/Double/Vector.cs index 9098a52d..b873887e 100644 --- a/src/Numerics/LinearAlgebra/Double/Vector.cs +++ b/src/Numerics/LinearAlgebra/Double/Vector.cs @@ -30,6 +30,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double using Distributions; using Generic; using Properties; + using Storage; using Threading; /// @@ -40,15 +41,9 @@ namespace MathNet.Numerics.LinearAlgebra.Double { /// /// Initializes a new instance of the Vector class. - /// Constructs a Vector with the given size. /// - /// - /// The size of the Vector to construct. - /// - /// - /// If is less than one. - /// - protected Vector(int size) : base(size) + protected Vector(VectorStorage storage) + : base(storage) { } diff --git a/src/Numerics/LinearAlgebra/Generic/Vector.cs b/src/Numerics/LinearAlgebra/Generic/Vector.cs index 91948081..370793c6 100644 --- a/src/Numerics/LinearAlgebra/Generic/Vector.cs +++ b/src/Numerics/LinearAlgebra/Generic/Vector.cs @@ -34,6 +34,7 @@ namespace MathNet.Numerics.LinearAlgebra.Generic using Distributions; using Numerics; using Properties; + using Storage; using Threading; /// @@ -60,60 +61,152 @@ namespace MathNet.Numerics.LinearAlgebra.Generic private static readonly T One = Common.SetOne(); /// - /// Initializes a new instance of the Vector class. - /// Constructs a Vector with the given size. + /// Initializes a new instance of the Vector class. /// - /// - /// The size of the Vector to construct. - /// - /// - /// If is less than one. - /// - protected Vector(int size) + protected Vector(VectorStorage storage) { - if (size < 1) - { - throw new ArgumentOutOfRangeException("size", Resources.ArgumentMustBePositive); - } - - Count = size; + Storage = storage; + Count = storage.Length; } /// - /// Gets he number of elements in the vector. + /// Gets the raw vector data storage. /// - public int Count - { - get; - private set; - } + public VectorStorage Storage { get; private set; } + + /// + /// Gets the number of items. + /// + public int Count { get; private set; } /// Gets or sets the value at the given . /// The index of the value to get or set. /// The value of the vector at the given . /// If is negative or /// greater than the size of the vector. - public virtual T this[int index] + public T this[int index] { - get + get { return Storage[index]; } + set { Storage[index] = value;} + } + + /// Gets the value at the given without range checking.. + /// The index of the value to get or set. + /// The value of the vector at the given . + public T At(int index) + { + return Storage.At(index); + } + + /// Sets the at the given without range checking.. + /// The index of the value to get or set. + /// The value to set. + public void At(int index, T value) + { + Storage.At(index, value); + } + + /// + /// Resets all values to zero. + /// + public void Clear() + { + Storage.Clear(); + } + + /// + /// Sets all values of a subvector to zero. + /// + public void ClearSubVector(int index, int count) + { + if (count < 1) { - if (index < 0 || index >= Count) - { - throw new ArgumentOutOfRangeException("index"); - } + throw new ArgumentOutOfRangeException("count", Resources.ArgumentMustBePositive); + } - return At(index); + if (index + count > Count || index < 0) + { + throw new ArgumentOutOfRangeException("index"); } - set + Storage.Clear(index, count); + } + + /// + /// Returns a deep-copy clone of the vector. + /// + /// + /// A deep-copy clone of the vector. + /// + public Vector Clone() + { + var result = CreateVector(Count); + Storage.CopyTo(result.Storage, skipClearing: true); + return result; + } + + /// + /// Copies the values of this vector into the target vector. + /// + /// + /// The vector to copy elements into. + /// + /// + /// If is . + /// + /// + /// If is not the same size as this vector. + /// + public void CopyTo(Vector target) + { + if (target == null) { - if (index < 0 || index >= Count) - { - throw new ArgumentOutOfRangeException("index"); - } + throw new ArgumentNullException("target"); + } + + if (ReferenceEquals(this, target) || ReferenceEquals(Storage, target.Storage)) + { + return; + } - At(index, value); + if (Count != target.Count) + { + throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); + } + + Storage.CopyTo(target.Storage); + } + + /// + /// Copies the requested elements from this vector to another. + /// + /// + /// The vector to copy the elements to. + /// + /// + /// The element to start copying from. + /// + /// + /// The element to start copying to. + /// + /// + /// The number of elements to copy. + /// + public void CopySubVectorTo(Vector destination, int sourceIndex, int targetIndex, int count) + { + if (destination == null) + { + throw new ArgumentNullException("destination"); } + + // TODO: refactor range checks + Storage.CopySubVectorTo(destination.Storage, sourceIndex, targetIndex, count); + } + + [Obsolete("Use CopySubVectorTo instead.")] + public void CopyTo(Vector destination, int sourceIndex, int targetIndex, int count) + { + CopySubVectorTo(destination, sourceIndex, targetIndex, count); } /// @@ -1166,115 +1259,6 @@ namespace MathNet.Numerics.LinearAlgebra.Generic #region Copying and Conversion - /// - /// Returns a deep-copy clone of the vector. - /// - /// - /// A deep-copy clone of the vector. - /// - public Vector Clone() - { - var retrunVector = CreateVector(Count); - CopyTo(retrunVector); - return retrunVector; - } - - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public virtual void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - CommonParallel.For( - 0, - Count, - index => target[index] = this[index]); - } - - /// - /// Copies the requested elements from this vector to another. - /// - /// - /// The vector to copy the elements to. - /// - /// - /// The element to start copying from. - /// - /// - /// The element to start copying to. - /// - /// - /// The number of elements to copy. - /// - public virtual void CopyTo(Vector destination, int offset, int destinationOffset, int count) - { - if (destination == null) - { - throw new ArgumentNullException("destination"); - } - - if (offset >= Count) - { - throw new ArgumentOutOfRangeException("offset"); - } - - if (offset + count > Count) - { - throw new ArgumentOutOfRangeException("count"); - } - - if (destinationOffset >= destination.Count) - { - throw new ArgumentOutOfRangeException("destinationOffset"); - } - - if (destinationOffset + count > destination.Count) - { - throw new ArgumentOutOfRangeException("count"); - } - - if (ReferenceEquals(this, destination)) - { - var tmpVector = destination.CreateVector(destination.Count); - CopyTo(tmpVector); - - CommonParallel.For( - 0, - count, - index => destination[destinationOffset + index] = tmpVector[offset + index]); - } - else - { - CommonParallel.For( - 0, - count, - index => destination[destinationOffset + index] = this[offset + index]); - } - } - /// /// Returns the data contained in the vector as an array. /// @@ -1605,23 +1589,5 @@ namespace MathNet.Numerics.LinearAlgebra.Generic } #endregion - - /// - /// Resets all values to zero. - /// - public virtual void Clear() - { - CommonParallel.For(0, Count, index => this[index] = default(T)); - } - - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - internal protected abstract T At(int index); - - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - internal protected abstract void At(int index, T value); } } diff --git a/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs b/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs index a76560d8..726d742f 100644 --- a/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/DenseMatrix.cs @@ -39,8 +39,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single [Serializable] public class DenseMatrix : Matrix { - readonly DenseColumnMajorMatrixStorage _storage; - /// /// Number of rows. /// @@ -59,15 +57,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// Gets the matrix's data. /// /// The matrix's data. - readonly float[] _data; + readonly float[] _values; internal DenseMatrix(DenseColumnMajorMatrixStorage storage) : base(storage) { - _storage = storage; - _rowCount = _storage.RowCount; - _columnCount = _storage.ColumnCount; - _data = _storage.Data; + _rowCount = storage.RowCount; + _columnCount = storage.ColumnCount; + _values = storage.Data; } /// @@ -109,9 +106,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single public DenseMatrix(int rows, int columns, float value) : this(rows, columns) { - for (var i = 0; i < _data.Length; i++) + for (var i = 0; i < _values.Length; i++) { - _data[i] = value; + _values[i] = value; } } @@ -139,7 +136,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single { for (var j = 0; j < _columnCount; j++) { - _data[(j * _rowCount) + i] = array[i, j]; + _values[(j * _rowCount) + i] = array[i, j]; } } } @@ -161,7 +158,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The matrix's data. public float[] Data { - get { return _data; } + get { return _values; } } /// @@ -203,7 +200,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single var index = j * _rowCount; for (var i = 0; i < _rowCount; i++) { - ret._data[(i * _columnCount) + j] = _data[index + i]; + ret._values[(i * _columnCount) + j] = _values[index + i]; } } @@ -214,21 +211,21 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The L1 norm of the matrix. public override float L1Norm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.OneNorm, _rowCount, _columnCount, _values); } /// Calculates the Frobenius norm of this matrix. /// The Frobenius norm of this matrix. public override float FrobeniusNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.FrobeniusNorm, _rowCount, _columnCount, _values); } /// Calculates the infinity norm of this matrix. /// The infinity norm of this matrix. public override float InfinityNorm() { - return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _data); + return Control.LinearAlgebraProvider.MatrixNorm(Norm.InfinityNorm, _rowCount, _columnCount, _values); } #region Static constructors for special matrices. @@ -246,7 +243,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single var m = new DenseMatrix(order); for (var i = 0; i < order; i++) { - m._data[(i * order) + i] = 1.0f; + m._values[(i * order) + i] = 1.0f; } return m; @@ -271,7 +268,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.AddArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.AddArrays(_values, denseOther._values, denseResult._values); } } @@ -290,7 +287,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.SubtractArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.SubtractArrays(_values, denseOther._values, denseResult._values); } } @@ -308,7 +305,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult._values); } } @@ -332,14 +329,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0f, - denseResult.Data); + denseResult.Values); } } @@ -363,14 +360,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0f, - denseResult._data); + denseResult._values); } } @@ -394,14 +391,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single Algorithms.LinearAlgebra.Transpose.DontTranspose, Algorithms.LinearAlgebra.Transpose.Transpose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0f, - denseResult._data); + denseResult._values); } } @@ -425,14 +422,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseRight.Data, + denseRight.Values, denseRight.Count, 1, 0.0f, - denseResult.Data); + denseResult.Values); } } @@ -456,14 +453,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single Algorithms.LinearAlgebra.Transpose.Transpose, Algorithms.LinearAlgebra.Transpose.DontTranspose, 1.0f, - _data, + _values, _rowCount, _columnCount, - denseOther._data, + denseOther._values, denseOther._rowCount, denseOther._columnCount, 0.0f, - denseResult._data); + denseResult._values); } } @@ -481,7 +478,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.ScaleArray(-1, _data, denseResult._data); + Control.LinearAlgebraProvider.ScaleArray(-1, _values, denseResult._values); } } @@ -501,7 +498,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseMultiplyArrays(_values, denseOther._values, denseResult._values); } } @@ -521,7 +518,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.PointWiseDivideArrays(_data, denseOther._data, denseResult._data); + Control.LinearAlgebraProvider.PointWiseDivideArrays(_values, denseOther._values, denseResult._values); } } @@ -547,8 +544,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single CommonParallel.For( 0, - _data.Length, - index => denseResult._data[index] %= divisor); + _values.Length, + index => denseResult._values[index] %= divisor); } } @@ -567,7 +564,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single var sum = 0.0f; for (var i = 0; i < _rowCount; i++) { - sum += _data[(i * _rowCount) + i]; + sum += _values[(i * _rowCount) + i]; } return sum; diff --git a/src/Numerics/LinearAlgebra/Single/DenseVector.cs b/src/Numerics/LinearAlgebra/Single/DenseVector.cs index 51b7a4cd..ac1bd409 100644 --- a/src/Numerics/LinearAlgebra/Single/DenseVector.cs +++ b/src/Numerics/LinearAlgebra/Single/DenseVector.cs @@ -33,6 +33,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single using Generic; using NumberTheory; using Properties; + using Storage; using Threading; /// @@ -41,6 +42,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single [Serializable] public class DenseVector : Vector { + /// + /// Number of elements + /// + readonly int _length; + + /// + /// Gets the vector's data. + /// + readonly float[] _values; + + internal DenseVector(DenseVectorStorage storage) + : base(storage) + { + _length = storage.Length; + _values = storage.Data; + } + /// /// Initializes a new instance of the class with a given size. /// @@ -51,9 +69,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// If is less than one. /// public DenseVector(int size) - : base(size) + : this(new DenseVectorStorage(size)) { - Data = new float[size]; } /// @@ -72,9 +89,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single public DenseVector(int size, float value) : this(size) { - for (var index = 0; index < Data.Length; index++) + for (var index = 0; index < _values.Length; index++) { - Data[index] = value; + _values[index] = value; } } @@ -88,31 +105,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single public DenseVector(Vector other) : this(other.Count) { - var vector = other as DenseVector; - if (vector == null) - { - CommonParallel.For( - 0, - Data.Length, - index => this[index] = other[index]); - } - else - { - Buffer.BlockCopy(vector.Data, 0, Data, 0, Data.Length * Constants.SizeOfFloat); - } - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public DenseVector(DenseVector other) - : this(other.Count) - { - Buffer.BlockCopy(other.Data, 0, Data, 0, Data.Length * Constants.SizeOfFloat); + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -121,20 +114,18 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The array to create this vector from. /// The vector does not copy the array, but keeps a reference to it. Any /// changes to the vector will also change the array. - public DenseVector(float[] array) : base(array.Length) + public DenseVector(float[] array) + : this(new DenseVectorStorage(array.Length, array)) { - Data = array; } /// - /// Gets the vector's internal data. + /// Gets the vector's data. /// - /// The vector's internal data. - /// Changing values in the array also changes the corresponding value in vector. Use with care. - internal float[] Data + /// The vector's data. + public float[] Values { - get; - private set; + get { return _values; } } /// @@ -145,14 +136,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// /// A reference to the internal date of the given vector. /// - public static implicit operator float[](DenseVector vector) + public static explicit operator float[](DenseVector vector) { if (vector == null) { throw new ArgumentNullException(); } - return vector.Data; + return vector.Values; } /// @@ -178,10 +169,10 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// This vector as a column matrix. public override Matrix ToColumnMatrix() { - var matrix = new DenseMatrix(Count, 1); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(_length, 1); + for (var i = 0; i < _values.Length; i++) { - matrix.At(i, 0, Data[i]); + matrix.At(i, 0, _values[i]); } return matrix; @@ -193,33 +184,15 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// This vector as a row matrix. public override Matrix ToRowMatrix() { - var matrix = new DenseMatrix(1, Count); - for (var i = 0; i < Data.Length; i++) + var matrix = new DenseMatrix(1, _length); + for (var i = 0; i < _values.Length; i++) { - matrix.At(0, i, Data[i]); + matrix.At(0, i, _values[i]); } return matrix; } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - /// If is negative or - /// greater than the size of the vector. - public override float this[int index] - { - get - { - return Data[index]; - } - - set - { - Data[index] = value; - } - } - /// /// Creates a matrix with the given dimensions using the same storage type /// as this vector. @@ -253,49 +226,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single return new DenseVector(size); } - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as DenseVector; - if (otherVector == null) - { - CommonParallel.For( - 0, - Data.Length, - index => target[index] = Data[index]); - } - else - { - Buffer.BlockCopy(Data, 0, otherVector.Data, 0, Data.Length * Constants.SizeOfFloat); - } - } - /// /// Adds a scalar to each element of the vector and stores the result in the result vector. /// @@ -312,8 +242,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] + scalar); + _values.Length, + index => dense._values[index] = _values[index] + scalar); } } @@ -328,7 +258,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, 1.0f, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, 1.0f, odense.Values, rdense.Values); } else { @@ -397,8 +327,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] - scalar); + _values.Length, + index => dense._values[index] = _values[index] - scalar); } } @@ -413,7 +343,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single var odense = other as DenseVector; if (rdense != null && odense != null) { - Control.LinearAlgebraProvider.AddVectorToScaledVector(Data, -1.0f, odense.Data, rdense.Data); + Control.LinearAlgebraProvider.AddVectorToScaledVector(_values, -1.0f, odense.Values, rdense.Values); } else { @@ -472,11 +402,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new DenseVector(Count); + var result = new DenseVector(_length); CommonParallel.For( 0, - Data.Length, - index => result[index] = -Data[index]); + _values.Length, + index => result[index] = -_values[index]); return result; } @@ -496,7 +426,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - Control.LinearAlgebraProvider.ScaleArray(scalar, Data, denseResult.Data); + Control.LinearAlgebraProvider.ScaleArray(scalar, _values, denseResult.Values); } } @@ -510,7 +440,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single { var denseVector = other as DenseVector; - return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(Data, denseVector.Data); + return denseVector == null ? base.DoDotProduct(other) : Control.LinearAlgebraProvider.DotProduct(_values, denseVector.Values); } /// @@ -572,7 +502,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide"); } - return Control.LinearAlgebraProvider.DotProduct(leftSide.Data, rightSide.Data); + return Control.LinearAlgebraProvider.DotProduct(leftSide.Values, rightSide.Values); } /// @@ -602,16 +532,16 @@ namespace MathNet.Numerics.LinearAlgebra.Single var denseResult = result as DenseVector; if (denseResult == null) { - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - result.At(index, Data[index] % divisor); + result.At(index, _values[index] % divisor); } } else { - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - denseResult.Data[index] = Data[index] % divisor; + denseResult._values[index] = _values[index] % divisor; } } } @@ -640,10 +570,10 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override int AbsoluteMinimumIndex() { var index = 0; - var min = Math.Abs(Data[index]); - for (var i = 1; i < Count; i++) + var min = Math.Abs(_values[index]); + for (var i = 1; i < _length; i++) { - var test = Math.Abs(Data[i]); + var test = Math.Abs(_values[i]); if (test < min) { index = i; @@ -660,7 +590,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The value of the absolute minimum element. public override float AbsoluteMinimum() { - return Math.Abs(Data[AbsoluteMinimumIndex()]); + return Math.Abs(_values[AbsoluteMinimumIndex()]); } /// @@ -669,7 +599,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The value of the absolute maximum element. public override float AbsoluteMaximum() { - return Math.Abs(Data[AbsoluteMaximumIndex()]); + return Math.Abs(_values[AbsoluteMaximumIndex()]); } /// @@ -679,10 +609,10 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override int AbsoluteMaximumIndex() { var index = 0; - var max = Math.Abs(Data[index]); - for (var i = 1; i < Count; i++) + var max = Math.Abs(_values[index]); + for (var i = 1; i < _length; i++) { - var test = Math.Abs(Data[i]); + var test = Math.Abs(_values[i]); if (test > max) { index = i; @@ -706,7 +636,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// If is not positive. public override Vector SubVector(int index, int length) { - if (index < 0 || index >= Count) + if (index < 0 || index >= _length) { throw new ArgumentOutOfRangeException("index"); } @@ -716,7 +646,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single throw new ArgumentOutOfRangeException("length"); } - if (index + length > Count) + if (index + length > _length) { throw new ArgumentOutOfRangeException("length"); } @@ -725,8 +655,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single CommonParallel.For( index, - index + length, - i => result.Data[i - index] = Data[i]); + index + length, + i => result._values[i - index] = _values[i]); return result; } @@ -743,7 +673,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single throw new ArgumentNullException("values"); } - if (values.Length != Count) + if (values.Length != _length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "values"); } @@ -751,7 +681,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single CommonParallel.For( 0, values.Length, - i => Data[i] = values[i]); + i => _values[i] = values[i]); } /// @@ -761,13 +691,13 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override int MaximumIndex() { var index = 0; - var max = Data[0]; - for (var i = 1; i < Count; i++) + var max = _values[0]; + for (var i = 1; i < _length; i++) { - if (max < Data[i]) + if (max < _values[i]) { index = i; - max = Data[i]; + max = _values[i]; } } @@ -781,13 +711,13 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override int MinimumIndex() { var index = 0; - var min = Data[0]; - for (var i = 1; i < Count; i++) + var min = _values[0]; + for (var i = 1; i < _length; i++) { - if (min > Data[i]) + if (min > _values[i]) { index = i; - min = Data[i]; + min = _values[i]; } } @@ -802,9 +732,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single { var sum = 0.0f; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Data[i]; + sum += _values[i]; } return sum; @@ -818,9 +748,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single { var sum = 0.0f; - for (var i = 0; i < Count; i++) + for (var i = 0; i < _length; i++) { - sum += Math.Abs(Data[i]); + sum += Math.Abs(_values[i]); } return sum; @@ -842,8 +772,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] * other[index]); + _values.Length, + index => dense._values[index] = _values[index] * other[index]); } } @@ -864,8 +794,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single { CommonParallel.For( 0, - Data.Length, - index => dense.Data[index] = Data[index] / other[index]); + _values.Length, + index => dense._values[index] = _values[index] / other[index]); } } @@ -897,7 +827,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single { for (var j = 0; j < v.Count; j++) { - matrix.At(i, j, u.Data[i] * v.Data[j]); + matrix.At(i, j, u._values[i] * v._values[j]); } }); return matrix; @@ -937,19 +867,19 @@ namespace MathNet.Numerics.LinearAlgebra.Single if (2.0 == p) { - return Data.Aggregate(0f, SpecialFunctions.Hypotenuse); + return _values.Aggregate(0f, SpecialFunctions.Hypotenuse); } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(Data, (i, v) => Math.Abs(v), Math.Max, 0f); + return CommonParallel.Aggregate(_values, (i, v) => Math.Abs(v), Math.Max, 0f); } var sum = 0.0; - for (var index = 0; index < Count; index++) + for (var index = 0; index < _length; index++) { - sum += Math.Pow(Math.Abs(Data[index]), p); + sum += Math.Pow(Math.Abs(_values[index]), p); } return (float)Math.Pow(sum, 1.0 / p); @@ -1115,29 +1045,5 @@ namespace MathNet.Numerics.LinearAlgebra.Single } #endregion - - /// - /// Resets all values to zero. - /// - public override void Clear() - { - Array.Clear(Data, 0, Data.Length); - } - - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - internal protected override float At(int index) - { - return Data[index]; - } - - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - internal protected override void At(int index, float value) - { - Data[index] = value; - } } } diff --git a/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs b/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs index 376df7d6..41fdc738 100644 --- a/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/DiagonalMatrix.cs @@ -377,12 +377,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single return; } - if (_data.Length != denseSource.Data.Length) + if (_data.Length != denseSource.Values.Length) { throw new ArgumentException(Resources.ArgumentVectorsSameLength, "source"); } - Buffer.BlockCopy(denseSource.Data, 0, _data, 0, denseSource.Data.Length * Constants.SizeOfFloat); + Buffer.BlockCopy(denseSource.Values, 0, _data, 0, denseSource.Values.Length * Constants.SizeOfFloat); } /// diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseCholesky.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseCholesky.cs index e5c29589..bbef476f 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseCholesky.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseCholesky.cs @@ -167,11 +167,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } // Copy the contents of input to result. - Buffer.BlockCopy(dinput.Data, 0, dresult.Data, 0, dinput.Data.Length * Constants.SizeOfFloat); + Buffer.BlockCopy(dinput.Values, 0, dresult.Values, 0, dinput.Values.Length * Constants.SizeOfFloat); // Cholesky solve by overwriting result. var dfactor = (DenseMatrix)CholeskyFactor; - Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Data, 1); + Control.LinearAlgebraProvider.CholeskySolveFactored(dfactor.Data, dfactor.RowCount, dresult.Values, 1); } } } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs index 9460a9d7..321574fa 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseGramSchmidt.cs @@ -216,7 +216,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do GramSchmidt factorization for dense vectors at the moment."); } - _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Data, 1, dresult.Data); + _provider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixQ.RowCount, MatrixQ.ColumnCount, null, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseLU.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseLU.cs index 03eb9371..e307bb38 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseLU.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseLU.cs @@ -168,11 +168,11 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization } // Copy the contents of input to result. - Buffer.BlockCopy(dinput.Data, 0, dresult.Data, 0, dinput.Data.Length * Constants.SizeOfFloat); + Buffer.BlockCopy(dinput.Values, 0, dresult.Values, 0, dinput.Values.Length * Constants.SizeOfFloat); // LU solve by overwriting result. var dfactors = (DenseMatrix)Factors; - Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Data); + Control.LinearAlgebraProvider.LUSolveFactored(1, dfactors.Data, dfactors.RowCount, Pivots, dresult.Values); } /// diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs index 9157b412..46ae3c1b 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseQR.cs @@ -187,7 +187,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do QR factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.QRSolveFactored(((DenseMatrix)MatrixQ).Data, ((DenseMatrix)MatrixR).Data, MatrixR.RowCount, MatrixR.ColumnCount, Tau, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Single/Factorization/DenseSvd.cs b/src/Numerics/LinearAlgebra/Single/Factorization/DenseSvd.cs index c210034b..19eed0fd 100644 --- a/src/Numerics/LinearAlgebra/Single/Factorization/DenseSvd.cs +++ b/src/Numerics/LinearAlgebra/Single/Factorization/DenseSvd.cs @@ -69,7 +69,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization VectorS = new DenseVector(nm); MatrixU = new DenseMatrix(matrix.RowCount); MatrixVT = new DenseMatrix(matrix.ColumnCount); - Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); + Control.LinearAlgebraProvider.SingularValueDecomposition(computeVectors, ((DenseMatrix)matrix.Clone()).Data, matrix.RowCount, matrix.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data); } /// @@ -125,7 +125,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do SVD factorization for dense matrices at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, input.ColumnCount, dresult.Data); } /// @@ -175,7 +175,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single.Factorization throw new NotSupportedException("Can only do SVD factorization for dense vectors at the moment."); } - Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Data, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Data, 1, dresult.Data); + Control.LinearAlgebraProvider.SvdSolveFactored(MatrixU.RowCount, MatrixVT.ColumnCount, ((DenseVector)VectorS).Values, ((DenseMatrix)MatrixU).Data, ((DenseMatrix)MatrixVT).Data, dinput.Values, 1, dresult.Values); } } } diff --git a/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs b/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs index 271ec7a6..52b3596d 100644 --- a/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs +++ b/src/Numerics/LinearAlgebra/Single/SparseMatrix.cs @@ -97,6 +97,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The number of columns. /// /// The value which we assign to each element of the matrix. + [Obsolete("Use a dense matrix instead.")] public SparseMatrix(int rows, int columns, float value) : this(rows, columns) { diff --git a/src/Numerics/LinearAlgebra/Single/SparseVector.cs b/src/Numerics/LinearAlgebra/Single/SparseVector.cs index 8fa7cd46..f0c1adc9 100644 --- a/src/Numerics/LinearAlgebra/Single/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Single/SparseVector.cs @@ -37,6 +37,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single using Generic; using NumberTheory; using Properties; + using Storage; using Threading; /// @@ -46,24 +47,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single [Serializable] public class SparseVector : Vector { - /// - /// Gets the vector's internal data. The array containing the actual values; only the non-zero values are stored. - /// - private float[] _nonZeroValues = new float[0]; + readonly SparseVectorStorage _storage; - /// - /// The indices of the non-zero entries. - /// - private int[] _nonZeroIndices = new int[0]; - - /// - /// Gets the number of non zero elements in the vector. - /// - /// The number of non zero elements. - public int NonZerosCount + internal SparseVector(SparseVectorStorage storage) + : base(storage) { - get; - private set; + _storage = storage; } #region Constructors @@ -77,7 +66,8 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// /// If is less than one. /// - public SparseVector(int size) : base(size) + public SparseVector(int size) + : this(new SparseVectorStorage(size)) { } @@ -94,55 +84,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// /// If is less than one. /// - public SparseVector(int size, float value) : this(size) + [Obsolete("Use a dense vector instead.")] + public SparseVector(int size, float value) + : this(new SparseVectorStorage(size)) { if (value == 0.0) { - // Skip adding values return; } - // We already know that this vector is "full", let's allocate all needed memory - _nonZeroValues = new float[size]; - _nonZeroIndices = new int[size]; - NonZerosCount = size; + var valueCount = _storage.ValueCount = size; + var indices = _storage.Indices = new int[valueCount]; + var values = _storage.Values = new float[valueCount]; - CommonParallel.For( - 0, - Count, - index => - { - _nonZeroValues[index] = value; - _nonZeroIndices[index] = index; - }); - } - - /// - /// Initializes a new instance of the class by - /// copying the values from another. - /// - /// - /// The vector to create the new vector from. - /// - public SparseVector(Vector other) : this(other.Count) - { - var vector = other as SparseVector; - if (vector == null) + for (int i = 0; i < values.Length; i++) { - for (var i = 0; i < other.Count; i++) - { - this[i] = other.At(i); - } - } - else - { - _nonZeroValues = new float[vector.NonZerosCount]; - _nonZeroIndices = new int[vector.NonZerosCount]; - NonZerosCount = vector.NonZerosCount; - - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - Buffer.BlockCopy(vector._nonZeroValues, 0, _nonZeroValues, 0, vector.NonZerosCount * Constants.SizeOfFloat); - Buffer.BlockCopy(vector._nonZeroIndices, 0, _nonZeroIndices, 0, vector.NonZerosCount * Constants.SizeOfInt); + values[i] = value; + indices[i] = i; } } @@ -153,15 +111,10 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// /// The vector to create the new vector from. /// - public SparseVector(SparseVector other) : this(other.Count) + public SparseVector(Vector other) + : this(new SparseVectorStorage(other.Count)) { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - _nonZeroValues = new float[other.NonZerosCount]; - _nonZeroIndices = new int[other.NonZerosCount]; - NonZerosCount = other.NonZerosCount; - - Buffer.BlockCopy(other._nonZeroValues, 0, _nonZeroValues, 0, other.NonZerosCount * Constants.SizeOfFloat); - Buffer.BlockCopy(other._nonZeroIndices, 0, _nonZeroIndices, 0, other.NonZerosCount * Constants.SizeOfInt); + other.Storage.CopyTo(Storage, skipClearing: true); } /// @@ -169,11 +122,12 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// /// The array to create this vector from. /// The vector copy the array. Any changes to the vector will NOT change the array. - public SparseVector(IList array) : this(array.Count) + public SparseVector(IList array) + : this(new SparseVectorStorage(array.Count)) { for (var i = 0; i < array.Count; i++) { - this[i] = array[i]; + Storage.At(i, array[i]); } } @@ -186,9 +140,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override Matrix ToColumnMatrix() { var matrix = new SparseMatrix(Count, 1); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(_nonZeroIndices[i], 0, _nonZeroValues[i]); + matrix.At(_storage.Indices[i], 0, _storage.Values[i]); } return matrix; @@ -201,44 +155,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override Matrix ToRowMatrix() { var matrix = new SparseMatrix(1, Count); - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - matrix.At(0, _nonZeroIndices[i], _nonZeroValues[i]); + matrix.At(0, _storage.Indices[i], _storage.Values[i]); } return matrix; } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - /// If is negative or - /// greater than the size of the vector. - public override float this[int index] - { - get - { - // If index is out of bounds - if ((index < 0) || (index >= Count)) - { - throw new IndexOutOfRangeException(); - } - - return At(index); - } - - set - { - // If index is out of bounds - if ((index < 0) || (index >= Count)) - { - throw new IndexOutOfRangeException(); - } - - At(index, value); - } - } - /// /// Creates a matrix with the given dimensions using the same storage type /// as this vector. @@ -272,68 +196,6 @@ namespace MathNet.Numerics.LinearAlgebra.Single return new SparseVector(size); } - /// - /// Resets all values to zero. - /// - public override void Clear() - { - NonZerosCount = 0; - } - - /// - /// Copies the values of this vector into the target vector. - /// - /// - /// The vector to copy elements into. - /// - /// - /// If is . - /// - /// - /// If is not the same size as this vector. - /// - public override void CopyTo(Vector target) - { - if (target == null) - { - throw new ArgumentNullException("target"); - } - - if (Count != target.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); - } - - if (ReferenceEquals(this, target)) - { - return; - } - - var otherVector = target as SparseVector; - if (otherVector == null) - { - target.Clear(); - - for (var index = 0; index < NonZerosCount; index++) - { - target.At(_nonZeroIndices[index], _nonZeroValues[index]); - } - } - else - { - // Lets copy only needed data. Portion of needed data is determined by NonZerosCount value - otherVector._nonZeroValues = new float[NonZerosCount]; - otherVector._nonZeroIndices = new int[NonZerosCount]; - otherVector.NonZerosCount = NonZerosCount; - - if (NonZerosCount != 0) - { - Buffer.BlockCopy(_nonZeroValues, 0, otherVector._nonZeroValues, 0, NonZerosCount*Constants.SizeOfFloat); - Buffer.BlockCopy(_nonZeroIndices, 0, otherVector._nonZeroIndices, 0, NonZerosCount*Constants.SizeOfInt); - } - } - } - #region Operators and supplementary functions /// @@ -371,15 +233,15 @@ namespace MathNet.Numerics.LinearAlgebra.Single } //populate the non zero values from this - for (int j = 0; j < NonZerosCount; j++) + for (int j = 0; j < _storage.ValueCount; j++) { - vnonZeroValues[_nonZeroIndices[j]] = _nonZeroValues[j] + scalar; + vnonZeroValues[_storage.Indices[j]] = _storage.Values[j] + scalar; } //assign this vectors arrary to the new arrays. - _nonZeroValues = vnonZeroValues; - _nonZeroIndices = vnonZeroIndices; - NonZerosCount = Count; + _storage.Values = vnonZeroValues; + _storage.Indices = vnonZeroIndices; + _storage.ValueCount = Count; } else @@ -421,22 +283,22 @@ namespace MathNet.Numerics.LinearAlgebra.Single if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] += otherSparse._nonZeroValues[j++]; + _storage.Values[i++] += otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != 0.0) { - InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], otherValue); + InsertAtUnchecked(i++, otherSparse._storage.Indices[j], otherValue); } j++; } @@ -446,25 +308,25 @@ namespace MathNet.Numerics.LinearAlgebra.Single { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] + otherSparse.At(next)); + result.At(next, _storage.Values[i] + otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) + otherSparse._nonZeroValues[j]); + result.At(next, At(next) + otherSparse._storage.Values[j]); } j++; } @@ -567,22 +429,22 @@ namespace MathNet.Numerics.LinearAlgebra.Single if (ReferenceEquals(this, resultSparse)) { int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (i < NonZerosCount && j < otherSparse.NonZerosCount && _nonZeroIndices[i] == otherSparse._nonZeroIndices[j]) + if (i < _storage.ValueCount && j < otherSparse._storage.ValueCount && _storage.Indices[i] == otherSparse._storage.Indices[j]) { - _nonZeroValues[i++] -= otherSparse._nonZeroValues[j++]; + _storage.Values[i++] -= otherSparse._storage.Values[j++]; } - else if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + else if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { i++; } else { - var otherValue = otherSparse._nonZeroValues[j]; + var otherValue = otherSparse._storage.Values[j]; if (otherValue != 0.0) { - InsertAtUnchecked(i++, otherSparse._nonZeroIndices[j], -otherValue); + InsertAtUnchecked(i++, otherSparse._storage.Indices[j], -otherValue); } j++; } @@ -592,25 +454,25 @@ namespace MathNet.Numerics.LinearAlgebra.Single { result.Clear(); int i = 0, j = 0, last = -1; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] <= otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] <= otherSparse._storage.Indices[j]) { - var next = _nonZeroIndices[i]; + var next = _storage.Indices[i]; if (next != last) { last = next; - result.At(next, _nonZeroValues[i] - otherSparse.At(next)); + result.At(next, _storage.Values[i] - otherSparse.At(next)); } i++; } else { - var next = otherSparse._nonZeroIndices[j]; + var next = otherSparse._storage.Indices[j]; if (next != last) { last = next; - result.At(next, At(next) - otherSparse._nonZeroValues[j]); + result.At(next, At(next) - otherSparse._storage.Values[j]); } j++; } @@ -669,23 +531,18 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// Added as an alternative to the unary negation operator. public override Vector Negate() { - var result = new SparseVector(Count) - { - _nonZeroValues = new float[NonZerosCount], - _nonZeroIndices = new int[NonZerosCount], - NonZerosCount = NonZerosCount - }; + var result = new SparseVectorStorage(Count); + var valueCount = result.ValueCount = _storage.ValueCount; + var indices = result.Indices = new int[valueCount]; + var values = result.Values = new float[valueCount]; - if (NonZerosCount != 0) + if (valueCount != 0) { - CommonParallel.For( - 0, - NonZerosCount, - index => result._nonZeroValues[index] = -_nonZeroValues[index]); - Buffer.BlockCopy(_nonZeroIndices, 0, result._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); + CommonParallel.For(0, valueCount, index => values[index] = -_storage.Values[index]); + Buffer.BlockCopy(_storage.Indices, 0, indices, 0, valueCount * Constants.SizeOfInt); } - return result; + return new SparseVector(result); } /// @@ -706,7 +563,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } var copy = new SparseVector(this); - Control.LinearAlgebraProvider.ScaleArray(scalar, copy._nonZeroValues, copy._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, copy._storage.Values, copy._storage.Values); return copy; } @@ -741,23 +598,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single if (sparseResult == null) { result.Clear(); - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - result.At(_nonZeroIndices[index], scalar * _nonZeroValues[index]); + result.At(_storage.Indices[index], scalar * _storage.Values[index]); } } else { if (!ReferenceEquals(this, result)) { - sparseResult.NonZerosCount = NonZerosCount; - sparseResult._nonZeroIndices = new int[NonZerosCount]; - Buffer.BlockCopy(_nonZeroIndices, 0, sparseResult._nonZeroIndices, 0, NonZerosCount * Constants.SizeOfInt); - sparseResult._nonZeroValues = new float[NonZerosCount]; - Buffer.BlockCopy(_nonZeroValues, 0, sparseResult._nonZeroValues, 0, NonZerosCount * Constants.SizeOfFloat); + sparseResult._storage.ValueCount = _storage.ValueCount; + sparseResult._storage.Indices = new int[_storage.ValueCount]; + Buffer.BlockCopy(_storage.Indices, 0, sparseResult._storage.Indices, 0, _storage.ValueCount * Constants.SizeOfInt); + sparseResult._storage.Values = new float[_storage.ValueCount]; + Buffer.BlockCopy(_storage.Values, 0, sparseResult._storage.Values, 0, _storage.ValueCount * Constants.SizeOfFloat); } - Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._nonZeroValues, sparseResult._nonZeroValues); + Control.LinearAlgebraProvider.ScaleArray(scalar, sparseResult._storage.Values, sparseResult._storage.Values); } } @@ -776,16 +633,16 @@ namespace MathNet.Numerics.LinearAlgebra.Single if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * _nonZeroValues[i]; + result += _storage.Values[i] * _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i] * other.At(_nonZeroIndices[i]); + result += _storage.Values[i] * other.At(_storage.Indices[i]); } } @@ -880,17 +737,17 @@ namespace MathNet.Numerics.LinearAlgebra.Single { if (ReferenceEquals(this, result)) { - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - _nonZeroValues[index] %= divisor; + _storage.Values[index] %= divisor; } } else { result.Clear(); - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - result.At(_nonZeroIndices[index], _nonZeroValues[index] % divisor); + result.At(_storage.Indices[index], _storage.Values[index] % divisor); } } } @@ -918,17 +775,17 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The index of absolute minimum element. public override int AbsoluteMinimumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { // No non-zero elements. Return 0 return 0; } var index = 0; - var min = Math.Abs(_nonZeroValues[index]); - for (var i = 1; i < NonZerosCount; i++) + var min = Math.Abs(_storage.Values[index]); + for (var i = 1; i < _storage.ValueCount; i++) { - var test = Math.Abs(_nonZeroValues[i]); + var test = Math.Abs(_storage.Values[i]); if (test < min) { index = i; @@ -936,7 +793,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -1006,23 +863,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The index of absolute maximum element. public override int MaximumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0; } var index = 0; - var max = _nonZeroValues[0]; - for (var i = 1; i < NonZerosCount; i++) + var max = _storage.Values[0]; + for (var i = 1; i < _storage.ValueCount; i++) { - if (max < _nonZeroValues[i]) + if (max < _storage.Values[i]) { index = i; - max = _nonZeroValues[i]; + max = _storage.Values[i]; } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -1031,23 +888,23 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The index of minimum element. public override int MinimumIndex() { - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0; } var index = 0; - var min = _nonZeroValues[0]; - for (var i = 1; i < NonZerosCount; i++) + var min = _storage.Values[0]; + for (var i = 1; i < _storage.ValueCount; i++) { - if (min > _nonZeroValues[i]) + if (min > _storage.Values[i]) { index = i; - min = _nonZeroValues[i]; + min = _storage.Values[i]; } } - return _nonZeroIndices[index]; + return _storage.Indices[index]; } /// @@ -1057,9 +914,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override float Sum() { var result = 0.0f; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += _nonZeroValues[i]; + result += _storage.Values[i]; } return result; @@ -1072,9 +929,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override float SumMagnitudes() { var result = 0.0f; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - result += Math.Abs(_nonZeroValues[i]); + result += Math.Abs(_storage.Values[i]); } return result; @@ -1089,17 +946,17 @@ namespace MathNet.Numerics.LinearAlgebra.Single { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] *= _nonZeroValues[i]; + _storage.Values[i] *= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, other.At(index) * _nonZeroValues[i]); + var index = _storage.Indices[i]; + result.At(index, other.At(index) * _storage.Values[i]); } } } @@ -1113,17 +970,17 @@ namespace MathNet.Numerics.LinearAlgebra.Single { if (ReferenceEquals(this, other)) { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - _nonZeroValues[i] /= _nonZeroValues[i]; + _storage.Values[i] /= _storage.Values[i]; } } else { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - var index = _nonZeroIndices[i]; - result.At(index, _nonZeroValues[i] / other.At(index)); + var index = _storage.Indices[i]; + result.At(index, _storage.Values[i] / other.At(index)); } } } @@ -1149,13 +1006,13 @@ namespace MathNet.Numerics.LinearAlgebra.Single } var matrix = new SparseMatrix(u.Count, v.Count); - for (var i = 0; i < u.NonZerosCount; i++) + for (var i = 0; i < u._storage.ValueCount; i++) { - for (var j = 0; j < v.NonZerosCount; j++) + for (var j = 0; j < v._storage.ValueCount; j++) { - if (u._nonZeroIndices[i] == v._nonZeroIndices[j]) + if (u._storage.Indices[i] == v._storage.Indices[j]) { - matrix.At(i, j, u._nonZeroValues[i] * v._nonZeroValues[j]); + matrix.At(i, j, u._storage.Values[i] * v._storage.Values[j]); } } } @@ -1191,25 +1048,25 @@ namespace MathNet.Numerics.LinearAlgebra.Single throw new ArgumentOutOfRangeException("p"); } - if (NonZerosCount == 0) + if (_storage.ValueCount == 0) { return 0.0f; } if (2.0 == p) { - return _nonZeroValues.Aggregate(0f, SpecialFunctions.Hypotenuse); + return _storage.Values.Aggregate(0f, SpecialFunctions.Hypotenuse); } if (Double.IsPositiveInfinity(p)) { - return CommonParallel.Aggregate(0, NonZerosCount, i => Math.Abs(_nonZeroValues[i]), Math.Max, 0f); + return CommonParallel.Aggregate(0, _storage.ValueCount, i => Math.Abs(_storage.Values[i]), Math.Max, 0f); } var sum = 0.0; - for (var index = 0; index < NonZerosCount; index++) + for (var index = 0; index < _storage.ValueCount; index++) { - sum += Math.Pow(Math.Abs(_nonZeroValues[index]), p); + sum += Math.Pow(Math.Abs(_storage.Values[index]), p); } return (float)Math.Pow(sum, 1.0 / p); @@ -1376,97 +1233,33 @@ namespace MathNet.Numerics.LinearAlgebra.Single #endregion - /// - /// Gets the value at the given index. - /// - /// Value real index in array - /// The value at the given index. - internal protected override float At(int index) - { - // Search if item idex exists in NonZeroIndices array in range "0 - real nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - return itemIndex >= 0 ? _nonZeroValues[itemIndex] : 0.0f; - } - - /// - /// Delete, Add or Update the value in NonZeroValues and NonZeroIndices - /// - /// Value real index in array - /// The value to set. - /// This method assume that index is between 0 and Array Size - internal protected override void At(int index, float value) - { - // Search if "index" already exists in range "0 - real nonzero values count" - var itemIndex = Array.BinarySearch(_nonZeroIndices, 0, NonZerosCount, index); - - if (itemIndex >= 0) - { - // Item already exist at itemIndex - if (value == 0.0) - { - RemoveAtUnchecked(itemIndex); - } - else - { - _nonZeroValues[itemIndex] = value; - } - } - else - { - if (value != 0.0) - { - InsertAtUnchecked(~itemIndex, index, value); - } - } - } - private void InsertAtUnchecked(int itemIndex, int index, float value) { // Check if the storage needs to be increased - if ((NonZerosCount == _nonZeroValues.Length) && (NonZerosCount < Count)) + if ((_storage.ValueCount == _storage.Values.Length) && (_storage.ValueCount < Count)) { // Value and Indices arrays are completely full so we increase the size - var size = Math.Min(_nonZeroValues.Length + GrowthSize(), Count); - Array.Resize(ref _nonZeroValues, size); - Array.Resize(ref _nonZeroIndices, size); + var size = Math.Min(_storage.Values.Length + GrowthSize(), Count); + Array.Resize(ref _storage.Values, size); + Array.Resize(ref _storage.Indices, size); } // Move all values (with a position larger than index) in the value array // to the next position // Move all values (with a position larger than index) in the columIndices // array to the next position - for (var i = NonZerosCount - 1; i > itemIndex - 1; i--) + for (var i = _storage.ValueCount - 1; i > itemIndex - 1; i--) { - _nonZeroValues[i + 1] = _nonZeroValues[i]; - _nonZeroIndices[i + 1] = _nonZeroIndices[i]; + _storage.Values[i + 1] = _storage.Values[i]; + _storage.Indices[i + 1] = _storage.Indices[i]; } // Add the value and the column index - _nonZeroValues[itemIndex] = value; - _nonZeroIndices[itemIndex] = index; + _storage.Values[itemIndex] = value; + _storage.Indices[itemIndex] = index; // increase the number of non-zero numbers by one - NonZerosCount += 1; - } - - private void RemoveAtUnchecked(int itemIndex) - { - // Value is zero. Let's delete it from Values and Indices array - for (var i = itemIndex + 1; i < NonZerosCount; i++) - { - _nonZeroValues[i - 1] = _nonZeroValues[i]; - _nonZeroIndices[i - 1] = _nonZeroIndices[i]; - } - - NonZerosCount -= 1; - - // Check whether we need to shrink the arrays. This is reasonable to do if - // there are a lot of non-zero elements and storage is two times bigger - if ((NonZerosCount > 1024) && (NonZerosCount < _nonZeroIndices.Length / 2)) - { - Array.Resize(ref _nonZeroValues, NonZerosCount); - Array.Resize(ref _nonZeroIndices, NonZerosCount); - } + _storage.ValueCount += 1; } /// @@ -1477,19 +1270,19 @@ namespace MathNet.Numerics.LinearAlgebra.Single private int GrowthSize() { int delta; - if (_nonZeroValues.Length > 1024) + if (_storage.Values.Length > 1024) { - delta = _nonZeroValues.Length / 4; + delta = _storage.Values.Length / 4; } else { - if (_nonZeroValues.Length > 256) + if (_storage.Values.Length > 256) { delta = 512; } else { - delta = _nonZeroValues.Length > 64 ? 128 : 32; + delta = _storage.Values.Length > 64 ? 128 : 32; } } @@ -1502,7 +1295,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single { if (Count > 20) { - return String.Format("SparseVectorOfSingle({0},{1},{2})", Count, NonZerosCount, GetHashCode()); + return String.Format("SparseVectorOfSingle({0},{1},{2})", Count, _storage.ValueCount, GetHashCode()); } return base.ToString(format, formatProvider); @@ -1516,14 +1309,14 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// public override int GetHashCode() { - var hashNum = Math.Min(NonZerosCount, 20); + var hashNum = Math.Min(_storage.ValueCount, 20); long hash = 0; for (var i = 0; i < hashNum; i++) { #if PORTABLE - hash ^= Precision.DoubleToInt64Bits(this._nonZeroValues[i].GetHashCode()); + hash ^= Precision.DoubleToInt64Bits(this._storage.Values[i].GetHashCode()); #else - hash ^= BitConverter.DoubleToInt64Bits(_nonZeroValues[i].GetHashCode()); + hash ^= BitConverter.DoubleToInt64Bits(_storage.Values[i].GetHashCode()); #endif } @@ -1567,27 +1360,27 @@ namespace MathNet.Numerics.LinearAlgebra.Single } int i = 0, j = 0; - while (i < NonZerosCount || j < otherSparse.NonZerosCount) + while (i < _storage.ValueCount || j < otherSparse._storage.ValueCount) { - if (j >= otherSparse.NonZerosCount || i < NonZerosCount && _nonZeroIndices[i] < otherSparse._nonZeroIndices[j]) + if (j >= otherSparse._storage.ValueCount || i < _storage.ValueCount && _storage.Indices[i] < otherSparse._storage.Indices[j]) { - if (_nonZeroValues[i++] != 0f) + if (_storage.Values[i++] != 0f) { return false; } continue; } - if (i >= NonZerosCount || j < otherSparse.NonZerosCount && otherSparse._nonZeroIndices[j] < _nonZeroIndices[i]) + if (i >= _storage.ValueCount || j < otherSparse._storage.ValueCount && otherSparse._storage.Indices[j] < _storage.Indices[i]) { - if (otherSparse._nonZeroValues[j++] != 0f) + if (otherSparse._storage.Values[j++] != 0f) { return false; } continue; } - if (!_nonZeroValues[i].AlmostEqual(otherSparse._nonZeroValues[j])) + if (!_storage.Values[i].AlmostEqual(otherSparse._storage.Values[j])) { return false; } @@ -1615,9 +1408,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// public override IEnumerable> GetIndexedEnumerator() { - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - yield return new Tuple(_nonZeroIndices[i], _nonZeroValues[i]); + yield return new Tuple(_storage.Indices[i], _storage.Values[i]); } } @@ -1630,9 +1423,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single public override float[] ToArray() { var ret = new float[Count]; - for (var i = 0; i < NonZerosCount; i++) + for (var i = 0; i < _storage.ValueCount; i++) { - ret[_nonZeroIndices[i]] = _nonZeroValues[i]; + ret[_storage.Indices[i]] = _storage.Values[i]; } return ret; diff --git a/src/Numerics/LinearAlgebra/Single/Vector.cs b/src/Numerics/LinearAlgebra/Single/Vector.cs index f5655156..8cc9e8e7 100644 --- a/src/Numerics/LinearAlgebra/Single/Vector.cs +++ b/src/Numerics/LinearAlgebra/Single/Vector.cs @@ -30,6 +30,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single using Distributions; using Generic; using Properties; + using Storage; using Threading; /// @@ -40,16 +41,9 @@ namespace MathNet.Numerics.LinearAlgebra.Single { /// /// Initializes a new instance of the Vector class. - /// Constructs a Vector with the given size. /// - /// - /// The size of the Vector to construct. - /// - /// - /// If is less than one. - /// - protected Vector(int size) - : base(size) + protected Vector(VectorStorage storage) + : base(storage) { } diff --git a/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs index 92b08994..b6c2afb5 100644 --- a/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs @@ -90,7 +90,6 @@ namespace MathNet.Numerics.LinearAlgebra.Storage throw new ArgumentException(message, "target"); } - //Buffer.BlockCopy(Data, 0, target.Data, 0, Data.Length * System.Runtime.InteropServices.Marshal.SizeOf(typeof(T))); Array.Copy(Data, 0, target.Data, 0, Data.Length); } @@ -113,14 +112,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage void CopySubVectorTo(DenseVectorStorage target, int sourceIndex, int targetIndex, int count) { - if (ReferenceEquals(this, target)) - { - throw new NotSupportedException(); - } - ValidateSubVectorRange(target, sourceIndex, targetIndex, count); - - //Buffer.BlockCopy(Data, sourceIndex, target.Data, targetIndex, count * System.Runtime.InteropServices.Marshal.SizeOf(typeof(T))); Array.Copy(Data, sourceIndex, target.Data, targetIndex, count); } } diff --git a/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs new file mode 100644 index 00000000..a4743aae --- /dev/null +++ b/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs @@ -0,0 +1,333 @@ +using System; +using MathNet.Numerics.Properties; + +namespace MathNet.Numerics.LinearAlgebra.Storage +{ + [Serializable] + public class SparseVectorStorage : VectorStorage + where T : struct, IEquatable, IFormattable + { + // [ruegg] public fields are OK here + + readonly T _zero; + + /// + /// Array that contains the indices of the non-zero values. + /// + public int[] Indices; + + /// + /// Array that contains the non-zero elements of the vector. + /// + public T[] Values; + + /// + /// Gets the number of non-zero elements in the vector. + /// + public int ValueCount; + + internal SparseVectorStorage(int length, T zero = default(T)) + : base(length) + { + _zero = zero; + Indices = new int[0]; + Values = new T[0]; + ValueCount = 0; + } + + /// + /// Retrieves the requested element without range checking. + /// + public override T At(int index) + { + // Search if item idex exists in NonZeroIndices array in range "0 - nonzero values count" + var itemIndex = Array.BinarySearch(Indices, 0, ValueCount, index); + return itemIndex >= 0 ? Values[itemIndex] : _zero; + } + + /// + /// Sets the element without range checking. + /// + public override void At(int index, T value) + { + // Search if "index" already exists in range "0 - nonzero values count" + var itemIndex = Array.BinarySearch(Indices, 0, ValueCount, index); + if (itemIndex >= 0) + { + // Non-zero item found in matrix + if (_zero.Equals(value)) + { + // Delete existing item + RemoveAtIndexUnchecked(itemIndex); + } + else + { + // Update item + Values[itemIndex] = value; + } + } + else + { + // Item not found. Add new value + if (!_zero.Equals(value)) + { + InsertAtIndexUnchecked(~itemIndex, index, value); + } + } + } + + void InsertAtIndexUnchecked(int itemIndex, int index, T value) + { + // Check if the storage needs to be increased + if ((ValueCount == Values.Length) && (ValueCount < Length)) + { + // Value and Indices arrays are completely full so we increase the size + var size = Math.Min(Values.Length + GrowthSize(), Length); + Array.Resize(ref Values, size); + Array.Resize(ref Indices, size); + } + + // Move all values (with a position larger than index) in the value array to the next position + // Move all values (with a position larger than index) in the columIndices array to the next position + Array.Copy(Values, itemIndex, Values, itemIndex + 1, ValueCount - itemIndex); + Array.Copy(Indices, itemIndex, Indices, itemIndex + 1, ValueCount - itemIndex); + + // Add the value and the column index + Values[itemIndex] = value; + Indices[itemIndex] = index; + + // increase the number of non-zero numbers by one + ValueCount += 1; + } + + void RemoveAtIndexUnchecked(int itemIndex) + { + // Value is zero. Let's delete it from Values and Indices array + Array.Copy(Values, itemIndex + 1, Values, itemIndex, ValueCount - itemIndex - 1); + Array.Copy(Indices, itemIndex + 1, Indices, itemIndex, ValueCount - itemIndex - 1); + + ValueCount -= 1; + + // Check whether we need to shrink the arrays. This is reasonable to do if + // there are a lot of non-zero elements and storage is two times bigger + if ((ValueCount > 1024) && (ValueCount < Indices.Length / 2)) + { + Array.Resize(ref Values, ValueCount); + Array.Resize(ref Indices, ValueCount); + } + } + + /// + /// Calculates the amount with which to grow the storage array's if they need to be + /// increased in size. + /// + /// The amount grown. + int GrowthSize() + { + int delta; + if (Values.Length > 1024) + { + delta = Values.Length / 4; + } + else + { + if (Values.Length > 256) + { + delta = 512; + } + else + { + delta = Values.Length > 64 ? 128 : 32; + } + } + + return delta; + } + + public override void Clear() + { + ValueCount = 0; + } + + public override void Clear(int index, int count) + { + if (index == 0 && count == Length) + { + Clear(); + return; + } + + var first = Array.BinarySearch(Indices, 0, ValueCount, index); + var last = Array.BinarySearch(Indices, 0, ValueCount, index + count - 1); + if (first < 0) first = ~first; + if (last < 0) last = ~last - 1; + int itemCount = last - first + 1; + + if (itemCount > 0) + { + Array.Copy(Values, first + count, Values, first, ValueCount - first - count); + Array.Copy(Indices, first + count, Indices, first, ValueCount - first - count); + + ValueCount -= count; + } + + // Check whether we need to shrink the arrays. This is reasonable to do if + // there are a lot of non-zero elements and storage is two times bigger + if ((ValueCount > 1024) && (ValueCount < Indices.Length / 2)) + { + Array.Resize(ref Values, ValueCount); + Array.Resize(ref Indices, ValueCount); + } + } + + /// Parameters assumed to be validated already. + public override void CopyTo(VectorStorage target, bool skipClearing = false) + { + var sparseTarget = target as SparseVectorStorage; + if (sparseTarget != null) + { + CopyTo(sparseTarget); + return; + } + + // FALL BACK + + if (!skipClearing) + { + target.Clear(); + } + + if (ValueCount != 0) + { + for (int i = 0; i < ValueCount; i++) + { + target.At(Indices[i], Values[i]); + } + } + } + + void CopyTo(SparseVectorStorage target) + { + if (ReferenceEquals(this, target)) + { + return; + } + + if (Length != target.Length) + { + var message = string.Format(Resources.ArgumentMatrixDimensions2, Length, target.Length); + throw new ArgumentException(message, "target"); + } + + target.ValueCount = ValueCount; + target.Values = new T[ValueCount]; + target.Indices = new int[ValueCount]; + + if (ValueCount != 0) + { + Array.Copy(Values, target.Values, ValueCount); + Buffer.BlockCopy(Indices, 0, target.Indices, 0, ValueCount * Constants.SizeOfInt); + } + } + + public override void CopySubVectorTo(VectorStorage target, + int sourceIndex, int targetIndex, int count, + bool skipClearing = false) + { + ValidateSubVectorRange(target, sourceIndex, targetIndex, count); + + var sparseTarget = target as SparseVectorStorage; + if (sparseTarget != null) + { + CopySubVectorTo(sparseTarget, sourceIndex, targetIndex, count, skipClearing); + return; + } + + // FALL BACK + + var offset = targetIndex - sourceIndex; + + var sourceFirst = Array.BinarySearch(Indices, 0, ValueCount, sourceIndex); + var sourceLast = Array.BinarySearch(Indices, 0, ValueCount, sourceIndex + count - 1); + if (sourceFirst < 0) sourceFirst = ~sourceFirst; + if (sourceLast < 0) sourceLast = ~sourceLast - 1; + + if (!skipClearing) + { + target.Clear(targetIndex, count); + } + + for (int i = sourceFirst; i <= sourceLast; i++) + { + target.At(Indices[i] + offset, Values[i]); + } + } + + void CopySubVectorTo(SparseVectorStorage target, + int sourceIndex, int targetIndex, int count, + bool skipClearing) + { + var offset = targetIndex - sourceIndex; + + var sourceFirst = Array.BinarySearch(Indices, 0, ValueCount, sourceIndex); + var sourceLast = Array.BinarySearch(Indices, 0, ValueCount, sourceIndex + count - 1); + if (sourceFirst < 0) sourceFirst = ~sourceFirst; + if (sourceLast < 0) sourceLast = ~sourceLast - 1; + int sourceCount = sourceLast - sourceFirst + 1; + + // special case when copying to itself + if (ReferenceEquals(this, target)) + { + var values = new T[sourceCount]; + var indices = new int[sourceCount]; + + Array.Copy(Values, sourceFirst, values, 0, sourceCount); + for (int i = 0; i < indices.Length; i++) + { + indices[i] = Indices[i + sourceFirst]; + } + + if (!skipClearing) + { + Clear(targetIndex, count); + } + + for (int i = sourceFirst; i <= sourceLast; i++) + { + At(indices[i] + offset, values[i]); + } + + return; + } + + // special case for empty target - much faster + if (target.ValueCount == 0) + { + var values = new T[sourceCount]; + var indices = new int[sourceCount]; + + Array.Copy(Values, sourceFirst, values, 0, sourceCount); + for (int i = 0; i < indices.Length; i++) + { + indices[i] = Indices[i + sourceFirst] + offset; + } + + target.ValueCount = sourceCount; + target.Values = values; + target.Indices = indices; + + return; + } + + if (!skipClearing) + { + target.Clear(targetIndex, count); + } + + for (int i = sourceFirst; i <= sourceLast; i++) + { + target.At(Indices[i] + offset, Values[i]); + } + } + } +} diff --git a/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs index 949573c2..2993a0b9 100644 --- a/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs @@ -184,12 +184,22 @@ namespace MathNet.Numerics.LinearAlgebra.Storage throw new ArgumentNullException("target"); } + ValidateSubVectorRange(target, sourceIndex, targetIndex, count); + if (ReferenceEquals(this, target)) { - throw new NotSupportedException(); - } + var tmp = new T[count]; + for (int i = 0; i < tmp.Length; i++) + { + tmp[i] = At(i + sourceIndex); + } + for (int i = 0; i < tmp.Length; i++) + { + At(i + targetIndex, tmp[i]); + } - ValidateSubVectorRange(target, sourceIndex, targetIndex, count); + return; + } for (int i = sourceIndex, ii = targetIndex; i < sourceIndex + count; i++, ii++) { diff --git a/src/Numerics/Numerics.csproj b/src/Numerics/Numerics.csproj index 4eda2e3b..8a6a9de7 100644 --- a/src/Numerics/Numerics.csproj +++ b/src/Numerics/Numerics.csproj @@ -124,6 +124,7 @@ + diff --git a/src/Portable/Portable.csproj b/src/Portable/Portable.csproj index 5030f2a2..cb907d41 100644 --- a/src/Portable/Portable.csproj +++ b/src/Portable/Portable.csproj @@ -906,6 +906,9 @@ LinearAlgebra\Storage\SparseCompressedRowMatrixStorage.cs + + LinearAlgebra\Storage\SparseVectorStorage.cs + LinearAlgebra\Storage\VectorStorage.cs diff --git a/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs index bae3cc3b..32bc2c74 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs @@ -35,6 +35,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex using System.Numerics; using LinearAlgebra.Complex; using LinearAlgebra.Generic; + using LinearAlgebra.Storage; using NUnit.Framework; /// @@ -273,40 +274,41 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex public void CheckSparseMechanismBySettingValues() { var vector = new SparseVector(10000); + var storage = (SparseVectorStorage)vector.Storage; // Add non-zero elements vector[200] = new Complex(1.5, 1); Assert.AreEqual(new Complex(1.5, 1), vector[200]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[500] = new Complex(3.5, 1); Assert.AreEqual(new Complex(3.5, 1), vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = new Complex(5.5, 1); Assert.AreEqual(new Complex(5.5, 1), vector[800]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[0] = new Complex(7.5, 1); Assert.AreEqual(new Complex(7.5, 1), vector[0]); - Assert.AreEqual(4, vector.NonZerosCount); + Assert.AreEqual(4, storage.ValueCount); // Remove non-zero elements vector[200] = Complex.Zero; Assert.AreEqual(Complex.Zero, vector[200]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[500] = Complex.Zero; Assert.AreEqual(Complex.Zero, vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = Complex.Zero; Assert.AreEqual(Complex.Zero, vector[800]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[0] = Complex.Zero; Assert.AreEqual(Complex.Zero, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -325,11 +327,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex // Multiply by 0 vector *= 0; + + var storage = (SparseVectorStorage)vector.Storage; Assert.AreEqual(Complex.Zero, vector[200]); Assert.AreEqual(Complex.Zero, vector[500]); Assert.AreEqual(Complex.Zero, vector[800]); Assert.AreEqual(Complex.Zero, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -371,7 +375,8 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex Assert.AreEqual(Data[i] * zeroArray[i], result[i]); } - Assert.AreEqual(2, result.NonZerosCount); + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(2, resultStorage.ValueCount); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex/UserDefinedVector.cs b/src/UnitTests/LinearAlgebraTests/Complex/UserDefinedVector.cs index 8f341ee8..c168d1e3 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex/UserDefinedVector.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex/UserDefinedVector.cs @@ -29,65 +29,56 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex using System.Numerics; using LinearAlgebra.Complex; using LinearAlgebra.Generic; + using LinearAlgebra.Storage; /// /// User-defined vector implementation (internal class for testing purposes) /// internal class UserDefinedVector : Vector { - /// - /// Values storage - /// - private readonly Complex[] _data; - - /// - /// Initializes a new instance of the class with a given size. - /// - /// The size of the vector. - public UserDefinedVector(int size) : base(size) + class UserDefinedVectorStorage : VectorStorage { - _data = new Complex[size]; - } + public readonly Complex[] Data; - /// - /// Initializes a new instance of the class for an array. - /// - /// The array to create this vector from. - public UserDefinedVector(Complex[] data) : base(data.Length) - { - _data = (Complex[])data.Clone(); - } + public UserDefinedVectorStorage(int size) + : base(size) + { + Data = new Complex[size]; + } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - public override Complex this[int index] - { - get + public UserDefinedVectorStorage(int size, Complex[] data) + : base(size) { - return _data[index]; + Data = data; } - set + public override Complex At(int index) { - _data[index] = value; + return Data[index]; + } + + public override void At(int index, Complex value) + { + Data[index] = value; } } - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - protected internal override Complex At(int index) + /// + /// Initializes a new instance of the class with a given size. + /// + /// The size of the vector. + public UserDefinedVector(int size) + : base(new UserDefinedVectorStorage(size)) { - return this[index]; } - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - protected internal override void At(int index, Complex value) + /// + /// Initializes a new instance of the class for an array. + /// + /// The array to create this vector from. + public UserDefinedVector(Complex[] data) + : base(new UserDefinedVectorStorage(data.Length, (Complex[])data.Clone())) { - this[index] = value; } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex/VectorTests.cs b/src/UnitTests/LinearAlgebraTests/Complex/VectorTests.cs index f0188b44..b5c962a6 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex/VectorTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex/VectorTests.cs @@ -97,7 +97,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex var vector = CreateVector(Data); var other = CreateVector(Data.Length); - vector.CopyTo(other, 2, 2, 2); + vector.CopySubVectorTo(other, 2, 2, 2); AssertHelpers.AreEqual(Complex.Zero, other[0]); AssertHelpers.AreEqual(Complex.Zero, other[1]); @@ -113,7 +113,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex public void CanCopyPartialVectorToSelf() { var vector = CreateVector(Data); - vector.CopyTo(vector, 0, 2, 2); + vector.CopySubVectorTo(vector, 0, 2, 2); AssertHelpers.AreEqual(new Complex(1, 1), vector[0]); AssertHelpers.AreEqual(new Complex(2, 1), vector[1]); diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Complex32/SparseVectorTest.cs index b0746766..eb7a40f4 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex32/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex32/SparseVectorTest.cs @@ -34,6 +34,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 using System.Collections.Generic; using LinearAlgebra.Complex32; using LinearAlgebra.Generic; + using LinearAlgebra.Storage; using NUnit.Framework; using Complex32 = Numerics.Complex32; @@ -273,40 +274,41 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 public void CheckSparseMechanismBySettingValues() { var vector = new SparseVector(10000); + var storage = (SparseVectorStorage)vector.Storage; // Add non-zero elements vector[200] = new Complex32(1.5f, 1); Assert.AreEqual(new Complex32(1.5f, 1), vector[200]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[500] = new Complex32(3.5f, 1); Assert.AreEqual(new Complex32(3.5f, 1), vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = new Complex32(5.5f, 1); Assert.AreEqual(new Complex32(5.5f, 1), vector[800]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[0] = new Complex32(7.5f, 1); Assert.AreEqual(new Complex32(7.5f, 1), vector[0]); - Assert.AreEqual(4, vector.NonZerosCount); + Assert.AreEqual(4, storage.ValueCount); // Remove non-zero elements vector[200] = Complex32.Zero; Assert.AreEqual(Complex32.Zero, vector[200]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[500] = Complex32.Zero; Assert.AreEqual(Complex32.Zero, vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = Complex32.Zero; Assert.AreEqual(Complex32.Zero, vector[800]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[0] = Complex32.Zero; Assert.AreEqual(Complex32.Zero, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -325,11 +327,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 // Multiply by 0 vector *= 0; + + var storage = (SparseVectorStorage)vector.Storage; Assert.AreEqual(Complex32.Zero, vector[200]); Assert.AreEqual(Complex32.Zero, vector[500]); Assert.AreEqual(Complex32.Zero, vector[800]); Assert.AreEqual(Complex32.Zero, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -371,7 +375,8 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 Assert.AreEqual(Data[i] * zeroArray[i], result[i]); } - Assert.AreEqual(2, result.NonZerosCount); + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(2, resultStorage.ValueCount); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/UserDefinedVector.cs b/src/UnitTests/LinearAlgebraTests/Complex32/UserDefinedVector.cs index 9b2fb9f3..f405ea1f 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex32/UserDefinedVector.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex32/UserDefinedVector.cs @@ -28,6 +28,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 { using LinearAlgebra.Complex32; using LinearAlgebra.Generic; + using LinearAlgebra.Storage; using Complex32 = Numerics.Complex32; /// @@ -35,59 +36,49 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 /// internal class UserDefinedVector : Vector { - /// - /// Values storage - /// - private readonly Complex32[] _data; - - /// - /// Initializes a new instance of the class with a given size. - /// - /// The size of the vector. - public UserDefinedVector(int size) : base(size) + class UserDefinedVectorStorage : VectorStorage { - _data = new Complex32[size]; - } + public readonly Complex32[] Data; - /// - /// Initializes a new instance of the class for an array. - /// - /// The array to create this vector from. - public UserDefinedVector(Complex32[] data) : base(data.Length) - { - _data = (Complex32[])data.Clone(); - } + public UserDefinedVectorStorage(int size) + : base(size) + { + Data = new Complex32[size]; + } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - public override Complex32 this[int index] - { - get + public UserDefinedVectorStorage(int size, Complex32[] data) + : base(size) { - return _data[index]; + Data = data; } - set + public override Complex32 At(int index) { - _data[index] = value; + return Data[index]; + } + + public override void At(int index, Complex32 value) + { + Data[index] = value; } } - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - protected internal override Complex32 At(int index) + /// + /// Initializes a new instance of the class with a given size. + /// + /// The size of the vector. + public UserDefinedVector(int size) + : base(new UserDefinedVectorStorage(size)) { - return this[index]; } - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - protected internal override void At(int index, Complex32 value) + /// + /// Initializes a new instance of the class for an array. + /// + /// The array to create this vector from. + public UserDefinedVector(Complex32[] data) + : base(new UserDefinedVectorStorage(data.Length, (Complex32[])data.Clone())) { - this[index] = value; } /// diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/VectorTests.cs b/src/UnitTests/LinearAlgebraTests/Complex32/VectorTests.cs index 0a478666..46f6fb57 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex32/VectorTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex32/VectorTests.cs @@ -97,7 +97,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 var vector = CreateVector(Data); var other = CreateVector(Data.Length); - vector.CopyTo(other, 2, 2, 2); + vector.CopySubVectorTo(other, 2, 2, 2); AssertHelpers.AreEqual(Complex32.Zero, other[0]); AssertHelpers.AreEqual(Complex32.Zero, other[1]); @@ -113,7 +113,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32 public void CanCopyPartialVectorToSelf() { var vector = CreateVector(Data); - vector.CopyTo(vector, 0, 2, 2); + vector.CopySubVectorTo(vector, 0, 2, 2); AssertHelpers.AreEqual(new Complex32(1, 1), vector[0]); AssertHelpers.AreEqual(new Complex32(2, 1), vector[1]); diff --git a/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs index 52b8bcf5..406a12e0 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs @@ -34,6 +34,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double using System.Collections.Generic; using LinearAlgebra.Double; using LinearAlgebra.Generic; + using LinearAlgebra.Storage; using NUnit.Framework; /// @@ -272,40 +273,41 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double public void CheckSparseMechanismBySettingValues() { var vector = new SparseVector(10000); + var storage = (SparseVectorStorage)vector.Storage; // Add non-zero elements vector[200] = 1.5; Assert.AreEqual(1.5, vector[200]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[500] = 3.5; Assert.AreEqual(3.5, vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = 5.5; Assert.AreEqual(5.5, vector[800]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[0] = 7.5; Assert.AreEqual(7.5, vector[0]); - Assert.AreEqual(4, vector.NonZerosCount); + Assert.AreEqual(4, storage.ValueCount); // Remove non-zero elements vector[200] = 0; Assert.AreEqual(0, vector[200]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[500] = 0; Assert.AreEqual(0, vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = 0; Assert.AreEqual(0, vector[800]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[0] = 0; Assert.AreEqual(0, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -324,11 +326,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double // Multiply by 0 vector *= 0; + + var storage = (SparseVectorStorage)vector.Storage; Assert.AreEqual(0, vector[200]); Assert.AreEqual(0, vector[500]); Assert.AreEqual(0, vector[800]); Assert.AreEqual(0, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -370,7 +374,8 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double Assert.AreEqual(Data[i] * zeroArray[i], result[i]); } - Assert.AreEqual(2, result.NonZerosCount); + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(2, resultStorage.ValueCount); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Double/UserDefinedVector.cs b/src/UnitTests/LinearAlgebraTests/Double/UserDefinedVector.cs index 03a7f618..ace40e08 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/UserDefinedVector.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/UserDefinedVector.cs @@ -28,65 +28,56 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double { using LinearAlgebra.Double; using LinearAlgebra.Generic; + using LinearAlgebra.Storage; /// /// User-defined vector implementation (internal class for testing purposes) /// internal class UserDefinedVector : Vector { - /// - /// Values storage - /// - private readonly double[] _data; - - /// - /// Initializes a new instance of the class with a given size. - /// - /// The size of the vector. - public UserDefinedVector(int size) : base(size) + class UserDefinedVectorStorage : VectorStorage { - _data = new double[size]; - } + public readonly double[] Data; - /// - /// Initializes a new instance of the class for an array. - /// - /// The array to create this vector from. - public UserDefinedVector(double[] data) : base(data.Length) - { - _data = (double[])data.Clone(); - } + public UserDefinedVectorStorage(int size) + : base(size) + { + Data = new double[size]; + } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - public override double this[int index] - { - get + public UserDefinedVectorStorage(int size, double[] data) + : base(size) { - return _data[index]; + Data = data; } - set + public override double At(int index) { - _data[index] = value; + return Data[index]; + } + + public override void At(int index, double value) + { + Data[index] = value; } } - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - protected internal override double At(int index) + /// + /// Initializes a new instance of the class with a given size. + /// + /// The size of the vector. + public UserDefinedVector(int size) + : base(new UserDefinedVectorStorage(size)) { - return this[index]; } - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - protected internal override void At(int index, double value) + /// + /// Initializes a new instance of the class for an array. + /// + /// The array to create this vector from. + public UserDefinedVector(double[] data) + : base(new UserDefinedVectorStorage(data.Length, (double[])data.Clone())) { - this[index] = value; } /// diff --git a/src/UnitTests/LinearAlgebraTests/Double/VectorTests.cs b/src/UnitTests/LinearAlgebraTests/Double/VectorTests.cs index 75fc6035..3c128ea3 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/VectorTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/VectorTests.cs @@ -95,7 +95,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double var vector = CreateVector(Data); var other = CreateVector(Data.Length); - vector.CopyTo(other, 2, 2, 2); + vector.CopySubVectorTo(other, 2, 2, 2); Assert.AreEqual(0.0, other[0]); Assert.AreEqual(0.0, other[1]); @@ -111,7 +111,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double public void CanCopyPartialVectorToSelf() { var vector = CreateVector(Data); - vector.CopyTo(vector, 0, 2, 2); + vector.CopySubVectorTo(vector, 0, 2, 2); Assert.AreEqual(1.0, vector[0]); Assert.AreEqual(2.0, vector[1]); diff --git a/src/UnitTests/LinearAlgebraTests/Single/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Single/SparseVectorTest.cs index 435632e0..24a6146f 100644 --- a/src/UnitTests/LinearAlgebraTests/Single/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Single/SparseVectorTest.cs @@ -34,7 +34,9 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single using System.Collections.Generic; using LinearAlgebra.Generic; using LinearAlgebra.Single; + using LinearAlgebra.Storage; using NUnit.Framework; + using Complex32 = Numerics.Complex32; /// /// Sparse vector tests. @@ -272,40 +274,41 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single public void CheckSparseMechanismBySettingValues() { var vector = new SparseVector(10000); + var storage = (SparseVectorStorage)vector.Storage; // Add non-zero elements vector[200] = 1.5f; Assert.AreEqual(1.5f, vector[200]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[500] = 3.5f; Assert.AreEqual(3.5f, vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = 5.5f; Assert.AreEqual(5.5f, vector[800]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[0] = 7.5f; Assert.AreEqual(7.5f, vector[0]); - Assert.AreEqual(4, vector.NonZerosCount); + Assert.AreEqual(4, storage.ValueCount); // Remove non-zero elements vector[200] = 0; Assert.AreEqual(0, vector[200]); - Assert.AreEqual(3, vector.NonZerosCount); + Assert.AreEqual(3, storage.ValueCount); vector[500] = 0; Assert.AreEqual(0, vector[500]); - Assert.AreEqual(2, vector.NonZerosCount); + Assert.AreEqual(2, storage.ValueCount); vector[800] = 0; Assert.AreEqual(0, vector[800]); - Assert.AreEqual(1, vector.NonZerosCount); + Assert.AreEqual(1, storage.ValueCount); vector[0] = 0; Assert.AreEqual(0, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -324,11 +327,13 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single // Multiply by 0 vector *= 0; + + var storage = (SparseVectorStorage)vector.Storage; Assert.AreEqual(0, vector[200]); Assert.AreEqual(0, vector[500]); Assert.AreEqual(0, vector[800]); Assert.AreEqual(0, vector[0]); - Assert.AreEqual(0, vector.NonZerosCount); + Assert.AreEqual(0, storage.ValueCount); } /// @@ -370,7 +375,8 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single Assert.AreEqual(Data[i] * zeroArray[i], result[i]); } - Assert.AreEqual(2, result.NonZerosCount); + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(2, resultStorage.ValueCount); } /// diff --git a/src/UnitTests/LinearAlgebraTests/Single/UserDefinedVector.cs b/src/UnitTests/LinearAlgebraTests/Single/UserDefinedVector.cs index 37f0a9e9..6700b4a2 100644 --- a/src/UnitTests/LinearAlgebraTests/Single/UserDefinedVector.cs +++ b/src/UnitTests/LinearAlgebraTests/Single/UserDefinedVector.cs @@ -28,65 +28,56 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single { using LinearAlgebra.Generic; using LinearAlgebra.Single; + using LinearAlgebra.Storage; /// /// User-defined vector implementation (internal class for testing purposes) /// internal class UserDefinedVector : Vector { - /// - /// Values storage - /// - private readonly float[] _data; - - /// - /// Initializes a new instance of the class with a given size. - /// - /// The size of the vector. - public UserDefinedVector(int size) : base(size) + class UserDefinedVectorStorage : VectorStorage { - _data = new float[size]; - } + public readonly float[] Data; - /// - /// Initializes a new instance of the class for an array. - /// - /// The array to create this vector from. - public UserDefinedVector(float[] data) : base(data.Length) - { - _data = (float[])data.Clone(); - } + public UserDefinedVectorStorage(int size) + : base(size) + { + Data = new float[size]; + } - /// Gets or sets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - public override float this[int index] - { - get + public UserDefinedVectorStorage(int size, float[] data) + : base(size) { - return _data[index]; + Data = data; } - set + public override float At(int index) { - _data[index] = value; + return Data[index]; + } + + public override void At(int index, float value) + { + Data[index] = value; } } - /// Gets the value at the given . - /// The index of the value to get or set. - /// The value of the vector at the given . - protected internal override float At(int index) + /// + /// Initializes a new instance of the class with a given size. + /// + /// The size of the vector. + public UserDefinedVector(int size) + : base(new UserDefinedVectorStorage(size)) { - return this[index]; } - /// Sets the at the given . - /// The index of the value to get or set. - /// The value to set. - protected internal override void At(int index, float value) + /// + /// Initializes a new instance of the class for an array. + /// + /// The array to create this vector from. + public UserDefinedVector(float[] data) + : base(new UserDefinedVectorStorage(data.Length, (float[])data.Clone())) { - this[index] = value; } /// diff --git a/src/UnitTests/LinearAlgebraTests/Single/VectorTests.cs b/src/UnitTests/LinearAlgebraTests/Single/VectorTests.cs index 9fff4251..a84ef969 100644 --- a/src/UnitTests/LinearAlgebraTests/Single/VectorTests.cs +++ b/src/UnitTests/LinearAlgebraTests/Single/VectorTests.cs @@ -95,7 +95,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single var vector = CreateVector(Data); var other = CreateVector(Data.Length); - vector.CopyTo(other, 2, 2, 2); + vector.CopySubVectorTo(other, 2, 2, 2); Assert.AreEqual(0.0f, other[0]); Assert.AreEqual(0.0f, other[1]); @@ -111,7 +111,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single public void CanCopyPartialVectorToSelf() { var vector = CreateVector(Data); - vector.CopyTo(vector, 0, 2, 2); + vector.CopySubVectorTo(vector, 0, 2, 2); Assert.AreEqual(1.0f, vector[0]); Assert.AreEqual(2.0f, vector[1]);