BCPOL
Minimizes a function of N variables subject to bounds on the variables using a direct search complex algorithm.
Required Arguments
FCN — User-supplied subroutine to evaluate the function to be minimized. The usage is CALL FCN (N, X, F), where
N – Length of X. (Input)
X – Vector of length N at which point the function is evaluated. (Input)
X should not be changed by FCN.
F – The computed function value at the point X. (Output)
FCN must be declared EXTERNAL in the calling program.
IBTYPE — Scalar indicating the types of bounds on variables. (Input)
|
IBTYPE |
Action |
|
0 |
User will supply all the bounds. |
|
1 |
All variables are nonnegative. |
|
2 |
All variables are nonpositive. |
|
3 |
User supplies only the bounds on first variable. All other variables will have the same bounds. |
XLB — Vector of length N containing the lower bounds on the variables. (Input, if
IBTYPE = 0; output, if IBTYPE = 1 or 2; input/output, if IBTYPE = 3)
XUB — Vector of length N containing the upper bounds on the variables. (Input, if
IBTYPE = 0; output, if IBTYPE = 1 or 2; input/output, if IBTYPE = 3)
X — Real vector of length N containing the best estimate of the minimum found. (Output)
Optional Arguments
N — The number of variables. (Input)
Default: N = SIZE (X,1).
XGUESS — Real vector of length N that contains an initial guess to the minimum. (Input)
Default: XGUESS = 0.0.
FTOL — The convergence tolerance. (Input)
First convergence criterion: The algorithm stops when the relative error in the function values is less than FTOL, i.e. when F(worst) ‑ F(best) < FTOL * (1 + ABS(F(best))) where F(worst) and F(best) are the function values of the current worst and best point, respectively.
Second convergence criterion: The algorithm stops when the standard deviation of the function values at the 2 * N current points is less than FTOL. If the subroutine terminates prematurely, try again with a smaller value FTOL.
Default: FTOL = 1.0e-4 for single and 1.0d-8 for double precision.
MAXFCN — The number of function evaluations. (Input/Output)
On input, maximum allowed number of function evaluations. On output, actual number of function evaluations needed.
Default: MAXFCN = 300.
FVALUE — Function value at the computed solution. (Output)
XCPLX — Real N x 2*N matrix containing the 2 * N initial complex points. (Input)
If the XCPLX optional argument is used, then the initial guess must be specified in XCPLX(:,1) and there is no need to provide optional argument XGUESS. Thus, if both optional arguments XCPLX and XGUESS are supplied, XGUESS will be ignored.
REFLCOEF — Real scalar containing the reflection coefficient. (Input)
REFLCOEF must be greater than 0.
Default: REFLCOEF = 1.0.
EXPNCOEF — Real scalar containing the expansion coefficient. (Input)
EXPNCOEF must be greater than 1.
Default: EXPNCOEF = 2.0.
CNTRCOEF — Real scalar containing the contraction coefficient. (Input)
CNTRCOEF must be greater than 0 and less than 1.
Default: CNTRCOEF = 0.5.
FORTRAN 90 Interface
Generic: CALL BCPOL (FCN, IBTYPE, XLB, XUB, X [, …])
Specific: The specific interface names are S_BCPOL and D_BCPOL.
FORTRAN 77 Interface
Single: CALL BCPOL (FCN, N, XGUESS, IBTYPE, XLB, XUB, FTOL, MAXFCN, X, FVALUE)
Double: The double precision name is DBCPOL.
Note: To maintain backward compatibility, the argument list for the Fortran 77 interface does not include the optional arguments XCPLX, REFLCOEF, EXPNCOEF, or CNTRCOEF, and therefore cannot be used to specify these arguments.
Description
The routine BCPOL uses the complex method to find a minimum point of a function of n variables. The method is based on function comparison; no smoothness is assumed. It starts with an initial complex of 2n points x1, x2, …, x2n which is either user‑specified (using optional argument XCPLX, see above) or is otherwise, by default, randomly initialized. At each iteration, a new point is generated to replace the “worst” point xj, which has the largest function value among these 2n points, as described below.
Each iteration begins by determining the best and two worst points in the present complex, and then constructing a new “reflection” point xr by the formula
xr = c + α(c ‑ xj)
where
is the centroid of the 2n - 1 best points and α(α > 0) is the reflection coefficient. (See optional argument REFLCOEF above). Depending on how the new reflection point xr compares with the existing complex points, the iteration proceeds as follows:
If xr is neither better than the best point nor worse than the second worst point, then xr replaces the worst point xj and, if neither of the stopping criteria (see below) are satisfied, a new iteration begins.
If xr is a best point, that is, if f(xr) ≤ f(xi) for i = 1, …, 2n, an expansion point xe is computed to see if an even better point can be obtained by moving further in the reflection direction:
xe = c + β(xr ‑ c)
where β(β > 1) is called the expansion coefficient (see optional argument EXPNCOEF above), and worst point xj is replaced by the better of xe and xr and, if neither of the stopping criteria are satisfied, a new iteration begins.
If xr and xj are the two worst points, then the complex is contracted to try to get a better new contraction point xc :
where 
(0 < < 1) is called the contraction coefficient. (See optional argument CNTRCOEF above.) If the contraction step is successful (i.e. if min( f(xr), f(xj) ) > f(xc)), then worst point xj is replaced by xc. If the contraction step is unsuccessful, then the complex is shrunk by moving the 2n ‑ 1 worst points halfway towards the current best point. Following the contraction step, if neither of the stopping criteria are satisfied, a new iteration begins.
Whenever the new point generated is beyond the bound, it will be set to the bound. If, at the end of an iteration, one of the following stopping criteria is satisfied, then the process ends with the best point returned as the optimum.
Criterion 1:
fworst‑ fbest ≤ ɛf (1. + |fbest|)
Criterion 2:
where fi = f(xi), fj = f(xj), and ɛf is a given tolerance. For a complete description, see Nelder and Mead (1965) or Gill et al. (1981).
Comments
1. Workspace may be explicitly provided, if desired, by use of B2POL/DB2POL. The reference is:
CALL B2POL (FCN, N, XGUESS, IBTYPE, XLB, XUB, FTOL, MAXFCN, X, FVALUE, WK)
The additional argument is:
WK — Real work vector of length 2 * N**2 + 5 * N
Note: To maintain backward compatibility, the argument list for B2POL/DB2POL does not include the optional arguments XCPLX, REFLCOEF, EXPNCOEF, or CNTRCOEF, and therefore cannot be used to specify these arguments.
2. Informational error
|
Type |
Code |
Description |
|
3 |
1 |
The maximum number of function evaluations is exceeded. |
3. Since BCPOL uses only function-value information at each step to determine a new approximate minimum, it could be quite inefficient on smooth problems compared to other methods such as those implemented in routine BCONF, which takes into account derivative information at each iteration. Hence, routine BCPOL should only be used as a last resort. Briefly, a set of 2 * N points in an N-dimensional space is called a complex. The minimization process iterates by replacing the point with the largest function value by a new point with a smaller function value. The iteration continues until all the points cluster sufficiently close to a minimum.
Example 1
The problem
is solved with an initial guess (1.2, 1.0), and the solution is printed.
USE BCPOL_INT
USE UMACH_INT
IMPLICIT NONE
! Variable declarations
INTEGER N
PARAMETER (N=2)
!
INTEGER IBTYPE, K, NOUT
REAL FTOL, FVALUE, X(N), XGUESS(N), XLB(N), XUB(N)
EXTERNAL FCN
!
! Initializations
! XGUESS = (-1.2, 1.0)
! XLB = (-2.0, -1.0)
! XUB = ( 0.5, 2.0)
DATA XGUESS/-1.2, 1.0/, XLB/-2.0E0, -1.0E0/, XUB/0.5E0, 2.0E0/
!
FTOL = 1.0E-5
IBTYPE = 0
!
CALL BCPOL (FCN, IBTYPE, XLB, XUB, X, xguess=xguess, ftol=ftol, &
fvalue=fvalue)
!
CALL UMACH (2, NOUT)
WRITE (NOUT,99999) (X(K),K=1,N), FVALUE
99999 FORMAT (' The best estimate for the minimum value of the', /, &
' function is X = (', 2(2X,F4.2), ')', /, ' with ', &
'function value FVALUE = ', E12.6)
!
END
! External function to be minimized
SUBROUTINE FCN (N, X, F)
INTEGER N
REAL X(N), F
!
F = 100.0*(X(2)-X(1)*X(1))**2 + (1.0-X(1))**2
RETURN
END
The best estimate for the minimum value of the
function is X = ( 0.50 0.25)
with function value FVALUE = 0.250002E+00
Example 2
This example is intended to illustrate the use of optional arguments available for BCPOL, and how their use can affect the number of function calls needed to complete the optimization process. The same problem as in Example 1 is approached in five ways, including the use of a function penalty in an attempt to constrain the solution space. In each case the number of function evaluations required is output.
Solution 1 uses XLB and XUB to impose bounds, as in Example 1. Solutions 2 through 5 use a function penalty to impose constraints, and to varying degrees make use of optional arguments XCPLX, REFLCOEF, EXPNCOEF, and CNTRCOEF.
Note, the actual number of function calls required to complete the minimization process may vary, depending on the computing platform and precision used.
!
! BCPOL_EX_2.F90
USE BCPOL_INT
USE UMACH_INT
IMPLICIT NONE
!
integer :: IBTYPE, MAXFCN, N
integer, parameter :: NNN=2
real(kind(1d0)) :: FTOL
integer :: K, NOUT
real(kind(1d0)) :: FVALUE
integer :: ICTYPE
integer :: I, II, JJ, IARGTST
real(kind(1d0)) :: X(NNN), XLB(NNN), XUB(NNN)
real(kind(1d0)) :: XCPLX0(NNN,2*NNN), XCPLX(NNN,2*NNN)
real(kind(1d0)) :: ALFA, BETA, DGAMMA
!
! SPECIFICATIONS FOR SUBROUTINES
external FCN, FCN2
real(kind(1d0)), dimension(2) :: XGUESS = (/ -1.2d0, 1.0d0 /)
real(kind(1d0)), dimension(2) :: XCN = (/ -2.0d0, -1.0d0 /), XCX = (/0.5d0, 2.0d0/)
call umach (2, nout)
!
IBTYPE = 0
FTOL = 1.0d-15
N = 2
MAXFCN = 500
! Constraints imposed with Bounds
XLB(1) = -2.0d0
XLB(2) = -1.0d0
XUB(1) = 0.5d0
XUB(2) = 2.0d0
!
WRITE (NOUT,'(/" Solution 1:")')
WRITE (NOUT,'(/" Bounds: XLB = (-2.0, -1.0); XUB = (0.5, 2.0)")')
WRITE (NOUT,'( " Constraints imposed with Bounds")')
WRITE (NOUT,'( " using no optional complex or step-size arguments:")')
!
CALL BCPOL (FCN, IBTYPE, XLB, XUB, X, n=N, xguess=XGUESS, &
ftol=FTOL, maxfcn=MAXFCN, fvalue=FVALUE)
!
WRITE (NOUT,99999) (X(K),K=1,N), FVALUE, MAXFCN
!
! Lower and Upper Bounds:
! -2. <= X(1) <= 2. ; -2. <= X(2) <= 2.;
! Constraints:
! -2. <= X(1) <= 0.5; -1. <= X(2) <= 2.;
!
! Constraints imposed with function penalties
! User-supplied initial complex XCPLX
! chosen to be within constraint boundaries
XLB(1) = -2.0d0
XUB(1) = 2.0d0
XLB(2) = -2.0d0
XUB(2) = 2.0d0
!
XCPLX0(1,1) = XGUESS(1)
XCPLX0(2,1) = XGUESS(2)
XCPLX0(1,2) = XCX(1)
XCPLX0(2,2) = XCX(2)
XCPLX0(1,3) = XCN(1)
XCPLX0(2,3) = XCN(2)
XCPLX0(1,4) = XCX(1)
XCPLX0(2,4) = XCN(2)
!
IBTYPE = 0
ALFA = 1.0d0
!
WRITE (NOUT,'(/"----------------------------------------------"/)')
WRITE (NOUT,'(" Solutions 2 - 5:"/)')
WRITE (NOUT,'(" Constraints: -2. <= x(1) <= 0.5; -1. <= x(2) <= 2.")')
WRITE (NOUT,'(" Bounds: XLB = (-2.0, -2.0); XUB = (2.0, 2.0)")')
WRITE (NOUT,'(" Constraints imposed with function penalties")')
WRITE (NOUT,'(" using optional complex and step-size arguments:")')
!
!-----------------------------------------------------------------------
!
WRITE (NOUT,'(/"----------------------------------------------"/)')
WRITE (NOUT,'( " Solutions 2 - 3:"/)')
WRITE (NOUT,'( " User Specified Initial Complex:")')
WRITE (NOUT,'(/"----------------------------------------------"/)')
WRITE (NOUT,'( " Solution 2: step-size coefficients with default values:"/)')
!
DO JJ=1,2*N
DO II=1,N
XCPLX(II,JJ) = XCPLX0(II,JJ)
END DO
END DO
MAXFCN = 10000
BETA = 2.0D+01
DGAMMA = 0.5D+00
!
WRITE (NOUT,99998) ALFA, BETA, DGAMMA
!
CALL BCPOL (FCN2, IBTYPE, XLB, XUB, X, &
xguess=XGUESS, ftol=FTOL, maxfcn=MAXFCN, &
fvalue=FVALUE, xcplx=XCPLX)
!
WRITE (NOUT,99999) (X(K),K=1,N), FVALUE, MAXFCN
!
!-----------------------------------------------------------------------
!
WRITE (NOUT, '(/"----------------------------------------------"/)')
WRITE (NOUT,'(" Solution 3: step-size coefficients with user-specified values:"/)')
!
DO JJ=1,2*N
DO II=1,N
XCPLX(II,JJ) = XCPLX0(II,JJ)
END DO
END DO
MAXFCN = 10000
!
BETA = 0.31841776469083554D+01
DGAMMA = 0.33464404002126491D+00
!
WRITE (NOUT,99998) ALFA, BETA, DGAMMA
!
CALL BCPOL (FCN2, IBTYPE, XLB, XUB, X, &
xguess=XGUESS, ftol=FTOL, maxfcn=MAXFCN, &
fvalue=FVALUE, xcplx=XCPLX, &
reflcoef=ALFA, expncoef=BETA, cntrcoef=DGAMMA)
!
WRITE (NOUT,99999) (X(K),K=1,N), FVALUE, MAXFCN
!
!-----------------------------------------------------------------------
!
WRITE (NOUT,'(/"----------------------------------------------"/)')
WRITE (NOUT,'( " Solutions 4 - 5:"/)')
WRITE (NOUT,'( " Randomly Generated Initial Complex:")')
WRITE (NOUT,'(/"----------------------------------------------"/)')
WRITE (NOUT,'( " Solution 4: step-size coefficients with default values:"/)')
!
MAXFCN = 10000
BETA = 2.0D+01
DGAMMA = 0.5D+00
!
WRITE (NOUT,99998) ALFA, BETA, DGAMMA
!
CALL BCPOL (FCN2, IBTYPE, XLB, XUB, X, &
xguess=XGUESS, ftol=FTOL, maxfcn=MAXFCN, &
fvalue=FVALUE)
!
WRITE (NOUT,99999) (X(K),K=1,N), FVALUE, MAXFCN
!
!-----------------------------------------------------------------------
!
WRITE (NOUT, '(/"----------------------------------------------"/)')
WRITE (NOUT,'(" Solution 5: step-size coefficients with user-specified values:"/)')
!
MAXFCN = 10000
!
BETA = 0.18204845270362373D+02
DGAMMA = 0.31542073037934792D+00
!
WRITE (NOUT,99998) ALFA, BETA, DGAMMA
!
CALL BCPOL (FCN2, IBTYPE, XLB, XUB, X, &
xguess=XGUESS, ftol=FTOL, maxfcn=MAXFCN, &
fvalue=FVALUE, &
reflcoef=ALFA, expncoef=BETA, cntrcoef=DGAMMA)
!
WRITE (NOUT,99999) (X(K),K=1,N), FVALUE, MAXFCN
!
!-----------------------------------------------------------------------
!
99998 FORMAT (' REFLCOEF = ', D25.17 / &
' EXPNCOEF = ', D25.17 / &
' CNTRCOEF = ', D25.17 )
!
99999 FORMAT (/' The best estimate for the minimum value of the' / &
' function is X = ( ', d14.7,' , ', d14.7, ' )' / &
' with function value FVALUE = ', d14.7 / &
' and # fcn calls MAXFCN = ', I7 )
end
!
!-----------------------------------------------------------------------
!
SUBROUTINE FCN (N, X, F)
integer :: N
real(kind(1d0)) :: X(N), F
!
F = 100.0*(X(2)-X(1)*X(1))**2 + (1.0-X(1))**2
!
RETURN
END
!
!-----------------------------------------------------------------------
!
SUBROUTINE FCN2 (N, X, F)
integer :: N
real(kind(1d0)) :: X(N), F
!
IF ( ((-2.0 .LE. X(1)) .AND. (X(1) .LE. 0.5)) .AND. &
((-1.0 .LE. X(2)) .AND. (X(2) .LE. 2.0))) THEN
F = 100.0*(X(2)-X(1)*X(1))**2 + (1.0-X(1))**2
ELSE
F = 1000.
END IF
!
RETURN
END
Solution 1:
Bounds: XLB = (-2.0, -1.0); XUB = (0.5, 2.0)
Constraints imposed with Bounds
using no optional complex or step-size arguments:
The best estimate for the minimum value of the
function is X = ( 0.5000000D+00 , 0.2500000D+00 )
with function value FVALUE = 0.2500000D+00
and # fcn calls MAXFCN = 226
----------------------------------------------
Solutions 2 - 5:
Constraints: -2. <= x(1) <= 0.5; -1. <= x(2) <= 2.
Bounds: XLB = (-2.0, -2.0); XUB = (2.0, 2.0)
Constraints imposed with function penalties
using optional complex and step-size arguments:
----------------------------------------------
Solutions 2 - 3:
User Specified Initial Complex:
----------------------------------------------
Solution 2: step-size coefficients with default values:
REFLCOEF = 0.10000000000000000D+01
EXPNCOEF = 0.20000000000000000D+02
CNTRCOEF = 0.50000000000000000D+00
The best estimate for the minimum value of the
function is X = ( 0.5000000D+00 , 0.2500000D+00 )
with function value FVALUE = 0.2500000D+00
and # fcn calls MAXFCN = 379
----------------------------------------------
Solution 3: step-size coefficients with user-specified values:
REFLCOEF = 0.10000000000000000D+01
EXPNCOEF = 0.31841776469083554D+01
CNTRCOEF = 0.33464404002126491D+00
The best estimate for the minimum value of the
function is X = ( 0.5000000D+00 , 0.2500000D+00 )
with function value FVALUE = 0.2500000D+00
and # fcn calls MAXFCN = 323
----------------------------------------------
Solutions 4 - 5:
Randomly Generated Initial Complex:
----------------------------------------------
Solution 4: step-size coefficients with default values:
REFLCOEF = 0.10000000000000000D+01
EXPNCOEF = 0.20000000000000000D+02
CNTRCOEF = 0.50000000000000000D+00
The best estimate for the minimum value of the
function is X = ( 0.5000000D+00 , 0.2500000D+00 )
with function value FVALUE = 0.2500000D+00
and # fcn calls MAXFCN = 430
----------------------------------------------
Solution 5: step-size coefficients with user-specified values:
REFLCOEF = 0.10000000000000000D+01
EXPNCOEF = 0.18204845270362373D+02
CNTRCOEF = 0.31542073037934792D+00
The best estimate for the minimum value of the
function is X = ( 0.5000000D+00 , 0.2500000D+00 )
with function value FVALUE = 0.2500000D+00
and # fcn calls MAXFCN = 294
