PROGRAM Slow_Fourier

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

CALL fourier_transform (u_trans, u_fct)

print*, "zu Aufgabenteil a)"
print*, "u_tranformiert"

DO n=1, nmax
  print*, "n=", n, u_trans(n)
END DO

CALL inverse_transform (u_t_inv, u_trans)

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

DO n=1, nmax
  print*, "n=", n, u_t_inv(n)
END DO

END PROGRAM



SUBROUTINE the_function(u_init) !!!!!!!!!!!HIER VLL FEHLER WEGEN n=?????
  IMPLICIT NONE
  INTEGER, PARAMETER :: nmax=20
  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(k)/DBLE(nmax))+A2*sin(2*PI*b2*DBLE(k)/DBLE(nmax)))
    END DO

END SUBROUTINE the_function


SUBROUTINE fourier_transform (u,u_in)
  IMPLICIT NONE

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

  DO n=1, nmax
    Uk_SUM = CMPLX (0.d0,0.d0)
    DO k=0,nmax-1
  !  u_in=
     Uk_SUM = Uk_SUM + CMPLX(u_in(n)*cos(2.d0*PI*DBLE(k)*DBLE(n)/DBLE(nmax)),u_in(n)*sin(2*PI*DBLE(k)*DBLE(n)/DBLE(nmax)))
    END DO
    u(n)=(DBLE(Uk_SUM))
  END DO

END SUBROUTINE fourier_transform


SUBROUTINE inverse_transform (inv,u_in)
  IMPLICIT NONE

  COMPLEX(KIND=4) :: Uk_SUM
  INTEGER, PARAMETER :: nmax=20
  DOUBLE PRECISION :: PI=ACOS(-1.D0)
  INTEGER :: k, n
  DOUBLE PRECISION, 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
      Uk_SUM = Uk_SUM + CMPLX(u_in(n)*cos(2.d0*PI*DBLE(k)*DBLE(n)/DBLE(nmax)),u_in(n)*sin(-2*PI*DBLE(k)*DBLE(n)/DBLE(nmax)))
    END DO
    inv(n)=(DBLE(Uk_SUM))/dble(nmax)
  END DO

END SUBROUTINE inverse_transform