Showing results for 
Search instead for 
Did you mean: 
New Contributor II

Fast Fourier Stuff


How do you make a complex array allocatable as a double if all the other real elements are the Real(kind=dp)?

Why am I using real(kind=dp) if all it really is is a double?

Or am I overthinking this?

I am using your mkl routines for FFT in the FFT subroutine
  !  Wulf.f90
    !  Wulf - Entry point of console application.

    !  PROGRAM: Wulf
    !  PURPOSE:  Entry point for the console application.

    program Wulf
    use BASE
    use FastFourierTransform
    use S

    implicit none
    integer i,ISTAT,count, dummy
    REAL (KIND=dp) time
    Double Complex :: X(NF)
    REAL (KIND=dp), ALLOCATABLE :: aXYZ(:,:)
    !REAL (KIND=dp) aXYZ(count1, 3)

    open(sm,FILE="m.csv",status = "unknown")
    open(1, file="", STATUS = 'old')

    ALLOCATE (aXYZ(count1,count),STAT=istat)
    IF (istat.NE.0) THEN
        WRITE (*, *) '*** Could not allocate some arrays in LINSOLVE'
    END IF
    time = 0.1d0

    do i =1,NF
        if(i .eq. 1) then
            x(i) = (0.70154d0,0.0d0)
        elseif(i .eq. 2) then
            x(i) =(-2.5018d0,0.0d0)
        elseif(i .eq. 3) then
            x(i) = (-0.35385d0,0.0d0)
        elseif(i .eq. 4) then
            x(i) = (-0.82359d0,0.0d0)
        elseif(i .eq. 5) then
            x(i) = (-1.5771d0,0.0d0)
        elseif(i .eq. 6) then
            x(i) = (0.50797d0,0.0d0)
        elseif(i .eq. 7) then
            x(i) = (0.28198d0,0.0d0)
        elseif(i .eq. 8) then
            x(i) = (0.03348d0, 0.0d0)
        end if
    end do

    ! Variables

    ! Body of Wulf
    print *, 'Hello World'
    call pardiso_unsym(aXYZ,count,1)
    call FFT(sm, NF, time, X)
    end program Wulf


0 Kudos
2 Replies

I am not familiar with the MKL routines. Your question confuses me as I see no complex variables in the code.

I have to assume that "dp" is a PARAMETER constant with a value of 8. (Hopefully from a call to SELECTED_REAL_KIND.) Given that, real(kind-dp) is double precision real, and complex(kind=dp) is double precision complex. Note that this is a case where the nonstandard *n syntax differs from our kind numbers (COMPLEX(8) is COMPLEX*16).

0 Kudos
New Contributor II


mecej4 - showed me how to use the following code about 2 years ago when I was developing the Water Supply Analysis program. I now have it in a base module that sits below everything, makes it nice and easy to control variables.

I presume the 15 is number of places and the 307 is the maximum exponent - never really looked it up.

Module Base

    INTEGER, PARAMETER :: dp = selected_real_kind(15, 307)


The MKL routines are merely fairly standard FFT routines, although the manual is somewhat obscure. Double Complex :: X(NF)is in the program - as shown in the first post - but I could not find the complex(kind=dp) call - so thank you.

I tried it and it worked. 

I was moving from NAG to MKL for the FFT, the FFT is used for the data from a Newmark-Beta output for a vibrating beam. Unfortunately MKL does not have a Newmark Beta routine - which is not a big problem. I used an algorithm from Gavin at Duke, he has the best notes I came across - clear concise and set out well. Newmark's original paper is not that clear.


0 Kudos