量子アニーリング（HOBO）でTakuzuパズルを解いてみた

from hobotan import *
import numpy as np
import matplotlib.pyplot as plt

#問題（空マスは5にした）
question = [[1,5,5,5,5,5],
            [5,5,5,5,5,0],
            [1,5,0,5,0,5],
            [5,5,5,5,0,5],
            [5,5,5,5,5,5],
            [1,1,5,5,5,5]]


#量子ビットの用意
q = symbols_list([6, 6], 'q{}_{}')

#式を用意
Hconst = 0

#問題を反映
for i in range(6):
    for j in range(6):
        if question[i][j] != 5:
            Hconst += (q[i, j] - question[i][j])**2

#よこ方向の白は3個
for i in range(6):
    Hconst += (np.sum(q[i, :]) - 3)**2

#たて方向の白は3個
for j in range(6):
    Hconst += (np.sum(q[:, j]) - 3)**2

#直線の3マスが連続で1になったら臨時ペナルティ
for i in range(6):
    for j in range(6 - 3 + 1):
        Hconst += np.prod(q[i, j:j+3])
for j in range(6):
    for i in range(6 - 3 + 1):
        Hconst += np.prod(q[i:i+3, j])

#直線の3マスが連続で0になったら臨時ペナルティ
for i in range(6):
    for j in range(6 - 3 + 1):
        Hconst += np.prod((1 - q[i, j:j+3]))
for j in range(6):
    for i in range(6 - 3 + 1):
        Hconst += np.prod((1 - q[i:i+3, j]))

#式の合体
H = Hconst

#HOBOテンソルにコンパイル
hobo, offset = Compile(H).get_hobo()
print(f'offset\n{offset}')

#サンプラー選択
solver = sampler.MIKASAmpler()

#サンプリング
result = solver.run(hobo, shots=1000, T_num=10000)

#上位3件
for r in result[:3]:
    print(f'Energy {r[1]}, Occurrence {r[2]}')

    #さくっと配列に
    arr, subs = Auto_array(r[0]).get_ndarray('q{}_{}')
    
    #さくっと画像に
    img, subs = Auto_array(r[0]).get_image('q{}_{}')
    plt.figure(figsize=(3, 3))
    plt.imshow(img)
    plt.show()

tensor order = 2
----- compile -----
MODE: GPU
DEVICE: cuda:0
-------------------
offset
160.0
===== sampler =====
MODE: GPU
DEVICE: cuda:0
===================
.......... 1000/10000 min=-156.0 mean=-149.4350128173828
.......... 2000/10000 min=-156.0 mean=-150.55599975585938
.......... 3000/10000 min=-158.0 mean=-151.4910125732422
（中略）
.......... 8000/10000 min=-160.0 mean=-156.65200805664062
.......... 9000/10000 min=-160.0 mean=-156.7930145263672
.......... 10000/10000 min=-160.0 mean=-156.9390106201172
Energy -160.0, Occurrence 6
（画像）
Energy -160.0, Occurrence 3
（画像）
Energy -159.0, Occurrence 1
（画像）

from itertools import combinations

#よこ配列が重複しない
for a, b in combinations(range(6), 2):
    Hconst += (1 - (q[a, 0]-q[b, 0])**2) * (1 - (q[a, 1]-q[b, 1])**2) * (1 - (q[a, 2]-q[b, 2])**2)\
            * (1 - (q[a, 3]-q[b, 3])**2) * (1 - (q[a, 4]-q[b, 4])**2) * (1 - (q[a, 5]-q[b, 5])**2)

#たて配列が重複しない
for a, b in combinations(range(6), 2):
    Hconst += (1 - (q[0, a]-q[0, b])**2) * (1 - (q[1, a]-q[1, b])**2) * (1 - (q[2, a]-q[2, b])**2)\
            * (1 - (q[3, a]-q[3, b])**2) * (1 - (q[4, a]-q[4, b])**2) * (1 - (q[5, a]-q[5, b])**2)

tensor order = 12
----- compile -----
MODE: GPU
DEVICE: cuda:0
-------------------
offset
190.0
===== sampler =====
MODE: GPU
DEVICE: cuda:0
===================
.......... 1000/10000 min=-184.0 mean=-179.20901489257812
.......... 2000/10000 min=-186.0 mean=-180.26100158691406
.......... 3000/10000 min=-186.0 mean=-181.1070098876953
（中略）
.......... 8000/10000 min=-190.0 mean=-186.239013671875
.......... 9000/10000 min=-190.0 mean=-186.4110107421875
.......... 10000/10000 min=-190.0 mean=-186.58700561523438
Energy -190.0, Occurrence 7
（画像）
Energy -189.0, Occurrence 3
（画像）
Energy -189.0, Occurrence 1
（画像）