PROGRAM Slow_Fourier

IMPLICIT NONE
  INTEGER, PARAMETER :: nmax=20
!from subroutines
  DOUBLE PRECISION, DIMENSION(nmax) :: u_fct
  COMPLEX(KIND=4), DIMENSION(nmax) :: u_trans
  DOUBLE PRECISION, DIMENSION(nmax) :: u_t_inv
!index
  INTEGER :: n
!::::::::::Main::::::::::!
CALL the_function (u_fct, nmax)

CALL fourier_transform (u_trans, u_fct, nmax)

print*,"======================================================================================"
print*, "  zu Aufgabenteil a)"
print*, "     u_tranformiert"
print*, "      Realt.  Imagin„rt."

DO n=1, nmax
  WRITE(*,FMT="(A2,I3,2(F8.3))"), "n= ", n, u_trans(n)
END DO
print*,"======================================================================================"



CALL inverse_transform (u_t_inv, u_trans, nmax)

print*,"======================================================================================"
print*, "  zu Aufgabenteil d)"
print*, "     u_t_inverse, u_fct"

DO n=1, nmax
  WRITE(*,FMT="(A3,I3,2(F8.3))"), "n=", n, u_t_inv(n), u_fct(n)
END DO

print*,"======================================================================================"
END PROGRAM



SUBROUTINE the_function(u_init, nmax) !!!!!!!!!!!HIER VLL FEHLER WEGEN n=?????
  IMPLICIT NONE
  INTEGER, INTENT(IN ) :: nmax
  DOUBLE PRECISION, DIMENSION(nmax) :: u_init
  DOUBLE PRECISION, PARAMETER :: b1=4.D0
  DOUBLE PRECISION, PARAMETER :: b2=2.D0
  DOUBLE PRECISION, PARAMETER :: A1=8.D0
  DOUBLE PRECISION, PARAMETER :: A2=5.D0
  INTEGER :: k, n
  DOUBLE PRECISION :: PI=ACOS(-1.D0)

    DO n=1,nmax
     u_init(n)=(A1*COS(2*PI*b1*(DBLE(n)-1.d0)/DBLE(nmax))+A2*SIN(2*PI*b2*(DBLE(n)-1.d0)/DBLE(nmax)))
    END DO

END SUBROUTINE the_function


SUBROUTINE fourier_transform (u,u_in, nmax)
  IMPLICIT NONE

  COMPLEX(KIND=4) :: Uk_SUM
  INTEGER, INTENT(IN) :: nmax
  DOUBLE PRECISION, DIMENSION(nmax), intent(in) :: u_in
  DOUBLE PRECISION :: PI=ACOS(-1.D0)
  INTEGER :: k, n
  COMPLEX, DIMENSION(nmax), INTENT (out) :: U

  !cms DO n=1, nmax hier war die Aufgabenstellung irrefuehrend
   DO n=0,nmax-1
    Uk_SUM = CMPLX (0.d0,0.d0)
    DO k=0,nmax-1
  !  u_in=
     Uk_SUM = Uk_SUM + CMPLX(u_in(k+1)*COS(2.d0*PI*DBLE(k)*DBLE(n)/DBLE(nmax)),u_in(k+1)*SIN(2.D0*PI*DBLE(k)*DBLE(n)/DBLE(nmax)))
    END DO
     U(n) = Uk_SUM
  END DO

END SUBROUTINE fourier_transform


SUBROUTINE inverse_transform (inv,u_in, nmax)
  IMPLICIT NONE

  COMPLEX(KIND=4) :: Uk_SUM
  INTEGER, INTENT(IN) :: nmax
  DOUBLE PRECISION :: PI=ACOS(-1.D0)
  INTEGER :: k, n
  COMPLEX(KIND=4), DIMENSION(nmax), INTENT(IN)  :: u_in
  DOUBLE PRECISION, DIMENSION(nmax), intent (out) :: inv

   DO n=1, nmax
    Uk_SUM = CMPLX (0.d0,0.d0)
    DO k=0,nmax-1
     !cms Uk_SUM = Uk_SUM + CMPLX(-u_in(k+1)*cos(2.d0*PI*DBLE(k)*DBLE(n)/DBLE(nmax)),u_in(k+1)*sin(2.D0*PI*DBLE(k)*DBLE(n)/DBLE(nmax))) Achtung: Minus Zeichen, n-1
      Uk_SUM = Uk_SUM + u_in(k)*CMPLX(COS(-2.d0*PI*DBLE(k)*DBLE(n-1)/DBLE(nmax)),SIN(-2.D0*PI*DBLE(k)*DBLE(n-1)/DBLE(nmax)))
    END DO
    inv(n)=(DBLE(Uk_SUM))/DBLE(nmax)
  END DO

END SUBROUTINE inverse_transform