Allocatable Arrays and Dynamic Memory

real(dp), allocatable :: u(:, :)        ! rank-2 array: shape unknown until allocated
real(dp), allocatable :: t0             ! allocatable scalar
character(:), allocatable :: label      ! deferred-length character scalar
type(grid_t), allocatable :: grids(:)   ! array of a derived type
class(solver_t), allocatable :: s       ! polymorphic scalar: dynamic type chosen at allocation

allocate ([type-spec ::] allocation-list [, alloc-opt-list])
! allocation is   allocate-object [(allocate-shape-spec-list)] [lbracket allocate-coarray-spec rbracket]
! allocate-shape-spec is   [lower-bound-expr :] upper-bound-expr
! alloc-opt is one of   stat=stat-variable   errmsg=errmsg-variable   source=source-expr   mold=source-expr

program allocation_basics
  use, intrinsic :: iso_fortran_env, only: real64, int32
  implicit none
  integer, parameter :: dp = real64
  real(dp), allocatable :: u(:, :)      ! deferred-shape rank-2 array, unallocated at start
  real(dp), allocatable :: t0           ! allocatable scalar
  integer(int32) :: nx, ny

  nx = 5_int32
  ny = 4_int32
  print '(a, l1)', 'allocated(u) at start:      ', allocated(u)

  allocate (u(nx, ny))                  ! bounds 1:nx, 1:ny
  u = 0.0_dp
  print '(a, l1, a, 2i3, a, 2i3, a, 2i3)', 'allocated(u) after allocate: ', allocated(u), &
    '  shape=', shape(u), '  lbound=', lbound(u), '  ubound=', ubound(u)
  deallocate (u)
  print '(a, l1)', 'allocated(u) after deallocate: ', allocated(u)

  allocate (u(0:nx - 1, -1:ny - 2))     ! explicit lower bounds
  print '(a, 2i3, a, 2i3)', 'lbound=', lbound(u), '  ubound=', ubound(u)
  deallocate (u)

  allocate (u(3, 0))                    ! upper bound < lower bound: zero-sized but allocated
  print '(a, l1, a, i0, a, 2i3)', 'zero-sized: allocated=', allocated(u), '  size=', size(u), '  shape=', shape(u)

  allocate (t0, source=300.0_dp)        ! allocatable scalar, value from source=
  print '(a, f6.1)', 't0=', t0
end program allocation_basics

allocated(u) at start:      F
allocated(u) after allocate: T  shape=  5  4  lbound=  1  1  ubound=  5  4
allocated(u) after deallocate: F
lbound=  0 -1  ubound=  4  2
zero-sized: allocated=T  size=0  shape=  3  0
t0= 300.0

program source_and_mold
  use, intrinsic :: iso_fortran_env, only: real64
  implicit none
  integer, parameter :: dp = real64
  integer :: i
  real(dp) :: field(0:3, 2)
  real(dp), allocatable :: copy(:, :), shape_only(:, :), typed(:), zeros(:)
  character(:), allocatable :: label, trimmed
  character(len=12) :: name = 'heatflow'

  field = reshape([real(dp) :: (real(i, dp), i=1, 8)], [4, 2])

  ! source=: same shape, lower bounds and values as the source expression
  allocate (copy, source=field)
  print '(a, 2i3, a, 2i3, a, f5.1)', 'source: lbound=', lbound(copy), ' ubound=', ubound(copy), ' copy(3,2)=', copy(3, 2)

  ! mold=: shape, type and type parameters only; the value is NOT copied
  allocate (shape_only, mold=field)
  print '(a, 2i3, a, 2i3)', 'mold:   lbound=', lbound(shape_only), ' ubound=', ubound(shape_only)

  ! a scalar source= initializes every element of an array allocation
  allocate (zeros(6), source=-1.0_dp)
  print '(a, i0, a, f5.1)', 'scalar source: size=', size(zeros), ' sum=', sum(zeros)

  ! type-spec: required for a deferred-length character unless source= or mold= is given
  allocate (character(len=len_trim(name)) :: label)
  print '(a, i0)', 'after allocate: len(label)=', len(label)
  label = name                         ! intrinsic assignment reallocates to len(name) = 12
  print '(a, i0, 3a)', 'after label = name: len(label)=', len(label), ' [', label, ']'
  trimmed = trim(name)                 ! allocated on assignment with the trimmed length
  print '(a, i0, 3a)', 'trimmed: len=', len(trimmed), ' [', trimmed, ']'

  allocate (real(dp) :: typed(5))      ! the type-spec form also works for non-character types
  print '(a, i0)', 'typed size=', size(typed)
end program source_and_mold

source: lbound=  0  1 ubound=  3  2 copy(3,2)=  8.0
mold:   lbound=  0  1 ubound=  3  2
scalar source: size=6 sum= -6.0
after allocate: len(label)=8
after label = name: len(label)=12 [heatflow    ]
trimmed: len=8 [heatflow]
typed size=5

if (.not. allocated(self%work)) allocate (self%work(nx, ny))

module field_store
  implicit none
  private
  public :: tracker, scope_demo, keep_demo, make_tracked, reset_tracked

  type :: tracker                        ! prints a line when an object of this type is finalized
    character(len=8) :: tag = '?'
  contains
    final :: tracker_final
  end type tracker

contains

  subroutine tracker_final(self)
    type(tracker), intent(inout) :: self
    print '(2a)', '    final subroutine runs for: ', trim(self%tag)
  end subroutine tracker_final

  subroutine scope_demo()
    type(tracker), allocatable :: local            ! unsaved local: unallocated on every call
    type(tracker), allocatable, save :: kept       ! saved: keeps its status between calls
    print '(a, l1, a, l1)', '  on entry: allocated(local)=', allocated(local), ' allocated(kept)=', allocated(kept)
    allocate (local, source=tracker('local'))
    if (.not. allocated(kept)) allocate (kept, source=tracker('kept'))
    print '(a)', '  leaving scope_demo'
  end subroutine scope_demo

  subroutine keep_demo()
    block
      type(tracker), allocatable :: in_block
      allocate (in_block, source=tracker('in_block'))
      print '(a)', '  end of the block construct'
    end block
    print '(a)', '  after the block'
  end subroutine keep_demo

  function make_tracked(tag) result(r)             ! the result keeps its allocation
    character(len=*), intent(in) :: tag
    type(tracker), allocatable :: r
    allocate (r, source=tracker(tag))
  end function make_tracked

  subroutine reset_tracked(t)
    type(tracker), allocatable, intent(out) :: t   ! an allocated actual argument is deallocated on entry
    print '(a, l1)', '    inside reset_tracked: allocated(t)=', allocated(t)
  end subroutine reset_tracked
end module field_store

program auto_dealloc
  use field_store
  implicit none
  type(tracker), allocatable :: t

  print '(a)', 'first call of scope_demo:'
  call scope_demo()
  print '(a)', 'second call of scope_demo:'
  call scope_demo()

  print '(a)', 'block construct:'
  call keep_demo()

  print '(a)', 'function result, then intent(out):'
  t = make_tracked('result')            ! the function result is a temporary, deallocated after the statement
  print '(a, l1)', '  t allocated: ', allocated(t)
  call reset_tracked(t)

  print '(a)', 'end of the main program'
end program auto_dealloc

first call of scope_demo:
  on entry: allocated(local)=F allocated(kept)=F
  leaving scope_demo
    final subroutine runs for: local
second call of scope_demo:
  on entry: allocated(local)=F allocated(kept)=T
  leaving scope_demo
    final subroutine runs for: local
block construct:
  end of the block construct
    final subroutine runs for: in_block
  after the block
function result, then intent(out):
    final subroutine runs for: result
  t allocated: T
    final subroutine runs for: result
    inside reset_tracked: allocated(t)=F
end of the main program

program reallocation_on_assignment
  use, intrinsic :: iso_fortran_env, only: real64
  implicit none
  integer, parameter :: dp = real64
  real(dp), allocatable :: x(:), history(:)
  real(dp) :: y(0:4), r
  integer :: step

  y = [10.0_dp, 20.0_dp, 30.0_dp, 40.0_dp, 50.0_dp]

  ! 1. assignment to an unallocated array allocates it with the shape of the expression
  x = y(1:3)
  print '(a, i0, a, i0, a, i0)', 'x = y(1:3):       size=', size(x), ' lbound=', lbound(x, 1), ' ubound=', ubound(x, 1)

  ! 2. a different shape reallocates; the new lower bound is lbound(expr)
  x = y
  print '(a, i0, a, i0, a, i0)', 'x = y:            size=', size(x), ' lbound=', lbound(x, 1), ' ubound=', ubound(x, 1)
  x = y(2:4)                            ! a section always has lower bound 1
  print '(a, i0, a, i0, a, i0)', 'x = y(2:4):       size=', size(x), ' lbound=', lbound(x, 1), ' ubound=', ubound(x, 1)

  ! 3. an array constructor has lower bound 1, a scalar keeps the current shape
  x = [x, 99.0_dp]                      ! the right-hand side is evaluated first, then x is reallocated
  print '(a, i0, a, f5.1)', 'x = [x, 99.0]:    size=', size(x), ' last=', x(size(x))
  x = 0.0_dp
  print '(a, i0, a, f5.1)', 'x = 0.0:          size=', size(x), ' x(1)=', x(1)

  ! 4. the "growing history" idiom: simple, but every step copies the whole array
  allocate (history(0))
  do step = 1, 5
    r = 1.0_dp/real(step, dp)
    history = [history, r]
  end do
  print '(a, i0, a, f6.4)', 'history: size=', size(history), ' last=', history(size(history))

  ! 5. x(:) = ... turns reallocation off: the shapes must already agree
  x(:) = 1.0_dp
  print '(a, i0)', 'x(:) = 1.0 keeps size ', size(x)
end program reallocation_on_assignment

x = y(1:3):       size=3 lbound=1 ubound=3
x = y:            size=5 lbound=0 ubound=4
x = y(2:4):       size=3 lbound=1 ubound=3
x = [x, 99.0]:    size=4 last= 99.0
x = 0.0:          size=4 x(1)=  0.0
history: size=5 last=0.2000
x(:) = 1.0 keeps size 4

! Deliberately wrong: a(:) = ... does not reallocate, so the shapes must already agree.
program section_assignment_mismatch
  implicit none
  real, allocatable :: a(:)
  a = [1.0, 2.0, 3.0]
  a(:) = [1.0, 2.0, 3.0, 4.0, 5.0]
  print *, a
end program section_assignment_mismatch

At line 6 of file section_mismatch.f90
Fortran runtime error: Array bound mismatch for dimension 1 of array 'a' (3/5)

fatal Fortran runtime error(section_mismatch.f90:6): Assign: mismatching element counts in array assignment (to 3, from 5)
Aborted

module residual_log
  use, intrinsic :: iso_fortran_env, only: real64
  implicit none
  private
  public :: dp, log_t

  integer, parameter :: dp = real64

  type :: log_t                          ! growable list of residuals with spare capacity
    real(dp), allocatable :: values(:)
    integer :: n = 0
    integer :: grow_count = 0
  contains
    procedure :: append => log_append
  end type log_t

contains

  subroutine log_append(self, r)
    class(log_t), intent(inout) :: self
    real(dp), intent(in) :: r
    real(dp), allocatable :: bigger(:)
    if (.not. allocated(self%values)) allocate (self%values(4))
    if (self%n == size(self%values)) then
      allocate (bigger(2*size(self%values)))
      bigger(1:self%n) = self%values(1:self%n)
      call move_alloc(from=bigger, to=self%values)   ! hands the new block over without a second copy
      self%grow_count = self%grow_count + 1
    end if
    self%n = self%n + 1
    self%values(self%n) = r
  end subroutine log_append

end module residual_log

program move_alloc_demo
  use residual_log
  implicit none
  type(log_t) :: history
  real(dp), allocatable, target :: a(:), b(:)
  real(dp), pointer :: p(:)
  integer :: k, stat
  character(len=100) :: msg

  do k = 1, 100
    call history%append(1.0_dp/real(k, dp))
  end do
  print '(a, i0, a, i0, a, i0)', 'n=', history%n, ' capacity=', size(history%values), ' reallocations=', history%grow_count

  ! what move_alloc does to the two variables
  allocate (a(0:2))
  a = [1.0_dp, 2.0_dp, 3.0_dp]
  allocate (b(10))
  p => a
  call move_alloc(a, b, stat=stat, errmsg=msg)
  print '(a, l1, a, l1, a, i0)', 'a allocated=', allocated(a), '  b allocated=', allocated(b), '  stat=', stat
  print '(a, 2i3, a, 3f5.1)', 'b took the bounds of a: lbound/ubound=', lbound(b, 1), ubound(b, 1), '  values', b
  print '(a, l1, a, f5.1)', 'associated(p, b)=', associated(p, b), '  p(0)=', p(0)   ! b has TARGET, so p follows the move

  ! moving from an unallocated variable leaves the destination unallocated
  call move_alloc(a, b)
  print '(a, l1)', 'after moving the unallocated a: allocated(b)=', allocated(b)
end program move_alloc_demo

n=100 capacity=128 reallocations=5
a allocated=F  b allocated=T  stat=0
b took the bounds of a: lbound/ubound=  0  2  values  1.0  2.0  3.0
associated(p, b)=T  p(0)=  1.0
after moving the unallocated a: allocated(b)=F

program stack_versus_heap
  use, intrinsic :: iso_fortran_env, only: real64
  implicit none
  integer, parameter :: dp = real64
  integer, parameter :: n = 50000000          ! 50 million reals of 8 bytes = 400 MB
  character(len=16) :: mode

  call get_command_argument(1, mode)
  select case (trim(mode))
  case ('automatic')
    call with_automatic(n)
  case ('allocatable')
    call with_allocatable(n)
  case default
    print '(a)', 'usage: a.out automatic|allocatable'
  end select

contains

  subroutine with_automatic(m)                 ! automatic array: size is a dummy argument
    integer, intent(in) :: m
    real(dp) :: work(m)
    work = 1.0_dp
    print '(a, es10.3)', 'automatic:   sum = ', sum(work)
  end subroutine with_automatic

  subroutine with_allocatable(m)
    integer, intent(in) :: m
    real(dp), allocatable :: work(:)
    integer :: stat
    allocate (work(m), stat=stat)
    if (stat /= 0) then
      print '(a)', 'allocatable: allocation failed, the program can report it and go on'
      return
    end if
    work = 1.0_dp
    print '(a, es10.3)', 'allocatable: sum = ', sum(work)
  end subroutine with_allocatable

end program stack_versus_heap

module heatflow_scratch
  implicit none
contains
  subroutine smooth(u)
    real, intent(inout) :: u(:)
    real :: scratch(100000)          ! fixed-size local array of 400 000 bytes
    integer :: n
    n = min(size(u), size(scratch))
    scratch(1:n) = u(1:n)
    u(2:n - 1) = 0.5*(scratch(1:n - 2) + scratch(3:n))
  end subroutine smooth
end module heatflow_scratch

program big_local
  use heatflow_scratch
  implicit none
  real :: u(5) = [0.0, 1.0, 2.0, 3.0, 4.0]
  call smooth(u)
  print '(5f5.1)', u
end program big_local

big_local.f90:6:13:

    6 |     real :: scratch(100000)          ! fixed-size local array of 400 000 bytes
      |             1~~~~~~~~~~~~~~
Warning: Array 'scratch' at (1) is larger than limit set by '-fmax-stack-var-size=', moved from stack to static storage. This makes the procedure unsafe when called recursively, or concurrently from multiple threads. Consider increasing the '-fmax-stack-var-size=' limit (or use '-frecursive', which implies unlimited '-fmax-stack-var-size') - or change the code to use an ALLOCATABLE array. If the variable is never accessed concurrently, this warning can be ignored, and the variable could also be declared with the SAVE attribute. [-Wsurprising]
  0.0  1.0  2.0  3.0  4.0

program leaks
  use, intrinsic :: iso_fortran_env, only: real64
  implicit none
  integer, parameter :: dp = real64
  character(len=16) :: which

  call get_command_argument(1, which)
  select case (trim(which))
  case ('allocatable')
    call with_allocatable()
  case ('pointer')
    call with_pointer()
  case default
    print '(a)', 'usage: a.out allocatable|pointer'
  end select

contains

  subroutine with_allocatable()
    real(dp), allocatable :: work(:)
    integer :: k
    do k = 1, 3
      allocate (work(1000))        ! allocating an allocated variable is an error, so every pass must free it
      work = real(k, dp)
      deallocate (work)
    end do
    allocate (work(1000))          ! never deallocated explicitly: freed automatically on return
    work = 0.0_dp
  end subroutine with_allocatable

  subroutine with_pointer()
    real(dp), pointer :: work(:)
    integer :: k
    do k = 1, 3
      allocate (work(1000))        ! the previous block is silently lost: pointers may be re-allocated
      work = real(k, dp)
    end do
  end subroutine with_pointer

end program leaks

$ valgrind --leak-check=full ./a.out allocatable
==12==   in use at exit: 0 bytes in 0 blocks
==12==   total heap usage: 22 allocs, 22 frees, 37,575 bytes allocated
==12== All heap blocks were freed -- no leaks are possible

$ valgrind --leak-check=full ./a.out pointer
==14==   in use at exit: 24,000 bytes in 3 blocks
==14==   total heap usage: 21 allocs, 18 frees, 29,571 bytes allocated
==14== 24,000 bytes in 3 blocks are definitely lost in loss record 1 of 1
==14== LEAK SUMMARY:
==14==    definitely lost: 24,000 bytes in 3 blocks
==14==    indirectly lost: 0 bytes in 0 blocks

! Deliberately wrong: each case triggers one runtime error; the case is chosen on the command line.
program allocation_runtime_errors
  implicit none
  real, allocatable :: a(:)
  character(len=16) :: which

  call get_command_argument(1, which)
  select case (trim(which))
  case ('reference')            ! reading an unallocated array
    print *, a(1)
  case ('allocate-twice')       ! allocate without stat= on an allocated array
    allocate (a(3))
    allocate (a(4))
  case ('deallocate-twice')     ! deallocate without stat= on an unallocated array
    allocate (a(3))
    deallocate (a)
    deallocate (a)
  case ('bounds')               ! subscript beyond the allocated extent
    allocate (a(3))
    a(4) = 1.0
  case ('whole')                ! whole-array use of an unallocated array
    print *, sum(a)
  case default
    print '(a)', 'usage: a.out reference|allocate-twice|deallocate-twice|bounds|whole'
  end select
end program allocation_runtime_errors

$ ./a.out reference
At line 10 of file runtime_errors.f90
Fortran runtime error: Index '1' of dimension 1 of array 'a' below lower bound of 281474053255328
$ ./a.out allocate-twice
At line 13 of file runtime_errors.f90
Fortran runtime error: Attempting to allocate already allocated variable 'a'
$ ./a.out deallocate-twice
At line 17 of file runtime_errors.f90
Fortran runtime error: Attempt to DEALLOCATE unallocated 'a'
$ ./a.out bounds
At line 20 of file runtime_errors.f90
Fortran runtime error: Index '4' of dimension 1 of array 'a' above upper bound of 3
$ ./a.out whole
Program received signal SIGSEGV: Segmentation fault - invalid memory reference.

$ ./a.out reference
Segmentation fault
$ ./a.out allocate-twice
(no output, exit status 0)
$ ./a.out deallocate-twice
(no output, exit status 0)
$ ./a.out bounds
(no output, exit status 0)
$ ./a.out whole
 0.