# R/chebyshev.R In pracma: Practical Numerical Math Functions

#### Documented in chebApproxchebCoeffchebPoly

```##
##  c h e b y s h e v . R
##

chebPoly <- function(n, x = NULL) {
stopifnot(is.numeric(n), length(n) == 1,
floor(n) == ceiling(n), n >= 0)

N <- max(2, n+1)
T <- matrix(0, N, N)

# Preset degree 0 and 1
T[1, 1] <- 1
T[2, 2] <- 1

# Use recursion formula
if (n >= 2)
for (i in 3:N)
T[i, ] <- 2 * c(0, T[i-1, 1:(N-1)]) - T[i-2, ]

# Reverse with highest coefficient first
T <- T[ ,ncol(T):1]

if  (is.null(x)) {
return(T)
} else if (is.numeric(x)) {
return(polyval(T[n+1, ], x))
} else
stop("Argument 'x' must be a numeric vector.")
}

chebCoeff <- function(fun, a, b, n) {
N <- n+1
# Map interval [a, b] to [-1, 1]
c1 <- 0.5*(b-a)
c2 <- 0.5*(b+a)

# Evaluate function at Chebyshev points (in [a, b])
k  <- c(0:(N-1))
y  <- cos(pi*(k+0.5)/N)
fy <- fun(c1*y + c2)

# Now compute the Chebyshev coefficients
c0 <- 2.0 / N
K  <- matrix(k+0.5, N, 1)
cheb <- c0 * fy %*% cos(pi/N * K %*% k)

# Remove too small coefficients and return
eps <- 2 * .Machine\$double.eps    # machine precision
cheb <- ifelse(abs(cheb) < eps, 0, cheb)
return(drop(cheb))
}

chebApprox <- function(x, fun, a, b, n) {
# Compute the Chebyshev polynomials of fun in [a, b]
cP <- chebPoly(n)
cC <- chebCoeff(fun, a, b, n)
p  <- drop(cC %*% cP)
c0 <- cC[1]

# Map x into [-1, 1] and evaluate the Chebyshev polynomial
xx <- (2*x - (b+a))/(b-a)
yy <- polyval(p, xx) - c0/2
return(yy)
}
```

## Try the pracma package in your browser

Any scripts or data that you put into this service are public.

pracma documentation built on Nov. 10, 2023, 1:14 a.m.