4

在数值分析中,我们的学生有义务在 R 中实现代码,给定函数 f(x) 找到其傅里叶插值 tN(x) 并计算插值误差

$||f(x)-t_{N}(x)||=\int_{0}^{2\pi}$ $|f(x)-t_{N}(x)|^2$ 

或各种不同的 $N$ 我首先尝试根据这个公式计算 d 系数:

$d = \frac 1N  M  y$

其中 M 表示 DFT 矩阵,y 表示一系列等距函数值

$y_j = f(x_j)$ and 
$x_j = e^{\frac{2*pi*i}N*j}$ 
for $j = 1,..,N-1$. 

我的目标是得出一个可以描述为的总和:

$t_{N}(x) = \Sigma_{k=0}^{N-1} d_k * e^{i*k*x}$

稍后将更容易集成到随后的加法符号中。

f <- function(x) 3/(6+4*cos(x)) #first function to compare with
g <- function(x) sin(32*x) #second one
xj <- function(x,n) 2*pi*x/n

M <- function(n){
   w = exp(-2*pi*1i/n)
   m = outer(0:(n-1),0:(n-1))
   return(w^m)
}

y <- function(n){
   f(xj(0:(n-1),n))
} 
transformFunction <- function(n, f){
   d = 1/n * t(M(n)) %*% f(xj(0:(n-1),n))
   script <- paste(d[1])
   for(i in 2:n)
   script <- paste0(script,paste0("+",d[i],"*exp(1i*x*",i,")"))
   #trans <- sum(d[1:n] * exp(1i*x*(0:(n-1))))
   return(script)
 } 

最初,变换函数的主要目的是返回一个函数——或者更确切地说:一个数学表达式——然后可以使用它来声明我的傅里叶插值函数。问题是,根据我相当有限的知识,我无法集成仍然嵌套了总和的函数(这就是我在代码中注释相应行的原因)。出于绝对的绝望,我随后尝试以文本形式粘贴每个总和,只是为了再次将它们解析为表达式。所以剩下的主要问题是:我如何以允许我将它们用作函数并稍后集成它们的方式返回数学表达式?对于任何误解或混淆,以及我看似业余的编码,我深表歉意。提前致谢!

4

1 回答 1

1

R 中的函数可以返回任何类,特别是 class 的对象function。因此,您可以创建trans一个函数x并返回它。

由于该integrate函数需要矢量化函数,因此我们Vectorize在输出之前使用。

transformFunction <- function(n, f){
    d = 1/n * t(M(n)) %*% f(xj(0:(n-1),n))

    ## Output function
    trans <- function(x) sum(d[1:n] * exp(1i*x*(0:(n-1))))
    ## Vectorize output for the integrate function
    Vectorize(trans)
} 

要集成,现在只需使用以下输出创建一个新变量transformFunction

myint <- transformFunction(n = 10,f = f)

测试:(integrate只能处理实值函数)

integrate(function(x) Re(myint(x)),0,2)$value + 
    1i*integrate(function(x) Im(myint(x)),0,2)$value
# [1] 1.091337-0.271636i
于 2016-06-14T07:47:17.663 回答