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 10397 for NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src – NEMO

Ignore:
Timestamp:
2018-12-14T17:27:24+01:00 (5 years ago)
Author:
smasson
Message:

dev_r10164_HPC09_ESIWACE_PREP_MERGE: action 5 introduce mpp_delay_sum, see #2133

Location:
NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src
Files:
11 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/ICE/icedyn_adv_umx.F90

    r10386 r10397  
    7373      INTEGER  ::   ipl                     ! third dimention of tracer array 
    7474 
    75       REAL(wp) ::   zcfl(2), zusnit, zdt      !   -      - 
     75      REAL(wp) ::   zusnit, zdt 
     76      REAL(wp), DIMENSION(1)                  ::   zcflprv, zcflnow   ! send zcflnow and receive zcflprv 
    7677      REAL(wp), ALLOCATABLE, DIMENSION(:,:)   ::   zudy, zvdx, zcu_box, zcv_box 
    7778      REAL(wp), ALLOCATABLE, DIMENSION(:,:,:) ::   zpato 
     
    8687      !     When needed, the advection split is applied at the next time-step in order to avoid blocking global comm. 
    8788      !     ...this should not affect too much the stability... Was ok on the tests we did... 
    88       zcfl(1) =               MAXVAL( ABS( pu_ice(:,:) ) * rdt_ice * r1_e1u(:,:) ) 
    89       zcfl(1) = MAX( zcfl(1), MAXVAL( ABS( pv_ice(:,:) ) * rdt_ice * r1_e2v(:,:) ) ) 
     89      zcflnow(1) =                  MAXVAL( ABS( pu_ice(:,:) ) * rdt_ice * r1_e1u(:,:) ) 
     90      zcflnow(1) = MAX( zcflnow(1), MAXVAL( ABS( pv_ice(:,:) ) * rdt_ice * r1_e2v(:,:) ) ) 
    9091       
    91       ! non-blocking global communication send zcfl(1) and receive zcfl(2) 
    92       IF( lk_mpp )   CALL mpp_delay_max( 'icedyn_adv_umx', 'advumx_delay', zcfl, kt == nitend - nn_fsbc + 1 ) 
    93  
    94       IF( zcfl(2) > .5 ) THEN   ;   initad = 2   ;   zusnit = 0.5_wp   ! zcfl(2) corresponding to zcfl(1) of the previous time-step 
    95       ELSE                      ;   initad = 1   ;   zusnit = 1.0_wp 
     92      ! non-blocking global communication send zcflnow and receive zcflprv 
     93      CALL mpp_delay_max( 'icedyn_adv_umx', 'cflice', zcflnow(:), zcflprv(:), kt == nitend - nn_fsbc + 1 ) 
     94 
     95      IF( zcflprv(1) > .5 ) THEN   ;   initad = 2   ;   zusnit = 0.5_wp 
     96      ELSE                         ;   initad = 1   ;   zusnit = 1.0_wp 
    9697      ENDIF 
    9798 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/ICE/icerst.F90

    r10380 r10397  
    9595      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    9696      !! 
    97       INTEGER ::   ji, jk ,jl   ! dummy loop indices 
     97      INTEGER ::   jk    ! dummy loop indices 
    9898      INTEGER ::   iter 
    9999      CHARACTER(len=25) ::   znam 
     
    115115      CALL iom_rstput( iter, nitrst, numriw, 'nn_fsbc', REAL( nn_fsbc, wp ) )      ! time-step  
    116116      CALL iom_rstput( iter, nitrst, numriw, 'kt_ice' , REAL( iter   , wp ) )      ! date 
    117       DO ji = 1, nbdelay 
    118          CALL iom_rstput( iter, nitrst, numriw, c_delaylist(ji), todelay(ji) ) 
    119       END DO 
     117      CALL iom_delay_rst( 'WRITE', 'ICE', numriw )   ! save only ice delayed global communication variables 
    120118 
    121119      ! Prognostic variables 
     
    169167      !! ** purpose  :   read restart file 
    170168      !!---------------------------------------------------------------------- 
    171       INTEGER           ::   ji, jk, jl 
     169      INTEGER           ::   jk 
    172170      LOGICAL           ::   llok 
    173171      INTEGER           ::   id1            ! local integer 
    174172      CHARACTER(len=25) ::   znam 
    175173      CHARACTER(len=2)  ::   zchar, zchar1 
    176       REAL(wp)          ::   zfice, ziter, zcfl 
     174      REAL(wp)          ::   zfice, ziter 
    177175      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   z3d   ! 3D workspace 
    178176      !!---------------------------------------------------------------------- 
     
    236234      CALL iom_get( numrir, jpdom_autoglo, 'v_ice', v_ice ) 
    237235 
    238       ! fields needed for ice_dyn_adv_umx 
    239       zcfl =            MAXVAL( ABS( u_ice(:,:) ) * rdt_ice * r1_e1u(:,:) ) 
    240       zcfl = MAX( zcfl, MAXVAL( ABS( v_ice(:,:) ) * rdt_ice * r1_e2v(:,:) ) ) 
    241  
    242       DO ji = 1, nbdelay 
    243          IF( iom_varid( numrir, c_delaylist(ji), ldstop = .FALSE. ) > 0 )  THEN 
    244             CALL iom_get( numrir, c_delaylist(ji), todelay(ji) ) 
    245             ndelayid(ji) = 0   ! set to 0 to speficy that the value was read in the restart 
    246          ENDIF 
    247       END DO 
     236      CALL iom_delay_rst( 'READ', 'ICE', numrir )   ! read only ice delayed global communication variables 
    248237 
    249238      ! fields needed for Met Office (Jules) coupling 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/IOM/iom.F90

    r10380 r10397  
    5858#endif 
    5959   PUBLIC iom_init, iom_swap, iom_open, iom_close, iom_setkt, iom_varid, iom_get 
    60    PUBLIC iom_getatt, iom_putatt, iom_gettime, iom_rstput, iom_put 
     60   PUBLIC iom_chkatt, iom_getatt, iom_putatt, iom_gettime, iom_rstput, iom_delay_rst, iom_put 
    6161   PUBLIC iom_use, iom_context_finalize 
    6262 
     
    14011401   END SUBROUTINE iom_gettime 
    14021402 
     1403   !!---------------------------------------------------------------------- 
     1404   !!                   INTERFACE iom_chkatt 
     1405   !!---------------------------------------------------------------------- 
     1406   SUBROUTINE iom_chkatt( kiomid, cdatt, llok, ksize, cdvar ) 
     1407      INTEGER         , INTENT(in   )                 ::   kiomid    ! Identifier of the file 
     1408      CHARACTER(len=*), INTENT(in   )                 ::   cdatt     ! Name of the attribute 
     1409      LOGICAL         , INTENT(  out)                 ::   llok      ! Error code 
     1410      INTEGER         , INTENT(  out), OPTIONAL       ::   ksize     ! Size of the attribute array 
     1411      CHARACTER(len=*), INTENT(in   ), OPTIONAL       ::   cdvar     ! Name of the variable 
     1412      ! 
     1413      IF( kiomid > 0 ) THEN 
     1414         IF( iom_file(kiomid)%nfid > 0 )   CALL iom_nf90_chkatt( kiomid, cdatt, llok, ksize=ksize, cdvar=cdvar ) 
     1415      ENDIF 
     1416      ! 
     1417   END SUBROUTINE iom_chkatt 
    14031418 
    14041419   !!---------------------------------------------------------------------- 
     
    16431658   END SUBROUTINE iom_rp3d 
    16441659 
     1660 
     1661  SUBROUTINE iom_delay_rst( cdaction, cdcpnt, kncid ) 
     1662      !!--------------------------------------------------------------------- 
     1663      !!   Routine iom_delay_rst: used read/write restart related to mpp_delay 
     1664      !! 
     1665      !!--------------------------------------------------------------------- 
     1666      CHARACTER(len=*), INTENT(in   ) ::   cdaction        ! 
     1667      CHARACTER(len=*), INTENT(in   ) ::   cdcpnt 
     1668      INTEGER         , INTENT(in   ) ::   kncid 
     1669      ! 
     1670      INTEGER  :: ji 
     1671      INTEGER  :: indim 
     1672      LOGICAL  :: llattexist 
     1673      REAL(wp), ALLOCATABLE, DIMENSION(:) ::   zreal1d 
     1674      !!--------------------------------------------------------------------- 
     1675      ! 
     1676      !                                      =================================== 
     1677      IF( TRIM(cdaction) == 'READ' ) THEN   ! read restart related to mpp_delay ! 
     1678         !                                   =================================== 
     1679         DO ji = 1, nbdelay 
     1680            IF ( c_delaycpnt(ji) == cdcpnt ) THEN 
     1681               CALL iom_chkatt( kncid, 'DELAY_'//c_delaylist(ji), llattexist, indim ) 
     1682               IF( llattexist )  THEN 
     1683                  ALLOCATE( todelay(ji)%z1d(indim) ) 
     1684                  CALL iom_getatt( kncid, 'DELAY_'//c_delaylist(ji), todelay(ji)%z1d(:) ) 
     1685                  ndelayid(ji) = 0   ! set to 0 to specify that the value was read in the restart 
     1686               ENDIF 
     1687           ENDIF 
     1688         END DO 
     1689         !                                   ==================================== 
     1690      ELSE                                  ! write restart related to mpp_delay ! 
     1691         !                                   ==================================== 
     1692         DO ji = 1, nbdelay   ! save only ocean delayed global communication variables 
     1693            IF ( c_delaycpnt(ji) == cdcpnt ) THEN 
     1694               IF( ASSOCIATED(todelay(ji)%z1d) ) THEN 
     1695                  CALL mpp_delay_rcv(ji)   ! make sure %z1d is received 
     1696                  CALL iom_putatt( kncid, 'DELAY_'//c_delaylist(ji), todelay(ji)%z1d(:) ) 
     1697               ENDIF 
     1698            ENDIF 
     1699         END DO 
     1700         ! 
     1701      ENDIF 
     1702       
     1703   END SUBROUTINE iom_delay_rst 
     1704   
     1705    
    16451706 
    16461707   !!---------------------------------------------------------------------- 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/IOM/iom_nf90.F90

    r10380 r10397  
    3030 
    3131   PUBLIC iom_nf90_open  , iom_nf90_close, iom_nf90_varid, iom_nf90_get, iom_nf90_gettime, iom_nf90_rstput 
    32    PUBLIC iom_nf90_getatt, iom_nf90_putatt 
     32   PUBLIC iom_nf90_chkatt, iom_nf90_getatt, iom_nf90_putatt 
    3333 
    3434   INTERFACE iom_nf90_get 
     
    321321 
    322322 
     323   SUBROUTINE iom_nf90_chkatt( kiomid, cdatt, llok, ksize, cdvar ) 
     324      !!----------------------------------------------------------------------- 
     325      !!                  ***  ROUTINE  iom_nf90_chkatt  *** 
     326      !! 
     327      !! ** Purpose : check existence of attribute with NF90 
     328      !!              (either a global attribute (default) or a variable 
     329      !!               attribute if optional variable name is supplied (cdvar)) 
     330      !!----------------------------------------------------------------------- 
     331      INTEGER         , INTENT(in   ) ::   kiomid   ! Identifier of the file 
     332      CHARACTER(len=*), INTENT(in   ) ::   cdatt    ! attribute name 
     333      LOGICAL         , INTENT(  out) ::   llok     ! error code 
     334      INTEGER         , INTENT(  out), OPTIONAL     & 
     335                      &               ::   ksize    ! attribute size 
     336      CHARACTER(len=*), INTENT(in   ), OPTIONAL     & 
     337                      &               ::   cdvar    ! name of the variable 
     338      ! 
     339      INTEGER                         ::   if90id   ! temporary integer 
     340      INTEGER                         ::   isize    ! temporary integer 
     341      INTEGER                         ::   ivarid   ! NetCDF variable Id 
     342      !--------------------------------------------------------------------- 
     343      ! 
     344      if90id = iom_file(kiomid)%nfid 
     345      IF( PRESENT(cdvar) ) THEN 
     346         ! check the variable exists in the file 
     347         llok = NF90_INQ_VARID( if90id, TRIM(cdvar), ivarid ) == nf90_noerr 
     348         IF( llok ) & 
     349            ! check the variable has the attribute required 
     350            llok = NF90_Inquire_attribute(if90id, ivarid, cdatt, len=isize ) == nf90_noerr 
     351      ELSE 
     352         llok = NF90_Inquire_attribute(if90id, NF90_GLOBAL, cdatt, len=isize ) == nf90_noerr 
     353      ENDIF 
     354      ! 
     355      IF( PRESENT(ksize) ) ksize = isize 
     356      ! 
     357      IF( .not. llok) & 
     358         CALL ctl_warn('iom_nf90_chkatt: no attribute '//cdatt//' found') 
     359      ! 
     360   END SUBROUTINE iom_nf90_chkatt 
     361 
     362 
     363   !!---------------------------------------------------------------------- 
     364   !!                   INTERFACE iom_nf90_getatt 
     365   !!---------------------------------------------------------------------- 
     366 
    323367   SUBROUTINE iom_nf90_getatt( kiomid, cdatt, katt0d, katt1d, patt0d, patt1d, cdatt0d, cdvar) 
    324368      !!----------------------------------------------------------------------- 
     
    397441      INTEGER                         ::   if90id   ! temporary integer 
    398442      INTEGER                         ::   ivarid   ! NetCDF variable Id 
     443      INTEGER                         ::   isize    ! Attribute size 
     444      INTEGER                         ::   itype    ! Attribute type 
    399445      LOGICAL                         ::   llok     ! temporary logical 
    400       LOGICAL                         ::   lenddef  ! temporary logical 
     446      LOGICAL                         ::   llatt     ! temporary logical 
     447      LOGICAL                         ::   lldata   ! temporary logical 
    401448      CHARACTER(LEN=100)              ::   clinfo   ! info character 
    402449      !--------------------------------------------------------------------- 
    403450      ! 
    404451      if90id = iom_file(kiomid)%nfid 
    405       lenddef = .false. 
    406452      IF( PRESENT(cdvar) ) THEN 
    407          ! check the variable exists in the file 
    408          llok = NF90_INQ_VARID( if90id, TRIM(cdvar), ivarid ) == nf90_noerr 
     453         llok = NF90_INQ_VARID( if90id, TRIM(cdvar), ivarid ) == nf90_noerr   ! is the variable in the file? 
    409454         IF( .NOT. llok ) THEN 
    410             CALL ctl_warn('iom_nf90_putatt: no variable '//TRIM(cdvar)//' found') 
     455            CALL ctl_warn('iom_nf90_putatt: no variable '//TRIM(cdvar)//' found'   & 
     456               &        , '                 no attribute '//cdatt//' written' ) 
     457            RETURN 
    411458         ENDIF 
    412459      ELSE 
    413          llok = .true. 
    414460         ivarid = NF90_GLOBAL 
    415461      ENDIF 
    416       ! 
    417       IF( llok) THEN 
    418          clinfo = 'iom_nf90_putatt, file: '//TRIM(iom_file(kiomid)%name)//', att: '//TRIM(cdatt) 
    419          IF( iom_file(kiomid)%irec /= -1 ) THEN    
    420             ! trick: irec used to know if the file is in define mode or not 
    421             ! if it is not then temporarily put it into define mode 
    422             CALL iom_nf90_check(NF90_REDEF( if90id ), clinfo) 
    423             lenddef = .true. 
    424          ENDIF 
     462      llatt = NF90_Inquire_attribute(if90id, ivarid, cdatt, len = isize, xtype = itype ) == nf90_noerr 
     463      ! 
     464      ! trick: irec used to know if the file is in define mode or not 
     465      lldata = iom_file(kiomid)%irec /= -1   ! default: go back in define mode if in data mode 
     466      IF( lldata .AND. llatt ) THEN          ! attribute already there. Do we really need to go back in define mode? 
     467         ! do we have the appropriate type? 
     468         IF(PRESENT( katt0d) .OR. PRESENT( katt1d))   llok = itype == NF90_INT 
     469         IF(PRESENT( patt0d) .OR. PRESENT( patt1d))   llok = itype == NF90_DOUBLE 
     470         IF(PRESENT(cdatt0d)                      )   llok = itype == NF90_CHAR 
     471         ! and do we have the appropriate size? 
     472         IF(PRESENT( katt0d))   llok = llok .AND. isize == 1 
     473         IF(PRESENT( katt1d))   llok = llok .AND. isize == SIZE(katt1d) 
     474         IF(PRESENT( patt0d))   llok = llok .AND. isize == 1 
     475         IF(PRESENT( patt1d))   llok = llok .AND. isize == SIZE(patt1d) 
     476         IF(PRESENT(cdatt0d))   llok = llok .AND. isize == LEN_TRIM(cdatt0d) 
    425477         ! 
    426          IF(PRESENT( katt0d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =  katt0d), clinfo) 
    427          IF(PRESENT( katt1d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =  katt1d), clinfo) 
    428          IF(PRESENT( patt0d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =  patt0d), clinfo) 
    429          IF(PRESENT( patt1d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =  patt1d), clinfo) 
    430          IF(PRESENT(cdatt0d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values = cdatt0d), clinfo) 
    431          ! 
    432          IF( lenddef ) THEN   ! file was in data mode on entry; put it back in that mode 
    433             CALL iom_nf90_check(NF90_ENDDEF( if90id ), clinfo) 
    434          ENDIF 
    435       ELSE 
    436          CALL ctl_warn('iom_nf90_putatt: no attribute '//TRIM(cdatt)//' written') 
    437       ENDIF 
     478         lldata = .NOT. llok 
     479      ENDIF 
     480      ! 
     481      clinfo = 'iom_nf90_putatt, file: '//TRIM(iom_file(kiomid)%name)//', att: '//TRIM(cdatt) 
     482      IF(lldata)   CALL iom_nf90_check(NF90_REDEF( if90id ), clinfo)   ! leave data mode to define mode 
     483      ! 
     484      IF(PRESENT( katt0d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =       katt0d) , clinfo) 
     485      IF(PRESENT( katt1d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =       katt1d) , clinfo) 
     486      IF(PRESENT( patt0d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =       patt0d) , clinfo) 
     487      IF(PRESENT( patt1d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values =       patt1d) , clinfo) 
     488      IF(PRESENT(cdatt0d))   CALL iom_nf90_check(NF90_PUT_ATT(if90id, ivarid, cdatt, values = trim(cdatt0d)), clinfo) 
     489      ! 
     490      IF(lldata)   CALL iom_nf90_check(NF90_ENDDEF( if90id ), clinfo)   ! leave define mode to data mode 
    438491      ! 
    439492   END SUBROUTINE iom_nf90_putatt 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/IOM/restart.F90

    r10380 r10397  
    2828   USE iom             ! I/O module 
    2929   USE diurnal_bulk 
     30   USE lib_mpp         ! distribued memory computing library 
    3031 
    3132   IMPLICIT NONE 
     
    140141      !!---------------------------------------------------------------------- 
    141142      INTEGER, INTENT(in) ::   kt   ! ocean time-step 
    142       INTEGER             :: ji 
    143143      !!---------------------------------------------------------------------- 
    144144                     IF(lwxios) CALL iom_swap(      cwxios_context          ) 
    145145                     CALL iom_rstput( kt, nitrst, numrow, 'rdt'    , rdt       , ldxios = lwxios)   ! dynamics time step 
    146                      DO ji = 1, nbdelay 
    147                         CALL iom_rstput( kt, nitrst, numrow, c_delaylist(ji), todelay(ji), ldxios = lwxios) 
    148                      END DO 
     146                     CALL iom_delay_rst( 'WRITE', 'OCE', numrow )   ! save only ocean delayed global communication variables 
    149147 
    150148      IF ( .NOT. ln_diurnal_only ) THEN 
     
    249247      !!---------------------------------------------------------------------- 
    250248      REAL(wp) ::   zrdt 
    251       INTEGER  ::   ji, jk 
     249      INTEGER  ::   jk 
    252250      REAL(wp), DIMENSION(jpi, jpj, jpk) :: w3d 
    253251      !!---------------------------------------------------------------------- 
     
    261259      ENDIF 
    262260 
    263       DO ji = 1, nbdelay 
    264          IF( iom_varid( numror, c_delaylist(ji), ldstop = .FALSE. ) > 0 )  THEN 
    265             CALL iom_get( numror, c_delaylist(ji), todelay(ji), ldxios = lrxios ) 
    266             ndelayid(ji) = 0   ! set to 0 to speficy that the value was read in the restart 
    267          ENDIF 
    268       END DO 
     261      CALL iom_delay_rst( 'READ', 'OCE', numror )   ! read only ocean delayed global communication variables 
    269262       
    270263      ! Diurnal DSST  
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/LBC/lib_mpp.F90

    r10386 r10397  
    8383   PUBLIC   mpp_lbc_north_icb 
    8484   PUBLIC   mpp_min, mpp_max, mpp_sum, mpp_minloc, mpp_maxloc 
    85    PUBLIC   mpp_delay_max 
     85   PUBLIC   mpp_delay_max, mpp_delay_sum, mpp_delay_rcv 
    8686   PUBLIC   mppscatter, mppgather 
    8787   PUBLIC   mpp_ini_znl 
     
    164164   INTEGER, PUBLIC                               ::   numcom = -1                  !: logical unit for communicaton report 
    165165   LOGICAL, PUBLIC                               ::   l_full_nf_update = .TRUE.    !: logical for a full (2lines) update of bc at North fold report 
    166    INTEGER,                    PARAMETER, PUBLIC ::   nbdelay = 1       !: number of delayed operations 
    167 !$AGRIF_DO_NOT_TREAT 
    168    CHARACTER(len=32), DIMENSION(nbdelay), PUBLIC ::   c_delaylist       !: name (used as id) of allreduce-delayed operations 
    169    DATA c_delaylist(1) / 'advumx_delay' / 
    170 !$AGRIF_END_DO_NOT_TREAT 
    171    REAL(wp),          DIMENSION(nbdelay), PUBLIC ::   todelay           !: current value of the delayed operations 
    172    INTEGER,           DIMENSION(nbdelay), PUBLIC ::   ndelayid = -1     !: mpi request id of the delayed operations 
     166   INTEGER,                    PARAMETER, PUBLIC ::   nbdelay = 2       !: number of delayed operations 
     167   !: name (used as id) of allreduce-delayed operations 
     168   CHARACTER(len=32), DIMENSION(nbdelay), PUBLIC ::   c_delaylist = (/ 'cflice', 'fwb' /) 
     169   !: component name where the allreduce-delayed operation is performed 
     170   CHARACTER(len=3),  DIMENSION(nbdelay), PUBLIC ::   c_delaycpnt = (/ 'ICE'   , 'OCE' /) 
     171   TYPE, PUBLIC ::   DELAYARR 
     172      REAL(   wp), POINTER, DIMENSION(:) ::  z1d => NULL() 
     173      COMPLEX(wp), POINTER, DIMENSION(:) ::  y1d => NULL() 
     174   END TYPE DELAYARR 
     175   TYPE( DELAYARR ), DIMENSION(nbdelay), PUBLIC  ::   todelay               
     176   INTEGER,          DIMENSION(nbdelay), PUBLIC  ::   ndelayid = -1     !: mpi request id of the delayed operations 
    173177 
    174178   ! timing summary report 
     
    573577   END SUBROUTINE mppscatter 
    574578 
    575    
    576    SUBROUTINE mpp_delay_max( cdname, cdelay, pinout, ldlast, kcom ) 
    577       !!---------------------------------------------------------------------- 
    578       !!                   ***  routine mpp_delay_max  *** 
    579       !! 
    580       !! ** Purpose :   performed delayed mpp_max, the result is received on next call 
     579    
     580   SUBROUTINE mpp_delay_sum( cdname, cdelay, y_in, pout, ldlast, kcom ) 
     581     !!---------------------------------------------------------------------- 
     582      !!                   ***  routine mpp_delay_sum  *** 
     583      !! 
     584      !! ** Purpose :   performed delayed mpp_sum, the result is received on next call 
    581585      !! 
    582586      !!---------------------------------------------------------------------- 
    583587      CHARACTER(len=*), INTENT(in   )               ::   cdname  ! name of the calling subroutine 
    584588      CHARACTER(len=*), INTENT(in   )               ::   cdelay  ! name (used as id) of the delayed operation 
    585       REAL(wp),         INTENT(inout), DIMENSION(2) ::   pinout  ! pinout(1): in data, pinout(2): out data 
     589      COMPLEX(wp),      INTENT(in   ), DIMENSION(:) ::   y_in 
     590      REAL(wp),         INTENT(  out), DIMENSION(:) ::   pout 
    586591      LOGICAL,          INTENT(in   )               ::   ldlast  ! true if this is the last time we call this routine 
    587592      INTEGER,          INTENT(in   ), OPTIONAL     ::   kcom 
    588593      !! 
    589       INTEGER ::   ji 
     594      INTEGER ::   ji, isz 
    590595      INTEGER ::   idvar 
    591       INTEGER ::   ierror, ilocalcomm 
     596      INTEGER ::   ierr, ilocalcomm 
     597      COMPLEX(wp), ALLOCATABLE, DIMENSION(:) ::   ytmp 
    592598      !!---------------------------------------------------------------------- 
    593599      ilocalcomm = mpi_comm_oce 
    594600      IF( PRESENT(kcom) )   ilocalcomm = kcom 
    595601 
     602      isz = SIZE(y_in) 
     603       
    596604      IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, ld_glb = .TRUE. ) 
    597605 
    598606      idvar = -1 
    599607      DO ji = 1, nbdelay 
    600          IF( TRIM(cdelay) == trim(c_delaylist(ji)) ) idvar = ji 
     608         IF( TRIM(cdelay) == TRIM(c_delaylist(ji)) ) idvar = ji 
     609      END DO 
     610      IF ( idvar == -1 )   CALL ctl_stop( 'STOP',' mpp_delay_sum : please add a new delayed exchange for '//TRIM(cdname) ) 
     611 
     612      IF ( ndelayid(idvar) == 0 ) THEN         ! first call    with restart: %z1d defined in iom_delay_rst 
     613         !                                       -------------------------- 
     614         IF ( SIZE(todelay(idvar)%z1d) /= isz ) THEN                  ! Check dimension coherence 
     615            IF(lwp) WRITE(numout,*) ' WARNING: the nb of delayed variables in restart file is not the model one' 
     616            DEALLOCATE(todelay(idvar)%z1d) 
     617            ndelayid(idvar) = -1                                      ! do as if we had no restart 
     618         ELSE 
     619            ALLOCATE(todelay(idvar)%y1d(isz)) 
     620            todelay(idvar)%y1d(:) = CMPLX(todelay(idvar)%z1d(:), 0., wp)   ! create %y1d, complex variable needed by mpi_sumdd 
     621         END IF 
     622      ENDIF 
     623       
     624      IF( ndelayid(idvar) == -1 ) THEN         ! first call without restart: define %y1d and %z1d from y_in with blocking allreduce 
     625         !                                       -------------------------- 
     626         ALLOCATE(todelay(idvar)%z1d(isz), todelay(idvar)%y1d(isz))   ! allocate also %z1d as used for the restart 
     627         CALL mpi_allreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ierr )   ! get %y1d 
     628         todelay(idvar)%z1d(:) = REAL(todelay(idvar)%y1d(:), wp)      ! define %z1d from %y1d 
     629      ENDIF 
     630 
     631      IF( ndelayid(idvar) > 0 )   CALL mpp_delay_rcv( idvar )         ! make sure %z1d is received 
     632 
     633      ! send back pout from todelay(idvar)%z1d defined at previous call 
     634      pout(:) = todelay(idvar)%z1d(:) 
     635 
     636      ! send y_in into todelay(idvar)%y1d with a non-blocking communication 
     637      CALL mpi_iallreduce( y_in(:), todelay(idvar)%y1d(:), isz, MPI_DOUBLE_COMPLEX, mpi_sumdd, ilocalcomm, ndelayid(idvar), ierr ) 
     638 
     639   END SUBROUTINE mpp_delay_sum 
     640 
     641    
     642   SUBROUTINE mpp_delay_max( cdname, cdelay, p_in, pout, ldlast, kcom ) 
     643      !!---------------------------------------------------------------------- 
     644      !!                   ***  routine mpp_delay_max  *** 
     645      !! 
     646      !! ** Purpose :   performed delayed mpp_max, the result is received on next call 
     647      !! 
     648      !!---------------------------------------------------------------------- 
     649      CHARACTER(len=*), INTENT(in   )                 ::   cdname  ! name of the calling subroutine 
     650      CHARACTER(len=*), INTENT(in   )                 ::   cdelay  ! name (used as id) of the delayed operation 
     651      REAL(wp),         INTENT(in   ), DIMENSION(:)   ::   p_in    !  
     652      REAL(wp),         INTENT(  out), DIMENSION(:)   ::   pout    !  
     653      LOGICAL,          INTENT(in   )                 ::   ldlast  ! true if this is the last time we call this routine 
     654      INTEGER,          INTENT(in   ), OPTIONAL       ::   kcom 
     655      !! 
     656      INTEGER ::   ji, isz 
     657      INTEGER ::   idvar 
     658      INTEGER ::   ierr, ilocalcomm 
     659      !!---------------------------------------------------------------------- 
     660      ilocalcomm = mpi_comm_oce 
     661      IF( PRESENT(kcom) )   ilocalcomm = kcom 
     662 
     663      isz = SIZE(p_in) 
     664 
     665      IF( narea == 1 .AND. numcom == -1 ) CALL mpp_report( cdname, ld_glb = .TRUE. ) 
     666 
     667      idvar = -1 
     668      DO ji = 1, nbdelay 
     669         IF( TRIM(cdelay) == TRIM(c_delaylist(ji)) ) idvar = ji 
    601670      END DO 
    602671      IF ( idvar == -1 )   CALL ctl_stop( 'STOP',' mpp_delay_max : please add a new delayed exchange for '//TRIM(cdname) ) 
    603672 
    604       IF( ndelayid(idvar) == -1 ) THEN        ! first call without restart: get pinout(2) from pinout(1) with a blocking allreduce 
    605          CALL mpi_allreduce(      pinout(1), pinout(2), 1, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ierror ) 
    606       ELSE IF ( ndelayid(idvar) == 0 ) THEN   ! first call    with restart: get pinout(2) from todelay with a blocking allreduce 
    607          CALL mpi_allreduce( todelay(idvar), pinout(2), 1, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ierror ) 
    608       ELSE                                    ! from the second call, get value from previous time step 
     673      IF ( ndelayid(idvar) == 0 ) THEN         ! first call    with restart: %z1d defined in iom_delay_rst 
     674         !                                       -------------------------- 
     675         IF ( SIZE(todelay(idvar)%z1d) /= isz ) THEN                  ! Check dimension coherence 
     676            IF(lwp) WRITE(numout,*) ' WARNING: the nb of delayed variables in restart file is not the model one' 
     677            DEALLOCATE(todelay(idvar)%z1d) 
     678            ndelayid(idvar) = -1                                      ! do as if we had no restart 
     679         END IF 
     680      ENDIF 
     681 
     682      IF( ndelayid(idvar) == -1 ) THEN         ! first call without restart: define %z1d from p_in with a blocking allreduce 
     683         !                                       -------------------------- 
     684         ALLOCATE(todelay(idvar)%z1d(isz)) 
     685         CALL mpi_allreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ierr )   ! get %z1d 
     686      ENDIF 
     687 
     688      IF( ndelayid(idvar) > 0 )   CALL mpp_delay_rcv( idvar )         ! make sure %z1d is received 
     689 
     690      ! send back pout from todelay(idvar)%z1d defined at previous call 
     691      pout(:) = todelay(idvar)%z1d(:) 
     692 
     693      ! send p_in into todelay(idvar)%z1d with a non-blocking communication 
     694      CALL mpi_iallreduce( p_in(:), todelay(idvar)%z1d(:), isz, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ndelayid(idvar), ierr ) 
     695 
     696   END SUBROUTINE mpp_delay_max 
     697 
     698    
     699   SUBROUTINE mpp_delay_rcv( kid ) 
     700      !!---------------------------------------------------------------------- 
     701      !!                   ***  routine mpp_delay_rcv  *** 
     702      !! 
     703      !! ** Purpose :  force barrier for delayed mpp (needed for restart)  
     704      !! 
     705      !!---------------------------------------------------------------------- 
     706      INTEGER,INTENT(in   )      ::  kid  
     707      INTEGER ::   ierr 
     708      !!---------------------------------------------------------------------- 
     709      IF( ndelayid(kid) /= -2 ) THEN   
    609710         IF( ln_timing ) CALL tic_tac( .TRUE., ld_global = .TRUE.) 
    610          CALL mpi_wait(ndelayid(idvar), MPI_STATUS_IGNORE, ierror )   ! make sure todelay(idvar) is received 
     711         CALL mpi_wait( ndelayid(kid), MPI_STATUS_IGNORE, ierr )                        ! make sure todelay(kid) is received 
    611712         IF( ln_timing ) CALL tic_tac(.FALSE., ld_global = .TRUE.) 
    612          pinout(2) = todelay(idvar)           ! send back pinout(2) from todelay(idvar) defined at previous call 
    613       ENDIF 
    614  
    615       IF( ldlast ) THEN      ! last call: put pinout(1) in todelay that will be stored in the restart files 
    616          todelay(idvar) = pinout(1) 
    617       ELSE                   ! send pinout(1) to todelay(idvar), to be received at next call 
    618          CALL mpi_iallreduce( pinout(1), todelay(idvar), 1, MPI_DOUBLE_PRECISION, mpi_max, ilocalcomm, ndelayid(idvar), ierror ) 
    619       ENDIF 
    620  
    621    END SUBROUTINE mpp_delay_max 
     713         IF( ASSOCIATED(todelay(kid)%y1d) )   todelay(kid)%z1d(:) = REAL(todelay(kid)%y1d(:), wp)  ! define %z1d from %y1d 
     714         ndelayid(kid) = -2   ! add flag to know that mpi_wait was already called on kid 
     715      ENDIF 
     716   END SUBROUTINE mpp_delay_rcv 
     717 
    622718    
    623719   !!---------------------------------------------------------------------- 
     
    13761472      ! 
    13771473      ! find the smallest common frequency: default = frequency product, if multiple, choose the larger of the 2 frequency 
    1378       ncom_freq = ncom_fsbc * ncom_dttrc 
    1379       IF( MOD( MAX(ncom_fsbc,ncom_dttrc), MIN(ncom_fsbc,ncom_dttrc) ) == 0 ) ncom_freq = MAX(ncom_fsbc,ncom_dttrc) 
     1474      IF( ncom_dttrc /= 1 )   CALL ctl_stop( 'STOP', 'mpp_report, ncom_dttrc /= 1 not coded...' )  
     1475      ncom_freq = ncom_fsbc 
    13801476      ! 
    13811477      IF ( ncom_stp == nit000+ncom_freq ) THEN   ! avoid to count extra communications in potential initializations at nit000 
     
    13841480            IF( .NOT. ALLOCATED(    crname_lbc) ) ALLOCATE(     crname_lbc(ncom_rec_max  ) ) 
    13851481            n_sequence_lbc = n_sequence_lbc + 1 
    1386             IF( n_sequence_lbc > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lnk_generic, increase ncom_rec_max' )   ! deadlock 
     1482            IF( n_sequence_lbc > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' )   ! deadlock 
    13871483            crname_lbc(n_sequence_lbc) = cdname     ! keep the name of the calling routine 
    13881484            ncomm_sequence(n_sequence_lbc,1) = kpk*kpl   ! size of 3rd and 4th dimensions 
     
    13921488            IF( .NOT. ALLOCATED(crname_glb) ) ALLOCATE( crname_glb(ncom_rec_max) ) 
    13931489            n_sequence_glb = n_sequence_glb + 1 
    1394             IF( n_sequence_glb > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lnk_generic, increase ncom_rec_max' )   ! deadlock 
     1490            IF( n_sequence_glb > ncom_rec_max ) CALL ctl_stop( 'STOP', 'lib_mpp, increase ncom_rec_max' )   ! deadlock 
    13951491            crname_glb(n_sequence_glb) = cdname     ! keep the name of the calling routine 
    13961492         ENDIF 
     
    15041600   INTEGER, PARAMETER, PUBLIC               ::   nbdelay = 0   ! make sure we don't enter loops: DO ji = 1, nbdelay 
    15051601   CHARACTER(len=32), DIMENSION(1), PUBLIC  ::   c_delaylist = 'empty' 
    1506    REAL(wp), PUBLIC, DIMENSION(1)           ::   todelay 
     1602   CHARACTER(len=32), DIMENSION(1), PUBLIC  ::   c_delaycpnt = 'empty' 
     1603   TYPE ::   DELAYARR 
     1604      REAL(   wp), POINTER, DIMENSION(:) ::  z1d => NULL() 
     1605      COMPLEX(wp), POINTER, DIMENSION(:) ::  y1d => NULL() 
     1606   END TYPE DELAYARR 
     1607   TYPE( DELAYARR ), DIMENSION(1), PUBLIC  ::   todelay               
    15071608   INTEGER,  PUBLIC, DIMENSION(1)           ::   ndelayid = -1 
    15081609   !!---------------------------------------------------------------------- 
     
    16731774#  undef OPERATION_MAXLOC 
    16741775 
    1675    SUBROUTINE mpp_delay_max( cdname, cdelay, pinout, ldlast, kcom ) 
     1776   SUBROUTINE mpp_delay_sum( cdname, cdelay, y_in, pout, ldlast, kcom ) 
    16761777      CHARACTER(len=*), INTENT(in   )               ::   cdname  ! name of the calling subroutine 
    16771778      CHARACTER(len=*), INTENT(in   )               ::   cdelay  ! name (used as id) of the delayed operation 
    1678       REAL(wp),         INTENT(inout), DIMENSION(2) ::   pinout  ! pinout(1): in data, pinout(2): out data 
     1779      COMPLEX(wp),      INTENT(in   ), DIMENSION(:) ::   y_in 
     1780      REAL(wp),         INTENT(  out), DIMENSION(:) ::   pout 
    16791781      LOGICAL,          INTENT(in   )               ::   ldlast  ! true if this is the last time we call this routine 
    16801782      INTEGER,          INTENT(in   ), OPTIONAL     ::   kcom 
    1681       WRITE(*,*) 'mpp_delay_max: You should not have seen this print! error?', cdname 
     1783      ! 
     1784      pout(:) = REAL(y_in(:), wp) 
     1785   END SUBROUTINE mpp_delay_sum 
     1786 
     1787   SUBROUTINE mpp_delay_max( cdname, cdelay, p_in, pout, ldlast, kcom ) 
     1788      CHARACTER(len=*), INTENT(in   )               ::   cdname  ! name of the calling subroutine 
     1789      CHARACTER(len=*), INTENT(in   )               ::   cdelay  ! name (used as id) of the delayed operation 
     1790      REAL(wp),         INTENT(in   ), DIMENSION(:) ::   p_in 
     1791      REAL(wp),         INTENT(  out), DIMENSION(:) ::   pout 
     1792      LOGICAL,          INTENT(in   )               ::   ldlast  ! true if this is the last time we call this routine 
     1793      INTEGER,          INTENT(in   ), OPTIONAL     ::   kcom 
     1794      ! 
     1795      pout(:) = p_in(:) 
    16821796   END SUBROUTINE mpp_delay_max 
    16831797 
     1798   SUBROUTINE mpp_delay_rcv( kid ) 
     1799      INTEGER,INTENT(in   )      ::  kid  
     1800      WRITE(*,*) 'mpp_delay_rcv: You should not have seen this print! error?', kid 
     1801   END SUBROUTINE mpp_delay_rcv 
     1802    
    16841803   SUBROUTINE mppstop( ldfinal, ld_force_abort ) 
    16851804      LOGICAL, OPTIONAL, INTENT(in) :: ldfinal    ! source process number 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/LBC/mppini.F90

    r10357 r10397  
    182182      ! then we calculate them here now that we have our communicator size 
    183183      IF( jpni < 1 .OR. jpnj < 1 )   CALL mpp_init_bestpartition( mppsize, jpni, jpnj ) 
     184!!$      IF( jpni < 1 .OR. jpnj < 1 ) THEN 
     185!!$         CALL mpp_init_bestpartition( mppsize, jpni, jpnj ) 
     186!!$      ELSE 
     187!!$         CALL mpp_init_bestpartition( mppsize, inbi, inbj ) 
     188!!$         CALL mpp_basic_decomposition( jpni, jpnj, jpimax, jpjmax ) 
     189!!$         CALL mpp_basic_decomposition( inbi, inbj,  iimax,  ijmax ) 
     190!!$         IF( iimax*ijmax < jpimax*jpjmax ) THEN 
     191!!$            WRITE(ctmp1,*) '   The chossen domain decomposition ', jpni, jpnj 
     192!!$            WRITE(ctmp2,*) '   exceeds the maximum number of ocean subdomains = ', inijmin 
     193!!$            WRITE(ctmp3,*) '   we suppressed ', jpni*jpnj - mppsize, ' land subdomains ' 
     194!!$            WRITE(ctmp4,*) '   BUT we had to keep ', mppsize - inijmin, ' land subdomains that are useless...' 
     195!!$            CALL ctl_warn( 'mpp_init:', '~~~~~~~~ ', ctmp1, ctmp2, ctmp3, ctmp4 ) 
     196!!$         ENDIF 
     197!!$      ENDIF 
    184198       
    185199      ! look for land mpi subdomains... 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/SBC/sbcfwb.F90

    r10365 r10397  
    7171      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   ztmsk_neg, ztmsk_pos, z_wgt ! 2D workspaces 
    7272      REAL(wp), ALLOCATABLE, DIMENSION(:,:) ::   ztmsk_tospread, zerp_cor    !   -      - 
     73      REAL(wp)   ,DIMENSION(1) ::   z_fwfprv   
     74      COMPLEX(wp),DIMENSION(1) ::   y_fwfnow   
    7375      !!---------------------------------------------------------------------- 
    7476      ! 
     
    102104         ! 
    103105         IF( MOD( kt-1, kn_fsbc ) == 0 ) THEN 
    104             z_fwf = glob_sum( 'sbcfwb', e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) - snwice_fmass(:,:) ) ) / area   ! sum over the global domain 
    105             zcoef = z_fwf * rcp 
    106             emp(:,:) = emp(:,:) - z_fwf              * tmask(:,:,1) 
     106            y_fwfnow(1) = local_sum( e1e2t(:,:) * ( emp(:,:) - rnf(:,:) + fwfisf(:,:) - snwice_fmass(:,:) ) ) 
     107            CALL mpp_delay_sum( 'sbcfwb', 'fwb', y_fwfnow(:), z_fwfprv(:), kt == nitend - nn_fsbc + 1 ) 
     108            z_fwfprv(1) = z_fwfprv(1) / area 
     109            zcoef = z_fwfprv(1) * rcp 
     110            emp(:,:) = emp(:,:) - z_fwfprv(1)        * tmask(:,:,1) 
    107111            qns(:,:) = qns(:,:) + zcoef * sst_m(:,:) * tmask(:,:,1) ! account for change to the heat budget due to fw correction 
    108112         ENDIF 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/lib_fortran.F90

    r10314 r10397  
    2727   PUBLIC   glob_sum      ! used in many places (masked with tmask_i) 
    2828   PUBLIC   glob_sum_full ! used in many places (masked with tmask_h, ie only over the halos) 
     29   PUBLIC   local_sum     ! used in trcrad, local operation before glob_sum_delay 
    2930   PUBLIC   DDPDD         ! also used in closea module 
    3031   PUBLIC   glob_min, glob_max 
     
    3940      MODULE PROCEDURE glob_sum_full_2d, glob_sum_full_3d 
    4041   END INTERFACE 
     42   INTERFACE local_sum 
     43      MODULE PROCEDURE local_sum_2d, local_sum_3d 
     44   END INTERFACE 
    4145   INTERFACE glob_min 
    4246      MODULE PROCEDURE glob_min_2d, glob_min_3d 
     
    126130#     undef DIM_3d 
    127131#  undef GLOBMINMAX_CODE 
     132 
     133!                          ! FUNCTION local_sum ! 
     134 
     135   FUNCTION local_sum_2d( ptab ) 
     136      !!---------------------------------------------------------------------- 
     137      REAL(wp),  INTENT(in   ) ::   ptab(:,:) ! array on which operation is applied 
     138      COMPLEX(wp)              ::  local_sum_2d 
     139      ! 
     140      !!----------------------------------------------------------------------- 
     141      ! 
     142      COMPLEX(wp)::   ctmp 
     143      REAL(wp)   ::   ztmp 
     144      INTEGER    ::   ji, jj    ! dummy loop indices 
     145      INTEGER    ::   ipi, ipj  ! dimensions 
     146      !!----------------------------------------------------------------------- 
     147      ! 
     148      ipi = SIZE(ptab,1)   ! 1st dimension 
     149      ipj = SIZE(ptab,2)   ! 2nd dimension 
     150      ! 
     151      ctmp = CMPLX( 0.e0, 0.e0, wp )   ! warning ctmp is cumulated 
     152 
     153      DO jj = 1, ipj 
     154         DO ji = 1, ipi 
     155            ztmp =  ptab(ji,jj) * tmask_i(ji,jj) 
     156            CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
     157         END DO 
     158      END DO 
     159      ! 
     160      local_sum_2d = ctmp 
     161        
     162   END FUNCTION local_sum_2d 
     163 
     164   FUNCTION local_sum_3d( ptab ) 
     165      !!---------------------------------------------------------------------- 
     166      REAL(wp),  INTENT(in   ) ::   ptab(:,:,:) ! array on which operation is applied 
     167      COMPLEX(wp)              ::  local_sum_3d 
     168      ! 
     169      !!----------------------------------------------------------------------- 
     170      ! 
     171      COMPLEX(wp)::   ctmp 
     172      REAL(wp)   ::   ztmp 
     173      INTEGER    ::   ji, jj, jk   ! dummy loop indices 
     174      INTEGER    ::   ipi, ipj, ipk    ! dimensions 
     175      !!----------------------------------------------------------------------- 
     176      ! 
     177      ipi = SIZE(ptab,1)   ! 1st dimension 
     178      ipj = SIZE(ptab,2)   ! 2nd dimension 
     179      ipk = SIZE(ptab,3)   ! 3rd dimension 
     180      ! 
     181      ctmp = CMPLX( 0.e0, 0.e0, wp )   ! warning ctmp is cumulated 
     182 
     183      DO jk = 1, ipk 
     184        DO jj = 1, ipj 
     185          DO ji = 1, ipi 
     186             ztmp =  ptab(ji,jj,jk) * tmask_i(ji,jj) 
     187             CALL DDPDD( CMPLX( ztmp, 0.e0, wp ), ctmp ) 
     188          END DO 
     189        END DO 
     190      END DO 
     191      ! 
     192      local_sum_3d = ctmp 
     193        
     194   END FUNCTION local_sum_3d 
    128195 
    129196 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/OCE/lib_fortran_generic.h90

    r10314 r10397  
    5050      ipk = K_SIZE(ptab)   ! 3rd dimension 
    5151      ! 
    52       ztmp = 0.e0 
    53       ctmp = CMPLX( 0.e0, 0.e0, wp ) 
     52      ctmp = CMPLX( 0.e0, 0.e0, wp )   ! warning ctmp is cumulated 
    5453    
    5554      DO jk = 1, ipk 
  • NEMO/branches/2018/dev_r10164_HPC09_ESIWACE_PREP_MERGE/src/TOP/trcrst.F90

    r10380 r10397  
    2323   USE iom 
    2424   USE daymod 
     25   USE lib_mpp 
    2526    
    2627   IMPLICIT NONE 
     
    116117      END DO 
    117118      ! 
     119      CALL iom_delay_rst( 'READ', 'TOP', numrtr )   ! read only TOP delayed global communication variables 
     120       
    118121   END SUBROUTINE trc_rst_read 
    119122 
     
    127130      !! 
    128131      INTEGER  :: jn 
    129       REAL(wp) :: zarak0 
    130132      !!---------------------------------------------------------------------- 
    131133      ! 
     
    141143      END DO 
    142144      ! 
     145      CALL iom_delay_rst( 'WRITE', 'TOP', numrtw )   ! save only TOP delayed global communication variables 
     146     
    143147      IF( kt == nitrst ) THEN 
    144148          CALL trc_rst_stat            ! statistics 
Note: See TracChangeset for help on using the changeset viewer.