diff --git a/src/Numerics/LinearAlgebra/Double/SparseVector.cs b/src/Numerics/LinearAlgebra/Double/SparseVector.cs index 4b58d572..9ec18d65 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseVector.cs @@ -153,7 +153,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double if (ReferenceEquals(this, result)) { - //populate a new vector with the scalar + // populate a new vector with the scalar var vnonZeroValues = new double[Count]; var vnonZeroIndices = new int[Count]; for (int index = 0; index < Count; index++) @@ -162,7 +162,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double vnonZeroValues[index] = scalar; } - //populate the non zero values from this + // populate the non zero values from this var indices = _storage.Indices; var values = _storage.Values; for (int j = 0; j < _storage.ValueCount; j++) @@ -170,7 +170,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double vnonZeroValues[indices[j]] = values[j] + scalar; } - //assign this vectors arrary to the new arrays. + // assign this vectors array to the new arrays. _storage.Values = vnonZeroValues; _storage.Indices = vnonZeroIndices; _storage.ValueCount = Count; diff --git a/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs index a83a5f3f..e92adce9 100644 --- a/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/DenseVectorStorage.cs @@ -172,6 +172,29 @@ namespace MathNet.Numerics.LinearAlgebra.Storage { Array.Copy(Data, 0, denseTarget.Data, 0, Data.Length); } + + return; + } + + var sparseTarget = target as SparseVectorStorage; + if (sparseTarget != null) + { + var indices = new List(); + var values = new List(); + + for (int i = 0; i < Data.Length; i++) + { + var item = Data[i]; + if (!Zero.Equals(item)) + { + values.Add(item); + indices.Add(i); + } + } + + sparseTarget.Indices = indices.ToArray(); + sparseTarget.Values = values.ToArray(); + sparseTarget.ValueCount = values.Count; return; } @@ -364,5 +387,60 @@ namespace MathNet.Numerics.LinearAlgebra.Storage target.At(i, f(i, Data[i])); } } + + internal override void Map2ToUnchecked(VectorStorage target, VectorStorage other, Func f, Zeros zeros = Zeros.AllowSkip, ExistingData existingData = ExistingData.Clear) + { + if (target is SparseVectorStorage) + { + // Recursive to dense target at first, since the operation is + // effectively dense anyway because at least one operand is dense + var intermediate = new DenseVectorStorage(target.Length); + Map2ToUnchecked(intermediate, other, f, zeros, ExistingData.AssumeZeros); + intermediate.CopyTo(target, existingData); + return; + } + + var denseTarget = target as DenseVectorStorage; + var denseOther = other as DenseVectorStorage; + if (denseTarget != null && denseOther != null) + { + CommonParallel.For(0, Data.Length, 4096, (a, b) => + { + for (int i = a; i < b; i++) + { + denseTarget.Data[i] = f(Data[i], denseOther.Data[i]); + } + }); + + return; + } + + var sparseOther = other as SparseVectorStorage; + if (denseTarget != null && sparseOther != null) + { + T[] targetData = denseTarget.Data; + int[] otherIndices = sparseOther.Indices; + T[] otherValues = sparseOther.Values; + int otherValueCount = sparseOther.ValueCount; + + int k = 0; + for (int i = 0; i < Data.Length; i++) + { + if (k < otherValueCount && otherIndices[k] == i) + { + targetData[i] = f(Data[i], otherValues[k]); + k++; + } + else + { + targetData[i] = f(Data[i], Zero); + } + } + + return; + } + + base.Map2ToUnchecked(target, other, f, zeros, existingData); + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs index 115c3584..eb95135f 100644 --- a/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs @@ -30,6 +30,7 @@ using System; using System.Collections.Generic; +using System.Diagnostics; using System.Linq; using MathNet.Numerics.Properties; using MathNet.Numerics.Threading; @@ -768,5 +769,156 @@ namespace MathNet.Numerics.LinearAlgebra.Storage base.MapIndexedToUnchecked(target, f, zeros, existingData); } + + internal override void Map2ToUnchecked(VectorStorage target, VectorStorage other, Func f, Zeros zeros = Zeros.AllowSkip, ExistingData existingData = ExistingData.Clear) + { + var processZeros = zeros == Zeros.Include || !Zero.Equals(f(Zero, Zero)); + + var denseTarget = target as DenseVectorStorage; + var denseOther = other as DenseVectorStorage; + + if (denseTarget == null && (denseOther != null || processZeros)) + { + // The handling is effectively dense but we're supposed to push + // to a sparse target. Let's use a dense target instead, + // then copy it normalized back to the sparse target. + var intermediate = new DenseVectorStorage(target.Length); + Map2ToUnchecked(intermediate, other, f, zeros, ExistingData.AssumeZeros); + intermediate.CopyTo(target, existingData); + return; + } + + if (denseOther != null) + { + T[] targetData = denseTarget.Data; + T[] otherData = denseOther.Data; + + int k = 0; + for (int i = 0; i < otherData.Length; i++) + { + if (k < ValueCount && Indices[k] == i) + { + targetData[i] = f(Values[k], otherData[i]); + k++; + } + else + { + targetData[i] = f(Zero, otherData[i]); + } + } + + return; + } + + var sparseOther = other as SparseVectorStorage; + if (sparseOther != null && denseTarget != null) + { + T[] targetData = denseTarget.Data; + int[] otherIndices = sparseOther.Indices; + T[] otherValues = sparseOther.Values; + int otherValueCount = sparseOther.ValueCount; + + if (processZeros) + { + int p = 0, q = 0; + for (int i = 0; i < targetData.Length; i++) + { + var left = p < ValueCount && Indices[p] == i ? Values[p++] : Zero; + var right = q < otherValueCount && otherIndices[q] == i ? otherValues[q++] : Zero; + targetData[i] = f(left, right); + } + } + else + { + if (existingData == ExistingData.Clear) + { + denseTarget.Clear(); + } + + int p = 0, q = 0; + while (p < ValueCount || q < otherValueCount) + { + if (q >= otherValueCount || p < ValueCount && Indices[p] < otherIndices[q]) + { + targetData[Indices[p]] = f(Values[p], Zero); + p++; + } + else if (p >= ValueCount || q < otherValueCount && Indices[p] > otherIndices[q]) + { + targetData[otherIndices[q]] = f(Zero, otherValues[q]); + q++; + } + else + { + Debug.Assert(Indices[p] == otherIndices[q]); + targetData[Indices[p]] = f(Values[p], otherValues[q]); + p++; + q++; + } + } + } + + return; + } + + var sparseTarget = target as SparseVectorStorage; + if (sparseOther != null && sparseTarget != null) + { + var indices = new List(); + var values = new List(); + int[] otherIndices = sparseOther.Indices; + T[] otherValues = sparseOther.Values; + int otherValueCount = sparseOther.ValueCount; + + int p = 0, q = 0; + while (p < ValueCount || q < otherValueCount) + { + if (q >= otherValueCount || p < ValueCount && Indices[p] < otherIndices[q]) + { + var value = f(Values[p], Zero); + if (!Zero.Equals(value)) + { + indices.Add(Indices[p]); + values.Add(value); + } + + p++; + } + else if (p >= ValueCount || q < otherValueCount && Indices[p] > otherIndices[q]) + { + var value = f(Zero, otherValues[q]); + if (!Zero.Equals(value)) + { + indices.Add(otherIndices[q]); + values.Add(value); + } + + q++; + } + else + { + Debug.Assert(Indices[p] == otherIndices[q]); + var value = f(Values[p], otherValues[q]); + if (!Zero.Equals(value)) + { + indices.Add(Indices[p]); + values.Add(value); + } + + p++; + q++; + } + } + + sparseTarget.Indices = indices.ToArray(); + sparseTarget.Values = values.ToArray(); + sparseTarget.ValueCount = values.Count; + return; + } + + // FALL BACK + + base.Map2ToUnchecked(target, other, f, zeros, existingData); + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs index b44fb7cf..bb1edeb0 100644 --- a/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/VectorStorage.cs @@ -469,5 +469,40 @@ namespace MathNet.Numerics.LinearAlgebra.Storage target.At(i, f(i, At(i))); } } + + public void Map2To(VectorStorage target, VectorStorage other, Func f, + Zeros zeros = Zeros.AllowSkip, ExistingData existingData = ExistingData.Clear) + { + if (target == null) + { + throw new ArgumentNullException("target"); + } + + if (other == null) + { + throw new ArgumentNullException("other"); + } + + if (Length != target.Length) + { + throw new ArgumentException(Resources.ArgumentVectorsSameLength, "target"); + } + + if (Length != other.Length) + { + throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other"); + } + + Map2ToUnchecked(target, other, f, zeros, existingData); + } + + internal virtual void Map2ToUnchecked(VectorStorage target, VectorStorage other, Func f, + Zeros zeros = Zeros.AllowSkip, ExistingData existingData = ExistingData.Clear) + { + for (int i = 0; i < Length; i++) + { + target.At(i, f(At(i), other.At(i))); + } + } } } diff --git a/src/Numerics/LinearAlgebra/Vector.cs b/src/Numerics/LinearAlgebra/Vector.cs index 04b168fd..4056ff89 100644 --- a/src/Numerics/LinearAlgebra/Vector.cs +++ b/src/Numerics/LinearAlgebra/Vector.cs @@ -457,5 +457,23 @@ namespace MathNet.Numerics.LinearAlgebra Storage.MapIndexedToUnchecked(result.Storage, f, zeros, ExistingData.AssumeZeros); return result; } + + /// + /// Applies a function to each value pair of two vectors and replaces the value in the result vector. + /// + public void Map2(Func f, Vector other, Vector result, Zeros zeros = Zeros.AllowSkip) + { + Storage.Map2To(result.Storage, other.Storage, f, zeros, ExistingData.Clear); + } + + /// + /// Applies a function to each value pair of two vectors and returns the results as a new vector. + /// + public Vector Map2(Func f, Vector other, Zeros zeros = Zeros.AllowSkip) + { + var result = Build.SameAs(this); + Storage.Map2To(result.Storage, other.Storage, f, zeros, ExistingData.AssumeZeros); + return result; + } } } diff --git a/src/UnitTests/LinearAlgebraTests/VectorStorageTests.cs b/src/UnitTests/LinearAlgebraTests/VectorStorageTests.cs new file mode 100644 index 00000000..87e6637d --- /dev/null +++ b/src/UnitTests/LinearAlgebraTests/VectorStorageTests.cs @@ -0,0 +1,184 @@ +// +// Math.NET Numerics, part of the Math.NET Project +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com +// +// Copyright (c) 2009-2014 Math.NET +// +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +using MathNet.Numerics.LinearAlgebra; +using MathNet.Numerics.LinearAlgebra.Storage; +using NUnit.Framework; + +namespace MathNet.Numerics.UnitTests.LinearAlgebraTests +{ + [TestFixture, Category("LA")] + public class VectorStorageTests + { + [Test] + public void Map2SkipZeros() + { + double[] a = { 1.0, 2.0, 0.0, 4.0, 0.0, 6.0 }; + double[] b = { 11.0, 12.0, 13.0, 0.0, 0.0, 16.0 }; + var adense = DenseVectorStorage.OfEnumerable(a); + var asparse = SparseVectorStorage.OfEnumerable(a); + var bdense = DenseVectorStorage.OfEnumerable(b); + var bsparse = SparseVectorStorage.OfEnumerable(b); + + var rdense = new DenseVectorStorage(a.Length); + var rsparse = new SparseVectorStorage(a.Length); + + var expected = new DenseVectorStorage(6, new[] { 12.0, 14.0, 13.0, 4.0, 0.0, 22.0 }); + + rdense.Clear(); + adense.Map2To(rdense, bdense, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "dense*dense->dense"); + + rsparse.Clear(); + adense.Map2To(rsparse, bdense, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "dense*dense->sparse"); + + rdense.Clear(); + adense.Map2To(rdense, bsparse, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "dense*sparse->dense"); + + rsparse.Clear(); + adense.Map2To(rsparse, bsparse, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "dense*sparse->sparse"); + + rdense.Clear(); + asparse.Map2To(rdense, bdense, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "sparse*dense->dense"); + + rsparse.Clear(); + asparse.Map2To(rsparse, bdense, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "sparse*dense->sparse"); + + rdense.Clear(); + asparse.Map2To(rdense, bsparse, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "sparse*sparse->dense"); + + rsparse.Clear(); + asparse.Map2To(rsparse, bsparse, (u, v) => u+v, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "sparse*sparse->sparse"); + } + + [Test] + public void Map2ForceIncludeZeros() + { + double[] a = { 1.0, 2.0, 0.0, 4.0, 0.0, 6.0 }; + double[] b = { 11.0, 12.0, 13.0, 0.0, 0.0, 16.0 }; + var adense = DenseVectorStorage.OfEnumerable(a); + var asparse = SparseVectorStorage.OfEnumerable(a); + var bdense = DenseVectorStorage.OfEnumerable(b); + var bsparse = SparseVectorStorage.OfEnumerable(b); + + var rdense = new DenseVectorStorage(a.Length); + var rsparse = new SparseVectorStorage(a.Length); + + var expected = new DenseVectorStorage(6, new[] { 13.0, 15.0, 14.0, 5.0, 1.0, 23.0 }); + + rdense.Clear(); + adense.Map2To(rdense, bdense, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rdense.Equals(expected), "dense*dense->dense"); + + rsparse.Clear(); + adense.Map2To(rsparse, bdense, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rsparse.Equals(expected), "dense*dense->sparse"); + + rdense.Clear(); + adense.Map2To(rdense, bsparse, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rdense.Equals(expected), "dense*sparse->dense"); + + rsparse.Clear(); + adense.Map2To(rsparse, bsparse, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rsparse.Equals(expected), "dense*sparse->sparse"); + + rdense.Clear(); + asparse.Map2To(rdense, bdense, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rdense.Equals(expected), "sparse*dense->dense"); + + rsparse.Clear(); + asparse.Map2To(rsparse, bdense, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rsparse.Equals(expected), "sparse*dense->sparse"); + + rdense.Clear(); + asparse.Map2To(rdense, bsparse, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rdense.Equals(expected), "sparse*sparse->dense"); + + rsparse.Clear(); + asparse.Map2To(rsparse, bsparse, (u, v) => u+v+1.0, Zeros.Include); + Assert.That(rsparse.Equals(expected), "sparse*sparse->sparse"); + } + + [Test] + public void Map2AutoIncludeZeros() + { + double[] a = { 1.0, 2.0, 0.0, 4.0, 0.0, 6.0 }; + double[] b = { 11.0, 12.0, 13.0, 0.0, 0.0, 16.0 }; + var adense = DenseVectorStorage.OfEnumerable(a); + var asparse = SparseVectorStorage.OfEnumerable(a); + var bdense = DenseVectorStorage.OfEnumerable(b); + var bsparse = SparseVectorStorage.OfEnumerable(b); + + var rdense = new DenseVectorStorage(a.Length); + var rsparse = new SparseVectorStorage(a.Length); + + var expected = new DenseVectorStorage(6, new[] { 13.0, 15.0, 14.0, 5.0, 1.0, 23.0 }); + + rdense.Clear(); + adense.Map2To(rdense, bdense, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "dense*dense->dense"); + + rsparse.Clear(); + adense.Map2To(rsparse, bdense, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "dense*dense->sparse"); + + rdense.Clear(); + adense.Map2To(rdense, bsparse, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "dense*sparse->dense"); + + rsparse.Clear(); + adense.Map2To(rsparse, bsparse, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "dense*sparse->sparse"); + + rdense.Clear(); + asparse.Map2To(rdense, bdense, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "sparse*dense->dense"); + + rsparse.Clear(); + asparse.Map2To(rsparse, bdense, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "sparse*dense->sparse"); + + rdense.Clear(); + asparse.Map2To(rdense, bsparse, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rdense.Equals(expected), "sparse*sparse->dense"); + + rsparse.Clear(); + asparse.Map2To(rsparse, bsparse, (u, v) => u+v+1.0, Zeros.AllowSkip); + Assert.That(rsparse.Equals(expected), "sparse*sparse->sparse"); + } + } +} diff --git a/src/UnitTests/UnitTests.csproj b/src/UnitTests/UnitTests.csproj index fe2b9319..31dbc0b7 100644 --- a/src/UnitTests/UnitTests.csproj +++ b/src/UnitTests/UnitTests.csproj @@ -350,6 +350,7 @@ +