weightsGJ {rjacobi}R Documentation

Weights of Gauss-Jacobi Quadrature

Description

Calculates the weights of the Q point Gauss-Jacobi quadrature.

Usage

z = weightsGJ(z, a=0.0, b=0.0)
z = weightsGLJ(z, a=0.0, b=0.0)
z = weightsGRJM(z, a=0.0, b=0.0)
z = weightsGRJP(z, a=0.0, b=0.0)

Arguments

z Quadrature node calculated with zerosGJ according to quadrature type
a α parameter of Jacobi polynomial. Defaults to 0
b β parameter of Jacobi polynomial. Defaults to 0

Details

Different types of quadrature are possible:

These functions find the quadrature nodes that will be used in subsequent calculations.

With these weights, integrals can be calculated according to the following expression:

int_{-1}^{1} f(x)dx approx sum_{i=1}{Q}w_i f(x_i)

Value

A vector containing the weights of the quadrature

References

Abramowitz, Milton and Stegun, Irene (editors); "Handbook of Mathematical functions", Dover Publications 1965.

Karniadakis, George Em and Sherwin, Spencer; "Spectral/hp Element Methods for Computational Fluid Dynamics", Oxford Science Publications, 2nd edition, 2005.

See Also

zerosGJ diffmatGJ

Examples

### This example will shoud the exponential convergence when integrating
### smooth functions.

Q <- 2:15
f <- function(x)cos(5*x)
n <- length(Q)
integr <- double(n)

for (i in 1:n){
  z <- zerosGJ(Q[i])
  w <- weightsGJ(z)
  integr[i] <- sum(w * f(z))
}
err <- abs(integr - 2*sin(5)/5)

plot(Q, err, ty='b', log='y', main="Error of integral of cos(5*x)",
xlab="No Quadrature points", ylab="Error")




[Package rjacobi version 0.9.2 Index]