12
$\begingroup$

I am solving a physical problem using implicit numerical scheme. This leads me to solving a linear equation with tridiagonal matrix. I've coded this algorithm from Wikipedia. I wonder if there is an efficient library which allows to solve this type of equation in optimised way. An important note is that matrix itself changes only when system parameters are changing, so I had an opportunity to precalculate some algorithm steps for a nice peformance bonus. I am using C++.

$\endgroup$
  • $\begingroup$ How big of a system is it, does it need to be parallel? $\endgroup$ – aterrel Nov 30 '11 at 22:59
  • 1
    $\begingroup$ Size depends on required accuracy (from hundred to tens of thousands values). Now I am coding on one-core computer, but it's possible to get access to university supercomputer with many cpus available, so parallelism support would be nice. $\endgroup$ – gmk Nov 30 '11 at 23:17
15
$\begingroup$

You should probably start with the LAPACK implementation, ?gtsv, e.g., dgtsv. If you want a distributed-memory version, then you might want to start with ScaLAPACK's p?gtsv.

EDIT: Since your matrix does not change very often, you can avoid redundantly factoring the tridiagonal matrix by breaking up the LAPACK routine ?gtsv into the factorization step, ?gttrf, and the solve stage, ?gttrs. Similarly named routines exist in ScaLAPACK that serve the same purpose.

$\endgroup$
  • $\begingroup$ Thank you, it looks like what I need. I'll try now to run this routines from my code. $\endgroup$ – gmk Nov 30 '11 at 23:59
  • 1
    $\begingroup$ Since you are calling it from C++, make sure to declare the prototype inside of an extern "C" { } block. Depending on your system, you may need to append an underscore to the routine name. $\endgroup$ – Jack Poulson Dec 1 '11 at 0:06
2
$\begingroup$

For distributed parallel systems: I have not tried ScaLAPACK, which has a parallel tridiagonal solver, for which there are examples available online. I have tried with some success a method proposed by David Moulton in a LANL publication. Coding this up might be more than you want to do, but by using LAPACK, it is strait forward.

$\endgroup$
1
$\begingroup$

I found an interesting recursive algorithm here at page 975. It looks promising, I wonder what more experienced people say about it.

$\endgroup$
  • $\begingroup$ Numerical Recipes has some errors in it. In terms of a source of codes to use, it isn't the best, although some consider it a classic. I'd be surprised if ScaLAPACK didn't implement an algorithm at least as efficient as recursive cyclic reduction. $\endgroup$ – Geoff Oxberry Mar 6 '12 at 10:49

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.