From 941f0ca790b070490f06d07ddf829bdb719470cd Mon Sep 17 00:00:00 2001 From: Jurgen Van Gael Date: Mon, 26 Jul 2010 05:10:30 +0800 Subject: [PATCH] Added sparse interface to the F# front-end. Fixed unit tests for new functional friendly vector interface. --- src/FSharp/DenseMatrix.fs | 20 +-- src/FSharp/DenseVector.fs | 4 +- src/FSharp/FSharp.fsproj | 4 +- src/FSharp/Main.fs | 4 +- src/FSharp/Matrix.fs | 4 +- src/FSharp/SparseMatrix.fs | 75 +++++++++ src/FSharp/SparseVector.fs | 48 ++++++ src/FSharp/Vector.fs | 10 +- src/FSharpUnitTests/Program.fs | 52 ++++-- .../LinearAlgebra/Double/SparseVector.cs | 152 ++---------------- .../Double/SparseVectorTest.cs | 5 +- 11 files changed, 202 insertions(+), 176 deletions(-) create mode 100644 src/FSharp/SparseMatrix.fs create mode 100644 src/FSharp/SparseVector.fs diff --git a/src/FSharp/DenseMatrix.fs b/src/FSharp/DenseMatrix.fs index b7638275..a7d81416 100644 --- a/src/FSharp/DenseMatrix.fs +++ b/src/FSharp/DenseMatrix.fs @@ -1,6 +1,9 @@ // // Math.NET Numerics, part of the Math.NET Project // http://mathnet.opensourcedotnet.info +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com // // Copyright (c) 2009 Math.NET // @@ -40,13 +43,6 @@ module DenseMatrix = for j=0 to m-1 do A.[i,j] <- f i j A - - /// Create the identity matrix. - let inline identity (n:int) = - let A = new DenseMatrix(n, n, 0.0) - for i=0 to n-1 do - A.[i,i] <- 1.0 - A /// Create a matrix from a list of float lists. Every list in the master list specifies a row. let inline ofList (fll: float list list) = @@ -95,17 +91,11 @@ module DenseMatrix = /// Initialize a matrix by calling a construction function for every row. let inline initRow (n: int) (m: int) (f: int -> #Vector) = let A = new DenseMatrix(n,m) - for i=0 to n-1 do - let row = f i - if row.Count <> m then failwith "Row generator does not create rows of the appropriate size." - for j=0 to m-1 do A.[i, j] <- row.[j] + for i=0 to n-1 do A.SetRow(i, f i) A /// Initialize a matrix by calling a construction function for every column. let inline initCol (n: int) (m: int) (f: int -> #Vector) = let A = new DenseMatrix(n,m) - for i=0 to m-1 do - let col = f i - if col.Count <> n then failwith "Column generator does not create columns of the appropriate size." - for j=0 to n-1 do A.[j, i] <- col.[j] + for i=0 to m-1 do A.SetColumn(i, f i) A \ No newline at end of file diff --git a/src/FSharp/DenseVector.fs b/src/FSharp/DenseVector.fs index ddf594da..3e56484e 100644 --- a/src/FSharp/DenseVector.fs +++ b/src/FSharp/DenseVector.fs @@ -1,6 +1,8 @@ // // Math.NET Numerics, part of the Math.NET Project -// http://mathnet.opensourcedotnet.info +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com // // Copyright (c) 2009 Math.NET // diff --git a/src/FSharp/FSharp.fsproj b/src/FSharp/FSharp.fsproj index b1df4e62..0da44018 100644 --- a/src/FSharp/FSharp.fsproj +++ b/src/FSharp/FSharp.fsproj @@ -52,9 +52,11 @@ + + - + diff --git a/src/FSharp/Main.fs b/src/FSharp/Main.fs index ede4b47d..8956b227 100644 --- a/src/FSharp/Main.fs +++ b/src/FSharp/Main.fs @@ -1,6 +1,8 @@ // // Math.NET Numerics, part of the Math.NET Project -// http://mathnet.opensourcedotnet.info +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com // // Copyright (c) 2009 Math.NET // diff --git a/src/FSharp/Matrix.fs b/src/FSharp/Matrix.fs index 977aafa2..472c4b42 100644 --- a/src/FSharp/Matrix.fs +++ b/src/FSharp/Matrix.fs @@ -1,6 +1,8 @@ // // Math.NET Numerics, part of the Math.NET Project -// http://mathnet.opensourcedotnet.info +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com // // Copyright (c) 2009 Math.NET // diff --git a/src/FSharp/SparseMatrix.fs b/src/FSharp/SparseMatrix.fs new file mode 100644 index 00000000..d43d64eb --- /dev/null +++ b/src/FSharp/SparseMatrix.fs @@ -0,0 +1,75 @@ +// +// 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 Math.NET +// +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double + +open MathNet.Numerics.LinearAlgebra + +/// A module which implements functional sparse vector operations. +module SparseMatrix = + + /// Create a matrix from a list of float lists. Every list in the master list specifies a row. + let inline ofList (rows: int) (cols: int) (fll: list) = + let A = new Double.SparseMatrix(rows, cols) + fll |> List.iter (fun (i, j, x) -> A.[i,j] <- x) + A + + /// Create a matrix from a list of sequences. Every sequence in the master sequence specifies a row. + let inline ofSeq (rows: int) (cols: int) (fss: #seq) = + let A = new Double.SparseMatrix(rows, cols) + fss |> Seq.iter (fun (i, j, x) -> A.[i,j] <- x) + A + + /// Create a square matrix with constant diagonal entries. + let inline constDiag (n: int) (f: float) = + let A = new Double.SparseMatrix(n,n) + for i=0 to n-1 do + A.[i,i] <- f + A + + /// Create a square matrix with the vector elements on the diagonal. + let inline diag (v: #Vector) = + let n = v.Count + let A = new Double.SparseMatrix(n,n) + for i=0 to n-1 do + A.[i,i] <- v.Item(i) + A + + /// Initialize a matrix by calling a construction function for every row. + let inline initRow (n: int) (m: int) (f: int -> #Vector) = + let A = new Double.SparseMatrix(n,m) + for i=0 to n-1 do A.SetRow(i, f i) + A + + /// Initialize a matrix by calling a construction function for every column. + let inline initCol (n: int) (m: int) (f: int -> #Vector) = + let A = new Double.SparseMatrix(n,m) + for i=0 to n-1 do A.SetColumn(i, f i) + A \ No newline at end of file diff --git a/src/FSharp/SparseVector.fs b/src/FSharp/SparseVector.fs new file mode 100644 index 00000000..e5c748bc --- /dev/null +++ b/src/FSharp/SparseVector.fs @@ -0,0 +1,48 @@ +// +// 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 Math.NET +// +// Permission is hereby granted, free of charge, to any person +// obtaining a copy of this software and associated documentation +// files (the "Software"), to deal in the Software without +// restriction, including without limitation the rights to use, +// copy, modify, merge, publish, distribute, sublicense, and/or sell +// copies of the Software, and to permit persons to whom the +// Software is furnished to do so, subject to the following +// conditions: +// +// The above copyright notice and this permission notice shall be +// included in all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, +// EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES +// OF MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND +// NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT +// HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, +// WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR +// OTHER DEALINGS IN THE SOFTWARE. +// + +namespace MathNet.Numerics.LinearAlgebra.Double + +open MathNet.Numerics.LinearAlgebra + +/// A module which implements functional sparse vector operations. +module SparseVector = + + /// Create a sparse vector with a given dimension from a list of entry, value pairs. + let inline ofList (dim: int) (fl: list) = + let v = new Double.SparseVector(dim) + fl |> List.iter (fun (i, f) -> v.[i] <- f) + v + + /// Create a sparse vector with a given dimension from a sequence of entry, value pairs. + let inline ofSeq (dim: int) (fs: #seq) = + let v = new Double.SparseVector(dim) + fs |> Seq.iter (fun (i, f) -> v.[i] <- f) + v \ No newline at end of file diff --git a/src/FSharp/Vector.fs b/src/FSharp/Vector.fs index 0955eb3f..fd042361 100644 --- a/src/FSharp/Vector.fs +++ b/src/FSharp/Vector.fs @@ -1,6 +1,8 @@ // // Math.NET Numerics, part of the Math.NET Project -// http://mathnet.opensourcedotnet.info +// http://numerics.mathdotnet.com +// http://github.com/mathnet/mathnet-numerics +// http://mathnetnumerics.codeplex.com // // Copyright (c) 2009 Math.NET // @@ -56,12 +58,10 @@ module Vector = () /// In-place vector addition. - let inline addInPlace (v: #Vector) (w: #Vector) = - v.Add w + let inline addInPlace (v: #Vector) (w: #Vector) = v.Add(w, v) /// In place vector subtraction. - let inline subInPlace (v: #Vector) (w: #Vector) = - v.Subtract w + let inline subInPlace (v: #Vector) (w: #Vector) = v.Subtract(w, v) /// Functional map operator for vectors. /// diff --git a/src/FSharpUnitTests/Program.fs b/src/FSharpUnitTests/Program.fs index b27567df..a653e975 100644 --- a/src/FSharpUnitTests/Program.fs +++ b/src/FSharpUnitTests/Program.fs @@ -25,6 +25,20 @@ let DenseVectorTests = ] +/// Unit tests for the sparse vector type. +let SparseVectorTests = + + /// A small uniform vector. + let smallv = new DenseVector( [|0.0;0.3;0.0;0.0;0.0|] ) :> Vector + + specs "SparseVector" [ + spec "SparseVector.ofList" + ((SparseVector.ofList 5 [ (1,0.3) ] :> Vector) |> should equal smallv) + spec "SparseVector.ofSeq" + ((SparseVector.ofSeq 5 (List.toSeq [ (1,0.3) ]) :> Vector) |> should equal smallv) + ] + + /// Unit tests for the vector type. let VectorTests = @@ -106,7 +120,7 @@ let MatrixTests = (Matrix.foralli (fun i j x -> x = float i * 100.0 + float j) largeM |> should equal true) spec "Matrix.existsi" (Matrix.existsi (fun i j x -> x = float i * 100.0 + float j) largeM |> should equal true) - (*spec "Matrix.map" + spec "Matrix.map" (Matrix.map (fun x -> 2.0 * x) smallM |> should equal (2.0 * smallM)) spec "Matrix.mapi" (Matrix.mapi (fun i j x -> float i * 100.0 + float j + x) largeM |> should equal (2.0 * largeM)) @@ -117,7 +131,7 @@ let MatrixTests = spec "Matrix.inplaceMapi" ( let N = largeM.Clone() Matrix.inplaceMapi (fun i j x -> 2.0 * (float i * 100.0 + float j) + x) N - N |> should equal (3.0 * largeM))*) + N |> should equal (3.0 * largeM)) spec "Matrix.nonZeroEntries" (Seq.length (Matrix.nonZeroEntries smallM) |> should equal 4) spec "Matrix.sum" @@ -145,8 +159,6 @@ let DenseMatrixTests = specs "DenseMatrix" [ spec "DenseMatrix.init" (DenseMatrix.init 100 100 (fun i j -> float i * 100.0 + float j) |> should equal largeM) - spec "DenseMatrix.identity" - (DenseMatrix.identity 10 |> should equal (DenseMatrix.init 10 10 (fun i j -> if i = j then 1.0 else 0.0))) spec "DenseMatrix.ofList" (DenseMatrix.ofList [[0.3;0.3];[0.3;0.3]] |> should equal smallM) spec "DenseMatrix.ofSeq" @@ -156,18 +168,32 @@ let DenseMatrixTests = spec "DenseMatrix.initDense" (DenseMatrix.initDense 100 100 (seq { for i in 0 .. 99 do for j in 0 .. 99 -> (i,j, float i * 100.0 + float j)}) |> should equal largeM) - (*spec "DenseMatrix.constDiag" - (DenseMatrix.constDiag 100 2.0 |> should equal (2.0 * (DenseMatrix.identity 100))) + spec "DenseMatrix.constDiag" + (DenseMatrix.constDiag 100 2.0 |> should equal (2.0 * (DenseMatrix.Identity 100))) spec "DenseMatrix.diag" - (DenseMatrix.diag (new DenseVector(100, 2.0)) |> should equal (2.0 * (DenseMatrix.identity 100))) + (DenseMatrix.diag (new DenseVector(100, 2.0)) |> should equal (2.0 * (DenseMatrix.Identity 100))) spec "DenseMatrix.init_row" - (DenseMatrix.init_row 100 100 (fun i -> (DenseVector.init 100 (fun j -> float i * 100.0 + float j))) |> should equal largeM) + (DenseMatrix.initRow 100 100 (fun i -> (DenseVector.init 100 (fun j -> float i * 100.0 + float j))) |> should equal largeM) spec "DenseMatrix.init_col" - (DenseMatrix.init_col 100 100 (fun j -> (DenseVector.init 100 (fun i -> float i * 100.0 + float j))) |> should equal largeM) - spec "DenseMatrix.of_rowvector" - (DenseMatrix.of_rowvector (new DenseVector(10,3.0)) |> should equal ((new DenseMatrix(1,10,3.0)))) - spec "DenseMatrix.of_vector" - (DenseMatrix.of_vector (new DenseVector(10,3.0)) |> should equal ((new DenseMatrix(10,1,3.0))))*) + (DenseMatrix.initCol 100 100 (fun j -> (DenseVector.init 100 (fun i -> float i * 100.0 + float j))) |> should equal largeM) + ] + + +/// Unit tests for the sparse matrix type. +let SparseMatrixTests = + + /// A small uniform vector. + let smallM = DenseMatrix.init 4 4 (fun i j -> if i = 1 && j = 2 then 1.0 else 0.0) :> Matrix + + specs "SparseMatrix" [ + spec "SparseMatrix.ofList" + ((SparseMatrix.ofList 4 4 [(1,2,1.0)] :> Matrix) |> should equal smallM) + spec "SparseMatrix.ofSeq" + ((SparseMatrix.ofSeq 4 4 (Seq.ofList [(1,2,1.0)]) :> Matrix) |> should equal smallM) + spec "SparseMatrix.constDiag" + (SparseMatrix.constDiag 100 2.0 |> should equal (2.0 * (SparseMatrix.Identity 100))) + spec "SparseMatrix.diag" + (SparseMatrix.diag (new DenseVector(100, 2.0)) |> should equal (2.0 * (SparseMatrix.Identity 100))) ] diff --git a/src/Numerics/LinearAlgebra/Double/SparseVector.cs b/src/Numerics/LinearAlgebra/Double/SparseVector.cs index 7621ec29..228929be 100644 --- a/src/Numerics/LinearAlgebra/Double/SparseVector.cs +++ b/src/Numerics/LinearAlgebra/Double/SparseVector.cs @@ -381,8 +381,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentException(Resources.ArgumentVectorsSameLength, "result"); } - CopyTo(result); - result.Add(scalar); + base.Add(scalar, result); } /// @@ -582,9 +581,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide"); } - var ret = leftSide.Clone(); - ret.Add(rightSide); - return ret; + return leftSide.Add(rightSide); } /// @@ -626,8 +623,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentException(Resources.ArgumentVectorsSameLength, "result"); } - CopyTo(result); - result.Subtract(scalar); + base.Subtract(scalar, result); } /// @@ -757,9 +753,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentException(Resources.ArgumentVectorsSameLength, "rightSide"); } - var ret = leftSide.Clone(); - ret.Subtract(rightSide); - return ret; + return leftSide.Subtract(rightSide); } /// @@ -856,9 +850,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentNullException("leftSide"); } - var ret = (SparseVector)leftSide.Clone(); - ret.Multiply(rightSide); - return ret; + return (SparseVector)leftSide.Multiply(rightSide); } /// @@ -875,9 +867,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentNullException("rightSide"); } - var ret = (SparseVector)rightSide.Clone(); - ret.Multiply(leftSide); - return ret; + return (SparseVector)rightSide.Multiply(leftSide); } /// @@ -922,9 +912,7 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentNullException("leftSide"); } - var ret = (SparseVector)leftSide.Clone(); - ret.Multiply(1.0 / rightSide); - return ret; + return (SparseVector) leftSide.Multiply(1.0 / rightSide); } /// @@ -1114,131 +1102,19 @@ namespace MathNet.Numerics.LinearAlgebra.Double throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other"); } - // We cannot iterate using NonZeroCount because the value may be changed (if multiply by 0) - var copy = (SparseVector)this.Clone(); - for (var i = 0; i < Count; i++) - { - copy[i] *= other[i]; - } - return copy; - } - - /// - /// 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. - /// If the other vector is . - /// If the result vector is . - /// If this vector and are not the same size. - /// If this vector and are not the same size. - public override void PointwiseMultiply(Vector other, Vector result) - { - if (result == null) - { - throw new ArgumentNullException("result"); - } - - if (other == null) - { - throw new ArgumentNullException("other"); - } - - if (Count != other.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other"); - } - - if (Count != result.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "result"); - } - - if (ReferenceEquals(this, result) || ReferenceEquals(other, result)) + var copy = new SparseVector(Count); + for (int i = 0; i < this._nonZeroIndices.Length; i++) { - var tmp = result.CreateVector(result.Count); - PointwiseMultiply(other, tmp); - tmp.CopyTo(result); - } - else - { - CopyTo(result); - result.PointwiseMultiply(other); - } - } - - /// - /// Pointwise divide this vector with another vector. - /// - /// The vector to pointwise divide this one by. - /// A new vector which is the pointwise division of the two vectors. - /// If the other vector is . - /// If this vector and are not the same size. - public override Vector PointwiseDivide(Vector other) - { - if (other == null) - { - throw new ArgumentNullException("other"); - } - - if (Count != other.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other"); + var d = this._nonZeroValues[i] * other[this._nonZeroIndices[i]]; + if (d != 0.0) + { + copy[this._nonZeroIndices[i]] = d; + } } - // base implementation iterates though all elements, but we need only take non-zeros - var copy = (SparseVector)this.Clone(); - for (var i = 0; i < NonZerosCount; i++) - { - copy[_nonZeroIndices[i]] /= other[_nonZeroIndices[i]]; - } return copy; } - /// - /// Pointwise divide this vector with another vector and stores the result into the result vector. - /// - /// The vector to pointwise divide this one by. - /// The vector to store the result of the pointwise division. - /// If the other vector is . - /// If the result vector is . - /// If this vector and are not the same size. - /// If this vector and are not the same size. - public override void PointwiseDivide(Vector other, Vector result) - { - if (result == null) - { - throw new ArgumentNullException("result"); - } - - if (other == null) - { - throw new ArgumentNullException("other"); - } - - if (Count != other.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "other"); - } - - if (Count != result.Count) - { - throw new ArgumentException(Resources.ArgumentVectorsSameLength, "result"); - } - - if (ReferenceEquals(this, result) || ReferenceEquals(other, result)) - { - var tmp = result.CreateVector(result.Count); - PointwiseDivide(other, tmp); - tmp.CopyTo(result); - } - else - { - CopyTo(result); - result.PointwiseDivide(other); - } - } - /// /// Outer product of two vectors /// diff --git a/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs b/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs index 28d12226..c08baaaa 100644 --- a/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs +++ b/src/UnitTests/LinearAlgebraTests/Double/SparseVectorTest.cs @@ -152,7 +152,10 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double { var vector = new SparseVector(_data); var other = +vector; - Assert.AreSame(vector, other, "Should be the same vector"); + for (var i = 0; i < _data.Length; i++) + { + Assert.AreEqual(vector[i], other[i]); + } } [Test]