Merge pull request #5 from timholy/sl/spec-builder

Added sparse_factors function

src/WoodburyMatrices.jl

@@ -98,4 +98,54 @@ function A_ldiv_B!(W::Woodbury, B::AbstractVector)
     B
 end
 
+"""
+Given a type `T`, an integer `n` and `m` tuples of the form (i, j, v), build
+sparse matrices `rows`, `vals`, `cols` such that the product
+`out = rows * vals * cols` is equivalent to:
+
+```julia
+out = zeros(T, n, n)
+
+for (i, j, v) in args
+    out[i, j] = v
+end
+```
+
+The first two components (`i` and `j`) of each tuple should be integers
+whereas the third component should be of type `T`
+
+Example:
+
+
+```
+julia> r, v, c = WoodburyMatrices.sparse_factors(Float64, 3,
+                                                 (1, 1, 2.0),
+                                                 (2, 2, 3.0),
+                                                 (3, 3, 4.0));
+
+julia> r*c*v - Diagonal([2.0, 3.0, 4.0])
+3x3 sparse matrix with 0 Float64 entries:
+```
+
+"""
+function sparse_factors{T}(::Type{T}, n::Int, args::Tuple{Int, Int, Any}...)
+    m = length(args)
+    rows = spzeros(T, n, m)
+    cols = spzeros(T, m, n)
+    val_diag = zeros(m)
+
+    ix = 1
+    for (i, (row, col, val)) in enumerate(args)
+        rows[row, ix] = 1
+        cols[ix, col] = 1
+        val_diag[ix] = val
+        ix += 1
+    end
+
+    vals = spdiagm(val_diag)
+
+    rows, vals, cols
+end
+
+
 end

test/runtests.jl

@@ -96,3 +96,19 @@ Wc = copy(W)
 # Mismatched sizes
 @test_throws DimensionMismatch Woodbury(rand(5,5),rand(5,2),rand(2,3),rand(3,5))
 @test_throws DimensionMismatch Woodbury(rand(5,5),rand(5,2),rand(3,3),rand(2,5))
+
+# spec builder
+n = 10
+num_nz = 3
+for i in 1:5  # repeat 5 times
+    args = [(rand(1:n), rand(1:n), rand()) for j in 1:num_nz]
+    r, v, c = WoodburyMatrices.sparse_factors(Float64, n, args...)
+    out = r*v*c
+
+    by_hand = zeros(n, n)
+    for (i, j, v) in args
+        by_hand[i, j] = v
+    end
+
+    @test maxabs(out - by_hand) == 0.0
+end