Performing split step fourier using MKL dfti

Nadia_A_ — Tue, 16 May 2017 06:11:54 GMT

Hi,

I am trying to perform split step fourier method with the following code. There seems to be a problem with the Fourier portion because when I compare my data with the values obtained with another program, they match perfectly just until I begin the Fourier transform. Here is my code:

program ss

Use MKL_DFTI

implicit none

real*8 , parameter :: distance = 100.0d0
real*8 , parameter ::beta2 = 1.0d0
real*8, parameter :: N = 1.0d0
real*8, parameter :: nt = 1024.0d0
real*8, parameter :: Tmax = 32.0d0
real*8 :: step_num
real*8 :: deltaz
real*8 :: dtau
real*8, dimension (1024) :: tau, omega, uu
real*8, parameter :: pi = 3.141592653590d0
type(DFTI_DESCRIPTOR), POINTER :: handle
complex*16 :: hhz, C= (0,1.0d0)
complex*16, dimension (1024):: dispersion, xy,temp,temp1,temp2
integer :: i, status

step_num = 100000
deltaz = distance/step_num
dtau = (2.0d0*Tmax)/nt

do i = 1,1024
tau (i) = dble(i-513)*dtau

if (i<=512) then
omega (i) = dble(i-1)*pi/Tmax
else
omega (i) = dble(i-1025)*pi/Tmax
end if

end do

uu = 1/dcosh(tau)

dispersion = exp(0.50d0*C*beta2*(omega**2)*deltaz)

hhz = 1*(0,1.0d0)*N**2*deltaz

temp = uu*exp(abs(uu)**2*hhz/2)

status = DftiCreateDescriptor(handle, DFTI_DOUBLE,DFTI_COMPLEX,1,1024)
status = DftiSetValue(handle,DFTI_FORWARD_SCALE,1.0d0)
status = DftiSetValue(handle,DFTI_BACKWARD_SCALE,1/dble(1024))
status = DftiCommitDescriptor(handle)

do j=1,step_num

status = DftiComputeBackward(handle,temp)

do i= 1,1024

temp1(i) = temp(i)*dispersion(i)

end do

Status = DftiComputeForward(handle,temp1)

do i= 1,1024

temp2(i) = temp1(i)*exp(abs(temp1(i))**2*hhz)

end do

end do

status = DftiFreeDescriptor(handle)

xy = temp2*exp(-abs(temp1)**2*hhz/2)

do i=1,1024

open (10, file = 'spec.dat')
write (10,*) tau(i), abs(xy(i)**2)

end do

end program ss

Hi Nadia,

Zhen_Z_Intel — Wed, 17 May 2017 03:03:00 GMT

Hi Nadia,

You compute temp backward for 100000 times. I am afraid the result for 100000 later would be small enough which be seen as 0 value. Thus, temp1 and temp2 also be 0. I tried to print the original temp value, and result for calculating once, the result is fine. And I also found the result after 30 times computing already extend the boundary of smallest double floating number (around 10E-60). Please reduce step_num to see the result, Thanks.

Best regards,
Fiona

topic Performing split step fourier using MKL dfti in Intel® oneAPI Math Kernel Library

Performing split step fourier using MKL dfti

Hi Nadia,