Coverage for /usr/lib/python3/dist-packages/scipy/optimize/_direct_py.py: 19%

58 statements  

« prev     ^ index     » next       coverage.py v7.9.1, created at 2025-06-14 15:55 +0200

1from __future__ import annotations 

2from typing import ( # noqa: UP035 

3 Any, Callable, Iterable, TYPE_CHECKING 

4) 

5 

6import numpy as np 

7from scipy.optimize import OptimizeResult 

8from ._constraints import old_bound_to_new, Bounds 

9from ._direct import direct as _direct # type: ignore 

10 

11if TYPE_CHECKING: 

12 import numpy.typing as npt 

13 

14__all__ = ['direct'] 

15 

16ERROR_MESSAGES = ( 

17 "Number of function evaluations done is larger than maxfun={}", 

18 "Number of iterations is larger than maxiter={}", 

19 "u[i] < l[i] for some i", 

20 "maxfun is too large", 

21 "Initialization failed", 

22 "There was an error in the creation of the sample points", 

23 "An error occurred while the function was sampled", 

24 "Maximum number of levels has been reached.", 

25 "Forced stop", 

26 "Invalid arguments", 

27 "Out of memory", 

28) 

29 

30SUCCESS_MESSAGES = ( 

31 ("The best function value found is within a relative error={} " 

32 "of the (known) global optimum f_min"), 

33 ("The volume of the hyperrectangle containing the lowest function value " 

34 "found is below vol_tol={}"), 

35 ("The side length measure of the hyperrectangle containing the lowest " 

36 "function value found is below len_tol={}"), 

37) 

38 

39 

40def direct( 

41 func: Callable[[npt.ArrayLike, tuple[Any]], float], 

42 bounds: Iterable | Bounds, 

43 *, 

44 args: tuple = (), 

45 eps: float = 1e-4, 

46 maxfun: int | None = None, 

47 maxiter: int = 1000, 

48 locally_biased: bool = True, 

49 f_min: float = -np.inf, 

50 f_min_rtol: float = 1e-4, 

51 vol_tol: float = 1e-16, 

52 len_tol: float = 1e-6, 

53 callback: Callable[[npt.ArrayLike], None] | None = None 

54) -> OptimizeResult: 

55 """ 

56 Finds the global minimum of a function using the 

57 DIRECT algorithm. 

58 

59 Parameters 

60 ---------- 

61 func : callable 

62 The objective function to be minimized. 

63 ``func(x, *args) -> float`` 

64 where ``x`` is an 1-D array with shape (n,) and ``args`` is a tuple of 

65 the fixed parameters needed to completely specify the function. 

66 bounds : sequence or `Bounds` 

67 Bounds for variables. There are two ways to specify the bounds: 

68 

69 1. Instance of `Bounds` class. 

70 2. ``(min, max)`` pairs for each element in ``x``. 

71 

72 args : tuple, optional 

73 Any additional fixed parameters needed to 

74 completely specify the objective function. 

75 eps : float, optional 

76 Minimal required difference of the objective function values 

77 between the current best hyperrectangle and the next potentially 

78 optimal hyperrectangle to be divided. In consequence, `eps` serves as a 

79 tradeoff between local and global search: the smaller, the more local 

80 the search becomes. Default is 1e-4. 

81 maxfun : int or None, optional 

82 Approximate upper bound on objective function evaluations. 

83 If `None`, will be automatically set to ``1000 * N`` where ``N`` 

84 represents the number of dimensions. Will be capped if necessary to 

85 limit DIRECT's RAM usage to app. 1GiB. This will only occur for very 

86 high dimensional problems and excessive `max_fun`. Default is `None`. 

87 maxiter : int, optional 

88 Maximum number of iterations. Default is 1000. 

89 locally_biased : bool, optional 

90 If `True` (default), use the locally biased variant of the 

91 algorithm known as DIRECT_L. If `False`, use the original unbiased 

92 DIRECT algorithm. For hard problems with many local minima, 

93 `False` is recommended. 

94 f_min : float, optional 

95 Function value of the global optimum. Set this value only if the 

96 global optimum is known. Default is ``-np.inf``, so that this 

97 termination criterion is deactivated. 

98 f_min_rtol : float, optional 

99 Terminate the optimization once the relative error between the 

100 current best minimum `f` and the supplied global minimum `f_min` 

101 is smaller than `f_min_rtol`. This parameter is only used if 

102 `f_min` is also set. Must lie between 0 and 1. Default is 1e-4. 

103 vol_tol : float, optional 

104 Terminate the optimization once the volume of the hyperrectangle 

105 containing the lowest function value is smaller than `vol_tol` 

106 of the complete search space. Must lie between 0 and 1. 

107 Default is 1e-16. 

108 len_tol : float, optional 

109 If `locally_biased=True`, terminate the optimization once half of 

110 the normalized maximal side length of the hyperrectangle containing 

111 the lowest function value is smaller than `len_tol`. 

112 If `locally_biased=False`, terminate the optimization once half of 

113 the normalized diagonal of the hyperrectangle containing the lowest 

114 function value is smaller than `len_tol`. Must lie between 0 and 1. 

115 Default is 1e-6. 

116 callback : callable, optional 

117 A callback function with signature ``callback(xk)`` where ``xk`` 

118 represents the best function value found so far. 

119 

120 Returns 

121 ------- 

122 res : OptimizeResult 

123 The optimization result represented as a ``OptimizeResult`` object. 

124 Important attributes are: ``x`` the solution array, ``success`` a 

125 Boolean flag indicating if the optimizer exited successfully and 

126 ``message`` which describes the cause of the termination. See 

127 `OptimizeResult` for a description of other attributes. 

128 

129 Notes 

130 ----- 

131 DIviding RECTangles (DIRECT) is a deterministic global 

132 optimization algorithm capable of minimizing a black box function with 

133 its variables subject to lower and upper bound constraints by sampling 

134 potential solutions in the search space [1]_. The algorithm starts by 

135 normalising the search space to an n-dimensional unit hypercube. 

136 It samples the function at the center of this hypercube and at 2n 

137 (n is the number of variables) more points, 2 in each coordinate 

138 direction. Using these function values, DIRECT then divides the 

139 domain into hyperrectangles, each having exactly one of the sampling 

140 points as its center. In each iteration, DIRECT chooses, using the `eps` 

141 parameter which defaults to 1e-4, some of the existing hyperrectangles 

142 to be further divided. This division process continues until either the 

143 maximum number of iterations or maximum function evaluations allowed 

144 are exceeded, or the hyperrectangle containing the minimal value found 

145 so far becomes small enough. If `f_min` is specified, the optimization 

146 will stop once this function value is reached within a relative tolerance. 

147 The locally biased variant of DIRECT (originally called DIRECT_L) [2]_ is 

148 used by default. It makes the search more locally biased and more 

149 efficient for cases with only a few local minima. 

150 

151 A note about termination criteria: `vol_tol` refers to the volume of the 

152 hyperrectangle containing the lowest function value found so far. This 

153 volume decreases exponentially with increasing dimensionality of the 

154 problem. Therefore `vol_tol` should be decreased to avoid premature 

155 termination of the algorithm for higher dimensions. This does not hold 

156 for `len_tol`: it refers either to half of the maximal side length 

157 (for ``locally_biased=True``) or half of the diagonal of the 

158 hyperrectangle (for ``locally_biased=False``). 

159 

160 This code is based on the DIRECT 2.0.4 Fortran code by Gablonsky et al. at 

161 https://ctk.math.ncsu.edu/SOFTWARE/DIRECTv204.tar.gz . 

162 This original version was initially converted via f2c and then cleaned up 

163 and reorganized by Steven G. Johnson, August 2007, for the NLopt project. 

164 The `direct` function wraps the C implementation. 

165 

166 .. versionadded:: 1.9.0 

167 

168 References 

169 ---------- 

170 .. [1] Jones, D.R., Perttunen, C.D. & Stuckman, B.E. Lipschitzian 

171 optimization without the Lipschitz constant. J Optim Theory Appl 

172 79, 157-181 (1993). 

173 .. [2] Gablonsky, J., Kelley, C. A Locally-Biased form of the DIRECT 

174 Algorithm. Journal of Global Optimization 21, 27-37 (2001). 

175 

176 Examples 

177 -------- 

178 The following example is a 2-D problem with four local minima: minimizing 

179 the Styblinski-Tang function 

180 (https://en.wikipedia.org/wiki/Test_functions_for_optimization). 

181 

182 >>> from scipy.optimize import direct, Bounds 

183 >>> def styblinski_tang(pos): 

184 ... x, y = pos 

185 ... return 0.5 * (x**4 - 16*x**2 + 5*x + y**4 - 16*y**2 + 5*y) 

186 >>> bounds = Bounds([-4., -4.], [4., 4.]) 

187 >>> result = direct(styblinski_tang, bounds) 

188 >>> result.x, result.fun, result.nfev 

189 array([-2.90321597, -2.90321597]), -78.3323279095383, 2011 

190 

191 The correct global minimum was found but with a huge number of function 

192 evaluations (2011). Loosening the termination tolerances `vol_tol` and 

193 `len_tol` can be used to stop DIRECT earlier. 

194 

195 >>> result = direct(styblinski_tang, bounds, len_tol=1e-3) 

196 >>> result.x, result.fun, result.nfev 

197 array([-2.9044353, -2.9044353]), -78.33230330754142, 207 

198 

199 """ 

200 # convert bounds to new Bounds class if necessary 

201 if not isinstance(bounds, Bounds): 

202 if isinstance(bounds, list) or isinstance(bounds, tuple): 

203 lb, ub = old_bound_to_new(bounds) 

204 bounds = Bounds(lb, ub) 

205 else: 

206 message = ("bounds must be a sequence or " 

207 "instance of Bounds class") 

208 raise ValueError(message) 

209 

210 lb = np.ascontiguousarray(bounds.lb, dtype=np.float64) 

211 ub = np.ascontiguousarray(bounds.ub, dtype=np.float64) 

212 

213 # validate bounds 

214 # check that lower bounds are smaller than upper bounds 

215 if not np.all(lb < ub): 

216 raise ValueError('Bounds are not consistent min < max') 

217 # check for infs 

218 if (np.any(np.isinf(lb)) or np.any(np.isinf(ub))): 

219 raise ValueError("Bounds must not be inf.") 

220 

221 # validate tolerances 

222 if (vol_tol < 0 or vol_tol > 1): 

223 raise ValueError("vol_tol must be between 0 and 1.") 

224 if (len_tol < 0 or len_tol > 1): 

225 raise ValueError("len_tol must be between 0 and 1.") 

226 if (f_min_rtol < 0 or f_min_rtol > 1): 

227 raise ValueError("f_min_rtol must be between 0 and 1.") 

228 

229 # validate maxfun and maxiter 

230 if maxfun is None: 

231 maxfun = 1000 * lb.shape[0] 

232 if not isinstance(maxfun, int): 

233 raise ValueError("maxfun must be of type int.") 

234 if maxfun < 0: 

235 raise ValueError("maxfun must be > 0.") 

236 if not isinstance(maxiter, int): 

237 raise ValueError("maxiter must be of type int.") 

238 if maxiter < 0: 

239 raise ValueError("maxiter must be > 0.") 

240 

241 # validate boolean parameters 

242 if not isinstance(locally_biased, bool): 

243 raise ValueError("locally_biased must be True or False.") 

244 

245 def _func_wrap(x, args=None): 

246 x = np.asarray(x) 

247 if args is None: 

248 f = func(x) 

249 else: 

250 f = func(x, *args) 

251 # always return a float 

252 return np.asarray(f).item() 

253 

254 # TODO: fix disp argument 

255 x, fun, ret_code, nfev, nit = _direct( 

256 _func_wrap, 

257 np.asarray(lb), np.asarray(ub), 

258 args, 

259 False, eps, maxfun, maxiter, 

260 locally_biased, 

261 f_min, f_min_rtol, 

262 vol_tol, len_tol, callback 

263 ) 

264 

265 format_val = (maxfun, maxiter, f_min_rtol, vol_tol, len_tol) 

266 if ret_code > 2: 

267 message = SUCCESS_MESSAGES[ret_code - 3].format( 

268 format_val[ret_code - 1]) 

269 elif 0 < ret_code <= 2: 

270 message = ERROR_MESSAGES[ret_code - 1].format(format_val[ret_code - 1]) 

271 elif 0 > ret_code > -100: 

272 message = ERROR_MESSAGES[abs(ret_code) + 1] 

273 else: 

274 message = ERROR_MESSAGES[ret_code + 99] 

275 

276 return OptimizeResult(x=np.asarray(x), fun=fun, status=ret_code, 

277 success=ret_code > 2, message=message, 

278 nfev=nfev, nit=nit)