!--------------------------------------------------------------------------------
! Copyright (c) 2016 Peter Grünberg Institut, Forschungszentrum Jülich, Germany
! This file is part of FLEUR and available as free software under the conditions
! of the MIT license as expressed in the LICENSE file in more detail.
!--------------------------------------------------------------------------------

      MODULE m_ifft
      use m_juDFT
      CONTAINS
      INTEGER FUNCTION  ifft235 (iofile,ksfft,n,gmaxp)
C^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
C   this function checks wether n can be expressed as :
C   n =  (2**P) * (3**Q) * (5**R) to match withe the MFFT - routines
C   used in this program. If close to n there is a number n' which
c   is solely expressable as n'= 2**p then this number sis choosen
c   and is outputted by the function. If n is not expressable as
C   n =  (2**P) * (3**Q) * (5**R) a number n' is choosen close to ns
c   is selected which fullfilles these requirements.
c
c
c                            Stefan Bl"ugel , kfa, Oct. 1993
C^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^^
      IMPLICIT NONE
c
c----> declaration part
c
      INTEGER iofile,ksfft,n
      REAL    gmaxp
c
c----> local variabels
c
      INTEGER np,nn,nl,iexp,is,itrial
      REAL rl
      REAL fac_1,two,fac_2,fac_3,one,fac_4
      PARAMETER ( two=2.0e0,fac_1=1.001e0,fac_2=1.7e0,fac_3=1.95e0 )
      PARAMETER ( one=1.0e0,fac_4=0.03e0 )
      INTRINSIC abs,log,mod,real
c
c
      ifft235=0
      IF ( n.eq.0 ) THEN
         WRITE (iofile,*) 'n should not be zero'
          CALL juDFT_error("n should not be zero",calledby="ifft235")
      ENDIF
c
c====> RADIX 2 CALCULATION ONLY
c
      IF ( ksfft .EQ. 0 ) THEN
          rl   = log(real(n))/log(two)
          nl   = rl
          IF ( abs( n - 2 ** nl ) .GT. abs(n - 2 ** (nl+1) ) ) THEN
               np = nl + 1
          ELSE
               np = nl
          ENDIF
          IF ( ( np .eq. nl ) .AND. ( gmaxp * real(nl)/rl .lt. fac_1 ) )
     >         np = nl + 1
          ifft235 = 2 ** np
          RETURN
       ENDIF
c
c====> RADIX 2 , 3, 5 CALCULATION
c
c
c----> since gmaxp is large enough, try also smaller n
c
      IF ( gmaxp .GT. fac_2 .AND. gmaxp .LT. fac_3 ) THEN
           is = -1
      ELSE
           is =  1
      ENDIF
      np   = n
      nn   = n
c
      DO 20 itrial = 1 , 200
c
c----> check whether there is a number close by, which is 2**NL
c      ( FFT very fast )
c
          rl   = log(real(nn))/log(two)
          nl   = rl
          IF (( abs( real(nl)/rl - one ) .LT. fac_4 ).AND.
     +        ( 2**nl.GE.n ) ) THEN
             IF (gmaxp * real(nl)/rl .GT. one) ifft235 = 2 ** nl
             RETURN
          ELSE IF (( abs( real(nl+1)/rl - one ) .lt. fac_4 ).AND.
     +             ( 2**(nl+1).GE.n ) ) THEN
             ifft235 = 2 ** (nl + 1)
             RETURN
          ELSE
c
c----> no , no binary number is arround, check wether number can be
c      divided by 2,3 or 5
c
             DO 10 iexp = 1 , nl+1
                IF (mod(nn,2) .EQ. 0) THEN
                   nn = nn / 2
                ELSE IF (mod(nn,3) .EQ. 0 ) THEN
                   nn = nn / 3
                ELSE IF (mod(nn,5) .EQ. 0 ) THEN
                   nn = nn / 5
                ENDIF
                IF ( nn.eq.2 .OR. nn.eq.3 .OR. nn.eq.5 ) THEN
c---> o.k.
                   ifft235 = np
                   RETURN
                ENDIF
  10         ENDDO
c
c----> no , n  cannot be expressed as (2**P) * (3**Q) * (5**R)
c      change number
             nn = n + (is ** (itrial)) * ((itrial + 1)/2)
             np = n + (is ** (itrial)) * ((itrial + 1)/2)
          ENDIF
  20    ENDDO
c
      WRITE (iofile,'('' to few trials '', i3)') itrial
      STOP
c
      END FUNCTION ifft235
C-----------------------------------------------------------------------
      INTEGER FUNCTION i2357(ii)

!
! simple setup to determine fft-length for ESSL calls
!
      IMPLICIT NONE
      INTEGER, INTENT (IN)  :: ii
      INTEGER i,j,k,h,m,n,nn

      nn = 12582912
      DO m = 0,1 
        DO k = 0,1 
          DO j = 0,1 
            DO i = 0,1 
              DO h = 1,25
                n = 2**h * 3**i * 5**j * 7**k * 11**m
                IF ( (n.GT.ii).AND.(n.LT.nn) ) nn = n
              ENDDO
            ENDDO
          ENDDO
        ENDDO
      ENDDO
      i2357 = nn

      END FUNCTION i2357
 
      END MODULE m_ifft