Nota

Esta página foi gerada a partir do ‘tutoriais/otimizacao/5_otimizador_admm.ipynb’__.

Execute interativamente no IBM Quantum lab.

Otimizador ADMM

Introdução

O Otimizador ADMM pode resolver classes de problemas de otimização mistos-binários restritos, daqui em diante (MBCO), que geralmente aparecem em pesquisas logísticas, financeiras e de operações. Em particular, o Otimizador ADMM aqui projetado pode resolver o seguinte problema de otimização \((P)\):

\[\min_{x \in \mathcal{X},u\in\mathcal{U} \subseteq \mathbb{R}^l } \quad q(x) + \varphi(u),\]

sujeito às restrições:

\[\mathrm{s.t.:~} \quad G x = b, \quad g(x) \leq 0, \quad \ell(x, u) \leq 0,\]

com as premissas funcionais correspondentes.

1. A função \(q: \mathbb{R}^n \to \mathbb{R}\) é quadrática, p.ex., \(q(x) = x^{\intercal} Q x + a^{\intercal} x\) para uma determinada matriz quadrada simétrica \(Q \in \mathbb{R}^n \times \mathbb{R}^n, Q = Q^{\intercal}\), e vetor \(a \in \mathbb{R}^n\);

2. O conjunto \(\mathcal{X} = \{0,1\}^n = \{x_{(i)} (1-x_{(i)}) = 0, \forall i\}\) impõe as restrições binárias;

3. Matriz \(G\in\mathbb{R}^n \times \mathbb{R}^{n'}\), vector \(b \in \mathbb{R}^{n'}\), e função \(g: \mathbb{R}^n \to \mathbb{R}\) são convexos;

4. Função \(\varphi: \mathbb{R}^l \to \mathbb{R}\) é convexa e \(\mathcal{U}\) é um conjunto convexo;

5. Função \(\ell: \mathbb{R}^n\times \mathbb{R}^l \to \mathbb{R}\) é conjuntamente convexa em \(x, u\).

Para resolver problemas com MBO, [1] Heurística proposta para \((P)\) baseado no Método Alternativo de Direção dos Multiplicadores (ADMM) [2]. ADMM é um operador que divide algoritmos com um longo histórico na otimização convexa, e é sabido que tem propriedades de convergência residual, objetiva e dual de variáveis, desde que os pressupostos de convexidade sejam válidos.

O método [1] (referido como 3-ADMM-H) aproveita o procedimento de divisão do operador ADMM para conceber uma decomposição para certas classes de MBOs em: - um subproblema QUBO para ser resolvido no dispositivo quântico através de algoritmos variacionais, tais como VQE ou QAOA; - subproblema contínuo convexo restritos, que pode ser resolvido eficientemente com os solucionadores de otimização clássica.

O algoritmo 3-ADMM-H funciona da seguinte forma:

0. Fase de inicialização (defina os parâmetros, o QUBO e os solucionadores convexos);

1. For each ADMM iterations ($k = 1, 2, \ldots, $) until termination:

   • Resolva um subproblema QUBO devidamente definido (com um solucionador clássico ou quântico);

   • Resolva problemas convexos devidamente definidos (com um solucionador clássico);

   • Atualize as variáveis duais.

2. Retornar otimizadores e custo.

Uma discussão abrangente sobre as condições para a convergência, viabilidade e otimização do algoritmo podem ser encontrados em [1]. Uma variante com 2 blocos ADMM, isto é, um subproblema da QUBO e um subproblema contínuo convexo restrito, também são introduzidos em [1].

Referências

[1] C. Gambella and A. Simonetto, Multi-block ADMM heuristics for mixed-binary optimization, on classical and quantum computers, arXiv preprint arXiv:2001.02069 (2020).

[2] S. Boyd, N. Parikh, E. Chu, B. Peleato, and J. Eckstein, Distributed optimization and statistical learning via the alternating direction method of multipliers, Foundations and Trends in Machine learning, 3, 1–122 (2011).

Inicialização

Em primeiro lugar, carregamos todos os pacotes de que precisamos.

import time
from typing import List, Optional, Any
import numpy as np
import matplotlib.pyplot as plt

from docplex.mp.model import Model

from qiskit import BasicAer
from qiskit.aqua.algorithms import QAOA, NumPyMinimumEigensolver
from qiskit.optimization.algorithms import CobylaOptimizer, MinimumEigenOptimizer
from qiskit.optimization.problems import QuadraticProgram
from qiskit.optimization.algorithms.admm_optimizer import ADMMParameters, ADMMOptimizer

# If CPLEX is installed, you can uncomment this line to import the CplexOptimizer.
# CPLEX can be used in this tutorial to solve the convex continuous problem,
# but also as a reference to solve the QUBO, or even the full problem.
#
# from qiskit.optimization.algorithms import CplexOptimizer

Primeiramente, inicializamos todos os algoritmos que planejamos usar durante este tutorial.

Para resolver os problemas QUBO podemos escolher entre - MinimumEigenOptimizer usando diferentes MinimumEigensolver, tais como VQE, QAOA ou NumpyMinimumEigensolver (clássico) - GroverOptimizer - CplexOptimizer (clássico, se CPLEX estiver instalado)

e para resolver os problemas contínuos convexos, podemos escolher entre os seguintes solucionadores clássicos: - CplexOptimizer (se o CPLEX estiver instalado) - CobylaOptimizer

Caso o CPLEX não esteja disponível, o CobylaOptimizer (para problemas contínuos convexos) e o MinimumEigenOptimizer usando o NumpyMinimumEigensolver (para QUBOs) podem ser usados como alternativas clássicas ao CPLEX para teste, validação e benchmarking.

# define COBYLA optimizer to handle convex continuous problems.
cobyla = CobylaOptimizer()

# define QAOA via the minimum eigen optimizer
qaoa = MinimumEigenOptimizer(QAOA(quantum_instance=BasicAer.get_backend('statevector_simulator')))

# exact QUBO solver as classical benchmark
exact = MinimumEigenOptimizer(NumPyMinimumEigensolver()) # to solve QUBOs

# in case CPLEX is installed it can also be used for the convex problems, the QUBO,
# or as a benchmark for the full problem.
#
# cplex = CplexOptimizer()

Exemplo

Testamos o algoritmo 3-ADMM-H em um problema Quadrático Mixed-Binary simples com restrições de igualdade e desigualdade (Exemplo 6 relatado em [1]). Primeiramente, construímos um problema docplex e, em seguida, carregamos em um QuadraticProgram.

# construct model using docplex
mdl = Model('ex6')

v = mdl.binary_var(name='v')
w = mdl.binary_var(name='w')
t = mdl.binary_var(name='t')
u = mdl.continuous_var(name='u')

mdl.minimize(v + w + t + 5 * (u-2)**2)
mdl.add_constraint(v + 2 * w + t + u <= 3, "cons1")
mdl.add_constraint(v + w + t >= 1, "cons2")
mdl.add_constraint(v + w == 1, "cons3")

# load quadratic program from docplex model
qp = QuadraticProgram()
qp.from_docplex(mdl)
print(qp.export_as_lp_string())

\ This file has been generated by DOcplex
\ ENCODING=ISO-8859-1
\Problem name: ex6

Minimize
 obj: v + w + t - 20 u + [ 10 u^2 ]/2 + 20
Subject To
 cons1: v + 2 w + t + u <= 3
 cons2: v + w + t >= 1
 cons3: v + w = 1

Bounds
 0 <= v <= 1
 0 <= w <= 1
 0 <= t <= 1

Binaries
 v w t
End

Solução clássica

3-ADMM-H precisa de um otimizador QUBO para resolver o subproblema do QUBO e um otimizador contínuo para resolver o subproblema contínuo convexo restrito. Primeiramente, resolvemos o problema de forma clássica: usamos o MinimumEigenOptimizer com o NumPyMinimumEigenSolver como um solucionador QUBO clássico e exato e usamos o CobylaOptimizer como um solucionador contínuo convexo. 3-ADMM-H suporta qualquer outro solucionador adequado disponível no Qiskit. Por exemplo, VQE, QAOA e GroverOptimizer podem ser chamados como solucionadores quânticos, como demonstrado mais tarde. Se o CPLEX estiver instalado, o CplexOptimizer também pode ser usado tanto como um solucionador QUBO ou convexo.

Parâmetros

Os 3-ADMM-H são envolvidos na classe ADMMParameters. Valores de parâmetros customizados podem ser definidos como argumentos da classe. Neste exemplo, os parâmetros \(\rho, \beta\) são inicializados para \(1001\) e \(1000\), respectivamente. A penalização factor_c das restrições de igualdade \(Gx = b\) são definidos para \(900\). A tolerância tol para convergência residual primitiva está ajustada para 1.e-6. Neste caso, a implementação de 3 blocos está garantida para convergir para o Theorem 4 de [1], porque a restrição de desigualdade com a variável contínua é sempre ativa. A implementação de dois blocos pode ser executada definindo three_block=False, e praticamente converge para uma solução viável e não otimizada.

admm_params = ADMMParameters(
                            rho_initial=1001,
                            beta=1000,
                            factor_c=900,
                            maxiter=100,
                            three_block=True, tol=1.e-6
                        )

Chamando algoritmo 3-ADMM-H

Para invocar o algoritmo 3-ADMM-H, uma instância da classe ADMMOptimizer precisa ser criada. Isto gera parâmetros específicos de ADMM e otimizadores de subproblemas separadamente para o construtor. A solução retornada é uma instância da classe OptimizationResult.

# define QUBO optimizer
qubo_optimizer = exact
# qubo_optimizer = cplex  # uncomment to use CPLEX instead

# define classical optimizer
convex_optimizer = cobyla
# convex_optimizer = cplex  # uncomment to use CPLEX instead

# initialize ADMM with classical QUBO and convex optimizer
admm = ADMMOptimizer(params=admm_params,
                     qubo_optimizer=qubo_optimizer,
                     continuous_optimizer=convex_optimizer)

# run ADMM to solve problem
result = admm.solve(qp)

Resultado do Solver Clássico

A solução 3-ADMM-H pode ser impressa e visualizada. O atributo x da solução contém respectivamente, os valores das variáveis de decisão binária e os valores das variáveis de decisão contínua. O fval é o valor objetivo da solução.

print("x={}".format(result.x))
print("fval={:.2f}".format(result.fval))

x=[1. 0. 0. 2.]
fval=1.00

As estatísticas de solução podem ser acessadas e visualizadas no campo state. Aqui mostramos a convergência de 3-ADMM-H, em termos de resíduos primitivos.

plt.plot(result.state.residuals)
plt.xlabel("Iterations")
plt.ylabel("Residuals")
plt.show()

Solução Quântica

Agora resolvemos o mesmo problema de otimização com QAOA que o otimizador de QUBO, rodando em dispositivos quânticos simulados. Primeiro, é necessário selecionar o otimizador clássico do eigensolver QAOA. Em seguida, o backend da simulação é definido. Finalmente, o eigensolver é envolvido na classe MinimumEigenOptimizer. Uma nova instância de ADMMOptimizer é preenchida com QAOA como otimizador de QUBO.

# define QUBO optimizer
qubo_optimizer = qaoa

# define classical optimizer
convex_optimizer = cobyla
# convex_optimizer = cplex  # uncomment to use CPLEX instead

# initialize ADMM with quantum QUBO optimizer and classical convex optimizer
admm_q = ADMMOptimizer(params=admm_params,
                       qubo_optimizer=qubo_optimizer,
                       continuous_optimizer=convex_optimizer)

# run ADMM to solve problem
result_q = admm_q.solve(qp)

Resultados de Solver Quânticos

Apresentamos aqui os resultados obtidos a partir da solução quântica. Como no exemplo acima de x significa a solução, o fval é o valor objetivo.

print("x={}".format(result_q.x))
print("fval={:.2f}".format(result_q.fval))

x=[1. 0. 0. 2.]
fval=1.00

plt.clf()
plt.plot(result_q.state.residuals)
plt.xlabel("Iterations")
plt.ylabel("Residuals")
plt.show()

import qiskit.tools.jupyter
%qiskit_version_table
%qiskit_copyright

Version Information

Qiskit Software	Version
Qiskit	None
Terra	0.15.1
Aer	0.6.1
Ignis	0.4.0
Aqua	0.7.5
IBM Q Provider	0.8.0
System information
Python	3.7.4 (default, Aug 13 2019, 15:17:50) [Clang 4.0.1 (tags/RELEASE_401/final)]
OS	Darwin
CPUs	2
Memory (Gb)	12.0
Fri Aug 28 09:24:00 2020 EDT

This code is a part of Qiskit

© Copyright IBM 2017, 2020.

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.