Converting sigma model Hamiltonian into a circuit

PennyLane Help
1

Hi all, I would like to ask some specific questions (I’m sorry its very specific, forgive me) about a paper I’m reading. In 1903.06577, authors are playing with the sigma model hamiltonian, and they are showing how they constructed the time evolution circuit with respect to this hamiltonian. They are using four-dimensional Hilbert space (per site), so each site is defined via two-qubit. I was wondering how I can compute the expectation value of this hamiltonian on a quantum circuit instead of time evolution.

So for the interaction terms (please correct me if I’m wrong), they are giving the j_{1,2,3} definition (eq. 8) so I believe it is correct to use that definition to form the hamiltonian i.e. j_1 = \frac{\mathbf{1}\otimes\sigma_2}{\sqrt{3}} = qml.Hamiltonian([1/np.sqrt(3.)], [qml.Identity(0) @ qml.PauliY(1)]). However, I’m not sure about the kinetic term, they are converting everything into this new T basis they are using but without the conversion, (e.g. 1812.00944 eq 2), the kinetic term for the sigma model is \frac{1}{2\beta}\sum_k^3 J^2_k where J_k is the angular momentum operator. So would it be correct to write it as qml.Hamiltonian([(1/2*betta)*(hbar/2)**2]*3, [qml.PauliX(0) @ qml.PauliX(1), qml.PauliY(0) @ qml.PauliY(1), qml.PauliZ(0) @ qml.PauliZ(1)])? If I write the entire sum I believe this will give me an identity matrix with a phase. So if I write a simple circuit to calculate the expectation value for only the terms k=0 in kinetic term (for simplicity) and first interaction term

import pennylane as qml

dev = qml.device("default.qubit", wires = 2)

@qml.qnode(dev)
def circuit(): 
    # prepare an initial state (not important)
    for wire in range(2):
        qml.RY(np.pi/2., wire)

    return qml.expval(qml.Hamiltonian([(1/2)**2, 1/np.sqrt(3.)], 
                                      [qml.PauliX(0) @ qml.PauliX(1), qml.Identity(0) @ qml.PauliY(1)]))

Would this be correct? Any insight and/or reference is highly appreciated.

Thanks for your time

2

Hi @jackaraz, thank you for asking this question here. We are taking a look at it and should be back soon with an answer.

3

Hi @CatalinaAlbornoz, thank you so very much!!

In the meantime, I was looking into other sources and I believe for the interaction terms I can use raising and lowering operators i.e. S^{\pm} = \sigma_x \pm i\sigma_y and S_z so the full summation should be \sum_n \left( S^+_nS^-_{n+1} + S^-_nS^+_{n+1} + S^z_n S^z_{n+1} \right) where n indicates qubit here. So the first interaction term in the Hamiltonian, I believe, can be written as qml.Hamiltonian([1, -1j, 1j, 1], [X(0)@X(1), X(0)@Y(1), Y(0)@X(1), Y(0)@Y(1)]) where X = lambda wire: qml.PauliX(wire) and the same for Y. Additionally, since J^2 should give identity, I believe it can be written as [qml.Identity(n)@qml.Identity(n+1) for n in range(0, Nqubit, 2)] inside the hamiltonian with the proper constant in front. Tho, I’m not sure if I can use qubit information as such. For instance, if I use 6D Hilbert space instead of 4, should I just extend the qubit structure, i.e. X(0)@X(1)@X(2), or should I change the entire formulation of my hamiltonian to have my raising and lowering operators in different basis?

4

Hi @jackaraz . Thanks for the question. This seems really interesting.

As far as I can tell (though I may be misinterpreting the notation), the kinetic term is projectors:

\mathcal{H}^{0} = g^2 \left( | T_1\rangle\langle T_1 | + | T_2 \rangle\langle T_2| + |T_3\rangle\langle T_3 | \right)

We can see from this we recover the relationship in the paper:

h_{1,1}^{0} = g^2 \langle T_1| T_1 \rangle \langle T_1 | T_1 \rangle = g^2

The question then is: What does this look like in terms of qubits on a quantum computer?

If we make the mapping:

|T_0 \rangle \rightarrow |00\rangle
| T_1 \rangle \rightarrow | 0 1\rangle
| T_2 \rangle \rightarrow |1 0 \rangle
| T_3 \rangle \rightarrow | 11 \rangle

We could express the kinetic hamiltonian as a sum of the projectors over those states:

proj1 = qml.Projector([0,1], wires=(0,1))
proj2 = qml.Projector([1,0], wires=(0,1))
proj3 = qml.Projector([1,1], wires=(0,1))

h = qml.Hamiltonian([g**2, g**2, g**2], [proj1, proj2, proj3])

This seems to be consistent with eq. 11 in the paper.

The expectation value of this Hamiltonian would be equivalent to measuring the probabilities and then summing over all but the first term:

@qml.qnode(dev)
def circ():
    return qml.probs(wires=(0,1))

g**2 * sum(circ()[1:])

I’m still trying to figure out the interaction hamiltonian component. Hope this helps a little.

1 Like
5

Hi @christina, thanks a lot for your answer! I have couple questions about your answer; I’m assuming in your circ function you assumed an initial state otherthan |00\rangle right? i.e.

@qml.qnode(dev)
def circ():
    qml.RY(np.random.uniform(0, np.pi, 1)[0], 0)
    qml.RY(np.random.uniform(0, np.pi, 1)[0], 1)
    return qml.probs(wires=(0,1))

Also I believe 2 qubit operators are not supported yet so the expectation value of the Hamiltonian that you wrote raises a following error;

ValueError: Can only compute sparse representation for tensors whose operations act on consecutive wires; got Projector([0, 1], wires=[0, 1]).

And I can see that in the code this feature listed as todo:

File ~/pennylane/operation.py:2017, in Tensor.sparse_matrix(self, wires)
   2014 for o in self.obs:
   2015     if len(o.wires) > 1:
   2016         # todo: deal with multi-qubit operations that do not act on consecutive qubits
-> 2017         raise ValueError(
   2018             f"Can only compute sparse representation for tensors whose operations "
   2019             f"act on consecutive wires; got {o}."
   2020         )
   2021     # store the single-qubit ops according to the order of their wires
   2022     idx = wires.index(o.wires)

So I believe a workaround for now would be

proj1 = qml.Projector([0], wires=(0)) @ qml.Projector([1], wires=(1))
proj2 = qml.Projector([1], wires=(0)) @ qml.Projector([0], wires=(1))
proj3 = qml.Projector([1], wires=(0)) @ qml.Projector([1], wires=(1))

h = qml.Hamiltonian([g**2, g**2, g**2], [proj1, proj2, proj3])

@qml.qnode(dev)
def circ():
    qml.RY(np.random.uniform(-np.pi, np.pi, 1)[0], 0)
    qml.RY(np.random.uniform(-np.pi, np.pi, 1)[0], 1)
    return qml.expval(h)

Thanks a lot this was very helpful!

6

Glad you found a workaround.

And yes, in my example I was just trying to show the expval part and left out state prep.