Semidefinite programming

Toolkit to solve the semi-definite program defined by the optimization problem and its set of constraints.

---

source

picos2np

 picos2np (variable)

Converts picos variable (even sparse) to numpy matrix.

SdP formalism

In order to formulate and solve the SdP there are two main needed items: * an object containing the problem details in terms of picos.Constant with a to_sdp() method, such as a Hamiltonian with Hamiltonian.to_sdp(). * a layout based in np.ndarrays, e.g., L = [np.array([0, 1]), np.array([0, 1, 2]) determining the problem constraints.

---

source

SdPSolver

 SdPSolver (solver='cvxopt')

Initialize self. See help(type(self)) for accurate signature.

The SdPSolver must be tailored to the problem instances that we aim to solve. Defining the problem and the solver should go hand-to-hand and they are handled by the SdPEnvironment during the optimization. Essentially, the environment takes the result of problem.to_sdp() and feeds it into the solver.

SdPSolvers have two main methods: solve and ojimetro. The first one, solves a problem instance given a layout, such as solver.solve(problem.to_sdp(), layout). The second one estimates the cost associated to solving the SdP given a layout, with the syntax solver.ojimetro(layout). We can choose any metric for the cost, for example, the amount of free variables in the resulting SdP.

Ground state energy approximation

Finding the ground state energy of a quantum Hamiltonian is a fundamental problem in quantum physics. It constitutes global optimization problem over all possible quantum states \(\rho\) of the form \[E_0=\min_\rho \text{Tr}\left[\rho H\right]\,,\] where \(E_0\) denotes the ground state energy for the Hamiltonian \(H\). Nonetheless, the complexity of such optimization task grows exponentially with the physical system’s size, so we often rely on approximations. Here, we consider a relaxation of the problem based on an SdP.

---

source

SdPEnergySolver

 SdPEnergySolver (solver='cvxopt')

Solver to compute lower bounds of the ground state energy of local Hamiltonians. It follows the method described in [1] (see refs below) to formulate the objective and enforce the problem constraints.

The SdPEnergySolver solves the problem of finding lower bounds to the ground state energy of local Hamiltonians. Hence, the SdPEnvironment should receive a Hamiltonian as problem instance together with this solver. We use this case to illustrate the general behaviour of the SdPSolvers and how to use them.

N = 6
chain = Chain1D(N)
B, J = 1., [i%3 for i in range(N)]
H = XXHamiltonian(chain, B, J)

Solving a problem instance

Now that we have defined a problem, we use the solve method to define the objective function, which corresponds to the expected energy, and the problem constraints, given by a layout.

---

source

SdPEnergySolver.solve

 SdPEnergySolver.solve (interactions, layout)

Creates and solves the SdP associated to the given layout and Hamiltonian interactions. The result is a lower bound of the Hamiltonian’s ground state energy.

The method takes the Hamiltonian.interactions, provided by the Hamiltonian.to_sdp() method, and a layout with which we define the set of compatibility constraints imposed in the SdP, as we explain in [1].

Given the Hamiltonian, let’s define a layout in order to compute an energy bound.

simple_layout = fill_layout([np.array([0, N-1])], N)
simple_layout

[array([0, 5]), array([1]), array([2]), array([3]), array([4])]

And now we can solve the associated SdP to the Hamiltonian with the given set of constraints.

solver = SdPEnergySolver()
solver.solve(H.to_sdp(), simple_layout)

-15.999999999012134

Tightening the constraints, we can obtain a better energy bound. For instance, adding a 3-body element to the previous set we will obtain a higher energy.

stronger_layout = fill_layout([np.array([1, 2, 3]), np.array([0, N-1])], N)
stronger_layout

[array([1, 2, 3]), array([0, 5]), array([4])]

solver.solve(H.to_sdp(), stronger_layout)

-12.472135949820899

Cost estimation

A key element of any optimization problem is its computational cost. We can take any consistent measure of the problem’s cost as a valid estimation for the posterior optimization with the reinforcement learning framework. For example, since the actual complexity of solving an SdP depends on the algorithm that we use, we take the number of free parameters as a metric for the associated cost.

In order to get a cost estimation, we rely on the ojimetro method.

---

source

SdPEnergySolver.ojimetro

 SdPEnergySolver.ojimetro (layout)

Estimates the amount of free parameters involved in the SdP associated to the layout.

Following our example, the first SdP was a wild relaxation of the original problem.

solver.ojimetro(simple_layout)

31

However, the second one was a tighter one.

solver.ojimetro(stronger_layout)

83

With the first set of constraints, the resulting SdP had 31 free variables to optimize, while the second SdP had to deal with 83. Tighter energy bounds usually come at the cost of higher computational costs.

assert solver.ojimetro(simple_layout) < solver.ojimetro(stronger_layout)

Alternative methods to lower bound the ground state energy

We implement two other methods to obtain lower bounds of many-body Hamiltonians. We rewrite the Anderson bound [2] and the method introduced by Uskov and Lichovskiy in [3] in SdP form to compare with the proposed methodology above. The SdP formulation ensures that the obtained result is the actual global minimum.

The main limitation of these methods is not taking into account compatibility constraints between the reduced density matrices spanning the system. Furthermore, they rely on symmetries (such as translational invariance) and can, thus, not be used to solve inhomogeneous systems. Conversely, the approach we introduced above can be nautrally complemented by introducing additional constratints stemming from any previously known symmetry from the system.

---

source

SdPEnergyAndersonSolver

 SdPEnergyAndersonSolver (solver='cvxopt')

Solver to compute bounds to the ground state energy of local Hamiltonians implementing the methods described in [2] (see references below). Finds the so-called Anderson bound. Assumes the system is one-dimensional.

---

source

SdPEnergyUskovLichkovskiySolver

 SdPEnergyUskovLichkovskiySolver (solver='cvxopt')

Solver to compute ground state energy bounds of local Hamiltonians. It follows the method introduced in [3] (see references below) to improve over the Anderson bound. The method assumes the reduced states have certain symmetries, such as translational invariance.

We can quickly reproduce the results from [3] (see Table 2) by properly adjusting the cluster size and using the same Hamiltonian.

N = 8
chain = Chain1D(N)
B, J = 0, 1
H = XYZHamiltonian(chain, B, J)

solver_anderson = SdPEnergyAndersonSolver()
solver_uskov = SdPEnergyUskovLichkovskiySolver()

cluster_size = 5
layout = [np.sort(np.arange(i, i+cluster_size)%N) for i in np.arange(0, N, cluster_size-3)]

solver_anderson.solve(H.to_sdp(), cluster_size)/N, solver_uskov.solve(H.to_sdp(), cluster_size)/N, solver.solve(H.to_sdp(), layout)/N

(-1.9278862525095595, -1.8685170845449277, -1.868517086761974)

Entanglement witnessing

We take an energy-based approach to build entanglement witnesses here. The main principle is, given a Hamiltonian, we compute the minimum possible energy obtainable within the set of separable states. Then, if we find a state that yields a lower energy, it means that it is entangled [4].

Formally, \[\Delta E = \langle H\rangle - E_{\text{sep}}\] denotes the energy gap between the expected energy of a state, \(\langle H\rangle = \text{Tr}\left[\rho H\right]\), and the minimum energy obtainable by separable states, \(E_{\text{sep}}\). If \(\Delta E < 0\), \(\rho\) is an entangled state, as it is beyond the set of separable states. The larger the difference, \(|\Delta E|\), the more robust the entanglement. Actually, \(|\Delta E|\) quantifies the energy the system must receive to become separable.

---

source

SdPWitnessSolver

 SdPWitnessSolver (solver='cvxopt')

Solver to detect entanglement from energy bounds. It follows the formalism described in [1] to obtain a lower bound of the minimum energy for separable states. If we can find a state with a lower energy than that, the state is entangled.

The SdPWitnessSolver provides a lower bound for \(E_{\text{sep}}\) given a Hamiltonian. We build upon the SdPEnergySolver imposing additional constraints to enforce separability on our states. However, imposing such conditions to their full extent is extremely costly. Instead, we impose that our states are positive under partial transposition (PPT), which is a relaxation that is exact in the one-dimensional limit. Hence, the bound on \(E_{\text{sep}}\).

Let’s see an example! First of all, we need to build our problem.

N = 6
chain = Chain1D(N)
B, J = 1, 2
H = XXHamiltonian(chain, B, J)

Now we can instantiate our solver.

solver_witness = SdPWitnessSolver()

We define a layout of 2-body terms.

layout = fill_layout([np.sort([i, (i+1)%N]) for i in range(N)], N)
layout

[array([0, 1]),
 array([1, 2]),
 array([2, 3]),
 array([3, 4]),
 array([4, 5]),
 array([0, 5])]

And now we’re ready to solve it.

solver_witness.solve(H.to_sdp(), layout)

-12.749999991956518

In [4] we find analytic expressions of \(E_{\text{sep}}\) for some Hamiltonians, including the XX model we consider here.

def exact_XX_sep_bound(N, B, J, d):
    b = B/J/d
    if b <= 2:
        return -d*N*J*(1+b**2/4)
    return -d*N*J*b

exact_XX_sep_bound(N, B, J, 1)

-12.75

We see that, since we’ve chosen to implement the Hamiltonian in a one-dimensional chain, the result we obtained is exact. Now we can use the SdPEnergySolver to obtain the ground state energy and compare it to the separable threshold.

solver_energy = SdPEnergySolver()
solver_energy.solve(H.to_sdp(), layout)

-23.99999993375356

The ground state is entangled! :)

If we now change the Hamiltonian parameters, we can land on a separable ground state.

B, J = 5, 1
H = XXHamiltonian(chain, B, J)

solver_witness.solve(H.to_sdp(), layout), solver_energy.solve(H.to_sdp(), layout)

(-29.9999999937129, -29.999999800727085)

In this case, both energies match. This confirms that the ground state is, indeed, a separable state.

References

[1] B. Requena, G. Muñoz-Gil, M. Lewenstein, V. Dunjko, J. Tura. Certificates of quantum many-body properties assisted by machine learning. Phys. Rev. Research 5, 013097 (2023)

[2] P. W. Anderson. Limits on the Energy of the Antiferromagnetic Ground State. Physical Review 83, 1260 (1951)

[3] F. Uskov and O. Lychkovskiy. A variational lower bound on the ground state of a many-body system and the squaring parametrization of density matrices. Journal of Physics: Conference Series 1163 012057 (2019)

[4] G. Tóth. Entanglement Witnesses in Spin Models. Physical Review A 71 010301 (2005).