tnc.py 16 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. from __future__ import division, print_function, absolute_import
  29. from scipy.optimize import moduleTNC, approx_fprime
  30. from .optimize import MemoizeJac, OptimizeResult, _check_unknown_options
  31. from numpy import inf, array, zeros, asfarray
  32. __all__ = ['fmin_tnc']
  33. MSG_NONE = 0 # No messages
  34. MSG_ITER = 1 # One line per iteration
  35. MSG_INFO = 2 # Informational messages
  36. MSG_VERS = 4 # Version info
  37. MSG_EXIT = 8 # Exit reasons
  38. MSG_ALL = MSG_ITER + MSG_INFO + MSG_VERS + MSG_EXIT
  39. MSGS = {
  40. MSG_NONE: "No messages",
  41. MSG_ITER: "One line per iteration",
  42. MSG_INFO: "Informational messages",
  43. MSG_VERS: "Version info",
  44. MSG_EXIT: "Exit reasons",
  45. MSG_ALL: "All messages"
  46. }
  47. INFEASIBLE = -1 # Infeasible (lower bound > upper bound)
  48. LOCALMINIMUM = 0 # Local minimum reached (|pg| ~= 0)
  49. FCONVERGED = 1 # Converged (|f_n-f_(n-1)| ~= 0)
  50. XCONVERGED = 2 # Converged (|x_n-x_(n-1)| ~= 0)
  51. MAXFUN = 3 # Max. number of function evaluations reached
  52. LSFAIL = 4 # Linear search failed
  53. CONSTANT = 5 # All lower bounds are equal to the upper bounds
  54. NOPROGRESS = 6 # Unable to progress
  55. USERABORT = 7 # User requested end of minimization
  56. RCSTRINGS = {
  57. INFEASIBLE: "Infeasible (lower bound > upper bound)",
  58. LOCALMINIMUM: "Local minimum reached (|pg| ~= 0)",
  59. FCONVERGED: "Converged (|f_n-f_(n-1)| ~= 0)",
  60. XCONVERGED: "Converged (|x_n-x_(n-1)| ~= 0)",
  61. MAXFUN: "Max. number of function evaluations reached",
  62. LSFAIL: "Linear search failed",
  63. CONSTANT: "All lower bounds are equal to the upper bounds",
  64. NOPROGRESS: "Unable to progress",
  65. USERABORT: "User requested end of minimization"
  66. }
  67. # Changes to interface made by Travis Oliphant, Apr. 2004 for inclusion in
  68. # SciPy
  69. def fmin_tnc(func, x0, fprime=None, args=(), approx_grad=0,
  70. bounds=None, epsilon=1e-8, scale=None, offset=None,
  71. messages=MSG_ALL, maxCGit=-1, maxfun=None, eta=-1,
  72. stepmx=0, accuracy=0, fmin=0, ftol=-1, xtol=-1, pgtol=-1,
  73. rescale=-1, disp=None, callback=None):
  74. """
  75. Minimize a function with variables subject to bounds, using
  76. gradient information in a truncated Newton algorithm. This
  77. method wraps a C implementation of the algorithm.
  78. Parameters
  79. ----------
  80. func : callable ``func(x, *args)``
  81. Function to minimize. Must do one of:
  82. 1. Return f and g, where f is the value of the function and g its
  83. gradient (a list of floats).
  84. 2. Return the function value but supply gradient function
  85. separately as `fprime`.
  86. 3. Return the function value and set ``approx_grad=True``.
  87. If the function returns None, the minimization
  88. is aborted.
  89. x0 : array_like
  90. Initial estimate of minimum.
  91. fprime : callable ``fprime(x, *args)``, optional
  92. Gradient of `func`. If None, then either `func` must return the
  93. function value and the gradient (``f,g = func(x, *args)``)
  94. or `approx_grad` must be True.
  95. args : tuple, optional
  96. Arguments to pass to function.
  97. approx_grad : bool, optional
  98. If true, approximate the gradient numerically.
  99. bounds : list, optional
  100. (min, max) pairs for each element in x0, defining the
  101. bounds on that parameter. Use None or +/-inf for one of
  102. min or max when there is no bound in that direction.
  103. epsilon : float, optional
  104. Used if approx_grad is True. The stepsize in a finite
  105. difference approximation for fprime.
  106. scale : array_like, optional
  107. Scaling factors to apply to each variable. If None, the
  108. factors are up-low for interval bounded variables and
  109. 1+|x| for the others. Defaults to None.
  110. offset : array_like, optional
  111. Value to subtract from each variable. If None, the
  112. offsets are (up+low)/2 for interval bounded variables
  113. and x for the others.
  114. messages : int, optional
  115. Bit mask used to select messages display during
  116. minimization values defined in the MSGS dict. Defaults to
  117. MGS_ALL.
  118. disp : int, optional
  119. Integer interface to messages. 0 = no message, 5 = all messages
  120. maxCGit : int, optional
  121. Maximum number of hessian*vector evaluations per main
  122. iteration. If maxCGit == 0, the direction chosen is
  123. -gradient if maxCGit < 0, maxCGit is set to
  124. max(1,min(50,n/2)). Defaults to -1.
  125. maxfun : int, optional
  126. Maximum number of function evaluation. if None, maxfun is
  127. set to max(100, 10*len(x0)). Defaults to None.
  128. eta : float, optional
  129. Severity of the line search. if < 0 or > 1, set to 0.25.
  130. Defaults to -1.
  131. stepmx : float, optional
  132. Maximum step for the line search. May be increased during
  133. call. If too small, it will be set to 10.0. Defaults to 0.
  134. accuracy : float, optional
  135. Relative precision for finite difference calculations. If
  136. <= machine_precision, set to sqrt(machine_precision).
  137. Defaults to 0.
  138. fmin : float, optional
  139. Minimum function value estimate. Defaults to 0.
  140. ftol : float, optional
  141. Precision goal for the value of f in the stopping criterion.
  142. If ftol < 0.0, ftol is set to 0.0 defaults to -1.
  143. xtol : float, optional
  144. Precision goal for the value of x in the stopping
  145. criterion (after applying x scaling factors). If xtol <
  146. 0.0, xtol is set to sqrt(machine_precision). Defaults to
  147. -1.
  148. pgtol : float, optional
  149. Precision goal for the value of the projected gradient in
  150. the stopping criterion (after applying x scaling factors).
  151. If pgtol < 0.0, pgtol is set to 1e-2 * sqrt(accuracy).
  152. Setting it to 0.0 is not recommended. Defaults to -1.
  153. rescale : float, optional
  154. Scaling factor (in log10) used to trigger f value
  155. rescaling. If 0, rescale at each iteration. If a large
  156. value, never rescale. If < 0, rescale is set to 1.3.
  157. callback : callable, optional
  158. Called after each iteration, as callback(xk), where xk is the
  159. current parameter vector.
  160. Returns
  161. -------
  162. x : ndarray
  163. The solution.
  164. nfeval : int
  165. The number of function evaluations.
  166. rc : int
  167. Return code, see below
  168. See also
  169. --------
  170. minimize: Interface to minimization algorithms for multivariate
  171. functions. See the 'TNC' `method` in particular.
  172. Notes
  173. -----
  174. The underlying algorithm is truncated Newton, also called
  175. Newton Conjugate-Gradient. This method differs from
  176. scipy.optimize.fmin_ncg in that
  177. 1. It wraps a C implementation of the algorithm
  178. 2. It allows each variable to be given an upper and lower bound.
  179. The algorithm incorporates the bound constraints by determining
  180. the descent direction as in an unconstrained truncated Newton,
  181. but never taking a step-size large enough to leave the space
  182. of feasible x's. The algorithm keeps track of a set of
  183. currently active constraints, and ignores them when computing
  184. the minimum allowable step size. (The x's associated with the
  185. active constraint are kept fixed.) If the maximum allowable
  186. step size is zero then a new constraint is added. At the end
  187. of each iteration one of the constraints may be deemed no
  188. longer active and removed. A constraint is considered
  189. no longer active is if it is currently active
  190. but the gradient for that variable points inward from the
  191. constraint. The specific constraint removed is the one
  192. associated with the variable of largest index whose
  193. constraint is no longer active.
  194. Return codes are defined as follows::
  195. -1 : Infeasible (lower bound > upper bound)
  196. 0 : Local minimum reached (|pg| ~= 0)
  197. 1 : Converged (|f_n-f_(n-1)| ~= 0)
  198. 2 : Converged (|x_n-x_(n-1)| ~= 0)
  199. 3 : Max. number of function evaluations reached
  200. 4 : Linear search failed
  201. 5 : All lower bounds are equal to the upper bounds
  202. 6 : Unable to progress
  203. 7 : User requested end of minimization
  204. References
  205. ----------
  206. Wright S., Nocedal J. (2006), 'Numerical Optimization'
  207. Nash S.G. (1984), "Newton-Type Minimization Via the Lanczos Method",
  208. SIAM Journal of Numerical Analysis 21, pp. 770-778
  209. """
  210. # handle fprime/approx_grad
  211. if approx_grad:
  212. fun = func
  213. jac = None
  214. elif fprime is None:
  215. fun = MemoizeJac(func)
  216. jac = fun.derivative
  217. else:
  218. fun = func
  219. jac = fprime
  220. if disp is not None: # disp takes precedence over messages
  221. mesg_num = disp
  222. else:
  223. mesg_num = {0:MSG_NONE, 1:MSG_ITER, 2:MSG_INFO, 3:MSG_VERS,
  224. 4:MSG_EXIT, 5:MSG_ALL}.get(messages, MSG_ALL)
  225. # build options
  226. opts = {'eps': epsilon,
  227. 'scale': scale,
  228. 'offset': offset,
  229. 'mesg_num': mesg_num,
  230. 'maxCGit': maxCGit,
  231. 'maxiter': maxfun,
  232. 'eta': eta,
  233. 'stepmx': stepmx,
  234. 'accuracy': accuracy,
  235. 'minfev': fmin,
  236. 'ftol': ftol,
  237. 'xtol': xtol,
  238. 'gtol': pgtol,
  239. 'rescale': rescale,
  240. 'disp': False}
  241. res = _minimize_tnc(fun, x0, args, jac, bounds, callback=callback, **opts)
  242. return res['x'], res['nfev'], res['status']
  243. def _minimize_tnc(fun, x0, args=(), jac=None, bounds=None,
  244. eps=1e-8, scale=None, offset=None, mesg_num=None,
  245. maxCGit=-1, maxiter=None, eta=-1, stepmx=0, accuracy=0,
  246. minfev=0, ftol=-1, xtol=-1, gtol=-1, rescale=-1, disp=False,
  247. callback=None, **unknown_options):
  248. """
  249. Minimize a scalar function of one or more variables using a truncated
  250. Newton (TNC) algorithm.
  251. Options
  252. -------
  253. eps : float
  254. Step size used for numerical approximation of the jacobian.
  255. scale : list of floats
  256. Scaling factors to apply to each variable. If None, the
  257. factors are up-low for interval bounded variables and
  258. 1+|x] fo the others. Defaults to None
  259. offset : float
  260. Value to subtract from each variable. If None, the
  261. offsets are (up+low)/2 for interval bounded variables
  262. and x for the others.
  263. disp : bool
  264. Set to True to print convergence messages.
  265. maxCGit : int
  266. Maximum number of hessian*vector evaluations per main
  267. iteration. If maxCGit == 0, the direction chosen is
  268. -gradient if maxCGit < 0, maxCGit is set to
  269. max(1,min(50,n/2)). Defaults to -1.
  270. maxiter : int
  271. Maximum number of function evaluation. if None, `maxiter` is
  272. set to max(100, 10*len(x0)). Defaults to None.
  273. eta : float
  274. Severity of the line search. if < 0 or > 1, set to 0.25.
  275. Defaults to -1.
  276. stepmx : float
  277. Maximum step for the line search. May be increased during
  278. call. If too small, it will be set to 10.0. Defaults to 0.
  279. accuracy : float
  280. Relative precision for finite difference calculations. If
  281. <= machine_precision, set to sqrt(machine_precision).
  282. Defaults to 0.
  283. minfev : float
  284. Minimum function value estimate. Defaults to 0.
  285. ftol : float
  286. Precision goal for the value of f in the stopping criterion.
  287. If ftol < 0.0, ftol is set to 0.0 defaults to -1.
  288. xtol : float
  289. Precision goal for the value of x in the stopping
  290. criterion (after applying x scaling factors). If xtol <
  291. 0.0, xtol is set to sqrt(machine_precision). Defaults to
  292. -1.
  293. gtol : float
  294. Precision goal for the value of the projected gradient in
  295. the stopping criterion (after applying x scaling factors).
  296. If gtol < 0.0, gtol is set to 1e-2 * sqrt(accuracy).
  297. Setting it to 0.0 is not recommended. Defaults to -1.
  298. rescale : float
  299. Scaling factor (in log10) used to trigger f value
  300. rescaling. If 0, rescale at each iteration. If a large
  301. value, never rescale. If < 0, rescale is set to 1.3.
  302. """
  303. _check_unknown_options(unknown_options)
  304. epsilon = eps
  305. maxfun = maxiter
  306. fmin = minfev
  307. pgtol = gtol
  308. x0 = asfarray(x0).flatten()
  309. n = len(x0)
  310. if bounds is None:
  311. bounds = [(None,None)] * n
  312. if len(bounds) != n:
  313. raise ValueError('length of x0 != length of bounds')
  314. if mesg_num is not None:
  315. messages = {0:MSG_NONE, 1:MSG_ITER, 2:MSG_INFO, 3:MSG_VERS,
  316. 4:MSG_EXIT, 5:MSG_ALL}.get(mesg_num, MSG_ALL)
  317. elif disp:
  318. messages = MSG_ALL
  319. else:
  320. messages = MSG_NONE
  321. if jac is None:
  322. def func_and_grad(x):
  323. f = fun(x, *args)
  324. g = approx_fprime(x, fun, epsilon, *args)
  325. return f, g
  326. else:
  327. def func_and_grad(x):
  328. f = fun(x, *args)
  329. g = jac(x, *args)
  330. return f, g
  331. """
  332. low, up : the bounds (lists of floats)
  333. if low is None, the lower bounds are removed.
  334. if up is None, the upper bounds are removed.
  335. low and up defaults to None
  336. """
  337. low = zeros(n)
  338. up = zeros(n)
  339. for i in range(n):
  340. if bounds[i] is None:
  341. l, u = -inf, inf
  342. else:
  343. l,u = bounds[i]
  344. if l is None:
  345. low[i] = -inf
  346. else:
  347. low[i] = l
  348. if u is None:
  349. up[i] = inf
  350. else:
  351. up[i] = u
  352. if scale is None:
  353. scale = array([])
  354. if offset is None:
  355. offset = array([])
  356. if maxfun is None:
  357. maxfun = max(100, 10*len(x0))
  358. rc, nf, nit, x = moduleTNC.minimize(func_and_grad, x0, low, up, scale,
  359. offset, messages, maxCGit, maxfun,
  360. eta, stepmx, accuracy, fmin, ftol,
  361. xtol, pgtol, rescale, callback)
  362. funv, jacv = func_and_grad(x)
  363. return OptimizeResult(x=x, fun=funv, jac=jacv, nfev=nf, nit=nit, status=rc,
  364. message=RCSTRINGS[rc], success=(-1 < rc < 3))
  365. if __name__ == '__main__':
  366. # Examples for TNC
  367. def example():
  368. print("Example")
  369. # A function to minimize
  370. def function(x):
  371. f = pow(x[0],2.0)+pow(abs(x[1]),3.0)
  372. g = [0,0]
  373. g[0] = 2.0*x[0]
  374. g[1] = 3.0*pow(abs(x[1]),2.0)
  375. if x[1] < 0:
  376. g[1] = -g[1]
  377. return f, g
  378. # Optimizer call
  379. x, nf, rc = fmin_tnc(function, [-7, 3], bounds=([-10, 1], [10, 10]))
  380. print("After", nf, "function evaluations, TNC returned:", RCSTRINGS[rc])
  381. print("x =", x)
  382. print("exact value = [0, 1]")
  383. print()
  384. example()