Performance, C vs. Fortran

Here’s what I eventually ended up with. Maybe you can see what you think. It has close to the same number of iterations by using the new suggested epsilon and putting fabs in the C code.

module rhofunc
implicit none
public
integer, parameter :: dp = kind(0.d0)
contains
pure real(dp) function rho(x, y)
   real(dp), intent(in) :: x, y
   if(x > 0.6_dp .and. x < 0.8_dp .and. y > 0.6_dp .and. y < 0.8) then
      rho = 1.0_dp
   else if(x > 0.2_dp .and. x < 0.4_dp .and. y > 0.2_dp .and. y < 0.4_dp) then
      rho = -1.0_dp
   else
      rho = 0.0_dp
   endif
end function

end module

program poisson
use rhofunc, only: rho
implicit none
integer, parameter :: dp = kind(0.d0), M = 100
integer :: i, j
real(dp), parameter :: epsilon0 = 3.85e-8_dp, target = 1e-6, a = 0.01
real(dp) :: delta, start, end, phiprime(0:M+1, 0:M+1), phi(0:M+1, 0:M+1), a2
integer :: iters
real(dp) :: rho_mat(M, M)

iters = 0
phi = 0.0_dp

delta = 1.0_dp

call cpu_time(start)

do i = 1, M
   do j = 1, M
      rho_mat(i, j) = rho(i*a, j*a)
   enddo
enddo

a2 = a**2.0_dp

do while(.true.)
   phiprime(1:M, 1:M) = (phi(2:M+1, 1:M) + phi(0:M-1, 1:M) + phi(1:M, 2:M+1) + phi(1:M, 0:M-1))/4.0_dp &
                        + a2/4.0_dp/epsilon0*rho_mat
   delta = maxval(abs(phiprime - phi))
   iters = iters + 1
   if(delta > target) then
      phi = phiprime
   else
      exit
   endif
enddo

call cpu_time(end)
print *, iters
print *, end - start
end program
1 Like