21-6. 量子アニーリング（D-wave）で野口英世を並び替えて北里柴三郎を作ろうとした

SolverFailureError: Problem not accepted because user *** has no remaining quota in project DEV

pip install -U dwave-ocean-sdk --user

import cv2
import numpy as np
import matplotlib.pyplot as plt

from pyqubo import Binary
from dwave.system.composites import EmbeddingComposite
from dwave.system.samplers import DWaveSampler
from dwave_qbsolv import QBSolv

#============================
# パラメータ
#============================
#画像
noguchi_file = 'vccup2_img1.png'
kitazato_file = 'vccup2_img2.png'

#たてよこ分割数
num = 2 #2, 4, 8, 10, 20など200の約数を指定可能

#1ブロックの幅
w = 200 // num

#============================
# 関数
#============================
#ndarray形式の画像を表示
def show(img):
    plt.figure(figsize=(8, 8))
    plt.imshow(img, vmin = 0, vmax = 255)
    plt.gray()
    plt.show()
    plt.close()
    print()

#ndarray形式のパーツを90°回転（反時計回りにn回）
def rotate(img, n):
    return np.rot90(img, n)

#パーツセットから1枚の画像を再構成
def marge(ims):
    tmp = [cv2.hconcat([ims[i*num + j] for j in range(num)]) for i in range(num)] #横長の短冊にして
    img = cv2.vconcat([tmp[i] for i in range(num)]) #縦に並べる
    return img


#============================
# 画像の読み込み
#============================
#野口英世をグレースケール（1ch）で読み込む（cv2で読むとndarrayになっている）
noguchi_img = cv2.imread(noguchi_file, cv2.IMREAD_GRAYSCALE)
#サイズの確認
print('noguchi_img.shape:{}'.format(noguchi_img.shape))
#表示
show(noguchi_img)

#北里柴三郎も読み込む
kitazato_img = cv2.imread(kitazato_file, cv2.IMREAD_GRAYSCALE)
print('kitazato_img.shape:{}'.format(kitazato_img.shape))
#表示
show(kitazato_img)


#============================
# 画像をパーツに分解する
#============================
#野口
noguchi_ims = []
for tanzaku in np.vsplit(noguchi_img, num): #横長の短冊にして
    noguchi_ims += np.hsplit(tanzaku, num) #正方形に切る
noguchi_ims = np.array(noguchi_ims)
print('noguchi_ims.shape:{}'.format(noguchi_ims.shape))

#北里
kitazato_ims = []
for tanzaku in np.vsplit(kitazato_img, num): #横長の短冊にして
    kitazato_ims += np.hsplit(tanzaku, num) #正方形に切る
kitazato_ims = np.array(kitazato_ims)
print('kitazato_ims.shape:{}'.format(kitazato_ims.shape))

#条件を記述する
H = (q0 - 1 + k)**2

#============================
# どの行と列を見ても一つだけ「１」の状態にするための設定
#============================
#量子ビットを用意する
for i in range(num**2**2):
    command = 'q{:03} = Binary(\'q{:03}\')'.format(i, i)
    print(command)
    exec(command)

#条件を記述する関数（ｎ個の量子ビットからｍ個だけ１にする）
def H_add(qbits, num):
    command = '(' + ' + '.join(['q{:03}'.format(i) for i in qbits]) + ' - {})**2'.format(num)
    print(command)
    return eval(command)

#条件を記述していく
H = 0

#ブロックについてのfor
for i in range(num**2):
    #あるブロックには野口パーツから一つしか採用されない
    qbits = [i*(num**2) + j for j in range(num**2)]
    H += H_add(qbits, 1)
    
#野口パーツについてのfor
for i in range(num**2):
    #ある野口パーツは一つのブロックでしか採用されない
    qbits = [i + j*(num**2) for j in range(num**2)]
    H += H_add(qbits, 1)

#条件を記述する関数２（ｎ番目の量子ビットを１にしにくくする）
def H_add2(qbit, num):
    command = '(q{:03} - 1 + {})**2'.format(qbit, num)
    print(command)
    return eval(command)

H += H_add2(量子ビット番号, 1へのなりにくさ)

#============================
# 北里柴三郎に近くするための設定
#============================
#条件を記述する関数２（ｎ番目の量子ビットを１にしにくくする）
def H_add2(qbit, num):
    command = '(q{:03} - 1 + {})**2'.format(qbit, num)
    print(command)
    return eval(command)

#距離計算、ブロックと野口パーツの組み合わせで調べる
scores = []
for i in range(num**2):
    for j in range(num**2):
        #ブロック[i]に対する野口パーツ[j]の距離を計算
        score = np.mean((noguchi_ims[i] - noguchi_ims[j])**2) #np.mean((kitazato_ims[i] - noguchi_ims[j])**2)
        scores.append(score**2)

#scoresが0-1になるよう規格化
scores /= np.max(scores)
print(np.min(scores))
print(np.max(scores))

#スコアに応じて重みを調整
for i in range(num**2**2):
    #スコアは距離なので、小さいほど結果を1にしやすくしたい
    H += H_add2(i, 1.0*scores[i])

#============================
# D-waveでサンプリング
#============================
#コンパイル
model = H.compile()
qubo, offset = model.to_qubo()
print('qubo\n{}'.format(qubo))

#D-waveの設定
sampler = EmbeddingComposite(DWaveSampler(endpoint='https://cloud.dwavesys.com/sapi/',
                                          token='your_token',
                                          solver='DW_2000Q_6'))

#問題を分割してサンプリングしまくる
response = QBSolv().sample_qubo(qubo, solver=sampler, num_repeats=1000, timeout=5.0)

#レスポンスの確認
count = 0
print('response')
for (sample, energy, occurrences) in response.data():
    print('Sample:{} Energy:{} Occurrences:{}'.format(list(sample.values()), energy, occurrences))
    
    #１番目を保存しておく
    if count == 0:
        best_sample = np.array(list(sample.values()))
    
    count += 1
    if count > 20:
        break

#最頻解の確認
print('best_sample\n{}'.format(best_sample))

#============================
# 野口パーツを当てはめる
#============================
#新しいパーツ配列を用意
new_ims = np.zeros(kitazato_ims.shape, int)

#ブロックについてのfor
for i in range(num**2):
    #あるブロックに対応する量子ビット
    candidate = best_sample[i*(num**2) : (i+1)*(num**2)]
    #採用された野口パーツ
    new_ims[i] = noguchi_ims[np.argmax(candidate)]

#画像を再構成して表示
new_img = marge(new_ims)
show(new_img)

import cv2
import numpy as np
import matplotlib.pyplot as plt

from pyqubo import Binary
import dimod
from dwave.system import LeapHybridSampler

#============================
# D-waveでサンプリング
#============================
#コンパイル
model = H.compile()
qubo, offset = model.to_qubo()
bqm = dimod.BQM.from_numpy_matrix(qubo)

#設定とサンプリング
response = LeapHybridSampler(token='your_token').sample(bqm)