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 10702 for NEMO/trunk/src/OCE/ICB/icbutl.F90 – NEMO

Ignore:
Timestamp:
2019-02-20T10:44:07+01:00 (5 years ago)
Author:
mathiot
Message:

merge branch fix_ticket2238_solution1 into the trunk (ticket #2238)

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/trunk/src/OCE/ICB/icbutl.F90

    r10691 r10702  
    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      ! 
     
    163152      pui = icb_utl_bilin_h( ui_e , pi, pj, 'U', .false. )    ! sea-ice velocities 
    164153      pvi = icb_utl_bilin_h( vi_e , pi, pj, 'V', .false. ) 
    165       phi = icb_utl_bilin_h( hicth, pi, pj, 'T', .true.  )    ! ice thickness 
     154      phi = icb_utl_bilin_h( hi_e , pi, pj, 'T', .true.  )    ! ice thickness 
    166155#else 
    167156      pui = 0._wp 
     
    206195         ! since we're looking for four T points containing quadrant we're in of  
    207196         ! current T cell 
    208          ii = MAX(1, INT( pi     )) 
    209          ij = MAX(1, INT( pj     ))    ! T-point 
     197         ii = MAX(0, INT( pi     )) 
     198         ij = MAX(0, INT( pj     ))    ! T-point 
    210199         zi = pi - REAL(ii,wp) 
    211200         zj = pj - REAL(ij,wp) 
    212201      CASE ( 'U' ) 
    213          ii = MAX(1, INT( pi-0.5 )) 
    214          ij = MAX(1, INT( pj     ))    ! U-point 
    215          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) 
    216205         zj = pj - REAL(ij,wp) 
    217206      CASE ( 'V' ) 
    218          ii = MAX(1, INT( pi     )) 
    219          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 
    220209         zi = pi - REAL(ii,wp) 
    221          zj = pj - 0.5 - REAL(ij,wp) 
     210         zj = pj - 0.5_wp - REAL(ij,wp) 
    222211      CASE ( 'F' ) 
    223          ii = MAX(1, INT( pi-0.5 )) 
    224          ij = MAX(1, INT( pj-0.5 ))    ! F-point 
    225          zi = pi - 0.5 - REAL(ii,wp) 
    226          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) 
    227216      END SELECT 
    228217      ! 
    229218      ! find position in this processor. Prevent near edge problems (see #1389) 
    230       ! 
    231       IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1   ;  
    232       ELSEIF( ii > mig(jpi) ) THEN   ;   ii = jpi ; 
    233       ELSE                           ;   ii = mi1(ii) 
    234       ENDIF 
    235       IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1   ; 
    236       ELSEIF( ij > mjg(jpj) ) THEN   ;   ij = jpj ; 
    237       ELSE                           ;   ij  = mj1(ij) 
    238       ENDIF 
    239       ! 
    240       IF( ii == jpi ) ii = ii-1 
    241       IF( ij == jpj ) ij = ij-1 
     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 
    242229      ! 
    243230      ! define mask array  
     
    402389      REAL(wp)                , INTENT(in) ::   pi, pj               ! targeted coordinates in (i,j) referential 
    403390      ! 
    404       INTEGER  ::   ii, ij, icase   ! local integer 
     391      INTEGER  ::   ii, ij, icase, ierr   ! local integer 
    405392      ! 
    406393      ! weights corresponding to corner points of a T cell quadrant 
     
    424411 
    425412      ! find position in this processor. Prevent near edge problems (see #1389) 
    426       IF    ( ii < mig( 1 ) ) THEN   ;   ii = 1 
    427       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 
    428417      ELSE                           ;   ii = mi1(ii) 
    429418      ENDIF 
    430       IF    ( ij < mjg( 1 ) ) THEN   ;   ij = 1 
    431       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 
    432421      ELSE                           ;   ij  = mj1(ij) 
    433422      ENDIF 
    434423      ! 
    435       IF( ii == jpi )   ii = ii-1       
    436       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)') 
    437428      ! 
    438429      IF(    0.0_wp <= zi .AND. zi < 0.5_wp   ) THEN 
     
    466457      ENDIF 
    467458      ! 
    468       icb_utl_bilin_e = ( ze01 * (1.-zi) + ze11 * zi ) *     zj    & 
    469          &            + ( 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) 
    470461      ! 
    471462   END FUNCTION icb_utl_bilin_e 
Note: See TracChangeset for help on using the changeset viewer.