Deploy and run a Qiskit Function template for AQC + Trotter Hamiltonian dynamics

Overview

This is an experiment-agnostic Qiskit Function template for Hamiltonian dynamics. Given a 1D nearest-neighbor Pauli Hamiltonian, a prepared initial state (optional), and a set of observables, it runs Trotter time-evolution, approximate quantum compilation (AQC) circuit compression, and mitigated execution, then returns each observable's time series. Swap the setup (PRE) and the analysis (POST) and the same core drives a different experiment:

The template is published in the Qiskit Function templates repository, alongside the other application templates. This notebook deploys it to your own Qiskit Serverless account. Run it once, and any notebook can then call the function with serverless.load("aqc-dynamics-function").

For a worked scientific example, see Simulate neutron scattering with an AQC + Trotter dynamics Serverless workflow, which calls this function to compute the dynamical structure factor of KCuF$_3$. This notebook covers deployment and the input contract instead.

Requirements

Before starting, be sure you have the following in this notebook's kernel environment:

• Qiskit SDK v2.0 or later (pip install qiskit).
• The Qiskit IBM Catalog client (pip install qiskit-ibm-catalog), which deploys and runs workloads on Qiskit Serverless.

The function's own scientific dependencies (qiskit-addon-aqc-tensor, cotengrust, qiskit-aer) do not need to be installed locally.

Get the template source files

The function is a small Python package that Qiskit Serverless runs in the cloud, so its source has to exist as local files that are uploaded at deploy time. The package is published in the Qiskit Function templates repository.

Download source_files

The download is a single zip, named after the full path of the directory in the repository:

qiskit-community qiskit-function-templates main physics aqc_trotter source_files.zip

1. Unzip it into the directory that holds this notebook.
2. Rename the extracted folder from that long name to source_files.

Your working directory then looks like this:

your-working-directory/
├── function-template-aqc-trotter.ipynb    <- this notebook
└── source_files/                         <- the renamed folder
    ├── __init__.py
    ├── program.py
    └── source/
        ├── __init__.py
        ├── _serverless.py
        ├── app_function.py
        ├── aqc.py
        ├── build.py
        ├── execute.py
        └── hamiltonian.py

The name has to be exactly source_files, because that is the working_dir Step 3 uploads.

program.py is the entry point the gateway invokes. Everything under source/ is the implementation, split by stage: Hamiltonian and Trotter synthesis, AQC compression, and execution. None of it needs editing to run the examples that follow. Step 3 uploads the whole directory, so repeat that step whenever you change a file.

1. Authentication

Use qiskit-ibm-catalog to authenticate to QiskitServerless with your API key (token) and CRN (instance), which you can find on the IBM Quantum® Platform dashboard. With these credentials you can instantiate the serverless client locally to upload or run the selected function:

from qiskit_ibm_catalog import QiskitServerless
serverless = QiskitServerless(channel="ibm_quantum_platform", token="MY_TOKEN", instance="MY_CRN")

You can optionally use save_account() to save your credentials in your local environment (see the Set up your IBM Cloud® account guide). Note that this writes your credentials to the same file as QiskitRuntimeService.save_account():

QiskitServerless.save_account(channel="ibm_quantum_platform", token="MY_TOKEN", instance="MY_CRN")

If the account is saved, there is no need to provide the token to authenticate:

from qiskit_ibm_catalog import QiskitServerless

# Authenticate to the remote cluster
# In this case, loading a saved account
serverless = QiskitServerless()

# REPLACE WITH YOUR OWN CREDENTIALS or SAVED ACCOUNT
# serverless = QiskitServerless(channel="ibm_quantum_platform", token="MY_TOKEN", instance="MY_CRN")

2. Declare dependencies

Packages the function needs on top of the managed base serverless image.

DEPENDENCIES = [
    "qiskit-addon-aqc-tensor[quimb-jax]==0.3.1",
    "qiskit-aer==0.17.2",
    "cotengrust==0.2.0",
]

3. Define and upload the function

from qiskit_ibm_catalog import QiskitFunction

fn = QiskitFunction(
    title="aqc-dynamics-function",
    entrypoint="program.py",
    working_dir="source_files/",
    dependencies=DEPENDENCIES,
)
serverless.upload(fn)

Output:

QiskitFunction(aqc-dynamics-function)

4. Verify it registered

next(p for p in serverless.list() if p.title == "aqc-dynamics-function")

Output:

QiskitFunction(aqc-dynamics-function)

Function reference

This is a brief introduction. Every field is documented in full in the AQC Dynamics Template README: the complete inputs table with its validation rules, the output fields, the execution backends, and further worked examples. What follows is the short version, enough to read the examples that follow.

Inputs

Every run is a single fn.run(...) call. Only the first three inputs in the table are required: hamiltonian, t_steps, and aqc_segments. Everything after them is optional and falls back to the default shown, so a minimal call is three arguments and the rest of the table is the functionality you can opt into. The Hamiltonian's num_qubits sets the chain length, so there is no separate size input.

Execution backends

All three paths share the same code and the same mitigation settings. They differ only in where the circuits run.

Both simulator paths still call the deployed function, so they need a saved Serverless account even though they use no QPU time. The two examples that follow run the same workload on statevector first, then on runtime.

Output

job.result() returns a plain dictionary:

{
    "times": [...],                  # length t_steps + 1, t_k = k * dt (t=0 is the prepared state)
    "expectation_values": [[...]],   # shape (n_times, n_observables)
    "observable_labels": [...],      # for example: ["Z_0", "ZZ_0_1"]
    "metadata": {
        "n", "t_steps", "dt", "tier",
        "aqc_compressed_steps": 5,   # total compressed steps (= sum of segment n_steps)
        "aqc_segments": [            # per segment: the plan plus its own results
            {"n_steps": 3, "ansatz_steps": 1, "steps": [1, 2, 3], "n_params": 133,
             "fidelities": {"1": ..., "2": ..., "3": ...}},
            {"n_steps": 2, "ansatz_steps": 2, "steps": [4, 5], "n_params": 245,
             "fidelities": {"4": ..., "5": ...}},
        ],
        "execution_backend",
        "aqc_fidelities": {"1": ..., "2": ...},  # flat per-step fidelity, all compressed steps
        "circuit_stats": {           # per-step 2q depth and gate count, full Trotter vs AQC
            "1": {"full_trotter": {"depth_2q": ..., "num_2q_gates": ...},
                "aqc_trotter":  {"depth_2q": ..., "num_2q_gates": ...}},
            "2": {...},
        },
        "warnings": [...],           # non-fatal notices; for example, a cotengrust fallback
        "resource_usage": {          # per stage; QPU_TIME is the charged QPU time
            "RUNNING: OPTIMIZING_FOR_HARDWARE": {"CPU_TIME": ...},
            "RUNNING: WAITING_FOR_QPU": {"CPU_TIME": ...},
            "RUNNING: EXECUTING_QPU": {"QPU_TIME": ...},
        },
    },
}

aqc_fidelities and circuit_stats are the two to read first: together they tell you whether the compression stayed faithful and whether it actually saved depth. On runtime, resource_usage reports the queue wait separately from the QPU time you are charged for. A rejected input fails fast as a structured ServerlessError (code 4615).

Simulator example

Run the function on the exact statevector backend first. It spends no QPU time and validates the deployment end to end. The model here is an eight-qubit transverse-field Ising chain, and observables is omitted so the function measures the default per-site $Z$.

The compression plan is the input worth understanding. Each segment {"n_steps": k, "ansatz_steps": m} compresses k consecutive Trotter steps into an ansatz built from an m-step Trotter target, and any steps beyond sum(n_steps) run as plain Trotter. Early, low-entanglement steps compress well into a shallow single-layer ansatz; later, more-entangled steps need a deeper one.

from qiskit.quantum_info import SparsePauliOp

fn = serverless.load("aqc-dynamics-function")

n = 8
H = SparsePauliOp.from_sparse_list(
    [("ZZ", [i, i + 1], 1.0) for i in range(n - 1)]
    + [("X", [i], 0.8) for i in range(n)],
    num_qubits=n,
)

job = fn.run(
    t_steps=8,
    aqc_segments=[
        {
            "n_steps": 4,
            "ansatz_steps": 1,
        },  # early steps -> shallow 1-layer ansatz
        {
            "n_steps": 2,
            "ansatz_steps": 2,
        },  # later steps -> deeper 2-layer ansatz
    ],
    hamiltonian=H,
    aqc_options={"max_bond": 32},
    backend="statevector",
)
print("job ID:", job.job_id)

Output:

job ID: ee1f3793-e995-427d-81d1-5924549beb38

Follow the run and read the result

status() reports both the coarse job lifecycle and the per-stage sub-status the function publishes as it runs. The same stages apply to the hardware run later in this guide:

QUEUED -> INITIALIZING -> RUNNING: OPTIMIZING_FOR_HARDWARE -> RUNNING: WAITING_FOR_QPU -> RUNNING: EXECUTING_QPU -> RUNNING: POST_PROCESSING -> DONE

Terminal states are DONE, ERROR, and CANCELED. This statevector run has no QPU queue, so it skips RUNNING: WAITING_FOR_QPU. Use job.logs() at any point to see the per-stage logs, including the AQC fidelity reached at each step.

print(job.status())  # re-run until this reports DONE

Output:

DONE

import numpy as np

result = job.result()
ev = np.array(result["expectation_values"])

print("observables:", result["observable_labels"])
print("shape:", ev.shape, "-> (n_times, n_observables)")
print("first row (t = 0, the prepared state):", np.round(ev[0], 4))
print("last row (t = t_steps * dt):", np.round(ev[-1], 4))
print(
    "AQC fidelities:",
    {k: round(v, 4) for k, v in result["metadata"]["aqc_fidelities"].items()},
)

# What the compression bought: 2-qubit depth at the final time step.
stats = result["metadata"]["circuit_stats"][
    str(result["metadata"]["t_steps"])
]
print(
    "2q depth at the final step:",
    stats["full_trotter"]["depth_2q"],
    "(full Trotter) ->",
    stats["aqc_trotter"]["depth_2q"],
    "(AQC + Trotter)",
)

Output:

observables: ['Z_0', 'Z_1', 'Z_2', 'Z_3', 'Z_4', 'Z_5', 'Z_6', 'Z_7']
shape: (9, 8) -> (n_times, n_observables)
first row (t = 0, the prepared state): [1. 1. 1. 1. 1. 1. 1. 1.]
last row (t = t_steps * dt): [0.1442 0.2956 0.4686 0.4877 0.4869 0.4686 0.2963 0.1441]
AQC fidelities: {'1': 1.0, '2': 1.0, '3': 1.0, '4': 1.0, '5': 1.0, '6': 0.9999}
2q depth at the final step: 210 (full Trotter) -> 79 (AQC + Trotter)

Hardware example

A function call with backend="runtime" transpiles and executes on a real IBM Quantum processor, with the function's built-in error mitigation: dynamical decoupling (XY4), gate twirling, and twirled readout error extinction (TREX). backend_name selects the device; omit it and the function takes the least busy one.

Nothing about the science code changes. What differs from the simulator example is the chain length, the number of Trotter steps, the compression plan, the backend, and the explicit mitigation settings covered in the following section.

Sizing the job for the control hardware

estimator_options is the input worth setting deliberately. Gate twirling builds num_randomizations separate randomized circuits for every PUB, and the whole job, every PUB with all of its randomizations, has to fit in the instruction memory of the QPU's classical control system. The function defaults to 1000 randomizations, so a 10-step evolution submits 11 PUBs of 1000 circuits each: roughly 11,000 circuit instances in a single job.

Exceed what the control system holds and the job fails with error 6073. Job limits gives the thresholds and how to count against them, the main one being 26.8 million control-system instructions per qubit, applied per job rather than per PUB. Dynamical decoupling adds gates that count toward it.

Two inputs control the size:

• estimator_options sets the shot budget. Total shots is num_randomizations * shots_per_randomization, so you can trade randomizations against shots per randomization, keep the statistics, and still shrink the program. The following cell uses 100 randomizations at 200 shots each, which is 20,000 shots per observable and about a tenth of the circuit instances the defaults would submit. See TwirlingOptions and Estimator options for the full set of fields.
• batches splits the PUBs across that many separate runtime jobs, which is the remedy error 6073 itself suggests and why the per-job framing matters. Setting batches=4 sends roughly three PUBs per job instead of eleven at once, and the jobs go out together in one batch so the group queues once rather than each job queueing separately.

Remember that a supplied estimator_options replaces the function's defaults wholesale rather than merging into them, so dynamical decoupling and TREX are restated in the following cell to keep them switched on.

from qiskit.quantum_info import SparsePauliOp

fn = serverless.load("aqc-dynamics-function")

n = 10
H = SparsePauliOp.from_sparse_list(
    [("ZZ", [i, i + 1], 1.0) for i in range(n - 1)]
    + [("X", [i], 0.8) for i in range(n)],
    num_qubits=n,
)

job = fn.run(
    t_steps=10,
    aqc_segments=[
        {"n_steps": 3, "ansatz_steps": 1},
        {"n_steps": 3, "ansatz_steps": 2},
    ],
    hamiltonian=H,
    aqc_options={"max_bond": 32},
    backend="runtime",
    backend_name="ibm_marrakesh",
    # The function defaults to 1000 twirling randomizations, which was too large
    # for this device. Total shots is num_randomizations *
    # shots_per_randomization, so this is 20,000 shots per observable.
    estimator_options={
        "dynamical_decoupling": {"enable": True, "sequence_type": "XY4"},
        "twirling": {
            "enable_gates": True,
            "num_randomizations": 100,
            "shots_per_randomization": 200,
        },
        "resilience": {"measure_mitigation": True},
    },
)
print("job ID (save this to reconnect later):", job.job_id)

Output:

job ID (save this to reconnect later): 7229a8bf-9f83-4785-8dd4-489844abc2d9

# Reconnect to a previously submitted job by its ID. Only needed in a NEW kernel
# session; if you are still in the session where you submitted, the `job` handle
# from the preceding cell is already live, so skip this cell. Replace the ID that follows with your own.
job = serverless.get_job_by_id("<your job ID>")

# Re-run this until it reports DONE, then fetch the result in the following cell.
print(job.status())

Output:

DONE

# Run this only once the preceding status cell reports DONE. result() blocks until
# the job finishes, so calling it earlier just waits.
result = job.result()
ev = np.array(result["expectation_values"])

print("backend:", result["metadata"]["execution_backend"])
print("shape:", ev.shape, "-> (n_times, n_observables)")
print("last row (t = t_steps * dt):", np.round(ev[-1], 4))
print(
    "AQC fidelities:",
    {k: round(v, 4) for k, v in result["metadata"]["aqc_fidelities"].items()},
)

# What the compression bought: 2-qubit depth at the final time step.
stats = result["metadata"]["circuit_stats"][
    str(result["metadata"]["t_steps"])
]
print(
    "2q depth at the final step:",
    stats["full_trotter"]["depth_2q"],
    "(full Trotter) ->",
    stats["aqc_trotter"]["depth_2q"],
    "(AQC + Trotter)",
)

Output:

backend: runtime
shape: (11, 10) -> (n_times, n_observables)
last row (t = t_steps * dt): [0.1504 0.1361 0.218  0.2144 0.2275 0.1783 0.1749 0.1599 0.0915 0.0922]
AQC fidelities: {'1': 1.0, '2': 1.0, '3': 1.0, '4': 1.0, '5': 0.9999, '6': 0.9999}
2q depth at the final step: 342 (full Trotter) -> 171 (AQC + Trotter)

Next steps