gauss.quad | R Documentation |
Calculate nodes and weights for Gaussian quadrature.
gauss.quad(n, kind = "legendre", alpha = 0, beta = 0)
n |
number of nodes and weights |
kind |
kind of Gaussian quadrature, one of |
alpha |
parameter for Jacobi or Laguerre quadrature, must be greater than -1 |
beta |
parameter for Jacobi quadrature, must be greater than -1 |
The integral from a
to b
of w(x)*f(x)
is approximated by sum(w*f(x))
where x
is the vector of nodes and w
is the vector of weights. The approximation is exact if f(x)
is a polynomial of order no more than 2n-1
.
The possible choices for w(x)
, a
and b
are as follows:
Legendre quadrature: w(x)=1
on (-1,1)
.
Chebyshev quadrature of the 1st kind: w(x)=1/sqrt(1-x^2)
on (-1,1)
.
Chebyshev quadrature of the 2nd kind: w(x)=sqrt(1-x^2)
on (-1,1)
.
Hermite quadrature: w(x)=exp(-x^2)
on (-Inf,Inf)
.
Jacobi quadrature: w(x)=(1-x)^alpha*(1+x)^beta
on (-1,1)
. Note that Chebyshev quadrature is a special case of this.
Laguerre quadrature: w(x)=x^alpha*exp(-x)
on (0,Inf)
.
The algorithm used to generated the nodes and weights is explained in Golub and Welsch (1969).
A list containing the components
nodes |
vector of values at which to evaluate the function |
weights |
vector of weights to give the function values |
Gordon Smyth, using Netlib Fortran code https://netlib.org/go/gaussq.f, and including a suggestion from Stephane Laurent
Golub, G. H., and Welsch, J. H. (1969). Calculation of Gaussian quadrature rules. Mathematics of Computation 23, 221-230.
Golub, G. H. (1973). Some modified matrix eigenvalue problems. Siam Review 15, 318-334.
Smyth, G. K. (1998). Numerical integration. In: Encyclopedia of Biostatistics, P. Armitage and T. Colton (eds.), Wiley, London, pages 3088-3095. http://www.statsci.org/smyth/pubs/NumericalIntegration-Preprint.pdf
Smyth, G. K. (1998). Polynomial approximation. In: Encyclopedia of Biostatistics, P. Armitage and T. Colton (eds.), Wiley, London, pages 3425-3429. http://www.statsci.org/smyth/pubs/PolyApprox-Preprint.pdf
Stroud, AH, and Secrest, D (1966). Gaussian Quadrature Formulas. Prentice-Hall, Englewood Cliffs, N.J.
gauss.quad.prob
, integrate
# mean of gamma distribution with alpha=6 out <- gauss.quad(10,"laguerre",alpha=5) sum(out$weights * out$nodes) / gamma(6)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.