Utility Function and the Penalties
def penalise(solution):
sol_dists = distances[solution[0:-1], solution[1:]]
sol_dists = numpy.append(sol_dists, (distances[solution[0], solution[-1]],))
sol_pens = penalties[solution[0:-1], solution[1:]]
sol_pens = numpy.append(sol_pens, (penalties[solution[0], solution[-1]],))
utility = sol_dists/(sol_pens+1)
index_a = numpy.argmax(utility)
index_b = index_a+1 if index_a != count-1 else 0
a = solution[index_a]
b = solution[index_b]
penalties[a,b] += 1.0
penalties[b,a] += 1.0