Solveur de programme quadratique (QP) qui ne dépend que de NumPy/SciPy?

Je voudrais que les étudiants résolvent un programme quadratique dans une tâche sans qu'ils aient à installer un logiciel supplémentaire comme cvxopt, etc. Existe-t-il une implémentation Python disponible qui ne dépend que de NumPy/SciPy?

27
demandé sur ali_m 2013-06-09 16:50:55

4 réponses

J'ai rencontré une bonne solution et je voulais la sortir. Il existe une implémentation python de LOQO dans la boîte à outils d'apprentissage automatique ELEFANT de NICTA (http://elefant.forge.nicta.com.au à partir de cette publication). Jetez un oeil à l'optimisation.intpointsolver. Cela a été codé par Alex Smola, et j'ai utilisé une version C du même code avec beaucoup de succès.

5
répondu Tom Vacek 2013-11-04 14:37:02

Je ne suis pas très familier avec la programmation quadratique, mais je pense que vous pouvez résoudre ce genre de problème en utilisant simplement les algorithmes de minimisation contraints de scipy.optimize. Voici un exemple:

import numpy as np
from scipy import optimize
from matplotlib import pyplot as plt
from mpl_toolkits.mplot3d.axes3d import Axes3D

# minimize
#     F = x[1]^2 + 4x[2]^2 -32x[2] + 64

# subject to:
#      x[1] + x[2] <= 7
#     -x[1] + 2x[2] <= 4
#      x[1] >= 0
#      x[2] >= 0
#      x[2] <= 4

# in matrix notation:
#     F = (1/2)*x.T*H*x + c*x + c0

# subject to:
#     Ax <= b

# where:
#     H = [[2, 0],
#          [0, 8]]

#     c = [0, -32]

#     c0 = 64

#     A = [[ 1, 1],
#          [-1, 2],
#          [-1, 0],
#          [0, -1],
#          [0,  1]]

#     b = [7,4,0,0,4]

H = np.array([[2., 0.],
              [0., 8.]])

c = np.array([0, -32])

c0 = 64

A = np.array([[ 1., 1.],
              [-1., 2.],
              [-1., 0.],
              [0., -1.],
              [0.,  1.]])

b = np.array([7., 4., 0., 0., 4.])

x0 = np.random.randn(2)

def loss(x, sign=1.):
    return sign * (0.5 * np.dot(x.T, np.dot(H, x))+ np.dot(c, x) + c0)

def jac(x, sign=1.):
    return sign * (np.dot(x.T, H) + c)

cons = {'type':'ineq',
        'fun':lambda x: b - np.dot(A,x),
        'jac':lambda x: -A}

opt = {'disp':False}

def solve():

    res_cons = optimize.minimize(loss, x0, jac=jac,constraints=cons,
                                 method='SLSQP', options=opt)

    res_uncons = optimize.minimize(loss, x0, jac=jac, method='SLSQP',
                                   options=opt)

    print '\nConstrained:'
    print res_cons

    print '\nUnconstrained:'
    print res_uncons

    x1, x2 = res_cons['x']
    f = res_cons['fun']

    x1_unc, x2_unc = res_uncons['x']
    f_unc = res_uncons['fun']

    # plotting
    xgrid = np.mgrid[-2:4:0.1, 1.5:5.5:0.1]
    xvec = xgrid.reshape(2, -1).T
    F = np.vstack([loss(xi) for xi in xvec]).reshape(xgrid.shape[1:])

    ax = plt.axes(projection='3d')
    ax.hold(True)
    ax.plot_surface(xgrid[0], xgrid[1], F, rstride=1, cstride=1,
                    cmap=plt.cm.jet, shade=True, alpha=0.9, linewidth=0)
    ax.plot3D([x1], [x2], [f], 'og', mec='w', label='Constrained minimum')
    ax.plot3D([x1_unc], [x2_unc], [f_unc], 'oy', mec='w',
              label='Unconstrained minimum')
    ax.legend(fancybox=True, numpoints=1)
    ax.set_xlabel('x1')
    ax.set_ylabel('x2')
    ax.set_zlabel('F')

Sortie:

Constrained:
  status: 0
 success: True
    njev: 4
    nfev: 4
     fun: 7.9999999999997584
       x: array([ 2.,  3.])
 message: 'Optimization terminated successfully.'
     jac: array([ 4., -8.,  0.])
     nit: 4

Unconstrained:
  status: 0
 success: True
    njev: 3
    nfev: 5
     fun: 0.0
       x: array([ -2.66453526e-15,   4.00000000e+00])
 message: 'Optimization terminated successfully.'
     jac: array([ -5.32907052e-15,  -3.55271368e-15,   0.00000000e+00])
     nit: 3

entrez la description de l'image ici

33
répondu ali_m 2015-02-26 09:51:25

Cela peut être une réponse tardive, mais j'ai trouvé CVXOPT - http://cvxopt.org/ - comme couramment utilisés libre bibliothèque python pour Quadratic Programming. Cependant, il n'est pas facile à installer, car il nécessite l'installation d'autres dépendances.

9
répondu Curious 2016-10-06 23:06:37

mystic fournit une implémentation Python pure d'algorithmes d'optimisation non linéaires/non convexes avec des fonctionnalités de contraintes avancées qui ne se trouvent généralement que dans les solveurs QP. mystic fournit en fait des contraintes plus robustes que la plupart des solveurs QP. Cependant, si vous recherchez une vitesse algorithmique d'optimisation, ce qui suit n'est pas pour vous. mystic n'est pas lent, mais il s'agit de python pur par opposition aux liaisons python à C. Si vous recherchez des fonctionnalités de flexibilité et de contraintes QP dans un solveur non-linéaire, alors vous pourriez être intéressé.

"""
Maximize: f = 2*x[0]*x[1] + 2*x[0] - x[0]**2 - 2*x[1]**2

Subject to: -2*x[0] + 2*x[1] <= -2
             2*x[0] - 4*x[1] <= 0
               x[0]**3 -x[1] == 0

where: 0 <= x[0] <= inf
       1 <= x[1] <= inf
"""
import numpy as np
import mystic.symbolic as ms
import mystic.solvers as my
import mystic.math as mm

# generate constraints and penalty for a nonlinear system of equations 
ieqn = '''
   -2*x0 + 2*x1 <= -2
    2*x0 - 4*x1 <= 0'''
eqn = '''
     x0**3 - x1 == 0'''
cons = ms.generate_constraint(ms.generate_solvers(ms.simplify(eqn,target='x1')))
pens = ms.generate_penalty(ms.generate_conditions(ieqn), k=1e3)
bounds = [(0., None), (1., None)]

# get the objective
def objective(x, sign=1):
  x = np.asarray(x)
  return sign * (2*x[0]*x[1] + 2*x[0] - x[0]**2 - 2*x[1]**2)

# solve    
x0 = np.random.rand(2)
sol = my.fmin_powell(objective, x0, constraint=cons, penalty=pens, disp=True,
                     bounds=bounds, gtol=3, ftol=1e-6, full_output=True,
                     args=(-1,))

print 'x* = %s; f(x*) = %s' % (sol[0], -sol[1])

Il convient de noter que {[2] } peut appliquer de manière générique des contraintes D'égalité et d'inégalité LP, QP et d'ordre supérieur à n'importe quel optimiseur donné, pas seulement à un solveur QP spécial. Deuxièmement, {[2] } peut digérer les mathématiques symboliques, donc la facilité de définir/entrer les contraintes est un peu plus agréable que de travailler avec les matrices et les dérivées de fonctions. mystic dépend de numpy, et d'utiliser scipy s'il est installé (toutefois, scipy n'est pas requis). mystic utilise sympy pour gérer les contraintes symboliques, mais ce n'est pas non plus nécessaire pour l'optimisation en général.

Sortie:

Optimization terminated successfully.
         Current function value: -2.000000
         Iterations: 3
         Function evaluations: 103

x* = [ 2.  1.]; f(x*) = 2.0

Obtenez mystic ici: https://github.com/uqfoundation

3
répondu Mike McKerns 2015-10-02 19:21:16