Python >> Tutoriel Python >  >> Python Tag >> SciPy

Minimiser la fonction quadratique soumise à des contraintes d'égalité linéaire avec SciPy

Voici comment ce problème pourrait être résolu en utilisant nlopt qui est une bibliothèque pour l'optimisation non linéaire qui m'a beaucoup impressionné.

Tout d'abord, la fonction objectif et le gradient sont tous deux définis à l'aide de la même fonction :

def obj_func(x, grad):
    if grad.size > 0:
        grad[:] = obj_jac(x)
    return ( ( ( x/x0 - 1 )) ** 2 ).sum()

def obj_jac(x):
    return 2. * ( x - x0 ) / x0 ** 2

def constr_func(x, grad):
    if grad.size > 0:
        grad[:] = constr_jac(x)
    return x.sum() - target

def constr_jac(x):
    return np.ones(n)

Ensuite, pour exécuter la minimisation à l'aide de Nelder-Mead et SLSQP :

opt = nlopt.opt(nlopt.LN_NELDERMEAD,len(x0)-1)
opt.set_min_objective(unconstr_func)
opt.set_ftol_abs(1e-15)
xopt = opt.optimize(x0[1:].copy())
xopt = np.hstack([target - xopt.sum(), xopt])
fval = opt.last_optimum_value()
print_res(xopt,fval,"Nelder-Mead");

opt = nlopt.opt(nlopt.LD_SLSQP,len(x0))
opt.set_min_objective(obj_func)
opt.add_equality_constraint(constr_func)
opt.set_ftol_abs(1e-15)
xopt = opt.optimize(x0.copy())
fval = opt.last_optimum_value()
print_res(xopt,fval,"SLSQP w/ jacobian");

Et voici les résultats :

 *****  Nelder-Mead  ***** 

obj func value at solution 0.00454545454546
result:  3
starting values:  [ 10000.  20000.  30000.  40000.  50000.]
ending values:    [10090 20363 30818 41454 52272]
% diff [0 1 2 3 4]
target achieved? 155000.0 155000.0


 *****  SLSQP w/ jacobian  ***** 

obj func value at solution 0.00454545454545
result:  3
starting values:  [ 10000.  20000.  30000.  40000.  50000.]
ending values:    [10090 20363 30818 41454 52272]
% diff [0 1 2 3 4]
target achieved? 155000.0 155000.0

En testant cela, je pense avoir découvert quel était le problème avec la tentative initiale. Si je règle la tolérance absolue sur la fonction à 1e-8 c'est ce que j'obtiens par défaut des fonctions scipy :

 *****  Nelder-Mead  ***** 

obj func value at solution 0.0045454580693
result:  3
starting values:  [ 10000.  20000.  30000.  40000.  50000.]
ending values:    [10090 20363 30816 41454 52274]
% diff [0 1 2 3 4]
target achieved? 155000.0 155000.0


 *****  SLSQP w/ jacobian  ***** 

obj func value at solution 0.0146361108503
result:  3
starting values:  [ 10000.  20000.  30000.  40000.  50000.]
ending values:    [10999 21000 31000 41000 51000]
% diff [9 5 3 2 2]
target achieved? 155000.0 155000.0

c'est exactement ce que vous voyiez. Donc, je suppose que le minimiseur se retrouve quelque part dans l'espace de vraisemblance pendant SLSQP où le saut suivant est inférieur à 1e-8 depuis la dernière place.