MATLAB自定义逆DFT函数存在微小误差，求技术分析

我因实际需求需在实空间中使用傅里叶模式基的子集，遂尝试在MATLAB中编写自定义逆DFT函数。但提取傅里叶系数后，将自定义逆DFT程序与MATLAB自带ifft函数及原始波形对比时，发现存在微小误差。不确定该问题是浮点运算问题还是频率bins设置错误，这是我首次使用fft，恳请提供技术见解。

附上测试代码：

% 1D fft and inverse fft
clear; clc; close all;
L=1; % total domain size
n=100;
x=linspace(0,L,n); %define x values
dx=x(1)-x(2);
F=1*cos(2*pi*5*x); %declare the function
figure(1)
plot(x,F,'linewidth',2,'Displayname','Original');
Fhat=fftshift(fft(F));%compute fft and shift to center
if(mod(n,2)==0) % even split
    dk=2*pi/L.*[-n/2:1:1:n/2-1];
else %odd split
    dk=2*pi/L*[-(n-1)/2:1:(n-1)/2];
end
figure % plot FFt results
plot(dk/(2*pi),abs(Fhat),'displayname','Frequency')
% Now here is my inverse fft
f=zeros(1,n);
for i=1:n
    f=f+Fhat(i).*exp(1j*dk(i).*x)
end
f=f*1/n;
figure(1); hold on;
plot(x,real(f),'b-.','linewidth',2,'displayname','From my own invDFT')
plot(x,ifft(ifftshift(Fhat)),'M--','linewidth',2,'displayname','From matlab ifft')
xlabel('x','fontsize',16)
ylabel('y','fontsize',16)
legend()

波形对比图：

问题根源分析 #

1. 频率bins语法错误：
   代码中偶数n时的dk定义存在语法错误：-n/2:1:1:n/2-1是无效的MATLAB冒号表达式（冒号仅支持start:step:end格式，此处多了一个冗余的1）。正确写法应为：

   if(mod(n,2)==0) % even split
       dk=2*pi/L * (-n/2 : n/2 - 1);
   else %odd split
       dk=2*pi/L * (-(n-1)/2 : (n-1)/2);
   end
   

   这个错误会导致频率点的数量和取值完全偏离预期，是自定义逆DFT出现明显误差的核心原因。

2. 采样间隔与FFT周期假设不匹配：
   你用linspace(0,L,n)生成x，采样间隔为L/(n-1)，但MATLAB FFT隐含假设信号是周期为L的周期信号，采样间隔应为L/n。若要严格匹配FFT的频率定义，建议将x改为：

   x = linspace(0, L*(1-1/n), n);
   

   这样频率轴dk的取值就能完全对应FFT的频率分量，进一步缩小误差。

3. 浮点运算误差：
   修正频率bins后，自定义逆DFT与MATLAB ifft的差异会缩小到浮点运算精度级别（通常为1e-15量级），这属于正常现象——MATLAB的ifft内部经过优化，数值稳定性更好，但本质上和正确的自定义逆DFT结果一致。

验证效果 #

修正dk的定义后，重新运行代码，自定义逆DFT的结果会与MATLAB ifft几乎完全重合，误差仅为浮点精度级别的噪声。

内容的提问来源于stack exchange，提问作者Mart