Синтез гейта Evolution у бозонному просторі

Часову еволюцію $e^{-itH}$ гамільтоніана, який є сумою некомутівних частин, не можна реалізувати точно гейт за гейтом. Формула добутку наближає її послідовністю коротких еволюцій під дією окремих частин, а похибка наближення називається похибкою Троттера.

Формулу добутку можна застосувати у двох місцях, і вони не взаємозамінні:

1. У бозонному просторі, до того як обрано будь-яке обрізання чи кодування. Формула бачить групи оператора (див. Групування членів оператора: використання структури оператора) і може їх чергувати.

2. У кубітному просторі, всередині PauliEvolutionGate з Qiskit, після відображення оператора. Формула бачить один оператор Паулі на групу і може лише розбивати члени всередині нього.

Лише перший спосіб може зменшити похибку між групами, що не комутують, а вона зазвичай домінує. Цей посібник показує, як ним користуватися.

Важливо

Обидва кроки необов'язкові. Гейт Evolution, який ніколи не декомпонується в бозонному просторі, потрапляє на етап «бозон у кубіт» цілим, і там відображається та синтезується. Це не вносить жодної похибки на бозонному рівні, лише ту, що вносить формула кубітного рівня.

Налаштування

Двовузловий ланцюжок Бозе-Хаббарда — найменший приклад, де групи не комутують (повну модель див. у Симуляція динаміки Бозе-Хаббарда з вибором кодування). Його згруповано в одну ермітову групу для зв'язку та одну діагональну групу для внутрішньовузлової взаємодії:

>>> import numpy as np
>>> from qiskit.quantum_info import SparsePauliOp, Statevector
>>> from qiskit.transpiler import PassManager
>>> from scipy.sparse import csc_matrix
>>> from scipy.sparse.linalg import expm_multiply
>>> from qiskit_bosons.circuit import BosonicCircuit
>>> from qiskit_bosons.circuit.library import Evolution, InitializeModes
>>> from qiskit_bosons.mappers import LevelSpec
>>> from qiskit_bosons.mappers.library import occupation_encoding
>>> from qiskit_bosons.operators import BosonOperator, ann, cre
>>> from qiskit_bosons.transpiler.presets import generate_preset_occupation_pass_manager
>>>
>>> num_sites, max_occupation, total_time = 2, 2, 0.8
>>> bonds = [(0, 1)]
>>> terms = {}
>>> for left, right in bonds:
...     terms[(cre(left), ann(right))] = -1.0
...     terms[(cre(right), ann(left))] = -1.0
>>> for site in range(num_sites):
...     terms[(cre(site), cre(site), ann(site), ann(site))] = 1.0
>>> hamiltonian = BosonOperator.from_dict(terms)
>>> hamiltonian.groups = [
...     len(bonds)
...     if len({mode for _, mode in term}) == 1
...     else bonds.index(tuple(sorted({mode for _, mode in term})))
...     for term, _ in hamiltonian.iter_terms()
... ]
>>> hamiltonian.num_groups()
2

Щоб виміряти похибку, потрібен еталон. Правильним еталоном є точна еволюція гамільтоніана у відображенні з тим самим обрізанням, адже обрізання — окреме наближення:

>>> spec = LevelSpec(max_occupation, "unary_projected")
>>> mapped = occupation_encoding(hamiltonian, spec).simplify()
>>> matrix = csc_matrix(SparsePauliOp.from_sparse_observable(mapped).to_matrix(sparse=True))
>>> pm = generate_preset_occupation_pass_manager(spec)
>>> pm.qubit = PassManager()  # keep each group as one exact PauliEvolutionGate
>>>
>>> def build(**evolution_kwargs):
...     circuit = BosonicCircuit(num_sites)
...     circuit.append(InitializeModes([1, 1]), circuit.modes)
...     circuit.append(
...         Evolution(num_sites, hamiltonian, time=total_time, **evolution_kwargs), circuit.modes
...     )
...     return circuit
>>>
>>> # Under the unary encodings, occupation ``n`` of mode ``m`` sets qubit ``m * d + n``, with
>>> # ``d = max_occupation + 1`` levels per mode. That places |1, 1> on the basis state below.
>>> levels = max_occupation + 1
>>> index = sum(1 << (mode * levels + n) for mode, n in enumerate([1, 1]))
>>> initial = Statevector.from_int(index, 2**(num_sites * levels)).data
>>> exact = expm_multiply(-1j * total_time * matrix, initial)

>>> def error(circuit):
...     """Distance to the exact evolution, of the circuit's product of group exponentials."""
...     lowered = pm.run(circuit)
...     state = Statevector.from_int(0, 2**lowered.num_qubits).evolve(lowered).data
...     return float(np.max(np.abs(state - exact)))

Примітка

error навмисно симулює знижене коло без декомпозиції його екземплярів PauliEvolutionGate. Statevector застосовує кожен як точну експоненту його групи, тож залишається лише добуток експонент груп, який обрав бозонний синтез.

Типова поведінка

Без аргументу synthesis гейт Evolution декомпонується за допомогою BosonicLieTrotter, формули першого порядку: один множник на групу, кожен на весь час.

>>> dict(build().decompose().count_ops())
{'Evolution': 2, 'InitializeModes': 1}

Без декомпозиції весь гейт передається на кубітний етап, і бозонної похибки Троттера немає взагалі:

>>> error(build()) < 1e-6
True

Після декомпозиції з'являється справжня похибка, бо дві групи не комутують:

>>> error(build().decompose()) > 0.1
True

Зменшення похибки

Передай гейту інший метод синтезу. Вищий порядок компонує множники симетрично, так що провідні члени похибки скорочуються; більше повторень ділять час на коротші кроки.

>>> from qiskit_bosons.circuit.library.synthesis import BosonicLieTrotter, BosonicSuzukiTrotter
>>>
>>> for label, synthesis in [
...     ("Lie-Trotter", BosonicLieTrotter()),
...     ("Suzuki, order 2", BosonicSuzukiTrotter(order=2)),
...     ("Suzuki, order 4", BosonicSuzukiTrotter(order=4)),
...     ("Suzuki, order 4*", BosonicSuzukiTrotter(order=4, merge_adjacent=False)),
... ]:
...     circuit = build(synthesis=synthesis).decompose()
...     print(f"{label:16s} {circuit.count_ops()['Evolution']:3d} factors")
Lie-Trotter        2 factors
Suzuki, order 2    3 factors
Suzuki, order 4   11 factors
Suzuki, order 4*  15 factors

Кожне збільшення порядку на два компонує п'ять проходів попереднього порядку (зростання $5^{(k-2)/2}$ у рекурсії Сузукі, застосованій до $2n - 1$ множників порядку 2), і саме це число наведено в останньому рядку, позначеному зірочкою. Проходи сходяться на одному множнику, тож сусідні гейти еволюціонують під дією одного оператора і зливаються в один гейт із сумарним часом. Це заощаджує один гейт на шов і є точним. Передай merge_adjacent=False, щоб зберегти підручникову послідовність. Так чи інакше, вищий порядок дає меншу похибку:

>>> lie = error(build(synthesis=BosonicLieTrotter()).decompose())
>>> second = error(build(synthesis=BosonicSuzukiTrotter(order=2)).decompose())
>>> fourth = error(build(synthesis=BosonicSuzukiTrotter(order=4)).decompose())
>>> lie > second > fourth
True

Повторення обмінюють глибину на точність за іншою ставкою. Похибка першого порядку спадає лінійно за розміром кроку, другого — квадратично:

>>> lie_errors = [error(build(synthesis=BosonicLieTrotter(reps=r)).decompose()) for r in (1, 2, 4, 8)]
>>> bool(all(a > b for a, b in zip(lie_errors, lie_errors[1:])))
True
>>> bool(lie_errors[0] / lie_errors[-1] > 6)  # about 8, i.e. linear in 1 / reps
True
>>> suzuki_errors = [
...     error(build(synthesis=BosonicSuzukiTrotter(order=2, reps=r)).decompose()) for r in (1, 2, 4)
... ]
>>> bool(suzuki_errors[0] / suzuki_errors[-1] > 10)  # about 16, i.e. quadratic in 1 / reps
True

Яка угода краща, залежить від оператора. Зауваж, що оператор з однією групою точний за будь-якого порядку, бо єдиний множник комутує сам із собою, і тоді формула згортається в один гейт, а не розбиває його знову.

Вимірювання вартості

Менша похибка коштує більше множників, а кожен множник коштує заплутувальних гейтів після того, як його синтезував кубітний етап. Скільки саме, залежить від кодування, тож порахуй їх через пайплайн пресета з фіксованим базисом і без оптимізації на кубітному рівні (optimization_level=0), щоб побачити ціну кожного методу:

>>> basis = ["rz", "sx", "x", "cx"]
>>>
>>> def cx_count(circuit, encoding):
...     spec = LevelSpec(max_occupation, encoding)
...     lowered = generate_preset_occupation_pass_manager(
...         spec, optimization_level=0, basis_gates=basis
...     ).run(circuit)
...     return lowered.count_ops().get("cx", 0)
>>>
>>> methods = [
...     ("Lie-Trotter", BosonicLieTrotter()),
...     ("Suzuki, order 2", BosonicSuzukiTrotter(order=2)),
...     ("Suzuki, order 4", BosonicSuzukiTrotter(order=4)),
... ]
>>> encodings = ("gray", "unary_unprojected")
>>> counts = {
...     encoding: [
...         cx_count(build(synthesis=method).decompose(), encoding) for _, method in methods
...     ]
...     for encoding in encodings
... }
>>> # the order buys accuracy at a price that grows with it
>>> bool(counts["gray"][0] < counts["gray"][1] < counts["gray"][2])
True

Точні кількості залежать від того, як Qiskit синтезує гейти еволюції, а це змінюється між випусками, тому графік показує їх замість таблиці:

Кожна група стовпчиків — це метод: Лі-Троттер, Сузукі порядку 2 і Сузукі порядку 4. Порядок купує точність за ціну, що зростає з множником рекурсії, а унарне кодування платить більше за код Грея для кожного методу, бо кожен його множник відображається на більше заплутувальних гейтів. Обирай метод за потрібною точністю, а кодування — під апаратуру, і тримай ці два вибори окремо: бозонне коло вище ніколи не змінюється.

Застосування одного методу до всього кола

Щоб задати synthesis для кожного гейта, треба протягнути вибір крізь усе, що його будує. Прохід B2BSynthesis робить вибір один раз для всього пайплайна. Його config називає метод для кожного типу гейта, і він одночасно обирає метод та розгортає гейти:

>>> from qiskit_bosons.transpiler import BosonicPassManager
>>> from qiskit_bosons.transpiler.passes import B2BSynthesis
>>>
>>> trotterize = BosonicPassManager([B2BSynthesis({"Evolution": BosonicSuzukiTrotter(order=2)})])
>>> dict(trotterize.run(build()).count_ops())
{'Evolution': 3, 'InitializeModes': 1}

Метод також можна назвати, із аргументами для його __init__, замість того щоб будувати: B2BSynthesis({"Evolution": ("SuzukiTrotter", (), {"order": 2})}). Доступні імена — це імена з B2BSynthesisPluginManager, а метод, не зареєстрований ніде, можна передати як екземпляр, як вище.

З apply=False прохід лише обирає метод і залишає розгортання до моменту, коли воно потрібне, наприклад щоб спершу оглянути коло. Атомарний гейт не чіпається ніколи, навіть не позначається:

>>> selected = BosonicPassManager(
...     [B2BSynthesis({"Evolution": ("SuzukiTrotter", (), {"order": 2})}, apply=False)]
... ).run(build())
>>> selected.count_ops()["Evolution"]  # selected, not yet expanded
1
>>> selected.decompose().count_ops()["Evolution"]  # expanded on demand
3

Нерозгорнутий гейт потрапляє на етап «бозон у кубіт» цілим, де він відображається без звернення до його synthesis, тож усі методи тоді дають однаковий результат.

Необов'язковий filter обмежує прохід гейтами, яким це корисно. Еволюції, групи яких усі комутують, як-от гамільтоніан лише з операторів числа, синтезуються точно за будь-якого порядку, тож вищий порядок лише додав би там глибини.

Обмеження проходу гейтами, яким це корисно

Суто діагональний оператор унаочнює це, бо оператори числа комутують між собою. Його групи точні вже на першому порядку, а другий порядок марно витрачає заплутувальні гейти:

>>> # one group per term, each Hermitian; the order is fixed because order two evolves the first
>>> # group twice, so which group comes first changes the cost
>>> diagonal = BosonOperator.from_terms_with_groups(
...     [([cre(0), ann(0), cre(1), ann(1)], 1.0, 0), ([cre(0), ann(0)], 0.5, 1)]
... )
>>>
>>> def diagonal_cost(synthesis):
...     circuit = BosonicCircuit(num_sites)
...     circuit.append(Evolution(num_sites, diagonal, time=total_time, synthesis=synthesis), circuit.modes)
...     return cx_count(circuit.decompose(), "gray")
>>>
>>> bool(diagonal_cost(BosonicSuzukiTrotter(order=2)) > diagonal_cost(BosonicLieTrotter()))
True

Обидва кола точні, тож друге витрачає більше заплутувальних гейтів без користі. (Симетрична формула другого порядку застосовує першу групу двічі з половинним часом, тож вартість залежить від того, яка група йде першою.) filter приймає предикат над DAGOpNode і залишає відхилені гейти недоторканими, тож один прохід може підвищити порядок для гейтів, яким це корисно, а решту лишити на першому порядку. Предикат нижче приймає гейт лише тоді, коли його оператор переносить бозони між модами, тобто коли групи не комутують:

>>> def is_diagonal(operator):
...     """Whether every term only counts bosons, rather than moving them between modes."""
...     return all(
...         sorted(mode for created, mode in actions if created)
...         == sorted(mode for created, mode in actions if not created)
...         for actions, _ in operator.iter_terms()
...     )
>>>
>>> def mixed_circuit():
...     circuit = BosonicCircuit(num_sites)
...     circuit.append(Evolution(num_sites, hamiltonian, time=total_time), circuit.modes)
...     circuit.append(Evolution(num_sites, diagonal, time=total_time), circuit.modes)
...     return circuit
>>>
>>> order_2 = BosonicSuzukiTrotter(order=2)
>>> everywhere = cx_count(
...     BosonicPassManager([B2BSynthesis({"Evolution": order_2})]).run(mixed_circuit()), "gray"
... )
>>> restricted = cx_count(
...     BosonicPassManager(
...         [
...             B2BSynthesis(
...                 {"Evolution": order_2},
...                 filter=lambda node: not is_diagonal(node.op.operator),
...             )
...         ]
...     ).run(mixed_circuit()),
...     "gray",
... )
>>> bool(everywhere > restricted)  # the filter saves entangling gates
True

Якщо запустити прохід знову, нічого не зміниться: множники, які він випускає, є атомарними.

Атомарні множники

Evolution, позначений як atomic, не має визначення, тож decompose() залишає його на місці, скільки б разів його не застосовували. Синтез позначає так множники, які випускає, і так само робить BosonicQDrift для членів, які вибирає.

Це більше, ніж просто охайність. Подальше розбиття ермітової групи по членах дає множники, які окремо не є ермітовими, а експонента неермітового оператора не є унітарною:

>>> circuit = build()
>>> once = circuit.decompose()
>>> dict(circuit.decompose(reps=5).count_ops()) == dict(once.count_ops())
True

Обережно

Кожен множник сам по собі має бути ермітовим. Ермітів оператор, розбитий на групи, може дати неермітові, наприклад член стрибка, відокремлений від спряженого. Перевір за допомогою groups_are_hermitian(). Оператор без груп розбивається по членах, і тоді кубітний етап відхиляє множники з Operator contains complex coefficients. Симетризація парного порядку частково скорочує похибку неермітового множника, тож неправильно згрупований оператор може виглядати краще на порядку 2, ніж на порядку 1, залишаючись хибним.

Повторення короткого кроку

Як альтернатива на рівні кола, побудуй один короткий крок Троттера і повтори його. repeat() розгортає повторення в одне коло і залишає межі суміжними, щоб проходи могли діяти крізь них:

>>> step = BosonicCircuit(num_sites)
>>> step.append(Evolution(num_sites, hamiltonian, time=total_time / 4), step.modes)
>>> dict(step.repeat(4).count_ops())
{'Evolution': 4}

Передай insert_barriers=True, щоб розділити повторення. Бар'єр відображається на кожен кубіт регістрів, яких торкається, бо кодуванню потрібно кілька кубітів на моду, а мода не має власного покубітного образу:

>>> [inst.operation.name for inst in step.repeat(3, insert_barriers=True)._inner.data]
['Evolution', 'barrier', 'Evolution', 'barrier', 'Evolution']

Примітка

repeat не пов'язаний із reps методу синтезу. Синтез ділить одну еволюцію на більше коротших множників. Повторення кола, що містить одну еволюцію, лише множить її повний час.

Поєднання з формулою кубітного рівня

Ці два рівні компонуються, і де обидва наближають, точність результату визначається слабшим із них. Підвищення бозонного порядку, коли формула з боку кубітів лишається першого порядку, дає мало. Про формулу з боку кубітів див. MapperFnEvolutionSynthesis. Вона зберігає порядок членів Паулі, який обрав відображувач, тож group_wise() може випускати члени кожної групи поруч, а simplify() може натомість зафіксувати канонічний порядок.

Наступні кроки

• Транспіляція бозонних кіл пояснює етапи, що йдуть після бозонного синтезу, і де серед них стоїть B2BSynthesis.

• Групування членів оператора: використання структури оператора пояснює, як призначати й перевіряти групи, які чергує формула добутку.

• Симуляція динаміки Бозе-Хаббарда з вибором кодування застосовує схему повторень до повного прикладу.