autosymptr(f, B::AbstractMatrix, [syms=nothing]; atol=0, rtol=sqrt(eps()), maxevals=typemax(Int64), buffer=nothing)

Computes the integral of f to within the specified tolerances and returns a tuple (I, E, numevals, rules) containing the estimated integral, the estimated error, the number of integrand evaluations. The integral on the symmetrized domain needs to be mapped back to the full domain by the caller, since for array-valued integrands this depends on the representation of the integrand under the action of the symmetries.

Note that if a vector buffer is provided to store integrand evaluations, the integrand evaluations will be parallelized and it will be assumed that the integrand is threadsafe.

Convergence depends on periodicity

If the routine takes a long time to return, double check that the period of the function f along the basis vectors in the columns of B is consistent.


These methods are internal and may change between AutoSymPTR.jl releases

PTR{d}(x::AbstractVector{T}) where {d,T}

Stores a d dimensional Cartesian product grid of SVector{d,T}. Similar to Iterators.product(ntuple(n->x, d)...). Uses the same number of grid points per dimension.

alloc_cache(::Type{T}, ::Val{d}, rule)

Initialize an empty buffer of PTR rules evaluated from rule(T,Val(d)) where T is the domain type and d is the number of dimensions.

For developers

Providing a special rule, (e.g. PTR)

AffineQuad(rule, A, b, c, vol)

A wrapper to an iterable and indexable quadrature rule that applies an affine coordinate transformation to the nodes of the form A*(x+c)+b. While the weights should be rescaled by the Jacobian determinant vol=abs(det(A)), the value is stored in the vol field of the struct so that the caller can minimize multiplication by applying the rescaling after computing the rule.

While the quadrature rule should be unitless, the affine transform may have units.

symptr_rule(npt, ::Val{d}, syms, T::Integer=UInt8) where {d}

Returns wsym, a npt^d array containing the T::Integer weights of the symmetrized PTR quadrature, aka the Monkhorst Pack weights, flag a tuple of boolean arrays that indicate whether a slice of wsym has nonzero entries, and nsym the number of non-zero weights. (see Algorithm 3. in Kaye et al.). The algorithm is serial and has time and memory requirements scaling as O(npt^d).

The arithmetic precision for the calculation is inferred from the element type of syms and the integer size needs to be increased proportionally to the size of syms. The default value of UInt8 is sufficient for the 48 possible symmetries in 3d.

InplaceIntegrand(f!, result::AbstractArray)

Constructor for a InplaceIntegrand accepting an integrand of the form f!(y,x). The caller also provides an output array needed to store the result of the quadrature. Intermediate y arrays are allocated during the calculation, and the final result is may or may not be written to result, so use the IntegralSolution immediately after the calculation to read the result, and don't expect it to persist if the same integrand is used for another calculation.