*DECK DINTYD
      SUBROUTINE DINTYD (T, K, YH, NYH, DKY, IFLAG)
C***BEGIN PROLOGUE  DINTYD
C***SUBSIDIARY
C***PURPOSE  Subsidiary to DDEBDF
C***LIBRARY   SLATEC
C***TYPE      DOUBLE PRECISION (INTYD-S, DINTYD-D)
C***AUTHOR  Watts, H. A., (SNLA)
C***DESCRIPTION
C
C   DINTYD approximates the solution and derivatives at T by polynomial
C   interpolation. Must be used in conjunction with the integrator
C   package DDEBDF.
C ----------------------------------------------------------------------
C DINTYD computes interpolated values of the K-th derivative of the
C dependent variable vector Y, and stores it in DKY.
C This routine is called by DDEBDF with K = 0,1 and T = TOUT, but may
C also be called by the user for any K up to the current order.
C (see detailed instructions in LSODE usage documentation.)
C ----------------------------------------------------------------------
C The computed values in DKY are gotten by interpolation using the
C Nordsieck history array YH.  This array corresponds uniquely to a
C vector-valued polynomial of degree NQCUR or less, and DKY is set
C to the K-th derivative of this polynomial at T.
C The formula for DKY is..
C              Q
C  DKY(I)  =  Sum  C(J,K) * (T - TN)**(J-K) * H**(-J) * YH(I,J+1)
C             J=K
C where  C(J,K) = J*(J-1)*...*(J-K+1), Q = NQCUR, TN = TCUR, H = HCUR.
C The quantities  NQ = NQCUR, L = NQ+1, N = NEQ, TN, and H are
C communicated by common.  The above sum is done in reverse order.
C IFLAG is returned negative if either K or T is out of bounds.
C ----------------------------------------------------------------------
C
C***SEE ALSO  DDEBDF
C***ROUTINES CALLED  (NONE)
C***COMMON BLOCKS    DDEBD1
C***REVISION HISTORY  (YYMMDD)
C   820301  DATE WRITTEN
C   890911  Removed unnecessary intrinsics.  (WRB)
C   891214  Prologue converted to Version 4.0 format.  (BAB)
C   900328  Added TYPE section.  (WRB)
C   910722  Updated AUTHOR section.  (ALS)
C***END PROLOGUE  DINTYD
C
      INTEGER I, IC, IER, IFLAG, IOWND, IOWNS, J, JB, JB2, JJ, JJ1,
     1      JP1, JSTART, K, KFLAG, L, MAXORD, METH, MITER, N, NFE,
     2      NJE, NQ, NQU, NST, NYH
      DOUBLE PRECISION C, DKY, EL0, H, HMIN, HMXI, HU, R, ROWND,
     1      ROWNS, S, T, TN, TP, UROUND, YH
      DIMENSION YH(NYH,*),DKY(*)
      COMMON /DDEBD1/ ROWND,ROWNS(210),EL0,H,HMIN,HMXI,HU,TN,UROUND,
     1                IOWND(14),IOWNS(6),IER,JSTART,KFLAG,L,METH,MITER,
     2                MAXORD,N,NQ,NST,NFE,NJE,NQU
C
C     BEGIN BLOCK PERMITTING ...EXITS TO 130
C***FIRST EXECUTABLE STATEMENT  DINTYD
         IFLAG = 0
         IF (K .LT. 0 .OR. K .GT. NQ) GO TO 110
            TP = TN - HU*(1.0D0 + 100.0D0*UROUND)
            IF ((T - TP)*(T - TN) .LE. 0.0D0) GO TO 10
               IFLAG = -2
C     .........EXIT
               GO TO 130
   10       CONTINUE
C
            S = (T - TN)/H
            IC = 1
            IF (K .EQ. 0) GO TO 30
               JJ1 = L - K
               DO 20 JJ = JJ1, NQ
                  IC = IC*JJ
   20          CONTINUE
   30       CONTINUE
            C = IC
            DO 40 I = 1, N
               DKY(I) = C*YH(I,L)
   40       CONTINUE
            IF (K .EQ. NQ) GO TO 90
               JB2 = NQ - K
               DO 80 JB = 1, JB2
                  J = NQ - JB
                  JP1 = J + 1
                  IC = 1
                  IF (K .EQ. 0) GO TO 60
                     JJ1 = JP1 - K
                     DO 50 JJ = JJ1, J
                        IC = IC*JJ
   50                CONTINUE
   60             CONTINUE
                  C = IC
                  DO 70 I = 1, N
                     DKY(I) = C*YH(I,JP1) + S*DKY(I)
   70             CONTINUE
   80          CONTINUE
C     .........EXIT
               IF (K .EQ. 0) GO TO 130
   90       CONTINUE
            R = H**(-K)
            DO 100 I = 1, N
               DKY(I) = R*DKY(I)
  100       CONTINUE
         GO TO 120
  110    CONTINUE
C
            IFLAG = -1
  120    CONTINUE
  130 CONTINUE
      RETURN
C     ----------------------- END OF SUBROUTINE DINTYD
C     -----------------------
      END