メインコンテンツへスキップ

QUICK-PDE: ColibriTDによるQiskit Function

備考

Qiskit Functionsは、IBM Quantum® Premium Plan、Flex Plan、およびOn-Prem（IBM Quantum Platform API経由）Planユーザーが利用できる実験的な機能です。プレビューリリースの状態であり、変更される可能性があります。

概要​

ここで紹介する偏微分方程式（PDE）ソルバーは、Quantum Innovative Computing Kit（QUICK）プラットフォーム（QUICK-PDE）の一部であり、Qiskit Functionとしてパッケージ化されています。QUICK-PDE functionを使用すると、IBM Quantum QPU上でドメイン固有の偏微分方程式を解くことができます。この functionは、ColibriTDのH-DES説明論文に記載されているアルゴリズムに基づいています。このアルゴリズムは、計算流体力学（CFD）や材料変形（MD）を皮切りに、複雑なマルチフィジックス問題を解くことができ、他のユースケースも近日公開予定です。

微分方程式を解くために、試行解は直交関数の線形結合（通常はチェビシェフ多項式、より具体的には関数をエンコードするQubit数をnnとした場合の2n2^n個）としてエンコードされ、可変量子Circuit（VQC）の角度によってパラメータ化されます。ansatzは関数をエンコードする状態を生成し、その組み合わせによってすべての点で関数を評価できるオブザーバブルによって評価されます。そのあと、微分方程式がエンコードされた損失関数を評価し、以下に示すようなハイブリッドループで角度を微調整します。試行解は徐々に実際の解に近づき、満足のいく結果が得られるまで繰り返されます。

QUICK-PDE functionのワークフロー

このハイブリッドループに加えて、異なるオプティマイザーを連鎖させることもできます。これは、グローバルオプティマイザーで良い角度のセットを見つけ、その後、より細かく調整されたオプティマイザーで勾配に沿って最良の隣接角度のセットを探索したい場合に便利です。計算流体力学（CFD）の場合、デフォルトの最適化シーケンスが最良の結果をもたらします。一方、材料変形（MD）の場合、デフォルトでも良い結果が得られますが、問題固有の利点のためにさらに設定をカスタマイズすることもできます。

関数の各変数に対してQubit数を指定します（変更可能です）。10個の同一Circuitを重ねて、1つの大きなCircuit全体の異なるQubitで10個の同一オブザーバブルを評価することで、CMA最適化プロセス内でノイズ軽減を行い、ノイズ学習法を活用して必要なショット数を大幅に削減できます。

入出力​

計算流体力学​

非粘性バーガーズ方程式は、非粘性流体の流れを次のようにモデル化します：

∂u∂t+u∂u∂x=0,\frac{\partial u}{\partial t} + u\frac{\partial u}{\partial x} = 0,

uuは流体速度場を表します。このユースケースには時間的境界条件があります：初期条件を選択し、システムを緩和させることができます。現在、受け入れられる初期条件は線形関数のみです：ax+bax + b。

CFDの微分方程式の引数は、次のように固定グリッド上にあります：

  • ttは0から0.95の間で30サンプル点。xxは0から0.95の間でステップサイズ0.2375。

材料変形​

このユースケースは、空間に固定された棒材がもう一方の端で引っ張られる一次元引張試験における弾塑性変形に焦点を当てています。問題は次のように記述されます：

u′−σ3K−23ϵ0(σ′σ03)n=0u' - \frac{\sigma}{3K} - \frac{2}{\sqrt{3}}\epsilon_0\left(\frac{\sigma'}{\sigma_0\sqrt{3}}\right)^n = 0

σ′−b=0,\sigma' - b = 0,

KKは引き伸ばされる材料の体積弾性率、nnはべき乗則の指数、bbは単位質量あたりの力、ϵ0\epsilon_0は比例応力限界、σ0\sigma_0は比例ひずみ限界、uuは応力関数、σ\sigmaはひずみ関数を表します。

対象となる棒材は単位長さです。このユースケースには、表面応力tt（棒材を引き伸ばすのに必要な仕事量）の境界条件があります。

MDの微分方程式の引数は、次のように固定グリッド上にあります：

  • xxは0から1の間でステップサイズ0.04。

入力​

QUICK-PDE Qiskit Functionを実行するために、次のパラメーターを調整できます：

名前型説明ユースケース固有例
use_caseLiteral["MD", "CFD"]解く微分方程式のシステムを選択するトグルいいえ"CFD"
parametersdict[str, Any]微分方程式のパラメーター（詳細は次の表を参照）いいえ{"a": 1.0, "b": 1.0}
nb_qubitsOptional[dict[str, dict[str, int]]]関数および変数ごとのQubit数。最適値はfunctionによって選択されますが、より良い組み合わせを探したい場合はデフォルト値を上書きできますいいえ{"u": {"x": 1, "t":3}}
depthOptional[dict[str, int]]関数ごとのansatzの深さ。最適値はfunctionによって選択されますが、より良い組み合わせを探したい場合はデフォルト値を上書きできますいいえ{"u": 4}
optimizerOptional[list[str]]使用するオプティマイザー。cma Pythonライブラリの"CMAES"か、scipyのオプティマイザーのいずれかMD"SLSQP"
shotsOptional[list[int]]各Circuitの実行に使用するショット数。複数の最適化ステップが必要なため、リストの長さは使用するオプティマイザーの数（CFDの場合は4）と等しくなければなりません。MDはデフォルト[50_000] * nb_optimizers、CFDはデフォルト[5_00, 2_000, 5_000, 10_000]いいえ[15_000, 30_000]
optimizer_optionsOptional[dict[str, Any]]オプティマイザーに渡すオプション。この入力の詳細は使用するオプティマイザーによって異なります。詳細については、使用するオプティマイザーのドキュメントを参照してくださいMD{"maxiter": 50 }
initializationOptional[Literal["RANDOM", "PHYSICALLY_INFORMED"]]ランダムな角度か、賢く選ばれた角度のどちらから開始するかを設定します。賢く選ばれた角度は100%のケースで機能しない可能性があることに注意してください。デフォルトは"PHYSICALLY_INFORMED"です。いいえ"RANDOM"
backend_nameOptional[str]使用するBackend名。いいえ"ibm_torino"
modeOptional[Literal["job", "session", "batch"]]使用する実行モード。デフォルトは"job"です。いいえ"job"

微分方程式のパラメーター（物理パラメーターおよび境界条件）は、次の形式に従う必要があります：

ユースケースキー値の型説明例
CFDafloatuuの初期値の係数1.0
CFDbfloatuuの初期値のオフセット1.0
MDtfloat表面応力12.0
MDKfloat体積弾性率100.0
MDnintべき乗則4.0
MDbfloat単位質量あたりの力10.0
MDepsilon_0float比例応力限界0.1
MDsigma_0float比例ひずみ限界5.0

出力​

出力はサンプル点のリストと、それらの各点での関数の値を含む辞書です：

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit-ibm-catalog
from numpy import array
solution = {
    "functions": {
        "u": array(
            [
                [0.01, 0.1, 1],
                [0.02, 0.2, 2],
                [0.03, 0.3, 3],
                [0.04, 0.4, 4],
            ]
        ),
    },
    "samples": {
        "t": array([0.1, 0.2, 0.3, 0.4]),
        "x": array([0.5, 0.6, 0.7]),
    },
}

解の配列の形状は変数のサンプルに依存します：

assert len(solution["functions"]["u"].shape) == len(solution["samples"])
for col_size, samples in zip(
    solution["functions"]["u"].shape, solution["samples"].values()
):
    assert col_size == len(samples)

関数変数のサンプル点と解の配列の次元とのマッピングは、変数名の英数字順で行われます。たとえば、変数が"t"と"x"の場合、solution["functions"]["u"]の行は固定された"t"での解の値を表し、solution["functions"]["u"]の列は固定された"x"での解の値を表します。

特定の座標セットにおける関数の値を取得する方法の例を以下に示します：

# u(t=0.2, x=0.7) == 2
assert solution["samples"]["t"][1] == 0.2
assert solution["samples"]["x"][2] == 0.7
assert solution["functions"]["u"][1, 2] == 2

ベンチマーク​

次の表は、このfunctionの様々な実行に関する統計を示しています。

例Qubit数初期化エラー合計時間（分）ランタイム使用量（分）
非粘性バーガーズ方程式50PHYSICALLY_INFORMED10−210^{-2}6625
弾塑性一次元引張試験18RANDOM10−210^{-2}123100

はじめに​

QUICK-PDE functionへのアクセスリクエストフォームに記入してください。そのあと、すでにお使いのローカル環境にアカウントを保存済みであると仮定して、次のようにfunctionを選択します：

from qiskit_ibm_catalog import QiskitFunctionsCatalog

catalog = QiskitFunctionsCatalog(channel="ibm_quantum_platform")

quick = catalog.load("colibritd/quick-pde")

例​

はじめに、次の例のいずれかを試してみてください：

計算流体力学（CFD）​

初期条件をu(0,x)=xu(0,x) = xに設定した場合、結果は次のようになります：

# launch the simulation with initial conditions u(0,x) = a*x + b
job = quick.run(use_case="cfd", physical_parameters={"a": 1.0, "b": 0.0})

Qiskit Function ワークロードのステータスを確認するか、次のように結果を取得します：

print(job.status())
solution = job.result()
'QUEUED'
import numpy as np
import matplotlib.pyplot as plt

_ = plt.figure()
ax = plt.axes(projection="3d")

# plot the solution using the 3d plotting capabilities of pyplot
t, x = np.meshgrid(solution["samples"]["t"], solution["samples"]["x"])
ax.plot_surface(
    t,
    x,
    solution["functions"]["u"],
    edgecolor="royalblue",
    lw=0.25,
    rstride=26,
    cstride=26,
    alpha=0.3,
)
ax.scatter(t, x, solution["functions"]["u"], marker=".")
ax.set(xlabel="t", ylabel="x", zlabel="u(t,x)")

plt.show()

前のコードセルの出力

材料変形​

材料変形のユースケースでは、使用する材料の物理パラメーターと印加する力が必要です。次のように指定します：

import matplotlib.pyplot as plt

# select the properties of your material
job = quick.run(
    use_case="md",
    physical_parameters={
        "t": 12.0,
        "K": 100.0,
        "n": 4.0,
        "b": 10.0,
        "epsilon_0": 0.1,
        "sigma_0": 5.0,
    },
)

# plot the result
solution = job.result()

_ = plt.figure()
stress_plot = plt.subplot(211)
plt.plot(solution["samples"]["x"], solution["functions"]["u"])
strain_plot = plt.subplot(212)
plt.plot(solution["samples"]["x"], solution["functions"]["sigma"])

plt.show()

前のコードセルの出力

エラーメッセージの取得​

ワークロードのステータスがERRORの場合、job.error_message()を使用してエラーメッセージを取得し、デバッグに役立ててください。次のように実行します：

job = quick.run(use_case="mdf", physical_params={})

print(job.error_message())

# or write a wrapper around it for a more human readable version
def pprint_error(job):
    print("".join(eval(job.error_message())["error"]))

print("___")
pprint_error(job)
{"error": ["qiskit.exceptions.QiskitError: 'Unknown argument \"physical_params\", did you make a typo? -- https://docs.quantum.ibm.com/errors#1804'\n"]}
___
qiskit.exceptions.QiskitError: 'Unknown argument "physical_params", did you make a typo? -- https://docs.quantum.ibm.com/errors#1804'

サポートを受ける​

サポートについては、qiskit-function-support@colibritd.com にお問い合わせください。

次のステップ​

推奨事項