I have followed this tutorial and based on it, I've written the following function in qiskit, which can explicitly calculate the ground states of a transverse-field Ising Hamiltonian.
from qiskit import * import numpy as np def Hamiltonian(n,h): pow_n=2**n qc = np.empty(2*n-1, dtype=object) #Creating the quantum circuits that are used in the calculation of the Hamiltonian based on the number of qubits for i in range(0, 2*n-1): #2n-1 is the number of factors on the n-site Hamiltonian qr = QuantumRegister(n) qc[i] = QuantumCircuit(qr) #create quantum circuits for each factor of the Hamiltonian #print(i) if (i<=n-2): #for the first sum of the Hamiltonian qc[i].z(i) #value of current spin qc[i].z(i+1) #and value of its neighboring spin else: #for the second sum of the Hamiltonian qc[i].x(2*n-2-i) #2*n-2 gives the proper index since counting starts at 0 #Run each circuit in the simulator simulator = Aer.get_backend('unitary_simulator') result = np.empty(2*n-1, dtype=object) unitary = np.empty(2*n-1, dtype=object) Hamiltonian_Matrix=0 #Get the results for each circuit in unitary form for i in range(0, 2*n-1): result[i] = execute(qc[i], backend=simulator).result() unitary[i] = result[i].get_unitary() #print(unitary[i]) #And calculate the Hamiltonian matrix according to the formula if (i<=n-2): Hamiltonian_Matrix=np.add(Hamiltonian_Matrix,-unitary[i]) else: Hamiltonian_Matrix=np.add(Hamiltonian_Matrix,-h*unitary[i]) print("The",pow_n,"x",pow_n, "Hamiltonian Matrix is:") print(Hamiltonian_Matrix) #Now that we have the Hamiltonian #find the eigenvalues and eigenvectors w, v = np.linalg.eig(Hamiltonian_Matrix) print("Eigenvectors") print(v) print("Eigenvalues") print(w) minimum=w min_spot=0 for i in range(1, pow_n): if w[i]<minimum: min_spot=i minimum=w[i] print(min_spot) groundstate = v[:,min_spot] #the probability to measure each basic state of n qubits probability = np.square(groundstate).real print("The probability for each of the",pow_n,"base states is:") print(probability) print("The probabilities for each of the",pow_n,"base states add up to:") print ("%.2f" % np.sum(probability))
My problem with this piece of code I've written is that it can only run on a unitary simulator. To my understanding (which may lack some of the underlying physics), the Hamiltonian itself is not a "purely" quantum calculation, since there are additions to be made which cannot be expressed with a quantum (unitary) gate, and this is why the resulting Hamiltonian matrix is also not unitary. For example, if you run Hamiltonian(3, 1), the Hamiltonian matrix is:
[[-2.+0.j -1.+0.j -1.+0.j 0.+0.j -1.+0.j 0.+0.j 0.+0.j 0.+0.j] [-1.+0.j 0.+0.j 0.+0.j -1.+0.j 0.+0.j -1.+0.j 0.+0.j 0.+0.j] [-1.+0.j 0.+0.j 2.+0.j -1.+0.j 0.+0.j 0.+0.j -1.+0.j 0.+0.j] [ 0.+0.j -1.+0.j -1.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j -1.+0.j] [-1.+0.j 0.+0.j 0.+0.j 0.+0.j 0.+0.j -1.+0.j -1.+0.j 0.+0.j] [ 0.+0.j -1.+0.j 0.+0.j 0.+0.j -1.+0.j 2.+0.j 0.+0.j -1.+0.j] [ 0.+0.j 0.+0.j -1.+0.j 0.+0.j -1.+0.j 0.+0.j 0.+0.j -1.+0.j] [ 0.+0.j 0.+0.j 0.+0.j -1.+0.j 0.+0.j -1.+0.j -1.+0.j -2.+0.j]]
Does this mean that there is no way for this approach to run on a real quantum computer where all you can do is measurements on the qubits? I've seen different approaches online such as QAOA or the use of transformations, but I thought if it's so easy to do it with unitaries and some additions, there should be a way to do it with measurements as well.