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 10288 for NEMO/branches/2018/dev_r9866_HPC_03_globcom/src/OCE/BDY/bdyice.F90 – NEMO

Ignore:
Timestamp:
2018-11-07T18:25:49+01:00 (5 years ago)
Author:
francesca
Message:

reduce global communications, see #2010

File:
1 edited

Legend:

Unmodified
Added
Removed
  • NEMO/branches/2018/dev_r9866_HPC_03_globcom/src/OCE/BDY/bdyice.F90

    • Property svn:keywords set to Id
    r9657 r10288  
    1818   USE ice             ! sea-ice: variables 
    1919   USE icevar          ! sea-ice: operations 
    20    USE iceitd          ! sea-ice: rebining 
     20   USE icecor          ! sea-ice: corrections 
    2121   USE icectl          ! sea-ice: control prints 
    2222   USE phycst          ! physical constant 
     
    4141   !!---------------------------------------------------------------------- 
    4242   !! NEMO/OCE 4.0 , NEMO Consortium (2018) 
    43    !! $Id: bdyice.F90 8306 2017-07-10 10:18:03Z clem $ 
    44    !! Software governed by the CeCILL licence (./LICENSE) 
     43   !! $Id$ 
     44   !! Software governed by the CeCILL license (see ./LICENSE) 
    4545   !!---------------------------------------------------------------------- 
    4646CONTAINS 
     
    5050      !!                  ***  SUBROUTINE bdy_ice  *** 
    5151      !! 
    52       !! ** Purpose : - Apply open boundary conditions for ice (SI3) 
     52      !! ** Purpose : Apply open boundary conditions for sea ice 
    5353      !! 
    5454      !!---------------------------------------------------------------------- 
    5555      INTEGER, INTENT(in) ::   kt   ! Main time step counter 
    5656      ! 
    57       INTEGER ::   ib_bdy   ! Loop index 
     57      INTEGER ::   jbdy   ! BDY set index 
    5858      !!---------------------------------------------------------------------- 
    5959      ! 
    60       IF( ln_timing )   CALL timing_start('bdy_ice') 
     60      IF( ln_timing )   CALL timing_start('bdy_ice_thd') 
    6161      ! 
    6262      CALL ice_var_glo2eqv 
    6363      ! 
    64       DO ib_bdy = 1, nb_bdy 
    65          ! 
    66          SELECT CASE( cn_ice(ib_bdy) ) 
     64      DO jbdy = 1, nb_bdy 
     65         ! 
     66         SELECT CASE( cn_ice(jbdy) ) 
    6767         CASE('none')   ;   CYCLE 
    68          CASE('frs' )   ;   CALL bdy_ice_frs( idx_bdy(ib_bdy), dta_bdy(ib_bdy), kt, ib_bdy ) 
     68         CASE('frs' )   ;   CALL bdy_ice_frs( idx_bdy(jbdy), dta_bdy(jbdy), kt, jbdy ) 
    6969         CASE DEFAULT 
    7070            CALL ctl_stop( 'bdy_ice : unrecognised option for open boundaries for ice fields' ) 
     
    7373      END DO 
    7474      ! 
    75       CALL ice_var_zapsmall 
     75      CALL ice_cor( kt , 0 )      ! -- In case categories are out of bounds, do a remapping 
     76      !                           !    i.e. inputs have not the same ice thickness distribution (set by rn_himean) 
     77      !                           !         than the regional simulation 
    7678      CALL ice_var_agg(1) 
    7779      ! 
    7880      IF( ln_icectl )   CALL ice_prt( kt, iiceprt, jiceprt, 1, ' - ice thermo bdy - ' ) 
    79       IF( ln_timing )   CALL timing_stop('bdy_ice') 
     81      IF( ln_timing )   CALL timing_stop('bdy_ice_thd') 
    8082      ! 
    8183   END SUBROUTINE bdy_ice 
    8284 
    8385 
    84    SUBROUTINE bdy_ice_frs( idx, dta, kt, ib_bdy ) 
     86   SUBROUTINE bdy_ice_frs( idx, dta, kt, jbdy ) 
    8587      !!------------------------------------------------------------------------------ 
    8688      !!                 ***  SUBROUTINE bdy_ice_frs  *** 
    8789      !!                     
    88       !! ** Purpose : Apply the Flow Relaxation Scheme for sea-ice fields in the case  
    89       !!              of unstructured open boundaries. 
     90      !! ** Purpose : Apply the Flow Relaxation Scheme for sea-ice fields 
    9091      !!  
    9192      !! Reference : Engedahl H., 1995: Use of the flow relaxation scheme in a three- 
     
    9596      TYPE(OBC_DATA),  INTENT(in) ::   dta     ! OBC external data 
    9697      INTEGER,         INTENT(in) ::   kt      ! main time-step counter 
    97       INTEGER,         INTENT(in) ::   ib_bdy  ! BDY set index 
     98      INTEGER,         INTENT(in) ::   jbdy    ! BDY set index 
    9899      ! 
    99100      INTEGER  ::   jpbound            ! 0 = incoming ice 
    100101      !                                ! 1 = outgoing ice 
    101       INTEGER  ::   jb, jk, jgrd, jl   ! dummy loop indices 
    102       INTEGER  ::   ji, jj, ii, ij     ! local scalar 
     102      INTEGER  ::   i_bdy, jgrd        ! dummy loop indices 
     103      INTEGER  ::   ji, jj, jk, jl, ib, jb 
    103104      REAL(wp) ::   zwgt, zwgt1        ! local scalar 
    104105      REAL(wp) ::   ztmelts, zdh 
     
    108109      ! 
    109110      DO jl = 1, jpl 
    110          DO jb = 1, idx%nblenrim(jgrd) 
    111             ji    = idx%nbi(jb,jgrd) 
    112             jj    = idx%nbj(jb,jgrd) 
    113             zwgt  = idx%nbw(jb,jgrd) 
    114             zwgt1 = 1.e0 - idx%nbw(jb,jgrd) 
    115             a_i(ji,jj,jl) = ( a_i(ji,jj,jl) * zwgt1 + dta%a_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Leads fraction  
    116             h_i(ji,jj,jl) = ( h_i(ji,jj,jl) * zwgt1 + dta%h_i(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Ice depth  
    117             h_s(ji,jj,jl) = ( h_s(ji,jj,jl) * zwgt1 + dta%h_s(jb,jl) * zwgt ) * tmask(ji,jj,1)  ! Snow depth 
     111         DO i_bdy = 1, idx%nblenrim(jgrd) 
     112            ji    = idx%nbi(i_bdy,jgrd) 
     113            jj    = idx%nbj(i_bdy,jgrd) 
     114            zwgt  = idx%nbw(i_bdy,jgrd) 
     115            zwgt1 = 1.e0 - idx%nbw(i_bdy,jgrd) 
     116            a_i(ji,jj,jl) = ( a_i(ji,jj,jl) * zwgt1 + dta%a_i(i_bdy,jl) * zwgt ) * tmask(ji,jj,1)  ! Leads fraction  
     117            h_i(ji,jj,jl) = ( h_i(ji,jj,jl) * zwgt1 + dta%h_i(i_bdy,jl) * zwgt ) * tmask(ji,jj,1)  ! Ice depth  
     118            h_s(ji,jj,jl) = ( h_s(ji,jj,jl) * zwgt1 + dta%h_s(i_bdy,jl) * zwgt ) * tmask(ji,jj,1)  ! Snow depth 
    118119 
    119120            ! ----------------- 
     
    124125 
    125126            ! Then, a) transfer the snow excess into the ice (different from icethd_dh) 
    126             zdh = MAX( 0._wp, ( rhosn * h_s(ji,jj,jl) + ( rhoic - rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 ) 
     127            zdh = MAX( 0._wp, ( rhos * h_s(ji,jj,jl) + ( rhoi - rau0 ) * h_i(ji,jj,jl) ) * r1_rau0 ) 
    127128            ! Or, b) transfer all the snow into ice (if incoming ice is likely to melt as it comes into a warmer environment) 
    128             !zdh = MAX( 0._wp, h_s(ji,jj,jl) * rhosn / rhoic ) 
     129            !zdh = MAX( 0._wp, h_s(ji,jj,jl) * rhos / rhoi ) 
    129130 
    130131            ! recompute h_i, h_s 
    131132            h_i(ji,jj,jl) = MIN( hi_max(jl), h_i(ji,jj,jl) + zdh ) 
    132             h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoic / rhosn )  
     133            h_s(ji,jj,jl) = MAX( 0._wp, h_s(ji,jj,jl) - zdh * rhoi / rhos )  
    133134 
    134135         ENDDO 
    135          CALL lbc_bdy_lnk( a_i(:,:,jl), 'T', 1., ib_bdy ) 
    136          CALL lbc_bdy_lnk( h_i(:,:,jl), 'T', 1., ib_bdy ) 
    137          CALL lbc_bdy_lnk( h_s(:,:,jl), 'T', 1., ib_bdy ) 
    138136      ENDDO 
    139       ! retrieve at_i 
    140       at_i(:,:) = 0._wp 
     137      CALL lbc_bdy_lnk( a_i(:,:,:), 'T', 1., jbdy ) 
     138      CALL lbc_bdy_lnk( h_i(:,:,:), 'T', 1., jbdy ) 
     139      CALL lbc_bdy_lnk( h_s(:,:,:), 'T', 1., jbdy ) 
     140 
    141141      DO jl = 1, jpl 
    142          at_i(:,:) = a_i(:,:,jl) + at_i(:,:) 
    143       END DO 
    144  
    145       DO jl = 1, jpl 
    146          DO jb = 1, idx%nblenrim(jgrd) 
    147             ji    = idx%nbi(jb,jgrd) 
    148             jj    = idx%nbj(jb,jgrd) 
     142         DO i_bdy = 1, idx%nblenrim(jgrd) 
     143            ji = idx%nbi(i_bdy,jgrd) 
     144            jj = idx%nbj(i_bdy,jgrd) 
    149145 
    150146            ! condition on ice thickness depends on the ice velocity 
    151147            ! if velocity is outward (strictly), then ice thickness, volume... must be equal to adjacent values 
    152             jpbound = 0   ;   ii = ji   ;   ij = jj 
    153             ! 
    154             IF( u_ice(ji+1,jj  ) < 0. .AND. umask(ji-1,jj  ,1) == 0. ) jpbound = 1; ii = ji+1; ij = jj 
    155             IF( u_ice(ji-1,jj  ) > 0. .AND. umask(ji+1,jj  ,1) == 0. ) jpbound = 1; ii = ji-1; ij = jj 
    156             IF( v_ice(ji  ,jj+1) < 0. .AND. vmask(ji  ,jj-1,1) == 0. ) jpbound = 1; ii = ji  ; ij = jj+1 
    157             IF( v_ice(ji  ,jj-1) > 0. .AND. vmask(ji  ,jj+1,1) == 0. ) jpbound = 1; ii = ji  ; ij = jj-1 
    158             ! 
    159             IF( nn_ice_dta(ib_bdy) == 0 ) jpbound = 0; ii = ji; ij = jj   ! case ice boundaries = initial conditions 
    160             !                                                             !      do not make state variables dependent on velocity 
    161             ! 
    162             rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ii,ij) - 0.01 ) ) ! 0 if no ice 
    163             ! 
    164             ! concentration and thickness 
    165             a_i(ji,jj,jl) = a_i(ii,ij,jl) * rswitch 
    166             h_i(ji,jj,jl) = h_i(ii,ij,jl) * rswitch 
    167             h_s(ji,jj,jl) = h_s(ii,ij,jl) * rswitch 
    168             ! 
    169             ! Ice and snow volumes 
    170             v_i(ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl) 
    171             v_s(ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl) 
    172             ! 
    173             SELECT CASE( jpbound ) 
    174             ! 
    175             CASE( 0 )   ! velocity is inward 
    176                ! 
    177                ! Ice salinity, age, temperature 
    178                s_i (ji,jj,jl)   = rswitch * rn_ice_sal(ib_bdy)  + ( 1.0 - rswitch ) * rn_simin 
    179                oa_i(ji,jj,jl)   = rswitch * rn_ice_age(ib_bdy) * a_i(ji,jj,jl) 
    180                t_su(ji,jj,jl)   = rswitch * rn_ice_tem(ib_bdy)  + ( 1.0 - rswitch ) * rn_ice_tem(ib_bdy) 
     148            jpbound = 0   ;   ib = ji   ;   jb = jj 
     149            ! 
     150            IF( u_ice(ji+1,jj  ) < 0. .AND. umask(ji-1,jj  ,1) == 0. )   jpbound = 1 ; ib = ji+1 ; jb = jj 
     151            IF( u_ice(ji-1,jj  ) > 0. .AND. umask(ji+1,jj  ,1) == 0. )   jpbound = 1 ; ib = ji-1 ; jb = jj 
     152            IF( v_ice(ji  ,jj+1) < 0. .AND. vmask(ji  ,jj-1,1) == 0. )   jpbound = 1 ; ib = ji   ; jb = jj+1 
     153            IF( v_ice(ji  ,jj-1) > 0. .AND. vmask(ji  ,jj+1,1) == 0. )   jpbound = 1 ; ib = ji   ; jb = jj-1 
     154            ! 
     155            IF( nn_ice_dta(jbdy) == 0 )   jpbound = 0 ; ib = ji ; jb = jj   ! case ice boundaries = initial conditions 
     156            !                                                               !      do not make state variables dependent on velocity 
     157            ! 
     158            IF( a_i(ib,jb,jl) > 0._wp ) THEN   ! there is ice at the boundary 
     159               ! 
     160               a_i(ji,jj,jl) = a_i(ib,jb,jl) ! concentration 
     161               h_i(ji,jj,jl) = h_i(ib,jb,jl) ! thickness ice 
     162               h_s(ji,jj,jl) = h_s(ib,jb,jl) ! thickness snw 
     163               ! 
     164               SELECT CASE( jpbound ) 
     165                  ! 
     166               CASE( 0 )   ! velocity is inward 
     167                  ! 
     168                  oa_i(ji,jj,  jl) = rn_ice_age(jbdy) * a_i(ji,jj,jl) ! age 
     169                  a_ip(ji,jj,  jl) = 0._wp                            ! pond concentration 
     170                  v_ip(ji,jj,  jl) = 0._wp                            ! pond volume 
     171                  t_su(ji,jj,  jl) = rn_ice_tem(jbdy)                 ! temperature surface 
     172                  t_s (ji,jj,:,jl) = rn_ice_tem(jbdy)                 ! temperature snw 
     173                  t_i (ji,jj,:,jl) = rn_ice_tem(jbdy)                 ! temperature ice 
     174                  s_i (ji,jj,  jl) = rn_ice_sal(jbdy)                 ! salinity 
     175                  sz_i(ji,jj,:,jl) = rn_ice_sal(jbdy)                 ! salinity profile 
     176                  ! 
     177               CASE( 1 )   ! velocity is outward 
     178                  ! 
     179                  oa_i(ji,jj,  jl) = oa_i(ib,jb,  jl) ! age 
     180                  a_ip(ji,jj,  jl) = a_ip(ib,jb,  jl) ! pond concentration 
     181                  v_ip(ji,jj,  jl) = v_ip(ib,jb,  jl) ! pond volume 
     182                  t_su(ji,jj,  jl) = t_su(ib,jb,  jl) ! temperature surface 
     183                  t_s (ji,jj,:,jl) = t_s (ib,jb,:,jl) ! temperature snw 
     184                  t_i (ji,jj,:,jl) = t_i (ib,jb,:,jl) ! temperature ice 
     185                  s_i (ji,jj,  jl) = s_i (ib,jb,  jl) ! salinity 
     186                  sz_i(ji,jj,:,jl) = sz_i(ib,jb,:,jl) ! salinity profile 
     187                  ! 
     188               END SELECT 
     189               ! 
     190               IF( nn_icesal == 1 ) THEN     ! if constant salinity 
     191                  s_i (ji,jj  ,jl) = rn_icesal 
     192                  sz_i(ji,jj,:,jl) = rn_icesal 
     193               ENDIF 
     194               ! 
     195               ! global fields 
     196               v_i (ji,jj,jl) = h_i(ji,jj,jl) * a_i(ji,jj,jl)                       ! volume ice 
     197               v_s (ji,jj,jl) = h_s(ji,jj,jl) * a_i(ji,jj,jl)                       ! volume snw 
     198               sv_i(ji,jj,jl) = MIN( s_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content 
    181199               DO jk = 1, nlay_s 
    182                   t_s(ji,jj,jk,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rt0 
     200                  e_s(ji,jj,jk,jl) = rhos * ( rcpi * ( rt0 - t_s(ji,jj,jk,jl) ) + rLfus )   ! enthalpy in J/m3 
     201                  e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s           ! enthalpy in J/m2 
     202               END DO                
     203               DO jk = 1, nlay_i 
     204                  ztmelts          = - rTmlt  * sz_i(ji,jj,jk,jl)             ! Melting temperature in C 
     205                  t_i(ji,jj,jk,jl) = MIN( t_i(ji,jj,jk,jl), ztmelts + rt0 )   ! Force t_i to be lower than melting point => likely conservation issue 
     206                  ! 
     207                  e_i(ji,jj,jk,jl) = rhoi * ( rcpi  * ( ztmelts - ( t_i(ji,jj,jk,jl) - rt0 ) )           &   ! enthalpy in J/m3 
     208                     &                      + rLfus * ( 1._wp - ztmelts / ( t_i(ji,jj,jk,jl) - rt0 ) )   & 
     209                     &                      - rcp   *   ztmelts )                   
     210                  e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * v_i(ji,jj,jl) * r1_nlay_i                            ! enthalpy in J/m2 
    183211               END DO 
    184                DO jk = 1, nlay_i 
    185                   t_i (ji,jj,jk,jl) = rswitch * rn_ice_tem(ib_bdy) + ( 1.0 - rswitch ) * rt0  
    186                   sz_i(ji,jj,jk,jl) = rswitch * rn_ice_sal(ib_bdy) + ( 1.0 - rswitch ) * rn_simin 
    187                END DO 
    188                ! 
    189             CASE( 1 )   ! velocity is outward 
    190                ! 
    191                ! Ice salinity, age, temperature 
    192                s_i (ji,jj,jl)   = rswitch * s_i (ii,ij,jl)  + ( 1.0 - rswitch ) * rn_simin 
    193                oa_i(ji,jj,jl)   = rswitch * oa_i(ii,ij,jl) 
    194                t_su(ji,jj,jl)   = rswitch * t_su(ii,ij,jl)  + ( 1.0 - rswitch ) * rt0 
    195                DO jk = 1, nlay_s 
    196                   t_s(ji,jj,jk,jl) = rswitch * t_s(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rt0 
    197                END DO 
    198                DO jk = 1, nlay_i 
    199                   t_i (ji,jj,jk,jl) = rswitch * t_i (ii,ij,jk,jl) + ( 1.0 - rswitch ) * rt0 
    200                   sz_i(ji,jj,jk,jl) = rswitch * sz_i(ii,ij,jk,jl) + ( 1.0 - rswitch ) * rn_simin 
    201                END DO 
    202                ! 
    203             END SELECT 
    204             ! 
    205             IF( nn_icesal == 1 ) THEN     ! constant salinity : overwrite rn_icesal 
    206                s_i (ji,jj  ,jl) = rn_icesal 
    207                sz_i(ji,jj,:,jl) = rn_icesal 
     212               ! 
     213            ELSE   ! no ice at the boundary 
     214               ! 
     215               a_i (ji,jj,  jl) = 0._wp 
     216               h_i (ji,jj,  jl) = 0._wp 
     217               h_s (ji,jj,  jl) = 0._wp 
     218               oa_i(ji,jj,  jl) = 0._wp 
     219               a_ip(ji,jj,  jl) = 0._wp 
     220               v_ip(ji,jj,  jl) = 0._wp 
     221               t_su(ji,jj,  jl) = rt0 
     222               t_s (ji,jj,:,jl) = rt0 
     223               t_i (ji,jj,:,jl) = rt0  
     224                
     225               IF( nn_icesal == 1 ) THEN     ! if constant salinity 
     226                  s_i (ji,jj  ,jl) = rn_icesal 
     227                  sz_i(ji,jj,:,jl) = rn_icesal 
     228               ELSE                          ! if variable salinity 
     229                  s_i (ji,jj,jl)   = rn_simin 
     230                  sz_i(ji,jj,:,jl) = rn_simin 
     231               ENDIF 
     232               ! 
     233               ! global fields 
     234               v_i (ji,jj,  jl) = 0._wp 
     235               v_s (ji,jj,  jl) = 0._wp 
     236               sv_i(ji,jj,  jl) = 0._wp 
     237               e_s (ji,jj,:,jl) = 0._wp 
     238               e_i (ji,jj,:,jl) = 0._wp 
     239 
    208240            ENDIF 
    209             ! 
    210             ! contents 
    211             sv_i(ji,jj,jl)  = MIN( s_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) 
    212             DO jk = 1, nlay_s 
    213                ! Snow energy of melting 
    214                e_s(ji,jj,jk,jl) = rswitch * rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus ) 
    215                ! Multiply by volume, so that heat content in J/m2 
    216                e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s 
    217             END DO 
    218             DO jk = 1, nlay_i 
    219                ztmelts          = - tmut * sz_i(ji,jj,jk,jl) + rt0 !Melting temperature in K                   
    220                ! heat content per unit volume 
    221                e_i(ji,jj,jk,jl) = rswitch * rhoic * & 
    222                   (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) & 
    223                   +   lfus    * ( 1.0 - (ztmelts-rt0) / MIN((t_i(ji,jj,jk,jl)-rt0),-epsi20) ) & 
    224                   - rcp      * ( ztmelts - rt0 ) ) 
    225                ! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2 
    226                e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * a_i(ji,jj,jl) * h_i(ji,jj,jl) * r1_nlay_i 
    227             END DO 
    228             ! 
     241                         
    229242         END DO 
    230243         ! 
    231          CALL lbc_bdy_lnk( a_i(:,:,jl), 'T', 1., ib_bdy ) 
    232          CALL lbc_bdy_lnk( h_i(:,:,jl), 'T', 1., ib_bdy ) 
    233          CALL lbc_bdy_lnk( h_s(:,:,jl), 'T', 1., ib_bdy ) 
    234          CALL lbc_bdy_lnk( v_i(:,:,jl), 'T', 1., ib_bdy ) 
    235          CALL lbc_bdy_lnk( v_s(:,:,jl), 'T', 1., ib_bdy ) 
    236          ! 
    237          CALL lbc_bdy_lnk( sv_i(:,:,jl), 'T', 1., ib_bdy ) 
    238          CALL lbc_bdy_lnk(  s_i(:,:,jl), 'T', 1., ib_bdy ) 
    239          CALL lbc_bdy_lnk( oa_i(:,:,jl), 'T', 1., ib_bdy ) 
    240          CALL lbc_bdy_lnk( t_su(:,:,jl), 'T', 1., ib_bdy ) 
    241          DO jk = 1, nlay_s 
    242             CALL lbc_bdy_lnk(t_s(:,:,jk,jl), 'T', 1., ib_bdy ) 
    243             CALL lbc_bdy_lnk(e_s(:,:,jk,jl), 'T', 1., ib_bdy ) 
    244          END DO 
    245          DO jk = 1, nlay_i 
    246             CALL lbc_bdy_lnk(t_i(:,:,jk,jl), 'T', 1., ib_bdy ) 
    247             CALL lbc_bdy_lnk(e_i(:,:,jk,jl), 'T', 1., ib_bdy ) 
    248          END DO 
    249          ! 
    250       END DO !jl 
    251       ! 
    252       ! --- In case categories are out of bounds, do a remapping --- ! 
    253       !     i.e. inputs have not the same ice thickness distribution  
    254       !          (set by rn_himean) than the regional simulation 
    255       IF( jpl > 1 )   CALL ice_itd_reb( kt ) 
     244      END DO ! jl 
     245 
     246      CALL lbc_bdy_lnk( a_i (:,:,:)  , 'T', 1., jbdy ) 
     247      CALL lbc_bdy_lnk( h_i (:,:,:)  , 'T', 1., jbdy ) 
     248      CALL lbc_bdy_lnk( h_s (:,:,:)  , 'T', 1., jbdy ) 
     249      CALL lbc_bdy_lnk( oa_i(:,:,:)  , 'T', 1., jbdy ) 
     250      CALL lbc_bdy_lnk( a_ip(:,:,:)  , 'T', 1., jbdy ) 
     251      CALL lbc_bdy_lnk( v_ip(:,:,:)  , 'T', 1., jbdy ) 
     252      CALL lbc_bdy_lnk( s_i (:,:,:)  , 'T', 1., jbdy ) 
     253      CALL lbc_bdy_lnk( t_su(:,:,:)  , 'T', 1., jbdy ) 
     254      CALL lbc_bdy_lnk( v_i (:,:,:)  , 'T', 1., jbdy ) 
     255      CALL lbc_bdy_lnk( v_s (:,:,:)  , 'T', 1., jbdy ) 
     256      CALL lbc_bdy_lnk( sv_i(:,:,:)  , 'T', 1., jbdy ) 
     257      CALL lbc_bdy_lnk( t_s (:,:,:,:), 'T', 1., jbdy ) 
     258      CALL lbc_bdy_lnk( e_s (:,:,:,:), 'T', 1., jbdy ) 
     259      CALL lbc_bdy_lnk( t_i (:,:,:,:), 'T', 1., jbdy ) 
     260      CALL lbc_bdy_lnk( e_i (:,:,:,:), 'T', 1., jbdy ) 
    256261      !       
    257262   END SUBROUTINE bdy_ice_frs 
     
    262267      !!                 ***  SUBROUTINE bdy_ice_dyn  *** 
    263268      !!                     
    264       !! ** Purpose : Apply dynamics boundary conditions for sea-ice in the cas of unstructured open boundaries. 
    265       !!              u_ice and v_ice are equal to the value of the adjacent grid point if this latter is not ice free 
    266       !!              if adjacent grid point is ice free, then u_ice and v_ice are equal to ocean velocities 
     269      !! ** Purpose : Apply dynamics boundary conditions for sea-ice. 
    267270      !! 
    268       !! 2013-06 : C. Rousset 
     271      !! ** Method :  if this adjacent grid point is not ice free, then u_ice and v_ice take its value 
     272      !!              if                          is     ice free, then u_ice and v_ice are unchanged by BDY 
     273      !!                                                           they keep values calculated in rheology 
     274      !! 
    269275      !!------------------------------------------------------------------------------ 
    270276      CHARACTER(len=1), INTENT(in)  ::   cd_type   ! nature of velocity grid-points 
    271277      ! 
    272       INTEGER  ::   jb, jgrd           ! dummy loop indices 
    273       INTEGER  ::   ji, jj             ! local scalar 
    274       INTEGER  ::   ib_bdy             ! Loop index 
     278      INTEGER  ::   i_bdy, jgrd      ! dummy loop indices 
     279      INTEGER  ::   ji, jj           ! local scalar 
     280      INTEGER  ::   jbdy             ! BDY set index 
    275281      REAL(wp) ::   zmsk1, zmsk2, zflag 
    276282      !!------------------------------------------------------------------------------ 
    277       ! 
    278       DO ib_bdy=1, nb_bdy 
    279          ! 
    280          SELECT CASE( cn_ice(ib_bdy) ) 
     283      IF( ln_timing )   CALL timing_start('bdy_ice_dyn') 
     284      ! 
     285      DO jbdy=1, nb_bdy 
     286         ! 
     287         SELECT CASE( cn_ice(jbdy) ) 
    281288         ! 
    282289         CASE('none') 
     
    285292         CASE('frs') 
    286293            ! 
    287             IF( nn_ice_dta(ib_bdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions  
    288             !                                              !      do not change ice velocity (it is only computed by rheology) 
     294            IF( nn_ice_dta(jbdy) == 0 ) CYCLE            ! case ice boundaries = initial conditions  
     295            !                                            !      do not change ice velocity (it is only computed by rheology) 
    289296            SELECT CASE ( cd_type ) 
    290297            !      
    291298            CASE ( 'U' )   
    292299               jgrd = 2      ! u velocity 
    293                DO jb = 1, idx_bdy(ib_bdy)%nblenrim(jgrd) 
    294                   ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd) 
    295                   jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd) 
    296                   zflag = idx_bdy(ib_bdy)%flagu(jb,jgrd) 
     300               DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd) 
     301                  ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 
     302                  jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 
     303                  zflag = idx_bdy(jbdy)%flagu(i_bdy,jgrd) 
    297304                  ! 
    298305                  IF ( ABS( zflag ) == 1. ) THEN  ! eastern and western boundaries 
     
    301308                     zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji-1,jj) ) ) ! 0 if no ice 
    302309                     !   
    303                      ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce) 
     310                     ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then do not change u_ice) 
    304311                     u_ice (ji,jj) = u_ice(ji+1,jj) * 0.5_wp * ABS( zflag + 1._wp ) * zmsk1 + & 
    305312                        &            u_ice(ji-1,jj) * 0.5_wp * ABS( zflag - 1._wp ) * zmsk2 + & 
    306                         &            u_oce(ji  ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 
     313                        &            u_ice(ji  ,jj) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 
    307314                  ELSE                             ! everywhere else 
    308                      !u_ice(ji,jj) = u_oce(ji,jj) 
    309315                     u_ice(ji,jj) = 0._wp 
    310316                  ENDIF 
    311                   ! mask ice velocities 
    312                   rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ji,jj) - 0.01_wp ) ) ! 0 if no ice 
    313                   u_ice(ji,jj) = rswitch * u_ice(ji,jj) 
    314317                  ! 
    315318               END DO 
    316                CALL lbc_bdy_lnk( u_ice(:,:), 'U', -1., ib_bdy ) 
     319               CALL lbc_bdy_lnk( u_ice(:,:), 'U', -1., jbdy ) 
    317320               ! 
    318321            CASE ( 'V' ) 
    319322               jgrd = 3      ! v velocity 
    320                DO jb = 1, idx_bdy(ib_bdy)%nblenrim(jgrd) 
    321                   ji    = idx_bdy(ib_bdy)%nbi(jb,jgrd) 
    322                   jj    = idx_bdy(ib_bdy)%nbj(jb,jgrd) 
    323                   zflag = idx_bdy(ib_bdy)%flagv(jb,jgrd) 
     323               DO i_bdy = 1, idx_bdy(jbdy)%nblenrim(jgrd) 
     324                  ji    = idx_bdy(jbdy)%nbi(i_bdy,jgrd) 
     325                  jj    = idx_bdy(jbdy)%nbj(i_bdy,jgrd) 
     326                  zflag = idx_bdy(jbdy)%flagv(i_bdy,jgrd) 
    324327                  ! 
    325328                  IF ( ABS( zflag ) == 1. ) THEN  ! northern and southern boundaries 
     
    328331                     zmsk2 = 1._wp - MAX( 0.0_wp, SIGN ( 1.0_wp , - vt_i(ji,jj-1) ) ) ! 0 if no ice 
    329332                     !   
    330                      ! u_ice = u_ice of the adjacent grid point except if this grid point is ice-free (then u_ice = u_oce) 
     333                     ! v_ice = v_ice of the adjacent grid point except if this grid point is ice-free (then do not change v_ice) 
    331334                     v_ice (ji,jj) = v_ice(ji,jj+1) * 0.5_wp * ABS( zflag + 1._wp ) * zmsk1 + & 
    332335                        &            v_ice(ji,jj-1) * 0.5_wp * ABS( zflag - 1._wp ) * zmsk2 + & 
    333                         &            v_oce(ji,jj  ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 
     336                        &            v_ice(ji,jj  ) * ( 1._wp - MIN( 1._wp, zmsk1 + zmsk2 ) ) 
    334337                  ELSE                             ! everywhere else 
    335                      !v_ice(ji,jj) = v_oce(ji,jj) 
    336338                     v_ice(ji,jj) = 0._wp 
    337339                  ENDIF 
    338                   ! mask ice velocities 
    339                   rswitch = MAX( 0.0_wp , SIGN ( 1.0_wp , at_i(ji,jj) - 0.01 ) ) ! 0 if no ice 
    340                   v_ice(ji,jj) = rswitch * v_ice(ji,jj) 
    341340                  ! 
    342341               END DO 
    343                CALL lbc_bdy_lnk( v_ice(:,:), 'V', -1., ib_bdy ) 
     342               CALL lbc_bdy_lnk( v_ice(:,:), 'V', -1., jbdy ) 
    344343               ! 
    345344            END SELECT 
     
    350349         ! 
    351350      END DO 
     351      ! 
     352      IF( ln_timing )   CALL timing_stop('bdy_ice_dyn') 
    352353      ! 
    353354    END SUBROUTINE bdy_ice_dyn 
     
    359360CONTAINS 
    360361   SUBROUTINE bdy_ice( kt )      ! Empty routine 
     362      IMPLICIT NONE 
     363      INTEGER, INTENT( in ) :: kt 
    361364      WRITE(*,*) 'bdy_ice: You should not have seen this print! error?', kt 
    362365   END SUBROUTINE bdy_ice 
Note: See TracChangeset for help on using the changeset viewer.