diff --git a/examples/QAOA_LABS_optimization.ipynb b/examples/QAOA_LABS_optimization.ipynb index a658870d6..ba527b652 100644 --- a/examples/QAOA_LABS_optimization.ipynb +++ b/examples/QAOA_LABS_optimization.ipynb @@ -117,8 +117,8 @@ "" ], "text/plain": [ - " N p overlap gamma \n", - "35 10 1 0.127219 [0.0788139605] \\\n", + " N p overlap gamma \\\n", + "35 10 1 0.127219 [0.0788139605] \n", "36 10 2 0.196232 [0.06896497310000001, 0.1511922366] \n", "37 10 3 0.258441 [0.063816785, 0.1396294913, 0.1538820941] \n", "38 10 4 0.318326 [0.0655247286, 0.1319443875, 0.1430007569, 0.1... \n", @@ -139,7 +139,7 @@ ], "source": [ "N = 10\n", - "parameters = parameter_utils.get_LABS_wrt_best_overlap(N)\n", + "parameters = parameter_utils.get_best_known_parameters_for_LABS_wrt_overlap(N)\n", "known_p = parameters.p.max()\n", "print(f\"Maximum p available for N={N} is {known_p}\")\n", "parameters.head()" @@ -159,7 +159,7 @@ "outputs": [], "source": [ "p = known_p + 1\n", - "gamma, beta = parameter_utils.get_LABS_wrt_best_overlap_for_p(N, known_p)" + "gamma, beta = parameter_utils.get_best_known_parameters_for_LABS_wrt_overlap_for_p(N, known_p)" ] }, { @@ -188,20 +188,6 @@ "execution_count": 5, "metadata": {}, "outputs": [ - { - "name": "stderr", - "output_type": "stream", - "text": [ - "/usr/lib/python3.11/site-packages/numba/cuda/dispatcher.py:539: NumbaPerformanceWarning: \u001b[1mGrid size 1 will likely result in GPU under-utilization due to low occupancy.\u001b[0m\n", - " warn(NumbaPerformanceWarning(msg))\n", - "/usr/lib/python3.11/site-packages/numba/cuda/dispatcher.py:539: NumbaPerformanceWarning: \u001b[1mGrid size 1 will likely result in GPU under-utilization due to low occupancy.\u001b[0m\n", - " warn(NumbaPerformanceWarning(msg))\n", - "/usr/lib/python3.11/site-packages/numba/cuda/dispatcher.py:539: NumbaPerformanceWarning: \u001b[1mGrid size 1 will likely result in GPU under-utilization due to low occupancy.\u001b[0m\n", - " warn(NumbaPerformanceWarning(msg))\n", - "/usr/lib/python3.11/site-packages/numba/cuda/dispatcher.py:539: NumbaPerformanceWarning: \u001b[1mGrid size 1 will likely result in GPU under-utilization due to low occupancy.\u001b[0m\n", - " warn(NumbaPerformanceWarning(msg))\n" - ] - }, { "name": "stdout", "output_type": "stream", @@ -241,7 +227,7 @@ "name": "stdout", "output_type": "stream", "text": [ - "Success probability at p=9 after optimization is 0.5242988366332573\n" + "Success probability at p=9 after optimization is 0.5242896565265064\n" ] } ], @@ -264,9 +250,9 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3 (ipykernel)", + "display_name": "py39", "language": "python", - "name": "python3" + "name": "py39" }, "language_info": { "codemirror_mode": { @@ -278,9 +264,9 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.11.3" + "version": "3.9.4" } }, "nbformat": 4, - "nbformat_minor": 2 + "nbformat_minor": 4 } diff --git a/examples/QAOA_portfolio_optimization.ipynb b/examples/QAOA_portfolio_optimization.ipynb index eb3bc4f4a..9986743c2 100644 --- a/examples/QAOA_portfolio_optimization.ipynb +++ b/examples/QAOA_portfolio_optimization.ipynb @@ -36,7 +36,7 @@ }, { "cell_type": "code", - "execution_count": 1, + "execution_count": null, "id": "c3de09e7-5315-4170-8263-6041a7819fac", "metadata": { "tags": [] @@ -58,7 +58,7 @@ }, { "cell_type": "code", - "execution_count": 2, + "execution_count": null, "id": "57544f68-6538-4dda-a47b-254252caa54a", "metadata": { "tags": [] @@ -72,7 +72,7 @@ }, { "cell_type": "code", - "execution_count": 3, + "execution_count": null, "id": "e49f9152-e49c-4cbd-a2b3-20894a927a9a", "metadata": { "tags": [] @@ -84,7 +84,7 @@ }, { "cell_type": "code", - "execution_count": 4, + "execution_count": null, "id": "a800d2c6-d7a9-449a-befc-e1620676847e", "metadata": { "tags": [] @@ -96,7 +96,7 @@ }, { "cell_type": "code", - "execution_count": 5, + "execution_count": null, "id": "60942736-20e4-4519-a584-65f613b7bddf", "metadata": { "tags": [] @@ -118,7 +118,7 @@ }, { "cell_type": "code", - "execution_count": 6, + "execution_count": null, "id": "c4cb1e0a-f0eb-48e4-a1dc-8d342462c156", "metadata": { "tags": [] @@ -131,7 +131,7 @@ }, { "cell_type": "code", - "execution_count": 7, + "execution_count": null, "id": "1fa46296-89cb-40e7-aaf8-09b708eda12d", "metadata": { "tags": [] @@ -145,7 +145,7 @@ }, { "cell_type": "code", - "execution_count": 8, + "execution_count": null, "id": "4bb70cb8-db79-45b3-b808-96b7d6647b68", "metadata": { "tags": [] @@ -158,20 +158,12 @@ }, { "cell_type": "code", - "execution_count": 9, + "execution_count": null, "id": "238479eb-db75-4e19-a159-36a2bed9c2ae", "metadata": { "tags": [] }, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "energy = 0.004962484517944697, Approximation ratio = 0.7564585339722991\n" - ] - } - ], + "outputs": [], "source": [ "po_energy = qaoa_obj(x0).real\n", "po_ar = (po_energy-best_portfolio[1])/(best_portfolio[0]-best_portfolio[1])\n", @@ -185,12 +177,14 @@ "source": [ "# Optimize QAOA parameters\n", "\n", - "Note that we are using NLopt for optimization as it supports better-performing BOBYQA optimizer. Run `pip install nlopt` to install this dependency." + "Note that we are using NLopt for optimization as it supports better-performing BOBYQA optimizer. \n", + "\n", + "Run `pip install nlopt` to install this dependency." ] }, { "cell_type": "code", - "execution_count": 10, + "execution_count": null, "id": "9d730369-3d9a-4e0e-bd2b-a81d92baa7b7", "metadata": { "tags": [] @@ -203,7 +197,7 @@ }, { "cell_type": "code", - "execution_count": 11, + "execution_count": null, "id": "55a761e4-2122-4b00-bc1f-39730d6a3dd7", "metadata": { "tags": [] @@ -230,20 +224,12 @@ }, { "cell_type": "code", - "execution_count": 12, + "execution_count": null, "id": "86b7ff7f-43f4-4242-95b9-82784f0632a5", "metadata": { "tags": [] }, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "energy = 0.004135567715736568, Approximation ratio = 0.8050399835268952\n" - ] - } - ], + "outputs": [], "source": [ "_, opt_energy = minimize_nlopt(qaoa_obj, x0, p=1, rhobeg=0.01/1)\n", "opt_ar = (opt_energy-best_portfolio[1])/(best_portfolio[0]-best_portfolio[1])\n", @@ -252,20 +238,12 @@ }, { "cell_type": "code", - "execution_count": 13, + "execution_count": null, "id": "90c5c6af-003c-4724-a943-c747bd24f3c6", "metadata": { "tags": [] }, - "outputs": [ - { - "name": "stdout", - "output_type": "stream", - "text": [ - "energy = 0.004135567784198798, Approximation ratio = 0.8050399795047319\n" - ] - } - ], + "outputs": [], "source": [ "res = minimize(qaoa_obj, x0, method='COBYLA', options={'rhobeg':0.001})\n", "print(f\"energy = {res.fun}, Approximation ratio = {(res.fun-best_portfolio[1])/(best_portfolio[0]-best_portfolio[1])}\")" @@ -281,7 +259,7 @@ }, { "cell_type": "code", - "execution_count": 14, + "execution_count": null, "id": "0c6e7780-0392-4a9f-b0a7-4b23c0af84dc", "metadata": { "tags": [] @@ -293,7 +271,7 @@ }, { "cell_type": "code", - "execution_count": 15, + "execution_count": null, "id": "d5fec9f4-09e0-409f-a613-07f23194c94d", "metadata": { "tags": [] @@ -322,7 +300,7 @@ }, { "cell_type": "code", - "execution_count": 16, + "execution_count": null, "id": "bc26bf47-72ce-4944-9a9f-d217e4c77894", "metadata": { "tags": [] @@ -336,7 +314,7 @@ }, { "cell_type": "code", - "execution_count": 17, + "execution_count": null, "id": "0c4f35dc-9402-48ce-89ef-c5a80883aaf8", "metadata": { "tags": [] @@ -351,7 +329,7 @@ }, { "cell_type": "code", - "execution_count": 18, + "execution_count": null, "id": "8581ac1f-2675-4ccb-b311-ee7e658b2782", "metadata": { "tags": [] @@ -367,7 +345,7 @@ }, { "cell_type": "code", - "execution_count": 19, + "execution_count": null, "id": "2ceed19a-9391-4eed-95b9-dbfc3b8d6a9d", "metadata": { "tags": [] @@ -391,9 +369,9 @@ ], "metadata": { "kernelspec": { - "display_name": "Python 3 (ipykernel)", + "display_name": "qokit", "language": "python", - "name": "python3" + "name": "qokit" }, "language_info": { "codemirror_mode": { @@ -405,7 +383,7 @@ "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", - "version": "3.8.11" + "version": "3.9.4" } }, "nbformat": 4, diff --git a/qokit/assets/transferred_Dan_mean_0317.json b/qokit/assets/transferred_Dan_mean_0317.json index 1e80bf2f9..3b18a4a2f 100644 --- a/qokit/assets/transferred_Dan_mean_0317.json +++ b/qokit/assets/transferred_Dan_mean_0317.json @@ -16677,7 +16677,6 @@ 1.0, 1.0 ] -<<<<<<< HEAD:assets/transferred_Dan_mean_0317.json }, "278":{ "overlap transferred":0.0000138523, @@ -17218,7 +17217,5 @@ 1.0, 1.0 ] -======= ->>>>>>> main:qokit/assets/transferred_Dan_mean_0317.json } -} \ No newline at end of file +} diff --git a/qokit/parameter_utils.py b/qokit/parameter_utils.py index ea428c0b9..e0923907f 100644 --- a/qokit/parameter_utils.py +++ b/qokit/parameter_utils.py @@ -142,23 +142,6 @@ def convert_to_gamma_beta(*args, parameterization: QAOAParameterization | str): return gamma, beta -def set_parameterized_func( - parameterization: QAOAParameterization | str, probabilities_from_gamma_beta: Callable, compute_objective_from_probabilities: Callable -): - """ - parameterization: the approach to parametrize the QAOA parameters - probabilities_from_gamma_beta: callable function to estimate state vector from QAOA parameters - compute_objective_from_probabilities: callable function to computate objective from state vector - """ - - def f(*args): - gamma, beta = convert_to_gamma_beta(*args, parameterization=parameterization) - probabilities = probabilities_from_gamma_beta(gamma, beta) - return compute_objective_from_probabilities(probabilities) - - return f - - def get_sk_gamma_beta(p, parameterization: QAOAParameterization | str = "gamma beta"): """ Load the look-up table for initial points from diff --git a/qokit/portfolio_optimization.py b/qokit/portfolio_optimization.py index 18d3bd4d8..eb67c9db4 100644 --- a/qokit/portfolio_optimization.py +++ b/qokit/portfolio_optimization.py @@ -1,14 +1,12 @@ -############################################################################### -# // SPDX-License-Identifier: Apache-2.0 -# // Copyright : JP Morgan Chase & Co -############################################################################### """ Helper functions for the portfolio optimization problem """ +from __future__ import annotations +from collections.abc import Sequence import numpy as np from functools import partial import itertools -from typing import Optional +from typing import Any from qokit.parameter_utils import get_sk_gamma_beta @@ -57,7 +55,7 @@ def kbits(N, K): yield np.array(s) -def portfolio_brute_force(po_problem: dict, return_bitstring=False) -> (float, float, float): +def portfolio_brute_force(po_problem: dict, return_bitstring=False) -> tuple[float, float, float] | tuple[float, float, float, float]: N = po_problem["N"] K = po_problem["K"] min_constrained = float("inf") @@ -82,7 +80,7 @@ def portfolio_brute_force(po_problem: dict, return_bitstring=False) -> (float, f return min_constrained, min_x, max_constrained, max_x, mean_constrained -def get_data(N, seed=1, real=False) -> (float, float): +def get_data(N, seed=1, real=False) -> tuple[float, float]: """ load portofolio data from qiskit-finance (Yahoo) https://github.com/Qiskit/qiskit-finance/blob/main/docs/tutorials/11_time_series.ipynb @@ -151,7 +149,7 @@ def get_data(N, seed=1, real=False) -> (float, float): return means, cov -def get_problem(N, K, q, seed=1, pre=False) -> dict: +def get_problem(N, K, q, seed=1, pre=False) -> dict[str, Any]: """generate the portofolio optimziation problem dict""" po_problem = {} po_problem["N"] = N @@ -284,7 +282,7 @@ def get_sk_ini(p: int): return X0 -def alignment_para_to_qokit_scale(gammas: Optional[list], betas: Optional[list]): +def alignment_para_to_qokit_scale(gammas: Sequence[float] | None, betas: Sequence[float] | None): """Converts from format in alignment project into the scale used in qokit """ diff --git a/qokit/qaoa_objective.py b/qokit/qaoa_objective.py index f7b2e4199..33bea38a8 100644 --- a/qokit/qaoa_objective.py +++ b/qokit/qaoa_objective.py @@ -12,7 +12,7 @@ from .fur import choose_simulator, choose_simulator_xyring, QAOAFastSimulatorBase import typing -from .parameter_utils import from_fourier_basis, set_parameterized_func +from .parameter_utils import from_fourier_basis import qokit.parameter_utils from qokit.parameter_utils import QAOAParameterization from .qaoa_circuit_portfolio import measure_circuit @@ -38,6 +38,7 @@ def _get_qiskit_objective( objective: str = "expectation", terms=None, parameterization: str | QAOAParameterization = "theta", + mixer: str = "x", ): N = parameterized_circuit.num_qubits if objective == "expectation": @@ -85,15 +86,29 @@ def compute_objective_from_probabilities(probabilities): return en, 1 - overlap else: - raise ValueError(f"Unknown objective passed to get_qaoa_labs_objective: {objective}, allowed ['expectation', 'overlap', 'expectation and overlap']") + raise ValueError(f"Unknown objective passed to get_qaoa_objective: {objective}, allowed ['expectation', 'overlap', 'expectation and overlap']") - backend = Aer.get_backend("aer_simulator_statevector") + if mixer == "x": + backend = Aer.get_backend("aer_simulator_statevector") - def g(gamma, beta): - qc = parameterized_circuit.bind_parameters(list(np.hstack([beta, gamma]))) - sv = np.asarray(backend.run(qc).result().get_statevector()) - probs = np.abs(sv) ** 2 - return compute_objective_from_probabilities(probs) + def g(gamma, beta): + qc = parameterized_circuit.bind_parameters(list(np.hstack([beta, gamma]))) + sv = np.asarray(backend.run(qc).result().get_statevector()) + probs = np.abs(sv) ** 2 + return compute_objective_from_probabilities(probs) + + elif mixer == "xy": + backend = Aer.get_backend("statevector_simulator") + + def g(gamma, beta): + qc = parameterized_circuit.bind_parameters(np.hstack([np.asarray(beta) / 2, np.asarray(gamma) / 2])) + result = execute(qc, backend).result() + sv = reverse_array_index_bit_order(result.get_statevector()) + probs = np.abs(sv) ** 2 + return compute_objective_from_probabilities(probs) + + else: + raise ValueError(f"Unknown mixer type passed to get_qaoa_objective: {mixer}, allowed ['x', 'xy']") return g @@ -109,6 +124,9 @@ def get_qaoa_objective( objective: str = "expectation", parameterized_circuit=None, simulator: str = "auto", + mixer: str = "x", + initial_state: np.ndarray | None = None, + n_trotters: int = 1, ) -> typing.Callable: """Return QAOA objective to be minimized @@ -132,6 +150,13 @@ def get_qaoa_objective( If simulator == 'auto', implementation is chosen automatically (either the fastest CPU simulator or a GPU simulator if CUDA is available) If simulator == 'qiskit', implementation in qaoa_qiskit is used + mixer : str + If mixer == 'x', then uses the default Pauli X as the mixer + If mixer == 'xy', then uses the ring-XY as the mixer + initial_state : np.ndarray + The initial state for QAOA, default is the uniform superposition state (corresponding to the X mixer) + n_trotters : int + Number of Trotter steps in each mixer layer for the xy mixer Returns ------- @@ -141,7 +166,7 @@ def get_qaoa_objective( # -- Qiskit edge case if simulator == "qiskit": - g = _get_qiskit_objective(parameterized_circuit, precomputed_objectives, precomputed_optimal_bitstrings, objective, terms, parameterization) + g = _get_qiskit_objective(parameterized_circuit, precomputed_objectives, precomputed_optimal_bitstrings, objective, terms, parameterization, mixer) def fq(*args): gamma, beta = qokit.parameter_utils.convert_to_gamma_beta(*args, parameterization=parameterization) @@ -150,7 +175,13 @@ def fq(*args): return fq # -- - simulator_cls = choose_simulator(name=simulator) + if mixer == "x": + simulator_cls = choose_simulator(name=simulator) + elif mixer == "xy": + simulator_cls = choose_simulator_xyring(name=simulator) + else: + raise ValueError(f"Unknown mixer type passed to get_qaoa_objective: {mixer}, allowed ['x', 'xy']") + sim = simulator_cls(N, terms=terms, costs=precomputed_diagonal_hamiltonian) if precomputed_objectives is None: precomputed_objectives = sim.get_cost_diagonal() @@ -162,7 +193,7 @@ def fq(*args): # -- Final function def f(*args): gamma, beta = qokit.parameter_utils.convert_to_gamma_beta(*args, parameterization=parameterization) - result = sim.simulate_qaoa(gamma, beta) + result = sim.simulate_qaoa(gamma, beta, initial_state, n_trotters=n_trotters) if objective == "expectation": return sim.get_expectation(result, costs=precomputed_objectives, preserve_state=False) elif objective == "overlap": @@ -217,83 +248,3 @@ def get_qaoa_labs_overlap_fourier(N: int, p: int, parameterization: str = "freq" """ return get_qaoa_labs_objective(N, p, objective="overlap", parameterization=parameterization, **kwargs) - - -def get_qaoa_xy_objective( - N: int, - p: int | None = None, - sv0: np.ndarray | None = None, - mixer: str | None = None, - T: int | None = None, - precomputed_diagonal_hamiltonian=None, - precomputed_objectives=None, - precomputed_optimal_bitstrings=None, - parameterized_circuit=None, - parameterization: str = "theta", - objective: str = "expectation", - simulator: str = "auto", -): - """Return QAOA objective to be minimized - - Parameters - ---------- - N : int - Number of qubits - p : int - Number of QAOA layers (number of parameters will be 2*p) - sv0: - Initial state of the QAOA, choice = ["dicke"] - mixer: str - Chosen Mixer, choice: ["trotter_ring"] - parameterization : str - If parameterization == 'theta', then f takes one parameter (gamma and beta concatenated) - If parameterization == 'gamma beta', then f takes two parameters (gamma and beta) - For below Fourier parameters, q=p - If parameterization == 'freq', then f takes one parameter (fourier parameters u and v concatenated) - If parameterization == 'u v', then f takes two parameters (fourier parameters u and v) - objective : str - If objective == 'expectation', then returns f(theta) = - < theta | C_{PO} | theta > (minus for minimization) - If simulator == 'auto', implementation is chosen automatically - (either the fastest CPU simulator or a GPU simulator if CUDA is available) - If simulator == 'qiskit', implementation in qaoa_qiskit is used - - Returns - ------- - f : callable - Function returning the negative of expected value of QAOA with parameters theta - """ - # Prepare function that computes objective from precomputed energies / bitstrings - if objective == "expectation": - if precomputed_objectives is None: - raise ValueError(f"precomputed_objectives are required when using the {objective} objective") - - def compute_objective_from_probabilities(probabilities): - return precomputed_objectives.dot(probabilities) - - else: - raise ValueError(f"Unknown objective passed to get_qaoa_xy_objective: {objective}, allowed ['expectation']") - - if simulator == "qiskit": - # Prepare parameterized circuit - backend = Aer.get_backend("statevector_simulator") - - def probabilities_from_gamma_beta(gamma, beta): - qc = parameterized_circuit.bind_parameters(np.hstack([np.asarray(beta) / 2, np.asarray(gamma) / 2])) - result = execute(qc, backend).result() - sv = reverse_array_index_bit_order(result.get_statevector()) - return np.abs(sv) ** 2 - - else: - if simulator == "auto": - sim_class = choose_simulator_xyring() - else: - sim_class = choose_simulator_xyring(name=simulator) - sim = sim_class(N, precomputed_diagonal_hamiltonian) - - def probabilities_from_gamma_beta(gamma, beta): - res = sim.simulate_qaoa(gamma, beta, sv0=sv0, n_trotters=T) - return sim.get_probabilities(res) - - f = set_parameterized_func(parameterization, probabilities_from_gamma_beta, compute_objective_from_probabilities) - - return f diff --git a/qokit/qaoa_objective_portfolio.py b/qokit/qaoa_objective_portfolio.py index bd3993166..e6c18b6ae 100644 --- a/qokit/qaoa_objective_portfolio.py +++ b/qokit/qaoa_objective_portfolio.py @@ -1,13 +1,9 @@ -############################################################################### -# // SPDX-License-Identifier: Apache-2.0 -# // Copyright : JP Morgan Chase & Co -############################################################################### from __future__ import annotations import numpy as np from .utils import precompute_energies, reverse_array_index_bit_order -from .portfolio_optimization import get_configuration_cost_kw, po_obj_func +from .portfolio_optimization import get_configuration_cost_kw, po_obj_func, portfolio_brute_force from qokit.qaoa_circuit_portfolio import generate_dicke_state_fast, get_parameterized_qaoa_circuit -from .qaoa_objective import get_qaoa_xy_objective +from .qaoa_objective import get_qaoa_objective def get_qaoa_portfolio_objective( @@ -71,14 +67,16 @@ def get_qaoa_portfolio_objective( if mixer == "trotter_ring": pass else: - raise ValueError(f"Unknown mixer passed to get_qaoa_portfolio_objective: {ini}, allowed ['trotter_ring']") + raise ValueError(f"Unknown mixer passed to get_qaoa_portfolio_objective: {mixer}, allowed ['trotter_ring']") - return get_qaoa_xy_objective( + if objective in ["overlap", "expectation and overlap"] and precomputed_optimal_bitstrings is None: + bf_result = portfolio_brute_force(po_problem, return_bitstring=True) + precomputed_optimal_bitstrings = bf_result[1].reshape(1, -1) + assert precomputed_optimal_bitstrings.shape[1] == N # only one optimal bitstring + + return get_qaoa_objective( N=N, p=p, - sv0=sv0, - mixer=mixer, - T=T, precomputed_diagonal_hamiltonian=po_problem["scale"] * precomputed_energies, precomputed_objectives=precomputed_energies, precomputed_optimal_bitstrings=precomputed_optimal_bitstrings, @@ -86,4 +84,7 @@ def get_qaoa_portfolio_objective( parameterization=parameterization, objective=objective, simulator=simulator, + mixer="xy", + initial_state=sv0, + n_trotters=T, ) diff --git a/tests/test_portfolio_optimization.py b/tests/test_portfolio_optimization.py index 82658c080..93c2eb14e 100644 --- a/tests/test_portfolio_optimization.py +++ b/tests/test_portfolio_optimization.py @@ -1,8 +1,4 @@ -############################################################################### -# // SPDX-License-Identifier: Apache-2.0 -# // Copyright : JP Morgan Chase & Co -############################################################################### -from qokit.portfolio_optimization import get_problem, get_problem_H, get_problem_H_bf, po_obj_func +from qokit.portfolio_optimization import get_problem, get_problem_H, get_problem_H_bf, po_obj_func, portfolio_brute_force, get_sk_ini from qokit.utils import precompute_energies, reverse_array_index_bit_order from qokit.qaoa_circuit_portfolio import ( circuit_measurement_function, @@ -18,6 +14,7 @@ import numpy as np import pytest from qiskit import execute, Aer +from functools import reduce simulators_to_run = get_available_simulator_names("xyring") @@ -131,3 +128,28 @@ def test_get_problem_H(): H_problem = get_problem_H(po_problem) H_problem_bf = get_problem_H_bf(po_problem) assert np.allclose(H_problem, H_problem_bf) + + +def test_portfolio_AR(): + po_problem = get_problem(N=6, K=3, q=0.5, seed=1, pre="rule") + p = 1 + qaoa_obj = get_qaoa_portfolio_objective(po_problem=po_problem, p=p, ini="dicke", mixer="trotter_ring", T=1, simulator="python") + best_portfolio = portfolio_brute_force(po_problem, return_bitstring=False) + x0 = get_sk_ini(p=p) + po_energy = qaoa_obj(x0).real + po_ar = (po_energy - best_portfolio[1]) / (best_portfolio[0] - best_portfolio[1]) + # a problem with known AR > 0.7564 + assert po_ar > 0.75 + + +def test_best_bitstring(): + N = 6 + po_problem = get_problem(N=N, K=3, q=0.5, seed=1, pre="rule") + bf_result = portfolio_brute_force(po_problem, return_bitstring=True) + precomputed_optimal_bitstrings = bf_result[1].reshape(1, -1) + assert precomputed_optimal_bitstrings.shape[1] == N + bitstring_loc = np.array([reduce(lambda a, b: 2 * a + b, x) for x in precomputed_optimal_bitstrings]) + assert len(bitstring_loc) == 1 + po_obj = po_obj_func(po_problem) + precomputed_energies = reverse_array_index_bit_order(precompute_energies(po_obj, N)).real + assert np.isclose(precomputed_energies[bitstring_loc[0]], bf_result[0])