In my first article, I posted I tried to solve TSP with QAOA. However, this used a so-called modularized version of qiskit aqua, which is possible without having to design the circuit yourself. So this time I would like to use QAOA without using qiskit aqua. First, use max cut as the first step. To be honest, I think that other problems can be solved just by changing the Hamiltonian, so I hope that it will be a reference for those who are doing quantum programming with qiskit. There are several sites that explain quantum programming in Japanese, but since they all use different quantum programming languages, I think that my article written in qiskit is also meaningful.
In this article as well, I would like to explain the qiskit code, assuming that the explanation of QAOA itself will be left to other people.
Reference: Quantum Native Dojo
code
First of all, the whole code
python
# coding: utf-8
from qiskit import QuantumCircuit, QuantumRegister, ClassicalRegister, execute
from qiskit.aqua.utils import tensorproduct
from qiskit import BasicAer
from qiskit.quantum_info.analysis import average_data
from scipy.optimize import minimize
import numpy as np
import time
def classical_minimize(cost_func, initial_params, options, method='powell'):
print('classical minimize is starting now... ')
start_time = time.time()
result = minimize(cost_func, initial_params, options=options, method=method)
print('running time: {}'.format(time.time() - start_time))
print('opt_cost: {}'.format(result.fun))
print('opt_params: {}'.format(result.x))
return result.x
class MAXCUT:
#Initial value setting
def __init__(self, n_qubits, weight_matrix, p=1, num_steps=1):
self.n_qubits = n_qubits
self.weight_matrix = weight_matrix
self.P = p
self.num_steps = num_steps
self.edge_list = []
self._make_edge_list()
#Creating an edge list from weight matrix
def _make_edge_list(self):
for i in range(self.weight_matrix.shape[0]):
for j in range(i+1, self.weight_matrix.shape[1]):
if self.weight_matrix[i][j] != 0:
self.edge_list.append([i, j])
'------------------------------------------------------------------------------------------------------------------'
def Zi_Zj(self, q1, q2):
I_mat = np.array([[1, 0], [0, 1]])
Z_mat = np.array([[1, 0], [0, -1]])
if q1 == 0 or q2 == 0:
tensor = Z_mat
else:
tensor = I_mat
for i in range(1, self.n_qubits):
if i == q1 or i == q2:
tensor = tensorproduct(tensor, Z_mat)
else:
tensor = tensorproduct(tensor, I_mat)
return tensor
def observable(self):
obs = np.zeros((2**self.n_qubits, 2**self.n_qubits))
for a_edge in self.edge_list:
q1, q2 = a_edge
obs = obs + 0.5 * self.Zi_Zj(q1, q2)
return obs
'------------------------------------------------------------------------------------------------------------------'
# U_C(gamma)
def add_U_C(self, qc, gamma):
for a_edge in self.edge_list:
q1, q2 = a_edge
qc.cx(q1, q2)
qc.rz(-2**gamma, q2)
qc.cx(q1, q2)
return qc
# U_X(beta)
def add_U_X(self, qc, beta):
for i in range(self.n_qubits):
qc.rx(-2*beta, i)
return qc
def QAOA_output_onelayer(self, params, run=False):
beta, gamma = params
qr = QuantumRegister(self.n_qubits)
cr = ClassicalRegister(self.n_qubits)
qc = QuantumCircuit(qr, cr)
qc.h(range(self.n_qubits))
qc = self.add_U_C(qc, gamma)
qc = self.add_U_X(qc, beta)
qc.measure(range(self.n_qubits), range(self.n_qubits))
NUM_SHOTS = 10000
seed = 1234
backend = BasicAer.get_backend('qasm_simulator')
results = execute(qc, backend, shots=NUM_SHOTS, seed_simulator=seed).result()
counts = results.get_counts(qc)
expectation = average_data(counts, self.observable())
return expectation
'------------------------------------------------------------------------------------------------------------------'
def minimize(self):
initial_params = np.array([0.1, 0.1])
opt_params = classical_minimize(self.QAOA_output_onelayer, initial_params,
options={'maxiter':500}, method='powell')
return opt_params
def run(self):
opt_params = self.minimize()
beta_opt, gamma_opt = opt_params
qr = QuantumRegister(self.n_qubits)
cr = ClassicalRegister(self.n_qubits)
qc = QuantumCircuit(qr, cr)
qc.h(range(self.n_qubits))
qc = self.add_U_C(qc, gamma_opt)
qc = self.add_U_X(qc, beta_opt)
qc.measure(range(self.n_qubits), range(self.n_qubits))
NUM_SHOTS = 10000
seed = 1234
backend = BasicAer.get_backend('qasm_simulator')
results = execute(qc, backend, shots=NUM_SHOTS, seed_simulator=seed).result()
print(results.get_counts())
if __name__ == '__main__':
weight_matrix = np.array([[0, 1, 0, 1],
[1, 0, 1, 0],
[0, 1, 0, 1],
[1, 0, 1, 0]])
cut = MAXCUT(4, weight_matrix, p=1)
cut.run()
I will pick up and explain only the functions that interest me.
python
def Zi_Zj(self, q1, q2):
I_mat = np.array([[1, 0], [0, 1]])
Z_mat = np.array([[1, 0], [0, -1]])
if q1 == 0 or q2 == 0:
tensor = Z_mat
else:
tensor = I_mat
for i in range(1, self.n_qubits):
if i == q1 or i == q2:
tensor = tensorproduct(tensor, Z_mat)
else:
tensor = tensorproduct(tensor, I_mat)
return tensor
def observable(self):
obs = np.zeros((2**self.n_qubits, 2**self.n_qubits))
for a_edge in self.edge_list:
q1, q2 = a_edge
obs = obs + 0.5 * self.Zi_Zj(q1, q2)
return obs
python
<\beta,\gamma|C(Z)|\beta,\gamma>
The C (Z) required to calculate is created here. At this time, it is necessary to calculate the tensor product, but it can be calculated with the tensor product of qiskit.aqua.utils.
To calculate the expected value with qiskit
python
results = execute(qc, backend, shots=NUM_SHOTS, seed_simulator=seed).result()
counts = results.get_counts(qc)
expectation = average_data(counts, self.observable())
You can do it by the method of.
The execution result is as follows.
{'0101': 2578, '1110': 171, '1101': 146, '1001': 793, '1111': 141, '0011': 815, '0111': 122, '0010': 159, '0000': 161, '0100': 170, '0001': 151, '0110': 802, '1000': 160, '1100': 811, '1010': 2682, '1011': 138}
Larger values, 0101 and 1010, can be retrieved as solutions.
Recommended Posts