The R Project SVN R

Rev

Rev 12778 | Show entire file | Ignore whitespace | Details | Blame | Last modification | View Log | RSS feed

Rev 12778 Rev 86632
Line 48... Line 48...
48 
c     cleve moler, university of new mexico, argonne national lab.
 48 
c     cleve moler, university of new mexico, argonne national lab.
49 
c
 49 
c
50 
c     subroutines and functions
 50 
c     subroutines and functions
51 
c
 51 
c
52 
c     blas daxpy,ddot
 52 
c     blas daxpy,ddot
53 
c     fortran min0
 53 
c     fortran min
54 
c
 54 
c
55 
c     internal variables
 55 
c     internal variables
56 
c
 56 
c
57 
      double precision ddot,t
 57 
      double precision ddot,t
58 
      integer k,kb,la,lb,lm
 58 
      integer k,kb,la,lb,lm
59 
c
 59 
c
60 
c     solve trans(r)*y = b
 60 
c     solve trans(r)*y = b
61 
c
 61 
c
62 
      do 10 k = 1, n
 62 
      do 10 k = 1, n
63 
         lm = min0(k-1,m)
 63 
         lm = min(k-1,m)
64 
         la = m + 1 - lm
 64 
         la = m + 1 - lm
65 
         lb = k - lm
 65 
         lb = k - lm
66 
         t = ddot(lm,abd(la,k),1,b(lb),1)
 66 
         t = ddot(lm,abd(la,k),1,b(lb),1)
67 
         b(k) = (b(k) - t)/abd(m+1,k)
 67 
         b(k) = (b(k) - t)/abd(m+1,k)
68 
   10 continue
 68 
   10 continue
69 
c
 69 
c
70 
c     solve r*x = y
 70 
c     solve r*x = y
71 
c
 71 
c
72 
      do 20 kb = 1, n
 72 
      do 20 kb = 1, n
73 
         k = n + 1 - kb
 73 
         k = n + 1 - kb
74 
         lm = min0(k-1,m)
 74 
         lm = min(k-1,m)
75 
         la = m + 1 - lm
 75 
         la = m + 1 - lm
76 
         lb = k - lm
 76 
         lb = k - lm
77 
         b(k) = b(k)/abd(m+1,k)
 77 
         b(k) = b(k)/abd(m+1,k)
78 
         t = -b(k)
 78 
         t = -b(k)
79 
         call daxpy(lm,t,abd(la,k),1,b(lb),1)
 79 
         call daxpy(lm,t,abd(la,k),1,b(lb),1)