Sparse Linear Algebra

ADCME augments TensorFlow APIs by adding sparse linear algebra support. In ADCME, sparse matrices are represented by SparseTensor. This data structure stores indices, rows and cols of the sparse matrices and keep track of relevant information such as whether it is diagonal for performance consideration. The default is row major (due to TensorFlow backend).

When evaluating SparseTensor, the output will be SparseMatrixCSC, the native Julia representation of sparse matrices

A = run(sess, s) # A has type SparseMatrixCSC{Float64,Int64}

Sparse Matrix Construction

By passing columns (Int64), rows (Int64) and values (Float64) arrays

ii = [1;2;3;4]
jj = [1;2;3;4]
vv = [1.0;1.0;1.0;1.0]
s = SparseTensor(ii, jj, vv, 4, 4)

By passing a SparseMatrixCSC

using SparseArrays
s = SparseTensor(sprand(10,10,0.3))

By passing a dense array (tensor or numerical array)

D = Array(sprand(10,10,0.3)) # a dense array
d = constant(D)
s = dense_to_sparse(d)

There are also special constructors.

Description	Code
Diagonal matrix with diagonal `v`	`spdiag(v)`
Empty matrix with size `m`, `n`	`spzero(m, n)`
Identity matrix with size `m`	`spdiag(m)`

Matrix Traits

Size of the matrices
```
size(s) # (10,20)
size(s,1) # 10
```
Return row, col, val arrays (also known as COO arrays)
```
ii,jj,vv = find(s)
```

Arithmetic Operations

Add Subtract
```
s = s1 + s2
s = s1 - s2
```
Scalar Product
```
s = 2.0 * s1
s = s1 / 2.0
```
Sparse Product
```
s = s1 * s2
```
Transposition
```
s = s1'
```

Sparse Solvers

Solve a linear system (s is a square matrix)
```
sol = s\rhs
```
Solve a least square system (s is a tall matrix)
```
sol = s\rhs
```

Note

The least square solvers are implemented using Eigen sparse linear packages, and the gradients are also implemented. Thus, the following codes will work as expected (the gradients functions will correctly compute the gradients):

ii = [1;2;3;4]
jj = [1;2;3;4]
vv = constant([1.0;1.0;1.0;1.0])
rhs = constant(rand(4))
s = SparseTensor(ii, jj, vv, 4, 4)
sol = s\rhs
run(sess, sol)
run(sess, gradients(sum(sol), rhs))
run(sess, gradients(sum(sol), vv))

Assembling Sparse Matrix

In many applications, we want to accumulate row, col and val to assemble a sparse matrix in iterations. For this purpose, we provide the SparseAssembler utilities.

ADCME.SparseAssembler — Function

SparseAssembler(handle::Union{PyObject, <:Integer}, n::Union{PyObject, <:Integer}, tol::Union{PyObject, <:Real}=0.0)

Creates a SparseAssembler for accumulating row, col, val for sparse matrices.

handle: an integer handle for creating a sparse matrix. If the handle already exists, SparseAssembler return the existing sparse matrix handle. If you are creating different sparse matrices, the handles should be different.
n: Number of rows of the sparse matrix.
tol (optional): Tolerance. SparseAssembler will treats any values less than tol as zero.

Example 1

handle = SparseAssembler(100, 5, 1e-8)
op1 = accumulate(handle, 1, [1;2;3], [1.0;2.0;3.0])
op2 = accumulate(handle, 2, [1;2;3], [1.0;2.0;3.0])
J = assemble(5, 5, [op1;op2])

J will be a SparseTensor object.

Example 2

handle = SparseAssembler(0, 5)
op1 = accumulate(handle, 1, [1;2;3], ones(3))
op2 = accumulate(handle, 1, [3], [1.])
op3 = accumulate(handle, 2, [1;3], ones(2))
J = assemble(5, 5, [op1;op2;op3]) # op1, op2, op3 are parallel
Array(run(sess, J))≈[1.0  1.0  2.0  0.0  0.0
                1.0  0.0  1.0  0.0  0.0
                0.0  0.0  0.0  0.0  0.0
                0.0  0.0  0.0  0.0  0.0
                0.0  0.0  0.0  0.0  0.0]

source

Base.accumulate — Function

accumulate(handle::PyObject, row::Union{PyObject, <:Integer}, cols::Union{PyObject, Array{<:Integer}}, vals::Union{PyObject, Array{<:Real}})

Accumulates row-th row. It adds the value to the sparse matrix

for k = 1:length(cols)
    A[row, cols[k]] += vals[k]
end

handle is the handle created by SparseAssembler.

See SparseAssembler for an example.

Note

The function accumulate returns a op::PyObject. Only when op is executed, the nonzero values are populated into the sparse matrix.

source

ADCME.assemble — Function

assemble(m::Union{PyObject, <:Integer}, n::Union{PyObject, <:Integer}, ops::PyObject)

Assembles the sparse matrix from the ops created by accumulate. ops is either a single output from accumulate, or concated from several ops

op1 = accumulate(handle, 1, [1;2;3], [1.0;2.0;3.0])
op2 = accumulate(handle, 2, [1;2;3], [1.0;2.0;3.0])
op = [op1;op2] # equivalent to `vcat([op1, op2]...)`

m and n are rows and columns of the sparse matrix.

See SparseAssembler for an example.

source