Ошибка сегментации FFTW 3

У меня возникает ошибка сегментации при использовании FFTW с фортраном, и я не смог понять, почему. Не повезло и в онлайн-поиске. Вот мой код:

 integer, parameter :: Nx=128, Ny=Nx, Nz=Nx
 integer, parameter :: NORM=Nx*Ny*Nz,Ntop=MAX(Nx,Ny,Nz)/2, LX1=Nx/2+1
 double precision, parameter :: Etot=0.5
 double complex, dimension(LX1,Ny,Nz) :: fbx, fby, fbz
 double precision, dimension(Nx,Ny,Nz) :: bx, by, bz
 real
 ! ...
 ! Assign values to fbx, fby, fbz 
 ! ...
 call dfftw_plan_dft_c2r_3d(plan,Nx,Ny,Nz,fbx,bx,FFTW_ESTIMATE)
 call dfftw_execute_dft_c2r(plan, fbx, bx)
 call dfftw_execute_dft_c2r(plan, fby, by)
 call dfftw_execute_dft_c2r(plan, fbz, bz)
 call dfftw_destroy_plan(plan)
 ! et = sum(abs(fbx)**2+abs(fby)**2+abs(fbz)**2)
 ! bx=bx*sqrt(Etot/et); by=by*sqrt(Etot/et); bz=bz*sqrt(Etot/et)

 ! Write bx, by, bz to files

Код работает нормально таким образом. Проблема возникает, когда я раскомментирую строки, вычисляющие et и масштабирующие bx. Когда я раскомментирую эти строки, я получаю сообщение об ошибке сегментации, говорящее «недопустимая ссылка на память». Ошибка seg возникает при первом выполнении dfftw_execute_dft_c2r().

Я пробовал как fftw 3.2, так и fftw 3.4 с ifort и gfortran на двух разных машинах (со всеми возможными комбинациями), но безуспешно.

Это сводит меня с ума!! Мне нужна помощь, пожалуйста!! Почему это происходит именно так?

Спасибо!


person toylas    schedule 10.06.2014    source источник
comment
Вы cimoile с проверкой ошибок? Вы пробовали отладчик?   -  person Vladimir F    schedule 10.06.2014
comment
Я скомпилировал с проверкой ошибок, и отладчик сообщает мне только, где ошибка (как я упоминал в своем исходном сообщении).   -  person toylas    schedule 12.06.2014


Ответы (1)


Вы используете один и тот же план для разных массивов. FFTW работает только в том случае, если они выровнены в памяти - см., например. последняя точка здесь, а также соответствующие части в документах почему это трудно гарантировать.

Два предложения, которые я видел в документах или где-либо еще, либо

  1. Создайте временные массивы, создайте план для этих временных массивов и скопируйте данные во временные массивы и из них до/после вызова execute_dft. Я использовал это, чтобы избежать подобных ошибок в прошлом, но это не лучший вариант для производительности. Тем не менее, это простое изменение, которое позволяет легко проверить, является ли это причиной проблемы.
  2. Создайте несколько планов.

Кроме того, похоже, что вы используете ссылку устаревший интерфейс Fortran. Если возможно, лучше использовать "Современный" интерфейс, но мой опыт в прошлом подсказывал, что некоторые компиляторы не очень хорошо с этим справляются.

person Ian    schedule 10.06.2014
comment
Я попробовал оба ваших предложения, и они не сработали. Насколько я понимаю, вы можете создать план один раз и использовать его снова и снова, если вы используете его на массивах именно того размера, для которого вы создали план. Более того, если бы это была проблема с планом, это выдало бы ошибку, несмотря ни на что. Сумасшествие заключается в том, что он отлично работает, если строки et и bx закомментированы, и перестает работать на более раннем этапе (первое выполнение плана), как только я раскомментирую эти строки. - person toylas; 12.06.2014
comment
Это странно. Примером файла, где я использовал его в прошлом, будет этот файл на github, но я не вижу ничего очевидного, чего бы не было в этом предложении. Если FFTW топает по всей памяти, закомментированные строки в вашем коде находятся там, где я ожидаю, что он потерпит неудачу. - person Ian; 13.06.2014