API
Spectral to Correlation Function
SmoQySynthAC.spectral_to_imaginary_time_correlation_function
— Functionspectral_to_imaginary_time_correlation_function(;
# KEYWORD ARGUMENTS
τ::AbstractVector{T},
β::T,
spectral_function::Function,
kernel_function::Function,
tol::T = 1e-10,
bounds = (-Inf, +Inf),
normalize::Bool = false
) where {T<:AbstractFloat}
Calculate and return the imaginary-time correlation function
\[C(\tau) = \int_{-\infty}^{\infty} d\omega \ K_\beta(\omega, \tau) \ A(\omega).\]
on a grid of $\tau$ (τ
) values, given a spectral function $A(\omega)$ (spectral_function
) and kernel function $K_\beta(\omega,\tau)$ (kernel_function
). This integral is evaluated within a specified tolerance tol
.
Arguments
τ::AbstractVector{T}
: Vector of imaginary time such thatτ[end] = β
equal the inverse temperature.spectral_function::Function
: The spectral function $A(\omega)$ that takes a single argument.kernel_function::Function
: The kernel function $K_\beta(\omega,\tau)$ that takes three arguments as shown.tol::T = 1e-10
: Specified precision with which $C(\tau)$ is evaluatedbounds = (-Inf, +Inf)
: Bounds on integration domain.normalize::Bool = false
: Whether the spectral function needs to be normalized to unity.
SmoQySynthAC.spectral_to_matsubara_correlation_function
— Functionspectral_to_matsubara_correlation_function(;
# KEYWORD ARGUMENTS
n::AbstractVector{Int},
β::T,
spectral_function::Function,
kernel_function::Function,
tol::T = 1e-10,
bounds = (-Inf, 0.0, +Inf),
normalize::Bool = false
) where {T<:AbstractFloat}
Calculate and return the Matsubara correlation function
\[C({\rm i} \omega_n) = \int_{-\infty}^{\infty} d\omega \ K_\beta(\omega, \omega_n) \ A(\omega).\]
for a vector of $n \in \mathbb{Z}$ values, given a spectral function $A(\omega)$ (spectral_function
) and kernel function $K_\beta(\omega,\omega_n)$ (kernel_function
). This integral is evaluated within a specified tolerance tol
.
Note that the kernel function should be called as kernel_function(ω, n, β)
where n
specifies the Matsubara frequency, which is evaluated internally as either $\omega_n = (2n+1)\pi/\beta$ or $\omega_n = 2n\pi/\beta$ depending on whether the kernel function is fermionic of bosonic respectively.
Arguments
n::AbstractVector{Int}
: Vector of integers specifying Matsubara frequencies for which $C({\rm i}\omega_n)$ will be evaluated.spectral_function::Function
: The spectral function $A(\omega)$ that takes a single argument.kernel_function::Function
: The kernel function $K_\beta(\omega,\omega_n)$ that takes three arguments as shown.tol::T = 1e-10
: Specified precision with which $C({\rm i}\omega_n)$ is evaluated.bounds = (-Inf, +Inf)
: Bounds on integration domain.normalize::Bool = false
: Whether the spectral function needs to be normalized to unity.
Add Noise to an Imaginary Time Correlation Function
SmoQySynthAC.add_noise
— Functionadd_noise(;
# KEYWORD ARGUMENTS
Cτ_exact::AbstractVector{T},
τ::AbstractVector{T},
σ::T,
ξ::T,
sum_rule::Function = C0 -> 1 - C0,
noise_type::Symbol = :TruncatedNormal
) where {T<:AbstractFloat}
add_noise(
# ARGUMENTS
N_samples::Int;
# KEYWORD ARGUMENTS
Cτ_exact::AbstractVector{T},
τ::AbstractVector{T},
σ::T,
ξ::T,
sum_rule::Function = C0 -> 1 - C0,
noise_type::Symbol = :TruncatedNormal
) where {T<:AbstractFloat}
Add noise to an imaginary time correlation function $C(\tau)$ that is exponentially correlated in imaginary time.
Arguments
N_samples
(optional): Number of random samples to generate. If passed this function returns a matrix where the columsn correspond to the different samples.
Keyword Arguments
Cτ_exact::AbstractVector{T}
: Vector containing the exact values for $C(\tau)$.τ::AbstractVector{T}
: Vector specifying the imaginary time $\tau$ grid that $C(\tau)$ is evaluated on. Assumes that the last element equals the inverse temperature, i.e.τ[end] = β
.σ::T
: Standard deviation of the noise; controls the typical amplitude of the error.ξ::T
: Correlation length associated with the noise in imaginary time.sum_rule::Function = C0 -> 1 - C0
: Enforces sum rule, or bounday condition in imaginay time. Default behavior assumes a fermionic correlation function, enforcing that $C(\beta) = 1 - C(0)$.noise_type::Symbol = :TruncatedNormal
: Distribution that the noise is sampled from prior to correlations being introduced in imaginary time. Available options arenoise_type ∈ (:TruncatedNormal, :Gamma, :Normal)
.
Additional Information
If the correlation function $C(\tau_i)$ the corresponding noisy correlation function is given by
\[C_{\rm noisy}(\tau_i) = C(\tau_i) + \frac{\sum_j e^{-|\tau_j-\tau_i|/\xi} R_j}{\sum_j e^{-2|\tau_j-\tau_i|/\xi}},\]
where the sum is performed assuming periodic boundary conditions. If noise_type = :Normal
then the random numbers of sampled according to $R_j \sim {\rm Normal}(0,\sigma)$. Otherwise,
\[R_j \sim {\rm sign}(C(\tau_j)) \cdot P(|C(\tau_j)|, \sigma) - C(\tau_j)\]
where $P(|C(\tau_j)|, \sigma)$ is either a Gamma or Truncated Normal distribution with mean and standard deviation given by $|C(\tau_j)|$ and $\sigma$ respectively.
The support for the Truncated Normal and Gamma distributions is $[0,\infty)$. Note that in the case that a Normal distribution is used it is possible for $C_{\rm noisy}(\tau_i)$ to have a different sign that $C(\tau_i)$, whereas this is not possible with the other two distributions.
SmoQySynthAC.add_noise!
— Functionadd_noise!(
# ARGUMENTS
Cτ_noisy::AbstractVector{T};
# KEYWORD ARGUMENTS
Cτ_exact::AbstractVector{T},
τ::AbstractVector{T},
σ::T,
ξ::T,
sum_rule::Function = C0 -> 1 - C0,
noise_type::Symbol = :TruncatedNormal
) where {T<:AbstractFloat}
add_noise!(
# ARGUMENTS
Cτ_noisy::AbstractMatrix{T};
# KEYWORD ARGUMENTS
Cτ_exact::AbstractVector{T},
τ::AbstractVector{T},
σ::T,
ξ::T,
sum_rule::Function = C0 -> 1 - C0,
noise_type::Symbol = :TruncatedNormal
) where {T<:AbstractFloat}
Add noise to an imaginary time correlation function $C(\tau)$ that is exponentially correlated in imaginary time. See add_noise
for more information.
Arguments
Cτ_noisy
: The array to which the noisy $C(\tau)$ data will be written. If a matrix and not a vector then the columns correspond to samples and the rows to the imaginary time slices $\tau$.
Keyword Arguments
Cτ_exact::AbstractVector{T}
: Vector containing the exact values for $C(\tau)$.τ::AbstractVector{T}
: Vector specifying the imaginary time $\tau$ grid that $C(\tau)$ is evaluated on. Assumes that the last element equals the inverse temperature, i.e.τ[end] = β
.σ::T
: Standard deviation of the noise; controls the typical amplitude of the error.ξ::T
: Correlation length associated with the noise in imaginary time.sum_rule::Function = C0 -> 1 - C0
: Enforces sum rule, or bounday condition in imaginay time. Default behavior assumes a fermionic correlation function, enforcing that $C(\beta) = 1 - C(0)$.noise_type::Symbol = :TruncatedNormal
: Distribution that the noise is sampled from prior to correlations being introduced in imaginary time. Available options arenoise_type ∈ (:TruncatedNormal, :Gamma, :Normal)
.
Kernel Functions
kernel_mat
kernel_tau_fermi
kernel_mat_fermi
kernel_tau_bose
kernel_mat_bose
kernel_tau_bose_alt
kernel_mat_bose_alt
kernel_tau_sym_bose
kernel_mat_sym_bose
kernel_tau_sym_bose_alt
kernel_mat_sym_bose_alt
SmoQySynthAC.kernel_mat
— Functionkernel_mat(ω::T, ω_n::T) where {T<:AbstractFloat}
Generic Matsubara frequency kernel
\[K_\beta(\omega, \omega_n) = \frac{1}{\omega - {\rm i}\omega_n},\]
where $\omega_n = (2n+1)\pi/\beta \ (= 2n\pi/\beta)$ for fermions (bosons) with $n \in \mathbb{Z}$. This kernel assumes the sign convention
\[C(\tau) = \langle \hat{a}(\tau) \hat{a}^\dagger(0) \rangle\]
for a the Green's function, where it is assumed that $\tau \in [0,\beta)$.
Fermionic Kernel Functions
SmoQySynthAC.kernel_tau_fermi
— Functionkernel_tau_fermi(ω::T, τ::T, β::T) where {T<:AbstractFloat}
The imaginary time fermionic kernel
\[\begin{align*} K_\beta(\omega,\tau) & = \overbrace{\left(\frac{e^{-\tau\omega}}{1+e^{-\beta\omega}}\right)}^{\text{numerically unstable}} \\ & = \underbrace{\left( e^{\tau\omega} + e^{(\tau-\beta)\omega} \right)^{-1}}_{\text{numerically stable}}, \end{align*}\]
where it is assumed that $\tau \in [0,\beta)$.
SmoQySynthAC.kernel_mat_fermi
— Functionkernel_mat_fermi(ω::T, n::Int, β::T) where {T<:AbstractFloat}
The fermionic matsubara frequency kernel
\[K_\beta(\omega, \omega_n) = \frac{1}{{\omega - {\rm i}\omega_n}},\]
where $\omega_n = (2n+1)\pi/\beta$ for fermions with $n \in \mathbb{Z}$.
Bosonic Kernel Functions
SmoQySynthAC.kernel_tau_bose
— Functionkernel_tau_bose(ω::T, τ::T, β::T) where {T<:AbstractFloat}
The imaginary time bosonic kernel
\[\begin{align*} K_\beta(\omega,\tau) & = \overbrace{\left(\frac{e^{-\tau\omega}}{1-e^{-\beta\omega}}\right)}^{\text{numerically unstable}} \\ & = \underbrace{\left( e^{\tau\omega} - e^{(\tau-\beta)\omega} \right)^{-1}}_{\text{numerically stable}}, \end{align*}\]
where it is assumed that $\tau \in [0,\beta)$.
SmoQySynthAC.kernel_mat_bose
— Functionkernel_mat_bose(ω::T, n::Int, β::T) where {T<:AbstractFloat}
The bosonic matsubara frequency kernel
\[K_\beta(\omega, \omega_n) = \frac{1}{\omega - {\rm i}\omega_n},\]
where $\omega_n = 2n\pi/\beta$ for bosons with $n \in \mathbb{Z}$.
Modified Bosonic Kernel Functions
SmoQySynthAC.kernel_tau_bose_alt
— Functionkernel_tau_bose_alt(ω::T, τ::T, β::T) where {T<:AbstractFloat}
The imaginary time bosonic kernel
\[\begin{align*} K_\beta(\omega,\tau) & = \overbrace{\left(\frac{\omega e^{-\tau\omega}}{1-e^{-\beta\omega}}\right)}^{\text{numerically unstable}} \\ & = \underbrace{\omega \left( e^{\tau\omega} - e^{(\tau-\beta)\omega} \right)^{-1}}_{\text{numerically stable}}, \end{align*}\]
where it is assumed that $\tau \in [0,\beta)$ and $K_\beta(0,\tau) = \beta^{-1}$.
SmoQySynthAC.kernel_mat_bose_alt
— Functionkernel_mat_bose_alt(ω::T, n::Int, β::T) where {T<:AbstractFloat}
The bosonic matsubara frequency kernel
\[K_\beta(\omega, \omega_n) = \frac{\omega}{\omega - {\rm i}\omega_n},\]
where $\omega_n = 2n\pi/\beta$ for bosons with $n \in \mathbb{Z}$ and $K_\beta(0,0) = -1$.
Symmetric Bosonic Kernel Function
SmoQySynthAC.kernel_tau_sym_bose
— Functionkernel_tau_sym_bose(ω::T, τ::T, β::T) where {T<:AbstractFloat}
The imaginary time symmetrized bosonic kernel
\[\begin{align*} K_\beta(\omega,\tau) & = \overbrace{\left(\frac{e^{-\tau\omega} + e^{-(\beta-\tau)\omega}}{1-e^{-\beta\omega}}\right)}^\text{numerically unstable} \\ & = \underbrace{\left( e^{\tau\omega} - e^{(\tau-\beta)\omega} \right)^{-1} - \left(e^{-\tau\omega} - e^{-(\tau-\beta)\omega} \right)^{-1}}_{\text{numerically stable}}, \end{align*}\]
where it is assumed that $\tau \in [0,\beta)$.
SmoQySynthAC.kernel_mat_sym_bose
— Functionkernel_mat_sym_bose(ω::T, n::Int, β::T) where {T<:AbstractFloat}
The symmetrized bosonic matsubara frequency kernel
\[\begin{align*} K_β(\omega, \omega_n) & = \frac{2\omega}{\omega_n^2 + \omega^2}, \end{align*}\]
where $\omega_n = 2n\pi/\beta$ for bosons with $n \in \mathbb{Z}$.
Modified Symmetric Bosonic Kernel Function
SmoQySynthAC.kernel_tau_sym_bose_alt
— Functionkernel_tau_sym_bose_alt(ω::T, τ::T, β::T) where {T<:AbstractFloat}
The imaginary time symmetrized bosonic kernel
\[\begin{align*} K_\beta(\omega,\tau) & = \overbrace{\omega \left(\frac{e^{-\tau\omega} + e^{-(\beta-\tau)\omega}}{1-e^{-\beta\omega}}\right)}^\text{numerically unstable} \\ & = \underbrace{\omega\left[\left( e^{\tau\omega} - e^{(\tau-\beta)\omega} \right)^{-1} - \left(e^{-\tau\omega} - e^{-(\tau-\beta)\omega} \right)^{-1}\right]}_{\text{numerically stable}}, \end{align*}\]
where it is assumed that $\tau \in [0,\beta)$ and $K_\beta(0,\tau) = 2/\beta$.
SmoQySynthAC.kernel_mat_sym_bose_alt
— Functionkernel_mat_sym_bose_alt(ω::T, n::Int, β::T) where {T<:AbstractFloat}
The symmetrized bosonic matsubara frequency kernel
\[\begin{align} K_β(\omega, \omega_n) & = \frac{2\omega^2}{\omega_n^2 + \omega^2}, \end{align}\]
where $\omega_n = 2n\pi/\beta$ for bosons with $n \in \mathbb{Z}$ and $K_\beta(0,0) = 2$.
Quantum Statistics Functions
SmoQySynthAC.fermi
— Functionfermi(ϵ::T, β::T) where {T<:AbstractFloat}
The Fermi-Dirac function
\[\begin{align*} f_\beta(\epsilon) & = \overbrace{\left( \frac{1}{e^{\beta\epsilon} + 1}\right)}^{\text{numerically unstable}} \\ & = \underbrace{\frac{1}{2}\left( 1 - \tanh\left(\frac{\beta\epsilon}{2}\right) \right)}_{\text{numerically stable}}, \end{align*}\]
where $\epsilon$ is energy and $\beta$ is inverse temperature.
SmoQySynthAC.bose
— Functionbose(ϵ::T, β::T) where {T<:AbstractFloat}
The Bose-Einstein function
\[\begin{align*} n_\beta(\epsilon) & = \overbrace{\left( \frac{1}{e^{\beta\epsilon} - 1}\right)}^{\text{numerically unstable}} \\ & = \underbrace{\frac{1}{2}\left(\coth\left(\frac{\beta\epsilon}{2}\right) - 1 \right)}_{\text{numerically stable}}, \end{align*}\]
where $\epsilon$ is energy and $\beta$ is inverse temperature.