Help us shape the future of Go! Take the 2026 Go Developer Survey

sparse

package module
v0.1.0 Latest Latest
Warning

This package is not in the latest version of its module.

Go to latest
Published: Oct 20, 2026 License: BSD-3-Clause Imports: 6 Imported by: 0

README

sparse

Sparse matrices, fill-reducing orderings and iterative solvers in pure Go (CGO_ENABLED=0, every architecture the fleet's CI tests). Module github.com/go-sparse/sparse, BSD-3-Clause.

a := gallery.Laplacian2D(100, 100)           // *sparse.CSR[float64]
b := make([]float64, a.Rows)
ic, _ := iter.NewIC0(a)                      // incomplete Cholesky
res, _ := iter.CG(a, b, &iter.Options{RTol: 1e-8, M: ic})
fmt.Println(res.Info, res.Iterations)        // 0 = converged, as in SciPy

What is provided

package provides SciPy MATLAB
sparse COO, CSR, CSC, DOK over float64 and complex128; conversions (duplicates summed, indices sorted); O(1) Transpose, T, H; MulVec/MulVecT (SpMV), MulDense (SpMM), Mul (SpGEMM), Add/Sub/Axpby; Scale, ScaleRows, ScaleCols; Diagonal(k), Tril, Triu, Permute, Prune, norms; Eye, Diags; LinearOperator scipy.sparse.{coo,csr,csc,dok}_array, diags, eye, LinearOperator sparse, spdiags, speye, nnz, tril, triu
nd *ndarray.Array (go-ndarray v0.12.0) in and out: vectors, dense matrices, SpMV and SpMM – –
mtx Matrix Market read/write: coordinate and array, real/integer/pattern/complex, general/symmetric/skew-symmetric/hermitian, gzip input scipy.io.mmread, mmwrite mmread.m
order RCM, AMD, COLAMD; Bandwidth, Profile, CholeskyNNZ, QRFillNNZ reverse_cuthill_mckee symrcm, amd, symamd, colamd
iter CG, GMRES(m), BiCGSTAB, MINRES, LSQR, LSMR; preconditioners Jacobi, SSOR, IC0, ILU0, ILUT cg, gmres, bicgstab, minres, lsqr, lsmr pcg, gmres, bicgstab, minres, lsqr, ichol, ilu
gallery Laplacians (1-D/2-D/3-D), convection-diffusion, seeded random sparse matrices – gallery

Solvers follow scipy.sparse.linalg: the stopping rule is ‖b − A x‖ ≤ max(rtol·‖b‖, atol), rtol defaults to 1e-5, maxiter to 10n (5n for MINRES), info is 0 on success, the iteration limit when it was reached, and negative on breakdown; GMRES uses SciPy's left preconditioning and inner-tolerance control. The solvers take any sparse.LinearOperator.

Not provided yet: sparse direct factorisations (LU, Cholesky, QR), nested dissection, eigensolvers, algebraic multigrid, Rutherford-Boeing I/O.

Orderings: provenance

AMD and COLAMD are written from the papers (Amestoy, Davis and Duff, SIAM J. Matrix Anal. Appl. 17(4), 1996; Davis, Gilbert, Larimore and Ng, ACM TOMS 30(3), 2004) on one quotient-graph engine. SuiteSparse AMD and COLAMD (BSD-3-Clause) are the reference algorithms; CHOLMOD, UMFPACK and CSparse (GPL/LGPL) were not consulted. Tie-breaking differs, so permutations differ; quality is compared by fill instead.

Verification

docs/verification.md compares every solver, LSQR/LSMR and the orderings with SciPy 1.18.1 / NumPy 2.5.4 (and SuiteSparse AMD through cvxopt as a binary oracle). A Go harness (cmd/genfixtures) drives SciPy over a pipe, using only the minimal SciPy call per request, and records the answers in internal/conformance/testdata/oracle.json; the regression tests re-run the Go side from that file, so they need neither Python nor network. cmd/conformance -corpus also compares live on small Harwell-Boeing matrices from NIST Matrix Market, downloaded into a cache outside the repository (their licence is unstated) and checked against pinned SHA-256 sums.

Results on the committed suite: CG, GMRES, MINRES, LSQR and LSMR iteration counts equal SciPy's; BiCGSTAB agrees within a few iterations. AMD Cholesky fill is within 1% of SuiteSparse AMD; COLAMD LU fill is within 4% of SuperLU's COLAMD; RCM bandwidth equals SciPy's.

Performance

go run ./cmd/perfbench measures SpMV, CG and GMRES(30) on 2-D Laplacians of 1e5 and 1e6 nonzeros against SciPy on one core, interleaved, one at a time, and refuses to measure while the host's load average is above 8. Each batch lasts about 200 ms on both sides, after one untimed call: shorter batches, each run after waiting on the other process, measured macOS waking the process on an efficiency core more than the kernel. The table is docs/perf.md, measured on 2026-10-20 at load 5 (Go time / SciPy time):

work SpMV CG GMRES(30)
2-D Laplacian, 1.1e5 nonzeros 0.97 0.70 0.56
2-D Laplacian, 1.0e6 nonzeros 0.98 0.89 0.82

SpMV sums four terms per loop step, still one chain in storage order, so the result is bit-identical to SciPy's csr_matvec order.

Column indices are int32, as in scipy.sparse: SpMV streams one index per value, so their width is a third of its memory traffic. Dot products and the vector updates use the fleet's SIMD kernels (go-simd/floats v0.4.0 Dot, Axpy, Xpby, AxpyDot, generated by go-asmgen; the update kernels were added there for this package). The Krylov solvers fuse each vector update with the dot product that follows it (CG: r −= αq with the next r·r, one AxpyDot; GMRES: each Gram–Schmidt step with the next coefficient), one pass instead of two; without a preconditioner no copy is made.

Untrusted input

mtx.Read is meant to be safe on hostile files. mtx.ReadLimited takes Limits (defaults: dimensions 2^27, 2^28 entries, 64 KiB lines, 16 GiB read after decompression, 4096 comment lines); declared sizes are checked before anything is allocated from them, the stored-value count is computed in int64 (so it cannot overflow a 32-bit int), memory grows only as entries are actually read, and a gzip bomb is cut off at the byte limit. NewCSR/NewCSC validate hand-built arrays, COO.Validate does the same for a hand-built COO, and FromDense checks rows*cols without overflow. mtx.FuzzRead has a committed corpus under mtx/testdata/fuzz and its seeds run in every go test.

Quality gates

gofmt -l, go vet, go test -race with 100% statement coverage over all packages (-coverpkg=./...), 386 and the other architectures in .github/workflows/ci.yml, govulncheck.

Tooling rule

Generators, harnesses and comparisons are Go programs. The only Python is the oracle's own SciPy call (internal/scipy/script.go).

Documentation

Overview

Package sparse provides sparse matrices (COO, CSR, CSC, DOK), conversions, products, sums and scaling in pure Go, plus the LinearOperator abstraction used by the iterative solvers in package iter.

Element types are float64 and complex128 (the Scalar constraint). Indices are Go ints, so the package works unchanged on 32-bit targets.

Equivalents: scipy.sparse.{coo,csr,csc,dok}_array, MATLAB sparse, find, nnz, spdiags, speye.

Index

Constants

View Source
const MaxDim = math.MaxInt32

MaxDim is the largest number of rows or columns of a CSR or CSC matrix: minor indices are stored as int32.

Variables

View Source
var (
	// ErrShape reports incompatible dimensions.
	ErrShape = errors.New("sparse: incompatible shapes")
	// ErrIndex reports an index outside the matrix.
	ErrIndex = errors.New("sparse: index out of range")
	// ErrFormat reports a malformed compressed structure.
	ErrFormat = errors.New("sparse: malformed structure")
)

Errors returned (or wrapped) by this package.

Functions

This section is empty.

Types

type COO

type COO[T Scalar] struct {
	Rows, Cols int
	I, J       []int
	V          []T
}

COO is a coordinate-format matrix: parallel slices of row index, column index and value. Duplicate (i, j) entries are allowed and are summed on conversion, as in scipy.sparse.coo_array.

func NewCOO

func NewCOO[T Scalar](rows, cols int) *COO[T]

NewCOO returns an empty rows x cols COO matrix.

func (*COO[T]) Append

func (a *COO[T]) Append(i, j int, v T) error

Append adds entry (i, j) = v.

func (*COO[T]) Dims

func (a *COO[T]) Dims() (int, int)

Dims returns the number of rows and columns.

func (*COO[T]) NNZ

func (a *COO[T]) NNZ() int

NNZ returns the number of stored entries (duplicates counted).

func (*COO[T]) ToCSC

func (a *COO[T]) ToCSC() *CSC[T]

ToCSC converts to CSC, summing duplicates and sorting row indices.

func (*COO[T]) ToCSR

func (a *COO[T]) ToCSR() *CSR[T]

ToCSR converts to CSR, summing duplicates and sorting column indices.

func (*COO[T]) Transpose

func (a *COO[T]) Transpose() *COO[T]

Transpose returns the transposed COO matrix (shares no storage).

func (*COO[T]) Validate

func (a *COO[T]) Validate() error

Validate checks a COO matrix assembled by hand (for example from untrusted arrays): equal slice lengths, non-negative dimensions and every index inside the matrix. Matrices built with Append, Transpose or mtx.Read are valid by construction; ToCSR and ToCSC assume validity.

type CSC

type CSC[T Scalar] struct {
	Rows, Cols int
	Indptr     []int
	Indices    []int32
	Data       []T
}

CSC is a compressed-sparse-column matrix. Column j occupies Indices[Indptr[j]:Indptr[j+1]] (row indices, strictly increasing).

func NewCSC

func NewCSC[T Scalar](rows, cols int, indptr []int, indices []int32, data []T) (*CSC[T], error)

NewCSC wraps the given arrays (without copying) after validation.

func (*CSC[T]) At

func (a *CSC[T]) At(i, j int) T

At returns A[i, j].

func (*CSC[T]) Clone

func (a *CSC[T]) Clone() *CSC[T]

Clone returns a deep copy.

func (*CSC[T]) Dense

func (a *CSC[T]) Dense() []T

Dense returns the matrix as a row-major dense slice.

func (*CSC[T]) Diagonal

func (a *CSC[T]) Diagonal(k int) []T

Diagonal returns the k-th diagonal.

func (*CSC[T]) Dims

func (a *CSC[T]) Dims() (int, int)

Dims returns the number of rows and columns.

func (*CSC[T]) MulVec

func (a *CSC[T]) MulVec(dst, x []T)

MulVec computes dst = A x (column-oriented scatter).

func (*CSC[T]) MulVecT

func (a *CSC[T]) MulVecT(dst, x []T)

MulVecT computes dst = Aᵀ x.

func (*CSC[T]) NNZ

func (a *CSC[T]) NNZ() int

NNZ returns the number of stored entries.

func (*CSC[T]) Scale

func (a *CSC[T]) Scale(s T)

Scale multiplies every entry by s in place.

func (*CSC[T]) ToCOO

func (a *CSC[T]) ToCOO() *COO[T]

ToCOO converts to coordinate format.

func (*CSC[T]) ToCSR

func (a *CSC[T]) ToCSR() *CSR[T]

ToCSR converts to CSR (a copy).

func (*CSC[T]) Transpose

func (a *CSC[T]) Transpose() *CSR[T]

Transpose returns Aᵀ as a CSR matrix sharing A's arrays (O(1)).

type CSR

type CSR[T Scalar] struct {
	Rows, Cols int
	Indptr     []int
	Indices    []int32
	Data       []T
}

CSR is a compressed-sparse-row matrix. Row i occupies Indices[Indptr[i]:Indptr[i+1]] / Data[...] with strictly increasing column indices (canonical form, as produced by every constructor here). Column indices are int32, as scipy.sparse stores them whenever they fit: the products stream one index per stored value, so their width is a third of the memory traffic (dimensions are limited to MaxDim).

func Add

func Add[T Scalar](a, b *CSR[T]) (*CSR[T], error)

Add returns A + B.

func Axpby

func Axpby[T Scalar](alpha T, a *CSR[T], beta T, b *CSR[T]) (*CSR[T], error)

Axpby returns alpha*A + beta*B (sparse + sparse).

func Diags

func Diags[T Scalar](diagonals [][]T, offsets []int, rows, cols int) (*CSR[T], error)

Diags builds a rows x cols matrix from diagonals placed at the given offsets (scipy.sparse.diags / MATLAB spdiags). diagonals[t] has the length of diagonal offsets[t].

func Eye

func Eye[T Scalar](n int) *CSR[T]

Eye returns the n x n identity.

func FromDense

func FromDense[T Scalar](rows, cols int, d []T) (*CSR[T], error)

FromDense builds a CSR matrix from a row-major dense slice, dropping zeros.

func Mul

func Mul[T Scalar](a, b *CSR[T]) (*CSR[T], error)

Mul returns the sparse product A B (Gustavson's row-by-row algorithm).

func NewCSR

func NewCSR[T Scalar](rows, cols int, indptr []int, indices []int32, data []T) (*CSR[T], error)

NewCSR wraps the given arrays (without copying) after validating that they form a canonical CSR matrix.

func Sub

func Sub[T Scalar](a, b *CSR[T]) (*CSR[T], error)

Sub returns A - B.

func (*CSR[T]) At

func (a *CSR[T]) At(i, j int) T

At returns A[i, j] (zero when not stored) by binary search.

func (*CSR[T]) Clone

func (a *CSR[T]) Clone() *CSR[T]

Clone returns a deep copy.

func (*CSR[T]) Dense

func (a *CSR[T]) Dense() []T

Dense returns the matrix as a row-major dense slice of length Rows*Cols.

func (*CSR[T]) Diagonal

func (a *CSR[T]) Diagonal(k int) []T

Diagonal returns the k-th diagonal (k > 0 above the main one).

func (*CSR[T]) Dims

func (a *CSR[T]) Dims() (int, int)

Dims returns the number of rows and columns.

func (*CSR[T]) H

func (a *CSR[T]) H() *CSR[T]

H returns the conjugate transpose Aᴴ as a new CSR matrix.

func (*CSR[T]) IsSymmetric

func (a *CSR[T]) IsSymmetric(tol float64) bool

IsSymmetric reports whether A equals Aᵀ within |a_ij - a_ji| <= tol.

func (*CSR[T]) MulDense

func (a *CSR[T]) MulDense(dst, b []T, k int)

MulDense computes dst = A B for row-major dense B (Cols x k) and dst (Rows x k): the SpMM kernel (scipy A @ B for a dense B).

func (*CSR[T]) MulVec

func (a *CSR[T]) MulVec(dst, x []T)

MulVec computes dst = A x.

func (*CSR[T]) MulVecT

func (a *CSR[T]) MulVecT(dst, x []T)

MulVecT computes dst = Aᵀ x (plain transpose, no conjugation).

func (*CSR[T]) NNZ

func (a *CSR[T]) NNZ() int

NNZ returns the number of stored entries.

func (*CSR[T]) Norm1

func (a *CSR[T]) Norm1() float64

Norm1 is the maximum absolute column sum.

func (*CSR[T]) NormFro

func (a *CSR[T]) NormFro() float64

NormFro is the Frobenius norm.

func (*CSR[T]) NormInf

func (a *CSR[T]) NormInf() float64

NormInf is the maximum absolute row sum.

func (*CSR[T]) Permute

func (a *CSR[T]) Permute(p, q []int) (*CSR[T], error)

Permute returns B with B[i, j] = A[p[i], q[j]]. A nil permutation is the identity.

func (*CSR[T]) Prune

func (a *CSR[T]) Prune(tol float64)

Prune removes stored entries with |value| <= tol (tol = 0 removes explicit zeros) in place.

func (*CSR[T]) Scale

func (a *CSR[T]) Scale(s T)

Scale multiplies every entry by s in place.

func (*CSR[T]) ScaleCols

func (a *CSR[T]) ScaleCols(d []T) error

ScaleCols computes A * diag(d) in place.

func (*CSR[T]) ScaleRows

func (a *CSR[T]) ScaleRows(d []T) error

ScaleRows computes diag(d) * A in place.

func (*CSR[T]) T

func (a *CSR[T]) T() *CSR[T]

T returns Aᵀ as a new CSR matrix.

func (*CSR[T]) ToCOO

func (a *CSR[T]) ToCOO() *COO[T]

ToCOO converts to coordinate format.

func (*CSR[T]) ToCSC

func (a *CSR[T]) ToCSC() *CSC[T]

ToCSC converts to CSC (a copy).

func (*CSR[T]) Transpose

func (a *CSR[T]) Transpose() *CSC[T]

Transpose returns Aᵀ as a CSC matrix sharing A's arrays (O(1)).

func (*CSR[T]) Tril

func (a *CSR[T]) Tril(k int) *CSR[T]

Tril returns the lower triangle (entries with j - i <= k).

func (*CSR[T]) Triu

func (a *CSR[T]) Triu(k int) *CSR[T]

Triu returns the upper triangle (entries with j - i >= k).

func (*CSR[T]) Validate

func (a *CSR[T]) Validate() error

Validate checks the canonical-form invariants.

type DOK

type DOK[T Scalar] struct {
	Rows, Cols int
	// contains filtered or unexported fields
}

DOK is a dictionary-of-keys matrix, convenient for incremental construction (scipy.sparse.dok_array).

func NewDOK

func NewDOK[T Scalar](rows, cols int) *DOK[T]

NewDOK returns an empty rows x cols DOK matrix.

func (*DOK[T]) At

func (d *DOK[T]) At(i, j int) T

At returns the value at (i, j).

func (*DOK[T]) Dims

func (d *DOK[T]) Dims() (int, int)

Dims returns the number of rows and columns.

func (*DOK[T]) Len

func (d *DOK[T]) Len() int

Len returns the number of stored entries.

func (*DOK[T]) Set

func (d *DOK[T]) Set(i, j int, v T) error

Set stores v at (i, j); storing zero deletes the entry.

func (*DOK[T]) ToCOO

func (d *DOK[T]) ToCOO() *COO[T]

ToCOO returns the entries as a COO matrix (order unspecified).

func (*DOK[T]) ToCSR

func (d *DOK[T]) ToCSR() *CSR[T]

ToCSR converts to CSR.

type FuncOperator

type FuncOperator struct {
	Rows, Cols int
	Mv, MvT    func(dst, x []float64)
}

FuncOperator builds an operator from functions. MvT may be nil when no transpose product is available.

func (FuncOperator) Dims

func (f FuncOperator) Dims() (int, int)

Dims implements LinearOperator.

func (FuncOperator) MulVec

func (f FuncOperator) MulVec(dst, x []float64)

MulVec implements LinearOperator.

func (FuncOperator) MulVecT

func (f FuncOperator) MulVecT(dst, x []float64)

MulVecT implements TransposeOperator; it panics when MvT is nil.

type LinearOperator

type LinearOperator interface {
	// Dims returns the number of rows and columns of A.
	Dims() (rows, cols int)
	// MulVec stores A x in dst; len(x) == cols, len(dst) == rows.
	MulVec(dst, x []float64)
}

LinearOperator is a real matrix-free linear map x -> A x, the abstraction the iterative solvers work on (scipy.sparse.linalg.LinearOperator). *CSR[float64] and *CSC[float64] implement it.

func Shifted

func Shifted(a LinearOperator, shift float64) LinearOperator

Shifted returns the operator A + shift*I.

type Scalar

type Scalar interface{ float64 | complex128 }

Scalar is the element type of a sparse matrix.

type TransposeOperator

type TransposeOperator interface {
	LinearOperator
	// MulVecT stores Aᵀ x in dst; len(x) == rows, len(dst) == cols.
	MulVecT(dst, x []float64)
}

TransposeOperator is a LinearOperator that can also apply Aᵀ, which LSQR requires. *CSR[float64] and *CSC[float64] implement it.

Directories

Path Synopsis
cmd
conformance command
Command conformance writes docs/verification.md: the committed fixture compared with the Go implementation, and with -corpus a live comparison against SciPy on the NIST Harwell-Boeing matrices, downloaded into a cache directory outside the repository.
Command conformance writes docs/verification.md: the committed fixture compared with the Go implementation, and with -corpus a live comparison against SciPy on the NIST Harwell-Boeing matrices, downloaded into a cache directory outside the repository.
genfixtures command
Command genfixtures runs the conformance suite against a live SciPy oracle and writes the fixture file the regression tests read.
Command genfixtures runs the conformance suite against a live SciPy oracle and writes the fixture file the regression tests read.
perfbench command
Command perfbench compares this module with SciPy on one core: SpMV, CG and GMRES on 2-D Laplacians of about 1e5 and 1e6 nonzeros.
Command perfbench compares this module with SciPy on one core: SpMV, CG and GMRES on 2-D Laplacians of about 1e5 and 1e6 nonzeros.
Package gallery generates deterministic test matrices: finite-difference Laplacians, a convection-diffusion operator, and seeded random sparse matrices.
Package gallery generates deterministic test matrices: finite-difference Laplacians, a convection-diffusion operator, and seeded random sparse matrices.
internal
conformance
Package conformance compares this module with SciPy.
Package conformance compares this module with SciPy.
scipy
Package scipy drives SciPy as a verification oracle: one Python process per Session, spoken to in JSON lines over a pipe.
Package scipy drives SciPy as a verification oracle: one Python process per Session, spoken to in JSON lines over a pipe.
Package iter provides Krylov-subspace solvers for sparse linear systems (CG, GMRES(m), BiCGSTAB, MINRES, LSQR) and preconditioners (Jacobi, SSOR, IC(0), ILU(0), ILUT).
Package iter provides Krylov-subspace solvers for sparse linear systems (CG, GMRES(m), BiCGSTAB, MINRES, LSQR) and preconditioners (Jacobi, SSOR, IC(0), ILU(0), ILUT).
Package mtx reads and writes the Matrix Market exchange format (https://math.nist.gov/MatrixMarket/formats.html): coordinate and array layouts; real, integer, pattern and complex fields; general, symmetric, skew-symmetric and hermitian symmetry.
Package mtx reads and writes the Matrix Market exchange format (https://math.nist.gov/MatrixMarket/formats.html): coordinate and array layouts; real, integer, pattern and complex fields; general, symmetric, skew-symmetric and hermitian symmetry.
Package nd connects sparse matrices to go-ndarray arrays: dense vectors and matrices go in and out of the solvers as *ndarray.Array.
Package nd connects sparse matrices to go-ndarray arrays: dense vectors and matrices go in and out of the solvers as *ndarray.Array.
Package order computes sparse-matrix orderings and measures their quality.
Package order computes sparse-matrix orderings and measures their quality.

Jump to

Keyboard shortcuts

? : This menu
/ : Search site
f or F : Jump to
y or Y : Canonical URL