_tnc.py 17 KB

123456789101112131415161718192021222324252627282930313233343536373839404142434445464748495051525354555657585960616263646566676869707172737475767778798081828384858687888990919293949596979899100101102103104105106107108109110111112113114115116117118119120121122123124125126127128129130131132133134135136137138139140141142143144145146147148149150151152153154155156157158159160161162163164165166167168169170171172173174175176177178179180181182183184185186187188189190191192193194195196197198199200201202203204205206207208209210211212213214215216217218219220221222223224225226227228229230231232233234235236237238239240241242243244245246247248249250251252253254255256257258259260261262263264265266267268269270271272273274275276277278279280281282283284285286287288289290291292293294295296297298299300301302303304305306307308309310311312313314315316317318319320321322323324325326327328329330331332333334335336337338339340341342343344345346347348349350351352353354355356357358359360361362363364365366367368369370371372373374375376377378379380381382383384385386387388389390391392393394395396397398399400401402403404405406407408409410411412413414415416417418419420421422423424425426427428429430431432433434435436437438439440441
  1. # TNC Python interface
  2. # @(#) $Jeannot: tnc.py,v 1.11 2005/01/28 18:27:31 js Exp $
  3. # Copyright (c) 2004-2005, Jean-Sebastien Roy (js@jeannot.org)
  4. # Permission is hereby granted, free of charge, to any person obtaining a
  5. # copy of this software and associated documentation files (the
  6. # "Software"), to deal in the Software without restriction, including
  7. # without limitation the rights to use, copy, modify, merge, publish,
  8. # distribute, sublicense, and/or sell copies of the Software, and to
  9. # permit persons to whom the Software is furnished to do so, subject to
  10. # the following conditions:
  11. # The above copyright notice and this permission notice shall be included
  12. # in all copies or substantial portions of the Software.
  13. # THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS
  14. # OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
  15. # MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT.
  16. # IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY
  17. # CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION OF CONTRACT,
  18. # TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION WITH THE
  19. # SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
  20. """
  21. TNC: A Python interface to the TNC non-linear optimizer
  22. TNC is a non-linear optimizer. To use it, you must provide a function to
  23. minimize. The function must take one argument: the list of coordinates where to
  24. evaluate the function; and it must return either a tuple, whose first element is the
  25. value of the function, and whose second argument is the gradient of the function
  26. (as a list of values); or None, to abort the minimization.
  27. """
  28. import warnings
  29. from scipy.optimize import _moduleTNC as moduleTNC
  30. from ._optimize import (MemoizeJac, OptimizeResult, _check_unknown_options,
  31. _prepare_scalar_function)
  32. from ._constraints import old_bound_to_new
  33. from numpy import inf, array, zeros, asfarray
  34. __all__ = ['fmin_tnc']
  35. MSG_NONE = 0 # No messages
  36. MSG_ITER = 1 # One line per iteration
  37. MSG_INFO = 2 # Informational messages
  38. MSG_VERS = 4 # Version info
  39. MSG_EXIT = 8 # Exit reasons
  40. MSG_ALL = MSG_ITER + MSG_INFO + MSG_VERS + MSG_EXIT
  41. MSGS = {
  42. MSG_NONE: "No messages",
  43. MSG_ITER: "One line per iteration",
  44. MSG_INFO: "Informational messages",
  45. MSG_VERS: "Version info",
  46. MSG_EXIT: "Exit reasons",
  47. MSG_ALL: "All messages"
  48. }
  49. INFEASIBLE = -1 # Infeasible (lower bound > upper bound)
  50. LOCALMINIMUM = 0 # Local minimum reached (|pg| ~= 0)
  51. FCONVERGED = 1 # Converged (|f_n-f_(n-1)| ~= 0)
  52. XCONVERGED = 2 # Converged (|x_n-x_(n-1)| ~= 0)
  53. MAXFUN = 3 # Max. number of function evaluations reached
  54. LSFAIL = 4 # Linear search failed
  55. CONSTANT = 5 # All lower bounds are equal to the upper bounds
  56. NOPROGRESS = 6 # Unable to progress
  57. USERABORT = 7 # User requested end of minimization
  58. RCSTRINGS = {
  59. INFEASIBLE: "Infeasible (lower bound > upper bound)",
  60. LOCALMINIMUM: "Local minimum reached (|pg| ~= 0)",
  61. FCONVERGED: "Converged (|f_n-f_(n-1)| ~= 0)",
  62. XCONVERGED: "Converged (|x_n-x_(n-1)| ~= 0)",
  63. MAXFUN: "Max. number of function evaluations reached",
  64. LSFAIL: "Linear search failed",
  65. CONSTANT: "All lower bounds are equal to the upper bounds",
  66. NOPROGRESS: "Unable to progress",
  67. USERABORT: "User requested end of minimization"
  68. }
  69. # Changes to interface made by Travis Oliphant, Apr. 2004 for inclusion in
  70. # SciPy
  71. def fmin_tnc(func, x0, fprime=None, args=(), approx_grad=0,
  72. bounds=None, epsilon=1e-8, scale=None, offset=None,
  73. messages=MSG_ALL, maxCGit=-1, maxfun=None, eta=-1,
  74. stepmx=0, accuracy=0, fmin=0, ftol=-1, xtol=-1, pgtol=-1,
  75. rescale=-1, disp=None, callback=None):
  76. """
  77. Minimize a function with variables subject to bounds, using
  78. gradient information in a truncated Newton algorithm. This
  79. method wraps a C implementation of the algorithm.
  80. Parameters
  81. ----------
  82. func : callable ``func(x, *args)``
  83. Function to minimize. Must do one of:
  84. 1. Return f and g, where f is the value of the function and g its
  85. gradient (a list of floats).
  86. 2. Return the function value but supply gradient function
  87. separately as `fprime`.
  88. 3. Return the function value and set ``approx_grad=True``.
  89. If the function returns None, the minimization
  90. is aborted.
  91. x0 : array_like
  92. Initial estimate of minimum.
  93. fprime : callable ``fprime(x, *args)``, optional
  94. Gradient of `func`. If None, then either `func` must return the
  95. function value and the gradient (``f,g = func(x, *args)``)
  96. or `approx_grad` must be True.
  97. args : tuple, optional
  98. Arguments to pass to function.
  99. approx_grad : bool, optional
  100. If true, approximate the gradient numerically.
  101. bounds : list, optional
  102. (min, max) pairs for each element in x0, defining the
  103. bounds on that parameter. Use None or +/-inf for one of
  104. min or max when there is no bound in that direction.
  105. epsilon : float, optional
  106. Used if approx_grad is True. The stepsize in a finite
  107. difference approximation for fprime.
  108. scale : array_like, optional
  109. Scaling factors to apply to each variable. If None, the
  110. factors are up-low for interval bounded variables and
  111. 1+|x| for the others. Defaults to None.
  112. offset : array_like, optional
  113. Value to subtract from each variable. If None, the
  114. offsets are (up+low)/2 for interval bounded variables
  115. and x for the others.
  116. messages : int, optional
  117. Bit mask used to select messages display during
  118. minimization values defined in the MSGS dict. Defaults to
  119. MGS_ALL.
  120. disp : int, optional
  121. Integer interface to messages. 0 = no message, 5 = all messages
  122. maxCGit : int, optional
  123. Maximum number of hessian*vector evaluations per main
  124. iteration. If maxCGit == 0, the direction chosen is
  125. -gradient if maxCGit < 0, maxCGit is set to
  126. max(1,min(50,n/2)). Defaults to -1.
  127. maxfun : int, optional
  128. Maximum number of function evaluation. If None, maxfun is
  129. set to max(100, 10*len(x0)). Defaults to None. Note that this function
  130. may violate the limit because of evaluating gradients by numerical
  131. differentiation.
  132. eta : float, optional
  133. Severity of the line search. If < 0 or > 1, set to 0.25.
  134. Defaults to -1.
  135. stepmx : float, optional
  136. Maximum step for the line search. May be increased during
  137. call. If too small, it will be set to 10.0. Defaults to 0.
  138. accuracy : float, optional
  139. Relative precision for finite difference calculations. If
  140. <= machine_precision, set to sqrt(machine_precision).
  141. Defaults to 0.
  142. fmin : float, optional
  143. Minimum function value estimate. Defaults to 0.
  144. ftol : float, optional
  145. Precision goal for the value of f in the stopping criterion.
  146. If ftol < 0.0, ftol is set to 0.0 defaults to -1.
  147. xtol : float, optional
  148. Precision goal for the value of x in the stopping
  149. criterion (after applying x scaling factors). If xtol <
  150. 0.0, xtol is set to sqrt(machine_precision). Defaults to
  151. -1.
  152. pgtol : float, optional
  153. Precision goal for the value of the projected gradient in
  154. the stopping criterion (after applying x scaling factors).
  155. If pgtol < 0.0, pgtol is set to 1e-2 * sqrt(accuracy).
  156. Setting it to 0.0 is not recommended. Defaults to -1.
  157. rescale : float, optional
  158. Scaling factor (in log10) used to trigger f value
  159. rescaling. If 0, rescale at each iteration. If a large
  160. value, never rescale. If < 0, rescale is set to 1.3.
  161. callback : callable, optional
  162. Called after each iteration, as callback(xk), where xk is the
  163. current parameter vector.
  164. Returns
  165. -------
  166. x : ndarray
  167. The solution.
  168. nfeval : int
  169. The number of function evaluations.
  170. rc : int
  171. Return code, see below
  172. See also
  173. --------
  174. minimize: Interface to minimization algorithms for multivariate
  175. functions. See the 'TNC' `method` in particular.
  176. Notes
  177. -----
  178. The underlying algorithm is truncated Newton, also called
  179. Newton Conjugate-Gradient. This method differs from
  180. scipy.optimize.fmin_ncg in that
  181. 1. it wraps a C implementation of the algorithm
  182. 2. it allows each variable to be given an upper and lower bound.
  183. The algorithm incorporates the bound constraints by determining
  184. the descent direction as in an unconstrained truncated Newton,
  185. but never taking a step-size large enough to leave the space
  186. of feasible x's. The algorithm keeps track of a set of
  187. currently active constraints, and ignores them when computing
  188. the minimum allowable step size. (The x's associated with the
  189. active constraint are kept fixed.) If the maximum allowable
  190. step size is zero then a new constraint is added. At the end
  191. of each iteration one of the constraints may be deemed no
  192. longer active and removed. A constraint is considered
  193. no longer active is if it is currently active
  194. but the gradient for that variable points inward from the
  195. constraint. The specific constraint removed is the one
  196. associated with the variable of largest index whose
  197. constraint is no longer active.
  198. Return codes are defined as follows::
  199. -1 : Infeasible (lower bound > upper bound)
  200. 0 : Local minimum reached (|pg| ~= 0)
  201. 1 : Converged (|f_n-f_(n-1)| ~= 0)
  202. 2 : Converged (|x_n-x_(n-1)| ~= 0)
  203. 3 : Max. number of function evaluations reached
  204. 4 : Linear search failed
  205. 5 : All lower bounds are equal to the upper bounds
  206. 6 : Unable to progress
  207. 7 : User requested end of minimization
  208. References
  209. ----------
  210. Wright S., Nocedal J. (2006), 'Numerical Optimization'
  211. Nash S.G. (1984), "Newton-Type Minimization Via the Lanczos Method",
  212. SIAM Journal of Numerical Analysis 21, pp. 770-778
  213. """
  214. # handle fprime/approx_grad
  215. if approx_grad:
  216. fun = func
  217. jac = None
  218. elif fprime is None:
  219. fun = MemoizeJac(func)
  220. jac = fun.derivative
  221. else:
  222. fun = func
  223. jac = fprime
  224. if disp is not None: # disp takes precedence over messages
  225. mesg_num = disp
  226. else:
  227. mesg_num = {0:MSG_NONE, 1:MSG_ITER, 2:MSG_INFO, 3:MSG_VERS,
  228. 4:MSG_EXIT, 5:MSG_ALL}.get(messages, MSG_ALL)
  229. # build options
  230. opts = {'eps': epsilon,
  231. 'scale': scale,
  232. 'offset': offset,
  233. 'mesg_num': mesg_num,
  234. 'maxCGit': maxCGit,
  235. 'maxfun': maxfun,
  236. 'eta': eta,
  237. 'stepmx': stepmx,
  238. 'accuracy': accuracy,
  239. 'minfev': fmin,
  240. 'ftol': ftol,
  241. 'xtol': xtol,
  242. 'gtol': pgtol,
  243. 'rescale': rescale,
  244. 'disp': False}
  245. res = _minimize_tnc(fun, x0, args, jac, bounds, callback=callback, **opts)
  246. return res['x'], res['nfev'], res['status']
  247. def _minimize_tnc(fun, x0, args=(), jac=None, bounds=None,
  248. eps=1e-8, scale=None, offset=None, mesg_num=None,
  249. maxCGit=-1, maxiter=None, eta=-1, stepmx=0, accuracy=0,
  250. minfev=0, ftol=-1, xtol=-1, gtol=-1, rescale=-1, disp=False,
  251. callback=None, finite_diff_rel_step=None, maxfun=None,
  252. **unknown_options):
  253. """
  254. Minimize a scalar function of one or more variables using a truncated
  255. Newton (TNC) algorithm.
  256. Options
  257. -------
  258. eps : float or ndarray
  259. If `jac is None` the absolute step size used for numerical
  260. approximation of the jacobian via forward differences.
  261. scale : list of floats
  262. Scaling factors to apply to each variable. If None, the
  263. factors are up-low for interval bounded variables and
  264. 1+|x] fo the others. Defaults to None.
  265. offset : float
  266. Value to subtract from each variable. If None, the
  267. offsets are (up+low)/2 for interval bounded variables
  268. and x for the others.
  269. disp : bool
  270. Set to True to print convergence messages.
  271. maxCGit : int
  272. Maximum number of hessian*vector evaluations per main
  273. iteration. If maxCGit == 0, the direction chosen is
  274. -gradient if maxCGit < 0, maxCGit is set to
  275. max(1,min(50,n/2)). Defaults to -1.
  276. maxiter : int, optional
  277. Maximum number of function evaluations. If `maxfun` is also provided
  278. then `maxiter` is ignored.
  279. Default is None.
  280. .. deprecated :: 1.9.0
  281. `maxiter` is deprecated in favor of `maxfun` and will removed in
  282. SciPy 1.11.0.
  283. eta : float
  284. Severity of the line search. If < 0 or > 1, set to 0.25.
  285. Defaults to -1.
  286. stepmx : float
  287. Maximum step for the line search. May be increased during
  288. call. If too small, it will be set to 10.0. Defaults to 0.
  289. accuracy : float
  290. Relative precision for finite difference calculations. If
  291. <= machine_precision, set to sqrt(machine_precision).
  292. Defaults to 0.
  293. minfev : float
  294. Minimum function value estimate. Defaults to 0.
  295. ftol : float
  296. Precision goal for the value of f in the stopping criterion.
  297. If ftol < 0.0, ftol is set to 0.0 defaults to -1.
  298. xtol : float
  299. Precision goal for the value of x in the stopping
  300. criterion (after applying x scaling factors). If xtol <
  301. 0.0, xtol is set to sqrt(machine_precision). Defaults to
  302. -1.
  303. gtol : float
  304. Precision goal for the value of the projected gradient in
  305. the stopping criterion (after applying x scaling factors).
  306. If gtol < 0.0, gtol is set to 1e-2 * sqrt(accuracy).
  307. Setting it to 0.0 is not recommended. Defaults to -1.
  308. rescale : float
  309. Scaling factor (in log10) used to trigger f value
  310. rescaling. If 0, rescale at each iteration. If a large
  311. value, never rescale. If < 0, rescale is set to 1.3.
  312. finite_diff_rel_step : None or array_like, optional
  313. If `jac in ['2-point', '3-point', 'cs']` the relative step size to
  314. use for numerical approximation of the jacobian. The absolute step
  315. size is computed as ``h = rel_step * sign(x) * max(1, abs(x))``,
  316. possibly adjusted to fit into the bounds. For ``method='3-point'``
  317. the sign of `h` is ignored. If None (default) then step is selected
  318. automatically.
  319. maxfun : int
  320. Maximum number of function evaluations. If None, `maxfun` is
  321. set to max(100, 10*len(x0)). Defaults to None.
  322. """
  323. _check_unknown_options(unknown_options)
  324. fmin = minfev
  325. pgtol = gtol
  326. x0 = asfarray(x0).flatten()
  327. n = len(x0)
  328. if bounds is None:
  329. bounds = [(None,None)] * n
  330. if len(bounds) != n:
  331. raise ValueError('length of x0 != length of bounds')
  332. new_bounds = old_bound_to_new(bounds)
  333. if mesg_num is not None:
  334. messages = {0:MSG_NONE, 1:MSG_ITER, 2:MSG_INFO, 3:MSG_VERS,
  335. 4:MSG_EXIT, 5:MSG_ALL}.get(mesg_num, MSG_ALL)
  336. elif disp:
  337. messages = MSG_ALL
  338. else:
  339. messages = MSG_NONE
  340. sf = _prepare_scalar_function(fun, x0, jac=jac, args=args, epsilon=eps,
  341. finite_diff_rel_step=finite_diff_rel_step,
  342. bounds=new_bounds)
  343. func_and_grad = sf.fun_and_grad
  344. """
  345. low, up : the bounds (lists of floats)
  346. if low is None, the lower bounds are removed.
  347. if up is None, the upper bounds are removed.
  348. low and up defaults to None
  349. """
  350. low = zeros(n)
  351. up = zeros(n)
  352. for i in range(n):
  353. if bounds[i] is None:
  354. l, u = -inf, inf
  355. else:
  356. l,u = bounds[i]
  357. if l is None:
  358. low[i] = -inf
  359. else:
  360. low[i] = l
  361. if u is None:
  362. up[i] = inf
  363. else:
  364. up[i] = u
  365. if scale is None:
  366. scale = array([])
  367. if offset is None:
  368. offset = array([])
  369. if maxfun is None:
  370. if maxiter is not None:
  371. warnings.warn(
  372. "'maxiter' has been deprecated in favor of 'maxfun'"
  373. " and will be removed in SciPy 1.11.0.",
  374. DeprecationWarning, stacklevel=3
  375. )
  376. maxfun = maxiter
  377. else:
  378. maxfun = max(100, 10*len(x0))
  379. rc, nf, nit, x, funv, jacv = moduleTNC.tnc_minimize(
  380. func_and_grad, x0, low, up, scale,
  381. offset, messages, maxCGit, maxfun,
  382. eta, stepmx, accuracy, fmin, ftol,
  383. xtol, pgtol, rescale, callback
  384. )
  385. # the TNC documentation states: "On output, x, f and g may be very
  386. # slightly out of sync because of scaling". Therefore re-evaluate
  387. # func_and_grad so they are synced.
  388. funv, jacv = func_and_grad(x)
  389. return OptimizeResult(x=x, fun=funv, jac=jacv, nfev=sf.nfev,
  390. nit=nit, status=rc, message=RCSTRINGS[rc],
  391. success=(-1 < rc < 3))