I need to find solutions to a tridiagonal system. By this I mean a set of linear equations Ax = d where A is a square matrix containing elements A[i,i-1], A[i,i] and A[i,i+1] for i in 1:nrow, and zero elsewhere. R is probably not the ideal way to do this, but this is part of a larger problem that requires R.
In my application it is much easier (and much faster) to generate the diagonal and off-diagonal elements of A as vectors, i.e. a = A[i,i-1], b = A[i,i] and c = A[i,i+1]. So I have three vectors that define A, along with a solution vector d. The conventional method of solving such systems is to use the so-called "Thomas algorithm", see e.g. <http://www.enseeiht.fr/hmf/travaux/CD0001/travaux/optmfn/hi/01pa/hyb74/node24.html>. This is very easy to code, but much more difficult to "vectorize". Is anyone aware of a library that contains a fast implementation of this algorithm? Another alternative is to use backsolve. I can easily eliminate the lower diagonal a, but I'm still left with b and c, whereas backsolve requires a matrix. Again, I can write a function to read b and c into a matrix, but this requires loops, and is too slow. Is there a vectorized way of doing it? Of course, the diag command works for b, but what about c? In Octave, diag allows for an offset, but R apparently does not. I would appreciate any and all assistance you experts can offer. Thanks in advance. Will Harvey ______________________________________________ [EMAIL PROTECTED] mailing list https://www.stat.math.ethz.ch/mailman/listinfo/r-help
