steady_state

Graph.steady_state(hamiltonian, lindblad_terms, method='QR', *, name=None)

Find the steady state of a time-independent open quantum system.

The Hamiltonian and Lindblad operators that you provide have to give rise to a unique steady state.

Parameters:

hamiltonian (Tensor or spmatrix) – A 2D array of shape (D, D) representing the time-independent Hamiltonian of the system, \(H_{\rm s}\).
lindblad_terms (list[tuple[float, np.ndarray or Tensor or scipy.sparse.spmatrix]]) – A list of pairs, \((\gamma_j, L_j)\), representing the positive decay rate \(\gamma_j\) and the Lindblad operator \(L_j\) for each coupling channel \(j\). You must provide at least one Lindblad term.
method (str, optional) – The method used to find the steady state. Must be one of “QR”, “EIGEN_SPARSE”, or “EIGEN_DENSE”. The “QR” method obtains the steady state through a QR decomposition and is suitable for small quantum systems with dense representation. The “EIGEN_SPARSE” and “EIGEN_DENSE” methods find the steady state as the eigenvector whose eigenvalue is closest to zero, using either a sparse or a dense representation of the generator. Defaults to “QR”.
name (str or None, optional) – The name of the node.

Returns:

The density matrix representing the steady state of the system.

Return type:

Tensor

Warning

This function currently does not support calculating the gradient with respect to its inputs. Therefore, it cannot be used in a graph for a calculate_optimization or calculate_stochastic_optimization call, which will raise a RuntimeError. Please use gradient-free optimization if you want to perform an optimization task with this function. You can learn more about it in the How to optimize controls using gradient-free optimization user guide.

See also

density_matrix_evolution_pwc: State evolution of an open quantum system.

Notes

Under the Markovian approximation, the dynamics of an open quantum system can be described by the GKS–Lindblad master equation [1] [2],

\[\frac{{\rm d}\rho_{\rm s}(t)}{{\rm d}t} = {\mathcal L} (\rho_{\rm s}(t)) ,\]

where the Lindblad superoperator \({\mathcal L}\) is defined as

\[{\mathcal L} (\rho_{\rm s}(t)) = -i [H_{\rm s}(t), \rho_{\rm s}(t)] + \sum_j \gamma_j {\mathcal D}[L_j] \rho_{\rm s}(t) ,\]

where \({\mathcal D}\) is a superoperator describing the decoherent process in the system evolution and defined as

\[{\mathcal D}[X]\rho := X \rho X^\dagger - \frac{1}{2}\left( X^\dagger X \rho + \rho X^\dagger X \right)\]

for any system operator \(X\).

This function computes the steady state of \({\mathcal L}\) by solving

\[\frac{{\rm d}\rho_{\rm s}(t)}{{\rm d}t} = 0 .\]

The function assumes that \(H_{\rm s}\) is time independent and that the dynamics generated by \({\mathcal L}\) give rise to a unique steady state. That is, the generated quantum dynamical map has to be ergodic [3].

References

Examples

Compute the steady state of the single qubit open system dynamics according to the Hamiltonian \(H=\omega\sigma_z\) and the single Lindblad operator \(L=\sigma_-\).

>>> omega = 0.8
>>> gamma = 0.5
>>> hamiltonian = omega * graph.pauli_matrix("Z")
>>> lindblad_terms = [(gamma, graph.pauli_matrix("M"))]
>>> graph.steady_state(hamiltonian, lindblad_terms, name="steady_state")
<Tensor: name="steady_state", operation_name="steady_state", shape=(2, 2)>
>>> result = qctrl.functions.calculate_graph(graph=graph, output_node_names=["steady_state"])
>>> result.output["steady_state"]["value"]
array([[1.+0.j 0.-0.j]
       [0.-0.j 0.-0.j]])