Digitized Counterdiabatic Quantum Optimization#

class DCQOProblem(Q, H_init, H_prob, A_lam, agp_coeffs, lam_func, H_control=None, qarg_prep=None)[source]#

General structure to formulate Digitized Counterdiabatic Quantum Optimization problems.

This class is used to solve DCQO problems with the algorithms COLD (counterdiabatic optimized local driving) or LCD (local counterdiabatic driving). To run the COLD algorithm on the problem, you need to specify the control Hamiltonian H_control and the inverse scheduling function g_func. These are not needed for the LCD algorithm. To learn more about counterdiabatic driving, make sure to check out the tutorial.

Parameters:

Qnp.array: The QUBO matrix. Must be symmetric: the cost is read as \(y = x^T Q x\) with off-diagonal entries counted twice, and the Ising coefficients \(J_{ij} = \frac{1}{2} Q_{ij}\), \(h_i = -\frac{1}{2} \sum_j Q_{ij}\) derived from it are only valid for symmetric \(Q\). Symmetrize an upper-triangular matrix with Q = 0.5 * (Q + Q.T) before passing it in.
H_initQubitOperator: Hamiltonian, the system is at the time t=0.
H_probQubitOperator: Hamiltonian, the system evolves to for t=T.
A_lamQubitOperator: Operator holding an appoximation for the adiabatic gauge potential (AGP).
agp_coeffscallable: The parameters for the adiabatic gauge potential (AGP). If the COLD method is being used, they must depend on the optimization pulses in H_control.
lam_funccallable: A function \(\lambda(t, T)\) mapping \(t \in [0, T]\) to \(\lambda \in [0, 1]\). This function needs to return a sympy expression with \(t\) and \(T\) as sympy.Symbols.
H_controlQubitOperator, optional: Hamiltonian specifying the control pulses for the COLD method. If not given, the LCD method is used automatically.
qarg_prepcallable, optional: A function receiving a QuantumVariable for preparing the inital state. By default, the groundstate of the x-operator \(\ket{-}^n\) is prepared.

Examples

For a quick demonstration we build a DCQO problem instance for a 4x4 QUBO. We choose a first order AGP ansatz with uniform coefficients and solve it with LCD.

import numpy as np
import sympy as sp
from qrisp.operators.qubit import QubitOperator, X, Y, Z
from qrisp.algorithms.cold import DCQOProblem
from qrisp import QuantumVariable

Q = np.array([
    [-1.2,  0.40, 0.0,  0.0],
    [ 0.40,  0.30, 0.20, 0.0],
    [ 0.0,   0.20,-1.1,  0.30],
    [ 0.0,   0.0,  0.30,-0.80]
])
N = Q.shape[0]

# Define QUBO problem hamiltonian
h = -0.5 * np.sum(Q, axis=1)
J = 0.5 * Q

H_init = 1 * sum([X(i) for i in range(N)])

# QubitOperator.sum merges all terms in a single pass -- for a dense
# QUBO (O(N^2) coupling terms), Python's built-in sum() would instead
# fold pairwise via repeated __add__ calls, each copying the full
# running total, costing O(N^4) instead of O(N^2).
H_prob = QubitOperator.sum(
    J[i][j] * Z(i) * Z(j) for i in range(N) for j in range(i)
) + QubitOperator.sum(h[i] * Z(i) for i in range(N))

# Create AGP
A_lam = sum([Y(i) for i in range(N)]) # uniform

# Function for uniform AGP coefficients. This minimizes the action
# S = Tr[G^2], G = d_lambda H + i[A_lam, H], for the ansatz A_lam = sum_i Y_i.
def alpha(lam):
    J_sq = np.sum(J**2) - np.sum(np.diag(J)**2)   # sum_{i != j} J_ij**2
    nom = np.sum(h)
    denom = 2 * (N * (1 - lam)**2 + np.sum((lam * h)**2) + lam**2 * J_sq)
    return [-nom/denom]*N

# Simple scheduling function 0 -> 1
def lam():
    t, T = sp.symbols("t T", real=True)
    lam_expr = t/T
    return lam_expr

# Create problem instance
lcd_problem = DCQOProblem(Q, H_init, H_prob, A_lam, alpha, lam)

# Run problem with LCD algorithm
qarg = QuantumVariable(N)
res = lcd_problem.run(qarg, N_steps=4, T=12, method="LCD")
print(res)

{
    '1011': [0.5666, np.float64(-2.5)],
    '1111': [0.0986, np.float64(-1.0)],
    '1001': [0.087, np.float64(-2.0)],
    '1010': [0.0846, np.float64(-2.3)],
    '0111': [0.0586, np.float64(-0.6000000000000001)],
    '1000': [0.031, np.float64(-1.2)],
    '1110': [0.0218, np.float64(-0.8000000000000003)],
    '0011': [0.0216, np.float64(-1.3)],
    '0110': [0.0074, np.float64(-0.40000000000000013)],
    '0001': [0.0056, np.float64(-0.8)],
    '1101': [0.0052, np.float64(-0.9)],
    '0101': [0.005, np.float64(-0.5)],
    '0100': [0.0026, np.float64(0.3)],
    '0000': [0.0022, np.float64(0.0)],
    '0010': [0.0018, np.float64(-1.1)],
    '1100': [0.0004, np.float64(-0.09999999999999998)]
}

We get a dictionary where the key is the quantum state and the values are lists of [probability, cost]. So our most likely result is ‘1011’ with probability around 0.57 and the QUBO cost \(x^T Q x = -2.5\).

Methods#

`DCQOProblem.run`(qarg, N_steps, T, method[, ...])	Run the DCQO problem instance with given quantum arguments, number of timesteps, evolution time and method.
`DCQOProblem.apply_lcd_hamiltonian`(qarg, ...)	Simulate the local counterdiabatic driving (LCD) Hamiltonian on a quantum argument via trotterization.
`DCQOProblem.apply_cold_hamiltonian`(qarg, ...)	Simulate counterdiabatic optimized local driving (COLD) Hamiltonian on a quantumvariable via trotterization.
`DCQOProblem.compile_U_cold`(qarg, N_opt, ...)	Compile the circuit created by the COLD Hamiltonian application step.
`DCQOProblem.optimization_routine`(qarg, ...)	Subroutine for the optimization method used in COLD.
`DCQOProblem.QUBO_cost`(res)	Returns the cost y = x^T Q x for a given binary array x.

Quick run#

Quickly run LCD/COLD with our predefined operators and functions in QUBO.py:

solve_QUBO(Q: array, problem_args: dict, run_args: dict)[source]#

Solves a QUBO Matrix using counterdiabatic driving.

This method uses the pre-defined COLD/LCD operators (hamiltonian, scheduling function, AGP parameters) as described in the tutorial. To define your own operators, create a DCQO instance and use the run method.

Parameters:

Qnp.array

QUBO Matrix to solve. Must be symmetric: the cost is read as \(y = x^T Q x\) with off-diagonal entries counted twice, and the Ising coefficients \(J_{ij} = \frac{1}{2} Q_{ij}\), \(h_i = -\frac{1}{2} \sum_j Q_{ij}\) derived from it are only valid for symmetric \(Q\). Symmetrize an upper-triangular matrix with Q = 0.5 * (Q + Q.T) before passing it in.

problem_argsdict

Holds arguments for DCQO problem creation:

method : str – “COLD” or “LCD”.
uniform : bool.
agp_type : str, optional – “local” (default) or “nc”. Applies to both methods. nc approximates the AGP with first-order nested commutators and is the better approximation, especially at short evolution times. For COLD, nc requires "uniform": True.

run_argsdict

Holds arguments for running the DCQO instance (N_steps, T, N_opt, CRAB, objective, precision, backend, exp_value_backend). All options are also listed here: DCQOProblem.run().

Returns:

resultdict: The dictionary holding the QUBO vector results, with their probabilitites and cost. They are ordered from most to least likely and the dictionary entries are {“state”: [prob, cost]}.

Examples

import numpy as np
from qrisp.algorithms.cold import solve_QUBO

Q = np.array([[-0.6,  0.2, -0.5, -0.4, -0.6,  0.0],
            [ 0.2, -1.0,  0.0,  0.0,  0.0,  0.0],
            [-0.5,  0.0, -1.2,  0.5, -0.5,  0.0],
            [-0.4,  0.0,  0.5, -0.8,  0.0,  0.1],
            [-0.6,  0.0, -0.5,  0.0, -1.2,  0.0],
            [ 0.0,  0.0,  0.0,  0.1,  0.0,  0.3]])

problem_args = {"method": "COLD", "uniform": False}

run_args = {"N_steps": 20, "T": 1, "N_opt": 1, "CRAB": False, "bounds": (-3, 3)}

result = solve_QUBO(Q, problem_args, run_args)

print(result)

{
    '111110': [0.9816, np.float64(-7.3999999999999995)],
    '111111': [0.0058, np.float64(-6.8999999999999995)],
    '111100': [0.0048, np.float64(-4.0)],
    '011110': [0.0032, np.float64(-4.2)],
    '110110': [0.0022, np.float64(-5.199999999999999)],
    '111010': [0.0012, np.float64(-6.8)],
    '001110': [0.0004, np.float64(-3.2)],
    '011100': [0.0002, np.float64(-2.0)],
    '110010': [0.0002, np.float64(-3.5999999999999996)],
    '011010': [0.0002, np.float64(-4.4)],
    '101110': [0.0002, np.float64(-6.8)]
}