quadrature_least_squares, a C++ code which computes weights for "sub-interpolatory" quadrature rules.
A large class of quadrature rules may be computed by specifying a set of N abscissas, or sample points, X(1:N), determining the Lagrange interpolation basis functions L(1:N), and then setting a weight vector W by
W(i) = I(L(i))
after which, the integral of any function f(x) is estimated by
I(f) \approx Q(f) = sum ( 1 <= i <= N ) W(i) * f(X(i))
We call this an interpolatory rule because the function f(x) has first been interpolated by
f(x) \approx sum ( 1 <= i <= N ) L(i) * f(X(i))
after which, we apply the integration operator:
I(f) \approx I(sum ( 1 <= i <= N ) L(i) * f(X(i)))
= sum ( 1 <= i <= N ) I(L(i)) * f(X(i))
= sum ( 1 <= i <= N ) W(i) * f(X(i)).
For badly chosen sets of X, or high values of N, or unruly functions f(x), interpolation may be a bad way to approximate the function. An alternative is to seek a polynomial interpolant of degree D < N-1, and then integrate that. We might call this a "sub-interpolatory" rule.
As it turns out, a natural way to seek such a rule is to write out the N by D+1 Vandermonde matrix and use a least squares solver. Even though the N by N Vandermonde matrix is ill-conditioned for Gauss elimination, a least squares approach can produce usable solutions from the N by D+1 matrix.
The outline of this procedure was devised by Professor Mac Hyman of Tulane University.
The computer code and data files described and made available on this web page are distributed under the MIT license
quadrature_least_squares is available in a C version and a C++ version and a FORTRAN90 version and a MATLAB version.
CLENSHAW_CURTIS_RULE, a C++ code which defines a multiple dimension Clenshaw Curtis quadrature rule.
QR_SOLVE, a C++ code which computes the least squares solution of a rectangular linear system A*x=b.
QUADMOM, a C++ code which computes a Gaussian quadrature rule for a weight function rho(x) based on the Golub-Welsch procedure that only requires knowledge of the moments of rho(x).
QUADRATURE_GOLUB_WELSCH, a C++ code which computes the points and weights of a Gaussian quadrature rule using the Golub-Welsch procedure, assuming that the points have been specified.
QUADRATURE_WEIGHTS_VANDERMONDE, a C++ code which computes the weights of a quadrature rule using the Vandermonde matrix, assuming that the points have been specified.
QUADRULE, a C++ code which defines quadrature rules for approximating an integral over a 1D domain.
TEST_INT, a C++ code which defines test integrands for 1D quadrature rules.
TRUNCATED_NORMAL_RULE, a C++ code which computes a quadrature rule for a normal probability density function (PDF), also called a Gaussian distribution, that has been truncated to [A,+oo), (-oo,B] or [A,B].
VANDERMONDE, a C++ code which carries out certain operations associated with the Vandermonde matrix.