محاكاة ديناميكيات Bose-Hubbard مع اختيار الترميز

يبني هذا الدليل دائرة لـالديناميكيات الزمنية لنموذج Bose-Hubbard أحادي البعد، بدءًا من كتابة الهاملتوني وحتى مقارنة كثافات المواقع بعد التطور مع مرجع دقيق. ما يجعل هذا المسار مثيرًا للاهتمام هو فكرة واحدة:

لا تختر ترميز البوزونات إلى كيوبتات عند بناء الدائرة. ابنِ الدائرة على أنماط بوزونية، ودع المترجم (transpiler) يقرر كيف تتحول هذه الأنماط إلى كيوبتات.

في حالة البوزونات، تحتاج دائمًا تقريبًا إلى تأجيل الترميز. يستطيع النمط البوزوني حمل عدد غير محدود من الجسيمات، لذا فإن فضاء Fock له لا نهائي الأبعاد ولا يوجد ترميز إلا بعد اختيار اقتطاع $d$ أولًا (عدد مستويات الإشغال التي تحتفظ بها). في الشيفرة يُعطى الاقتطاع على هيئة أكبر إشغال يمكن أن يحمله النمط، max_occupation $= d - 1$، فالنمط الذي يحتفظ بالمستويات $0, 1, 2$ له max_occupation = 2 و$d = 3$.

وهذا يجعل الاقتطاع معاملًا للترميز وليس للفيزياء، وهو معامل تريد تغييره: إن كان صغيرًا جدًا كانت الإجابة خاطئة، وإن كان كبيرًا جدًا صارت الدائرة عريضة بلا داعٍ. وإذا أجّلت الترميز إلى وقت الترجمة (transpile)، يمكنك مسح قيم الاقتطاع دون تغيير الدائرة أو الهاملتوني أو التجميع.

انظر أيضًا

دليل المعاملات لبنية بيانات المعامل، ودليل المُحوِّلات للترميزات ووصفة كتابة محوِّل مخصص، ودليل الترجمة لمراحل المترجم (transpiler) المشار إليها في كل مكان.

الاقتطاع وكلفته

عندما تقتطع نمطًا إلى $d$ من المستويات، فإنك تكسر المبدّل القانوني عند أعلى مستوى محتفظ به [1]، ولذلك لا يوجد ترميز يمثّل تماثلًا بين الجبور: فالترميز البوزوني يحوّل الحالات، والمعاملات التي تعمل داخل الفضاء الجزئي المقتطع. وهناك نتيجتان تهمان في كل ما يلي:

1. يجب أن يحدّ $d$ الإشغال الذي تستطيع الديناميكيات بلوغه، لا إشغال الحالة الابتدائية فحسب.

2. التحويل لا يتبادل مع الترتيب الطبيعي، لذا قرّر مرة واحدة هل الهاملتوني لديك معرّف بصيغة مرتبة طبيعيًا، والتزم بذلك.

يشتق دليل المُحوِّلات العيب ويُظهر النتيجتين على شيفرة قيد التشغيل.

1. الهاملتوني

نموذج Bose-Hubbard على شبكة ذات روابط $\langle i j \rangle$ هو

$$H = -t \sum_{\langle i j \rangle} \left(b^\dagger_i b_j + b^\dagger_j b_i\right) + \frac{U}{2} \sum_i n_i (n_i - 1) \, ,$$

مع $n_i = b^\dagger_i b_i$. يعدّ حدّ التفاعل أزواج البوزونات على موقع، لذا يتطلب أن يحمل الموقع جسيمين في آن واحد. وبعد الترتيب الطبيعي يصبح $n_i(n_i-1) = b^\dagger_i b^\dagger_i b_i b_i$، وهي الصيغة المستخدمة فيما يلي.

تدخل الهندسة فقط عبر قائمة الحواف، لذا لا شيء هنا خاص بالسلسلة:

>>> from qiskit_bosons.operators import BosonOperator, ann, cre
>>>
>>> def bose_hubbard(edges, num_sites, tunneling, interaction):
...     """Build the normal-ordered Bose-Hubbard Hamiltonian on an arbitrary lattice."""
...     terms = {}
...     for i, j in edges:
...         terms[(cre(i), ann(j))] = -tunneling
...         terms[(cre(j), ann(i))] = -tunneling
...     for i in range(num_sites):
...         # U/2 * n_i (n_i - 1), already normal-ordered
...         terms[(cre(i), cre(i), ann(i), ann(i))] = 0.5 * interaction
...     return BosonOperator.from_dict(terms)
>>>
>>> num_sites = 4
>>> edges = [(i, i + 1) for i in range(num_sites - 1)]
>>> hamiltonian = bose_hubbard(edges, num_sites, tunneling=1.0, interaction=2.0)
>>> print(f"Number of terms: {len(hamiltonian)}")
Number of terms: 10

هناك خاصيتان تستحقان الفحص أولًا. يجب أن يكون الهاملتوني هرميتيًا كي يكون $e^{-itH}$ وحداويًا، وهو يحفظ العدد الكلي للبوزونات:

>>> hamiltonian.is_hermitian()
True
>>> hamiltonian.conserves_particle_number()
True

والآن السؤال البوزوني: كم مستوى يتطلبه هذا المعامل؟

>>> hamiltonian.max_occupation()
2

يطبّق حد التفاعل معاملَي إنشاء على نمط واحد، لذا يجب أن يحمل النمط بوزونين كي يعمل الحد أصلًا. وهذا يعني إشغالًا أقصى قدره 2، وهو الرقم الذي يأخذه LevelSpec على أنه max_occupation. النمط الذي يحتفظ بالإشغالات $0, 1, 2$ له $d = \texttt{max\_occupation} + 1 = 3$ مستويات، لذا فإن max_occupation >= 2 هو المتطلب نفسه $d \geq 3$.

هذا حد أدنى فقط: فهو ناتج عن تطبيق واحد للمعامل، والتطبيق المتكرر عبر تطور زمني قد يصل أبعد. وهو يستبعد القيم الحدّية الصغيرة جزمًا، وتحسم الخطوة 6 حجم القيمة الحدّية المطلوبة لهذا النموذج وهذه الحالة. تُبلغ mode_occupations() عن الحد لكل نمط على حدة، وهذا ما تُقارَن به القيمة الحدّية لكل نمط.

2. التجميع قبل التحويل

يُجرى التطور بطريقة Trotter، وطريقة تقسيم الهاملتوني تحدد ما يُسلَّم إلى التوليف. وهناك اعتباران يدفعان في الاتجاه نفسه.

الأول فيزيائي: حد القفز الفردي $b^\dagger_i b_j$ ليس هرميتيًا. الزوج $b^\dagger_i b_j + b^\dagger_j b_i$ فقط هو الهرميتي. وإذا طوّرت تحت حد واحد غير هرميتي، فلن تكون العملية وحداوية، وسيحمل المعامل المحوَّل معاملات مركّبة يرفضها PauliEvolutionGate. لذا يجب أن يبقى حدّا القفز لكل رابطة معًا.

الثاني متعلق بالكفاءة: جميع حدود التفاعل الموقعي قطرية في أساس الإشغال ومتبادلة فيما بينها، لذا يمكنها أن تشترك في مجموعة واحدة وتُوَلَّف معًا. والأمر نفسه ينطبق على الروابط التي لا تمس موقعًا مشتركًا. في السلسلة، كل رابطة ثانية تشكّل مجموعة كهذه: الرابطتان 0 و2 تعملان على زوجين منفصلين من الأنماط، فتتبادل حدود القفز فيهما ويمكنهما أن تشتركا في مجموعة، وتشكّل الرابطة 1 مجموعة خاصة بها. فإذا جمّعت الروابط الزوجية والفردية، فأنت تحتاج إلى ثلاث مجموعات، بينما تحتاج مجموعة لكل رابطة إلى أربع. المجموعات الأقل تعني عوامل أقل، ودائرة أقل عمقًا أيضًا: فمع رمز Gray عند max_occupation = 2 يكون عمق البوابات ثنائية الكيوبت 95 مقابل 139. تقرر الشبكة أي الروابط يمكن أن تشترك في مجموعة، ولهذا يُختار التقسيم هنا.

تنقل السمة groups هذا القرار إلى الدائرة: تتحلل بوابة Evolution على معامل مجمَّع إلى تطور فرعي واحد لكل مجموعة. يُختار التقسيم على المستوى البوزوني، بينما لا تزال بنية الشبكة مرئية، ويبقى سليمًا عبر التحويل. يغطي دليل التجميع الفحوص التي ينبغي أن يجتازها التجميع.

>>> def group_of(actions, edges):
...     """Assign each term to the even bonds, the odd bonds, or the on-site group."""
...     creations = sorted(mode for action, mode in actions if action)
...     annihilations = sorted(mode for action, mode in actions if not action)
...     if creations == annihilations:
...         return 2  # diagonal: the on-site interaction group
...     bond = tuple(sorted(set(creations) | set(annihilations)))
...     return edges.index(bond) % 2  # bonds 0 and 2 share a group, bond 1 has its own
>>>
>>> hamiltonian.groups = [group_of(actions, edges) for actions, _ in hamiltonian.iter_terms()]
>>> hamiltonian.num_groups()
3

يجب أن تكون كل مجموعة هرميتية بحد ذاتها، وهذا ما يوفره الإقران السابق. تفحص groups_are_hermitian() كل المجموعات في استدعاء واحد:

>>> from qiskit_bosons.operators.terms import groups_are_hermitian
>>> groups_are_hermitian(hamiltonian)
[True, True, True]

لمعرفة ما تحتويه كل مجموعة، افصل المجموعات:

>>> for index, group in enumerate(hamiltonian.split_out_groups()):
...     print(f"group {index}: {len(group)} terms, hermitian={group.is_hermitian()}")
group 0: 4 terms, hermitian=True
group 1: 2 terms, hermitian=True
group 2: 4 terms, hermitian=True

التجميع الزوجي والفردي مجرد اختيار، والبديل سهل القياس. إذا أعطيت كل رابطة مجموعة خاصة بها، فأنت تحتاج إلى أربع مجموعات، وتكلّف الدائرة نفسها حينئذٍ عمقًا أكبر للبوابات ثنائية الكيوبت. يُحتسب العمق بعد تخفيض الدائرة باستخدام رمز Gray عند max_occupation = 2، وهو ما توضحه الخطوة 4:

>>> from qiskit.transpiler import PassManager
>>> from qiskit_bosons.circuit import BosonicCircuit
>>> from qiskit_bosons.circuit.library import Evolution
>>> from qiskit_bosons.mappers import LevelSpec
>>> from qiskit_bosons.transpiler.presets import generate_preset_occupation_pass_manager
>>>
>>> def two_qubit_depth(operator):
...     """The two-qubit depth of one evolution under ``operator``, lowered with the Gray code."""
...     probe = BosonicCircuit(num_sites)
...     probe.append(Evolution(num_sites, operator, time=0.4), probe.modes)
...     pass_manager = generate_preset_occupation_pass_manager(LevelSpec(2, "gray"))
...     pass_manager.qubit = PassManager()
...     lowered = pass_manager.run(probe.decompose())
...     return lowered.decompose(reps=6).depth(lambda instruction: len(instruction.qubits) == 2)
>>>
>>> per_bond = bose_hubbard(edges, num_sites, tunneling=1.0, interaction=2.0)
>>> per_bond.groups = [
...     len(edges) if sorted(m for a, m in actions if a) == sorted(m for a, m in actions if not a)
...     else edges.index(tuple(sorted({m for _, m in actions})))
...     for actions, _ in per_bond.iter_terms()
... ]
>>> per_bond.num_groups(), two_qubit_depth(per_bond)
(4, 139)
>>> hamiltonian.num_groups(), two_qubit_depth(hamiltonian)
(3, 95)

3. الدائرة البوزونية

لاحظ ما هو غائب عن الدائرة. لا يوجد اقتطاع ولا ترميز ولا عدد كيوبتات. إنها تعمل على أربعة أنماط بوزونية:

>>> from qiskit_bosons.circuit import BosonicCircuit
>>> from qiskit_bosons.circuit.library import Evolution, InitializeModes
>>>
>>> total_time = 0.4
>>> occupation = [1, 1, 0, 0]  # one boson on each of the first two sites
>>>
>>> circuit = BosonicCircuit(num_sites)
>>> circuit.append(InitializeModes(occupation), circuit.modes)
>>> circuit.append(Evolution(num_sites, hamiltonian, time=total_time), circuit.modes)
>>> dict(circuit.count_ops())
{'InitializeModes': 1, 'Evolution': 1}

الأسلاك هنا أنماط بوزونية، وليست كيوبتات. لاحظ أن InitializeModes يأخذ أعدادًا صحيحة وليس أعلامًا: فالنمط يمكنه حمل أي عدد من البوزونات، لذا فحالة Fock هي قائمة من أعداد الإشغال.

يحوّل استدعاء واحد لـdecompose() التجميع من الخطوة 2 إلى بوابة واحدة لكل مجموعة:

>>> decomposed = circuit.decompose()
>>> dict(decomposed.count_ops())
{'Evolution': 3, 'InitializeModes': 1}

تحتوي المجموعة الأولى على الروابط بين النمطين 0 و1 وبين النمطين 2 و3، لذا تمتد بوابتها على الأنماط الأربعة جميعها، وكذلك مجموعة الحدود الموقعية. وللرابطة الوسطى بوابة خاصة بها على النمطين 1 و2.

إذا تخطيت هذا التفكيك فسيتم تجاهل التجميع ببساطة. إن الطريقة التي تُجري التقسيم هي synthesis الخاصة بالبوابة، وتكون BosonicLieTrotter هي الافتراضية؛ ويعرض توليف بوابة Evolution في الفضاء البوزوني البدائل.

4. اختر الترميز وقت الترجمة (transpile)

هنا تظهر الفائدة. تُخفَّض الدائرة نفسها أربع مرات، تحت أربع قيم مختلفة لـLevelEncoding، بتغيير معامل واحد. يقرن LevelSpec الترميز باقتطاع كل نمط، ويشتق الإعداد المسبق التالي التخطيط وتحضير الحالة وتحويل المعامل من هذه المواصفة الواحدة (انظر وصف الاقتطاع باستخدام LevelSpec). تُقبل سلسلة نصية مثل "gray" أينما يُقبل LevelEncoding.

>>> from qiskit.transpiler import PassManager
>>> from qiskit_bosons.mappers import LevelEncoding, LevelSpec
>>> from qiskit_bosons.transpiler.presets import generate_preset_occupation_pass_manager
>>>
>>> max_occupation = 2  # the largest occupation per mode, so d = max_occupation + 1 = 3 levels
>>>
>>> encodings = (
...     LevelEncoding.UNARY_PROJECTED,
...     LevelEncoding.UNARY_UNPROJECTED,
...     LevelEncoding.STANDARD_BINARY,
...     LevelEncoding.GRAY,
... )
>>> for encoding in encodings:
...     pass_manager = generate_preset_occupation_pass_manager(
...         LevelSpec(max_occupation, encoding)
...     )
...     pass_manager.qubit = PassManager()  # keep the circuit as the pipeline synthesized it
...     qubit_circuit = pass_manager.run(decomposed)
...     basis = qubit_circuit.decompose(reps=6)
...     two_qubit = basis.depth(lambda instruction: len(instruction.qubits) == 2)
...     print(f"{encoding!s:18s} qubits={qubit_circuit.num_qubits:2d}  2q-depth={two_qubit:4d}")
unary_projected    qubits=12  2q-depth= 645
unary_unprojected  qubits=12  2q-depth= 383
standard_binary    qubits= 8  2q-depth= 218
gray               qubits= 8  2q-depth=  95

هذه الأرقام تناقض التخمين الطبيعي. تستخدم الترميزات الأحادية $d$ من الكيوبتات لكل نمط مقابل $\lceil \log_2 d \rceil$ للرموز المدمجة، وعند $d = 3$ تتفوق الرموز المدمجة على المحورين معًا: كيوبتات أقل وبوابات تشابك أقل بكثير. يختلف النوعان الأحاديان في المعامل الذي يُصدرانه، لا في تخطيطهما: يحمل LevelEncoding.UNARY_PROJECTED مسقطًا على كل كيوبت ثانٍ من النمط، مما يجعله الأغلى بين الأربعة في التطوير، بينما يترك LevelEncoding.UNARY_UNPROJECTED المساقط جانبًا ويتفق معه على الفضاء الجزئي أحادي الإشغال (one-hot). يعرض دليل المُحوِّلات المعاملين جنبًا إلى جنب، ويشرح لماذا يتفوق رمز Gray على الثنائي القياسي (SB) في حدود القفز.

أفاد Sawaya وآخرون بالترتيب نفسه [2]:

مع أن الرمز الأحادي يميل إلى طلب كيوبتات أكثر لكن عمليات أقل، فإن الرمز SB أو Gray غالبًا ما يكون أكثر كفاءة من حيث عدد الكيوبتات وعدد العمليات معًا

استنتاجهم أن أفضل رمز يعتمد على الهاملتوني وعلى $d$، وهذه أقوى حجة لعدم تثبيت الاختيار داخل الدائرة.

ملاحظة

في ظل ترميز المستويات يمتد النمط على عدة كيوبتات، لذا يكون سجل الكيوبتات أعرض من سجل الأنماط. يستخدم الإعداد المسبق OccupationB2QLayout، الذي يشتق عرض السجل من LevelSpec. يستعرض دليل الترجمة المراحل التي يشغّلها الإعداد المسبق.

5. التحقق مقابل الديناميكيات الدقيقة

لا شيء من الخطوات السابقة يثبت أن الدائرة صحيحة. لذلك قارن كثافات المواقع $\langle n_m \rangle$ بعد التطور مع مرجع دقيق محسوب على الهاملتوني المرمَّز.

يجب حساب المرجع في الفضاء المقتطع نفسه الذي تعمل فيه الدائرة، لأن المقارنة بحساب غير مقتطع ستخلط خطأ Trotter بخطأ الاقتطاع، وهما شيئان مختلفان يستحقان الفصل (انظر اقتطاع فضاء Fock).

>>> import numpy as np
>>> from qiskit.quantum_info import SparsePauliOp, Statevector
>>> from scipy.sparse.linalg import expm_multiply
>>> from qiskit_bosons.mappers.library import occupation_encoding
>>>
>>> def sparse_matrix(operator, encoding):
...     """Map a bosonic operator and hand back a sparse qubit matrix."""
...     spec = LevelSpec(max_occupation, encoding).resolve(num_sites)
...     observable = occupation_encoding(operator, spec)
...     return SparsePauliOp.from_sparse_observable(observable).to_matrix(sparse=True)
>>>
>>> def fock_index(occupation, encoding):
...     """The computational basis index encoding a product Fock state."""
...     spec = LevelSpec(max_occupation, encoding).resolve(num_sites)
...     return sum(
...         encoding.level_to_bits(n) << offset for n, offset in zip(occupation, spec.offsets)
...     )
>>>
>>> def densities(state, encoding):
...     """The per-site boson number <n_m> of an encoded state vector."""
...     out = []
...     for mode in range(num_sites):
...         number = BosonOperator.from_dict({(cre(mode), ann(mode)): 1.0})
...         matrix = sparse_matrix(number, encoding)
...         out.append(float((state.conj() @ (matrix @ state)).real))
...     return out

وبعد توفر هذه، إليك التطور الدقيق في الفضاء المقتطع:

>>> encoding = LevelEncoding.GRAY
>>> num_qubits = LevelSpec(max_occupation, encoding).resolve(num_sites).total_width
>>>
>>> initial = np.zeros(2**num_qubits, dtype=complex)
>>> initial[fock_index(occupation, encoding)] = 1.0
>>>
>>> evolved_exact = expm_multiply(
...     -1j * total_time * sparse_matrix(hamiltonian, encoding), initial
... )
>>> reference = densities(evolved_exact, encoding)
>>> for site, value in enumerate(reference):
...     print(f"site {site}: {value:.6f}")
site 0: 0.993834
site 1: 0.858853
site 2: 0.141309
site 3: 0.006004

تنتشر البوزونات من الطرف المشغول للسلسلة نحو الطرف الفارغ، ويبقى البوزونان محفوظين:

>>> bool(np.isclose(sum(reference), 2.0))
True

والآن الدائرة. لاحظ أن المرجع السابق بُني من initial مباشرة، لذا يجب أن تحتوي الدائرة المقارَنة معه على التطور فقط، لأن البوابة InitializeModes ستعيد تحضير حالة محضَّرة بالفعل.

>>> from qiskit_bosons.circuit.library.synthesis import BosonicLieTrotter
>>>
>>> def trotter_densities(steps, encoding):
...     """Site densities after `steps` first-order Trotter steps of the grouped evolution."""
...     circuit = BosonicCircuit(num_sites)
...     circuit.append(
...         Evolution(
...             num_sites, hamiltonian, time=total_time, synthesis=BosonicLieTrotter(reps=steps)
...         ),
...         circuit.modes,
...     )
...     pass_manager = generate_preset_occupation_pass_manager(
...         LevelSpec(max_occupation, encoding)
...     )
...     qubit_circuit = pass_manager.run(circuit.decompose())
...     # decompose BEFORE evolving: see the warning that follows
...     basis = qubit_circuit.decompose(reps=6)
...     state = np.asarray(Statevector(initial).evolve(basis).data)
...     return densities(state, encoding)
>>>
>>> for steps in (1, 2, 4, 8, 16):
...     error = max(
...         abs(a - b) for a, b in zip(trotter_densities(steps, encoding), reference)
...     )
...     print(f"steps={steps:2d}  max density error = {error:.4f}")
steps= 1  max density error = 0.0843
steps= 2  max density error = 0.0389
steps= 4  max density error = 0.0183
steps= 8  max density error = 0.0088
steps=16  max density error = 0.0043

ينخفض الخطأ إلى النصف كلما تضاعف عدد الخطوات، وهذه هي بصمة $O(t^2/\text{steps})$ لصيغة جداء من الرتبة الأولى. وهذا دليل على أن الفرق هو خطأ Trotter وأن الترميز أمين.

مهم

تأتي التكرارات من synthesis الخاصة بالبوابة، التي تطبّق Trotter في الفضاء البوزوني حيث لا يزال التجميع موجودًا: تطبّق كل خطوة كل مجموعة مرة واحدة، وتتبع الخطوات بعضها بعضًا. وهذا يختلف عن تمرير LieTrotter مع reps > 1 إلى إضافة التوليف على مستوى الكيوبت. هناك، تكون كل مجموعة بالفعل PauliEvolutionGate مستقلة، فلا تقسّم الصيغة إلا الحدود داخل المجموعة ولا تتداخل المجموعات أبدًا. المجموعات غير متبادلة فيما بينها، لذا يبقى ذلك خطوة Trotter خارجية واحدة لا ينكمش خطؤها مهما كبر reps. يقارن توليف بوابة Evolution في الفضاء البوزوني بين الاثنين ويغطي الرتب الأعلى.

تحذير

استدعاء decompose() السابق أساسي. تستخدم evolve() التعريف الدقيق لـPauliEvolutionGate وتتجاوز صيغة الجداء تمامًا. وإذا تحققت مقابل دائرة غير مفكَّكة، فسترى تطابقًا شبه تام حتى مع توليف خاطئ عمدًا. كما أن الدائرة غير المفكَّكة تفرض بناء معامل كثيف بحجم $2^{n} \times 2^{n}$، بينما تطبّق الدائرة المفكَّكة البوابات واحدة تلو الأخرى. قد تعامل المحاكيات الأخرى البوابة بشكل مختلف، لذا يعتمد هذا الفحص على المحاكي، بينما تحتاج الأجهزة الفعلية دائمًا إلى الصيغة المفكَّكة.

أخيرًا، الفحص المحوري: يجب ألا تعتمد الفيزياء على الترميز الذي وقع عليه الاختيار.

>>> def exact_densities(encoding):
...     """Exact truncated-space densities, computed under a given encoding."""
...     width = LevelSpec(max_occupation, encoding).resolve(num_sites).total_width
...     start = np.zeros(2**width, dtype=complex)
...     start[fock_index(occupation, encoding)] = 1.0
...     evolved = expm_multiply(
...         -1j * total_time * sparse_matrix(hamiltonian, encoding), start
...     )
...     return densities(evolved, encoding)
>>>
>>> all(
...     bool(np.allclose(exact_densities(name), reference, atol=1e-10))
...     for name in encodings
... )
True

تختلف أعداد الكيوبتات وبُنى Pauli وأعماق الدوائر جميعها، وتتفق الفيزياء حتى عشر خانات عشرية. الترميز خيار تنفيذي، ولهذا ينتمي إلى المترجم (transpiler). يتفق النوعان الأحاديان هنا لأن الحالة المحاكاة تبقى على الفضاء الجزئي أحادي الإشغال (one-hot)؛ وخارجه سيختلفان.

6. الاقتطاع معامل تقارب

الخطوة الأخيرة خاصة بالبوزونات. كل ما سبق ثبّت max_occupation = 2، أي $d = 3$. لكن الاقتطاع ليس اختيارًا حرًا: إنه معامل تقارب، ولا يمكن الوثوق بالإجابة إلا بعد أن تتوقف عن التحرك (انظر اقتطاع فضاء Fock).

لا تحمل الدائرة أي اقتطاع، لذا لا تتغير الدائرة نفسها عندما تتغير القيمة الحدّية: يُعاد فقط إجراء الترجمة (transpilation). ما يتغير هو حجم الناتج، ولا يكون المسح رخيصًا إلى هذا الحد إلا لأن النظام صغير جدًا. المرجع التالي محاكاة كلاسيكية دقيقة للديناميكيات المقتطعة، وهو ممكن لأربعة مواقع لكنه ينمو أُسّيًا مع عدد الكيوبتات. وفي مسألة تتجاوز هذا المدى لا يوجد مرجع دقيق للمقارنة، ويجب الحكم على التقارب بوسائل أخرى.

>>> def densities_at_cutoff(max_occ, encoding=LevelEncoding.GRAY):
...     """Exact truncated-space densities at a given cutoff."""
...     spec = LevelSpec(max_occ, encoding).resolve(num_sites)
...
...     start = np.zeros(2**spec.total_width, dtype=complex)
...     start[
...         sum(encoding.level_to_bits(n) << offset for n, offset in zip(occupation, spec.offsets))
...     ] = 1.0
...
...     def matrix(operator):
...         observable = occupation_encoding(operator, spec)
...         return SparsePauliOp.from_sparse_observable(observable).to_matrix(sparse=True)
...
...     evolved = expm_multiply(-1j * total_time * matrix(hamiltonian), start)
...     out = []
...     for mode in range(num_sites):
...         number = BosonOperator.from_dict({(cre(mode), ann(mode)): 1.0})
...         out.append(float((evolved.conj() @ (matrix(number) @ evolved)).real))
...     return out
>>>
>>> previous = None
>>> for max_occ in (1, 2, 3, 4):
...     current = densities_at_cutoff(max_occ)
...     # Report the drift bucketed rather than to full precision: after it reaches machine
...     # epsilon, the exact digits are arithmetic noise and vary between BLAS builds.
...     if previous is None:
...         drift = "--"
...     else:
...         gap = max(abs(a - b) for a, b in zip(current, previous))
...         drift = f"{gap:.0e}" if gap > 1e-12 else "< 1e-12"
...     print(f"d={max_occ + 1}  "
...           + " ".join(f"{x:.6f}" for x in current)
...           + f"   drift={drift}")
...     previous = current
d=2  0.993986 0.858173 0.141827 0.006014   drift=--
d=3  0.993834 0.858853 0.141309 0.006004   drift=7e-04
d=4  0.993834 0.858853 0.141309 0.006004   drift=< 1e-12
d=5  0.993834 0.858853 0.141309 0.006004   drift=< 1e-12

عند max_occupation = 1 ($d = 2$) يحمل كل موقع بوزونًا واحدًا على الأكثر، لذا لا يستطيع حد التفاعل، الذي يحتاج إلى اثنين على موقع، أن يعمل إطلاقًا، وتنحرف الكثافات بحوالي $7 \times 10^{-4}$. وعند إضافة المستوى الثالث تلتقط النتيجة التفاعل، ومن max_occupation = 3 فصاعدًا تصبح النتيجة ثابتة حتى دقة الآلة: تقاربت.

سبب توقفها عند max_occupation = 2 هو عدد الجسيمات. الهاملتوني يحفظ عدد البوزونات، وهذه الحالة تحمل اثنين، فلا يمكن لأي موقع أن يحمل أكثر من اثنين أبدًا، مهما كان زمن التطور وسواء كان التفاعل تنافريًا أو تجاذبيًا. وإذا طوّرت حتى $t = 6$ بدلًا من $t = 0.4$، أو قلبت إشارة $U$، تبقى الكثافات عند max_occupation = 2 وما فوق دون تغيير حتى $10^{-13}$.

ما يرفع المتطلب فعلًا هو عدد البوزونات. ثلاثة بوزونات على السلسلة نفسها تتغير بحوالي $2 \times 10^{-3}$ بين max_occupation = 2 و3 بعد التطور القصير نفسه، وبحوالي $0.3$ بعد $t = 6$. وهذا كلام عن هذا النموذج وكمية محفوظة، وهو حد فقط للمعاملات التي تحفظ عدد الجسيمات؛ وليس قاعدة عامة، لذا يجب مسح القيمة الحدّية لا تخمينها.

للمسح ثمن على جانب الدائرة أيضًا، وهو ما لا يُظهره المرجع الدقيق. الدائرة البوزونية نفسها، مخفَّضة برمز Gray عند كل قيمة حدّية:

>>> for max_occ in (1, 2, 3, 4):
...     pass_manager = generate_preset_occupation_pass_manager(
...         LevelSpec(max_occ, LevelEncoding.GRAY)
...     )
...     pass_manager.qubit = PassManager()
...     lowered = pass_manager.run(decomposed)
...     depth = lowered.decompose(reps=6).depth(lambda instruction: len(instruction.qubits) == 2)
...     print(f"d={max_occ + 1}  qubits={lowered.num_qubits:2d}  2q-depth={depth:4d}")
d=2  qubits= 4  2q-depth=   8
d=3  qubits= 8  2q-depth=  95
d=4  qubits= 8  2q-depth= 217
d=5  qubits=12  2q-depth=1176

كل مستوى إضافي يكلّف عمقًا، ويقفز عدد الكيوبتات كلما تجاوز $d$ قوة للعدد اثنين. وعندما تتقارب في $d$، فأنت توازن بين الدقة وحجم الدائرة، لذا فالمسح أكثر من مجرد فحص عددي.

لاحظ أيضًا أن كون $d = 2$ صحيحًا تقريبًا مصادفة تتعلق بالمُشاهَد: فـكثافات المواقع بالكاد تتحرك، بينما الحالة الكامنة خاطئة نوعيًا (ليس فيها مكوّن مزدوج الإشغال). إذا تقارب مُشاهَد واحد، فهذا لا يثبت أن غيره قد تقارب.

ما ينبغي أخذه من هذا

• الترميز مسألة تخص الـ Transpiler. الدائرة البوزونية في الخطوة 3 لا تذكر أي اقتطاع ولا ترميز ولا عدد كيوبتات. أربعة ترميزات حوّلتها إلى أربع دوائر مختلفة، وكلها تعيد إنتاج نفس الفيزياء. ولم يكن لازم نغيّر أي شيء قبل الـ Transpiler.

• بالنسبة للبوزونات، هذا مش مجرد راحة. الاقتطاع معامل تقارب لازم تمسحه على قيم مختلفة، وهو موجود داخل الترميز. الدائرة المبنية على ترميز ثابت لازم تنبني من جديد لكل قيمة من $d$.

• جمّع قبل التعيين، وجمّع من أجل الهرميتية. حد القفز الواحد ليس هرميتيًا ولا يمكن رفعه إلى أُس؛ أما زوج الرابطة فهو هرميتي وبالتالي يمكن رفعه إلى أُس. حدّد التقسيم على المستوى البوزوني، حيث لسا الشبكة ظاهرة. هذا ضروري فيزيائيًا، والخاصية groups موجودة لتنقل هذا القرار إلى الأمام.

• لا تفترض أن الترميز ذو الوزن المنخفض هو الأرخص. عند $d = 3$ الشيفرات المدمجة تتفوق على النوعين الأحاديين (unary) في عدد الكيوبتات والبوابات معًا، والنوع المُسقَط هو الأغلى بينها كلها. قِس بنفسك؛ نقطة التقاطع تعتمد على $d$.

• الاقتطاع يكسر الجبر، بطريقة يمكن قياسها. فقط عند أعلى مستوى محتفَظ به، بمقدار $d$. وهذا بالضبط سبب إن التعيين لا يتبادل مع الترتيب الطبيعي.

الخطوات التالية

في اتجاهين يمكن توسيع هذا المثال. إذا استبدلت قائمة حواف السلسلة بروابط شبكة مربعة، ما بتحتاج أي تغيير ثاني، لأن الدالة bose_hubbard المعرّفة قبل مستقلة أصلًا عن الهندسة. وكمان الترميز المخصص من عندك (مثلًا شيفرة أحادية بالكتل تتوسط بين الأحادي والمدمج) يتركّب عبر نفس الوصلة MapperFnEvolutionSynthesis التي يستخدمها الإعداد المسبق.

• مبادئ تصميم تمثيلات المؤثرات يشرح بنية بيانات المؤثر وراء الهاملتوني.

• تعيين تمثيلات المؤثرات يغطي الترميزات وكيف تكتب Mapper خاص فيك.

• عمل Transpile للدوائر البوزونية يمشي معك عبر المراحل التي يشغّلها الإعداد المسبق.

• توليف بوابة Evolution في الفضاء البوزوني يوضح كيف تقلل خطأ Trotter في الفضاء البوزوني.

مهم

كل عمق مذكور محسوب على سجل مجرد ومتصل بالكامل. التعيين على اتصالية العتاد الحقيقي مرحلة منفصلة، وكتل الكيوبتات لكل نمط في ترميز الإشغال تعطي الموجِّه بنية مختلفة تمامًا عن كيوبت واحد لكل نمط.

[1]

R. Somma, G. Ortiz, E. Knill, and J. Gubernatis, Quantum Simulations of Physics Problems (2003), arXiv:quant-ph/0304063, Sec. 4, Eq. 22.

[2]

N. P. D. Sawaya, T. Menke, T. H. Kyaw, S. Johri, A. Aspuru-Guzik, and G. G. Guerreschi, Resource-efficient digital quantum simulation of d-level systems for photonic, vibrational, and spin-s Hamiltonians, npj Quantum Information 6, 49 (2020), arXiv:1909.12847.