dq.smesolve
smesolve(
H: ArrayLike  TimeArray,
jump_ops: list[ArrayLike  TimeArray],
etas: ArrayLike,
rho0: ArrayLike,
tsave: ArrayLike,
*,
tmeas: ArrayLike  None = None,
ntrajs: int = 10,
exp_ops: list[ArrayLike]  None = None,
solver: Solver  None = None,
gradient: Gradient  None = None,
options: Options = Options()
) > Result
Solve the diffusive stochastic master equation (SME).
Warning
This function has not been ported to JAX yet. The following documentation is a draft API, copied from the old PyTorch version of the library.
This function computes the evolution of the density matrix \(\rho(t)\) at time \(t\), starting from an initial state \(\rho_0\), according to the diffusive SME in ItΓ΄ form (with \(\hbar=1\) and where time is implicit(1)) $$ \begin{split} \dd\rho =&~ i[H, \rho]\,\dt + \sum_{k=1}^N \left( L_k \rho L_k^\dag  \frac{1}{2} L_k^\dag L_k \rho  \frac{1}{2} \rho L_k^\dag L_k \right)\dt \\ &+ \sum_{k=1}^N \sqrt{\eta_k} \left( L_k \rho + \rho L_k^\dag  \tr{(L_k+L_k^\dag)\rho}\rho \right)\dd W_k, \end{split} $$ where \(H\) is the system's Hamiltonian, \(\{L_k\}\) is a collection of jump operators, each continuously measured with efficiency \(0\leq\eta_k\leq1\) (\(\eta_k=0\) for purely dissipative loss channels) and \(\dd W_k\) are independent Wiener processes.
 With explicit time dependence:
 \(\rho\to\rho(t)\)
 \(H\to H(t)\)
 \(L_k\to L_k(t)\)
 \(\dd W_k\to \dd W_k(t)\)
Diffusive vs. jump SME
In quantum optics the diffusive SME corresponds to homodyne or heterodyne detection schemes, as opposed to the jump SME which corresponds to photon counting schemes. No solver for the jump SME is provided yet, if this is needed don't hesitate to open an issue on GitHub.
The measured signals \(I_k=\dd y_k/\dt\) verifies (again time is implicit): $$ \dd y_k =\sqrt{\eta_k} \tr{(L_k + L_k^\dag) \rho} \dt + \dd W_k. $$
Signal normalisation
Sometimes the signals are defined with a different but equivalent normalisation \(\dd y_k' = \dd y_k/(2\sqrt{\eta_k})\).
The signals \(I_k\) are singular quantities, the solver returns the averaged signals
\(J_k\) defined for a time interval \([t_0, t_1)\) by:
$$
J_k([t_0, t_1)) = \frac{1}{t_1t_0}\int_{t_0}^{t_1} I_k(t)\, \dt
= \frac{1}{t_1t_0}\int_{t_0}^{t_1} \dd y_k(t).
$$
The time intervals for integration are defined by the argument tmeas
, which
defines len(tmeas)  1
intervals. By default, tmeas = tsave
, so the signals
are averaged between the times at which the states are saved.
Defining a timedependent Hamiltonian or jump operator
If the Hamiltonian or the jump operators depend on time, they can be converted
to timearrays using dq.pwc()
,
dq.modulated()
, or
dq.timecallable()
. See the
Timedependent operators
tutorial for more details.
Running multiple simulations concurrently
The Hamiltonian H
, the jump operators jump_ops
and the initial density
matrix rho0
can be batched to solve multiple SMEs concurrently. All other
arguments are common to every batch. See the
Batching simulations
tutorial for
more details.
Parameters

H
(arraylike or timearray of shape (bH?, n, n))
β
Hamiltonian.

jump_ops
(list of arraylike or timearray, of shape (nL, n, n))
β
List of jump operators.

etas
(arraylike of shape (nL,))
β
Measurement efficiencies, must be of the same length as
jump_ops
with values between 0 and 1. For a purely dissipative loss channel, set the corresponding efficiency to 0. No measurement signal will be returned for such channels. 
rho0
(arraylike of shape (brho?, n, 1) or (brho?, n, n))
β
Initial state.

tsave
(arraylike of shape (ntsave,))
β
Times at which the states and expectation values are saved. The equation is solved from
tsave[0]
totsave[1]
, or fromt0
totsave[1]
ift0
is specified inoptions
. 
tmeas
(arraylike of shape (ntmeas,), optional)
β
Times between which measurement signals are averaged and saved. Defaults to
tsave
. 
ntrajs
β
Number of stochastic trajectories to solve concurrently.

exp_ops
(list of arraylike, of shape (nE, n, n), optional)
β
List of operators for which the expectation value is computed.

solver
β
Solver for the integration.

gradient
β
Algorithm used to compute the gradient.

options
β
Generic options, see
dq.Options
.