share/lbfgs/mcstep.f

   1       SUBROUTINE MCSTEP(STX,FX,DX,STY,FY,DY,STP,FP,DP,BRACKT,
   2      *                 STPMIN,STPMAX,INFO)
   3       INTEGER INFO
   4       DOUBLE PRECISION STX,FX,DX,STY,FY,DY,STP,FP,DP,STPMIN,STPMAX
   5       LOGICAL BRACKT,BOUND
   6 C
   7 C     SUBROUTINE MCSTEP
   8 C
   9 C     THE PURPOSE OF MCSTEP IS TO COMPUTE A SAFEGUARDED STEP FOR
  10 C     A LINESEARCH AND TO UPDATE AN INTERVAL OF UNCERTAINTY FOR
  11 C     A MINIMIZER OF THE FUNCTION.
  12 C
  13 C     THE PARAMETER STX CONTAINS THE STEP WITH THE LEAST FUNCTION
  14 C     VALUE. THE PARAMETER STP CONTAINS THE CURRENT STEP. IT IS
  15 C     ASSUMED THAT THE DERIVATIVE AT STX IS NEGATIVE IN THE
  16 C     DIRECTION OF THE STEP. IF BRACKT IS SET TRUE THEN A
  17 C     MINIMIZER HAS BEEN BRACKETED IN AN INTERVAL OF UNCERTAINTY
  18 C     WITH ENDPOINTS STX AND STY.
  19 C
  20 C     THE SUBROUTINE STATEMENT IS
  21 C
  22 C       SUBROUTINE MCSTEP(STX,FX,DX,STY,FY,DY,STP,FP,DP,BRACKT,
  23 C                        STPMIN,STPMAX,INFO)
  24 C
  25 C     WHERE
  26 C
  27 C       STX, FX, AND DX ARE VARIABLES WHICH SPECIFY THE STEP,
  28 C         THE FUNCTION, AND THE DERIVATIVE AT THE BEST STEP OBTAINED
  29 C         SO FAR. THE DERIVATIVE MUST BE NEGATIVE IN THE DIRECTION
  30 C         OF THE STEP, THAT IS, DX AND STP-STX MUST HAVE OPPOSITE
  31 C         SIGNS. ON OUTPUT THESE PARAMETERS ARE UPDATED APPROPRIATELY.
  32 C
  33 C       STY, FY, AND DY ARE VARIABLES WHICH SPECIFY THE STEP,
  34 C         THE FUNCTION, AND THE DERIVATIVE AT THE OTHER ENDPOINT OF
  35 C         THE INTERVAL OF UNCERTAINTY. ON OUTPUT THESE PARAMETERS ARE
  36 C         UPDATED APPROPRIATELY.
  37 C
  38 C       STP, FP, AND DP ARE VARIABLES WHICH SPECIFY THE STEP,
  39 C         THE FUNCTION, AND THE DERIVATIVE AT THE CURRENT STEP.
  40 C         IF BRACKT IS SET TRUE THEN ON INPUT STP MUST BE
  41 C         BETWEEN STX AND STY. ON OUTPUT STP IS SET TO THE NEW STEP.
  42 C
  43 C       BRACKT IS A LOGICAL VARIABLE WHICH SPECIFIES IF A MINIMIZER
  44 C         HAS BEEN BRACKETED. IF THE MINIMIZER HAS NOT BEEN BRACKETED
  45 C         THEN ON INPUT BRACKT MUST BE SET FALSE. IF THE MINIMIZER
  46 C         IS BRACKETED THEN ON OUTPUT BRACKT IS SET TRUE.
  47 C
  48 C       STPMIN AND STPMAX ARE INPUT VARIABLES WHICH SPECIFY LOWER
  49 C         AND UPPER BOUNDS FOR THE STEP.
  50 C
  51 C       INFO IS AN INTEGER OUTPUT VARIABLE SET AS FOLLOWS:
  52 C         IF INFO = 1,2,3,4,5, THEN THE STEP HAS BEEN COMPUTED
  53 C         ACCORDING TO ONE OF THE FIVE CASES BELOW. OTHERWISE
  54 C         INFO = 0, AND THIS INDICATES IMPROPER INPUT PARAMETERS.
  55 C
  56 C     SUBPROGRAMS CALLED
  57 C
  58 C       FORTRAN-SUPPLIED ... ABS,MAX,MIN,SQRT
  59 C
  60 C     ARGONNE NATIONAL LABORATORY. MINPACK PROJECT. JUNE 1983
  61 C     JORGE J. MORE', DAVID J. THUENTE
  62 C
  63       DOUBLE PRECISION GAMMA,P,Q,R,S,SGND,STPC,STPF,STPQ,THETA
  64       INFO = 0
  65 C
  66 C     CHECK THE INPUT PARAMETERS FOR ERRORS.
  67 C
  68       IF ((BRACKT .AND. (STP .LE. MIN(STX,STY) .OR.
  69      *    STP .GE. MAX(STX,STY))) .OR.
  70      *    DX*(STP-STX) .GE. 0.0D0 .OR. STPMAX .LT. STPMIN) RETURN
  71 C
  72 C     DETERMINE IF THE DERIVATIVES HAVE OPPOSITE SIGN.
  73 C
  74       SGND = DP*(DX/DABS(DX))
  75 C
  76 C     FIRST CASE. A HIGHER FUNCTION VALUE.
  77 C     THE MINIMUM IS BRACKETED. IF THE CUBIC STEP IS CLOSER
  78 C     TO STX THAN THE QUADRATIC STEP, THE CUBIC STEP IS TAKEN,
  79 C     ELSE THE AVERAGE OF THE CUBIC AND QUADRATIC STEPS IS TAKEN.
  80 C
  81       IF (FP .GT. FX) THEN
  82          INFO = 1
  83          BOUND = .TRUE.
  84          THETA = 3*(FX - FP)/(STP - STX) + DX + DP
  85          S = MAX(DABS(THETA),DABS(DX),DABS(DP))
  86          GAMMA = S*DSQRT((THETA/S)**2 - (DX/S)*(DP/S))
  87          IF (STP .LT. STX) GAMMA = -GAMMA
  88          P = (GAMMA - DX) + THETA
  89          Q = ((GAMMA - DX) + GAMMA) + DP
  90          R = P/Q
  91          STPC = STX + R*(STP - STX)
  92          STPQ = STX + ((DX/((FX-FP)/(STP-STX)+DX))/2)*(STP - STX)
  93          IF (DABS(STPC-STX) .LT. DABS(STPQ-STX)) THEN
  94             STPF = STPC
  95          ELSE
  96            STPF = STPC + (STPQ - STPC)/2
  97            END IF
  98          BRACKT = .TRUE.
  99 C
 100 C     SECOND CASE. A LOWER FUNCTION VALUE AND DERIVATIVES OF
 101 C     OPPOSITE SIGN. THE MINIMUM IS BRACKETED. IF THE CUBIC
 102 C     STEP IS CLOSER TO STX THAN THE QUADRATIC (SECANT) STEP,
 103 C     THE CUBIC STEP IS TAKEN, ELSE THE QUADRATIC STEP IS TAKEN.
 104 C
 105       ELSE IF (SGND .LT. 0.0D0) THEN
 106          INFO = 2
 107          BOUND = .FALSE.
 108          THETA = 3*(FX - FP)/(STP - STX) + DX + DP
 109          S = MAX(DABS(THETA),DABS(DX),DABS(DP))
 110          GAMMA = S*DSQRT((THETA/S)**2 - (DX/S)*(DP/S))
 111          IF (STP .GT. STX) GAMMA = -GAMMA
 112          P = (GAMMA - DP) + THETA
 113          Q = ((GAMMA - DP) + GAMMA) + DX
 114          R = P/Q
 115          STPC = STP + R*(STX - STP)
 116          STPQ = STP + (DP/(DP-DX))*(STX - STP)
 117          IF (DABS(STPC-STP) .GT. DABS(STPQ-STP)) THEN
 118             STPF = STPC
 119          ELSE
 120             STPF = STPQ
 121             END IF
 122          BRACKT = .TRUE.
 123 C
 124 C     THIRD CASE. A LOWER FUNCTION VALUE, DERIVATIVES OF THE
 125 C     SAME SIGN, AND THE MAGNITUDE OF THE DERIVATIVE DECREASES.
 126 C     THE CUBIC STEP IS ONLY USED IF THE CUBIC TENDS TO INFINITY
 127 C     IN THE DIRECTION OF THE STEP OR IF THE MINIMUM OF THE CUBIC
 128 C     IS BEYOND STP. OTHERWISE THE CUBIC STEP IS DEFINED TO BE
 129 C     EITHER STPMIN OR STPMAX. THE QUADRATIC (SECANT) STEP IS ALSO
 130 C     COMPUTED AND IF THE MINIMUM IS BRACKETED THEN THE THE STEP
 131 C     CLOSEST TO STX IS TAKEN, ELSE THE STEP FARTHEST AWAY IS TAKEN.
 132 C
 133       ELSE IF (DABS(DP) .LT. DABS(DX)) THEN
 134          INFO = 3
 135          BOUND = .TRUE.
 136          THETA = 3*(FX - FP)/(STP - STX) + DX + DP
 137          S = MAX(DABS(THETA),DABS(DX),DABS(DP))
 138 C
 139 C        THE CASE GAMMA = 0 ONLY ARISES IF THE CUBIC DOES NOT TEND
 140 C        TO INFINITY IN THE DIRECTION OF THE STEP.
 141 C
 142          GAMMA = S*DSQRT(MAX(0.0D0,(THETA/S)**2 - (DX/S)*(DP/S)))
 143          IF (STP .GT. STX) GAMMA = -GAMMA
 144          P = (GAMMA - DP) + THETA
 145          Q = (GAMMA + (DX - DP)) + GAMMA
 146          R = P/Q
 147          IF (R .LT. 0.0D0 .AND. GAMMA .NE. 0.0D0) THEN
 148             STPC = STP + R*(STX - STP)
 149          ELSE IF (STP .GT. STX) THEN
 150             STPC = STPMAX
 151          ELSE
 152             STPC = STPMIN
 153             END IF
 154          STPQ = STP + (DP/(DP-DX))*(STX - STP)
 155          IF (BRACKT) THEN
 156             IF (DABS(STP-STPC) .LT. DABS(STP-STPQ)) THEN
 157                STPF = STPC
 158             ELSE
 159                STPF = STPQ
 160                END IF
 161          ELSE
 162             IF (DABS(STP-STPC) .GT. DABS(STP-STPQ)) THEN
 163                STPF = STPC
 164             ELSE
 165                STPF = STPQ
 166                END IF
 167             END IF
 168 C
 169 C     FOURTH CASE. A LOWER FUNCTION VALUE, DERIVATIVES OF THE
 170 C     SAME SIGN, AND THE MAGNITUDE OF THE DERIVATIVE DOES
 171 C     NOT DECREASE. IF THE MINIMUM IS NOT BRACKETED, THE STEP
 172 C     IS EITHER STPMIN OR STPMAX, ELSE THE CUBIC STEP IS TAKEN.
 173 C
 174       ELSE
 175          INFO = 4
 176          BOUND = .FALSE.
 177          IF (BRACKT) THEN
 178             THETA = 3*(FP - FY)/(STY - STP) + DY + DP
 179             S = MAX(DABS(THETA),DABS(DY),DABS(DP))
 180             GAMMA = S*DSQRT((THETA/S)**2 - (DY/S)*(DP/S))
 181             IF (STP .GT. STY) GAMMA = -GAMMA
 182             P = (GAMMA - DP) + THETA
 183             Q = ((GAMMA - DP) + GAMMA) + DY
 184             R = P/Q
 185             STPC = STP + R*(STY - STP)
 186             STPF = STPC
 187          ELSE IF (STP .GT. STX) THEN
 188             STPF = STPMAX
 189          ELSE
 190             STPF = STPMIN
 191             END IF
 192          END IF
 193 C
 194 C     UPDATE THE INTERVAL OF UNCERTAINTY. THIS UPDATE DOES NOT
 195 C     DEPEND ON THE NEW STEP OR THE CASE ANALYSIS ABOVE.
 196 C
 197       IF (FP .GT. FX) THEN
 198          STY = STP
 199          FY = FP
 200          DY = DP
 201       ELSE
 202          IF (SGND .LT. 0.0D0) THEN
 203             STY = STX
 204             FY = FX
 205             DY = DX
 206             END IF
 207          STX = STP
 208          FX = FP
 209          DX = DP
 210          END IF
 211 C
 212 C     COMPUTE THE NEW STEP AND SAFEGUARD IT.
 213 C
 214       STPF = MIN(STPMAX,STPF)
 215       STPF = MAX(STPMIN,STPF)
 216       STP = STPF
 217       IF (BRACKT .AND. BOUND) THEN
 218          IF (STY .GT. STX) THEN
 219             STP = MIN(STX+0.66D0*(STY-STX),STP)
 220          ELSE
 221             STP = MAX(STX+0.66D0*(STY-STX),STP)
 222             END IF
 223          END IF
 224       RETURN
 225 C
 226 C     LAST LINE OF SUBROUTINE MCSTEP.
 227 C
 228       END