SqDRIFT回路を生成する

サンプルベース・クォンタム対角化(SQD)では、ビット列をサンプリングする元となるアンザッツを選ぶ必要があります。SqDRIFTという変種は、代わりに、対象のハミルトニアンから直接構築した時間発展回路のアンサンブルを使います。これは、ハミルトニアンの係数に基づいて、より小さな時間発展演算子をサブサンプリングすることで実現され、qDRIFTトロッター化法として知られています。

この入門ガイドでは、そのようなランダム化回路のアンサンブルを生成する方法を示します。

1. ハミルトニアンのセットアップ

このガイドでは、FCIDUMPファイルからN2の電子構造ハミルトニアンを読み込みます。FermionOperatorを構築する方法は他にもあります。そのドキュメントや、qiskit_fermions.operators.libraryもぜひ参照してください。

PYTHON

>>> from qiskit_fermions.operators.library import FCIDump
>>> from qiskit_fermions.operators import FermionOperator
>>>
>>> fcidump = FCIDump.from_file("docs/guides/n2.fcidump")
>>> num_modes = 2 * fcidump.norb
>>> hamil = FermionOperator.from_fcidump(fcidump)

C

#include <qiskit_fermions.h>

QfFCIDump* fcidump = NULL;
qf_fcidump_from_file("docs/guides/n2.fcidump", &fcidump);
QfFermionOperator* hamil = qf_ferm_op_from_fcidump(fcidump);
uint32_t num_modes = 2 * qf_fcidump_norb(fcidump);

2. ハミルトニアンの項をグループ化する

同一の係数を持つ関連する項をグループ化することで、電子構造ハミルトニアンに存在する多くの対称性を活用します。この操作は、qDRIFTプロトコルがサンプリングする演算子係数の分布を変えますが、収束保証には影響しません。重要な点として、対称性で関連する項をグループ化すると、Pauli項が有利に相殺され、状態をそれらの作用の下で時間発展させる際の回路の深さが全体として短くなります。

qiskit_fermions.operators.terms.groupingモジュールは、演算子の項をグループ化するための便利な関数を提供します。詳細はこのガイドで説明しています。

注意

group_terms_by_electronic_structure()の実装は、ハミルトニアンの項が正規順序化されていることを前提としています。

PYTHON

>>> from qiskit_fermions.operators.terms.grouping import group_terms_by_electronic_structure
>>> from qiskit_fermions.operators.terms.ordering import canonical_order
>>>
>>> canon = canonical_order(hamil.normal_ordered().simplify(atol=1e-16))
>>> exit_code = group_terms_by_electronic_structure(canon, num_modes, two_body_physicist_order=False)
>>> assert exit_code is None
>>> print(canon.groups)  # the groups attribute now contains some list of group indices
[0, ...]

C

QfFermionOperator* normal = qf_ferm_op_normal_ordered(hamil, NULL);
QfFermionOperator* simplified = qf_ferm_op_simplify(normal, 1e-16);
QfFermionOperator* canon = qf_ferm_op_canonical_order(simplified);
QfExitCode exit = qf_ferm_op_group_terms_by_electronic_structure(canon, num_modes, false);
assert(exit == QfExitCode_Success);

qf_ferm_op_free(normal);
qf_ferm_op_free(simplified);

ヒント

完全な電子構造ハミルトニアンには、時間発展回路に含めても観測されるビット列に影響を与えず、したがってサンプリングのオーバーヘッドを増やすだけの項が含まれています。そのため、次のステップでEvolutionゲートを構築する前のこの時点で、そのような項をハミルトニアンから除外することをお勧めします。

これに該当するのは、占有数基底で対角な項、つまり数演算子の積($a^\dagger_i a_i$)です。これには、時間発展が回路にグローバル位相をもたらすだけの定数エネルギーオフセット、時間発展が単一量子ビットのZ回転に相当する個々の数演算子、および$n_i n_j$のような高次の積が含まれます。これらはいずれもモードの占有を変えられないため、除外しても励起の内容は失われません。ステップ5での描画ごとの棄却とは異なり、この除外はハミルトニアン自体に対して一度だけ行われるため、プロトコルがサンプリングする係数分布は、描画の途中で再正規化されるのではなく、除外後の演算子から一貫して再構築されます。

filter_diagonal_terms()関数は、そのような項を演算子からインプレースで取り除きます。

PYTHON

>>> from qiskit_fermions.operators.terms.filtering import filter_diagonal_terms
>>>
>>> filter_diagonal_terms(canon)

C

qf_ferm_op_filter_diagonal_terms(canon);

ここで一度だけ除外するほうが、繰り返し除外するよりもはるかに低コストです。FermionicQDriftはトランスパイルされる回路ごとに一度実行されるため、呼び出しのたびではなく最初にハミルトニアンを除外しておけば、そこから生成されるすべての回路で作業をやり直す必要がなくなります。

ヒント

同じ理由で、ここでハミルトニアンの項をグループのインデックス順に並べておくと効果的です。グループのインデックスは項ごとのタグであり、その項がどこに位置するかを示すものではないため、グループの項は一般に演算子全体に散らばっています。上で適用したグループ化は、最初に出現した順にインデックスを割り当てるため、グループの2番目の項は最初の項から離れた位置になります。

FermionicQDriftはグループ全体をサンプリングし、描画された各グループをsplit_out_groups()で検索します。散らばった項の中からグループを見つけるにはすべてを走査する必要があるため、その検索コストは、描画されるグループが1つでも40個でも同じです。group_order()は各グループを連続した1つの区間にまとめ、その後の検索はグループ境界に対する二分探索となり、そのコストは保持している項ではなく描画されたグループに応じたものになります。検索自体のコストは劇的に下がります。ただし、走査がなくなると、マッピングと合成が支配的になるため、特定のトランスパイルでどれだけ節約できるかは、パイプラインの他の処理に依存します。

PYTHON

>>> from qiskit_fermions.operators.terms.ordering import group_order
>>>
>>> canon = group_order(canon)
>>> print(canon.groups[:4])
[0, 0, ...]

C

QfFermionOperator* grouped = qf_ferm_op_group_order(canon);

これは演算子の値を変えずに項を並べ替えるもので、そこから何回回路のランダム化を描画しても、一度行うだけで済みます。

3. 時間発展回路を準備する

ハミルトニアンの時間発展回路と、サンプルを描画する元となる基本回路を準備します。qiskit_fermions.circuit.libraryには、Qiskitの規約に沿ってこれを行うために必要なすべてのコンポーネントが含まれています。

qDRIFTプロトコルは合成法FermionicQDriftであり、Evolutionゲートはそれをsynthesisとして保持します。この手法には、生成する回路をカスタマイズするためのいくつかのパラメーターがあります。

• 励起グループ数で表した各回路の長さ: num_groups

• 任意で、乱数生成器のシード: rng

• 任意で、トランスパイル自体がボトルネックになる場合に使う、事前計算したサンプリング重み: weights(下記の(任意)サンプリング重みを事前計算するを参照)

手法はゲートが保持するため、重みなど、ゲートの演算子を記述するものはその演算子のすぐ隣に置かれます。

PYTHON

>>> from qiskit_fermions.circuit import FermionicCircuit
>>> from qiskit_fermions.circuit.library import Evolution
>>> from qiskit_fermions.circuit.library.synthesis import FermionicQDrift
>>>
>>> time = 1.0  # you can choose a desired scaling factor here
>>> num_groups = 10
>>> qdrift = FermionicQDrift(num_groups, rng=19)
>>> evo_gate = Evolution(num_modes, canon, time, synthesis=qdrift)
>>>
>>> circ = FermionicCircuit(num_modes)
>>> circ.append(evo_gate, circ.modes)

C

// WARNING: Qiskit's C API does not yet allow us to implement circuits
// with custom gate definitions.

注

この例では、フェルミオンモードに粒子を初期化することも、その最終状態を測定することもしていません。

4. 回路をトランスパイルする

qiskit_fermions.transpilerモジュールはQiskitのトランスパイル・パイプラインと直接統合されており、上で構築したFermionicCircuitをQuantumCircuitへ直接トランスパイルできます。

jordan_wigner()のフェルミオン-量子ビット・マッピングを使って、フェルミオンで表現されたハミルトニアンをPauli文字列で表現し直します。これは、F2QSynthesisトランスパイル・パスを使うことで、トランスパイル処理の一部として直接行えます。generate_preset_jw_pass_manager()を使ってFermionicStagedPassManagerを構築すると、すべての回路命令でJordan-Wigner符号化が一貫して使われることが保証されます。

重要な点として、QiskitのDecomposeパスをトランスパイル・パイプラインのoptimizationステージに追加します。これは、上で準備したqDRIFT法である、ゲートが保持する合成法を使ってEvolutionゲートを分解します。これにより、ハミルトニアン全体の時間発展が回路で使われなくなります。その深さは、現在利用可能な量子コンピューティング・ハードウェアの能力を超えてしまうためです。

代わりに、回路がトランスパイルされるたびに、回路ごとにハミルトニアンの項のgroupsを固定数だけサブサンプリングします。これにより、トランスパイル・パイプラインを繰り返し実行することで、SqDRIFTアルゴリズムが必要とする複数の回路ランダム化を生成できます。生成する回路の数はnum_sqdrift_randomizationsです。

PYTHON

>>> from qiskit.transpiler.passes import Decompose
>>> from qiskit_fermions.transpiler import FermionicPassManager
>>> from qiskit_fermions.transpiler.presets import generate_preset_jw_pass_manager
>>>
>>> pm = generate_preset_jw_pass_manager()
>>> pm.optimization = FermionicPassManager([Decompose(["Evolution"])])
>>>
>>> num_sqdrift_randomizations = 10
>>> sqdrift_circuits = [
...     pm.run(circ) for _ in range(num_sqdrift_randomizations)
... ]

C

// WARNING: Qiskit's C API does not yet allow us to implement circuits
// with custom gate definitions, which we therefore also cannot transpile
// via this API.

注

前の例では、ランダム化の列を再現できるように、FermionicQDrift合成法のrngシードを固定しています。

ヒント

同じパイプラインは、他の2か所でもqDRIFT法を受け付けます。F2FSynthesisパスは、すべてのEvolutionゲートに対して1つの手法を名前で受け取ります(F2FSynthesis({"Evolution": FermionicQDrift(num_groups)}))。これは、どのゲートも自身の演算子に固有のものを何も保持していない場合に必要なすべてです。次のステップでは、3つ目であるQDriftTrotterizationパスを使います。

5. 描画された励起を占有数でフィルタリングする

前のステップで除外した対角項に加えて、描画された励起が、占有がすでに固定されている(すべて占有またはすべて非占有の)モードの集合の内部だけで作用する場合があり、その場合は粒子を一方の側からもう一方へ移動させることができません。そのような励起はサンプリングされるビット列について何も教えてくれないため、それに割り当てたスロットは無駄になります。同じFermionicQDrift法をラップするQDriftTrotterizationパスは、そのような励起を描画した時点で棄却して代わりを描画し直すため、結果の回路のnum_groups個のスロットはそれぞれ、両側を結合する励起を持つことになります。

警告

このフィルタリングは無料ではなく、そのため別個のパスになっています。これが保つのは回路の深さではなく、サンプリングされるnum_groups個のスロットという予算です。その予算は固定なので、棄却された励起は削除されるのではなく置き換えられ、安価なものは、合成により高いコストがかかる結合励起に取って代わられます。同じシードから描画した場合、フィルタリングされた回路はフィルタリングされていない回路より深くなると考えてください。

描画を受理された励起に対して再正規化すると、qDRIFTプロトコルがサンプリングする分布が変わります。また、ステップ2のグループ化とは異なり、プロトコルの収束保証も損なわれます。サンプリングされた積は、渡したハミルトニアンの下での発展に平均して一致しなくなります。保持されたすべての励起は受理確率の逆数で重み付けされ、棄却されたものは何も寄与せず、係数の大きさの合計も、ゲートごとの発展時間も、それを補償するようには調整されません。

受理規則は過大近似でもあります。モードのサポートだけを比較するため、振幅の情報は持たず、追跡される集合は増える一方です。占有集合に閉じた励起は積状態では確かに作用しませんが、先に受理された励起によって状態がもつれ、ほとんどのモードが「不確定」と記されるようになると、後の描画ではこの論理は弱まります。

このフィルタリングは、このガイドが説明するビット列サンプリングのワークフロー、つまり回路がSQDの後処理に入力され、サンプリングされたビット列だけが重要となる場合に使ってください。これらの回路から期待値を推定する場合や、他のいかなる形でもqDRIFTの誤差上界に頼る場合には使わないでください。その結果には、このパスが補正しないバイアスが含まれます。

このフィルタリングでは、どのモードが最初に占有されているかを知る必要があります。そのため、回路内でEvolutionゲートの前にInitializeModesゲートが必要です。InitializeModes.from_hartree_fock()は、これを構築する便利な方法です。ここでN2のHartree-Fock参照状態(14個の空間軌道にalpha電子7個とbeta電子7個)に対して1つ追加し、フィルタリングありとなしでサンプリングされた励起を比較します。

PYTHON

>>> from qiskit_fermions.circuit.library import InitializeModes
>>> from qiskit_fermions.transpiler.passes import F2FSynthesis, QDriftTrotterization
>>>
>>> init = InitializeModes.from_hartree_fock(fcidump.norb, (7, 7))
>>>
>>> hf_circ = FermionicCircuit(num_modes)
>>> hf_circ.append(init, hf_circ.modes)
>>> hf_circ.append(Evolution(num_modes, canon, time), hf_circ.modes)
>>>
>>> num_groups = 5
>>> qdrift_unfiltered = F2FSynthesis(
...     {"Evolution": FermionicQDrift(num_groups, rng=3480)}
... )
>>> qdrift_trivial = QDriftTrotterization(
...     FermionicQDrift(num_groups, rng=3480)
... )
>>>
>>> for instruction in FermionicPassManager(qdrift_unfiltered).run(hf_circ)._inner.data:
...     if instruction.operation.name == "Evolution":
...         print(sorted(instruction.operation.operator.get_support()))
[2, 4]
[41, 45, 52]
[15, 45, 55]
[41, 52, 53]
[10, 16, 37, 38]
>>>
>>> for instruction in FermionicPassManager(qdrift_trivial).run(hf_circ)._inner.data:
...     if instruction.operation.name == "Evolution":
...         print(sorted(instruction.operation.operator.get_support()))
[0, 1, 6, 7]
[0, 1, 28, 29]
[4, 13, 55]
[13, 20, 40, 41]
[0, 1, 28, 29]

C

// WARNING: Qiskit's C API does not yet allow us to implement circuits
// with custom gate definitions, which we therefore also cannot transpile
// via this API.

フィルタリングなしで描画された励起は、占有集合(0-6と28-34)にまったく触れておらず、したがって、そのいずれも占有モードと非占有モードの間で粒子を移動させることができません。フィルタリングありでは、5つすべてが棄却され、占有モードと非占有モードを実際に結合する励起に置き換えられます。たとえば、最初に受理された励起[0, 1, 6, 7]は、占有モード0、1、6と非占有モード7の間で粒子を移動させます。

励起が受理されると、そのサポート内のすべてのモードは「不確定」になり、したがって後のサンプルでどちらの役割も担える資格を得るため、受理される励起が増えるにつれて、占有モード集合と非占有モード集合は増え続けます。これが、2番目の励起[0, 1, 28, 29]が受理可能になる理由です。その4つのモードはすべて元々占有されていたものなので、元々非占有だったどのモードとも結合しません。モード0と1が、最初の励起がそれらに触れたことで不確定になった(したがって結合の「非占有」側として資格を得た)ためだけに、受理されます。

この増加によって、この規則は検定ではなく過大近似になります。不確定になるモードが増えるにつれて棄却される描画は減り、フィルタリングがサンプリング分布に最も大きく作用するのは、最初のほうの描画においてです。

フィルタリングが分布をどれだけ動かしたかを判断するには、パスが回路のメタデータに記録する棄却カウントを読みます。フィルタリングされたEvolutionゲートごとに1つのエントリーが保存され、受理数と総数の比が受理確率を推定します。これは、保持された励起が過大に重み付けされた係数です。

PYTHON

>>> # a freshly seeded pass, so the counts do not depend on earlier draws
>>> qdrift_counted = QDriftTrotterization(
...     FermionicQDrift(num_groups, rng=3480)
... )
>>> filtered = FermionicPassManager(qdrift_counted).run(hf_circ)
>>> discarded = filtered.metadata.get("filter_trivial.discarded")
>>> emitted = filtered.metadata.get("filter_trivial.emitted")
>>> print(f"{emitted[0]} kept, {discarded[0]} discarded")
5 kept, 7 discarded
>>> print(f"acceptance probability: {emitted[0] / (emitted[0] + discarded[0]):.2f}")
acceptance probability: 0.42

C

// WARNING: Qiskit's C API does not yet allow us to implement circuits
// with custom gate definitions, which we therefore also cannot transpile
// via this API.

注

フィルタリングには、照合の基準となる占有情報が必要です。PrepareSlaterDeterminantゲートもInitializeModesと同様にこれを与え、OrbitalRotationは、作用するすべてのモードを「不確定」と記します。これらがいずれもない場合や、与えられたモードがすべて片側(すべて占有、またはすべて非占有)に集中している場合、パスはUserWarningを発し、そのEvolutionゲートのサンプリングをフィルタリングせずに残します。残りの項がもはや追跡される2つの集合を結合できないハミルトニアンでは、代わりに、QDriftTrotterization.MAX_SAMPLE_RETRIES回連続で棄却された後、パスはRuntimeErrorを発生させます。

(任意)サンプリング重みを事前計算する

これは、トランスパイル自体がボトルネックになっている場合にのみ使ってください。得られる回路には何の変化もありません。

FermionicQDriftは、各グループをその係数の大きさに比例する確率で描画し、その大きさを呼び出しのたびにハミルトニアンから導出します。つまり、ハミルトニアンの項ごとに1つの値を求め、それをグループごとに1つへと縮約します。この手法はステートレスであるため、結果が毎回同じであるにもかかわらず、ランダム化のたびにその縮約を繰り返します。weights引数を使うと、group_coeff_means()で一度だけ計算して結果を渡すことができ、コストをループの外に移しつつ、手法がキャッシュされたハミルトニアン依存の状態を持たない状態を保てます。

これは、ステップ2の項の並べ替えに対応するものです。あちらはグループの検索をループの外に出し、こちらは重みの縮約を外に出します。それぞれは、もう一方がない場合にのみ目立つため、両方を併せて適用する価値があります。canonはすでに上で並べ替え済みなので、この2つはそのまま組み合わせられます。

重みは特定の演算子のグループを記述するため、その演算子のEvolutionゲートに属します。そのゲート用に手法を構築し、ステップ3とまったく同様に、ゲートのsynthesisとして渡します。

PYTHON

>>> from qiskit_fermions.operators.terms.grouping import group_coeff_means
>>>
>>> weights = group_coeff_means(canon)  # `canon` is group-ordered, from step 2
>>> qdrift_weighted = FermionicQDrift(10, rng=19, weights=weights)
>>> weighted_gate = Evolution(num_modes, canon, time, synthesis=qdrift_weighted)
>>>
>>> weighted_circ = FermionicCircuit(num_modes)
>>> weighted_circ.append(weighted_gate, weighted_circ.modes)
>>>
>>> weighted_circuits = [
...     pm.run(weighted_circ) for _ in range(num_sqdrift_randomizations)
... ]
>>> len(weighted_circuits)
10

C

// WARNING: Qiskit's C API does not yet allow us to implement circuits
// with custom gate definitions, which we therefore also cannot transpile
// via this API.

これらはまさにパスが計算していたはずの重みなので、回路は以前と同じです。繰り返される縮約がなくなっただけです。

演算子がグループを持つ場合、グループ化された演算子はグループ単位でサンプリングされるため、項ごとではなくグループごとに1つのエントリーが想定されます。エントリーは非負でなければなりません。重みは大きさであり、qDRIFTプロトコルは$H = \sum_j h_j H_j$を正の$h_j$で分解するためで、係数の符号は$H_j$に属し、パスは演算子から直接それを読み取ります。

注意

重みは相対的な好みではなく、絶対的な大きさです。その合計は発展時間も決定するため、すべてのエントリーを$\gamma$倍すると、サンプリングは同一のまま、発展時間は$t$ではなく$\gamma t$になります。したがって、ハミルトニアン自身の係数の大きさ以外を渡すと、アンサンブルが近似する発展が変わってしまいます。これが、これがサンプリングを調整する手段ではなく性能向上のためのツールである理由です。

重みは、記述する演算子を持つゲートに属するため、1つの回路が複数のEvolutionゲートを持ち、それぞれが自身の演算子の重みを持つこともあります。ステップ5のQDriftTrotterizationパスも同じ方法でそれらを尊重します。ゲートが保持するFermionicQDriftでそのゲートを描画し、何も保持していないゲートには、パス自体の構築時に渡されたものにフォールバックします。ゲート上の他の種類の合成法はqDRIFTに関する記述ではないため、パスは無視します。F2FSynthesis({"Evolution": ...})のように、パスに手法をゲート名で渡すと、1つの手法がすべてのEvolutionゲートに適用されるため、そこで渡した重みが意味を持つのは、ゲートが1つだけの回路に限られます。重みの長さは、適用される各ゲートに対して検査されます。

(任意)フェルミオンモードのインデックス付けを最適化する

フェルミオンモードのインデックスを付け替えることでフェルミオン励起のスパンの距離を最小化する追加の最適化ステップを、トランスパイル・パイプラインに加えることができます。この最適化はSqDRIFTの論文で導入されたもので、build_excitation_span_minimization_model()によって実装されています。RelabelModesパスを使うことで、トランスパイラー・パイプラインに簡単に挿入できます。

PYTHON

>>> from pyomo.environ import SolverFactory
>>> from qiskit_fermions.transpiler.passes import RelabelModes
>>>
>>> solver = SolverFactory("appsi_highs")
>>> solver.options["time_limit"] = 10
>>>
>>> qdrift = F2FSynthesis({"Evolution": FermionicQDrift(5, rng=19)})
>>> relabel = RelabelModes(solver=solver)
>>>
>>> pm.optimization = FermionicPassManager([qdrift, relabel])
>>>
>>> relabeled_circ = pm.run(circ)
>>> # if the automatic mode relabeling was successful, the circuit's
>>> # metadata will contain the mode `permutation` information

C

// WARNING: This feature is not available via the C API.

注

RelabelModes内の自動最適化(build_excitation_span_minimization_model()を利用します)を使うには、HAS_PYOMOで管理されるオプションの依存関係が必要です。

重要

ラベルの付け替えは、ある量子ビットがどのモードを担うかを変えるため、部分空間対角化は、それを考慮に入れた場合にのみ正しいビット列を得られます。その考慮をどこで行うかは、どこで測定するかによって決まります。

トランスパイルの前に、FermionicCircuit.measure_all()でFermionicCircuitに測定を追加すれば、何もする必要はありません。各測定は読み取るモードを名前で指定しているため、ラベルの付け替えがそれを一緒に引き継ぎ、サンプリングされたビット列は元のモードでインデックス付けされます。

代わりに、トランスパイルされたQuantumCircuitを測定する場合は、回路のメタデータに記録されたpermutationを使ってビット列を後処理する必要があります。手順はRelabelModesのドキュメントを参照してください。

次のステップ

これで回路のアンサンブルを生成できたので、そこからビット列をサンプリングできます。そのためには、回路をハードウェアで実行する必要があります。詳しい手順はQiskitドキュメントを参照してください。

ビット列のサンプルが得られたら、それらをハミルトニアンの係数と組み合わせて、SQDアドオンのチュートリアルで説明されているとおりにSQDの後処理を行えます。