solsy.f

Go to the documentation of this file.

      SUBROUTINE solsy (WM, IWM, X, TEM)
CLLL. OPTIMIZE
      INTEGER iwm
      INTEGER iownd, iowns,
     1   icf, ierpj, iersl, jcur, jstart, kflag, l, meth, miter,
     2   maxord, maxcor, msbp, mxncf, n, nq, nst, nfe, nje, nqu
      INTEGER i, meband, ml, mu
      DOUBLE PRECISION wm, x, tem
      DOUBLE PRECISION rowns, 
     1   ccmax, el0, h, hmin, hmxi, hu, rc, tn, uround
      DOUBLE PRECISION di, hl0, phl0, r 
      dimension wm(*), iwm(*), x(*), tem(*)
      COMMON /ls0001/ rowns(209),
     2   ccmax, el0, h, hmin, hmxi, hu, rc, tn, uround,
     3   iownd(14), iowns(6), 
     4   icf, ierpj, iersl, jcur, jstart, kflag, l, meth, miter,
     5   maxord, maxcor, msbp, mxncf, n, nq, nst, nfe, nje, nqu
C-----------------------------------------------------------------------
C THIS ROUTINE MANAGES THE SOLUTION OF THE LINEAR SYSTEM ARISING FROM 
C A CHORD ITERATION.  IT IS CALLED IF MITER .NE. 0.
C IF MITER IS 1 OR 2, IT CALLS DGETRS TO ACCOMPLISH THIS.
C IF MITER = 3 IT UPDATES THE COEFFICIENT H*EL0 IN THE DIAGONAL
C MATRIX, AND THEN COMPUTES THE SOLUTION.
C IF MITER IS 4 OR 5, IT CALLS DGBTRS.
C COMMUNICATION WITH SOLSY USES THE FOLLOWING VARIABLES..
C WM    = REAL WORK SPACE CONTAINING THE INVERSE DIAGONAL MATRIX IF
C         MITER = 3 AND THE LU DECOMPOSITION OF THE MATRIX OTHERWISE. 
C         STORAGE OF MATRIX ELEMENTS STARTS AT WM(3).
C         WM ALSO CONTAINS THE FOLLOWING MATRIX-RELATED DATA..
C         WM(1) = SQRT(UROUND) (NOT USED HERE),
C         WM(2) = HL0, THE PREVIOUS VALUE OF H*EL0, USED IF MITER = 3.
C IWM   = INTEGER WORK SPACE CONTAINING PIVOT INFORMATION, STARTING AT
C         IWM(21), IF MITER IS 1, 2, 4, OR 5.  IWM ALSO CONTAINS BAND 
C         PARAMETERS ML = IWM(1) AND MU = IWM(2) IF MITER IS 4 OR 5.
C X     = THE RIGHT-HAND SIDE VECTOR ON INPUT, AND THE SOLUTION VECTOR
C         ON OUTPUT, OF LENGTH N.
C TEM   = VECTOR OF WORK SPACE OF LENGTH N, NOT USED IN THIS VERSION. 
C IERSL = OUTPUT FLAG (IN COMMON).  IERSL = 0 IF NO TROUBLE OCCURRED. 
C         IERSL = 1 IF A SINGULAR MATRIX AROSE WITH MITER = 3.
C THIS ROUTINE ALSO USES THE COMMON VARIABLES EL0, H, MITER, AND N.
C-----------------------------------------------------------------------
      iersl = 0
      go to(100, 100, 300, 400, 400), miter
 100  CALL dgetrs( 'N', n, 1, wm(3), n, iwm(21), x, n, inlpck)
      RETURN
C
 300  phl0 = wm(2)
      hl0 = h*el0
      wm(2) = hl0
      IF (hl0 .EQ. phl0) go to 330
      r = hl0/phl0
      DO 320 i = 1,n
        di = 1.0d0 - r*(1.0d0 - 1.0d0/wm(i+2))
        IF (dabs(di) .EQ. 0.0d0) go to 390
 320    wm(i+2) = 1.0d0/di
 330  DO 340 i = 1,n
 340    x(i) = wm(i+2)*x(i)
      RETURN
 390  iersl = 1
      RETURN
C
 400  ml = iwm(1)
      mu = iwm(2)
      meband = 2*ml + mu + 1
      CALL dgbtrs( 'N', n, ml, mu, 1, wm(3), meband, iwm(21), x, n, 
     * inlpck)
      RETURN
C----------------------- END OF SUBROUTINE SOLSY -----------------------
      END