FT Report

FT
Exec Time (%)	Codelet Name	Error (%)
5.7	appft_appft__40	5.9
68.0	appft_appft__1	1.3
28.6	fft3d_fftxyz__152	0.1
70.3	fft3d_swarztrauber__27	4.2
5.1	auxfnct_evolve__167	0.6
28.6	fft3d_fftxyz__112	1.4
70.3	fft3d_swarztrauber__27	4.2
27.4	fft3d_fftxyz__137	0.0
70.3	fft3d_swarztrauber__27	4.2
2.2	auxfnct_compute_initial_conditions__148	25.2

Exec Time : 1.257131e+10 cycles

Matching : 95.4%

                c---------------------------------------------------------------------   
      subroutine Swarztrauber(is,m,vlen,n,x,xd1,exponent)

      implicit none
      include 'global.h'
c---------------------------------------------------------------------
c   Computes NY N-point complex-to-complex FFTs of X using an algorithm due
c   to Swarztrauber.  X is both the input and the output array, while Y is a 
c   scratch array.  It is assumed that N = 2^M.  Before calling 
c   Swarztrauber to 
c   perform FFTs
c---------------------------------------------------------------------
      integer is,m,vlen,n,xd1
      double complex x(xd1,n), exponent(n)      

      integer i,j,l
      double complex u1,x11,x21
      integer k, n1,li,lj,lk,ku,i11,i12,i21,i22
      
      if (timers_enabled) call timer_start(4)
c---------------------------------------------------------------------
c   Perform one variant of the Stockham FFT.
c---------------------------------------------------------------------
      n1 = n / 2
      lj = 1
      li = 2 ** m
      do l = 1, m, 2
        lk = lj
        lj = 2 * lk
	li = li / 2
        ku = li + 1

        do i = 0, li - 1
          i11 = i * lk + 1
          i12 = i11 + n1
          i21 = i * lj + 1
          i22 = i21 + lk
        
          if (is .ge. 1) then
            u1 = exponent(ku+i)
          else
            u1 = dconjg (exponent(ku+i))
          endif
          do k = 0, lk - 1
            do j = 1, vlen
              x11 = x(j,i11+k)
              x21 = x(j,i12+k)
              scr(j,i21+k) = x11 + x21
              scr(j,i22+k) = u1 * (x11 - x21)
            end do
          end do
        end do

        if (l .eq. m) then
          do k = 1, n
            do j = 1, vlen
              x(j,k) = scr(j,k)
            enddo
          enddo
        else
          lk = lj
          lj = 2 * lk
          li = li / 2
          ku = li + 1

          do i = 0, li - 1
            i11 = i * lk + 1
            i12 = i11 + n1
            i21 = i * lj + 1
            i22 = i21 + lk
        
            if (is .ge. 1) then
              u1 = exponent(ku+i)
            else
              u1 = dconjg (exponent(ku+i))
            endif
            do k = 0, lk - 1
              do j = 1, vlen
                x11 = scr(j,i11+k)
                x21 = scr(j,i12+k)
                x(j,i21+k) = x11 + x21
                x(j,i22+k) = u1 * (x11 - x21)
              end do
            end do
          end do
        endif
      end do
      if (timers_enabled) call timer_stop(4)

      return
      end

c---------------------------------------------------------------------
      subroutine fftXYZ(sign,x,xout,exp1,exp2,exp3,n1,n2,n3)
        include 'global.h'
        integer sign,n1,n2,n3
        double complex x(n1+1,n2,n3)
        double complex xout((n1+1)*n2*n3)
        double complex exp1(n1), exp2(n2), exp3(n3)
        integer i, j, k, log
        integer bls,ble
	integer len
        integer blkp

        if (timers_enabled) call timer_start(3)	   

        fftblock=CacheSize/n1
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n1)
        if (timers_enabled) call timer_start(7) 
        do k = 1, n3
          do bls = 1, n2, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n2) ble = n2
	    len=ble-bls+1
            do j = bls, ble
            do i = 1, n1
              plane(j-bls+1+blkp*(i-1)) = x(i,j,k)
            end do
            end do
            call Swarztrauber(sign,log,len,n1,plane,blkp,exp1)     
            do j = bls, ble
            do i = 1, n1
              x(i,j,k)=plane(j-bls+1+blkp*(i-1))
            end do
            end do
          end do
        end do
        if (timers_enabled) call timer_stop(7)

        fftblock=CacheSize/n2
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n2 )
        if (timers_enabled) call timer_start(8)	   
        do k = 1, n3
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            call Swarztrauber(sign,log,len,n2,x(bls,1,k),n1+1,exp2)
	  enddo
        end do
        if (timers_enabled) call timer_stop(8)

        fftblock=CacheSize/n3
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
       	log = ilog2(n3)
        if (timers_enabled) call timer_start(9)
        do k = 1, n2
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            do i = 1,n3
            do j = bls, ble
              plane(j-bls+1+blkp*(i-1)) = x(j,k,i)
            end do
            end do
            call Swarztrauber(sign,log,len,n3,plane,blkp,exp3)
            do i = 0,n3-1
            do j = bls, ble
              xout(j+(n1+1)*(k-1+n2*i)) = plane(j-bls+1+blkp*i)
            end do		   
            end do
          end do
        end do  	       
        if (timers_enabled) call timer_stop(9)
        if (timers_enabled) call timer_stop(3)
        return
      end

invivofft3d_fftxyz__112

Real time: 3.605644e+09 -- Predicted time: 3.657419e+09

Callcount: 8

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	2.00380837394855	4.510172e+08	4.589938e+08	1.7
2	### 2	5.98581043570881	4.600324e+08	4.498596e+08	2.2

                c---------------------------------------------------------------------
c compute the roots-of-unity array that will be used for subsequent FFTs. 
c---------------------------------------------------------------------
      subroutine CompExp (n, exponent)

      implicit none
      integer n
      double complex exponent(n) 
      integer ilog2
      external ilog2      
     
      integer m,nu,ku,i,j,ln
      double precision t, ti, pi 
      data pi /3.141592653589793238d0/

      nu = n
      m = ilog2(n)
      exponent(1) = m
      ku = 2
      ln = 1
      do j = 1, m
         t = pi / ln
         do i = 0, ln - 1
            ti = i * t
            exponent(i+ku) = dcmplx(cos(ti),sin(ti))
         enddo        
         ku = ku + ln
         ln = 2 * ln
      enddo
            
      return
      end

c---------------------------------------------------------------------
c---------------------------------------------------------------------
c---------------------------------------------------------------------
      integer function ilog2(n)
      implicit none
      integer n
c---------------------------------------------------------------------
c---------------------------------------------------------------------  
      integer nn, lg
      if (n .eq. 1) then
         ilog2=0
         return
      endif
      lg = 1
      nn = 2
      do while (nn .lt. n)
         nn = nn*2
         lg = lg+1
      end do
      ilog2 = lg
      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------

      subroutine ipow46(a, exponent, result)
c---------------------------------------------------------------------
c compute a^exponent mod 2^46
c---------------------------------------------------------------------

      implicit none
      double precision a, result, dummy, q, r
      integer exponent, n, n2
      external randlc
      double precision randlc
c---------------------------------------------------------------------
c Use
c   a^n = a^(n/2)*a^(n/2) if n even else
c   a^n = a*a^(n-1)       if n odd
c---------------------------------------------------------------------
      result = 1
      if (exponent .eq. 0) return
      q = a
      r = 1
      n = exponent

      do while (n .gt. 1)
         n2 = n/2
         if (n2 * 2 .eq. n) then
            dummy = randlc(q, q) 
            n = n2
         else
            dummy = randlc(r, q)
            n = n-1
         endif
      end do
      dummy = randlc(r, q)
      result = r
      return
      end
c---------------------------------------------------------------------
      subroutine CalculateChecksum(csum,iterN,u,d1,d2,d3)
        implicit none
        integer iterN
        integer d1,d2,d3
        double complex csum
        double complex u(d1+1,d2,d3)
        integer i, i1, ii, ji, ki
        csum = dcmplx (0.0, 0.0)
        do i = 1, 1024
          i1 = i
          ii = mod (i1, d1) + 1
          ji = mod (3 * i1, d2) + 1
          ki = mod (5 * i1, d3) + 1
          csum = csum + u(ii,ji,ki)
        end do
        csum = csum/dble(d1*d2*d3)
        write(*,30) iterN, csum
 30     format (' T =',I5,5X,'Checksum =',1P2D22.12)
      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------
      subroutine compute_initial_conditions(u0,d1,d2,d3)

      implicit none
      include 'npbparams.h'
      integer d1,d2,d3
      double complex u0(d1+1,d2,d3), tmp(maxdim)
      double precision x0, start, an, dummy
      double precision RanStarts(maxdim)

      integer i,j,k
      double precision seed, a
      parameter (seed = 314159265.d0, a = 1220703125.d0)
      external randlc
      double precision randlc
      
      start = seed                                    
c---------------------------------------------------------------------
c Jump to the starting element for our first plane.
c---------------------------------------------------------------------
      call ipow46(a, 0, an)
      dummy = randlc(start, an)
      call ipow46(a, 2*d1*d2, an)
c---------------------------------------------------------------------
c Go through by z planes filling in one square at a time.
c---------------------------------------------------------------------
      RanStarts(1) = start
      do k = 2, d3 
         dummy = randlc(start, an)
         RanStarts(k) = start
      end do
      
      do k = 1, d3 
         x0 = RanStarts(k)
         do j = 1, d2 
           call vranlc(2*d1, x0, a, tmp)
           do i = 1, d1 
             u0(i,j,k)=tmp(i)
           end do
         end do
      end do

      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------
      subroutine evolve(x,y,twiddle,nx,ny,nz)
      double complex x(nx+1,ny,nz),y(nx+1,ny,nz)
      real*8 twiddle(nx+1,ny,nz)
      integer nx,ny,nz
      integer i,j,k
           do i = 1, nz
             do k = 1, ny
               do j = 1, nx
                   y(j,k,i)=y(j,k,i)*twiddle(j,k,i)
                   x(j,k,i)=y(j,k,i)
                 end do
              end do
           end do
      
      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------

invivoauxfnct_evolve__167

Real time: 6.506289e+08 -- Predicted time: 6.546663e+08

Callcount: 6

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	6.01663993807673	1.088093e+08	1.092862e+08	0.4

                c---------------------------------------------------------------------   
      subroutine Swarztrauber(is,m,vlen,n,x,xd1,exponent)

      implicit none
      include 'global.h'
c---------------------------------------------------------------------
c   Computes NY N-point complex-to-complex FFTs of X using an algorithm due
c   to Swarztrauber.  X is both the input and the output array, while Y is a 
c   scratch array.  It is assumed that N = 2^M.  Before calling 
c   Swarztrauber to 
c   perform FFTs
c---------------------------------------------------------------------
      integer is,m,vlen,n,xd1
      double complex x(xd1,n), exponent(n)      

      integer i,j,l
      double complex u1,x11,x21
      integer k, n1,li,lj,lk,ku,i11,i12,i21,i22
      
      if (timers_enabled) call timer_start(4)
c---------------------------------------------------------------------
c   Perform one variant of the Stockham FFT.
c---------------------------------------------------------------------
      n1 = n / 2
      lj = 1
      li = 2 ** m
      do l = 1, m, 2
        lk = lj
        lj = 2 * lk
	li = li / 2
        ku = li + 1

        do i = 0, li - 1
          i11 = i * lk + 1
          i12 = i11 + n1
          i21 = i * lj + 1
          i22 = i21 + lk
        
          if (is .ge. 1) then
            u1 = exponent(ku+i)
          else
            u1 = dconjg (exponent(ku+i))
          endif
          do k = 0, lk - 1
            do j = 1, vlen
              x11 = x(j,i11+k)
              x21 = x(j,i12+k)
              scr(j,i21+k) = x11 + x21
              scr(j,i22+k) = u1 * (x11 - x21)
            end do
          end do
        end do

        if (l .eq. m) then
          do k = 1, n
            do j = 1, vlen
              x(j,k) = scr(j,k)
            enddo
          enddo
        else
          lk = lj
          lj = 2 * lk
          li = li / 2
          ku = li + 1

          do i = 0, li - 1
            i11 = i * lk + 1
            i12 = i11 + n1
            i21 = i * lj + 1
            i22 = i21 + lk
        
            if (is .ge. 1) then
              u1 = exponent(ku+i)
            else
              u1 = dconjg (exponent(ku+i))
            endif
            do k = 0, lk - 1
              do j = 1, vlen
                x11 = scr(j,i11+k)
                x21 = scr(j,i12+k)
                x(j,i21+k) = x11 + x21
                x(j,i22+k) = u1 * (x11 - x21)
              end do
            end do
          end do
        endif
      end do
      if (timers_enabled) call timer_stop(4)

      return
      end

c---------------------------------------------------------------------
      subroutine fftXYZ(sign,x,xout,exp1,exp2,exp3,n1,n2,n3)
        include 'global.h'
        integer sign,n1,n2,n3
        double complex x(n1+1,n2,n3)
        double complex xout((n1+1)*n2*n3)
        double complex exp1(n1), exp2(n2), exp3(n3)
        integer i, j, k, log
        integer bls,ble
	integer len
        integer blkp

        if (timers_enabled) call timer_start(3)	   

        fftblock=CacheSize/n1
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n1)
        if (timers_enabled) call timer_start(7) 
        do k = 1, n3
          do bls = 1, n2, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n2) ble = n2
	    len=ble-bls+1
            do j = bls, ble
            do i = 1, n1
              plane(j-bls+1+blkp*(i-1)) = x(i,j,k)
            end do
            end do
            call Swarztrauber(sign,log,len,n1,plane,blkp,exp1)     
            do j = bls, ble
            do i = 1, n1
              x(i,j,k)=plane(j-bls+1+blkp*(i-1))
            end do
            end do
          end do
        end do
        if (timers_enabled) call timer_stop(7)

        fftblock=CacheSize/n2
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n2 )
        if (timers_enabled) call timer_start(8)	   
        do k = 1, n3
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            call Swarztrauber(sign,log,len,n2,x(bls,1,k),n1+1,exp2)
	  enddo
        end do
        if (timers_enabled) call timer_stop(8)

        fftblock=CacheSize/n3
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
       	log = ilog2(n3)
        if (timers_enabled) call timer_start(9)
        do k = 1, n2
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            do i = 1,n3
            do j = bls, ble
              plane(j-bls+1+blkp*(i-1)) = x(j,k,i)
            end do
            end do
            call Swarztrauber(sign,log,len,n3,plane,blkp,exp3)
            do i = 0,n3-1
            do j = bls, ble
              xout(j+(n1+1)*(k-1+n2*i)) = plane(j-bls+1+blkp*i)
            end do		   
            end do
          end do
        end do  	       
        if (timers_enabled) call timer_stop(9)
        if (timers_enabled) call timer_stop(3)
        return
      end

invivofft3d_fftxyz__152

Real time: 3.601005e+09 -- Predicted time: 3.595193e+09

Callcount: 8

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	1	4.745844e+08	5.245534e+08	9.5
2	### 2	1	4.718219e+08	4.620311e+08	2.0
3	### 3	5.9979291150855	4.416169e+08	4.329994e+08	1.9

                         subroutine appft (niter, total_time, verified)
         include  'global.h'
         integer niter
	 double precision total_time
	 logical verified
!
! Local variables
!
         integer i, j, k, kt, n12, n22, n32, ii, jj, kk, ii2, ik2
         double precision ap
         double precision twiddle(nx+1,ny,nz)
        
         double complex xnt(nx+1,ny,nz),y(nx+1,ny,nz),
     &                  pad1(128),pad2(128)
         common /mainarrays/ xnt,pad1,y,pad2,twiddle

         double complex exp1(nx), exp2(ny), exp3(nz)

         do i=1,15
           call timer_clear(i)
	 end do         

         call timer_start(2)         
         call compute_initial_conditions(xnt,nx,ny,nz)
         
         call CompExp( nx, exp1 )
         call CompExp( ny, exp2 )
         call CompExp( nz, exp3 )           
         call fftXYZ(1,xnt,y,exp1,exp2,exp3,nx,ny,nz)
         call timer_stop(2)      

         call timer_start(1)
         if (timers_enabled) call timer_start(13)
	 
         n12 = nx/2
         n22 = ny/2
         n32 = nz/2
         ap = - 4.d0 * alpha * pi ** 2
         do i = 1, nz
	   ii = i-1-((i-1)/n32)*nz
	   ii2 = ii*ii
           do k = 1, ny
	     kk = k-1-((k-1)/n22)*ny
	     ik2 = ii2 + kk*kk
             do j = 1, nx
	         jj = j-1-((j-1)/n12)*nx
                 twiddle(j,k,i) = exp(ap*dble(jj*jj + ik2))
               end do
            end do
         end do
         if (timers_enabled) call timer_stop(13)      
         
         if (timers_enabled) call timer_start(12)
         call compute_initial_conditions(xnt,nx,ny,nz)             
         if (timers_enabled) call timer_stop(12)      
         if (timers_enabled) call timer_start(15)      
         call fftXYZ(1,xnt,y,exp1,exp2,exp3,nx,ny,nz)
         if (timers_enabled) call timer_stop(15)      

         do kt = 1, niter
	   if (timers_enabled) call timer_start(11)      
	   call evolve(xnt,y,twiddle,nx,ny,nz)
           if (timers_enabled) call timer_stop(11)      
           if (timers_enabled) call timer_start(15)      
           call fftXYZ(-1,xnt,xnt,exp1,exp2,exp3,nx,ny,nz)
           if (timers_enabled) call timer_stop(15)      
           if (timers_enabled) call timer_start(10)      
           call CalculateChecksum(sums(kt),kt,xnt,nx,ny,nz)           
           if (timers_enabled) call timer_stop(10)      
         end do
!
! Verification test.
!
         if (timers_enabled) call timer_start(14)      
         call verify(nx, ny, nz, niter, sums, verified)
         if (timers_enabled) call timer_stop(14)      
         call timer_stop(1)

	 total_time = timer_read(1)
	 if (.not.timers_enabled) return

	 print*,'FT subroutine timers '    
         write(*,40) 'FT total                  ', timer_read(1)
         write(*,40) 'WarmUp time               ', timer_read(2)
         write(*,40) 'fftXYZ body               ', timer_read(3)
         write(*,40) 'Swarztrauber              ', timer_read(4)
         write(*,40) 'X time                    ', timer_read(7)
         write(*,40) 'Y time                    ', timer_read(8)
         write(*,40) 'Z time                    ', timer_read(9)
         write(*,40) 'CalculateChecksum         ', timer_read(10)
         write(*,40) 'evolve                    ', timer_read(11)
         write(*,40) 'compute_initial_conditions', timer_read(12)
         write(*,40) 'twiddle                   ', timer_read(13)
         write(*,40) 'verify                    ', timer_read(14)
         write(*,40) 'fftXYZ                    ', timer_read(15)
         write(*,40) 'Benchmark time            ', total_time
   40    format(' ',A26,' =',F9.4)

         return
      end

invivoappft_appft__1

Real time: 8.548762e+09 -- Predicted time: 8.432556e+09

Callcount: 1

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	1	8.432556e+09	8.561607e+09	1.5

                c---------------------------------------------------------------------   
      subroutine Swarztrauber(is,m,vlen,n,x,xd1,exponent)

      implicit none
      include 'global.h'
c---------------------------------------------------------------------
c   Computes NY N-point complex-to-complex FFTs of X using an algorithm due
c   to Swarztrauber.  X is both the input and the output array, while Y is a 
c   scratch array.  It is assumed that N = 2^M.  Before calling 
c   Swarztrauber to 
c   perform FFTs
c---------------------------------------------------------------------
      integer is,m,vlen,n,xd1
      double complex x(xd1,n), exponent(n)      

      integer i,j,l
      double complex u1,x11,x21
      integer k, n1,li,lj,lk,ku,i11,i12,i21,i22
      
      if (timers_enabled) call timer_start(4)
c---------------------------------------------------------------------
c   Perform one variant of the Stockham FFT.
c---------------------------------------------------------------------
      n1 = n / 2
      lj = 1
      li = 2 ** m
      do l = 1, m, 2
        lk = lj
        lj = 2 * lk
	li = li / 2
        ku = li + 1

        do i = 0, li - 1
          i11 = i * lk + 1
          i12 = i11 + n1
          i21 = i * lj + 1
          i22 = i21 + lk
        
          if (is .ge. 1) then
            u1 = exponent(ku+i)
          else
            u1 = dconjg (exponent(ku+i))
          endif
          do k = 0, lk - 1
            do j = 1, vlen
              x11 = x(j,i11+k)
              x21 = x(j,i12+k)
              scr(j,i21+k) = x11 + x21
              scr(j,i22+k) = u1 * (x11 - x21)
            end do
          end do
        end do

        if (l .eq. m) then
          do k = 1, n
            do j = 1, vlen
              x(j,k) = scr(j,k)
            enddo
          enddo
        else
          lk = lj
          lj = 2 * lk
          li = li / 2
          ku = li + 1

          do i = 0, li - 1
            i11 = i * lk + 1
            i12 = i11 + n1
            i21 = i * lj + 1
            i22 = i21 + lk
        
            if (is .ge. 1) then
              u1 = exponent(ku+i)
            else
              u1 = dconjg (exponent(ku+i))
            endif
            do k = 0, lk - 1
              do j = 1, vlen
                x11 = scr(j,i11+k)
                x21 = scr(j,i12+k)
                x(j,i21+k) = x11 + x21
                x(j,i22+k) = u1 * (x11 - x21)
              end do
            end do
          end do
        endif
      end do
      if (timers_enabled) call timer_stop(4)

      return
      end

c---------------------------------------------------------------------
      subroutine fftXYZ(sign,x,xout,exp1,exp2,exp3,n1,n2,n3)
        include 'global.h'
        integer sign,n1,n2,n3
        double complex x(n1+1,n2,n3)
        double complex xout((n1+1)*n2*n3)
        double complex exp1(n1), exp2(n2), exp3(n3)
        integer i, j, k, log
        integer bls,ble
	integer len
        integer blkp

        if (timers_enabled) call timer_start(3)	   

        fftblock=CacheSize/n1
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n1)
        if (timers_enabled) call timer_start(7) 
        do k = 1, n3
          do bls = 1, n2, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n2) ble = n2
	    len=ble-bls+1
            do j = bls, ble
            do i = 1, n1
              plane(j-bls+1+blkp*(i-1)) = x(i,j,k)
            end do
            end do
            call Swarztrauber(sign,log,len,n1,plane,blkp,exp1)     
            do j = bls, ble
            do i = 1, n1
              x(i,j,k)=plane(j-bls+1+blkp*(i-1))
            end do
            end do
          end do
        end do
        if (timers_enabled) call timer_stop(7)

        fftblock=CacheSize/n2
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n2 )
        if (timers_enabled) call timer_start(8)	   
        do k = 1, n3
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            call Swarztrauber(sign,log,len,n2,x(bls,1,k),n1+1,exp2)
	  enddo
        end do
        if (timers_enabled) call timer_stop(8)

        fftblock=CacheSize/n3
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
       	log = ilog2(n3)
        if (timers_enabled) call timer_start(9)
        do k = 1, n2
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            do i = 1,n3
            do j = bls, ble
              plane(j-bls+1+blkp*(i-1)) = x(j,k,i)
            end do
            end do
            call Swarztrauber(sign,log,len,n3,plane,blkp,exp3)
            do i = 0,n3-1
            do j = bls, ble
              xout(j+(n1+1)*(k-1+n2*i)) = plane(j-bls+1+blkp*i)
            end do		   
            end do
          end do
        end do  	       
        if (timers_enabled) call timer_stop(9)
        if (timers_enabled) call timer_stop(3)
        return
      end

invivofft3d_fftxyz__137

Real time: 3.447745e+09 -- Predicted time: 3.444628e+09

Callcount: 8

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	1	4.305089e+08	4.429517e+08	2.8
2	### 2	6.99305129404454	4.310163e+08	4.335512e+08	0.5

                         subroutine appft (niter, total_time, verified)
         include  'global.h'
         integer niter
	 double precision total_time
	 logical verified
!
! Local variables
!
         integer i, j, k, kt, n12, n22, n32, ii, jj, kk, ii2, ik2
         double precision ap
         double precision twiddle(nx+1,ny,nz)
        
         double complex xnt(nx+1,ny,nz),y(nx+1,ny,nz),
     &                  pad1(128),pad2(128)
         common /mainarrays/ xnt,pad1,y,pad2,twiddle

         double complex exp1(nx), exp2(ny), exp3(nz)

         do i=1,15
           call timer_clear(i)
	 end do         

         call timer_start(2)         
         call compute_initial_conditions(xnt,nx,ny,nz)
         
         call CompExp( nx, exp1 )
         call CompExp( ny, exp2 )
         call CompExp( nz, exp3 )           
         call fftXYZ(1,xnt,y,exp1,exp2,exp3,nx,ny,nz)
         call timer_stop(2)      

         call timer_start(1)
         if (timers_enabled) call timer_start(13)
	 
         n12 = nx/2
         n22 = ny/2
         n32 = nz/2
         ap = - 4.d0 * alpha * pi ** 2
         do i = 1, nz
	   ii = i-1-((i-1)/n32)*nz
	   ii2 = ii*ii
           do k = 1, ny
	     kk = k-1-((k-1)/n22)*ny
	     ik2 = ii2 + kk*kk
             do j = 1, nx
	         jj = j-1-((j-1)/n12)*nx
                 twiddle(j,k,i) = exp(ap*dble(jj*jj + ik2))
               end do
            end do
         end do
         if (timers_enabled) call timer_stop(13)      
         
         if (timers_enabled) call timer_start(12)
         call compute_initial_conditions(xnt,nx,ny,nz)             
         if (timers_enabled) call timer_stop(12)      
         if (timers_enabled) call timer_start(15)      
         call fftXYZ(1,xnt,y,exp1,exp2,exp3,nx,ny,nz)
         if (timers_enabled) call timer_stop(15)      

         do kt = 1, niter
	   if (timers_enabled) call timer_start(11)      
	   call evolve(xnt,y,twiddle,nx,ny,nz)
           if (timers_enabled) call timer_stop(11)      
           if (timers_enabled) call timer_start(15)      
           call fftXYZ(-1,xnt,xnt,exp1,exp2,exp3,nx,ny,nz)
           if (timers_enabled) call timer_stop(15)      
           if (timers_enabled) call timer_start(10)      
           call CalculateChecksum(sums(kt),kt,xnt,nx,ny,nz)           
           if (timers_enabled) call timer_stop(10)      
         end do
!
! Verification test.
!
         if (timers_enabled) call timer_start(14)      
         call verify(nx, ny, nz, niter, sums, verified)
         if (timers_enabled) call timer_stop(14)      
         call timer_stop(1)

	 total_time = timer_read(1)
	 if (.not.timers_enabled) return

	 print*,'FT subroutine timers '    
         write(*,40) 'FT total                  ', timer_read(1)
         write(*,40) 'WarmUp time               ', timer_read(2)
         write(*,40) 'fftXYZ body               ', timer_read(3)
         write(*,40) 'Swarztrauber              ', timer_read(4)
         write(*,40) 'X time                    ', timer_read(7)
         write(*,40) 'Y time                    ', timer_read(8)
         write(*,40) 'Z time                    ', timer_read(9)
         write(*,40) 'CalculateChecksum         ', timer_read(10)
         write(*,40) 'evolve                    ', timer_read(11)
         write(*,40) 'compute_initial_conditions', timer_read(12)
         write(*,40) 'twiddle                   ', timer_read(13)
         write(*,40) 'verify                    ', timer_read(14)
         write(*,40) 'fftXYZ                    ', timer_read(15)
         write(*,40) 'Benchmark time            ', total_time
   40    format(' ',A26,' =',F9.4)

         return
      end

invivoappft_appft__40

Real time: 7.234248e+08 -- Predicted time: 6.806448e+08

Callcount: 1

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	1	6.806448e+08	7.224093e+08	5.7

                c---------------------------------------------------------------------   
      subroutine Swarztrauber(is,m,vlen,n,x,xd1,exponent)

      implicit none
      include 'global.h'
c---------------------------------------------------------------------
c   Computes NY N-point complex-to-complex FFTs of X using an algorithm due
c   to Swarztrauber.  X is both the input and the output array, while Y is a 
c   scratch array.  It is assumed that N = 2^M.  Before calling 
c   Swarztrauber to 
c   perform FFTs
c---------------------------------------------------------------------
      integer is,m,vlen,n,xd1
      double complex x(xd1,n), exponent(n)      

      integer i,j,l
      double complex u1,x11,x21
      integer k, n1,li,lj,lk,ku,i11,i12,i21,i22
      
      if (timers_enabled) call timer_start(4)
c---------------------------------------------------------------------
c   Perform one variant of the Stockham FFT.
c---------------------------------------------------------------------
      n1 = n / 2
      lj = 1
      li = 2 ** m
      do l = 1, m, 2
        lk = lj
        lj = 2 * lk
	li = li / 2
        ku = li + 1

        do i = 0, li - 1
          i11 = i * lk + 1
          i12 = i11 + n1
          i21 = i * lj + 1
          i22 = i21 + lk
        
          if (is .ge. 1) then
            u1 = exponent(ku+i)
          else
            u1 = dconjg (exponent(ku+i))
          endif
          do k = 0, lk - 1
            do j = 1, vlen
              x11 = x(j,i11+k)
              x21 = x(j,i12+k)
              scr(j,i21+k) = x11 + x21
              scr(j,i22+k) = u1 * (x11 - x21)
            end do
          end do
        end do

        if (l .eq. m) then
          do k = 1, n
            do j = 1, vlen
              x(j,k) = scr(j,k)
            enddo
          enddo
        else
          lk = lj
          lj = 2 * lk
          li = li / 2
          ku = li + 1

          do i = 0, li - 1
            i11 = i * lk + 1
            i12 = i11 + n1
            i21 = i * lj + 1
            i22 = i21 + lk
        
            if (is .ge. 1) then
              u1 = exponent(ku+i)
            else
              u1 = dconjg (exponent(ku+i))
            endif
            do k = 0, lk - 1
              do j = 1, vlen
                x11 = scr(j,i11+k)
                x21 = scr(j,i12+k)
                x(j,i21+k) = x11 + x21
                x(j,i22+k) = u1 * (x11 - x21)
              end do
            end do
          end do
        endif
      end do
      if (timers_enabled) call timer_stop(4)

      return
      end

c---------------------------------------------------------------------
      subroutine fftXYZ(sign,x,xout,exp1,exp2,exp3,n1,n2,n3)
        include 'global.h'
        integer sign,n1,n2,n3
        double complex x(n1+1,n2,n3)
        double complex xout((n1+1)*n2*n3)
        double complex exp1(n1), exp2(n2), exp3(n3)
        integer i, j, k, log
        integer bls,ble
	integer len
        integer blkp

        if (timers_enabled) call timer_start(3)	   

        fftblock=CacheSize/n1
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n1)
        if (timers_enabled) call timer_start(7) 
        do k = 1, n3
          do bls = 1, n2, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n2) ble = n2
	    len=ble-bls+1
            do j = bls, ble
            do i = 1, n1
              plane(j-bls+1+blkp*(i-1)) = x(i,j,k)
            end do
            end do
            call Swarztrauber(sign,log,len,n1,plane,blkp,exp1)     
            do j = bls, ble
            do i = 1, n1
              x(i,j,k)=plane(j-bls+1+blkp*(i-1))
            end do
            end do
          end do
        end do
        if (timers_enabled) call timer_stop(7)

        fftblock=CacheSize/n2
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
	log = ilog2( n2 )
        if (timers_enabled) call timer_start(8)	   
        do k = 1, n3
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            call Swarztrauber(sign,log,len,n2,x(bls,1,k),n1+1,exp2)
	  enddo
        end do
        if (timers_enabled) call timer_stop(8)

        fftblock=CacheSize/n3
	if(fftblock.ge.BlockMax) fftblock=BlockMax
	blkp=fftblock+1
       	log = ilog2(n3)
        if (timers_enabled) call timer_start(9)
        do k = 1, n2
          do bls = 1, n1, fftblock
	    ble = bls + fftblock - 1
	    if ( ble .gt. n1) ble = n1
	    len=ble-bls+1
            do i = 1,n3
            do j = bls, ble
              plane(j-bls+1+blkp*(i-1)) = x(j,k,i)
            end do
            end do
            call Swarztrauber(sign,log,len,n3,plane,blkp,exp3)
            do i = 0,n3-1
            do j = bls, ble
              xout(j+(n1+1)*(k-1+n2*i)) = plane(j-bls+1+blkp*i)
            end do		   
            end do
          end do
        end do  	       
        if (timers_enabled) call timer_stop(9)
        if (timers_enabled) call timer_stop(3)
        return
      end

invivofft3d_swarztrauber__27

Real time: 8.847104e+09 -- Predicted time: 9.241310e+09

Callcount: 32768

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
14276	### 1	7156.8835542227	4.192960e+05	4.134800e+05	1.3
4180	### 2	8195.06987453204	3.568270e+05	3.344280e+05	6.2
17561	### 3	1029.32478321909	5.200480e+05	4.843600e+05	6.8
2234	### 4	16399.5501721471	1.695740e+05	1.603280e+05	5.4

                c---------------------------------------------------------------------
c compute the roots-of-unity array that will be used for subsequent FFTs. 
c---------------------------------------------------------------------
      subroutine CompExp (n, exponent)

      implicit none
      integer n
      double complex exponent(n) 
      integer ilog2
      external ilog2      
     
      integer m,nu,ku,i,j,ln
      double precision t, ti, pi 
      data pi /3.141592653589793238d0/

      nu = n
      m = ilog2(n)
      exponent(1) = m
      ku = 2
      ln = 1
      do j = 1, m
         t = pi / ln
         do i = 0, ln - 1
            ti = i * t
            exponent(i+ku) = dcmplx(cos(ti),sin(ti))
         enddo        
         ku = ku + ln
         ln = 2 * ln
      enddo
            
      return
      end

c---------------------------------------------------------------------
c---------------------------------------------------------------------
c---------------------------------------------------------------------
      integer function ilog2(n)
      implicit none
      integer n
c---------------------------------------------------------------------
c---------------------------------------------------------------------  
      integer nn, lg
      if (n .eq. 1) then
         ilog2=0
         return
      endif
      lg = 1
      nn = 2
      do while (nn .lt. n)
         nn = nn*2
         lg = lg+1
      end do
      ilog2 = lg
      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------

      subroutine ipow46(a, exponent, result)
c---------------------------------------------------------------------
c compute a^exponent mod 2^46
c---------------------------------------------------------------------

      implicit none
      double precision a, result, dummy, q, r
      integer exponent, n, n2
      external randlc
      double precision randlc
c---------------------------------------------------------------------
c Use
c   a^n = a^(n/2)*a^(n/2) if n even else
c   a^n = a*a^(n-1)       if n odd
c---------------------------------------------------------------------
      result = 1
      if (exponent .eq. 0) return
      q = a
      r = 1
      n = exponent

      do while (n .gt. 1)
         n2 = n/2
         if (n2 * 2 .eq. n) then
            dummy = randlc(q, q) 
            n = n2
         else
            dummy = randlc(r, q)
            n = n-1
         endif
      end do
      dummy = randlc(r, q)
      result = r
      return
      end
c---------------------------------------------------------------------
      subroutine CalculateChecksum(csum,iterN,u,d1,d2,d3)
        implicit none
        integer iterN
        integer d1,d2,d3
        double complex csum
        double complex u(d1+1,d2,d3)
        integer i, i1, ii, ji, ki
        csum = dcmplx (0.0, 0.0)
        do i = 1, 1024
          i1 = i
          ii = mod (i1, d1) + 1
          ji = mod (3 * i1, d2) + 1
          ki = mod (5 * i1, d3) + 1
          csum = csum + u(ii,ji,ki)
        end do
        csum = csum/dble(d1*d2*d3)
        write(*,30) iterN, csum
 30     format (' T =',I5,5X,'Checksum =',1P2D22.12)
      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------
      subroutine compute_initial_conditions(u0,d1,d2,d3)

      implicit none
      include 'npbparams.h'
      integer d1,d2,d3
      double complex u0(d1+1,d2,d3), tmp(maxdim)
      double precision x0, start, an, dummy
      double precision RanStarts(maxdim)

      integer i,j,k
      double precision seed, a
      parameter (seed = 314159265.d0, a = 1220703125.d0)
      external randlc
      double precision randlc
      
      start = seed                                    
c---------------------------------------------------------------------
c Jump to the starting element for our first plane.
c---------------------------------------------------------------------
      call ipow46(a, 0, an)
      dummy = randlc(start, an)
      call ipow46(a, 2*d1*d2, an)
c---------------------------------------------------------------------
c Go through by z planes filling in one square at a time.
c---------------------------------------------------------------------
      RanStarts(1) = start
      do k = 2, d3 
         dummy = randlc(start, an)
         RanStarts(k) = start
      end do
      
      do k = 1, d3 
         x0 = RanStarts(k)
         do j = 1, d2 
           call vranlc(2*d1, x0, a, tmp)
           do i = 1, d1 
             u0(i,j,k)=tmp(i)
           end do
         end do
      end do

      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------
      subroutine evolve(x,y,twiddle,nx,ny,nz)
      double complex x(nx+1,ny,nz),y(nx+1,ny,nz)
      real*8 twiddle(nx+1,ny,nz)
      integer nx,ny,nz
      integer i,j,k
           do i = 1, nz
             do k = 1, ny
               do j = 1, nx
                   y(j,k,i)=y(j,k,i)*twiddle(j,k,i)
                   x(j,k,i)=y(j,k,i)
                 end do
              end do
           end do
      
      return
      end
c---------------------------------------------------------------------
c---------------------------------------------------------------------

invivoauxfnct_compute_initial_conditions__148

Real time: 2.814867e+08 -- Predicted time: 2.104654e+08

Callcount: 2

Invocation	Cluster	Part	Invitro (cycles)	Invivo (cycles)	Error (%)
1	### 1	1	1.054505e+08	1.795510e+08	41.2
2	### 2	1	1.050149e+08	1.054767e+08	0.4

FT

Exec Time : 1.257131e+10 cycles

Matching : 95.4%

__invivo__fft3d_fftxyz__112

Real time: 3.605644e+09 -- Predicted time: 3.657419e+09

Callcount: 8

__invivo__auxfnct_evolve__167

Real time: 6.506289e+08 -- Predicted time: 6.546663e+08

Callcount: 6

__invivo__fft3d_fftxyz__152

Real time: 3.601005e+09 -- Predicted time: 3.595193e+09

Callcount: 8

__invivo__appft_appft__1

Real time: 8.548762e+09 -- Predicted time: 8.432556e+09

Callcount: 1

__invivo__fft3d_fftxyz__137

Real time: 3.447745e+09 -- Predicted time: 3.444628e+09

Callcount: 8

__invivo__appft_appft__40

Real time: 7.234248e+08 -- Predicted time: 6.806448e+08

Callcount: 1

__invivo__fft3d_swarztrauber__27

Real time: 8.847104e+09 -- Predicted time: 9.241310e+09

Callcount: 32768

__invivo__auxfnct_compute_initial_conditions__148

Real time: 2.814867e+08 -- Predicted time: 2.104654e+08

Callcount: 2

invivofft3d_fftxyz__112

invivoauxfnct_evolve__167

invivofft3d_fftxyz__152

invivoappft_appft__1

invivofft3d_fftxyz__137

invivoappft_appft__40

invivofft3d_swarztrauber__27

invivoauxfnct_compute_initial_conditions__148