sparsegrad.impl.sparse package

Submodules

sparsegrad.impl.sparse.sparse module

This module contains implementation details sparse matrix operations

class sparsegrad.impl.sparse.sparse.sdcsr(mshape, s=array(1), diag=array(1), M=None)

Bases: object

Scaled matrix, which is stored as

s \cdot diag( \mathbf{diag} ) \cdot \mathbf{M}

where s is scalar, diag is row scaling vector (scalar and vector allowed), and M is general part (None is allowed to indicate diagonal matrix).

mshape stores matrix shape. None for mshape[0] denotes differentation of scalar. None for mshape[1] denotes differentiation with repsect to scalar.

No copies of M, diag are made, therefore they must be constant objects.

broadcast(output)

Return broadcast matrix \mathbf{B_{output}} for broadcasting x to output, this matrix being Jacobian of x

chain(output, x)

Apply chain rule for elementwise operation

Jacobian of elementwise operation is diag(\mathbf{x}). Return diag(\mathbf{B_{output}} \cdot \mathbf{x})\cdot\mathbf{B_{output}}\cdot\mathbf{self}

classmethod fma(output, *terms)

Apply chain rule to elementwise functions

Returns sum(d.chain(output,x) for x,d in terms)

classmethod fma2(output, *terms)

Apply chain rule to elementwise functions

Returns sum(d.chain(output,x) for x,d in terms)

getitem_arrayp(output, idx)

Generate Jacobian matrix for operation output=x[idx], this matrix being Jacobian of x. idx is array with all entries positive.

getitem_general(output, idx)

Generate Jacobian matrix for operation output=x[idx], this matrix being Jacobian of x. General version.

classmethod new(mshape, diag=array(1), M=None)

Alternative constructor, which checks dimension of diag and assigns to scalar/vector part properly

rdot(y, other)

Return Jacobian of \mathbf{y} = \mathbf{other} \cdot \mathbf{self}, with \cdot denoting matrix multiplication.

sum()

Return Jacobian of y=sum(x), this matrix being Jacobian of x

tovalue()

Return this matrix as standard CSR matrix. The result is cached.

vstack(output, parts)

Return Jacobian of output=hstack(parts)

zero(output)

Return empty Jacobian, which would result from output=0*x, this matrix being Jacobian of x.

class sparsegrad.impl.sparse.sparse.sparsity_csr(mshape, s=None, diag=None, M=None)

Bases: sparsegrad.impl.sparse.sparse.sdcsr

This is a variant of matrix only propagating sparsity information

chain(output, x)

Apply chain rule for elementwise operation

Jacobian of elementwise operation is diag(\mathbf{x}). Return diag(\mathbf{B_{output}} \cdot \mathbf{x})\cdot\mathbf{B_{output}}\cdot\mathbf{self}

classmethod fma(output, *terms)

Apply chain rule to elementwise functions

Returns sum(d.chain(output,x) for x,d in terms)

classmethod fma2(output, *terms)

Apply chain rule to elementwise functions

Returns sum(d.chain(output,x) for x,d in terms)

rdot(y, other)

Return Jacobian of \mathbf{y} = \mathbf{other} \cdot \mathbf{self}, with \cdot denoting matrix multiplication.

sparsegrad.impl.sparse.sparse.sample_csr_rows(csr, rows)

return (indptr,ix) such that csr[rows]=csr_matrix((csr.data[ix],csr.indices[ix],indptr))

sparsegrad.impl.sparse.sparse.csr_matrix

alias of sparsegrad.impl.sparse.sparse.csr_matrix_nochecking

sparsegrad.impl.sparse.sparse.csc_matrix

alias of sparsegrad.impl.sparse.sparse.csc_matrix_unchecked

Module contents