Solves linear SPD systems

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

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

Havard Rue hrue@r-inla.org

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")
#> Error in system(paste(shQuote(inla.call.no.remote()), "-s -m qsolve",     "-r", reordering, "-S", smtp, Qfile, Xfile, Bfile, method),     intern = TRUE): error in running command
print(paste("err", sum(abs( Q %*% X - B))))
#> Error in h(simpleError(msg, call)): error in evaluating the argument 'x' in selecting a method for function 'print': object 'X' not found

L = t(chol(Q))
X = inla.qsolve(Q, B, method = "forward")
#> Error in system(paste(shQuote(inla.call.no.remote()), "-s -m qsolve",     "-r", reordering, "-S", smtp, Qfile, Xfile, Bfile, method),     intern = TRUE): error in running command
print(paste("err", sum(abs( L %*% X - B))))
#> Error in h(simpleError(msg, call)): error in evaluating the argument 'x' in selecting a method for function 'print': object 'X' not found

X = inla.qsolve(Q, B, method = "backward")
#> Error in system(paste(shQuote(inla.call.no.remote()), "-s -m qsolve",     "-r", reordering, "-S", smtp, Qfile, Xfile, Bfile, method),     intern = TRUE): error in running command
print(paste("err", sum(abs( t(L) %*% X - B))))
#> Error in h(simpleError(msg, call)): error in evaluating the argument 'x' in selecting a method for function 'print': object 'X' not found

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")
#> Error in system(paste(shQuote(inla.call.no.remote()), "-s -m qsolve",     "-r", reordering, "-S", smtp, Qfile, Xfile, Bfile, method),     intern = TRUE): error in running command
print(paste("err", sum(abs( t(L) %*% X - B))))
#> Error in h(simpleError(msg, call)): error in evaluating the argument 'x' in selecting a method for function 'print': object 'X' not found
unlink(Q.file)
unlink(B.file)