quimb.evo ========= .. py:module:: quimb.evo .. autoapi-nested-parse:: Easy and efficient time evolutions. Contains an evolution class, Evolution to easily and efficiently manage time evolution of quantum states according to the Schrodinger equation, and related functions. Attributes ---------- .. autoapisummary:: quimb.evo.CALLABLE_TIME_INDEP_CLASSES Classes ------- .. autoapisummary:: quimb.evo.Try2Then3Args quimb.evo.Evolution Functions --------- .. autoapisummary:: quimb.evo.schrodinger_eq_ket quimb.evo.schrodinger_eq_ket_timedep quimb.evo.schrodinger_eq_dop quimb.evo.schrodinger_eq_dop_timedep quimb.evo.schrodinger_eq_dop_vectorized quimb.evo.lindblad_eq quimb.evo.lindblad_eq_vectorized quimb.evo._calc_evo_eq Module Contents --------------- .. py:data:: CALLABLE_TIME_INDEP_CLASSES .. py:function:: schrodinger_eq_ket(ham) Wavefunction schrodinger equation. :param ham: Time-independant Hamiltonian governing evolution. :type ham: operator :returns: **psi_dot(t, y)** -- Function to calculate psi_dot(t) at psi(t). :rtype: callable .. py:function:: schrodinger_eq_ket_timedep(ham) Wavefunction time dependent schrodinger equation. :param ham: Time-dependant Hamiltonian governing evolution, such that ``ham(t)`` returns an operator representation of the Hamiltonian at time ``t``. :type ham: callable :returns: **psi_dot(t, y)** -- Function to calculate psi_dot(t) at psi(t). :rtype: callable .. py:function:: schrodinger_eq_dop(ham) Density operator schrodinger equation, but with flattened input/output. Note that this assumes both `ham` and `rho` are hermitian in order to speed up the commutator, non-hermitian hamiltonians as used to model loss should be treated explicilty or with `schrodinger_eq_dop_vectorized`. :param ham: Time-independant Hamiltonian governing evolution. :type ham: operator :returns: **rho_dot(t, y)** -- Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form). :rtype: callable .. py:function:: schrodinger_eq_dop_timedep(ham) Time dependent density operator schrodinger equation, but with flattened input/output. Note that this assumes both `ham(t)` and `rho` are hermitian in order to speed up the commutator, non-hermitian hamiltonians as used to model loss should be treated explicilty or with `schrodinger_eq_dop_vectorized`. :param ham: Time-dependant Hamiltonian governing evolution, such that ``ham(t)`` returns an operator representation of the Hamiltonian at time ``t``. :type ham: callable :returns: **rho_dot(t, y)** -- Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form). :rtype: callable .. py:function:: schrodinger_eq_dop_vectorized(ham) Density operator schrodinger equation, but with flattened input/output and vectorised superoperator mode (no reshaping required). Note that this is probably only more efficient for sparse Hamiltonians. :param ham: :type ham: time-independant hamiltonian governing evolution :returns: **rho_dot(t, y)** -- Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form). :rtype: callable .. py:function:: lindblad_eq(ham, ls, gamma) Lindblad equation, but with flattened input/output. :param ham: Time-independant hamiltonian governing evolution. :type ham: operator :param ls: Lindblad operators. :type ls: sequence of matrices :param gamma: Dampening strength. :type gamma: float :returns: **rho_dot(t, y)** -- Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form). :rtype: callable .. py:function:: lindblad_eq_vectorized(ham, ls, gamma, sparse=False) Lindblad equation, but with flattened input/output and vectorised superoperation mode (no reshaping required). :param ham: Time-independant hamiltonian governing evolution. :type ham: operator :param ls: Lindblad operators. :type ls: sequence of matrices :param gamma: Dampening strength. :type gamma: float :returns: **rho_dot(t, y)** -- Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form). :rtype: callable .. py:function:: _calc_evo_eq(isdop, issparse, isopen=False, timedep=False) Choose an appropirate dynamical equation to evolve with. .. py:class:: Try2Then3Args(fn) .. py:attribute:: fn .. py:attribute:: num_args :value: None .. py:method:: first_call(t, p, H) .. py:method:: __call__(t, p, H) .. py:class:: Evolution(p0, ham, t0=0, compute=None, int_stop=None, method='integrate', int_small_step=False, expm_backend='AUTO', expm_opts=None, progbar=False) Bases: :py:obj:`object` A class for evolving quantum systems according to Schrodinger equation. The evolution can be performed in a number of ways: - diagonalise the Hamiltonian (or use already diagonalised system). - integrate the complex ODE, that is, the Schrodinger equation, using scipy. Here either a mid- or high-order Dormand-Prince adaptive time stepping scheme is used (see :class:`scipy.integrate.complex_ode`). :param p0: Inital state, either vector or operator. If vector, converted to ket. :type p0: quantum state :param ham: Governing Hamiltonian, if tuple then assumed to contain ``(eigvals, eigvecs)`` of presolved system. If callable (but not a SciPy ``LinearOperator``), assume a time-dependent hamiltonian such that ``ham(t)`` is the Hamiltonian at time ``t``. In this case, the latest call to ``ham`` will be cached (and made immutable) in case it is needed by callbacks passed to ``compute``. :type ham: operator, tuple (1d array, operator), or callable :param t0: Initial time (i.e. time of state ``p0``), defaults to zero. :type t0: float, optional :param compute: Function(s) to compute on the state at each time step. Function(s) should take args (t, pt) or (t, pt, ham) if the Hamiltonian is required. If ham is required, it will be passed in to the function exactly as given to this ``Evolution`` instance, except if ``method`` is ``'solve'``, in which case it will be passed in as the solved system ``(eigvals, eigvecs)``. If supplied with: - single callable : ``Evolution.results`` will contain the results as a list, - dict of callables : ``Evolution.results`` will contain the results as a dict of lists with corresponding keys to those given in ``compute``. :type compute: callable, or dict of callable, optional :param int_stop: A condition to terminate the integration early if ``method`` is ``'integrate'``. This callable is called at every successful integration step and should take args (t, pt) or (t, pt, ham) similar to the function(s) in the ``compute`` argument. It should return ``-1`` to stop the integration, otherwise it should return ``None`` or ``0``. :type int_stop: callable, optional :param method: How to evolve the system: - ``'integrate'``: use definite integration. Get system at each time step, only need action of Hamiltonian on state. Generally efficient. - ``'solve'``: diagonalise dense hamiltonian. Best for small systems and allows arbitrary time steps without loss of precision. - ``'expm'``: compute the evolved state using the action of the operator exponential in a 'single shot' style. Only needs action of Hamiltonian, for very large systems can use distributed MPI. :type method: {'integrate', 'solve', 'expm'} :param int_small_step: If ``method='integrate'``, whether to use a low or high order integrator to give naturally small or large steps. :type int_small_step: bool, optional :param expm_backend: How to perform the expm_multiply function if ``method='expm'``. Can further specifiy ``'slepc-krylov'``, or ``'slepc-expokit'``. :type expm_backend: {'auto', 'scipy', 'slepc'} :param expm_opts: Supplied to :func:`~quimb.linalg.base_linalg.expm_multiply` function if ``method='expm'``. :type expm_opts: dict :param progbar: Whether to show a progress bar when calling ``at_times`` or integrating with the ``update_to`` method. :type progbar: bool, optional .. py:attribute:: _p0 .. py:attribute:: _isdop .. py:attribute:: _d .. py:attribute:: _progbar :value: False .. py:attribute:: _timedep .. py:attribute:: _method :value: 'integrate' .. py:method:: _setup_callback(fn, int_stop) Setup callbacks in the correct place to compute into _results .. py:method:: _setup_solved_ham() Solve the hamiltonian if needed and find the initial state in the energy eigenbasis for quick evolution later. .. py:method:: _start_integrator(ham, small_step) Initialize a stepping integrator. .. py:method:: _update_to_expm_ket(t) Update the simulation to time ``t``, without explicitly computing the operator exponential itself. .. py:method:: _update_to_solved_ket(t) Update simulation consisting of a solved hamiltonian and a wavefunction to time `t`. .. py:method:: _update_to_solved_dop(t) Update simulation consisting of a solved hamiltonian and a density operator to time `t`. .. py:method:: _update_to_integrate(t) Update simulation consisting of unsolved hamiltonian. .. py:method:: update_to(t) Update the simulation to time ``t`` using relevant method. :param t: Time to update the evolution to. :type t: float .. py:method:: at_times(ts) Generator expression to yield state af list of times. :param ts: Times at which to evolve to, then yield the state. :type ts: sequence of floats :Yields: **pt** (*quantum state*) -- Quantum state of evolution at next time in ``ts``. .. rubric:: Notes If integrating, currently any compute callbacks will be called at every *integration* step, not just the times `ts` -- i.e. in general len(Evolution.results) != len(ts) and if the adaptive step times are needed they should be added as a callback, e.g. ``compute['t'] = lambda t, _: return t``. .. py:property:: t Current time of simulation. :type: float .. py:property:: pt State of the system at the current time (t). :type: quantum state .. py:property:: results Results of the compute callback(s) for each time step. :type: list, or dict of lists, optional