# Block Tridiagonal and Block Bidiagonal solver

 Register Blogs Members List Search Today's Posts Mark Forums Read

 September 16, 2003, 21:38 Block Tridiagonal and Block Bidiagonal solver #1 Amit Guest   Posts: n/a Hi, Can anybody tell me where I can find an effecient Block Tridiagonal and a Block Bidiagonal Solver in C. I need this to solve implicit Euler equations. Thanks in advance. Amit.

 September 17, 2003, 20:01 Re: Block Tridiagonal and Block Bidiagonal solver #2 carlos Guest   Posts: n/a Have you tried the Numerical Recipes?... Look for it at the IMSL Fortran Library (Fortran 90) Hope that helps you.

 August 2, 2011, 15:02 #3 Member   Join Date: Jul 2011 Location: US Posts: 39 Rep Power: 7 There is a Block Tridiagonal Solver available with full source code here that I have used for many CFD projects. This is mostly applicable to 1D CFD problems, however.

August 2, 2011, 17:34
#4
Senior Member

Arjun
Join Date: Mar 2009
Location: Nurenberg, Germany
Posts: 482
Rep Power: 13
Quote:
 Originally Posted by Docfreezzzz There is a Block Tridiagonal Solver available with full source code here that I have used for many CFD projects. This is mostly applicable to 1D CFD problems, however.

This is based on Thomas' algorithm so in case of block matrices it is NOT efficient. If your case is large lets say block of 2000 * 2000 , it will take lot of time.

If possible one should try Block Cyclic reduction methods, they are very efficient for block matrices.

 August 3, 2011, 10:24 #5 Member   Join Date: Jul 2011 Location: US Posts: 39 Rep Power: 7 You are correct. However, that method is much more complex and generally harder to follow for the inexperienced. Also, may I suggest that for most tridiagonal problems the additional complexity is not warranted for such a small time savings. Any code that I have ever needed the capability to solve a block tridiagonal system is so fast to begin with that an increase in speed would be almost imperceptible. __________________ CFD engineering resource

August 3, 2011, 17:46
#6
Senior Member

Arjun
Join Date: Mar 2009
Location: Nurenberg, Germany
Posts: 482
Rep Power: 13
Quote:
 Originally Posted by Docfreezzzz You are correct. However, that method is much more complex and generally harder to follow for the inexperienced. Also, may I suggest that for most tridiagonal problems the additional complexity is not warranted for such a small time savings. Any code that I have ever needed the capability to solve a block tridiagonal system is so fast to begin with that an increase in speed would be almost imperceptible.
If take say mesh of 1000 x 1000 ( or say n x m )and discretesize Poisson problem on it.

The operation count in block cyclic reduction is roughly nmlog(m) .

But if you use the method of Thomas, you are looking at inverting m x m system n times . And have to do it 2 times because of forward substitution and then backward substitution. Now with LU factorization you operation count for inverting the system would be say 3m^2 So total cost is

(3m^2 ) x n x 2 + other operations like matrix matrix multiplications etc.

For 1000x 1000 mesh you are looking at

3 x 1000 x 1000 x 1000 x 2

compared to 1000 x 1000 x log2(1000) of block cyclic reduction.

The ratio by which block cyclic reduction is faster (roughly) : 6000 / log2(1000) .

Which I would not call close. In fact many iterative methods would be much faster than what Thomas algorithm would achieve in this case.

Quote:
 Originally Posted by Docfreezzzz You are correct. However, that method is much more complex and generally harder to follow for the inexperienced.
Yes, writing that code is difficult, so easy way out is to download fishpack and use the routine from there.
It is really difficult to code it yourself. I did it but would not recommend anyone to do it.

PS: There are approaximate block cyclic reduction algorithms exists, they ,instead of
nmlog2(m) , can do in nklog2(m) where k ~ 10 or 15.

That means for uniform mesh I could reduce the timing from 1000x1000xlog2(1000) to 15x1000xlog2(1000) . All this is different story though.

 August 3, 2011, 18:36 #7 Member   Join Date: Jul 2011 Location: US Posts: 39 Rep Power: 7 The cost of solving the discrete Euler equations on a 1000 node 1D grid (which results in a block tridiagonal system) with the algorithm I posted is more like 3*1000*(3^2), since the block sizes are only 3x3. This makes the situation far less dire than you propose with the Poisson equation. However, if you read my post carefully I was not disagreeing with you in the slightest. I was merely stating what a feel is a near universal truth. "Premature optimization is pure evil". That is, if you don't know what the costs are to begin with, in terms of wall clock time, then optimizing a routine based on theoretical operations counts will only get you into trouble. Often this sort of activity results in an algorithm that is only moderately faster in practice due to cache trashing, etc. and a whole boatload of undecipherable code for the next in line to deal with. Also, most 1D Euler codes which is what I assume the OP is writing are often not of much interest past a few thousand nodes/cvs. In this case, I would propose that the additional effort is simply not worth the end result. As for iterative methods, these are the only shot most of us get when dealing with "real-world" CFD problems due to the expense of direct methods. In fact, were the OP developing a more complex code I would suggest dropping direct methods altogether as they are completely unnecessary and often unusable on large problems due to memory limitations. Kudos. Thanks for the education on block cyclic methods. I had never had a need for these but I'll keep that in my toolbox should it ever arise. __________________ CFD engineering resource

August 3, 2011, 19:04
#8
Senior Member

Arjun
Join Date: Mar 2009
Location: Nurenberg, Germany
Posts: 482
Rep Power: 13
I give you that inverting 3x3 matrix is not very expensive. But it seems you missed my point. Which was Thomas algorithm is not very efficient for general block tridiagonal system. It could be used for smaller blocks where inverting them is not an issue.

This is why I said in my first post:

Quote:
 Originally Posted by arjun This is based on Thomas' algorithm so in case of block matrices it is NOT efficient. If your case is large lets say block of 2000 * 2000 , it will take lot of time. If possible one should try Block Cyclic reduction methods, they are very efficient for block matrices.

February 18, 2015, 17:53
#9
New Member

AJ
Join Date: Feb 2012
Location: Austin, TX
Posts: 10
Rep Power: 6
Quote:
 Originally Posted by Docfreezzzz There is a Block Tridiagonal Solver available with full source code here that I have used for many CFD projects. This is mostly applicable to 1D CFD problems, however.
I want to use this code. Can you please update the link?

 Thread Tools Display Modes Linear Mode

 Posting Rules You may not post new threads You may not post replies You may not post attachments You may not edit your posts BB code is On Smilies are On [IMG] code is On HTML code is OffTrackbacks are On Pingbacks are On Refbacks are On Forum Rules

 Similar Threads Thread Thread Starter Forum Replies Last Post Abullah Main CFD Forum 2 November 22, 2007 01:21 Abdullah Main CFD Forum 8 February 23, 2007 14:15 Abdulhafid M. Elfaghi Main CFD Forum 2 December 23, 2006 13:20 Grace Main CFD Forum 0 February 23, 2002 15:18 John Main CFD Forum 2 July 30, 2001 22:10

All times are GMT -4. The time now is 22:16.