CPS-615 Fall 1996 -- Final Project -- The HPF Source Code

HPF Distribution Methods in a Model CFD Problem

Michael McMahon

!To compile: f90 -fixed -o executable_name -wsf #of_processors vent.f90
!This program solves laplace's equation for Psi over a grid that
!is shaped into an air vent.  The grid is shaped using the logical
!variable Hide.  The equation is solved using a red black SOR method
!with the relaxation parameter set to the best value I could find,
!omega=1.92.  The execution time and number of iterations required
!are output upon completion.  The critical line is the HPF DISTRIBUTE
!line.  This line allows the method of distribution to be altered
!and the run time to be reduced to a minimum.

        Program Invis
        Parameter(n=100)
        Parameter(m=200)
        Parameter(z=4)
        Real,dimension(1:n,1:m)  ::  Psi,Psiold
        Real omega,tol,tstart,tend
        Logical,dimension(1:n,1:m) :: hide,even
        Integer i,j
        Integer iteration

!****************HPF Directives*************************

!HPF$ DISTRIBUTE Psi(*,cyclic(10))
!HPF$ ALIGN (:,:) with psi :: Psiold,hide,even

!*******************************************************

!Start to shape grid and begin timing
!Logical variable even is used for red-black algorithm
        hide = .false.
        even = .false.

        call timer(tstart,0.0)
        Forall(i=1:n,j=1:m,mod((i+j),2) .eq. 0) even(i,j)=.true.
        Forall(i=1:n-1,j=1:m/4,(i-j) > 0) hide(i,j) = .true.
        Forall(i=n/2+1:n-1,j=m/4:3*m/4) hide(i,j) = .true.
        Forall(i=1:n-1,j=3*m/4:m,(i+j-m) > 0) hide(i,j) = .true.

!Set boundary conditions

        Forall(i=1:n/2,j=1:m/4,(i-j).eq.0) Psi(i,j) = 50.0
        Forall(j=m/4:3*m/4) Psi(n/2,j)=50.0
        Forall(i=1:n/2,j=3*m/4:m/4,(i+j-200).eq.0) Psi(i,j)=50.0
        Forall(j=1:m) Psi(n,j) = 0.0

!Entrance and exit to vary linearly between top and bottom edges

        Forall(i=1:n) Psi(i,1)=50.0-50.0*real(i/n)
        Forall(i=1:n) Psi(i,m)=50.0-50.0*real(i/n)

!Preliminaries
        iteration=0
        Psiold=Psi
        omega=1.92
        tol=1.0
!on your marks, get set, go

        Do While(tol .ge. 1.0e-4)

!-------------red (even) sites----------------------
        Where(hide .eq. .true. .and. even .eq. .true.)
                psi=omega/4.0*(cshift(psi,dim=1,shift=1) 
     $       + cshift(psi,dim=1,shift= -1) 
     $       + cshift(psi,dim=2,shift=1)
     $       + cshift(psi,dim=2,shift= -1)) +(1.0-omega)*psi
        end where
!-------------------black (odd) sites------------------
        Where((hide .eq. .true.) .and. even .eq. .false.)
                psi=omega/4.0*(cshift(psi,dim=1,shift=1) 
     $       + cshift(psi,dim=1,shift= -1) 
     $       + cshift(psi,dim=2,shift=1)
     $       + cshift(psi,dim=2,shift= -1)) +(1.0-omega)*psi
        end where
!-------------------black (odd) sites------------------
        Where((hide .eq. .true.) .and. even .eq. .false.)
                psi=omega/4.0*(cshift(psi,dim=1,shift=1) 
     $       + cshift(psi,dim=1,shift= -1) 
     $       + cshift(psi,dim=2,shift=1)
     $       + cshift(psi,dim=2,shift= -1))+(1.0-omega)*psi
        end where

        iteration=iteration+1
        tol=maxval(abs(psiold-psi))

        psiold=psi

        end do
        call timer(tend,tstart)

        write(6,*) 'This many interations:',iteration
        write(6,*) 'And it took this long:',tend,'seconds'
        end program

!++++++++++++++++
         Subroutine timer(return_time,initial_time)
         real, intent(in)  :: initial_time
         real,  intent(out) :: return_time
         integer finish,rate
          CALL system_clock( COUNT=finish,COUNT_RATE=rate)
          return_time = FLOAT(finish) / FLOAT(rate) - initial_time
          end