Skip to content

Rust crate for Gauss-Legendre quadrature with the method of I. Bogaert. This makes the computation of node-weight pairs O(1) and easily parallelized

Notifications You must be signed in to change notification settings

JSorngard/fast_gauleg

Folders and files

NameName
Last commit message
Last commit date

Latest commit

 
 
 
 
 
 
 
 
 

Repository files navigation

fast_gauleg

A crate for numerical integration with Gauss-Legendre quadrature. This crate uses the method described in the paper "Iteration-Free Computation of Gauss-Legendre Quadrature Nodes and Weights" by I. Bogaert which enables computation of node-weight pairs in O(1) time complexity (and optionally in parallel) while maintaining an accuracy of a few ulps in the nodes and weights (see paper for details).

The integrals are computed with the formula

$$\int_a^b\!f(x)\;\mathrm{d}x\approx\frac{b-a}{2}\sum_{k=1}^nw_n(k)f\left(\frac{b-a}{2}x_n(k)+\frac{a+b}{2}\right)$$

where $n$ is a chosen integer, $a$ and $b$ are finite numbers, and $w_n$ and $x_n$ are specific numbers derived from Legendre polynomials.

The serial parts of this code has been merged into the gauss_quad crate in this PR.

Examples

Integrate a degree 5 polynomial while only evaluating it at three points:

use fast_gauleg::glq_integrate;
use approx::assert_relative_eq;

assert_relative_eq!(
    glq_integrate(-1.0, 1.0, |x| x.powf(5.0), 3.try_into().unwrap()),
    0.0,
);

Very large quadrature rules are feasible, and can be evaluated in parallel:

assert_relative_eq!(
    par_glq_integrate(0.0, 1.0, |x| x.ln(), 100_000_000.try_into().unwrap()),
    -1.0
);

This is because the computation of individual node-weight pairs is O(1)

use fast_gauleg::GlqPair;
// This:
let a = GlqPair::new(1_000_000_000, 500_000_000);
// Is as fast as this:
let b = GlqPair::new(24, 12);

About

Rust crate for Gauss-Legendre quadrature with the method of I. Bogaert. This makes the computation of node-weight pairs O(1) and easily parallelized

Topics

Resources

Stars

Watchers

Forks

Releases

No releases published

Packages

No packages published

Languages