setupAggregate(agg; groundEnergy=true, verbose=false)

Generate all basic data from the Aggregate. Returns aggInds, vibindices, aggIndLen, basis, FCFact, FCProd, Ham, Ham_0, Ham_I.

evolutionApproximate!(op_array, op0, tspan, Hamiltonian)

Calculate approximate time evolution of the op0 inplace based on $U(t)$ that is calculated for $t_0$ and $t_\text{step}$. See evolutionOperator. This method returns Vector of Arrays.

evolutionApproximate!(ket_array, ket0, tspan, Hamiltonian)

Calculate approximate time evolution of the ket0 inplace based on $U(t)$ that is calculated for $t_0$ and $t_\text{step}$. See evolutionOperator. Argument ket_array is Vector of Arrays.

evolutionApproximate(ket0, tspan, Hamiltonian)

Calculate approximate time evolution of the ket0 based on $U(t)$ that is calculated for $t_0$ and $t_\text{step}$. See evolutionOperator. This method returns Vector of Ket states.

evolutionApproximate(op0, tspan, Hamiltonian)

Calculate approximate time evolution of the op0 based on $U(t)$ that is calculated for $t_0$ and $t_\text{step}$. See evolutionOperator. This method returns Vector of Operators.

evolutionExact!(op_array, op0, tspan, Hamiltonian; 
	diagonalize = true, approximate = false)

Calculate exact time evolution of the op0 state inplace see evolutionOperatorIterator. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. The approximate option asusmes that tspan is made of equidistant points, therefore $U(t)$ has to be calculated for $t_0$ and $t_\text{step}$.

evolutionExact!(ket_array, ket0, tspan, Hamiltonian; 
	diagonalize = true, approximate = false)

Calculate exact time evolution of the ket0 state inplace see evolutionOperatorIterator. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. The approximate option asusmes that tspan is made of equidistant points, therefore $U(t)$ has to be calculated for $t_0$ and $t_\text{step}$.

evolutionExact(ket0, tspan, Hamiltonian; diagonalize = true, approximate = false)

Calculate exact time evolution of the ket0 state see evolutionOperatorIterator. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. The approximate option asusmes that tspan is made of equidistant points, therefore $U(t)$ has to be calculated for $t_0$ and $t_\text{step}$.

evolutionExact(op0, tspan, Hamiltonian; diagonalize = true, approximate = false)

Calculate exact time evolution of the op0 state see evolutionOperatorIterator. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. The approximate option asusmes that tspan is made of equidistant points, therefore $U(t)$ has to be calculated for $t_0$ and $t_\text{step}$.

evolutionOperator(Hamiltonian, t)

Get evolution operator as Operator using the definition

$U(t) = e^{-i H t / \hbar}, \quad \hbar = 1$.

evolutionOperatorA(H_lambda, S, Sinv, t)

Get evolution operator as Array using the definition

$U(t) = S e^{-i H_\lambda t / \hbar} S^{-1},$

where $\hbar = 1$, $H_{\lambda,ii} = \lambda_i$ are eigenvalues of $H$, so $H$ has to be non-singular, otherwise $H_{\lambda,ij} = 0, i \neq j$. $S$ is obtained from eigendecomposition of $H$, for example

H_lambda, S = eigen(Hamiltonian.data)
Sinv = inv(S)
H_lambda = diagm(H_lambda)

and arguments have to be Arrays.

evolutionOperatorArray(Hamiltonian, tspan)

Get evolution operators as Vector or Operators using the definition

$U(t) = e^{-i H t / \hbar}, \quad \hbar = 1$.

evolutionOperatorIterator(Hamiltonian, tspan; diagonalize = true, approximate = false)

Resumable function that returns evolution operator as Operator type at the time t from tspan. See evolutionOperator. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. The approximate option asusmes that tspan is made of equidistant points, therefore $U(t)$ has to be calculated for $t_0$ and $t_\text{step}$.

evolutionSuperOperator(Hamiltonian, t)

Get evolution operator as SuperOperator using the definition

$\mathcal{U}(t) \:\cdot\: = U(t) \:\cdot\: U^\dagger(t) = e^{-i H t / \hbar} \:\cdot\: e^{i H t / \hbar}, \quad \hbar = 1$.

evolutionOperatorArray(Hamiltonian, tspan)

Get evolution superoperators as Vector or SuperOperators using the definition

$\mathcal{U}(t) \:\cdot\: = U(t) \:\cdot\: U^\dagger(t) = e^{-i H t / \hbar} \:\cdot\: e^{i H t / \hbar}, \quad \hbar = 1$.

evolutionSuperOperatorIterator(Hamiltonian, tspan; 
	diagonalize = true, approximate = false)

Resumable function that returns evolution operator as Operator type at the time t from tspan. See evolutionSuperOperator. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. The approximate option asusmes that tspan is made of equidistant points, therefore $U(t)$ has to be calculated for $t_0$ and $t_\text{step}$.

evolution_exact(rho0, tspan, Hamiltonian; diagonalize = false)

Calculate approximate time evolution of the rho0 based on $U(t)$ that is calculated for $t_0$ and $t_\text{step}$. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. See evolutionOperator. This method returns tspan and Vector of Operators.

evolution_exact(rho0, tspan, Hamiltonian; diagonalize = false)

Calculate exact time evolution of the rho0 inplace based on $U(t)$. The diagonalize argument decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition. See evolutionOperator. This method returns tspan and Vector of Operators.

thermal_state(T, mu_array, Ham, aggIndices; 
	boltzmann_const = 0.69503476, diagonalize = false, diagonal = false)

Get initial state as thermal state excited with ultra-fast laser pulse. In this version we suppose that after the thermal state is excited with laser pulse, the whole population of ground state is distributed over electric states in mu_array. We assume Condon approximation.

$\rho_\text{thermal} = \exp( -\frac{i}{\hbar} H ), \quad \hbar = 1$`.

Arguments

• T: Temperature of the initial thermal state.
• mu_array: Vector of electric states in local basis, see electronicIndices. The first index is for ground state of the aggregate the rest are first excited states in local basis.
• Ham: Arbitrary operator specifying the Hamiltonian.
• aggIndices: Aggregate indices, see getIndices.
• boltzmann_const: Boltzmann const in $\mathrm{cm^{-1}}$.
• diagonalize: Decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition.
• diagonal: Return only the diagonal part of the excited density matrix.

thermal_state_composite(T, mu_weighted, Ham, aggIndices; 
	boltzmann_const::Float64 = 0.69503476, diagonalize::Bool=false, diagonal=false)

Functionality of this method is similar to thermal_state, but the final state is constructed from partial thermal_states with weight specified in mu_weighted. For example

julia thermal_state_composite(T, [0.0, 0.8, 0.2], ...) = 0.8 * thermal_state(T, [1, 2, 1], ...) + 0.2 * thermal_state(T, [1, 1, 2], ...)`

$\rho_\text{thermal} = \exp( -\frac{i}{\hbar} H ), \quad \hbar = 1$.

Arguments

• T: Temperature of the initial thermal state.
• mu_weighted: Vector of weights of electric states in local basis, see electronicIndices.
• Ham: Arbitrary operator specifying the Hamiltonian.
• aggIndices: Aggregate indices, see getIndices.
• boltzmann_const: Boltzmann const in $\mathrm{cm^{-1}}$.
• diagonalize: Decompose Hamiltonian into $\lambda_i, S, S^{-1}$ and calculate the exponential using eigenvalue decomposition.
• diagonal: Return only the diagonal part of the excited density matrix.

liouvilleVonNeumann(rho0, tspan, H; 
	reltol=1.0e-12, abstol=1.0e-12, fout=nothing, alg=OrdinaryDiffEq.Tsit5())

Integrate Liouville-von Neumann equation to evolve states or compute propagators

$\frac{d}{d t} \rho(t) = - \frac{i}{\hbar} [ \hat H, \rho(t) ], \quad \hbar = 1$.

Arguments

• rho0: Initial state vector (can be a bra or a ket) or initial propagator.
• tspan: Vector specifying the points of time for which output should be displayed.
• H: Arbitrary operator specifying the Hamiltonian.
• reltol: Relative tolerance for OrdinaryDiffEq solver and its inner states.
• abstol: Absolute tolerance for OrdinaryDiffEq solver and its inner states.
• fout=nothing: If given, this function fout(t, rho) is called every time an output should be displayed. ATTENTION: The state rho is neither normalized nor permanent! It is still in use by the ode solver and therefore must not be changed.
• alg: Algorithm with which OrdinaryDiffEq will solve LvN equation.

MemoryKernel_1_traced(H_II_t, H_II_tau, W_bath, agg, FCProd, aggIndices, indicesMap)

Calculate the first part of Memory Kernel with the definition

$\mathcal{M}_1(t, \tau) = \operatorname{tr}_B \{ \hat{H}_I^{(I)}(t) \hat{H}_I^{(I)}(\tau) W_\text{bath} \}$.

Arguments

• H_II_t: Interaction Hamiltonian in interaction picutre at the time t, $\hat{H}_I^{(I)}(t)$.
• H_II_tau: Interaction Hamiltonian in interaction picutre at the time tau, $\hat{H}_I^{(I)}(\tau)$.
• W_bath: Density matrix representing bath part of the density matrix, see get_rho_bath.
• agg: Aggregate of molecules, see Aggregate.
• aggIndices: Aggregate indices, see getIndices.
• indicesMap: Aggregate vibrational indices, see getVibIndices.

MemoryKernel_2_traced(H_II_t, H_II_tau, W_bath, agg, FCProd, aggIndices, indicesMap)

Calculate the second part of Memory Kernel with the definition

$\mathcal{M}_1(t, \tau) = \operatorname{tr}_B \{ \hat{H}_I^{(I)}(t) W_\text{bath} \hat{H}_I^{(I)}(\tau) \}$.

Arguments

• H_II_t: Interaction Hamiltonian in interaction picutre at the time t, $\hat{H}_I^{(I)}(t)$.
• H_II_tau: Interaction Hamiltonian in interaction picutre at the time tau, $\hat{H}_I^{(I)}(\tau)$.
• W_bath: Density matrix representing bath part of the density matrix, see get_rho_bath.
• agg: Aggregate of molecules, see Aggregate.
• aggIndices: Aggregate indices, see getIndices.
• indicesMap: Aggregate vibrational indices, see getVibIndices.

MemoryKernel_3_traced(H_II_t, H_II_tau, W_bath, agg, FCProd, aggIndices, indicesMap)

Calculate the third part of Memory Kernel with the definition

$\mathcal{M}_1(t, \tau) = \operatorname{tr}_B \{ \hat{H}_I^{(I)}(\tau) W_\text{bath} \hat{H}_I^{(I)}(t) \}$.

Arguments

• H_II_t: Interaction Hamiltonian in interaction picutre at the time t, $\hat{H}_I^{(I)}(t)$.
• H_II_tau: Interaction Hamiltonian in interaction picutre at the time tau, $\hat{H}_I^{(I)}(\tau)$.
• W_bath: Density matrix representing bath part of the density matrix, see get_rho_bath.
• agg: Aggregate of molecules, see Aggregate.
• aggIndices: Aggregate indices, see getIndices.
• indicesMap: Aggregate vibrational indices, see getVibIndices.

MemoryKernel_4_traced(H_II_t, H_II_tau, W_bath, agg, FCProd, aggIndices, indicesMap)

Calculate the fourth part of Memory Kernel with the definition

$\mathcal{M}_1(t, \tau) = \operatorname{tr}_B \{ W_\text{bath} \hat{H}_I^{(I)}(t) \hat{H}_I^{(I)}(\tau) \}$.

Arguments

• H_II_t: Interaction Hamiltonian in interaction picutre at the time t, $\hat{H}_I^{(I)}(t)$.
• H_II_tau: Interaction Hamiltonian in interaction picutre at the time tau, $\hat{H}_I^{(I)}(\tau)$.
• W_bath: Density matrix representing bath part of the density matrix, see get_rho_bath.
• agg: Aggregate of molecules, see Aggregate.
• aggIndices: Aggregate indices, see getIndices.
• indicesMap: Aggregate vibrational indices, see getVibIndices.

MemoryKernel_traced(H_II_t, H_II_tau, W_bath, agg, FCProd, aggIndices, indicesMap)

Calculate Memory Kernel with the definition

$\mathcal{M}(t, \tau) = \operatorname{tr}_B \{ [ \hat{H}_I^{(I)}(t), [ \hat{H}_I^{(I)}(\tau), W_\text{bath} ]]\}$.

Arguments

• H_II_t: Interaction Hamiltonian in interaction picutre at the time t, $\hat{H}_I^{(I)}(t)$.
• H_II_tau: Interaction Hamiltonian in interaction picutre at the time tau, $\hat{H}_I^{(I)}(\tau)$.
• W_bath: Density matrix representing bath part of the density matrix, see get_rho_bath.
• agg: Aggregate of molecules, see Aggregate.
• aggIndices: Aggregate indices, see getIndices.
• indicesMap: Aggregate vibrational indices, see getVibIndices.

tracedistance(rho, sigma)

Trace distance between rho and sigma. It is defined as

It calls tracenorm which in turn either uses tracenorm_h or tracenorm_nh depending if $ρ-σ$ is hermitian or not.

tracedistance_h(rho, sigma)

Trace distance between rho and sigma. It uses the identity

where $λ_i$ are the eigenvalues of rho - sigma.

tracedistance_nh(rho, sigma)

Trace distance between rho and sigma. Note that in this case rho and sigma don't have to be represented by square matrices (i.e. they can have different left-hand and right-hand bases). It uses the identity

where $σ_i$ are the singular values of rho - sigma.

tracenorm(rho)

Trace norm of rho. It is defined as

Depending if rho is hermitian either tracenorm_h or tracenorm_nh is called.

tracenorm_h(rho)

Trace norm of rho. It uses the identity

where $λ_i$ are the eigenvalues of rho.

tracenorm_nh(rho)

Trace norm of rho. Note that in this case rho doesn't have to be represented by a square matrix (i.e. it can have different left-hand and right-hand bases). It uses the identity

where $σ_i$ are the singular values of rho.

Mode{T}(omega, shift)

Stuct which purpose is to model vibrational LHO mode in Molecule.

$V(q) = \hbar \omega (q - q_0)^2, \quad \hbar = 1$

Arguments

• omega: The frequency of LHO ($\omega$).
• shift: The shift of the coordinate of LHO ($q_0$).

Molecule{T,C1,C2}(modes, Nvib, E)

Stuct which purpose is to model a molecule in Aggregate.

Arguments

• modes: Vector of modes (Mode).
• Nvib: Maximum number of vibrational states for all modes.
• E: Energy of ground and excited state of molecule (HOMO, LUMO).

electronicIndices(molCount)

Get the electric indices for all states on Aggregate.

Arguments

• molCount: Number of molecules in Aggregate.

franckCondonFactors(size, shift)

Get Franck-Condon factors for LHO mode calculated using ShiftOperator.

getMolFrequencies(mol)

Get shifts of every mode on the Molecule state.

Arguments

• mol: Instance of Molecule.

getMolShifts(mol)

Get shifts of every mode on the Molecule state.

Arguments

• mol: Instance of Molecule.

getMolStateEnergy(mol, molElState, molVibState)

Get the energy of the Molecule state.

Arguments

• mol: Instance of Molecule.
• molElState: Electric state of the molecule in local basis (i.e. 1 or 2, where 1 is ground state).
• molVibState: Vibrational state of the molecule (e.g. [1, 5, 2, 2]).

getMolStateFC(mol, molElState1, molVibState1, molElState2, molVibState2)

Get the energy of the Molecule state.

Arguments

• mol: Instance of Molecule.
• molElState1, molElState2: Electric state of the molecule in local basis (i.e. 1 or 2, where 1 is ground state).
• molVibState1, molVibState2: Vibrational state of the molecule (e.g. [1, 5, 2, 2]).

vibrationalIndices(maxInds)

Get the vibrational indices for all states on Molecule.

Arguments

• maxInds: Vector of maximum number of vibrational states.

AnnihilationOperator{BL,BR}(basis_l, basis_r)

Dense annihilation operator as a mutable struct.

CreationOperator{BL,BR}(basis_l, basis_r)

Dense creation operator as a mutable struct.

MomentumOperator{BL,BR}(basisl, basisr)

Dense momentum operator as a mutable struct without mass and frequency.

PositionOperator{BL,BR}(basis_l, basis_r)

Dense position operator as a mutable struct without mass and frequency.

ShiftOperator{BL,BR}(basis_l, basis_r, shift)

Dense shift operator as a mutable struct using the definition

$D(\alpha) = \exp(\alpha a^\dagger - \alpha^* a)$.

Arguments

• basis_l: Bra basis.
• basis_r: Ket basis.
• shift: Shift or $\alpha$ parameter, can be complex number.

OneDenseOperator(basis_l, basis_r)

DenseOperator with ones on the diagonal.

schroedinger(psi0, tspan, H; 
	reltol=1.0e-12, abstol=1.0e-12, fout=nothing, alg=OrdinaryDiffEq.Tsit5())

Integrate Schroedinger equation to evolve states or compute propagators

$\frac{d}{d t}\vert\psi(t)\rangle = - \frac{i}{\hbar} \hat H\vert\psi(t)\rangle, \quad \hbar = 1$.

Arguments

• psi0: Initial state vector (can be a bra or a ket) or initial propagator.
• tspan: Vector specifying the points of time for which output should be displayed.
• H: Arbitrary operator specifying the Hamiltonian.
• reltol: Relative tolerance for OrdinaryDiffEq solver and its inner states.
• abstol: Absolute tolerance for OrdinaryDiffEq solver and its inner states.
• fout=nothing: If given, this function fout(t, psi) is called every time an output should be displayed. ATTENTION: The state psi is neither normalized nor permanent! It is still in use by the ode solver and therefore must not be changed.
• alg: Algorithm with which OrdinaryDiffEq will solve Schroedinger equation.

schroedinger_dynamic(psi0, tspan, f; 
	reltol=1.0e-12, abstol=1.0e-12, fout=nothing, alg=OrdinaryDiffEq.Tsit5())

Integrate time-dependent Schroedinger equation to evolve states or compute propagators

$\frac{d}{d t}\vert\psi(t)\rangle = - \frac{i}{\hbar} \hat H(t)\vert\rho(t)\rangle, \quad \hbar = 1$.

Arguments

• psi0: Initial state vector (can be a bra or a ket) or initial propagator.
• tspan: Vector specifying the points of time for which output should be displayed.
• f: Function f(t, psi) -> H returning the time and or state dependent Hamiltonian.
• reltol: Relative tolerance for OrdinaryDiffEq solver and its inner states.
• abstol: Absolute tolerance for OrdinaryDiffEq solver and its inner states.
• fout=nothing: If given, this function fout(t, psi) is called every time an output should be displayed. ATTENTION: The state psi is neither normalized nor permanent! It is still in use by the ode solver and therefore must not be changed.
• alg: Algorithm with which OrdinaryDiffEq will solve Schroedinger equation.

Commutator(A)

Create commutator in a form of superoperator from a given operator.

$\text{Commutator}(A) \:\cdot\: = [A, \:\cdot\:]$

integrate(tspan, df::Function, x0::Vector{ComplexF64},
        state::T, dstate::T, fout::Function; kwargs...)

Integrate using OrdinaryDiffEq

ad(rho, rho_bath, agg, FCProd, aggIndices, vibindices)

This is the inverse operation to the trace over bath trace_bath and get_rho_bath defined as follows

$\rho = \operatorname{ad}\{\rho_\text{tr}, \rho_\text{bath} \}$

correlation_function(t, rho0_bath, Ham_0, Ham_I, agg, FCProd, aggInds, vibindices)

Get time dependent correlation function for a specified time t using following definition

$C(t) = \operatorname{tr}_B \{ \hat{H}_I^{(I)}(t) \hat{H}_I \rho_\text{bath}(0) \}$

get_rho_bath(rho, agg, FCProd, aggIndices, vibindices; justCopy=false)

This method will return the bath part of rho knowing the result of trace_bath defined as follows

$\rho_\text{bath} = \operatorname{tr}_S \{\rho\}$

$\rho_{\text{bath}, ab} = \rho_{ab} / \langle a \vert \operatorname{tr}_B \{ \rho \}\vert b \rangle$

trace_bath(rho, a, b, agg, FCProd, aggIndices, vibindices)

Trace out bath degrees of freedom from rho without the product of Franck-Condon factors. The trace will be done only on the Hilber space for electric bra part a and ket part b. Input density matrix rho is for the whole Hilber space. This method returns number.

trace_bath_part(rho, a, b, agg, FCProd, aggIndices, vibindices)

Trace out bath degrees of freedom from rho without the product of Franck-Condon factors. The trace will be done only on the Hilber space for electric bra part a and ket part b. Input density matrix rho is only for the subspace. This method returns number.

trace_bath_slow(rho, agg, FCFact, aggIndices, vibindices)

Trace out bath degrees of freedom from rho without the product of Franck-Condon factors.

trace_bath(rho, agg, FCProd, aggIndices, vibindices)

Trace out bath degrees of freedom from rho

$\rho_\text{tr} = \operatorname{tr}_B \{\rho\} = \sum_{k} \langle k \vert \left( \sum_{ab} \rho_{am, bn} \vert am \rangle \langle bn \vert \right)\vert k \rangle$