rodinia-srad: fortan: Initial port

src/benchmarks/rodinia-srad/fortran/meson.build

+cc = meson.get_compiler('c')
+
+sources = [
+	'srad.F08',
+]
+
+dependencies = [
+	dependency('openmp'),
+]
+
+# Proxy the RAND_MAX value of libc to the fortran code.
+rand_max = cc.get_define('RAND_MAX', prefix: '#include <stdlib.h>')
+
+arguments = [
+	'-DRAND_MAX=@0@'.format(rand_max),
+]
+
+options = [
+	'fortran_std=f2008',
+]
+
+executable(
+	'rodinia-srad-fortran',
+	sources,
+	fortran_args: arguments,
+	dependencies: dependencies,
+	override_options: options,
+)

src/benchmarks/rodinia-srad/fortran/srad.F08

+program srad
+   use, intrinsic :: iso_c_binding
+   use, intrinsic :: iso_fortran_env
+
+   implicit none
+
+   ! Use the RNG from libc so that all programs produce the same output
+   interface
+      subroutine libc_srand(seed) bind(C, name="srand")
+         import :: c_int
+         integer(c_int), intent(in), value :: seed
+      end subroutine
+
+      integer(c_int) function libc_rand() bind(C, name="rand")
+         import :: c_int
+      end function
+   end interface
+
+   call main()
+
+contains
+   subroutine usage()
+      character(len=1024) :: program
+
+      call get_command_argument(0, program)
+
+      write (error_unit, *) "Usage: "//trim(program)//" <rows> <cols> <y1> <y2> <x1> <x2> <lambda> <no. of iter>"
+      write (error_unit, *) "<rows>        - number of rows"
+      write (error_unit, *) "<cols>        - number of columns"
+      write (error_unit, *) "<y1>          - y1 value of the speckle"
+      write (error_unit, *) "<y2>          - y2 value of the speckle"
+      write (error_unit, *) "<x1>          - x1 value of the speckle"
+      write (error_unit, *) "<x2>          - x2 value of the speckle"
+      write (error_unit, *) "<lambda>      - lambda (0, 1)"
+      write (error_unit, *) "<no. of iter> - number of iterations"
+
+      stop
+   end subroutine
+
+   subroutine random_matrix(I)
+      real, dimension(:, :), intent(inout) :: I
+
+      integer :: x
+      integer :: y
+
+      call libc_srand(7)
+
+      do x = 1, size(I, 1)
+         do y = 1, size(I, 2)
+            I(y, x) = libc_rand()/real(RAND_MAX)
+         end do
+      end do
+   end subroutine
+
+   subroutine main()
+      integer :: rows, cols
+      integer :: r1, r2
+      integer :: c1, c2
+
+      real :: lambda
+      integer :: niter
+
+      character(len=32) :: arg
+
+      real, allocatable, dimension(:, :) :: I, J, c
+      integer, allocatable, dimension(:) :: iN, iS, jW, jE
+      real, allocatable, dimension(:, :) :: dN, dS, dW, dE
+
+      integer :: x, y
+      integer :: iter
+
+      real :: size_R
+
+      real :: tmp, sum, sum2
+      real :: meanROI, varROI, q0sqr
+
+      real :: Jc, G2, L, num, den, qsqr
+      real :: cN, cS, cW, cE, D
+
+      if (command_argument_count() /= 8) then
+         call usage()
+      end if
+
+      call get_command_argument(1, arg)
+      read (arg, *) rows
+
+      call get_command_argument(2, arg)
+      read (arg, *) cols
+
+      if (mod(rows, 16) /= 0 .or. mod(cols, 16) /= 0) then
+         error stop "rows and cols must be multiples of 16"
+      end if
+
+      call get_command_argument(3, arg)
+      read (arg, *) r1
+
+      call get_command_argument(4, arg)
+      read (arg, *) r2
+
+      call get_command_argument(5, arg)
+      read (arg, *) c1
+
+      call get_command_argument(6, arg)
+      read (arg, *) c2
+
+      call get_command_argument(7, arg)
+      read (arg, *) lambda
+
+      call get_command_argument(8, arg)
+      read (arg, *) niter
+
+      size_R = (r2 - r1 + 1)*(c2 - c1 + 1)
+
+      allocate (I(cols, rows))
+      allocate (J(cols, rows))
+      allocate (c(cols, rows))
+
+      allocate (iN(rows))
+      allocate (iS(rows))
+      allocate (jW(cols))
+      allocate (jE(cols))
+
+      allocate (dN(cols, rows))
+      allocate (dS(cols, rows))
+      allocate (dW(cols, rows))
+      allocate (dE(cols, rows))
+
+      do x = 1, rows
+         iN(x) = x - 1
+         iS(x) = x + 1
+      end do
+
+      do y = 1, cols
+         jW(y) = y - 1
+         jE(y) = y + 1
+      end do
+
+      iN(1) = 1
+      iS(rows) = rows
+
+      jW(1) = 1
+      jE(cols) = cols
+
+      print '(a)', "Randomizing the input matrix"
+      call random_matrix(I)
+
+      J = exp(I)
+
+      print '(a)', "Start the SRAD main loop"
+
+      do iter = 1, niter
+         sum = 0
+         sum2 = 0
+
+         do x = r1, r2
+            do y = c1, c2
+               tmp = J(y + 1, x + 1)
+
+               sum = sum + tmp
+               sum2 = sum2 + tmp*tmp
+            end do
+         end do
+
+         meanROI = sum/size_R
+         varROI = (sum2/size_R) - meanROI*meanROI
+         q0sqr = varROI/(meanROI*meanROI)
+
+         !$omp parallel do shared(J, dN, dS, dW, dE, c, iN, iS, jW, jE) private(x, y, Jc, G2, L, num, den, qsqr)
+         do x = 1, size(J, 2)
+            do y = 1, size(J, 1)
+               Jc = J(y, x)
+
+               ! directional derivates
+               dN(y, x) = J(y, iN(x)) - Jc
+               dS(y, x) = J(y, iS(x)) - Jc
+               dW(y, x) = J(jW(y), x) - Jc
+               dE(y, x) = J(jE(y), x) - Jc
+
+               G2 = ((dN(y, x)*dN(y, x)) + &
+                    (dS(y, x)*dS(y, x)) + &
+                    (dW(y, x)*dW(y, x)) + &
+                    (dE(y, x)*dE(y, x)))/(Jc*Jc)
+
+               L = (dN(y, x) + dS(y, x) + dW(y, x) + dE(y, x))/Jc
+
+               num = (0.5*G2) - ((1.0/16.0)*L*L)
+               den = 1 + (0.25*L)
+
+               qsqr = num/(den*den)
+
+               ! diffusion coefficent (equ 33)
+               den = (qsqr - q0sqr)/(q0sqr*(1 + q0sqr))
+               c(y, x) = 1.0/(1.0 + den)
+
+               ! saturate diffusion coefficent
+               if (c(y, x) < 0) then
+                  c(y, x) = 0
+               else if (c(y, x) > 1) then
+                  c(y, x) = 1
+               end if
+            end do
+         end do
+         !$omp end parallel do
+
+         !$omp parallel do shared(J, c, lambda) private(x, y, D, cS, cN, cW, cE)
+         do x = 1, size(J, 2)
+            do y = 1, size(J, 1)
+               ! diffusion coefficent
+               cN = c(y, x)
+               cS = c(y, iS(x))
+               cW = c(y, x)
+               cE = c(jE(y), x)
+
+               ! divergence (equ 58)
+               D = cN*dN(y, x) + cS*dS(y, x) + cW*dW(y, x) + cE*dE(y, x)
+
+               ! image update(equ 61)
+               J(y, x) = J(y, x) + 0.25*lambda*D;
+            end do
+         end do
+         !$omp end parallel do
+      end do
+
+      do x = 1, rows
+         do y = 1, cols
+            write (*, '(f0.5)', advance="no") J(y, x)
+            write (*, '(a)', advance="no") " "
+         end do
+
+         print '(a)', ""
+      end do
+
+      print '(a)', "Computation Done"
+
+      deallocate (I)
+      deallocate (J)
+      deallocate (iN)
+      deallocate (iS)
+      deallocate (jW)
+      deallocate (jE)
+      deallocate (dN)
+      deallocate (dS)
+      deallocate (dW)
+      deallocate (dE)
+
+      deallocate (c)
+
+   end subroutine
+end program srad

src/benchmarks/rodinia-srad/meson.build

 subdir('c')
 subdir('cpp')
+subdir('fortran')