subroutine ma28id(n, nz, aorg, irnorg, icnorg, licn, a, icn,
     * ikeep, rhs, x, r, w, mtype, prec, iflag)
c this subroutine uses the factors from an earlier call to ma28a/ad
c     or ma28b/bd to solve the system of equations with iterative
c     refinement.
c
c the parameters are...
c
c n is equal to the order of the matrix. it is not altered by the
c     subroutine.
c nz is equal to the number of entries in the original matrix.  it is
c     not altered by the subroutine.
c for this entry the original matrix must have been saved in
c     aorg,irnorg,icnorg where entry aorg(k) is in row irnorg(k) and
c     column icnorg(k), k=1,...nz.  information about the factors of a
c     is communicated to this subroutine via the parameters licn, a, icn
c     and ikeep where:
c aorg is an array of length nz.  not altered by ma28i/id.
c irnorg is an array of length nz.  not altered by ma28i/id.
c icnorg is an array of length nz.  not altered by ma28i/id.
c licn is equal to the length of arrays a and icn. it is not altered by
c     the subroutine.
c a is an array of length licn. it must be unchanged since the last call
c     to ma28a/ad or ma28b/bd. it is not altered by the subroutine.
c icn, ikeep are the arrays (of lengths licn and 5*n, respectively) of
c     the same names as in the previous all to ma28a/ad. they should be
c     unchanged since this earlier call and    they are not altered by
c     ma28i/id.
c the  other parameters are as follows:
c rhs is an array of length n. the user must set rhs(i) to contain the
c     value of the i th component of the right hand side. it is not
c     altered by ma28i/id.
c x is an array of length n. if an initial guess of the solution is
c     given (istart equal to 1), then the user must set x(i) to contain
c     the value of the i th component of the estimated solution.  on
c     exit, x(i) contains the i th component of the solution vector.
c r is an array of length n. it need not be set on entry.  on exit, r(i)
c     contains the i th component of an estimate of the error if maxit
c     is greater than 0.
c w is an array of length n. it is used as workspace by ma28i/id.
c mtype must be set to determine whether ma28i/id will solve a*x=rhs
c     (mtype equal to 1) or at*x=rhs (mtype ne 1, zero say). it is not
c     altered by ma28i/id.
c prec should be set by the user to the relative accuracy required. the
c     iterative refinement will terminate if the magnitude of the
c     largest component of the estimated error relative to the largest
c     component in the solution is less than prec.  it is not altered by
c     ma28i/id.
c iflag is a diagnostic flag which will be set to zero on successful
c     exit from ma28i/id, otherwise it will have a non-zero value. the
c     non-zero value iflag can have on exit from ma28i/id are ...
c     -16    indicating that more than maxit iteartions are required.
c     -17    indicating that more convergence was too slow.
c
      integer n, nz, licn, mtype, iflag
      integer icnorg(nz), irnorg(nz), ikeep(n,5), icn(licn)
      double precision a(licn), aorg(nz), rhs(n), r(n), x(n), w(n), prec
c
c private and common variables.
c see block data for comments on variables in common.
c see comments in code for use of private variables.
c
      logical lblock, grow, lbig
      integer lp, mp, ndrop, maxit, noiter, istart
      double precision tol, themax, big, dxmax, errmax, dres,
     * cgce, d, dd, conver, zero
      common /ma28ed/ lp, mp, lblock, grow
      common /ma28hd/ tol, themax, big, dxmax, errmax, dres, cgce,
     * ndrop, maxit, noiter, nsrch, istart, lbig
c
      data zero /0.0d0/
c
c   initialization of noiter, errmax and iflag.
c
      noiter = 0
      errmax = zero
      iflag = 0
c
c  jump if a starting vector has been supplied by the user.
c
      if (istart.eq.1) go to 20
c
c  make a copy of the right-hand side vector.
c
      do 10 i=1,n
        x(i) = rhs(i)
   10 continue
c
c  find the first solution.
c
      call ma28cd(n, a, licn, icn, ikeep, x, w, mtype)
c
c  stop the computations if   maxit=0.
c
   20 if (maxit.eq.0) go to 160
c
c  calculate the max-norm of the first solution.
c
      dd = 0.0
      do 30 i=1,n
        dd = dmax1(dd,dabs(x(i)))
   30 continue
      dxmax = dd
c
c  begin the iterative process.
c
      do 120 iterat=1,maxit
        d = dd
c
c  calculate the residual vector.
c
        do 40 i=1,n
          r(i) = rhs(i)
   40   continue
        if (mtype.eq.1) go to 60
        do 50 i=1,nz
          nrow = irnorg(i)
          ncol = icnorg(i)
          r(ncol) = r(ncol) - aorg(i)*x(nrow)
   50   continue
        go to 80
c mtype=1.
   60   do 70 i=1,nz
          nrow = irnorg(i)
          ncol = icnorg(i)
          r(nrow) = r(nrow) - aorg(i)*x(ncol)
   70   continue
   80   dres = 0.0
c
c  find the max-norm of the residual vector.
c
        do 90 i=1,n
          dres = dmax1(dres,dabs(r(i)))
   90   continue
c
c  stop the calculations if the max-norm of
c  the residual vector is zero.
c
        if (dres.eq.0.0) go to 150
c
c  calculate the correction vector.
c
        noiter = noiter + 1
        call ma28cd(n, a, licn, icn, ikeep, r, w, mtype)
c
c  find the max-norm of the correction vector.
c
        dd = 0.0
        do 100 i=1,n
          dd = dmax1(dd,dabs(r(i)))
  100   continue
c
c  check the convergence.
c
        if (dd.gt.d*cgce .and. iterat.ge.2) go to 130
        if (dxmax*10.0+dd.eq.dxmax*10.0) go to 140
c
c  attempt to improve the solution.
c
        dxmax = 0.0
        do 110 i=1,n
          x(i) = x(i) + r(i)
          dxmax = dmax1(dxmax,dabs(x(i)))
  110   continue
c
c  check the stopping criterion.
c
        if (dd.lt.prec*dxmax) go to 140
  120 continue
c more than maxit iterations required.
      iflag = -16
      write (lp,99999) iflag, maxit
      go to 140
c convergence rate unacceptably slow.
  130 iflag = -17
      conver = dd/d
      write (lp,99998) iflag, conver, cgce
c
c  the iterative process is terminated.
c
  140 errmax = dd
  150 continue
  160 return
99999 format (41h error return from ma28i/id with iflag = , i3/6h more ,
     * 4hthan, i5, 20h iterations required)
99998 format (38h error return from ma28i with iflag = , i3/9h converge,
     * 12hnce rate of , 1pe9.2, 9h too slow/24h maximum acceptable rate,
     * 8h set to , 1pe9.2)
      end