Changeset 10726


Ignore:
Timestamp:
2019-02-27T16:06:35+01:00 (13 months ago)
Author:
mathiot
Message:

changes related to bug fixes described in tickets #2228, #2229, #2238 and #1595 in NEMO 4.0

Location:
NEMO/releases/release-4.0/src/OCE
Files:
5 edited

Legend:

Unmodified
Added
Removed
  • NEMO/releases/release-4.0/src/OCE/ICB/icb_oce.F90

    r10425 r10726  
    8686   ! particularly for MPP when iceberg can lie inside T grid but outside U, V, or f grid 
    8787   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   uo_e, vo_e 
    88    REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ff_e, tt_e, fr_e, hicth 
     88   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ff_e, tt_e, fr_e 
    8989   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ua_e, va_e 
    9090   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ssh_e 
     91   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   tmask_e, umask_e, vmask_e 
    9192#if defined key_si3 || defined key_cice 
    92    REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   ui_e, vi_e 
     93   REAL(wp), PUBLIC, DIMENSION(:,:), ALLOCATABLE ::   hi_e, ui_e, vi_e 
    9394#endif 
    9495 
     
    169170      ! 
    170171      ! expanded arrays for bilinear interpolation 
    171       ALLOCATE( uo_e(0:jpi+1,0:jpj+1) , ua_e(0:jpi+1,0:jpj+1) ,   & 
    172          &      vo_e(0:jpi+1,0:jpj+1) , va_e(0:jpi+1,0:jpj+1) ,   & 
     172      ALLOCATE( uo_e(0:jpi+1,0:jpj+1) , ua_e(0:jpi+1,0:jpj+1) ,    & 
     173         &      vo_e(0:jpi+1,0:jpj+1) , va_e(0:jpi+1,0:jpj+1) ,    & 
    173174#if defined key_si3 || defined key_cice 
    174175         &      ui_e(0:jpi+1,0:jpj+1) ,                            & 
    175176         &      vi_e(0:jpi+1,0:jpj+1) ,                            & 
     177         &      hi_e(0:jpi+1,0:jpj+1) ,                            & 
    176178#endif 
    177179         &      ff_e(0:jpi+1,0:jpj+1) , fr_e(0:jpi+1,0:jpj+1)  ,   & 
    178180         &      tt_e(0:jpi+1,0:jpj+1) , ssh_e(0:jpi+1,0:jpj+1) ,   & 
    179          &      hicth(0:jpi+1,0:jpj+1),                            & 
    180181         &      first_width(nclasses) , first_length(nclasses) ,   & 
    181182         &      src_calving (jpi,jpj) ,                            & 
    182183         &      src_calving_hflx(jpi,jpj) , STAT=ill) 
     184      icb_alloc = icb_alloc + ill 
     185 
     186      ALLOCATE( tmask_e(0:jpi+1,0:jpj+1), umask_e(0:jpi+1,0:jpj+1), vmask_e(0:jpi+1,0:jpj+1), & 
     187         &      STAT=ill) 
    183188      icb_alloc = icb_alloc + ill 
    184189 
  • NEMO/releases/release-4.0/src/OCE/ICB/icbclv.F90

    r10425 r10726  
    2424   USE icbdia         ! iceberg diagnostics 
    2525   USE icbutl         ! iceberg utility routines 
     26   USE icb_oce        ! iceberg parameters  
    2627 
    2728   IMPLICIT NONE 
     
    5758      ! this assumes that input is given as equivalent water flux so that pure water density is appropriate 
    5859 
    59       zfact = ( (1000._wp)**3 / ( NINT(rday) * nyear_len(1) ) ) * 850._wp 
    60       berg_grid%calving(:,:) = src_calving(:,:) * tmask_i(:,:) * zfact 
     60      zfact = ( (1000._wp)**3 / ( NINT(rday) * nyear_len(1) ) ) * rn_rho_bergs 
     61      berg_grid%calving(:,:) = src_calving(:,:) * zfact * tmask_i(:,:) * tmask(:,:,1) 
    6162 
    6263      ! Heat in units of W/m2, and mask (just in case) 
    63       berg_grid%calving_hflx(:,:) = src_calving_hflx(:,:) * tmask_i(:,:) 
     64      berg_grid%calving_hflx(:,:) = src_calving_hflx(:,:) * tmask_i(:,:) * tmask(:,:,1) 
    6465 
    6566      IF( ll_first_call .AND. .NOT. l_restarted_bergs ) THEN      ! This is a hack to simplify initialization 
  • NEMO/releases/release-4.0/src/OCE/ICB/icbini.F90

    r10570 r10726  
    7474      !                          ! allocate gridded fields 
    7575      IF( icb_alloc() /= 0 )   CALL ctl_stop( 'STOP', 'icb_alloc : unable to allocate arrays' ) 
    76  
     76      ! 
     77      !                          ! initialised variable with extra haloes to zero 
     78      uo_e(:,:) = 0._wp   ;   vo_e(:,:) = 0._wp   ; 
     79      ua_e(:,:) = 0._wp   ;   va_e(:,:) = 0._wp   ; 
     80      ff_e(:,:) = 0._wp   ;   tt_e(:,:) = 0._wp   ; 
     81      fr_e(:,:) = 0._wp   ; 
     82#if defined key_si3 
     83      hi_e(:,:) = 0._wp   ; 
     84      ui_e(:,:) = 0._wp   ;   vi_e(:,:) = 0._wp   ; 
     85#endif 
     86      ssh_e(:,:) = 0._wp  ;  
     87      ! 
    7788      !                          ! open ascii output file or files for iceberg status information 
    7889      !                          ! note that we choose to do this on all processors since we cannot 
     
    224235      src_calving_hflx(:,:) = 0._wp 
    225236 
     237      ! definition of extended surface masked needed by icb_bilin_h 
     238      tmask_e(:,:) = 0._wp   ;   tmask_e(1:jpi,1:jpj) = tmask(:,:,1) 
     239      umask_e(:,:) = 0._wp   ;   umask_e(1:jpi,1:jpj) = umask(:,:,1) 
     240      vmask_e(:,:) = 0._wp   ;   vmask_e(1:jpi,1:jpj) = vmask(:,:,1) 
     241      CALL lbc_lnk_icb( 'icbini', tmask_e, 'T', +1._wp, 1, 1 ) 
     242      CALL lbc_lnk_icb( 'icbini', umask_e, 'T', +1._wp, 1, 1 ) 
     243      CALL lbc_lnk_icb( 'icbini', vmask_e, 'T', +1._wp, 1, 1 ) 
     244      ! 
    226245      ! assign each new iceberg with a unique number constructed from the processor number 
    227246      ! and incremented by the total number of processors 
  • NEMO/releases/release-4.0/src/OCE/ICB/icbutl.F90

    r10570 r10726  
    7070      ! and ssh which is used to calculate gradients 
    7171 
    72       uo_e(:,:) = 0._wp   ;   uo_e(1:jpi,1:jpj) = ssu_m(:,:) * umask(:,:,1) 
    73       vo_e(:,:) = 0._wp   ;   vo_e(1:jpi,1:jpj) = ssv_m(:,:) * vmask(:,:,1) 
    74       ff_e(:,:) = 0._wp   ;   ff_e(1:jpi,1:jpj) = ff_f (:,:)  
    75       tt_e(:,:) = 0._wp   ;   tt_e(1:jpi,1:jpj) = sst_m(:,:) 
    76       fr_e(:,:) = 0._wp   ;   fr_e(1:jpi,1:jpj) = fr_i (:,:) 
    77       ua_e(:,:) = 0._wp   ;   ua_e(1:jpi,1:jpj) = utau (:,:) * umask(:,:,1) ! maybe mask useless because mask applied in sbcblk 
    78       va_e(:,:) = 0._wp   ;   va_e(1:jpi,1:jpj) = vtau (:,:) * vmask(:,:,1) ! maybe mask useless because mask applied in sbcblk 
     72      uo_e(1:jpi,1:jpj) = ssu_m(:,:) * umask(:,:,1) 
     73      vo_e(1:jpi,1:jpj) = ssv_m(:,:) * vmask(:,:,1) 
     74      ff_e(1:jpi,1:jpj) = ff_f (:,:)  
     75      tt_e(1:jpi,1:jpj) = sst_m(:,:) 
     76      fr_e(1:jpi,1:jpj) = fr_i (:,:) 
     77      ua_e(1:jpi,1:jpj) = utau (:,:) * umask(:,:,1) ! maybe mask useless because mask applied in sbcblk 
     78      va_e(1:jpi,1:jpj) = vtau (:,:) * vmask(:,:,1) ! maybe mask useless because mask applied in sbcblk 
    7979      ! 
    8080      CALL lbc_lnk_icb( 'icbutl', uo_e, 'U', -1._wp, 1, 1 ) 
     
    8686      CALL lbc_lnk_icb( 'icbutl', tt_e, 'T', +1._wp, 1, 1 ) 
    8787#if defined key_si3 
    88       hicth(:,:) = 0._wp ;  hicth(1:jpi,1:jpj) = hm_i (:,:)   
    89       ui_e(:,:) = 0._wp ;   ui_e(1:jpi, 1:jpj) = u_ice(:,:) 
    90       vi_e(:,:) = 0._wp ;   vi_e(1:jpi, 1:jpj) = v_ice(:,:) 
     88      hi_e(1:jpi, 1:jpj) = hm_i (:,:)   
     89      ui_e(1:jpi, 1:jpj) = u_ice(:,:) 
     90      vi_e(1:jpi, 1:jpj) = v_ice(:,:) 
    9191      !       
    9292      ! compute ssh slope using ssh_lead if embedded 
    9393      zssh_lead_m(:,:) = ice_var_sshdyn(ssh_m, snwice_mass, snwice_mass_b) 
    94       ssh_e(:,:) = 0._wp ;  ssh_e(1:jpi, 1:jpj) = zssh_lead_m(:,:) * tmask(:,:,1) 
    95       ! 
    96       CALL lbc_lnk_icb( 'icbutl', hicth, 'T', +1._wp, 1, 1 ) 
     94      ssh_e(1:jpi, 1:jpj) = zssh_lead_m(:,:) * tmask(:,:,1) 
     95      ! 
     96      CALL lbc_lnk_icb( 'icbutl', hi_e , 'T', +1._wp, 1, 1 ) 
    9797      CALL lbc_lnk_icb( 'icbutl', ui_e , 'U', -1._wp, 1, 1 ) 
    9898      CALL lbc_lnk_icb( 'icbutl', vi_e , 'V', -1._wp, 1, 1 ) 
    9999#else 
    100       ssh_e(:,:) = 0._wp ;  ssh_e(1:jpi, 1:jpj) = ssh_m(:,:) * tmask(:,:,1) 
     100      ssh_e(1:jpi, 1:jpj) = ssh_m(:,:) * tmask(:,:,1) 
    101101#endif 
    102  
    103       !! special for ssh which is used to calculate slope 
    104       !! so fudge some numbers all the way around the boundary 
    105       ssh_e(0    ,    :) = ssh_e(1  ,  :) 
    106       ssh_e(jpi+1,    :) = ssh_e(jpi,  :) 
    107       ssh_e(:    ,    0) = ssh_e(:  ,  1) 
    108       ssh_e(:    ,jpj+1) = ssh_e(:  ,jpj) 
    109       ssh_e(0,0)         = ssh_e(1,1) 
    110       ssh_e(jpi+1,0)     = ssh_e(jpi,1) 
    111       ssh_e(0,jpj+1)     = ssh_e(1,jpj) 
    112       ssh_e(jpi+1,jpj+1) = ssh_e(jpi,jpj) 
    113102      CALL lbc_lnk_icb( 'icbutl', ssh_e, 'T', +1._wp, 1, 1 ) 
    114103      ! 
     
    131120      !!             is half the off shore value, wile the normal-to-the-coast value is zero. 
    132121      !!             This is OK as a starting point. 
     122      !!       !!pm  HARD CODED: - rho_air now computed in sbcblk (what are the effect ?) 
     123      !!                         - drag coefficient (should it be namelist parameter ?) 
    133124      !! 
    134125      !!---------------------------------------------------------------------- 
     
    142133      !!---------------------------------------------------------------------- 
    143134 
    144       pe1 = icb_utl_bilin_e( e1t, e1u, e1v, e1f, pi, pj )     ! scale factors 
     135      pe1 = icb_utl_bilin_e( e1t, e1u, e1v, e1f, pi, pj )      ! scale factors 
    145136      pe2 = icb_utl_bilin_e( e2t, e2u, e2v, e2f, pi, pj ) 
    146137      ! 
    147       puo  = icb_utl_bilin_h( uo_e, pi, pj, 'U' )             ! ocean velocities 
    148       pvo  = icb_utl_bilin_h( vo_e, pi, pj, 'V' ) 
    149       psst = icb_utl_bilin_h( tt_e, pi, pj, 'T' )             ! SST 
    150       pcn  = icb_utl_bilin_h( fr_e , pi, pj, 'T' )            ! ice concentration 
    151       pff  = icb_utl_bilin_h( ff_e , pi, pj, 'F' )            ! Coriolis parameter 
    152       ! 
    153       pua  = icb_utl_bilin_h( ua_e , pi, pj, 'U' )            ! 10m wind 
    154       pva  = icb_utl_bilin_h( va_e , pi, pj, 'V' )            ! here (ua,va) are stress => rough conversion from stress to speed 
    155       zcd  = 1.22_wp * 1.5e-3_wp                              ! air density * drag coefficient 
     138      puo  = icb_utl_bilin_h( uo_e, pi, pj, 'U', .false.  )    ! ocean velocities 
     139      pvo  = icb_utl_bilin_h( vo_e, pi, pj, 'V', .false. ) 
     140      psst = icb_utl_bilin_h( tt_e, pi, pj, 'T', .true.   )    ! SST 
     141      pcn  = icb_utl_bilin_h( fr_e, pi, pj, 'T', .true.   )    ! ice concentration 
     142      pff  = icb_utl_bilin_h( ff_e, pi, pj, 'F', .false.  )    ! Coriolis parameter 
     143      ! 
     144      pua  = icb_utl_bilin_h( ua_e, pi, pj, 'U', .true.   )    ! 10m wind 
     145      pva  = icb_utl_bilin_h( va_e, pi, pj, 'V', .true.   )    ! here (ua,va) are stress => rough conversion from stress to speed 
     146      zcd  = 1.22_wp * 1.5e-3_wp                               ! air density * drag coefficient  
    156147      zmod = 1._wp / MAX(  1.e-20, SQRT(  zcd * SQRT( pua*pua + pva*pva)  )  ) 
    157148      pua  = pua * zmod                                       ! note: stress module=0 necessarly implies ua=va=0 
     
    159150 
    160151#if defined key_si3 
    161       pui = icb_utl_bilin_h( ui_e , pi, pj, 'U' )              ! sea-ice velocities 
    162       pvi = icb_utl_bilin_h( vi_e , pi, pj, 'V' ) 
    163       phi = icb_utl_bilin_h( hicth, pi, pj, 'T' )              ! ice thickness 
     152      pui = icb_utl_bilin_h( ui_e , pi, pj, 'U', .false. )    ! sea-ice velocities 
     153      pvi = icb_utl_bilin_h( vi_e , pi, pj, 'V', .false. ) 
     154      phi = icb_utl_bilin_h( hi_e , pi, pj, 'T', .true.  )    ! ice thickness 
    164155#else 
    165156      pui = 0._wp 
     
    169160 
    170161      ! Estimate SSH gradient in i- and j-direction (centred evaluation) 
    171       pssh_i = ( icb_utl_bilin_h( ssh_e, pi+0.1_wp, pj, 'T' ) -   & 
    172          &       icb_utl_bilin_h( ssh_e, pi-0.1_wp, pj, 'T' )  ) / ( 0.2_wp * pe1 ) 
    173       pssh_j = ( icb_utl_bilin_h( ssh_e, pi, pj+0.1_wp, 'T' ) -   & 
    174          &       icb_utl_bilin_h( ssh_e, pi, pj-0.1_wp, 'T' )  ) / ( 0.2_wp * pe2 ) 
     162      pssh_i = ( icb_utl_bilin_h( ssh_e, pi+0.1_wp, pj, 'T', .true. ) -   & 
     163         &       icb_utl_bilin_h( ssh_e, pi-0.1_wp, pj, 'T', .true. )  ) / ( 0.2_wp * pe1 ) 
     164      pssh_j = ( icb_utl_bilin_h( ssh_e, pi, pj+0.1_wp, 'T', .true. ) -   & 
     165         &       icb_utl_bilin_h( ssh_e, pi, pj-0.1_wp, 'T', .true. )  ) / ( 0.2_wp * pe2 ) 
    175166      ! 
    176167   END SUBROUTINE icb_utl_interp 
    177168 
    178169 
    179    REAL(wp) FUNCTION icb_utl_bilin_h( pfld, pi, pj, cd_type ) 
     170   REAL(wp) FUNCTION icb_utl_bilin_h( pfld, pi, pj, cd_type, plmask ) 
    180171      !!---------------------------------------------------------------------- 
    181172      !!                  ***  FUNCTION icb_utl_bilin  *** 
     
    191182      REAL(wp)                            , INTENT(in) ::   pi, pj    ! targeted coordinates in (i,j) referential 
    192183      CHARACTER(len=1)                    , INTENT(in) ::   cd_type   ! type of pfld array grid-points: = T , U , V or F points 
     184      LOGICAL                             , INTENT(in) ::   plmask    ! special treatment of mask point 
    193185      ! 
    194186      INTEGER  ::   ii, ij   ! local integer 
    195187      REAL(wp) ::   zi, zj   ! local real 
     188      REAL(wp) :: zw1, zw2, zw3, zw4 
     189      REAL(wp), DIMENSION(4) :: zmask 
    196190      !!---------------------------------------------------------------------- 
    197191      ! 
     
    201195         ! since we're looking for four T points containing quadrant we're in of  
    202196         ! current T cell 
    203          ii = MAX(1, INT( pi     )) 
    204          ij = MAX(1, INT( pj     ))    ! T-point 
     197         ii = MAX(0, INT( pi     )) 
     198         ij = MAX(0, INT( pj     ))    ! T-point 
    205199         zi = pi - REAL(ii,wp) 
    206200         zj = pj - REAL(ij,wp) 
    207201      CASE ( 'U' ) 
    208          ii = MAX(1, INT( pi-0.5 )) 
    209          ij = MAX(1, INT( pj     ))    ! U-point 
    210          zi = pi - 0.5 - REAL(ii,wp) 
     202         ii = MAX(0, INT( pi-0.5_wp )) 
     203         ij = MAX(0, INT( pj     ))    ! U-point 
     204         zi = pi - 0.5_wp - REAL(ii,wp) 
    211205         zj = pj - REAL(ij,wp) 
    212206      CASE ( 'V' ) 
    213          ii = MAX(1, INT( pi     )) 
    214          ij = MAX(1, INT( pj-0.5 ))    ! V-point 
     207         ii = MAX(0, INT( pi     )) 
     208         ij = MAX(0, INT( pj-0.5_wp ))    ! V-point 
    215209         zi = pi - REAL(ii,wp) 
    216          zj = pj - 0.5 - REAL(ij,wp) 
     210         zj = pj - 0.5_wp - REAL(ij,wp) 
    217211      CASE ( 'F' ) 
    218          ii = MAX(1, INT( pi-0.5 )) 
    219          ij = MAX(1, INT( pj-0.5 ))    ! F-point 
    220          zi = pi - 0.5 - REAL(ii,wp) 
    221          zj = pj - 0.5 - REAL(ij,wp) 
     212         ii = MAX(0, INT( pi-0.5_wp )) 
     213         ij = MAX(0, INT( pj-0.5_wp ))    ! F-point 
     214         zi = pi - 0.5_wp - REAL(ii,wp) 
     215         zj = pj - 0.5_wp - REAL(ij,wp) 
    222216      END SELECT 
    223217      ! 
    224218      ! find position in this processor. Prevent near edge problems (see #1389) 
    225       ! 
    226       IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1 
    227       ELSEIF( ii > mig(jpi) ) THEN   ;   ii = jpi 
    228       ELSE                           ;   ii = mi1(ii) 
    229       ENDIF 
    230       IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1 
    231       ELSEIF( ij > mjg(jpj) ) THEN   ;   ij = jpj 
    232       ELSE                           ;   ij  = mj1(ij) 
    233       ENDIF 
    234       ! 
    235       IF( ii == jpi )   ii = ii-1       
    236       IF( ij == jpj )   ij = ij-1 
    237       ! 
    238       icb_utl_bilin_h = ( pfld(ii,ij  ) * (1.-zi) + pfld(ii+1,ij  ) * zi ) * (1.-zj)   & 
    239          &            + ( pfld(ii,ij+1) * (1.-zi) + pfld(ii+1,ij+1) * zi ) *     zj 
     219      ! (PM) will be useless if extra halo is used in NEMO 
     220      ! 
     221      IF    ( ii <= mig(1)-1 ) THEN   ;   ii = 0 
     222      ELSEIF( ii  > mig(jpi) ) THEN   ;   ii = jpi 
     223      ELSE                            ;   ii = mi1(ii) 
     224      ENDIF 
     225      IF    ( ij <= mjg(1)-1 ) THEN   ;   ij = 0 
     226      ELSEIF( ij  > mjg(jpj) ) THEN   ;   ij = jpj 
     227      ELSE                            ;   ij = mj1(ij) 
     228      ENDIF 
     229      ! 
     230      ! define mask array  
     231      IF (plmask) THEN 
     232         ! land value is not used in the interpolation 
     233         SELECT CASE ( cd_type ) 
     234         CASE ( 'T' ) 
     235            zmask = (/tmask_e(ii,ij), tmask_e(ii+1,ij), tmask_e(ii,ij+1), tmask_e(ii+1,ij+1)/) 
     236         CASE ( 'U' ) 
     237            zmask = (/umask_e(ii,ij), umask_e(ii+1,ij), umask_e(ii,ij+1), umask_e(ii+1,ij+1)/) 
     238         CASE ( 'V' ) 
     239            zmask = (/vmask_e(ii,ij), vmask_e(ii+1,ij), vmask_e(ii,ij+1), vmask_e(ii+1,ij+1)/) 
     240         CASE ( 'F' ) 
     241            ! F case only used for coriolis, ff_f is not mask so zmask = 1 
     242            zmask = 1. 
     243         END SELECT 
     244      ELSE 
     245         ! land value is used during interpolation 
     246         zmask = 1. 
     247      END iF 
     248      ! 
     249      ! compute weight 
     250      zw1 = zmask(1) * (1._wp-zi) * (1._wp-zj) 
     251      zw2 = zmask(2) *        zi  * (1._wp-zj) 
     252      zw3 = zmask(3) * (1._wp-zi) *        zj 
     253      zw4 = zmask(4) *        zi  *        zj 
     254      ! 
     255      ! compute interpolated value 
     256      icb_utl_bilin_h = ( pfld(ii,ij)*zw1 + pfld(ii+1,ij)*zw2 + pfld(ii,ij+1)*zw3 + pfld(ii+1,ij+1)*zw4 ) / MAX(1.e-20, zw1+zw2+zw3+zw4)  
    240257      ! 
    241258   END FUNCTION icb_utl_bilin_h 
     
    372389      REAL(wp)                , INTENT(in) ::   pi, pj               ! targeted coordinates in (i,j) referential 
    373390      ! 
    374       INTEGER  ::   ii, ij, icase   ! local integer 
     391      INTEGER  ::   ii, ij, icase, ierr   ! local integer 
    375392      ! 
    376393      ! weights corresponding to corner points of a T cell quadrant 
     
    394411 
    395412      ! find position in this processor. Prevent near edge problems (see #1389) 
    396       IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1 
    397       ELSEIF( ii > mig(jpi) ) THEN   ;   ii = jpi 
     413      ! 
     414      ierr = 0 
     415      IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1       ; ierr = ierr + 1 
     416      ELSEIF( ii > mig(jpi) ) THEN   ;   ii = jpi     ; ierr = ierr + 1 
    398417      ELSE                           ;   ii = mi1(ii) 
    399418      ENDIF 
    400       IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1 
    401       ELSEIF( ij > mjg(jpj) ) THEN   ;   ij = jpj 
     419      IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1       ; ierr = ierr + 1 
     420      ELSEIF( ij > mjg(jpj) ) THEN   ;   ij = jpj     ; ierr = ierr + 1 
    402421      ELSE                           ;   ij  = mj1(ij) 
    403422      ENDIF 
    404423      ! 
    405       IF( ii == jpi )   ii = ii-1       
    406       IF( ij == jpj )   ij = ij-1 
     424      IF( ii == jpi ) THEN ; ii = ii-1 ; ierr = ierr + 1 ; END IF      
     425      IF( ij == jpj ) THEN ; ij = ij-1 ; ierr = ierr + 1 ; END IF 
     426      ! 
     427      IF ( ierr > 0 ) CALL ctl_stop('STOP','icb_utl_bilin_e: an icebergs coordinates is out of valid range (out of bound error)') 
    407428      ! 
    408429      IF(    0.0_wp <= zi .AND. zi < 0.5_wp   ) THEN 
     
    436457      ENDIF 
    437458      ! 
    438       icb_utl_bilin_e = ( ze01 * (1.-zi) + ze11 * zi ) *     zj    & 
    439          &            + ( ze00 * (1.-zi) + ze10 * zi ) * (1.-zj) 
     459      icb_utl_bilin_e = ( ze01 * (1._wp-zi) + ze11 * zi ) *        zj    & 
     460         &            + ( ze00 * (1._wp-zi) + ze10 * zi ) * (1._wp-zj) 
    440461      ! 
    441462   END FUNCTION icb_utl_bilin_e 
  • NEMO/releases/release-4.0/src/OCE/LBC/mpp_loc_generic.h90

    r10425 r10726  
    1717#      define MPI_OPERATION mpi_maxloc 
    1818#      define LOC_OPERATION MAXLOC 
     19#      define ERRVAL -HUGE 
    1920#   endif 
    2021#   if defined OPERATION_MINLOC 
    2122#      define MPI_OPERATION mpi_minloc 
    2223#      define LOC_OPERATION MINLOC 
     24#      define ERRVAL HUGE 
    2325#   endif 
    2426 
     
    4244      ! 
    4345      idim = SIZE(kindex) 
    44       ALLOCATE ( ilocs(idim) ) 
    4546      ! 
    46       ilocs = LOC_OPERATION( ARRAY_IN(:,:,:) , mask= MASK_IN(:,:,:) == 1._wp ) 
    47       zmin  = ARRAY_IN(ilocs(1),ilocs(2),ilocs(3)) 
    48       ! 
    49       kindex(1) = ilocs(1) + nimpp - 1 
     47      IF ( ALL(MASK_IN(:,:,:) /= 1._wp) ) THEN 
     48         ! special case for land processors 
     49         zmin = ERRVAL(zmin) 
     50         index0 = 0 
     51      ELSE 
     52         ALLOCATE ( ilocs(idim) ) 
     53         ! 
     54         ilocs = LOC_OPERATION( ARRAY_IN(:,:,:) , mask= MASK_IN(:,:,:) == 1._wp ) 
     55         zmin  = ARRAY_IN(ilocs(1),ilocs(2),ilocs(3)) 
     56         ! 
     57         kindex(1) = mig( ilocs(1) ) 
    5058#  if defined DIM_2d || defined DIM_3d    /* avoid warning when kindex has 1 element */ 
    51       kindex(2) = ilocs(2) + njmpp - 1 
     59         kindex(2) = mjg( ilocs(2) ) 
    5260#  endif 
    5361#  if defined DIM_3d                      /* avoid warning when kindex has 2 elements */ 
    54       kindex(3) = ilocs(3) 
     62         kindex(3) = ilocs(3) 
    5563#  endif 
    56       !  
    57       DEALLOCATE (ilocs) 
    58       ! 
    59       index0 = kindex(1)-1   ! 1d index starting at 0 
     64         !  
     65         DEALLOCATE (ilocs) 
     66         ! 
     67         index0 = kindex(1)-1   ! 1d index starting at 0 
    6068#  if defined DIM_2d || defined DIM_3d   /* avoid warning when kindex has 1 element */ 
    61       index0 = index0 + jpiglo * (kindex(2)-1) 
     69         index0 = index0 + jpiglo * (kindex(2)-1) 
    6270#  endif 
    6371#  if defined DIM_3d                     /* avoid warning when kindex has 2 elements */ 
    64       index0 = index0 + jpiglo * jpjglo * (kindex(3)-1) 
     72         index0 = index0 + jpiglo * jpjglo * (kindex(3)-1) 
    6573#  endif 
     74      END IF 
    6675      zain(1,:) = zmin 
    6776      zain(2,:) = REAL(index0, wp) 
     
    98107#undef LOC_OPERATION 
    99108#undef INDEX_TYPE 
     109#undef ERRVAL 
Note: See TracChangeset for help on using the changeset viewer.