如何修复 FFT 微分程序中的偏移和缩放错误？
How to fix shift and scaling errors in an FFT-powered differentiation program?
所以我正在尝试，就像我之前的许多其他人一样，筛选通过傅里叶变换工作进行微分所必需的相移和归一化系数的所有天体排列。
我正在尝试使用尽可能少的代码，完全使用 numpy 功能来保持代码干净，同时使用其数组操作，据我所知，这比显式 python 循环更快。
来自 函数关系的 table 的第 106 行，一维，在关于傅立叶变换的维基百科文章的 Important Fourier Transforms 部分，我知道傅立叶变换有这个 属性，其中实数微分 space 等于频率乘法 space 使得 d^nf/dx^n = ifft[F*(i*w)^n]，其中 F 是傅立叶f.
的变换
寻找更可靠的参考来定义 w 在 d^nf/dx^n = ifft[F*(i*w)^n] 中的含义，我找到了一篇关于使用 DFT 进行微分的数学的 PDF 论文。在论文中，作者讨论了给定函数有无限多个所谓的三角插值，并展示了一种获得 "less wiggly" 解的方法，这是唯一的。
所以我使用了他的定义并编写了一小段 python 代码来进行微分。我想出的代码是这个（注意最后两行，因为它们会改变）
N = 51
a = 0
b = 2*np.pi
X = np.linspace(a,b, N)
L = abs(b-a)
TwoPiOverL = 2.0*np.pi/L
freqs = np.fft.fftfreq(N,1./N)
# THESE COME FROM [1] #
D1 = TwoPiOverL*freqs*1j
D2 = -(TwoPiOverL*freqs)**2
#=========================#
S = np.sin(X)
fftS=np.fft.fft(S)
C = np.cos(X)
# COMPUTING THE DERIVATIVES #
Y1 = np.fft.ifft(D1*fftS).real
Y2 = np.fft.ifft(D2*fftS).real
在生成的图形中，符号~d/dx sin表示根据傅里叶变换乘法计算的正弦函数的导数，而不是解析表达式。
测试 1
第一次测试，仅使用 PDF 中的定义：
所以第一个测试是一团糟。导数在域的中间很好地插值，但端点很疯狂。在 PDF 本身的第 6 节开头，作者说：
All of the above applies to spectral differentiation of periodic functions y(x). If the function is non-periodic,
artifacts will appear in the derivatives due to the implicit discontinuities at the endpoints.
听起来不错，这意味着我在某处搞砸了，因为如果整篇论文都在讨论周期信号，而如果我的信号恰好落在非周期 class 中，那么如果我可以使用周期性信号，还有改进的余地。
但我正在使用 DFT，这意味着我正在对数组应用数值过程。不包含有关周期性信息的离散数据。于是我想：
那么，我如何告诉我的傅立叶微分器信号是周期性的？在代码中这意味着什么？在代码中，我必须更改哪里才能将该信息传达给我的程序？
Whosebug 里好像有答案。在 this post 关于使用 FFT 的数值微分（是的，多余的，但我会到达那里）中，@dkato 给出了我解释为更改代码中最后两行的说明
Y1 = np.fft.ifft(D1*fftS).real
Y2 = np.fft.ifft(D2*fftS).real
至
Y1 = np.fft.ifft(np.exp(D1)*fftS).real
Y2 = np.fft.ifft(np.exp(D2)*fftS).real
所以我这样做了，然后测试了我会得到什么。
测试 2
第二次测试，使用@dkato 更正：
从上图可以看出，通过@dkato 提供的校正，端点摆动消失了，但导数的估计仍然很远。我正在使用 51 个数据点，所以我猜测这不是由于扩展不准确造成的。看来我计算导数的方式不正确。
分析一阶导数与其傅里叶变换估计之间精度损失最大的来源似乎是向右移动，而二阶导数不仅移动而且缩放不当。然而，有趣的是，当我们使用
Y2 = np.fft.ifft(-np.exp(D2)*fftS).real
计算二阶导数，我们看到 all of the error due to shifting goes away, leaving only the scaling。我不知道这是 属性 的一般 FFT 驱动的二阶导数，还是仅仅是因为我使用正弦函数作为我的 y(x)，特别是。
最后的重点
我很高兴端点处奇怪、丑陋的摆动消失了，这意味着我必须以某种方式告诉我的傅立叶微分器我正在处理周期信号。但是我的差异估计的这些奇怪的移动和缩放仍然很麻烦，而且我终究无法指出哪里出了问题。正如您可能已经猜到的那样，部分原因是我对 DFT 中所有运动部件的操作经验不多。
以这种身份，我有两个问题要问社区，如果你愿意帮助我解决这个问题。
1) 如果我确实以某种方式告诉我的傅立叶微分器我正在使用周期函数，那么我究竟是怎么做到的？有效代码和无效代码背后的原因是什么？
和
2) 解析导数和估计导数之间这些奇怪的 scale/shift 不一致的根源是什么？我该如何修复它们？
先回答问题的症结：
What is the source of these weird scale/shift incongruities between the analytical and estimated derivatives? How can I fix them?
问题从以下行开始：
X = np.linspace(a,b, N)
默认情况下 numpy.linspace 包括端点，这会在计算的周期函数扩展中引入不连续性，因为端点是重复的。
为了说明这可能是一个问题，您可以查看 N=4 的计算函数 sin(X):
sin(0)
sin(2*np.pi/3)
sin(4*np.pi/3)
sin(2*np.pi)
其相位在每个后续值处递增 2*np.pi/3。重复这些值以获得周期为 N=4 的周期性扩展，将产生
sin(0)
sin(2*np.pi/3)
sin(4*np.pi/3)
sin(2*np.pi)
sin(0) # == sin(2*np.pi)
sin(2*np.pi/3) # == sin(2*np.pi + 2*np.pi/3)
sin(4*np.pi/3) # == sin(2*np.pi + 4*np.pi/3)
sin(2*np.pi) # == sin(2*np.pi + 2*np.pi)
你可以看到相位在每一步增加 2*np.pi/3 除了在块边界处它保持不变。这种小的不连续性使得利用三角插值基础获得良好的近似变得更加困难。由此产生的插值包含显着的振荡，推导运算进一步放大了振荡。
要解决此问题，您应该排除端点（否则保持初始解决方案不变）：
X = np.linspace(a,b, N, endpoint=False)
就@dkato's post中的解决方案而言，它看起来只是实现了一个时移操作。对于一个看起来还不错，但与微分无关的正弦信号。
最后一点，回答部分
If it is true that I have somehow told my fourier differentiator that I am using a periodic function, how exactly did I do it?
是的，这是真的，您在开始使用离散傅立叶变换（及其 FFT 实现）的那一刻就做到了。
所以我正在尝试，就像我之前的许多其他人一样，筛选通过傅里叶变换工作进行微分所必需的相移和归一化系数的所有天体排列。
我正在尝试使用尽可能少的代码，完全使用 numpy 功能来保持代码干净，同时使用其数组操作，据我所知，这比显式 python 循环更快。
来自 函数关系的 table 的第 106 行，一维，在关于傅立叶变换的维基百科文章的 Important Fourier Transforms 部分，我知道傅立叶变换有这个 属性，其中实数微分 space 等于频率乘法 space 使得 d^nf/dx^n = ifft[F*(i*w)^n]，其中 F 是傅立叶f.
寻找更可靠的参考来定义 w 在 d^nf/dx^n = ifft[F*(i*w)^n] 中的含义，我找到了一篇关于使用 DFT 进行微分的数学的 PDF 论文。在论文中，作者讨论了给定函数有无限多个所谓的三角插值，并展示了一种获得 "less wiggly" 解的方法，这是唯一的。
所以我使用了他的定义并编写了一小段 python 代码来进行微分。我想出的代码是这个（注意最后两行，因为它们会改变）
N = 51
a = 0
b = 2*np.pi
X = np.linspace(a,b, N)
L = abs(b-a)
TwoPiOverL = 2.0*np.pi/L
freqs = np.fft.fftfreq(N,1./N)
# THESE COME FROM [1] #
D1 = TwoPiOverL*freqs*1j
D2 = -(TwoPiOverL*freqs)**2
#=========================#
S = np.sin(X)
fftS=np.fft.fft(S)
C = np.cos(X)
# COMPUTING THE DERIVATIVES #
Y1 = np.fft.ifft(D1*fftS).real
Y2 = np.fft.ifft(D2*fftS).real
在生成的图形中，符号~d/dx sin表示根据傅里叶变换乘法计算的正弦函数的导数，而不是解析表达式。
测试 1
第一次测试，仅使用 PDF 中的定义：
所以第一个测试是一团糟。导数在域的中间很好地插值，但端点很疯狂。在 PDF 本身的第 6 节开头，作者说：
All of the above applies to spectral differentiation of periodic functions y(x). If the function is non-periodic, artifacts will appear in the derivatives due to the implicit discontinuities at the endpoints.
听起来不错，这意味着我在某处搞砸了，因为如果整篇论文都在讨论周期信号，而如果我的信号恰好落在非周期 class 中，那么如果我可以使用周期性信号，还有改进的余地。
但我正在使用 DFT，这意味着我正在对数组应用数值过程。不包含有关周期性信息的离散数据。于是我想：
那么，我如何告诉我的傅立叶微分器信号是周期性的？在代码中这意味着什么？在代码中，我必须更改哪里才能将该信息传达给我的程序？
Whosebug 里好像有答案。在 this post 关于使用 FFT 的数值微分（是的，多余的，但我会到达那里）中，@dkato 给出了我解释为更改代码中最后两行的说明
Y1 = np.fft.ifft(D1*fftS).real
Y2 = np.fft.ifft(D2*fftS).real
至
Y1 = np.fft.ifft(np.exp(D1)*fftS).real
Y2 = np.fft.ifft(np.exp(D2)*fftS).real
所以我这样做了，然后测试了我会得到什么。
测试 2
第二次测试，使用@dkato 更正：
从上图可以看出，通过@dkato 提供的校正，端点摆动消失了，但导数的估计仍然很远。我正在使用 51 个数据点，所以我猜测这不是由于扩展不准确造成的。看来我计算导数的方式不正确。
分析一阶导数与其傅里叶变换估计之间精度损失最大的来源似乎是向右移动，而二阶导数不仅移动而且缩放不当。然而，有趣的是，当我们使用
Y2 = np.fft.ifft(-np.exp(D2)*fftS).real
计算二阶导数，我们看到 all of the error due to shifting goes away, leaving only the scaling。我不知道这是 属性 的一般 FFT 驱动的二阶导数，还是仅仅是因为我使用正弦函数作为我的 y(x)，特别是。
最后的重点
我很高兴端点处奇怪、丑陋的摆动消失了，这意味着我必须以某种方式告诉我的傅立叶微分器我正在处理周期信号。但是我的差异估计的这些奇怪的移动和缩放仍然很麻烦，而且我终究无法指出哪里出了问题。正如您可能已经猜到的那样，部分原因是我对 DFT 中所有运动部件的操作经验不多。
以这种身份，我有两个问题要问社区，如果你愿意帮助我解决这个问题。
1) 如果我确实以某种方式告诉我的傅立叶微分器我正在使用周期函数，那么我究竟是怎么做到的？有效代码和无效代码背后的原因是什么？
和
2) 解析导数和估计导数之间这些奇怪的 scale/shift 不一致的根源是什么？我该如何修复它们？
先回答问题的症结：
What is the source of these weird scale/shift incongruities between the analytical and estimated derivatives? How can I fix them?
问题从以下行开始：
X = np.linspace(a,b, N)
默认情况下 numpy.linspace 包括端点，这会在计算的周期函数扩展中引入不连续性，因为端点是重复的。
为了说明这可能是一个问题，您可以查看 N=4 的计算函数 sin(X):
sin(0)
sin(2*np.pi/3)
sin(4*np.pi/3)
sin(2*np.pi)
其相位在每个后续值处递增 2*np.pi/3。重复这些值以获得周期为 N=4 的周期性扩展，将产生
sin(0)
sin(2*np.pi/3)
sin(4*np.pi/3)
sin(2*np.pi)
sin(0) # == sin(2*np.pi)
sin(2*np.pi/3) # == sin(2*np.pi + 2*np.pi/3)
sin(4*np.pi/3) # == sin(2*np.pi + 4*np.pi/3)
sin(2*np.pi) # == sin(2*np.pi + 2*np.pi)
你可以看到相位在每一步增加 2*np.pi/3 除了在块边界处它保持不变。这种小的不连续性使得利用三角插值基础获得良好的近似变得更加困难。由此产生的插值包含显着的振荡，推导运算进一步放大了振荡。
要解决此问题，您应该排除端点（否则保持初始解决方案不变）：
X = np.linspace(a,b, N, endpoint=False)
就@dkato's post中的解决方案而言，它看起来只是实现了一个时移操作。对于一个看起来还不错，但与微分无关的正弦信号。
最后一点，回答部分
If it is true that I have somehow told my fourier differentiator that I am using a periodic function, how exactly did I do it?
是的，这是真的，您在开始使用离散傅立叶变换（及其 FFT 实现）的那一刻就做到了。