share/odepack/fortran/nnsc.f

   1       subroutine nnsc
   2      *     (n, r, c, il, jl, ijl, l, d, iu, ju, iju, u, z, b, tmp)
   3 c*** subroutine nnsc
   4 c*** numerical solution of sparse nonsymmetric system of linear
   5 c      equations given ldu-factorization (compressed pointer storage)
   6 c
   7 c
   8 c       input variables..  n, r, c, il, jl, ijl, l, d, iu, ju, iju, u, b
   9 c       output variables.. z
  10 c
  11 c       parameters used internally..
  12 c fia   - tmp   - temporary vector which gets result of solving  ly = b.
  13 c       -           size = n.
  14 c
  15 c  internal variables..
  16 c    jmin, jmax - indices of the first and last positions in a row of
  17 c      u or l  to be used.
  18 c
  19       integer r(*), c(*), il(*), jl(*), ijl(*), iu(*), ju(*), iju(*)
  20 c     real l(*), d(*), u(*), b(*), z(*), tmp(*), tmpk, sum
  21       double precision  l(*), d(*), u(*), b(*), z(*), tmp(*), tmpk,sum
  22 c
  23 c  ******  set tmp to reordered b  *************************************
  24       do 1 k=1,n
  25    1    tmp(k) = b(r(k))
  26 c  ******  solve  ly = b  by forward substitution  *********************
  27       do 3 k=1,n
  28         jmin = il(k)
  29         jmax = il(k+1) - 1
  30         tmpk = -d(k) * tmp(k)
  31         tmp(k) = -tmpk
  32         if (jmin .gt. jmax) go to 3
  33         ml = ijl(k) - jmin
  34         do 2 j=jmin,jmax
  35    2      tmp(jl(ml+j)) = tmp(jl(ml+j)) + tmpk * l(j)
  36    3    continue
  37 c  ******  solve  ux = y  by back substitution  ************************
  38       k = n
  39       do 6 i=1,n
  40         sum = -tmp(k)
  41         jmin = iu(k)
  42         jmax = iu(k+1) - 1
  43         if (jmin .gt. jmax) go to 5
  44         mu = iju(k) - jmin
  45         do 4 j=jmin,jmax
  46    4      sum = sum + u(j) * tmp(ju(mu+j))
  47    5    tmp(k) = -sum
  48         z(c(k)) = -sum
  49         k = k - 1
  50    6    continue
  51       return
  52       end