quimb.evo#

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.

Functions

lindblad_eq(ham, ls, gamma)

Lindblad equation, but with flattened input/output.

lindblad_eq_vectorized(ham, ls, gamma[, sparse])

Lindblad equation, but with flattened input/output and vectorised superoperation mode (no reshaping required).

schrodinger_eq_dop(ham)

Density operator schrodinger equation, but with flattened input/output.

schrodinger_eq_dop_timedep(ham)

Time dependent density operator schrodinger equation, but with flattened input/output.

schrodinger_eq_dop_vectorized(ham)

Density operator schrodinger equation, but with flattened input/output and vectorised superoperator mode (no reshaping required).

schrodinger_eq_ket(ham)

Wavefunction schrodinger equation.

schrodinger_eq_ket_timedep(ham)

Wavefunction time dependent schrodinger equation.

Classes

Evolution(p0, ham[, t0, compute, int_stop, ...])

A class for evolving quantum systems according to Schrodinger equation.

Try2Then3Args(fn)

class quimb.evo.Evolution(p0, ham, t0=0, compute=None, int_stop=None, method='integrate', int_small_step=False, expm_backend='AUTO', expm_opts=None, progbar=False)[source]#

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 scipy.integrate.complex_ode).

Parameters
  • p0 (quantum state) – Inital state, either vector or operator. If vector, converted to ket.

  • ham (operator, tuple (1d array, operator), or callable) – 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.

  • t0 (float, optional) – Initial time (i.e. time of state p0), defaults to zero.

  • compute (callable, or dict of callable, optional) –

    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.

  • int_stop (callable, optional) – 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.

  • method ({'integrate', 'solve', 'expm'}) –

    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.

  • int_small_step (bool, optional) – If method='integrate', whether to use a low or high order integrator to give naturally small or large steps.

  • expm_backend ({'auto', 'scipy', 'slepc'}) – How to perform the expm_multiply function if method='expm'. Can further specifiy 'slepc-krylov', or 'slepc-expokit'.

  • expm_opts (dict) – Supplied to expm_multiply() function if method='expm'.

  • progbar (bool, optional) – Whether to show a progress bar when calling at_times or integrating with the update_to method.

at_times(ts)[source]#

Generator expression to yield state af list of times.

Parameters

ts (sequence of floats) – Times at which to evolve to, then yield the state.

Yields

pt (quantum state) – Quantum state of evolution at next time in ts.

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.

property pt#

State of the system at the current time (t).

Type

quantum state

property results#

Results of the compute callback(s) for each time step.

Type

list, or dict of lists, optional

property t#

Current time of simulation.

Type

float

update_to(t)[source]#

Update the simulation to time t using relevant method.

Parameters

t (float) – Time to update the evolution to.

quimb.evo.lindblad_eq(ham, ls, gamma)[source]#

Lindblad equation, but with flattened input/output.

Parameters
  • ham (operator) – Time-independant hamiltonian governing evolution.

  • ls (sequence of matrices) – Lindblad operators.

  • gamma (float) – Dampening strength.

Returns

rho_dot(t, y) – Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form).

Return type

callable

quimb.evo.lindblad_eq_vectorized(ham, ls, gamma, sparse=False)[source]#

Lindblad equation, but with flattened input/output and vectorised superoperation mode (no reshaping required).

Parameters
  • ham (operator) – Time-independant hamiltonian governing evolution.

  • ls (sequence of matrices) – Lindblad operators.

  • gamma (float) – Dampening strength.

Returns

rho_dot(t, y) – Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form).

Return type

callable

quimb.evo.schrodinger_eq_dop(ham)[source]#

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.

Parameters

ham (operator) – 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).

Return type

callable

quimb.evo.schrodinger_eq_dop_timedep(ham)[source]#

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.

Parameters

ham (callable) – Time-dependant Hamiltonian governing evolution, such that ham(t) returns an operator representation of the Hamiltonian at time t.

Returns

rho_dot(t, y) – Function to calculate rho_dot(t) at rho(t), input and output both in ravelled (1D form).

Return type

callable

quimb.evo.schrodinger_eq_dop_vectorized(ham)[source]#

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.

Parameters

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).

Return type

callable

quimb.evo.schrodinger_eq_ket(ham)[source]#

Wavefunction schrodinger equation.

Parameters

ham (operator) – Time-independant Hamiltonian governing evolution.

Returns

psi_dot(t, y) – Function to calculate psi_dot(t) at psi(t).

Return type

callable

quimb.evo.schrodinger_eq_ket_timedep(ham)[source]#

Wavefunction time dependent schrodinger equation.

Parameters

ham (callable) – Time-dependant Hamiltonian governing evolution, such that ham(t) returns an operator representation of the Hamiltonian at time t.

Returns

psi_dot(t, y) – Function to calculate psi_dot(t) at psi(t).

Return type

callable