qsolve: Solves linear SPD systems

In inbo/INLA: Full Bayesian Analysis of Latent Gaussian Models using Integrated Nested Laplace Approximations

Description

This routine use the GMRFLib implementation to solve linear systems with a SPD matrix.

Usage

     inla.qsolve(Q, B, reordering, method = c("solve", "forward", "backward"))
 

Arguments

Q	A SPD matrix, either as a (dense) matrix, sparse-matrix or a filename containing the matrix (in the fmesher-format).
B	The right hand side matrix, either as a (dense) matrix, sparse-matrix or a filename containing the matrix (in the fmesher-format). (Must be a matrix or sparse-matrix even if ncol(B) is 1.)
reordering	The type of reordering algorithm to be used for TAUCS; either one of the names listed in inla.reorderings() or the output from inla.qreordering(Q). The default is "auto" which try several reordering algorithm and use the best one for this particular matrix.
method	The system to solve, one of "solve", "forward" or "backward". Let Q = L L^T, where L is lower triangular (the Cholesky triangle), then method="solve" solves L L^T X = B or equivalently Q X = B, method="forward" solves L X = B, and method="backward" solves L^T X = B.

Value

inla.qsolve returns a matrix X, which is the solution of Q X = B, L X = B or L^T X = B depending on the value of method.

Author(s)

Havard Rue [e-mail]

Examples

n = 10
QQ = matrix(runif(n^2), n, n)
Q = inla.as.dgTMatrix(QQ %*% t(QQ))
B = matrix(runif(n^2-n), n, n-1)

X = inla.qsolve(Q, B, method = "solve")
print(paste("err", sum(abs( Q %*% X - B))))

L = t(chol(Q))
X = inla.qsolve(Q, B, method = "forward")
print(paste("err", sum(abs( L %*% X - B))))

X = inla.qsolve(Q, B, method = "backward")
print(paste("err", sum(abs( t(L) %*% X - B))))

Q.file = INLA:::inla.write.fmesher.file(Q)
B.file = INLA:::inla.write.fmesher.file(B)
X = inla.qsolve(Q.file, B.file, method = "backward")
print(paste("err", sum(abs( t(L) %*% X - B))))
unlink(Q.file)
unlink(B.file)

inbo/INLA documentation built on Dec. 6, 2019, 9:51 a.m.