在没有内存泄漏的情况下为 Fortran FFT (FFTW) 使用对齐内存
Using aligned memory for Fortran FFTs (FFTW) without memory leaks
我想使用 FFTW 的现代 Fortran 界面,但使用的方式允许简单的函数调用,例如 ifftshift(fft_c2c(vec)*exp(vec)) 等等。这是我对如何做到这一点的理解(我也明白每次调用都做一个新计划并不是最有效的事情)。目前这段代码是有效的(returns 正确的结果);但是,存在内存泄漏,导致重复调用导致丢失。我不太确定在哪里!我曾希望 return 变量 `fft' 与唯一未释放内存的关联不会导致泄漏,但这显然不是真的。我错过了什么,我怎样才能更好地构建我想用适当的现代 fortran 做什么?谢谢!
function fft_c2c(x) result(fft)
integer :: N
type(C_PTR) :: plan
complex(C_DOUBLE_COMPLEX), pointer :: fft(:)
complex(C_DOUBLE_COMPLEX), dimension(:), intent(in) :: x
! Use an auxiliary array that is allocated with fftw_alloc_complex
! to ensure memory alignment for performance, see FFTW docs
complex(C_DOUBLE_COMPLEX), pointer :: x_align(:)
type(C_PTR) :: p
N = size(x)
p = fftw_alloc_complex(int(N, C_SIZE_T))
call c_f_pointer(p, fft, [N]);
p = fftw_alloc_complex(int(N, C_SIZE_T))
call c_f_pointer(p, x_align, [N]);
plan = fftw_plan_dft_1d(N, x_align, fft, FFTW_FORWARD, FFTW_MEASURE);
! FFTW overwrites x_align and fft during planning process, so assign
! data here
x_align = x
call fftw_execute_dft(plan, x_align, fft);
call fftw_free(p);
end function fft_c2c
你不能轻易做到这一点。您在 Fortran 上强制使用 "modern"="everything is a function",这里它不太适合(或根本不适合)。
对于内存泄漏,规则很简单 - 释放所有指针。将它们用于结果变量是内存泄漏的保证。如果需要local allocated aligned memory,需要在本地分配,把数据拷贝到那里,拷贝出来再deallocate。
Fortran 中的每个指针都需要显式释放,没有引用计数或垃圾收集来为您释放它们。
你考虑只使用带有适当标志的非对齐内存并测量差异,反正你似乎并不关心最高性能。
最后,在每次转换之前做FFTW_MEASURE
不仅仅是 "not the most efficient thing",它是 绝对的性能灾难。您至少应该使用 FFTW_ESTIMATE
来缓解它。
我想使用 FFTW 的现代 Fortran 界面,但使用的方式允许简单的函数调用,例如 ifftshift(fft_c2c(vec)*exp(vec)) 等等。这是我对如何做到这一点的理解(我也明白每次调用都做一个新计划并不是最有效的事情)。目前这段代码是有效的(returns 正确的结果);但是,存在内存泄漏,导致重复调用导致丢失。我不太确定在哪里!我曾希望 return 变量 `fft' 与唯一未释放内存的关联不会导致泄漏,但这显然不是真的。我错过了什么,我怎样才能更好地构建我想用适当的现代 fortran 做什么?谢谢!
function fft_c2c(x) result(fft)
integer :: N
type(C_PTR) :: plan
complex(C_DOUBLE_COMPLEX), pointer :: fft(:)
complex(C_DOUBLE_COMPLEX), dimension(:), intent(in) :: x
! Use an auxiliary array that is allocated with fftw_alloc_complex
! to ensure memory alignment for performance, see FFTW docs
complex(C_DOUBLE_COMPLEX), pointer :: x_align(:)
type(C_PTR) :: p
N = size(x)
p = fftw_alloc_complex(int(N, C_SIZE_T))
call c_f_pointer(p, fft, [N]);
p = fftw_alloc_complex(int(N, C_SIZE_T))
call c_f_pointer(p, x_align, [N]);
plan = fftw_plan_dft_1d(N, x_align, fft, FFTW_FORWARD, FFTW_MEASURE);
! FFTW overwrites x_align and fft during planning process, so assign
! data here
x_align = x
call fftw_execute_dft(plan, x_align, fft);
call fftw_free(p);
end function fft_c2c
你不能轻易做到这一点。您在 Fortran 上强制使用 "modern"="everything is a function",这里它不太适合(或根本不适合)。
对于内存泄漏,规则很简单 - 释放所有指针。将它们用于结果变量是内存泄漏的保证。如果需要local allocated aligned memory,需要在本地分配,把数据拷贝到那里,拷贝出来再deallocate。
Fortran 中的每个指针都需要显式释放,没有引用计数或垃圾收集来为您释放它们。
你考虑只使用带有适当标志的非对齐内存并测量差异,反正你似乎并不关心最高性能。
最后,在每次转换之前做FFTW_MEASURE
不仅仅是 "not the most efficient thing",它是 绝对的性能灾难。您至少应该使用 FFTW_ESTIMATE
来缓解它。