import numpy as np
from scipy.optimize import minimize
from quafu import QuantumCircuit, Task
import matplotlib.pyplot as plt
# 1. 定义五个点的坐标
points = np.array([
[0.0, 0.0], # 城市0
[0.3, 0.2], # 城市1
[0.7, 0.6] # 城市2
])
# 2. 计算距离矩阵
def compute_distance_matrix(points):
n = len(points)
dist_matrix = np.zeros((n, n))
for i in range(n):
for j in range(n):
if i != j:
dist_matrix[i, j] = np.linalg.norm(points[i] - points[j])
return dist_matrix
distance_matrix = compute_distance_matrix(points)
print("距离矩阵:\n", distance_matrix)
# 3. QAOA参数设置
num_cities = 3
num_qubits = num_cities**2 # 5x5=25个量子比特
p = 3 # QAOA层数
lambda_penalty = np.max(distance_matrix) * 2 # 惩罚系数
# 4. 构建QAOA量子电路
def build_qaoa_circuit(params, shots=1000):
"""
构建QAOA量子电路
:param params: 优化参数 [γ1, γ2, γ3, β1, β2, β3]
:param shots: 测量次数
:return: 量子电路和任务对象
"""
gamma = params[:p]
beta = params[p:]
qc = QuantumCircuit(num_qubits)
# 初始态制备:均匀叠加态
for q in range(num_qubits):
qc.h(q)
# QAOA层
for layer in range(p):
# 问题哈密顿量演化 (约束和目标)
add_constraints(qc, gamma[layer], lambda_penalty)
add_objective(qc, gamma[layer], distance_matrix)
# 混合哈密顿量演化
for q in range(num_qubits):
qc.rx(q, 2 * beta[layer])
qc.measure(list(range(num_qubits)))
# 创建任务
task = Task()
task.config(backend="ScQ-Sim10", shots=shots) # 使用Quafu模拟器
#task.load(qc) # 显式加载量子电路
res = task.send(qc)
print(res.counts) #counts
print(res.probabilities) #probabilities
#res.plot_probabilities()
return qc, task
def add_constraints(qc, gamma, penalty):
"""添加约束条件演化算子"""
# 每个时间步只有一个城市被访问
for t in range(num_cities):
qubits = [t * num_cities + i for i in range(num_cities)]
add_zz_terms(qc, gamma * penalty, qubits)
# 每个城市只被访问一次
for i in range(num_cities):
qubits = [t * num_cities + i for t in range(num_cities)]
add_zz_terms(qc, gamma * penalty, qubits)
def add_objective(qc, gamma, dist_matrix):
"""添加目标函数演化算子"""
for i in range(num_cities):
for j in range(num_cities):
if i != j:
for t in range(num_cities):
# 连接城市i在时间t和城市j在时间t+1
q1 = t * num_cities + i
q2 = ((t + 1) % num_cities) * num_cities + j
add_zz_terms(qc, gamma * dist_matrix[i, j], [q1, q2])
def add_zz_terms(qc, angle, qubits):
"""添加ZZ相互作用项"""
if len(qubits) < 2:
return
for i in range(len(qubits)):
qc.h(qubits[i])
for i in range(len(qubits) - 1):
qc.cnot(qubits[i], qubits[i + 1])
qc.rz(qubits[-1], 2 * angle)
for i in range(len(qubits) - 2, -1, -1):
qc.cnot(qubits[i], qubits[i + 1])
for i in range(len(qubits)):
qc.h(qubits[i])
# 5. 计算期望能量
def calculate_energy(counts, dist_matrix):
total_energy = 0
total_shots = sum(counts.values())
for bitstring, count in counts.items():
path = decode_path(bitstring)
if path is None:
# 无效路径赋予高惩罚(如最大距离乘以城市数)
energy = lambda_penalty * num_cities
else:
path_length = sum(dist_matrix[path[t], path[(t+1)%num_cities]]
for t in range(num_cities))
penalty = calculate_penalty(bitstring)
energy = path_length + lambda_penalty * penalty
total_energy += energy * count
return total_energy / total_shots
def decode_path(bitstring):
"""改进版路径解码,无效时返回默认路径"""
path = [-1] * num_cities
# 检查比特串长度是否匹配
if len(bitstring) != num_qubits:
return list(range(num_cities)) # 返回默认路径
for t in range(num_cities):
city_count = 0
for i in range(num_cities):
pos = t * num_cities + i
if pos < len(bitstring) and bitstring[pos] == '1':
city_count += 1
path[t] = i
if city_count != 1:
return list(range(num_cities)) # 返回默认路径
if len(set(path)) != num_cities:
return list(range(num_cities)) # 返回默认路径
return path
def calculate_penalty(bitstring):
"""计算约束违反惩罚"""
penalty = 0
# 检查每个时间步是否只有一个城市
for t in range(num_cities):
count = 0
for i in range(num_cities):
pos = t * num_cities + i
if bitstring[pos] == '1':
count += 1
if count != 1:
penalty += (count - 1)**2
# 检查每个城市是否只访问一次
for i in range(num_cities):
count = 0
for t in range(num_cities):
pos = t * num_cities + i
if bitstring[pos] == '1':
count += 1
if count != 1:
penalty += (count - 1)**2
return penalty
# 6. 优化函数
def run_qaoa(shots=1000, max_iter=50):
"""运行QAOA优化过程"""
# 初始化参数
init_params = np.random.uniform(0, np.pi, 2 * p)
# 存储优化过程
history = {'params': [], 'energy': []}
def cost_function(params):
# 构建并运行量子电路
qc, task = build_qaoa_circuit(params, shots)
#async_result = task.submit() # 返回AsyncResult对象
#res = async_result.wait() # 同步等待结果
res = task.send(qc, wait=True) # 同步提交并获取结果
# 错误处理(统一4空格缩进)
if not hasattr(res, 'counts') or len(res.counts) == 0:
print("任务失败:结果无效或为空")
return float('inf')
# 计算能量(与if语句同级)
energy = calculate_energy(res.counts, distance_matrix)
# 记录历史
history['params'].append(params.copy())
history['energy'].append(energy)
print(f"迭代 {len(history['energy'])}: 能量 = {energy:.4f}")
return energy
# 使用COBYLA优化器
result = minimize(cost_function, init_params, method='COBYLA',
options={'maxiter': max_iter})
return result, history
# 7. 主程序
#import matplotlib.pyplot as plt
#plt.rcParams['font.sans-serif'] = ['SimHei', 'Microsoft YaHei', 'WenQuanYi Micro Hei'] # 多字体备选
#plt.rcParams['axes.unicode_minus'] = False # 解决负号显示问题
if __name__ == "__main__":
# 运行QAOA优化
result, history = run_qaoa(shots=1000, max_iter=30)
# 输出优化结果
print("\n优化完成!")
# print(f"最优参数: γ = {result.x[:p]}, β = {result.x[p:]}")
# print(f"最小能量: {result.fun:.4f}")
# 获取最优路径
qc_opt, task_opt = build_qaoa_circuit(result.x, shots=5000)
res_opt = task_opt.send(qc_opt, wait=True)
best_bitstring = max(res_opt.counts, key=res_opt.counts.get)
best_path = decode_path(best_bitstring)
if best_path is not None: # 显式检查
print("\n最优路径:", best_path)
print("路径长度:", sum(distance_matrix[best_path[i], best_path[(i+1)%num_cities]]
for i in range(num_cities))) # 使用动态num_cities
else:
print("未找到有效路径,请调整QAOA参数或增加shots")
print("\n最优路径:", best_path)
print("路径长度:", sum(distance_matrix[best_path[i], best_path[(i+1)%5]]
for i in range(5)))
# 可视化优化过程
plt.plot(history['energy'], '-o')
plt.xlabel('iteration times')
plt.ylabel('Energy')
plt.title('QAOA optimization process')
plt.grid(True)
plt.savefig('qaoa_optimization.png')
plt.close()
plt.show()
# 可视化路径o
plt.figure(figsize=(8, 6))
for i, (x, y) in enumerate(points):
plt.scatter(x, y, s=200, label=f'City{i}')
plt.text(x, y, f' {i}', fontsize=12)
# 绘制路径
for i in range(num_cities):
start = points[best_path[i]]
end = points[best_path[(i+1)%num_cities]]
plt.plot([start[0], end[0]], [start[1], end[1]], '-b')
plt.title('Optimal travel route')
plt.xlabel('X')
plt.ylabel('Y')
plt.grid(True)
plt.legend()
plt.savefig('optimal_path.png')
plt.close()
plt.show()WARNING:matplotlib.font_manager:findfont: Generic family 'sans-serif' not found because none of the following families were found: SimHei, Microsoft YaHei, WenQuanYi Micro Hei
WARNING:matplotlib.font_manager:findfont: Generic family 'sans-serif' not found because none of the following families were found: SimHei, Microsoft YaHei, WenQuanYi Micro Hei
WARNING:matplotlib.font_manager:findfont: Generic family 'sans-serif' not found because none of the following families were found: SimHei, Microsoft YaHei, WenQuanYi Micro Hei
WARNING:matplotlib.font_manager:findfont: Generic family 'sans-serif' not found because none of the following families were found: SimHei, Microsoft YaHei, WenQuanYi Micro Hei
WARNING:matplotlib.font_manager:findfont: Generic family 'sans-serif' not found because none of the following families were found: SimHei, Microsoft YaHei, WenQuanYi Micro Hei
最新发布