Representation Comparison

YAQS supports multiple state representations for analog evolution. Each path targets a different scaling regime; the table below summarizes when each is appropriate.

For how to set representation on State, see Initializing Quantum States.

Choosing a representation

Path

When to use

Notes

"mps" (default)

Larger systems and tensor-network-friendly Hamiltonians

TJM trajectories; tune num_traj, max_bond_dim, and accuracy presets

"vector"

MCWF / state-vector quantum trajectories

Exponential memory in qubits; single-trajectory wavefunction dynamics

"density_matrix"

Lindblad master-equation evolution

Exponential memory; deterministic ensemble average without trajectory sampling

Practical guidance:

  • Start with preset="balanced" (or "fast" while exploring) on AnalogSimParams and increase num_traj until observables stabilize.

  • Tighten max_bond_dim / svd_threshold when entanglement growth demands it.

  • For trade-offs between unravellings and trajectory cost, see [1] (references).

The sections below run the same noisy benchmark on all three paths so you can validate agreement on small systems.

1. Noisy open-system benchmark

 1import matplotlib.pyplot as plt
 2
 3from mqt.yaqs import AnalogSimParams, Hamiltonian, NoiseModel, Observable, Simulator, State
 4
 5sim = Simulator(show_progress=False)
 6
 7L = 3
 8H = Hamiltonian.ising(L, J=1.0, g=0.5)
 9noise = NoiseModel([{"name": "pauli_z", "sites": [i], "strength": 0.2} for i in range(L)])
10obs = Observable("x", sites=[0])
11
12t_max = 5.0
13dt = 0.05
14
15print("Running density_matrix...")
16params_rho = AnalogSimParams(observables=[obs], elapsed_time=t_max, dt=dt)
17result_rho = sim.run(State(L, initial="zeros", representation="density_matrix"), H, params_rho, noise)
18res_rho = result_rho.expectation_values[0].flatten()
19times = params_rho.times
20
21print("Running vector...")
22params_vector = AnalogSimParams(observables=[obs], elapsed_time=t_max, dt=dt, num_traj=500)
23result_vector = sim.run(State(L, initial="zeros", representation="vector"), H, params_vector, noise)
24res_vector = result_vector.expectation_values[0].flatten()
25
26print("Running mps...")
27params_mps = AnalogSimParams(observables=[obs], elapsed_time=t_max, dt=dt, num_traj=500, max_bond_dim=16)
28result_mps = sim.run(State(L, initial="zeros", representation="mps"), H, params_mps, noise)
29res_mps = result_mps.expectation_values[0].flatten()
Running density_matrix...
Running vector...
Running mps...
 1fig, ax = plt.subplots(figsize=(6, 3.5))
 2ax.plot(times, res_rho, label="density_matrix (exact)", linewidth=2, color="black")
 3ax.plot(times, res_vector, label="vector (500 traj)", linestyle="--")
 4ax.plot(times, res_mps, label="mps (500 traj)", linestyle=":")
 5ax.set_xlabel("Time")
 6ax.set_ylabel(r"$\langle X_0 \rangle$")
 7ax.legend()
 8ax.set_title("Open-system evolution across representations")
 9ax.grid(alpha=0.3)
10plt.tight_layout()
11plt.show()

Note

vector and mps curves are Monte Carlo means over 500 trajectories; statistical error scales as \(1/\sqrt{N_{\mathrm{traj}}}\). The density_matrix path returns the deterministic ensemble average directly.

2. Noiseless cross-check

With noise_model=None, all three representations should agree on unitary observables (single trajectory for mps and vector):

1obs_z = Observable("z", sites=[0])
2params_mps_u = AnalogSimParams(observables=[obs_z], elapsed_time=1.0, dt=0.1, max_bond_dim=16)
3params_rho_u = AnalogSimParams(observables=[obs_z], elapsed_time=1.0, dt=0.1)
4
5z_mps = sim.run(State(L, initial="zeros", representation="mps"), H, params_mps_u, None).expectation_values[0][-1]
6z_vec = sim.run(State(L, initial="zeros", representation="vector"), H, params_mps_u, None).expectation_values[0][-1]
7z_rho = sim.run(State(L, initial="zeros", representation="density_matrix"), H, params_rho_u, None).expectation_values[0][-1]
8print(f"Noiseless ⟨Z₀⟩ at t=1: mps={z_mps:.6f}, vector={z_vec:.6f}, density_matrix={z_rho:.6f}")
Noiseless ⟨Z₀⟩ at t=1: mps=0.672316, vector=0.672315, density_matrix=0.672315