Operatordarstellungen abbilden

Ein Mapper wandelt einen BosonOperator in einen Qubit-Operator um, hier eine SparseObservable. Der Operator selbst weiß nichts von Qubits: Er speichert Leiteroperatoren auf Modi, und ein bosonischer Modus fasst beliebig viele Teilchen. Alles, was den Operator endlich macht (wie viele Besetzungsniveaus du behältst und wie diese Niveaus auf Qubits angeordnet werden), wird beim Mapping festgelegt.

Der Bibliotheks-Mapper occupation_encoding() arbeitet in der Besetzungszahlbasis (Fock-Basis): Das Niveau $n$ eines Modus ist der Zustand mit $n$ Bosonen. Er nimmt den Operator und eine LevelSpec entgegen, die die Abschneidung und die Kodierung enthält. Diese Anleitung erklärt, was diese beiden Entscheidungen bedeuten und was jede kostet. Für einen durchgängigen Weg, der sie verwendet, siehe Bose-Hubbard-Dynamik mit wählbarer Kodierung simulieren; wie der Transpiler dieselben Entscheidungen auf einen ganzen Circuit anwendet, steht unter Bosonische Circuits transpilieren.

Den Fock-Raum abschneiden

Ein Modus, der auf $d = n_\text{max} + 1$ Niveaus beschränkt wird, wobei $n_\text{max}$ sein max_occupation ist, genügt der kanonischen Kommutatorrelation nicht mehr. Schreibt man $\bar{b}$ für den auf diese Niveaus abgeschnittenen Vernichtungsoperator,

$$[\bar{b}, \bar{b}^\dagger] = 1 - d \, \lvert d-1 \rangle\langle d-1 \rvert \, ,$$

gilt die Relation auf jedem Niveau außer dem höchsten [1]. Du kannst das sehen, indem du den Kommutator der beiden Leiteroperatoren abbildest und seine Diagonale auf den Niveaus eines Modus liest. Die folgende Hilfsfunktion baut die dichte Matrix eines abgebildeten Operators für einen einzelnen Modus, was nur für sehr kleine Register sinnvoll ist:

>>> import numpy as np
>>> from qiskit.quantum_info import SparsePauliOp
>>> from qiskit_bosons.mappers import LevelEncoding, LevelSpec
>>> from qiskit_bosons.mappers.library import occupation_encoding
>>> from qiskit_bosons.operators import BosonOperator, ann, cre
>>>
>>> def level_diagonal(operator, max_occupation, encoding):
...     """The diagonal of a mapped one-mode operator on the codewords of the levels."""
...     spec = LevelSpec(max_occupation, encoding)
...     matrix = SparsePauliOp.from_sparse_observable(occupation_encoding(operator, spec)).to_matrix()
...     codewords = [spec.encoding.level_to_bits(level) for level in range(max_occupation + 1)]
...     return [round(float(matrix[c, c].real), 6) for c in codewords]
>>>
>>> b = BosonOperator.from_dict({(ann(0),): 1.0})
>>> b_dagger = BosonOperator.from_dict({(cre(0),): 1.0})
>>> commutator = (b @ b_dagger) - (b_dagger @ b)
>>> level_diagonal(commutator, 2, "unary_unprojected")  # d = 3
[1.0, 1.0, -2.0]

Der letzte Eintrag ist $1 - d$. Weil das bedeutet, dass der Aufspann von $\{\bar{b}^\dagger, \bar{b}\}$ unter dem Kommutator nicht abgeschlossen ist, ist keine Kodierung ein Isomorphismus von Algebren. Eine Kodierung bildet Zustände ab, zusammen mit den Operatoren, die innerhalb des abgeschnittenen Unterraums wirken. Die folgenden zwei Konsequenzen sind in der Praxis wichtig.

Die Abschneidung muss die erreichbare Besetzung begrenzen. Es genügt nicht, die Besetzung zu begrenzen, mit der du startest. Ein Strahlteiler, der auf $\lvert 1, 1 \rangle$ wirkt, bevölkert $\lvert 2, 0 \rangle$ und $\lvert 0, 2 \rangle$, sodass eine Abschneidung bei eins stillschweigend einen falschen Zustand erzeugt:

>>> from scipy.linalg import expm
>>>
>>> beamsplitter = BosonOperator.from_dict({(cre(0), ann(1)): 1.0, (cre(1), ann(0)): 1.0})
>>>
>>> def stays_in_11(max_occupation):
...     """Probability to remain in |1, 1> after a 50:50 beamsplitter."""
...     spec = LevelSpec(max_occupation, "gray").resolve(2)
...     observable = occupation_encoding(beamsplitter, spec, spec.total_width)
...     matrix = SparsePauliOp.from_sparse_observable(observable).to_matrix()
...     unitary = expm(-1j * (np.pi / 4) * matrix)
...     index = sum(spec.encoding.level_to_bits(1) << offset for offset in spec.offsets)
...     return round(float(abs(unitary[index, index]) ** 2), 6)
>>>
>>> stays_in_11(1), stays_in_11(2)
(1.0, 0.0)

Bei einer Abschneidung von eins ist die Wahrscheinlichkeit, in $\lvert 1, 1 \rangle$ zu bleiben, eins, weil die beiden Zustände, die die Dynamik bevölkern würde, nicht existieren. Bei einer Abschneidung von zwei verschwindet sie, was der Hong-Ou-Mandel-Effekt ist. Deshalb ist max_occupation ein Pflichtargument ohne Standardwert: Es ist eine Modellierungsentscheidung und ein Konvergenzparameter, den du durchfahren sollst.

Das Mapping vertauscht nicht mit der Normalordnung. Die Operatoren $b b^\dagger$ und $b^\dagger b + 1$ sind vor der Abschneidung gleich, werden aber auf unterschiedliche Qubit-Operatoren abgebildet:

>>> product = BosonOperator.from_dict({(ann(0), cre(0)): 1.0})
>>> ordered = BosonOperator.from_dict({(cre(0), ann(0)): 1.0}) + BosonOperator.one()
>>> bool(product.equiv(ordered))
True
>>> level_diagonal(product, 3, "unary_unprojected")  # d = 4
[1.0, 2.0, 3.0, 0.0]
>>> level_diagonal(ordered, 3, "unary_unprojected")
[1.0, 2.0, 3.0, 4.0]

Die beiden stimmen auf jedem Niveau außer dem höchsten überein. Auf das oberste Niveau angewandt, hat $b^\dagger$ keinen Platz mehr, und die Amplitude wird abgeschnitten, sodass $b b^\dagger$ dort null ist, während $b^\dagger b + 1$ den Wert $d$ ergibt. Entscheide einmal, ob dein Hamiltonoperator in normalgeordneter Form definiert ist, und bilde diese Form konsequent ab.

Drei Methoden helfen dir bei der Wahl einer Abschneidung, und alle liefern untere Schranken. Der Operator meldet die größte Anzahl von Erzeugungsoperatoren, die er in einem einzelnen Term auf jeden Modus anwendet, und eine Spezifikation vergleicht diese mit ihren Abschneidungen:

>>> quartic = BosonOperator.from_dict({(cre(0), cre(0), ann(0), ann(0)): 1.0})
>>> quartic.mode_occupations(), quartic.max_occupation()
([2], 2)
>>> LevelSpec(1, "gray").truncation(quartic)  # one level short
[1]
>>> LevelSpec(2, "gray").truncation(quartic)
[0]

Ein Eintrag ungleich null bedeutet, dass die Abschneidung unter der Schranke liegt, und der Wert gibt an, um wie viel. Diese Hilfsfunktionen kennen die Dynamik nicht, daher macht eine bestandene Prüfung eine Abschneidung nicht hinreichend; der Strahlteiler oben hat eine Schranke von eins, braucht aber eine Abschneidung von zwei. Sie können auch ein Defizit für einen Term melden, der auf den beibehaltenen Niveaus null ist, wie es $b^\dagger b^\dagger b b$ bei einer Abschneidung von eins ist.

Die Abschneidung mit einer LevelSpec beschreiben

Eine LevelSpec bündelt die Abschneidung jedes Modus mit der Kodierung und leitet daraus ab, wie viele Qubits jeder Modus belegt und wo seine Qubits beginnen. Alles, was darüber übereinstimmen muss (der Mapper, das Transpiler-Layout, die Zustandspräparation), erhält dieselbe Spezifikation. Die Kodierung ist ein Pflichtargument, ebenso die Abschneidung: Beide haben keinen Standardwert.

Die Abschneidung ist entweder eine einzelne Ganzzahl, die für alle Modi gilt, oder eine Ganzzahl pro Modus:

>>> shared = LevelSpec(2, "gray")
>>> per_mode = LevelSpec([4, 2, 2], "gray")
>>> shared.num_modes is None, per_mode.num_modes
(True, 3)

Eine aus einer einzelnen Ganzzahl aufgebaute Spezifikation lässt die Anzahl der Modi offen, daher musst du sie auflösen, bevor du nach Breiten fragst. Eine Sequenz wird nie gebroadcastet, egal wie lang sie ist, und sie legt die Anzahl der Modi fest:

>>> per_mode.widths, per_mode.offsets, per_mode.total_width
([3, 2, 2], [0, 3, 5], 7)
>>> shared.resolve(3).widths
[2, 2, 2]
>>> try:
...     shared.widths
... except ValueError as error:
...     print(error)
this spec was built from a single cutoff, so its number of modes is open; call `resolve(num_modes)` first

Modus $m$ belegt die Qubits offsets[m] bis offsets[m] + widths[m] - 1, sodass Modi mit unterschiedlichen Abschneidungen unterschiedlich viele Qubits belegen können und die Offsets die laufende Summe der Breiten sind. Du musst dem Mapper die Anzahl der Modi nicht übergeben. Er entnimmt sie der Spezifikation oder, bei einer Spezifikation mit einer einzelnen Abschneidung, dem Operator über min_num_modes(). Das Register kann breiter sein, als die Modi es brauchen; die überzähligen Qubits sind untätig und stehen am Ende:

>>> number = BosonOperator.from_dict({(cre(0), ann(0)): 1.0})
>>> number.min_num_modes()
1
>>> occupation_encoding(number, LevelSpec(1, "gray")).num_qubits
1
>>> occupation_encoding(number, LevelSpec(1, "gray"), 4).num_qubits
4

Ein Operator, der auf einen Modus jenseits des letzten einer Spezifikation pro Modus wirkt, wird abgelehnt:

>>> far = BosonOperator.from_dict({(cre(5), ann(5)): 1.0})
>>> try:
...     occupation_encoding(far, LevelSpec([1, 1], "gray"))
... except ValueError as error:
...     print(error)
the operator acts on 6 modes, but the spec covers only 2

Einige Grenzen ergeben sich daraus, wie die Kodierungen Qubits zählen. Bei einer kompakten Kodierung belegt ein Modus mit nur einem Niveau (max_occupation = 0) überhaupt kein Qubit, weil sein einziges Niveau nichts braucht, um sich von einem anderen zu unterscheiden. Ein solcher Modus ist erlaubt und verschiebt die Qubits der nachfolgenden Modi nicht, aber eine Spezifikation, in der jeder Modus kein Qubit belegen würde, wird abgelehnt. Die unären Kodierungen fassen höchstens 128 Niveaus pro Modus, und eine größere Abschneidung löst einen Fehler aus, der eine kompakte Kodierung vorschlägt:

>>> LevelSpec([0, 1], "gray").resolve(2).widths
[0, 1]
>>> try:
...     LevelSpec(0, "gray")
... except ValueError as error:
...     print(str(error).split(":")[0])
the spec would span no qubits
>>> try:
...     LevelSpec(128, "unary_projected")
... except ValueError as error:
...     print(str(error).split(" (")[0])
max_occupation = 128 is too large for the unary encodings, which hold at most 128 levels per mode

Eine Kodierung wählen

Eine LevelEncoding bestimmt, welcher Basiszustand eines Qubit-Blocks jedes Niveau darstellt. Es gibt vier, und keine ist im Allgemeinen die beste [2]. Die Anzahl der Qubits pro Modus und das Bitmuster jedes Niveaus ergeben sich aus der Kodierung:

>>> encodings = (
...     LevelEncoding.UNARY_PROJECTED,
...     LevelEncoding.UNARY_UNPROJECTED,
...     LevelEncoding.STANDARD_BINARY,
...     LevelEncoding.GRAY,
... )
>>> rows = ["encoding           qubits for d = 1, 2, 3, 4, 5, 8    level 0 to 4 as bit patterns"]
>>> for encoding in encodings:
...     widths = [encoding.num_qubits(levels) for levels in (1, 2, 3, 4, 5, 8)]
...     patterns = [f"{encoding.level_to_bits(level):05b}" for level in range(5)]
...     rows.append(f"{encoding!s:18s} {widths!s:34s} {patterns}")
>>> print("\n".join(rows))
encoding           qubits for d = 1, 2, 3, 4, 5, 8    level 0 to 4 as bit patterns
unary_projected    [1, 2, 3, 4, 5, 8]                 ['00001', '00010', '00100', '01000', '10000']
unary_unprojected  [1, 2, 3, 4, 5, 8]                 ['00001', '00010', '00100', '01000', '10000']
standard_binary    [0, 1, 2, 2, 3, 3]                 ['00000', '00001', '00010', '00011', '00100']
gray               [0, 1, 2, 2, 3, 3]                 ['00000', '00001', '00011', '00010', '00110']

Die Bitmuster sind Ganzzahlen, in denen Bit $q$ der Wert von Qubit $q$ ist. Die unären Kodierungen verbrauchen ein Qubit pro Niveau und setzen eines davon. Die kompakten Kodierungen verbrauchen $\lceil \log_2 d \rceil$ Qubits und speichern die Niveaunummer, binär oder als Gray-Code, in dem sich benachbarte Niveaus in einem Bit unterscheiden. Der Zeichenkettenname einer Kodierung wird überall akzeptiert, wo die Kodierung selbst akzeptiert wird.

Der Kompromiss besteht zwischen Qubits und Operatorgewicht. Die kompakten Kodierungen brauchen exponentiell weniger Qubits, aber ein Leiteroperator verbindet Niveaus, deren Muster sich in mehreren Bits unterscheiden können, sodass seine Terme auf mehr Qubits wirken. Die folgenden Abschnitte betrachten, was das kostet, beginnend mit den beiden unären Kodierungen, die identische Layouts haben und sich nur im Operator unterscheiden.

Die beiden unären Kodierungen vergleichen

In der unären Kodierung verschiebt ein Leiteroperator eine Anregung zwischen zwei benachbarten Qubits. Um das als Operator auf den gesamten Modus zu schreiben, musst du festlegen, was mit den übrigen Qubits geschieht. Die beiden Varianten entscheiden das unterschiedlich:

• unary_projected setzt einen Projektor $\lvert 0 \rangle\langle 0 \rvert$ auf jedes andere Qubit des Modus. Das erzwingt die One-Hot-Bedingung im Operator selbst: Der Term wirkt nur, wenn der Rest des Modus leer ist.

• unary_unprojected lässt diese Qubits unberührt, sodass der Term auf ihnen als Identität wirkt.

Betrachte ein einzelnes $b^\dagger$ bei $d = 4$ unter beiden Varianten. Der Mapper gibt eine SparseObservable zurück, deren Terme Zeichenketten aus Einzel-Qubit-Bit-Termen sind. Hier sind die ersten vier Terme jeder Variante mit den Qubits, auf die sie wirken:

>>> for encoding in ("unary_projected", "unary_unprojected"):
...     observable = occupation_encoding(b_dagger, LevelSpec(3, encoding))
...     terms = observable.to_sparse_list()
...     print(encoding, "->", len(terms), "terms on", observable.num_qubits, "qubits")
...     for label, qubits, coefficient in terms[:4]:
...         print(f"   {label:5s} on qubits {qubits}   coefficient {complex(coefficient):.4f}")
unary_projected -> 12 terms on 4 qubits
   XX00  on qubits [0, 1, 2, 3]   coefficient 0.2500+0.0000j
   XY00  on qubits [0, 1, 2, 3]   coefficient 0.0000-0.2500j
   YX00  on qubits [0, 1, 2, 3]   coefficient 0.0000+0.2500j
   YY00  on qubits [0, 1, 2, 3]   coefficient 0.2500+0.0000j
unary_unprojected -> 12 terms on 4 qubits
   XX    on qubits [0, 1]   coefficient 0.2500+0.0000j
   XY    on qubits [0, 1]   coefficient 0.0000-0.2500j
   YX    on qubits [0, 1]   coefficient 0.0000+0.2500j
   YY    on qubits [0, 1]   coefficient 0.2500+0.0000j

Beide haben zwölf Terme mit denselben Koeffizienten. Die projizierten Terme tragen einen Bit-Term 0 auf den beiden Qubits, die die nicht projizierten Terme nicht erwähnen, sodass jeder projizierte Term alle vier Qubits berührt und jeder nicht projizierte Term zwei. Allgemein berührt ein projizierter Term alle $d$ Qubits des Modus, weil er zwei Faktoren trägt, jeweils $X$ oder $Y$, und einen Projektor auf jedem anderen Qubit.

Die beiden Operatoren sind auf den Zuständen gleich, für die die Kodierung gedacht ist, den One-Hot-Basiszuständen (den Codewörtern), und unterscheiden sich außerhalb davon. Zur Überprüfung nimm einen kleinen Operator auf zwei Modi bei $d = 3$ und vergleiche die Matrizen auf den neun Codewörtern:

>>> hamiltonian = BosonOperator.from_dict({
...     (cre(0), ann(1)): 1.0, (cre(1), ann(0)): 1.0,
...     (cre(0), ann(0)): 0.5, (cre(1), cre(1), ann(1), ann(1)): 0.25,
... })
>>> matrices = {}
>>> for encoding in ("unary_projected", "unary_unprojected"):
...     observable = occupation_encoding(hamiltonian, LevelSpec(2, encoding))
...     matrices[encoding] = SparsePauliOp.from_sparse_observable(observable).to_matrix()
>>> unary = LevelEncoding.UNARY_PROJECTED
>>> codewords = sorted(
...     unary.level_to_bits(first) | (unary.level_to_bits(second) << 3)
...     for first in range(3) for second in range(3)
... )
>>> projected, unprojected = matrices["unary_projected"], matrices["unary_unprojected"]
>>> bool(np.allclose(projected[np.ix_(codewords, codewords)], unprojected[np.ix_(codewords, codewords)]))
True
>>> bool(np.allclose(projected, unprojected))
False

Beide lassen außerdem den Codewort-Unterraum invariant, sodass ein auf einem Codewort präparierter Zustand auf Codewörtern bleibt. Die Wahl spielt nur für Zustände eine Rolle, die diesen Unterraum verlassen, etwa durch Rauschen. Dort unterdrückt der projizierte Operator die Entwicklung, weil seine Projektoren versagen, während der nicht projizierte weiter wirkt. Keines der beiden Verhalten ist eine physikalische Anforderung. Der nicht projizierte Operator ist nur auf dem One-Hot-Unterraum definiert, daher ist seine Wahl eine Modellierungsentscheidung und keine freie Optimierung.

Was die Projektoren kosten, sind Gates. Entwickle ein hermitesches $b^\dagger + b$ auf einem Modus für kurze Zeit und zähle die Zwei-Qubit-Gates nach der Transpilation durch Qiskit:

>>> from qiskit import QuantumCircuit, transpile
>>> from qiskit.circuit.library import PauliEvolutionGate
>>>
>>> def cx_count(observable):
...     circuit = QuantumCircuit(observable.num_qubits)
...     circuit.append(PauliEvolutionGate(observable, time=0.1), range(observable.num_qubits))
...     basis = ["cx", "rz", "sx", "x"]
...     return transpile(circuit, basis_gates=basis, optimization_level=1).count_ops().get("cx", 0)
>>>
>>> hermitian = b + b_dagger
>>> levels = (3, 4, 5, 6)
>>> counts = {
...     str(encoding): [
...         cx_count(occupation_encoding(hermitian, LevelSpec(d - 1, encoding))) for d in levels
...     ]
...     for encoding in encodings
... }

Die genauen Zahlen hängen davon ab, wie Qiskit den Evolution-Gate synthetisiert, was sich zwischen Releases ändert, daher stelle sie grafisch dar:

Die Anzahl ohne Projektoren wächst linear mit $d$, weil jeder der $d - 1$ Übergänge eine Zwei-Qubit-Rotation ist. Die Anzahl mit Projektoren wächst deutlich schneller, weil jeder Term auf alle $d$ Qubits wirkt. Die kompakten Kodierungen brauchen weniger Qubits, aber das Gewicht ihrer Terme hängt davon ab, wie sich die Bitmuster benachbarter Niveaus unterscheiden: Der Gray-Code, in dem sie sich in einem Bit unterscheiden, ist die günstigere der beiden kompakten Kodierungen, und die Standard-Binärkodierung springt, wenn ein Übertrag mehrere Bits auf einmal kippt ($d = 5$). Ob eine kompakte Kodierung auch die nicht projizierte unäre schlägt, hängt von $d$ und davon ab, wie Qiskit die Zeitentwicklung synthetisiert; lies daher den Trend ab und nicht die Rangfolge bei einer einzelnen Größe.

Die Projektoren spielen auch eine Rolle, wenn du einen abgebildeten Operator in Paulis entwickelst, was from_sparse_observable() tut, indem es jeden Projektor durch eine Summe zweier Pauli-Operatoren ersetzt. Ein Term mit $k$ Projektoren wird zu $2^k$ Pauli-Termen. Hier ist ein einzelnes $b^\dagger$ bei kleinem $d$, der einzigen Größe, bei der du das versuchen solltest:

>>> rows = [f"{'d':>2s} {'stored (both)':>14s} {'expanded, projected':>20s} {'expanded, unprojected':>22s}"]
>>> for levels in (3, 4, 5):
...     stored = len(occupation_encoding(b_dagger, LevelSpec(levels - 1, "unary_projected")))
...     expanded = [
...         len(SparsePauliOp.from_sparse_observable(
...             occupation_encoding(b_dagger, LevelSpec(levels - 1, encoding))).simplify())
...         for encoding in ("unary_projected", "unary_unprojected")
...     ]
...     rows.append(f"{levels:2d} {stored:14d} {expanded[0]:20d} {expanded[1]:22d}")
>>> print("\n".join(rows))
 d  stored (both)  expanded, projected  expanded, unprojected
 3              8                   16                      8
 4             12                   48                     12
 5             16                  128                     16

Beide Varianten speichern dieselbe Anzahl von Termen, aber die Entwicklung der projizierten wächst exponentiell in $d$. Lies bei allem, was größer als ein Spielbeispiel ist, die Terme der SparseObservable direkt, wie es die obige Ausgabe tut, statt sie umzuwandeln.

Die kompakten Kodierungen können ihre Projektoren nicht weglassen. Dort sind die Null-Bits Teil der Niveaunummer, sodass das Entfernen eines $\lvert 0 \rangle\langle 0 \rvert$ ändert, welche Niveaus ein Term verbindet, und der Operator auf den Codewörtern nicht mehr derselbe ist. Deshalb gibt es nur die unäre Kodierung in zwei Varianten.

Ein Produkt von Leitern wird zuerst abgeschnitten

Eine weitere Entwurfsentscheidung zeigt sich bei den Kosten diagonaler Terme. Der Mapper bildet nicht jeden Leiteroperator einzeln ab und multipliziert die Ergebnisse. Er sammelt die Leiteroperatoren, die ein Term auf einen Modus anwendet, multipliziert sie als $d \times d$-Matrix in der Fock-Basis, schneidet dieses Produkt einmal ab und kodiert das Ergebnis [2]. Ein Produkt wie der Teilchenzahloperator $b^\dagger b$ bleibt daher diagonal:

>>> number_observable = occupation_encoding(number, LevelSpec(3, "unary_unprojected"))
>>> [(label, qubits) for label, qubits, _ in number_observable.to_sparse_list()]
[('1', [1]), ('1', [2]), ('1', [3])]
>>> cx_count(number_observable)
0

Seine Terme sind die Projektoren $\lvert 1 \rangle\langle 1 \rvert$ auf dem Qubit jedes Niveaus über null, gewichtet mit der Niveaunummer (Niveau null hat Gewicht null und somit keinen Term). Eine Zeitentwicklung davon ist eine Menge von Ein-Qubit-Rotationen und braucht kein Zwei-Qubit-Gate. Multipliziert man stattdessen die beiden abgebildeten Leiteroperatoren, erhält man auf den Codewörtern denselben Operator, außerhalb davon aber einen anderen, mit nichtdiagonalen Zeichenketten, die Gates kosten:

>>> unprojected = LevelSpec(3, "unary_unprojected")
>>> as_pauli = lambda operator: SparsePauliOp.from_sparse_observable(
...     occupation_encoding(operator, unprojected, 4)).simplify()
>>> composed = as_pauli(b_dagger).dot(as_pauli(b)).simplify()
>>> direct = as_pauli(number)
>>> len(direct), len(composed)
(4, 16)
>>> cx_count(direct)
0
>>> bool(cx_count(composed) > 0)
True
>>> codewords = [LevelEncoding.UNARY_UNPROJECTED.level_to_bits(level) for level in range(4)]
>>> bool(np.allclose(direct.to_matrix()[np.ix_(codewords, codewords)],
...                  composed.to_matrix()[np.ix_(codewords, codewords)]))
True
>>> bool(np.allclose(direct.to_matrix(), composed.to_matrix()))
False

Diagonale Terme wie die On-Site-Wechselwirkung des Bose-Hubbard-Modells sind solche Produkte, und dort spielt diese Entscheidung eine Rolle.

Einen eigenen Mapper schreiben

Wenn du eine Kodierung willst, die die Bibliothek nicht bereitstellt, baue sie aus dem Bild eines einzelnen Leiteroperators auf. map_boson_action_generators() nimmt eine Funktion, die eine einzelne BosonAction abbildet, und eine Funktion, die die Identität zurückgibt, und baut daraus den Mapper für ganze Operatoren. Er multipliziert die Bilder der Aktionen jedes Terms der Reihe nach, skaliert mit dem Koeffizienten und summiert die Terme. Die Ausgabe kann ein beliebiger Typ sein, der skalare Multiplikation und Komposition unterstützt.

Nimm als Beispiel die Kodierung für $d = 2$, die einen Modus in einem Qubit speichert. Ein Erzeugungsoperator ist $\lvert 1 \rangle\langle 0 \rvert$, was in der Konvention von Qiskit $(X - iY)/2$ ist, und ein Vernichtungsoperator ist sein Adjungierter:

>>> from qiskit_bosons.mappers import map_boson_action_generators
>>> from qiskit_bosons.operators import BosonAction
>>>
>>> num_qubits = 2
>>>
>>> def map_action(action: BosonAction) -> SparsePauliOp:
...     is_creation, mode = action
...     sign = -1j if is_creation else 1j
...     return SparsePauliOp.from_sparse_list(
...         [("X", [mode], 0.5), ("Y", [mode], 0.5 * sign)], num_qubits=num_qubits
...     )
>>>
>>> def identity() -> SparsePauliOp:
...     return SparsePauliOp.from_sparse_list([("", [], 1.0)], num_qubits)
>>>
>>> model = BosonOperator.from_dict({
...     (cre(0), ann(1)): -1.0, (cre(1), ann(0)): -1.0, (ann(0), cre(0)): 0.7, (cre(1), ann(1)): 0.3,
... })
>>> custom = map_boson_action_generators(model, map_action, identity)
>>> library = SparsePauliOp.from_sparse_observable(occupation_encoding(model, LevelSpec(1, "gray")))
>>> bool(np.allclose(custom.to_matrix(), library.to_matrix()))
True

Dieser handgeschriebene Mapper multipliziert die abgebildeten Leiteroperatoren und stimmt hier mit dem Bibliotheks-Mapper überein, weil jeder Modus zwei Niveaus hat und ein Qubit den ganzen Modus fasst. Bei größerem $d$ sind, wie der vorige Abschnitt gezeigt hat, das Zusammensetzen der abgebildeten Leitern und das vorherige Abschneiden des Produkts außerhalb der Codewörter nicht derselbe Operator.

Um einen solchen Mapper in einem Circuit zu verwenden, verpacke ihn in eine Funktion mit der Signatur (operator, num_qubits) -> SparseObservable und übergib sie an MapperFnEvolutionSynthesis, das damit den Operator jedes Evolution-Gates abbildet. Die Transpiler-Pipeline wird unter Bosonische Circuits transpilieren beschrieben.

>>> from qiskit.quantum_info import SparseObservable
>>> from qiskit_bosons.transpiler.passes import MapperFnEvolutionSynthesis
>>>
>>> def mapper_fn(operator, num_qubits):
...     return SparseObservable.from_sparse_pauli_op(
...         map_boson_action_generators(operator, map_action, identity)
...     )
>>>
>>> synthesis = MapperFnEvolutionSynthesis(mapper_fn)

Zwei Details sind wissenswert. Erstens muss der Qubit-Operator reelle Koeffizienten haben, weil PauliEvolutionGate einen Imaginärteil ablehnt. Ein hermitescher Operator wird auf einen reellen abgebildet, aber erst, nachdem sich aufhebende Terme zusammengeführt wurden. Dieser handgeschriebene Mapper führt sie nicht zusammen, daher hinterlassen ein Hopping-Term und sein Adjungierter imaginäre Koeffizienten, bis du das Ergebnis vereinfachst. Der Bibliotheks-Mapper führt sein eigenes Ergebnis zusammen; ein eigener muss das selbst tun, zum Beispiel indem du ihn in simplify() einwickelst:

>>> hopping = BosonOperator.from_dict({(cre(0), ann(1)): 1.0})
>>> hermitian_hopping = hopping + hopping.adjoint()
>>> raw = mapper_fn(hermitian_hopping, 2)
>>> bool(np.allclose(np.asarray(raw.coeffs).imag, 0))
False
>>> merged = raw.simplify()
>>> bool(np.allclose(np.asarray(merged.coeffs).imag, 0))
True
>>> sorted(label for label, _, _ in merged.to_sparse_list())
['XX', 'YY']

Zweitens synthetisiert eine Produktformel die Terme in der Reihenfolge, in der sie sie erhält, und die Reihenfolge, die ein Mapper erzeugt, ist nicht garantiert stabil. Die Hilfsfunktionen simplify() und group_wise() umhüllen eine Mapper-Funktion, um sie zu steuern: simplify führt doppelte Zeichenketten zusammen und sortiert die Terme in eine kanonische Reihenfolge, während group_wise jeweils eine Gruppe des Operators abbildet, sodass die Terme einer Gruppe benachbart bleiben. Kombiniere die beiden nicht, weil das Sortieren die Gruppierung verwirft. Die Anleitung zur Gruppierung zeigt, wie du Gruppen zuweist, und die Anleitung zur Synthese zeigt, wie die Reihenfolge den Circuit beeinflusst.

Nächste Schritte

• Folge Bose-Hubbard-Dynamik mit wählbarer Kodierung simulieren für eine vollständige Simulation, die die Abschneidung durchfährt und die Kodierungen mit der exakten Dynamik vergleicht.

• Siehe Bosonische Circuits transpilieren, wie dieselbe LevelSpec auf einen bosonischen Circuit angewendet wird.

• Lies Entwurfsprinzipien der Operatordarstellungen, wie ein BosonOperator seine Terme speichert und ordnet, einschließlich der Kontraktionsregel hinter dem obigen Normalordnungsbeispiel.

• Lies Bosonische Modi rotieren für den Gate, der Bosonen zwischen Modi verschiebt und bei dem die Warnung zur erreichbaren Besetzung am wichtigsten ist.

• Konsultiere die API-Referenz von qiskit_bosons.mappers und qiskit_bosons.mappers.library.

[1]

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

[2] (1,2)

N. P. D. Sawaya, T. Menke, T. H. Kyaw, S. Johri, A. Aspuru-Guzik und 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.