Symmetric, Tridiagonal, Positive Definite System Solver
(from Golub-Van Loan)

Given an n-by-n symmetric, positive definite matrix A and a vector b in R^n,
the following algorithm overwrites b with the solution to Ax=b. It is assumed
that the diagonal of A is stored in d(1:n) and the superdiagonal in e(1:n-1).

for k=2:n
   t=e(k-1);
   e(k-1)=t/d(k-1);
   d(k)=d(k)-te(k-1);
end

for k=2:n
   b(k)=b(k)-e(k-1)b(k-1)
end

b(n)=b(n)/d(n)

for k=n-1:-1:1
   b(k)=b(k)/d(k)-e(k)b(k+1)
end