「世界初、来店不要のフルオーダーメイド靴」のプレスリリースを行いました。

21-4. D-waveマシンでテント・アンド・ツリーパズルを解く(大規模編)

やること

20-2はGAでテント・アンド・ツリーの大規模な問題を解きました。

21-2はD-waveで同じく小規模な問題を解きました。

21-3ではD-waveで大規模な問題を解く準備を行いました。

さっそくD-waveで大規模なテント・アンド・ツリーに挑戦してみましょう

実行環境

WinPython3.6をおすすめしています。

WinPython - Browse /WinPython_3.6/3.6.7.0 at SourceForge.net
Portable Scientific Python 2/3 32/64bit Distribution for Windows

パズル問題

必要なパッケージをimportします。

import numpy as np
from scipy import signal
from dwave.system.composites import EmbeddingComposite
from dwave.system.samplers import DWaveSampler
from pyqubo import Binary
from dwave_qbsolv import QBSolv

パズル問題は、まずは21-2と同じ5*5サイズのものを記述します。

#============================
# パズル
#============================
pazzle = ['-----',
          '-t-t-',
          '-----',
          'tt---',
          '----t']
horizontal = [2,0,1,1,1]
vertical = [1,1,0,2,1]


#============================
# パズルの変換とフィルタの用意
#============================
tree = []
for line in pazzle:
    line = line.replace('-', '0')
    line = line.replace('t', '1')
    tree.append(list(line))
tree = np.array(tree, dtype=int)
print('tree\n{}'.format(tree))

tent = np.zeros(tree.shape, dtype=int)
print('tent\n{}'.format(tent))

horizontal = np.array(horizontal)
print('horizontal\n{}'.format(horizontal))
vertical = np.array(vertical)
print('vertical\n{}'.format(vertical))

#畳み込み用のフィルタ
filter_4 = np.array([[0, 1, 0],
                     [1, 0, 1],
                     [0, 1, 0]], dtype=int)
filter_8 = np.array([[1, 1, 1],
                     [1, 0, 1],
                     [1, 1, 1]], dtype=int)
tree
[[0 0 0 0 0]
 [0 1 0 1 0]
 [0 0 0 0 0]
 [1 1 0 0 0]
 [0 0 0 0 1]]
tent
[[0 0 0 0 0]
 [0 0 0 0 0]
 [0 0 0 0 0]
 [0 0 0 0 0]
 [0 0 0 0 0]]
horizontal
[2 0 1 1 1]
vertical
[1 1 0 2 1]

量子ビットの用意

詳しい説明は割愛します。

#============================
# 量子ビットの用意
#============================
#木の4近傍のマス(=テント候補地)
mask_candidate = (signal.convolve2d(tree, filter_4, mode='same') > 0) * (tree == 0)
print('mask_candidate\n{}'.format(mask_candidate))

#量子ビットの数
num_qbits = np.sum(mask_candidate)
print('num_qbits\n{}'.format(num_qbits))

#テント候補地にインデックスを振る(テント候補地以外は-1)
candidate_index = np.ones(tree.shape, dtype=int) * -1
candidate_index[mask_candidate==True] = range(num_qbits)
print('candidate_index\n{}'.format(candidate_index))

#量子ビットを用意する
for i in range(num_qbits):
    command = 'q{:03} = Binary(\'q{:03}\')'.format(i, i)
    print(command)
    exec(command)
mask_candidate
[[False  True False  True False]
 [ True False  True False  True]
 [ True  True False  True False]
 [False False  True False  True]
 [ True  True False  True False]]
num_qbits
13
candidate_index
[[-1  0 -1  1 -1]
 [ 2 -1  3 -1  4]
 [ 5  6 -1  7 -1]
 [-1 -1  8 -1  9]
 [10 11 -1 12 -1]]
q000 = Binary('q000')
q001 = Binary('q001')
q002 = Binary('q002')
q003 = Binary('q003')
q004 = Binary('q004')
q005 = Binary('q005')
q006 = Binary('q006')
q007 = Binary('q007')
q008 = Binary('q008')
q009 = Binary('q009')
q010 = Binary('q010')
q011 = Binary('q011')
q012 = Binary('q012')

自動的に量子ビットが13個定義されました。

条件1|各木の4近傍のマスは1つか2つだけ1になる

面倒くさくてコードはあまり洗練していませんが、条件が自動的に設定されます。

#============================
# 条件を記述する
#============================
#条件を記述する関数
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

#条件1、各木の4近傍のマスについて、1つか2つだけ1になる
#各木について
for i, j in np.array(np.where(tree==1)).T:
    #注目する木の座標
    print(i, j)
    
    #一つの木
    mask_single_tree = np.zeros(tree.shape, dtype=int)
    mask_single_tree[i, j] = 1
    
    #その周囲のテント候補地
    mask_around_single_tree = signal.convolve2d(mask_single_tree, filter_4, mode='same') * mask_candidate
    #print('mask_around_single_tree\n{}'.format(mask_around_single_tree))
    
    #そのインデックス
    qbits_index = candidate_index[mask_around_single_tree>0]
    #print('qbits_index\n{}'.format(qbits_index))
    
    #これらは1つか2つだけ1になる
    H += H_add(qbits_index, 1.5)
1 1
(q000+q002+q003+q006-1.5)**2
1 3
(q001+q003+q004+q007-1.5)**2
3 0
(q005+q010-1.5)**2
3 1
(q006+q008+q011-1.5)**2
4 4
(q009+q012-1.5)**2

条件2|8近傍で隣接する2つのマスは0つまたは1つだけ1になる

#条件2、8近傍で隣接する2つのマスのうち、0つまたは1つだけ1になる
#各テント候補地について
for i, j in np.array(np.where(mask_candidate==1)).T:
    #他のテント候補地が
    for k, l in np.array(np.where(mask_candidate==1)).T:
        #8近傍の関係であれば
        if (k, l) == (i, j+1) or (k, l) == (i+1, j-1) or (k, l) == (i+1, j) or (k, l) == (i+1, j+1):
            #注目する2つのテント候補地
            print((i, j), (k, l))
            
            #これらのインデックス
            qbits_index = [candidate_index[i, j], candidate_index[k, l]]
            #print(qbits_index)
            
            #これらは0つまたは1つだけ1になる
            H += H_add(qbits_index, 0.5)
(0, 1) (1, 0)
(q000+q002-0.5)**2
(0, 1) (1, 2)
(q000+q003-0.5)**2
(0, 3) (1, 2)
(q001+q003-0.5)**2
(0, 3) (1, 4)
(q001+q004-0.5)**2
(1, 0) (2, 0)
(q002+q005-0.5)**2
(1, 0) (2, 1)
(q002+q006-0.5)**2
(1, 2) (2, 1)
(q003+q006-0.5)**2
(1, 2) (2, 3)
(q003+q007-0.5)**2
(1, 4) (2, 3)
(q004+q007-0.5)**2
(2, 0) (2, 1)
(q005+q006-0.5)**2
(2, 1) (3, 2)
(q006+q008-0.5)**2
(2, 3) (3, 2)
(q007+q008-0.5)**2
(2, 3) (3, 4)
(q007+q009-0.5)**2
(3, 2) (4, 1)
(q008+q011-0.5)**2
(3, 2) (4, 3)
(q008+q012-0.5)**2
(3, 4) (4, 3)
(q009+q012-0.5)**2
(4, 0) (4, 1)
(q010+q011-0.5)**2

17個設定されました。

条件3|各行、各列の “テントの数” と “数字” が一致する

#条件3、各行、各列の "テントの数" と "数字" が一致する
#vertical
for i in range(len(vertical)):
    #この列のインデックス
    qbits_index = candidate_index[i, :][candidate_index[i, :] > -1]
    print(qbits_index)
    
    #これらは指定された数だけ1になる
    H += H_add(qbits_index, vertical[i])

#horizontal
for i in range(len(horizontal)):
    #この列のインデックス
    qbits_index = candidate_index[:, i][candidate_index[:, i] > -1]
    print(qbits_index)
    
    #これらは指定された数だけ1になる
    H += H_add(qbits_index, horizontal[i])
[0 1]
(q000+q001-1)**2
[2 3 4]
(q002+q003+q004-1)**2
[5 6 7]
(q005+q006+q007-0)**2
[8 9]
(q008+q009-2)**2
[10 11 12]
(q010+q011+q012-1)**2
[ 2  5 10]
(q002+q005+q010-2)**2
[ 0  6 11]
(q000+q006+q011-0)**2
[3 8]
(q003+q008-1)**2
[ 1  7 12]
(q001+q007+q012-1)**2
[4 9]
(q004+q009-1)**2

5*5サイズの問題で動作を確認してみる

最後のパーツを追加して実行します。

#コンパイル
model = H.compile()
qubo, offset = model.to_qubo()
print('qubo\n{}'.format(qubo))

#D-waveの設定
D_OBJECT = EmbeddingComposite(DWaveSampler(endpoint='https://cloud.dwavesys.com/sapi',
                                           token='your_token',
                                           solver='DW_2000Q_2_1')) #ここでは実機を指定

#問題を分割してサンプリングしまくる
response = QBSolv().sample_qubo(qubo, solver=D_OBJECT, solver_limit=200)

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

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


#============================
# テントの位置を確認
#============================
tent[mask_candidate] = best_sample
print('tent\n{}'.format(tent))
response
Sample:[0, 1, 1, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0] Energy:-24.0 Occurrences:48
Sample:[1, 0, 0, 0, 1, 0, 0, 0, 1, 1, 1, 0, 0] Energy:-20.0 Occurrences:2
Sample:[1, 1, 0, 0, 0, 0, 0, 0, 1, 1, 1, 0, 0] Energy:-20.0 Occurrences:1
Sample:[0, 0, 0, 1, 0, 0, 0, 0, 1, 1, 1, 0, 0] Energy:-20.0 Occurrences:3
best_sample
[0 1 1 0 0 0 0 0 1 1 1 0 0]
tent
[[0 0 0 1 0]
 [1 0 0 0 0]
 [0 0 0 0 0]
 [0 0 1 0 1]
 [1 0 0 0 0]]

21-2と同じ結果になりました。きちんと動作しています。

18*18サイズの問題

ここからが本番です。問題を書き換えるだけで、すべて自動的にやってくれます(感無量)。

#============================
# パズル
#============================
pazzle = ['---tt-----t-------',
          't-------tt---t--t-',
          't-----------t-t---',
          '-----t---t--------',
          '--t-t--------tt---',
          '-----t----------t-',
          '--t--t----t--t---t',
          '------t-t-t-----t-',
          '-tt---t-----t-----',
          '---t----t--t--t---',
          '---------t------t-',
          '------t--t------t-',
          '-t------------t---',
          '--t-t---t-----t---',
          '-----------t------',
          '---tt--t---------t',
          't-----------t-----',
          '----t-tt-t---tt-t-']
horizontal = [3,4,2,5,2,6,2,2,6,2,4,3,3,4,1,5,2,5]
vertical = [5,3,4,2,4,3,4,3,4,2,5,0,6,1,2,5,0,8]
...
best_sample
[1 0 1 0 1 1 1 0 0 1 0 1 0 0 0 0 1 1 1 0 0 1 0 1 0 0 0 0 0 1 0 0 1 0 1 1 1
 0 0 1 0 0 0 0 1 1 0 0 1 1 0 1 1 1 0 0 0 0 0 0 0 0 0 0 0 1 0 1 1 0 1 1 0 1
 1 0 0 0 0 1 0 0 0 0 0 0 0 0 1 0 1 1 1 0 1 0 0 1 0 0 0 0 0 0 0 0 1 1 1 0 0
 1 1 0 1 0 0 0 1 0 0 0 0 0 0 0 0 0 1 0 1 0 0 1 1 0 1 0 1 1 0 0 0 0 0 0 0 0
 0 0 0 1 1 1 1 1 1 1 1]
tent
[[1 0 0 0 0 1 0 0 0 1 0 1 0 1 0 0 0 0]
 [0 0 0 1 0 0 0 1 0 0 0 0 0 0 0 0 0 1]
 [0 1 0 0 0 1 0 0 0 0 0 1 0 0 0 1 0 0]
 [0 0 0 0 0 0 0 0 1 0 0 0 0 1 0 0 0 0]
 [0 1 0 1 0 1 0 0 0 0 0 0 0 0 0 1 0 0]
 [0 0 0 0 0 0 0 0 0 0 1 0 0 1 0 0 0 1]
 [0 1 0 0 1 0 1 0 1 0 0 0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0 0 0 0 0 0 1 0 0 1 0 1]
 [0 0 0 1 0 1 0 0 1 0 1 0 0 0 0 0 0 0]
 [0 1 0 0 0 0 0 0 0 0 0 0 0 0 0 1 0 0]
 [0 0 0 1 0 0 1 0 1 0 0 1 0 0 0 0 0 1]
 [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
 [1 0 1 0 1 0 0 0 0 1 0 0 0 1 0 0 1 0]
 [0 0 0 0 0 0 0 1 0 0 0 0 0 0 0 0 0 0]
 [0 0 0 0 0 0 0 0 0 0 1 0 0 0 1 0 0 0]
 [0 0 1 0 0 1 0 0 1 0 0 0 1 0 0 0 1 0]
 [0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0 0]
 [1 0 0 1 0 1 0 0 1 0 1 0 1 0 0 1 0 1]]

解けました!

タイトルとURLをコピーしました