Unwanted sine modes appearing in convolution using FFT
还没有人认领这个 Issue。
评估
调研方向
首先运行 issue 中的 Julia MWE,重点检查 fft 频率构造、weight_hat 的计算以及 real(ifft(I_hat)) 的结果。将 FFT 卷积与手动计算的积分进行比较,并确定剩余振荡是来自卷积设置还是 FastTransforms.jl;记录其原因以及可复现的修正方法或限制。
由索引模型根据 Issue 内容生成。
描述
Hi All, I don't know if this is an issue with FastTransforms.jl specifically or something that I've done, but any guidance that could be provide will be appreciated.
To give a little bit of background, I'm trying to perform the following convolution (nV, r and r' are vectors):
n_\mathbf{V}(r) = \int d\mathbf{r}' \rho (r') \delta(|r-r'|-R) \frac{r-r'}{|r-r'|}
Where rho(r) is a profile which looks like (I've made sure it's periodic):
I have obtained an analytical expression for the fourier transform of the last two terms in the integral as:
\hat{w} = 4\pi i \frac{k}{|k|^3} (\sin(|k|R)-|k|R\cos(|k|R))
If I use just regular Float64, the resulting nV profile looks fine:
Except when I zoom in:
While these values are small, in a later part of my code, I need to divide this profile by a small number which makes this 'noise' blow up and affects calculations downstream.
I tried using FastTransforms.jl since it would let me go to higher precision. However, even when I use BigFloat, although it removes the noise, some unphysical oscillations remain:
Im no longer certain that these remaining oscillations are due to floating point errors. I also call these oscillations unphysical as, aside from the nature of this convolution integral, if I perform this convolution integral manually, the oscillations disappear completely (and I'm able to obtain the integral to higher precision).
Am I doing something wrong in the way I'm performing the convolution integral? Thanks!
I've attached an MWE below:
using FFTW, FastTransforms
# Generating the density profile
tanh_prof(x,start,stop,shift,coef) = 1/2*(start-stop)*tanh((x-shift)*coef)+1/2*(start+stop)
ub = 10
lb = -10
z = LinRange(lb,ub,4000)
rho1 = 1e3
rho2 = 1e-12
rho = @. tanh_prof(z,rho1,rho2,(ub/4+3*lb/4),20)*(z<=0) +
tanh_prof(z,rho2,rho1,(3*ub/4+lb/4),20)*(z>0)
# Obtaining the fourier transform of w_hat
freq = fftfreq(4000,4000/20)
R = 1/2*2pi
weight_hat = @. 0.0 - 4π*im*freq/abs(freq)^3*(sin(abs(freq)*R)-R*abs(freq)*cos(abs(freq)*R)) *(freq != 0.0)
# Performing the convolution integral
rho_hat = fft(rho)
I_hat = @. rho_hat*weight_hat
nV = real(ifft(I_hat))
- 主要语言
- Julia
- 星标
- 282
- 派生
- 27
- 平均合并
- 53 分钟
- 30 天内合并 PR
- 1
环境准备
这个项目没有提供开发容器、Dockerfile 或贡献指南,环境需要你自己搭建:先看它的 README,通用步骤见我们的新手贡献指南。
从这里开始
- 先读完整个 Issue,再读项目的贡献指南。
- 在 Issue 下留言说明你要接手 —— 这能避免两个人做同样的事。
- Fork 仓库,在一个分支上完成修改。
- 提交 Pull Request,并在描述里引用这个 Issue 编号。
JuliaApproximation/FastTransforms.jl 的其他 Issue
-
Loading FastTransforms.jl can make FFTs via FFTW.jl 100x slower due to threading conflicts可能已有人在做 @dlfivefifty 于 91 天前认领。 未关闭
JuliaApproximation/FastTransforms.jl#267 · 2 个 reaction · 已指派 2 人 ·
-
难度 3/5 1-2 天 新手友好度 38/100
JuliaApproximation/FastTransforms.jl#266 · 1 条评论 ·
-
难度 4/5 3-5 天 新手友好度 42/100
JuliaApproximation/FastTransforms.jl#263 · 2 条评论 ·
-
Allocating lmul!未关闭
难度 4/5 3-5 天 新手友好度 25/100
JuliaApproximation/FastTransforms.jl#253 · 7 条评论 ·
-
难度 2/5 1-3 小时 新手友好度 45/100
查看 JuliaApproximation/FastTransforms.jl 的全部 Issue
相似的 Issue
-
documentation
难度 2/5 1-3 小时 新手友好度 62/100
ohno/Antique.jl#165 ·
-
难度 2/5 1-3 小时 新手友好度 62/100
-
难度 2/5 1-3 小时 新手友好度 68/100
grame-cncm/faust#1344 · 1 条评论 ·
维护者通常 1 天内回复
-
难度 2/5 1-3 小时 新手友好度 70/100
SciML/DiffEqNoiseProcess.jl#342 ·
-
难度 1/5 1 小时以内 新手友好度 88/100