[R] Solving a tridiagonal system

[Gabor Grothendieck]

I just noticed that you defined a to be the lower diagonal
whereas I had it as the upper diagonal so the previous email
should be as follows to correspond to your notation:


# define test vectors for a, b and c
va <- -(1:4); vb <- 11:15; vc <- 1:4

# diag.num is a matrix whose ith diagonal equals i
diag.num <- -outer(seq(vb),seq(vb),"-")

x <- diag(vb)
x[diag.num == -1] <- va
x[diag.num == 1] <- vc


From: Gabor Grothendieck <ggrothendieck at myway.com>
To: <will_harvey03 at yahoo.com>, <r-help at stat.math.ethz.ch> 
Subject: RE: [R] Solving a tridiagonal system 

 
 

The following will create a matrix x with given sub diagonal,
diagonal and super diagonal.

# define test vectors for a, b and c
va <- -(1:4); vb <- 11:15; vc <- 1:4

# diag.num is a matrix whose ith super diagonal equals i and 
# sub diagonal equals -i
diag.num <- -outer(seq(vb),seq(vb),"-")

x <- diag(vb)
x[diag.num == 1] <- va
x[diag.num == -1] <- vc

---

Date: Wed, 1 Oct 2003 05:56:34 -0700 (PDT) 
From: Will Harvey <will_harvey03 at yahoo.com>
To: <r-help at stat.math.ethz.ch> 
Subject: [R] Solving a tridiagonal system 

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

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help




message or move to Select a Folder----- Folders ------InboxDraftsSentTrashBulk Mail---- My Folders ----R
Back to Inbox | <Prev Next> 

As AttachmentAs Inline Text 












Make My Way Your Home Page | Spread the Word

My Settings: Overview | Search | Email | Chat | Portfolio | Calendar | Groups | Profile 


IMPORTANT: We do not present our users with pop-ups, banners or any other non-text advertising. Nor do we send email
to our users. If you see or receive one of these items, it is coming from an outside source, either as a result of something you
have previously downloaded or as an "exit" pop-up from the site you just visited. It is not coming from our site.

Privacy Policy Terms of Service Partner with Us Our Mission Advertise with Us Sign In Sign Out Help Center

© 2002-2003 My Way









_______________________________________________
No banners. No pop-ups. No kidding.
Introducing My Way - http://www.myway.com

______________________________________________
R-help at stat.math.ethz.ch mailing list
https://www.stat.math.ethz.ch/mailman/listinfo/r-help

 
 
 
 message or move to Select a Folder----- Folders ------InboxDraftsSentTrashBulk Mail---- My Folders ----R
 Back to Inbox | <Prev Next> 
 
   As AttachmentAs Inline Text 
 
 



    
 
 




Make My Way Your Home Page  |  Spread the Word

My Settings: Overview | Search | Email | Chat | Portfolio | Calendar | Groups | Profile 
 

IMPORTANT: We do not present our users with pop-ups, banners or any other non-text advertising. Nor do we send email
to our users. If you see or receive one of these items, it is coming from an outside source, either as a result of something you
have previously downloaded or as an "exit" pop-up from the site you just visited. It is not coming from our site.

Privacy Policy   Terms of Service   Partner with Us   Our Mission   Advertise with Us   Sign In   Sign Out   Help Center

© 2002-2003 My Way

 
 

  
 

        

_______________________________________________
No banners. No pop-ups. No kidding.
Introducing My Way - http://www.myway.com