py-tedopa Documentation

Module contents

tedopa.tedopa module

This module comprises functions for simulating the time-evolution of one- and two-site open quantum systems in a bosonic environment via tedopa1() and tedopa2() and helper functions.

Todo

Decide if better to convert the arguments of tedopa1 and tedopa2 to kwargs.

tedopa.tedopa._get_annihilation_op(dim)[source]

Creates the annihilation operator

Parameters:dim (int) – Dimension of the site it should act on
Returns:The annihilation operator
Return type:numpy.ndarray
tedopa.tedopa._get_parameters(n, j, domain, g, ncap)[source]

Calculate the parameters needed for mapping the Hamiltonian to a 1D chain

Parameters:
  • n (int) – Number of recurrence coefficients required (rc.recurrenceCoefficients() actually returns one more and the system site does not need one, so the argument n-2 is passed)
  • j (types.LambdaType) – Spectral function \(J(\omega)\) as defined in Chin et al. doi: 10.1063 / 1.3490188
  • domain (list[float]) – Domain on which \(J(\omega)\) is defined, for example [0, np.inf]
  • g (float) – Constant \(g\), assuming that for \(J(\omega)\) it is \(g(\omega)=g\omega\).
  • ncap (int) – Number internally used by py-orthpol to determine accuracy of the returned recurrence coefficients. Must be <= 60000, the higher the longer the calculation of the recurrence coefficients takes and the more accurate it becomes.
Returns:

omegas, ts, c0 as defined in the paper
Return type:

list[list[float], list[float], float]
tedopa.tedopa._get_singlesite_ops(h_loc, params, bs, b_daggers)[source]

List of the operators acting on every single site after chain mapping

Parameters:
  • h_loc (numpy.ndarray) – Local Hamiltonian
  • params (list) – Parameters as returned by _get_parameters()
  • bs (list) – The list of annihilation operators acting on each site of the chain
  • b_daggers (list) – The list of creation operators acting on each site of the chain
Returns:

List of operators acting on every single site
Return type:

list
tedopa.tedopa._get_twosite_ops(a, params, bs, b_daggers)[source]

List of the operators acting on every two adjacent sites after chain mapping

Parameters:
  • a (numpy.ndarray) – Interaction operator provided by the user
  • params (list) – Parameters as returned by _get_parameters()
  • bs (list) – The list of annihilation operators acting on each site of the chain
  • b_daggers (list) – The list of creation operators acting on each site of the chain
Returns:

List of operators acting on every two adjacent sites
Return type:

list
tedopa.tedopa.get_times(ts_full, ts_system, len_state, sys_position, sys_length)[source]

This is a function specifically designed for TEDOPA systems. It calculates the proper ‘ts’ and ‘subsystems’ input lists for tmps.evolve() from a list of times where the full state shall be returned and a list of times where only the reduced state of the system in question shall be returned. ts then basically is a concatenation of ts_full and ts_system, while subsystems will indicate that at the respective time in ts either the full state or only a reduced density matrix should be returned.

Parameters:
  • ts_full (list[float]) – List of times where the full state including environment chain should be returned
  • ts_system (list[float]) – List of times where only the reduced density matrix of the system should be returned
  • len_state (int) – The length of the state
  • sys_position (int) – The position of the system (first site would be 0)
  • sys_length (int) – Length of the system, i.e. number of sites the system is comprised of
Returns:

Times and subsystems in the form that has to be provided to tmps.evolve()
Return type:

tuple(list[float], list[list[int]])
tedopa.tedopa.map(h_loc, a, state_shape, j, domain, g, ncap)[source]

Map the Hamiltonian of one site coupled to bosonic environmentself.

This function calculates the operators acting on every single site of the resulting chain and the operators acting on every two adjacent sites in the chain from the local Hamiltonian and the spectral density.

The mapping is based on an algorithm introduced by Chin et al. in Journal of Mathematical Physics 51, 092109 (2010); doi: 10.1063/1.3490188.

Parameters:
  • h_loc (numpy.ndarray) – Matrix representation of the local Hamiltonian of the one-site system.
  • a (numpy.ndarray) – Interaction operator. This is the site-part of the tensor product that comprises the interaction Hamiltonian and is defined as \(\hat{A}\) in Chin et al. doi: 10.1063 / 1.3490188
  • state_shape (list[list[int]]) – The shape of the chain on which the Hamiltonian is to be applied. For example [[3, 3], [2, 2], [2, 2]] for a system comprised of 3 sites, where the first one has 2 physical legs each of dimension 3, the second has 2 physical legs each of dimension 2 and so on. This is the typical MPA shape used by mpnum (state.shape if state is an mpnum.MPArray).
  • j (types.LambdaType) – Spectral function \(J(\omega)\) as defined in Chin et al. doi: 10.1063 / 1.3490188
  • domain (list[float]) – Domain on which \(J(\omega)\) is defined, for example [0, np.inf]
  • g (float) – Constant \(g\), assuming that for \(J(\omega)\) it is \(g(\omega)=g\omega\).
  • ncap (int) – Number internally used by py-orthpol.
Returns:

Terms of the effective Hamiltonian acting on the chain after the chain mapping. Two lists, one with the single-site operators and the other with adjacent-site operators that act on two sites. See the input parameters of tmps.evolve()
Return type:

list[list[numpy.ndarray]]
tedopa.tedopa.tedopa1(h_loc, a, state, method, j, domain, ts_full, ts_system, trotter_compr=None, compr=None, g=1, trotter_order=2, num_trotter_slices=100, ncap=20000, v=0)[source]

TEDOPA for a single site coupled to a bosonic bath.

This function proceeds in two steps: (i) Map the Hamiltonian from a system composed of one site that is linearly coupled to a reservoir of bosonic modes with a given spectral function to the same site coupled to a 1D chain of bosonic modes and (ii) perform time evolution.

The performed mapping is based on an algorithm introduced by Chin et al. in Journal of Mathematical Physics 51, 092109 (2010); doi: 10.1063/1.3490188.

The inputs to this function include the initial state, the local Hamiltonian, the spectral density, time of evolution, query times and some performance parameters. The output provides the MPAs of the evolved state at the requested times.

Parameters:
  • h_loc (numpy.ndarray) – Matrix representation of the local Hamiltonian of the one-site system.
  • a (numpy.ndarray) – Interaction operator. This is the site-part of the tensor product that comprises the interaction Hamiltonian and is defined as \(\hat{A}\) in Chin et al.
  • state (mpnum.MPArray) – The initial state of the system which is to be evolved.
  • method (str) – The parameterization of the intial state. Either ‘mps’, ‘mpo’ or ‘pmps’. Determines which method is used in the simulation.
  • j (types.LambdaType) – Spectral function \(J(\omega)\) as defined in Chin et al.
  • domain (list[float]) – Domain on which \(J(\omega)\) is defined, for example [0, np.inf]
  • ts_full (list[float]) – The times for which the evolution should be computed and the whole state chain returned.
  • ts_system (list[float]) – The times for which the evolution should be computed and the reduced density matrix of only the system should be returned.
  • trotter_compr (dict) – Compression parameters used in the iterations of Trotter (in the form required by mpnum.MPArray.compress(). If unsure, look at https://github.com/dseuss/mpnum/blob/master/examples/mpnum_intro.ipynb .) If omitted, some default compression will be used that will probably work but might lead to problems. See tmps._set_compr_params() for more information.
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation, where trotter_compr is used. compr = dict(method=’svd’, rank=10) would for example ensure that the ranks of any MPA never exceed 10 during all of the calculations. An accepted relative error for the compression can be provided in addition to or instead of ranks, which would lead to e.g. compr = dict(method=’svd’, rank=10, relerr=1e-12). If omitted, some default compression will be used that will probably work but might lead to problems. See tmps._set_compr_params() for more information.
  • g (float) – Cutoff \(g\), assuming that for \(J(\omega)\) it is \(g(\omega)=g\omega\).
  • trotter_order (int) – Order of Trotter - Suzuki decomposition to be used. Currently only 2 and 4 are implemented
  • num_trotter_slices (int) – Number of Trotter slices to be used for the largest t in ts_full or ts_system. If ts_system=[10, 25, 30] and num_trotter_slices=100, then the program would use 100/30*10=33, 100/30*25=83 and 100/30*30=100 Trotter slices to calculate the time evolution for the three times.
  • ncap (int) – Number internally used by py-orthpol to determine accuracy of the returned recurrence coefficients. Must be <= 60000, the higher the longer the calculation of the recurrence coefficients takes and the more accurate it becomes.
  • v (int) – Level of verbose output. 0 means no output, 1 means that some basic output showing the progress of calculations is produced. 2 will in addition show the bond dimensions of the state after every couple of iterations, 3 will show bond dimensions after every Trotter iteration.
Returns:

The first list is an array of times for which the states are actually computed (might differ slightly from the times in ts_full and ts_system, since the ones in times have to be multiples of tau). The second list is an array of the corresponding evolved states.
Return type:

list[list[float], list[mpnum.MPArray]]
tedopa.tedopa.tedopa2(h_loc, a_twosite, state, method, sys_position, js, domains, ts_full, ts_system, trotter_compr=None, compr=None, gs=(1, 1), trotter_order=2, num_trotter_slices=100, ncap=20000, v=0)[source]

TEDOPA for two sites coupled to each other and each individually to its own bosonic bath.

This function proceeds in two steps: (i) Map the Hamiltonian from a system composed of two coupled sites that are each linearly coupled to their possibly distinct reservoir of bosonic modes with a given spectral function to a 1D chain representing the whole setup and (ii) perform time evolution.

The performed mapping is based on an algorithm introduced by Chin et al. in Journal of Mathematical Physics 51, 092109 (2010); doi: 10.1063/1.3490188.

The inputs to this function include the initial state, the local Hamiltonian, the spectral densities, time of evolution, query times and some performance parameters. The output provides the MPAs of the evolved state at the requested times.

The first elements in the lists js, domains, etc. always refer to the first (left) site and the second elements in the lists refer to the second (right) site of the system

Parameters:
  • h_loc (numpy.ndarray) – Matrix representation of the local Hamiltonian of the two-site system.
  • a_twosite (list[numpy.ndarray]) – List of two matrices, each of which represents the site-part of the tensor product interaction Hamiltonian for the two sites.
  • state (mpnum.MPArray) – The initial state of the system which is to be evolved.
  • method (str) – The parameterization of the intial state. Either ‘mps’, ‘mpo’ or ‘pmps’. Determines which method is used in the simulation.
  • sys_position (int) – Which index, in the chain representing the state, is the position of the first site of the system (starting at 0). E.g. sys_position = 2 if the chain-length is is 6 and sites are of the form env-env-sys-sys-env-env.
  • js (list[types.LambdaType]) – Spectral functions \(J(\omega)\) for the two environments as defined by Chin et al.
  • domains (list[list[float]]) – Domains on which the \(J(\omega)\) are defined. Can be different for the two sites, for example, [[0, np.inf], [0,1]]
  • ts_full (list[float]) – The times for which the evolution should be computed and the whole state chain returned.
  • ts_system (list[float]) – The times for which the evolution should be computed and the reduced density matrix of only the system should be returned.
  • trotter_compr (dict) – Compression parameters used in the iterations of Trotter (in the form required by mpnum.MPArray.compress(). If unsure, look at https://github.com/dseuss/mpnum/blob/master/examples/mpnum_intro.ipynb .) If omitted, some default compression will be used that will probably work but might lead to problems. See tmps._set_compr_params() for more information.
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation, where trotter_compr is used. compr = dict(method=’svd’, rank=10) would for example ensure that the ranks of any MPA never exceed 10 during all of the calculations. An accepted relative error for the compression can be provided in addition to or instead of ranks, which would lead to e.g. compr = dict(method=’svd’, rank=10, relerr=1e-12). If omitted, some default compression will be used that will probably work but might lead to problems. See tmps._set_compr_params() for more information.
  • gs (list[float]) – List of cutoffs \(g\), assuming that for \(J(\omega)\) it is \(g(\omega)=g\omega\).
  • trotter_order (int) – Order of Trotter-Suzuki decomposition to be used. Currently only 2 and 4 are implemented
  • num_trotter_slices (int) – Number of Trotter slices to be used for the largest t in ts_full or ts_system. If ts_system=[10, 25, 30] and num_trotter_slices=100, then the program would use 100/30*10=33, 100/30*25=83 and 100/30*30=100 Trotter slices to calculate the time evolution for the three times.
  • ncap (int) – Number internally used by py-orthpol to determine accuracy of the returned recurrence coefficients. Must be <= 60000, the higher the longer the calculation of the recurrence coefficients takes and the more accurate it becomes.
  • v (int) – Level of verbose output. 0 means no output, 1 means that some basic output showing the progress of calculations is produced. 2 will in addition show the bond dimensions of the state after every couple of iterations, 3 will show bond dimensions after every Trotter iteration.
Returns:

The first list is an array of times for which the states are actually computed (might differ slightly from the times in ts_full and ts_system, since the ones in times have to be multiples of tau). The second list is an array of the corresponding evolved states.
Return type:

list[list[float], list[mpnum.MPArray]]

tedopa.tedopa_models module

This module contains functions which enable the user to run TEDOPA for certain settings without having to know mpnum. This is a wrapper for the tedopa/tedopa.py module.

TEDOPA proceeds in two steps: (i) Map the Hamiltonian from a system, composed of one or two sites that are each linearly coupled to a distinct reservoir of bosonic modes with a given spectral function, to the same site/sites coupled to a 1D chain of oscillators (one chain per system site) and (ii) perform time evolution.

The performed mapping is based on an algorithm introduced by Chin et al. in Journal of Mathematical Physics 51, 092109 (2010); doi: 10.1063/1.3490188.

If something is unclear in the functions here, please look at tedopa/tedopa.py or of course the example notebooks for deeper understanding.

Todo

Might be good to have a thermal state of a given Hamiltonian. => make a function create_bosonic_thermal_state()

tedopa.tedopa_models.calculate_expectation_values(states, observable)[source]

Calculates the expectation values \(\langle M \rangle_i\) of the the observable \(M\) with respect to the states \(\{\rho\}\)

\[\langle M \rangle_i = \text{tr}(\rho_i M)\]

For this function to work, the observable has to have the same dimension as the density matrix which represents the state, i.e. all states must have the same dimensions. This function is meant to work with a list of evolved states of the same system at different times.

Parameters:
  • states (list[mpnum.MPArray]) – List of states \(\{\rho\}\). They are assumed to be MPOs and already normalized
  • observable (numpy.ndarray) – The matrix representing the observable \(M\) in global form (as opposed to local form. Global form is just the usual form to write matrices in QM. For more information, see the mpnum documentation.)
Returns:

List of expectation values for the states
Return type:

list[mpnum.MPArray.dtype]
tedopa.tedopa_models.create_bosonic_vacuum_state(system_site_states, len_chain, dim_oscillators)[source]

This function will create the state of a system, comprised of one or two sites, each connected to a (possibly distinct) environment of bosonic modes in vacuum state. It is assumed that the state is not entangled and the single sites of the system are each in a pure state. It will be represented as an MPS rather than an MPO.

Parameters:
  • system_site_states (list[numpy.ndarray]) – List of the vectors describing the state of the system site or sites, i.e. [\(|\psi_1 \rangle\)] or [\(|\psi_1 \rangle\), \(|\psi_2 \rangle\)].
  • len_chain (list[int]) – List specifying the length of the chains representing the environment of the respective system site.
  • dim_oscillators (list[int]) – List specifying the dimensions of the oscillators in the respective chain.
Returns:

The vacuum state in MPS form.
Return type:

mpnum.MPArray
tedopa.tedopa_models.tedopa1_for_bosonic_vacuum_state(system_site_state, len_chain, dim_oscillators, h_loc, a, j, domain, ts, observable, g=1, trotter_order=2, num_trotter_slices=100, v=1)[source]

This function will take information about a system of only one site, the coupling to its environment and some simulation parameters, and performs TEDOPA on it. The initial state of the environment will be set to a bosonic vacuum state. The function also takes an observable and will return the expectation values for that observable with respect to the reduced density matrices of the evolved one-site system.

Parameters:
  • system_site_state (numpy.ndarray) – The vector \(|\psi \rangle\) describing the state of the one-site system.
  • len_chain (int) – The length of the chain representing the environment of the system site.
  • dim_oscillators (int) – The dimension of the oscillators in the chain. When choosing e.g. 6, each of the len_chain number of oscillators in the chain will have a 6x6 density matrix
  • h_loc (numpy.ndarray) – Matrix representation of the local Hamiltonian of the one-site system.
  • a (numpy.ndarray) – Interaction operator. This is the site-part of the tensor product that comprises the interaction Hamiltonian and is defined as \(\hat{A}\) in Chin et al.
  • j (types.LambdaType) – Spectral function \(J(\omega)\) as defined in Chin et al. doi: 10.1063 / 1.3490188
  • domain (list[float]) – Domain on which \(J(\omega)\) is defined, for example [0, np.inf]
  • ts (list[float]) – The times for which the evolved states, and based on that the expectation values of the observable, should be computed.
  • observable (numpy.ndarray) – The observable used to determine the expectation values at the times in ts for the reduced density matrices of the evolved one-site system. It must be of the same dimension as the density matrix of the system in question.
  • g (float) – Cutoff \(g\), assuming that for \(J(\omega)\) it is \(g(\omega)=g\omega\).
  • trotter_order (int) – Order of Trotter - Suzuki decomposition to be used. Currently only 2 and 4 are implemented
  • num_trotter_slices (int) – Number of Trotter slices to be used for the largest t in ts. If ts=[10, 25, 30] and num_trotter_slices=100, then the program would use 100/30*10=33, 100/30*25=83 and 100/30*30=100 Trotter slices to calculate the time evolution for the three times.
  • v (int) – Level of verbose output. 0 means no output, 1 will show the progress of the calculations. 2 and 3 will also indicate bond dimensions in the evolved state. For more information see tedopa.tedopa1()
Returns:

A list of times corresponding to the respective expectation values and a list of the expectation values of the provided observable for evolved states
Return type:

tuple[list[float], list[mpnum.MPArray.dtype]]
tedopa.tedopa_models.tedopa2_for_bosonic_vacuum_state(system_site_state, len_chain, dim_oscillators, h_loc, a_twosite, js, domains, ts, observable, gs=1, trotter_order=2, num_trotter_slices=100, v=0)[source]

This function will take information about a system comprised of two not entangled sites, the coupling to their two environments and some simulation parameters, and performs TEDOPA on it. The initial state of both environments will be set to a bosonic vacuum state. The function also takes an observable and will return the expectation values for that observable with respect to the reduced density matrices of the evolved two-site system.

Parameters:
  • system_site_state (list[numpy.ndarray]) – List of the vectors describing the states of the system sites, i.e. [\(|\psi_1 \rangle\), \(|\psi_2 \rangle\)].
  • len_chain (list[int]) – The lengths of the two chains representing the environments of the two system sites.
  • dim_oscillators (list[int]) – The dimensions of the oscillators in the two chains. When choosing e.g. [6, 8], each of the len_chain[0] number of oscillators in the first chain will have a 6x6 density matrix, and each of the oscillators in the second chain 8x8 density matrices.
  • h_loc (numpy.ndarray) – Matrix representation of the local Hamiltonian of the two-site system.
  • a_twosite (list[numpy.ndarray]) – List of two matrices, each of which represents the site-part of the tensor product interaction Hamiltonian for the two sites. See \(\hat{A}\) in Chin et al. doi: 10.1063 / 1.3490188
  • js (list[types.LambdaType]) – Spectral functions \(J(\omega)\) for the two environments as defined by Chin et al.
  • domains (list[list[float]]) – Domains on which the \(J(\omega)\) are defined. Can be different for the two sites, for example, [[0, np.inf], [0,1]]
  • ts (list[float]) – The times for which the evolved states, and based on that the expectation values of the observable, should be computed.
  • observable (numpy.ndarray) – The observable used to determine the expectation values at the times in ts for the reduced density matrices of the evolved two-site system. It must be of the same dimension as the density matrix of the system in question.
  • gs (list[float]) – List of cutoffs \(g\), assuming that for \(J(\omega)\) it is \(g(\omega)=g\omega\).
  • trotter_order (int) – Order of Trotter - Suzuki decomposition to be used. Currently only 2 and 4 are implemented
  • num_trotter_slices (int) – Number of Trotter slices to be used for the largest t in ts. If ts=[10, 25, 30] and num_trotter_slices=100, then the program would use 100/30*10=33, 100/30*25=83 and 100/30*30=100 Trotter slices to calculate the time evolution for the three times.
  • v (int) – Level of verbose output. 0 means no output, 1 will show the progress of the calculations. 2 and 3 will also indicate bond dimensions in the evolved state. For more information see tedopa.tedopa2()
Returns:

A list of times corresponding to the respective expectation values and a list of the expectation values of the provided observable for evolved states
Return type:

tuple[list[float], list[mpnum.MPArray.dtype]]

tedopa._recurrence_coefficients module

This module contains functions to calculate recurrence coefficients, based on the package py-orthpol.

tedopa._recurrence_coefficients._j_to_hsquared(func, lb, rb, g)[source]

Transform spectral density J to square of coupling function h(x)

This transforms the given spectral density \(J(\omega)\) to the square \(h^2(x)\) of the coupling function assuming a linear dispersion relation \(g(x)=gx\).

Parameters:
  • func (types.LambdaType) – J(omega)
  • lb (float) – Left boundary
  • rb (float) – Right boundary
  • g (float) – Constant slope of the linear dispersion relations
Returns:

lb, rb, \(h^2\), where lb and rb are the new left and right boundaries for \(h^2\).
Return type:

list[float, float, types.LambdaType]
tedopa._recurrence_coefficients.recurrenceCoefficients(n, lb, rb, j, g, ncap=60000)[source]

Calculate recurrence coefficients for given spectral density

Recurrence coeffcients for an arbitrary measure are defined as follows. Given some measure \(d\mu(x)\), which defines the set \(\{\pi_n(x)\in \mathbb{P}_n,n=0,1,2,\ldots\}\) of monic orthogonal polynomials with respect to the measure, the following recurrence relation holds.

\[\pi_{k+1}(x)=(x- \alpha_k)\pi_k(x)-\beta_k\pi_{k-1}(x), \quad k=0,1,2...\]

where \(\pi_{-1}(x)\equiv 0\), and the recurrence coefficients are:

\[\alpha_k=\frac{\langle x \pi_k,\pi_k\rangle}{\langle \pi_k,\pi_k\rangle}, \quad \beta_k=\frac{\langle \pi_k,\pi_k\rangle}{\langle \pi_{k-1},\pi_{k-1}\rangle}.\]

The TEDOPA mapping for a given spectral density relies on calculating the recurrence coefficients with respect to the measure \(d\mu(x)=h^2(x)dx\), where \(J(\omega)=\pi h^2[g^{-1}(\omega)] \frac{dg^{-1}(\omega)}{d\omega}\) and \(g(x) = gx\). Thus, this function first calculates the function \(h^2(x)\) from \(J(\omega)\) and then calls py-orthpol package to find the recurrence coefficients.

For more details, see Journal of Mathematical Physics 51, 092109 (2010); doi: 10.1063/1.3490188 and ACM Trans. Math. Softw. 20, 21-62 (1994); doi: 10.1145/174603.174605

Note that the input j must be a python lambda function representing the spectral density \(J(\omega)\)

Parameters:
  • n (int) – Number of recurrence coefficients required.
  • g (float) – Constant g, assuming that for J(omega) it is g(omega)=g*omega.
  • lb (float) – Left bound of interval on which J is defined.
  • rb (float) – Right bound of interval on which J is defined.
  • j (types.LambdaType) – \(J(\omega)\) defined on the interval (lb, rb)
  • ncap (int) – Number internally used by py-orthpol to determine the accuracy with which the recurrence coefficients are calculated. Must be >n and <=60000. Between 10000 and 60000 recommended, the higher the number the higher the accuracy and the longer the execution time. (Default value = 60000)
Returns:

A list of two lists, each with the respective recurrence coefficients \(\{\alpha_i: i = 1,2,\dots,n\}\) and \(\{\beta_i: i = 1,2,\dots,n\}\) defined above.
Return type:

list[list[float], list[float]]

tedopa.tmps module

This module provides functions for time evolution of a state given as MPS, MPO or PMPS via the tMPS algorithm.

This is based on functions which calculate the time evolution of an operator in MPS, MPO or PMPS form from Hamiltonians acting on every single and every two adjacent sites.

tMPS is a method to evolve a one dimensional quantum state, represented as an MPS, MPO or PMPS, in time. It requires that the Hamiltonian is only comprised of terms acting on single sites or pairs of adjacent sites of the state. This allows the Hamiltonian \(H\) to be written as

\[H = \sum_j h_{j, j+1},\]

where \(j\) is the index of the respective site in the state. These components can be grouped into those acting on even and those acting on odd sites, leading to a time evolution operator

\[U(\tau) = \text{e}^{\mathrm{i}(H_{\text{even}}+H_{\text{ odd}})\tau},\]

with

\[H_{\text{even}} = \sum_{j\text{ even}} h_{j, j+1}\]
\[H_{\text{odd}} = \sum_{j\text{ odd}} h_{j, j+1}\]

This allows to perform Trotter-Suzuki decompositions of \(U\), for example of second order:

\[U(\tau) = \text{e}^{\mathrm{i} H_{\text{odd}} \tau/2} \text{e}^{\mathrm{i} H_{\text{even}} \tau} \text{e}^{\mathrm{i} H_{\text{odd}} \tau/2} + \mathcal{O}(\tau^3).\]

These decompositions provide the advantage that \(U\) does not need to be calculated as a whole matrix, which could potentially become way too big. Since the elements within \(H_{\text{even}}\) and those within \(H_{\text{odd}}\) commute, \(U\) can be broken up into smaller pieces which a computer can handle even for very large systems.

For more information, see chapter 7 in Schollwöck’s paper Annals of Physics 326, 96-192 (2011); doi: 10.1016/j.aop.2010.09.012

In this file, evolve() is the main function to be called to evolve a state in time. It will itself call _trotter_slice() which will call _trotter_two() or _trotter_four() to calculate the \(U( \tau)\) representing one Trotter slice. When that is done, evolve() will take it and pass it on to _time_evolution() which will then go through the Trotter iterations, thus actually evolving the state in time, and store the requested results on the way.

tedopa.tmps._append(times, states, compr_errors, trot_errors, tau, i, j, step_numbers, subsystems, state, accumulated_overlap, accumulated_trotter_error, method)[source]

Function to append time evolved state etc to output of _time_evolution()

Parameters:
  • times (list[float]) – List containing the times to which the states are evolved
  • states (list[mpnum.MPArray]) – List containing the evolved states
  • compr_errors (list[float]) – List containing the respective compression errors
  • trot_errors (list[float]) – List containing the respective Trotter errors
  • tau (float) – The time of one Trotter slice
  • i (int) – Number indicating which is the current Trotter slice
  • j (int) – Number indicating how many times a state related to the current i has been appended already
  • step_numbers (list[int]) – List containing the time steps
  • subsystems (list[list[int]]) – List of sites for which the subsystem should be returned at the respective time
  • state (mpnum.MPArray) – The current state
  • accumulated_overlap (float) – The accumulated overlap error
  • accumulated_trotter_error (float) – The accumulated Trotter error
  • method (str) – Method to use as defined in evolve()
Returns:

Nothing, changes happen in place
Return type:

None
tedopa.tmps._get_h_list(hamiltonians, num_sites)[source]

Convert given list of Hamiltonians into form suitable for exponentiation.

If only one Hamiltonian acting on every single site and one acting on every two adjacent sites is given, transform it into the form returned. If not, check whether the lengths of the lists match the number of sites.

Parameters:
  • hamiltonians (list) – Hamiltonians as in evolve()
  • num_sites (int) – Number of sites of the state to be evolved
Returns:

A list with two items: The first is a list of Hamiltonians acting on the single sites, like [h1, h2, h3, …] and the second is a list of Hamiltonians acting on each two adjacent sites, like [h12, h23, h34, …]
Return type:

list[list[numpy.ndarray], list[numpy.ndarray]]
tedopa.tmps._get_subsystems_list(subsystems, len_step_numbers)[source]

This function just brings subsystems, which indicates which subsystem should be returned at the respective step number, in the right form.

Parameters:
  • subsystems (list) – Same as that described in evolve()
  • len_step_numbers (int) – The length of the array containing the step numbers for which the evolved state is to be stored.
Returns:

Sites for which the subsystem should be returned at the respective time,
Return type:

list[list[int]]
tedopa.tmps._get_u_list_even(dims, h_single, h_adjacent, tau)[source]

Calculates individual operator exponentials of adjacent even-odd sites, i.e. transforms \(\{h_{j,j+1} : j \text{ even}\}\) into \(\{ \text{e}^{\mathrm{i} h_{j,j+1} \tau} : j \text{ even}\}\)

Parameters:
  • dims (list) – The dimensions of the single sites of the state U will be applied to
  • h_single (list) – The Hamiltonians acting on every single site
  • h_adjacent (list) – The Hamiltonians acting on every two adjacent sites
  • tau (float) – The time step for the time evolution of U
Returns:

List of operators acting on even adjacent sites, like [u23, u45, …]
Return type:

list[numpy.ndarray]
tedopa.tmps._get_u_list_odd(dims, h_single, h_adjacent, tau)[source]

Calculates individual operator exponentials of adjacent odd-even sites, i.e. transforms \(\{h_{j, j+1} : j \text{ odd}\}\) into \(\{ \text{e}^{\mathrm{i} h_{j,j+1} \tau} : j \text{ odd}\}\)

>>> dims = [2, 2, 2, 2]
>>> sx = np.array([[0, 1], [1, 0]])
>>> sz = np.array([[1, 0], [0, -1]])
>>> tau = 1
>>> actual_result = _get_u_list_odd(dims, [sx] * 4, [np.kron(sz, sz)] * 3, tau)
>>> expected_result = [expm(-1j * tau * (np.kron(sz,sz) + .5 * (np.kron(sx, np.identity(2)) + np.kron(np.identity(2), sx))))] * 2
>>> print(np.array_equal(expected_result, actual_result))
True
Parameters:
  • dims (list) – The dimensions of the single sites of the state U will be applied to
  • h_single (list) – The Hamiltonians acting on every single site
  • h_adjacent (list) – The Hamiltonians acting on every two adjacent sites
  • tau (float) – The time step for the time evolution of U
Returns:

List of operators acting on odd adjacent sites, like [u12, u34, …]
Return type:

list[numpy.ndarray]
tedopa.tmps._set_compr_params()[source]

A function to set default compression parameters if none were provided. They will suffice in many cases, but might well lead to problems described in the Pitfalls section of the introduction notebook. If that is the case, try to find and provide your own suitable compression parameters, which is also recommended to have more control over the calculations and their precision. For more information on this, read the introduction notebook and make use of the verbose output option to monitor bond dimensions during calculations.

Returns:Some default compression and Trotter compression parameters
Return type:list[dict]
tedopa.tmps._time_evolution(state, us, step_numbers, subsystems, tau, method, trotter_compr, v)[source]

Implements time-evolution via Trotter-Suzuki decomposition

Parameters:
  • state (mpnum.MPArray) – The state to be evolved in time
  • us (list[mpnum.MPArray]) – List of ordered operator exponentials for a single Trotter slice
  • step_numbers (list[int]) – List of time steps as generated by _times_to_steps()
  • subsystems (list[list[int]]) – Sites for which the subsystem states should be returned at the respective times
  • tau (float) – Duration of one Trotter slice. As defined in _times_to_steps()
  • method (str) – Which method to use. Either ‘mps’, ‘mpo’ or ‘pmps’.
  • trotter_compr (dict) – Compression parameters used in the iterations of Trotter-Suzuki decomposition.
  • v (int) – Level of verbose output. 0 means no output, 1 means that some basic output showing the progress of calculations is produced. 2 will in addition show the bond dimensions of the state after every couple of iterations, 3 will show bond dimensions after every Trotter iteration.
Returns:

A list with five items: (i) The list of times for which the density matrices have been computed (ii) The list indicating which subsystems of the system are returned at the respective time of the first list (iii) The list of density matrices as MPO or PMPS as mpnum.MPArray, depending on the input “method”. If that was MPS, the full states will still be MPSs, the reduced ones will be MPOs.
Return type:

list[list[float], list[list[int]], list[mpnum.MPArray]
tedopa.tmps._times_to_steps(ts, num_trotter_slices)[source]

Based on the requested times ts, calculate Trotter step numbers at which (subsystems of) evolved states need to be saved.

>>> _times_to_steps([10, 25, 30], 100)
([33, 83, 100], 0.3)
>>> _times_to_steps([8, 26, 19], 80)
([25, 80, 58], 0.325)
Parameters:
  • ts (list[float]) – The times for which the evolution should be computed and the state of the full system or a subsystem returned (i.e. it’s reduced density matrix (for now only works with method=’mpo’. If the method is not mpo, omit subsystems)). The algorithm will calculate the evolution using the given number of Trotter steps for the largest number in ts. On the way there it will store the evolved states for smaller times. NB: Beware of memory overload since len(t) mpnum.MPArrays will be stored
  • num_trotter_slices (int) – Number of Trotter slices to be used for the largest t in ts.
Returns:

step numbers, tau = maximal t / num_trotter_slices
Return type:

tuple[list[int], float]
tedopa.tmps._trotter_four(hamiltonians, tau, num_sites, compr)[source]

Get a list of ordered operator exponentials for one fourth-order Trotter slice.

Based on the description in the documentation of _trotter_slice() and on the paper by Schollwöck, \(N\) = 11, with

\[\tau_1 = \tau_{11} = \frac{\tau}{2(4 - 4^{1/3})},\]
\[\tau_2 = \tau_3 = \tau_4 = \tau_8 = \tau_9 = \tau_{10} = \frac{\tau}{4 - 4^{1 / 3}},\]
\[\tau_5 = \tau_7 = \frac{\tau (1 - 3)}{2(4 - 4^{1 / 3})},\]

and

\[\tau_6 = \frac{\tau (1 - 4)}{4 - 4^{1 / 3}}.\]
Parameters:
  • hamiltonians (list) – List of two lists of Hamiltonians, the Hamiltonians in the first acting on every single site, the Hamiltonians in the second acting on every pair of two adjacent sites
  • tau (float) – As defined in _times_to_steps()
  • num_sites (int) – Number of sites of the state to be evolved
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation where trotter_compr is used
Returns:

The time evolution operator parts, which, applied one after another, give one Trotter slice
Return type:

list[mpnum.MPArray]
tedopa.tmps._trotter_slice(hamiltonians, tau, num_sites, trotter_order, compr)[source]

Get a list of ordered operator exponentials for one Trotter slice.

The Trotter-Suzuki formula approximates the time-evolution during a single Trotter slice

\[U(\tau) = \text{e}^{\mathrm{i}\sum_{j=1}^m H_j \tau},\]

with

\[{U}^\prime(\tau) =\prod_{p=1}^N U_p,\]

which is a product of \(N\) operator exponentials

\[U_p := {\text{e}^{H_{j_p}\tau_p}}\]

of \(H_j\). Here \(\{\tau_p\}\) is a sequence of real numbers such that \(\sum_p \tau_p = \tau\). The \(H_{j_p}\) for a certain \(p\) are all elements of the Hamiltonian acting either on even or on odd pairs of adjacent sites. This ensures that within one \(U_p\) all terms in the exponential commute.

This function returns the list of operators \(U_p\) as MPOs.

For more information on Trotter-Suzuki, see chapter 7 in Schollwöck’s paper Annals of Physics 326, 96-192 (2011); doi: 10.1016/j.aop.2010.09.012.

Parameters:
  • hamiltonians (list) – List of two lists of Hamiltonians, the Hamiltonians in the first acting on every single site, the Hamiltonians in the second acting on every pair of two adjacent sites
  • tau (float) – As defined in _times_to_steps()
  • num_sites (int) – Number of sites of the state to be evolved
  • trotter_order (int) – Order of Trotter-Suzuki decomposition to be used
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation where trotter_compr is used
Returns:

The time evolution operator parts, which, applied one after another, give one Trotter slice
Return type:

list[mpnum.MPArray]
tedopa.tmps._trotter_two(hamiltonians, tau, num_sites, compr)[source]

Get a list of ordered operator exponentials for one second-order Trotter slice.

Based on the description in the documentation of _trotter_slice() and on the paper by Schollwöck, \(N\) = 3, with \(\tau_1 =\tau_3 = \tau/2\) and \(\tau_2=\tau\).

Parameters:
  • hamiltonians (list) – List of two lists of Hamiltonians, the Hamiltonians in the first acting on every single site, the Hamiltonians in the second acting on every pair of two adjacent sites
  • tau (float) – As defined in _times_to_steps()
  • num_sites (int) – Number of sites of the state to be evolved
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation where trotter_compr is used
Returns:

The time evolution operator parts, which, applied one after another, give one Trotter slice
Return type:

list[mpnum.MPArray]
tedopa.tmps._u_list_to_mpo_even(dims, u_even, compr)[source]

Transforms list of matrices on even-odd sites to MPO acting on full state. So the list of u_even \(\{u_{j,j+1} : j \text{ even}\}\), which are numpy.ndarrays, is transformed into \(\bigotimes_j u_{j,j+1} : j \text{ even}\) of the type mpnum.MPArray.

Parameters:
  • dims (list) – List of dimensions of each site
  • u_even (list) – List of time evolution operators acting on even adjacent sites
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation where trotter_compr is used
Returns:

The MPO for the full state acting on even-odd adjacent sites
Return type:

mpnum.MPArray
tedopa.tmps._u_list_to_mpo_odd(dims, u_odd, compr)[source]

Transforms list of matrices on odd-even sites to MPO acting on full state. So the list of u_odd \(\{u_{j,j+1} : j \text{ odd}\}\), which are numpy.ndarrays, is transformed into \(\bigotimes_j u_{j,j+1} : j \text{ odd}\) of the type mpnum.MPArray.

Parameters:
  • dims (list) – List of dimensions of each site
  • u_odd (list) – List of time evolution operators acting on odd adjacent sites
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation where trotter_compr is used
Returns:

The MPO for the full state acting on odd-even adjacent sites
Return type:

mpnum.MPArray
tedopa.tmps.evolve(state, hamiltonians, num_trotter_slices, method, trotter_order, ts, trotter_compr=None, compr=None, subsystems=None, v=0)[source]

Evolve a one dimensional MPS, MPO or PMPS state using tMPS as described in the module’s documentation.

The initial state, Hamiltonians and certain parameters are required. The output is a list of times and a list of the evolved states at these times. Those states might be subsystems of the whole evolved system, which allows for the user to keep memory consumption small by focusing on the subsystems of interest.

Todo

Raise exception if hamiltonians are not of the right dimension

Todo

Implement tracking of compression errors.

Todo

Get variable compression to work (might involve changing mpnum).

Parameters:
  • state (mpnum.MPArray) – The state to be evolved in time. The state has to be an MPS, MPO or PMPS, depending on which method is chosen
  • hamiltonians (list) – Either a list containing the Hamiltonian acting on every single site and the Hamiltonian acting on every two adjacents sites, like [H_i, H_ij], or a list containing a list of Hamiltonians acting on the single sites and a list of Hamiltonians acting on each two adjacent sites, like [[h1, h2, h3, ...], [h12, h23, h34, ...]]
  • num_trotter_slices (int) – Number of Trotter slices to be used for evolution over time equal to the largest t in ts.
  • method (str) – Which method to use. Either ‘mps’, ‘mpo’ or ‘pmps’.
  • trotter_order (int) – Order of Trotter-Suzuki decomposition to be used. Currently only 2 and 4 are implemented
  • ts (list[float]) – The times for which the evolution should be computed and the state of the full system or a subsystem returned (i.e. it’s reduced density matrix). The algorithm will calculate the evolution using the given number of Trotter steps for the largest number in ts. On the way there it will store the evolved states for smaller times. NB: Beware of memory overload since len(t) number of mpnum.MPArrays will be stored
  • trotter_compr (dict) – Compression parameters used in the iterations of Trotter (in the form required by mpnum.MPArray.compress(). If unsure, look at https://github.com/dseuss/mpnum/blob/master/examples/mpnum_intro.ipynb .) If omitted, some default compression will be used that will probably work but might lead to problems. See _set_compr_params() for more information.
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation, where trotter_compr is used. compr = dict(method=’svd’, rank=10) would for example ensure that the ranks of any MPA never exceed 10 during all of the calculations. An accepted relative error for the compression can be provided in addition to or instead of ranks, which would lead to e.g. compr = dict(method=’svd’, rank=10, relerr=1e-12). If omitted, some default compression will be used that will probably work but might lead to problems. See _set_compr_params() for more information.
  • subsystems (list) – A list defining for which subsystem the reduced density matrix or whether the full state should be returned for a time in ts. This can be a list of the length same as that of ts looking like [[a1, b1], [a2, b2], ...] or just a list like [a, b]. In the first case the respective subsystem for every entry in ts will be returned, in the second case the same subsystem will be returned for all entries in ts. [a, b] will lead to a return of the reduced density matrix of the sites from a up to, but not including, b. For example [0, 3] if the reduced density matrix of the first three sites should be returned. A time can occur twice in ts and then different subsystems to be returned can be defined for that same time. If this parameter is omitted, the full system will be returned for every time in ts.
  • v (int) – Level of verbose output. 0 means no output, 1 means that some basic output showing the progress of calculations is produced. 2 will in addition show the bond dimensions of the state after every couple of iterations, 3 will show bond dimensions after every Trotter iteration.
Returns:

A list with five items: (i) The list of times for which the density matrices have been computed (ii) The list indicating which subsystems of the system are returned at the respective time of the first list (iii) The list of density matrices as MPO or PMPS as mpnum.MPArray, depending on the input “method”. If that was MPS, the full states will still be MPSs, the reduced ones will be MPOs.
Return type:

list[list[float], list[list[int]], list[mpnum.MPArray]]
tedopa.tmps.matrix_to_mpo(matrix, shape, compr=None)[source]

Convert matrix to MPO

Converts given \(M \times N\) matrix in global form into an MPO with the given shape. The number of legs per site must be the same for all sites.

>>> matrix = np.array([[1, 0, 0], [0, 0, 0], [0, 0, 0]])
>>> mpo = matrix_to_mpo(matrix, [[3, 3]])
>>> print(mpo.to_array_global())
[[1 0 0]
 [0 0 0]
 [0 0 0]]
Parameters:
  • matrix (numpy.ndarray) – The matrix to be transformed to an MPO
  • shape (list) – The shape the single sites of the resulting MPO should have, as used in mpnum. For example three sites with two legs each might look like [[3, 3], [2, 2], [2, 2]]. Format same as numpy.ndarray.shape
  • compr (dict) – Parameters for the compression which is executed on every MPA during the calculations, except for the Trotter calculation where trotter_compr is used
Returns:

The MPO with shape shape representing the matrix
Return type:

mpnum.MPArray
tedopa.tmps.normalize(state, method)[source]

Normalize a state.

Parameters:
  • state (mpnum.MPArray) – The state to be normalized
  • method (str) – Whether it is a MPS, MPO or PMPS state
Returns:

The normalized state
Return type:

mpnum.MPArray