Example usage and detailed explanation
We will explain the inner workings of SparseIR.jl by means of an example use case, adapted from the sparse-ir paper.
Problem statement
Let us perform self-consistent second-order perturbation theory for the single impurity Anderson model at finite temperature. Its Hamiltonian is given by
\[ H = U c^\dagger_\uparrow c^\dagger_\downarrow c_\downarrow c_\uparrow + \sum_{p\sigma} \big(V_{p\sigma} f_{p\sigma}^\dagger c_\sigma + V_{p\sigma}^* c_\sigma^\dagger c_\sigma^\dagger\big) + \sum_{p\sigma} \epsilon_{p} f_{p\sigma}^\dagger f_{p\sigma}\]
where $U$ is the electron interaction strength, $c_\sigma$ annihilates an electron on the impurity, $f_{p\sigma}$ annihilates an electron in the bath, $\dagger$ denotes the Hermitian conjugate, $p\in\mathbb R$ is bath momentum, and $\sigma\in\{\uparrow, \downarrow\}$ is spin. The hybridization strength $V_{p\sigma}$ and bath energies $\epsilon_p$ are chosen such that the non-interacting density of states is semi-elliptic with a half-bandwidth of one, $\rho_0(\omega) = \frac2\pi\sqrt{1-\omega^2}$, $U=1.2$, $\beta=10$, and the system is assumed to be half-filled.
Treatment
We first import SparseIR and construct an appropriate basis ($\omega_\mathrm{max} = 8$ should be more than enough for this example):
julia> using SparseIR
+Guide · SparseIR.jl Example usage and detailed explanation
We will explain the inner workings of SparseIR.jl by means of an example use case, adapted from the sparse-ir paper.
Problem statement
Let us perform self-consistent second-order perturbation theory for the single impurity Anderson model at finite temperature. Its Hamiltonian is given by
\[ H = U c^\dagger_\uparrow c^\dagger_\downarrow c_\downarrow c_\uparrow + \sum_{p\sigma} \big(V_{p\sigma} f_{p\sigma}^\dagger c_\sigma + V_{p\sigma}^* c_\sigma^\dagger c_\sigma^\dagger\big) + \sum_{p\sigma} \epsilon_{p} f_{p\sigma}^\dagger f_{p\sigma}\]
where $U$ is the electron interaction strength, $c_\sigma$ annihilates an electron on the impurity, $f_{p\sigma}$ annihilates an electron in the bath, $\dagger$ denotes the Hermitian conjugate, $p\in\mathbb R$ is bath momentum, and $\sigma\in\{\uparrow, \downarrow\}$ is spin. The hybridization strength $V_{p\sigma}$ and bath energies $\epsilon_p$ are chosen such that the non-interacting density of states is semi-elliptic with a half-bandwidth of one, $\rho_0(\omega) = \frac2\pi\sqrt{1-\omega^2}$, $U=1.2$, $\beta=10$, and the system is assumed to be half-filled.
Treatment
We first import SparseIR and construct an appropriate basis ($\omega_\mathrm{max} = 8$ should be more than enough for this example):
julia> using SparseIR
julia> basis = FiniteTempBasis(fermion, 10, 8)
FiniteTempBasis{LogisticKernel, Float64}(fermion, 10.0, 8.0)
There's quite a lot happening behind the scenes in this first innocuous-looking statement, so let's break it down: Because we did not specify otherwise, the constructor chose the analytic continuation kernel for fermions, LogisticKernel(80.0), defined by
\[\begin{equation}
@@ -49,4 +49,4 @@
S_\ell \vec v_\ell &\approx K^\mathrm{T} \vec u_\ell.
\end{aligned}
\right\}
-\end{equation}\]
Together with the property $\vec u_\ell^\mathrm{T} \vec u_{\ell'} \approx \delta_{\ell\ell'} \approx \vec v_\ell^\mathrm{T} \vec v_{\ell'}$ we have successfully translated the original SVE problem into an SVD, because
\[ K = \sum_\ell S_\ell \vec u_\ell \vec v_\ell^\mathrm{T}.\]
The next step is calling the matrices function which computes the matrix $K$ derived in the previous step.
Note The function is named in the plural because in the centrosymmetric case it actually returns two matrices $K_+$ and $K_-$, one for the even and one for the odd kernel. These matrices' SVDs are later concatenated, so for simplicity, we will refer to $K$ from here on out.
Info Special care is taken here to avoid FP-arithmetic cancellation around $x = -1$ and $x = +1$.
Note that in the plot, the matrices are rotated 90 degrees to the left to make the connection with the (subregion $[0, 1] × [0, 1]$ of the) previous figure more obvious. Thus we can see how the choice of sampling points has magnified and brought to the matrices' centers the regions of interest. Furthermore, elements with absolute values smaller than $10\%$ of the maximum have been omitted to emphasize the structure; this should however not be taken to mean that there is any sparsity to speak of we could exploit in the next step.
Take the truncated singular value decompostion (TSVD) of $K$, or rather, of $K_+$ and $K_-$. We use here a custom TSVD routine written by Markus Wallerberger which combines a homemade rank-revealing QR decomposition with GenericLinearAlgebra.svd!. This is necessary because there is currently no TSVD for arbitrary types available.
Via the function truncate, we throw away superfluous terms in our expansion. More specifically, we choose $L$ in \eqref{SVE} such that $S_\ell / S_0 > \varepsilon$ for all $\ell \leq L$. Here $\varepsilon$ is our selected precision, in our case it's equal to the double precision machine epsilon, $2^{-52} \approx 2.22 \times 10^{-16}$.
Settings
This document was generated with Documenter.jl version 0.27.25 on Tuesday 12 September 2023. Using Julia version 1.9.3.
Together with the property $\vec u_\ell^\mathrm{T} \vec u_{\ell'} \approx \delta_{\ell\ell'} \approx \vec v_\ell^\mathrm{T} \vec v_{\ell'}$ we have successfully translated the original SVE problem into an SVD, because
\[ K = \sum_\ell S_\ell \vec u_\ell \vec v_\ell^\mathrm{T}.\]
The next step is calling the matrices function which computes the matrix $K$ derived in the previous step.
Note The function is named in the plural because in the centrosymmetric case it actually returns two matrices $K_+$ and $K_-$, one for the even and one for the odd kernel. These matrices' SVDs are later concatenated, so for simplicity, we will refer to $K$ from here on out.
Info Special care is taken here to avoid FP-arithmetic cancellation around $x = -1$ and $x = +1$.
Note that in the plot, the matrices are rotated 90 degrees to the left to make the connection with the (subregion $[0, 1] × [0, 1]$ of the) previous figure more obvious. Thus we can see how the choice of sampling points has magnified and brought to the matrices' centers the regions of interest. Furthermore, elements with absolute values smaller than $10\%$ of the maximum have been omitted to emphasize the structure; this should however not be taken to mean that there is any sparsity to speak of we could exploit in the next step.
Take the truncated singular value decompostion (TSVD) of $K$, or rather, of $K_+$ and $K_-$. We use here a custom TSVD routine written by Markus Wallerberger which combines a homemade rank-revealing QR decomposition with GenericLinearAlgebra.svd!. This is necessary because there is currently no TSVD for arbitrary types available.
Via the function truncate, we throw away superfluous terms in our expansion. More specifically, we choose $L$ in \eqref{SVE} such that $S_\ell / S_0 > \varepsilon$ for all $\ell \leq L$. Here $\varepsilon$ is our selected precision, in our case it's equal to the double precision machine epsilon, $2^{-52} \approx 2.22 \times 10^{-16}$.
![\"Logistic](\"../assets/img/kernel.png\")
","category":"page"},{"location":"guide/#SVE","page":"Guide","title":"SVE","text":"","category":"section"},{"location":"guide/","page":"Guide","title":"Guide","text":"Central is the singular value expansion's (SVE) computation, which is handled by the function SVEResult: Its purpose is constructing the decomposition","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"beginequationlabelSVE\n K(x y) approx sum_ell = 0^L U_ell(x) S_ell V_ell(y)\nendequation","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"where U_ell(x) and V_ell(y) are called K's left and right singular functions respectively and S_ell are its singular values. The singular functions are form an orthonormal basis by construction, i.e.","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"beginequation\n int ddx U_ell(x) U_ell(x) = delta_ellell = int ddy V_ell(y) V_ell(y)\nendequation","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"and thus","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"beginequation labelcoeff1\nleft\nbeginaligned\n S_ell U_ell(x) = int ddy K(x y) V_ell(y) \n S_ell V_ell(y) = int ddx K(x y) U_ell(x)\nendaligned\nright\nendequation","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"Here and in what follows, unless otherwise indicated, integrals are taken to be over the interval -11.","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"The function first calls the choose_accuracy helper and thereby sets the appropriate working precision. Because we did not specify a working accuracy varepsilon, it chooses for us varepsilon approx 22 times 10^-16 and working type Float64x2 - a 128 bits floating point type provided by the MultiFloats.jl package - because in computing the SVD we incur a precision loss of about half our input bits, leaving us with full double accuracy results only if we use quad precision during the computation.\nThen - by calling out to the CentrosymmSVE constructor - a support grid x_i times y_j the kernel will later be evaluated on is built. Along with these support points weights w_i and z_j are computed. These points and weights consist of repeated scaled Gauss integration rules, such that\nbeginequation labelintrules\n int ddx f(x) approx sum_i f(x_i) w_i\n quadtextandquad\n int ddy g(y) approx sum_j g(y_j) z_j\nendequation\nTo get an idea regarding the distribution of these sampling points, refer to following figure, which shows x_i times y_j for Lambda = 80:\n(Image: Sampling point distribution)\nnote: Note\nThe points do not cover -1 1 -1 1 but only 0 1 0 1. This is actually a special case as we exploit the kernel's centrosymmetry, i.e. K(x y) = K(-x -y). It is straightforward to show that the left/right singular vectors then can be chosen as either odd or even functions.Consequentially, they are singular functions of a reduced kernel K^mathrmred_pm on 0 1 0 1 that is given as either:beginequation\n K^mathrmred_pm(x y) = K(x y) pm K(x -y)\nendequationIt is these reduced kernels we will actually sample from, gaining a 4-fold speedup in constructing the SVE. (Image: abc)\nUsing the integration rules \\eqref{intrules} allows us to approximate \\eqref{coeff1} by\nbeginequation labelcoeff2\nleft\nbeginaligned\n S_ell U_ell(x_i) approx sum_j K(x_i y_j) V_ell(y_j) z_j forall i \n S_ell V_ell(y_j) approx sum_i K(x_i y_j) U_ell(x_i) w_i forall j\nendaligned\nright\nendequation\nwhich we now multiply by sqrtw_i and sqrtz_j respectively, yielding\nbeginequation labelcoeff3\nleft\nbeginaligned\n S_ell sqrtw_i U_ell(x_i) approx sum_j sqrtw_i K(x_i y_j) sqrtz_j sqrtz_j V_ell(y_j) \n S_ell sqrtz_j V_ell(y_j) approx sum_i sqrtw_i K(x_i y_j) sqrtz_j sqrtw_i U_ell(x_i)\nendaligned\nright\nendequation\nIf we now define vectors vec u_ell, vec v_ell and a matrix K with entries u_ell i equiv sqrtw_i U_ell(x_i), v_ell j equiv sqrtz_j V_ell(y_j) and K_ij equiv sqrtw_i K(x_i y_j) sqrtz_j, then\nbeginequation labelcoeff4\nleft\nbeginaligned\n S_ell u_ell i approx sum_j K_ij v_ell j \n S_ell v_ell j approx sum_i K_ij u_ell i\nendaligned\nright\nendequation\nor\nbeginequation labelcoeff5\nleft\nbeginaligned\n S_ell vec u_ell approx K^phantommathrmT vec v_ell \n S_ell vec v_ell approx K^mathrmT vec u_ell\nendaligned\nright\nendequation\nTogether with the property vec u_ell^mathrmT vec u_ell approx delta_ellell approx vec v_ell^mathrmT vec v_ell we have successfully translated the original SVE problem into an SVD, because\n K = sum_ell S_ell vec u_ell vec v_ell^mathrmT\nThe next step is calling the matrices function which computes the matrix K derived in the previous step.\nnote: Note\nThe function is named in the plural because in the centrosymmetric case it actually returns two matrices K_+ and K_-, one for the even and one for the odd kernel. These matrices' SVDs are later concatenated, so for simplicity, we will refer to K from here on out.\ninfo: Info\nSpecial care is taken here to avoid FP-arithmetic cancellation around x = -1 and x = +1.\n(Image: Kernel matrices) Note that in the plot, the matrices are rotated 90 degrees to the left to make the connection with the (subregion 0 1 0 1 of the) previous figure more obvious. Thus we can see how the choice of sampling points has magnified and brought to the matrices' centers the regions of interest. Furthermore, elements with absolute values smaller than 10 of the maximum have been omitted to emphasize the structure; this should however not be taken to mean that there is any sparsity to speak of we could exploit in the next step.\nTake the truncated singular value decompostion (TSVD) of K, or rather, of K_+ and K_-. We use here a custom TSVD routine written by Markus Wallerberger which combines a homemade rank-revealing QR decomposition with GenericLinearAlgebra.svd!. This is necessary because there is currently no TSVD for arbitrary types available.\nVia the function truncate, we throw away superfluous terms in our expansion. More specifically, we choose L in \\eqref{SVE} such that S_ell S_0 varepsilon for all ell leq L. Here varepsilon is our selected precision, in our case it's equal to the double precision machine epsilon, 2^-52 approx 222 times 10^-16.","category":"page"},{"location":"private/","page":"Private","title":"Private","text":"CurrentModule = SparseIR","category":"page"},{"location":"private/#Private-names-index","page":"Private","title":"Private names index","text":"","category":"section"},{"location":"private/","page":"Private","title":"Private","text":"Modules = [SparseIR]\nPrivate = true\nPublic = false","category":"page"},{"location":"private/#Core.Int-Tuple{MatsubaraFreq}","page":"Private","title":"Core.Int","text":"Get prefactor n for the Matsubara frequency ω = n*π/β\n\n\n\n\n\n","category":"method"},{"location":"private/#Core.Integer-Tuple{MatsubaraFreq}","page":"Private","title":"Core.Integer","text":"Get prefactor n for the Matsubara frequency ω = n*π/β\n\n\n\n\n\n","category":"method"},{"location":"private/#Core.Union-Union{Tuple{MatsubaraFreq{S}}, Tuple{S}} where S","page":"Private","title":"Core.Union","text":"(polyFT::PiecewiseLegendreFT)(ω)\n\nObtain Fourier transform of polynomial for given MatsubaraFreq ω.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.AbstractAugmentation","page":"Private","title":"SparseIR.AbstractAugmentation","text":"AbstractAugmentation\n\nScalar function in imaginary time/frequency.\n\nThis represents a single function in imaginary time and frequency, together with some auxiliary methods that make it suitable for augmenting a basis.\n\nSee also: AugmentedBasis\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.AbstractBasis","page":"Private","title":"SparseIR.AbstractBasis","text":"AbstractBasis\n\nAbstract base class for bases on the imaginary-time axis.\n\nLet basis be an abstract basis. Then we can expand a two-point propagator G(τ), where τ is imaginary time, into a set of basis functions:\n\nG(τ) == sum(basis.u[l](τ) * g[l] for l in 1:length(basis)) + ϵ(τ),\n\nwhere basis.u[l] is the l-th basis function, g[l] is the associated expansion coefficient and ϵ(τ) is an error term. Similarly, the Fourier transform Ĝ(n), where n is now a Matsubara frequency, can be expanded as follows:\n\nĜ(n) == sum(basis.uhat[l](n) * g[l] for l in 1:length(basis)) + ϵ(n),\n\nwhere basis.uhat[l] is now the Fourier transform of the basis function.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.AbstractKernel","page":"Private","title":"SparseIR.AbstractKernel","text":"(kernel::AbstractKernel)(x, y[, x₊, x₋])\n\nEvaluate kernel at point (x, y).\n\nThe parameters x₊ and x₋, if given, shall contain the values of x - xₘᵢₙ and xₘₐₓ - x, respectively. This is useful if either difference is to be formed and cancellation expected.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.AbstractKernel-2","page":"Private","title":"SparseIR.AbstractKernel","text":"AbstractKernel\n\nIntegral kernel K(x, y).\n\nAbstract base type for an integral kernel, i.e. a AbstractFloat binary function K(x y) used in a Fredhold integral equation of the first kind:\n\n u(x) = K(x y) v(y) dy\n\nwhere x x_mathrmmin x_mathrmmax and y y_mathrmmin y_mathrmmax. For its SVE to exist, the kernel must be square-integrable, for its singular values to decay exponentially, it must be smooth.\n\nIn general, the kernel is applied to a scaled spectral function ρ(y) as:\n\n K(x y) ρ(y) dy\n\nwhere ρ(y) = w(y) ρ(y).\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.AbstractSVEHints","page":"Private","title":"SparseIR.AbstractSVEHints","text":"AbstractSVEHints\n\nDiscretization hints for singular value expansion of a given kernel.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.AbstractSampling","page":"Private","title":"SparseIR.AbstractSampling","text":"AbstractSampling\n\nAbstract type for sparse sampling.\n\nEncodes the \"basis transformation\" of a propagator from the truncated IR basis coefficients G_ir[l] to time/frequency sampled on sparse points G(x[i]) together with its inverse, a least squares fit:\n\n ________________ ___________________\n | | evaluate | |\n | Basis |---------------->| Value on |\n | coefficients |<----------------| sampling points |\n |________________| fit |___________________|\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.CentrosymmSVE","page":"Private","title":"SparseIR.CentrosymmSVE","text":"CentrosymmSVE <: AbstractSVE\n\nSVE of centrosymmetric kernel in block-diagonal (even/odd) basis.\n\nFor a centrosymmetric kernel K, i.e., a kernel satisfying: K(x, y) == K(-x, -y), one can make the following ansatz for the singular functions:\n\nu[l](x) = ured[l](x) + sign[l] * ured[l](-x)\nv[l](y) = vred[l](y) + sign[l] * ured[l](-y)\n\nwhere sign[l] is either +1 or -1. This means that the singular value expansion can be block-diagonalized into an even and an odd part by (anti-)symmetrizing the kernel:\n\nK_even = K(x, y) + K(x, -y)\nK_odd = K(x, y) - K(x, -y)\n\nThe lth basis function, restricted to the positive interval, is then the singular function of one of these kernels. If the kernel generates a Chebyshev system [1], then even and odd basis functions alternate.\n\n[1]: A. Karlin, Total Positivity (1968).\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.LogisticKernelOdd","page":"Private","title":"SparseIR.LogisticKernelOdd","text":"LogisticKernelOdd <: AbstractReducedKernel\n\nFermionic analytical continuation kernel, odd.\n\nIn dimensionless variables x = 2τβ - 1, y = βωΛ, the fermionic integral kernel is a function on -1 1 -1 1:\n\n K(x y) = -fracsinh(Λ x y 2)cosh(Λ y 2)\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.PiecewiseLegendreFT","page":"Private","title":"SparseIR.PiecewiseLegendreFT","text":"PiecewiseLegendreFT <: Function\n\nFourier transform of a piecewise Legendre polynomial.\n\nFor a given frequency index n, the Fourier transform of the Legendre function is defined as:\n\n p̂(n) == ∫ dx exp(im * π * n * x / (xmax - xmin)) p(x)\n\nThe polynomial is continued either periodically (freq=:even), in which case n must be even, or antiperiodically (freq=:odd), in which case n must be odd.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.PiecewiseLegendrePoly","page":"Private","title":"SparseIR.PiecewiseLegendrePoly","text":"PiecewiseLegendrePoly <: Function\n\nPiecewise Legendre polynomial.\n\nModels a function on the interval xmin xmax as a set of segments on the intervals Si = ai ai+1, where on each interval the function is expanded in scaled Legendre polynomials.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.PiecewiseLegendrePolyVector","page":"Private","title":"SparseIR.PiecewiseLegendrePolyVector","text":"PiecewiseLegendrePolyVector\n\nAlias for Vector{PiecewiseLegendrePoly}.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.PowerModel","page":"Private","title":"SparseIR.PowerModel","text":"PowerModel\n\nModel from a high-frequency series expansion::\n\nA(iω) == sum(A[n] / (iω)^(n+1) for n in 1:N)\n\nwhere iω == i * π2 * wn is a reduced imaginary frequency, i.e., wn is an odd/even number for fermionic/bosonic frequencies.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.ReducedKernel","page":"Private","title":"SparseIR.ReducedKernel","text":"ReducedKernel\n\nRestriction of centrosymmetric kernel to positive interval.\n\nFor a kernel K on -1 1 -1 1 that is centrosymmetric, i.e. K(x y) = K(-x -y), it is straight-forward to show that the left/right singular vectors can be chosen as either odd or even functions.\n\nConsequentially, they are singular functions of a reduced kernel K_mathrmred on 0 1 0 1 that is given as either:\n\n K_mathrmred(x y) = K(x y) pm K(x -y)\n\nThis kernel is what this type represents. The full singular functions can be reconstructed by (anti-)symmetrically continuing them to the negative axis.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.RegularizedBoseKernelOdd","page":"Private","title":"SparseIR.RegularizedBoseKernelOdd","text":"RegularizedBoseKernelOdd <: AbstractReducedKernel\n\nBosonic analytical continuation kernel, odd.\n\nIn dimensionless variables x = 2 τ β - 1, y = β ω Λ, the fermionic integral kernel is a function on -1 1 -1 1:\n\n K(x y) = -y fracsinh(Λ x y 2)sinh(Λ y 2)\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.Rule","page":"Private","title":"SparseIR.Rule","text":"Rule{T<:AbstractFloat}\n\nQuadrature rule.\n\nApproximation of an integral over [a, b] by a sum over discrete points x with weights w:\n\n f(x) ω(x) dx _i f(x_i) w_i\n\nwhere we generally have superexponential convergence for smooth f(x) in the number of quadrature points.\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.SVEResult-Tuple{SparseIR.AbstractKernel}","page":"Private","title":"SparseIR.SVEResult","text":"SVEResult(kernel::AbstractKernel;\n Twork=nothing, ε=nothing, lmax=typemax(Int),\n n_gauss=nothing, svd_strat=:auto,\n sve_strat=iscentrosymmetric(kernel) ? CentrosymmSVE : SamplingSVE\n)\n\nPerform truncated singular value expansion of a kernel.\n\nPerform a truncated singular value expansion (SVE) of an integral kernel kernel : [xmin, xmax] x [ymin, ymax] -> ℝ:\n\nkernel(x, y) == sum(s[l] * u[l](x) * v[l](y) for l in (1, 2, 3, ...)),\n\nwhere s[l] are the singular values, which are ordered in non-increasing fashion, u[l](x) are the left singular functions, which form an orthonormal system on [xmin, xmax], and v[l](y) are the right singular functions, which form an orthonormal system on [ymin, ymax].\n\nThe SVE is mapped onto the singular value decomposition (SVD) of a matrix by expanding the kernel in piecewise Legendre polynomials (by default by using a collocation).\n\nArguments\n\nK::AbstractKernel: Integral kernel to take SVE from.\nε::Real: Accuracy target for the basis: attempt to have singular values down to a relative magnitude of ε, and have each singular value and singular vector be accurate to ε. A Twork with a machine epsilon of ε^2 or lower is required to satisfy this. Defaults to 2.2e-16 if xprec is available, and 1.5e-8 otherwise.\ncutoff::Real: Relative cutoff for the singular values. A Twork with machine epsilon of cutoff is required to satisfy this. Defaults to a small multiple of the machine epsilon.\nNote that cutoff and ε serve distinct purposes. cutoff reprsents the accuracy to which the kernel is reproduced, whereas ε is the accuracy to which the singular values and vectors are guaranteed.\nlmax::Integer: Maximum basis size. If given, only at most the lmax most significant singular values and associated singular functions are returned.\n`n_gauss (int): Order of Legendre polynomials. Defaults to kernel hinted value.\nTwork: Working data type. Defaults to a data type with machine epsilon of at mostε^2and at mostcutoff`, or otherwise most accurate data type available.\nsve_strat::AbstractSVE: SVE to SVD translation strategy. Defaults to SamplingSVE, optionally wrapped inside of a CentrosymmSVE if the kernel is centrosymmetric.\nsvd_strat ('fast' or 'default' or 'accurate'): SVD solver. Defaults to fast (ID/RRQR) based solution when accuracy goals are moderate, and more accurate Jacobi-based algorithm otherwise.\n\nReturns: An SVEResult containing the truncated singular value expansion.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.SamplingSVE","page":"Private","title":"SparseIR.SamplingSVE","text":"SamplingSVE <: AbstractSVE\n\nSVE to SVD translation by sampling technique [1].\n\nMaps the singular value expansion (SVE) of a kernel kernel onto the singular value decomposition of a matrix A. This is achieved by choosing two sets of Gauss quadrature rules: (x, wx) and (y, wy) and approximating the integrals in the SVE equations by finite sums. This implies that the singular values of the SVE are well-approximated by the singular values of the following matrix:\n\nA[i, j] = √(wx[i]) * K(x[i], y[j]) * √(wy[j])\n\nand the values of the singular functions at the Gauss sampling points can be reconstructed from the singular vectors u and v as follows:\n\nu[l,i] ≈ √(wx[i]) u[l](x[i])\nv[l,j] ≈ √(wy[j]) u[l](y[j])\n\n[1] P. Hansen, Discrete Inverse Problems, Ch. 3.1\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.Statistics","page":"Private","title":"SparseIR.Statistics","text":"Statistics(zeta)\n\nAbstract type for quantum statistics (fermionic/bosonic/etc.)\n\n\n\n\n\n","category":"type"},{"location":"private/#SparseIR.accuracy","page":"Private","title":"SparseIR.accuracy","text":"accuracy(basis::AbstractBasis)\n\nAccuracy of the basis.\n\nUpper bound to the relative error of reprensenting a propagator with the given number of basis functions (number between 0 and 1).\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.canonicalize!-Tuple{Any, Any}","page":"Private","title":"SparseIR.canonicalize!","text":"canonicalize!(u, v)\n\nCanonicalize basis.\n\nEach SVD (u[l], v[l]) pair is unique only up to a global phase, which may differ from implementation to implementation and also platform. We fix that gauge by demanding u[l](1) > 0. This ensures a diffeomorphic connection to the Legendre polynomials as Λ → 0.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.choose_accuracy-Tuple{Any, Any, Any}","page":"Private","title":"SparseIR.choose_accuracy","text":"choose_accuracy(ε, Twork[, svd_strat])\n\nChoose work type and accuracy based on specs and defaults\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.compute_unl_inner-Tuple{SparseIR.PiecewiseLegendrePoly, Any}","page":"Private","title":"SparseIR.compute_unl_inner","text":"compute_unl_inner(poly, wn)\n\nCompute piecewise Legendre to Matsubara transform.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.conv_radius","page":"Private","title":"SparseIR.conv_radius","text":"conv_radius(kernel)\n\nConvergence radius of the Matsubara basis asymptotic model.\n\nFor improved relative numerical accuracy, the IR basis functions on the Matsubara axis uhat(basis, n) can be evaluated from an asymptotic expression for abs(n) > conv_radius. If isinf(conv_radius), then the asymptotics are unused (the default).\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.default_matsubara_sampling_points","page":"Private","title":"SparseIR.default_matsubara_sampling_points","text":"default_matsubara_sampling_points(basis::AbstractBasis; positive_only=false)\n\nDefault sampling points on the imaginary frequency axis.\n\nArguments\n\npositive_only::Bool: Only return non-negative frequencies. This is useful if the object to be fitted is symmetric in Matsubura frequency, ĝ(ω) == conj(ĝ(-ω)), or, equivalently, real in imaginary time.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.default_tau_sampling_points","page":"Private","title":"SparseIR.default_tau_sampling_points","text":"default_tau_sampling_points(basis::AbstractBasis)\n\nDefault sampling points on the imaginary time/x axis.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.deriv-Union{Tuple{SparseIR.PiecewiseLegendrePoly}, Tuple{n}, Tuple{SparseIR.PiecewiseLegendrePoly, Val{n}}} where n","page":"Private","title":"SparseIR.deriv","text":"deriv(poly[, ::Val{n}=Val(1)])\n\nGet polynomial for the nth derivative.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.eval_matrix","page":"Private","title":"SparseIR.eval_matrix","text":"eval_matrix(T, basis, x)\n\nReturn evaluation matrix from coefficients to sampling points. T <: AbstractSampling.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.find_extrema-Tuple{SparseIR.PiecewiseLegendreFT}","page":"Private","title":"SparseIR.find_extrema","text":"find_extrema(polyFT::PiecewiseLegendreFT; part=nothing, grid=DEFAULT_GRID)\n\nObtain extrema of Fourier-transformed polynomial.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.finite_temp_bases","page":"Private","title":"SparseIR.finite_temp_bases","text":"finite_temp_bases(β, ωmax, ε, sve_result=SVEResult(LogisticKernel(β * ωmax); ε))\n\nConstruct FiniteTempBasis objects for fermion and bosons using the same LogisticKernel instance.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.from_IR","page":"Private","title":"SparseIR.from_IR","text":"from_IR(dlr::DiscreteLehmannRepresentation, gl::AbstractArray, dims=1)\n\nFrom IR to DLR. gl`: Expansion coefficients in IR.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.get_symmetrized-Tuple{SparseIR.AbstractKernel, Any}","page":"Private","title":"SparseIR.get_symmetrized","text":"get_symmetrized(kernel, sign)\n\nConstruct a symmetrized version of kernel, i.e. kernel(x, y) + sign * kernel(x, -y).\n\nwarning: Beware!\nBy default, this returns a simple wrapper over the current instance which naively performs the sum. You may want to override this to avoid cancellation.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.get_tnl-Tuple{Any, Any}","page":"Private","title":"SparseIR.get_tnl","text":"get_tnl(l, w)\n\nFourier integral of the l-th Legendre polynomial::\n\nTₗ(ω) == ∫ dx exp(iωx) Pₗ(x)\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.giw-Tuple{Any, Integer}","page":"Private","title":"SparseIR.giw","text":"giw(polyFT, wn)\n\nReturn model Green's function for reduced frequencies\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.iscentrosymmetric","page":"Private","title":"SparseIR.iscentrosymmetric","text":"iscentrosymmetric(kernel)\n\nReturn true if kernel(x, y) == kernel(-x, -y) for all values of x and y in range. This allows the kernel to be block-diagonalized, speeding up the singular value expansion by a factor of 4. Defaults to false.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.iswellconditioned-Tuple{SparseIR.AbstractBasis}","page":"Private","title":"SparseIR.iswellconditioned","text":"iswellconditioned(basis::AbstractBasis)\n\nReturns true if the sampling is expected to be well-conditioned.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.joinrules-Union{Tuple{AbstractArray{SparseIR.Rule{T}, 1}}, Tuple{T}} where T","page":"Private","title":"SparseIR.joinrules","text":"joinrules(rules)\n\nJoin multiple Gauss quadratures together.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.legder-Union{Tuple{AbstractMatrix{T}}, Tuple{T}, Tuple{AbstractMatrix{T}, Any}} where T","page":"Private","title":"SparseIR.legder","text":"legder\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.legendre-Union{Tuple{Any}, Tuple{T}, Tuple{Any, Type{T}}} where T","page":"Private","title":"SparseIR.legendre","text":"legendre(n[, T])\n\nGauss-Legendre quadrature with n points on [-1, 1].\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.legendre_collocation","page":"Private","title":"SparseIR.legendre_collocation","text":"legendre_collocation(rule, n=length(rule.x))\n\nGenerate collocation matrix from Gauss-Legendre rule.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.legvander-Union{Tuple{T}, Tuple{AbstractVector{T}, Integer}} where T","page":"Private","title":"SparseIR.legvander","text":"legvander(x, deg)\n\nPseudo-Vandermonde matrix of degree deg.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.matop!-Union{Tuple{N}, Tuple{T}, Tuple{S}, Tuple{AbstractArray{S, N}, Any, AbstractArray{T, N}, Any, Any}} where {S, T, N}","page":"Private","title":"SparseIR.matop!","text":"matop!(buffer, mat, arr::AbstractArray, op, dim)\n\nApply the operator op to the matrix mat and to the array arr along the first dimension (dim=1) or the last dimension (dim=N).\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.matop_along_dim!-Union{Tuple{N}, Tuple{T}, Tuple{Any, Any, AbstractArray{T, N}, Any, Any}} where {T, N}","page":"Private","title":"SparseIR.matop_along_dim!","text":"matop_along_dim!(buffer, mat, arr::AbstractArray, dim::Integer, op)\n\nApply the operator op to the matrix mat and to the array arr along the dimension dim, writing the result to buffer.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.matrices-Tuple{SparseIR.SamplingSVE}","page":"Private","title":"SparseIR.matrices","text":"matrices(sve::AbstractSVE)\n\nSVD problems underlying the SVE.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.matrix_from_gauss-Union{Tuple{T}, Tuple{Any, SparseIR.Rule{T}, SparseIR.Rule{T}}} where T","page":"Private","title":"SparseIR.matrix_from_gauss","text":"matrix_from_gauss(kernel, gauss_x, gauss_y)\n\nCompute matrix for kernel from Gauss rules.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.movedim-Union{Tuple{N}, Tuple{T}, Tuple{AbstractArray{T, N}, Pair}} where {T, N}","page":"Private","title":"SparseIR.movedim","text":"movedim(arr::AbstractArray, src => dst)\n\nMove arr's dimension at src to dst while keeping the order of the remaining dimensions unchanged.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.ngauss","page":"Private","title":"SparseIR.ngauss","text":"ngauss(hints)\n\nGauss-Legendre order to use to guarantee accuracy.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.nsvals-Tuple{SparseIR.SVEHintsLogistic}","page":"Private","title":"SparseIR.nsvals","text":"nsvals(hints)\n\nUpper bound for number of singular values.\n\nUpper bound on the number of singular values above the given threshold, i.e. where s[l] ≥ ε * first(s).\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.phase_stable-Tuple{Any, Integer}","page":"Private","title":"SparseIR.phase_stable","text":"phase_stable(poly, wn)\n\nPhase factor for the piecewise Legendre to Matsubara transform.\n\nCompute the following phase factor in a stable way:\n\nexp.(iπ/2 * wn * cumsum(poly.Δx))\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.piecewise-Tuple{Any, Vector}","page":"Private","title":"SparseIR.piecewise","text":"piecewise(rule, edges)\n\nPiecewise quadrature with the same quadrature rule, but scaled.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.postprocess-Tuple{SparseIR.SamplingSVE, Any, Any, Any}","page":"Private","title":"SparseIR.postprocess","text":"postprocess(sve::AbstractSVE, u, s, v)\n\nConstruct the SVE result from the SVD.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.rescale-Tuple{FiniteTempBasis, Any}","page":"Private","title":"SparseIR.rescale","text":"rescale(basis::FiniteTempBasis, new_β)\n\nReturn a basis for different temperature.\n\nUses the same kernel with the same ε, but a different temperature. Note that this implies a different UV cutoff ωmax, since Λ == β * ωmax stays constant.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.reseat-Tuple{SparseIR.Rule, Any, Any}","page":"Private","title":"SparseIR.reseat","text":"reseat(rule, a, b)\n\nReseat quadrature rule to new domain.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.roots-Tuple{SparseIR.PiecewiseLegendrePoly}","page":"Private","title":"SparseIR.roots","text":"roots(poly)\n\nFind all roots of the piecewise polynomial poly.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.scale-Tuple{Any, Any}","page":"Private","title":"SparseIR.scale","text":"scale(rule, factor)\n\nScale weights by factor.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.segments_x-Union{Tuple{SparseIR.SVEHintsLogistic}, Tuple{T}, Tuple{SparseIR.SVEHintsLogistic, Type{T}}} where T","page":"Private","title":"SparseIR.segments_x","text":"segments_x(sve_hints::AbstractSVEHints[, T])\n\nSegments for piecewise polynomials on the x axis.\n\nList of segments on the x axis for the associated piecewise polynomial. Should reflect the approximate position of roots of a high-order singular function in x.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.segments_y-Union{Tuple{SparseIR.SVEHintsLogistic}, Tuple{T}, Tuple{SparseIR.SVEHintsLogistic, Type{T}}} where T","page":"Private","title":"SparseIR.segments_y","text":"segments_y(sve_hints::AbstractSVEHints[, T])\n\nSegments for piecewise polynomials on the y axis.\n\nList of segments on the y axis for the associated piecewise polynomial. Should reflect the approximate position of roots of a high-order singular function in y.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.shift_xmid-Tuple{Any, Any}","page":"Private","title":"SparseIR.shift_xmid","text":"shift_xmid(knots, Δx)\n\nReturn midpoint relative to the nearest integer plus a shift.\n\nReturn the midpoints xmid of the segments, as pair (diff, shift), where shift is in (0, 1, -1) and diff is a float such that xmid == shift + diff to floating point accuracy.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.significance","page":"Private","title":"SparseIR.significance","text":"significance(basis::AbstractBasis)\n\nReturn vector σ, where 0 ≤ σ[i] ≤ 1 is the significance level of the i-th basis function. If ϵ is the desired accuracy to which to represent a propagator, then any basis function where σ[i] < ϵ can be neglected.\n\nFor the IR basis, we simply have that σ[i] = s[i] / first(s).\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.split-Tuple{Any, Real}","page":"Private","title":"SparseIR.split","text":"split(poly, x)\n\nSplit segment.\n\nFind segment of poly's domain that covers x.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.statistics-Union{Tuple{SparseIR.AbstractBasis{S}}, Tuple{S}} where S<:SparseIR.Statistics","page":"Private","title":"SparseIR.statistics","text":"statistics(basis::AbstractBasis)\n\nQuantum statistic (Statistics instance, Fermionic() or Bosonic()).\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.sve_hints","page":"Private","title":"SparseIR.sve_hints","text":"sve_hints(kernel, ε)\n\nProvide discretisation hints for the SVE routines.\n\nAdvises the SVE routines of discretisation parameters suitable in tranforming the (infinite) SVE into an (finite) SVD problem.\n\nSee also AbstractSVEHints.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.to_IR","page":"Private","title":"SparseIR.to_IR","text":"to_IR(dlr::DiscreteLehmannRepresentation, g_dlr::AbstractArray, dims=1)\n\nFrom DLR to IR. g_dlr`: Expansion coefficients in DLR.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.truncate-Tuple{Any, Any, Any}","page":"Private","title":"SparseIR.truncate","text":"truncate(u, s, v; rtol=0.0, lmax=typemax(Int))\n\nTruncate singular value expansion.\n\nArguments\n\n- `u`, `s`, `v`: Thin singular value expansion\n- `rtol`: Only singular values satisfying `s[l]/s[1] > rtol` are retained.\n- `lmax`: At most the `lmax` most significant singular values are retained.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.value-Tuple{MatsubaraFreq, Real}","page":"Private","title":"SparseIR.value","text":"Get value of the Matsubara frequency ω = n*π/β\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.valueim-Tuple{MatsubaraFreq, Real}","page":"Private","title":"SparseIR.valueim","text":"Get complex value of the Matsubara frequency iω = iπ/β * n\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.weight_func","page":"Private","title":"SparseIR.weight_func","text":"weight_func(kernel, statistics::Statistics)\n\nReturn the weight function for the given statistics.\n\nFermion: w(x) == 1\nBoson: w(y) == 1/tanh(Λ*y/2)\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.workarrlength-Tuple{SparseIR.AbstractSampling, AbstractArray}","page":"Private","title":"SparseIR.workarrlength","text":"workarrlength(smpl::AbstractSampling, al; dim=1)\n\nReturn length of workarr for fit!.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.xrange","page":"Private","title":"SparseIR.xrange","text":"xrange(kernel)\n\nReturn a tuple (x_mathrmmin x_mathrmmax) delimiting the range of allowed x values.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.ypower","page":"Private","title":"SparseIR.ypower","text":"ypower(kernel)\n\nPower with which the y coordinate scales.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.yrange","page":"Private","title":"SparseIR.yrange","text":"yrange(kernel)\n\nReturn a tuple (y_mathrmmin y_mathrmmax) delimiting the range of allowed y values.\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.zeta-Tuple{MatsubaraFreq}","page":"Private","title":"SparseIR.zeta","text":"Get statistics ζ for Matsubara frequency ω = (2*m+ζ)*π/β\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.Λ","page":"Private","title":"SparseIR.Λ","text":"Λ(basis::AbstractBasis)\nlambda(basis::AbstractBasis)\n\nBasis cutoff parameter, Λ = β * ωmax, or None if not present\n\n\n\n\n\n","category":"function"},{"location":"private/#SparseIR.β-Tuple{SparseIR.AbstractBasis}","page":"Private","title":"SparseIR.β","text":"β(basis::AbstractBasis)\nbeta(basis::AbstractBasis)\n\nInverse temperature or nothing if unscaled basis.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR.ωmax","page":"Private","title":"SparseIR.ωmax","text":"ωmax(basis::AbstractBasis)\nwmax(basis::AbstractBasis)\n\nReal frequency cutoff or nothing if unscaled basis.\n\n\n\n\n\n","category":"function"},{"location":"private/","page":"Private","title":"Private","text":"Modules = [SparseIR._LinAlg]\nPrivate = true\nPublic = true","category":"page"},{"location":"private/#SparseIR._LinAlg.givens_lmul-Union{Tuple{T}, Tuple{Tuple{T, T}, Any}} where T","page":"Private","title":"SparseIR._LinAlg.givens_lmul","text":"Apply Givens rotation to vector:\n\n [ a ] = [ c s ] [ x ]\n [ b ] [ -s c ] [ y ]\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.givens_params-Union{Tuple{T}, Tuple{T, T}} where T<:AbstractFloat","page":"Private","title":"SparseIR._LinAlg.givens_params","text":"Compute Givens rotation R matrix that satisfies:\n\n[ c s ] [ f ] [ r ]\n[ -s c ] [ g ] = [ 0 ]\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.rrqr!-Union{Tuple{AbstractMatrix{T}}, Tuple{T}} where T<:AbstractFloat","page":"Private","title":"SparseIR._LinAlg.rrqr!","text":"Truncated rank-revealing QR decomposition with full column pivoting.\n\nDecomposes a (m, n) matrix A into the product:\n\nA[:,piv] == Q * R\n\nwhere Q is an (m, k) isometric matrix, R is a (k, n) upper triangular matrix, piv is a permutation vector, and k is chosen such that the relative tolerance tol is met in the equality above.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.rrqr-Union{Tuple{AbstractMatrix{T}}, Tuple{T}} where T<:AbstractFloat","page":"Private","title":"SparseIR._LinAlg.rrqr","text":"Truncated rank-revealing QR decomposition with full column pivoting.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.svd2x2-Union{Tuple{T}, NTuple{4, T}} where T","page":"Private","title":"SparseIR._LinAlg.svd2x2","text":"Perform the SVD of an arbitrary two-by-two matrix:\n\n [ a11 a12 ] = [ cu -su ] [ smax 0 ] [ cv sv ]\n [ a21 a22 ] [ su cu ] [ 0 smin ] [ -sv cv ]\n\nNote that smax and smin can be negative.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.svd2x2-Union{Tuple{T}, Tuple{T, T, T}} where T<:AbstractFloat","page":"Private","title":"SparseIR._LinAlg.svd2x2","text":"Perform the SVD of upper triangular two-by-two matrix:\n\n [ f g ] = [ cu -su ] [ smax 0 ] [ cv sv ]\n [ 0 h ] [ su cu ] [ 0 smin ] [ -sv cv ]\n\nNote that smax and smin can be negative.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.svd_jacobi!-Union{Tuple{AbstractMatrix{T}}, Tuple{T}} where T","page":"Private","title":"SparseIR._LinAlg.svd_jacobi!","text":"Singular value decomposition using Jacobi rotations.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.svd_jacobi-Union{Tuple{AbstractMatrix{T}}, Tuple{T}} where T","page":"Private","title":"SparseIR._LinAlg.svd_jacobi","text":"Singular value decomposition using Jacobi rotations.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.truncate_qr_result-Union{Tuple{T}, Tuple{LinearAlgebra.QRPivoted{T, S, C} where {S<:AbstractMatrix{T}, C<:AbstractVector{T}}, Integer}} where T","page":"Private","title":"SparseIR._LinAlg.truncate_qr_result","text":"Truncate RRQR result low-rank\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.tsvd!-Union{Tuple{AbstractMatrix{T}}, Tuple{T}} where T<:AbstractFloat","page":"Private","title":"SparseIR._LinAlg.tsvd!","text":"Truncated singular value decomposition.\n\nDecomposes an (m, n) matrix A into the product:\n\nA == U * (s .* VT)\n\nwhere U is a (m, k) matrix with orthogonal columns, VT is a (k, n) matrix with orthogonal rows and s are the singular values, a set of k nonnegative numbers in non-ascending order. The SVD is truncated in the sense that singular values below tol are discarded.\n\n\n\n\n\n","category":"method"},{"location":"private/#SparseIR._LinAlg.tsvd-Union{Tuple{AbstractMatrix{T}}, Tuple{T}} where T<:AbstractFloat","page":"Private","title":"SparseIR._LinAlg.tsvd","text":"Truncated singular value decomposition.\n\n\n\n\n\n","category":"method"},{"location":"public/","page":"Public","title":"Public","text":"CurrentModule = SparseIR","category":"page"},{"location":"public/#Public-names-index","page":"Public","title":"Public names index","text":"","category":"section"},{"location":"public/","page":"Public","title":"Public","text":"Modules = [SparseIR]\nPrivate = false\nPublic = true","category":"page"},{"location":"public/#SparseIR.SparseIR","page":"Public","title":"SparseIR.SparseIR","text":"Intermediate representation (IR) for many-body propagators.\n\n\n\n\n\n","category":"module"},{"location":"public/#SparseIR.AugmentedBasis","page":"Public","title":"SparseIR.AugmentedBasis","text":"AugmentedBasis <: AbstractBasis\n\nAugmented basis on the imaginary-time/frequency axis.\n\nGroups a set of additional functions, augmentations, with a given basis. The augmented functions then form the first basis functions, while the rest is provided by the regular basis, i.e.:\n\nu[l](x) == l < naug ? augmentations[l](x) : basis.u[l-naug](x),\n\nwhere naug = length(augmentations) is the number of added basis functions through augmentation. Similar expressions hold for Matsubara frequencies.\n\nAugmentation is useful in constructing bases for vertex-like quantities such as self-energies [wallerberger2021] and when constructing a two-point kernel that serves as a base for multi-point functions [shinaoka2018].\n\nwarning: Warning\nBases augmented with TauConst and TauLinear tend to be poorly conditioned. Care must be taken while fitting and compactness should be enforced if possible to regularize the problem.While vertex bases, i.e. bases augmented with MatsubaraConst, stay reasonably well-conditioned, it is still good practice to treat the Hartree–Fock term separately rather than including it in the basis, if possible.\n\nSee also: MatsubaraConst for vertex basis [wallerberger2021], TauConst, TauLinear for multi-point [shinaoka2018]\n\n[wallerberger2021]: https://doi.org/10.1103/PhysRevResearch.3.033168\n\n[shinaoka2018]: https://doi.org/10.1103/PhysRevB.97.205111\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.Bosonic","page":"Public","title":"SparseIR.Bosonic","text":"Bosonic statistics.\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.DiscreteLehmannRepresentation","page":"Public","title":"SparseIR.DiscreteLehmannRepresentation","text":"DiscreteLehmannRepresentation <: AbstractBasis\n\nDiscrete Lehmann representation (DLR) with poles selected according to extrema of IR.\n\nThis class implements a variant of the discrete Lehmann representation (DLR) 1. Instead of a truncated singular value expansion of the analytic continuation kernel K like the IR, the discrete Lehmann representation is based on a \"sketching\" of K. The resulting basis is a linear combination of discrete set of poles on the real-frequency axis, continued to the imaginary-frequency axis:\n\n G(iv) == sum(a[i] / (iv - w[i]) for i in range(L))\n\nWarning The poles on the real-frequency axis selected for the DLR are based on a rank-revealing decomposition, which offers accuracy guarantees. Here, we instead select the pole locations based on the zeros of the IR basis functions on the real axis, which is a heuristic. We do not expect that difference to matter, but please don't blame the DLR authors if we were wrong :-)\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.Fermionic","page":"Public","title":"SparseIR.Fermionic","text":"Fermionic statistics.\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.FiniteTempBasis","page":"Public","title":"SparseIR.FiniteTempBasis","text":"FiniteTempBasis <: AbstractBasis\n\nIntermediate representation (IR) basis for given temperature.\n\nFor a continuation kernel K from real frequencies, ω ∈ [-ωmax, ωmax], to imaginary time, τ ∈ [0, β], this type stores the truncated singular value expansion or IR basis:\n\nK(τ, ω) ≈ sum(u[l](τ) * s[l] * v[l](ω) for l in 1:L)\n\nThis basis is inferred from a reduced form by appropriate scaling of the variables.\n\nFields\n\nu::PiecewiseLegendrePolyVector: Set of IR basis functions on the imaginary time (tau) axis. These functions are stored as piecewise Legendre polynomials.\nTo obtain the value of all basis functions at a point or a array of points x, you can call the function u(x). To obtain a single basis function, a slice or a subset l, you can use u[l].\nuhat::PiecewiseLegendreFT: Set of IR basis functions on the Matsubara frequency (wn) axis. These objects are stored as a set of Bessel functions.\nTo obtain the value of all basis functions at a Matsubara frequency or a array of points wn, you can call the function uhat(wn). Note that we expect reduced frequencies, which are simply even/odd numbers for bosonic/fermionic objects. To obtain a single basis function, a slice or a subset l, you can use uhat[l].\ns: Vector of singular values of the continuation kernel\nv::PiecewiseLegendrePoly: Set of IR basis functions on the real frequency (w) axis. These functions are stored as piecewise Legendre polynomials.\nTo obtain the value of all basis functions at a point or a array of points w, you can call the function v(w). To obtain a single basis function, a slice or a subset l, you can use v[l].\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.FiniteTempBasis-Union{Tuple{S}, Tuple{Real, Real}, Tuple{Real, Real, Any}} where S","page":"Public","title":"SparseIR.FiniteTempBasis","text":"FiniteTempBasis{S}(β, ωmax, ε=nothing; max_size=nothing, args...)\n\nConstruct a finite temperature basis suitable for the given S (Fermionic or Bosonic) and cutoffs β and ωmax.\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.FiniteTempBasisSet","page":"Public","title":"SparseIR.FiniteTempBasisSet","text":"FiniteTempBasisSet\n\nType for holding IR bases and sparse-sampling objects.\n\nAn object of this type holds IR bases for fermions and bosons and associated sparse-sampling objects.\n\nFields\n\nbasis_f::FiniteTempBasis: Fermion basis\nbasis_b::FiniteTempBasis: Boson basis\ntau::Vector{Float64}: Sampling points in the imaginary-time domain\nwn_f::Vector{Int}: Sampling fermionic frequencies\nwn_b::Vector{Int}: Sampling bosonic frequencies\nsmpltauf::TauSampling: Sparse sampling for tau & fermion\nsmpltaub::TauSampling: Sparse sampling for tau & boson\nsmplwnf::MatsubaraSampling: Sparse sampling for Matsubara frequency & fermion\nsmplwnb::MatsubaraSampling: Sparse sampling for Matsubara frequency & boson\nsve_result::Tuple{PiecewiseLegendrePoly,Vector{Float64},PiecewiseLegendrePoly}: Results of SVE\n\nGetters\n\nbeta::Float64: Inverse temperature\nωmax::Float64: Cut-off frequency\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.LogisticKernel","page":"Public","title":"SparseIR.LogisticKernel","text":"LogisticKernel <: AbstractKernel\n\nFermionic/bosonic analytical continuation kernel.\n\nIn dimensionless variables x = 2 τβ - 1, y = β ωΛ, the integral kernel is a function on -1 1 -1 1:\n\n K(x y) = frace^-Λ y (x + 1) 21 + e^-Λ y\n\nLogisticKernel is a fermionic analytic continuation kernel. Nevertheless, one can model the τ dependence of a bosonic correlation function as follows:\n\n frace^-Λ y (x + 1) 21 - e^-Λ y ρ(y) dy = K(x y) ρ(y) dy\n\nwith\n\n ρ(y) = w(y) ρ(y)\n\nwhere the weight function is given by\n\n w(y) = frac1tanh(Λ y2)\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.MatsubaraFreq","page":"Public","title":"SparseIR.MatsubaraFreq","text":"MatsubaraFreq(n)\n\nPrefactor n of the Matsubara frequency ω = n*π/β\n\nStruct representing the Matsubara frequency ω entering the Fourier transform of a propagator G(τ) on imaginary time τ to its Matsubara equivalent Ĝ(iω) on the imaginary-frequency axis:\n\n β\nĜ(iω) = ∫ dτ exp(iωτ) G(τ) with ω = n π/β,\n 0\n\nwhere β is inverse temperature and by convention we include the imaginary unit in the frequency argument, i.e, Ĝ(iω). The frequencies depend on the statistics of the propagator, i.e., we have that:\n\nG(τ - β) = ± G(τ)\n\nwhere + is for bosons and - is for fermions. The frequencies are restricted accordingly.\n\nBosonic frequency (S == Fermionic): n even (periodic in β)\nFermionic frequency (S == Bosonic): n odd (anti-periodic in β)\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.MatsubaraSampling","page":"Public","title":"SparseIR.MatsubaraSampling","text":"MatsubaraSampling <: AbstractSampling\n\nSparse sampling in Matsubara frequencies.\n\nAllows the transformation between the IR basis and a set of sampling points in (scaled/unscaled) imaginary frequencies.\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.MatsubaraSampling-Tuple{SparseIR.AbstractBasis}","page":"Public","title":"SparseIR.MatsubaraSampling","text":"MatsubaraSampling(basis; positive_only=false,\n sampling_points=default_matsubara_sampling_points(basis; positive_only))\n\nConstruct a MatsubaraSampling object. If not given, the sampling_points are chosen as the (discrete) extrema of the highest-order basis function in Matsubara. This turns out to be close to optimal with respect to conditioning for this size (within a few percent).\n\nBy setting positive_only=true, one assumes that functions to be fitted are symmetric in Matsubara frequency, i.e.:\n\n G(iν) = conj(G(-iν))\n\nor equivalently, that they are purely real in imaginary time. In this case, sparse sampling is performed over non-negative frequencies only, cutting away half of the necessary sampling space.\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.RegularizedBoseKernel","page":"Public","title":"SparseIR.RegularizedBoseKernel","text":"RegularizedBoseKernel <: AbstractKernel\n\nRegularized bosonic analytical continuation kernel.\n\nIn dimensionless variables x = 2 τβ - 1, y = β ωΛ, the fermionic integral kernel is a function on -1 1 -1 1:\n\n K(x y) = y frace^-Λ y (x + 1) 2e^-Λ y - 1\n\nCare has to be taken in evaluating this expression around y = 0.\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.TauSampling","page":"Public","title":"SparseIR.TauSampling","text":"TauSampling <: AbstractSampling\n\nSparse sampling in imaginary time.\n\nAllows the transformation between the IR basis and a set of sampling points in (scaled/unscaled) imaginary time.\n\n\n\n\n\n","category":"type"},{"location":"public/#SparseIR.TauSampling-Tuple{SparseIR.AbstractBasis}","page":"Public","title":"SparseIR.TauSampling","text":"TauSampling(basis[; sampling_points])\n\nConstruct a TauSampling object. If not given, the sampling_points are chosen as the extrema of the highest-order basis function in imaginary time. This turns out to be close to optimal with respect to conditioning for this size (within a few percent).\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.evaluate!-Union{Tuple{N}, Tuple{T}, Tuple{S}, Tuple{AbstractArray{T, N}, SparseIR.AbstractSampling, AbstractArray{S, N}}} where {S, T, N}","page":"Public","title":"SparseIR.evaluate!","text":"evaluate!(buffer::AbstractArray{T,N}, sampling, al; dim=1) where {T,N}\n\nLike evaluate, but write the result to buffer. Please use dim = 1 or N to avoid allocating large temporary arrays internally.\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.evaluate-Union{Tuple{N}, Tuple{T}, Tuple{Tmat}, Tuple{S}, Tuple{SparseIR.AbstractSampling{S, Tmat}, AbstractArray{T, N}}} where {S, Tmat, T, N}","page":"Public","title":"SparseIR.evaluate","text":"evaluate(sampling, al; dim=1)\n\nEvaluate the basis coefficients al at the sparse sampling points.\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.fit!-Union{Tuple{N}, Tuple{T}, Tuple{S}, Tuple{Array{S, N}, SparseIR.AbstractSampling, Array{T, N}}} where {S, T, N}","page":"Public","title":"SparseIR.fit!","text":"fit!(buffer::Array{S,N}, smpl::AbstractSampling, al::Array{T,N}; \n dim=1, workarr::Vector{S}) where {S,T,N}\n\nLike fit, but write the result to buffer. Use dim = 1 or dim = N to avoid allocating large temporary arrays internally. The length of workarr cannot be smaller than SparseIR.workarrlength(smpl, al).\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.fit-Union{Tuple{N}, Tuple{T}, Tuple{Tmat}, Tuple{S}, Tuple{SparseIR.AbstractSampling{S, Tmat}, AbstractArray{T, N}}} where {S, Tmat, T, N}","page":"Public","title":"SparseIR.fit","text":"fit(sampling, al::AbstractArray{T,N}; dim=1)\n\nFit basis coefficients from the sparse sampling points Please use dim = 1 or N to avoid allocating large temporary arrays internally.\n\n\n\n\n\n","category":"method"},{"location":"public/#SparseIR.overlap-Union{Tuple{F}, Tuple{SparseIR.PiecewiseLegendrePoly, F}} where F","page":"Public","title":"SparseIR.overlap","text":"overlap(poly::PiecewiseLegendrePoly, f; \n rtol=eps(T), return_error=false, maxevals=10^4, points=T[])\n\nEvaluate overlap integral of poly with arbitrary function f.\n\nGiven the function f, evaluate the integral\n\n∫ dx f(x) poly(x)\n\nusing adaptive Gauss-Legendre quadrature.\n\npoints is a sequence of break points in the integration interval where local difficulties of the integrand may occur (e.g. singularities, discontinuities).\n\n\n\n\n\n","category":"method"},{"location":"","page":"Home","title":"Home","text":"CurrentModule = SparseIR","category":"page"},{"location":"#SparseIR.jl","page":"Home","title":"SparseIR.jl","text":"","category":"section"},{"location":"","page":"Home","title":"Home","text":"Documentation for SparseIR.jl.","category":"page"},{"location":"","page":"Home","title":"Home","text":"There is a guide available which details SparseIR.jl's inner workings by means of a worked example.","category":"page"},{"location":"","page":"Home","title":"Home","text":"For listings of all documented names, see Public names index and the Private names index.","category":"page"}]
+[{"location":"guide/#guide","page":"Guide","title":"Example usage and detailed explanation","text":"","category":"section"},{"location":"guide/","page":"Guide","title":"Guide","text":"We will explain the inner workings of SparseIR.jl by means of an example use case, adapted from the sparse-ir paper.","category":"page"},{"location":"guide/#Problem-statement","page":"Guide","title":"Problem statement","text":"","category":"section"},{"location":"guide/","page":"Guide","title":"Guide","text":"Let us perform self-consistent second-order perturbation theory for the single impurity Anderson model at finite temperature. Its Hamiltonian is given by H = U c^dagger_uparrow c^dagger_downarrow c_downarrow c_uparrow + sum_psigma big(V_psigma f_psigma^dagger c_sigma + V_psigma^* c_sigma^dagger c_sigma^daggerbig) + sum_psigma epsilon_p f_psigma^dagger f_psigmawhere U is the electron interaction strength, c_sigma annihilates an electron on the impurity, f_psigma annihilates an electron in the bath, dagger denotes the Hermitian conjugate, pinmathbb R is bath momentum, and sigmainuparrow downarrow is spin. The hybridization strength V_psigma and bath energies epsilon_p are chosen such that the non-interacting density of states is semi-elliptic with a half-bandwidth of one, rho_0(omega) = frac2pisqrt1-omega^2, U=12, beta=10, and the system is assumed to be half-filled.","category":"page"},{"location":"guide/#Treatment","page":"Guide","title":"Treatment","text":"","category":"section"},{"location":"guide/","page":"Guide","title":"Guide","text":"We first import SparseIR and construct an appropriate basis (omega_mathrmmax = 8 should be more than enough for this example):","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"julia> using SparseIR\n\njulia> basis = FiniteTempBasis(fermion, 10, 8)\nFiniteTempBasis{LogisticKernel, Float64}(fermion, 10.0, 8.0)","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"There's quite a lot happening behind the scenes in this first innocuous-looking statement, so let's break it down: Because we did not specify otherwise, the constructor chose the analytic continuation kernel for fermions, LogisticKernel(80.0), defined by","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"beginequation\n K(x y) = frace^-Λ y (x + 1) 21 + e^-Λ y\nendequation","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"for us, where 80.0 is the value of the scale parameter Lambda = betaomega_mathrmmax, shown below.","category":"page"},{"location":"guide/","page":"Guide","title":"Guide","text":"