dq.jssesolve

jssesolve(
    H: QArrayLike | TimeQArray,
    jump_ops: list[QArrayLike | TimeQArray],
    psi0: QArrayLike,
    tsave: ArrayLike,
    keys: PRNGKeyArray,
    *,
    exp_ops: list[QArrayLike] | None = None,
    method: Method = Event(),
    gradient: Gradient | None = None,
    options: Options = Options()
) -> JSSESolveResult

Solve the jump stochastic Schrödinger equation (SSE).

The jump SSE describes the evolution of a quantum system measured by an ideal jump detector (for example photodetection in quantum optics). This function computes the evolution of the state vector $\ket{\psi(t)}$ at time $t$, starting from an initial state $\ket{\psi_0}$, according to the jump SSE ($\hbar=1$, time is implicit(1)) $$ \dd\!\ket\psi = \left[ -iH \dt - \frac12 \sum_{k=1}^N \left( L_k^\dag L_k - \braket{L_k^\dag L_k} \right) \dt + \sum_{k=1}^N \left( \frac{L_k}{\sqrt{\braket{L_k^\dag L_k}}} - 1 \right) \dd N_k \right] \!\ket\psi $$ where $H$ is the system's Hamiltonian, $\{L_k\}$ is a collection of jump operators, each continuously measured with perfect efficiency, and $\dd N_k$ are independent point processes with law $$ \begin{split} \mathbb{P}[\dd N_k = 0] &= 1 - \mathbb{P}[\dd N_k = 1], \\ \mathbb{P}[\dd N_k = 1] &= \braket{L_k^\dag L_k} \dt. \end{split} $$

With explicit time dependence:
- $\ket\psi\to\ket{\psi(t)}$
- $H\to H(t)$
- $L_k\to L_k(t)$
- $\dd N_k\to \dd N_k(t)$

The continuous-time measurements are defined by the point processes $\dd N_k$. The solver returns the times at which the detector clicked, $I_k = \{t \in [t_0, t_\text{end}[ \,|\, \dd N_k(t)=1\}$.

Parameters

H (qarray-like or time-qarray of shape (...H, n, n)) –
Hamiltonian.
jump_ops (list of qarray-like or time-qarray, each of shape (...Lk, n, n)) –
List of jump operators.
psi0 (qarray-like of shape (...psi0, n, 1)) –
Initial state.
tsave (array-like of shape (ntsave,)) –
Times at which the states and expectation values are saved. The equation is solved from tsave[0] to tsave[-1].
keys (list of PRNG keys) –
PRNG keys used to sample the point processes. The number of elements defines the number of sampled stochastic trajectories.
exp_ops (list of array-like, each of shape (n, n), optional) –
List of operators for which the expectation value is computed.
method –
Method for the integration. Defaults to dq.method.Event (supported: Event).
gradient –
Algorithm used to compute the gradient. The default is method-dependent, refer to the documentation of the chosen method for more details.
options –
Generic options (supported: save_states, cartesian_batching, t0, save_extra, nmaxclick).
Detailed options API
```
dq.Options(
    save_states: bool = True,
    cartesian_batching: bool = True,
    t0: ScalarLike | None = None,
    save_extra: callable[[Array], PyTree] | None = None,
    nmaxclick: int = 10_000,
)
```
Parameters
- save_states - If True, the state is saved at every time in tsave, otherwise only the final state is returned.
- cartesian_batching - If True, batched arguments are treated as separated batch dimensions, otherwise the batching is performed over a single shared batched dimension.
- t0 - Initial time. If None, defaults to the first time in tsave.
- save_extra (function, optional) - A function with signature f(QArray) -> PyTree that takes a state as input and returns a PyTree. This can be used to save additional arbitrary data during the integration, accessible in result.extra.
- nmaxclick - Maximum buffer size for result.clicktimes, should be set higher than the expected maximum number of clicks.

Returns

dq.JSSESolveResult object holding the result of the jump SSE integration. Use result.states to access the saved states, result.expects to access the saved expectation values and result.clicktimes to access the detector click times.

Detailed result API

dq.JSSESolveResult

For the shape indications we define ntrajs as the number of trajectories (ntrajs = len(keys)).

Attributes

states (qarray of shape (..., ntrajs, nsave, n, 1)) - Saved states with nsave = ntsave, or nsave = 1 if options.save_states=False.
final_state (qarray of shape (..., ntrajs, n, 1)) - Saved final state.
expects (array of shape (..., ntrajs, len(exp_ops), ntsave) or None) - Saved expectation values, if specified by exp_ops.
clicktimes (array of shape (..., ntrajs, len(jump_ops), nmaxclick)) - Times at which the detectors clicked. Variable-length array padded with None up to nmaxclick.
nclicks (array of shape (..., ntrajs, len(jump_ops)) - Number of clicks for each jump operator.
noclick_states (array of shape (..., nsave, n, 1)) - Saved states for the no-click trajectory. Only for the Event method with smart_sampling=True.
noclick_prob (array of shape (..., nsave)) - Probability of the no-click trajectory. Only for the Event method with smart_sampling=True.
extra (PyTree or None) - Extra data saved with save_extra() if specified in options.
keys (PRNG key array of shape (ntrajs,)) - PRNG keys used to sample the point processes.
infos (PyTree or None) - Method-dependent information on the resolution.
tsave (array of shape (ntsave,)) - Times for which results were saved.
method (Method) - Method used.
gradient (Gradient) - Gradient used.
options (Options) - Options used.

Advanced use-cases

Defining a time-dependent Hamiltonian or jump operator

If the Hamiltonian or the jump operators depend on time, they can be converted to time-arrays using dq.pwc(), dq.modulated(), or dq.timecallable(). See the Time-dependent operators tutorial for more details.

Running multiple simulations concurrently

The Hamiltonian H, the jump operators jump_ops and the initial state psi0 can be batched to solve multiple SSEs concurrently. All other arguments (including the PRNG key) are common to every batch. The resulting states, measurements and expectation values are batched according to the leading dimensions of H, jump_ops and psi0. The behaviour depends on the value of the cartesian_batching option.

If cartesian_batching = True (default value)If cartesian_batching = False

The results leading dimensions are

... = ...H, ...L0, ...L1, (...), ...psi0

For example if:

H has shape (2, 3, n, n),
jump_ops = [L0, L1] has shape [(4, 5, n, n), (6, n, n)],
psi0 has shape (7, n, 1),

then result.states has shape (2, 3, 4, 5, 6, 7, ntrajs, ntsave, n, 1).

The results leading dimensions are

... = ...H = ...L0 = ...L1 = (...) = ...psi0  # (once broadcasted)

For example if:

H has shape (2, 3, n, n),
jump_ops = [L0, L1] has shape [(3, n, n), (2, 1, n, n)],
psi0 has shape (3, n, 1),

then result.states has shape (2, 3, ntrajs, ntsave, n, 1).

See the Batching simulations tutorial for more details.