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.

Attributes

eye

Alias for identity().

qu

Alias of quimbify().

eigh

CALLABLE_TIME_INDEP_CLASSES

Classes

qarray

Thin subclass of numpy.ndarray with some convenient quantum

Lazy

A simple class representing an unconstructed matrix. This can be passed

continuous_progbar

A continuous version of tqdm, so that it can be updated with a float

Try2Then3Args

Evolution

A class for evolving quantum systems according to Schrodinger equation.

Functions

dag(qob)

Conjugate transpose.

dot(a, b)

Matrix multiplication, dispatched to dense or sparse functions.

explt(l, t)

Complex exponenital as used in solution to schrodinger equation.

isop(qob)

Checks if qob is an operator.

issparse(qob)

Checks if qob is explicitly sparse.

ldmul(diag, mat)

Accelerated left diagonal multiplication. Equivalent to

make_immutable(mat)

Make array read only, in-place.

rdmul(mat, diag)

Accelerated left diagonal multiplication.

norm(A[, ntype])

Operator norms.

expm_multiply(mat, vec[, backend])

Compute the action of expm(mat) on vec.

norm_fro_approx(A, **kwargs)

Calculate the approximate frobenius norm of any hermitian linear

ensure_dict(x)

Make sure x is a dict, creating an empty one if x is None.

schrodinger_eq_ket(ham)

Wavefunction schrodinger equation.

schrodinger_eq_ket_timedep(ham)

Wavefunction time dependent schrodinger equation.

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

schrodinger_eq_dop_vectorized(ham)

Density operator schrodinger equation, but with flattened input/output

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

_calc_evo_eq(isdop, issparse[, isopen, timedep])

Choose an appropirate dynamical equation to evolve with.

Module Contents

quimb.evo.dag(qob)[source]

Conjugate transpose.

quimb.evo.dot(a, b)[source]

Matrix multiplication, dispatched to dense or sparse functions.

Parameters:
  • a (dense or sparse operator) – First array.

  • b (dense or sparse operator) – Second array.

Returns:

Dot product of a and b.

Return type:

dense or sparse operator

quimb.evo.explt(l, t)[source]

Complex exponenital as used in solution to schrodinger equation.

quimb.evo.eye[source]

Alias for identity().

quimb.evo.isop(qob)[source]

Checks if qob is an operator.

quimb.evo.issparse(qob)[source]

Checks if qob is explicitly sparse.

quimb.evo.ldmul(diag, mat)[source]

Accelerated left diagonal multiplication. Equivalent to numpy.diag(diag) @ mat, but faster than numpy.

Parameters:
  • diag (vector or 1d-array) – Vector representing the diagonal of a matrix.

  • mat (dense or sparse matrix) – A normal (non-diagonal) matrix.

Returns:

Dot product of the matrix whose diagonal is diag and mat.

Return type:

dense or sparse matrix

quimb.evo.make_immutable(mat)[source]

Make array read only, in-place.

Parameters:

mat (sparse or dense array) – Matrix to make immutable.

class quimb.evo.qarray(shape, dtype=float, buffer=None, offset=0, strides=None, order=None)[source]

Bases: numpy.ndarray

Thin subclass of numpy.ndarray with some convenient quantum linear algebra related methods and attributes (.H, &, etc.), and matrix-like preservation of at least 2-dimensions so as to distiguish kets and bras.

property H
property A
__array__()[source]
__and__(other)[source]
normalize(inplace=True)[source]
nmlz(inplace=True)[source]
chop(inplace=True)[source]
tr()[source]
partial_trace(dims, keep)[source]
ptr(dims, keep)[source]
__str__()[source]

Return str(self).

__repr__()[source]

Return repr(self).

quimb.evo.qu[source]

Alias of quimbify().

quimb.evo.rdmul(mat, diag)[source]

Accelerated left diagonal multiplication.

Equivalent to mat @ numpy.diag(diag), but faster.

Parameters:
  • mat (dense or sparse matrix) – A normal (non-diagonal) matrix.

  • diag (vector or 1d-array) – Vector representing the diagonal of a matrix.

Returns:

Dot product of mat and the matrix whose diagonal is diag.

Return type:

dense or sparse matrix

quimb.evo.eigh
quimb.evo.norm(A, ntype=2, **kwargs)[source]

Operator norms.

Parameters:
  • A (operator) – The operator to find norm of.

  • ntype (str) –

    Norm to calculate, if any of:

    • {2, ‘2’, ‘spectral’}: largest singular value

    • {‘f’, ‘fro’}: frobenius norm

    • {‘t’, ‘nuc’, ‘tr’, ‘trace’}: sum of singular values

Returns:

x – The operator norm.

Return type:

float

quimb.evo.expm_multiply(mat, vec, backend='AUTO', **kwargs)[source]

Compute the action of expm(mat) on vec.

Parameters:
  • mat (operator) – Operator with which to act with exponential on vec.

  • vec (vector-like) – Vector to act with exponential of operator on.

  • backend ({'AUTO', 'SCIPY', 'SLEPC', 'SLEPC-KRYLOV', 'SLEPC-EXPOKIT'}) – Which backend to use.

  • kwargs – Supplied to backend function.

Returns:

Result of expm(mat) @ vec.

Return type:

vector

class quimb.evo.Lazy(fn, *args, shape=None, factor=None, **kwargs)[source]

A simple class representing an unconstructed matrix. This can be passed to, for example, MPI workers, who can then construct the matrix themselves. The main function fn should ideally take an ownership keyword to avoid forming every row.

This is essentially like using functools.partial and assigning the shape attribute.

Parameters:
  • fn (callable) – A function that constructs an operator.

  • shape – Shape of the constructed operator.

  • args – Supplied to fn.

  • kwargs – Supplied to fn.

Returns:

Lazy

Return type:

callable

Examples

Setup the lazy operator:

>>> H_lazy = Lazy(ham_heis, n=10, shape=(2**10, 2**10), sparse=True)
>>> H_lazy
<Lazy(ham_heis, shape=(1024, 1024), dtype=None)>

Build a matrix slice (usually done automatically by e.g. eigs):

>>> H_lazy(ownership=(256, 512))
<256x1024 sparse matrix of type '<class 'numpy.float64'>'
        with 1664 stored elements in Compressed Sparse Row format>
__imul__(x)[source]
__mul__(x)[source]
__rmul__(x)[source]
__call__(**kwargs)[source]
__repr__()[source]

Return repr(self).

quimb.evo.norm_fro_approx(A, **kwargs)[source]

Calculate the approximate frobenius norm of any hermitian linear operator:

\[\mathrm{Tr} \left[ A^{\dagger} A \right]\]
Parameters:
  • A (linear operator like) – Operator with a dot method, assumed to be hermitian, to estimate the frobenius norm of.

  • kwargs – Supplied to approx_spectral_function().

Return type:

float

class quimb.evo.continuous_progbar(*args, total=100, **kwargs)[source]

Bases: tqdm.tqdm

A continuous version of tqdm, so that it can be updated with a float within some pre-given range, rather than a number of steps.

Parameters:
  • args ((stop) or (start, stop)) – Stopping point (and starting point if len(args) == 2) of window within which to evaluate progress.

  • total (int) – The number of steps to represent the continuous progress with.

  • kwargs – Supplied to tqdm.tqdm

cupdate(x)[source]

‘Continuous’ update of progress bar.

Parameters:

x (float) – Current position within the range [self.start, self.stop].

quimb.evo.ensure_dict(x)[source]

Make sure x is a dict, creating an empty one if x is None.

quimb.evo.CALLABLE_TIME_INDEP_CLASSES
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

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.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._calc_evo_eq(isdop, issparse, isopen=False, timedep=False)[source]

Choose an appropirate dynamical equation to evolve with.

class quimb.evo.Try2Then3Args(fn)[source]
first_call(t, p, H)[source]
__call__(t, p, H)[source]
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]

Bases: 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 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.

_setup_callback(fn, int_stop)[source]

Setup callbacks in the correct place to compute into _results

_setup_solved_ham()[source]

Solve the hamiltonian if needed and find the initial state in the energy eigenbasis for quick evolution later.

_start_integrator(ham, small_step)[source]

Initialize a stepping integrator.

_update_to_expm_ket(t)[source]

Update the simulation to time t, without explicitly computing the operator exponential itself.

_update_to_solved_ket(t)[source]

Update simulation consisting of a solved hamiltonian and a wavefunction to time t.

_update_to_solved_dop(t)[source]

Update simulation consisting of a solved hamiltonian and a density operator to time t.

_update_to_integrate(t)[source]

Update simulation consisting of unsolved hamiltonian.

update_to(t)[source]

Update the simulation to time t using relevant method.

Parameters:

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

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 t
Current time of simulation.
Type:

float

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