SUBROUTINE  G7ITB(B, D, G, IV, LIV, LV, P, PS, V, X, Y)
C
C  ***  CARRY OUT NL2SOL-LIKE ITERATIONS FOR GENERALIZED LINEAR   ***
C  ***  REGRESSION PROBLEMS (AND OTHERS OF SIMILAR STRUCTURE)     ***
C  ***  HAVING SIMPLE BOUNDS ON THE PARAMETERS BEING ESTIMATED.   ***
C
C  ***  PARAMETER DECLARATIONS  ***
C
      INTEGER LIV, LV, P, PS
      INTEGER IV(LIV)
      REAL B(2,P), D(P), G(P), V(LV), X(P), Y(P)
C
C--------------------------  PARAMETER USAGE  --------------------------
C
C B.... VECTOR OF LOWER AND UPPER BOUNDS ON X.
C D.... SCALE VECTOR.
C IV... INTEGER VALUE ARRAY.
C LIV.. LENGTH OF IV.  MUST BE AT LEAST 80.
C LH... LENGTH OF H = P*(P+1)/2.
C LV... LENGTH OF V.  MUST BE AT LEAST P*(3*P + 19)/2 + 7.
C G.... GRADIENT AT X (WHEN IV(1) = 2).
C HC... GAUSS-NEWTON HESSIAN AT X (WHEN IV(1) = 2).
C P.... NUMBER OF PARAMETERS (COMPONENTS IN X).
C PS... NUMBER OF NONZERO ROWS AND COLUMNS IN S.
C V.... FLOATING-POINT VALUE ARRAY.
C X.... PARAMETER VECTOR.
C Y.... PART OF YIELD VECTOR (WHEN IV(1)= 2, SCRATCH OTHERWISE).
C
C  ***  DISCUSSION  ***
C
C         G7ITB IS SIMILAR TO G7LIT, EXCEPT FOR THE EXTRA PARAMETER B
C     --  G7ITB ENFORCES THE BOUNDS  B(1,I) .LE. X(I) .LE. B(2,I),
C     I = 1(1)P.
C         G7ITB PERFORMS NL2SOL-LIKE ITERATIONS FOR A VARIETY OF
C     REGRESSION PROBLEMS THAT ARE SIMILAR TO NONLINEAR LEAST-SQUARES
C     IN THAT THE HESSIAN IS THE SUM OF TWO TERMS, A READILY-COMPUTED
C     FIRST-ORDER TERM AND A SECOND-ORDER TERM.  THE CALLER SUPPLIES
C     THE FIRST-ORDER TERM OF THE HESSIAN IN HC (LOWER TRIANGLE, STORED
C     COMPACTLY BY ROWS), AND  G7ITB BUILDS AN APPROXIMATION, S, TO THE
C     SECOND-ORDER TERM.  THE CALLER ALSO PROVIDES THE FUNCTION VALUE,
C     GRADIENT, AND PART OF THE YIELD VECTOR USED IN UPDATING S.
C      G7ITB DECIDES DYNAMICALLY WHETHER OR NOT TO _USE_ S WHEN CHOOSING
C     THE NEXT STEP TO TRY...  THE HESSIAN APPROXIMATION USED IS EITHER
C     HC ALONE (GAUSS-NEWTON MODEL) OR HC + S (AUGMENTED MODEL).
C     IF PS .LT. P, THEN ROWS AND COLUMNS PS+1...P OF S ARE KEPT
C     CONSTANT.  THEY WILL BE ZERO UNLESS THE CALLER SETS IV(INITS) TO
C     1 OR 2 AND SUPPLIES NONZERO VALUES FOR THEM, OR THE CALLER SETS
C     IV(INITS) TO 3 OR 4 AND THE FINITE-DIFFERENCE INITIAL S THEN
C     COMPUTED HAS NONZERO VALUES IN THESE ROWS.
C
C        IF IV(INITS) IS 3 OR 4, THEN THE INITIAL S IS COMPUTED BY
C     FINITE DIFFERENCES.  3 MEANS _USE_ FUNCTION DIFFERENCES, 4 MEANS
C     _USE_ GRADIENT DIFFERENCES.  FINITE DIFFERENCING IS DONE THE SAME
C     WAY AS IN COMPUTING A COVARIANCE MATRIX (WITH IV(COVREQ) = -1, -2,
C     1, OR 2).
C
C        FOR UPDATING S,  G7ITB ASSUMES THAT THE GRADIENT HAS THE FORM
C     OF A SUM OVER I OF RHO(I,X)*GRAD(R(I,X)), WHERE GRAD DENOTES THE
C     GRADIENT WITH RESPECT TO X.  THE TRUE SECOND-ORDER TERM THEN IS
C     THE SUM OVER I OF RHO(I,X)*HESSIAN(R(I,X)).  IF X = X0 + STEP,
C     THEN WE WISH TO UPDATE S SO THAT S*STEP IS THE SUM OVER I OF
C     RHO(I,X)*(GRAD(R(I,X)) - GRAD(R(I,X0))).  THE CALLER MUST SUPPLY
C     PART OF THIS IN Y, NAMELY THE SUM OVER I OF
C     RHO(I,X)*GRAD(R(I,X0)), WHEN CALLING  G7ITB WITH IV(1) = 2 AND
C     IV(MODE) = 0 (WHERE MODE = 38).  G THEN CONTANS THE OTHER PART,
C     SO THAT THE DESIRED YIELD VECTOR IS G - Y.  IF PS .LT. P, THEN
C     THE ABOVE DISCUSSION APPLIES ONLY TO THE FIRST PS COMPONENTS OF
C     GRAD(R(I,X)), STEP, AND Y.
C
C        PARAMETERS IV, P, V, AND X ARE THE SAME AS THE CORRESPONDING
C     ONES TO   N2GB (AND NL2SOL), EXCEPT THAT V CAN BE SHORTER
C     (SINCE THE PART OF V THAT   N2GB USES FOR STORING D, J, AND R IS
C     NOT NEEDED).  MOREOVER, COMPARED WITH   N2GB (AND NL2SOL), IV(1)
C     MAY HAVE THE TWO ADDITIONAL OUTPUT VALUES 1 AND 2, WHICH ARE
C     EXPLAINED BELOW, AS IS THE _USE_ OF IV(TOOBIG) AND IV(NFGCAL).
C     THE VALUES IV(D), IV(J), AND IV(R), WHICH ARE OUTPUT VALUES FROM
C       N2GB (AND   N2FB), ARE NOT REFERENCED BY  G7ITB OR THE
C     SUBROUTINES IT CALLS.
C
C        WHEN  G7ITB IS FIRST CALLED, I.E., WHEN  G7ITB IS CALLED WITH
C     IV(1) = 0 OR 12, V(F), G, AND HC NEED NOT BE INITIALIZED.  TO
C     OBTAIN THESE STARTING VALUES,  G7ITB RETURNS FIRST WITH IV(1) = 1,
C     THEN WITH IV(1) = 2, WITH IV(MODE) = -1 IN BOTH CASES.  ON
C     SUBSEQUENT RETURNS WITH IV(1) = 2, IV(MODE) = 0 IMPLIES THAT
C     Y MUST ALSO BE SUPPLIED.  (NOTE THAT Y IS USED FOR SCRATCH -- ITS
C     INPUT CONTENTS ARE LOST.  BY CONTRAST, HC IS NEVER CHANGED.)
C     ONCE CONVERGENCE HAS BEEN OBTAINED, IV(RDREQ) AND IV(COVREQ) MAY
C     IMPLY THAT A FINITE-DIFFERENCE HESSIAN SHOULD BE COMPUTED FOR USE
C     IN COMPUTING A COVARIANCE MATRIX.  IN THIS CASE  G7ITB WILL MAKE
C     A NUMBER OF RETURNS WITH IV(1) = 1 OR 2 AND IV(MODE) POSITIVE.
C     WHEN IV(MODE) IS POSITIVE, Y SHOULD NOT BE CHANGED.
C
C IV(1) = 1 MEANS THE CALLER SHOULD SET V(F) (I.E., V(10)) TO F(X), THE
C             FUNCTION VALUE AT X, AND CALL  G7ITB AGAIN, HAVING CHANGED
C             NONE OF THE OTHER PARAMETERS.  AN EXCEPTION OCCURS IF F(X)
C             CANNOT BE EVALUATED (E.G. IF OVERFLOW WOULD OCCUR), WHICH
C             MAY HAPPEN BECAUSE OF AN OVERSIZED STEP.  IN THIS CASE
C             THE CALLER SHOULD SET IV(TOOBIG) = IV(2) TO 1, WHICH WILL
C             CAUSE  G7ITB TO IGNORE V(F) AND TRY A SMALLER STEP.  NOTE
C             THAT THE CURRENT FUNCTION EVALUATION COUNT IS AVAILABLE
C             IN IV(NFCALL) = IV(6).  THIS MAY BE USED TO IDENTIFY
C             WHICH COPY OF SAVED INFORMATION SHOULD BE USED IN COM-
C             PUTING G, HC, AND Y THE NEXT TIME  G7ITB RETURNS WITH
C             IV(1) = 2.  SEE MLPIT FOR AN EXAMPLE OF THIS.
C IV(1) = 2 MEANS THE CALLER SHOULD SET G TO G(X), THE GRADIENT OF F AT
C             X.  THE CALLER SHOULD ALSO SET HC TO THE GAUSS-NEWTON
C             HESSIAN AT X.  IF IV(MODE) = 0, THEN THE CALLER SHOULD
C             ALSO COMPUTE THE PART OF THE YIELD VECTOR DESCRIBED ABOVE.
C             THE CALLER SHOULD THEN CALL  G7ITB AGAIN (WITH IV(1) = 2).
C             THE CALLER MAY ALSO CHANGE D AT THIS TIME, BUT SHOULD NOT
C             CHANGE X.  NOTE THAT IV(NFGCAL) = IV(7) CONTAINS THE
C             VALUE THAT IV(NFCALL) HAD DURING THE RETURN WITH
C             IV(1) = 1 IN WHICH X HAD THE SAME VALUE AS IT NOW HAS.
C             IV(NFGCAL) IS EITHER IV(NFCALL) OR IV(NFCALL) - 1.  MLPIT
C             IS AN EXAMPLE WHERE THIS INFORMATION IS USED.  IF G OR HC
C             CANNOT BE EVALUATED AT X, THEN THE CALLER MAY SET
C             IV(NFGCAL) TO 0, IN WHICH CASE  G7ITB WILL RETURN WITH
C             IV(1) = 15.
C
C  ***  GENERAL  ***
C
C     CODED BY DAVID M. GAY.
C
C        (SEE NL2SOL FOR REFERENCES.)
C
C+++++++++++++++++++++++++++  DECLARATIONS  ++++++++++++++++++++++++++++
C
C  ***  LOCAL VARIABLES  ***
C
      LOGICAL HAVQTR, HAVRM
      INTEGER DUMMY, DIG1, G01, H1, HC1, I, I1, IPI, IPIV0, IPIV1,
     1        IPIV2, IPN, J, K, L, LMAT1, LSTGST, P1, P1LEN, PP1, PP1O2,
     2        QTR1, RMAT1, RSTRST, STEP1, STPMOD, S1, TD1, TEMP1, TEMP2,
     3        TG1, W1, WLM1, X01
      REAL E, GI, STTSST, T, T1, XI
C
C     ***  CONSTANTS  ***
C
      REAL HALF, NEGONE, ONE, ONEP2, ZERO
C
C  ***  EXTERNAL FUNCTIONS AND SUBROUTINES  ***
C
      LOGICAL STOPX
      REAL  D7TPR,  RLDST,  V2NRM
      EXTERNAL A7SST,  D7TPR,  F7DHB,  G7QSB,I7COPY, I7PNVR, I7SHFT,
     1         ITSUM,  L7MSB,  L7SQR,  L7TVM, L7VML, PARCK,  Q7RSH,
     2          RLDST,  S7DMP,  S7IPR,  S7LUP,  S7LVM, STOPX,  V2NRM,
     3         V2AXY, V7CPY,  V7IPR,  V7SCP,  V7VMP
C
C A7SST.... ASSESSES CANDIDATE STEP.
C  D7TPR... RETURNS INNER PRODUCT OF TWO VECTORS.
C  F7DHB... COMPUTE FINITE-DIFFERENCE HESSIAN (FOR INIT. S MATRIX).
C  G7QSB... COMPUTES GOLDFELD-QUANDT-TROTTER STEP (AUGMENTED MODEL).
C I7COPY.... COPIES ONE INTEGER VECTOR TO ANOTHER.
C I7PNVR... INVERTS PERMUTATION ARRAY.
C I7SHFT... SHIFTS AN INTEGER VECTOR.
C ITSUM.... PRINTS ITERATION SUMMARY AND INFO ON INITIAL AND FINAL X.
C  L7MSB... COMPUTES LEVENBERG-MARQUARDT STEP (GAUSS-NEWTON MODEL).
C  L7SQR... COMPUTES L * L**T FROM LOWER TRIANGULAR MATRIX L.
C  L7TVM... COMPUTES L**T * V, V = VECTOR, L = LOWER TRIANGULAR MATRIX.
C L7VML.... COMPUTES L * V, V = VECTOR, L = LOWER TRIANGULAR MATRIX.
C PARCK.... CHECK VALIDITY OF IV AND V INPUT COMPONENTS.
C  Q7RSH... SHIFTS A QR FACTORIZATION.
C  RLDST... COMPUTES V(RELDX) = RELATIVE STEP SIZE.
C  S7DMP... MULTIPLIES A SYM. MATRIX FORE AND AFT BY A DIAG. MATRIX.
C  S7IPR... APPLIES PERMUTATION TO (LOWER TRIANG. OF) SYM. MATRIX.
C  S7LUP... PERFORMS QUASI-NEWTON UPDATE ON COMPACTLY STORED LOWER TRI-
C             ANGLE OF A SYMMETRIC MATRIX.
C  S7LVM... MULTIPLIES COMPACTLY STORED SYM. MATRIX TIMES VECTOR.
C STOPX.... RETURNS .TRUE. IF THE BREAK KEY HAS BEEN PRESSED.
C  V2NRM... RETURNS THE 2-NORM OF A VECTOR.
C V2AXY.... COMPUTES SCALAR TIMES ONE VECTOR PLUS ANOTHER.
C V7CPY.... COPIES ONE VECTOR TO ANOTHER.
C  V7IPR... APPLIES A PERMUTATION TO A VECTOR.
C  V7SCP... SETS ALL ELEMENTS OF A VECTOR TO A SCALAR.
C  V7VMP... MULTIPLIES (DIVIDES) VECTORS COMPONENTWISE.
C
C  ***  SUBSCRIPTS FOR IV AND V  ***
C
      INTEGER CNVCOD, COSMIN, COVMAT, COVREQ, DGNORM, DIG,
     1        DSTNRM, F, FDH, FDIF, FUZZ, F0, GTSTEP, H, HC, IERR,
     2        INCFAC, INITS, IPIVOT, IRC, IVNEED, KAGQT, KALM, LMAT,
     3        LMAX0, LMAXS, MODE, MODEL, MXFCAL, MXITER, NEXTIV, NEXTV,
     4        NFCALL, NFGCAL, NFCOV, NGCOV, NGCALL, NITER, NVSAVE, P0,
     5        PC, PERM, PHMXFC, PREDUC, QTR, RADFAC, RADINC, RADIUS,
     6        RAD0, RDREQ, REGD, RELDX, RESTOR, RMAT, S, SIZE, STEP,
     7        STGLIM, STPPAR, SUSED, SWITCH, TOOBIG, TUNER4, TUNER5,
     8        VNEED, VSAVE, W, WSCALE, XIRC, X0
C
C  ***  IV SUBSCRIPT VALUES  ***
C
C  ***  (NOTE THAT P0 AND PC ARE STORED IN IV(G0) AND IV(STLSTG) RESP.)
C
C/6
C     DATA CNVCOD/55/, COVMAT/26/, COVREQ/15/, DIG/37/, FDH/74/, H/56/,
C    1     HC/71/, IERR/75/, INITS/25/, IPIVOT/76/, IRC/29/, IVNEED/3/,
C    2     KAGQT/33/, KALM/34/, LMAT/42/, MODE/35/, MODEL/5/,
C    3     MXFCAL/17/, MXITER/18/, NEXTIV/46/, NEXTV/47/, NFCALL/6/,
C    4     NFGCAL/7/, NFCOV/52/, NGCOV/53/, NGCALL/30/, NITER/31/,
C    5     P0/48/, PC/41/, PERM/58/, QTR/77/, RADINC/8/, RDREQ/57/,
C    6     REGD/67/, RESTOR/9/, RMAT/78/, S/62/, STEP/40/, STGLIM/11/,
C    7     SUSED/64/, SWITCH/12/, TOOBIG/2/, VNEED/4/, VSAVE/60/, W/65/,
C    8     XIRC/13/, X0/43/
C/7
      PARAMETER (CNVCOD=55, COVMAT=26, COVREQ=15, DIG=37, FDH=74, H=56,
     1           HC=71, IERR=75, INITS=25, IPIVOT=76, IRC=29, IVNEED=3,
     2           KAGQT=33, KALM=34, LMAT=42, MODE=35, MODEL=5,
     3           MXFCAL=17, MXITER=18, NEXTIV=46, NEXTV=47, NFCALL=6,
     4           NFGCAL=7, NFCOV=52, NGCOV=53, NGCALL=30, NITER=31,
     5           P0=48, PC=41, PERM=58, QTR=77, RADINC=8, RDREQ=57,
     6           REGD=67, RESTOR=9, RMAT=78, S=62, STEP=40, STGLIM=11,
     7           SUSED=64, SWITCH=12, TOOBIG=2, VNEED=4, VSAVE=60, W=65,
     8           XIRC=13, X0=43)
C/
C
C  ***  V SUBSCRIPT VALUES  ***
C
C/6
C     DATA COSMIN/47/, DGNORM/1/, DSTNRM/2/, F/10/, FDIF/11/, FUZZ/45/,
C    1     F0/13/, GTSTEP/4/, INCFAC/23/, LMAX0/35/, LMAXS/36/,
C    2     NVSAVE/9/, PHMXFC/21/, PREDUC/7/, RADFAC/16/, RADIUS/8/,
C    3     RAD0/9/, RELDX/17/, SIZE/55/, STPPAR/5/, TUNER4/29/,
C    4     TUNER5/30/, WSCALE/56/
C/7
      PARAMETER (COSMIN=47, DGNORM=1, DSTNRM=2, F=10, FDIF=11, FUZZ=45,
     1           F0=13, GTSTEP=4, INCFAC=23, LMAX0=35, LMAXS=36,
     2           NVSAVE=9, PHMXFC=21, PREDUC=7, RADFAC=16, RADIUS=8,
     3           RAD0=9, RELDX=17, SIZE=55, STPPAR=5, TUNER4=29,
     4           TUNER5=30, WSCALE=56)
C/
C
C
C/6
C     DATA HALF/0.5E+0/, NEGONE/-1.E+0/, ONE/1.E+0/, ONEP2/1.2E+0/,
C    1     ZERO/0.E+0/
C/7
      PARAMETER (HALF=0.5E+0, NEGONE=-1.E+0, ONE=1.E+0, ONEP2=1.2E+0,
     1           ZERO=0.E+0)
C/
C
C+++++++++++++++++++++++++++++++  BODY  ++++++++++++++++++++++++++++++++
C
      I = IV(1)
      IF (I .EQ. 1) GO TO 50
      IF (I .EQ. 2) GO TO 60
C
      IF (I .LT. 12) GO TO 10
      IF (I .GT. 13) GO TO 10
         IV(VNEED) = IV(VNEED) + P*(3*P + 25)/2 + 7
         IV(IVNEED) = IV(IVNEED) + 4*P
 10   CALL PARCK(1, D, IV, LIV, LV, P, V)
      I = IV(1) - 2
      IF (I .GT. 12) GO TO 999
      GO TO (360, 360, 360, 360, 360, 360, 240, 190, 240, 20, 20, 30), I
C
C  ***  STORAGE ALLOCATION  ***
C
 20   PP1O2 = P * (P + 1) / 2
      IV(S) = IV(LMAT) + PP1O2
      IV(X0) = IV(S) + PP1O2
      IV(STEP) = IV(X0) + 2*P
      IV(DIG) = IV(STEP) + 3*P
      IV(W) = IV(DIG) + 2*P
      IV(H) = IV(W) + 4*P + 7
      IV(NEXTV) = IV(H) + PP1O2
      IV(IPIVOT) = IV(PERM) + 3*P
      IV(NEXTIV) = IV(IPIVOT) + P
      IF (IV(1) .NE. 13) GO TO 30
         IV(1) = 14
         GO TO 999
C
C  ***  INITIALIZATION  ***
C
 30   IV(NITER) = 0
      IV(NFCALL) = 1
      IV(NGCALL) = 1
      IV(NFGCAL) = 1
      IV(MODE) = -1
      IV(STGLIM) = 2
      IV(TOOBIG) = 0
      IV(CNVCOD) = 0
      IV(COVMAT) = 0
      IV(NFCOV) = 0
      IV(NGCOV) = 0
      IV(RADINC) = 0
      IV(PC) = P
      V(RAD0) = ZERO
      V(STPPAR) = ZERO
      V(RADIUS) = V(LMAX0) / (ONE + V(PHMXFC))
C
C  ***  CHECK CONSISTENCY OF B AND INITIALIZE IP ARRAY  ***
C
      IPI = IV(IPIVOT)
      DO 40 I = 1, P
         IV(IPI) = I
         IPI = IPI + 1
         IF (B(1,I) .GT. B(2,I)) GO TO 680
 40      CONTINUE
C
C  ***  SET INITIAL MODEL AND S MATRIX  ***
C
      IV(MODEL) = 1
      IV(1) = 1
      IF (IV(S) .LT. 0) GO TO 710
      IF (IV(INITS) .GT. 1) IV(MODEL) = 2
      S1 = IV(S)
      IF (IV(INITS) .EQ. 0 .OR. IV(INITS) .GT. 2)
     1   CALL  V7SCP(P*(P+1)/2, V(S1), ZERO)
      GO TO 710
C
C  ***  NEW FUNCTION VALUE  ***
C
 50   IF (IV(MODE) .EQ. 0) GO TO 360
      IF (IV(MODE) .GT. 0) GO TO 590
C
      IF (IV(TOOBIG) .EQ. 0) GO TO 690
         IV(1) = 63
         GO TO 999
C
C  ***  MAKE SURE GRADIENT COULD BE COMPUTED  ***
C
 60   IF (IV(TOOBIG) .EQ. 0) GO TO 70
         IV(1) = 65
         GO TO 999
C
C  ***  NEW GRADIENT  ***
C
 70   IV(KALM) = -1
      IV(KAGQT) = -1
      IV(FDH) = 0
      IF (IV(MODE) .GT. 0) GO TO 590
      IF (IV(HC) .LE. 0 .AND. IV(RMAT) .LE. 0) GO TO 670
C
C  ***  CHOOSE INITIAL PERMUTATION  ***
C
      IPI = IV(IPIVOT)
      IPN = IPI + P - 1
      IPIV2 = IV(PERM) - 1
      K = IV(PC)
      P1 = P
      PP1 = P + 1
      RMAT1 = IV(RMAT)
      HAVRM = RMAT1 .GT. 0
      QTR1 = IV(QTR)
      HAVQTR = QTR1 .GT. 0
C     *** MAKE SURE V(QTR1) IS LEGAL (EVEN WHEN NOT REFERENCED) ***
      W1 = IV(W)
      IF (.NOT. HAVQTR) QTR1 = W1 + P
C
      DO 100 I = 1, P
         I1 = IV(IPN)
         IPN = IPN - 1
         IF (B(1,I1) .GE. B(2,I1)) GO TO 80
         XI = X(I1)
         GI = G(I1)
         IF (XI .LE. B(1,I1) .AND. GI .GT. ZERO) GO TO 80
         IF (XI .GE. B(2,I1) .AND. GI .LT. ZERO) GO TO 80
C           *** DISALLOW CONVERGENCE IF X(I1) HAS JUST BEEN FREED ***
            J = IPIV2 + I1
            IF (IV(J) .GT. K) IV(CNVCOD) = 0
            GO TO 100
 80      IF (I1 .GE. P1) GO TO 90
            I1 = PP1 - I
            CALL I7SHFT(P1, I1, IV(IPI))
            IF (HAVRM)
     1          CALL  Q7RSH(I1, P1, HAVQTR, V(QTR1), V(RMAT1), V(W1))
 90      P1 = P1 - 1
 100     CONTINUE
      IV(PC) = P1
C
C  ***  COMPUTE V(DGNORM) (AN OUTPUT VALUE IF WE STOP NOW)  ***
C
      V(DGNORM) = ZERO
      IF (P1 .LE. 0) GO TO 110
      DIG1 = IV(DIG)
      CALL  V7VMP(P, V(DIG1), G, D, -1)
      CALL  V7IPR(P, IV(IPI), V(DIG1))
      V(DGNORM) =  V2NRM(P1, V(DIG1))
 110  IF (IV(CNVCOD) .NE. 0) GO TO 580
      IF (IV(MODE) .EQ. 0) GO TO 510
      IV(MODE) = 0
      V(F0) = V(F)
      IF (IV(INITS) .LE. 2) GO TO 170
C
C  ***  ARRANGE FOR FINITE-DIFFERENCE INITIAL S  ***
C
      IV(XIRC) = IV(COVREQ)
      IV(COVREQ) = -1
      IF (IV(INITS) .GT. 3) IV(COVREQ) = 1
      IV(CNVCOD) = 70
      GO TO 600
C
C  ***  COME TO NEXT STMT AFTER COMPUTING F.D. HESSIAN FOR INIT. S  ***
C
 120  H1 = IV(FDH)
      IF (H1 .LE. 0) GO TO 660
      IV(CNVCOD) = 0
      IV(MODE) = 0
      IV(NFCOV) = 0
      IV(NGCOV) = 0
      IV(COVREQ) = IV(XIRC)
      S1 = IV(S)
      PP1O2 = PS * (PS + 1) / 2
      HC1 = IV(HC)
      IF (HC1 .LE. 0) GO TO 130
         CALL V2AXY(PP1O2, V(S1), NEGONE, V(HC1), V(H1))
         GO TO 140
 130  RMAT1 = IV(RMAT)
      LMAT1 = IV(LMAT)
      CALL  L7SQR(P, V(LMAT1), V(RMAT1))
      IPI = IV(IPIVOT)
      IPIV1 = IV(PERM) + P
      CALL I7PNVR(P, IV(IPIV1), IV(IPI))
      CALL  S7IPR(P, IV(IPIV1), V(LMAT1))
      CALL V2AXY(PP1O2, V(S1), NEGONE, V(LMAT1), V(H1))
C
C     *** ZERO PORTION OF S CORRESPONDING TO FIXED X COMPONENTS ***
C
 140  DO 160 I = 1, P
         IF (B(1,I) .LT. B(2,I)) GO TO 160
         K = S1 + I*(I-1)/2
         CALL  V7SCP(I, V(K), ZERO)
         IF (I .GE. P) GO TO 170
         K = K + 2*I - 1
         I1 = I + 1
         DO 150 J = I1, P
            V(K) = ZERO
            K = K + J
 150        CONTINUE
 160     CONTINUE
C
 170  IV(1) = 2
C
C
C-----------------------------  MAIN LOOP  -----------------------------
C
C
C  ***  PRINT ITERATION SUMMARY, CHECK ITERATION LIMIT  ***
C
 180  CALL ITSUM(D, G, IV, LIV, LV, P, V, X)
 190  K = IV(NITER)
      IF (K .LT. IV(MXITER)) GO TO 200
         IV(1) = 10
         GO TO 999
 200  IV(NITER) = K + 1
C
C  ***  UPDATE RADIUS  ***
C
      IF (K .EQ. 0) GO TO 220
      STEP1 = IV(STEP)
      DO 210 I = 1, P
         V(STEP1) = D(I) * V(STEP1)
         STEP1 = STEP1 + 1
 210     CONTINUE
      STEP1 = IV(STEP)
      T = V(RADFAC) *  V2NRM(P, V(STEP1))
      IF (V(RADFAC) .LT. ONE .OR. T .GT. V(RADIUS)) V(RADIUS) = T
C
C  ***  INITIALIZE FOR START OF NEXT ITERATION  ***
C
 220  X01 = IV(X0)
      V(F0) = V(F)
      IV(IRC) = 4
      IV(H) = -IABS(IV(H))
      IV(SUSED) = IV(MODEL)
C
C     ***  COPY X TO X0  ***
C
      CALL V7CPY(P, V(X01), X)
C
C  ***  CHECK STOPX AND FUNCTION EVALUATION LIMIT  ***
C
 230  IF (.NOT. STOPX(DUMMY)) GO TO 250
         IV(1) = 11
         GO TO 260
C
C     ***  COME HERE WHEN RESTARTING AFTER FUNC. EVAL. LIMIT OR STOPX.
C
 240  IF (V(F) .GE. V(F0)) GO TO 250
         V(RADFAC) = ONE
         K = IV(NITER)
         GO TO 200
C
 250  IF (IV(NFCALL) .LT. IV(MXFCAL) + IV(NFCOV)) GO TO 270
         IV(1) = 9
 260     IF (V(F) .GE. V(F0)) GO TO 999
C
C        ***  IN CASE OF STOPX OR FUNCTION EVALUATION LIMIT WITH
C        ***  IMPROVED V(F), EVALUATE THE GRADIENT AT X.
C
              IV(CNVCOD) = IV(1)
              GO TO 500
C
C. . . . . . . . . . . . .  COMPUTE CANDIDATE STEP  . . . . . . . . . .
C
 270  STEP1 = IV(STEP)
      TG1 = IV(DIG)
      TD1 = TG1 + P
      X01 = IV(X0)
      W1 = IV(W)
      H1 = IV(H)
      P1 = IV(PC)
      IPI = IV(PERM)
      IPIV1 = IPI + P
      IPIV2 = IPIV1 + P
      IPIV0 = IV(IPIVOT)
      IF (IV(MODEL) .EQ. 2) GO TO 280
C
C        ***  COMPUTE LEVENBERG-MARQUARDT STEP IF POSSIBLE...
C
         RMAT1 = IV(RMAT)
         IF (RMAT1 .LE. 0) GO TO 280
         QTR1 = IV(QTR)
         IF (QTR1 .LE. 0) GO TO 280
         LMAT1 = IV(LMAT)
         WLM1 = W1 + P
         CALL  L7MSB(B, D, G, IV(IERR), IV(IPIV0), IV(IPIV1),
     1               IV(IPIV2), IV(KALM), V(LMAT1), LV, P, IV(P0),
     2               IV(PC), V(QTR1), V(RMAT1), V(STEP1), V(TD1),
     3               V(TG1), V, V(W1), V(WLM1), X, V(X01))
C        *** H IS STORED IN THE END OF W AND HAS JUST BEEN OVERWRITTEN,
C        *** SO WE MARK IT INVALID...
         IV(H) = -IABS(H1)
C        *** EVEN IF H WERE STORED ELSEWHERE, IT WOULD BE NECESSARY TO
C        *** MARK INVALID THE INFORMATION G7QTS MAY HAVE STORED IN V...
         IV(KAGQT) = -1
         GO TO 330
C
 280  IF (H1 .GT. 0) GO TO 320
C
C     ***  SET H TO  D**-1 * (HC + T1*S) * D**-1.  ***
C
         P1LEN = P1*(P1+1)/2
         H1 = -H1
         IV(H) = H1
         IV(FDH) = 0
         IF (P1 .LE. 0) GO TO 320
C        *** MAKE TEMPORARY PERMUTATION ARRAY ***
         CALL I7COPY(P, IV(IPI), IV(IPIV0))
         J = IV(HC)
         IF (J .GT. 0) GO TO 290
            J = H1
            RMAT1 = IV(RMAT)
            CALL  L7SQR(P1, V(H1), V(RMAT1))
            GO TO 300
 290     CALL V7CPY(P*(P+1)/2, V(H1), V(J))
         CALL  S7IPR(P, IV(IPI), V(H1))
 300     IF (IV(MODEL) .EQ. 1) GO TO 310
            LMAT1 = IV(LMAT)
            S1 = IV(S)
            CALL V7CPY(P*(P+1)/2, V(LMAT1), V(S1))
            CALL  S7IPR(P, IV(IPI), V(LMAT1))
            CALL V2AXY(P1LEN, V(H1), ONE, V(LMAT1), V(H1))
 310     CALL V7CPY(P, V(TD1), D)
         CALL  V7IPR(P, IV(IPI), V(TD1))
         CALL  S7DMP(P1, V(H1), V(H1), V(TD1), -1)
         IV(KAGQT) = -1
C
C  ***  COMPUTE ACTUAL GOLDFELD-QUANDT-TROTTER STEP  ***
C
 320  LMAT1 = IV(LMAT)
      CALL  G7QSB(B, D, V(H1), G, IV(IPI), IV(IPIV1), IV(IPIV2),
     1            IV(KAGQT), V(LMAT1), LV, P, IV(P0), P1, V(STEP1),
     2            V(TD1), V(TG1), V, V(W1), X, V(X01))
      IF (IV(KALM) .GT. 0) IV(KALM) = 0
C
 330  IF (IV(IRC) .NE. 6) GO TO 340
         IF (IV(RESTOR) .NE. 2) GO TO 360
         RSTRST = 2
         GO TO 370
C
C  ***  CHECK WHETHER EVALUATING F(X0 + STEP) LOOKS WORTHWHILE  ***
C
 340  IV(TOOBIG) = 0
      IF (V(DSTNRM) .LE. ZERO) GO TO 360
      IF (IV(IRC) .NE. 5) GO TO 350
      IF (V(RADFAC) .LE. ONE) GO TO 350
      IF (V(PREDUC) .GT. ONEP2 * V(FDIF)) GO TO 350
         STEP1 = IV(STEP)
         X01 = IV(X0)
         CALL V2AXY(P, V(STEP1), NEGONE, V(X01), X)
         IF (IV(RESTOR) .NE. 2) GO TO 360
         RSTRST = 0
         GO TO 370
C
C  ***  COMPUTE F(X0 + STEP)  ***
C
 350  X01 = IV(X0)
      STEP1 = IV(STEP)
      CALL V2AXY(P, X, ONE, V(STEP1), V(X01))
      IV(NFCALL) = IV(NFCALL) + 1
      IV(1) = 1
      GO TO 710
C
C. . . . . . . . . . . . .  ASSESS CANDIDATE STEP  . . . . . . . . . . .
C
 360  RSTRST = 3
 370  X01 = IV(X0)
      V(RELDX) =  RLDST(P, D, X, V(X01))
      CALL A7SST(IV, LIV, LV, V)
      STEP1 = IV(STEP)
      LSTGST = X01 + P
      I = IV(RESTOR) + 1
      GO TO (410, 380, 390, 400), I
 380  CALL V7CPY(P, X, V(X01))
      GO TO 410
 390   CALL V7CPY(P, V(LSTGST), V(STEP1))
       GO TO 410
 400     CALL V7CPY(P, V(STEP1), V(LSTGST))
         CALL V2AXY(P, X, ONE, V(STEP1), V(X01))
         V(RELDX) =  RLDST(P, D, X, V(X01))
         IV(RESTOR) = RSTRST
C
C  ***  IF NECESSARY, SWITCH MODELS  ***
C
 410  IF (IV(SWITCH) .EQ. 0) GO TO 420
         IV(H) = -IABS(IV(H))
         IV(SUSED) = IV(SUSED) + 2
         L = IV(VSAVE)
         CALL V7CPY(NVSAVE, V, V(L))
 420  L = IV(IRC) - 4
      STPMOD = IV(MODEL)
      IF (L .GT. 0) GO TO (440,450,460,460,460,460,460,460,570,510), L
C
C  ***  DECIDE WHETHER TO CHANGE MODELS  ***
C
      E = V(PREDUC) - V(FDIF)
      S1 = IV(S)
      CALL  S7LVM(PS, Y, V(S1), V(STEP1))
      STTSST = HALF *  D7TPR(PS, V(STEP1), Y)
      IF (IV(MODEL) .EQ. 1) STTSST = -STTSST
      IF ( ABS(E + STTSST) * V(FUZZ) .GE.  ABS(E)) GO TO 430
C
C     ***  SWITCH MODELS  ***
C
         IV(MODEL) = 3 - IV(MODEL)
         IF (-2 .LT. L) GO TO 470
              IV(H) = -IABS(IV(H))
              IV(SUSED) = IV(SUSED) + 2
              L = IV(VSAVE)
              CALL V7CPY(NVSAVE, V(L), V)
              GO TO 230
C
 430  IF (-3 .LT. L) GO TO 470
C
C     ***  RECOMPUTE STEP WITH DIFFERENT RADIUS  ***
C
 440  V(RADIUS) = V(RADFAC) * V(DSTNRM)
      GO TO 230
C
C  ***  COMPUTE STEP OF LENGTH V(LMAXS) FOR SINGULAR CONVERGENCE TEST
C
 450  V(RADIUS) = V(LMAXS)
      GO TO 270
C
C  ***  CONVERGENCE OR FALSE CONVERGENCE  ***
C
 460  IV(CNVCOD) = L
      IF (V(F) .GE. V(F0)) GO TO 580
         IF (IV(XIRC) .EQ. 14) GO TO 580
              IV(XIRC) = 14
C
C. . . . . . . . . . . .  PROCESS ACCEPTABLE STEP  . . . . . . . . . . .
C
 470  IV(COVMAT) = 0
      IV(REGD) = 0
C
C  ***  SEE WHETHER TO SET V(RADFAC) BY GRADIENT TESTS  ***
C
      IF (IV(IRC) .NE. 3) GO TO 500
         STEP1 = IV(STEP)
         TEMP1 = STEP1 + P
         TEMP2 = IV(X0)
C
C     ***  SET  TEMP1 = HESSIAN * STEP  FOR _USE_ IN GRADIENT TESTS  ***
C
         HC1 = IV(HC)
         IF (HC1 .LE. 0) GO TO 480
              CALL  S7LVM(P, V(TEMP1), V(HC1), V(STEP1))
              GO TO 490
 480     RMAT1 = IV(RMAT)
         IPIV0 = IV(IPIVOT)
         CALL V7CPY(P, V(TEMP1), V(STEP1))
         CALL  V7IPR(P, IV(IPIV0), V(TEMP1))
         CALL  L7TVM(P, V(TEMP1), V(RMAT1), V(TEMP1))
         CALL L7VML(P, V(TEMP1), V(RMAT1), V(TEMP1))
         IPIV1 = IV(PERM) + P
         CALL I7PNVR(P, IV(IPIV1), IV(IPIV0))
         CALL  V7IPR(P, IV(IPIV1), V(TEMP1))
C
 490     IF (STPMOD .EQ. 1) GO TO 500
              S1 = IV(S)
              CALL  S7LVM(PS, V(TEMP2), V(S1), V(STEP1))
              CALL V2AXY(PS, V(TEMP1), ONE, V(TEMP2), V(TEMP1))
C
C  ***  SAVE OLD GRADIENT AND COMPUTE NEW ONE  ***
C
 500  IV(NGCALL) = IV(NGCALL) + 1
      G01 = IV(W)
      CALL V7CPY(P, V(G01), G)
      GO TO 690
C
C  ***  INITIALIZATIONS -- G0 = G - G0, ETC.  ***
C
 510  G01 = IV(W)
      CALL V2AXY(P, V(G01), NEGONE, V(G01), G)
      STEP1 = IV(STEP)
      TEMP1 = STEP1 + P
      TEMP2 = IV(X0)
      IF (IV(IRC) .NE. 3) GO TO 540
C
C  ***  SET V(RADFAC) BY GRADIENT TESTS  ***
C
C     ***  SET  TEMP1 = D**-1 * (HESSIAN * STEP  +  (G(X0) - G(X)))  ***
C
         K = TEMP1
         L = G01
         DO 520 I = 1, P
              V(K) = (V(K) - V(L)) / D(I)
              K = K + 1
              L = L + 1
 520          CONTINUE
C
C        ***  DO GRADIENT TESTS  ***
C
         IF ( V2NRM(P, V(TEMP1)) .LE. V(DGNORM) * V(TUNER4))  GO TO 530
              IF ( D7TPR(P, G, V(STEP1))
     1                  .GE. V(GTSTEP) * V(TUNER5))  GO TO 540
 530               V(RADFAC) = V(INCFAC)
C
C  ***  COMPUTE Y VECTOR NEEDED FOR UPDATING S  ***
C
 540  CALL V2AXY(PS, Y, NEGONE, Y, G)
C
C  ***  DETERMINE SIZING FACTOR V(SIZE)  ***
C
C     ***  SET TEMP1 = S * STEP  ***
      S1 = IV(S)
      CALL  S7LVM(PS, V(TEMP1), V(S1), V(STEP1))
C
      T1 =  ABS( D7TPR(PS, V(STEP1), V(TEMP1)))
      T =  ABS( D7TPR(PS, V(STEP1), Y))
      V(SIZE) = ONE
      IF (T .LT. T1) V(SIZE) = T / T1
C
C  ***  SET G0 TO WCHMTD CHOICE OF FLETCHER AND AL-BAALI  ***
C
      HC1 = IV(HC)
      IF (HC1 .LE. 0) GO TO 550
         CALL  S7LVM(PS, V(G01), V(HC1), V(STEP1))
         GO TO 560
C
 550  RMAT1 = IV(RMAT)
      IPIV0 = IV(IPIVOT)
      CALL V7CPY(P, V(G01), V(STEP1))
      I = G01 + PS
      IF (PS .LT. P) CALL  V7SCP(P-PS, V(I), ZERO)
      CALL  V7IPR(P, IV(IPIV0), V(G01))
      CALL  L7TVM(P, V(G01), V(RMAT1), V(G01))
      CALL L7VML(P, V(G01), V(RMAT1), V(G01))
      IPIV1 = IV(PERM) + P
      CALL I7PNVR(P, IV(IPIV1), IV(IPIV0))
      CALL  V7IPR(P, IV(IPIV1), V(G01))
C
 560  CALL V2AXY(PS, V(G01), ONE, Y, V(G01))
C
C  ***  UPDATE S  ***
C
      CALL  S7LUP(V(S1), V(COSMIN), PS, V(SIZE), V(STEP1), V(TEMP1),
     1            V(TEMP2), V(G01), V(WSCALE), Y)
      IV(1) = 2
      GO TO 180
C
C. . . . . . . . . . . . . .  MISC. DETAILS  . . . . . . . . . . . . . .
C
C  ***  BAD PARAMETERS TO ASSESS  ***
C
 570  IV(1) = 64
      GO TO 999
C
C
C  ***  CONVERGENCE OBTAINED -- SEE WHETHER TO COMPUTE COVARIANCE  ***
C
 580  IF (IV(RDREQ) .EQ. 0) GO TO 660
      IF (IV(FDH) .NE. 0) GO TO 660
      IF (IV(CNVCOD) .GE. 7) GO TO 660
      IF (IV(REGD) .GT. 0) GO TO 660
      IF (IV(COVMAT) .GT. 0) GO TO 660
      IF (IABS(IV(COVREQ)) .GE. 3) GO TO 640
      IF (IV(RESTOR) .EQ. 0) IV(RESTOR) = 2
      GO TO 600
C
C  ***  COMPUTE FINITE-DIFFERENCE HESSIAN FOR COMPUTING COVARIANCE  ***
C
 590  IV(RESTOR) = 0
 600  CALL  F7DHB(B, D, G, I, IV, LIV, LV, P, V, X)
      GO TO (610, 620, 630), I
 610  IV(NFCOV) = IV(NFCOV) + 1
      IV(NFCALL) = IV(NFCALL) + 1
      IV(1) = 1
      GO TO 710
C
 620  IV(NGCOV) = IV(NGCOV) + 1
      IV(NGCALL) = IV(NGCALL) + 1
      IV(NFGCAL) = IV(NFCALL) + IV(NGCOV)
      GO TO 690
C
 630  IF (IV(CNVCOD) .EQ. 70) GO TO 120
      GO TO 660
C
 640  H1 = IABS(IV(H))
      IV(FDH) = H1
      IV(H) = -H1
      HC1 = IV(HC)
      IF (HC1 .LE. 0) GO TO 650
           CALL V7CPY(P*(P+1)/2, V(H1), V(HC1))
           GO TO 660
 650  RMAT1 = IV(RMAT)
      CALL  L7SQR(P, V(H1), V(RMAT1))
C
 660  IV(MODE) = 0
      IV(1) = IV(CNVCOD)
      IV(CNVCOD) = 0
      GO TO 999
C
C  ***  SPECIAL RETURN FOR MISSING HESSIAN INFORMATION -- BOTH
C  ***  IV(HC) .LE. 0 AND IV(RMAT) .LE. 0
C
 670  IV(1) = 1400
      GO TO 999
C
C  ***  INCONSISTENT B  ***
C
 680  IV(1) = 82
      GO TO 999
C
C  *** SAVE, THEN INITIALIZE IPIVOT ARRAY BEFORE COMPUTING G ***
C
 690  IV(1) = 2
      J = IV(IPIVOT)
      IPI = IV(PERM)
      CALL I7PNVR(P, IV(IPI), IV(J))
      DO 700 I = 1, P
         IV(J) = I
         J = J + 1
 700     CONTINUE
C
C  ***  PROJECT X INTO FEASIBLE REGION (PRIOR TO COMPUTING F OR G)  ***
C
 710  DO 720 I = 1, P
         IF (X(I) .LT. B(1,I)) X(I) = B(1,I)
         IF (X(I) .GT. B(2,I)) X(I) = B(2,I)
 720     CONTINUE
      IV(TOOBIG) = 0
C
 999  RETURN
C
C  ***  LAST LINE OF  G7ITB FOLLOWS  ***
      END