API

API

Quantum-objects

Types

Abstract base class for all specialized bases.

The Basis class is meant to specify a basis of the Hilbert space of the studied system. Besides basis specific information all subclasses must implement a shape variable which indicates the dimension of the used Hilbert space. For a spin-1/2 Hilbert space this would be the vector Int[2]. A system composed of two spins would then have a shape vector Int[2 2].

Composite systems can be defined with help of the CompositeBasis class.

source
GenericBasis(N)

A general purpose basis of dimension N.

Should only be used rarely since it defeats the purpose of checking that the bases of state vectors and operators are correct for algebraic operations. The preferred way is to specify special bases for different systems.

source
CompositeBasis(b1, b2...)

Basis for composite Hilbert spaces.

Stores the subbases in a vector and creates the shape vector directly from the shape vectors of these subbases. Instead of creating a CompositeBasis directly tensor(b1, b2...) or b1 ⊗ b2 ⊗ … can be used.

source

Abstract base class for Bra and Ket states.

The state vector class stores the coefficients of an abstract state in respect to a certain basis. These coefficients are stored in the data field and the basis is defined in the basis field.

source
Bra(b::Basis[, data])

Bra state defined by coefficients in respect to the basis.

source
Ket(b::Basis[, data])

Ket state defined by coefficients in respect to the given basis.

source

Abstract base class for all operators.

All deriving operator classes have to define the fields basis_l and basis_r defining the left and right side bases.

For fast time evolution also at least the function gemv!(alpha, op::AbstractOperator, x::Ket, beta, result::Ket) should be implemented. Many other generic multiplication functions can be defined in terms of this function and are provided automatically.

source

Abstract type for operators with a data field.

This is an abstract type for operators that have a direct matrix representation stored in their .data field.

source
DenseOperator(b1[, b2, data])

Dense array implementation of Operator.

The matrix consisting of complex floats is stored in the data field.

source
SparseOperator(b1[, b2, data])

Sparse array implementation of Operator.

The matrix is stored as the julia built-in type SparseMatrixCSC in the data field.

source
LazyTensor(b1[, b2], indices, operators[, factor=1])

Lazy implementation of a tensor product of operators.

The suboperators are stored in the operators field. The indices field specifies in which subsystem the corresponding operator lives. Additionally, a complex factor is stored in the factor field which allows for fast multiplication with numbers.

source
LazySum([factors,] operators)

Lazy evaluation of sums of operators.

All operators have to be given in respect to the same bases. The field factors accounts for an additional multiplicative factor for each operator stored in the operators field.

source
LazyProduct(operators[, factor=1])
LazyProduct(op1, op2...)

Lazy evaluation of products of operators.

The factors of the product are stored in the operators field. Additionally a complex factor is stored in the factor field which allows for fast multiplication with numbers.

source

Base class for all super operator classes.

Super operators are bijective mappings from operators given in one specific basis to operators, possibly given in respect to another, different basis. To embed super operators in an algebraic framework they are defined with a left hand basis basis_l and a right hand basis basis_r where each of them again consists of a left and right hand basis.

\[A_{bl_1,bl_2} = S_{(bl_1,bl_2) ↔ (br_1,br_2)} B_{br_1,br_2} \\ A_{br_1,br_2} = B_{bl_1,bl_2} S_{(bl_1,bl_2) ↔ (br_1,br_2)}\]
source
DenseSuperOperator(b1[, b2, data])

SuperOperator stored as dense matrix.

source
SparseSuperOperator(b1[, b2, data])

SuperOperator stored as sparse matrix.

source

Functions

basisstate(b, index)

Basis vector specified by index as ket state.

For a composite system index can be a vector which then creates a tensor product state $|i_1⟩⊗|i_2⟩⊗…⊗|i_n⟩$ of the corresponding basis states.

source
basisstate(b::ManyBodyBasis, occupation::Vector{Int})

Return a ket state where the system is in the state specified by the given occupation numbers.

source
identityoperator(a::Basis[, b::Basis])

Return an identityoperator in the given bases.

source
diagonaloperator(b::Basis)

Create a diagonal operator of type SparseOperator.

source
randoperator(b1[, b2])

Calculate a random unnormalized dense operator.

source
spre(op)

Create a super-operator equivalent for right side operator multiplication.

For operators $A$, $B$ the relation

\[ \mathrm{spre}(A) B = A B\]

holds. op can be a dense or a sparse operator.

source

Create a super-operator equivalent for left side operator multiplication.

For operators $A$, $B$ the relation

\[ \mathrm{spost}(A) B = B A\]

holds. op can be a dense or a sparse operator.

source
liouvillian(H, J; rates, Jdagger)

Create a super-operator equivalent to the master equation so that $\dot ρ = S ρ$.

The super-operator $S$ is defined by

\[S ρ = -\frac{i}{ħ} [H, ρ] + \sum_i J_i ρ J_i^† - \frac{1}{2} J_i^† J_i ρ - \frac{1}{2} ρ J_i^† J_i\]

Arguments

  • H: Hamiltonian.
  • J: Vector containing the jump operators.
  • rates: Vector or matrix specifying the coefficients for the jump operators.
  • Jdagger: Vector containing the hermitian conjugates of the jump operators. If they are not given they are calculated automatically.
source
samebases(a, b)

Test if two objects have the same bases.

source
check_samebases(a, b)

Throw an IncompatibleBases error if the objects don't have the same bases.

source
multiplicable(a, b)

Check if two objects are multiplicable.

source
check_multiplicable(a, b)

Throw an IncompatibleBases error if the objects are not multiplicable.

source
basis(a)

Return the basis of an object.

If it's ambiguous, e.g. if an operator has a different left and right basis, an IncompatibleBases error is thrown.

source
dagger(x)

Hermitian conjugate.

source
tensor(x::Ket, y::Bra)

Outer product $|x⟩⟨y|$ of the given states.

source
tensor(x::Ket, y::Ket, z::Ket...)

Tensor product $|x⟩⊗|y⟩⊗|z⟩⊗…$ of the given states.

source
tensor(x, y, z...)

Tensor product of the given objects. Alternatively, the unicode symbol ⊗ (\otimes) can be used.

source
tensor(x::Basis, y::Basis, z::Basis...)

Create a CompositeBasis from the given bases.

Any given CompositeBasis is expanded so that the resulting CompositeBasis never contains another CompositeBasis.

source
tensor(x::AbstractOperator, y::AbstractOperator, z::AbstractOperator...)

Tensor product $\hat{x}⊗\hat{y}⊗\hat{z}⊗…$ of the given operators.

source
projector(a::Ket, b::Bra)

Projection operator $|a⟩⟨b|$.

source
projector(a::Ket)

Projection operator $|a⟩⟨a|$.

source
projector(a::Bra)

Projection operator $|a⟩⟨a|$.

source
dm(a::StateVector)

Create density matrix $|a⟩⟨a|$. Same as projector(a).

source
LinearAlgebra.normMethod.
norm(x::StateVector)

Norm of the given bra or ket state.

source
LinearAlgebra.trFunction.
tr(x::AbstractOperator)

Trace of the given operator.

source
ptrace(a, indices)

Partial trace of the given basis, state or operator.

The indices argument, which can be a single integer or a vector of integers, specifies which subsystems are traced out. The number of indices has to be smaller than the number of subsystems, i.e. it is not allowed to perform a full trace.

source
normalize(x::StateVector)

Return the normalized state so that norm(x) is one.

source
normalize(op)

Return the normalized operator so that its tr(op) is one.

source
normalize!(x::StateVector)

In-place normalization of the given bra or ket so that norm(x) is one.

source
normalize!(op)

In-place normalization of the given operator so that its tr(x) is one.

source
expect(op, state)

Expectation value of the given operator op for the specified state.

state can either be a (density) operator or a ket.

source
expect(index, op, state)

If an index is given, it assumes that op is defined in the subsystem specified by this number.

source
variance(op, state)

Variance of the given operator op for the specified state.

state can either be a (density) operator or a ket.

source
variance(index, op, state)

If an index is given, it assumes that op is defined in the subsystem specified by this number

source
embed(basis1[, basis2], indices::Vector, operators::Vector)

Tensor product of operators where missing indices are filled up with identity operators.

source
embed(basis1[, basis2], operators::Dict)

operators is a dictionary Dict{Vector{Int}, AbstractOperator}. The integer vector specifies in which subsystems the corresponding operator is defined.

source
permutesystems(a, perm)

Change the ordering of the subsystems of the given object.

For a permutation vector [2,1,3] and a given object with basis [b1, b2, b3] this function results in [b2, b1, b3].

source
Base.expMethod.
exp(op::AbstractOperator)

Operator exponential.

source
Base.expMethod.
exp(op::DenseSuperOperator)

Operator exponential which can for example used to calculate time evolutions.

source
gemv!(alpha, a, b, beta, result)

Fast in-place multiplication of operators with state vectors. It implements the relation result = beta*result + alpha*a*b. Here, alpha and beta are complex numbers, while result and either a or b are state vectors while the other one can be of any operator type.

source
gemm!(alpha, a, b, beta, result)

Fast in-place multiplication of of operators with DenseOperators. It implements the relation result = beta*result + alpha*a*b. Here, alpha and beta are complex numbers, while result and either a or b are dense operators while the other one can be of any operator type.

source
dense(op::AbstractOperator)

Convert an arbitrary Operator into a DenseOperator.

source
sparse(op::AbstractOperator)

Convert an arbitrary operator into a SparseOperator.

source

Exceptions

Exception that should be raised for an illegal algebraic operation.

source

Quantum systems

Fock

FockBasis(N)

Basis for a Fock space where N specifies a cutoff, i.e. what the highest included fock state is. Note that the dimension of this basis then is N+1.

source
number(b::FockBasis)

Number operator for the specified Fock space.

source
destroy(b::FockBasis)

Annihilation operator for the specified Fock space.

source
create(b::FockBasis)

Creation operator for the specified Fock space.

source
displace(b::FockBasis, alpha)

Displacement operator $D(α)$ for the specified Fock space.

source
fockstate(b::FockBasis, n)

Fock state $|n⟩$ for the specified Fock space.

source
coherentstate(b::FockBasis, alpha)

Coherent state $|α⟩$ for the specified Fock space.

source

Phase space

qfunc(a, α)
qfunc(a, x, y)
qfunc(a, xvec, yvec)

Husimi Q representation $⟨α|ρ|α⟩/π$ for the given state or operator a. The function can either be evaluated on one point α or on a grid specified by the vectors xvec and yvec. Note that conversion from x and y to α is done via the relation $α = \frac{1}{\sqrt{2}}(x + i y)$.

source
wigner(a, α)
wigner(a, x, y)
wigner(a, xvec, yvec)

Wigner function for the given state or operator a. The function can either be evaluated on one point α or on a grid specified by the vectors xvec and yvec. Note that conversion from x and y to α is done via the relation $α = \frac{1}{\sqrt{2}}(x + i y)$.

source
coherentspinstate(b::SpinBasis, θ::Real, ϕ::Real)

A coherent spin state |θ,ϕ⟩ is analogous to the coherent state of the linear harmonic oscillator. Coherent spin states represent a collection of identical two-level systems and can be described by two angles θ and ϕ (although this parametrization is not unique), similarly to a qubit on the Bloch sphere.

source
qfuncsu2(ket,Ntheta;Nphi=2Ntheta)
qfuncsu2(rho,Ntheta;Nphi=2Ntheta)

Husimi Q SU(2) representation $⟨θ,ϕ|ρ|θ,ϕ⟩/π$ for the given state.

The function calculates the SU(2) Husimi representation of a state on the generalised bloch sphere (0 < θ < π and 0 < ϕ < 2 π) with a given resolution (Ntheta, Nphi).

qfuncsu2(rho,θ,ϕ)
qfuncsu2(ket,θ,ϕ)

This version calculates the Husimi Q SU(2) function at a position given by θ and ϕ.

source
wignersu2(ket,Ntheta;Nphi=2Ntheta)
wignersu2(rho,Ntheta;Nphi=2Ntheta)

Wigner SU(2) representation for the given state with a resolution (Ntheta, Nphi).

The function calculates the SU(2) Wigner representation of a state on the generalised bloch sphere (0 < θ < π and 0 < ϕ < 2 π) with a given resolution by decomposing the state into the basis of spherical harmonics.

wignersu2(rho,θ,ϕ)
wignersu2(ket,θ,ϕ)

This version calculates the Wigner SU(2) function at a position given by θ and ϕ

source
ylm(l::Integer,m::Integer,theta::Real,phi::Real)

Spherical harmonics Y(l,m)(θ,ϕ) where l ∈ N, m = -l,-l+1,...,l-1,l, θ ∈ [0,π], and ϕ ∈ [0,2π).

This function calculates the value of Y(l,m) spherical harmonic at position θ and ϕ.

source

N-level

NLevelBasis(N)

Basis for a system consisting of N states.

source
transition(b::NLevelBasis, to::Int, from::Int)

Transition operator $|\mathrm{to}⟩⟨\mathrm{from}|$.

source
nlevelstate(b::NLevelBasis, n::Int)

State where the system is completely in the n-th level.

source

Spin

SpinBasis(n)

Basis for spin-n particles.

The basis can be created for arbitrary spinnumbers by using a rational number, e.g. SpinBasis(3//2). The Pauli operators are defined for all possible spin numbers.

source
sigmax(b::SpinBasis)

Pauli $σ_x$ operator for the given Spin basis.

source
sigmay(b::SpinBasis)

Pauli $σ_y$ operator for the given Spin basis.

source
sigmaz(b::SpinBasis)

Pauli $σ_z$ operator for the given Spin basis.

source
sigmap(b::SpinBasis)

Raising operator $σ_+$ for the given Spin basis.

source
sigmam(b::SpinBasis)

Lowering operator $σ_-$ for the given Spin basis.

source
spinup(b::SpinBasis)

Spin up state for the given Spin basis.

source
spindown(b::SpinBasis)

Spin down state for the given Spin basis.

source

Particle

PositionBasis(xmin, xmax, Npoints)
PositionBasis(b::MomentumBasis)

Basis for a particle in real space.

For simplicity periodic boundaries are assumed which means that the rightmost point defined by xmax is not included in the basis but is defined to be the same as xmin.

When a MomentumBasis is given as argument the exact values of $x_{min}$ and $x_{max}$ are due to the periodic boundary conditions more or less arbitrary and are chosen to be $-\pi/dp$ and $\pi/dp$ with $dp=(p_{max}-p_{min})/N$.

source
MomentumBasis(pmin, pmax, Npoints)
MomentumBasis(b::PositionBasis)

Basis for a particle in momentum space.

For simplicity periodic boundaries are assumed which means that pmax is not included in the basis but is defined to be the same as pmin.

When a PositionBasis is given as argument the exact values of $p_{min}$ and $p_{max}$ are due to the periodic boundary conditions more or less arbitrary and are chosen to be $-\pi/dx$ and $\pi/dx$ with $dx=(x_{max}-x_{min})/N$.

source
spacing(b::PositionBasis)

Difference between two adjacent points of the real space basis.

source
spacing(b::MomentumBasis)

Momentum difference between two adjacent points of the momentum basis.

source
samplepoints(b::PositionBasis)

x values of the real space basis.

source
samplepoints(b::MomentumBasis)

p values of the momentum basis.

source
Base.positionMethod.
position(b::PositionBasis)

Position operator in real space.

source
Base.positionMethod.
position(b:MomentumBasis)

Position operator in momentum space.

source
momentum(b::PositionBasis)

Momentum operator in real space.

source
momentum(b:MomentumBasis)

Momentum operator in momentum space.

source
potentialoperator(b::PositionBasis, V(x))

Operator representing a potential $V(x)$ in real space.

source
potentialoperator(b::MomentumBasis, V(x))

Operator representing a potential $V(x)$ in momentum space.

source
potentialoperator(b::CompositeBasis, V(x, y, z, ...))

Operator representing a potential $V$ in more than one dimension.

Arguments

  • b: Composite basis consisting purely either of PositionBasis or MomentumBasis. Note, that calling this with a composite basis in momentum space might consume a large amount of memory.
  • V: Function describing the potential. ATTENTION: The number of arguments accepted by V must match the spatial dimension. Furthermore, the order of the arguments has to match that of the order of the tensor product of bases (e.g. if b=bx⊗by⊗bz, then V(x,y,z)).
source
gaussianstate(b::PositionBasis, x0, p0, sigma)
gaussianstate(b::MomentumBasis, x0, p0, sigma)

Create a Gaussian state around x0 andp0 with width sigma.

In real space the gaussian state is defined as

\[\Psi(x) = \frac{\sqrt{\Delta x}}{\pi^{1/4}\sqrt{\sigma}} e^{i p_0 (x-\frac{x_0}{2}) - \frac{(x-x_0)^2}{2 \sigma^2}}\]

and is connected to the momentum space definition

\[\Psi(p) = \frac{\sqrt{\sigma} \sqrt{\Delta x}}{\pi^{1/4}} e^{-i x_0 (p-\frac{p_0}{2}) - \frac{1}{2}(p-p_0)^2 \sigma^2}\]

via a Fourier-transformation

\[\Psi(p) = \frac{1}{\sqrt{2\pi}} \int_{-\infty}^{\infty} e^{-ipx}\Psi(x) \mathrm{d}x\]

The state has the properties

  • $⟨p⟩ = p_0$
  • $⟨x⟩ = x_0$
  • $\mathrm{Var}(x) = \frac{σ^2}{2}$
  • $\mathrm{Var}(p) = \frac{1}{2 σ^2}$

Due to the numerically necessary discretization additional scaling factora $\sqrt{Δx}$ and $\sqrt{Δp}$ are used so that $Ψx_i = \sqrt{Δ x} Ψ(x_i)$ and $Ψp_i = \sqrt{Δ x} Ψ(p_i)$ so that the resulting Ket state is normalized.

source
FFTOperator

Abstract type for all implementations of FFT operators.

source
FFTOperators

Operator performing a fast fourier transformation when multiplied with a state that is a Ket or an Operator.

source
FFTKets

Operator that can only perform fast fourier transformations on Kets. This is much more memory efficient when only working with Kets.

source
transform(b1::PositionBasis, b2::FockBasis; x0=1)
transform(b1::FockBasis, b2::PositionBasis; x0=1)

Transformation operator between position basis and fock basis.

The coefficients are connected via the relation

\[ψ(x_i) = \sum_{n=0}^N ⟨x_i|n⟩ ψ_n\]

where $⟨x_i|n⟩$ is the value of the n-th eigenstate of a particle in a harmonic trap potential at position $x$, i.e.:

\[⟨x_i|n⟩ = π^{-\frac{1}{4}} \frac{e^{-\frac{1}{2}\left(\frac{x}{x_0}\right)^2}}{\sqrt{x_0}} \frac{1}{\sqrt{2^n n!}} H_n\left(\frac{x}{x_0}\right)\]
source
transform(b1::MomentumBasis, b2::PositionBasis)
transform(b1::PositionBasis, b2::MomentumBasis)

Transformation operator between position basis and momentum basis.

source
transform(b1::CompositeBasis, b2::CompositeBasis)

Transformation operator between two composite bases. Each of the bases has to contain bases of type PositionBasis and the other one a corresponding MomentumBasis.

source

Subspace bases

SubspaceBasis(basisstates)

A basis describing a subspace embedded a higher dimensional Hilbert space.

source
orthonormalize(b::SubspaceBasis)

Orthonormalize the basis states of the given SubspaceBasis

A modified Gram-Schmidt process is used.

source
projector(b1, b2)

Projection operator between subspaces and superspaces or between two subspaces.

source

Many-body

ManyBodyBasis(b, occupations)

Basis for a many body system.

The basis has to know the associated one-body basis b and which occupation states should be included. The occupations_hash is used to speed up checking if two many-body bases are equal.

source
fermionstates(Nmodes, Nparticles)
fermionstates(b, Nparticles)

Generate all fermionic occupation states for N-particles in M-modes. Nparticles can be a vector to define a Hilbert space with variable particle number.

source
bosonstates(Nmodes, Nparticles)
bosonstates(b, Nparticles)

Generate all bosonic occupation states for N-particles in M-modes. Nparticles can be a vector to define a Hilbert space with variable particle number.

source
number(b::ManyBodyBasis, index)

Particle number operator for the i-th mode of the many-body basis b.

source
number(b::ManyBodyBasis)

Total particle number operator.

source
destroy(b::ManyBodyBasis, index)

Annihilation operator for the i-th mode of the many-body basis b.

source
create(b::ManyBodyBasis, index)

Creation operator for the i-th mode of the many-body basis b.

source
manybodyoperator(b::ManyBodyBasis, op)

Create the many-body operator from the given one-body operator op.

The given operator can either be a one-body operator or a two-body interaction. Higher order interactions are at the moment not implemented.

The mathematical formalism for the one-body case is described by

\[X = \sum_{ij} a_i^† a_j ⟨u_i| x | u_j⟩\]

and for the interaction case by

\[X = \sum_{ijkl} a_i^† a_j^† a_k a_l ⟨u_i|⟨u_j| x |u_k⟩|u_l⟩\]

where $X$ is the N-particle operator, $x$ is the one-body operator and $|u⟩$ are the one-body states associated to the different modes of the N-particle basis.

source
onebodyexpect(op, state)

Expectation value of the one-body operator op in respect to the many-body state.

source

Metrics

tracenorm(rho)

Trace norm of rho.

It is defined as

\[T(ρ) = Tr\{\sqrt{ρ^† ρ}\}.\]

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

source
tracenorm_h(rho)

Trace norm of rho.

It uses the identity

\[T(ρ) = Tr\{\sqrt{ρ^† ρ}\} = \sum_i |λ_i|\]

where $λ_i$ are the eigenvalues of rho.

source
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

\[ T(ρ) = Tr\{\sqrt{ρ^† ρ}\} = \sum_i σ_i\]

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

source
tracedistance(rho, sigma)

Trace distance between rho and sigma.

It is defined as

\[T(ρ,σ) = \frac{1}{2} Tr\{\sqrt{(ρ - σ)^† (ρ - σ)}\}.\]

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

source
tracedistance_h(rho, sigma)

Trace distance between rho and sigma.

It uses the identity

\[T(ρ,σ) = \frac{1}{2} Tr\{\sqrt{(ρ - σ)^† (ρ - σ)}\} = \frac{1}{2} \sum_i |λ_i|\]

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

source
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

\[ T(ρ,σ) = \frac{1}{2} Tr\{\sqrt{(ρ - σ)^† (ρ - σ)}\} = \frac{1}{2} \sum_i σ_i\]

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

source
entropy_vn(rho)

Von Neumann entropy of a density matrix.

The Von Neumann entropy of a density operator is defined as

\[S(ρ) = -Tr(ρ \log(ρ)) = -\sum_n λ_n\log(λ_n)\]

where $λ_n$ are the eigenvalues of the density matrix $ρ$, $\log$ is the natural logarithm and $0\log(0) ≡ 0$.

Arguments

  • rho: Density operator of which to calculate Von Neumann entropy.
  • tol=1e-15: Tolerance for rounding errors in the computed eigenvalues.
source
fidelity(rho, sigma)

Fidelity of two density operators.

The fidelity of two density operators $ρ$ and $σ$ is defined by

\[F(ρ, σ) = Tr\left(\sqrt{\sqrt{ρ}σ\sqrt{ρ}}\right),\]

where $\sqrt{ρ}=\sum_n\sqrt{λ_n}|ψ⟩⟨ψ|$.

source
ptranspose(rho, index)

Partial transpose of rho with respect to subspace specified by index.

source
PPT(rho, index)

Peres-Horodecki criterion of partial transpose.

source
negativity(rho, index)

Negativity of rho with respect to subsystem index.

The negativity of a density matrix ρ is defined as

\[N(ρ) = \|ρᵀ\|,\]

where ρᵀ is the partial transpose.

source
logarithmic_negativity(rho, index)

The logarithmic negativity of a density matrix ρ is defined as

\[N(ρ) = \log₂\|ρᵀ\|,\]

where ρᵀ is the partial transpose.

source

Time-evolution

Schroedinger

timeevolution.schroedinger(tspan, psi0, H; fout)

Integrate Schroedinger equation.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • psi0: Initial state vector (can be a bra or a ket).
  • H: Arbitrary operator specifying the Hamiltonian.
  • 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.
source
timeevolution.schroedinger_dynamic(tspan, psi0, f; fout)

Integrate time-dependent Schroedinger equation.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • psi0: Initial state vector (can be a bra or a ket).
  • f: Function f(t, psi) -> H returning the time and or state dependent Hamiltonian.
  • 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.
source

Master

timeevolution.master(tspan, rho0, H, J; <keyword arguments>)

Time-evolution according to a master equation.

There are two implementations for integrating the master equation:

  • master_h: Usual formulation of the master equation.
  • master_nh: Variant with non-hermitian Hamiltonian.

For dense arguments the master function calculates the non-hermitian Hamiltonian and then calls master_nh which is slightly faster.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • rho0: Initial density operator. Can also be a state vector which is automatically converted into a density operator.
  • H: Arbitrary operator specifying the Hamiltonian.
  • J: Vector containing all jump operators which can be of any arbitrary operator type.
  • rates=nothing: Vector or matrix specifying the coefficients (decay rates) for the jump operators. If nothing is specified all rates are assumed to be 1.
  • Jdagger=dagger.(J): Vector containing the hermitian conjugates of the jump operators. If they are not given they are calculated automatically.
  • fout=nothing: If given, this function fout(t, rho) is called every time an output should be displayed. ATTENTION: The given state rho is not permanent! It is still in use by the ode solver and therefore must not be changed.
  • kwargs...: Further arguments are passed on to the ode solver.
source
timeevolution.master_h(tspan, rho0, H, J; <keyword arguments>)

Integrate the master equation with dmaster_h as derivative function.

Further information can be found at master.

source
timeevolution.master_nh(tspan, rho0, H, J; <keyword arguments>)

Integrate the master equation with dmaster_nh as derivative function.

In this case the given Hamiltonian is assumed to be the non-hermitian version:

\[H_{nh} = H - \frac{i}{2} \sum_k J^†_k J_k\]

Further information can be found at master.

source
timeevolution.master_dynamic(tspan, rho0, f; <keyword arguments>)

Time-evolution according to a master equation with a dynamic Hamiltonian and J.

There are two implementations for integrating the master equation with dynamic operators:

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • rho0: Initial density operator. Can also be a state vector which is automatically converted into a density operator.
  • f: Function f(t, rho) -> (H, J, Jdagger) or f(t, rho) -> (H, J, Jdagger, rates)
  • rates=nothing: Vector or matrix specifying the coefficients (decay rates) for the jump operators. If nothing is specified all rates are assumed to be 1.
  • fout=nothing: If given, this function fout(t, rho) is called every time an output should be displayed. ATTENTION: The given state rho is not permanent! It is still in use by the ode solver and therefore must not be changed.
  • kwargs...: Further arguments are passed on to the ode solver.
source
timeevolution.master_dynamic(tspan, rho0, f; <keyword arguments>)

Time-evolution according to a master equation with a dynamic non-hermitian Hamiltonian and J.

In this case the given Hamiltonian is assumed to be the non-hermitian version.

\[H_{nh} = H - \frac{i}{2} \sum_k J^†_k J_k\]

The given function can either be of the form f(t, rho) -> (Hnh, Hnhdagger, J, Jdagger) or f(t, rho) -> (Hnh, Hnhdagger, J, Jdagger, rates) For further information look at master_dynamic.

source
@skiptimechecks

Macro to skip checks during time-dependent problems. Useful for timeevolution.master_dynamic and similar functions.

source

Monte Carlo wave function

mcwf(tspan, psi0, H, J; <keyword arguments>)

Integrate the master equation using the MCWF method.

There are two implementations for integrating the non-hermitian schroedinger equation:

  • mcwf_h: Usual formulation with Hamiltonian + jump operators

separately.

  • mcwf_nh: Variant with non-hermitian Hamiltonian.

The mcwf function takes a normal Hamiltonian, calculates the non-hermitian Hamiltonian and then calls mcwf_nh which is slightly faster.

Arguments

  • tspan: Vector specifying the points of time for which output should

be displayed.

  • psi0: Initial state vector.
  • H: Arbitrary Operator specifying the Hamiltonian.
  • J: Vector containing all jump operators which can be of any arbitrary

operator type.

  • seed=rand(): Seed used for the random number generator.
  • rates=ones(): Vector of decay rates.
  • fout: 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 solve and therefore must not be changed.

  • Jdagger=dagger.(J): Vector containing the hermitian conjugates of the jump

operators. If they are not given they are calculated automatically.

  • display_beforeevent=false: fout is called before every jump.
  • display_afterevent=false: fout is called after every jump.
  • kwargs...: Further arguments are passed on to the ode solver.
source
mcwf_h(tspan, rho0, Hnh, J; <keyword arguments>)

Calculate MCWF trajectory where the Hamiltonian is given in hermitian form.

For more information see: mcwf

source
mcwf_nh(tspan, rho0, Hnh, J; <keyword arguments>)

Calculate MCWF trajectory where the Hamiltonian is given in non-hermitian form.

\[H_{nh} = H - \frac{i}{2} \sum_k J^†_k J_k\]

For more information see: mcwf

source
mcwf_dynamic(tspan, psi0, f; <keyword arguments>)

Integrate the master equation using the MCWF method with dynamic Hamiltonian and Jump operators.

The mcwf function takes a normal Hamiltonian, calculates the non-hermitian Hamiltonian and then calls mcwf_nh which is slightly faster.

Arguments

  • tspan: Vector specifying the points of time for which output should

be displayed.

  • psi0: Initial state vector.
  • f: Function f(t, psi) -> (H, J, Jdagger) or f(t, psi) -> (H, J, Jdagger, rates) that returns the time-dependent Hamiltonian and Jump operators.
  • seed=rand(): Seed used for the random number generator.
  • rates=ones(): Vector of decay rates.
  • fout: 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 solve and therefore must not be changed.

  • display_beforeevent=false: fout is called before every jump.
  • display_afterevent=false: fout is called after every jump.
  • kwargs...: Further arguments are passed on to the ode solver.
source
mcwf_nh_dynamic(tspan, rho0, f; <keyword arguments>)

Calculate MCWF trajectory where the dynamic Hamiltonian is given in non-hermitian form.

For more information see: mcwf_dynamic

source
diagonaljumps(rates, J)

Diagonalize jump operators.

The given matrix rates of decay rates is diagonalized and the corresponding set of jump operators is calculated.

Arguments

  • rates: Matrix of decay rates.
  • J: Vector of jump operators.
source

Spectral analysis

eigenstates(op::AbstractOperator[, n::Int; warning=true])

Calculate the lowest n eigenvalues and their corresponding eigenstates.

This is just a thin wrapper around julia's eigen and eigs functions. Which of them is used depends on the type of the given operator. If more control about the way the calculation is done is needed, use the functions directly. More details can be found at [http://docs.julialang.org/en/stable/stdlib/linalg/].

NOTE: Especially for small systems full diagonalization with Julia's eigen function is often more desirable. You can convert a sparse operator A to a dense one using dense(A).

If the given operator is non-hermitian a warning is given. This behavior can be turned off using the keyword warning=false.

source

For sparse operators by default it only returns the 6 lowest eigenvalues.

source
eigenenergies(op::AbstractOperator[, n::Int; warning=true])

Calculate the lowest n eigenvalues.

This is just a thin wrapper around julia's eigvals. If more control about the way the calculation is done is needed, use the function directly. More details can be found at [http://docs.julialang.org/en/stable/stdlib/linalg/].

If the given operator is non-hermitian a warning is given. This behavior can be turned off using the keyword warning=false.

source

For sparse operators by default it only returns the 6 lowest eigenvalues.

source
simdiag(ops; atol, rtol)

Simultaneously diagonalize commuting Hermitian operators specified in ops.

This is done by diagonalizing the sum of the operators. The eigenvalues are computed by $a = ⟨ψ|A|ψ⟩$ and it is checked whether the eigenvectors fulfill the equation $A|ψ⟩ = a|ψ⟩$.

Arguments

  • ops: Vector of sparse or dense operators.
  • atol=1e-14: kwarg of Base.isapprox specifying the tolerance of the approximate check
  • rtol=1e-14: kwarg of Base.isapprox specifying the tolerance of the approximate check

Returns

  • evals_sorted: Vector containing all vectors of the eigenvalues sorted by the eigenvalues of the first operator.
  • v: Common eigenvectors.
source

Steady-states

steadystate.master(H, J; <keyword arguments>)

Calculate steady state using long time master equation evolution.

Arguments

  • H: Arbitrary operator specifying the Hamiltonian.
  • J: Vector containing all jump operators which can be of any arbitrary operator type.
  • rho0=dm(basisstate(b)): Initial density operator. If not given the $|0⟩⟨0|$ state in respect to the choosen basis is used.
  • tol=1e-3: Tracedistance used as termination criterion.
  • hmin=1e-7: Minimal time step used in the time evolution.
  • rates=ones(N): Vector or matrix specifying the coefficients for the jump operators.
  • Jdagger=dagger.(Jdagger): Vector containing the hermitian conjugates of the jump operators. If they are not given they are calculated automatically.
  • fout=nothing: If given this function fout(t, rho) is called every time an output should be displayed. To limit copying to a minimum the given density operator rho is further used and therefore must not be changed.
  • kwargs...: Further arguments are passed on to the ode solver.
source
steadystate.eigenvector(L)
steadystate.eigenvector(H, J)

Find steady state by calculating the eigenstate with eigenvalue 0 of the Liouvillian matrix L, if it exists.

Keyword arguments:

  • tol = 1e-9: Check abs(eigenvalue) < tol to determine zero eigenvalue.
  • nev = 2: Number of calculated eigenvalues. If nev > 1 it is checked if there is only one eigenvalue with real part 0. No checks for nev = 1: use if faster or for avoiding convergence errors of eigs. Changing nev thus only makes sense when using SparseSuperOperator.
  • which = :LR: Find eigenvalues with largest real part. Keyword for eigs function (ineffective for DenseSuperOperator).
  • kwargs...: Keyword arguments for the Julia eigen or eigs function.
source
steadystate.liouvillianspectrum(L)
steadystate.liouvillianspectrum(H, J)

Calculate eigenspectrum of the Liouvillian matrix L. The eigenvalues and -states are sorted according to the absolute value of the eigenvalues.

Keyword arguments:

  • nev = min(10, length(L.basis_r[1])*length(L.basis_r[2])): Number of eigenvalues.
  • which = :LR: Find eigenvalues with largest real part. Keyword for eigs function (ineffective for DenseSuperOperator).
  • kwargs...: Keyword arguments for the Julia eigen or eigens function.
source

Time correlations

timecorrelations.correlation([tspan, ]rho0, H, J, op1, op2; <keyword arguments>)

Calculate two time correlation values $⟨A(t)B(0)⟩$.

The calculation is done by multiplying the initial density operator with $B$ performing a time evolution according to a master equation and then calculating the expectation value $\mathrm{Tr} \{A ρ\}$

Without the tspan argument the points in time are chosen automatically from the ode solver and the final time is determined by the steady state termination criterion specified in steadystate.master.

Arguments

  • tspan: Points of time at which the correlation should be calculated.
  • rho0: Initial density operator.
  • H: Operator specifying the Hamiltonian.
  • J: Vector of jump operators.
  • op1: Operator at time t.
  • op2: Operator at time t=0.
  • rates=ones(N): Vector or matrix specifying the coefficients (decay rates) for the jump operators.
  • Jdagger=dagger.(J): Vector containing the hermitian conjugates of the jump
  • kwargs...: Further arguments are passed on to the ode solver.
source
timecorrelations.spectrum([omega_samplepoints,] H, J, op; <keyword arguments>)

Calculate spectrum as Fourier transform of a correlation function

This is done with the Wiener-Khinchin theorem

\[S(ω, t) = 2\Re\left\{\int_0^{∞} dτ e^{-iωτ}⟨A^†(t+τ)A(t)⟩\right\}\]

The argument omega_samplepoints gives the list of frequencies where $S(ω)$ is caclulated. A corresponding list of times is calculated internally by means of a inverse discrete frequency fourier transform. If not given, the steady-state is computed before calculating the auto-correlation function.

Without the omega_samplepoints arguments the frequencies are chosen automatically.

Arguments

  • omega_samplepoints: List of frequency points at which the spectrum is calculated.
  • H: Operator specifying the Hamiltonian.
  • J: Vector of jump operators.
  • op: Operator for which the auto-correlation function is calculated.
  • rho0: Initial density operator.
  • tol=1e-4: Tracedistance used as termination criterion.
  • rates=ones(N): Vector or matrix specifying the coefficients for the jump operators.
  • Jdagger=dagger.(J): Vector containing the hermitian conjugates of the jump operators. If they are not given they are calculated automatically.
  • kwargs...: Further arguments are passed on to the ode solver.
source
timecorrelations.correlation2spectrum(tspan, corr; normalize_spec)

Calculate spectrum as Fourier transform of a correlation function with a given correlation function.

Arguments

  • tspan: List of time points corresponding to the correlation function.
  • corr: Correlation function of which the Fourier transform is to be calculated.
  • normalize_spec: Specify if spectrum should be normalized to its maximum.
source

Semi-classical

Semi-classical state.

It consists of a quantum part, which is either a Ket or a DenseOperator and a classical part that is specified as a complex vector of arbitrary length.

source
semiclassical.schroedinger_dynamic(tspan, state0, fquantum, fclassical[; fout, ...])

Integrate time-dependent Schrödinger equation coupled to a classical system.

Arguments

  • tspan: Vector specifying the points of time for which the output should be displayed.
  • psi0: Initial semi-classical state semiclassical.State.
  • fquantum: Function f(t, psi, u) -> H returning the time and or state dependent Hamiltonian.
  • fclassical: Function f(t, psi, u, du) calculating the possibly time and state dependent derivative of the classical equations and storing it in the vector du.
  • fout=nothing: If given, this function fout(t, state) is called every time an output should be displayed. ATTENTION: The given state is neither normalized nor permanent!
  • kwargs...: Further arguments are passed on to the ode solver.
source
semiclassical.master_dynamic(tspan, state0, fquantum, fclassical; <keyword arguments>)

Integrate time-dependent master equation coupled to a classical system.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • rho0: Initial semi-classical state semiclassical.State.
  • fquantum: Function f(t, rho, u) -> (H, J, Jdagger) returning the time and/or state dependent Hamiltonian and Jump operators.
  • fclassical: Function f(t, rho, u, du) calculating the possibly time and state dependent derivative of the classical equations and storing it in the complex vector du.
  • fout=nothing: If given, this function fout(t, state) is called every time an output should be displayed. ATTENTION: The given state is not permanent!
  • kwargs...: Further arguments are passed on to the ode solver.
source

Stochastics

stochastic.schroedinger(tspan, state0, H, Hs[; fout, ...])

Integrate stochastic Schrödinger equation.

Arguments

  • tspan: Vector specifying the points of time for which the output should be displayed.
  • psi0: Initial state as Ket.
  • H: Deterministic part of the Hamiltonian.
  • Hs: Stochastic part(s) of the Hamiltonian (either an operator or a vector of operators).
  • fout=nothing: If given, this function fout(t, state) is called every time an output should be displayed. ATTENTION: The given state is neither normalized nor permanent!
  • normalize_state=false: Specify whether or not to normalize the state after each time step taken by the solver.
  • kwargs...: Further arguments are passed on to the ode solver.
source
stochastic.schroedinger_dynamic(tspan, state0, fdeterm, fstoch[; fout, ...])

Integrate stochastic Schrödinger equation with dynamic Hamiltonian.

Arguments

  • tspan: Vector specifying the points of time for which the output should be displayed.
  • psi0: Initial state.
  • fdeterm: Function f(t, psi, u) -> H returning the deterministic (time- or state-dependent) part of the Hamiltonian.
  • fstoch: Function f(t, psi, u, du) -> Hs returning a vector that contains the stochastic terms of the Hamiltonian.
  • fout=nothing: If given, this function fout(t, state) is called every time an output should be displayed. ATTENTION: The given state is neither normalized nor permanent!
  • noise_processes=0: Number of distinct white-noise processes in the equation. This number has to be equal to the total number of noise operators returned by fstoch. If unset, the number is calculated automatically from the function output. NOTE: Set this number if you want to avoid an initial calculation of the function output!
  • normalize_state=false: Specify whether or not to normalize the state after each time step taken by the solver.
  • kwargs...: Further arguments are passed on to the ode solver.
source
stochastic.homodyne_carmichael(H0, C, theta)

Helper function that defines the functions needed to compute homodyne detection trajectories according to Carmichael with stochastic.schroedinger_dynamic.

Arguments

  • H0: The deterministic, time-independent system Hamiltonian.
  • C: Collapse operator (or vector of operators) of the detected output channel(s).
  • theta: The phase difference between the local oscillator and the signal field. Defines the operator of the measured quadrature as $X_\theta = C e^{-i\theta} + C^\dagger e^{i\theta}$. Needs to be a vector of the same length as C if C is a vector.
  • normalize_expect=true: Specifiy whether or not to normalize the state vector when the expectation value in the nonlinear term is calculated. NOTE: should only be set to false if the state is guaranteed to be normalized, e.g. by setting normalize_state=true in stochastic.schroedinger_dynamic.

Returns (fdeterm, fstoch), where fdeterm(t, psi) -> H and fstoch(t, psi) -> Hs are functions returning the deterministic and stochastic part of the Hamiltonian required for calling stochastic.schroedinger_dynamic.

The deterministic and stochastic parts of the Hamiltonian are constructed as

\[H_{det} = H_0 + H_{nl},\]

where

\[H_{nl} = iCe^{-i\theta} \langle X_\theta \rangle - \frac{i}{2} C^\dagger C,\]

and

\[H_s = iCe^{-i\theta}.\]
source
stochastic.master(tspan, rho0, H, J, C; <keyword arguments>)

Time-evolution according to a stochastic master equation.

For dense arguments the master function calculates the non-hermitian Hamiltonian and then calls master_nh which is slightly faster.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • rho0: Initial density operator. Can also be a state vector which is automatically converted into a density operator.
  • H: Deterministic part of the Hamiltonian.
  • J: Vector containing all deterministic jump operators which can be of any arbitrary operator type.
  • C: Vector containing the stochastic operators for a superoperator of the form C[i]*rho + rho*Cdagger[i].
  • rates=nothing: Vector or matrix specifying the coefficients (decay rates) for the jump operators. If nothing is specified all rates are assumed to be 1.
  • Jdagger=dagger.(J): Vector containing the hermitian conjugates of the jump operators. If they are not given they are calculated automatically.
  • fout=nothing: If given, this function fout(t, rho) is called every time an output should be displayed. ATTENTION: The given state rho is not permanent! It is still in use by the ode solver and therefore must not be changed.
  • kwargs...: Further arguments are passed on to the ode solver.
source
stochastic.master_dynamic(tspan, rho0, fdeterm, fstoch; <keyword arguments>)

Time-evolution according to a stochastic master equation with a dynamic Hamiltonian and J.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • rho0: Initial density operator. Can also be a state vector which is automatically converted into a density operator.
  • fdeterm: Function f(t, rho) -> (H, J, Jdagger) or f(t, rho) -> (H, J, Jdagger, rates) giving the deterministic part of the master equation.
  • fstoch: Function f(t, rho) -> (C, Cdagger) giving the stochastic superoperator of the form C[i]*rho + rho*Cdagger[i].
  • rates=nothing: Vector or matrix specifying the coefficients (decay rates) for the jump operators. If nothing is specified all rates are assumed to be 1.
  • fout=nothing: If given, this function fout(t, rho) is called every time an output should be displayed. ATTENTION: The given state rho is not permanent! It is still in use by the ode solver and therefore must not be changed.
  • noise_processes=0: Number of distinct white-noise processes in the equation. This number has to be equal to the total number of noise operators returned by fstoch and all optional functions. If unset, the number is calculated automatically from the function outputs. NOTE: Set this number if you want to avoid an initial calculation of function outputs!
  • kwargs...: Further arguments are passed on to the ode solver.
source
semiclassical.schroedinger_stochastic(tspan, state0, fquantum, fclassical[; fout, ...])

Integrate time-dependent Schrödinger equation coupled to a classical system.

Arguments

  • tspan: Vector specifying the points of time for which the output should be displayed.
  • state0: Initial semi-classical state semiclassical.State.
  • fquantum: Function f(t, psi, u) -> H returning the time and or state dependent Hamiltonian.
  • fclassical: Function f(t, psi, u, du) calculating the possibly time and state dependent derivative of the classical equations and storing it in the vector du.
  • fstoch_quantum=nothing: Function f(t, psi, u) -> Hs that returns a vector of operators corresponding to the stochastic terms of the Hamiltonian. NOTE: Either this function or fstoch_classical has to be defined.
  • fstoch_classical=nothing: Function f(t, psi, u, du) that calculates the stochastic terms of the derivative du. NOTE: Either this function or fstoch_quantum has to be defined.
  • fout=nothing: If given, this function fout(t, state) is called every time an output should be displayed. ATTENTION: The given state is neither normalized nor permanent!
  • noise_processes=0: Number of distinct quantum noise processes in the equation. This number has to be equal to the total number of noise operators returned by fstoch. If unset, the number is calculated automatically from the function output. NOTE: Set this number if you want to avoid an initial calculation of the function output!
  • noise_prototype_classical=nothing: The equivalent of the optional argument noise_rate_prototype in StochasticDiffEq for the classical stochastic function fstoch_classical only. Must be set for non-diagonal classical noise or combinations of quantum and classical noise. See the documentation for details.
  • normalize_state=false: Specify whether or not to normalize the state after each time step taken by the solver.
  • kwargs...: Further arguments are passed on to the ode solver.
source
stochastic.master_semiclassical(tspan, rho0, H, Hs, J; <keyword arguments>)

Time-evolution according to a stochastic master equation.

For dense arguments the master function calculates the non-hermitian Hamiltonian and then calls master_nh which is slightly faster.

Arguments

  • tspan: Vector specifying the points of time for which output should be displayed.
  • rho0: Initial density operator. Can also be a state vector which is automatically converted into a density operator.
  • fquantum: Function f(t, rho, u) -> (H, J, Jdagger) or f(t, rho, u) -> (H, J, Jdagger, rates) giving the deterministic part of the master equation.
  • fclassical: Function f(t, rho, u, du) that calculates the classical derivatives du.
  • fstoch_quantum=nothing: Function f(t, rho, u) -> C, Cdagger that returns the stochastic operator for the superoperator of the form C[i]*rho + rho*Cdagger[i].
  • fstoch_classical=nothing: Function f(t, rho, u, du) that calculates the stochastic terms of the derivative du.
  • rates=nothing: Vector or matrix specifying the coefficients (decay rates) for the jump operators. If nothing is specified all rates are assumed to be 1.
  • fout=nothing: If given, this function fout(t, rho) is called every time an output should be displayed. ATTENTION: The given state rho is not permanent! It is still in use by the ode solver and therefore must not be changed.
  • noise_processes=0: Number of distinct quantum noise processes in the equation. This number has to be equal to the total number of noise operators returned by fstoch. If unset, the number is calculated automatically from the function output. NOTE: Set this number if you want to avoid an initial calculation of the function output!
  • noise_prototype_classical=nothing: The equivalent of the optional argument noise_rate_prototype in StochasticDiffEq for the classical stochastic function fstoch_classical only. Must be set for non-diagonal classical noise or combinations of quantum and classical noise. See the documentation for details.
  • kwargs...: Further arguments are passed on to the ode solver.
source

State definitions

randstate(basis)

Calculate a random normalized ket state.

source
thermalstate(H,T)

Thermal state $exp(-H/T)/Tr[exp(-H/T)]$.

source
coherentthermalstate(basis::FockBasis,H,T,alpha)

Coherent thermal state $D(α)exp(-H/T)/Tr[exp(-H/T)]D^†(α)$.

source
phase_average(rho)

Returns the phase-average of $ρ$ containing only the diagonal elements.

source
passive_state(rho,IncreasingEigenenergies::Bool=true)

Passive state $π$ of $ρ$. IncreasingEigenenergies=true means that higher indices correspond to higher energies.

source

Printing

QuantumOptics.set_printing(; standard_order, rounding_tol)

Set options for REPL output.

Arguments

  • standard_order=false: For performance reasons, the order of the tensor product is inverted, i.e. tensor(a, b)=kron(b, a). When changing this to true, the output shown in the REPL will exhibit the correct order.
  • rounding_tol=1e-17: Tolerance for floating point errors shown in the output.
source