21-7. 量子アニーリング（D-wave）でマインスイーパを解いてみた

H += (q0 + q1 + q2 - 2)**2

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

#条件を記述する関数
def H_add(qbits, num):
    command = '(' + '+'.join(['q{:03}'.format(i) for i in qbits]) + '-{})**2'.format(num)
    print(command)
    return eval(command)

    #----- GA -----
    #パラメータ範囲
    near_num = np.sum(mask_near)
    para_range = [[0, 1]] * near_num
    print(para_range)
    
    #多様性がなくなるまで（自動停止するまで）集団を最適化
    myopt = vcopt()
    para, score = myopt.dcGA(para_range,   #パラメータ範囲
                             get_score,    #評価関数
                             -1,           #目標値（到達不可能）
                             show_pool_func='print', #Noneで非表示
                             pool_num=None) #100とか指定しても良いけど
    
    #評価値0に達していなければ終了
    if score > 0:
        print('GA err')
        break
    
    #集団を取得
    pool = myopt.pool
    
    #集団のスコアを調べてスコア0の個体のみに絞る
    scores = np.array([get_score(para) for para in pool])
    print(scores)
    pool = pool[scores==0]
    print(pool)
    #-------------

    #----- 量子アニーリング -----
    #隣接マスの個数の量子ビットを定義
    near_num = np.sum(mask_near)
    for i in range(near_num):
        command = 'q{:03} = Binary(\'q{:03}\')'.format(i, i)
        print(command)
        exec(command)
    
    #量子ビット番号の配列（他は-1）
    f_qbit = np.where(mask_near==True, 0, -1)
    f_qbit[f_qbit==0] = range(near_num)
    print(f_qbit)
    
    #条件対象の数字部分マスク
    mask_near_expand = signal.convolve2d(mask_near, g, mode='same', boundary='fill', fillvalue=0) > 0
    mask_num_select = mask_near_expand * mask_num
    print(mask_num_select)
    
    #条件対象の数字の個数の条件を記述
    H = 0
    ij_list = np.array(np.where(mask_num_select==True)).T
    for k, (i, j) in enumerate(ij_list):
        #print('~~~~~~~~~~~~~~~')
        #print(i, j)
        #条件の数字部分だけのマスク
        mask_tmp = np.zeros(f.shape, bool)
        mask_tmp[i, j] = True
        
        #それと関係がある隣接マスク
        mask_tmp_expand = signal.convolve2d(mask_tmp, g, mode='same', boundary='fill', fillvalue=0) > 0
        mask_near_select = mask_tmp_expand * mask_near
        #print(mask_near_select)
        
        #関係がある量子ビット番号リスト
        qbits = f_qbit[mask_near_select]
        print(qbits)
        
        #隣接の旗の数
        flag_num = np.sum(mask_tmp_expand * mask_flag)
        print(-flag_num)
        
        #式
        H += H_add(qbits, f_num[i, j] - flag_num)
    
    #コンパイル
    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 = sampler.sample_qubo(qubo, num_reads=1000)
    
    #最小エネルギーをチェック
    print('response')
    best_energy = 999
    for (sample, energy, occurrences, _) in response.data():
        print('Sample:{} Energy:{} Occurrences:{}'.format(list(sample.values()), energy, occurrences))
        if energy < best_energy:
            best_energy = energy
            #best_sample = np.array(list(sample.values()))
    print(best_energy)
    
    #最小エネルギーの配列を集める（Occurrencesの数だけ複数入れる）
    pool = []
    for (sample, energy, occurrences, _) in response.data():
        if energy == best_energy:
            for i in range(occurrences):
                pool.append(list(sample.values()))
    pool = np.array(pool)
    print(pool)
    #-------------