SUBROUTINE DD7DOG(DIG, LV, N, NWTSTP, STEP, V) C C *** COMPUTE DOUBLE DOGLEG STEP *** C C *** PARAMETER DECLARATIONS *** C INTEGER LV, N DOUBLE PRECISION DIG(N), NWTSTP(N), STEP(N), V(LV) C C *** PURPOSE *** C C THIS SUBROUTINE COMPUTES A CANDIDATE STEP (FOR USE IN AN UNCON- C STRAINED MINIMIZATION CODE) BY THE DOUBLE DOGLEG ALGORITHM OF C DENNIS AND MEI (REF. 1), WHICH IS A VARIATION ON POWELL*S DOGLEG C SCHEME (REF. 2, P. 95). C C-------------------------- PARAMETER USAGE -------------------------- C C DIG (INPUT) DIAG(D)**-2 * G -- SEE ALGORITHM NOTES. C G (INPUT) THE CURRENT GRADIENT VECTOR. C LV (INPUT) LENGTH OF V. C N (INPUT) NUMBER OF COMPONENTS IN DIG, G, NWTSTP, AND STEP. C NWTSTP (INPUT) NEGATIVE NEWTON STEP -- SEE ALGORITHM NOTES. C STEP (OUTPUT) THE COMPUTED STEP. C V (I/O) VALUES ARRAY, THE FOLLOWING COMPONENTS OF WHICH ARE C USED HERE... C V(BIAS) (INPUT) BIAS FOR RELAXED NEWTON STEP, WHICH IS V(BIAS) OF C THE WAY FROM THE FULL NEWTON TO THE FULLY RELAXED NEWTON C STEP. RECOMMENDED VALUE = 0.8 . C V(DGNORM) (INPUT) 2-NORM OF DIAG(D)**-1 * G -- SEE ALGORITHM NOTES. C V(DSTNRM) (OUTPUT) 2-NORM OF DIAG(D) * STEP, WHICH IS V(RADIUS) C UNLESS V(STPPAR) = 0 -- SEE ALGORITHM NOTES. C V(DST0) (INPUT) 2-NORM OF DIAG(D) * NWTSTP -- SEE ALGORITHM NOTES. C V(GRDFAC) (OUTPUT) THE COEFFICIENT OF DIG IN THE STEP RETURNED -- C STEP(I) = V(GRDFAC)*DIG(I) + V(NWTFAC)*NWTSTP(I). C V(GTHG) (INPUT) SQUARE-ROOT OF (DIG**T) * (HESSIAN) * DIG -- SEE C ALGORITHM NOTES. C V(GTSTEP) (OUTPUT) INNER PRODUCT BETWEEN G AND STEP. C V(NREDUC) (OUTPUT) FUNCTION REDUCTION PREDICTED FOR THE FULL NEWTON C STEP. C V(NWTFAC) (OUTPUT) THE COEFFICIENT OF NWTSTP IN THE STEP RETURNED -- C SEE V(GRDFAC) ABOVE. C V(PREDUC) (OUTPUT) FUNCTION REDUCTION PREDICTED FOR THE STEP RETURNED. C V(RADIUS) (INPUT) THE TRUST REGION RADIUS. D TIMES THE STEP RETURNED C HAS 2-NORM V(RADIUS) UNLESS V(STPPAR) = 0. C V(STPPAR) (OUTPUT) CODE TELLING HOW STEP WAS COMPUTED... 0 MEANS A C FULL NEWTON STEP. BETWEEN 0 AND 1 MEANS V(STPPAR) OF THE C WAY FROM THE NEWTON TO THE RELAXED NEWTON STEP. BETWEEN C 1 AND 2 MEANS A TRUE DOUBLE DOGLEG STEP, V(STPPAR) - 1 OF C THE WAY FROM THE RELAXED NEWTON TO THE CAUCHY STEP. C GREATER THAN 2 MEANS 1 / (V(STPPAR) - 1) TIMES THE CAUCHY C STEP. C C------------------------------- NOTES ------------------------------- C C *** ALGORITHM NOTES *** C C LET G AND H BE THE CURRENT GRADIENT AND HESSIAN APPROXIMA- C TION RESPECTIVELY AND LET D BE THE CURRENT SCALE VECTOR. THIS C ROUTINE ASSUMES DIG = DIAG(D)**-2 * G AND NWTSTP = H**-1 * G. C THE STEP COMPUTED IS THE SAME ONE WOULD GET BY REPLACING G AND H C BY DIAG(D)**-1 * G AND DIAG(D)**-1 * H * DIAG(D)**-1, C COMPUTING STEP, AND TRANSLATING STEP BACK TO THE ORIGINAL C VARIABLES, I.E., PREMULTIPLYING IT BY DIAG(D)**-1. C C *** REFERENCES *** C C 1. DENNIS, J.E., AND MEI, H.H.W. (1979), TWO NEW UNCONSTRAINED OPTI- C MIZATION ALGORITHMS WHICH USE FUNCTION AND GRADIENT C VALUES, J. OPTIM. THEORY APPLIC. 28, PP. 453-482. C 2. POWELL, M.J.D. (1970), A HYBRID METHOD FOR NON-LINEAR EQUATIONS, C IN NUMERICAL METHODS FOR NON-LINEAR EQUATIONS, EDITED BY C P. RABINOWITZ, GORDON AND BREACH, LONDON. C C *** GENERAL *** C C CODED BY DAVID M. GAY. C THIS SUBROUTINE WAS WRITTEN IN CONNECTION WITH RESEARCH SUPPORTED C BY THE NATIONAL SCIENCE FOUNDATION UNDER GRANTS MCS-7600324 AND C MCS-7906671. C C------------------------ EXTERNAL QUANTITIES ------------------------ C C *** INTRINSIC FUNCTIONS *** C/+ DOUBLE PRECISION DSQRT C/ C-------------------------- LOCAL VARIABLES -------------------------- C INTEGER I DOUBLE PRECISION CFACT, CNORM, CTRNWT, GHINVG, FEMNSQ, GNORM, 1 NWTNRM, RELAX, RLAMBD, T, T1, T2 DOUBLE PRECISION HALF, ONE, TWO, ZERO C C *** V SUBSCRIPTS *** C INTEGER BIAS, DGNORM, DSTNRM, DST0, GRDFAC, GTHG, GTSTEP, 1 NREDUC, NWTFAC, PREDUC, RADIUS, STPPAR C C *** DATA INITIALIZATIONS *** C C/6 C DATA HALF/0.5D+0/, ONE/1.D+0/, TWO/2.D+0/, ZERO/0.D+0/ C/7 PARAMETER (HALF=0.5D+0, ONE=1.D+0, TWO=2.D+0, ZERO=0.D+0) C/ C C/6 C DATA BIAS/43/, DGNORM/1/, DSTNRM/2/, DST0/3/, GRDFAC/45/, C 1 GTHG/44/, GTSTEP/4/, NREDUC/6/, NWTFAC/46/, PREDUC/7/, C 2 RADIUS/8/, STPPAR/5/ C/7 PARAMETER (BIAS=43, DGNORM=1, DSTNRM=2, DST0=3, GRDFAC=45, 1 GTHG=44, GTSTEP=4, NREDUC=6, NWTFAC=46, PREDUC=7, 2 RADIUS=8, STPPAR=5) C/ C C+++++++++++++++++++++++++++++++ BODY ++++++++++++++++++++++++++++++++ C NWTNRM = V(DST0) RLAMBD = ONE IF (NWTNRM .GT. ZERO) RLAMBD = V(RADIUS) / NWTNRM GNORM = V(DGNORM) GHINVG = TWO * V(NREDUC) V(GRDFAC) = ZERO V(NWTFAC) = ZERO IF (RLAMBD .LT. ONE) GO TO 30 C C *** THE NEWTON STEP IS INSIDE THE TRUST REGION *** C V(STPPAR) = ZERO V(DSTNRM) = NWTNRM V(GTSTEP) = -GHINVG V(PREDUC) = V(NREDUC) V(NWTFAC) = -ONE DO 20 I = 1, N 20 STEP(I) = -NWTSTP(I) GO TO 999 C 30 V(DSTNRM) = V(RADIUS) CFACT = (GNORM / V(GTHG))**2 C *** CAUCHY STEP = -CFACT * G. CNORM = GNORM * CFACT RELAX = ONE - V(BIAS) * (ONE - GNORM*CNORM/GHINVG) IF (RLAMBD .LT. RELAX) GO TO 50 C C *** STEP IS BETWEEN RELAXED NEWTON AND FULL NEWTON STEPS *** C V(STPPAR) = ONE - (RLAMBD - RELAX) / (ONE - RELAX) T = -RLAMBD V(GTSTEP) = T * GHINVG V(PREDUC) = RLAMBD * (ONE - HALF*RLAMBD) * GHINVG V(NWTFAC) = T DO 40 I = 1, N 40 STEP(I) = T * NWTSTP(I) GO TO 999 C 50 IF (CNORM .LT. V(RADIUS)) GO TO 70 C C *** THE CAUCHY STEP LIES OUTSIDE THE TRUST REGION -- C *** STEP = SCALED CAUCHY STEP *** C T = -V(RADIUS) / GNORM V(GRDFAC) = T V(STPPAR) = ONE + CNORM / V(RADIUS) V(GTSTEP) = -V(RADIUS) * GNORM V(PREDUC) = V(RADIUS)*(GNORM - HALF*V(RADIUS)*(V(GTHG)/GNORM)**2) DO 60 I = 1, N 60 STEP(I) = T * DIG(I) GO TO 999 C C *** COMPUTE DOGLEG STEP BETWEEN CAUCHY AND RELAXED NEWTON *** C *** FEMUR = RELAXED NEWTON STEP MINUS CAUCHY STEP *** C 70 CTRNWT = CFACT * RELAX * GHINVG / GNORM C *** CTRNWT = INNER PROD. OF CAUCHY AND RELAXED NEWTON STEPS, C *** SCALED BY GNORM**-1. T1 = CTRNWT - GNORM*CFACT**2 C *** T1 = INNER PROD. OF FEMUR AND CAUCHY STEP, SCALED BY C *** GNORM**-1. T2 = V(RADIUS)*(V(RADIUS)/GNORM) - GNORM*CFACT**2 T = RELAX * NWTNRM FEMNSQ = (T/GNORM)*T - CTRNWT - T1 C *** FEMNSQ = SQUARE OF 2-NORM OF FEMUR, SCALED BY GNORM**-1. T = T2 / (T1 + DSQRT(T1**2 + FEMNSQ*T2)) C *** DOGLEG STEP = CAUCHY STEP + T * FEMUR. T1 = (T - ONE) * CFACT V(GRDFAC) = T1 T2 = -T * RELAX V(NWTFAC) = T2 V(STPPAR) = TWO - T V(GTSTEP) = T1*GNORM**2 + T2*GHINVG V(PREDUC) = -T1*GNORM * ((T2 + ONE)*GNORM) 1 - T2 * (ONE + HALF*T2)*GHINVG 2 - HALF * (V(GTHG)*T1)**2 DO 80 I = 1, N 80 STEP(I) = T1*DIG(I) + T2*NWTSTP(I) C 999 RETURN C *** LAST LINE OF DD7DOG FOLLOWS *** END