SUBROUTINE DRN2GB(B, D, DR, IV, LIV, LV, N, ND, N1, N2, P, R, 1 RD, V, X) C C *** REVISED ITERATION DRIVER FOR NL2SOL WITH SIMPLE BOUNDS *** C INTEGER LIV, LV, N, ND, N1, N2, P INTEGER IV(LIV) DOUBLE PRECISION B(2,P), D(P), DR(ND,P), R(ND), RD(ND), V(LV), 1 X(P) C C-------------------------- PARAMETER USAGE -------------------------- C C B........ BOUNDS ON X. C D........ SCALE VECTOR. C DR....... DERIVATIVES OF R AT X. C IV....... INTEGER VALUES ARRAY. C LIV...... LENGTH OF IV... LIV MUST BE AT LEAST 4*P + 82. C LV....... LENGTH OF V... LV MUST BE AT LEAST 105 + P*(2*P+20). C N........ TOTAL NUMBER OF RESIDUALS. C ND....... MAX. NO. OF RESIDUALS PASSED ON ONE CALL. C N1....... LOWEST ROW INDEX FOR RESIDUALS SUPPLIED THIS TIME. C N2....... HIGHEST ROW INDEX FOR RESIDUALS SUPPLIED THIS TIME. C P........ NUMBER OF PARAMETERS (COMPONENTS OF X) BEING ESTIMATED. C R........ RESIDUALS. C V........ FLOATING-POINT VALUES ARRAY. C X........ PARAMETER VECTOR BEING ESTIMATED (INPUT = INITIAL GUESS, C OUTPUT = BEST VALUE FOUND). C C *** DISCUSSION *** C C THIS ROUTINE CARRIES OUT ITERATIONS FOR SOLVING NONLINEAR C LEAST SQUARES PROBLEMS. IT IS SIMILAR TO DRN2G, EXCEPT THAT C THIS ROUTINE ENFORCES THE BOUNDS B(1,I) .LE. X(I) .LE. B(2,I), C I = 1(1)P. C C *** GENERAL *** C C CODED BY DAVID M. GAY. C C+++++++++++++++++++++++++++++ DECLARATIONS ++++++++++++++++++++++++++ C C *** EXTERNAL FUNCTIONS AND SUBROUTINES *** C DOUBLE PRECISION DD7TPR, DV2NRM EXTERNAL DIVSET, DD7TPR,DD7UPD, DG7ITB,DITSUM,DL7VML, DQ7APL, 1 DQ7RAD, DR7TVM,DV7CPY, DV7SCP, DV2NRM C C DIVSET.... PROVIDES DEFAULT IV AND V INPUT COMPONENTS. C DD7TPR... COMPUTES INNER PRODUCT OF TWO VECTORS. C DD7UPD... UPDATES SCALE VECTOR D. C DG7ITB... PERFORMS BASIC MINIMIZATION ALGORITHM. C DITSUM.... PRINTS ITERATION SUMMARY, INFO ABOUT INITIAL AND FINAL X. C DL7VML.... COMPUTES L * V, V = VECTOR, L = LOWER TRIANGULAR MATRIX. C DQ7APL... APPLIES QR TRANSFORMATIONS STORED BY DQ7RAD. C DQ7RAD.... ADDS A NEW BLOCK OF ROWS TO QR DECOMPOSITION. C DR7TVM... MULT. VECTOR BY TRANS. OF UPPER TRIANG. MATRIX FROM QR FACT. C DV7CPY.... COPIES ONE VECTOR TO ANOTHER. C DV7SCP... SETS ALL ELEMENTS OF A VECTOR TO A SCALAR. C DV2NRM... RETURNS THE 2-NORM OF A VECTOR. C C C *** LOCAL VARIABLES *** C INTEGER G1, GI, I, IV1, IVMODE, JTOL1, L, LH, NN, QTR1, 1 RD1, RMAT1, YI, Y1 DOUBLE PRECISION T C DOUBLE PRECISION HALF, ZERO C C *** SUBSCRIPTS FOR IV AND V *** C INTEGER DINIT, DTYPE, DTINIT, D0INIT, F, G, JCN, JTOL, MODE, 1 NEXTV, NF0, NF00, NF1, NFCALL, NFCOV, NFGCAL, QTR, RDREQ, 1 REGD, RESTOR, RLIMIT, RMAT, TOOBIG, VNEED C C *** IV SUBSCRIPT VALUES *** C C/6 C DATA DTYPE/16/, G/28/, JCN/66/, JTOL/59/, MODE/35/, NEXTV/47/, C 1 NF0/68/, NF00/81/, NF1/69/, NFCALL/6/, NFCOV/52/, NFGCAL/7/, C 2 QTR/77/, RDREQ/57/, RESTOR/9/, REGD/67/, RMAT/78/, TOOBIG/2/, C 3 VNEED/4/ C/7 PARAMETER (DTYPE=16, G=28, JCN=66, JTOL=59, MODE=35, NEXTV=47, 1 NF0=68, NF00=81, NF1=69, NFCALL=6, NFCOV=52, NFGCAL=7, 2 QTR=77, RDREQ=57, RESTOR=9, REGD=67, RMAT=78, TOOBIG=2, 3 VNEED=4) C/ C C *** V SUBSCRIPT VALUES *** C C/6 C DATA DINIT/38/, DTINIT/39/, D0INIT/40/, F/10/, RLIMIT/46/ C/7 PARAMETER (DINIT=38, DTINIT=39, D0INIT=40, F=10, RLIMIT=46) C/ C/6 C DATA HALF/0.5D+0/, ZERO/0.D+0/ C/7 PARAMETER (HALF=0.5D+0, ZERO=0.D+0) C/ C C+++++++++++++++++++++++++++++++ BODY ++++++++++++++++++++++++++++++++ C LH = P * (P+1) / 2 IF (IV(1) .EQ. 0) CALL DIVSET(1, IV, LIV, LV, V) IV1 = IV(1) IF (IV1 .GT. 2) GO TO 10 NN = N2 - N1 + 1 IV(RESTOR) = 0 I = IV1 + 4 IF (IV(TOOBIG) .EQ. 0) GO TO (150, 130, 150, 120, 120, 150), I IF (I .NE. 5) IV(1) = 2 GO TO 40 C C *** FRESH START OR RESTART -- CHECK INPUT INTEGERS *** C 10 IF (ND .LE. 0) GO TO 220 IF (P .LE. 0) GO TO 220 IF (N .LE. 0) GO TO 220 IF (IV1 .EQ. 14) GO TO 30 IF (IV1 .GT. 16) GO TO 270 IF (IV1 .LT. 12) GO TO 40 IF (IV1 .EQ. 12) IV(1) = 13 IF (IV(1) .NE. 13) GO TO 20 IV(VNEED) = IV(VNEED) + P*(P+15)/2 20 CALL DG7ITB(B, D, X, IV, LIV, LV, P, P, V, X, X) IF (IV(1) .NE. 14) GO TO 999 C C *** STORAGE ALLOCATION *** C IV(G) = IV(NEXTV) IV(JCN) = IV(G) + 2*P IV(RMAT) = IV(JCN) + P IV(QTR) = IV(RMAT) + LH IV(JTOL) = IV(QTR) + 2*P IV(NEXTV) = IV(JTOL) + 2*P C *** TURN OFF COVARIANCE COMPUTATION *** IV(RDREQ) = 0 IF (IV1 .EQ. 13) GO TO 999 C 30 JTOL1 = IV(JTOL) IF (V(DINIT) .GE. ZERO) CALL DV7SCP(P, D, V(DINIT)) IF (V(DTINIT) .GT. ZERO) CALL DV7SCP(P, V(JTOL1), V(DTINIT)) I = JTOL1 + P IF (V(D0INIT) .GT. ZERO) CALL DV7SCP(P, V(I), V(D0INIT)) IV(NF0) = 0 IV(NF1) = 0 IF (ND .GE. N) GO TO 40 C C *** SPECIAL CASE HANDLING OF FIRST FUNCTION AND GRADIENT EVALUATION C *** -- ASK FOR BOTH RESIDUAL AND JACOBIAN AT ONCE C G1 = IV(G) Y1 = G1 + P CALL DG7ITB(B, D, V(G1), IV, LIV, LV, P, P, V, X, V(Y1)) IF (IV(1) .NE. 1) GO TO 260 V(F) = ZERO CALL DV7SCP(P, V(G1), ZERO) IV(1) = -1 QTR1 = IV(QTR) CALL DV7SCP(P, V(QTR1), ZERO) IV(REGD) = 0 RMAT1 = IV(RMAT) GO TO 100 C 40 G1 = IV(G) Y1 = G1 + P CALL DG7ITB(B, D, V(G1), IV, LIV, LV, P, P, V, X, V(Y1)) IF (IV(1) - 2) 50, 60, 260 C 50 V(F) = ZERO IF (IV(NF1) .EQ. 0) GO TO 240 IF (IV(RESTOR) .NE. 2) GO TO 240 IV(NF0) = IV(NF1) CALL DV7CPY(N, RD, R) IV(REGD) = 0 GO TO 240 C 60 CALL DV7SCP(P, V(G1), ZERO) IF (IV(MODE) .GT. 0) GO TO 230 RMAT1 = IV(RMAT) QTR1 = IV(QTR) RD1 = QTR1 + P CALL DV7SCP(P, V(QTR1), ZERO) IV(REGD) = 0 IF (ND .LT. N) GO TO 90 IF (N1 .NE. 1) GO TO 90 IF (IV(MODE) .LT. 0) GO TO 100 IF (IV(NF1) .EQ. IV(NFGCAL)) GO TO 70 IF (IV(NF0) .NE. IV(NFGCAL)) GO TO 90 CALL DV7CPY(N, R, RD) GO TO 80 70 CALL DV7CPY(N, RD, R) 80 CALL DQ7APL(ND, N, P, DR, RD, 0) CALL DR7TVM(ND, MIN0(N,P), V(Y1), V(RD1), DR, RD) IV(REGD) = 0 GO TO 110 C 90 IV(1) = -2 IF (IV(MODE) .LT. 0) IV(1) = -3 100 CALL DV7SCP(P, V(Y1), ZERO) 110 CALL DV7SCP(LH, V(RMAT1), ZERO) GO TO 240 C C *** COMPUTE F(X) *** C 120 T = DV2NRM(NN, R) IF (T .GT. V(RLIMIT)) GO TO 210 V(F) = V(F) + HALF * T**2 IF (N2 .LT. N) GO TO 250 IF (N1 .EQ. 1) IV(NF1) = IV(NFCALL) GO TO 40 C C *** COMPUTE Y *** C 130 Y1 = IV(G) + P YI = Y1 DO 140 L = 1, P V(YI) = V(YI) + DD7TPR(NN, DR(1,L), R) YI = YI + 1 140 CONTINUE IF (N2 .LT. N) GO TO 250 IV(1) = 2 IF (N1 .GT. 1) IV(1) = -3 GO TO 240 C C *** COMPUTE GRADIENT INFORMATION *** C 150 G1 = IV(G) IVMODE = IV(MODE) IF (IVMODE .LT. 0) GO TO 170 IF (IVMODE .EQ. 0) GO TO 180 IV(1) = 2 C C *** COMPUTE GRADIENT ONLY (FOR USE IN COVARIANCE COMPUTATION) *** C GI = G1 DO 160 L = 1, P V(GI) = V(GI) + DD7TPR(NN, R, DR(1,L)) GI = GI + 1 160 CONTINUE GO TO 200 C C *** COMPUTE INITIAL FUNCTION VALUE WHEN ND .LT. N *** C 170 IF (N .LE. ND) GO TO 180 T = DV2NRM(NN, R) IF (T .GT. V(RLIMIT)) GO TO 210 V(F) = V(F) + HALF * T**2 C C *** UPDATE D IF DESIRED *** C 180 IF (IV(DTYPE) .GT. 0) 1 CALL DD7UPD(D, DR, IV, LIV, LV, N, ND, NN, N2, P, V) C C *** COMPUTE RMAT AND QTR *** C QTR1 = IV(QTR) RMAT1 = IV(RMAT) CALL DQ7RAD(NN, ND, P, V(QTR1), .TRUE., V(RMAT1), DR, R) IV(NF1) = 0 IF (N1 .GT. 1) GO TO 200 IF (N2 .LT. N) GO TO 250 C C *** SAVE DIAGONAL OF R FOR COMPUTING Y LATER *** C RD1 = QTR1 + P L = RMAT1 - 1 DO 190 I = 1, P L = L + I V(RD1) = V(L) RD1 = RD1 + 1 190 CONTINUE C 200 IF (N2 .LT. N) GO TO 250 IF (IVMODE .GT. 0) GO TO 40 IV(NF00) = IV(NFGCAL) C C *** COMPUTE G FROM RMAT AND QTR *** C CALL DL7VML(P, V(G1), V(RMAT1), V(QTR1)) IV(1) = 2 IF (IVMODE .EQ. 0) GO TO 40 IF (N .LE. ND) GO TO 40 C C *** FINISH SPECIAL CASE HANDLING OF FIRST FUNCTION AND GRADIENT C Y1 = G1 + P IV(1) = 1 CALL DG7ITB(B, D, V(G1), IV, LIV, LV, P, P, V, X, V(Y1)) IF (IV(1) .NE. 2) GO TO 260 GO TO 40 C C *** MISC. DETAILS *** C C *** X IS OUT OF RANGE (OVERSIZE STEP) *** C 210 IV(TOOBIG) = 1 GO TO 40 C C *** BAD N, ND, OR P *** C 220 IV(1) = 66 GO TO 270 C C *** RECORD EXTRA EVALUATIONS FOR FINITE-DIFFERENCE HESSIAN *** C 230 IV(NFCOV) = IV(NFCOV) + 1 IV(NFCALL) = IV(NFCALL) + 1 IV(NFGCAL) = IV(NFCALL) IV(1) = -1 C C *** RETURN FOR MORE FUNCTION OR GRADIENT INFORMATION *** C 240 N2 = 0 250 N1 = N2 + 1 N2 = N2 + ND IF (N2 .GT. N) N2 = N GO TO 999 C C *** PRINT SUMMARY OF FINAL ITERATION AND OTHER REQUESTED ITEMS *** C 260 G1 = IV(G) 270 CALL DITSUM(D, V(G1), IV, LIV, LV, P, V, X) C 999 RETURN C *** LAST CARD OF DRN2GB FOLLOWS *** END .