Solution of a tridiagonal system of linear equations

Share:

Description

Solves the linear system of equations

Ax=B

where A has to be square and tridiagonal, i.e with nonzero elements only on, one band above, and one band below the diagonal.

Usage

1
Solve.tridiag ( diam1, dia, diap1, B=rep(0,times=length(dia)))

Arguments

diam1

a vector with (nonzero) elements below the diagonal.

dia

a vector with (nonzero) elements on the diagonal.

diap1

a vector with (nonzero) elements above the diagonal.

B

Right-hand side of the equations, a vector with length = number of rows of A, or a matrix with number of rows = number of rows of A.

Details

If the length of the vector dia is equal to N, then the lengths of diam1 and diap1 should be equal to N-1

Value

matrix with the solution, X, of the tridiagonal system of equations Ax=B. The number of columns of this matrix equals the number of columns of B.

Author(s)

Karline Soetaert <karline.soetaert@nioz.nl>

See Also

Solve.banded, the function to solve a banded system of linear equations.

Solve.block, the function to solve a block diagonal system of linear equations.

Solve the generalised inverse solution,

solve the R default

Examples

 1
 2
 3
 4
 5
 6
 7
 8
 9
10
11
12
13
14
15
16
17
18
19
20
21
22
# create tridagonal system: bands on diagonal, above and below
nn   <- 20                          # nr rows and columns of A
aa   <- runif(nn)
bb   <- runif(nn)
cc   <- runif(nn)

# full matrix
A                        <- matrix(nrow = nn, ncol = nn, data = 0)
diag(A)                  <- bb
A[cbind(1:(nn-1), 2:nn)] <- cc[-nn]
A[cbind(2:nn, 1:(nn-1))] <- aa[-1]
B <- runif(nn)

# solve as full matrix
solve(A, B)                           

# same,  now using tridiagonal algorithm
as.vector(Solve.tridiag(aa[-1], bb, cc[-nn], B))

# same, now with 3 different right hand sides
B3 <- cbind(B, B*2, B*3)
Solve.tridiag(aa[-1], bb, cc[-nn], B3)