Uso de Modelos y Solucionadores de Optimización Clásicos con Qiskit Optimization#

Podemos utilizar solucionadores de optimización clásicos (CPLEX y Gurobi) con Qiskit Optimization. Docplex y Gurobipy son las API de Python para CPLEX y Gurobi, respectivamente. Podemos cargar y guardar un modelo de optimización de Docplex y Gurobipy y podemos aplicar CPLEX y Gurobi a un QuadraticProgram.

Si deseas utilizar el solucionador de CPLEX, debes instalar pip install 'qiskit-optimization[cplex]'. Docplex se instala automáticamente, como dependencia, cuando instalas Qiskit Optimization.

Si deseas utilizar Gurobi y Gurobipy, debes instalar pip install 'qiskit-optimization[gurobi]'.

Nota: estos solucionadores, que se instalan mediante pip, son versiones gratuitas y vienen con algunas limitaciones, como la cantidad de variables. Los siguientes enlaces proporcionan más detalles:

CplexSolver y GurobiSolver#

Qiskit Optimization admite los solucionadores clásicos de CPLEX y Gurobi como CplexSolver y GurobiSolver, respectivamente. Podemos resolver un QuadraticProgram con CplexSolver y GurobiSolver de la siguiente manera.

[1]:
from qiskit_optimization.problems import QuadraticProgram

# define a problem
qp = QuadraticProgram()
qp.binary_var("x")
qp.integer_var(name="y", lowerbound=-1, upperbound=4)
qp.maximize(quadratic={("x", "y"): 1})
qp.linear_constraint({"x": 1, "y": -1}, "<=", 0)
print(qp.prettyprint())
Problem name:

Maximize
  x*y

Subject to
  Linear constraints (1)
    x - y <= 0  'c0'

  Integer variables (1)
    -1 <= y <= 4

  Binary variables (1)
    x

[2]:
from qiskit_optimization.algorithms import CplexOptimizer, GurobiOptimizer

cplex_result = CplexOptimizer().solve(qp)
gurobi_result = GurobiOptimizer().solve(qp)

print("cplex")
print(cplex_result.prettyprint())
print()
print("gurobi")
print(gurobi_result.prettyprint())
Restricted license - for non-production use only - expires 2025-11-24
cplex
objective function value: 4.0
variable values: x=1.0, y=4.0
status: SUCCESS

gurobi
objective function value: 4.0
variable values: x=1.0, y=4.0
status: SUCCESS

Podemos establecer el parámetro de resolución de CPLEX de la siguiente manera. Podemos mostrar el mensaje del solucionador de CPLEX configurando disp=True. Consulta Parameters of CPLEX para obtener detalles sobre los parámetros de CPLEX.

[3]:
result = CplexOptimizer(disp=True, cplex_parameters={"threads": 1, "timelimit": 0.1}).solve(qp)
print(result.prettyprint())
Version identifier: 22.1.1.0 | 2023-02-11 | 22d6266e5
CPXPARAM_Read_DataCheck                          1
CPXPARAM_Threads                                 1
CPXPARAM_TimeLimit                               0.10000000000000001
Found incumbent of value 0.000000 after 0.00 sec. (0.00 ticks)
Found incumbent of value 4.000000 after 0.00 sec. (0.00 ticks)

Root node processing (before b&c):
  Real time             =    0.00 sec. (0.00 ticks)
Sequential b&c:
  Real time             =    0.00 sec. (0.00 ticks)
                          ------------
Total (root+branch&cut) =    0.00 sec. (0.00 ticks)
objective function value: 4.0
variable values: x=1.0, y=4.0
status: SUCCESS

Obtenemos la misma solución óptima con QAOA de la siguiente manera.

[4]:
from qiskit_optimization.algorithms import MinimumEigenOptimizer

from qiskit_algorithms import QAOA
from qiskit_algorithms.optimizers import COBYLA
from qiskit.primitives import Sampler

meo = MinimumEigenOptimizer(QAOA(sampler=Sampler(), optimizer=COBYLA(maxiter=100)))
result = meo.solve(qp)
print(result.prettyprint())
print("\ndisplay the best 5 solution samples")
for sample in result.samples[:5]:
    print(sample)
objective function value: 4.0
variable values: x=1.0, y=4.0
status: SUCCESS

display the best 5 solution samples
SolutionSample(x=array([1., 4.]), fval=4.0, probability=0.024179976472345502, status=<OptimizationResultStatus.SUCCESS: 0>)
SolutionSample(x=array([1., 3.]), fval=3.0, probability=0.05819514104929059, status=<OptimizationResultStatus.SUCCESS: 0>)
SolutionSample(x=array([1., 2.]), fval=2.0, probability=0.16262615426689528, status=<OptimizationResultStatus.SUCCESS: 0>)
SolutionSample(x=array([1., 1.]), fval=1.0, probability=0.1342595428408242, status=<OptimizationResultStatus.SUCCESS: 0>)
SolutionSample(x=array([0., 0.]), fval=0.0, probability=0.048555539627953996, status=<OptimizationResultStatus.SUCCESS: 0>)

Traductores entre QuadraticProgram y Docplex/Gurobipy#

Qiskit Optimization puede cargar un QuadraticProgram desde un modelo Docplex y un modelo Gurobipy.

Primero, definimos un problema de optimización por Docplex y Gurobipy.

[5]:
# docplex model
from docplex.mp.model import Model

docplex_model = Model("docplex")
x = docplex_model.binary_var("x")
y = docplex_model.integer_var(-1, 4, "y")
docplex_model.maximize(x * y)
docplex_model.add_constraint(x <= y)
docplex_model.prettyprint()
// This file has been generated by DOcplex
// model name is: docplex
// single vars section
dvar bool x;
dvar int y;

maximize
 [ x*y ];

subject to {
 x <= y;

}
[23]:
# gurobi model
from io import StringIO
from tempfile import NamedTemporaryFile
from os import remove

import gurobipy as gp


def gpy_display(mdl):
    """Convenience function to pretty-print a Gurobipy-Model."""
    with NamedTemporaryFile(suffix=".lp") as tmp_file:
        mdl.write(tmp_file.name)

        with open(tmp_file.name, "r") as f:
            print(f.read())


gurobipy_model = gp.Model("gurobi")
x = gurobipy_model.addVar(vtype=gp.GRB.BINARY, name="x")
y = gurobipy_model.addVar(vtype=gp.GRB.INTEGER, lb=-1, ub=4, name="y")
gurobipy_model.setObjective(x * y, gp.GRB.MAXIMIZE)
gurobipy_model.addConstr(x - y <= 0)
gurobipy_model.update()
gpy_display(gurobipy_model)
\ Model gurobi
\ LP format - for model browsing. Use MPS format to capture full model detail.
Maximize
 [ 2 x * y ] / 2
Subject To
 R0: x - y <= 0
Bounds
 -1 <= y <= 4
Binaries
 x
Generals
 y
End

Podemos generar un objeto QuadraticProgram a partir de los modelos Docplex y Gurobipy. Vemos que los dos objetos QuadraticProgram generados a partir de Docplex y Gurobipy son idénticos.

[7]:
from qiskit_optimization.translators import from_docplex_mp, from_gurobipy

qp = from_docplex_mp(docplex_model)
print("QuadraticProgram obtained from docpblex")
print(qp.prettyprint())
print("-------------")
print("QuadraticProgram obtained from gurobipy")
qp2 = from_gurobipy(gurobipy_model)
print(qp2.prettyprint())
QuadraticProgram obtained from docpblex
Problem name: docplex

Maximize
  x*y

Subject to
  Linear constraints (1)
    x - y <= 0  'c0'

  Integer variables (1)
    -1 <= y <= 4

  Binary variables (1)
    x

-------------
QuadraticProgram obtained from gurobipy
Problem name: gurobi

Maximize
  x*y

Subject to
  Linear constraints (1)
    x - y <= 0  'R0'

  Integer variables (1)
    -1 <= y <= 4

  Binary variables (1)
    x

También podemos generar un modelo Docplex y un modelo Gurobipy de un QuadraticProgram.

[25]:
from qiskit_optimization.translators import to_gurobipy, to_docplex_mp

gmod = to_gurobipy(from_docplex_mp(docplex_model))
print("convert docplex to gurobipy via QuadraticProgram")
gpy_display(gmod)

dmod = to_docplex_mp(from_gurobipy(gurobipy_model))
print("\nconvert gurobipy to docplex via QuadraticProgram")
print(dmod.export_as_lp_string())
convert docplex to gurobipy via QuadraticProgram
\ Model docplex
\ LP format - for model browsing. Use MPS format to capture full model detail.
Maximize
 [ 2 x * y ] / 2
Subject To
 c0: x - y <= 0
Bounds
 -1 <= y <= 4
Binaries
 x
Generals
 y
End


convert gurobipy to docplex via QuadraticProgram
\ This file has been generated by DOcplex
\ ENCODING=ISO-8859-1
\Problem name: gurobi

Maximize
 obj: [ 2 x*y ]/2
Subject To
 R0: x - y <= 0

Bounds
 0 <= x <= 1
 -1 <= y <= 4

Binaries
 x

Generals
 y
End

Restricciones de los indicadores de Docplex#

from_docplex_mp admite restricciones de indicador, por ejemplo, u = 0 => x + y <= z (u: variable binaria) cuando convertimos un modelo Docplex en un QuadraticProgram. Convierte las restricciones de los indicadores en restricciones lineales utilizando la formulación big-M.

[9]:
ind_mod = Model("docplex")
x = ind_mod.binary_var("x")
y = ind_mod.integer_var(-1, 2, "y")
z = ind_mod.integer_var(-1, 2, "z")
ind_mod.maximize(3 * x + y - z)
ind_mod.add_indicator(x, y >= z, 1)
print(ind_mod.export_as_lp_string())
\ This file has been generated by DOcplex
\ ENCODING=ISO-8859-1
\Problem name: docplex

Maximize
 obj: 3 x + y - z
Subject To
 lc1: x = 1 -> y - z >= 0

Bounds
 0 <= x <= 1
 -1 <= y <= 2
 -1 <= z <= 2

Binaries
 x

Generals
 y z
End

Comparemos las soluciones del modelo con una restricción de indicador por

  1. aplicando CPLEX directamente al modelo Docplex (sin traducirlo a un QuadraticProgram. El solucionador CPLEX admite de forma nativa las restricciones del indicador),

  2. aplicando QAOA a un QuadraticProgram obtenido por from_docplex_mp.

Vemos que las soluciones son las mismas.

[10]:
qp = from_docplex_mp(ind_mod)
result = meo.solve(qp)  # apply QAOA to QuadraticProgram
print("QAOA")
print(result.prettyprint())
print("-----\nCPLEX")
print(ind_mod.solve())  # apply CPLEX directly to the Docplex model
QAOA
objective function value: 6.0
variable values: x=1.0, y=2.0, z=-1.0
status: SUCCESS
-----
CPLEX
solution for: docplex
objective: 6
status: OPTIMAL_SOLUTION(2)
x=1
y=2
z=-1

[11]:
import qiskit.tools.jupyter

%qiskit_version_table
%qiskit_copyright

Version Information

SoftwareVersion
qiskit0.45.1
qiskit_optimization0.7.0
qiskit_algorithms0.2.1
System information
Python version3.10.13
Python compilerGCC 13.2.1 20230728 (Red Hat 13.2.1-1)
Python buildmain, Aug 28 2023 00:00:00
OSLinux
CPUs12
Memory (Gb)31.056411743164062
Wed Dec 06 10:30:30 2023 CET

This code is a part of Qiskit

© Copyright IBM 2017, 2023.

This code is licensed under the Apache License, Version 2.0. You may
obtain a copy of this license in the LICENSE.txt file in the root directory
of this source tree or at http://www.apache.org/licenses/LICENSE-2.0.

Any modifications or derivative works of this code must retain this
copyright notice, and modified files need to carry a notice indicating
that they have been altered from the originals.