# Integral algorithms

AutoBZCore.QuadratureFunctionType
QuadratureFunction(; fun=trapz, npt=50, nthreads=1)

Quadrature rule for the standard interval [-1,1] computed from a function x, w = fun(npt). The nodes and weights should be set so the integral of f on [-1,1] is sum(w .* f.(x)). The default quadrature rule is trapz, although other packages provide rules, e.g.

using FastGaussQuadrature
alg = QuadratureFunction(fun=gausslegendre, npt=100)

nthreads sets the numbers of threads used to parallelize the quadrature only when the integrand is a BatchIntegrand, in which case the user must parallelize the integrand evaluations. For no threading set nthreads=1.

AutoBZCore.AuxQuadGKJLType
AuxQuadGKJL(; order = 7, norm = norm)

Generalization of the QuadGKJL provided by Integrals.jl that allows for AuxValued integrands for auxiliary integration and multi-threaded evaluation with the batch argument to IntegralProblem

AutoBZCore.ContQuadGKJLType
ContQuadGKJL(; order = 7, norm = norm, rho = 1.0, rootmeth = IteratedIntegration.ContQuadGK.NewtonDeflation())

A 1d contour deformation quadrature scheme for scalar, complex-valued integrands. It defaults to regular quadgk behavior on the real axis, but if it finds a root of 1/f nearby, in the sense of Bernstein ellipse for the standard segment [-1,1] with semiaxes cosh(rho) and sinh(rho), on either the upper/lower half planes, then it dents the contour away from the presumable pole.

AutoBZCore.MeroQuadGKJLType
MeroQuadGKJL(; order = 7, norm = norm, rho = 1.0, rootmeth = IteratedIntegration.MeroQuadGK.NewtonDeflation())

A 1d pole subtraction quadrature scheme for scalar, complex-valued integrands that are meromorphic. It defaults to regular quadgk behavior on the real axis, but if it finds nearby roots of 1/f, in the sense of Bernstein ellipse for the standard segment [-1,1] with semiaxes cosh(rho) and sinh(rho), it attempts pole subtraction on that segment.

## Cubature

AutoBZCore.MonkhorstPackType
MonkhorstPack(; npt=50, syms=nothing, nthreads=1)

Periodic trapezoidal rule with a fixed number of k-points per dimension, npt, using the PTR rule from AutoSymPTR.jl. nthreads sets the numbers of threads used to parallelize the quadrature only when the integrand is a BatchIntegrand, in which case the user must parallelize the integrand evaluations. For no threading set nthreads=1. The caller should check that the integral is converged w.r.t. npt.

AutoBZCore.AutoSymPTRJLType
AutoSymPTRJL(; norm=norm, a=1.0, nmin=50, nmax=1000, n₀=6, Δn=log(10), keepmost=2, nthreads=1)

Periodic trapezoidal rule with automatic convergence to tolerances passed to the solver with respect to norm using the routine autosymptr from AutoSymPTR.jl. nthreads sets the numbers of threads used to parallelize the quadrature only when the integrand is a BatchIntegrand, in which case the user must parallelize the integrand evaluations. For no threading set nthreads=1. This algorithm is the most efficient for smooth integrands.

## Meta-algorithms

AutoBZCore.NestedQuadType
NestedQuad(alg::IntegralAlgorithm)
NestedQuad(algs::IntegralAlgorithm...)

Nested integration by repeating one quadrature algorithm or composing a list of algorithms. The domain of integration must be an AbstractIteratedLimits from the IteratedIntegration.jl package. Analogous to nested_quad from IteratedIntegration.jl. The integrand should expect SVector inputs. Do not use this for very high-dimensional integrals, since the compilation time scales very poorly with respect to dimensionality. In order to improve the compilation time, FunctionWrappers.jl is used to enforce type stability of the integrand, so you should always pick the widest integration limit type so that inference works properly. For example, if ContQuadGKJL is used as an algorithm in the nested scheme, then the limits of integration should be made complex.

AutoBZCore.EvalCounterType
EvalCounter(::IntegralAlgorithm)

An algorithm which counts the evaluations used by another algorithm. The count is stored in the sol.numevals field.

AutoBZCore.AbsoluteEstimateType
AbsoluteEstimate(est_alg, abs_alg; kws...)

Most algorithms are efficient when using absolute error tolerances, but how do you know the size of the integral? One option is to estimate it using second algorithm.

A multi-algorithm to estimate an integral using an est_alg to generate a rough estimate of the integral that is combined with a user's relative tolerance to re-calculate the integral to higher accuracy using the abs_alg. The keywords passed to the algorithm may include reltol, abstol and maxiters and are given to the est_alg solver. They should limit the amount of work of est_alg so as to only generate an order-of-magnitude estimate of the integral. The tolerances passed to abs_alg are abstol=max(abstol,reltol*norm(I)) and reltol=0.

# BZ-specific integral algorithms

In order to make algorithms domain-agnostic, the BZ loaded from load_bz can be called with the algorithms below, which are wrappers for algorithms above with the additional capability of mapping integrals over the IBZ to the FBZ.

AutoBZCore.AutoBZAlgorithmType
AutoBZAlgorithm

Abstract supertype for Brillouin zone integration algorithms. All integration problems on the BZ get rescaled to fractional coordinates so that the Brillouin zone becomes [0,1]^d, and integrands should have this periodicity. If the integrand depends on the Brillouin zone basis, then it may have to be transformed to the Cartesian coordinates as a post-processing step. These algorithms also use the symmetries of the Brillouin zone and the integrand.

AutoBZCore.IAIType
IAI(alg::IntegralAlgorithm=AuxQuadGKJL())
IAI(algs::IntegralAlgorithm...)

Iterated-adaptive integration using nested_quad from IteratedIntegration.jl. This algorithm is the most efficient for localized integrands.

AutoBZCore.TAIType
TAI(; norm=norm, initdivs=1)

Tree-adaptive integration using hcubature from HCubature.jl. This routine is limited to integration over hypercube domains and may not use all symmetries.

AutoBZCore.PTRType
PTR(; npt=50, nthreads=1)

Periodic trapezoidal rule with a fixed number of k-points per dimension, npt, using the routine ptr from AutoSymPTR.jl. The caller should check that the integral is converged w.r.t. npt.

AutoBZCore.AutoPTRType
AutoPTR(; norm=norm, a=1.0, nmin=50, nmax=1000, n₀=6, Δn=log(10), keepmost=2, nthreads=1)

Periodic trapezoidal rule with automatic convergence to tolerances passed to the solver with respect to norm using the routine autosymptr from AutoSymPTR.jl. This algorithm is the most efficient for smooth integrands.

AutoBZCore.PTR_IAIFunction
PTR_IAI(; ptr=PTR(), iai=IAI())

Multi-algorithm that returns an IAI calculation with an abstol determined from the given reltol and a PTR estimate, I, as reltol*norm(I). This addresses the issue that IAI does not currently use a globally-adaptive algorithm and may not have the expected scaling with localization length unless an abstol is used since computational effort may be wasted via a reltol with the naive nested_quadgk.

AutoBZCore.AutoPTR_IAIFunction
AutoPTR_IAI(; reltol=1.0, ptr=AutoPTR(), iai=IAI())

Multi-algorithm that returns an IAI calculation with an abstol determined from an AutoPTR estimate, I, computed to reltol precision, and the rtol given to the solver as abstol=rtol*norm(I). This addresses the issue that IAI does not currently use a globally-adaptive algorithm and may not have the expected scaling with localization length unless an abstol is used since computational effort may be wasted via a reltol with the naive nested_quadgk.