All Modules
oqupy.base_api
Module for base classes of API objects.
- class oqupy.base_api.BaseAPIClass(name: Optional[str] = None, description: Optional[str] = None)[source]
Base class for API objects
- Parameters
name (str) – An optional name for the object.
description (str) – An optional description of the object.
- property description
Detailed description of the object.
- property name
Name of the object.
oqupy.bath
Module on physical information on the bath and its coupling to the system.
- class oqupy.bath.Bath(coupling_operator: ndarray, correlations: BaseCorrelations, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Represents the bath degrees of freedom with a specific coupling operator (to the system) and a specific auto-correlation function.
- Parameters
coupling_operator (np.ndarray) – The system operator to which the bath couples.
correlations (BaseCorrelations) – The bath’s auto correlation function.
name (str) – An optional name for the bath.
description (str) – An optional description of the bath.
- property correlations: BaseCorrelations
The correlations of the bath.
- property coupling_operator: ndarray
The diagonalised system operator to which the bath couples.
- property dimension: ndarray
Hilbert space dimension of the coupling operator.
- property unitary_transform: ndarray
The unitary that makes the coupling operator diagonal.
oqupy.bath_dynamics
Module for calculating bath dynamics as outlined in [Gribben2021].
[Gribben2021] D. Gribben, A. Strathearn, G. E. Fux, P. Kirton, and B. W. Lovett, Using the Environment to Understand non-Markovian Open Quantum Systems, arXiv:2106.04212 [quant-ph] (2021).
- class oqupy.bath_dynamics.TwoTimeBathCorrelations(system: BaseSystem, bath: Bath, process_tensor: BaseProcessTensor, initial_state: Optional[ndarray] = None, system_correlations: Optional[ndarray] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Class to facilitate calculation of two-time bath correlations.
- Parameters
system (BaseSystem) – The system.
bath (Bath) – The bath object containing all coupling information and temperature.
process_tensor (ProcessTensor) – The corresponding process tensor calculated for the given bath.
initial_state (ndarray) – Initial state of the system.
system_correlations (ndarray) – Optional previously calculated system correlations. This must be an upper triangular array with all ordered correlations up to a certain time.
name (str) – An optional name for the bath dynamics object.
description (str) – An optional description of the bath dynamics object.
- correlation(freq_1: float, time_1: float, freq_2: Optional[float] = None, time_2: Optional[float] = None, dw: Optional[tuple] = (1.0, 1.0), dagg: Optional[tuple] = (1, 0), interaction_picture: Optional[bool] = False, change_only: Optional[bool] = False, progress_type: Optional[str] = None) complex [source]
Function to calculate two-time correlation function between two frequency bands of a bath.
The calculation consists of a double integral of the form:
\[\int_0^t \int_0^{t'} \left\{ \mathrm{Re} \langle O(t')O(t'') \rangle \, K_R(t',t'') + i \,\mathrm{Im} \langle O(t')O(t'') \rangle \, K_I(t',t'') \right\} dt'' dt'\]where \(O\) is the system operator coupled to the bath and \(K_R\) and \(K_I\) are generally piecewise kernels which depend on the exact bath correlation function desired.
- Parameters
freq_1 (float) – Frequency of the earlier time operator.
time_1 (float) – Time the earlier operator acts.
freq_2 (float) – Frequency of the later time operator. If set to None will default to freq_2=freq_1.
time_2 (float) – Time the later operator acts. If set to None will default to time_2=time_1.
dw (tuple) – Width of the the frequency bands. By default this method returns a correlation density by setting the frequency bands to dw=(1.0, 1.0).
dagg (tuple) – Determines whether each operator is daggered or not e.g. (1,0) would correspond to \(< a^\dagger a >\).
interaction_picture (bool) – Option whether to generate the result within the bath interaction picture.
change_only (bool) – Option to include the initial occupation in the result.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.
- Returns
correlation – Bath correlation function <a^{dagg[0]}_{freq_2} (time_2) a^{dagg[1]}_{freq_1} (time_1)>
- Return type
complex
- generate_system_correlations(final_time: float, progress_type: Optional[str] = None) None [source]
Function to generate all ordered system correlations up to a given time using the process tensor.
- Parameters
final_time (float) – The latest time appearing in the generated system correlation functions.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.
- property initial_state: ndarray
The initial system state.
- occupation(freq: float, dw: Optional[float] = 1.0, change_only: Optional[bool] = False, progress_type: Optional[str] = None) Tuple[ndarray, ndarray] [source]
Function to calculate the change in bath occupation in a particular bandwidth.
- Parameters
freq (float) – Central frequency of the frequency band.
dw (float) – Width of the the frequency band. By default this method returns a a density by setting the frequency band dw=1.0.
change_only (bool) – Option to include the initial occupation (density) in the result.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.
- Returns
times (ndarray) – Times of the occupation dynamics.
bath_occupation (ndarray) – Occupation (density) (difference) of the bath in the specified frequency band.
- property system: BaseSystem
The system.
oqupy.contractions
Module for various applications involving contractions of the process tensor.
- oqupy.contractions.compute_correlations(system: BaseSystem, process_tensor: BaseProcessTensor, operator_a: ndarray, operator_b: ndarray, times_a: Union[int, slice, List[Union[int, slice]], float, Tuple[float, float]], times_b: Union[int, slice, List[Union[int, slice]], float, Tuple[float, float]], time_order: Optional[str] = 'ordered', initial_state: Optional[ndarray] = None, start_time: Optional[float] = 0.0, dt: Optional[float] = None, progress_type: Optional[str] = None) Tuple[ndarray, ndarray, ndarray] [source]
Compute system correlations for a given system Hamiltonian.
Times may be specified with indices, a single float, or a pair of floats specifying the start and end time. Indices may be integers, slices, or lists of integers and slices.
- Parameters
system (BaseSystem) – Object containing the system Hamiltonian.
process_tensor (BaseProcessTensor) – A process tensor object.
operator_a (ndarray) – System operator \(\hat{A}\).
operator_b (ndarray) – System operator \(\hat{B}\).
times_a (Union[Indices, float, Tuple[float, float]]) – Time(s) \(t_A\).
times_b (Union[Indices, float, Tuple[float, float]]) – Time(s) \(t_B\).
time_order (str (default =
'ordered'
)) – Which two time correlations to compute. Types are: {'ordered'
,'anti'
,'full'
}.initial_state (ndarray (default = None)) – Initial system state.
start_time (float (default = 0.0)) – Initial time.
dt (float (default = None)) – Time step size.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
'silent'
,'simple'
,'bar'
}. If None then the default progress type is used.
- Returns
times_a (ndarray) – The \(N\) times \(t^A_n\).
times_b (ndarray) – The \(M\) times \(t^B_m\).
correlations (ndarray) – The \(N \times M\) correlations \(\langle B(t^B_m) A(t^A_n) \rangle\). Entries that are outside the scope specified in time_order are set to be NaN + NaN j.
- oqupy.contractions.compute_dynamics(system: Union[System, TimeDependentSystem], initial_state: Optional[ndarray] = None, dt: Optional[float] = None, num_steps: Optional[int] = None, start_time: Optional[float] = 0.0, process_tensor: Optional[Union[List[BaseProcessTensor], BaseProcessTensor]] = None, control: Optional[Control] = None, record_all: Optional[bool] = True, progress_type: Optional[str] = None) Dynamics [source]
Compute the system dynamics for a given system Hamiltonian.
- Parameters
system (Union[System, TimeDependentSystem]) – Object containing the system Hamiltonian information.
initial_state (ndarray) – Initial system state.
dt (float) – Length of a single time step.
num_steps (int) – Optional number of time steps to be computed.
start_time (float) – Optional start time offset.
process_tensor (Union[List[BaseProcessTensor],BaseProcessTensor]) – Optional process tensor object or list of process tensor objects.
control (Control) – Optional control operations.
record_all (bool) – If false function only computes the final state.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.
- Returns
dynamics – The system dynamics for the given system Hamiltonian (accounting for the interaction with the environment).
- Return type
- oqupy.contractions.compute_dynamics_with_field(system: TimeDependentSystemWithField, initial_field: complex, initial_state: Optional[ndarray] = None, dt: Optional[float] = None, num_steps: Optional[int] = None, start_time: Optional[float] = 0.0, process_tensor: Optional[Union[List[BaseProcessTensor], BaseProcessTensor]] = None, control: Optional[Control] = None, record_all: Optional[bool] = True, epsrel: Optional[float] = 1.4901161193847656e-08, subdiv_limit: Optional[int] = 256, progress_type: Optional[str] = None) DynamicsWithField [source]
Compute the system and field dynamics for a given system Hamiltonian and field equation of motion.
- Parameters
system (TimeDependentSystemWithField) – Object containing the system Hamiltonian and field equation of motion.
initial_field (complex) – Initial field value.
initial_state (ndarray) – Initial system state.
dt (float) – Length of a single time step.
start_time (float) – Optional start time offset.
num_steps (int) – Optional number of time steps to be computed.
process_tensor (Union[List[BaseProcessTensor],BaseProcessTensor]) – Optional process tensor object or list of process tensor objects.
control (Control) – Optional control operations.
record_all (bool) – If false function only computes the final state.
subdiv_limit (int (default = config.SUBDIV_LIMIT)) – The maximum number of subdivisions used during the adaptive algorithm when integrating the system Liouvillian. If None then the Liouvillian is not integrated but sampled twice to to construct the system propagators at each timestep.
epsrel (float (default = config.INTEGRATE_EPSREL)) – The relative error tolerance for the adaptive algorithm when integrating the system Liouvillian.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.
- Returns
dynamics_with_field – The system and field dynamics for the given system Hamiltonian and field equation of motion (accounting for the interaction with the environment).
- Return type
oqupy.control
Module for system ‘control operations’ as discussed in [Pollock2018].
[Pollock2018] F. A. Pollock, C. Rodriguez-Rosario, T. Frauenheim, M. Paternostro, and K. Modi, Non-Markovian quantumprocesses: Complete framework and efficient characterization, Phys. Rev. A 97, 012127 (2018).
- class oqupy.control.ChainControl(hilbert_space_dimensions: List[int], name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Control operations on a linear system chain.
- Parameters
hilbert_space_dimensions (List[int]) – Hilbert space dimension for each chain site.
name (str) – An optional name for the chain controls.
description (str) – An optional description of the chain controls.
- add_single_site_control(control: ndarray, site: int, step: int, post: Optional[bool] = False, name: Optional[str] = None) None [source]
Add a control operation at site site and time step step.
- Parameters
control (ndarray) – Control operation in Liouville space.
site (int) – Site index.
step (int) – Timestep to which the control should be applied.
post (bool) – True if the control should be applied after the measurement of this time step.
name (Text) – An optional name to recognize a control operation.
- get_single_site_controls(step: int, post: bool) List[ndarray] [source]
Get a list of single site controls for the time step step.
- Parameters
step (int) – The time step.
post (bool) – If True (False) the set of control superoperators that should be applied after (before) the measurement is returned.
- Returns
superoperators_list – List of single site control superoperators.
- Return type
list[ndarray]
- property hs_dims
Hilbert space dimensions.
- class oqupy.control.Control(dimension: int, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Represents a set of system control operations.
A control operation is a superoperator that acts on the system instantaneously at a particular time, as described in [Pollock2018].
- Parameters
dimension (int) – The Hilbert space dimension of the system.
name (str) – An optional name for the set of control operations.
description (str) – An optional description of the set of control operations.
- add_continuous(control_fct: Callable[[ndarray, float], numpy.ndarray], post: Optional[bool] = False) None [source]
ToDo
- add_single(time: Union[int, float], control_operation: ndarray, post: Optional[bool] = False) None [source]
Adds a single control operation at time time.
- Parameters
time (Union[int, float]) – The time at which the operation should be applied. If type(time) is int then time is understood as the timestep to which it shall be applied.
control_operation (ndarray) – The control operation super operator of shape \(d^2 \times d^2\), where \(d\) is the system Hilbert space dimension.
post (bool) – If True (False) the operator is applied at the corresponding time step after (before) a possible measurement of the state.
- property dimension
Hilbert space dimension of the controlled system.
- get_controls(step: int, dt: Optional[float] = None, start_time: Optional[float] = 0.0) Tuple[ndarray, ndarray] [source]
Get the pre and post measurement control operation for a specific time step.
- Parameters
step (int) – The time step.
dt (float) – The time step length.
start_time (float) – The initial time step off-set.
- Returns
pre (ndarray) – The control superoperator that should be applied before a state measurement.
post (ndarray) – The control superoperator that should be applied after a state measurement.
oqupy.correlations
Module for environment correlations.
- class oqupy.correlations.BaseCorrelations(name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Base class for environment auto-correlations.
- correlation(tau: Any, epsrel: Optional[float] = 1.4901161193847656e-08, subdiv_limit: Optional[int] = 256) Any [source]
Auto-correlation function.
\[C(\tau) = C(t, t-\tau) \ = \langle F(t) F(t-\tau) \rangle_\mathrm{env}\]where \(\tau\) is the time difference tau and \(F(t)\) is the the environment part of the coupling operator in Heisenberg picture with respect to the environment Hamiltonian.
- Parameters
tau (ndarray) – Time difference \(\tau\)
epsrel (float) – Relative error tolerance.
subdiv_limit (int) – Maximal number of interval subdivisions for numerical integration.
- Returns
correlation – The auto-correlation function \(C(\tau)\) at time \(\tau\).
- Return type
ndarray
- correlation_2d_integral(delta: float, time_1: float, time_2: Optional[float] = None, shape: Optional[str] = 'square', epsrel: Optional[float] = 1.4901161193847656e-08, subdiv_limit: Optional[int] = 256) complex [source]
2D integrals of the correlation function
\[ \begin{align}\begin{aligned}\eta_\mathrm{square} = \int_{t_1}^{t_1+\Delta} \int_{0}^{\Delta} C(t'-t'') dt'' dt'\\\eta_\mathrm{upper-triangle} = \int_{t_1}^{t_1+\Delta} \int_{0}^{t'-t_1} C(t'-t'') dt'' dt'\\\eta_\mathrm{lower-triangle} = \int_{t_1}^{t_1+\Delta} \int_{t'-t_1}^{\Delta} C(t'-t'') dt'' dt'\\\eta_\mathrm{rectangle} = \int_{t_1}^{t_2} \int_{0}^{\Delta} C(t'-t'') dt'' dt'\end{aligned}\end{align} \]for shape either
'square'
,'upper-triangle'
,'lower-triangle'
, or'rectangle'
.- Parameters
delta (float) – Length of integration intervals.
time_1 (float) – Lower bound of integration interval of \(dt'\).
time_2 (float) – Upper bound of integration interval of \(dt'\) for shape =
'rectangle'
.shape (str (default =
'square'
)) – The shape of the 2D integral. Shapes are: {'square'
,'upper-triangle'
,'lower-triangle'
,'lower-triangle'
}epsrel (float) – Relative error tolerance.
subdiv_limit (int) – Maximal number of interval subdivisions for numerical integration.
- Returns
integral – The numerical value for the two dimensional integral \(\eta_\mathrm{shape}\).
- Return type
float
- class oqupy.correlations.CustomCorrelations(correlation_function: Callable[float, float], name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseCorrelations
Encodes a custom auto-correlation function
\[C(\tau) = C(t, t-\tau) = \langle F(t) F(t-\tau) \rangle_\mathrm{env}\]with time difference tau \(\tau\) and \(F(t)\) is the the environment part of the coupling operator in Heisenberg picture with respect to the environment Hamiltonian. We assume that \(C(\tau) = 0\) for all \(\tau > \tau_\mathrm{max}\).
- Parameters
correlation_function (callable) – The correlation function \(C\).
name (str) – An optional name for the correlations.
description (str) – An optional description of the correlations.
- correlation(tau: Any, epsrel: Optional[float] = None, subdiv_limit: Optional[int] = None) Any [source]
Auto-correlation function.
\[C(\tau) = C(t, t-\tau) \ = \langle F(t) F(t-\tau) \rangle_\mathrm{env}\]with time difference tau \(\tau\) and \(F(t)\) is the the environment part of the coupling operator in Heisenberg picture with respect to the environment Hamiltonian.
- Parameters
tau (ndarray) – Time difference \(\tau\)
epsrel (float) – Relative error tolerance (has no effect here).
subdiv_limit (int) – Maximal number of interval subdivisions for numerical integration (has no effect here).
- Returns
correlation – The auto-correlation function \(C(\tau)\) at time \(\tau\).
- Return type
ndarray
- correlation_2d_integral(delta: float, time_1: float, time_2: Optional[float] = None, shape: Optional[str] = 'square', epsrel: Optional[float] = 1.4901161193847656e-08, subdiv_limit: Optional[int] = 256) complex [source]
2D integrals of the correlation function
\[ \begin{align}\begin{aligned}\eta_\mathrm{square} = \int_{t_1}^{t_1+\Delta} \int_{0}^{\Delta} C(t'-t'') dt'' dt'\\\eta_\mathrm{upper-triangle} = \int_{t_1}^{t_1+\Delta} \int_{0}^{t'-t_1} C(t'-t'') dt'' dt'\\\eta_\mathrm{lower-triangle} = \int_{t_1}^{t_1+\Delta} \int_{t'-t_1}^{\Delta} C(t'-t'') dt'' dt'\\\eta_\mathrm{rectangle} = \int_{t_1}^{t_2} \int_{0}^{\Delta} C(t'-t'') dt'' dt'\end{aligned}\end{align} \]for shape either
'square'
,'upper-triangle'
,'lower-triangle'
, or'rectangle'
.- Parameters
delta (float) – Length of integration intervals.
time_1 (float) – Lower bound of integration interval of \(dt'\).
time_2 (float) – Upper bound of integration interval of \(dt'\) for shape =
'rectangle'
.shape (str (default =
'square'
)) – The shape of the 2D integral. Shapes are: {'square'
,'upper-triangle'
,'lower-triangle'
,'lower-triangle'
}epsrel (float) – Relative error tolerance.
subdiv_limit (int) – Maximal number of interval subdivisions for numerical integration.
- Returns
integral – The numerical value for the two dimensional integral \(\eta_\mathrm{shape}\).
- Return type
float
- class oqupy.correlations.CustomSD(j_function: Callable[float, float], cutoff: float, cutoff_type: Optional[str] = 'exponential', temperature: Optional[float] = 0.0, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseCorrelations
Correlations corresponding to a custom spectral density. The resulting spectral density is
\[J(\omega) = j(\omega) X(\omega,\omega_c) ,\]with j_function \(j\), cutoff \(\omega_c\) and a cutoff type \(X\).
If cutoff_type is
'hard'
then \(X(\omega,\omega_c)=\Theta(\omega_c-\omega)\), where \(\Theta\) is the Heaviside step function.'exponential'
then \(X(\omega,\omega_c)=\exp(-\omega/\omega_c)\).'gaussian'
then \(X(\omega,\omega_c)=\exp(-\omega^2/\omega_c^2)\).
- Parameters
j_function (callable) – The spectral density \(j\) without the cutoff.
cutoff (float) – The cutoff frequency \(\omega_c\).
cutoff_type (str (default =
'exponential'
)) – The cutoff type. Types are: {'hard'
,'exponential'
,'gaussian'
}temperature (float) – The environment’s temperature.
name (str) – An optional name for the correlations.
description (str) – An optional description of the correlations.
- correlation(tau: Any, epsrel: Optional[float] = 1.4901161193847656e-08, subdiv_limit: Optional[int] = 256) Any [source]
Auto-correlation function associated to the spectral density at the given temperature \(T\)
\[C(\tau) = \int_0^{\infty} J(\omega) \ \left[ \cos(\omega \tau) \ \coth\left( \frac{\omega}{2 T}\right) \ - i \sin(\omega \tau) \right] \mathrm{d}\omega .\]with time difference tau \(\tau\).
- Parameters
tau (ndarray) – Time difference \(\tau\)
epsrel (float) – Relative error tolerance.
subdiv_limit (int) – Maximal number of interval subdivisions for numerical integration.
- Returns
correlation – The auto-correlation function \(C(\tau)\) at time \(\tau\).
- Return type
ndarray
- correlation_2d_integral(delta: float, time_1: float, time_2: Optional[float] = None, shape: Optional[str] = 'square', epsrel: Optional[float] = 1.4901161193847656e-08, subdiv_limit: Optional[int] = 256) complex [source]
2D integrals of the correlation function
\[ \begin{align}\begin{aligned}\eta_\mathrm{square} = \int_{t_1}^{t_1+\Delta} \int_{0}^{\Delta} C(t'-t'') dt'' dt'\\\eta_\mathrm{upper-triangle} = \int_{t_1}^{t_1+\Delta} \int_{0}^{t'-t_1} C(t'-t'') dt'' dt'\\\eta_\mathrm{lower-triangle} = \int_{t_1}^{t_1+\Delta} \int_{t'-t_1}^{\Delta} C(t'-t'') dt'' dt'\\\eta_\mathrm{rectangle} = \int_{t_1}^{t_2} \int_{0}^{\Delta} C(t'-t'') dt'' dt'\end{aligned}\end{align} \]for shape either
'square'
,'upper-triangle'
,'lower-triangle'
, or'rectangle'
.- Parameters
delta (float) – Length of integration intervals.
time_1 (float) – Lower bound of integration interval of \(dt'\).
time_2 (float) – Upper bound of integration interval of \(dt'\) for shape =
'rectangle'
.shape (str (default =
'square'
)) – The shape of the 2D integral. Shapes are: {'square'
,'upper-triangle'
,'lower-triangle'
,'lower-triangle'
}epsrel (float) – Relative error tolerance.
subdiv_limit (int) – Maximal number of interval subdivisions for numerical integration.
- Returns
integral – The numerical value for the two dimensional integral \(\eta_\mathrm{shape}\).
- Return type
float
- spectral_density(omega: Any) Any [source]
The resulting spectral density (including the cutoff).
- Parameters
omega (ndarray) – The frequency \(\omega\) for which we want to know the spectral density.
- Returns
spectral_density – The resulting spectral density \(J(\omega)\) at the frequency \(\omega\).
- Return type
ndarray
- class oqupy.correlations.PowerLawSD(alpha: float, zeta: float, cutoff: float, cutoff_type: str = 'exponential', temperature: Optional[float] = 0.0, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
CustomSD
Correlations corresponding to the spectral density of the standard form
\[J(\omega) = 2 \alpha \frac{\omega^\zeta}{\omega_c^{\zeta-1}} \ X(\omega,\omega_c)\]with alpha \(\alpha\), zeta \(\zeta\) and a cutoff type \(X\).
If cutoff_type is
'hard'
then \(X(\omega,\omega_c)=\Theta(\omega_c-\omega)\), where \(\Theta\) is the Heaviside step function.'exponential'
then \(X(\omega,\omega_c)=\exp(-\omega/\omega_c)\).'gaussian'
then \(X(\omega,\omega_c)=\exp(-\omega^2/\omega_c^2)\).
- Parameters
alpha (float) – The coupling strength \(\alpha\).
zeta (float) – The exponent \(\zeta\) (corresponds to the dimensionality of the environment). The environment is called ohmic if \(\zeta=1\), superohmic if \(\zeta>1\) and subohmic if \(\zeta<1\)
cutoff (float) – The cutoff frequency \(\omega_c\).
cutoff_type (str (default =
'exponential'
)) – The cutoff type. Types are: {'hard'
,'exponential'
,'gaussian'
}temperature (float) – The environment’s temperature.
name (str) – An optional name for the correlations.
description (str) – An optional description of the correlations.
oqupy.dynamics
Module on the discrete time evolution of a density matrix.
- class oqupy.dynamics.BaseDynamics(name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Base class for objects recording the dynamics of an open quantum system. Consists at least of a list of times at which the dynamics has been computed and a list of states describing the system density matrix at those times.
- Parameters
name (str) – An optional name for the dynamics.
description (str) – An optional description of the dynamics.
- expectations(operator: Optional[ndarray] = None, real: Optional[bool] = False) Tuple[ndarray, ndarray] [source]
Return the time evolution of the expectation value of specific operator. The expectation for \(t\) is
\[\langle \hat{O}(t) \rangle = \mathrm{Tr}\{ \hat{O} \rho(t) \}\]with operator \(\hat{O}\).
- Parameters
operator (ndarray (default = None)) – The operator \(\hat{O}\). If operator is None then the trace of \(\rho(t)\) is returned.
real (bool (default = False)) – If set True then only the real part of the expectation is returned.
- Returns
times (ndarray) – The points in time \(t\).
expectations (ndarray) – Expectation values \(\langle \hat{O}(t) \rangle\).
- property shape: ndarray
Numpy shape of the states.
- property states: ndarray
States of the dynamics.
- property times: ndarray
Times of the dynamics.
- class oqupy.dynamics.Dynamics(times: Optional[List[float]] = None, states: Optional[List[ndarray]] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseDynamics
Represents a specific time evolution of a density matrix.
- Parameters
times (List[float] (default = None)) – A list of points in time.
states (List[ndarray] (default = None)) – A list of states at the times times.
name (str) – An optional name for the dynamics.
description (str) – An optional description of the dynamics.
- class oqupy.dynamics.DynamicsWithField(times: Optional[List[float]] = None, states: Optional[List[ndarray]] = None, fields: Optional[List[complex]] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseDynamics
Represents a specific time evolution of a density matrix together with a coherent field.
- Parameters
times (List[float] (default = None)) – A list of points in time.
states (List[ndarray] (default = None)) – A list of states at the times times.
fields (List[complex] (default = None)) – A list of fields at the times times.
name (str) – An optional name for the dynamics.
description (str) – An optional description of the dynamics.
- add(time: float, state: ndarray, field: complex) None [source]
Append a state and field at a specific time to the time evolution.
- Parameters
time (float) – The point in time.
state (ndarray) – The state at the time time.
field (complex) – The field at the time time.
- field_expectations() Tuple[ndarray, ndarray] [source]
Return the time evolution of the coherent field.
- Returns
times (ndarray) – The points in time \(t\).
field_expectations (ndarray) – Values \(\langle a(t) \rangle\).
- property fields: ndarray
Fields of the dynamics.
oqupy.exceptions
Custom TEMPO Warnings and Errors.
oqupy.helpers
Handy helper functions.
- oqupy.helpers.plot_correlations_with_parameters(correlations: BaseCorrelations, parameters: TempoParameters, ax: Optional[Axes] = None) Axes [source]
Plot the correlation function on a grid that corresponds to some tempo parameters. For comparison, it also draws a solid line that is 10% longer and has two more sampling points per interval.
- Parameters
correlations (BaseCorrelations) – The correlation object we are interested in.
parameters (TempoParameters) – The tempo parameters that determine the grid.
oqupy.mps_mpo
Module for MPSs and MPOs.
- class oqupy.mps_mpo.AugmentedMPS(gammas: List[ndarray], lambdas: Optional[List[ndarray]] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
An augmented matrix product state (as introduced in the supplemental material of [Fux2022]).
The full gamma tensors (one for each site) have the following axis: (L = left bond leg, P = physical leg, R = right bond leg, A = augmented leg).
Depending on the rank of the input gamma tensor, it is completed with legs of dimension 1 according to the following interpretation of the input:
rank = 1: Product state vectorized density matrix, i.e. (1, P, 1, 1)
rank = 2: Product state density matrix, i.e. (1, p*p, 1, 1)
rank = 3: Canonical MPS gamma tensor, i.e. (L, P, R, 1)
rank = 4: Augmented MPS gamma tensor, i.e. (L, P, R, A)
If no lambdas are given, they are assumed to be identities. A single lambda can be None (identity), a vector (giving the diagonal values) or a matrix (which must be diagonal).
- Parameters
gammas (List[ndarray]) – The input gamma tensors.
lambdas (List[ndarray]) – The input lambda diagonal matrices.
name (str) – An optional name for the augmented MPS.
description (str) – An optional description of the augmented MPS.
- property gammas: ndarray
“The gamma tensors.
- property lambdas: ndarray
“The values of the lambda matrices diagonals.
- class oqupy.mps_mpo.Gate(sites: List[int], tensors: List[ndarray])[source]
Bases:
object
Class representing an n-site gate in MPO form.
The axes of the MPO tensors are (L = left bond leg, I = input leg, O = output leg, R = right bond leg):
for n=1: (I, O)
for n=2: (I, O, R), (L, I, O)
for n>2: (I, O, R), (L, I, O, R), …, (L, I, O, R), (L, I, O)
- Parameters
sites (List[int]) – The site numbers onto which the MPO Gate acts.
tensors (List[ndarray]) – The MPO tensors of the gate.
- property sites: List[int]
The sites onto which the gate acts.
- property span: int
The span of sites onto which the gate acts.
- property tensors: List[ndarray]
The tensors of the MPO gate.
- class oqupy.mps_mpo.GateLayer(parallel: bool, gates: List[Gate])[source]
Bases:
object
A layer of gates.
- Parameters
parallel (bool) – Whether of not the gates are suitable for a parallel application.
gates (List[Gates]) – List of gates.
- property parallel: bool
Whether of not the gates are suitable for a parallel application.
- class oqupy.mps_mpo.NnGate(site: int, tensors: Tuple[ndarray, ndarray])[source]
Bases:
Gate
An MPO gate acting on a pair of neighboring sites.
- Parameters
site (int) – The index of the left site.
tensor (ndarray) – The two MPO tensors of shape.
- class oqupy.mps_mpo.SiteGate(site: int, tensor: ndarray)[source]
Bases:
Gate
An MPO gate acting on a single site.
- Parameters
site (int) – The site onto which the MPO gate acts.
tensor (ndarray) – The single site MPO (which is a matrix).
- class oqupy.mps_mpo.TebdPropagator(gate_layers: List[GateLayer])[source]
Bases:
object
TEBD (Time Evolving Block Decimation) Propagators consist of a list of GateLayers.
- Parameters
gate_layers (List[GateLayer]) – The gate layers that make up a full time step propagation in a TEBD tensor network.
- oqupy.mps_mpo.compute_nn_gate(liouvillian: ndarray, site: int, hs_dim_l: int, hs_dim_r: int, dt: float, epsrel: float) NnGate [source]
Compute nearest neighbor gate from Liouvillian.
- Parameters
liouvillian (ndarray) – The two site Liouvillian.
site (int) – The index of the left site.
hs_dim_l (int) – The Hilbert space dimension of the left site.
hs_dim_r (int) – The Hilbert space dimension of the right site.
dt (float) – The time step length.
epsrel (float) – The relative singular value truncation tolerance.
- Returns
nn_gate – Nearest neighbor gate.
- Return type
- oqupy.mps_mpo.compute_tebd_propagator(system_chain: Any, time_step: float, epsrel: float, order: int) TebdPropagator [source]
Compute a TebdPropagator object for a given SystemChain.
- Parameters
system_chain (SystemChain) – A SystemChain object that encodes the nearest neighbor Liouvillians.
time_step (float) – The time step length of the full TEBD propagator.
epsrel (float) – The relative singular value truncation tolerance.
order (int) – The expansion order.
- Returns
tebd_propagator – The TEBD Propagator.
- Return type
- oqupy.mps_mpo.compute_trotter_layers(nn_full_liouvillians: List[ndarray], hs_dims: List[int], dt: float, epsrel: float) Tuple[GateLayer, GateLayer] [source]
Compute even and odd Trotter layers.
- Parameters
nn_full_liouvillians (List[ndarrays]) – Full list of nearest neighbor Liouvillians.
hs_dims (List[int]) – Hilbert space dimensions of the chain sites.
dt (float) – The time step length.
epsrel (float) – The relative singular value truncation tolerance.
- Returns
gate_layer_even (GateLayer) – Gate layer with nearest neighbor gates with left sites having even indices.
gate_layer_odd (GateLayer) – Gate layer with nearest neighbor gates with left sites having odd indices.
oqupy.operators
Shorthand for commonly used operators.
- oqupy.operators.acommutator(operator: ndarray) ndarray [source]
Construct anti-commutator superoperator from operator.
- oqupy.operators.commutator(operator: ndarray) ndarray [source]
Construct commutator superoperator from operator.
- oqupy.operators.create(n: int) ndarray [source]
Bosonic creation operator of dimension n x n.
- Parameters
n (int) – Dimension of the Hilbert space.
- Returns
create – Creation operator matrix of dimension n x n.
- Return type
ndarray
- oqupy.operators.cross_acommutator(operator_1: ndarray, operator_2: ndarray) ndarray [source]
Construct anit-commutator of cross term (acting on two Hilbert spaces).
- oqupy.operators.cross_commutator(operator_1: ndarray, operator_2: ndarray) ndarray [source]
Construct commutator of cross term (acting on two Hilbert spaces).
- oqupy.operators.cross_left_right_super(operator_1_l: ndarray, operator_1_r: ndarray, operator_2_l: ndarray, operator_2_r: ndarray) ndarray [source]
Construct anit-commutator of cross term (acting on two Hilbert spaces).
- oqupy.operators.destroy(n: int) ndarray [source]
Bosonic annihilation operator of dimension n x n.
- Parameters
n (int) – Dimension of the Hilbert space.
- Returns
create – Annihilation operator matrix of dimension n x n.
- Return type
ndarray
- oqupy.operators.identity(n: int) ndarray [source]
Identity matrix of dimension n x n.
- Parameters
n (int) – Dimension of the square matrix.
- Returns
identity – Identity matrix of dimension n x n.
- Return type
ndarray
- oqupy.operators.left_right_super(left_operator: ndarray, right_operator: ndarray) ndarray [source]
Construct left and right acting superoperator from operators.
- oqupy.operators.left_super(operator: ndarray) ndarray [source]
Construct left acting superoperator from operator.
- oqupy.operators.preparation(density_matrix: ndarray) ndarray [source]
Construct the super operator that prepares the state.
- oqupy.operators.right_super(operator: ndarray) ndarray [source]
Construct right acting superoperator from operator.
oqupy.process_tensor
Module for the process tensor (PT) object. This code is based on [Pollock2018].
[Pollock2018] F. A. Pollock, C. Rodriguez-Rosario, T. Frauenheim, M. Paternostro, and K. Modi, Non-Markovian quantumprocesses: Complete framework and efficient characterization, Phys. Rev. A 97, 012127 (2018).
- class oqupy.process_tensor.BaseProcessTensor(hilbert_space_dimension: int, dt: Optional[float] = None, transform_in: Optional[ndarray] = None, transform_out: Optional[ndarray] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
ToDo
- property dt
ToDo.
- get_bond_dimensions() ndarray [source]
Return the bond dimensions of the MPS form of the process tensor.
- property hilbert_space_dimension
ToDo.
- property max_step: Union[int, float]
ToDo
- property transform_in
ToDo.
- property transform_out
ToDo.
- class oqupy.process_tensor.FileProcessTensor(mode: str, filename: Optional[str] = None, hilbert_space_dimension: Optional[int] = None, dt: Optional[float] = None, transform_in: Optional[ndarray] = None, transform_out: Optional[ndarray] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseProcessTensor
ToDo
- class oqupy.process_tensor.SimpleProcessTensor(hilbert_space_dimension: int, dt: Optional[float] = None, transform_in: Optional[ndarray] = None, transform_out: Optional[ndarray] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseProcessTensor
ToDo
- class oqupy.process_tensor.TrivialProcessTensor(hilbert_space_dimension: Optional[int] = 1, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseProcessTensor
ToDo
- oqupy.process_tensor.import_process_tensor(filename: str, process_tensor_type: Optional[str] = None) BaseProcessTensor [source]
ToDo.
oqupy.pt_tebd
Module for the process tensor approach to time evolving block decimation. This module is based on [Fux2022].
[Fux2022] G. E. Fux, D. Kilda, B. W. Lovett, and J. Keeling, Thermalization of a spin chain strongly coupled to its environment, arXiv:2201.05529 (2022).
- class oqupy.pt_tebd.PtTebd(initial_augmented_mps: AugmentedMPS, system_chain: SystemChain, process_tensors: List[Optional[BaseProcessTensor]], parameters: PtTebdParameters, chain_control: Optional[ChainControl] = None, start_time: Optional[float] = 0.0, start_step: Optional[int] = 0, dynamics_sites: Optional[List[Union[int, tuple]]] = None, backend_config: Optional[Dict] = None)[source]
Process tensor time evolving block decimation (PT-TEBD).
Backend configuration backend_config may have the following options:
‘parallel’ : ‘multiprocess’ / ‘multithread’
- Parameters
initial_augmented_mps (AugmentedMPS) – Initial augmented MPS.
system_chain (SystemChain) – Object encoding the system chain Liouvillians.
process_tensors (List[BaseProcessTensor]) – List of process tensors, one for each site. If a process tensor is ‘None’ it is assumed to be a trivial process tensor.
parameters (PtTebdParameters) – PT-TEBD computation parameters.
chain_control (ChainControl) – Optional control operations.
start_time (float) – Optional starting time stamp.
start_step (int) – Optional starting time step
dynamics_sites (List[Union[int,tuple]]) – Optional list of single sites or multiple site dynamics to be recorded.
backend_config (dict) – Optional backend configuration dictionary.
- property chain_control: ChainControl
The chain control object.
- compute(end_step: int, progress_type: Optional[str] = None) Dict [source]
Perform the PT-TEBD propagation up to time step ‘end_step’.
- Parameters
end_step (int) – The time step to which the propagation should be carried out.
progress_type (Text) – The progress report type during the computation. Types are: {
'silent'
,'simple'
,'bar'
}. If None then the default progress type is used.
- get_augmented_mps() AugmentedMPS [source]
Returns a copy of the current AugmentedMPS.
- get_current_density_matrix(sites: Union[int, tuple])[source]
Get the current density matrix of the site(s) ‘sites’.
- Parameters
sites (Union[int, tuple]) – The site(s).
- property step: int
The current step in the PT-TEBD computation.
- class oqupy.pt_tebd.PtTebdParameters(dt: float, epsrel: Optional[float] = 1e-05, order: Optional[int] = 2, name: Optional[str] = None, description: Optional[str] = None)[source]
Parameters for the process tensor time evolving block decimation computation.
- Parameters
dt (float) – Length of a time step \(\delta t\). - It should be small enough such that a trotterisation between the system Hamiltonian and the environment it valid, and the environment auto-correlation function is reasonably well sampled.
order (int) – Time evoling block decimation (TEBD) Trotterization order.
epsrel (float) – The maximal relative error in the singular value truncation (done in the underlying TEBD tensor network algorithm). - It must be small enough such that the numerical compression does not truncate relevant correlations.
name (str (default = None)) – An optional name for the PT-TEBD parameters object.
description (str (default = None)) – An optional description of the PT-TEBD parameters object.
- property dt: float
Length of a time step.
- property epsrel: float
The maximal relative error in the singular value truncation.
- property order: float
Time evoling block decimation (TEBD) Trotterization order.
oqupy.system
Module on physical information of the system.
- class oqupy.system.BaseSystem(dimension: int, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Base class for systems.
- property dimension: ndarray
Hilbert space dimension of the system.
- class oqupy.system.System(hamiltonian: ndarray, gammas: Optional[List[float]] = None, lindblad_operators: Optional[List[ndarray]] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseSystem
Represents a system (without any coupling to a non-Markovian bath). It is possible to include Lindblad terms in the master equation. The equations of motion for a system density matrix (without any coupling to a non-Markovian bath) is then:
\[\begin{split}\frac{d}{dt}\rho(t) = &-i [\hat{H}, \rho(t)] \\ &+ \sum_n^N \gamma_n \left( \hat{A}_n \rho(t) \hat{A}_n^\dagger - \frac{1}{2} \hat{A}_n^\dagger \hat{A}_n \rho(t) - \frac{1}{2} \rho(t) \hat{A}_n^\dagger \hat{A}_n \right)\end{split}\]with hamiltionian \(\hat{H}\), the rates gammas \(\gamma_n\) and linblad_operators \(\hat{A}_n\).
- Parameters
hamiltonian (ndarray) – System-only Hamiltonian \(\hat{H}\).
gammas (List(float)) – The rates \(\gamma_n\).
lindblad_operators (list(ndarray)) – The Lindblad operators \(\hat{A}_n\).
name (str) – An optional name for the system.
description (str) – An optional description of the system.
- property gammas: List[float]
List of gammas.
- property hamiltonian: ndarray
The system Hamiltonian.
- property lindblad_operators: List[ndarray]
List of lindblad operators.
- liouvillian() ndarray [source]
Returns the Liouvillian super-operator \(\mathcal{L}\) with
\[\mathcal{L}\rho = -i [\hat{H}, \rho] + \sum_n^N \gamma_n \left( \hat{A}_n \rho \hat{A}_n^\dagger - \frac{1}{2} \hat{A}_n^\dagger \hat{A}_n \rho - \frac{1}{2} \rho \hat{A}_n^\dagger \hat{A}_n \right) .\]- Returns
liouvillian – Liouvillian \(\mathcal{L}\).
- Return type
ndarray
- class oqupy.system.SystemChain(hilbert_space_dimensions: List[int], name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Represents a 1D chain of systems with nearest neighbor interactions.
- Parameters
hilbert_space_dimensions (List[int]) – Hilbert space dimension for each chain site.
name (str) – An optional name for the system chain.
description (str) – An optional description of the system chain.
- add_nn_dissipation(site: int, lindblad_operator_l: ndarray, lindblad_operator_r: ndarray, gamma: Optional[float] = 1.0) None [source]
Add two site lindblad dissipator
\[\mathcal{L} \rho_{n,n+1} = \gamma \left( \hat{A} \rho_{n,n+1} \hat{A}^\dagger - \frac{1}{2} \hat{A}^\dagger \hat{A} \rho_{n,n+1} - \frac{1}{2} \rho_{n,n+1} \hat{A}^\dagger \hat{A} \right)\]where \(\hat{A}=\hat{A}_l\otimes\hat{A}_r\), with site \(n\), lindblad_operator_l \(\hat{A}_l\), lindblad_operator_r \(\hat{A}_r\), and gamma \(\gamma\).
- Parameters
site (int) – Index of the left site \(n\).
lindblad_operator_l (ndarray) – Lindblad dissipator acting on the left site \(n\).
lindblad_operator_r (ndarray) – Lindblad dissipator acting on the right site \(n+1\).
gamma (float) – Optional multiplicative factor \(\gamma\).
- add_nn_hamiltonian(site: int, hamiltonian_l: ndarray, hamiltonian_r: ndarray) None [source]
Add a hamiltonian term to the Liouvillian of two neighboring sites:
\[\mathcal{L} \rho_{n,n+1} = -i [\hat{H}_l \otimes \hat{H}_r, \rho_{n,n+1}]\]with site \(n\), hamiltonian_l \(\hat{H}_l\) and hamiltonian_r \(\hat{H}_r\).
- Parameters
site (int) – Index of the left site \(n\).
hamiltonian_l (ndarray) – Hamiltonian acting on the left site \(n\).
hamiltonian_r (ndarray) – Hamiltonian acting on the right site \(n+1\).
- add_nn_liouvillian(site: int, liouvillian_l_r: ndarray) None [source]
Add Liouvillian of for the two neighboring sites site and site +1.
- Parameters
site (int) – Index of the left site \(n\).
liouvillian_l_r (ndarray) – Liouvillian acting on sites \(n\) and \(n+1\).
- add_site_dissipation(site: int, lindblad_operator: ndarray, gamma: Optional[float] = 1.0) None [source]
Add single site lindblad dissipator
\[\mathcal{L} \rho_n = \gamma \left( \hat{A} \rho_n \hat{A}^\dagger - \frac{1}{2} \hat{A}^\dagger \hat{A} \rho_n - \frac{1}{2} \rho_n \hat{A}^\dagger \hat{A} \right)\]with site \(n\), lindblad_operator \(\hat{A}\), and gamma \(\gamma\).
- Parameters
site (int) – Index of the site.
lindblad_operator (ndarray) – Lindblad dissipator acting on the single site.
gamma (float) – Optional multiplicative factor \(\gamma\).
- add_site_hamiltonian(site: int, hamiltonian: ndarray) None [source]
Add a hamiltonian term to a single site Liouvillian
\[\mathcal{L} \rho_n = -i [\hat{H}, \rho_n]\]with site \(n\) and hamiltonian \(\hat{H}\).
- Parameters
site (int) – Index of the site.
hamiltonian (ndarray) – Hamiltonian acting on the single site.
- add_site_liouvillian(site: int, liouvillian: ndarray) None [source]
Add a single site Liouvillian.
- Parameters
site (int) – Index of the site.
liouvillian (ndarray) – Liouvillian acting on the single site.
- get_nn_full_liouvillians() List[ndarray] [source]
Return the list of nearest neighbor Liouvillians (incorporating single site terms).
- property hs_dims
Hilbert space dimension for each chain site.
- property nn_liouvillians
The nearest neighbor Liouvillians.
- property site_liouvillians
The single site Liouvillians.
- class oqupy.system.TimeDependentSystem(hamiltonian: Callable[float, ndarray], gammas: Optional[List[Callable[float, float]]] = None, lindblad_operators: Optional[List[Callable[float, ndarray]]] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseSystem
Represents an explicitly time dependent system (without any coupling to a non-Markovian bath). It is possible to include (also explicitly time dependent) Lindblad terms in the master equation. The equations of motion for a system density matrix (without any coupling to a non-Markovian bath) is then:
\[\begin{split}\frac{d}{dt}\rho(t) = &-i [\hat{H}(t), \rho(t)] \\ &+ \sum_n^N \gamma_n(t) \left( \hat{A}_n(t) \rho(t) \hat{A}_n(t)^\dagger - \frac{1}{2} \hat{A}_n^\dagger(t) \hat{A}_n(t) \rho(t) - \frac{1}{2} \rho(t) \hat{A}_n^\dagger(t) \hat{A}_n(t) \right)\end{split}\]with the time dependent hamiltionian \(\hat{H}(t)\), the time dependent rates gammas \(\gamma_n(t)\) and the time dependent linblad_operators \(\hat{A}_n(t)\).
- Parameters
hamiltonian (callable) – System-only Hamiltonian \(\hat{H}(t)\).
gammas (List(callable)) – The rates \(\gamma_n(t)\).
lindblad_operators (list(callable)) – The Lindblad operators \(\hat{A}_n(t)\).
name (str) – An optional name for the system.
description (str) – An optional description of the system.
- property gammas: List[Callable[float, float]]
List of gammas.
- property hamiltonian: Callable[float, ndarray]
The system Hamiltonian.
- property lindblad_operators: List[Callable[float, ndarray]]
List of lindblad operators.
- liouvillian(t: float) ndarray [source]
Returns the Liouvillian super-operator \(\mathcal{L}(t)\) with
\[\mathcal{L}(t)\rho = -i [\hat{H}(t), \rho] + \sum_n^N \gamma_n \left( \hat{A}_n(t) \rho \hat{A}_n^\dagger(t) - \frac{1}{2} \hat{A}_n^\dagger(t) \hat{A}_n(t) \rho - \frac{1}{2} \rho \hat{A}_n^\dagger(t) \hat{A}_n(t) \right),\]with time \(t\).
- Parameters
t (float (default = None)) – time \(t\).
- Returns
liouvillian – Liouvillian \(\mathcal{L}(t)\) at time \(t\).
- Return type
ndarray
- class oqupy.system.TimeDependentSystemWithField(hamiltonian: Callable[[float, complex], numpy.ndarray], field_eom: Callable[[float, ndarray, complex], complex], gammas: Optional[List[Callable[float, float]]] = None, lindblad_operators: Optional[List[Callable[float, ndarray]]] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseSystem
Represents an explicitly time dependent system (without any coupling to a non-Markovian bath) with a coherent field (a complex scalar) \(\langle a \rangle\) that evolves according to a specified equation of motion \(\partial_t\langle a \rangle\).
It is possible to include (also explicitly time dependent) Lindblad terms in the master equation. The equations of motion for a system density matrix (without any coupling to a non-Markovian bath) is then:
\[\begin{split}\frac{d}{dt}\rho(t) = &-i [\hat{H}(t, \langle a \rangle), \rho(t)] \\ &+ \sum_n^N \gamma_n(t) \left( \hat{A}_n(t) \rho(t) \hat{A}_n(t)^\dagger - \frac{1}{2} \hat{A}_n^\dagger(t) \hat{A}_n(t) \rho(t) - \frac{1}{2} \rho(t) \hat{A}_n^\dagger(t) \hat{A}_n(t) \right)\end{split}\]with the hamiltionian \(\hat{H}(t, \langle a \rangle)\) depending on both time \(t\) and field \(\langle a \rangle\) (in general), the time dependent rates gammas \(\gamma_n(t)\) and the time dependent linblad_operators \(\hat{A}_n(t)\).
- Parameters
hamiltonian (callable) – System-only Hamiltonian \(\hat{H}(t, \langle a \rangle)\) where \(\langle a \rangle\) is the field at time \(t\).
field_eom (callable) – Field equation of motion \(\partial_t \langle a \rangle(t, \rho, \langle a \rangle)\) where \(\rho\) and \(\langle a \rangle\) are the system density matrix (a square matrix) and field at time \(t\).
gammas (list(callable)) – The rates \(\gamma_n(t)\).
lindblad_operators (list(callable)) – The Lindblad operators \(\hat{A}_n(t)\).
name (str) – An optional name for the system.
description (str) – An optional description of the system.
- property field_eom: Callable[[float, float, ndarray, complex], complex]
The field equation of motion.
- property gammas: List[Callable[float, float]]
List of gammas.
- property hamiltonian: Callable[[float, complex], numpy.ndarray]
The system Hamiltonian.
- property lindblad_operators: List[Callable[float, ndarray]]
List of lindblad operators.
- liouvillian(t0: float, t: float, state: ndarray, field: complex) ndarray [source]
Returns the Liouvillian super-operator \(\mathcal{L}(t)\) such that
\[\mathcal{L}(t)\rho = -i [\hat{H}(t, \langle a \rangle), \rho] + \sum_n^N \gamma_n \left( \hat{A}_n(t) \rho \hat{A}_n^\dagger(t) - \frac{1}{2} \hat{A}_n^\dagger(t) \hat{A}_n(t) \rho - \frac{1}{2} \rho \hat{A}_n^\dagger(t) \hat{A}_n(t) \right),\]with time \(t\).
- Parameters
t0 (float) – Start time of the current step
t (float) – Current time \(t\)
state (ndarray) – Flattened system density matrix (vector) at time \(t\)
field (complex) – Field value at time \(t\) obtained from the linearisation of the field at \(t\) using the field equation of motion.
- Returns
liouvillian – Liouvillian \(\mathcal{L}(t)\) at time \(t\).
- Return type
ndarray
oqupy.pt_tempo
Module for the process tensor time evolving matrix product operator algorithm (PT-TEMPO). This module is based on [Strathearn2018], [Pollock2018], [Jorgensen2019], and [Fux2021].
[Strathearn2018] A. Strathearn, P. Kirton, D. Kilda, J. Keeling and B. W. Lovett, Efficient non-Markovian quantum dynamics using time-evolving matrix product operators, Nat. Commun. 9, 3322 (2018).
[Pollock2018] F. A. Pollock, C. Rodriguez-Rosario, T. Frauenheim, M. Paternostro, and K. Modi, Non-Markovian quantumprocesses: Complete framework and efficient characterization, Phys. Rev. A 97, 012127 (2018).
[Jorgensen2019] M. R. Jørgensen and F. A. Pollock, Exploiting the causal tensor network structure of quantum processes to efficiently simulate non-markovian path integrals, Phys. Rev. Lett. 123, 240602 (2019)
[Fux2021] G. E. Fux, E. Butler, P. R. Eastham, B. W. Lovett, and J. Keeling, Efficient exploration of Hamiltonian parameter space for optimal control of non-Markovian open quantum systems, Phys. Rev. Lett. 126, 200401 (2021).
- class oqupy.pt_tempo.PtTempo(bath: Bath, start_time: float, end_time: float, parameters: TempoParameters, process_tensor_file: Optional[Union[str, bool]] = None, overwrite: Optional[bool] = False, backend_config: Optional[Dict] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Class to facilitate a PT-TEMPO computation.
- Parameters
bath (Bath) – The Bath (includes the coupling operator to the system).
parameters (TempoParameters) – The parameters for the PT-TEMPO computation.
start_time (float) – The start time.
backend_config (dict (default = None)) – The configuration of the backend. If backend_config is
None
then the default backend configuration is used.name (str (default = None)) – An optional name for the tempo object.
description (str (default = None)) – An optional description of the tempo object.
- compute(progress_type: Optional[str] = None) None [source]
Propagate (or continue to propagate) the TEMPO tensor network to time end_time.
- Parameters
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.
- property dimension: ndarray
Hilbert space dimension.
- get_process_tensor(progress_type: Optional[str] = None) BaseProcessTensor [source]
Returns a the computed process tensor. It performs the computation if it hasn’t been already done.
- Parameters
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
silent
,simple
,bar
}. If None then the default progress type is used.- Returns
process_tensor – The computed process tensor.
- Return type
- oqupy.pt_tempo.pt_tempo_compute(bath: Bath, start_time: float, end_time: float, parameters: Optional[TempoParameters] = None, tolerance: Optional[float] = 0.0039, process_tensor_file: Optional[Union[str, bool]] = None, overwrite: Optional[bool] = False, backend_config: Optional[Dict] = None, progress_type: Optional[str] = None, name: Optional[str] = None, description: Optional[str] = None) BaseProcessTensor [source]
Shortcut for creating a process tensor by performing a PT-TEMPO computation.
- Parameters
bath (Bath) – The Bath (includes the coupling operator to the system).
start_time (float) – The start time.
end_time (float) – The time to which the PT-TEMPO should be computed.
parameters (TempoParameters) – The parameters for the PT-TEMPO computation.
tolerance (float) – Tolerance for the parameter estimation (only applicable if parameters is None).
backend_config (dict (default = None)) – The configuration of the backend. If backend_config is
None
then the default backend configuration is used.progress_type (str (default = None)) – The progress report type during the computation. Types are: {
'silent'
,'simple'
,'bar'
}. If None then the default progress type is used.name (str (default = None)) – An optional name for the tempo object.
description (str (default = None)) – An optional description of the tempo object.
oqupy.tempo
Module on for the original time evolving matrix product operator (TEMPO) algorithm. This module is based on [Strathearn2017] and [Strathearn2018].
[Strathearn2017] A. Strathearn, B.W. Lovett, and P. Kirton, Efficient real-time path integrals for non-Markovian spin-boson models. New Journal of Physics, 19(9), p.093009 (2017).
[Strathearn2018] A. Strathearn, P. Kirton, D. Kilda, J. Keeling and B. W. Lovett, Efficient non-Markovian quantum dynamics using time-evolving matrix product operators, Nat. Commun. 9, 3322 (2018).
- class oqupy.tempo.BaseTempo(bath: Bath, parameters: TempoParameters, initial_state: ndarray, start_time: float, backend_config: Optional[Dict] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Base class for all TEMPO objects.
- Parameters
bath (Bath) – The Bath (includes the coupling operator to the system).
parameters (TempoParameters) – The parameters for the TEMPO computation.
initial_state (ndarray) – The initial density matrix of the system.
start_time (float) – The start time.
backend_config (dict (default = None)) – The configuration of the backend. If backend_config is
None
then the default backend configuration is used.name (str (default = None)) – An optional name for the tempo object.
description (str (default = None)) – An optional description of the tempo object.
- property dimension: ndarray
Hilbert space dimension.
- class oqupy.tempo.Tempo(system: Union[System, TimeDependentSystem], bath: Bath, parameters: TempoParameters, initial_state: ndarray, start_time: float, backend_config: Optional[Dict] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseTempo
Class representing the entire TEMPO tensornetwork as introduced in [Strathearn2018].
- Parameters
system (System or TimeDependentSystem) – The system.
bath (Bath) – The Bath (includes the coupling operator to the system).
parameters (TempoParameters) – The parameters for the TEMPO computation.
initial_state (ndarray) – The initial density matrix of the system.
start_time (float) – The start time.
backend_config (dict (default = None)) – The configuration of the backend. If backend_config is
None
then the default backend configuration is used.name (str (default = None)) – An optional name for the tempo object.
description (str (default = None)) – An optional description of the tempo object.
- compute(end_time: float, progress_type: Optional[str] = None) Dynamics [source]
Propagate (or continue to propagate) the TEMPO tensor network to time end_time.
- Parameters
end_time (float) – The time to which the TEMPO should be computed.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
'silent'
,'simple'
,'bar'
}. If None then the default progress type is used.
- Returns
dynamics – The instance of Dynamics associated with the TEMPO object.
- Return type
- class oqupy.tempo.TempoParameters(dt: float, dkmax: int, epsrel: float, add_correlation_time: Optional[float] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseAPIClass
Parameters for the TEMPO computation.
- Parameters
dt (float) – Length of a time step \(\delta t\). - It should be small enough such that a trotterisation between the system Hamiltonian and the environment it valid, and the environment auto-correlation function is reasonably well sampled.
dkmax (int) – Number of time steps \(K\in\mathbb{N}\) that should be included in the non-Markovian memory. - It must be large enough such that \(\delta t \times K\) is larger than the necessary memory time \(\tau_\mathrm{cut}\).
epsrel (float) – The maximal relative error in the singular value truncation (done in the underlying tensor network algorithm). - It must be small enough such that the numerical compression (using tensor network algorithms) does not truncate relevant correlations.
add_correlation_time (float) – Additional correlation time to include in the last influence functional as explained in [Strathearn2017].
name (str (default = None)) – An optional name for the tempo parameters object.
description (str (default = None)) – An optional description of the tempo parameters object.
- property add_correlation_time: float
Additional correlation time to include in the last influence functional.
- property dkmax: float
Number of time steps that should be included in the non-Markovian memory.
- property dt: float
Length of a time step.
- property epsrel: float
The maximal relative error in the singular value truncation.
- class oqupy.tempo.TempoWithField(system: TimeDependentSystemWithField, bath: Bath, parameters: TempoParameters, initial_state: ndarray, initial_field: complex, start_time: float, subdiv_limit: Optional[int] = 256, epsrel: Optional[float] = 1.4901161193847656e-08, backend_config: Optional[Dict] = None, name: Optional[str] = None, description: Optional[str] = None)[source]
Bases:
BaseTempo
Class representing the TEMPO tensornetwork with coherent field evolution as introduced in [FowlerWright2021].
- Parameters
system (TimeDependentSystemWithField) – The (time-dependent) system with a coherent field.
bath (Bath) – The Bath (includes the coupling operator to the system).
parameters (TempoParameters) – The parameters for the TEMPO computation.
initial_state (ndarray) – The initial density matrix of the system.
initial_field (complex) – The initial field value.
start_time (float) – The start time.
subdiv_limit (int (default = config.SUBDIV_LIMIT)) – The maximum number of subdivisions used during the adaptive algorithm when integrating the system Liouvillian. If None then the Liouvillian is not integrated but sampled twice to to construct the system propagators at each timestep.
epsrel (float (default = config.INTEGRATE_EPSREL)) – The relative error tolerance for the adaptive algorithm when integrating the system Liouvillian.
backend_config (dict (default = None)) – The configuration of the backend. If backend_config is
None
then the default backend configuration is used.name (str (default = None)) – An optional name for the tempo object.
description (str (default = None)) – An optional description of the tempo object.
- compute(end_time: float, progress_type: Optional[str] = None) DynamicsWithField [source]
Propagate (or continue to propagate) the TEMPO tensor network and coherent field to time end_time.
- Parameters
end_time (float) – The time to which the TEMPO should be computed.
progress_type (str (default = None)) – The progress report type during the computation. Types are: {
'silent'
,'simple'
,'bar'
}. If None then the default progress type is used.
- Returns
dynamics – The instance of DynamicsWithField associated with the TempoWithField object.
- Return type
DynamicsWithFields
- get_dynamics() DynamicsWithField [source]
Returns DynamicsWithField instance associated with the Tempo object.
- oqupy.tempo.guess_tempo_parameters(bath: Bath, start_time: float, end_time: float, system: Optional[BaseSystem] = None, tolerance: Optional[float] = 0.0039) TempoParameters [source]
Function to roughly estimate appropriate parameters for a TEMPO computation.
Warning
No guarantee that resulting TEMPO calculation converges towards the correct dynamics! Please refer to the TEMPO documentation and check convergence by varying the parameters for TEMPO manually.
- Parameters
bath (Bath) – The bath.
start_time (float) – The start time.
end_time (float) – The time to which the TEMPO should be computed.
system (BaseSystem) – The system.
tolerance (float) – Tolerance for the parameter estimation.
- Returns
tempo_parameters – Estimate of appropriate tempo parameters.
- Return type
- oqupy.tempo.influence_matrix(dk: int, parameters: TempoParameters, correlations: BaseCorrelations, coupling_acomm: ndarray, coupling_comm: ndarray)[source]
Compute the influence functional matrix.
- oqupy.tempo.tempo_compute(system: BaseSystem, bath: Bath, initial_state: ndarray, start_time: float, end_time: float, parameters: Optional[TempoParameters] = None, tolerance: Optional[float] = 0.0039, backend_config: Optional[Dict] = None, progress_type: Optional[str] = None, name: Optional[str] = None, description: Optional[str] = None) Dynamics [source]
Shortcut for creating a Tempo object and running the computation.
- Parameters
system (BaseSystem) – The system.
bath (Bath) – The Bath (includes the coupling operator to the system).
initial_state (ndarray) – The initial density matrix of the system.
start_time (float) – The start time.
end_time (float) – The time to which the TEMPO should be computed.
parameters (TempoParameters) – The parameters for the TEMPO computation.
tolerance (float) – Tolerance for the parameter estimation (only applicable if parameters is None).
backend_config (dict (default = None)) – The configuration of the backend. If backend_config is
None
then the default backend configuration is used.progress_type (str (default = None)) – The progress report type during the computation. Types are: {
'silent'
,'simple'
,'bar'
}. If None then the default progress type is used.name (str (default = None)) – An optional name for the tempo object.
description (str (default = None)) – An optional description of the tempo object.
oqupy.util
Module for utilities.
- class oqupy.util.ProgressBar(max_value, title=None)[source]
Class to display the computation progress with a nice progress bar.
- class oqupy.util.ProgressSilent(max_value, title=None)[source]
Class NOT to display the computation progress.
- class oqupy.util.ProgressSimple(max_value, title=None)[source]
Class to display the computation progress step by step.
- oqupy.util.add_singleton(tensor: ndarray, index: Optional[int] = - 1, copy: Optional[bool] = True) ndarray [source]
Add a singleton to a numpy tensor.
- oqupy.util.check_convert(variable: Any, conv_type: Any, name: Optional[str] = None, msg: Optional[str] = None)[source]
Attempt to convert variable into a specific type.
- oqupy.util.check_isinstance(variable: Any, types: Any, name: Optional[str] = None, msg: Optional[str] = None)[source]
Check that a variable is an instance of one of the given types.
- oqupy.util.check_true(expr: bool, msg: Optional[str] = None)[source]
Check that an specific expression is true.
- oqupy.util.create_delta(tensor: ndarray, index_scrambling: List[int]) ndarray [source]
Creates deltas in numpy tensor.
Warning
This is a computationally inefficient method to perform the task.
Todo
Make it better.
- oqupy.util.get_progress(progress_type: Optional[str] = None) BaseProgress [source]
Get a progress class from the progress_type.