見出し画像

QAP.06:制約条件を仕込む2(ヘルパー関数)【量子コンピュータ/アニーリング@Python/Fixstars Amplify】

【はじめに】

前回は「制約条件を数式化」し、「目的関数」の後ろに「ペナルティ関数」を付与した。

(※復習:前回のコードの抜粋)

... ...

coef = [1, 2, -3]

gen = SymbolGenerator(BinaryPoly) # BinaryPolyを生成させる
q = gen.array(len(coef)) # BinaryPolyを指定の数だけ生成

# 目的関数
f = coef[0]*q[0] + coef[1]* q[1] + coef[2]*q[2]

# ペナルティ関数(制約条件に対する数式)
g = (q[0]+q[1]+q[2]-2)**2

# 最終的な数式(エネルギー式):重みは10
h = f + 10*g

... ...


「制約条件」は、「等式制約」の他にも「不等式制約」「論理制約」…など様々あり、『与えられた条件をどのように「ペナルティ関数」として「数式表現する」のかは難しい』

これに対し「Fixstars Amplify」では「制約条件に対する数式化」等、『制約条件にかかわる部分を裏でいい感じにやってくれる「ヘルパー関数」』というものを用意している。

【ヘルパー関数】


「Fixstars Amplify」
が提供している制約条件に対する「ヘルパー関数」には「equa_to()」「less_equal()」「greater_equal()」「clamp()」「one_hot()」「penalty()」といったものがある。

制約条件に対するヘルパー関数

詳細は以下参照。

例えば「前回の例:x + y + z = 2(制約条件)」に対しては「equal_to()」を使うことができる。

... ...
coef = [1, 2, -3]
gen = SymbolGenerator(BinaryPoly)
q = gen.array(len(coef))
... ...


from amplify.constraint import equal_to

# 「x + y + z = 2(制約条件)」の制約条件オブジェクトを生成
g = equal_to(q[0]+q[1]+q[2],2)

print(g) # q_0 + q_1 + q_2 == 2
print(g.penalty) # 2 q_0 q_1 + 2 q_0 q_2 + 2 q_1 q_2 - 3 q_0 - 3 q_1 - 3 q_2 + 4

【制約条件オブジェクト】

ヘルパー関数」を経由して作成された「ペナルティ関数」は「制約条件オブジェクト」になっている。ドキュメントによるとこの「制約条件オブジェクト」の「ペナルティにかかわる設定値」は

・重み:1(変更可能)
・制約違反時は「ペナルティ値:1」として裏で計算をしている

となっている模様。

詳細な処理や仕様はドキュメントに任せるとして、便利なこととしては、

「制約条件」を「制約条件オブジェクト」として「目的関数の後ろ」に付け足す形にすることで、量子アニーリング計算後に「制約条件」を違反している「変数の組み合わせ」は「解候補」から外してくれる

という点。

例えば、「前回の【おまけ】」の場合。

▲「重みの値」によっては、量子アニーリングの挙動としては正しいが、「ペナルティ関数」の影響が弱すぎて「制約条件」を満たしていない「組み合わせ」が出ていた。

制約条件オブジェクト」を使えば、「条件を満たしていない変数の組み合わせを外してくれる」ので、出力された解候補がどんなものであるかを
確認する手間が軽減される。

・・・以上を踏まえて、前回の例題を「ヘルパー関数」を使った書き方にしてみる。

【例題】x + 2y - 3zに制約条件ある場合(ヘルパー関数使用版)

次の「目的関数」の値をできるだけ小さくする「(x, y, z)の組み合わせ」を求めよう。(※x, y, zは「0または1」のどちらかの値をとる)
・目的関数:x + 2y - 3z
・制約条件:x + y + z = 2 (※3つのうち2つの変数が1になる)

今回も実行環境は「Fixstars Amplify」+「Google Colab」。

【1】ライブラリのインストール

! pip install amplify
! pip install --upgrade amplify

【2】クライアントオブジェクトの作成

from amplify.client import FixstarsClient

client = FixstarsClient()
client.parameters.timeout = 1000 # タイムアウトは1000ミリ秒(1秒)
client.parameters.outputs.duplicate = True # みつかった解が複数あってもOK
client.parameters.outputs.num_outputs = 0  # 0: 見つかった解を全て出力

# トークン設定
client.token = "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX" # 無料登録して取得したトークンを指定する

※トークンは掲載上「XXX... ...XXX」としているが、実際は「Fixstars Amplify」に無料登録して取得したトークンを設定する。

【3】変数を用意する(イジング変数/QUBO変数)

今回は「QUBO変数」を使用する。

from amplify import BinaryPoly
from amplify import SymbolGenerator

# 数式:x + 2y - 3zの係数部分相当
coef = [1, 2, -3]

gen = SymbolGenerator(BinaryPoly) # BinaryPolyを生成させる
q = gen.array(len(coef)) # BinaryPolyを指定の数だけ生成
print(q)

【ここまでの実行結果】
[q_0, q_1, q_2]

▲ q_0, q_1, q_2の3つのQUBO変数を生成した。

【4】目的関数fを作る

「目的関数(組み合わせを探す数式):x + 2y - 3z」を作る。
※今回は変数は「(x, y, z) = (q_0, q_1, q_2)」と対応させている。

f = coef[0]*q[0] + coef[1]* q[1] + coef[2]*q[2]

# forループでまとめて作る場合
#f = sum(coef[i]*q[i] for i in range(len(coef)))


# アインシュタイン縮約表記の場合
#from amplify import einsum
#f = einsum("i,i",q,coef)

print(f)

【ここまでの実行結果】
q_0 + 2 q_1 - 3 q_2

ここまでは前回までと同じ。次の「ペナルティ関数gを作る」部分から変わってくる

【5】ペナルティ関数gを作る(ヘルパー関数使用版)

制約条件:x+y+z=2」に対する「ペナルティ関数g」を作成する。ここでは「Fixstars Amplify」が提供する「equal_to()」を使っていく。

from amplify.constraint import equal_to

# ペナルティ関数g(制約条件オブジェクトの生成)
g = equal_to(q[0] + q[1] + q[2], 2)

print(g)
print(g.penalty)

【ここまでの実行結果】
q_0 + q_1 + q_2 == 2
2 q_0 q_1 + 2 q_0 q_2 + 2 q_1 q_2 - 3 q_0 - 3 q_1 - 3 q_2 + 4

▲「制約条件オブジェクト」が「制約条件」とそれに対応する「ペナルティ関数」を生成して保持している。

【6】重みを決めて数式(エネルギー式)を完成させる

後は同じ。最後に「ペナルティ関数g」に対する「重み」を決めて、「最終的な数式(エネルギー式)」を完成させる。


■重みの設定方法について

「ヘルパー関数」で作成した「制約条件オブジェクト」の「デフォルトの重みは1」である。(はじめは非表示になっている。)

重みを変更したい場合は「値を乗算」すれば、指定した重みの値を持つオブジェクトを返してくれる。(表示されるようになる)

【例】重み10を設定した「制約条件オブジェクト」をつくる

# 重み10を設定した制約条件オブジェクトを返す
g2 = 10 * g 
print( g2 )

【実行結果例】
(q_0 + q_1 + q_2 == 2, 10)

▲「重み10」が設定された「制約条件オブジェクト」となった。

「重みの値」が表示されるようになり、その「重みの値」を更に変更したい場合は「[ ]演算子」でアクセスすれば変更可能である。

【例】作成した「制約条件オブジェクト」の重みを10→5に変更する

#重みを5に変更する
g2[1] = 5
print(g2)

【実行結果例】
(q_0 + q_1 + q_2 == 2, 5)



今回は「重み:0.5(※)」を設定して「制約条件オブジェクト」が制約条件を満たさない解を外してくれることを確認する。
前回の【おまけ】で図示した「制約条件を満たしていない組み合わせが最小値を出してしまう場合」の重みの値)

■重みを設定した数式(エネルギー式)

# 最終的な数式(エネルギー式)
# ペナルティ関数に対し「重み0.5」を設定
h = f + 0.5*g

【7】Solverオブジェクトを生成して量子アニーリング計算をする

# Solverオブジェクト生成
from amplify import Solver
solver = Solver(client)


# エネルギー式を与えて問題を解かせる
result = solver.solve(h)

【8】計算結果を確認する

# 計算結果があるかを確認
if len(result) == 0:
 print("no answer")

問題によっては「複数の答え(条件に合う変数の組み合わせ)」がでることもある。そのため、計算結果はイテレータで順次取り出していく。

from amplify import decode_solution

for sol in result:
    energy = sol.energy
    values = sol.values


    res = decode_solution(q, values)

    print(f"energy = {energy}, {values}  Result is  {q} = {res}")

【実行結果】
energy = -2.0, {2: 1, 1: 0, 0: 1} Result is [q_0, q_1, q_2] = [1. 0. 1.]

▲(q_0, q_1, q_2)=(1, 0, 1)で値 -2.0が答えとして返ってきた(正解)
※制約条件を満たさない(q_0, q_1, q_2)=(0, 0, 1)「値:-2.5」は外されている。


■補足:解候補だったが、制約条件で除外された組み合わせを見る

「solverオブジェクト」の「filter_solution」
を変更することで「実行可能解ではない(少なくとも一つの制約条件が満たされない)」も含めることができる。

これを使った場合の結果は次のようになる。

.. ....

# ペナルティの重みは適当につける
h = f + 0.5*g

... ...

from amplify import Solver
solver = Solver(client)
solver.filter_solution = False


result = solver.solve(h)

if len(result) == 0:
  print("no answer")


from amplify import decode_solution

for sol in result:
    energy = sol.energy
    values = sol.values


    res = decode_solution(q, values)

    print(f"energy = {energy}, {values}  Result is  {q} = {res}")

【実行結果】
energy = -3.0, {2: 1, 1: 0, 0: 0} Result is [q_0, q_1, q_2] = [0. 0. 1.]
energy = -2.0, {2: 1, 1: 0, 0: 1} Result is [q_0, q_1, q_2] = [1. 0. 1.]

▲「solver.filter_solution = False」により除外されていた解も保持・出力された

▲「実行可能解ではない(少なくとも一つの制約条件が満たされない)」も出力された。確かに「制約条件オブジェクト」が不適の解を外してくれていた、というのがわかる。

※energyは-3.0と出力されているがまともに計算すると-2.5になる。
(内部的な計算でペナルティ値1を使っているのかもしれないが、制約条件を違反してる時点で捨てているのか、最終的な出力計算結果(energy値)は目的関数部分がとりうる最小値を出力しているっぽい)

【9】全体コード

最後に全体コードを示しておく。

【例題】
・x + 2y - 3zを最小にするx, y, zを探す。
ただし
(i)「x + y + z = 2(制約条件)」を満たすこと。
(ii)制約条件には、「ヘルパー関数」を使い「制約条件オブジェクト」を生成すること。

【実行環境】
・Fixstars Amplify + Google Colab

■事前準備:ライブラリのインストール

! pip install amplify
! pip install --upgrade amplify

■ここからが全体コード

from amplify import BinaryPoly
from amplify import SymbolGenerator
from amplify import Solver
from amplify import decode_solution
from amplify.constraint import equal_to

from amplify.client import FixstarsClient

client = FixstarsClient()
client.parameters.timeout = 1000 # タイムアウトは1000ミリ秒(1秒)
client.parameters.outputs.duplicate = True # みつかった解が複数あってもOK
client.parameters.outputs.num_outputs = 0  # 0: 見つかった解を全て出力

# トークン設定
client.token = "XXXXXXXXXXXXXXXXXXXXXXXXXXXXXXXX" # 無料登録して取得したトークンを指定する



# 係数
coef = [1, 2, -3]

# QUBO変数を3つ用意
gen = SymbolGenerator(BinaryPoly)
q = gen.array(len(coef))
print(q) # [q_0, q_1, q_2]


# 目的関数f
f = coef[0]*q[0] + coef[1]* q[1] + coef[2]*q[2]

# forループでまとめて作る場合
#f = sum(coef[i]*q[i] for i in range(len(coef)))

# アインシュタイン縮約表記の場合
#from amplify import einsum
#f = einsum("i,i",q,coef)

print(f) # q_0 + 2 q_1 - 3 q_2


# ペナルティ関数g(制約条件オブジェクトの生成)
g = equal_to(q[0] + q[1] + q[2], 2)

print(g)
print(g.penalty)


# 最終的な数式(エネルギー式)
# ペナルティ関数に対し「重み0.5」を設定
h = f + 0.5*g


# Solverオブジェクト生成
solver = Solver(client)


# エネルギー式を与えて問題を解かせる
result = solver.solve(h)


# 結果の確認(回答有無チェック)
if len(result) == 0:
 print("no answer")


# 答えを出力
for sol in result:
   energy = sol.energy
   values = sol.values
   res=decode_solution(q, values)
   
   print(f"energy = {energy}, {values}  Result is  {q} = {res}") 

【実行結果】
energy = -2.0, {2: 1, 1: 0, 0: 1} Result is [q_0, q_1, q_2] = [1. 0. 1.]

▲(x, y, z)= (q_0, q_1, q_2) = (1, 0, 1)の時に、値:-2.0という結果が返ってきた。(正解)

もっと応援したいなと思っていただけた場合、よろしければサポートをおねがいします。いただいたサポートは活動費に使わせていただきます。