Routine to solve one-dimensional search in unconstrained minimization using 2-point quadratic interpolation, 3-point cubic interpolation and 4-point cubic interpolation.
BY G. N. VANDERPLAATS, APRIL, 1972.
OBJ = INITIAL FUNCTION VALUE.
| Type | Intent | Optional | Attributes | Name | ||
|---|---|---|---|---|---|---|
| class(conmin_class), | intent(inout) | :: | me | |||
| real(kind=wp), | intent(inout) | :: | x(:) | |||
| real(kind=wp), | intent(inout) | :: | s(:) | |||
| real(kind=wp), | intent(inout) | :: | slope |
SLOPE = INITIAL FUNCTION SLOPE = S-TRANSPOSE TIMES DF. SLOPE MUST BE NEGATIVE. |
||
| real(kind=wp), | intent(inout) | :: | alp |
PROPOSED MOVE PARAMETER |
||
| real(kind=wp), | intent(inout) | :: | fff | |||
| real(kind=wp), | intent(inout) | :: | a1 | |||
| real(kind=wp), | intent(inout) | :: | a2 | |||
| real(kind=wp), | intent(inout) | :: | a3 | |||
| real(kind=wp), | intent(inout) | :: | a4 | |||
| real(kind=wp), | intent(inout) | :: | f1 | |||
| real(kind=wp), | intent(inout) | :: | f2 | |||
| real(kind=wp), | intent(inout) | :: | f3 | |||
| real(kind=wp), | intent(inout) | :: | f4 | |||
| real(kind=wp), | intent(inout) | :: | app | |||
| integer, | intent(inout) | :: | ncal(2) | |||
| integer, | intent(out) | :: | kount | |||
| integer, | intent(inout) | :: | jgoto |
subroutine cnmn03(me, x, s, slope, alp, fff, a1, a2, a3, a4, f1, f2, f3, f4, app, & ncal, kount, jgoto) !! Routine to solve one-dimensional search in unconstrained !! minimization using 2-point quadratic interpolation, 3-point !! cubic interpolation and 4-point cubic interpolation. !! !! BY G. N. VANDERPLAATS, APRIL, 1972. !! !! OBJ = INITIAL FUNCTION VALUE. class(conmin_class), intent(inout) :: me real(wp), intent(inout) :: x(:), s(:) real(wp), intent(inout) :: slope !! SLOPE = INITIAL FUNCTION SLOPE = S-TRANSPOSE TIMES DF. !! SLOPE MUST BE NEGATIVE. real(wp), intent(inout) :: alp !! PROPOSED MOVE PARAMETER real(wp), intent(inout) :: fff, a1, a2, a3, a4, f1, f2, f3, f4 real(wp), intent(inout) :: app integer, intent(inout) :: ncal(2) integer, intent(out) :: kount integer, intent(inout) :: jgoto real(wp) :: aa, ab, ab2, ab3, ap, ff integer :: i, ii real(wp), parameter :: zro = 0.0_wp if (jgoto /= 0) then select case (jgoto) case (1); go to 30 case (2); go to 50 case (3); go to 80 case (4); go to 110 case (5); go to 150 case (6); go to 190 case (7); go to 230 end select end if ! ------------------------------------------------------------------ ! INITIAL INFORMATION (ALPHA=0) ! ------------------------------------------------------------------ if (slope >= 0.0_wp) then alp = 0.0_wp return end if if (me%iprint > 4) write (me%iunit, 5000) fff = me%obj a1 = 0.0_wp f1 = me%obj a2 = alp a3 = 0.0_wp f3 = 0.0_wp ap = a2 kount = 0 ! ------------------------------------------------------------------ ! MOVE A DISTANCE AP*S AND UPDATE FUNCTION VALUE ! ------------------------------------------------------------------ 10 kount = kount + 1 do i = 1, me%ndv x(i) = x(i) + ap*s(i) end do if (me%iprint > 4) write (me%iunit, 5100) ap if (me%iprint > 4) write (me%iunit, 5200) x(1:me%ndv) ncal(1) = ncal(1) + 1 jgoto = 1 return 30 f2 = me%obj if (me%iprint > 4) write (me%iunit, 5300) f2 if (f2 < f1) go to 90 ! ------------------------------------------------------------------ ! CHECK FOR ILL-CONDITIONING ! ------------------------------------------------------------------ if (kount <= 5) then ff = 2.0_wp*abs(f1) if (f2 < ff) go to 60 ff = 5.0_wp*abs(f1) if (f2 >= ff) then a2 = 0.5_wp*a2 ap = -a2 alp = a2 go to 10 end if end if f3 = f2 a3 = a2 a2 = 0.5_wp*a2 ! ------------------------------------------------------------------ ! UPDATE DESIGN VECTOR AND FUNCTION VALUE ! ------------------------------------------------------------------ ap = a2 - alp alp = a2 do i = 1, me%ndv x(i) = x(i) + ap*s(i) end do if (me%iprint > 4) write (me%iunit, 5100) a2 if (me%iprint > 4) write (me%iunit, 5200) x(1:me%ndv) ncal(1) = ncal(1) + 1 jgoto = 2 return 50 f2 = me%obj if (me%iprint > 4) write (me%iunit, 5300) f2 ! PROCEED TO CUBIC INTERPOLATION. go to 130 ! ------------------------------------------------------------------ ! ********** 2-POINT QUADRATIC INTERPOLATION ********** ! ------------------------------------------------------------------ 60 ii = 1 call me%cnmn04(ii, app, zro, a1, f1, slope, a2, f2, zro, zro, zro, zro) if (app < zro .or. app > a2) go to 90 f3 = f2 a3 = a2 a2 = app ! ------------------------------------------------------------------ ! UPDATE DESIGN VECTOR AND FUNCTION VALUE ! ------------------------------------------------------------------ ap = a2 - alp alp = a2 do i = 1, me%ndv x(i) = x(i) + ap*s(i) end do if (me%iprint > 4) write (me%iunit, 5100) a2 if (me%iprint > 4) write (me%iunit, 5200) x(1:me%ndv) ncal(1) = ncal(1) + 1 jgoto = 3 return 80 f2 = me%obj if (me%iprint > 4) write (me%iunit, 5300) f2 go to 120 90 a3 = 2.*a2 ! ------------------------------------------------------------------ ! UPDATE DESIGN VECTOR AND FUNCTION VALUE ! ------------------------------------------------------------------ ap = a3 - alp alp = a3 do i = 1, me%ndv x(i) = x(i) + ap*s(i) end do if (me%iprint > 4) write (me%iunit, 5100) a3 if (me%iprint > 4) write (me%iunit, 5200) x(1:me%ndv) ncal(1) = ncal(1) + 1 jgoto = 4 return 110 f3 = me%obj if (me%iprint > 4) write (me%iunit, 5300) f3 120 if (f3 < f2) go to 170 ! ------------------------------------------------------------------ ! ********** 3-POINT CUBIC INTERPOLATION ********** ! ------------------------------------------------------------------ 130 ii = 3 call me%cnmn04(ii, app, zro, a1, f1, slope, a2, f2, a3, f3, zro, zro) if (app < zro .or. app > a3) go to 170 ! ------------------------------------------------------------------ ! UPDATE DESIGN VECTOR AND FUNCTION VALUE. ! ------------------------------------------------------------------ ap = app - alp alp = app x(1:me%ndv) = x(1:me%ndv) + ap*s(1:me%ndv) if (me%iprint > 4) write (me%iunit, 5100) alp if (me%iprint > 4) write (me%iunit, 5200) x(1:me%ndv) ncal(1) = ncal(1) + 1 jgoto = 5 return 150 if (me%iprint > 4) write (me%iunit, 5300) me%obj ! ------------------------------------------------------------------ ! CHECK CONVERGENCE ! ------------------------------------------------------------------ aa = 1.-app/a2 ab2 = abs(f2) ab3 = abs(me%obj) ab = ab2 if (ab3 > ab) ab = ab3 if (ab < 1.0e-15_wp) ab = 1.0e-15_wp ab = (ab2 - ab3)/ab if (abs(ab) < 1.0e-15_wp .and. abs(aa) < 0.001_wp) go to 260 a4 = a3 f4 = f3 a3 = app f3 = me%obj if (a3 > a2) go to 200 a3 = a2 f3 = f2 a2 = app f2 = me%obj go to 200 ! ------------------------------------------------------------------ ! ********** 4-POINT CUBIC INTERPOLATION ********** ! ------------------------------------------------------------------ 170 a4 = 2.*a3 ! UPDATE DESIGN VECTOR AND FUNCTION VALUE. ap = a4 - alp alp = a4 x(1:me%ndv) = x(1:me%ndv) + ap*s(1:me%ndv) if (me%iprint > 4) write (me%iunit, 5100) alp if (me%iprint > 4) write (me%iunit, 5200) (x(i), i=1, me%ndv) ncal(1) = ncal(1) + 1 jgoto = 6 return 190 f4 = me%obj if (me%iprint > 4) write (me%iunit, 5300) f4 if (f4 <= f3) then a1 = a2 f1 = f2 a2 = a3 f2 = f3 a3 = a4 f3 = f4 go to 170 end if 200 ii = 4 call me%cnmn04(ii, app, a1, a1, f1, slope, a2, f2, a3, f3, a4, f4) if (app <= a1) then ap = a1 - alp alp = a1 me%obj = f1 do i = 1, me%ndv x(i) = x(i) + ap*s(i) end do go to 240 end if ! ------------------------------------------------------------------ ! UPDATE DESIGN VECTOR AND FUNCTION VALUE ! ------------------------------------------------------------------ ap = app - alp alp = app x(1:me%ndv) = x(1:me%ndv) + ap*s(1:me%ndv) if (me%iprint > 4) write (me%iunit, 5100) alp if (me%iprint > 4) write (me%iunit, 5200) x(1:me%ndv) ncal(1) = ncal(1) + 1 jgoto = 7 return 230 if (me%iprint > 4) write (me%iunit, 5300) me%obj ! ------------------------------------------------------------------ ! CHECK FOR ILL-CONDITIONING ! ------------------------------------------------------------------ 240 if (me%obj <= f2 .and. me%obj <= f3) then if (me%obj <= f1) go to 260 ap = a1 - alp alp = a1 me%obj = f1 else if (f2 >= f3) then me%obj = f3 ap = a3 - alp alp = a3 else me%obj = f2 ap = a2 - alp alp = a2 end if end if ! ------------------------------------------------------------------ ! UPDATE DESIGN VECTOR ! ------------------------------------------------------------------ x(1:me%ndv) = x(1:me%ndv) + ap*s(1:me%ndv) ! ------------------------------------------------------------------ ! CHECK FOR MULTIPLE MINIMA ! ------------------------------------------------------------------ 260 if (me%obj > fff) then ! INITIAL FUNCTION IS MINIMUM. x(1:me%ndv) = x(1:me%ndv) - alp*s(1:me%ndv) alp = 0.0_wp me%obj = fff end if jgoto = 0 return ! ------------------------------------------------------------------ ! FORMATS ! ------------------------------------------------------------------ 5000 format(///t6, & '* * * UNCONSTRAINED ONE-DIMENSIONAL SEARCH INFORMATION * * *') 5100 format(/t6, 'ALPHA =', e14.5/t6, 'X-VECTOR') 5200 format(t6, 6e13.5) 5300 format(/t6, 'OBJ =', e14.5) end subroutine cnmn03