Gegenbauer | R Documentation |
The Gegenbauer polynomials
\{C_k^{(\lambda)}(x)\}_{k = 0}^\infty
form a family of orthogonal polynomials on the interval [-1, 1]
with respect to the weight function (1 - x^2)^{\lambda - 1/2}
,
for \lambda > -1/2
, \lambda \neq 0
. They usually appear
when dealing with functions defined on
S^{p-1} := \{{\bf x} \in R^p : ||{\bf x}|| = 1\}
with index \lambda = p / 2 - 1
.
The Gegenbauer polynomials are somehow simpler to evaluate for
x = \cos(\theta)
, with \theta \in [0, \pi]
. This simplifies
also the connection with the Chebyshev polynomials
\{T_k(x)\}_{k = 0}^\infty
, which admit
the explicit expression
T_k(\cos(\theta)) = \cos(k\theta)
. The Chebyshev polynomials
appear as the limit of the Gegenbauer polynomials
(divided by \lambda
) when \lambda
goes to 0
, so they
can be regarded as the extension by continuity of
\{C_k^{(p/2 - 1)}(x)\}_{k = 0}^\infty
to the case p = 2
.
For a reasonably smooth function
\psi
defined on [0, \pi]
,
\psi(\theta) = \sum_{k = 0}^\infty b_{k, p}
C_k^{(p/2 - 1)}(\cos(\theta)),
provided that the coefficients
b_{k, p} := \frac{1}{c_{k, p}} \int_0^\pi \psi(\theta)
C_k^{(p/2 - 1)}(\cos(\theta)) (\sin(\theta))^{p - 2}\,\mathrm{d}\theta
are finite, where the normalizing constants are
c_{k, p} := \int_0^\pi (C_k^{(p/2 - 1)}(\cos(\theta)))^2
(\sin(\theta))^{p - 2} \,\mathrm{d}\theta.
The (squared) "Gegenbauer norm" of \psi
is
\|\psi\|_{G, p}^2 := \int_0^\pi \psi(\theta)^2
C_k^{(p/2 - 1)}(\cos(\theta)) (\sin(\theta))^{p - 2}\,\mathrm{d}\theta.
The previous expansion can be generalized for a 2-dimensional function
\psi
defined on [0, \pi] \times [0, \pi]
:
\psi(\theta_1, \theta_2) = \sum_{k = 0}^\infty \sum_{m = 0}^\infty
b_{k, m, p} C_k^{(p/2 - 1)}(\cos(\theta_1))
C_k^{(p/2 - 1)}(\cos(\theta_2)),
with coefficients
b_{k, m, p} := \frac{1}{c_{k, p} c_{m, p}} \int_0^\pi\int_0^\pi
\psi(\theta_1, \theta_2) C_k^{(p/2 - 1)}(\cos(\theta_1))
C_k^{(p/2 - 1)}(\cos(\theta_2)) (\sin(\theta_1))^{p - 2}
(\sin(\theta_2))^{p - 2}\,\mathrm{d}\theta_1\,\mathrm{d}\theta_2.
The (squared) "Gegenbauer norm" of \psi
is
\|\psi\|_{G, p}^2 := \int_0^\pi\int_0^\pi \psi(\theta_1, \theta_2)^2
C_k^{(p/2 - 1)}(\cos(\theta_1)) C_k^{(p/2 - 1)}(\cos(\theta_2))
(\sin(\theta_1))^{p - 2} (\sin(\theta_2))^{p - 2}
\,\mathrm{d}\theta_1\,\mathrm{d}\theta_2.
Gegen_polyn(theta, k, p)
Gegen_coefs(k, p, psi, Gauss = TRUE, N = 320, normalize = TRUE,
only_const = FALSE, tol = 1e-06, ...)
Gegen_series(theta, coefs, k, p, normalize = TRUE)
Gegen_norm(coefs, k, p, normalize = TRUE, cumulative = FALSE)
Gegen_polyn_2d(theta_1, theta_2, k, m, p)
Gegen_coefs_2d(k, m, p, psi, Gauss = TRUE, N = 320, normalize = TRUE,
only_const = FALSE, tol = 1e-06, ...)
Gegen_series_2d(theta_1, theta_2, coefs, k, m, p, normalize = TRUE)
Gegen_norm_2d(coefs, k, m, p, normalize = TRUE)
theta , theta_1 , theta_2 |
vectors with values in |
k , m |
vectors with the orders of the Gegenbauer polynomials. Must
be integers larger or equal than |
p |
integer giving the dimension of the ambient space |
psi |
function defined in |
Gauss |
use a Gauss–Legendre quadrature rule of |
N |
number of points used in the
Gauss–Legendre quadrature for computing the Gegenbauer coefficients.
Defaults to |
normalize |
consider normalized coefficients (divided by
|
only_const |
return only the normalizing constants |
tol |
tolerance passed to |
... |
further arguments to be passed to |
coefs |
for |
cumulative |
return the cumulative norm for increasing truncation of
the series? Defaults to |
The Gegen_polyn
function is a wrapper to the functions
gegenpoly_n
and
gegenpoly_array
in the
gsl-package
, which they interface the functions
defined in the header file gsl_sf_gegenbauer.h
(documented
here) of the
GNU Scientific Library.
Note that the function Gegen_polyn
computes the regular
unnormalized Gegenbauer polynomials.
For the case p = 2
, the Chebyshev polynomials are considered.
Gegen_polyn
: a matrix of size
c(length(theta), length(k))
containing the evaluation of the
length(k)
Gegenbauer polynomials at theta
.
Gegen_coefs
: a vector of size length(k)
containing
the coefficients b_{k, p}
.
Gegen_series
: the evaluation of the truncated series
expansion, a vector of size length(theta)
.
Gegen_norm
: the Gegenbauer norm of the truncated series,
a scalar if cumulative = FALSE
, otherwise a vector of size
length(k)
.
Gegen_polyn_2d
: a 4-dimensional array of size
c(length(theta_1), length(theta_2), length(k), length(m))
containing the evaluation of the length(k) * length(m)
2-dimensional Gegenbauer polynomials at the bivariate grid
spanned by theta_1
and theta_2
.
Gegen_coefs_2d
: a matrix of size
c(length(k), length(m))
containing the coefficients
b_{k, m, p}
.
Gegen_series_2d
: the evaluation of the truncated series
expansion, a matrix of size c(length(theta_1), length(theta_2))
.
Gegen_norm_2d
: the 2-dimensional Gegenbauer norm of the
truncated series, a scalar.
Galassi, M., Davies, J., Theiler, J., Gough, B., Jungman, G., Alken, P., Booth, M., and Rossi, F. (2009) GNU Scientific Library Reference Manual. Network Theory Ltd. http://www.gnu.org/software/gsl/
NIST Digital Library of Mathematical Functions. Release 1.0.20 of 2018-09-15. F. W. J. Olver, A. B. Olde Daalhuis, D. W. Lozier, B. I. Schneider, R. F. Boisvert, C. W. Clark, B. R. Miller, and B. V. Saunders, eds. https://dlmf.nist.gov/
## Representation of Gegenbauer polynomials (Chebyshev polynomials for p = 2)
th <- seq(0, pi, l = 500)
k <- 0:3
old_par <- par(mfrow = c(2, 2))
for (p in 2:5) {
matplot(th, t(Gegen_polyn(theta = th, k = k, p = p)), lty = 1,
type = "l", main = substitute(p == d, list(d = p)),
axes = FALSE, xlab = expression(theta), ylab = "")
axis(1, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
axis(2); box()
mtext(text = expression({C[k]^{p/2 - 1}}(cos(theta))), side = 2,
line = 2, cex = 0.75)
legend("bottomleft", legend = paste("k =", k), lwd = 2, col = seq_along(k))
}
par(old_par)
## Coefficients and series in p = 2
# Function in [0, pi] to be projected in Chebyshev polynomials
psi <- function(th) -sin(th / 2)
# Coefficients
p <- 2
k <- 0:4
(coefs <- Gegen_coefs(k = k, p = p, psi = psi))
# Series
plot(th, psi(th), type = "l", axes = FALSE, xlab = expression(theta),
ylab = "", ylim = c(-1.25, 0))
axis(1, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
axis(2); box()
col <- viridisLite::viridis(length(coefs))
for (i in seq_along(coefs)) {
lines(th, Gegen_series(theta = th, coefs = coefs[1:(i + 1)], k = 0:i,
p = p), col = col[i])
}
lines(th, psi(th), lwd = 2)
## Coefficients and series in p = 3
# Function in [0, pi] to be projected in Gegenbauer polynomials
psi <- function(th) tan(th / 3)
# Coefficients
p <- 3
k <- 0:10
(coefs <- Gegen_coefs(k = k, p = p, psi = psi))
# Series
plot(th, psi(th), type = "l", axes = FALSE, xlab = expression(theta),
ylab = "", ylim = c(0, 2))
axis(1, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
axis(2); box()
col <- viridisLite::viridis(length(coefs))
for (i in seq_along(coefs)) {
lines(th, Gegen_series(theta = th, coefs = coefs[1:(i + 1)], k = 0:i,
p = p), col = col[i])
}
lines(th, psi(th), lwd = 2)
## Surface representation
# Surface in [0, pi]^2 to be projected in Gegenbauer polynomials
p <- 3
psi <- function(th_1, th_2) A_theta_x(theta = th_1, x = cos(th_2),
p = p, as_matrix = TRUE)
# Coefficients
k <- 0:20
m <- 0:10
coefs <- Gegen_coefs_2d(k = k, m = m, p = p, psi = psi)
# Series
th <- seq(0, pi, l = 100)
col <- viridisLite::viridis(20)
old_par <- par(mfrow = c(2, 2))
image(th, th, A_theta_x(theta = th, x = cos(th), p = p), axes = FALSE,
col = col, zlim = c(0, 1), xlab = expression(theta[1]),
ylab = expression(theta[2]), main = "Original")
axis(1, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
axis(2, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
box()
for(K in c(5, 10, 20)) {
A <- Gegen_series_2d(theta_1 = th, theta_2 = th,
coefs = coefs[1:(K + 1), ], k = 0:K, m = m, p = p)
image(th, th, A, axes = FALSE, col = col, zlim = c(0, 1),
xlab = expression(theta[1]), ylab = expression(theta[2]),
main = paste(K, "x", m[length(m)], "coefficients"))
axis(1, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
axis(2, at = c(0, pi / 4, pi / 2, 3 * pi / 4, pi),
labels = expression(0, pi / 4, pi / 2, 3 * pi / 4, pi))
box()
}
par(old_par)
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.