Search code examples
fortran

Optimal way to create once, then frequently access to a large array in Fortran


I have a scientific simulation in Fortran. In it, I need to calculate a large array (say 10^5 elements, a size known / that can be preset when writing the code) "once" around the beginning, and then this array will need to be shared / accessed / read into a lot of times (basically at each step, so possibly millions of times).

What is the recommended way of sharing such an array with other subroutines ? Given that I want the code creating the large array to be its own entity.

Say there is case 1, where :

  • The array needs to be generated at runtime, and actually will be generated once per run of the Monte-Carlo (so it will actually be generated more than once during the simulation). The size of the array can be considered fixed and known / preset when writing the code.
  • a few different subroutines will need access to the array at each timestep (say ~10^6 timesteps per run)
  • each of them will read only one element at a time
  • each of them will read into it almost sequentially (and almost at the same points)

There is case 2, where :

  • The array needs to be generated at runtime, and will actually be generated only once (it will be the same for all runs of the Monte-Carlo)
  • The array will actually be noticeably smaller. The size of the array can be considered fixed and known / preset when writing the code.
  • It will only be used in one subroutine, but it will need to be almost read entirely every time (ie at every timestep).

Finally there is case 3, which is exactly as case 2, except the array could be physically preset / written in code as a data statement.

What would be the best practice in each case, and to optimize speed / minimize possible runtime penalties ?

Edit : this is the solution I have in place for now. Have a module of the form

    module big_array
    parameter (nbigarray=10**5)
    
    real my_by_array(nbigarray)
    
    contains
    
      subroutine init_big_array(arglist)
    
      *[here some code populating / creating *my_big_array*]*

      return
      end
    
    end module

run init_big_array once at the start of the program and use the module wherever my_big_array is needed.


Solution

  • Here is an intro to dynamic memory management that allows one to expose array to subprograms.

    module foo
    
       implicit none  ! Always include this statement
       !
       ! Make everything private and only expose those entities
       ! that are needed.
       !
       private
       public array, init_array, destroy_array
    
       integer, allocatable :: array(:)
    
       contains
          !
          ! Destroy array
          !
          subroutine destroy_array
             if (allocated(array)) then
                deallocate(array)
             end if   
          end subroutine
    
          subroutine init_array(n)
             integer, intent(in) :: n
             if (n < 0) then
                !
                ! Could implement a fallback.  Perhaps,
                ! allocate a zero-sized array.
                !
                write(*,'(A)') 'init_array: memory allocation issue'
                stop
             end if
             !
             ! Check if array is allocated and the requested size
             ! differs from the current allocated size.
             !
             if (allocated(array)) then
                if (size(array) /= n) then
                   call destroy_array
                   allocate(array(n))
                end if
             else
                allocate(array(n))
             end if   
          end subroutine
    end module foo
    
    program main
       use foo, only : array, init_array
       implicit none
       integer m
       !
       ! Perhaps, "read(*,*) m" to allow user to size
       ! memory appropriately.
       !
       m = 100
    
       call init_array(m)
       call bar
    end program
    
    subroutine bar
       use foo, only : array
       implicit none
       if (allocated(array) then
          print *, size(array)
       else
          write(*,'(A)') 'whoops'
       end if
    end subroutine bar