---
jupytext:
  formats: md:myst
  text_representation:
    extension: .md
    format_name: myst
kernelspec:
  display_name: Python 3
  language: python
  name: python3
---

Section: Implicit vs. Explicit Loops¶

Adapted from: https://github.com/gjbex/Fortran-MOOC/tree/master/source_code/implicit_loops

This program demonstrates the speed of implicit vs. explicit loops in Fortran.¶

The program will take in a command line argument and construct two argument X argument sized matrices composed of random numbers. The program will then compute the following:

\Large Result = \sum \sqrt{X^2 + Y^2}

(1)

where:
$Result$ is the summation result
$X$ and $Y$ are command line argument X command line argument sized matrices

program real32_vs_real64
    use, intrinsic :: iso_fortran_env, only : error_unit, I8 => INT64, &
        SP => REAL32, DP => REAL64
    implicit none
    integer(kind=I8) :: nr_values
    real(kind=SP), dimension(:), allocatable :: x_sp, y_sp
    real(kind=SP) :: result
    real :: start_time, end_time
    integer :: status, i

    call get_argument(nr_values)

    ! implicit loops
    allocate (x_sp(nr_values), y_sp(nr_values), stat=status)
    if (status /= 0) then
        write (unit=error_unit, fmt='(A)') 'error: can not allocate x_sp/y_sp'
        stop 1
    end if
    call random_number(x_sp)
    call random_number(y_sp)

    call cpu_time(start_time)
    result = sum(sqrt(x_sp**2 + y_sp**2))
    call cpu_time(end_time)
    print '(A, F20.8)', 'Implicit loop time [s]: ', end_time - start_time
    print '(A, F20.8)', 'Summation result: ', result
    deallocate (x_sp, y_sp)

    ! explicit loops
    allocate (x_sp(nr_values), y_sp(nr_values), stat=status)
    if (status /= 0) then
        write (unit=error_unit, fmt='(A)') 'error: can not allocate x_sp/y_sp'
        stop 1
    end if
    call random_number(x_sp)
    call random_number(y_sp)

    call cpu_time(start_time)
    result = 0.0_SP
    do i = 1, size(x_sp)
        result = result + sqrt(x_sp(i)**2 + y_sp(i)**2)
    end do
    call cpu_time(end_time)
    print '(A, F20.8)', 'Explicit loop time [s]: ', end_time - start_time
    print '(A, F20.8)', 'Summation result: ', result
    deallocate (x_sp, y_sp)

contains

    subroutine get_argument(nr_values)
        implicit none
        integer(Kind=I8), intent(out) :: nr_values
        character(len=1024) :: buffer, msg
        integer :: status

        if (command_argument_count() /= 1) then
            write (unit=error_unit, fmt='(A)') 'error: expect number of values'
            stop 2
        end if
        call get_command_argument(1, buffer)
        read (buffer, fmt=*, iomsg=msg, iostat=status) nr_values
        if (status /= 0) then
            write (unit=error_unit, fmt='(2A)') 'error: ', trim(msg)
            stop 3
        end if
    end subroutine get_argument

end program real32_vs_real64

The above program is compiled and run using Fortran Package Manager (fpm):

Build the Program using FPM (Fortran Package Manager)¶

import os
root_dir = ""
root_dir = os.getcwd()

code_dir = root_dir + "/" + "Fortran_Code/Section_Implicit_vs_Explicit_Loops"

%%capture
%%writefile section_implicit_vs_explicit_loops.f90
program real32_vs_real64
    use, intrinsic :: iso_fortran_env, only : error_unit, I8 => INT64, &
        SP => REAL32, DP => REAL64
    implicit none
    integer(kind=I8) :: nr_values
    real(kind=SP), dimension(:), allocatable :: x_sp, y_sp
    real(kind=SP) :: result
    real :: start_time, end_time
    integer :: status, i

    call get_argument(nr_values)

    ! implicit loops
    allocate (x_sp(nr_values), y_sp(nr_values), stat=status)
    if (status /= 0) then
        write (unit=error_unit, fmt='(A)') 'error: can not allocate x_sp/y_sp'
        stop 1
    end if
    call random_number(x_sp)
    call random_number(y_sp)

    call cpu_time(start_time)
    result = sum(sqrt(x_sp**2 + y_sp**2))
    call cpu_time(end_time)
    print '(A, F20.8)', 'Implicit loop time [s]: ', end_time - start_time
    print '(A, F20.8)', 'Summation result: ', result
    deallocate (x_sp, y_sp)

    ! explicit loops
    allocate (x_sp(nr_values), y_sp(nr_values), stat=status)
    if (status /= 0) then
        write (unit=error_unit, fmt='(A)') 'error: can not allocate x_sp/y_sp'
        stop 1
    end if
    call random_number(x_sp)
    call random_number(y_sp)

    call cpu_time(start_time)
    result = 0.0_SP
    do i = 1, size(x_sp)
        result = result + sqrt(x_sp(i)**2 + y_sp(i)**2)
    end do
    call cpu_time(end_time)
    print '(A, F20.8)', 'Explicit loop time [s]: ', end_time - start_time
    print '(A, F20.8)', 'Summation result: ', result
    deallocate (x_sp, y_sp)

contains

    subroutine get_argument(nr_values)
        implicit none
        integer(Kind=I8), intent(out) :: nr_values
        character(len=1024) :: buffer, msg
        integer :: status

        if (command_argument_count() /= 1) then
            write (unit=error_unit, fmt='(A)') 'error: expect number of values'
            stop 2
        end if
        call get_command_argument(1, buffer)
        read (buffer, fmt=*, iomsg=msg, iostat=status) nr_values
        if (status /= 0) then
            write (unit=error_unit, fmt='(2A)') 'error: ', trim(msg)
            stop 3
        end if
    end subroutine get_argument

end program real32_vs_real64

!bat *.f90

───────┬────────────────────────────────────────────────────────────────────────
       │ File: section_implicit_vs_explicit_loops.f90
───────┼────────────────────────────────────────────────────────────────────────
   1   │ program real32_vs_real64
   2   │     use, intrinsic :: iso_fortran_env, only : error_unit, I8 => INT64, 
       │ &
   3   │         SP => REAL32, DP => REAL64
   4   │     implicit none
   5   │     integer(kind=I8) :: nr_values
   6   │     real(kind=SP), dimension(:), allocatable :: x_sp, y_sp
   7   │     real(kind=SP) :: result
   8   │     real :: start_time, end_time
   9   │     integer :: status, i
  10   │ 
  11   │     call get_argument(nr_values)
  12   │ 
  13   │     ! implicit loops
  14   │     allocate (x_sp(nr_values), y_sp(nr_values), stat=status)
  15   │     if (status /= 0) then
  16   │         write (unit=error_unit, fmt='(A)') 'error: can not allocate x_s
       │ p/y_sp'
  17   │         stop 1
  18   │     end if
  19   │     call random_number(x_sp)
  20   │     call random_number(y_sp)
  21   │ 
  22   │     call cpu_time(start_time)
  23   │     result = sum(sqrt(x_sp**2 + y_sp**2))
  24   │     call cpu_time(end_time)
  25   │     print '(A, F20.8)', 'Implicit loop time [s]: ', end_time - start_ti
       │ me
  26   │     print '(A, F20.8)', 'Summation result: ', result
  27   │     deallocate (x_sp, y_sp)
  28   │ 
  29   │     ! explicit loops
  30   │     allocate (x_sp(nr_values), y_sp(nr_values), stat=status)
  31   │     if (status /= 0) then
  32   │         write (unit=error_unit, fmt='(A)') 'error: can not allocate x_s
       │ p/y_sp'
  33   │         stop 1
  34   │     end if
  35   │     call random_number(x_sp)
  36   │     call random_number(y_sp)
  37   │ 
  38   │     call cpu_time(start_time)
  39   │     result = 0.0_SP
  40   │     do i = 1, size(x_sp)
  41   │         result = result + sqrt(x_sp(i)**2 + y_sp(i)**2)
  42   │     end do
  43   │     call cpu_time(end_time)
  44   │     print '(A, F20.8)', 'Explicit loop time [s]: ', end_time - start_ti
       │ me
  45   │     print '(A, F20.8)', 'Summation result: ', result
  46   │     deallocate (x_sp, y_sp)
  47   │ 
  48   │ contains
  49   │ 
  50   │     subroutine get_argument(nr_values)
  51   │         implicit none
  52   │         integer(Kind=I8), intent(out) :: nr_values
  53   │         character(len=1024) :: buffer, msg
  54   │         integer :: status
  55   │ 
  56   │         if (command_argument_count() /= 1) then
  57   │             write (unit=error_unit, fmt='(A)') 'error: expect number of
       │  values'
  58   │             stop 2
  59   │         end if
  60   │         call get_command_argument(1, buffer)
  61   │         read (buffer, fmt=*, iomsg=msg, iostat=status) nr_values
  62   │         if (status /= 0) then
  63   │             write (unit=error_unit, fmt='(2A)') 'error: ', trim(msg)
  64   │             stop 3
  65   │         end if
  66   │     end subroutine get_argument
  67   │ 
  68   │ end program real32_vs_real64
───────┴────────────────────────────────────────────────────────────────────────

os.chdir(code_dir)

build_status = os.system("fpm build 2>/dev/null")

Run the Program using FPM (Fortran Package Manager)¶

exec_status = os.system("fpm run -- 10 2>/dev/null")

Implicit loop time [s]:           0.00000200
Summation result:           9.57156563
Explicit loop time [s]:           0.00000100
Summation result:           8.58767700

exec_status = os.system("fpm run -- 100 2>/dev/null")

Implicit loop time [s]:           0.00000100
Summation result:          82.67105103
Explicit loop time [s]:           0.00000100
Summation result:          78.12156677

exec_status = os.system("fpm run -- 1000 2>/dev/null")

Implicit loop time [s]:           0.00000300
Summation result:         738.62194824
Explicit loop time [s]:           0.00000400
Summation result:         779.73486328

exec_status = os.system("fpm run -- 10000 2>/dev/null")

Implicit loop time [s]:           0.00001100
Summation result:        7668.53564453
Explicit loop time [s]:           0.00002200
Summation result:        7716.35302734

exec_status = os.system("fpm run -- 100000 2>/dev/null")

Implicit loop time [s]:           0.00010800
Summation result:       76268.64062500
Explicit loop time [s]:           0.00022100
Summation result:       76432.85156250

exec_status = os.system("fpm run -- 1000000 2>/dev/null")

Implicit loop time [s]:           0.00211100
Summation result:      765269.12500000
Explicit loop time [s]:           0.00220800
Summation result:      765385.31250000

exec_status = os.system("fpm run -- 10000000 2>/dev/null")

Implicit loop time [s]:           0.02267501
Summation result:     7671492.50000000
Explicit loop time [s]:           0.02569100
Summation result:     7669864.50000000

exec_status = os.system("fpm run -- 100000000 2>/dev/null")

Implicit loop time [s]:           0.10609502
Summation result:    33554432.00000000
Explicit loop time [s]:           0.23883396
Summation result:    33554432.00000000