tem_trackmem Subroutine

public subroutine tem_trackmem(memfile, iteration)

Write the current memory status into the memfile.

The file will be opened and closed, so this might be slow. Only the root process writes to this file, but the data is gathered from all processes. The output will be prepended by the current date. We will track the VmHWM: from /proc/self/status, MemFree:, Buffers: and Cached: from /proc/meminfo and pgmajfault from /proc/vmstat The min, max and average across all processes will be recorded.

Arguments

Type IntentOptional Attributes Name
character(len=*), intent(in) :: memfile
integer, intent(in) :: iteration

Source Code

  subroutine tem_trackmem(memfile, iteration)
    character(len=*), intent(in) :: memfile
    integer, intent(in)          :: iteration

    integer :: myMem(5), minMem(5), maxMem(5)
    integer :: nProcs
    integer :: iError
    integer :: funit
    integer :: myRank
    logical :: fexists
    character(len=8) :: fstat
    character(len=8) :: today
    character(len=10) :: now
    real :: myScale(5)
    real :: sumMem(5)
    real :: avgMem(5)

    if (trim(memfile) /= '') then
      call MPI_Comm_Size(MPI_COMM_WORLD, nProcs, iError)

      myMem(1) = my_status_int('VmHWM:')
      myMem(2:4) = my_status_int_vec( key = ['MemFree:', 'Buffers:', &
        &                                    'Cached: '],            &
        &                             info = '/proc/meminfo'         )
      myMem(5) = my_status_int(key = 'pgmajfault', info = '/proc/vmstat')
      call MPI_Reduce(myMem, minMem, 5, MPI_INTEGER, MPI_MIN, 0, &
        &             MPI_COMM_WORLD, iError)
      call MPI_Reduce( myMem, maxMem, 5, MPI_INTEGER, MPI_MAX, 0, &
        &              MPI_COMM_WORLD, iError )

      myScale(1:4) = real(myMem(1:4))/1024.0
      myScale(5) = real(myMem(5))
      call MPI_Reduce( myScale, sumMem, 5, MPI_REAL, MPI_SUM, 0, &
        &              MPI_COMM_WORLD, iError                    )

      avgMem = sumMem / real(nProcs)

      inquire(file = trim(memfile), exist = fexists)
      if (fexists) then
        fstat = 'old'
      else
        fstat = 'new'
      end if

      call MPI_Comm_rank(MPI_COMM_WORLD, myRank, iError)
      if (myRank == 0) then
        call tem_open( file     = trim(memfile), &
          &            newunit  = funit,         &
          &            status   = trim(fstat),   &
          &            action   = 'write',       &
          &            position = 'append',      &
          &            form     = 'formatted'    )
        call date_and_time(date = today, time = now)
        write(funit,'(a)') today // ': ' // now
        write(funit, '(a,I5)') '  iteration ', iteration
        write(funit,'(a,3(1x,f11.3))') '   VmHWM   (min/max/avg):', &
          &                                  minMem(1)/1024.0, &
          &                                  maxMem(1)/1024.0, avgMem(1)
        write(funit,'(a,3(1x,f11.3))') '   MemFree (min/max/avg):', &
          &                                  minMem(2)/1024.0, &
          &                                  maxMem(2)/1024.0, avgMem(2)
        write(funit,'(a,3(1x,f11.3))') '   Buffers (min/max/avg):', &
          &                                  minMem(3)/1024.0, &
          &                                  maxMem(3)/1024.0, avgMem(3)
        write(funit,'(a,3(1x,f11.3))') '   Cached  (min/max/avg):', &
          &                                  minMem(4)/1024.0, &
          &                                  maxMem(4)/1024.0, avgMem(4)
        write(funit,'(a,i0,1x,i0,1x,f11.3)') '   pgmajflt(min/max/avg):', &
          &                                  minMem(5), maxMem(5), avgMem(5)
        write(funit,'(a)') ''
        close(funit)
      end if

    end if

  end subroutine tem_trackmem