solve_cg {sanic} | R Documentation |
Solve Systems of Equations Using Iterative Methods
Description
Iterative solvers using the Conjugate Gradient method for sparse systems of
equations Ax = b
. Three different types are available: (1)
stabilized bi-conjugate gradient (BiCGSTAB) for square matrices, (2)
conjugate gradient for rectangular least-squares (LSCG), and (3) classic
conjugate gradient (CG) for symmetric positive definite matrices.
Usage
solve_cg(
a,
b,
x0,
type = c("BiCGSTAB", "LSCG", "CG"),
iter,
tol,
precond = 1L,
verbose = FALSE
)
Arguments
a |
Square numeric matrix with the coefficients of the linear system.
Dense and sparse matrices are supported, but the format must be sparse (see
|
b |
Numeric vector or matrix at the right-hand side of the linear system. If missing, 'b' is set to an identity matrix and 'a' is inverted. |
x0 |
Numeric vector or matrix with an initial guess. Must be of the same dimension as 'b'. |
type |
Character scalar. Whether to use the BiCGSTAB, least squares CG or classic CG method. |
iter |
Integer scalar with the maximum number of iterations. Defaults to the theoretical maximum, i.e. the number of columns in 'a'. |
tol |
Numeric scalar with the desired tolerance. Defaults to the machine precision. |
precond |
Integer scalar indicating the type of preconditioner to be used. Defaults to diagonal preconditioning. See the Details for further information. |
verbose |
Logical scalar. Whether to print iterations and tolerance. |
Details
Preconditioners can be set to 0
for no / identity preconditioning,
1
(default) for Jacobi / diagonal preconditioning, or 2
for
incomplete factorisation. Not all schemes are available for every type:
* type = "BiCGSTAB"
The default is precond = 1
for diagonal
preconditioning. Set precond = 0
for no preconditioning, or
precond = 2
for an incomplete LUT preconditioner.
* type = "LSCG"
The default is precond = 1
for diagonal least
squares preconditioning. Set precond = 0
for no preconditioning.
* type = "CG"
The default is precond = 1
for diagonal
preconditioning. Set precond = 0
for no preconditioning, or
precond = 2
for an incomplete Cholesky preconditioner.
Value
Solves for x
and returns a numeric matrix with the results.
Examples
set.seed(42)
x <- rnorm(3)
# Solve via BiCGSTAB for square matrices
A <- matrix(rnorm(9), nrow = 3, ncol = 3)
b <- A %*% x
norm(solve_cg(A, b, type = "B") - x)
# Solve via LSCG for rectangular matrices
A <- matrix(rnorm(12), nrow = 4, ncol = 3)
b <- A %*% x
norm(solve_cg(A, b, type = "LS") - x)
# Solve via classic CG for symmetric matrices
A <- crossprod(matrix(rnorm(9), nrow = 3, ncol = 3))
b <- A %*% x
norm(solve_cg(A, b, type = "CG") - x)
# The input matrix A should always be in sparse format
A <- sparsify(crossprod(matrix(rnorm(9), nrow = 3, ncol = 3)))
# The right-hand side should be a dense matrix
b <- as.matrix(A %*% x)
# We can check the speed of convergence and quality directly
solve_cg(A, b, verbose = TRUE)
# And provide guesses as starting value
solve_cg(A, b, x0 = x, verbose = TRUE)