Discrete Lehmann representation

function DLRGrid(Euv, β, rtol, isFermi::Bool; symmetry::Symbol = :none, rebuild = false, folder = nothing, algorithm = :functional, verbose = false) function DLRGrid(; isFermi::Bool, β = -1.0, beta = -1.0, Euv = 1.0, symmetry::Symbol = :none, rtol = 1e-14, rebuild = false, folder = nothing, algorithm = :functional, verbose = false)

Create DLR grids

#Arguments:

• Euv : the UV energy scale of the spectral density
• β or beta : inverse temeprature
• isFermi : bool is fermionic or bosonic
• symmetry : particle-hole symmetric :ph, or particle-hole asymmetric :pha, or :none
• rtol : tolerance absolute error
• rebuild : set false to load DLR basis from the file, set true to recalculate the DLR basis on the fly
• folder : if rebuild is true and folder is set, then dlrGrid will be rebuilt and saved to the specified folder if rebuild is false and folder is set, then dlrGrid will be loaded from the specified folder
• algorithm : if rebuild = true, then set :functional to use the functional algorithm to generate the DLR basis, or set :discrete to use the matrix algorithm.
• verbose : false not to print DLRGrid to terminal, or true to print

struct DLRGrid

DLR grids for imaginary-time/Matsubara frequency correlators

#Members:

• isFermi: bool is fermionic or bosonic
• symmetry: particle-hole symmetric :ph, or particle-hole asymmetric :pha, or :none
• Euv : the UV energy scale of the spectral density
• β or beta : inverse temeprature
• Λ or lambda: cutoff = UV Energy scale of the spectral density * inverse temperature
• rtol: tolerance absolute error
• size : number of DLR basis
• ω or omega : selected representative real-frequency grid
• n : selected representative Matsubara-frequency grid (integer)
• ωn or omegaN : (2n+1)π/β
• τ or tau : selected representative imaginary-time grid

Base.size(dlrGrid::DLRGrid) = (length(dlrGrid.ω),) Base.length(dlrGrid::DLRGrid) = length(dlrGrid.ω) rank(dlrGrid::DLRGrid) = length(dlrGrid.ω)

get the rank of the DLR grid, namely the number of the DLR coefficients.

function dlr2matfreq(dlrGrid::DLRGrid, dlrcoeff, nGrid = dlrGrid.n; axis = 1, verbose = true)

DLR representation to Matsubara-frequency representation

#Members:

• dlrGrid : DLRGrid
• dlrcoeff : DLR coefficients
• nGrid : expected fine Matsubara-freqeuncy grids (integer)
• axis : Matsubara-frequency axis in the data dlrcoeff
• verbose : true to print warning information

function dlr2tau(dlrGrid::DLRGrid, dlrcoeff, τGrid = dlrGrid.τ; axis = 1, verbose = true)

DLR representation to imaginary-time representation

#Members:

• dlrGrid : DLRGrid
• dlrcoeff : DLR coefficients
• τGrid : expected fine imaginary-time grids
• axis : imaginary-time axis in the data dlrcoeff
• verbose : true to print warning information

function matfreq2dlr(dlrGrid::DLRGrid, green, nGrid = dlrGrid.n; error = nothing, axis = 1, sumrule = nothing, verbose = true)

Matsubara-frequency representation to DLR representation

#Members:

• dlrGrid : DLRGrid struct.
• green : green's function in Matsubara-frequency domain
• nGrid : the n grid that Green's function is defined on.
• error : error the Green's function.
• axis : the Matsubara-frequency axis in the data green
• sumrule : enforce the sum rule
• verbose : true to print warning information

function matfreq2matfreq(dlrGrid, green, nNewGrid, nGrid = dlrGrid.n; error = nothing, axis = 1, sumrule = nothing, verbose = true)

Fourier transform from Matsubara-frequency to imaginary-time using the DLR representation

#Members:

• dlrGrid : DLRGrid
• green : green's function in Matsubara-freqeuncy repsentation
• nNewGrid : expected fine Matsubara-freqeuncy grids (integer)
• nGrid : the n grid that Green's function is defined on.
• error : error the Green's function.
• axis : Matsubara-frequency axis in the data green
• sumrule : enforce the sum rule
• verbose : true to print warning information

function matfreq2tau(dlrGrid, green, τNewGrid = dlrGrid.τ, nGrid = dlrGrid.n; error = nothing, axis = 1, sumrule = nothing, verbose = true)

Fourier transform from Matsubara-frequency to imaginary-time using the DLR representation

#Members:

• dlrGrid : DLRGrid
• green : green's function in Matsubara-freqeuncy repsentation
• τNewGrid : expected fine imaginary-time grids
• nGrid : the n grid that Green's function is defined on.
• error : error the Green's function.
• axis : Matsubara-frequency axis in the data green
• sumrule : enforce the sum rule
• verbose : true to print warning information

function tau2dlr(dlrGrid::DLRGrid, green, τGrid = dlrGrid.τ; error = nothing, axis = 1, sumrule = nothing, verbose = true)

imaginary-time domain to DLR representation

#Members:

• dlrGrid : DLRGrid struct.
• green : green's function in imaginary-time domain.
• τGrid : the imaginary-time grid that Green's function is defined on.
• error : error the Green's function.
• axis : the imaginary-time axis in the data green.
• sumrule : enforce the sum rule
• verbose : true to print warning information

function tau2matfreq(dlrGrid, green, nNewGrid = dlrGrid.n, τGrid = dlrGrid.τ; error = nothing, axis = 1, sumrule = nothing, verbose = true)

Fourier transform from imaginary-time to Matsubara-frequency using the DLR representation

#Members:

• dlrGrid : DLRGrid
• green : green's function in imaginary-time domain
• nNewGrid : expected fine Matsubara-freqeuncy grids (integer)
• τGrid : the imaginary-time grid that Green's function is defined on.
• error : error the Green's function.
• axis : the imaginary-time axis in the data green
• sumrule : enforce the sum rule
• verbose : true to print warning information

function tau2tau(dlrGrid, green, τNewGrid, τGrid = dlrGrid.τ; error = nothing, axis = 1, sumrule = nothing, verbose = true)

Interpolation from the old imaginary-time grid to a new grid using the DLR representation

#Members:

• dlrGrid : DLRGrid
• green : green's function in imaginary-time domain
• τNewGrid : expected fine imaginary-time grids
• τGrid : the imaginary-time grid that Green's function is defined on.
• error : error the Green's function.
• axis : the imaginary-time axis in the data green
• sumrule : enforce the sum rule
• verbose : true to print warning information

Spectral representation related functions

boseEinstein(ω)

Compute the Fermi Dirac function. Assume $k_B T/\hbar=1$

\[f(ω) = 1/(e^{ωβ}-1)\]

Arguments

• ω: frequency
• β: the inverse temperature

density(isFermi::Bool, ω, β)

Compute the imaginary-time kernel of different type. Assume $k_B T/\hbar=1$

Arguments

• isFermi: fermionic or bosonic
• ω: energy
• β: the inverse temperature

fermiDirac(ω)

Compute the Fermi Dirac function. Assume $k_B T/\hbar=1$

\[f(ω) = 1/(e^{ωβ}+1)\]

Arguments

• ω: frequency
• β: the inverse temperature

kernelBoseT(τ, ω, β)

Compute the imaginary-time bosonic kernel. Machine accuracy ~eps(g) is guaranteed``

\[g(τ>0) = e^{-ωτ}/(1-e^{-ωβ}), g(τ≤0) = -e^{-ωτ}/(1-e^{ωβ})\]

Arguments

• τ: the imaginary time, must be (-β, β]
• ω: frequency
• β: the inverse temperature

kernelBoseT_PH(τ, ω, β)

Compute the imaginary-time kernel for correlation function $⟨O(τ)O(0)⟩$. Machine accuracy ~eps(C) is guaranteed``

\[K(τ) = e^{-ω|τ|}+e^{-ω(β-|τ|)}\]

Arguments

• τ: the imaginary time, must be (-β, β]
• ω: frequency, ω>=0
• β: the inverse temperature

kernelBoseT_PHA(τ, ω, β)

Compute the imaginary-time kernel for correlation function $⟨O(τ)O(0)⟩$. Machine accuracy ~eps(C) is guaranteed``

\[K(τ) = e^{-ω|τ|}-e^{-ω(β-|τ|)}\]

Arguments

• τ: the imaginary time, must be (0, β]
• ω: frequency, ω>=0
• β: the inverse temperature

kernelBoseT_regular(τ, ω, β)

Compute the imaginary-time bosonic kernel with a regulator near ω=0. Machine accuracy ~eps(g) is guaranteed``

\[g(τ>0) = e^{-ωτ}/(1+e^{-ωβ}), g(τ≤0) = e^{-ωτ}/(1+e^{ωβ})\]

Arguments

• τ: the imaginary time, must be (-β, β]
• ω: frequency
• β: the inverse temperature

kernelBoseΩ(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the bosonic kernel with Matsubara frequency.

\[g(iω_n) = -1/(iω_n-ω),\]

where $ω_n=2nπ/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelBoseΩ_PH(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the Matsubara-frequency kernel for a correlator $⟨O(τ)O(0)⟩_{iω_n}$.

\[K(iω_n) = \frac{2ω}{ω^2+ω_n^2}(1-e^{-ωβ}),\]

where $ω_n=2nπ/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelBoseΩ_PHA(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the Matsubara-frequency kernel for a anormalous fermionic correlator with particle-hole symmetry.

\[K(iω_n) = -\frac{2iω_n}{ω^2+ω_n^2}(1-e^{-ωβ}),\]

where $ω_n=2nπ/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the fermionic Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelBoseΩ_regular(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the bosonic kernel in Matsubara frequency with a regulartor near ω=0

\[g(iω_n) = -(1-e^{-ωβ})/(1+e^{-ωβ})/(iω_n-ω),\]

where $ω_n=2nπ/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelFermiT(τ, ω, β)

Compute the imaginary-time fermionic kernel. Machine accuracy ~eps(g) is guaranteed``

\[g(τ>0) = e^{-ωτ}/(1+e^{-ωβ}), g(τ≤0) = -e^{-ωτ}/(1+e^{ωβ})\]

Arguments

• τ: the imaginary time, must be (-β, β]
• ω: frequency
• β: the inverse temperature

kernelFermiT_PH(τ, ω, β)

Compute the imaginary-time kernel for correlation function $⟨O(τ)O(0)⟩$. Machine accuracy ~eps(C) is guaranteed``

\[K(τ) = e^{-ω|τ|}+e^{-ω(β-|τ|)}\]

Arguments

• τ: the imaginary time, must be (-β, β]
• ω: frequency, ω>=0
• β: the inverse temperature

kernelFermiT_PHA(τ, ω, β)

Compute the imaginary-time kernel for correlation function $⟨O(τ)O(0)⟩$. Machine accuracy ~eps(C) is guaranteed``

\[K(τ) = e^{-ω|τ|}-e^{-ω(β-|τ|)}\]

Arguments

• τ: the imaginary time, must be (0, β]
• ω: frequency, ω>=0
• β: the inverse temperature

kernelFermiΩ(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the fermionic kernel with Matsubara frequency.

\[g(iω_n) = -1/(iω_n-ω),\]

where $ω_n=(2n+1)π/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelFermiΩ_PH(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the Matsubara-frequency kernel for a correlator $⟨O(τ)O(0)⟩_{iω_n}$.

\[K(iω_n) = -\frac{2iω_n}{ω^2+ω_n^2}(1+e^{-ωβ}),\]

where $ω_n=(2n+1)π/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelFermiΩ_PHA(n::Int, ω::T, β::T) where {T <: AbstractFloat}

Compute the Matsubara-frequency kernel for a anormalous fermionic correlator with particle-hole symmetry.

\[K(iω_n) = \frac{2ω}{ω^2+ω_n^2}(1+e^{-ωβ}),\]

where $ω_n=(2n+1)π/β$. The convention here is consist with the book "Quantum Many-particle Systems" by J. Negele and H. Orland, Page 95

Arguments

• n: index of the fermionic Matsubara frequency
• ω: energy
• β: the inverse temperature

kernelT(isFermi, symmetry, τGrid::AbstractVector{T}, ωGrid::AbstractVector{T}, β::T, regularized::Bool = false; type = T) where {T<:AbstractFloat}

Compute kernel with given τ and ω grids.

kernelT(::Val{isFermi}, ::Val{symmetry}, τ::T, ω::T, β::T, regularized::Bool = false) where {T<:AbstractFloat}

Compute the imaginary-time kernel of different type.

Arguments

• isFermi: fermionic or bosonic. It should be wrapped as Val(isFermi).
• symmetry: :ph, :pha, or :none. It should be wrapped as Val(symmetry).
• τ: the imaginary time, must be (-β, β]
• ω: energy
• β: the inverse temperature
• regularized: use regularized kernel or not

kernelΩ(isFermi, symmetry, nGrid::Vector{Int}, ωGrid::Vector{T}, β::T, regularized::Bool = false; type = Complex{T}) where {T<:AbstractFloat}

Compute kernel matrix with given ωn (integer!) and ω grids.

kernelΩ(::Val{isFermi}, ::Val{symmetry}, n::Int, ω::T, β::T, regularized::Bool = false) where {T<:AbstractFloat}

Compute the imaginary-time kernel of different type. Assume $k_B T/\hbar=1$

Arguments

• isFermi: fermionic or bosonic. It should be wrapped as Val(isFermi).
• symmetry: :ph, :pha, or :none. It should be wrapped as Val(symmetry).
• n: index of the Matsubara frequency
• ω: energy
• β: the inverse temperature
• regularized: use regularized kernel or not

MultiPole(β, isFermi::Bool, symmetry::Symbol, Grid, type::Symbol, poles, regularized::Bool = true)
MultiPole(dlr, type::Symbol, poles, Grid = dlrGrid(dlr, type); regularized::Bool = true)

Generate Green's function from a spectral density with delta peaks specified by the argument $poles$. Return the function on Grid and the systematic error.

#Arguments

• dlr: dlrGrid struct
• β : inverse temperature
• isFermi: is fermionic or bosonic
• symmetry: particle-hole symmetric :ph, particle-hole antisymmetric :pha, or :none
• Grid: grid to evalute on
• type: imaginary-time with :τ, or Matsubara-frequency with :n
• poles: a list of frequencies for the delta functions
• regularized: use regularized bosonic kernel if symmetry = :none

SemiCircle(Euv, β, isFermi::Bool, Grid, type::Symbol, symmetry::Symbol = :none; rtol = nothing, degree = 24, regularized::Bool = true)
SemiCircle(dlr, type::Symbol, Grid = dlrGrid(dlr, type); degree = 24, regularized::Bool = true)

Generate Green's function from a semicircle spectral density. Return the function on Grid and the systematic error.

#Arguments

• dlr: dlrGrid struct
• Euv : ultraviolet energy cutoff
• β : inverse temperature
• isFermi: is fermionic or bosonic
• Grid: grid to evalute on
• type: imaginary-time with :τ, or Matsubara-frequency with :n
• symmetry: particle-hole symmetric :ph, particle-hole antisymmetric :pha, or :none
• rtol: accuracy to achieve
• degree: polynomial degree for integral
• regularized: use regularized bosonic kernel if symmetry = :none

function build(dlrGrid, print::Bool = true) Construct discrete Lehmann representation

#Arguments:

• dlrGrid: struct that contains the information to construct the DLR grid. The following entries are required: Λ: the dimensionless scale β*Euv, rtol: the required relative accuracy, isFermi: fermionic or bosonic, symmetry: particle-hole symmetry/antisymmetry or none
• print: print the internal information or not

function preciseKernelT(dlrGrid, τ, ω, print::Bool = true)

Calculate the kernel matrix(τ, ω) for given τ, ω grids

Arguments

• τ: a CompositeChebyshevGrid struct or a simple one-dimensional array
• ω: a CompositeChebyshevGrid struct or a simple one-dimensional array
• print: print information or not

function build(dlrGrid, print::Bool = true) Construct discrete Lehmann representation

#Arguments:

• dlrGrid: struct that contains the information to construct the DLR grid. The following entries are required: Λ: the dimensionless scale β*Euv, rtol: the required relative accuracy, isFermi: fermionic or bosonic, symmetry: particle-hole symmetry/antisymmetry or none
• print: print the internal information or not

(1-exp(-Λ*x)/x

modified Gram-Schmidt process

<K(gi), K(gj)>

particle-hole asymmetric kernel: K(ω, τ)=e^{-ωτ}-e^{-ω(β-τ)}

KK=int_0^{Λ} dτ K(ω,t1)*K(ω2,t2)=(1-e^{t1+t2})/(t1+t2)+(1-e^{2β-t1-t2})/(2β-t1-t2)-(1-e^{β+t1-t2})/(β+t1-t2)-(1-e^{β-t1+t2})/(β-t1+t2)

particle=hole asymmetric kernel: K(ω, τ)=e^{-ωτ}-e^{-ω(β-τ)}

KK=int_0^{1/2} dτ K(ω1,τ)*K(ω2,τ)=(1-e^{ω1+ω2})/(ω1+ω2)-(e^{-ω2}-e^{-ω1})/(ω1-ω2)

particle-hole symmetric kernel: K(ω, τ)=e^{-ωτ}+e^{-ω(β-τ)}

KK=int_0^{Λ} dτ K(ω,t1)*K(ω2,t2)=(1-e^{t1+t2})/(t1+t2)+(1-e^{2β-t1-t2})/(2β-t1-t2)+(1-e^{β+t1-t2})/(β+t1-t2)+(1-e^{β-t1+t2})/(β-t1+t2)

particle-hole symmetric kernel: K(ω, τ)=e^{-ωτ}+e^{-ω(β-τ)}

KK=int_0^{1/2} dτ K(ω1,τ)*K(ω2,τ)=(1-e^{ω1+ω2})/(ω1+ω2)+(e^{-ω2}-e^{-ω1})/(ω1-ω2)

q1=sumj cj Kj q2=sumk dk Kk return <q1, q2> = sumjk cj*dk <Kj, K_k>

function barycheb(n, x, f, wc, xc)

Barycentric Lagrange interpolation at Chebyshev nodes Reference: Berrut, J.P. and Trefethen, L.N., 2004. Barycentric lagrange interpolation. SIAM review, 46(3), pp.501-517.

Arguments

• n: order of the Chebyshev interpolation
• x: coordinate to interpolate
• f: array of size n, function at the Chebyshev nodes
• wc: array of size n, Barycentric Lagrange interpolation weights
• xc: array of size n, coordinates of Chebyshev nodes

Returns

• Interpolation result

barychebinit(n)

Get Chebyshev nodes of first kind and corresponding barycentric Lagrange interpolation weights. Reference: Berrut, J.P. and Trefethen, L.N., 2004. Barycentric lagrange interpolation. SIAM review, 46(3), pp.501-517.

Arguments

• n: order of the Chebyshev interpolation

Returns

• Chebyshev nodes
• Barycentric Lagrange interpolation weights