# Solve the Optimization Problem¶

A previous section formulated the optimization problem. This section presents the core Classiq capabilities, generates a designated quantum solution, and executes the generated algorithm on a quantum backend.

This example uses the method for a common problem: the Max Independent Set (MIS) with networkx.star_graph(4), solved by a QAOA Penalty algorithm with three QAOA layers. See the problem library.

import networkx as nx
import pyomo.core as pyo

def mis(graph: nx.Graph) -> pyo.ConcreteModel:
problem = pyo.ConcreteModel("mis")
problem.x = pyo.Var(graph.nodes, domain=pyo.Binary)

@problem.Constraint(graph.edges)
def independent_rule(problem, node1, node2):
return problem.x[node1] + problem.x[node2] <= 1

problem.cost = pyo.Objective(expr=sum(list(problem.x.values())), sense=pyo.maximize)

return problem


The method builds a PYOMO problem, indicates the QAOA and optimizer configuration, and constructs a Classiq model that defines the quantum algorithm to optimize this problem. Synthesize the model to get a quantum program to execute (like any other model).

First, create the desired optimization problem as a PYOMO model. The following code snippet is a concise example of the application of the optimization engine.

import networkx as nx

graph = nx.star_graph(4)
mis_problem = mis(graph)


Construct a Classiq model for this optimization problem with the construct_combinatorial_optimization_model function.

from classiq.applications.combinatorial_optimization import (
QAOAConfig,
)
from classiq import construct_combinatorial_optimization_model

qaoa_config = QAOAConfig(num_layers=3)
mis_model = construct_combinatorial_optimization_model(
pyo_model=mis_problem, qaoa_config=qaoa_config
)

qmod
{
"functions": [
{
"name": "main",
"param_decls": {
"params_list": {
"kind": "array",
"element_type": {
"kind": "real"
},
"size": 6
}
},
"port_declarations": {
"target": {
"name": "target",
"size": {
"expr": "len(get_field(optimization_problem_to_hamiltonian(get_type(MyCombiProblem), 2.0)[0], 'pauli'))"
},
"direction": "output"
}
},
"body": [
{
"function": "allocate",
"positional_args": [
{
"expr": "len(target)"
},
{
"name": "target"
}
]
},
{
"function": "qaoa_penalty",
"params": {
"hamiltonian": {
"expr": "optimization_problem_to_hamiltonian(get_type(MyCombiProblem), 2.0)"
},
"params_list": {
"expr": "params_list"
},
"num_qubits": {
"expr": "len(target)"
},
"is_st": {
"expr": "True"
}
},
"inouts": {
"target": {
"name": "target"
}
}
}
]
}
],
"types": [
{
"name": "MyCombiProblem",
"variables": {
"x0": {
"kind": "int"
},
"x1": {
"kind": "int"
},
"x2": {
"kind": "int"
},
"x3": {
"kind": "int"
},
"x4": {
"kind": "int"
}
},
"constraints": [
{
"expr": "x0 + x1 <= 1"
},
{
"expr": "x0 + x2 <= 1"
},
{
"expr": "x0 + x3 <= 1"
},
{
"expr": "x0 + x4 <= 1"
}
],
"objective_type": "Max",
"objective_function": {
"expr": "x0 + x1 + x2 + x3 + x4"
}
}
],
"classical_execution_code": "\nvqe_result = vqe(\n    hamiltonian=optimization_problem_to_hamiltonian(get_type(MyCombiProblem), 2.0),\n    maximize=True,\n    initial_point=compute_qaoa_initial_point(optimization_problem_to_hamiltonian(get_type(MyCombiProblem), 2.0),3),\n    optimizer=Optimizer.COBYLA,\n    max_iteration=0,\n    tolerance=0.0,\n    step_size=0.0,\n    skip_compute_variance=False,\n    alpha_cvar=1.0\n)\nsolution = get_optimization_solution(get_type(MyCombiProblem), vqe_result, 2.0)\nhamiltonian = optimization_problem_to_hamiltonian(get_type(MyCombiProblem), 2.0)\nsave({'solution': solution, \"vqe_result\": vqe_result, \"hamiltonian\": hamiltonian})\n"
}


## Results¶

The quantum program is obtained by synthesizing the model. You can also view the quantum program interactively.

from classiq import synthesize, show

mis_quantum_program = synthesize(mis_model)
show(mis_quantum_program)


To get the results for the optimization problem, execute the received quantum program.

from classiq import execute

res = execute(mis_quantum_program).result()


The results show all obtained solutions, with the number of times they were encountered. Each solution contains the found objective function value and the variable assignment that gave it.

import pandas as pd

optimization_result = pd.DataFrame.from_records(res[0].value)

  probability   cost    solution    count
1   0.218   4.0 [0, 1, 1, 1, 1] 218
9   0.014   3.0 [0, 1, 1, 0, 1] 14
8   0.014   3.0 [0, 1, 1, 1, 0] 14
11  0.009   3.0 [0, 0, 1, 1, 1] 9
10  0.011   3.0 [0, 1, 0, 1, 1] 11


#### Convergence Graph¶

from classiq.execution import VQESolverResult

vqe_result = res[1].value
vqe_result.convergence_graph


### Histogram¶

optimization_result.hist("cost", weights=optimization_result["probability"])


### Best Sampled Solution¶

idx = optimization_result.cost.idxmax()
print(
"x =", optimization_result.solution[idx], ", cost =", optimization_result.cost[idx]
)

x = [0, 1, 1, 1, 1] , cost = 4.0


### Compare to a Classical Solver¶

from pyomo.opt import SolverFactory

solver = SolverFactory("couenne")
solver.solve(mis_problem)

mis_problem.display()

Model mis

Variables:
x : Size=5, Index=x_index
Key : Lower : Value                  : Upper : Fixed : Stale : Domain
0 :     0 : 1.9949147640946208e-09 :     1 : False : False : Binary
1 :     0 :                    1.0 :     1 : False : False : Binary
2 :     0 :                    1.0 :     1 : False : False : Binary
3 :     0 :                    1.0 :     1 : False : False : Binary
4 :     0 :                    1.0 :     1 : False : False : Binary

Objectives:
cost : Size=1, Index=None, Active=True
Key  : Active : Value
None :   True : 4.0000000019949145

Constraints:
independent_rule : Size=4
Key    : Lower : Body               : Upper
(0, 1) :  None : 1.0000000019949147 :   1.0
(0, 2) :  None : 1.0000000019949147 :   1.0
(0, 3) :  None : 1.0000000019949147 :   1.0
(0, 4) :  None : 1.0000000019949147 :   1.0