      SUBROUTINE TREGUP(NR,N,X,F,G,A,FCN,SC,SX,NWTAKE,STEPMX,STEPTL,
     +     DLT,IRETCD,XPLSP,FPLSP,XPLS,FPLS,MXTAKE,IPR,METHOD,UDIAG)
C
C PURPOSE
C -------
C DECIDE WHETHER TO ACCEPT XPLS=X+SC AS THE NEXT ITERATE AND UPDATE THE
C TRUST REGION DLT.
C
C PARAMETERS
C ----------
C NR           --> ROW DIMENSION OF MATRIX
C N            --> DIMENSION OF PROBLEM
C X(N)         --> OLD ITERATE X[K-1]
C F            --> FUNCTION VALUE AT OLD ITERATE, F(X)
C G(N)         --> GRADIENT AT OLD ITERATE, G(X), OR APPROXIMATE
C A(N,N)       --> CHOLESKY DECOMPOSITION OF HESSIAN IN
C                  LOWER TRIANGULAR PART AND DIAGONAL.
C                  HESSIAN OR APPROX IN UPPER TRIANGULAR PART
C FCN          --> NAME OF SUBROUTINE TO EVALUATE FUNCTION
C SC(N)        --> CURRENT STEP
C SX(N)        --> DIAGONAL SCALING MATRIX FOR X
C NWTAKE       --> BOOLEAN, =.TRUE. IF NEWTON STEP TAKEN
C STEPMX       --> MAXIMUM ALLOWABLE STEP SIZE
C STEPTL       --> RELATIVE STEP SIZE AT WHICH SUCCESSIVE ITERATES
C                  CONSIDERED CLOSE ENOUGH TO TERMINATE ALGORITHM
C DLT         <--> TRUST REGION RADIUS
C IRETCD      <--> RETURN CODE
C                    =0 XPLS ACCEPTED AS NEXT ITERATE;
C                       DLT TRUST REGION FOR NEXT ITERATION.
C                    =1 XPLS UNSATISFACTORY BUT ACCEPTED AS NEXT ITERATE
C                       BECAUSE XPLS-X .LT. SMALLEST ALLOWABLE
C                       STEP LENGTH.
C                    =2 F(XPLS) TOO LARGE.  CONTINUE CURRENT ITERATION
C                       WITH NEW REDUCED DLT.
C                    =3 F(XPLS) SUFFICIENTLY SMALL, BUT QUADRATIC MODEL
C                       PREDICTS F(XPLS) SUFFICIENTLY WELL TO CONTINUE
C                       CURRENT ITERATION WITH NEW DOUBLED DLT.
C XPLSP(N)    <--> WORKSPACE [VALUE NEEDS TO BE RETAINED BETWEEN
C                  SUCCESIVE CALLS OF K-TH GLOBAL STEP]
C FPLSP       <--> [RETAIN VALUE BETWEEN SUCCESSIVE CALLS]
C XPLS(N)     <--  NEW ITERATE X[K]
C FPLS        <--  FUNCTION VALUE AT NEW ITERATE, F(XPLS)
C MXTAKE      <--  BOOLEAN FLAG INDICATING STEP OF MAXIMUM LENGTH USED
C IPR          --> DEVICE TO WHICH TO SEND OUTPUT
C METHOD       --> ALGORITHM TO USE TO SOLVE MINIMIZATION PROBLEM
C                    =1 LINE SEARCH
C                    =2 DOUBLE DOGLEG
C                    =3 MORE-HEBDON
C UDIAG(N)     --> DIAGONAL OF HESSIAN IN A(.,.)
C
      DIMENSION X(N),XPLS(N),G(N)
      DIMENSION SX(N),SC(N),XPLSP(N)
      DIMENSION A(NR,1)
      LOGICAL NWTAKE,MXTAKE
      DIMENSION UDIAG(N)
C
      IPR=IPR
      MXTAKE=.FALSE.
      DO 100 I=1,N
        XPLS(I)=X(I)+SC(I)
  100 CONTINUE
      CALL FCN(N,XPLS,FPLS)
      DLTF=FPLS-F
C+++
      SLP=SDOT(N,G,1,SC,1)
C
C NEXT STATEMENT ADDED FOR CASE OF COMPILERS WHICH DO NOT OPTIMIZE
C EVALUATION OF NEXT "IF" STATEMENT (IN WHICH CASE FPLSP COULD BE
C UNDEFINED).
      IF(IRETCD.EQ.4) FPLSP=0.0
C$    WRITE(IPR,961) IRETCD,FPLS,FPLSP,DLTF,SLP
      IF(IRETCD.NE.3 .OR. (FPLS.LT.FPLSP .AND. DLTF.LE. 1.E-4*SLP))
     +                                                     GO TO 130
C     IF(IRETCD.EQ.3 .AND. (FPLS.GE.FPLSP .OR. DLTF.GT. 1.E-4*SLP))
C     THEN
C
C       RESET XPLS TO XPLSP AND TERMINATE GLOBAL STEP
C
        IRETCD=0
        DO 110 I=1,N
          XPLS(I)=XPLSP(I)
  110   CONTINUE
        FPLS=FPLSP
        DLT=.5*DLT
C$      WRITE(IPR,951)
        GO TO 230
C     ELSE
C
C       FPLS TOO LARGE
C
  130   IF(DLTF.LE. 1.E-4*SLP) GO TO 170
C       IF(DLTF.GT. 1.E-4*SLP)
C       THEN
C$        WRITE(IPR,952)
          RLN=0.
          DO 140 I=1,N
            RLN=MAX(RLN,ABS(SC(I))/MAX(ABS(XPLS(I)),1./SX(I)))
  140     CONTINUE
C$        WRITE(IPR,962) RLN
          IF(RLN.GE.STEPTL) GO TO 150
C         IF(RLN.LT.STEPTL)
C         THEN
C
C           CANNOT FIND SATISFACTORY XPLS SUFFICIENTLY DISTINCT FROM X
C
            IRETCD=1
C$          WRITE(IPR,954)
            GO TO 230
C         ELSE
C
C           REDUCE TRUST REGION AND CONTINUE GLOBAL STEP
C
  150       IRETCD=2
            DLTMP=-SLP*DLT/(2.*(DLTF-SLP))
C$          WRITE(IPR,963) DLTMP
            IF(DLTMP.GE. .1*DLT) GO TO 155
C           IF(DLTMP.LT. .1*DLT)
C           THEN
              DLT=.1*DLT
              GO TO 160
C           ELSE
  155         DLT=DLTMP
C           ENDIF
  160       CONTINUE
C$          WRITE(IPR,955)
            GO TO 230
C         ENDIF
C       ELSE
C
C         FPLS SUFFICIENTLY SMALL
C
  170     CONTINUE
C$        WRITE(IPR,958)
          DLTFP=0.
          IF (METHOD .EQ. 3) GO TO 180
C
C         IF (METHOD .EQ. 2)
C         THEN
C
          DO 177 I = 1, N
             TEMP = 0.0
             DO 173 J = I, N
                TEMP = TEMP + (A(J, I)*SC(J))
  173        CONTINUE
             DLTFP = DLTFP + TEMP*TEMP
  177     CONTINUE
          GO TO 190
C
C         ELSE
C
  180     DO 187 I = 1, N
             DLTFP = DLTFP + UDIAG(I)*SC(I)*SC(I)
             IF (I .EQ. N) GO TO 187
             TEMP = 0
             IP1 = I + 1
             DO 183 J = IP1, N
                TEMP = TEMP + A(I, J)*SC(I)*SC(J)
  183        CONTINUE
             DLTFP = DLTFP + 2.0*TEMP
  187     CONTINUE
C
C         END IF
C
  190     DLTFP = SLP + DLTFP/2.0
C$        WRITE(IPR,964) DLTFP,NWTAKE
          IF(IRETCD.EQ.2 .OR. (ABS(DLTFP-DLTF).GT. .1*ABS(DLTF))
     +         .OR. NWTAKE .OR. (DLT.GT. .99*STEPMX)) GO TO 210
C         IF(IRETCD.NE.2 .AND. (ABS(DLTFP-DLTF) .LE. .1*ABS(DLTF))
C    +         .AND. (.NOT.NWTAKE) .AND. (DLT.LE. .99*STEPMX))
C         THEN
C
C           DOUBLE TRUST REGION AND CONTINUE GLOBAL STEP
C
            IRETCD=3
            DO 200 I=1,N
              XPLSP(I)=XPLS(I)
  200       CONTINUE
            FPLSP=FPLS
            DLT=MIN(2.*DLT,STEPMX)
C$          WRITE(IPR,959)
            GO TO 230
C         ELSE
C
C           ACCEPT XPLS AS NEXT ITERATE.  CHOOSE NEW TRUST REGION.
C
  210       CONTINUE
C$          WRITE(IPR,960)
            IRETCD=0
            IF(DLT.GT. .99*STEPMX) MXTAKE=.TRUE.
            IF(DLTF.LT. .1*DLTFP) GO TO 220
C           IF(DLTF.GE. .1*DLTFP)
C           THEN
C
C             DECREASE TRUST REGION FOR NEXT ITERATION
C
              DLT=.5*DLT
              GO TO 230
C           ELSE
C
C             CHECK WHETHER TO INCREASE TRUST REGION FOR NEXT ITERATION
C
  220         IF(DLTF.LE. .75*DLTFP) DLT=MIN(2.*DLT,STEPMX)
C           ENDIF
C         ENDIF
C       ENDIF
C     ENDIF
  230 CONTINUE
C$    WRITE(IPR,953)
C$    WRITE(IPR,956) IRETCD,MXTAKE,DLT,FPLS
C$    WRITE(IPR,957)
C$    WRITE(IPR,965) (XPLS(I),I=1,N)
      RETURN
C
  951 FORMAT(55H TREGUP    RESET XPLS TO XPLSP. TERMINATION GLOBAL STEP)
  952 FORMAT(26H TREGUP    FPLS TOO LARGE.)
  953 FORMAT(38H0TREGUP    VALUES AFTER CALL TO TREGUP)
  954 FORMAT(54H TREGUP    CANNOT FIND SATISFACTORY XPLS DISTINCT FROM,
     +       27H X.  TERMINATE GLOBAL STEP.)
  955 FORMAT(53H TREGUP    REDUCE TRUST REGION. CONTINUE GLOBAL STEP.)
  956 FORMAT(21H TREGUP       IRETCD=,I3/
     +       21H TREGUP       MXTAKE=,L1/
     +       21H TREGUP       DLT   =,E20.13/
     +       21H TREGUP       FPLS  =,E20.13)
  957 FORMAT(32H TREGUP       NEW ITERATE (XPLS))
  958 FORMAT(35H TREGUP    FPLS SUFFICIENTLY SMALL.)
  959 FORMAT(54H TREGUP    DOUBLE TRUST REGION.  CONTINUE GLOBAL STEP.)
  960 FORMAT(50H TREGUP    ACCEPT XPLS AS NEW ITERATE.  CHOOSE NEW,
     +       38H TRUST REGION.  TERMINATE GLOBAL STEP.)
  961 FORMAT(18H TREGUP    IRETCD=,I5/
     +       18H TREGUP    FPLS  =,E20.13/
     +       18H TREGUP    FPLSP =,E20.13/
     +       18H TREGUP    DLTF  =,E20.13/
     +       18H TREGUP    SLP   =,E20.13)
  962 FORMAT(18H TREGUP    RLN   =,E20.13)
  963 FORMAT(18H TREGUP    DLTMP =,E20.13)
  964 FORMAT(18H TREGUP    DLTFP =,E20.13/
     +       18H TREGUP    NWTAKE=,L1)
  965 FORMAT(14H TREGUP       ,5(E20.13,3X))
      END
