lagrange {rjacobi}R Documentation

Lagrange Interpolants of Gauss-Jacobi Quadrature

Description

Calculates the ith Lagrange interpolation polinomial for Q point Gauss-Jacobi quadrature.

Usage

y = lagrangeGJ(x, i, z, a=0, b=0)
y = lagrangeGLJ(x, i, z, a=0, b=0)
y = lagrangeGRJM(x, i, z, a=0, b=0)
y = lagrangeGRJP(x, i, z, a=0, b=0)

Arguments

x Points where to calculate the Lagrange polynomial
i Which interpolant to calculate
z Quadrature nodes calculated by zerosGJ family of functions
a α parameter of Jacobi polynomial. Defaults to 0
b β parameter of Jacobi polynomial. Defaults to 0

Details

This function calculates

y = h_i(x)

where x is any point in the domain -1 < x < 1 and i is the corresponding Lagrangean interpolator with 1 <= i <= Q.

Different types of quadrature are possible:

This function is used internally (in C code) to calculate the interpolation matrix.

Value

The value of the Lagrangean interpolator

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

quadrature zerosGJ interpmatGJ

Examples


## Plots the Lagrange interpolator for GLJ Q=6

z <- zerosGLJ(6)

x <- seq(-1, 1, len=101)

plot(x, lagrangeGLJ(x, 1, z), ty='l', xlab="x", ylab="y",
main="Lagrange polynomials, GLJ, Q=6")
lines(x, lagrangeGLJ(x, 2, z))
lines(x, lagrangeGLJ(x, 3, z))
lines(x, lagrangeGLJ(x, 4, z))
lines(x, lagrangeGLJ(x, 5, z))
lines(x, lagrangeGLJ(x, 6, z))

abline(v=z, col='red')
abline(h=1, col='red')


[Package rjacobi version 0.9.2 Index]