From b3bb2d1e73ed0ea97a65fc37a58b0ef9c40ed84a Mon Sep 17 00:00:00 2001 From: Christoph Ruegg Date: Mon, 18 Apr 2016 20:20:13 +0200 Subject: [PATCH] Linear Algebra: fix bug in sparse vector pointwise multiply/divide to itself #390 --- .../LinearAlgebra/Complex/SparseVector.cs | 32 +--------- .../LinearAlgebra/Complex32/SparseVector.cs | 32 +--------- .../LinearAlgebra/Double/SparseVector.cs | 32 +--------- .../LinearAlgebra/Single/SparseVector.cs | 32 +--------- .../Complex/SparseVectorTest.cs | 2 +- .../Double/SparseVectorTest.cs | 64 +++++++++++++++++++ 6 files changed, 73 insertions(+), 121 deletions(-) diff --git a/src/Numerics/LinearAlgebra/Complex/SparseVector.cs b/src/Numerics/LinearAlgebra/Complex/SparseVector.cs index 2a2b3553..cbb04ced 100644 --- a/src/Numerics/LinearAlgebra/Complex/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex/SparseVector.cs @@ -778,7 +778,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex /// The vector to store the result of the pointwise multiplication. protected override void DoPointwiseMultiply(Vector other, Vector result) { - if (ReferenceEquals(this, other)) + if (ReferenceEquals(this, other) && ReferenceEquals(this, result)) { for (var i = 0; i < _storage.ValueCount; i++) { @@ -787,35 +787,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex } else { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, other.At(index) * _storage.Values[i]); - } - } - } - - /// - /// Pointwise multiplies this vector with another vector and stores the result into the result vector. - /// - /// The vector to pointwise multiply with this one. - /// The vector to store the result of the pointwise multiplication. - protected override void DoPointwiseDivide(Vector divisor, Vector result) - { - if (ReferenceEquals(this, divisor)) - { - for (var i = 0; i < _storage.ValueCount; i++) - { - _storage.Values[i] /= _storage.Values[i]; - } - } - else - { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, _storage.Values[i] / divisor.At(index)); - } + base.DoPointwiseMultiply(other, result); } } diff --git a/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs b/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs index bff2fe71..3a3193c9 100644 --- a/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Complex32/SparseVector.cs @@ -773,7 +773,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 /// The vector to store the result of the pointwise multiplication. protected override void DoPointwiseMultiply(Vector other, Vector result) { - if (ReferenceEquals(this, other)) + if (ReferenceEquals(this, other) && ReferenceEquals(this, result)) { for (var i = 0; i < _storage.ValueCount; i++) { @@ -782,35 +782,7 @@ namespace MathNet.Numerics.LinearAlgebra.Complex32 } else { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, other.At(index) * _storage.Values[i]); - } - } - } - - /// - /// Pointwise multiplies this vector with another vector and stores the result into the result vector. - /// - /// The vector to pointwise multiply with this one. - /// The vector to store the result of the pointwise multiplication. - protected override void DoPointwiseDivide(Vector divisor, Vector result) - { - if (ReferenceEquals(this, divisor)) - { - for (var i = 0; i < _storage.ValueCount; i++) - { - _storage.Values[i] /= _storage.Values[i]; - } - } - else - { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, _storage.Values[i] / divisor.At(index)); - } + base.DoPointwiseMultiply(other, result); } } diff --git a/src/Numerics/LinearAlgebra/Double/SparseVector.cs b/src/Numerics/LinearAlgebra/Double/SparseVector.cs index 27dcaab9..4704c857 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseVector.cs @@ -820,7 +820,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double /// The vector to store the result of the pointwise multiplication. protected override void DoPointwiseMultiply(Vector other, Vector result) { - if (ReferenceEquals(this, other)) + if (ReferenceEquals(this, other) && ReferenceEquals(this, result)) { for (var i = 0; i < _storage.ValueCount; i++) { @@ -829,35 +829,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double } else { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, other.At(index) * _storage.Values[i]); - } - } - } - - /// - /// Pointwise multiplies this vector with another vector and stores the result into the result vector. - /// - /// The vector to pointwise multiply with this one. - /// The vector to store the result of the pointwise multiplication. - protected override void DoPointwiseDivide(Vector divisor, Vector result) - { - if (ReferenceEquals(this, divisor)) - { - for (var i = 0; i < _storage.ValueCount; i++) - { - _storage.Values[i] /= _storage.Values[i]; - } - } - else - { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, _storage.Values[i] / divisor.At(index)); - } + base.DoPointwiseMultiply(other, result); } } diff --git a/src/Numerics/LinearAlgebra/Single/SparseVector.cs b/src/Numerics/LinearAlgebra/Single/SparseVector.cs index 3321caa0..0a88758d 100644 --- a/src/Numerics/LinearAlgebra/Single/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Single/SparseVector.cs @@ -821,7 +821,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single /// The vector to store the result of the pointwise multiplication. protected override void DoPointwiseMultiply(Vector other, Vector result) { - if (ReferenceEquals(this, other)) + if (ReferenceEquals(this, other) && ReferenceEquals(this, result)) { for (var i = 0; i < _storage.ValueCount; i++) { @@ -830,35 +830,7 @@ namespace MathNet.Numerics.LinearAlgebra.Single } else { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, other.At(index) * _storage.Values[i]); - } - } - } - - /// - /// Pointwise multiplies this vector with another vector and stores the result into the result vector. - /// - /// The vector to pointwise multiply with this one. - /// The vector to store the result of the pointwise multiplication. - protected override void DoPointwiseDivide(Vector divisor, Vector result) - { - if (ReferenceEquals(this, divisor)) - { - for (var i = 0; i < _storage.ValueCount; i++) - { - _storage.Values[i] /= _storage.Values[i]; - } - } - else - { - for (var i = 0; i < _storage.ValueCount; i++) - { - var index = _storage.Indices[i]; - result.At(index, _storage.Values[i] / divisor.At(index)); - } + base.DoPointwiseMultiply(other, result); } } diff --git a/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs index 77d14440..54d7fb7f 100644 --- a/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Complex/SparseVectorTest.cs @@ -362,7 +362,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex } /// - /// Test for issues #52. When setting previous non-zero values to zero, + /// Test for issue #52. When setting previous non-zero values to zero, /// DoMultiply would copy non-zero values to the result, but use the /// length of nonzerovalues instead of NonZerosCount. /// diff --git a/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs index 4ce7c23f..b5e73bb3 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs @@ -355,6 +355,70 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double Assert.AreEqual(2, resultStorage.ValueCount); } + /// + /// Can pointwise divide a sparse vector. + /// + [Test] + public void CanPointwiseDivideSparseVector() + { + var zeroArray = new[] { 0.0, 1.0, 0.0, 1.0, 0.0 }; + var vector1 = Vector.Build.SparseOfEnumerable(zeroArray); + var vector2 = Vector.Build.SparseOfEnumerable(Data); + var result = Vector.Build.Sparse(vector1.Count); + + vector1.PointwiseDivide(vector2, result); + + for (var i = 0; i < vector1.Count; i++) + { + Assert.AreEqual(zeroArray[i] / Data[i], result[i]); + } + + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(2, resultStorage.ValueCount); + } + + /// + /// Can pointwise multiple a sparse vector to itself. + /// + [Test] + public void CanPointwiseMultiplyToSelfSparseVector() + { + var zeroArray = new[] { 0.0, 2.0, 0.0, 2.0, 0.0 }; + var vector = Vector.Build.SparseOfEnumerable(zeroArray); + var result = Vector.Build.Sparse(vector.Count); + + vector.PointwiseMultiply(vector, result); + + for (var i = 0; i < vector.Count; i++) + { + Assert.AreEqual(zeroArray[i] * zeroArray[i], result[i]); + } + + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(2, resultStorage.ValueCount); + } + + /// + /// Can pointwise divide a sparse vector to itself. + /// + [Test] + public void CanPointwiseDivideToSelfSparseVector() + { + var zeroArray = new[] { 0.0, 2.0, 0.0, 2.0, 0.0 }; + var vector = Vector.Build.SparseOfEnumerable(zeroArray); + var result = Vector.Build.Sparse(vector.Count); + + vector.PointwiseDivide(vector, result); + + for (var i = 0; i < vector.Count; i++) + { + Assert.AreEqual(zeroArray[i] / zeroArray[i], result[i]); + } + + var resultStorage = (SparseVectorStorage)result.Storage; + Assert.AreEqual(5, resultStorage.ValueCount); + } + /// /// Can outer multiple two sparse vectors. Checking fix for workitem 5696. ///