Changes between Version 3 and Version 4 of 2011WP/2011Stream2/Dynamic Memory


Ignore:
Timestamp:
2010-11-15T14:08:39+01:00 (10 years ago)
Author:
trackstand
Comment:

Legend:

Unmodified
Added
Removed
Modified
  • 2011WP/2011Stream2/Dynamic Memory

    v3 v4  
    11= Discussion of coding approach/style for dynamic memory = 
    22 
    3 As a basis for the discussion, here's how I've currently coded NEMO (v.3.2) to use dynamic memory. 
     3As a basis for the discussion, here's how I've currently coded NEMO (v.3.2) to use dynamic memory. I've used !''...!'' to indicate places where I've missed-out chunks of code for clarity/brevity. 
    44 
    55Highest-level changes are in OPA_SRC/opa.F90 and the routine opa_init(): 
     
    146146   END SUBROUTINE opa_alloc 
    147147}}} 
     148 
     149Each of the modules must then be edited and the declaration of arrays with dimensions that depend upon jpi or jpj be changed. e.g. for OPA_SRC/DOM/dom_oce.F90: 
     150{{{ 
     151MODULE dom_oce 
     152   !!====================================================================== 
     153   !!                       ***  MODULE dom_oce  *** 
     154   !!        
     155   !! ** Purpose :   Define in memory all the ocean space domain variables 
     156   !!====================================================================== 
     157   !! History :  1.0  ! 2005-10  (A. Beckmann, G. Madec)  reactivate s-coordinate  
     158   !!---------------------------------------------------------------------- 
     159   !! NEMO/OPA 3.2 , LOCEAN-IPSL (2009)  
     160   !! $Id: dom_oce.F90 1601 2009-08-11 10:09:19Z ctlod $  
     161   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     162   !!---------------------------------------------------------------------- 
     163   USE par_oce      ! ocean parameters 
     164 
     165   IMPLICIT NONE 
     166... 
     167... 
     168#if defined key_mpp_dyndist 
     169   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: & !:  
     170#else 
     171   REAL(wp), PUBLIC,              DIMENSION(jpi,jpj,jpk) :: & !: 
     172#endif 
     173      tmask, umask,    &  !: land/ocean mask at T-, U-, V- and F-points 
     174      vmask, fmask        !: 
     175 
     176   REAL(wp), PUBLIC, DIMENSION(jpiglo) ::   tpol, fpol          !: north fold mask (nperio= 3 or 4) 
     177 
     178#if defined key_noslip_accurate 
     179   INTEGER, PUBLIC,                    DIMENSION(4,jpk) ::   npcoa   !: ??? 
     180#if defined key_mpp_dyndist 
     181   INTEGER, PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   nicoa, njcoa        !: ??? 
     182#else 
     183   INTEGER, PUBLIC,        DIMENSION(2*(jpi+jpj),4,jpk) ::   nicoa, njcoa        !: ??? 
     184#endif 
     185                                          
     186#endif 
     187 
     188   !!---------------------------------------------------------------------- 
     189   !! agrif domain 
     190   !!---------------------------------------------------------------------- 
     191#if defined key_agrif 
     192   LOGICAL, PUBLIC, PARAMETER ::   lk_agrif = .TRUE.    !: agrif flag 
     193#else 
     194   LOGICAL, PUBLIC, PARAMETER ::   lk_agrif = .FALSE.   !: agrif flag 
     195#endif 
     196 
     197 
     198#if defined key_mpp_dyndist 
     199   PUBLIC dom_oce_malloc 
     200#endif 
     201 
     202#if defined key_mpp_dyndist 
     203 
     204CONTAINS 
     205 
     206  FUNCTION dom_oce_malloc() 
     207    USE par_oce, Only: jpi, jpj, jpk, jpnij 
     208    IMPLICIT none 
     209    INTEGER :: dom_oce_malloc 
     210    INTEGER :: ierr, ierr_tmp 
     211     
     212    ierr = 0 
     213... 
     214... 
     215    ALLOCATE(tmask(jpi,jpj,jpk), umask(jpi,jpj,jpk),    &  
     216             vmask(jpi,jpj,jpk), fmask(jpi,jpj,jpk), Stat=ierr_tmp) 
     217    ierr = ierr + ierr_tmp 
     218 
     219#if defined key_noslip_accurate 
     220    ALLOCATE(nicoa(2*(jpi+jpj),4,jpk), njcoa(2*(jpi+jpj),4,jpk), & 
     221             Stat=ierr_tmp) 
     222    ierr = ierr + ierr_tmp 
     223#endif 
     224 
     225    IF(ierr != 0)THEN 
     226      WRITE() "Some error message detailing which module's ..._malloc() routine failed" 
     227    END IF 
     228 
     229    ! Error status passed back up 
     230    dom_oce_malloc = ierr 
     231 
     232  END FUNCTION dom_oce_malloc 
     233#endif 
     234}}} 
     235 
     236Finally (I think), there's the issue of work-space or local arrays. These are arrays local to subroutines that are declared with their dimensions in terms of jpi and jpj. Since jpi and jpj are no longer parameters, these arrays are now 'automatic' -- that is, they are allocated on the stack when the program enters the subroutine and free'd again when it leaves. Running out of stack memory will cause the program to abort abruptly and therefore it is more controlled to use allocate for big workspace arrays. (Although the amount of available stack can be set to 'unlimited' by using the shell's '''''ulimit''''' command, what limit this actually results in will be OS/set-up dependent.) 
     237I suggest therefore that we edit the code such that any large (3d or higher?) workspace arrays are made explicitly allocatable and then allocated the first time the subroutine is called. If we declare them with the SAVE attribute then the allocated memory will remain available to us in future calls, ''e.g.'': 
     238{{{ 
     239   SUBROUTINE dia_ptr( kt ) 
     240      !!---------------------------------------------------------------------- 
     241      !!                  ***  ROUTINE dia_ptr  *** 
     242      !!---------------------------------------------------------------------- 
     243      IMPLICIT none 
     244      INTEGER, INTENT(in) ::   kt   ! ocean time step index 
     245      !! 
     246      INTEGER  ::   jk, jj, ji   ! dummy loop 
     247      REAL(wp) ::   zsverdrup    ! conversion from m3/s to Sverdrup 
     248      REAL(wp) ::   zpwatt       ! conversion from W    to PW 
     249      REAL(wp) ::   zggram       ! conversion from g    to Pg 
     250#if defined key_mpp_dyndist 
     251      REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:) :: vt, vs   ! 3D workspace 
     252#else 
     253      REAL(wp),              DIMENSION(jpi,jpj,jpk) :: vt, vs   ! 3D workspace 
     254#endif 
     255      !!---------------------------------------------------------------------- 
     256 
     257      IF( kt == nit000 .OR. MOD( kt, nf_ptr ) == 0 )   THEN 
     258 
     259#if defined key_mpp_dyndist 
     260         IF(.not. ALLOCATED(vt))THEN 
     261            ALLOCATE(vt(jpi,jpj,jpk), vs(jpi,jpj,jpk), Stat=ji) 
     262            IF(ji .ne. 0)THEN 
     263               ! Need an interface to MPI_Abort() here really 
     264               STOP 'ERROR: Allocation of vt/vs failed in dia_ptr' 
     265            END IF 
     266         END IF 
     267#endif 
     268... 
     269... 
     270}}} 
     271If an allocate() is done in this way then handling any failure is difficult since there's no guarantee that the allocate() call will have failed for all MPI processes. Therefore, I think the best bet is to print out an error message and then do an MPI_Abort on MPI_COMM_WORLD.