Community
cancel
Showing results for 
Search instead for 
Did you mean: 
Pierre_L_4
Beginner
120 Views

Sherman Morrisson algorithm : almost tridiagonal matrices

Hello all,

Let M be a real tridiagonal matrix of size n greater than or equal to 2. We perfectly know how to solve systems of the form MX = Y with dgttrf followed by dgttrs. Now imagine that I am not interested in solving MX = Y anymore, but rather NX = Y where N = M + u.Tv where u,v are columns vectors of size n and where Tv means "transposed of v". Solving NX = Y is equivalent to solve MX = Y, MX = u and MX = v, and is "almost" not more expensive than solving MX = Y when the rank of the matrix u.Tv is small wrt n. (In my case this rank is inferior or equal to 2.)

This is due to the Sherman Morrisson formula, see here -->

http://en.wikipedia.org/wiki/Sherman%E2%80%93Morrison_formula

My question is : as far as I saw in MKL's documentation, there's no dedicated function doing this Sherman Morrisson trick : am I wrong or not ? Is there other function allowing me to solve system where the matrix defining the system is close to a tridiagonal matrix.

(If not, obviously I could then use dgttrf followed by a dgttrs with a three column matrix for the right hand side, obtained by putting X, u, v together.)

By close, I mean that there existe a tridiagonal matrix such that the difference between the matrix defining the system and this matrix is 1) a matrix of very small rank 2) and matrix with a very small number of non zero entries.

Thx

Pierre.

0 Kudos
5 Replies
120 Views

Hi Pierre,

Thanks a lot for a good question - you are correct, MKL doesn't support solution of system with tridiagonal+rank 1 matrices. But, as you proposed it can be solve by combination of dgttrf and dgttrs routine with different rhs. I didn't understand only one point from you proposal:

(If not, obviously I could then use dgttrf followed by a dgttrs with a three column matrix for the right hand side, obtained by putting X, u, v together.)

From my point of view you need to call dgttrs twice to implement numerator of fraction. The first one can be combine with computation of denominator be setting of several rhs but in any case you need to call dgttrs twice. Am I correct?

Thanks,

Alexander Kalinkin

 

Pierre_L_4
Beginner
120 Views

Hi Alexander, and thanks a lot !

I was inaaccurate I am sorry, I was meaning :

If not, obviously I could then use dgttrf (followed by dgttrs) with a three column matrix for the right hand side, obtained by putting X, u, v together

And you're right for calling dgttrs twice.

Kind Regards,

Pierre 

SergeyKostrov
Valued Contributor II
120 Views

Hi Pierre, You know that there are tens of thousand of different scientific algorithms and MKL has the most widelly used onces. I don't think it is possible at all to have all algorithms designed and implemented in scientific community during last 40 or 50 years in one magic library. For example, Intel IPP ( Integrated Performance Primitives ) library ( more than 10,000 different functions were at some point as far as I remember ) followed the path of having lots of different algorithms inside and the team is busy now with reducing the number of functions in order to reduce maintenance efforts.
Pierre_L_4
Beginner
120 Views

Yes I know. Speaking of which :

there is one thing I don't know but I would like to : I created some windows related topics in the intel mkl forum, and they are all visible, on the top (because I created them recently) :

why isn't it the case with the mac os x related topics I have created ?

Is mac os x intel mkl support still active ?

I am asking this because the two topics

http://software.intel.com/en-us/forums/topic/494329

and

http://software.intel.com/en-us/forums/topic/494365

are totally invisible in the intel mkl forum if you filter topic list to see most recent topics but where created in it.

Why so ? Only because it is mac related ?

PS : I really don't think "linking intel mkl is easy, especially on mac os x plateform.

Pierre_L_4
Beginner
120 Views

The two topics just appeared, by miracle. :-) Strange.

Reply