      SUBROUTINE NEWTON(NL, NR, NI)
C
C        ===============================================================
C
C **  THIS PROGRAM COMPUTES THE DIVIDED DIFFERENCES ARRAY AS FOLLOWS
C         NL COALESCED POINTS ON LEFT   - DERIV VALUES IN FLEFT
C         NR COALESCED POINTS ON RIGHT  -   -      -    - FRIGHT
C         NI DISTINCT  POINTS INBETWEEN - FNCTN    -    - FINTRP
C
C     THE POINTS ARE ORDERED XL = XLEFT (NSTACK)
C                            XR = XRIGHT(NSTACK)
C                            XINTRP ARRAY
C
C         LAYOUT OF THE DDTEMP DIVIDED DIFFERENCE ARRAY
C
C     NL=6    LLLLLL****II
C     NR=4    LLLLL****II     L = FIRST TRIANGLE
C     NI=2    LLLL****II
C             LLL****II       R = SECOND TRIANGLE
C             LL****II
C             L****II         * = FILL BETWEEN TRIANGLES
C             RRRRII
C             RRRII           I = COMPLETION FOR INTERPOLATION POINTS
C             RRII
C             RII       IDIF = HORIZONTAL COORD. = DIFFERENCE ORDER
C             II        IPT  = VERTICAL COORD. ASSOCIATED WITH POINTS
C             I
C
C
      DOUBLE PRECISION A, ACCUR, B, BLEFT, BRIGHT, CHARF, DDTEMP,
     * DSCTOL, ERROR, ERRORI, FACTOR, FINTRP, FLEFT, FRIGHT, NORM,
     * XBREAK, XDD, XINTRP, XLEFT, XRIGHT, BUFFER
      DIMENSION XBREAK(20), DBREAK(20), BLEFT(20), BRIGHT(20)
      DIMENSION XLEFT(50), XRIGHT(50)
      DIMENSION DDTEMP(20,20), FACTOR(20), FINTRP(18), FLEFT(10),
     * FRIGHT(10), XDD(20), XINTRP(18)
      INTEGER BOTH, BREAK, DBREAK, DEGREE, EDIST, RIGHT, RIGHTX, SMOOTH
      LOGICAL DISCRD
      DOUBLE PRECISION DIFFF, DIFFX
C
      COMMON /INPUTZ/ A, B, ACCUR, NORM, CHARF, XBREAK, BLEFT, BRIGHT,
     * DBREAK, DEGREE, SMOOTH, LEVEL, EDIST, NBREAK, KNTDIM, NPARDM
      COMMON /RESULZ/ ERROR, KNOTS
      COMMON /KONTRL/ DSCTOL, ERRORI, XLEFT, XRIGHT, BREAK, BOTH,
     * FACTOR, IBREAK, INTERP, LEFT, MAXAUX, MAXKNT, MAXPAR, MAXSTK,
     * NPAR, NSTACK, RIGHT, DISCRD, BUFFER
      COMMON /COMDIF/ DDTEMP, FINTRP, FLEFT, FRIGHT, XDD, XINTRP,
     * LEFTX, NINTRP, RIGHTX
C
C                 MAIN CALCULATION OF DIVIDED DIFFERENCES
C         DEFINE A FEW SHORT CONSTANTS
      NL1 = NL - 1
      NL2 = NL + 1
      NR1 = NR - 1
      NR2 = NR + 1
      NRL = NR + NL
C
C         PUT X-VALUES IN A SINGLE ARRAY WITH NDDX = NL+NR+NI POINTS
      DO 10 NDDX=1,NL
        XDD(NDDX) = XLEFT(NSTACK)
   10 CONTINUE
      NDDX = NL
      DO 20 K=1,NR
        NDDX = NDDX + 1
        XDD(NDDX) = XRIGHT(NSTACK)
   20 CONTINUE
C            CHECK IF THERE ARE ANY INTERPOLATION POINTS TO ADD TO XDD
      IF (NI.EQ.0) GO TO 40
      DO 30 K=1,NI
        NDDX = NDDX + 1
        XDD(NDDX) = XINTRP(K)
   30 CONTINUE
C
C           FILL BORDER OF FIRST TRIANGLE - SIZE NL.
   40 CONTINUE
C         TOP BORDER
      DO 50 IDIF=1,NL
        DDTEMP(IDIF,1) = FLEFT(IDIF)/FACTOR(IDIF)
   50 CONTINUE
      IF (NL1.EQ.0) GO TO 70
C                   BOTTOM BORDER
      DO 60 IDIF=1,NL1
        IPT = NL2 - IDIF
        DDTEMP(IDIF,IPT) = DDTEMP(IDIF,1)
   60 CONTINUE
C
C          FILL BORDER OF SECOND TRIANGLE - SIZE NR
   70 CONTINUE
C         TOP BORDER
      DO 80 IDIF=1,NR
        DDTEMP(IDIF,NL2) = FRIGHT(IDIF)/FACTOR(IDIF)
   80 CONTINUE
      IF (NRL.EQ.NL2) GO TO 100
C                   BOTTOM BORDER
      DO 90 IDIF=1,NR1
        IPT = NRL + 1 - IDIF
        DDTEMP(IDIF,IPT) = DDTEMP(IDIF,NL2)
   90 CONTINUE
C
C           FILL PARALLOGRAM BETWEEN THE TWO TRIANGLES JUST FILLED
C          FILL ENTRIES PARALLEL TO BOTTOM OF FIRST TRIANGLE
  100 CONTINUE
C
C         LOOP STEPPING ALONG TOP SIDE OF SECOND TRIANGLE
      DO 120 J=2,NR2
        IDIF = J
C             LOOP STEPPING PARALLEL TO BOTTOM SIDE OF FIRST TRIANGLE
        DO 110 K=2,NL2
          IPT = NL + 2 - K
          DIFFF = DDTEMP(IDIF-1,IPT+1) - DDTEMP(IDIF-1,IPT)
          IPT2 = IPT - 1 + IDIF
          DIFFX = XDD(IPT2) - XDD(IPT)
          DDTEMP(IDIF,IPT) = DIFFF/DIFFX
          IDIF = IDIF + 1
  110   CONTINUE
  120 CONTINUE
C
C         FILL IN BOTTOM DIAGONALS FOR INTERPOLATION POINTS, IF ANY
      IF (NI.EQ.0) GO TO 150
C         LOOP THROUGH THE INTERPOLATATION POINTS
      DO 140 J=1,NI
        IDIF = 2
        NRLJ = NRL + J
        DDTEMP(1,NRLJ) = FINTRP(J)
C         LOOP THROUGH THE DIFFERENCES (IDIF INDEX)
        NRLJ1 = NRLJ - 1
        DO 130 K=1,NRLJ1
          IPT = NRLJ - K
          DIFFF = DDTEMP(IDIF-1,IPT+1) - DDTEMP(IDIF-1,IPT)
          DIFFX = XDD(NRLJ) - XDD(IPT)
          DDTEMP(IDIF,IPT) = DIFFF/DIFFX
          IDIF = IDIF + 1
  130   CONTINUE
  140 CONTINUE
  150 CONTINUE
      RETURN
      END
