diff --git a/src/FSharp/LinearAlgebra.Double.Matrix.fs b/src/FSharp/LinearAlgebra.Double.Matrix.fs index 2d8d565d..acc236ba 100644 --- a/src/FSharp/LinearAlgebra.Double.Matrix.fs +++ b/src/FSharp/LinearAlgebra.Double.Matrix.fs @@ -36,6 +36,82 @@ open MathNet.Numerics.LinearAlgebra.Generic [] module Matrix = + /// Transform a vector into a 2D array. + let inline toArray2 (A: #Matrix) = A.ToArray() + + /// In-place map of every matrix element using a function. + let inline mapInPlace (f: float -> float) (A: #Matrix) = + A.MapInplace((fun x -> f x), true) + + /// In-place map of every matrix element using a position dependent function. + let inline mapiInPlace (f: int -> int -> float -> float) (A: #Matrix) = + A.MapIndexedInplace((fun i j x -> f i j x), true) + + /// In-place map of every matrix element using a function. + /// Zero-values may be skipped (relevant mostly for sparse matrices). + let inline mapnzInPlace (f: float -> float) (A: #Matrix) = + A.MapInplace((fun x -> f x), false) + + /// In-place map of every matrix element using a position dependent function. + /// Zero-values may be skipped (relevant mostly for sparse matrices). + let inline mapinzInPlace (f: int -> int -> float -> float) (A: #Matrix) = + A.MapIndexedInplace((fun i j x -> f i j x), false) + + /// In-place map every matrix column using the given position dependent function. + let inline mapColsInPlace (f: int -> Vector -> Vector) (A: #Matrix) = + for j = 0 to A.ColumnCount-1 do + A.SetColumn(j, f j (A.Column(j))) + + /// In-place map every matrix row using the given position dependent function. + let inline mapRowsInPlace (f: int -> Vector -> Vector) (A: #Matrix) = + for i = 0 to A.RowCount-1 do + A.SetRow(i, f i (A.Row(i))) + + [] + let inplaceMapi = mapiInPlace + [] + let inplaceMapCols = mapColsInPlace + [] + let inplaceMapRows = mapRowsInPlace + + /// Map every matrix element using the given function. + let inline map (f: float -> float) (A: #Matrix) = + let A = A.Clone() + A.MapInplace((fun x -> f x), true) + A + + /// Map every matrix element using the given function. + /// Zero-values may be skipped (relevant mostly for sparse matrices). + let inline mapnz (f: float -> float) (A: #Matrix) = + let A = A.Clone() + A.MapInplace((fun x -> f x), false) + A + + /// Map every matrix element using the given position dependent function. + let inline mapi (f: int -> int -> float -> float) (A: #Matrix) = + let A = A.Clone() + A.MapIndexedInplace((fun i j x -> f i j x), true) + A + + /// Map every matrix element using the given position dependent function. + /// Zero-values may be skipped (relevant mostly for sparse matrices). + let inline mapinz (f: int -> int -> float -> float) (A: #Matrix) = + let A = A.Clone() + A.MapIndexedInplace((fun i j x -> f i j x), false) + A + + /// Map every matrix column using the given position dependent function. + let inline mapCols (f: int -> Vector -> Vector) (A: #Matrix) = + let A = A.Clone() + mapColsInPlace f A + A + + /// Map every matrix row using the given position dependent function. + let inline mapRows (f: int -> Vector -> Vector) (A: #Matrix) = + let A = A.Clone() + mapRowsInPlace f A + A + /// Fold a function over all matrix elements. let inline fold (f: 'a -> float -> 'a) (acc0: 'a) (A: #Matrix) = let n = A.RowCount @@ -66,12 +142,6 @@ module Matrix = acc <- f i j acc (A.At(i,j)) acc - /// Create a 2D array from a matrix. - let inline toArray2 (A: #Matrix) = - let n = A.RowCount - let m = A.ColumnCount - Array2D.init n m (fun i j -> (A.Item(i,j))) - /// Checks whether a predicate holds for all elements of a matrix. let inline forall (p: float -> bool) (A: #Matrix) = let mutable b = true @@ -116,61 +186,9 @@ module Matrix = if j = A.ColumnCount then i <- i+1; j <- 0 b - /// Map every matrix element using the given function. - let inline map (f: float -> float) (A: #Matrix) = - let N = A.RowCount - let M = A.ColumnCount - let C = A.Clone() - for i=0 to N-1 do - for j=0 to M-1 do - C.At(i, j, f (C.At(i,j))) - C - - /// Map every matrix element using the given position dependent function. - let inline mapi (f: int -> int -> float -> float) (A: #Matrix) = - let N = A.RowCount - let M = A.ColumnCount - let C = A.Clone() - for i=0 to N-1 do - for j=0 to M-1 do - C.At(i, j, f i j (C.At(i,j))) - C - - /// In-place map every matrix column using the given position dependent function. - let inline inplaceMapCols (f: int -> Vector -> Vector) (A: #Matrix) = - for j = 0 to A.ColumnCount-1 do - A.SetColumn(j, f j (A.Column(j))) - () - - /// In-place map every matrix row using the given position dependent function. - let inline inplaceMapRows (f: int -> Vector -> Vector) (A: #Matrix) = - for i = 0 to A.RowCount-1 do - A.SetRow(i, f i (A.Row(i))) - () - - /// Map every matrix column using the given position dependent function. - let inline mapCols (f: int -> Vector -> Vector) (A: #Matrix) = - let A = A.Clone() - inplaceMapCols f A - A - - /// Map every matrix row using the given position dependent function. - let inline mapRows (f: int -> Vector -> Vector) (A: #Matrix) = - let A = A.Clone() - inplaceMapRows f A - A - /// In-place assignment. let inline inplaceAssign (f: int -> int -> float) (A: #Matrix) = - for i=0 to A.RowCount-1 do - for j=0 to A.ColumnCount-1 do - A.At(i, j, f i j) - - /// In-place map of every matrix element using a position dependent function. - let inline inplaceMapi (f: int -> int -> float -> float) (A: #Matrix) = - for i=0 to A.RowCount-1 do - for j=0 to A.ColumnCount-1 do - A.At(i, j, f i j (A.At(i,j))) + A.MapIndexedInplace((fun i j x -> f i j), true) /// Creates a sequence that iterates the non-zero entries in the matrix. let inline nonZeroEntries (A: #Matrix) = @@ -199,14 +217,12 @@ module Matrix = for i=0 to A.RowCount-1 do for j=0 to A.ColumnCount-1 do f (A.At(i,j)) - () /// Iterates over all elements of a matrix using the element indices. let inline iteri (f: int -> int -> float -> unit) (A: #Matrix) = for i=0 to A.RowCount-1 do for j=0 to A.ColumnCount-1 do f i j (A.At(i,j)) - () /// Fold one column. let inline foldCol (f: 'a -> float -> 'a) acc (A: #Matrix) k = @@ -298,10 +314,7 @@ module DenseMatrix = /// Create a matrix with the given entries. [] - let inline initDense (rows: int) (cols: int) (es: #seq) = - let A = new DenseMatrix(rows,cols) - Seq.iter (fun (i,j,f) -> A.At(i,j,f)) es - A + let inline initDense (rows: int) (cols: int) (es: #seq) = ofSeqi rows cols es /// Create a square matrix with constant diagonal entries. let inline constDiag (n: int) (f: float) = diff --git a/src/FSharp/LinearAlgebra.Double.Vector.fs b/src/FSharp/LinearAlgebra.Double.Vector.fs index 7857b3a5..a314ef58 100644 --- a/src/FSharp/LinearAlgebra.Double.Vector.fs +++ b/src/FSharp/LinearAlgebra.Double.Vector.fs @@ -45,24 +45,20 @@ module Vector = /// In-place mutation by applying a function to every element of the vector. let inline mapInPlace (f: float -> float) (v: #Vector) = v.MapInplace((fun x -> f x), true) - () /// In-place mutation by applying a function to every element of the vector. let inline mapiInPlace (f: int -> float -> float) (v: #Vector) = v.MapIndexedInplace((fun i x -> f i x), true) - () /// In-place mutation by applying a function to every element of the vector. /// Zero-values may be skipped (relevant mostly for sparse vectors). let inline mapnzInPlace (f: float -> float) (v: #Vector) = v.MapInplace((fun x -> f x), false) - () /// In-place mutation by applying a function to every element of the vector. /// Zero-values may be skipped (relevant mostly for sparse vectors). let inline mapinzInPlace (f: int -> float -> float) (v: #Vector) = v.MapIndexedInplace((fun i x -> f i x), false) - () /// Maps a vector to a new vector by applying a function to every element. let inline map f (v: #Vector) = diff --git a/src/FSharpUnitTests/DenseVectorTests.fs b/src/FSharpUnitTests/DenseVectorTests.fs index 8c8d0e7b..b99a52f2 100644 --- a/src/FSharpUnitTests/DenseVectorTests.fs +++ b/src/FSharpUnitTests/DenseVectorTests.fs @@ -4,6 +4,8 @@ open NUnit.Framework open FsUnit open MathNet.Numerics.LinearAlgebra.Generic open MathNet.Numerics.LinearAlgebra.Double +open MathNet.Numerics.Distributions +open MathNet.Numerics.Statistics /// Unit tests for the dense vector type. module DenseVectorTests = @@ -14,8 +16,23 @@ module DenseVectorTests = /// A large vector with increasingly large entries let largev = new DenseVector( Array.init 100 (fun i -> float i / 100.0) ) + [] + let ``DenseVector.zeroCreate`` () = + (DenseVector.zeroCreate 100) + largev |> should equal largev + + [] + let ``DenseVector.randomCreate`` () = + let m = DenseVector.randomCreate 100 (Normal.WithMeanStdDev(100.0,0.1)) + m.Values |> ArrayStatistics.Mean |> should (equalWithin 10.0) 100.0 + m.Count |> should equal 100 + + [] + let ``DenseVector.create`` () = + DenseVector.create 5 0.3 |> should equal smallv + [] let ``DenseVector.init`` () = + DenseVector.init 5 (fun i -> 0.3) |> should equal smallv DenseVector.init 100 (fun i -> float i / 100.0) |> should equal largev [] diff --git a/src/FSharpUnitTests/MatrixTests.fs b/src/FSharpUnitTests/MatrixTests.fs index c9b32d77..3a0c8c72 100644 --- a/src/FSharpUnitTests/MatrixTests.fs +++ b/src/FSharpUnitTests/MatrixTests.fs @@ -12,6 +12,9 @@ module MatrixTests = let smallM = DenseMatrix.OfArray( Array2D.create 2 2 0.3 ) let failingFoldBackM = DenseMatrix.init 2 3 (fun i j -> 1.0) + /// A small sparse matrix. + let sparseM = SparseMatrix.ofListi 2 3 [(1,0,0.3)] + /// A large vector with increasingly large entries let largeM = DenseMatrix.OfArray( Array2D.init 100 100 (fun i j -> float i * 100.0 + float j) ) @@ -41,6 +44,70 @@ module MatrixTests = m.[..1,2..] <- DenseMatrix(2,2,[|5.;6.;7.;8.|]); m |> should equal (DenseMatrix(4,4,[|0.;100.;200.;300.;1.;101.;201.;301.;5.;6.;202.;302.;7.;8.;203.;303.|])) + [] + let ``Matrix.toArray2`` () = + Matrix.toArray2 smallM |> should array2_equal (Array2D.create 2 2 0.3) + + [] + let ``Matrix.mapInPlace.Dense`` () = + let M = largeM.Clone() + M |> Matrix.mapInPlace (fun x -> 3.0 * x) + M |> should equal (3.0 * largeM) + + [] + let ``Matrix.mapInPlace.Sparse`` () = + let M = sparseM.Clone() + M |> Matrix.mapInPlace (fun x -> 3.0 * x) + M |> should equal (3.0 * sparseM) + + [] + let ``Matrix.mapnzInPlace.Sparse`` () = + let M = sparseM.Clone() + M |> Matrix.mapnzInPlace (fun x -> 3.0 * x) + M |> should equal (3.0 * sparseM) + + [] + let ``Matrix.mapiInPlace.Dense`` () = + let M = largeM.Clone() + M |> Matrix.mapiInPlace (fun i j x -> 2.0 * (float i * 100.0 + float j) + x) + M |> should equal (3.0 * largeM) + + [] + let ``Matrix.mapiInPlace.Sparse`` () = + let M = sparseM.Clone() + M |> Matrix.mapiInPlace (fun i j x -> if i=j then 2.0*x+1.0 else 2.0*x) + M |> should equal (2.0 * sparseM + SparseMatrix.init 2 3 (fun i j -> if i=j then 1.0 else 0.0)) + + [] + let ``Matrix.mapinzInPlace.Sparse`` () = + let M = sparseM.Clone() + M |> Matrix.mapinzInPlace (fun i j x -> 2.0*x) + M |> should equal (2.0 * sparseM) + + [] + let ``Matrix.map`` () = + Matrix.map (fun x -> 2.0 * x) smallM |> should equal (2.0 * smallM) + + [] + let ``Matrix.mapnz`` () = + Matrix.mapnz (fun x -> 2.0 * x) smallM |> should equal (2.0 * smallM) + + [] + let ``Matrix.mapi`` () = + Matrix.mapi (fun i j x -> float i * 100.0 + float j + x) largeM |> should equal (2.0 * largeM) + + [] + let ``Matrix.mapinz`` () = + Matrix.mapinz (fun i j x -> float i * 100.0 + float j + x) largeM |> should equal (2.0 * largeM) + + [] + let ``Matrix.mapCols`` () = + Matrix.mapCols (fun j col -> col.Add(float j)) smallM |> should (approximately_matrix_equal 14) (matrix [[0.3;1.3];[0.3;1.3]]) + + [] + let ``Matrix.mapRows`` () = + Matrix.mapRows (fun i row -> row.Add(float i)) smallM |> should (approximately_matrix_equal 14) (matrix [[0.3;0.3];[1.3;1.3]]) + [] let ``Matrix.fold`` () = Matrix.fold (fun a b -> a - b) 0.0 smallM |> should equal -1.2 @@ -57,10 +124,6 @@ module MatrixTests = let ``Matrix.foldi`` () = Matrix.foldi (fun i j acc x -> acc + x + float (i+j)) 0.0 smallM |> should equal 5.2 - [] - let ``Matrix.toArray2`` () = - Matrix.toArray2 smallM |> should array2_equal (Array2D.create 2 2 0.3) - [] let ``Matrix.forall`` () = Matrix.forall (fun x -> x = 0.3) smallM |> should equal true @@ -77,34 +140,12 @@ module MatrixTests = let ``Matrix.existsi`` () = Matrix.existsi (fun i j x -> x = float i * 100.0 + float j) largeM |> should equal true - [] - let ``Matrix.map`` () = - Matrix.map (fun x -> 2.0 * x) smallM |> should equal (2.0 * smallM) - - [] - let ``Matrix.mapi`` () = - Matrix.mapi (fun i j x -> float i * 100.0 + float j + x) largeM |> should equal (2.0 * largeM) - - [] - let ``Matrix.mapCols`` () = - Matrix.mapCols (fun j col -> col.Add(float j)) smallM |> should (approximately_matrix_equal 14) (matrix [[0.3;1.3];[0.3;1.3]]) - - [] - let ``Matrix.mapRows`` () = - Matrix.mapRows (fun i row -> row.Add(float i)) smallM |> should (approximately_matrix_equal 14) (matrix [[0.3;0.3];[1.3;1.3]]) - [] let ``Matrix.inplaceAssign`` () = let N = smallM.Clone() - Matrix.inplaceAssign (fun i j -> 0.0) N + N |> Matrix.inplaceAssign (fun i j -> 0.0) N |> should equal (0.0 * smallM) - [] - let ``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) - [] let ``Matrix.nonZeroEntries`` () = Seq.length (Matrix.nonZeroEntries smallM) |> should equal 4 diff --git a/src/FSharpUnitTests/VectorTests.fs b/src/FSharpUnitTests/VectorTests.fs index 569c2a6c..b4c50c0a 100644 --- a/src/FSharpUnitTests/VectorTests.fs +++ b/src/FSharpUnitTests/VectorTests.fs @@ -9,13 +9,13 @@ open MathNet.Numerics.LinearAlgebra.Double module VectorTests = /// A small uniform vector. - let smallv = new DenseVector([|0.3;0.3;0.3;0.3;0.3|]) :> Vector + let smallv = DenseVector([|0.3;0.3;0.3;0.3;0.3|]) - /// A small uniform vector. - let sparsev = SparseVector.OfIndexedEnumerable(5, [(1,0.3)]) :> Vector + /// A small sparse vector. + let sparsev = SparseVector.ofListi 5 [(1,0.3)] /// A large vector with increasingly large entries - let largev = new DenseVector(Array.init 100 (fun i -> float i / 100.0)) :> Vector + let largev = DenseVector(Array.init 100 (fun i -> float i / 100.0)) [] let ``Vector.GetSlice`` () = @@ -54,37 +54,37 @@ module VectorTests = [] let ``Vector.mapInPlace.Dense`` () = let w = smallv.Clone() - Vector.mapInPlace (fun x -> 2.0 * x) w + w |> Vector.mapInPlace (fun x -> 2.0 * x) w |> should equal (2.0 * smallv) [] let ``Vector.mapInPlace.Sparse`` () = let w = sparsev.Clone() - Vector.mapInPlace (fun x -> 2.0 * x) w + w |> Vector.mapInPlace (fun x -> 2.0 * x) w |> should equal (2.0 * sparsev) [] let ``Vector.mapnzInPlace.Sparse`` () = let w = sparsev.Clone() - Vector.mapnzInPlace (fun x -> 2.0 * x) w + w |> Vector.mapnzInPlace (fun x -> 2.0 * x) w |> should equal (2.0 * sparsev) [] let ``Vector.mapiInPlace.Dense`` () = let w = largev.Clone() - Vector.mapiInPlace (fun i x -> float i / 100.0) w + w |> Vector.mapiInPlace (fun i x -> float i / 100.0) w |> should equal (largev) - + [] let ``Vector.mapiInPlace.Sparse`` () = let w = sparsev.Clone() - Vector.mapiInPlace (fun i x -> 2.0 * float i * x) w + w |> Vector.mapiInPlace (fun i x -> 2.0 * float i * x) w |> should equal (2.0 * sparsev) [] let ``Vector.mapinzInPlace.Sparse`` () = let w = sparsev.Clone() - Vector.mapinzInPlace (fun i x -> 2.0 * float i * x) w + w |> Vector.mapinzInPlace (fun i x -> 2.0 * float i * x) w |> should equal (2.0 * sparsev) [] diff --git a/src/Numerics/LinearAlgebra/Generic/Matrix.cs b/src/Numerics/LinearAlgebra/Generic/Matrix.cs index f91221b0..322ea3c8 100644 --- a/src/Numerics/LinearAlgebra/Generic/Matrix.cs +++ b/src/Numerics/LinearAlgebra/Generic/Matrix.cs @@ -1438,5 +1438,26 @@ namespace MathNet.Numerics.LinearAlgebra.Generic { return Storage.ToRowMajorArray(); } + + /// + /// Applies a function to each value of this matrix and replaces the value with its result. + /// If forceMapZero is not set to true, zero values may or may not be skipped depending + /// on the actual data storage implementation (relevant mostly for sparse matrices). + /// + public void MapInplace(Func f, bool forceMapZeros = false) + { + Storage.MapInplace(f, forceMapZeros); + } + + /// + /// Applies a function to each value of this matrix and replaces the value with its result. + /// The row and column indices of each value (zero-based) are passed as first arguments to the function. + /// If forceMapZero is not set to true, zero values may or may not be skipped depending + /// on the actual data storage implementation (relevant mostly for sparse matrices). + /// + public void MapIndexedInplace(Func f, bool forceMapZeros = false) + { + Storage.MapIndexedInplace(f, forceMapZeros); + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs index eeec0191..da7f799d 100644 --- a/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/DenseColumnMajorMatrixStorage.cs @@ -364,5 +364,28 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } return ret; } + + // FUNCTIONAL COMBINATORS + + public override void MapInplace(Func f, bool forceMapZeros = false) + { + for (int i = 0; i < Data.Length; i++) + { + Data[i] = f(Data[i]); + } + } + + public override void MapIndexedInplace(Func f, bool forceMapZeros = false) + { + int index = 0; + for (int j = 0; j < ColumnCount; j++) + { + for (int i = 0; i < RowCount; i++) + { + Data[index] = f(i, j, Data[index]); + index++; + } + } + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/DiagonalMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/DiagonalMatrixStorage.cs index 5b5a7393..f2a08e94 100644 --- a/src/Numerics/LinearAlgebra/Storage/DiagonalMatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/DiagonalMatrixStorage.cs @@ -525,5 +525,27 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } return ret; } + + // FUNCTIONAL COMBINATORS + + public override void MapInplace(Func f, bool forceMapZeros = false) + { + // we deliberately ignore forceMapZeros since we would not actually + // support any non-zero results outside of the diagonal anyway + for (int i = 0; i < Data.Length; i++) + { + Data[i] = f(Data[i]); + } + } + + public override void MapIndexedInplace(Func f, bool forceMapZeros = false) + { + // we deliberately ignore forceMapZeros since we would not actually + // support any non-zero results outside of the diagonal anyway + for (int i = 0; i < Data.Length; i++) + { + Data[i] = f(i, i, Data[i]); + } + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/MatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/MatrixStorage.cs index 57eafc8d..ec01602e 100644 --- a/src/Numerics/LinearAlgebra/Storage/MatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/MatrixStorage.cs @@ -411,5 +411,29 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } return ret; } + + // FUNCTIONAL COMBINATORS + + public virtual void MapInplace(Func f, bool forceMapZeros = false) + { + for (int i = 0; i < RowCount; i++) + { + for (int j = 0; j < ColumnCount; j++) + { + At(i, j, f(At(i, j))); + } + } + } + + public virtual void MapIndexedInplace(Func f, bool forceMapZeros = false) + { + for (int i = 0; i < RowCount; i++) + { + for (int j = 0; j < ColumnCount; j++) + { + At(i, j, f(i, j, At(i, j))); + } + } + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs b/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs index 4d74db15..49c17df8 100644 --- a/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/SparseCompressedRowMatrixStorage.cs @@ -456,6 +456,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage rowPointers[row] = index; if (trows[row] != null) { + // TODO: Don't we need to sort here!? foreach (var item in trows[row]) { values.Add(item.Item2); @@ -907,5 +908,107 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } return ret; } + + // FUNCTIONAL COMBINATORS + + public override void MapInplace(Func f, bool forceMapZeros = false) + { + var newRowPointers = new int[RowCount]; + var newColumnIndices = new List(); + var newValues = new List(); + + if (forceMapZeros || !Zero.Equals(f(Zero))) + { + int k = 0; + for (int row = 0; row < RowCount; row++) + { + newRowPointers[row] = newValues.Count; + for (int col = 0; col < ColumnCount; col++) + { + var item = k < (row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount) && (ColumnIndices[k]) == col + ? f(Values[k++]) + : f(Zero); + if (!Zero.Equals(item)) + { + newValues.Add(item); + newColumnIndices.Add(col); + } + } + } + } + else + { + for (int row = 0; row < RowCount; row++) + { + newRowPointers[row] = newValues.Count; + var startIndex = RowPointers[row]; + var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + for (var j = startIndex; j < endIndex; j++) + { + var item = f(Values[j]); + if (!Zero.Equals(item)) + { + newValues.Add(item); + newColumnIndices.Add(ColumnIndices[j]); + } + } + } + } + + ColumnIndices = newColumnIndices.ToArray(); + Values = newValues.ToArray(); + ValueCount = newValues.Count; + Array.Copy(newRowPointers, RowPointers, RowCount); + } + + public override void MapIndexedInplace(Func f, bool forceMapZeros = false) + { + var newRowPointers = new int[RowCount]; + var newColumnIndices = new List(); + var newValues = new List(); + + if (forceMapZeros || !Zero.Equals(f(0,0,Zero))) + { + int k = 0; + for (int row = 0; row < RowCount; row++) + { + newRowPointers[row] = newValues.Count; + for (int col = 0; col < ColumnCount; col++) + { + var item = k < (row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount) && (ColumnIndices[k]) == col + ? f(row, col, Values[k++]) + : f(row, col, Zero); + if (!Zero.Equals(item)) + { + newValues.Add(item); + newColumnIndices.Add(col); + } + } + } + } + else + { + for (int row = 0; row < RowCount; row++) + { + newRowPointers[row] = newValues.Count; + var startIndex = RowPointers[row]; + var endIndex = row < RowPointers.Length - 1 ? RowPointers[row + 1] : ValueCount; + for (var j = startIndex; j < endIndex; j++) + { + var item = f(row, ColumnIndices[j], Values[j]); + if (!Zero.Equals(item)) + { + newValues.Add(item); + newColumnIndices.Add(ColumnIndices[j]); + } + } + } + } + + ColumnIndices = newColumnIndices.ToArray(); + Values = newValues.ToArray(); + ValueCount = newValues.Count; + Array.Copy(newRowPointers, RowPointers, RowCount); + } } } diff --git a/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs b/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs index 47d00ea1..ca6a606c 100644 --- a/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs +++ b/src/Numerics/LinearAlgebra/Storage/SparseVectorStorage.cs @@ -612,7 +612,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } else { - for (int i = 0; i < Values.Length; i++) + for (int i = 0; i < ValueCount; i++) { var item = f(Values[i]); if (!Zero.Equals(item)) @@ -646,7 +646,7 @@ namespace MathNet.Numerics.LinearAlgebra.Storage } else { - for (int i = 0; i < Values.Length; i++) + for (int i = 0; i < ValueCount; i++) { var item = f(Indices[i], Values[i]); if (!Zero.Equals(item))