Coverage for /usr/lib/python3/dist-packages/scipy/optimize/_cobyla_py.py: 19%
98 statements
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-14 15:55 +0200
« prev ^ index » next coverage.py v7.9.1, created at 2025-06-14 15:55 +0200
1"""
2Interface to Constrained Optimization By Linear Approximation
4Functions
5---------
6.. autosummary::
7 :toctree: generated/
9 fmin_cobyla
11"""
13import functools
14from threading import RLock
16import numpy as np
17from scipy.optimize import _cobyla as cobyla
18from ._optimize import (OptimizeResult, _check_unknown_options,
19 _prepare_scalar_function)
20try:
21 from itertools import izip
22except ImportError:
23 izip = zip
25__all__ = ['fmin_cobyla']
27# Workarund as _cobyla.minimize is not threadsafe
28# due to an unknown f2py bug and can segfault,
29# see gh-9658.
30_module_lock = RLock()
31def synchronized(func):
32 @functools.wraps(func)
33 def wrapper(*args, **kwargs):
34 with _module_lock:
35 return func(*args, **kwargs)
36 return wrapper
38@synchronized
39def fmin_cobyla(func, x0, cons, args=(), consargs=None, rhobeg=1.0,
40 rhoend=1e-4, maxfun=1000, disp=None, catol=2e-4,
41 *, callback=None):
42 """
43 Minimize a function using the Constrained Optimization By Linear
44 Approximation (COBYLA) method. This method wraps a FORTRAN
45 implementation of the algorithm.
47 Parameters
48 ----------
49 func : callable
50 Function to minimize. In the form func(x, \\*args).
51 x0 : ndarray
52 Initial guess.
53 cons : sequence
54 Constraint functions; must all be ``>=0`` (a single function
55 if only 1 constraint). Each function takes the parameters `x`
56 as its first argument, and it can return either a single number or
57 an array or list of numbers.
58 args : tuple, optional
59 Extra arguments to pass to function.
60 consargs : tuple, optional
61 Extra arguments to pass to constraint functions (default of None means
62 use same extra arguments as those passed to func).
63 Use ``()`` for no extra arguments.
64 rhobeg : float, optional
65 Reasonable initial changes to the variables.
66 rhoend : float, optional
67 Final accuracy in the optimization (not precisely guaranteed). This
68 is a lower bound on the size of the trust region.
69 disp : {0, 1, 2, 3}, optional
70 Controls the frequency of output; 0 implies no output.
71 maxfun : int, optional
72 Maximum number of function evaluations.
73 catol : float, optional
74 Absolute tolerance for constraint violations.
75 callback : callable, optional
76 Called after each iteration, as ``callback(x)``, where ``x`` is the
77 current parameter vector.
79 Returns
80 -------
81 x : ndarray
82 The argument that minimises `f`.
84 See also
85 --------
86 minimize: Interface to minimization algorithms for multivariate
87 functions. See the 'COBYLA' `method` in particular.
89 Notes
90 -----
91 This algorithm is based on linear approximations to the objective
92 function and each constraint. We briefly describe the algorithm.
94 Suppose the function is being minimized over k variables. At the
95 jth iteration the algorithm has k+1 points v_1, ..., v_(k+1),
96 an approximate solution x_j, and a radius RHO_j.
97 (i.e., linear plus a constant) approximations to the objective
98 function and constraint functions such that their function values
99 agree with the linear approximation on the k+1 points v_1,.., v_(k+1).
100 This gives a linear program to solve (where the linear approximations
101 of the constraint functions are constrained to be non-negative).
103 However, the linear approximations are likely only good
104 approximations near the current simplex, so the linear program is
105 given the further requirement that the solution, which
106 will become x_(j+1), must be within RHO_j from x_j. RHO_j only
107 decreases, never increases. The initial RHO_j is rhobeg and the
108 final RHO_j is rhoend. In this way COBYLA's iterations behave
109 like a trust region algorithm.
111 Additionally, the linear program may be inconsistent, or the
112 approximation may give poor improvement. For details about
113 how these issues are resolved, as well as how the points v_i are
114 updated, refer to the source code or the references below.
117 References
118 ----------
119 Powell M.J.D. (1994), "A direct search optimization method that models
120 the objective and constraint functions by linear interpolation.", in
121 Advances in Optimization and Numerical Analysis, eds. S. Gomez and
122 J-P Hennart, Kluwer Academic (Dordrecht), pp. 51-67
124 Powell M.J.D. (1998), "Direct search algorithms for optimization
125 calculations", Acta Numerica 7, 287-336
127 Powell M.J.D. (2007), "A view of algorithms for optimization without
128 derivatives", Cambridge University Technical Report DAMTP 2007/NA03
131 Examples
132 --------
133 Minimize the objective function f(x,y) = x*y subject
134 to the constraints x**2 + y**2 < 1 and y > 0::
136 >>> def objective(x):
137 ... return x[0]*x[1]
138 ...
139 >>> def constr1(x):
140 ... return 1 - (x[0]**2 + x[1]**2)
141 ...
142 >>> def constr2(x):
143 ... return x[1]
144 ...
145 >>> from scipy.optimize import fmin_cobyla
146 >>> fmin_cobyla(objective, [0.0, 0.1], [constr1, constr2], rhoend=1e-7)
147 array([-0.70710685, 0.70710671])
149 The exact solution is (-sqrt(2)/2, sqrt(2)/2).
153 """
154 err = "cons must be a sequence of callable functions or a single"\
155 " callable function."
156 try:
157 len(cons)
158 except TypeError as e:
159 if callable(cons):
160 cons = [cons]
161 else:
162 raise TypeError(err) from e
163 else:
164 for thisfunc in cons:
165 if not callable(thisfunc):
166 raise TypeError(err)
168 if consargs is None:
169 consargs = args
171 # build constraints
172 con = tuple({'type': 'ineq', 'fun': c, 'args': consargs} for c in cons)
174 # options
175 opts = {'rhobeg': rhobeg,
176 'tol': rhoend,
177 'disp': disp,
178 'maxiter': maxfun,
179 'catol': catol,
180 'callback': callback}
182 sol = _minimize_cobyla(func, x0, args, constraints=con,
183 **opts)
184 if disp and not sol['success']:
185 print(f"COBYLA failed to find a solution: {sol.message}")
186 return sol['x']
189@synchronized
190def _minimize_cobyla(fun, x0, args=(), constraints=(),
191 rhobeg=1.0, tol=1e-4, maxiter=1000,
192 disp=False, catol=2e-4, callback=None, bounds=None,
193 **unknown_options):
194 """
195 Minimize a scalar function of one or more variables using the
196 Constrained Optimization BY Linear Approximation (COBYLA) algorithm.
198 Options
199 -------
200 rhobeg : float
201 Reasonable initial changes to the variables.
202 tol : float
203 Final accuracy in the optimization (not precisely guaranteed).
204 This is a lower bound on the size of the trust region.
205 disp : bool
206 Set to True to print convergence messages. If False,
207 `verbosity` is ignored as set to 0.
208 maxiter : int
209 Maximum number of function evaluations.
210 catol : float
211 Tolerance (absolute) for constraint violations
213 """
214 _check_unknown_options(unknown_options)
215 maxfun = maxiter
216 rhoend = tol
217 iprint = int(bool(disp))
219 # check constraints
220 if isinstance(constraints, dict):
221 constraints = (constraints, )
223 if bounds:
224 i_lb = np.isfinite(bounds.lb)
225 if np.any(i_lb):
226 def lb_constraint(x, *args, **kwargs):
227 return x[i_lb] - bounds.lb[i_lb]
229 constraints.append({'type': 'ineq', 'fun': lb_constraint})
231 i_ub = np.isfinite(bounds.ub)
232 if np.any(i_ub):
233 def ub_constraint(x):
234 return bounds.ub[i_ub] - x[i_ub]
236 constraints.append({'type': 'ineq', 'fun': ub_constraint})
238 for ic, con in enumerate(constraints):
239 # check type
240 try:
241 ctype = con['type'].lower()
242 except KeyError as e:
243 raise KeyError('Constraint %d has no type defined.' % ic) from e
244 except TypeError as e:
245 raise TypeError('Constraints must be defined using a '
246 'dictionary.') from e
247 except AttributeError as e:
248 raise TypeError("Constraint's type must be a string.") from e
249 else:
250 if ctype != 'ineq':
251 raise ValueError("Constraints of type '%s' not handled by "
252 "COBYLA." % con['type'])
254 # check function
255 if 'fun' not in con:
256 raise KeyError('Constraint %d has no function defined.' % ic)
258 # check extra arguments
259 if 'args' not in con:
260 con['args'] = ()
262 # m is the total number of constraint values
263 # it takes into account that some constraints may be vector-valued
264 cons_lengths = []
265 for c in constraints:
266 f = c['fun'](x0, *c['args'])
267 try:
268 cons_length = len(f)
269 except TypeError:
270 cons_length = 1
271 cons_lengths.append(cons_length)
272 m = sum(cons_lengths)
274 # create the ScalarFunction, cobyla doesn't require derivative function
275 def _jac(x, *args):
276 return None
278 sf = _prepare_scalar_function(fun, x0, args=args, jac=_jac)
280 def calcfc(x, con):
281 f = sf.fun(x)
282 i = 0
283 for size, c in izip(cons_lengths, constraints):
284 con[i: i + size] = c['fun'](x, *c['args'])
285 i += size
286 return f
288 def wrapped_callback(x):
289 if callback is not None:
290 callback(np.copy(x))
292 info = np.zeros(4, np.float64)
293 xopt, info = cobyla.minimize(calcfc, m=m, x=np.copy(x0), rhobeg=rhobeg,
294 rhoend=rhoend, iprint=iprint, maxfun=maxfun,
295 dinfo=info, callback=wrapped_callback)
297 if info[3] > catol:
298 # Check constraint violation
299 info[0] = 4
301 return OptimizeResult(x=xopt,
302 status=int(info[0]),
303 success=info[0] == 1,
304 message={1: 'Optimization terminated successfully.',
305 2: 'Maximum number of function evaluations '
306 'has been exceeded.',
307 3: 'Rounding errors are becoming damaging '
308 'in COBYLA subroutine.',
309 4: 'Did not converge to a solution '
310 'satisfying the constraints. See '
311 '`maxcv` for magnitude of violation.',
312 5: 'NaN result encountered.'
313 }.get(info[0], 'Unknown exit status.'),
314 nfev=int(info[1]),
315 fun=info[2],
316 maxcv=info[3])