# qiskit: Traveling Salesman Problem using QAOA fails for more than 3 cities

I tried to implement the traveling salesman problem (TSP) using QAOA with qiskit. I worked with this qiskit QAOA tutorial and this qiskit minimum eigen optimizer tutorial, where they implement a TSP instance for three cities. However, when I tried to slightly change their example by considering four cities, the QAOA solution fails. I used the NumPyMinimumEigensolver classical solver to compare the QAOA solution with the classical obtained one. While the classical algorithm gives me the right solution, the QAOA solution even fails to visit each city once. I have tried to run QAOA with $$p=1$$ and $$p=2$$, both versions give me the same wrong result. I also used the qasm_simulatorinstead of the statevector_simulator, since for 16 qubits the latter is infeasable.

The weird thing to me is, that QAOA gives me a way too high cost value. I think this is due to the fact that QAOA gives me an invalid solution where the fourth city is never visited (QAOA's solution translates: city 1->2->3 and no city visited in fourth time step).

Does anyone know, why this slight change in the problem makes the algorithm fail? Is there anything wrong in my code?

My code:

from qiskit import Aer
from qiskit.optimization.applications.ising import tsp
from qiskit.aqua.algorithms import NumPyMinimumEigensolver, QAOA
from qiskit.aqua import aqua_globals, QuantumInstance
from qiskit.optimization.algorithms import MinimumEigenOptimizer

# Generating a TSP instance of n cities
n = 4
num_qubits = n ** 2
ins = tsp.random_tsp(n, seed=123)
print('distance\n', ins.w)

# Create a random TSP Hamiltonian
qubitOp, offset = tsp.get_operator(ins)
qp.from_ising(qubitOp, offset, linear=True)

aqua_globals.random_seed = 10598
quantum_instance = QuantumInstance(Aer.get_backend('qasm_simulator'),
seed_simulator=aqua_globals.random_seed,
seed_transpiler=aqua_globals.random_seed)
qaoa_mes = QAOA(quantum_instance=quantum_instance, initial_point=[0., 0.])
exact_mes = NumPyMinimumEigensolver()

# solving Quadratic Program using exact classical eigensolver
exact = MinimumEigenOptimizer(exact_mes)
exact_result = exact.solve(qp)
print("\nExact:\n", exact_result)

# solving the Problem using QAOA
qaoa = MinimumEigenOptimizer(qaoa_mes)
qaoa_result = qaoa.solve(qp)
print("\nQAOA:\n", qaoa_result)


As a result I get:

distance
[[ 0. 48. 91. 33.]
[48.  0. 63. 71.]
[91. 63.  0. 92.]
[33. 71. 92.  0.]]

Exact:
optimal function value: 236.0
optimal value: [0. 0. 0. 1. 1. 0. 0. 0. 0. 1. 0. 0. 0. 0. 1. 0.]
status: SUCCESS

QAOA:
optimal function value: 200111.0
optimal value: [1. 0. 0. 0. 0. 1. 0. 0. 0. 0. 1. 0. 0. 0. 0. 0.]
status: SUCCESS