New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 2648 – NEMO

Changeset 2648


Ignore:
Timestamp:
2011-03-03T17:13:18+01:00 (13 years ago)
Author:
cetlod
Message:

Changed OFF_SRC component to use dynamic memory

Location:
branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC
Files:
1 added
4 edited

Legend:

Unmodified
Added
Removed
  • branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/dommsk.F90

    r2528 r2648  
    1212   USE oce             ! ocean dynamics and tracers 
    1313   USE dom_oce         ! ocean space and time domain 
    14    USE in_out_manager  ! I/O manager 
     14   USE lib_mpp 
    1515 
    1616   IMPLICIT NONE 
     
    1919   PUBLIC   dom_msk    ! routine called by inidom.F90 
    2020 
    21 #if defined key_degrad 
    22    !! ------------------------------------------------ 
    23    !! Degradation method 
    24    !! -------------------------------------------------- 
    25    REAL(wp), PUBLIC, DIMENSION (jpi,jpj,jpk) ::   facvol  !! volume for degraded regions 
    26 #endif 
     21   REAL(wp), PUBLIC, ALLOCATABLE, SAVE, DIMENSION(:,:,:) ::   facvol  !! volume for degraded regions 
    2722 
    2823   !! * Substitutions 
     
    3429   !!---------------------------------------------------------------------- 
    3530CONTAINS 
     31 
    3632    
    3733   SUBROUTINE dom_msk 
     
    4945      !!               tpol     : ??? 
    5046      !!---------------------------------------------------------------------- 
     47      USE wrk_nemo, ONLY: iwrk_in_use, iwrk_not_released 
     48      USE wrk_nemo, ONLY: imsk => iwrk_2d_1 
    5149      INTEGER  ::   ji, jk                   ! dummy loop indices 
    5250      INTEGER  ::   iif, iil, ijf, ijl       ! local integers 
    53       INTEGER, DIMENSION(jpi,jpj) ::  imsk   ! 2D workspace 
    5451      !!--------------------------------------------------------------------- 
    5552      ! 
     53      IF( iwrk_in_use(2, 1) ) THEN 
     54         CALL ctl_stop('dom_msk: ERROR: requested workspace arrays unavailable.')  ;  RETURN 
     55      END IF 
     56      ! 
     57      CALL dom_msk_alloc 
     58 
    5659      ! Interior domain mask (used for global sum) 
    5760      ! -------------------- 
     
    9598      ENDIF 
    9699      ! 
     100      IF( iwrk_not_released(2, 1) ) CALL ctl_stop('dom_msk: failed to release workspace arrays.') 
     101      ! 
    97102   END SUBROUTINE dom_msk 
     103 
     104   SUBROUTINE dom_msk_alloc() 
     105      !!--------------------------------------------------------------------- 
     106      !!                 ***  ROUTINE dom_msk_alloc  *** 
     107      !!--------------------------------------------------------------------- 
     108#if defined key_degrad 
     109      INTEGER :: ierr  
     110 
     111      ALLOCATE( facvol(jpi,jpj,jpk), STAT=ierr ) 
     112      IF( ierr /= 0 )  & 
     113        &           CALL ctl_stop('STOP', 'dom_msk : unable to allocate facvol array') 
     114#endif 
     115 
     116   END SUBROUTINE dom_msk_alloc 
    98117 
    99118   !!====================================================================== 
  • branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/domrea.F90

    r2528 r2648  
    1616   USE dommsk          ! domain: masks 
    1717   USE lbclnk          ! lateral boundary condition - MPP exchanges 
    18    USE in_out_manager  ! I/O manager 
     18   USE lib_mpp  
     19   USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 
    1920 
    2021   IMPLICIT NONE 
     
    5354      !!---------------------------------------------------------------------- 
    5455      USE iom 
     56      USE wrk_nemo, ONLY: zmbk => wrk_2d_1, zprt => wrk_2d_2, zprw => wrk_2d_3 
    5557      !! 
    5658      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
    5759      INTEGER  ::   ik, inum0 , inum1 , inum2 , inum3 , inum4   ! local integers 
    5860      REAL(wp) ::   zrefdep         ! local real 
    59       REAL(wp), DIMENSION(jpi,jpj) ::   zmbk, zprt, zprw   ! 2D workspace 
    6061      !!---------------------------------------------------------------------- 
    6162 
     
    6364      IF(lwp) WRITE(numout,*) 'dom_rea : read NetCDF mesh and mask information file(s)' 
    6465      IF(lwp) WRITE(numout,*) '~~~~~~~' 
     66 
     67      IF( wrk_in_use(2, 1,2,3)  ) THEN 
     68         CALL ctl_stop('dom_rea: ERROR: requested workspace arrays unavailable.') ; RETURN 
     69      END IF 
    6570 
    6671      zmbk(:,:) = 0._wp 
     
    141146         CALL iom_get( inum3, jpdom_data, 'e2u', e2u ) 
    142147         CALL iom_get( inum3, jpdom_data, 'e2v', e2v ) 
     148 
     149         e1e2t(:,:) = e1t(:,:) * e2t(:,:)                              ! surface at T-points 
    143150 
    144151         CALL iom_get( inum3, jpdom_data, 'ff', ff ) 
     
    314321      END SELECT 
    315322      ! 
     323      IF( wrk_not_released(2, 1,2,3)  ) CALL ctl_stop('dom_rea:failed to release workspace arrays.') 
     324      ! 
    316325   END SUBROUTINE dom_rea 
    317326 
     
    327336      !! ** Action  : - update mbathy: level bathymetry (in level index) 
    328337      !!---------------------------------------------------------------------- 
     338      USE wrk_nemo, ONLY: zmbk => wrk_2d_4 
     339      ! 
    329340      INTEGER ::   ji, jj   ! dummy loop indices 
    330       REAL(wp), DIMENSION(jpi,jpj) ::   zmbk   ! 2D workspace  
    331       !!---------------------------------------------------------------------- 
     341      !!---------------------------------------------------------------------- 
     342 
    332343      ! 
    333344      IF(lwp) WRITE(numout,*) 
    334345      IF(lwp) WRITE(numout,*) '    zgr_bot_level : ocean bottom k-index of T-, U-, V- and W-levels ' 
    335346      IF(lwp) WRITE(numout,*) '    ~~~~~~~~~~~~~' 
     347      ! 
     348      IF( wrk_in_use(2, 4) ) THEN 
     349         CALL ctl_stop('dom_rea: ERROR: requested workspace arrays unavailable.')  ;  RETURN 
     350      END IF 
    336351      ! 
    337352      mbkt(:,:) = MAX( mbathy(:,:) , 1 )    ! bottom k-index of T-level (=1 over land) 
     
    347362      zmbk(:,:) = REAL( mbkv(:,:), wp )   ;   CALL lbc_lnk(zmbk,'V',1.)   ;   mbkv  (:,:) = MAX( INT( zmbk(:,:) ), 1 ) 
    348363      ! 
     364      IF( wrk_not_released(2, 4) ) CALL ctl_stop('dom_rea:failed to release workspace arrays.') 
     365      ! 
    349366   END SUBROUTINE zgr_bot_level 
    350367 
  • branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/dtadyn.F90

    r2559 r2648  
    6363   INTEGER ::   numfl_t, numfl_u, numfl_v, numfl_w 
    6464 
    65    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: tdta       ! temperature at two consecutive times 
    66    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: sdta       ! salinity at two consecutive times 
    67    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: udta       ! zonal velocity at two consecutive times 
    68    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: vdta       ! meridional velocity at two consecutive times 
    69    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wdta       ! vertical velocity at two consecutive times 
    70    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: avtdta     ! vertical diffusivity coefficient 
    71  
    72    REAL(wp), DIMENSION(jpi,jpj    ,2) :: hmlddta    ! mixed layer depth at two consecutive times 
    73    REAL(wp), DIMENSION(jpi,jpj    ,2) :: wspddta    ! wind speed at two consecutive times 
    74    REAL(wp), DIMENSION(jpi,jpj    ,2) :: frlddta    ! sea-ice fraction at two consecutive times 
    75    REAL(wp), DIMENSION(jpi,jpj    ,2) :: empdta     ! E-P at two consecutive times 
    76    REAL(wp), DIMENSION(jpi,jpj    ,2) :: qsrdta     ! short wave heat flux at two consecutive times 
    77    REAL(wp), DIMENSION(jpi,jpj    ,2) :: bblxdta    ! frequency of bbl in the x direction at 2 consecutive times  
    78    REAL(wp), DIMENSION(jpi,jpj    ,2) :: bblydta    ! frequency of bbl in the y direction at 2 consecutive times  
     65   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: tdta       ! temperature at two consecutive times 
     66   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: sdta       ! salinity at two consecutive times 
     67   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: udta       ! zonal velocity at two consecutive times 
     68   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: vdta       ! meridional velocity at two consecutive times 
     69   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wdta       ! vertical velocity at two consecutive times 
     70   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: avtdta     ! vertical diffusivity coefficient 
     71 
     72   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: hmlddta    ! mixed layer depth at two consecutive times 
     73   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: wspddta    ! wind speed at two consecutive times 
     74   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: frlddta    ! sea-ice fraction at two consecutive times 
     75   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: empdta     ! E-P at two consecutive times 
     76   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: qsrdta     ! short wave heat flux at two consecutive times 
     77   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: bblxdta    ! frequency of bbl in the x direction at 2 consecutive times  
     78   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: bblydta    ! frequency of bbl in the y direction at 2 consecutive times  
    7979   LOGICAL :: l_offbbl 
    8080#if defined key_ldfslp 
    81    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: uslpdta    ! zonal isopycnal slopes 
    82    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: vslpdta    ! meridional isopycnal slopes 
    83    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wslpidta   ! zonal diapycnal slopes 
    84    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: wslpjdta   ! meridional diapycnal slopes 
     81   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: uslpdta    ! zonal isopycnal slopes 
     82   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: vslpdta    ! meridional isopycnal slopes 
     83   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wslpidta   ! zonal diapycnal slopes 
     84   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: wslpjdta   ! meridional diapycnal slopes 
    8585#endif 
    8686#if ! defined key_degrad &&  defined key_traldf_c2d && defined key_traldf_eiv  
    87    REAL(wp), DIMENSION(jpi,jpj    ,2) :: aeiwdta    ! G&M coefficient 
     87   REAL(wp), ALLOCATABLE, SAVE,   DIMENSION(:,:,:) :: aeiwdta    ! G&M coefficient 
    8888#endif 
    8989#if defined key_degrad 
    90    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: ahtudta, ahtvdta, ahtwdta   ! Lateral diffusivity 
     90   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: ahtudta, ahtvdta, ahtwdta   ! Lateral diffusivity 
    9191# if defined key_traldf_eiv 
    92    REAL(wp), DIMENSION(jpi,jpj,jpk,2) :: aeiudta, aeivdta, aeiwdta   ! G&M coefficient 
     92   REAL(wp), ALLOCATABLE, SAVE, DIMENSION(:,:,:,:) :: aeiudta, aeivdta, aeiwdta   ! G&M coefficient 
    9393# endif 
    9494#endif 
     
    296296   END SUBROUTINE dta_dyn 
    297297 
     298   INTEGER FUNCTION dta_dyn_alloc() 
     299      !!--------------------------------------------------------------------- 
     300      !!                 ***  ROUTINE dta_dyn_alloc  *** 
     301      !!--------------------------------------------------------------------- 
     302 
     303      ALLOCATE( tdta    (jpi,jpj,jpk,2), sdta    (jpi,jpj,jpk,2),    & 
     304       &        udta    (jpi,jpj,jpk,2), vdta    (jpi,jpj,jpk,2),    & 
     305       &        wdta    (jpi,jpj,jpk,2), avtdta  (jpi,jpj,jpk,2),    & 
     306#if defined key_ldfslp  
     307       &        uslpdta (jpi,jpj,jpk,2), vslpdta (jpi,jpj,jpk,2),    & 
     308       &        wslpidta(jpi,jpj,jpk,2), wslpjdta(jpi,jpj,jpk,2),    & 
     309#endif 
     310#if defined key_degrad 
     311       &        ahtudta (jpi,jpj,jpk,2), ahtvdta (jpi,jpj,jpk,2),    & 
     312       &        ahtwdta (jpi,jpj,jpk,2),                             & 
     313# if defined key_traldf_eiv 
     314       &        aeiudta (jpi,jpj,jpk,2), aeivdta (jpi,jpj,jpk,2),    & 
     315       &        aeiwdta (jpi,jpj,jpk,2),                             & 
     316# endif 
     317#endif 
     318#if ! defined key_degrad &&  defined key_traldf_c2d && defined key_traldf_eiv 
     319       &        aeiwdta (jpi,jpj,    2),                             & 
     320#endif 
     321 
     322       &        hmlddta (jpi,jpj,    2), wspddta (jpi,jpj,    2),    & 
     323       &        frlddta (jpi,jpj,    2), qsrdta  (jpi,jpj,    2),    & 
     324       &        empdta  (jpi,jpj,    2),                         STAT=dta_dyn_alloc )  
     325 
     326      IF( dta_dyn_alloc /= 0 ) CALL ctl_warn('dta_dyn_alloc: failed to allocate facvol array.') 
     327 
     328   END FUNCTION dta_dyn_alloc 
    298329 
    299330   SUBROUTINE dynrea( kt, kenr ) 
     
    305336      !! ** Method : READ the kenr records of DATA and store in udta(...,2), ....   
    306337      !!---------------------------------------------------------------------- 
     338      USE wrk_nemo, ONLY: wrk_in_use, wrk_not_released 
     339      USE wrk_nemo, ONLY: zu    => wrk_3d_1 , zv    => wrk_3d_2 , zw    => wrk_3d_3 
     340      USE wrk_nemo, ONLY: zt    => wrk_3d_4 , zs    => wrk_3d_5 
     341      USE wrk_nemo, ONLY: zavt  => wrk_3d_6 , zhdiv => wrk_3d_7 
     342      USE wrk_nemo, ONLY: zahtu => wrk_3d_8 , zahtv => wrk_3d_9 , zahtw => wrk_3d_10 
     343      USE wrk_nemo, ONLY: zaeiu => wrk_3d_11, zaeiv => wrk_3d_12, zaeiw => wrk_3d_13 
     344      ! 
     345      USE wrk_nemo, ONLY: zemp  => wrk_2d_1 , zqsr  => wrk_2d_2 , zmld  => wrk_2d_3 
     346      USE wrk_nemo, ONLY: zice  => wrk_2d_4 , zwspd => wrk_2d_5  
     347      USE wrk_nemo, ONLY: ztaux => wrk_2d_6 , ztauy => wrk_2d_7 
     348      USE wrk_nemo, ONLY: zbblx => wrk_2d_8 , zbbly => wrk_2d_9 
     349      USE wrk_nemo, ONLY: zaeiw2d => wrk_2d_10 
     350      ! 
    307351      INTEGER, INTENT(in) ::   kt, kenr   ! time index 
    308352      !! 
    309353      INTEGER ::  jkenr 
    310       REAL(wp), DIMENSION(jpi,jpj,jpk) ::  zu, zv, zw, zt, zs, zavt , zhdiv              ! 3D workspace 
    311       REAL(wp), DIMENSION(jpi,jpj)     ::  zemp, zqsr, zmld, zice, zwspd, ztaux, ztauy   ! 2D workspace 
    312       REAL(wp), DIMENSION(jpi,jpj)     ::  zbblx, zbbly 
    313  
    314 #if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 
    315       REAL(wp), DIMENSION(jpi,jpj) :: zaeiw  
    316 #endif 
    317 #if defined key_degrad 
    318    REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zahtu, zahtv, zahtw  !  Lateral diffusivity 
    319 # if defined key_traldf_eiv 
    320    REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zaeiu, zaeiv, zaeiw  ! G&M coefficient 
    321 # endif 
    322 #endif 
    323       !!---------------------------------------------------------------------- 
    324  
    325       ! 0. Initialization 
     354      !!---------------------------------------------------------------------- 
     355 
     356      ! 0. Memory allocation 
     357      IF( wrk_in_use(3, 1,2,3,4,5,6,7,8,9,10,11,12,13) .OR. & 
     358          wrk_in_use(2, 1,2,3,4,5,6,7,8,9,10)               ) THEN 
     359         CALL ctl_stop('domrea/dta_dyn: requested workspace arrays unavailable')  ;  RETURN 
     360      END IF 
    326361       
    327362      ! cas d'un fichier non periodique : on utilise deux fois le premier et 
     
    390425 
    391426#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv  
    392       CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw (:,: ), jkenr ) 
     427      CALL iom_get( numfl_w, jpdom_data, 'soleaeiw', zaeiw2d(:,: ), jkenr ) 
    393428#endif 
    394429 
     
    413448 
    414449#if ! defined key_degrad && defined key_traldf_c2d && defined key_traldf_eiv 
    415       aeiwdta(:,:,2)  = zaeiw(:,:) * tmask(:,:,1) 
     450      aeiwdta(:,:,2)  = zaeiw2d(:,:) * tmask(:,:,1) 
    416451#endif 
    417452 
     
    451486      ENDIF 
    452487      !       
     488      IF( wrk_not_released(3, 1,2,3,4,5,6,7,8,9,10,11,12,13) .OR. & 
     489          wrk_not_released(2, 1,2,3,4,5,6,7,8,9,10)               ) THEN 
     490         CALL ctl_stop('domrea/dta_dyn: failed to release workspace arrays.') 
     491      END IF 
     492      ! 
    453493   END SUBROUTINE dynrea 
    454494 
     
    462502      !! ** Method : 
    463503      !!---------------------------------------------------------------------- 
    464       REAL(wp) ::   znspyr   !: number of time step per year 
     504      REAL(wp) :: znspyr   !: number of time step per year 
     505      INTEGER  :: ierr 
    465506      !! 
    466507      NAMELIST/namdyn/ ndtadyn, ndtatot, nsptint, lperdyn,  & 
    467508      &                cfile_grid_T, cfile_grid_U, cfile_grid_V, cfile_grid_W 
    468509      !!---------------------------------------------------------------------- 
     510 
     511      ierr = dta_dyn_alloc() 
     512      IF( ierr /= 0 )  CALL ctl_stop( 'STOP', 'dta_dyn_alloc : unable to allocate standard ocean arrays' ) 
    469513 
    470514      !  Define the dynamical input parameters 
  • branches/dev_r2586_dynamic_mem/NEMOGCM/NEMO/OFF_SRC/nemogcm.F90

    r2574 r2648  
    55   !!====================================================================== 
    66   !! History :  3.3  ! 2010-05  (C. Ethe)  Full reorganization of the off-line: phasing with the on-line 
     7   !!            4.0  ! 2011-01  (C. Ethe, A. R. Porter, STFC Daresbury) dynamical allocation 
    78   !!---------------------------------------------------------------------- 
    89 
     
    2627   USE trabbl          ! bottom boundary layer          (tra_bbl_init routine) 
    2728   USE zdfini          ! vertical physics: initialization 
     29   USE sbcmod          ! surface boundary condition       (sbc_init     routine) 
    2830   USE phycst          ! physical constant                  (par_cst routine) 
    2931   USE dtadyn          ! Lecture and Interpolation of the dynamical fields 
     
    137139      !                             !--------------------------------------------! 
    138140#if defined key_iomput 
    139       CALL init_ioclient( ilocal_comm )       ! nemo local communicator (used or not) given by the io_server 
    140       narea = mynode( cltxt, ilocal_comm )    ! Nodes selection 
     141      CALL  init_ioclient( ilocal_comm )                 ! exchange io_server nemo local communicator with the io_server 
     142      narea = mynode( cltxt, numnam, nstop, ilocal_comm )   ! Nodes selection 
    141143#else 
    142       narea = mynode( cltxt )                 ! Nodes selection (control print return in cltxt) 
     144      ilocal_comm = 0 
     145      narea = mynode( cltxt, numnam, nstop )                 ! Nodes selection (control print return in cltxt) 
    143146#endif 
     147 
    144148      narea = narea + 1                       ! mynode return the rank of proc (0 --> jpnij -1 ) 
    145149 
    146150      lwp = (narea == 1) .OR. ln_ctl          ! control of all listing output print 
     151 
     152      ! Decide on size of grid now that we have our communicator size 
     153      ! If we're not using dynamic memory then mpp_partition does nothing. 
     154 
     155#if   defined key_mpp_mpi   ||   defined key_mpp_shmem 
     156      CALL nemo_partition(mppsize) 
     157#else 
     158      jpni = 1 
     159      jpnj = 1 
     160      jpnij = jpni*jpnj 
     161#endif 
     162      ! Calculate domain dimensions given calculated jpni and jpnj 
     163      ! This used to be done in par_oce.F90 when they were parameters rather 
     164      ! than variables 
     165      jpi = ( jpiglo-2*jpreci + (jpni-1) ) / jpni + 2*jpreci !: first  dim. 
     166      jpj = ( jpjglo-2*jprecj + (jpnj-1) ) / jpnj + 2*jprecj !: second dim. 
     167      jpim1 = jpi-1                                          !: inner domain indices 
     168      jpjm1 = jpj-1                                          !:   "           " 
     169      jpkm1 = jpk-1                                          !:   "           " 
     170      jpij  = jpi*jpj                                        !:  jpi x j 
     171 
     172      ! Now we know the dimensions of the grid, allocate arrays 
     173      CALL nemo_alloc() 
    147174 
    148175      IF(lwp) THEN                            ! open listing units 
     
    182209 
    183210      !                                     ! Ocean physics 
     211                            CALL     sbc_init   ! Forcings : surface module 
    184212#if ! defined key_degrad 
    185213                            CALL ldf_tra_init   ! Lateral ocean tracer physics 
     
    307335   END SUBROUTINE nemo_closefile 
    308336 
     337   SUBROUTINE nemo_alloc 
     338     !!---------------------------------------------------------------------- 
     339     !!                     ***  ROUTINE nemo_alloc  *** 
     340     !! 
     341     !! ** Purpose :   Allocate all the dynamic arrays of the OPA modules 
     342     !! 
     343     !! ** Method  : 
     344     !!---------------------------------------------------------------------- 
     345     USE diawri,       ONLY: dia_wri_alloc 
     346     USE dom_oce,      ONLY: dom_oce_alloc 
     347     USE zdf_oce,      ONLY: zdf_oce_alloc 
     348     USE zdfmxl,       ONLY: zdf_mxl_alloc 
     349     USE ldftra_oce,   ONLY: ldftra_oce_alloc 
     350     USE trc_oce,      ONLY: trc_oce_alloc 
     351 
     352      USE wrk_nemo,    ONLY: wrk_alloc 
     353 
     354      INTEGER :: ierr 
     355      !!---------------------------------------------------------------------- 
     356 
     357      ierr =        oce_alloc       ()          ! ocean  
     358      ierr = ierr + dia_wri_alloc   () 
     359      ierr = ierr + dom_oce_alloc   ()          ! ocean domain 
     360      ierr = ierr + ldftra_oce_alloc()          ! ocean lateral  physics : tracers 
     361      ierr = ierr + zdf_oce_alloc   ()          ! ocean vertical physics 
     362      ierr = ierr + zdf_mxl_alloc   ()          ! ocean vertical physics 
     363      ! 
     364      ierr = ierr + lib_mpp_alloc   (numout)    ! mpp exchanges 
     365      ierr = ierr + trc_oce_alloc   ()          ! shared TRC / TRA arrays 
     366      ierr = ierr + wrk_alloc(numout, lwp) 
     367 
     368      IF( lk_mpp    )   CALL mpp_sum( ierr ) 
     369      IF( ierr /= 0 )   CALL ctl_stop( 'STOP', 'nemo_alloc : unable to allocate standard ocean arrays' ) 
     370      ! 
     371   END SUBROUTINE nemo_alloc 
     372 
     373   SUBROUTINE nemo_partition( num_pes ) 
     374      !!---------------------------------------------------------------------- 
     375      !!                 ***  ROUTINE nemo_partition  *** 
     376      !! 
     377      !! ** Purpose :    
     378      !! 
     379      !! ** Method  : 
     380      !!---------------------------------------------------------------------- 
     381      USE par_oce 
     382      INTEGER, INTENT(in) :: num_pes ! The number of MPI processes we have 
     383      ! Local variables 
     384      INTEGER, PARAMETER :: nfactmax = 20 
     385      INTEGER :: nfact ! The no. of factors returned 
     386      INTEGER :: ierr  ! Error flag 
     387      INTEGER :: i 
     388      INTEGER :: idiff, mindiff, imin ! For choosing pair of factors that are closest in value 
     389      INTEGER, DIMENSION(nfactmax) :: ifact ! Array of factors 
     390      !!---------------------------------------------------------------------- 
     391 
     392      ierr = 0 
     393 
     394      CALL factorise( ifact, nfactmax, nfact, num_pes, ierr ) 
     395 
     396      IF( nfact <= 1 ) THEN 
     397         WRITE (numout, *) 'WARNING: factorisation of number of PEs failed' 
     398         WRITE (numout, *) '       : using grid of ',num_pes,' x 1' 
     399         jpnj = 1 
     400         jpni = num_pes 
     401      ELSE 
     402         ! Search through factors for the pair that are closest in value 
     403         mindiff = 1000000 
     404         imin    = 1 
     405         DO i=1,nfact-1,2 
     406            idiff = ABS(ifact(i) - ifact(i+1)) 
     407            IF(idiff < mindiff)THEN 
     408               mindiff = idiff 
     409               imin = i 
     410            END IF 
     411         END DO 
     412         jpnj = ifact(imin) 
     413         jpni = ifact(imin + 1) 
     414      ENDIF 
     415      jpnij = jpni*jpnj 
     416 
     417      WRITE(*,*) 'ARPDBG: jpni = ',jpni,'jpnj = ',jpnj,'jpnij = ',jpnij 
     418      ! 
     419   END SUBROUTINE nemo_partition 
     420 
     421 
     422   SUBROUTINE factorise( ifax, maxfax, nfax, n, ierr ) 
     423      !!---------------------------------------------------------------------- 
     424      !!                     ***  ROUTINE factorise  *** 
     425      !! 
     426      !! ** Purpose :   return the prime factors of n. 
     427      !!                nfax factors are returned in array ifax which is of  
     428      !!                maximum dimension maxfax. 
     429      !! ** Method  : 
     430      !!---------------------------------------------------------------------- 
     431      INTEGER, INTENT(in)  :: n, maxfax 
     432      INTEGER, INTENT(Out) :: ierr, nfax 
     433      INTEGER, INTENT(out) :: ifax(maxfax) 
     434      ! Local variables. 
     435      INTEGER :: i, ifac, l, nu 
     436      INTEGER, PARAMETER :: ntest = 14 
     437      INTEGER :: lfax(ntest) 
     438 
     439      ! lfax contains the set of allowed factors. 
     440      data (lfax(i),i=1,ntest) / 16384, 8192, 4096, 2048, 1024, 512, 256,  & 
     441         &                         128,   64,   32,   16,    8,   4,   2  / 
     442      !!---------------------------------------------------------------------- 
     443 
     444      ! Clear the error flag and initialise output vars 
     445      ierr = 0 
     446      ifax = 1 
     447      nfax = 0 
     448 
     449      ! Find the factors of n. 
     450      IF( n == 1 ) GOTO 20 
     451 
     452      ! nu holds the unfactorised part of the number. 
     453      ! nfax holds the number of factors found. 
     454      ! l points to the allowed factor list. 
     455      ! ifac holds the current factor. 
     456 
     457      nu   = n 
     458      nfax = 0 
     459 
     460      DO l = ntest, 1, -1 
     461         ! 
     462         ifac = lfax(l) 
     463         IF(ifac > nu)CYCLE 
     464 
     465         ! Test whether the factor will divide. 
     466 
     467         IF( MOD(nu,ifac) == 0 ) THEN 
     468            ! 
     469            nfax = nfax+1            ! Add the factor to the list 
     470            IF( nfax > maxfax ) THEN 
     471               ierr = 6 
     472               write (*,*) 'FACTOR: insufficient space in factor array ',nfax 
     473               return 
     474            ENDIF 
     475            ifax(nfax) = ifac 
     476            ! Store the other factor that goes with this one 
     477            nfax = nfax + 1 
     478            ifax(nfax) = nu / ifac 
     479            !WRITE (*,*) 'ARPDBG, factors ',nfax-1,' & ',nfax,' are ', & 
     480            !            ifax(nfax-1),' and ',ifax(nfax) 
     481         ENDIF 
     482         ! 
     483      END DO 
     484 
     485   20 CONTINUE      ! Label 20 is the exit point from the factor search loop. 
     486      ! 
     487      RETURN 
     488      ! 
     489   END SUBROUTINE factorise 
    309490   !!====================================================================== 
    310491END MODULE nemogcm 
Note: See TracChangeset for help on using the changeset viewer.