Intel® oneAPI Math Kernel Library
Ask questions and share information with other developers who use Intel® Math Kernel Library.
6834 Discussions

## Performing split step fourier using MKL dfti Beginner
212 Views

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  