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 2528 for trunk/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90 – NEMO

Ignore:
Timestamp:
2010-12-27T18:33:53+01:00 (13 years ago)
Author:
rblod
Message:

Update NEMOGCM from branch nemo_v3_3_beta

File:
1 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

    r2477 r2528  
    44   !!           computation of the flux at the sea ice/ocean interface 
    55   !!====================================================================== 
    6    !! History :   -   !  2006-07 (M. Vancoppelle)  LIM3 original code 
    7    !!            3.0  !  2008-03 (C. Tallandier)  surface module 
    8    !!             -   !  2008-04 (C. Tallandier)  split in 2 + new ice-ocean coupling 
     6   !! History :   -   ! 2006-07 (M. Vancoppelle)  LIM3 original code 
     7   !!            3.0  ! 2008-03 (C. Tallandier)  surface module 
     8   !!             -   ! 2008-04 (C. Tallandier)  split in 2 + new ice-ocean coupling 
     9   !!            3.3  ! 2010-05 (G. Madec) decrease ocean & ice reference salinities in the Baltic sea 
     10   !!                 !                  + simplification of the ice-ocean stress calculation 
    911   !!---------------------------------------------------------------------- 
    1012#if defined key_lim3 
     
    1214   !!   'key_lim3'                                    LIM 3.0 sea-ice model 
    1315   !!---------------------------------------------------------------------- 
    14    !!   lim_sbc  : flux at the ice / ocean interface 
    15    !!---------------------------------------------------------------------- 
    16    USE oce 
     16   !!   lim_sbc_flx  : updates mass, heat and salt fluxes at the ocean surface 
     17   !!   lim_sbc_tau  : update i- and j-stresses, and its modulus at the ocean surface 
     18   !!---------------------------------------------------------------------- 
    1719   USE par_oce          ! ocean parameters 
    1820   USE par_ice          ! ice parameters 
     
    2123   USE sbc_oce          ! Surface boundary condition: ocean fields 
    2224   USE phycst           ! physical constants 
     25   USE albedo           ! albedo parameters 
    2326   USE ice              ! LIM sea-ice variables 
    24  
    2527   USE lbclnk           ! ocean lateral boundary condition 
    2628   USE in_out_manager   ! I/O manager 
    27    USE albedo           ! albedo parameters 
    2829   USE prtctl           ! Print control 
    2930 
     
    3435   PUBLIC   lim_sbc_tau   ! called by sbc_ice_lim 
    3536 
    36    REAL(wp)  ::   epsi16 = 1.e-16  ! constant values 
    37    REAL(wp)  ::   rzero  = 0.e0     
    38    REAL(wp)  ::   rone   = 1.e0 
    39  
    40    REAL(wp), DIMENSION(jpi,jpj) ::   utau_oce, vtau_oce   !: air-ocean surface i- & j-stress              [N/m2] 
    41    REAL(wp), DIMENSION(jpi,jpj) ::   tmod_io              !: modulus of the ice-ocean relative velocity   [m/s] 
    42    REAL(wp), DIMENSION(jpi,jpj) ::   ssu_mb  , ssv_mb     !: before mean ocean surface currents           [m/s] 
     37   REAL(wp)  ::   r1_rdtice            ! = 1. / rdt_ice  
     38   REAL(wp)  ::   epsi16 = 1.e-16_wp   ! constant values 
     39   REAL(wp)  ::   rzero  = 0._wp     
     40   REAL(wp)  ::   rone   = 1._wp 
     41 
     42   REAL(wp), DIMENSION(jpi,jpj) ::   utau_oce, vtau_oce   ! air-ocean surface i- & j-stress              [N/m2] 
     43   REAL(wp), DIMENSION(jpi,jpj) ::   tmod_io              ! modulus of the ice-ocean relative velocity   [m/s] 
     44 
     45   REAL(wp), DIMENSION(jpi,jpj) ::   soce_0, sice_0   ! constant SSS and ice salinity used in levitating sea-ice case 
    4346 
    4447   !! * Substitutions 
    4548#  include "vectopt_loop_substitute.h90" 
    4649   !!---------------------------------------------------------------------- 
    47    !! NEMO/LIM  3.2 ,  UCL-LOCEAN-IPSL (2009)  
     50   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4851   !! $Id$ 
    49    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    50    !!---------------------------------------------------------------------- 
    51  
     52   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     53   !!---------------------------------------------------------------------- 
    5254CONTAINS 
    53  
    54    SUBROUTINE lim_sbc_tau( kt, kcpl ) 
    55       !!------------------------------------------------------------------- 
    56       !!                ***  ROUTINE lim_sbc_tau *** 
    57       !!   
    58       !! ** Purpose : Update the ocean surface stresses due to the ice 
    59       !!          
    60       !! ** Action  : - compute the ice-ocean stress depending on kcpl: 
    61       !!              fluxes at the ice-ocean interface. 
    62       !!     Case 0  :  old LIM-3 way, computed at ice time-step only 
    63       !!     Case 1  :  at each ocean time step the stresses are computed 
    64       !!                using the current ocean velocity (now) 
    65       !!     Case 2  :  combination of half case 0 + half case 1 
    66       !!      
    67       !! ** Outputs : - utau   : sea surface i-stress (ocean referential) 
    68       !!              - vtau   : sea surface j-stress (ocean referential) 
    69       !! 
    70       !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90. 
    71       !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
    72       !!--------------------------------------------------------------------- 
    73       INTEGER ::   kt     ! number of ocean iteration 
    74       INTEGER ::   kcpl   ! ice-ocean coupling indicator: =0  LIM-3 old case 
    75       !                   !                               =1  stresses computed using now ocean velocity 
    76       !                   !                               =2  combination of 0 and 1 cases 
    77       !! 
    78       INTEGER  ::   ji, jj   ! dummy loop indices 
    79       REAL(wp) ::   zfrldu, zat_u, zu_ico, zutaui, zu_u, zu_ij, zu_im1j   ! temporary scalar 
    80       REAL(wp) ::   zfrldv, zat_v, zv_ico, zvtaui, zv_v, zv_ij, zv_ijm1   !    -         - 
    81       REAL(wp) ::   zsang, zztmp                                          !    -         - 
    82       REAL(wp), DIMENSION(jpi,jpj) ::   ztio_u, ztio_v   ! ocean stress below sea-ice 
    83 #if defined key_coupled     
    84       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb     ! albedo of ice under overcast sky 
    85       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalbp    ! albedo of ice under clear sky 
    86 #endif 
    87      !!--------------------------------------------------------------------- 
    88  
    89       IF( kt == nit000 ) THEN 
    90          IF(lwp) WRITE(numout,*) 
    91          IF(lwp) WRITE(numout,*) 'lim_sbc_tau : LIM 3.0 sea-ice - surface ocean momentum fluxes' 
    92          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
    93       ENDIF 
    94  
    95       IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    96 !CDIR NOVERRCHK 
    97          DO jj = 2, jpjm1                         !* modulus of the ice-ocean velocity 
    98 !CDIR NOVERRCHK 
    99             DO ji = fs_2, fs_jpim1 
    100                zu_ij   = u_ice(ji  ,jj) - ssu_m(ji  ,jj)               ! (i  ,j) 
    101                zu_im1j = u_ice(ji-1,jj) - ssu_m(ji-1,jj)               ! (i-1,j) 
    102                zv_ij   = v_ice(ji,jj  ) - ssv_m(ji,jj  )               ! (i,j  ) 
    103                zv_ijm1 = v_ice(ji,jj-1) - ssv_m(ji,jj-1)               ! (i,j-1) 
    104                zztmp =  0.5 * ( zu_ij * zu_ij + zu_im1j * zu_im1j   &  ! mean of the square values instead 
    105                   &           + zv_ij * zv_ij + zv_ijm1 * zv_ijm1 )    ! of the square of the mean values... 
    106                ! WARNING, here we make an approximation for case 1 and 2: taum is not computed at each time 
    107                ! step but only every nn_fsbc time step. This avoid mutiple avarage to pass from T -> U,V grids 
    108                ! and next from U,V grids to T grid. Taum is used in tke, which should not be too sensitive to 
    109                ! this approximaton... 
    110                taum(ji,jj) = ( 1. - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zztmp 
    111                tmod_io(ji,jj) = SQRT( zztmp )  
    112             END DO 
    113          END DO 
    114          CALL lbc_lnk( taum, 'T', 1. )   ;   CALL lbc_lnk( tmod_io, 'T', 1. ) 
    115       ENDIF 
    116  
    117       SELECT CASE( kcpl ) 
    118          !                                        !--------------------------------! 
    119       CASE( 0 )                                   !  LIM 3 old stress computation  !  (at ice timestep only) 
    120          !                                        !--------------------------------!  
    121          !                                             !* ice stress over ocean with a ice-ocean rotation angle 
    122          DO jj = 1, jpjm1 
    123             zsang  = SIGN( 1.e0, gphif(1,jj) ) * sangvg         ! change the sinus angle sign in the south hemisphere 
    124             DO ji = 1, fs_jpim1 
    125                zu_u = u_ice(ji,jj) - u_oce(ji,jj)               ! ice velocity relative to the ocean 
    126                zv_v = v_ice(ji,jj) - v_oce(ji,jj) 
    127                !                                                ! quadratic drag formulation with rotation 
    128 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1) 
    129                zutaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_u - zsang * zv_v )  
    130                zvtaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_v + zsang * zu_u )  
    131                !                                                ! bound for too large stress values 
    132                ! IMPORTANT: the exponential below prevents numerical oscillations in the ice-ocean stress 
    133                ! This is not physically based. A cleaner solution is offer in CASE kcpl=2 
    134                ztio_u(ji,jj) = zutaui * EXP( - ( tmod_io(ji,jj) + tmod_io(ji+1,jj) )  ) 
    135                ztio_v(ji,jj) = zvtaui * EXP( - ( tmod_io(ji,jj) + tmod_io(ji,jj+1) )  )  
    136             END DO 
    137          END DO 
    138          DO jj = 2, jpjm1                                       ! stress at the surface of the ocean 
    139             DO ji = fs_2, fs_jpim1   ! vertor opt. 
    140                zfrldu = 0.5 * ( ato_i(ji,jj) + ato_i(ji+1,jj  ) )   ! open-ocean fraction at U- & V-points (from T-point values) 
    141                zfrldv = 0.5 * ( ato_i(ji,jj) + ato_i(ji  ,jj+1) ) 
    142                !                                                    ! update surface ocean stress 
    143                utau(ji,jj) = zfrldu * utau(ji,jj) + ( 1. - zfrldu ) * ztio_u(ji,jj) 
    144                vtau(ji,jj) = zfrldv * vtau(ji,jj) + ( 1. - zfrldv ) * ztio_v(ji,jj) 
    145             END DO 
    146          END DO 
    147          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
    148  
    149          ! 
    150          !                                        !--------------------------------! 
    151       CASE( 1 )                                   !  Use the now velocity          !  (at each ocean timestep) 
    152          !                                        !--------------------------------!  
    153          IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    154             utau_oce(:,:) = utau(:,:)                 !* save the air-ocean stresses at ice time-step 
    155             vtau_oce(:,:) = vtau(:,:) 
    156          ENDIF 
    157          ! 
    158         DO jj = 2, jpjm1                              !* ice stress over ocean with a ice-ocean rotation angle 
    159             zsang  = SIGN(1.e0, gphif(1,jj-1) ) * sangvg        ! change the sinus angle sign in the south hemisphere 
    160             DO ji = fs_2, fs_jpim1 
    161                zat_u  = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5   ! ice area at u and V-points 
    162                zat_v  = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5 
    163                !                                                ! (u,v) ice-ocean velocity at (U,V)-point, resp. 
    164                zu_u   = u_ice(ji,jj) - ub(ji,jj,1) 
    165                zv_v   = v_ice(ji,jj) - vb(ji,jj,1) 
    166                !                                                ! quadratic drag formulation with rotation 
    167 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1) 
    168                zutaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_u - zsang * zv_v )  
    169                zvtaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_v + zsang * zu_u )  
    170                !                                                   ! stress at the ocean surface 
    171                utau(ji,jj) = ( 1.- zat_u ) * utau_oce(ji,jj) + zat_u * zutaui 
    172                vtau(ji,jj) = ( 1.- zat_v ) * vtau_oce(ji,jj) + zat_v * zvtaui 
    173             END DO 
    174          END DO 
    175          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
    176  
    177          ! 
    178          !                                        !--------------------------------! 
    179       CASE( 2 )                                   !  mixed 0 and 2 cases           !  (at each ocean timestep) 
    180          !                                        !--------------------------------!  
    181          IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    182             utau_oce(:,:) = utau (:,:)               !* save the air-ocean stresses at ice time-step 
    183             vtau_oce(:,:) = vtau (:,:) 
    184             ssu_mb  (:,:) = ssu_m(:,:)               !* save the ice-ocean velocity at ice time-step 
    185             ssv_mb  (:,:) = ssv_m(:,:) 
    186          ENDIF 
    187          DO jj = 2, jpjm1                            !* ice stress over ocean with a ice-ocean rotation angle 
    188             zsang  = SIGN(1.e0, gphif(1,jj-1) ) * sangvg 
    189             DO ji = fs_2, fs_jpim1 
    190                zat_u = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5     ! ice area at u and V-points 
    191                zat_v = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5  
    192                ! 
    193                zu_ico = u_ice(ji,jj) - 0.5 * ( ub(ji,jj,1) + ssu_mb(ji,jj) )   ! ice-oce velocity using ub and ssu_mb 
    194                zv_ico = v_ice(ji,jj) - 0.5 * ( vb(ji,jj,1) + ssv_mb(ji,jj) ) 
    195                !                                        ! quadratic drag formulation with rotation 
    196 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1) 
    197                zutaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_ico - zsang * zv_ico ) 
    198                zvtaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_ico + zsang * zu_ico ) 
    199                ! 
    200                utau(ji,jj) = ( 1.-zat_u ) * utau_oce(ji,jj) + zat_u * zutaui    ! stress at the ocean surface 
    201                vtau(ji,jj) = ( 1.-zat_v ) * vtau_oce(ji,jj) + zat_v * zvtaui 
    202             END DO 
    203          END DO 
    204          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
    205          ! 
    206       END SELECT 
    207  
    208       IF(ln_ctl)   CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
    209          &                       tab2d_2=vtau, clinfo2=' vtau    : '        , mask2=vmask ) 
    210       !   
    211    END SUBROUTINE lim_sbc_tau 
    212  
    21355 
    21456   SUBROUTINE lim_sbc_flx( kt ) 
     
    23476      !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
    23577      !!--------------------------------------------------------------------- 
    236       INTEGER ::   kt    ! number of iteration 
     78      INTEGER, INTENT(in) ::   kt    ! number of iteration 
    23779      !! 
    23880      INTEGER  ::   ji, jj           ! dummy loop indices 
     
    25395         IF(lwp) WRITE(numout,*) 'lim_sbc_flx : LIM 3.0 sea-ice - heat salt and mass ocean surface fluxes' 
    25496         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
     97         ! 
     98         r1_rdtice = 1. / rdt_ice 
     99         ! 
     100         soce_0(:,:) = soce 
     101         sice_0(:,:) = sice 
     102         ! 
     103         IF( cp_cfg == "orca" ) THEN           ! decrease ocean & ice reference salinities in the Baltic sea  
     104            WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.   & 
     105               &   54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp         )  
     106               soce_0(:,:) = 4._wp 
     107               sice_0(:,:) = 2._wp 
     108            END WHERE 
     109         ENDIF 
     110         ! 
    255111      ENDIF 
    256112 
     
    297153               ! fscmbq and ffltbif are obsolete 
    298154               !              &           + iflt * ffltbif(ji,jj) !!! only if one category is used 
    299                &           + ifral   * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) / rdt_ice   & 
    300                &           + ifrdv   * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) / rdt_ice                     & 
     155               &           + ifral   * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice   & 
     156               &           + ifrdv   * ( qfvbq(ji,jj) + qdtcn(ji,jj) )                   * r1_rdtice   & 
    301157               &           + fhmec(ji,jj)     & ! new contribution due to snow melt in ridging!! 
    302158               &           + fheat_rpo(ji,jj) & ! contribution from ridge formation 
     
    339195               WRITE(numout,*) ' qcmif     : ', qcmif(jiindx,jjindx) 
    340196               WRITE(numout,*) ' qldif     : ', qldif(jiindx,jjindx) 
    341                WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindx,jjindx) / rdt_ice 
    342                WRITE(numout,*) ' qldif / dt: ', qldif(jiindx,jjindx) / rdt_ice 
     197               WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindx,jjindx) * r1_rdtice 
     198               WRITE(numout,*) ' qldif / dt: ', qldif(jiindx,jjindx) * r1_rdtice 
    343199               WRITE(numout,*) ' ifrdv     : ', ifrdv 
    344200               WRITE(numout,*) ' qfvbq     : ', qfvbq(jiindx,jjindx) 
    345201               WRITE(numout,*) ' qdtcn     : ', qdtcn(jiindx,jjindx) 
    346                WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) / rdt_ice 
    347                WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) / rdt_ice 
     202               WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) * r1_rdtice 
     203               WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) * r1_rdtice 
    348204               WRITE(numout,*) ' ' 
    349205               WRITE(numout,*) ' fdtcn     : ', fdtcn(jiindx,jjindx) 
     
    374230 
    375231            !  computing freshwater exchanges at the ice/ocean interface 
    376             zpme = - emp(ji,jj)     * ( 1.0 - at_i(ji,jj) )  &   !  evaporation over oceanic fraction 
    377                &   + tprecip(ji,jj) *         at_i(ji,jj)    &   !  total precipitation 
    378                ! old fashioned way                
    379                !              &   - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) )  &   !  remov. snow precip over ice 
    380                &   - sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) )  &   !  remov. snow precip over ice 
    381                &   - rdmsnif(ji,jj) / rdt_ice                &   !  freshwaterflux due to snow melting  
    382                ! new contribution from snow falling when ridging 
    383                &   + fmmec(ji,jj) 
     232            zpme = - emp(ji,jj)     * ( 1.0 - at_i(ji,jj)          )  &   ! evaporation over oceanic fraction 
     233               &   + tprecip(ji,jj) *         at_i(ji,jj)             &   ! all precipitation reach the ocean 
     234               &   - sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) )  &   ! except solid precip intercepted by sea-ice 
     235               &   - rdmsnif(ji,jj) * r1_rdtice                       &   ! freshwaterflux due to snow melting  
     236               &   + fmmec(ji,jj)                                         ! snow falling when ridging 
     237 
    384238 
    385239            !  computing salt exchanges at the ice/ocean interface 
    386240            !  sice should be the same as computed with the ice model 
    387             zfons =  ( soce - sice ) * ( rdmicif(ji,jj) / rdt_ice )  
     241            zfons =  ( soce_0(ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice  
    388242            ! SOCE 
    389             zfons =  ( sss_m(ji,jj) - sice ) * ( rdmicif(ji,jj) / rdt_ice )  
     243            zfons =  ( sss_m (ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice 
    390244 
    391245            !CT useless            !  salt flux for constant salinity 
     
    445299   END SUBROUTINE lim_sbc_flx 
    446300 
     301 
     302   SUBROUTINE lim_sbc_tau( kt , pu_oce, pv_oce ) 
     303      !!------------------------------------------------------------------- 
     304      !!                ***  ROUTINE lim_sbc_tau *** 
     305      !!   
     306      !! ** Purpose : Update the ocean surface stresses due to the ice 
     307      !!          
     308      !! ** Action  : * at each ice time step (every nn_fsbc time step): 
     309      !!                - compute the modulus of ice-ocean relative velocity  
     310      !!                  (*rho*Cd) at T-point (C-grid) or I-point (B-grid) 
     311      !!                      tmod_io = rhoco * | U_ice-U_oce | 
     312      !!                - update the modulus of stress at ocean surface 
     313      !!                      taum = frld * taum + (1-frld) * tmod_io * | U_ice-U_oce | 
     314      !!              * at each ocean time step (every kt):  
     315      !!                  compute linearized ice-ocean stresses as 
     316      !!                      Utau = tmod_io * | U_ice - pU_oce | 
     317      !!                using instantaneous current ocean velocity (usually before) 
     318      !! 
     319      !!    NB: - ice-ocean rotation angle no more allowed 
     320      !!        - here we make an approximation: taum is only computed every ice time step 
     321      !!          This avoids mutiple average to pass from T -> U,V grids and next from U,V grids  
     322      !!          to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... 
     323      !! 
     324      !! ** Outputs : - utau, vtau   : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes 
     325      !!              - taum         : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes 
     326      !!--------------------------------------------------------------------- 
     327      INTEGER ,                     INTENT(in) ::   kt               ! ocean time-step index 
     328      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pu_oce, pv_oce   ! surface ocean currents 
     329      !! 
     330      INTEGER  ::   ji, jj   ! dummy loop indices 
     331      REAL(wp) ::   zat_u, zutau_ice, zu_t, zmodt   ! local scalar 
     332      REAL(wp) ::   zat_v, zvtau_ice, zv_t          !   -      - 
     333     !!--------------------------------------------------------------------- 
     334 
     335      IF( kt == nit000 ) THEN 
     336         IF(lwp) WRITE(numout,*) 
     337         IF(lwp) WRITE(numout,*) 'lim_sbc_tau : LIM-3 sea-ice - surface ocean momentum fluxes' 
     338         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
     339      ENDIF 
     340 
     341      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==!   (i.e. surface module time-step) 
     342!CDIR NOVERRCHK 
     343         DO jj = 2, jpjm1                             !* update the modulus of stress at ocean surface (T-point) 
     344!CDIR NOVERRCHK 
     345            DO ji = fs_2, fs_jpim1 
     346               !                                               ! 2*(U_ice-U_oce) at T-point 
     347               zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj)    
     348               zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1)  
     349               !                                              ! |U_ice-U_oce|^2 
     350               zmodt =  0.25_wp * (  zu_t * zu_t + zv_t * zv_t  ) 
     351               !                                               ! update the ocean stress modulus 
     352               taum(ji,jj) = ( 1._wp - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zmodt 
     353               tmod_io(ji,jj) = rhoco * SQRT( zmodt )          ! rhoco * |U_ice-U_oce| at T-point 
     354            END DO 
     355         END DO 
     356         CALL lbc_lnk( taum, 'T', 1. )   ;   CALL lbc_lnk( tmod_io, 'T', 1. ) 
     357         ! 
     358         utau_oce(:,:) = utau(:,:)                    !* save the air-ocean stresses at ice time-step 
     359         vtau_oce(:,:) = vtau(:,:) 
     360         ! 
     361      ENDIF 
     362         ! 
     363         !                                     !==  every ocean time-step  ==! 
     364         ! 
     365      DO jj = 2, jpjm1                                !* update the stress WITHOUT a ice-ocean rotation angle 
     366         DO ji = fs_2, fs_jpim1   ! Vect. Opt. 
     367            zat_u  = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5_wp   ! ice area at u and V-points 
     368            zat_v  = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5_wp 
     369            !                                                   ! linearized quadratic drag formulation 
     370            zutau_ice   = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) ) 
     371            zvtau_ice   = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) ) 
     372            !                                                   ! stresses at the ocean surface 
     373            utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_ice 
     374            vtau(ji,jj) = ( 1._wp - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice 
     375         END DO 
     376      END DO 
     377      CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
     378      ! 
     379      IF(ln_ctl)   CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
     380         &                       tab2d_2=vtau, clinfo2=' vtau    : '        , mask2=vmask ) 
     381      !   
     382   END SUBROUTINE lim_sbc_tau 
     383 
    447384#else 
    448385   !!---------------------------------------------------------------------- 
Note: See TracChangeset for help on using the changeset viewer.