Discussion of coding approach/style for dynamic memory

As a basis for the discussion, here's how I've currently coded NEMO (v.3.2) to use dynamic memory.

Highest-level changes are in OPA_SRC/opa.F90 and the routine opa_init():

#if defined key_mpp_dyndist
      ! Decide on size of grid now that we have our communicator size
      ! If we're not using dynamic memory then mpp_partition does nothing.

#if   defined key_mpp_mpi   ||   defined key_mpp_shmem
      CALL opa_partition()
#else
      jpni = 1
      jpnj = 1
      jpnij = jpni*jpnj
#endif
      ! Calculate domain dimensions given calculated jpi and jpj
      jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci !: first  dim.
      jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj !: second dim.
      jpim1 = jpi-1                                          !: inner domain indices
      jpjm1 = jpj-1                                          !:   "           "
      jpkm1 = jpk-1                                          !:   "           "
      jpij  = jpi*jpj                                        !:  jpi x j

      ! Now we know the dimensions of the grid, allocate arrays
      CALL opa_alloc()
#endif

The inclusion of the new code is currently controlled by the cpp key 'key_mpp_dyndist' but it would greatly improve the cleanliness of the code if we make a complete break from the static-memory version and thus can drop the use of this cpp key. My understanding of the conclusion of the meeting in Southampton was that this is what we're going to do.

This addition to opa_init() calls two further new routines, opa_partition() and opa_alloc():

   SUBROUTINE opa_partition
     USE par_oce
     IMPLICIT none
     INTEGER, PARAMETER :: nfactmax = 20
     INTEGER :: nfact ! The no. of factors returned
     INTEGER :: ierr  ! Error flag
     INTEGER :: i
     INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are
                                     ! closest in value
     INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors
     ierr = 0

#if ! defined key_mpp_dyndist
     RETURN ! If we aren't using dynamic memory then jpnj and jpni are set
            ! at compile time
#else

     CALL factorise(ifact, nfactmax, nfact, mppsize, ierr)

     IF(nfact <= 1)THEN
        WRITE (numout, *) 'WARNING: factorisation of number of PEs failed'
        WRITE (numout, *) '       : using grid of ',mppsize,' x 1'
        jpnj = 1
        jpni = mppsize
     ELSE
        ! Search through factors for the pair that are closest in value
        mindiff = 1000000
        imin    = 1
        DO i=1,nfact-1,2
           idiff = ABS(ifact(i) - ifact(i+1))
           IF(idiff < mindiff)THEN
              mindiff = idiff
              imin = i
           END IF
        END DO
        jpnj = ifact(imin)
        jpni = ifact(imin + 1)
     ENDIF
     jpnij = jpni*jpnj

     WRITE(*,*) 'ARPDBG: jpni = ',jpni,'jpnj = ',jpnj,'jpnij = ',jpnij

#endif

   END SUBROUTINE opa_partition

where factorise() returns the prime factors of its argument.

opa_alloc() oversees the allocation of all of the module-level arrays:

   SUBROUTINE opa_alloc
#if defined key_mpp_dyndist

#if defined key_lim2
     USE limwri_2,   ONLY: lim_wri_2_malloc
     USE limdmp_2,   ONLY: lim_dmp_2_malloc
     USE limhdf_2,   ONLY: lim_hdf_2_malloc
     USE ice_2,      ONLY: ice_2_malloc
     USE limsbc_2,   ONLY: lim_sbc_2_malloc
     USE thd_ice_2,  ONLY: thd_ice_2_malloc
     USE limdia_2,   ONLY: lim_dia_2_malloc
     USE dom_ice_2,  ONLY: dom_ice_2_malloc
#endif
     USE ldfslp,     ONLY: ldf_slp_malloc
     USE ldftra_oce, ONLY: ldftra_oce_malloc
     USE ldfdyn_oce, ONLY: ldfdyn_oce_malloc
#if defined key_vvl
     USE domvvl,     ONLY: dom_vvl_malloc
#endif
     USE dom_oce,    ONLY: dom_oce_malloc
...
...
     IMPLICIT none
     INTEGER, PARAMETER           :: NUMCALLS = 58
     INTEGER, DIMENSION(NUMCALLS) :: ierr
     INTEGER                      :: i, icount

     ierr = 0
     icount = 1

#if defined key_lim2
     ierr = ierr + lim_wri_2_malloc()
     ierr = ierr + lim_dmp_2_malloc()
     ierr = ierr + lim_hdf_2_malloc()
     ierr = ierr + ice_2_malloc()
     ierr = ierr + lim_sbc_2_malloc()
     ierr = ierr + thd_ice_2_malloc()
     ierr = ierr + lim_dia_2_malloc()
     ierr = ierr + dom_ice_2_malloc()
#endif
     ierr = ierr + ldf_slp_malloc()
     ierr = ierr + ldftra_oce_malloc()
     ierr = ierr + ldfdyn_oce_malloc()
#if defined key_vvl
     ierr = ierr + dom_vvl_malloc()
#endif
     ierr = ierr + dom_oce_malloc()
...
...
! Should do an MPI_SUM on ierr and then everyone can MPI_FINALIZE() if the
! memory allocation has failed on any one PE
     IF(ierr > 0)THEN
        WRITE(numout,*) 
        WRITE(numout,*) 'ERROR: Allocation of memory failed in opa_alloc'
        STOP 'ERROR: Allocation of memory failed in opa_alloc'
     END IF
#else
     RETURN ! Not using dynamic memory therefore do nothing
#endif

   END SUBROUTINE opa_alloc