Changeset 2370


Ignore:
Timestamp:
2010-11-10T08:48:54+01:00 (10 years ago)
Author:
gm
Message:

v3.3beta: ice-ocean stress at kt with VP & EVP (LIM-2 and -3)

Location:
branches/nemo_v3_3_beta/NEMOGCM
Files:
29 edited

Legend:

Unmodified
Added
Removed
  • branches/nemo_v3_3_beta/NEMOGCM/ARCH/arch-ifort_osx.fcm

    r2281 r2370  
    1515 
    1616 
    17 %NCDF_INC            -I/usr/local/pub/netcdf/4.0.1-ifort/include  
    18 %NCDF_LIB            -L /usr/local/pub/netcdf/4.0.1-ifort/lib -lnetcdf    
     17%NCDF_INC            -I/usr/local/include  
     18%NCDF_LIB            -L /usr/local/lib -lnetcdf    
    1919%FC                  mpif90 
    2020%FCFLAGS         -r8 -O3  -traceback  
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/GYRE/EXP00/namelist

    r2364 r2370  
    131131                           !  =1 use observed ice-cover      , 
    132132                           !  =2 ice-model used                             ("key_lim3" or "key_lim2) 
    133    nn_ico_cpl  = 0         !  ice-ocean coupling : =0 each nn_fsbc  
    134                            !                       =1 stresses recomputed each ocean time step ("key_lim3" only) 
    135                            !                       =2 combination of 0 and 1 cases             ("key_lim3" only) 
    136133   ln_dm2dc    = .false.   !  daily mean to diurnal cycle short wave (qsr) 
    137134   ln_rnf      = .false.   !  runoffs (T => fill namsbc_rnf) 
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/GYRE_LOBSTER/EXP00/namelist

    r2364 r2370  
    131131                           !  =1 use observed ice-cover      , 
    132132                           !  =2 ice-model used                             ("key_lim3" or "key_lim2) 
    133    nn_ico_cpl  = 0         !  ice-ocean coupling : =0 each nn_fsbc  
    134                            !                       =1 stresses recomputed each ocean time step ("key_lim3" only) 
    135                            !                       =2 combination of 0 and 1 cases             ("key_lim3" only) 
    136133   ln_dm2dc    = .false.   !  daily mean to diurnal cycle short wave (qsr) 
    137134   ln_rnf      = .false.   !  runoffs (T => fill namsbc_rnf) 
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/1_namelist

    r2367 r2370  
    131131                           !  =1 use observed ice-cover      , 
    132132                           !  =2 ice-model used                             ("key_lim3" or "key_lim2) 
    133    nn_ico_cpl  = 0         !  ice-ocean coupling : =0 each nn_fsbc  
    134                            !                       =1 stresses recomputed each ocean time step ("key_lim3" only) 
    135                            !                       =2 combination of 0 and 1 cases             ("key_lim3" only) 
    136133   ln_dm2dc    = .false.   !  daily mean to diurnal cycle short wave (qsr) 
    137134   ln_rnf      = .false.   !  runoffs (T => fill namsbc_rnf) 
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/EXP00/namelist

    r2367 r2370  
    131131                           !  =1 use observed ice-cover      , 
    132132                           !  =2 ice-model used                             ("key_lim3" or "key_lim2) 
    133    nn_ico_cpl  = 0         !  ice-ocean coupling : =0 each nn_fsbc  
    134                            !                       =1 stresses recomputed each ocean time step ("key_lim3" only) 
    135                            !                       =2 combination of 0 and 1 cases             ("key_lim3" only) 
    136133   ln_dm2dc    = .false.   !  daily mean to diurnal cycle short wave (qsr) 
    137134   ln_rnf      = .true.    !  runoffs (T => fill namsbc_rnf) 
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM/cpp_ORCA2_LIM.fcm

    r2369 r2370  
    1  bld::tool::fppkeys key_trabbl key_vectopt_loop key_orca_r2 key_lim2 key_dynspg_flt key_diaeiv key_ldfslp key_traldf_c2d key_traldf_eiv key_dynldf_c3d key_dtatem key_dtasal key_tradmp key_zdftke key_zdfddm key_iomput key_nproci=1 key_nprocj=1 
     1 bld::tool::fppkeys  key_trabbl key_vectopt_loop key_orca_r2 key_lim2 key_dynspg_flt key_diaeiv key_ldfslp key_traldf_c2d key_traldf_eiv key_dynldf_c3d key_dtatem key_dtasal key_tradmp key_zdftke key_zdfddm key_iomput key_nproci=1 key_nprocj=1 
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_LIM_PISCES/EXP00/namelist

    r2367 r2370  
    131131                           !  =1 use observed ice-cover      , 
    132132                           !  =2 ice-model used                             ("key_lim3" or "key_lim2) 
    133    nn_ico_cpl  = 0         !  ice-ocean coupling : =0 each nn_fsbc  
    134                            !                       =1 stresses recomputed each ocean time step ("key_lim3" only) 
    135                            !                       =2 combination of 0 and 1 cases             ("key_lim3" only) 
    136133   ln_dm2dc    = .false.   !  daily mean to diurnal cycle short wave (qsr) 
    137134   ln_rnf      = .true.    !  runoffs (T => fill namsbc_rnf) 
  • branches/nemo_v3_3_beta/NEMOGCM/CONFIG/ORCA2_OFF_PISCES/EXP00/namelist

    r2364 r2370  
    131131                           !  =1 use observed ice-cover      , 
    132132                           !  =2 ice-model used                             ("key_lim3" or "key_lim2) 
    133    nn_ico_cpl  = 0         !  ice-ocean coupling : =0 each nn_fsbc  
    134                            !                       =1 stresses recomputed each ocean time step ("key_lim3" only) 
    135                            !                       =2 combination of 0 and 1 cases             ("key_lim3" only) 
    136133   ln_dm2dc    = .false.   !  daily mean to diurnal cycle short wave (qsr) 
    137134   ln_rnf      = .true.    !  runoffs (T => fill namsbc_rnf) 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/ice_2.F90

    r2319 r2370  
    44   !! Sea Ice physics:  diagnostics variables of ice defined in memory 
    55   !!===================================================================== 
    6    !! History :  2.0  !  03-08  (C. Ethe)  F90: Free form and module 
     6   !! History :  2.0  !  2003-08  (C. Ethe)  F90: Free form and module 
    77   !!            3.3  !  2009-05  (G.Garric) addition of the lim2_evp cas 
    88   !!---------------------------------------------------------------------- 
     
    1111   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1212   !!---------------------------------------------------------------------- 
    13    !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 
    14    !! $Id$ 
    15    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    16    !!---------------------------------------------------------------------- 
    17    !! * Modules used 
    1813   USE par_ice_2          ! LIM sea-ice parameters 
    1914 
     
    6156   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ust2s         !: friction velocity 
    6257 
    63    !!* diagnostic quantities 
     58   !!* Ice Rheology 
    6459# if defined key_lim2_vp 
    6560   !                                                       !!* VP rheology * 
     61   LOGICAL , PUBLIC ::   lk_lim2_vp = .TRUE.               !: Visco-Plactic reology flag  
     62   ! 
    6663   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hsnm , hicm   !: mean snow and ice thicknesses 
    67    CHARACTER(len=1), PUBLIC             ::   cl_grid = 'B' !: type of grid used in ice dynamics, 'C' or 'B' 
    6864   ! 
    6965# else 
    7066   !                                                       !!* EVP rheology * 
    71    CHARACTER(len=1), PUBLIC             ::   cl_grid = 'C' !: type of grid used in ice dynamics, 'C' or 'B' 
     67   LOGICAL , PUBLIC::   lk_lim2_vp = .FALSE.               !: Visco-Plactic reology flag  
     68   ! 
    7269   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   stress1_i     !: first stress tensor element        
    7370   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   stress2_i     !: second stress tensor element 
     
    7673   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   divu_i        !: Divergence of the velocity field [s-1] -> limrhg.F90 
    7774   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   shear_i       !: Shear of the velocity field [s-1] -> limrhg.F90 
     75   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   at_i          !: ice fraction 
    7876   ! 
    79    REAL(wp), PUBLIC, DIMENSION(jpi,jpj)          :: at_i          !: 
    8077   REAL(wp), PUBLIC, DIMENSION(:,:)    , POINTER :: vt_s ,vt_i    !: mean snow and ice thicknesses 
    8178   REAL(wp), PUBLIC, DIMENSION(jpi,jpj), TARGET  :: hsnm , hicm   !: target vt_s,vt_i pointers  
    82    !   
    8379#endif 
    8480 
     
    131127#endif 
    132128 
     129   !!---------------------------------------------------------------------- 
     130   !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 
     131   !! $Id$ 
     132   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    133133   !!====================================================================== 
    134134END MODULE ice_2 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limdyn_2.F90

    r2332 r2370  
    44   !!   Sea-Ice dynamics :   
    55   !!====================================================================== 
    6    !! History :   1.0  !  01-04  (LIM)  Original code 
    7    !!             2.0  !  02-08  (C. Ethe, G. Madec)  F90, mpp 
    8    !!             2.0  !  03-08  (C. Ethe) add lim_dyn_init 
    9    !!             2.0  !  06-07  (G. Madec)  Surface module 
     6   !! History :  1.0  ! 2001-04  (LIM)  Original code 
     7   !!            2.0  ! 2002-08  (C. Ethe, G. Madec)  F90, mpp 
     8   !!            2.0  ! 2003-08  (C. Ethe) add lim_dyn_init 
     9   !!            2.0  ! 2006-07  (G. Madec)  Surface module 
     10   !!            3.3  ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 
    1011   !!--------------------------------------------------------------------- 
    1112#if defined key_lim2 
     
    1617   !!    lim_dyn_init_2 : initialization and namelist read 
    1718   !!---------------------------------------------------------------------- 
    18    USE dom_oce        ! ocean space and time domain 
    19    USE sbc_oce        ! 
    20    USE phycst         ! 
    21    USE ice_2          ! 
    22    USE dom_ice_2      ! 
    23    USE limistate_2    ! 
    24 #if defined key_lim2_vp 
    25    USE limrhg_2       ! ice rheology 
    26 #else 
    27    USE limrhg         ! LIM : EVP ice rheology 
    28 #endif 
    29    USE lbclnk         ! 
    30    USE lib_mpp        ! 
    31    USE in_out_manager ! I/O manager 
    32    USE prtctl         ! Print control 
     19   USE dom_oce          ! ocean space and time domain 
     20   USE sbc_oce          ! ocean surface boundary condition 
     21   USE phycst           ! physical constant 
     22   USE ice_2            ! LIM-2: ice variables 
     23   USE sbc_ice          ! Surface boundary condition: sea-ice fields 
     24   USE dom_ice_2        ! LIM-2: ice domain 
     25   USE limistate_2      ! LIM-2: initial state 
     26   USE limrhg_2         ! LIM-2: VP  ice rheology 
     27   USE limrhg           ! LIM  : EVP ice rheology 
     28   USE lbclnk           ! lateral boundary condition - MPP link 
     29   USE lib_mpp          ! MPP library 
     30   USE in_out_manager   ! I/O manager 
     31   USE prtctl           ! Print control 
    3332 
    3433   IMPLICIT NONE 
    3534   PRIVATE 
    3635 
    37    PUBLIC   lim_dyn_2 ! routine called by sbc_ice_lim 
    38  
    39    !! * Module variables 
    40    REAL(wp)  ::  rone    = 1.e0   ! constant value 
    41  
     36   PUBLIC   lim_dyn_2   ! routine called by sbc_ice_lim 
     37 
     38   !! * Substitutions 
    4239#  include "vectopt_loop_substitute.h90" 
    4340   !!---------------------------------------------------------------------- 
    4441   !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 
    4542   !! $Id$ 
    46    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    47    !!---------------------------------------------------------------------- 
    48  
     43   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     44   !!---------------------------------------------------------------------- 
    4945CONTAINS 
    5046 
     
    9086            i_jpj = jpj 
    9187            IF(ln_ctl)   CALL prt_ctl_info( 'lim_dyn  :    i_j1 = ', ivar1=i_j1, clinfo2=' ij_jpj = ', ivar2=i_jpj ) 
    92 #if defined key_lim2_vp 
    93             CALL lim_rhg_2( i_j1, i_jpj )             !  VP rheology 
    94 #else 
    95             CALL lim_rhg  ( i_j1, i_jpj )             ! EVP rheology 
    96 #endif 
     88            IF( lk_lim2_vp )   THEN   ;   CALL lim_rhg_2( i_j1, i_jpj )             !  VP rheology 
     89            ELSE                      ;   CALL lim_rhg  ( i_j1, i_jpj )             ! EVP rheology 
     90            ENDIF 
    9791            ! 
    9892         ELSE                                 ! optimization of the computational area 
     
    112106                  i_j1 = i_j1 + 1 
    113107               END DO 
    114 #if defined key_lim2_vp 
    115                i_j1 = MAX( 1, i_j1-1 ) 
    116                IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : NH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 
    117                !  
    118                CALL lim_rhg_2( i_j1, i_jpj ) 
    119                !  
    120 #else 
    121                i_j1 = MAX( 1, i_j1-2 ) 
     108               IF( lk_lim2_vp )   THEN             ! VP  rheology 
     109                  i_j1 = MAX( 1, i_j1-1 ) 
     110                  CALL lim_rhg_2( i_j1, i_jpj ) 
     111               ELSE                                ! EVP rheology 
     112                  i_j1 = MAX( 1, i_j1-2 ) 
     113                  CALL lim_rhg( i_j1, i_jpj ) 
     114               ENDIF 
    122115               IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : NH i_j1 = ', i_j1, 'ij_jpj = ', i_jpj 
    123                CALL lim_rhg( i_j1, i_jpj ) 
    124                !  
    125 #endif 
     116               ! 
    126117               ! Southern hemisphere 
    127118               i_j1  =  1  
     
    130121                  i_jpj = i_jpj - 1 
    131122               END DO 
    132 #if defined key_lim2_vp 
    133                i_jpj = MIN( jpj, i_jpj+2 ) 
    134                IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : SH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 
    135                !  
    136                CALL lim_rhg_2( i_j1, i_jpj ) 
    137                !  
    138 #else 
    139                i_jpj = MIN( jpj, i_jpj+1 ) 
     123               IF( lk_lim2_vp )   THEN             ! VP  rheology 
     124                  i_jpj = MIN( jpj, i_jpj+2 ) 
     125                  CALL lim_rhg_2( i_j1, i_jpj ) 
     126               ELSE                                ! EVP rheology 
     127                  i_jpj = MIN( jpj, i_jpj+1 ) 
     128                  CALL lim_rhg( i_j1, i_jpj ) 
     129               ENDIF 
    140130               IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : SH i_j1 = ', i_j1, 'ij_jpj = ', i_jpj 
    141                CALL lim_rhg( i_j1, i_jpj ) 
    142                !  
    143 #endif 
    144  
     131               ! 
    145132            ELSE                                 ! local domain extends over one hemisphere only 
    146133               !                                 ! Rheology is computed only over the ice cover 
     
    156143                  i_jpj = i_jpj - 1 
    157144               END DO 
    158                i_jpj = MIN( jpj, i_jpj+2) 
    159      
     145               i_jpj = MIN( jpj, i_jpj+2 ) 
     146               !  
     147               IF( lk_lim2_vp )   THEN             ! VP  rheology 
     148                  i_jpj = MIN( jpj, i_jpj+2 ) 
     149                  CALL lim_rhg_2( i_j1, i_jpj )                !  VP rheology 
     150               ELSE                                ! EVP rheology 
     151                  i_j1  = MAX( 1  , i_j1-2  ) 
     152                  i_jpj = MIN( jpj, i_jpj+1 ) 
     153                  CALL lim_rhg  ( i_j1, i_jpj )                ! EVP rheology 
     154               ENDIF 
    160155               IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : one hemisphere: i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 
    161                !  
    162 #if defined key_lim2_vp 
    163                i_jpj = MIN( jpj, i_jpj+2) 
    164                CALL lim_rhg_2( i_j1, i_jpj )                !  VP rheology 
    165 #else 
    166                i_j1 = MAX( 1, i_j1-2 ) 
    167                i_jpj = MIN( jpj, i_jpj+1) 
    168                CALL lim_rhg  ( i_j1, i_jpj )                ! EVP rheology 
    169 #endif 
    170156               ! 
    171157            ENDIF 
     
    177163         ! computation of friction velocity 
    178164         ! -------------------------------- 
    179          SELECT CASE( cl_grid ) 
    180          CASE( 'C' )       ! C-grid ice dynamics (EVP) 
    181                            ! ice-ocean & ice velocity at  ocean velocity points 
    182             zu_io(:,:) = u_ice(:,:) - ssu_m(:,:)   
     165         SELECT CASE( cp_ice_msh )           ! ice-ocean relative velocity at u- & v-pts 
     166         CASE( 'C' )                               ! EVP : C-grid ice dynamics 
     167            zu_io(:,:) = u_ice(:,:) - ssu_m(:,:)           ! ice-ocean & ice velocity at ocean velocity points 
    183168            zv_io(:,:) = v_ice(:,:) - ssv_m(:,:) 
    184             ! 
    185          CASE( 'B' )       ! B-grid ice dynamics (VP)  
    186                            ! ice-ocean velocity at U & V-points (u_ice v_ice at I-point ; ssu_m, ssv_m at U- & V-points) 
    187             DO jj = 1, jpjm1 
    188                DO ji = 1, jpim1   ! NO vector opt. 
    189                   zu_io(ji,jj) = 0.5 * ( u_ice(ji+1,jj+1) + u_ice(ji+1,jj  ) ) - ssu_m(ji,jj) 
    190                   zv_io(ji,jj) = 0.5 * ( v_ice(ji+1,jj+1) + v_ice(ji  ,jj+1) ) - ssv_m(ji,jj) 
     169         CASE( 'I' )                               ! VP  : B-grid ice dynamics (I-point)  
     170            DO jj = 1, jpjm1                               ! u_ice v_ice at I-point ; ssu_m, ssv_m at U- & V-points 
     171               DO ji = 1, jpim1   ! NO vector opt.         ! 
     172                  zu_io(ji,jj) = 0.5_wp * ( u_ice(ji+1,jj+1) + u_ice(ji+1,jj  ) ) - ssu_m(ji,jj) 
     173                  zv_io(ji,jj) = 0.5_wp * ( v_ice(ji+1,jj+1) + v_ice(ji  ,jj+1) ) - ssv_m(ji,jj) 
    191174               END DO 
    192175            END DO 
     
    194177 
    195178         ! frictional velocity at T-point 
     179         zcoef = 0.5_wp * cw 
    196180         DO jj = 2, jpjm1 
    197181            DO ji = 2, jpim1   ! NO vector opt. because of zu_io 
    198                ust2s(ji,jj) = 0.5 * cw                                                          & 
    199                   &         * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
    200                   &            + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
     182               ust2s(ji,jj) = zcoef * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
     183                  &                    + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
    201184            END DO 
    202185         END DO 
     
    207190         DO jj = 2, jpjm1 
    208191            DO ji = fs_2, fs_jpim1   ! vector opt. 
    209                ust2s(ji,jj) = zcoef * tms(ji,jj) * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
    210                   &                                     + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) 
     192               ust2s(ji,jj) = zcoef * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
     193                  &                        + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1)   ) * tms(ji,jj) 
    211194            END DO 
    212195         END DO 
     
    217200      ! 
    218201      IF(ln_ctl)   CALL prt_ctl(tab2d_1=ust2s , clinfo1=' lim_dyn  : ust2s :') 
    219  
     202      ! 
    220203   END SUBROUTINE lim_dyn_2 
    221204 
     
    265248         WRITE(numout,*) '       coefficient for the solution of int. stresses alphaevp = ', alphaevp 
    266249      ENDIF 
     250      ! 
     251      IF( angvg /= 0._wp .AND. .NOT.lk_lim2_vp ) THEN 
     252         CALL ctl_warn( 'lim_dyn_init_2: turning angle for oceanic stress not properly coded for EVP ',   & 
     253            &           '(see limsbc_2 module). We force  angvg = 0._wp'  ) 
     254         angvg = 0._wp 
     255      ENDIF 
    267256 
    268257      !  Initialization 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limmsh_2.F90

    r2319 r2370  
    44   !! LIM 2.0 ice model :   definition of the ice mesh parameters 
    55   !!====================================================================== 
     6   !! History :   -   ! 2001-04 (LIM) original code 
     7   !!            1.0  ! 2002-08 (C. Ethe, G. Madec) F90, module 
     8   !!            3.3  ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 
     9   !!---------------------------------------------------------------------- 
    610#if defined key_lim2 
    711   !!---------------------------------------------------------------------- 
     
    1014   !!   lim_msh_2   : definition of the ice mesh 
    1115   !!---------------------------------------------------------------------- 
    12    !! * Modules used 
    1316   USE phycst 
    1417   USE dom_oce 
     
    2023   PRIVATE 
    2124 
    22    !! * Accessibility 
    2325   PUBLIC lim_msh_2      ! routine called by ice_ini_2.F90 
    2426 
     
    2628   !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 
    2729   !! $Id$ 
    28    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    29    !!---------------------------------------------------------------------- 
    30  
     30   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     31   !!---------------------------------------------------------------------- 
    3132CONTAINS 
    3233 
     
    4344      !!  
    4445      !! ** Refer.  : Deleersnijder et al. Ocean Modelling 100, 7-10  
    45       !! 
    46       !! ** History : 
    47       !!         original    : 01-04 (LIM) 
    48       !!         addition    : 02-08 (C. Ethe, G. Madec) 
    4946      !!---------------------------------------------------------------------  
    50       !! * Local variables 
    5147      INTEGER :: ji, jj      ! dummy loop indices 
    5248      REAL(wp) ::   zusden   ! local scalars 
    53  
    5449#if defined key_lim2_vp 
    5550      REAL(wp) ::   zusden2           ! local scalars 
     
    267262         END DO 
    268263      END DO 
    269        
    270       !--lateral boundary conditions     
    271       CALL lbc_lnk( tmu(:,:), 'I', 1. ) 
     264      CALL lbc_lnk( tmu(:,:), 'I', 1. )      !--lateral boundary conditions     
    272265#else 
    273266      ! EVP rheology : ice velocity point are U- & V-points ; ice vorticity 
     
    278271      WHERE( fmask(:,:,1) == 1.e0 )   tmf(:,:) = 1.e0 
    279272#endif 
    280        
     273      ! 
    281274      ! unmasked and masked area of T-grid cell 
    282275      area(:,:) = e1t(:,:) * e2t(:,:) 
    283        
     276      ! 
    284277   END SUBROUTINE lim_msh_2 
    285278 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limrhg_2.F90

    r2319 r2370  
    77   !!            1.0  !  1994-12  (H. Goosse)  
    88   !!            2.0  !  2003-12  (C. Ethe, G. Madec)  F90, mpp 
    9    !!            " "  !  2006-08  (G. Madec)  surface module, ice-stress at I-point 
    10    !!            " "  !  2009-09  (G. Madec)  Huge verctor optimisation 
     9   !!             -   !  2006-08  (G. Madec)  surface module, ice-stress at I-point 
     10   !!             -   !  2009-09  (G. Madec)  Huge verctor optimisation 
    1111   !!            3.3  !  2009-05  (G.Garric, C. Bricaud) addition of the lim2_evp case 
    1212   !!---------------------------------------------------------------------- 
    1313#if defined   key_lim2   &&   defined key_lim2_vp 
    1414   !!---------------------------------------------------------------------- 
    15    !!   'key_lim2'                and                 LIM 2.0 sea-ice model 
     15   !!   'key_lim2'                AND                   LIM-2 sea-ice model 
    1616   !!   'key_lim2_vp'                                       VP ice rheology 
    17    !!---------------------------------------------------------------------- 
    1817   !!---------------------------------------------------------------------- 
    1918   !!   lim_rhg_2   : computes ice velocities 
     
    3635   PUBLIC   lim_rhg_2 ! routine called by lim_dyn 
    3736 
    38    REAL(wp) ::   rzero   = 0.e0   ! constant value: zero 
    39    REAL(wp) ::   rone    = 1.e0   !            and  one 
     37   REAL(wp) ::   rzero   = 0._wp   ! constant value: zero 
     38   REAL(wp) ::   rone    = 1._wp   !            and  one 
    4039 
    4140   !! * Substitutions 
     
    4443   !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 
    4544   !! $Id$ 
    46    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    47    !!---------------------------------------------------------------------- 
    48  
     45   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     46   !!---------------------------------------------------------------------- 
    4947CONTAINS 
    5048 
     
    6765      INTEGER ::   iter, jter          ! temporary integers 
    6866      CHARACTER (len=50) ::   charout 
    69       REAL(wp) ::   ze11  , ze12  , ze22  , ze21               ! temporary scalars 
    70       REAL(wp) ::   zt11  , zt12  , zt21  , zt22               !    "         " 
    71       REAL(wp) ::   zvis11, zvis21, zvis12, zvis22             !    "         " 
    72       REAL(wp) ::   zgphsx, ztagnx, zgsshx, zunw, zur, zusw    !    "         " 
    73       REAL(wp) ::   zgphsy, ztagny, zgsshy, zvnw, zvr          !    "         " 
     67      REAL(wp) ::   ze11  , ze12  , ze22  , ze21               ! local scalars 
     68      REAL(wp) ::   zt11  , zt12  , zt21  , zt22               !   -      - 
     69      REAL(wp) ::   zvis11, zvis21, zvis12, zvis22             !   -      - 
     70      REAL(wp) ::   zgphsx, ztagnx, zgsshx, zunw, zur, zusw    !   -      - 
     71      REAL(wp) ::   zgphsy, ztagny, zgsshy, zvnw, zvr          !   -      - 
    7472      REAL(wp) ::   zresm,  za, zac, zmod 
    7573      REAL(wp) ::   zmpzas, zstms, zindu, zusdtp, zmassdt, zcorlal 
     
    9189      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zzfrld, zztms 
    9290      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zi1, zi2, zmasst, zpresh 
    93  
    9491      !!------------------------------------------------------------------- 
    95  
    96 !!bug 
    97 !!    u_oce(:,:) = 0.e0 
    98 !!    v_oce(:,:) = 0.e0 
    99 !!    write(*,*) 'rhg min, max u & v', maxval(u_oce), minval(u_oce), maxval(v_oce), minval(v_oce) 
    100 !!bug 
    10192       
    10293      !  Store initial velocities 
    10394      !  ---------------- 
    104       zztms(:,0    ) = 0.e0       ;    zzfrld(:,0    ) = 0.e0 
    105       zztms(:,jpj+1) = 0.e0       ;    zzfrld(:,jpj+1) = 0.e0 
    106       zu0(:,0    ) = 0.e0         ;    zv0(:,0    ) = 0.e0 
    107       zu0(:,jpj+1) = 0.e0         ;    zv0(:,jpj+1) = 0.e0 
    108       zztms(:,1:jpj) = tms(:,:)   ;    zzfrld(:,1:jpj) = frld(:,:) 
    109       zu0(:,1:jpj) = u_ice(:,:)   ;    zv0(:,1:jpj) = v_ice(:,:) 
    110  
    111       zu_a(:,:)    = zu0(:,:)     ;   zv_a(:,:) = zv0(:,:) 
    112       zu_n(:,:)    = zu0(:,:)     ;   zv_n(:,:) = zv0(:,:) 
     95      zztms(:,0    ) = 0._wp        ;   zzfrld(:,0    ) = 0._wp 
     96      zztms(:,jpj+1) = 0._wp        ;   zzfrld(:,jpj+1) = 0._wp 
     97      zu0  (:,0    ) = 0._wp        ;   zv0   (:,0    ) = 0._wp 
     98      zu0  (:,jpj+1) = 0._wp        ;   zv0   (:,jpj+1) = 0._wp 
     99      zztms(:,1:jpj) = tms  (:,:)   ;   zzfrld(:,1:jpj) = frld (:,:) 
     100      zu0  (:,1:jpj) = u_ice(:,:)   ;   zv0   (:,1:jpj) = v_ice(:,:) 
     101      zu_a (:, :   ) = zu0  (:,:)   ;   zv_a  (:, :   ) = zv0  (:,:) 
     102      zu_n (:, :   ) = zu0  (:,:)   ;   zv_n  (:, :   ) = zv0  (:,:) 
    113103 
    114104!i 
    115       zi1   (:,:) = 0.e0 
    116       zi2   (:,:) = 0.e0 
    117       zpresh(:,:) = 0.e0 
    118       zmasst(:,:) = 0.e0 
     105      zi1   (:,:) = 0._wp 
     106      zi2   (:,:) = 0._wp 
     107      zpresh(:,:) = 0._wp 
     108      zmasst(:,:) = 0._wp 
    119109!i 
    120110!!gm violant 
    121       zfrld(:,:) =0.e0 
    122       zcorl(:,:) =0.e0 
    123       zmass(:,:) =0.e0 
    124       za1ct(:,:) =0.e0 
    125       za2ct(:,:) =0.e0 
     111      zfrld(:,:) =0._wp 
     112      zcorl(:,:) =0._wp 
     113      zmass(:,:) =0._wp 
     114      za1ct(:,:) =0._wp 
     115      za2ct(:,:) =0._wp 
    126116!!gm end 
    127117 
    128       zviszeta(:,:) = 0.e0 
    129       zviseta (:,:) = 0.e0 
    130  
    131 !i    zviszeta(:,0    ) = 0.e0    ;    zviseta(:,0    ) = 0.e0 
    132 !i    zviszeta(:,jpj  ) = 0.e0    ;    zviseta(:,jpj  ) = 0.e0 
    133 !i    zviszeta(:,jpj+1) = 0.e0    ;    zviseta(:,jpj+1) = 0.e0 
     118      zviszeta(:,:) = 0._wp 
     119      zviseta (:,:) = 0._wp 
     120 
     121!i    zviszeta(:,0    ) = 0._wp    ;    zviseta(:,0    ) = 0._wp 
     122!i    zviszeta(:,jpj  ) = 0._wp    ;    zviseta(:,jpj  ) = 0._wp 
     123!i    zviszeta(:,jpj+1) = 0._wp    ;    zviseta(:,jpj+1) = 0._wp 
    134124 
    135125 
     
    143133         DO ji = 1 , jpi 
    144134            ! only the sinus changes its sign with the hemisphere 
    145             zsang(ji,jj)  = SIGN( 1.e0, fcor(ji,jj) ) * sangvg   ! only the sinus changes its sign with the hemisphere 
     135            zsang(ji,jj)  = SIGN( 1._wp, fcor(ji,jj) ) * sangvg   ! only the sinus changes its sign with the hemisphere 
    146136            ! 
    147137            zmasst(ji,jj) = tms(ji,jj) * ( rhosn * hsnm(ji,jj) + rhoic * hicm(ji,jj) ) 
    148138            zpresh(ji,jj) = tms(ji,jj) *  pstarh * hicm(ji,jj) * EXP( -c_rhg * frld(ji,jj) ) 
    149139!!gm  :: stress given at I-point (F-point for the ocean) only compute the ponderation with the ice fraction (1-frld) 
    150             zi1(ji,jj)    = tms(ji,jj) * ( 1.0 - frld(ji,jj) ) 
    151             zi2(ji,jj)    = tms(ji,jj) * ( 1.0 - frld(ji,jj) ) 
     140            zi1(ji,jj)    = tms(ji,jj) * ( 1._wp - frld(ji,jj) ) 
     141            zi2(ji,jj)    = tms(ji,jj) * ( 1._wp - frld(ji,jj) ) 
    152142         END DO 
    153143      END DO 
     
    163153            zstms = zztms(ji,jj  ) * wght(ji,jj,2,2) + zztms(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    164154               &  + zztms(ji,jj-1) * wght(ji,jj,2,1) + zztms(ji-1,jj-1) * wght(ji,jj,1,1) 
    165             zusw  = 1.0 / MAX( zstms, epsd ) 
     155            zusw  = 1._wp / MAX( zstms, epsd ) 
    166156 
    167157            zt11 = zztms(ji  ,jj  ) * zzfrld(ji  ,jj  )  
     
    201191            ! Gradient of the sea surface height 
    202192            zgsshx =  (   (ssh_m(ji  ,jj  ) - ssh_m(ji-1,jj  ))/e1u(ji-1,jj  )   & 
    203                &       +  (ssh_m(ji  ,jj-1) - ssh_m(ji-1,jj-1))/e1u(ji-1,jj-1)   ) * 0.5 
     193               &       +  (ssh_m(ji  ,jj-1) - ssh_m(ji-1,jj-1))/e1u(ji-1,jj-1)   ) * 0.5_wp 
    204194            zgsshy =  (   (ssh_m(ji  ,jj  ) - ssh_m(ji  ,jj-1))/e2v(ji  ,jj-1)   & 
    205                &       +  (ssh_m(ji-1,jj  ) - ssh_m(ji-1,jj-1))/e2v(ji-1,jj-1)   ) * 0.5 
     195               &       +  (ssh_m(ji-1,jj  ) - ssh_m(ji-1,jj-1))/e2v(ji-1,jj-1)   ) * 0.5_wp 
    206196 
    207197            ! Computation of the velocity field taking into account the ice-ice interaction.                                  
     
    219209         !                                                ! ==================== !         
    220210         zindu = MOD( iter , 2 ) 
    221          zusdtp = ( zindu * 2.0 + ( 1.0 - zindu ) * 1.0 )  * REAL( nbiter ) / rdt_ice 
     211         zusdtp = ( zindu * 2._wp + ( 1._wp - zindu ) * 1._wp )  * REAL( nbiter ) / rdt_ice 
    222212 
    223213         ! Computation of free drift field for free slip boundary conditions. 
     
    241231               zdgi = zt12 + zt21 
    242232               ztrace2 = zdgp * zdgp  
    243                zdeter  = zt11 * zt22 - 0.25 * zdgi * zdgi 
     233               zdeter  = zt11 * zt22 - 0.25_wp * zdgi * zdgi 
    244234 
    245235               !  Creep limit depends on the size of the grid. 
    246                zdelta = MAX( SQRT( ztrace2 + ( ztrace2 - 4.0 * zdeter ) * usecc2 ),  creepl) 
     236               zdelta = MAX( SQRT( ztrace2 + ( ztrace2 - 4._wp * zdeter ) * usecc2 ),  creepl) 
    247237 
    248238               !-  Computation of viscosities. 
     
    256246            DO ji = 2, fs_jpim1   ! NO vector opt. 
    257247               !* zc1u , zc2v 
    258                zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
    259                zvis12 =       zviseta (ji-1,jj-1) + dm 
    260                zvis21 =       zviseta (ji-1,jj-1) 
    261                zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
     248               zvis11 = 2._wp * zviseta (ji-1,jj-1) + dm 
     249               zvis12 =         zviseta (ji-1,jj-1) + dm 
     250               zvis21 =         zviseta (ji-1,jj-1) 
     251               zvis22 =         zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
    262252               zdiag  = zvis22 * ( akappa(ji-1,jj-1,1,1) + akappa(ji-1,jj-1,2,1) ) 
    263253               zs11_11 =  zvis11 * akappa(ji-1,jj-1,1,1) + zdiag 
     
    266256               zs22_11 =  zvis11 * akappa(ji-1,jj-1,2,1) + zdiag 
    267257 
    268                zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
    269                zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
    270                zvis12 =       zviseta (ji,jj-1) + dm 
    271                zvis21 =       zviseta (ji,jj-1) 
     258               zvis11 = 2._wp * zviseta (ji,jj-1) + dm 
     259               zvis22 =         zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     260               zvis12 =         zviseta (ji,jj-1) + dm 
     261               zvis21 =         zviseta (ji,jj-1) 
    272262               zdiag = zvis22 * ( -akappa(ji,jj-1,1,1) + akappa(ji,jj-1,2,1) ) 
    273263               zs11_21 = -zvis11 * akappa(ji,jj-1,1,1) + zdiag 
     
    276266               zs21_21 = -zvis12 * akappa(ji,jj-1,1,2) + zvis21 * akappa(ji,jj-1,2,2) 
    277267 
    278                zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
    279                zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
    280                zvis12 =       zviseta (ji-1,jj) + dm 
    281                zvis21 =       zviseta (ji-1,jj) 
    282                zdiag = zvis22 * ( akappa(ji-1,jj,1,1) + akappa(ji-1,jj,2,1) ) 
     268               zvis11 = 2._wp * zviseta (ji-1,jj) + dm 
     269               zvis22 =         zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     270               zvis12 =         zviseta (ji-1,jj) + dm 
     271               zvis21 =         zviseta (ji-1,jj) 
     272               zdiag  = zvis22 * ( akappa(ji-1,jj,1,1) + akappa(ji-1,jj,2,1) ) 
    283273               zs11_12 =  zvis11 * akappa(ji-1,jj,1,1) + zdiag 
    284274               zs12_12 = -zvis12 * akappa(ji-1,jj,2,2) - zvis21 * akappa(ji-1,jj,1,2) 
     
    286276               zs21_12 = -zvis12 * akappa(ji-1,jj,1,2) - zvis21 * akappa(ji-1,jj,2,2) 
    287277 
    288                zvis11 = 2.0 * zviseta (ji,jj) + dm 
    289                zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
    290                zvis12 =       zviseta (ji,jj) + dm 
    291                zvis21 =       zviseta (ji,jj) 
     278               zvis11 = 2._wp * zviseta (ji,jj) + dm 
     279               zvis22 =         zviszeta(ji,jj) - zviseta(ji,jj) 
     280               zvis12 =         zviseta (ji,jj) + dm 
     281               zvis21 =         zviseta (ji,jj) 
    292282               zdiag = zvis22 * ( -akappa(ji,jj,1,1) + akappa(ji,jj,2,1) ) 
    293283               zs11_22 = -zvis11 * akappa(ji,jj,1,1) + zdiag 
     
    315305 
    316306               !* zc1v , zc2v. 
    317                zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
    318                zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
    319                zvis12 =       zviseta (ji-1,jj-1) + dm 
    320                zvis21 =       zviseta (ji-1,jj-1) 
     307               zvis11 = 2._wp * zviseta (ji-1,jj-1) + dm 
     308               zvis22 =         zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
     309               zvis12 =         zviseta (ji-1,jj-1) + dm 
     310               zvis21 =         zviseta (ji-1,jj-1) 
    321311               zdiag = zvis22 * ( akappa(ji-1,jj-1,1,2) + akappa(ji-1,jj-1,2,2) ) 
    322312               zs11_11 =  zvis11 * akappa(ji-1,jj-1,1,2) + zdiag 
     
    325315               zs21_11 =  zvis12 * akappa(ji-1,jj-1,1,1) - zvis21 * akappa(ji-1,jj-1,2,1) 
    326316  
    327                zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
    328                zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
    329                zvis12 =       zviseta (ji,jj-1) + dm 
    330                zvis21 =       zviseta (ji,jj-1) 
     317               zvis11 = 2._wp * zviseta (ji,jj-1) + dm 
     318               zvis22 =         zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     319               zvis12 =         zviseta (ji,jj-1) + dm 
     320               zvis21 =         zviseta (ji,jj-1) 
    331321               zdiag = zvis22 * ( akappa(ji,jj-1,1,2) + akappa(ji,jj-1,2,2) ) 
    332322               zs11_21 =  zvis11 * akappa(ji,jj-1,1,2) + zdiag 
     
    335325               zs21_21 = -zvis12 * akappa(ji,jj-1,1,1) - zvis21 * akappa(ji,jj-1,2,1) 
    336326 
    337                zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
    338                zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
    339                zvis12 =       zviseta (ji-1,jj) + dm 
    340                zvis21 =       zviseta (ji-1,jj) 
     327               zvis11 = 2._wp * zviseta (ji-1,jj) + dm 
     328               zvis22 =         zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     329               zvis12 =         zviseta (ji-1,jj) + dm 
     330               zvis21 =         zviseta (ji-1,jj) 
    341331               zdiag = zvis22 * ( akappa(ji-1,jj,1,2) - akappa(ji-1,jj,2,2) ) 
    342332               zs11_12 =  zvis11 * akappa(ji-1,jj,1,2) + zdiag 
     
    345335               zs21_12 =  zvis12 * akappa(ji-1,jj,1,1) - zvis21 * akappa(ji-1,jj,2,1) 
    346336 
    347                zvis11 = 2.0 * zviseta (ji,jj) + dm 
    348                zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
    349                zvis12 =       zviseta (ji,jj) + dm 
    350                zvis21 =       zviseta (ji,jj) 
     337               zvis11 = 2._wp * zviseta (ji,jj) + dm 
     338               zvis22 =         zviszeta(ji,jj) - zviseta(ji,jj) 
     339               zvis12 =         zviseta (ji,jj) + dm 
     340               zvis21 =         zviseta (ji,jj) 
    351341               zdiag = zvis22 * ( akappa(ji,jj,1,2) - akappa(ji,jj,2,2) ) 
    352342               zs11_22 =  zvis11 * akappa(ji,jj,1,2) + zdiag 
     
    388378                  ze22 = + akappa(ji,jj-1,2,2) * zv_a(ji+1,jj) + akappa(ji,jj-1,2,1) * zu_a(ji+1,jj) 
    389379                  ze21 =   akappa(ji,jj-1,1,1) * zv_a(ji+1,jj) - akappa(ji,jj-1,1,2) * zu_a(ji+1,jj) 
    390                   zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
    391                   zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
    392                   zvis12 =       zviseta (ji,jj-1) + dm 
    393                   zvis21 =       zviseta (ji,jj-1) 
     380                  zvis11 = 2._wp * zviseta (ji,jj-1) + dm 
     381                  zvis22 =         zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     382                  zvis12 =         zviseta (ji,jj-1) + dm 
     383                  zvis21 =         zviseta (ji,jj-1) 
    394384                  zdiag = zvis22 * ( ze11 + ze22 ) 
    395385                  zs11_21 =  zvis11 * ze11 + zdiag 
     
    406396                  ze21 =   akappa(ji-1,jj,1,1) * ( zv_a(ji  ,jj+1) - zv_a(ji-1,jj+1) )   & 
    407397                     &   - akappa(ji-1,jj,1,2) * ( zu_a(ji  ,jj+1) + zu_a(ji-1,jj+1) ) 
    408                   zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
    409                   zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
    410                   zvis12 =       zviseta (ji-1,jj) + dm 
    411                   zvis21 =       zviseta (ji-1,jj) 
     398                  zvis11 = 2._wp * zviseta (ji-1,jj) + dm 
     399                  zvis22 =         zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     400                  zvis12 =         zviseta (ji-1,jj) + dm 
     401                  zvis21 =         zviseta (ji-1,jj) 
    412402                  zdiag = zvis22 * ( ze11 + ze22 ) 
    413403                  zs11_12 =  zvis11 * ze11 + zdiag 
     
    424414                  ze21 =   akappa(ji,jj,1,1) * ( zv_a(ji+1,jj) + zv_a(ji+1,jj+1) - zv_a(ji  ,jj+1) )   & 
    425415                     &   - akappa(ji,jj,1,2) * ( zu_a(ji+1,jj) + zu_a(ji+1,jj+1) + zu_a(ji  ,jj+1) ) 
    426                   zvis11 = 2.0 * zviseta (ji,jj) + dm 
    427                   zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
    428                   zvis12 =       zviseta (ji,jj) + dm 
    429                   zvis21 =       zviseta (ji,jj) 
     416                  zvis11 = 2._wp * zviseta (ji,jj) + dm 
     417                  zvis22 =         zviszeta(ji,jj) - zviseta(ji,jj) 
     418                  zvis12 =         zviseta (ji,jj) + dm 
     419                  zvis21 =         zviseta (ji,jj) 
    430420                  zdiag = zvis22 * ( ze11 + ze22 ) 
    431421                  zs11_22 =  zvis11 * ze11 + zdiag 
     
    443433                  ze21 =   akappa(ji-1,jj-1,1,1) * ( zv_a(ji  ,jj-1) - zv_a(ji-1,jj-1) - zv_a(ji-1,jj) )   & 
    444434                     &  -  akappa(ji-1,jj-1,1,2) * ( zu_a(ji  ,jj-1) + zu_a(ji-1,jj-1) + zu_a(ji-1,jj) ) 
    445                   zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
    446                   zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
    447                   zvis12 =       zviseta (ji-1,jj-1) + dm 
    448                   zvis21 =       zviseta (ji-1,jj-1) 
     435                  zvis11 = 2._wp * zviseta (ji-1,jj-1) + dm 
     436                  zvis22 =         zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
     437                  zvis12 =         zviseta (ji-1,jj-1) + dm 
     438                  zvis21 =         zviseta (ji-1,jj-1) 
    449439                  zdiag = zvis22 * ( ze11 + ze22 ) 
    450440                  zs11_11 =  zvis11 * ze11 + zdiag 
     
    461451                  ze21 =   akappa(ji,jj-1,1,1) * ( zv_a(ji+1,jj-1) - zv_a(ji  ,jj-1) )   & 
    462452                     &   - akappa(ji,jj-1,1,2) * ( zu_a(ji+1,jj-1) + zu_a(ji  ,jj-1) ) 
    463                   zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
    464                   zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
    465                   zvis12 =       zviseta (ji,jj-1) + dm 
    466                   zvis21 =       zviseta (ji,jj-1) 
     453                  zvis11 = 2._wp * zviseta (ji,jj-1) + dm 
     454                  zvis22 =         zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     455                  zvis12 =         zviseta (ji,jj-1) + dm 
     456                  zvis21 =         zviseta (ji,jj-1) 
    467457                  zdiag = zvis22 * ( ze11 + ze22 ) 
    468458                  zs11_21 =  zs11_21 + zvis11 * ze11 + zdiag 
     
    475465                  ze22 = - akappa(ji-1,jj,2,2) * zv_a(ji-1,jj) + akappa(ji-1,jj,2,1) * zu_a(ji-1,jj) 
    476466                  ze21 = - akappa(ji-1,jj,1,1) * zv_a(ji-1,jj) - akappa(ji-1,jj,1,2) * zu_a(ji-1,jj) 
    477                   zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
    478                   zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
    479                   zvis12 =       zviseta (ji-1,jj) + dm 
    480                   zvis21 =       zviseta (ji-1,jj) 
     467                  zvis11 = 2._wp * zviseta (ji-1,jj) + dm 
     468                  zvis22 =         zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     469                  zvis12 =         zviseta (ji-1,jj) + dm 
     470                  zvis21 =         zviseta (ji-1,jj) 
    481471                  zdiag = zvis22 * ( ze11 + ze22 ) 
    482472                  zs11_12 =  zs11_12 + zvis11 * ze11 + zdiag 
     
    506496                  zvr     = zv_a(ji,jj) - v_oce(ji,jj) 
    507497!!!! 
    508                   zmod    = SQRT( zur*zur + zvr*zvr ) * ( 1.0 - zfrld(ji,jj) ) 
     498                  zmod    = SQRT( zur*zur + zvr*zvr ) * ( 1._wp - zfrld(ji,jj) ) 
    509499                  za      = rhoco * zmod 
    510500!!!! 
    511 !!gm chg resul    za      = rhoco * SQRT( zur*zur + zvr*zvr ) * ( 1.0 - zfrld(ji,jj) ) 
     501!!gm chg resul    za      = rhoco * SQRT( zur*zur + zvr*zvr ) * ( 1._wp - zfrld(ji,jj) ) 
    512502                  zac     = za * cangvg 
    513503                  zmpzas  = alpha * zcorl(ji,jj) + za * zsang(ji,jj) 
    514504                  zmassdt = zusdtp * zmass(ji,jj) 
    515                   zcorlal = ( 1.0 - alpha ) * zcorl(ji,jj) 
     505                  zcorlal = ( 1._wp - alpha ) * zcorl(ji,jj) 
    516506 
    517507                  za1 =  zmassdt * zu0(ji,jj) + zcorlal * zv0(ji,jj) + za1ct(ji,jj)   & 
     
    527517                  zunw   = (  ( za1 + zd1 ) * zc2 + ( za2 + zd2 ) * zc1 ) * zden 
    528518                  zvnw   = (  ( za2 + zd2 ) * zb1 - ( za1 + zd1 ) * zb2 ) * zden 
    529                   zmask  = ( 1.0 - MAX( rzero, SIGN( rone , 1.0 - zmass(ji,jj) ) ) ) * tmu(ji,jj) 
     519                  zmask  = ( 1._wp - MAX( rzero, SIGN( rone , 1._wp - zmass(ji,jj) ) ) ) * tmu(ji,jj) 
    530520 
    531521                  zu_n(ji,jj) = ( zu_a(ji,jj) + om * ( zunw - zu_a(ji,jj) ) * tmu(ji,jj) ) * zmask 
     
    582572#else 
    583573   !!---------------------------------------------------------------------- 
    584    !!   Default option          Dummy module       NO 2.0 LIM sea-ice model 
     574   !!   Default option        Dummy module      NO VP & LIM-2 sea-ice model 
    585575   !!---------------------------------------------------------------------- 
    586576CONTAINS 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limrst_2.F90

    r2319 r2370  
    121121      CALL iom_rstput( iter, nitrst, numriw, 'fsbbq'      , fsbbq (:,:)   ) 
    122122#if ! defined key_lim2_vp 
    123       CALL iom_rstput( iter, nitrst, numriw, 'stress1_i'  , stress1_i (:,:)   ) 
    124       CALL iom_rstput( iter, nitrst, numriw, 'stress2_i'  , stress2_i (:,:)   ) 
    125       CALL iom_rstput( iter, nitrst, numriw, 'stress12_i' , stress12_i(:,:)   ) 
     123      CALL iom_rstput( iter, nitrst, numriw, 'stress1_i'  , stress1_i (:,:) )    ! EVP rheology 
     124      CALL iom_rstput( iter, nitrst, numriw, 'stress2_i'  , stress2_i (:,:) ) 
     125      CALL iom_rstput( iter, nitrst, numriw, 'stress12_i' , stress12_i(:,:) ) 
    126126#endif 
    127127      CALL iom_rstput( iter, nitrst, numriw, 'sxice'      , sxice (:,:)   ) 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limsbc_2.F90

    r2319 r2370  
    22   !!====================================================================== 
    33   !!                       ***  MODULE limsbc_2   *** 
    4    !!           computation of the flux at the sea ice/ocean interface 
     4   !! LIM-2 :   updates the fluxes at the ocean surface with ice-ocean fluxes 
    55   !!====================================================================== 
    66   !! History :  LIM  ! 2000-01 (H. Goosse) Original code 
    77   !!            1.0  ! 2002-07 (C. Ethe, G. Madec) re-writing F90 
    88   !!            3.0  ! 2006-07 (G. Madec) surface module 
    9    !!            3.3  ! 2009-05 (G.Garric, C. Bricaud) addition of the lim2_evp case 
     9   !!            3.3  ! 2009-05 (G. Garric, C. Bricaud) addition of the lim2_evp case 
     10   !!             -   ! 2010-11 (G. Madec) ice-ocean stress computed at each ocean time-step 
    1011   !!---------------------------------------------------------------------- 
    1112#if defined key_lim2 
     
    1314   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    1415   !!---------------------------------------------------------------------- 
    15    !!---------------------------------------------------------------------- 
    16    !!   lim_sbc_2  : flux at the ice / ocean interface 
     16   !!   lim_sbc_flx_2  : update mass, heat and salt fluxes at the ocean surface 
     17   !!   lim_sbc_tau_2  : update i- and j-stresses, and its modulus at the ocean surface 
    1718   !!---------------------------------------------------------------------- 
    1819   USE par_oce          ! ocean parameters 
     20   USE phycst           ! physical constants 
    1921   USE dom_oce          ! ocean domain 
     22   USE dom_ice_2        ! LIM-2: ice domain 
     23   USE ice_2            ! LIM-2: ice variables 
    2024   USE sbc_ice          ! surface boundary condition: ice 
    2125   USE sbc_oce          ! surface boundary condition: ocean 
    22    USE phycst           ! physical constants 
    23    USE ice_2            ! LIM2: ice variables 
    2426 
    2527   USE lbclnk           ! ocean lateral boundary condition - MPP exchanges 
    2628   USE in_out_manager   ! I/O manager 
    2729   USE diaar5, ONLY :   lk_diaar5 
    28    USE iom              !  
     30   USE iom              ! I/O library 
    2931   USE albedo           ! albedo parameters 
    3032   USE prtctl           ! Print control 
     
    3436   PRIVATE 
    3537 
    36    PUBLIC lim_sbc_2     ! called by sbc_ice_lim_2 
    37  
    38    REAL(wp)  ::   r1_rdtice                    ! constant values 
    39    REAL(wp)  ::   epsi16 = 1.e-16              !     -      - 
    40    REAL(wp)  ::   rzero  = 0.e0                !     -      - 
    41    REAL(wp)  ::   rone   = 1.e0                !     -      - 
     38   PUBLIC   lim_sbc_flx_2   ! called by sbc_ice_lim_2 
     39   PUBLIC   lim_sbc_tau_2   ! called by sbc_ice_lim_2 
     40 
     41   REAL(wp)  ::   r1_rdtice            ! = 1. / rdt_ice  
     42   REAL(wp)  ::   epsi16 = 1.e-16_wp   ! constant values 
     43   REAL(wp)  ::   rzero  = 0._wp       !     -      - 
     44   REAL(wp)  ::   rone   = 1._wp       !     -      - 
    4245   ! 
    43    REAL(wp), DIMENSION(jpi,jpj) ::   soce_r, sice_r   ! constant SSS and ice salinity used in levitating sea-ice case 
     46   REAL(wp), DIMENSION(jpi,jpj) ::   soce_0, sice_0   ! constant SSS and ice salinity used in levitating sea-ice case 
     47 
     48   REAL(wp), DIMENSION(jpi,jpj) ::   utau_oce, vtau_oce   ! air-ocean surface i- & j-stress              [N/m2] 
     49   REAL(wp), DIMENSION(jpi,jpj) ::   tmod_io              ! modulus of the ice-ocean relative velocity   [m/s] 
    4450 
    4551   !! * Substitutions 
     
    4854   !! NEMO/LIM2 3.3 , UCL - NEMO Consortium (2010) 
    4955   !! $Id$ 
    50    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    51    !!---------------------------------------------------------------------- 
    52  
     56   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     57   !!---------------------------------------------------------------------- 
    5358CONTAINS 
    5459 
    55    SUBROUTINE lim_sbc_2( kt ) 
     60   SUBROUTINE lim_sbc_flx_2( kt ) 
    5661      !!------------------------------------------------------------------- 
    5762      !!                ***  ROUTINE lim_sbc_2 *** 
     
    7782      !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
    7883      !!--------------------------------------------------------------------- 
    79       INTEGER ::   kt    ! number of iteration 
     84      INTEGER, INTENT(in) ::   kt    ! number of iteration 
    8085      !! 
    81       INTEGER  ::   ji, jj           ! dummy loop indices 
     86      INTEGER  ::   ji, jj   ! dummy loop indices 
    8287      INTEGER  ::   ii0, ii1, ij0, ij1         ! local integers 
    8388      INTEGER  ::   ifvt, i1mfr, idfr, iflt    !   -       - 
    8489      INTEGER  ::   ial, iadv, ifral, ifrdv    !   -       - 
    85       REAL(wp) ::   zqsr, zqns, zsang, zmod, zfm   ! local scalars 
    86       REAL(wp) ::   zinda, zfons, zemp, zztmp      !   -      - 
    87       REAL(wp) ::   zfrldu, zutau, zu_io           !   -      - 
    88       REAL(wp) ::   zfrldv, zvtau, zv_io           !   -      - 
    89       REAL(wp), DIMENSION(jpi,jpj)   ::   ztio_u, ztio_v    ! 2D workspace 
    90       REAL(wp), DIMENSION(jpi,jpj)   ::   ztiomi, zqnsoce   !  -     - 
     90      REAL(wp) ::   zqsr, zqns, zfm            ! local scalars 
     91      REAL(wp) ::   zinda, zfons, zemp         !   -      - 
     92      REAL(wp), DIMENSION(jpi,jpj)   ::   zqnsoce       ! 2D workspace 
    9193      REAL(wp), DIMENSION(jpi,jpj,1) ::   zalb, zalbp   ! 2D/3D workspace 
    92  
    9394      !!--------------------------------------------------------------------- 
    9495      
    95        
    9696      IF( kt == nit000 ) THEN 
    9797         IF(lwp) WRITE(numout,*) 
    98 #if defined key_lim2_vp 
    99          IF(lwp) WRITE(numout,*) 'lim_sbc_2 : LIM 2.0 sea-ice (VP rheology)  - surface boundary condition' 
    100 #else 
    101          IF(lwp) WRITE(numout,*) 'lim_sbc_2 : LIM 2.0 sea-ice (EVP rheology) - surface boundary condition' 
    102 #endif 
    103          IF(lwp) WRITE(numout,*) '~~~~~~~~~   ' 
     98         IF(lwp) WRITE(numout,*) 'lim_sbc_flx_2 : LIM-2 sea-ice - surface boundary condition - Mass, heat & salt fluxes' 
     99         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   ' 
    104100         ! 
    105101         r1_rdtice = 1. / rdt_ice 
    106102         ! 
    107          soce_r(:,:) = soce 
    108          sice_r(:,:) = sice 
    109          ! 
    110          IF( cp_cfg == "orca" ) THEN 
    111            !   ocean/ice salinity in the Baltic sea  
    112            DO jj = 1, jpj 
    113               DO ji = 1, jpi 
    114                  IF( glamt(ji,jj) >= 14. .AND.  glamt(ji,jj) <= 32. .AND. gphit(ji,jj) >= 54. .AND. gphit(ji,jj) <= 66. ) THEN  
    115                    soce_r(ji,jj) = 4.e0  
    116                    sice_r(ji,jj) = 2.e0 
    117                  END IF 
    118               END DO 
    119            END DO 
    120            ! 
    121          END IF 
    122       END IF 
     103         soce_0(:,:) = soce                     ! constant SSS and ice salinity used in levitating sea-ice case 
     104         sice_0(:,:) = sice 
     105         ! 
     106         IF( cp_cfg == "orca" ) THEN            ! decrease ocean & ice reference salinities in the Baltic sea  
     107            WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.   & 
     108               &   54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp         )  
     109               soce_0(:,:) = 4._wp 
     110               sice_0(:,:) = 2._wp 
     111            END WHERE 
     112         ENDIF 
     113         ! 
     114      ENDIF 
    123115 
    124116      !------------------------------------------! 
     
    126118      !------------------------------------------! 
    127119 
    128 !!gm 
    129 !!gm CAUTION    
    130 !!gm re-verifies the non solar expression, especially over open ocen 
    131 !!gm 
    132120      zqnsoce(:,:) = qns(:,:) 
    133121      DO jj = 1, jpj 
     
    190178 
    191179            fsbbq(ji,jj) = ( 1.0 - ( ifvt + iflt ) ) * fscmbq(ji,jj)     ! ??? 
    192              
     180            ! 
    193181            qsr  (ji,jj) = zqsr                                          ! solar heat flux  
    194182            qns  (ji,jj) = zqns - fdtcn(ji,jj)                           ! non solar heat flux 
     
    203191      !      mass flux at the ocean surface      ! 
    204192      !------------------------------------------! 
    205  
    206 !!gm 
    207 !!gm CAUTION    
    208 !!gm re-verifies the emp & emps expression, especially the absence of 1-frld on zfm 
    209 !!gm 
    210193      DO jj = 1, jpj 
    211194         DO ji = 1, jpi 
    212              
     195            ! 
    213196#if defined key_coupled 
    214197            ! freshwater exchanges at the ice-atmosphere / ocean interface (coupled mode) 
     
    223206               !                                                   !  ice-covered fraction: 
    224207#endif             
    225  
     208            ! 
    226209            !  computing salt exchanges at the ice/ocean interface 
    227             zfons =  ( soce_r(ji,jj) - sice_r(ji,jj) ) * ( rdmicif(ji,jj) * r1_rdtice )  
    228              
     210            zfons =  ( soce_0(ji,jj) - sice_0(ji,jj) ) * ( rdmicif(ji,jj) * r1_rdtice )  
     211            ! 
    229212            !  converting the salt flux from ice to a freshwater flux from ocean 
    230213            zfm  = zfons / ( sss_m(ji,jj) + epsi16 ) 
    231              
     214            ! 
    232215            emps(ji,jj) = zemp + zfm      ! surface ocean concentration/dilution effect (use on SSS evolution) 
    233216            emp (ji,jj) = zemp            ! surface ocean volume flux (use on sea-surface height evolution) 
    234  
     217            ! 
    235218         END DO 
    236219      END DO 
    237220 
    238       IF( lk_diaar5 ) THEN 
     221      IF( lk_diaar5 ) THEN       ! AR5 diagnostics 
    239222         CALL iom_put( 'isnwmlt_cea'  ,                 rdmsnif(:,:) * r1_rdtice ) 
    240          CALL iom_put( 'fsal_virt_cea',   soce_r(:,:) * rdmicif(:,:) * r1_rdtice ) 
    241          CALL iom_put( 'fsal_real_cea', - sice_r(:,:) * rdmicif(:,:) * r1_rdtice ) 
    242       ENDIF 
    243  
    244       !------------------------------------------! 
    245       !    momentum flux at the ocean surface    ! 
    246       !------------------------------------------! 
    247  
    248       IF ( ln_limdyn ) THEN                        ! Update the stress over ice-over area (only in ice-dynamic case) 
    249          !                                         ! otherwise the atmosphere-ocean stress is used everywhere 
    250  
    251          ! ... ice stress over ocean with a ice-ocean rotation angle (at I-point) 
    252 !CDIR NOVERRCHK 
    253          DO jj = 1, jpj 
    254 !CDIR NOVERRCHK 
    255             DO ji = 1, jpi 
    256                ! ... change the cosinus angle sign in the south hemisphere 
    257                zsang  = SIGN(1.e0, gphif(ji,jj) ) * sangvg 
    258                ! ... ice velocity relative to the ocean at I-point 
    259                zu_io  = u_ice(ji,jj) - u_oce(ji,jj) 
    260                zv_io  = v_ice(ji,jj) - v_oce(ji,jj) 
    261                zmod   = SQRT( zu_io * zu_io + zv_io * zv_io ) 
    262                zztmp  = rhoco * zmod 
    263                ! ... components of ice stress over ocean with a ice-ocean rotation angle (at I-point) 
    264                ztio_u(ji,jj) = zztmp * ( cangvg * zu_io - zsang * zv_io ) 
    265                ztio_v(ji,jj) = zztmp * ( cangvg * zv_io + zsang * zu_io ) 
    266                ! ... module of ice stress over ocean (at I-point) 
    267                ztiomi(ji,jj) = zztmp * zmod 
    268                !  
    269             END DO 
    270          END DO 
    271  
    272          DO jj = 2, jpjm1 
    273             DO ji = 2, jpim1   ! NO vector opt. 
    274                ! ... components of ice-ocean stress at U and V-points  (from I-point values) 
    275 #if defined key_lim2_vp 
    276                zutau  = 0.5 * ( ztio_u(ji+1,jj) + ztio_u(ji+1,jj+1) )      ! VP rheology 
    277                zvtau  = 0.5 * ( ztio_v(ji,jj+1) + ztio_v(ji+1,jj+1) ) 
    278 #else 
    279                zutau  = ztio_u(ji,jj)                                      ! EVP rheology 
    280                zvtau  = ztio_v(ji,jj) 
    281 #endif 
    282                ! ... open-ocean (lead) fraction at U- & V-points (from T-point values) 
    283                zfrldu = 0.5 * ( frld(ji,jj) + frld(ji+1,jj  ) ) 
    284                zfrldv = 0.5 * ( frld(ji,jj) + frld(ji  ,jj+1) ) 
    285                ! ... update components of surface ocean stress (ice-cover wheighted) 
    286                utau(ji,jj) = zfrldu * utau(ji,jj) + ( 1. - zfrldu ) * zutau 
    287                vtau(ji,jj) = zfrldv * vtau(ji,jj) + ( 1. - zfrldv ) * zvtau 
    288                ! ... module of ice-ocean stress at T-points (from I-point values) 
    289                zztmp = 0.25 * ( ztiomi(ji,jj) + ztiomi(ji+1,jj) + ztiomi(ji,jj+1) + ztiomi(ji+1,jj+1) ) 
    290                ! ... update module of surface ocean stress (ice-cover wheighted) 
    291                taum(ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1. - frld(ji,jj) ) * zztmp 
    292                ! 
    293             END DO 
    294          END DO 
    295  
    296          ! boundary condition on the stress (utau,vtau,taum) 
    297          CALL lbc_lnk( utau, 'U', -1. ) 
    298          CALL lbc_lnk( vtau, 'V', -1. ) 
    299          CALL lbc_lnk( taum, 'T',  1. ) 
    300  
     223         CALL iom_put( 'fsal_virt_cea',   soce_0(:,:) * rdmicif(:,:) * r1_rdtice ) 
     224         CALL iom_put( 'fsal_real_cea', - sice_0(:,:) * rdmicif(:,:) * r1_rdtice ) 
    301225      ENDIF 
    302226 
     
    305229      !-----------------------------------------------! 
    306230 
    307       IF ( lk_cpl ) THEN            
     231      IF( lk_cpl ) THEN          ! coupled case 
    308232         ! Ice surface temperature  
    309233         tn_ice(:,:,1) = sist(:,:)          ! sea-ice surface temperature        
     
    314238      ENDIF 
    315239 
    316       IF(ln_ctl) THEN 
     240      IF(ln_ctl) THEN            ! control print 
    317241         CALL prt_ctl(tab2d_1=qsr   , clinfo1=' lim_sbc: qsr    : ', tab2d_2=qns   , clinfo2=' qns     : ') 
    318242         CALL prt_ctl(tab2d_1=emp   , clinfo1=' lim_sbc: emp    : ', tab2d_2=emps  , clinfo2=' emps    : ') 
     
    321245         CALL prt_ctl(tab2d_1=fr_i  , clinfo1=' lim_sbc: fr_i   : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice  : ') 
    322246      ENDIF  
    323     
    324     END SUBROUTINE lim_sbc_2 
     247      ! 
     248   END SUBROUTINE lim_sbc_flx_2 
     249 
     250 
     251   SUBROUTINE lim_sbc_tau_2( kt , pu_oce, pv_oce ) 
     252      !!------------------------------------------------------------------- 
     253      !!                ***  ROUTINE lim_sbc_tau *** 
     254      !!   
     255      !! ** Purpose : Update the ocean surface stresses due to the ice 
     256      !!          
     257      !! ** Action  : * at each ice time step (every nn_fsbc time step): 
     258      !!                - compute the modulus of ice-ocean relative velocity  
     259      !!                  at T-point (C-grid) or I-point (B-grid) 
     260      !!                      tmod_io = rhoco * | U_ice-U_oce | 
     261      !!                - update the modulus of stress at ocean surface 
     262      !!                      taum = frld * taum + (1-frld) * tmod_io * | U_ice-U_oce | 
     263      !!              * at each ocean time step (each kt):  
     264      !!                  compute linearized ice-ocean stresses as 
     265      !!                      Utau = tmod_io * | U_ice - pU_oce | 
     266      !!                using instantaneous current ocean velocity (usually before) 
     267      !! 
     268      !!    NB: - the averaging operator used depends on the ice dynamics grid (cp_ice_msh='I' or 'C') 
     269      !!        - ice-ocean rotation angle only allowed in cp_ice_msh='I' case 
     270      !!        - here we make an approximation: taum is only computed every ice time step 
     271      !!          This avoids mutiple average to pass from T -> U,V grids and next from U,V grids  
     272      !!          to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... 
     273      !! 
     274      !! ** Outputs : - utau, vtau   : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes 
     275      !!              - taum         : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes 
     276      !!--------------------------------------------------------------------- 
     277      INTEGER ,                     INTENT(in) ::   kt               ! ocean time-step index 
     278      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pu_oce, pv_oce   ! surface ocean currents 
     279      !! 
     280      INTEGER  ::   ji, jj   ! dummy loop indices 
     281      REAL(wp) ::   zfrldu, zat_u, zu_i, zutau_ice, zu_t, zmodt   ! local scalar 
     282      REAL(wp) ::   zfrldv, zat_v, zv_i, zvtau_ice, zv_t, zmodi   !   -      - 
     283      REAL(wp) ::   zsang, zumt                                          !    -         - 
     284      REAL(wp), DIMENSION(jpi,jpj) ::   ztio_u, ztio_v   ! ocean stress below sea-ice 
     285      !!--------------------------------------------------------------------- 
     286      ! 
     287      IF( kt == nit000 .AND. lwp ) THEN         ! control print 
     288         WRITE(numout,*) 
     289         WRITE(numout,*) 'lim_sbc_tau_2 : LIM 2.0 sea-ice - surface ocean momentum fluxes' 
     290         WRITE(numout,*) '~~~~~~~~~~~~~ ' 
     291         IF( lk_lim2_vp )   THEN   ;   WRITE(numout,*) '                VP  rheology - B-grid case' 
     292         ELSE                      ;   WRITE(numout,*) '                EVP rheology - C-grid case' 
     293         ENDIF 
     294      ENDIF 
     295      ! 
     296      SELECT CASE( cp_ice_msh )      
     297      !                             !-----------------------! 
     298      CASE( 'I' )                   !  B-grid ice dynamics  !   I-point (i.e. F-point with sea-ice indexation) 
     299         !                          !--=--------------------! 
     300         ! 
     301         IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==! (i.e. surface module time-step) 
     302!CDIR NOVERRCHK 
     303            DO jj = 1, jpj                               !* modulus of ice-ocean relative velocity at I-point 
     304!CDIR NOVERRCHK 
     305               DO ji = 1, jpi 
     306                  zu_i  = u_ice(ji,jj) - u_oce(ji,jj)                   ! ice-ocean relative velocity at I-point 
     307                  zv_i  = v_ice(ji,jj) - v_oce(ji,jj) 
     308                  tmod_io(ji,jj) = SQRT( zu_i * zu_i + zv_i * zv_i )    ! modulus of this velocity (at I-point) 
     309               END DO 
     310            END DO 
     311!CDIR NOVERRCHK 
     312            DO jj = 1, jpjm1                             !* update the modulus of stress at ocean surface (T-point) 
     313!CDIR NOVERRCHK 
     314               DO ji = 1, jpim1   ! NO vector opt. 
     315                  !                                               ! modulus of U_ice-U_oce at T-point 
     316                  zumt  = 0.25_wp * (  tmod_io(ji+1,jj) + tmod_io(ji+1,jj+1)    &    
     317                     &               + tmod_io(ji,jj+1) + tmod_io(ji+1,jj+1)  ) 
     318                  !                                               ! update the modulus of stress at ocean surface 
     319                  taum(ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1._wp - frld(ji,jj) ) * rhoco * zumt * zumt 
     320               END DO 
     321            END DO 
     322            CALL lbc_lnk( taum, 'T', 1. ) 
     323            ! 
     324            utau_oce(:,:) = utau(:,:)                    !* save the air-ocean stresses at ice time-step 
     325            vtau_oce(:,:) = vtau(:,:) 
     326            ! 
     327         ENDIF 
     328         ! 
     329         !                                        !==  at each ocean time-step  ==! 
     330         ! 
     331         !                                               !* ice/ocean stress WITH a ice-ocean rotation angle at I-point 
     332         DO jj = 2, jpj 
     333            zsang  = SIGN( 1._wp, gphif(1,jj) ) * sangvg          ! change the cosine angle sign in the SH  
     334            DO ji = 2, jpi    ! NO vect. opt. possible 
     335               ! ... ice-ocean relative velocity at I-point using instantaneous surface ocean current at u- & v-pts 
     336               zu_i = u_ice(ji,jj) - 0.5_wp * ( pu_oce(ji-1,jj  ) + pu_oce(ji-1,jj-1) ) * tmu(ji,jj) 
     337               zv_i = v_ice(ji,jj) - 0.5_wp * ( pv_oce(ji  ,jj-1) + pv_oce(ji-1,jj-1) ) * tmu(ji,jj) 
     338               ! ... components of stress with a ice-ocean rotation angle  
     339               zmodi = rhoco * tmod_io(ji,jj)                      
     340               ztio_u(ji,jj) = zmodi * ( cangvg * zu_i - zsang * zv_i ) 
     341               ztio_v(ji,jj) = zmodi * ( cangvg * zv_i + zsang * zu_i ) 
     342            END DO 
     343         END DO 
     344         !                                               !* surface ocean stresses at u- and v-points 
     345         DO jj = 2, jpjm1 
     346            DO ji = 2, jpim1   ! NO vector opt. 
     347               !                                   ! ice-ocean stress at U and V-points  (from I-point values) 
     348               zutau_ice  = 0.5_wp * ( ztio_u(ji+1,jj) + ztio_u(ji+1,jj+1) ) 
     349               zvtau_ice  = 0.5_wp * ( ztio_v(ji,jj+1) + ztio_v(ji+1,jj+1) ) 
     350               !                                   ! open-ocean (lead) fraction at U- & V-points (from T-point values) 
     351               zfrldu = 0.5_wp * ( frld(ji,jj) + frld(ji+1,jj  ) ) 
     352               zfrldv = 0.5_wp * ( frld(ji,jj) + frld(ji  ,jj+1) ) 
     353               !                                   ! update the surface ocean stress (ice-cover wheighted) 
     354               utau(ji,jj) = zfrldu * utau_oce(ji,jj) + ( 1._wp - zfrldu ) * zutau_ice 
     355               vtau(ji,jj) = zfrldv * vtau_oce(ji,jj) + ( 1._wp - zfrldv ) * zvtau_ice 
     356            END DO 
     357         END DO 
     358         CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )     ! lateral boundary condition 
     359         ! 
     360         ! 
     361         !                          !-----------------------! 
     362      CASE( 'C' )                   !  C-grid ice dynamics  !   U & V-points (same as in the ocean) 
     363         !                          !--=--------------------! 
     364         ! 
     365         IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==! (i.e. surface module time-step) 
     366!CDIR NOVERRCHK 
     367            DO jj = 2, jpjm1                          !* modulus of the ice-ocean velocity at T-point 
     368!CDIR NOVERRCHK 
     369               DO ji = fs_2, fs_jpim1 
     370                  zu_t  = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj)   ! 2*(U_ice-U_oce) at T-point 
     371                  zv_t  = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1)       
     372                  zmodt =  0.25_wp * (  zu_t * zu_t + zv_t * zv_t  )                      ! |U_ice-U_oce|^2 
     373                  !                                               ! update the modulus of stress at ocean surface 
     374                  taum   (ji,jj) = frld(ji,jj) * taum(ji,jj) + ( 1._wp - frld(ji,jj) ) * rhoco * zmodt 
     375                  tmod_io(ji,jj) = SQRT( zmodt ) * rhoco          ! rhoco*|Uice-Uoce| 
     376               END DO 
     377            END DO 
     378            CALL lbc_lnk( taum, 'T', 1. )   ;   CALL lbc_lnk( tmod_io, 'T', 1. ) 
     379            ! 
     380            utau_oce(:,:) = utau(:,:)                 !* save the air-ocean stresses at ice time-step 
     381            vtau_oce(:,:) = vtau(:,:) 
     382            ! 
     383         ENDIF 
     384         ! 
     385         !                                        !==  at each ocean time-step  ==! 
     386         ! 
     387         DO jj = 2, jpjm1                             !* ice stress over ocean WITHOUT a ice-ocean rotation angle 
     388            DO ji = fs_2, fs_jpim1 
     389               !                                            ! ocean area at u- & v-points 
     390               zfrldu  = 0.5_wp * ( frld(ji,jj) + frld(ji+1,jj  ) ) 
     391               zfrldv  = 0.5_wp * ( frld(ji,jj) + frld(ji  ,jj+1) ) 
     392               !                                            ! quadratic drag formulation without rotation 
     393               !                                            ! using instantaneous surface ocean current 
     394               zutau_ice = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) ) 
     395               zvtau_ice = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) ) 
     396               !                                            ! update the surface ocean stress (ice-cover wheighted) 
     397               utau(ji,jj) = zfrldu * utau_oce(ji,jj) + ( 1._wp - zfrldu ) * zutau_ice 
     398               vtau(ji,jj) = zfrldv * vtau_oce(ji,jj) + ( 1._wp - zfrldv ) * zvtau_ice 
     399            END DO 
     400         END DO 
     401         CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
     402         ! 
     403      END SELECT 
     404 
     405      IF(ln_ctl)   CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
     406         &                       tab2d_2=vtau, clinfo2=' vtau    : '        , mask2=vmask ) 
     407      !   
     408   END SUBROUTINE lim_sbc_tau_2 
    325409 
    326410#else 
    327411   !!---------------------------------------------------------------------- 
    328    !!   Default option :        Dummy module       NO LIM 2.0 sea-ice model 
    329    !!---------------------------------------------------------------------- 
    330 CONTAINS 
    331    SUBROUTINE lim_sbc_2         ! Dummy routine 
    332    END SUBROUTINE lim_sbc_2 
     412   !!   Default option         Empty module        NO LIM 2.0 sea-ice model 
     413   !!---------------------------------------------------------------------- 
    333414#endif  
    334415 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_2/limtrp_2.F90

    r2319 r2370  
    77   !!            2.0  !  2001-05 (G. Madec, R. Hordoir) opa norm 
    88   !!             -   !  2004-01 (G. Madec, C. Ethe)  F90, mpp 
    9    !!            3.3  !  2009-05  (G.Garric, C. Bricaud) addition of the lim2_evp case 
     9   !!            3.3  !  2009-05  (G. Garric, C. Bricaud) addition of the lim2_evp case 
    1010   !!---------------------------------------------------------------------- 
    1111#if defined key_lim2 
     
    4141   REAL(wp)  ::   rone   = 1.e0 
    4242 
    43  
    4443   !! * Substitution 
    4544#  include "vectopt_loop_substitute.h90" 
     
    8887         ! ice velocities at ocean U- and V-points (zui_u,zvi_v) 
    8988         ! --------------------------------------- 
    90          zvbord = 1.0 + ( 1.0 - bound )      ! zvbord=2 no-slip, =0 free slip boundary conditions         
    91 #if defined key_lim2_vp 
    92          DO jj = 1, jpjm1 
    93             DO ji = 1, jpim1   ! NO vector opt. 
    94                zui_u(ji,jj) = ( u_ice(ji+1,jj  ) + u_ice(ji+1,jj+1) ) / ( MAX( tmu(ji+1,jj  ) + tmu(ji+1,jj+1), zvbord ) ) 
    95                zvi_v(ji,jj) = ( v_ice(ji  ,jj+1) + v_ice(ji+1,jj+1) ) / ( MAX( tmu(ji  ,jj+1) + tmu(ji+1,jj+1), zvbord ) ) 
    96             END DO 
    97          END DO 
    98          CALL lbc_lnk( zui_u, 'U', -1. )   ;   CALL lbc_lnk( zvi_v, 'V', -1. )         ! Lateral boundary conditions 
    99 #else 
    100         zui_u(:,:) = u_ice(:,:)      ! EVP rheology: ice (u,v) at u- and v-points 
    101         zvi_v(:,:) = v_ice(:,:) 
    102 #endif 
     89         IF( lk_lim2_vp ) THEN      ! VP rheology : B-grid sea-ice dynamics (I-point ice velocity) 
     90            zvbord = 1._wp + ( 1._wp - bound )      ! zvbord=2 no-slip, =0 free slip boundary conditions         
     91            DO jj = 1, jpjm1 
     92               DO ji = 1, jpim1   ! NO vector opt. 
     93                  zui_u(ji,jj) = ( u_ice(ji+1,jj) + u_ice(ji+1,jj+1) ) / ( MAX( tmu(ji+1,jj)+tmu(ji+1,jj+1), zvbord ) ) 
     94                  zvi_v(ji,jj) = ( v_ice(ji,jj+1) + v_ice(ji+1,jj+1) ) / ( MAX( tmu(ji,jj+1)+tmu(ji+1,jj+1), zvbord ) ) 
     95               END DO 
     96            END DO 
     97            CALL lbc_lnk( zui_u, 'U', -1. )   ;   CALL lbc_lnk( zvi_v, 'V', -1. )      ! Lateral boundary conditions 
     98            ! 
     99         ELSE                       ! EVP rheology : C-grid sea-ice dynamics (u- & v-points ice velocity) 
     100            zui_u(:,:) = u_ice(:,:)      ! EVP rheology: ice (u,v) at u- and v-points 
     101            zvi_v(:,:) = v_ice(:,:) 
     102         ENDIF 
    103103 
    104104         ! CFL test for stability 
    105105         ! ---------------------- 
    106          zcfl  = 0.e0 
     106         zcfl  = 0._wp 
    107107         zcfl  = MAX( zcfl, MAXVAL( ABS( zui_u(1:jpim1, :     ) ) * rdt_ice / e1u(1:jpim1, :     ) ) ) 
    108108         zcfl  = MAX( zcfl, MAXVAL( ABS( zvi_v( :     ,1:jpjm1) ) * rdt_ice / e2v( :     ,1:jpjm1) ) ) 
     
    193193!        DO jj = 1, jpjm1          ! NB: has not to be defined on jpj line and jpi row 
    194194!           DO ji = 1 , fs_jpim1   ! vector opt. 
    195 !              IF( MIN( zs0a(ji,jj) , zs0a(ji+1,jj) ) == 0.e0 )   pahu(ji,jj) = 0.e0 
    196 !              IF( MIN( zs0a(ji,jj) , zs0a(ji,jj+1) ) == 0.e0 )   pahv(ji,jj) = 0.e0 
     195!              IF( MIN( zs0a(ji,jj) , zs0a(ji+1,jj) ) == 0.e0 )   pahu(ji,jj) = 0._wp 
     196!              IF( MIN( zs0a(ji,jj) , zs0a(ji,jj+1) ) == 0.e0 )   pahv(ji,jj) = 0._wp 
    197197!           END DO 
    198198!        END DO 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90

    r2287 r2370  
    44   !!   Sea-Ice dynamics :   
    55   !!====================================================================== 
     6   !! history :  1.0  ! 2002-08 (C. Ethe, G. Madec)  original VP code  
     7   !!            3.0  ! 2007-03 (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle)  LIM3: EVP-Cgrid 
     8   !!---------------------------------------------------------------------- 
    69#if defined key_lim3 
    710   !!---------------------------------------------------------------------- 
     
    1114   !!    lim_dyn_init : initialization and namelist read 
    1215   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    1416   USE phycst 
    1517   USE in_out_manager  ! I/O manager 
     
    3032   PRIVATE 
    3133 
    32    !! * Accessibility 
    33    PUBLIC lim_dyn  ! routine called by ice_step 
     34   PUBLIC   lim_dyn   ! routine called by ice_step 
    3435 
    3536   !! * Substitutions 
    3637#  include "vectopt_loop_substitute.h90" 
    37  
    38    !! * Module variables 
    39    REAL(wp)  ::  rone    = 1.e0   ! constant value 
    40  
    4138   !!---------------------------------------------------------------------- 
    4239   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4340   !! $Id$ 
    44    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    45    !!---------------------------------------------------------------------- 
    46  
     41   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     42   !!---------------------------------------------------------------------- 
    4743CONTAINS 
    4844 
     
    5955      !!              - computation of the stress at the ocean surface          
    6056      !!              - treatment of the case if no ice dynamic 
    61       !! History : 
    62       !!   1.0  !  01-04   (LIM)  Original code 
    63       !!   2.0  !  02-08   (C. Ethe, G. Madec)  F90, mpp 
    64       !!   3.0  !  2007-03 (M.A. Morales Maqueda, S. Bouillon, M. Vancoppenolle)  
    65       !!                   LIM3, EVP, C-grid 
    6657      !!------------------------------------------------------------------------------------ 
    6758      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    68       !! * Local variables 
     59      !! 
    6960      INTEGER  ::   ji, jj, jl, ja    ! dummy loop indices 
    7061      INTEGER  ::   i_j1, i_jpj       ! Starting/ending j-indices for rheology 
    71       REAL(wp) ::   zcoef             ! temporary scalar 
     62      REAL(wp) ::   zcoef             ! local scalar 
    7263      REAL(wp), DIMENSION(jpj)     ::   zind           ! i-averaged indicator of sea-ice 
    7364      REAL(wp), DIMENSION(jpj)     ::   zmsk           ! i-averaged of tmask 
     
    158149         zv_io(:,:) = v_ice(:,:) - ssv_m(:,:) 
    159150         ! frictional velocity at T-point 
     151         zcoef = 0.5 * cw 
    160152         DO jj = 2, jpjm1  
    161153            DO ji = fs_2, fs_jpim1   ! vector opt. 
    162                ust2s(ji,jj) = 0.5 * cw                                                          & 
    163                   &         * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
    164                   &            + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
     154               ust2s(ji,jj) = zcoef * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
     155                  &                    + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
    165156            END DO 
    166157         END DO 
     
    171162         DO jj = 2, jpjm1 
    172163            DO ji = fs_2, fs_jpim1   ! vector opt. 
    173                ust2s(ji,jj) = zcoef * tms(ji,jj) * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
    174                   &                                     + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) 
     164               ust2s(ji,jj) = zcoef * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
     165                  &                        + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1)   ) * tms(ji,jj) 
    175166            END DO 
    176167         END DO 
    177168         ! 
    178169      ENDIF 
    179  
    180170      CALL lbc_lnk( ust2s, 'T',  1. )   ! T-point 
    181171 
     
    218208         END DO 
    219209      ENDIF 
    220  
     210      ! 
    221211   END SUBROUTINE lim_dyn 
     212 
    222213 
    223214   SUBROUTINE lim_dyn_init 
     
    232223      !! 
    233224      !! ** input   :   Namelist namicedyn 
    234       !! 
    235       !! history : 
    236       !!  8.5  ! 03-08 (C. Ethe) original code 
    237       !!  9.0  ! 07-03 (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle) 
    238       !!               EVP-Cgrid-LIM3 
    239225      !!------------------------------------------------------------------- 
    240226      NAMELIST/namicedyn/ epsd, alpha,     & 
     
    244230      !!------------------------------------------------------------------- 
    245231 
    246       ! Define the initial parameters 
    247       ! ------------------------- 
    248  
    249       ! Read Namelist namicedyn 
    250       REWIND ( numnam_ice ) 
    251       READ   ( numnam_ice  , namicedyn ) 
    252       IF(lwp) THEN 
     232      REWIND( numnam_ice )                ! Read Namelist namicedyn 
     233      READ  ( numnam_ice  , namicedyn ) 
     234       
     235      IF(lwp) THEN                        ! control print 
    253236         WRITE(numout,*) 
    254237         WRITE(numout,*) 'lim_dyn_init : ice parameters for ice dynamics ' 
     
    272255         WRITE(numout,*) '   timescale for elastic waves                      telast = ', telast 
    273256         WRITE(numout,*) '   coefficient for the solution of int. stresses  alphaevp = ', alphaevp 
    274  
    275       ENDIF 
    276  
    277       usecc2 = 1.0 / ( ecc * ecc ) 
    278       rhoco  = rau0 * cw 
     257      ENDIF 
     258      ! 
     259      IF( angvg /= 0._wp ) THEN 
     260         CALL ctl_warn( 'lim_dyn_init: turning angle for oceanic stress not properly coded for EVP ',   & 
     261            &           '(see limsbc module). We force  angvg = 0._wp'  ) 
     262         angvg = 0._wp 
     263      ENDIF 
     264       
     265      usecc2 = 1._wp / ( ecc * ecc ) 
     266      rhoco  = rau0  * cw 
    279267      angvg  = angvg * rad 
    280268      sangvg = SIN( angvg ) 
    281269      cangvg = COS( angvg ) 
    282       pstarh = pstar / 2.0 
     270      pstarh = pstar * 0.5_wp 
    283271 
    284272      !  Diffusion coefficients. 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

    r2357 r2370  
    1111#if defined key_lim3 || (  defined key_lim2 && ! defined key_lim2_vp ) 
    1212   !!---------------------------------------------------------------------- 
    13    !!   'key_lim3'                                      LIM3 sea-ice model 
     13   !!   'key_lim3'               OR                     LIM-3 sea-ice model 
     14   !!   'key_lim2' AND NOT 'key_lim2_vp'             VP LIM-2 sea-ice model 
    1415   !!---------------------------------------------------------------------- 
    1516   !!   lim_rhg   : computes ice velocities 
    1617   !!---------------------------------------------------------------------- 
    17    !! * Modules used 
    18    USE phycst 
    19    USE par_oce 
    20    USE dom_oce 
    21    USE sbc_oce         ! Surface boundary condition: ocean fields 
    22    USE sbc_ice         ! Surface boundary condition: ice fields 
    23    USE lbclnk 
    24    USE lib_mpp 
    25    USE in_out_manager  ! I/O manager 
    26    USE prtctl          ! Print control 
     18   USE phycst           ! Physical constant 
     19   USE par_oce          ! Ocean parameters 
     20   USE dom_oce          ! Ocean domain 
     21   USE sbc_oce          ! Surface boundary condition: ocean fields 
     22   USE sbc_ice          ! Surface boundary condition: ice fields 
     23   USE lbclnk           ! Lateral Boundary Condition / MPP link 
     24   USE lib_mpp          ! MPP library 
     25   USE in_out_manager   ! I/O manager 
     26   USE prtctl           ! Print control 
    2727#if defined key_lim3 
    28    USE ice 
    29    USE dom_ice 
    30    USE iceini 
    31    USE limitd_me 
    32 #endif 
    33 #if defined key_lim2 && ! defined key_lim2_vp 
     28   USE ice              ! LIM-3: ice variables 
     29   USE dom_ice          ! LIM-3: ice domain 
     30   USE iceini           ! LIM-3: ice initialisation 
     31   USE limitd_me        ! LIM-3:  
     32#else 
    3433   USE ice_2            ! LIM2: ice variables 
    3534   USE dom_ice_2        ! LIM2: ice domain 
     
    3736#endif 
    3837 
    39  
    4038   IMPLICIT NONE 
    4139   PRIVATE 
    4240 
    43    !! * Routine accessibility 
    44    PUBLIC lim_rhg  ! routine called by lim_dyn 
    45  
     41   PUBLIC   lim_rhg   ! routine called by lim_dyn (or lim_dyn_2) 
     42 
     43   REAL(wp) ::   rzero   = 0._wp   ! constant values 
     44   REAL(wp) ::   rone    = 1._wp   ! constant values 
     45       
    4646   !! * Substitutions 
    4747#  include "vectopt_loop_substitute.h90" 
    48  
    49    !! * Module variables 
    50    REAL(wp)  ::           &  ! constant values 
    51       rzero   = 0.e0   ,  & 
    52       rone    = 1.e0 
    5348   !!---------------------------------------------------------------------- 
    5449   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    5550   !! $Id$ 
    56    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     51   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    5752   !!---------------------------------------------------------------------- 
    58  
    5953CONTAINS 
    6054 
    6155   SUBROUTINE lim_rhg( k_j1, k_jpj ) 
    62  
    6356      !!------------------------------------------------------------------- 
    6457      !!                 ***  SUBROUTINE lim_rhg  *** 
     
    107100      !!                 e.g. in the Canadian Archipelago 
    108101      !! 
    109       !! ** References : Hunke and Dukowicz, JPO97 
    110       !!                 Bouillon et al., 08, in prep (update this when 
    111       !!                 published) 
    112       !!                 Vancoppenolle et al., OM08 
     102      !! References : Hunke and Dukowicz, JPO97 
     103      !!              Bouillon et al., Ocean Modelling 2009 
     104      !!              Vancoppenolle et al., Ocean Modelling 2008 
     105      !!------------------------------------------------------------------- 
     106      INTEGER, INTENT(in) ::   k_j1    ! southern j-index for ice computation 
     107      INTEGER, INTENT(in) ::   k_jpj   ! northern j-index for ice computation 
    113108      !! 
    114       !!------------------------------------------------------------------- 
    115       ! * Arguments 
    116       ! 
    117       INTEGER, INTENT(in) :: & 
    118          k_j1 ,                      & !: southern j-index for ice computation 
    119          k_jpj                         !: northern j-index for ice computation 
    120  
    121       ! * Local variables 
    122       INTEGER ::   ji, jj              !: dummy loop indices 
    123  
    124       INTEGER  :: & 
    125          jter                          !: temporary integers 
    126  
     109      INTEGER ::   ji, jj   ! dummy loop indices 
     110      INTEGER ::   jter     ! local integers 
    127111      CHARACTER (len=50) ::   charout 
    128  
    129       REAL(wp) :: & 
    130          zt11, zt12, zt21, zt22,     & !: temporary scalars 
    131          ztagnx, ztagny,             & !: wind stress on U/V points                        
    132          delta                         ! 
    133  
    134       REAL(wp) :: & 
    135          za,                         & !: 
    136          zstms,                      & !: temporary scalar for ice strength 
    137          zsang,                      & !: temporary scalar for coriolis term 
    138          zmask                         !: mask for the computation of ice mass 
     112      REAL(wp) ::   zt11, zt12, zt21, zt22, ztagnx, ztagny, delta                         ! 
     113      REAL(wp) ::   za, zstms, zsang, zmask   ! local scalars 
    139114 
    140115      REAL(wp),DIMENSION(jpi,jpj) :: & 
     
    154129 
    155130      REAL(wp) :: & 
    156          dtevp,                      & !: time step for subcycling 
    157          dtotel,                     & !: 
    158          ecc2,                       & !: square of yield ellipse eccenticity 
    159          z0,                         & !: temporary scalar 
    160          zr,                         & !: temporary scalar 
    161          zcca, zccb,                 & !: temporary scalars 
    162          zu_ice2,                    & !:  
    163          zv_ice1,                    & !: 
    164          zddc, zdtc,                 & !: temporary array for delta on corners 
    165          zdst,                       & !: temporary array for delta on centre 
    166          zdsshx, zdsshy,             & !: term for the gradient of ocean surface 
    167          sigma1, sigma2                !: internal ice stress 
     131         dtevp,                      & ! time step for subcycling 
     132         dtotel,                     & ! 
     133         ecc2,                       & ! square of yield ellipse eccenticity 
     134         z0,                         & ! temporary scalar 
     135         zr,                         & ! temporary scalar 
     136         zcca, zccb,                 & ! temporary scalars 
     137         zu_ice2,                    & !  
     138         zv_ice1,                    & ! 
     139         zddc, zdtc,                 & ! temporary array for delta on corners 
     140         zdst,                       & ! temporary array for delta on centre 
     141         zdsshx, zdsshy,             & ! term for the gradient of ocean surface 
     142         sigma1, sigma2                ! internal ice stress 
     143 
     144      REAL(wp),DIMENSION(jpi,jpj) ::   zf1, zf2   ! arrays for internal stresses 
    168145 
    169146      REAL(wp),DIMENSION(jpi,jpj) :: & 
    170          zf1, zf2                      !: arrays for internal stresses 
    171  
    172       REAL(wp),DIMENSION(jpi,jpj) :: & 
    173          zdd, zdt,                   & !: Divergence and tension at centre of grid cells 
    174          zds,                        & !: Shear on northeast corner of grid cells 
    175          deltat,                     & !: Delta at centre of grid cells 
    176          deltac,                     & !: Delta on corners 
    177          zs1, zs2,                   & !: Diagonal stress tensor components zs1 and zs2  
    178          zs12                          !: Non-diagonal stress tensor component zs12 
     147         zdd, zdt,                   & ! Divergence and tension at centre of grid cells 
     148         zds,                        & ! Shear on northeast corner of grid cells 
     149         deltat,                     & ! Delta at centre of grid cells 
     150         deltac,                     & ! Delta on corners 
     151         zs1, zs2,                   & ! Diagonal stress tensor components zs1 and zs2  
     152         zs12                          ! Non-diagonal stress tensor component zs12 
    179153 
    180154      REAL(wp) :: & 
    181          zresm            ,          & !: Maximal error on ice velocity 
    182          zindb            ,          & !: ice (1) or not (0)       
    183          zdummy                        !: dummy argument 
    184  
    185       REAL(wp),DIMENSION(jpi,jpj) :: & 
    186          zu_ice           ,          & !: Ice velocity on previous time step 
    187          zv_ice           ,          & 
    188          zresr                         !: Local error on velocity 
     155         zresm            ,          & ! Maximal error on ice velocity 
     156         zindb            ,          & ! ice (1) or not (0)       
     157         zdummy                        ! dummy argument 
     158 
     159      REAL(wp),DIMENSION(jpi,jpj) ::   zu_ice, zv_ice, zresr   ! Local error on velocity 
    189160      !!------------------------------------------------------------------- 
    190  
    191161#if  defined key_lim2 && ! defined key_lim2_vp 
    192162# if defined key_agrif 
     
    205175      ! 
    206176      ! Put every vector to 0 
    207       zpresh(:,:) = 0.0 ; zc1(:,:) = 0.0 
    208       zpreshc(:,:) = 0.0 
    209       u_ice2(:,:)  = 0.0 ; v_ice1(:,:)  = 0.0 
    210       zdd(:,:)     = 0.0 ; zdt(:,:)     = 0.0 ; zds(:,:)     = 0.0 
     177      zpresh (:,:) = 0._wp   ;   zc1   (:,:) = 0._wp 
     178      zpreshc(:,:) = 0._wp 
     179      u_ice2 (:,:) = 0._wp   ;   v_ice1(:,:) = 0._wp 
     180      zdd    (:,:) = 0._wp   ;   zdt   (:,:) = 0._wp   ;   zds(:,:) = 0._wp 
    211181 
    212182#if defined key_lim3 
    213       ! Ice strength on T-points 
    214       CALL lim_itd_me_icestrength(ridge_scheme_swi) 
     183      CALL lim_itd_me_icestrength( ridge_scheme_swi )      ! LIM-3: Ice strength on T-points 
    215184#endif 
    216185 
    217       ! Ice mass and temp variables 
    218 !CDIR NOVERRCHK 
    219       DO jj = k_j1 , k_jpj 
     186!CDIR NOVERRCHK 
     187      DO jj = k_j1 , k_jpj       ! Ice mass and temp variables 
    220188!CDIR NOVERRCHK 
    221189         DO ji = 1 , jpi 
    222190            zc1(ji,jj)    = tms(ji,jj) * ( rhosn * vt_s(ji,jj) + rhoic * vt_i(ji,jj) ) 
    223191#if defined key_lim3 
    224             zpresh(ji,jj) = tms(ji,jj) *  strength(ji,jj) / 2. 
     192            zpresh(ji,jj) = tms(ji,jj) *  strength(ji,jj) * 0.5_wp 
    225193#endif 
    226194            ! tmi = 1 where there is ice or on land 
    227             tmi(ji,jj)    = 1.0 - ( 1.0 - MAX( 0.0 , SIGN ( 1.0 , vt_i(ji,jj) - & 
    228                epsd ) ) ) * tms(ji,jj) 
     195            tmi(ji,jj)    = 1._wp - ( 1._wp - MAX( 0._wp , SIGN ( 1._wp , vt_i(ji,jj) - epsd ) ) ) * tms(ji,jj) 
    229196         END DO 
    230197      END DO 
     
    271238         DO ji = fs_2, fs_jpim1 
    272239 
    273             zt11 = tms(ji,jj)*e1t(ji,jj) 
    274             zt12 = tms(ji+1,jj)*e1t(ji+1,jj) 
    275             zt21 = tms(ji,jj)*e2t(ji,jj) 
    276             zt22 = tms(ji,jj+1)*e2t(ji,jj+1) 
     240            zt11 = tms(ji  ,jj) * e1t(ji  ,jj) 
     241            zt12 = tms(ji+1,jj) * e1t(ji+1,jj) 
     242            zt21 = tms(ji,jj  ) * e2t(ji,jj  ) 
     243            zt22 = tms(ji,jj+1) * e2t(ji,jj+1) 
    277244 
    278245            ! Leads area. 
    279             zfrld1(ji,jj) = ( zt12 * ( 1.0 - at_i(ji,jj) ) + & 
    280                &                        zt11 * ( 1.0 - at_i(ji+1,jj) ) ) / ( zt11 + zt12 + epsd ) 
    281             zfrld2(ji,jj) = ( zt22 * ( 1.0 - at_i(ji,jj) ) + & 
    282                &                        zt21 * ( 1.0 - at_i(ji,jj+1) ) ) / ( zt21 + zt22 + epsd ) 
     246            zfrld1(ji,jj) = ( zt12 * ( 1.0 - at_i(ji,jj) ) + zt11 * ( 1.0 - at_i(ji+1,jj) ) ) / ( zt11 + zt12 + epsd ) 
     247            zfrld2(ji,jj) = ( zt22 * ( 1.0 - at_i(ji,jj) ) + zt21 * ( 1.0 - at_i(ji,jj+1) ) ) / ( zt21 + zt22 + epsd ) 
    283248 
    284249            ! Mass, coriolis coeff. and currents 
    285250            zmass1(ji,jj) = ( zt12*zc1(ji,jj) + zt11*zc1(ji+1,jj) ) / (zt11+zt12+epsd) 
    286251            zmass2(ji,jj) = ( zt22*zc1(ji,jj) + zt21*zc1(ji,jj+1) ) / (zt21+zt22+epsd) 
    287             zcorl1(ji,jj) = zmass1(ji,jj) * ( e1t(ji+1,jj)*fcor(ji,jj) + & 
    288                e1t(ji,jj)*fcor(ji+1,jj) ) & 
    289                / (e1t(ji,jj) + e1t(ji+1,jj) + epsd ) 
    290             zcorl2(ji,jj) = zmass2(ji,jj) * ( e2t(ji,jj+1)*fcor(ji,jj) + & 
    291                e2t(ji,jj)*fcor(ji,jj+1) ) & 
    292                / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 
     252            zcorl1(ji,jj) = zmass1(ji,jj) * ( e1t(ji+1,jj)*fcor(ji,jj) + e1t(ji,jj)*fcor(ji+1,jj) )   & 
     253               &                          / ( e1t(ji,jj) + e1t(ji+1,jj) + epsd ) 
     254            zcorl2(ji,jj) = zmass2(ji,jj) * ( e2t(ji,jj+1)*fcor(ji,jj) + e2t(ji,jj)*fcor(ji,jj+1) )   & 
     255               &                          / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 
    293256            ! 
    294257            u_oce1(ji,jj)  = u_oce(ji,jj) 
     
    314277            ! include it later 
    315278 
    316             zdsshx =  (ssh_m(ji+1,jj) - ssh_m(ji,jj))/e1u(ji,jj) 
    317             zdsshy =  (ssh_m(ji,jj+1) - ssh_m(ji,jj))/e2v(ji,jj) 
     279            zdsshx =  ( ssh_m(ji+1,jj) - ssh_m(ji,jj) ) / e1u(ji,jj) 
     280            zdsshy =  ( ssh_m(ji,jj+1) - ssh_m(ji,jj) ) / e2v(ji,jj) 
    318281 
    319282            za1ct(ji,jj) = ztagnx - zmass1(ji,jj) * grav * zdsshx 
     
    330293      ! Time step for subcycling 
    331294      dtevp  = rdt_ice / nevp 
    332       dtotel = dtevp / ( 2.0 * telast ) 
     295      dtotel = dtevp / ( 2._wp * telast ) 
    333296 
    334297      !-ecc2: square of yield ellipse eccenticrity (reminder: must become a namelist parameter) 
    335       ecc2 = ecc*ecc 
     298      ecc2 = ecc * ecc 
    336299 
    337300      !-Initialise stress tensor  
    338       zs1(:,:)  = stress1_i(:,:)  
    339       zs2(:,:)  = stress2_i(:,:) 
     301      zs1 (:,:) = stress1_i (:,:)  
     302      zs2 (:,:) = stress2_i (:,:) 
    340303      zs12(:,:) = stress12_i(:,:) 
    341304 
    342       !----------------------! 
     305      !                                               !----------------------! 
    343306      DO jter = 1 , nevp                              !    loop over jter    ! 
    344          !----------------------!         
     307         !                                            !----------------------!         
    345308         DO jj = k_j1, k_jpj-1 
    346             zu_ice(:,jj) = u_ice(:,jj) ! velocity at previous time step 
     309            zu_ice(:,jj) = u_ice(:,jj)    ! velocity at previous time step 
    347310            zv_ice(:,jj) = v_ice(:,jj) 
    348311         END DO 
     
    409372            END DO 
    410373         END DO 
    411          CALL lbc_lnk( v_ice1(:,:), 'U', -1. ) 
    412          CALL lbc_lnk( u_ice2(:,:), 'V', -1. ) 
     374         CALL lbc_lnk( v_ice1, 'U', -1. )   ;   CALL lbc_lnk( u_ice2, 'V', -1. )      ! lateral boundary cond. 
    413375 
    414376!CDIR NOVERRCHK 
     
    418380 
    419381               !- Calculate Delta at centre of grid cells 
    420                zdst       = (  e2u( ji  , jj   ) * v_ice1(ji,jj)            & 
    421                   &          - e2u( ji-1, jj   ) * v_ice1(ji-1,jj)          & 
    422                   &          + e1v( ji  , jj   ) * u_ice2(ji,jj)            & 
    423                   &          - e1v( ji  , jj-1 ) * u_ice2(ji,jj-1)          & 
    424                   &          )                                              & 
     382               zdst       = (  e2u(ji  , jj) * v_ice1(ji  ,jj)          & 
     383                  &          - e2u(ji-1, jj) * v_ice1(ji-1,jj)          & 
     384                  &          + e1v(ji, jj  ) * u_ice2(ji,jj  )          & 
     385                  &          - e1v(ji, jj-1) * u_ice2(ji,jj-1)          & 
     386                  &          )                                          & 
    425387                  &         / area(ji,jj) 
    426388 
    427                delta = SQRT( zdd(ji,jj)*zdd(ji,jj) +                        &  
    428                   &                       ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 )   
    429                deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 +                    & 
    430                   (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 
     389               delta = SQRT( zdd(ji,jj)*zdd(ji,jj) + ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 )   
     390               deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 + (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 
    431391 
    432392               !-Calculate stress tensor components zs1 and zs2  
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

    r2287 r2370  
    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 
     
    3537   PUBLIC   lim_sbc_tau   ! called by sbc_ice_lim 
    3638 
    37    REAL(wp)  ::   epsi16 = 1.e-16  ! constant values 
    38    REAL(wp)  ::   rzero  = 0.e0     
    39    REAL(wp)  ::   rone   = 1.e0 
    40  
    41    REAL(wp), DIMENSION(jpi,jpj) ::   utau_oce, vtau_oce   !: air-ocean surface i- & j-stress              [N/m2] 
    42    REAL(wp), DIMENSION(jpi,jpj) ::   tmod_io              !: modulus of the ice-ocean relative velocity   [m/s] 
    43    REAL(wp), DIMENSION(jpi,jpj) ::   ssu_mb  , ssv_mb     !: before mean ocean surface currents           [m/s] 
     39   REAL(wp)  ::   r1_rdtice            ! = 1. / rdt_ice  
     40   REAL(wp)  ::   epsi16 = 1.e-16_wp   ! constant values 
     41   REAL(wp)  ::   rzero  = 0._wp     
     42   REAL(wp)  ::   rone   = 1._wp 
     43 
     44   REAL(wp), DIMENSION(jpi,jpj) ::   utau_oce, vtau_oce   ! air-ocean surface i- & j-stress              [N/m2] 
     45   REAL(wp), DIMENSION(jpi,jpj) ::   tmod_io              ! modulus of the ice-ocean relative velocity   [m/s] 
     46 
     47   REAL(wp), DIMENSION(jpi,jpj) ::   soce_0, sice_0   ! constant SSS and ice salinity used in levitating sea-ice case 
    4448 
    4549   !! * Substitutions 
     
    4852   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4953   !! $Id$ 
    50    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    51    !!---------------------------------------------------------------------- 
    52  
     54   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     55   !!---------------------------------------------------------------------- 
    5356CONTAINS 
    54  
    55    SUBROUTINE lim_sbc_tau( kt, kcpl ) 
    56       !!------------------------------------------------------------------- 
    57       !!                ***  ROUTINE lim_sbc_tau *** 
    58       !!   
    59       !! ** Purpose : Update the ocean surface stresses due to the ice 
    60       !!          
    61       !! ** Action  : - compute the ice-ocean stress depending on kcpl: 
    62       !!              fluxes at the ice-ocean interface. 
    63       !!     Case 0  :  old LIM-3 way, computed at ice time-step only 
    64       !!     Case 1  :  at each ocean time step the stresses are computed 
    65       !!                using the current ocean velocity (now) 
    66       !!     Case 2  :  combination of half case 0 + half case 1 
    67       !!      
    68       !! ** Outputs : - utau   : sea surface i-stress (ocean referential) 
    69       !!              - vtau   : sea surface j-stress (ocean referential) 
    70       !! 
    71       !! References : Goosse, H. et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90. 
    72       !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
    73       !!--------------------------------------------------------------------- 
    74       INTEGER ::   kt     ! number of ocean iteration 
    75       INTEGER ::   kcpl   ! ice-ocean coupling indicator: =0  LIM-3 old case 
    76       !                   !                               =1  stresses computed using now ocean velocity 
    77       !                   !                               =2  combination of 0 and 1 cases 
    78       !! 
    79       INTEGER  ::   ji, jj   ! dummy loop indices 
    80       REAL(wp) ::   zfrldu, zat_u, zu_ico, zutaui, zu_u, zu_ij, zu_im1j   ! temporary scalar 
    81       REAL(wp) ::   zfrldv, zat_v, zv_ico, zvtaui, zv_v, zv_ij, zv_ijm1   !    -         - 
    82       REAL(wp) ::   zsang, zztmp                                          !    -         - 
    83       REAL(wp), DIMENSION(jpi,jpj) ::   ztio_u, ztio_v   ! ocean stress below sea-ice 
    84 #if defined key_coupled     
    85       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalb     ! albedo of ice under overcast sky 
    86       REAL(wp), DIMENSION(jpi,jpj,jpl) ::   zalbp    ! albedo of ice under clear sky 
    87 #endif 
    88      !!--------------------------------------------------------------------- 
    89  
    90       IF( kt == nit000 ) THEN 
    91          IF(lwp) WRITE(numout,*) 
    92          IF(lwp) WRITE(numout,*) 'lim_sbc_tau : LIM 3.0 sea-ice - surface ocean momentum fluxes' 
    93          IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
    94       ENDIF 
    95  
    96       IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    97 !CDIR NOVERRCHK 
    98          DO jj = 2, jpjm1                         !* modulus of the ice-ocean velocity 
    99 !CDIR NOVERRCHK 
    100             DO ji = fs_2, fs_jpim1 
    101                zu_ij   = u_ice(ji  ,jj) - ssu_m(ji  ,jj)               ! (i  ,j) 
    102                zu_im1j = u_ice(ji-1,jj) - ssu_m(ji-1,jj)               ! (i-1,j) 
    103                zv_ij   = v_ice(ji,jj  ) - ssv_m(ji,jj  )               ! (i,j  ) 
    104                zv_ijm1 = v_ice(ji,jj-1) - ssv_m(ji,jj-1)               ! (i,j-1) 
    105                zztmp =  0.5 * ( zu_ij * zu_ij + zu_im1j * zu_im1j   &  ! mean of the square values instead 
    106                   &           + zv_ij * zv_ij + zv_ijm1 * zv_ijm1 )    ! of the square of the mean values... 
    107                ! WARNING, here we make an approximation for case 1 and 2: taum is not computed at each time 
    108                ! step but only every nn_fsbc time step. This avoid mutiple avarage to pass from T -> U,V grids 
    109                ! and next from U,V grids to T grid. Taum is used in tke, which should not be too sensitive to 
    110                ! this approximaton... 
    111                taum(ji,jj) = ( 1. - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zztmp 
    112                tmod_io(ji,jj) = SQRT( zztmp )  
    113             END DO 
    114          END DO 
    115          CALL lbc_lnk( taum, 'T', 1. )   ;   CALL lbc_lnk( tmod_io, 'T', 1. ) 
    116       ENDIF 
    117  
    118       SELECT CASE( kcpl ) 
    119          !                                        !--------------------------------! 
    120       CASE( 0 )                                   !  LIM 3 old stress computation  !  (at ice timestep only) 
    121          !                                        !--------------------------------!  
    122          !                                             !* ice stress over ocean with a ice-ocean rotation angle 
    123          DO jj = 1, jpjm1 
    124             zsang  = SIGN( 1.e0, gphif(1,jj) ) * sangvg         ! change the sinus angle sign in the south hemisphere 
    125             DO ji = 1, fs_jpim1 
    126                zu_u = u_ice(ji,jj) - u_oce(ji,jj)               ! ice velocity relative to the ocean 
    127                zv_v = v_ice(ji,jj) - v_oce(ji,jj) 
    128                !                                                ! quadratic drag formulation with rotation 
    129 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1) 
    130                zutaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_u - zsang * zv_v )  
    131                zvtaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_v + zsang * zu_u )  
    132                !                                                ! bound for too large stress values 
    133                ! IMPORTANT: the exponential below prevents numerical oscillations in the ice-ocean stress 
    134                ! This is not physically based. A cleaner solution is offer in CASE kcpl=2 
    135                ztio_u(ji,jj) = zutaui * EXP( - ( tmod_io(ji,jj) + tmod_io(ji+1,jj) )  ) 
    136                ztio_v(ji,jj) = zvtaui * EXP( - ( tmod_io(ji,jj) + tmod_io(ji,jj+1) )  )  
    137             END DO 
    138          END DO 
    139          DO jj = 2, jpjm1                                       ! stress at the surface of the ocean 
    140             DO ji = fs_2, fs_jpim1   ! vertor opt. 
    141                zfrldu = 0.5 * ( ato_i(ji,jj) + ato_i(ji+1,jj  ) )   ! open-ocean fraction at U- & V-points (from T-point values) 
    142                zfrldv = 0.5 * ( ato_i(ji,jj) + ato_i(ji  ,jj+1) ) 
    143                !                                                    ! update surface ocean stress 
    144                utau(ji,jj) = zfrldu * utau(ji,jj) + ( 1. - zfrldu ) * ztio_u(ji,jj) 
    145                vtau(ji,jj) = zfrldv * vtau(ji,jj) + ( 1. - zfrldv ) * ztio_v(ji,jj) 
    146             END DO 
    147          END DO 
    148          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
    149  
    150          ! 
    151          !                                        !--------------------------------! 
    152       CASE( 1 )                                   !  Use the now velocity          !  (at each ocean timestep) 
    153          !                                        !--------------------------------!  
    154          IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    155             utau_oce(:,:) = utau(:,:)                 !* save the air-ocean stresses at ice time-step 
    156             vtau_oce(:,:) = vtau(:,:) 
    157          ENDIF 
    158          ! 
    159         DO jj = 2, jpjm1                              !* ice stress over ocean with a ice-ocean rotation angle 
    160             zsang  = SIGN(1.e0, gphif(1,jj-1) ) * sangvg        ! change the sinus angle sign in the south hemisphere 
    161             DO ji = fs_2, fs_jpim1 
    162                zat_u  = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5   ! ice area at u and V-points 
    163                zat_v  = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5 
    164                !                                                ! (u,v) ice-ocean velocity at (U,V)-point, resp. 
    165                zu_u   = u_ice(ji,jj) - ub(ji,jj,1) 
    166                zv_v   = v_ice(ji,jj) - vb(ji,jj,1) 
    167                !                                                ! quadratic drag formulation with rotation 
    168 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1) 
    169                zutaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_u - zsang * zv_v )  
    170                zvtaui   = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_v + zsang * zu_u )  
    171                !                                                   ! stress at the ocean surface 
    172                utau(ji,jj) = ( 1.- zat_u ) * utau_oce(ji,jj) + zat_u * zutaui 
    173                vtau(ji,jj) = ( 1.- zat_v ) * vtau_oce(ji,jj) + zat_v * zvtaui 
    174             END DO 
    175          END DO 
    176          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
    177  
    178          ! 
    179          !                                        !--------------------------------! 
    180       CASE( 2 )                                   !  mixed 0 and 2 cases           !  (at each ocean timestep) 
    181          !                                        !--------------------------------!  
    182          IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    183             utau_oce(:,:) = utau (:,:)               !* save the air-ocean stresses at ice time-step 
    184             vtau_oce(:,:) = vtau (:,:) 
    185             ssu_mb  (:,:) = ssu_m(:,:)               !* save the ice-ocean velocity at ice time-step 
    186             ssv_mb  (:,:) = ssv_m(:,:) 
    187          ENDIF 
    188          DO jj = 2, jpjm1                            !* ice stress over ocean with a ice-ocean rotation angle 
    189             zsang  = SIGN(1.e0, gphif(1,jj-1) ) * sangvg 
    190             DO ji = fs_2, fs_jpim1 
    191                zat_u = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5     ! ice area at u and V-points 
    192                zat_v = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5  
    193                ! 
    194                zu_ico = u_ice(ji,jj) - 0.5 * ( ub(ji,jj,1) + ssu_mb(ji,jj) )   ! ice-oce velocity using ub and ssu_mb 
    195                zv_ico = v_ice(ji,jj) - 0.5 * ( vb(ji,jj,1) + ssv_mb(ji,jj) ) 
    196                !                                        ! quadratic drag formulation with rotation 
    197 !!gm still an error in the rotation, but usually the angle is zero (zsang=0, cangvg=1) 
    198                zutaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * rhoco * ( cangvg * zu_ico - zsang * zv_ico ) 
    199                zvtaui = 0.5 * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * rhoco * ( cangvg * zv_ico + zsang * zu_ico ) 
    200                ! 
    201                utau(ji,jj) = ( 1.-zat_u ) * utau_oce(ji,jj) + zat_u * zutaui    ! stress at the ocean surface 
    202                vtau(ji,jj) = ( 1.-zat_v ) * vtau_oce(ji,jj) + zat_v * zvtaui 
    203             END DO 
    204          END DO 
    205          CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
    206          ! 
    207       END SELECT 
    208  
    209       IF(ln_ctl)   CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
    210          &                       tab2d_2=vtau, clinfo2=' vtau    : '        , mask2=vmask ) 
    211       !   
    212    END SUBROUTINE lim_sbc_tau 
    213  
    21457 
    21558   SUBROUTINE lim_sbc_flx( kt ) 
     
    23578      !!              Tartinville et al. 2001 Ocean Modelling, 3, 95-108. 
    23679      !!--------------------------------------------------------------------- 
    237       INTEGER ::   kt    ! number of iteration 
     80      INTEGER, INTENT(in) ::   kt    ! number of iteration 
    23881      !! 
    23982      INTEGER  ::   ji, jj           ! dummy loop indices 
     
    25497         IF(lwp) WRITE(numout,*) 'lim_sbc_flx : LIM 3.0 sea-ice - heat salt and mass ocean surface fluxes' 
    25598         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
     99         ! 
     100         r1_rdtice = 1. / rdt_ice 
     101         ! 
     102         soce_0(:,:) = soce 
     103         sice_0(:,:) = sice 
     104         ! 
     105         IF( cp_cfg == "orca" ) THEN           ! decrease ocean & ice reference salinities in the Baltic sea  
     106            WHERE( 14._wp <= glamt(:,:) .AND. glamt(:,:) <= 32._wp .AND.   & 
     107               &   54._wp <= gphit(:,:) .AND. gphit(:,:) <= 66._wp         )  
     108               soce_0(:,:) = 4._wp 
     109               sice_0(:,:) = 2._wp 
     110            END WHERE 
     111         ENDIF 
     112         ! 
    256113      ENDIF 
    257114 
     
    298155               ! fscmbq and ffltbif are obsolete 
    299156               !              &           + iflt * ffltbif(ji,jj) !!! only if one category is used 
    300                &           + ifral   * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) / rdt_ice   & 
    301                &           + ifrdv   * ( qfvbq(ji,jj) + qdtcn(ji,jj) ) / rdt_ice                     & 
     157               &           + ifral   * ( ial * qcmif(ji,jj) + (1 - ial) * qldif(ji,jj) ) * r1_rdtice   & 
     158               &           + ifrdv   * ( qfvbq(ji,jj) + qdtcn(ji,jj) )                   * r1_rdtice   & 
    302159               &           + fhmec(ji,jj)     & ! new contribution due to snow melt in ridging!! 
    303160               &           + fheat_rpo(ji,jj) & ! contribution from ridge formation 
     
    340197               WRITE(numout,*) ' qcmif     : ', qcmif(jiindx,jjindx) 
    341198               WRITE(numout,*) ' qldif     : ', qldif(jiindx,jjindx) 
    342                WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindx,jjindx) / rdt_ice 
    343                WRITE(numout,*) ' qldif / dt: ', qldif(jiindx,jjindx) / rdt_ice 
     199               WRITE(numout,*) ' qcmif / dt: ', qcmif(jiindx,jjindx) * r1_rdtice 
     200               WRITE(numout,*) ' qldif / dt: ', qldif(jiindx,jjindx) * r1_rdtice 
    344201               WRITE(numout,*) ' ifrdv     : ', ifrdv 
    345202               WRITE(numout,*) ' qfvbq     : ', qfvbq(jiindx,jjindx) 
    346203               WRITE(numout,*) ' qdtcn     : ', qdtcn(jiindx,jjindx) 
    347                WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) / rdt_ice 
    348                WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) / rdt_ice 
     204               WRITE(numout,*) ' qfvbq / dt: ', qfvbq(jiindx,jjindx) * r1_rdtice 
     205               WRITE(numout,*) ' qdtcn / dt: ', qdtcn(jiindx,jjindx) * r1_rdtice 
    349206               WRITE(numout,*) ' ' 
    350207               WRITE(numout,*) ' fdtcn     : ', fdtcn(jiindx,jjindx) 
     
    380237               !              &   - sprecip(ji,jj) * ( 1. - pfrld(ji,jj) )  &   !  remov. snow precip over ice 
    381238               &   - sprecip(ji,jj) * ( 1. - (pfrld(ji,jj)**betas) )  &   !  remov. snow precip over ice 
    382                &   - rdmsnif(ji,jj) / rdt_ice                &   !  freshwaterflux due to snow melting  
     239               &   - rdmsnif(ji,jj) * r1_rdtice                       &   !  freshwaterflux due to snow melting  
    383240               ! new contribution from snow falling when ridging 
    384241               &   + fmmec(ji,jj) 
     
    386243            !  computing salt exchanges at the ice/ocean interface 
    387244            !  sice should be the same as computed with the ice model 
    388             zfons =  ( soce - sice ) * ( rdmicif(ji,jj) / rdt_ice )  
     245            zfons =  ( soce_0(ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice  
    389246            ! SOCE 
    390             zfons =  ( sss_m(ji,jj) - sice ) * ( rdmicif(ji,jj) / rdt_ice )  
     247            zfons =  ( sss_m (ji,jj) - sice_0(ji,jj) ) * rdmicif(ji,jj) * r1_rdtice 
    391248 
    392249            !CT useless            !  salt flux for constant salinity 
     
    446303   END SUBROUTINE lim_sbc_flx 
    447304 
     305 
     306   SUBROUTINE lim_sbc_tau( kt , pu_oce, pv_oce ) 
     307      !!------------------------------------------------------------------- 
     308      !!                ***  ROUTINE lim_sbc_tau *** 
     309      !!   
     310      !! ** Purpose : Update the ocean surface stresses due to the ice 
     311      !!          
     312      !! ** Action  : * at each ice time step (every nn_fsbc time step): 
     313      !!                - compute the modulus of ice-ocean relative velocity  
     314      !!                  (*rho*Cd) at T-point (C-grid) or I-point (B-grid) 
     315      !!                      tmod_io = rhoco * | U_ice-U_oce | 
     316      !!                - update the modulus of stress at ocean surface 
     317      !!                      taum = frld * taum + (1-frld) * tmod_io * | U_ice-U_oce | 
     318      !!              * at each ocean time step (every kt):  
     319      !!                  compute linearized ice-ocean stresses as 
     320      !!                      Utau = tmod_io * | U_ice - pU_oce | 
     321      !!                using instantaneous current ocean velocity (usually before) 
     322      !! 
     323      !!    NB: - ice-ocean rotation angle no more allowed 
     324      !!        - here we make an approximation: taum is only computed every ice time step 
     325      !!          This avoids mutiple average to pass from T -> U,V grids and next from U,V grids  
     326      !!          to T grid. taum is used in TKE and GLS, which should not be too sensitive to this approximaton... 
     327      !! 
     328      !! ** Outputs : - utau, vtau   : surface ocean i- and j-stress (u- & v-pts) updated with ice-ocean fluxes 
     329      !!              - taum         : modulus of the surface ocean stress (T-point) updated with ice-ocean fluxes 
     330      !!--------------------------------------------------------------------- 
     331      INTEGER ,                     INTENT(in) ::   kt               ! ocean time-step index 
     332      REAL(wp), DIMENSION(jpi,jpj), INTENT(in) ::   pu_oce, pv_oce   ! surface ocean currents 
     333      !! 
     334      INTEGER  ::   ji, jj   ! dummy loop indices 
     335      REAL(wp) ::   zat_u, zutau_ice, zu_t, zmodt   ! local scalar 
     336      REAL(wp) ::   zat_v, zvtau_ice, zv_t          !   -      - 
     337     !!--------------------------------------------------------------------- 
     338 
     339      IF( kt == nit000 ) THEN 
     340         IF(lwp) WRITE(numout,*) 
     341         IF(lwp) WRITE(numout,*) 'lim_sbc_tau : LIM-3 sea-ice - surface ocean momentum fluxes' 
     342         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~ ' 
     343      ENDIF 
     344 
     345      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !==  Ice time-step only  ==!   (i.e. surface module time-step) 
     346!CDIR NOVERRCHK 
     347         DO jj = 2, jpjm1                             !* update the modulus of stress at ocean surface (T-point) 
     348!CDIR NOVERRCHK 
     349            DO ji = fs_2, fs_jpim1 
     350               !                                               ! 2*(U_ice-U_oce) at T-point 
     351               zu_t = u_ice(ji,jj) + u_ice(ji-1,jj) - u_oce(ji,jj) - u_oce(ji-1,jj)    
     352               zv_t = v_ice(ji,jj) + v_ice(ji,jj-1) - v_oce(ji,jj) - v_oce(ji,jj-1)  
     353               !                                              ! |U_ice-U_oce|^2 
     354               zmodt =  0.25_wp * (  zu_t * zu_t + zv_t * zv_t  ) 
     355               !                                               ! update the ocean stress modulus 
     356               taum(ji,jj) = ( 1._wp - at_i(ji,jj) ) * taum(ji,jj) + at_i(ji,jj) * rhoco * zmodt 
     357               tmod_io(ji,jj) = rhoco * SQRT( zmodt )          ! rhoco * |U_ice-U_oce| at T-point 
     358            END DO 
     359         END DO 
     360         CALL lbc_lnk( taum, 'T', 1. )   ;   CALL lbc_lnk( tmod_io, 'T', 1. ) 
     361         ! 
     362         utau_oce(:,:) = utau(:,:)                    !* save the air-ocean stresses at ice time-step 
     363         vtau_oce(:,:) = vtau(:,:) 
     364         ! 
     365      ENDIF 
     366         ! 
     367         !                                     !==  every ocean time-step  ==! 
     368         ! 
     369      DO jj = 2, jpjm1                                !* update the stress WITHOUT a ice-ocean rotation angle 
     370         DO ji = fs_2, fs_jpim1   ! Vect. Opt. 
     371            zat_u  = ( at_i(ji,jj) + at_i(ji+1,jj) ) * 0.5_wp   ! ice area at u and V-points 
     372            zat_v  = ( at_i(ji,jj) + at_i(ji,jj+1) ) * 0.5_wp 
     373            !                                                   ! linearized quadratic drag formulation 
     374            zutau_ice   = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji+1,jj) ) * ( u_ice(ji,jj) - pu_oce(ji,jj) ) 
     375            zvtau_ice   = 0.5_wp * ( tmod_io(ji,jj) + tmod_io(ji,jj+1) ) * ( v_ice(ji,jj) - pv_oce(ji,jj) ) 
     376            !                                                   ! stresses at the ocean surface 
     377            utau(ji,jj) = ( 1._wp - zat_u ) * utau_oce(ji,jj) + zat_u * zutau_ice 
     378            vtau(ji,jj) = ( 1._wp - zat_v ) * vtau_oce(ji,jj) + zat_v * zvtau_ice 
     379         END DO 
     380      END DO 
     381      CALL lbc_lnk( utau, 'U', -1. )   ;   CALL lbc_lnk( vtau, 'V', -1. )   ! lateral boundary condition 
     382      ! 
     383      IF(ln_ctl)   CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
     384         &                       tab2d_2=vtau, clinfo2=' vtau    : '        , mask2=vmask ) 
     385      !   
     386   END SUBROUTINE lim_sbc_tau 
     387 
    448388#else 
    449389   !!---------------------------------------------------------------------- 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_ice.F90

    r2287 r2370  
    99#if defined key_lim3 || defined key_lim2 
    1010   !!---------------------------------------------------------------------- 
    11    !!   'key_lim2' or 'key_lim3' :             LIM 2.0 or 3.0 sea-ice model 
     11   !!   'key_lim2' or 'key_lim3' :             LIM-2 or LIM-3 sea-ice model 
    1212   !!---------------------------------------------------------------------- 
    1313   USE par_oce          ! ocean parameters 
     
    2323 
    2424# if defined  key_lim2 
    25    LOGICAL         , PUBLIC, PARAMETER ::   lk_lim2        = .TRUE.    !: LIM-2 ice model 
    26    LOGICAL         , PUBLIC, PARAMETER ::   lk_lim3        = .FALSE.   !: no LIM-3 
    27    CHARACTER(len=1), PUBLIC            ::   cigr_type      = 'I'       !: 'I'-grid ice-velocity (B-grid lower left corner) 
     25   LOGICAL         , PUBLIC, PARAMETER ::   lk_lim2    = .TRUE.   !: LIM-2 ice model 
     26   LOGICAL         , PUBLIC, PARAMETER ::   lk_lim3    = .FALSE.  !: no LIM-3 
     27#  if defined key_lim2_vp 
     28   CHARACTER(len=1), PUBLIC, PARAMETER ::   cp_ice_msh = 'I'      !: VP : 'I'-grid ice-velocity (B-grid lower left corner) 
     29#  else 
     30   CHARACTER(len=1), PUBLIC, PARAMETER ::   cp_ice_msh = 'C'      !: EVP: 'C'-grid ice-velocity 
     31#  endif 
    2832# endif 
    2933# if defined  key_lim3 
    30    LOGICAL         , PUBLIC, PARAMETER ::   lk_lim2        = .FALSE.   !: no LIM-2 
    31    LOGICAL         , PUBLIC, PARAMETER ::   lk_lim3        = .TRUE.    !: LIM-3 ice model 
    32    CHARACTER(len=1), PUBLIC            ::   cigr_type      = 'C'       !: 'C'-grid ice-velocity 
     34   LOGICAL         , PUBLIC, PARAMETER ::   lk_lim2    = .FALSE.  !: no LIM-2 
     35   LOGICAL         , PUBLIC, PARAMETER ::   lk_lim3    = .TRUE.   !: LIM-3 ice model 
     36   CHARACTER(len=1), PUBLIC, PARAMETER ::   cp_ice_msh = 'C'      !: 'C'-grid ice-velocity 
    3337# endif 
    3438 
     
    4145   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpl) ::   alb_ice   !: albedo of ice 
    4246 
    43    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   utau_ice    !: u-stress over ice (I-point for LIM2 or U,V-point for LIM3)   [N/m2] 
    44    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   vtau_ice    !: v-stress over ice (I-point for LIM2 or U,V-point for LIM3)   [N/m2] 
    45    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fr1_i0      !: 1st fraction of sol. rad. which penetrate inside the ice cover 
    46    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fr2_i0      !: 2nd fraction of sol. rad. which penetrate inside the ice cover 
     47   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   utau_ice    !: u-stress over ice (I-pt for VP or U,V-pts for EVP)   [N/m2] 
     48   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   vtau_ice    !: v-stress over ice (I-pt for VP or U,V-pts for EVP)   [N/m2] 
     49   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fr1_i0      !: 1st fraction of Qsr which penetrates inside the ice cover 
     50   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fr2_i0      !: 2nd fraction of Qsr which penetrates inside the ice cover 
    4751   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   emp_ice     !: solid freshwater budget over ice: sublivation - snow 
    4852 
     
    5559   !!   Default option                      NO LIM 2.0 or 3.0 sea-ice model 
    5660   !!---------------------------------------------------------------------- 
    57    LOGICAL         , PUBLIC, PARAMETER ::   lk_lim2        = .FALSE.  !: no LIM-2 ice model 
    58    LOGICAL         , PUBLIC, PARAMETER ::   lk_lim3        = .FALSE.  !: no LIM-3 ice model 
    59    CHARACTER(len=1), PUBLIC            ::   cigr_type      = '-'      !: no grid ice-velocity 
     61   LOGICAL         , PUBLIC, PARAMETER ::   lk_lim2    = .FALSE.  !: no LIM-2 ice model 
     62   LOGICAL         , PUBLIC, PARAMETER ::   lk_lim3    = .FALSE.  !: no LIM-3 ice model 
     63   CHARACTER(len=1), PUBLIC, PARAMETER ::   cp_ice_msh = '-'      !: no grid ice-velocity 
    6064#endif 
    6165 
     
    6367   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    6468   !! $Id$  
    65    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    66    !!---------------------------------------------------------------------- 
    67  
     69   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    6870   !!====================================================================== 
    6971END MODULE sbc_ice 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbc_oce.F90

    r2287 r2370  
    44   !! Surface module :   variables defined in core memory  
    55   !!====================================================================== 
    6    !! History :  3.0   !  2006-06  (G. Madec)  Original code 
    7    !!             -    !  2008-08  (G. Madec)  namsbc moved from sbcmod 
    8    !!            3.3   !  2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps 
     6   !! History :  3.0  ! 2006-06  (G. Madec)  Original code 
     7   !!             -   ! 2008-08  (G. Madec)  namsbc moved from sbcmod 
     8   !!            3.3  ! 2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps 
     9   !!             -   ! 2010-11  (G. Madec) ice-ocean stress always computed at each ocean time-step 
    910   !!---------------------------------------------------------------------- 
    1011   USE par_oce          ! ocean parameters 
     
    3031   !                                             !:  = 1 global mean of e-p-r set to zero at each nn_fsbc time step 
    3132   !                                             !:  = 2 annual global mean of e-p-r set to zero 
    32    INTEGER , PUBLIC ::   nn_ico_cpl  = 0          !: ice-ocean coupling indicator 
    33    !                                             !:  = 0   LIM-3 old case 
    34    !                                             !:  = 1   stresses computed using now ocean velocity 
    35    !                                             !:  = 2   combination of 0 and 1 cases 
    3633 
    3734   !!---------------------------------------------------------------------- 
    3835   !!              Ocean Surface Boundary Condition fields 
    3936   !!---------------------------------------------------------------------- 
    40    LOGICAL , PUBLIC ::   lhftau = .FALSE.              !: HF tau used in TKE: mean(stress module) - module(mean stress) 
     37   LOGICAL , PUBLIC ::   lhftau = .FALSE.        !: HF tau used in TKE: mean(stress module) - module(mean stress) 
    4138   !!                                   !!   now    ! before   !! 
    4239   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   utau   , utau_b   !: sea surface i-stress (ocean referential)     [N/m2] 
     
    5350   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   emp_tot           !: total E-P over ocean and ice                 [Kg/m2/s] 
    5451   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rnf    , rnf_b    !: river runoff   [Kg/m2/s]   
    55    ! - ML - begin 
     52   !! 
    5653   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpts) ::  sbc_tsc, sbc_tsc_b  !: sbc content trend                      [K.m/s] 
    5754   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpk) ::   qsr_hc , qsr_hc_b   !: heat content trend due to qsr flux     [K.m/s] 
    58    ! - ML - end 
     55   !! 
    5956   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tprecip           !: total precipitation                          [Kg/m2/s] 
    6057   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sprecip           !: solid precipitation                          [Kg/m2/s] 
     
    6360   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   atm_co2           !: atmospheric pCO2                             [ppm] 
    6461#endif 
    65 !!$   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rrunoff        !: runoff 
    66 !!$   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   calving        !: calving 
    6762 
    6863   !!---------------------------------------------------------------------- 
     
    7974   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    8075   !! $Id$ 
    81    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     76   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    8277   !!====================================================================== 
    8378END MODULE sbc_oce 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_clio.F90

    r2287 r2370  
    8383   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    8484   !! $Id$  
    85    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     85   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    8686   !!---------------------------------------------------------------------- 
    87  
    8887CONTAINS 
    8988 
     
    137136 
    138137         ! (NB: frequency positive => hours, negative => months) 
    139          !            !    file     ! frequency !  variable  ! time intep !  clim   ! 'yearly' or ! weights  ! rotation  ! 
    140          !            !    name     !  (hours)  !   name     !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs      ! 
    141          sn_utau = FLD_N( 'utau'    ,    24     ,  'utau'    ,  .true.    , .false. ,   'yearly'  , ''       , ''         )  
    142          sn_vtau = FLD_N( 'vtau'    ,    24     ,  'vtau'    ,  .true.    , .false. ,   'yearly'  , ''       , ''         )  
    143          sn_wndm = FLD_N( 'mwnd10m' ,    24     ,  'm_10'    ,  .true.    , .false. ,   'yearly'  , ''       , ''         )  
    144          sn_tair = FLD_N( 'tair10m' ,    24     ,  't_10'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         )  
    145          sn_humi = FLD_N( 'humi10m' ,    24     ,  'q_10'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         )  
    146          sn_ccov = FLD_N( 'ccover'  ,    -1     ,  'cloud'   ,  .true.    , .false. ,   'yearly'  , ''       , ''         )  
    147          sn_prec = FLD_N( 'precip'  ,    -1     ,  'precip'  ,  .true.    , .false. ,   'yearly'  , ''       , ''         )  
     138         !            !    file    ! frequency ! variable ! time intep !  clim   ! 'yearly' or ! weights  ! rotation ! 
     139         !            !    name    !  (hours)  !  name    !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs    ! 
     140         sn_utau = FLD_N( 'utau'   ,    24     , 'utau'   ,  .true.    , .false. ,   'yearly'  , ''       , ''       )  
     141         sn_vtau = FLD_N( 'vtau'   ,    24     , 'vtau'   ,  .true.    , .false. ,   'yearly'  , ''       , ''       )  
     142         sn_wndm = FLD_N( 'mwnd10m',    24     , 'm_10'   ,  .true.    , .false. ,   'yearly'  , ''       , ''       )  
     143         sn_tair = FLD_N( 'tair10m',    24     , 't_10'   ,  .false.   , .false. ,   'yearly'  , ''       , ''       )  
     144         sn_humi = FLD_N( 'humi10m',    24     , 'q_10'   ,  .false.   , .false. ,   'yearly'  , ''       , ''       )  
     145         sn_ccov = FLD_N( 'ccover' ,    -1     , 'cloud'  ,  .true.    , .false. ,   'yearly'  , ''       , ''       )  
     146         sn_prec = FLD_N( 'precip' ,    -1     , 'precip' ,  .true.    , .false. ,   'yearly'  , ''       , ''       )  
    148147 
    149148         REWIND( numnam )                    ! ... read in namlist namsbc_clio 
     
    473472      INTEGER  ::   ijpl          ! number of ice categories (size of 3rd dim of input arrays) 
    474473      !! 
    475       REAL(wp) ::   zcoef, zmt1, zmt2, zmt3, ztatm3             ! temporary scalars 
     474      REAL(wp) ::   zcoef, zmt1, zmt2, zmt3, ztatm3     ! temporary scalars 
    476475      REAL(wp) ::   ztaevbk, zind1, zind2, zind3, ztamr         !    -         - 
    477476      REAL(wp) ::   zesi, zqsati, zdesidt                       !    -         - 
     
    496495      SELECT CASE( cd_grid ) 
    497496      CASE( 'C' )                          ! C-grid ice dynamics 
    498          ! Change from wind speed to wind stress over OCEAN (cao is used) 
    499          zcoef = cai / cao  
     497      zcoef  = cai / cao                         ! Change from air-sea stress to air-ice stress 
    500498!CDIR COLLAPSE 
    501499         DO jj = 1 , jpj 
     
    505503            END DO 
    506504         END DO 
    507       CASE( 'B' )                          ! B-grid ice dynamics 
    508          ! stress from ocean U- and V-points to ice U,V point 
    509 !CDIR COLLAPSE 
    510          DO jj = 2, jpj 
    511             DO ji = 2, jpi   ! B grid : no vector opt. 
    512                p_taui(ji,jj) = 0.5 * ( utau(ji-1,jj  ) + utau(ji-1,jj-1) ) 
    513                p_tauj(ji,jj) = 0.5 * ( vtau(ji  ,jj-1) + vtau(ji-1,jj-1) ) 
     505      CASE( 'I' )                          ! I-grid ice dynamics:  I-point (i.e. F-point lower-left corner) 
     506      zcoef  = 0.5_wp * cai / cao                ! Change from air-sea stress to air-ice stress!CDIR COLLAPSE 
     507         DO jj = 2, jpj         ! stress from ocean U- and V-points to ice U,V point 
     508            DO ji = 2, jpi   ! I-grid : no vector opt. 
     509               p_taui(ji,jj) = zcoef * ( utau(ji-1,jj  ) + utau(ji-1,jj-1) ) 
     510               p_tauj(ji,jj) = zcoef * ( vtau(ji  ,jj-1) + vtau(ji-1,jj-1) ) 
    514511            END DO 
    515512         END DO 
    516          CALL lbc_lnk( p_taui(:,:), 'I', -1. )   ! I-point (i.e. ice U-V point) 
    517          CALL lbc_lnk( p_tauj(:,:), 'I', -1. )   ! I-point (i.e. ice U-V point) 
     513         CALL lbc_lnk( p_taui(:,:), 'I', -1. )   ;   CALL lbc_lnk( p_tauj(:,:), 'I', -1. )   ! I-point 
    518514      END SELECT 
    519515 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcblk_core.F90

    r2291 r2370  
    7575   !! NEMO/OPA 3.3 , NEMO-consortium (2010)  
    7676   !! $Id$ 
    77    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     77   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    7878   !!---------------------------------------------------------------------- 
    7979CONTAINS 
     
    133133         ! 
    134134         ! (NB: frequency positive => hours, negative => months) 
    135          !            !    file     ! frequency !  variable  ! time intep !  clim   ! 'yearly' or ! weights  ! rotation  ! 
    136          !            !    name     !  (hours)  !   name     !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs      ! 
    137          sn_wndi = FLD_N( 'uwnd10m' ,    24     ,  'u_10'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         ) 
    138          sn_wndj = FLD_N( 'vwnd10m' ,    24     ,  'v_10'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         ) 
    139          sn_qsr  = FLD_N( 'qsw'     ,    24     ,  'qsw'     ,  .false.   , .false. ,   'yearly'  , ''       , ''         ) 
    140          sn_qlw  = FLD_N( 'qlw'     ,    24     ,  'qlw'     ,  .false.   , .false. ,   'yearly'  , ''       , ''         ) 
    141          sn_tair = FLD_N( 'tair10m' ,    24     ,  't_10'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         ) 
    142          sn_humi = FLD_N( 'humi10m' ,    24     ,  'q_10'    ,  .false.   , .false. ,   'yearly'  , ''       , ''         ) 
    143          sn_prec = FLD_N( 'precip'  ,    -1     ,  'precip'  ,  .true.    , .false. ,   'yearly'  , ''       , ''         ) 
    144          sn_snow = FLD_N( 'snow'    ,    -1     ,  'snow'    ,  .true.    , .false. ,   'yearly'  , ''       , ''         ) 
    145          sn_tdif = FLD_N( 'taudif'  ,    24     ,  'taudif'  ,  .true.    , .false. ,   'yearly'  , ''       , ''         ) 
     135         !            !    file    ! frequency ! variable ! time intep !  clim   ! 'yearly' or ! weights  ! rotation ! 
     136         !            !    name    !  (hours)  !  name    !   (T/F)    !  (T/F)  !  'monthly'  ! filename ! pairs    ! 
     137         sn_wndi = FLD_N( 'uwnd10m',    24     , 'u_10'   ,  .false.   , .false. ,   'yearly'  , ''       , ''       ) 
     138         sn_wndj = FLD_N( 'vwnd10m',    24     , 'v_10'   ,  .false.   , .false. ,   'yearly'  , ''       , ''       ) 
     139         sn_qsr  = FLD_N( 'qsw'    ,    24     , 'qsw'    ,  .false.   , .false. ,   'yearly'  , ''       , ''       ) 
     140         sn_qlw  = FLD_N( 'qlw'    ,    24     , 'qlw'    ,  .false.   , .false. ,   'yearly'  , ''       , ''       ) 
     141         sn_tair = FLD_N( 'tair10m',    24     , 't_10'   ,  .false.   , .false. ,   'yearly'  , ''       , ''       ) 
     142         sn_humi = FLD_N( 'humi10m',    24     , 'q_10'   ,  .false.   , .false. ,   'yearly'  , ''       , ''       ) 
     143         sn_prec = FLD_N( 'precip' ,    -1     , 'precip' ,  .true.    , .false. ,   'yearly'  , ''       , ''       ) 
     144         sn_snow = FLD_N( 'snow'   ,    -1     , 'snow'   ,  .true.    , .false. ,   'yearly'  , ''       , ''       ) 
     145         sn_tdif = FLD_N( 'taudif' ,    24     , 'taudif' ,  .true.    , .false. ,   'yearly'  , ''       , ''       ) 
    146146         ! 
    147147         REWIND( numnam )                          ! read in namlist namsbc_core 
     
    396396      !! caution : the net upward water flux has with mm/day unit 
    397397      !!--------------------------------------------------------------------- 
    398       REAL(wp), INTENT(in   ), DIMENSION(:,:,:)      ::   pst      ! ice surface temperature (>0, =rt0 over land) [Kelvin] 
    399       REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)    ::   pui      ! ice surface velocity (i- and i- components      [m/s] 
    400       REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj)    ::   pvi      !    at I-point (B-grid) or U & V-point (C-grid) 
    401       REAL(wp), INTENT(in   ), DIMENSION(:,:,:)      ::   palb     ! ice albedo (clear sky) (alb_ice_cs)               [%] 
    402       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj)    ::   p_taui   ! i- & j-components of surface ice stress        [N/m2] 
    403       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj)    ::   p_tauj   !   at I-point (B-grid) or U & V-point (C-grid) 
    404       REAL(wp), INTENT(  out), DIMENSION(:,:,:)      ::   p_qns    ! non solar heat flux over ice (T-point)         [W/m2] 
    405       REAL(wp), INTENT(  out), DIMENSION(:,:,:)      ::   p_qsr    !     solar heat flux over ice (T-point)         [W/m2] 
    406       REAL(wp), INTENT(  out), DIMENSION(:,:,:)      ::   p_qla    ! latent    heat flux over ice (T-point)         [W/m2] 
    407       REAL(wp), INTENT(  out), DIMENSION(:,:,:)      ::   p_dqns   ! non solar heat sensistivity  (T-point)         [W/m2] 
    408       REAL(wp), INTENT(  out), DIMENSION(:,:,:)      ::   p_dqla   ! latent    heat sensistivity  (T-point)         [W/m2] 
    409       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj)    ::   p_tpr    ! total precipitation          (T-point)      [Kg/m2/s] 
    410       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj)    ::   p_spr    ! solid precipitation          (T-point)      [Kg/m2/s] 
    411       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj)    ::   p_fr1    ! 1sr fraction of qsr penetration in ice (T-point)  [%] 
    412       REAL(wp), INTENT(  out), DIMENSION(jpi,jpj)    ::   p_fr2    ! 2nd fraction of qsr penetration in ice (T-point)  [%] 
    413       CHARACTER(len=1), INTENT(in   )                ::   cd_grid  ! ice grid ( C or B-grid) 
    414       INTEGER, INTENT(in   )                         ::   pdim     ! number of ice categories 
     398      REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   pst      ! ice surface temperature (>0, =rt0 over land) [Kelvin] 
     399      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj) ::   pui      ! ice surface velocity (i- and i- components      [m/s] 
     400      REAL(wp), INTENT(in   ), DIMENSION(jpi,jpj) ::   pvi      !    at I-point (B-grid) or U & V-point (C-grid) 
     401      REAL(wp), INTENT(in   ), DIMENSION(:,:,:)   ::   palb     ! ice albedo (clear sky) (alb_ice_cs)               [%] 
     402      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_taui   ! i- & j-components of surface ice stress        [N/m2] 
     403      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_tauj   !   at I-point (B-grid) or U & V-point (C-grid) 
     404      REAL(wp), INTENT(  out), DIMENSION(:,:,:)   ::   p_qns    ! non solar heat flux over ice (T-point)         [W/m2] 
     405      REAL(wp), INTENT(  out), DIMENSION(:,:,:)   ::   p_qsr    !     solar heat flux over ice (T-point)         [W/m2] 
     406      REAL(wp), INTENT(  out), DIMENSION(:,:,:)   ::   p_qla    ! latent    heat flux over ice (T-point)         [W/m2] 
     407      REAL(wp), INTENT(  out), DIMENSION(:,:,:)   ::   p_dqns   ! non solar heat sensistivity  (T-point)         [W/m2] 
     408      REAL(wp), INTENT(  out), DIMENSION(:,:,:)   ::   p_dqla   ! latent    heat sensistivity  (T-point)         [W/m2] 
     409      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_tpr    ! total precipitation          (T-point)      [Kg/m2/s] 
     410      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_spr    ! solid precipitation          (T-point)      [Kg/m2/s] 
     411      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_fr1    ! 1sr fraction of qsr penetration in ice (T-point)  [%] 
     412      REAL(wp), INTENT(  out), DIMENSION(jpi,jpj) ::   p_fr2    ! 2nd fraction of qsr penetration in ice (T-point)  [%] 
     413      CHARACTER(len=1), INTENT(in   )             ::   cd_grid  ! ice grid ( C or B-grid) 
     414      INTEGER, INTENT(in   )                      ::   pdim     ! number of ice categories 
    415415      !! 
    416416      INTEGER  ::   ji, jj, jl    ! dummy loop indices 
     
    449449      ! ----------------------------------------------------------------------------- ! 
    450450      SELECT CASE( cd_grid ) 
    451       CASE( 'B' )                  ! B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
     451      CASE( 'I' )                  ! B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
    452452         !                           and scalar wind at T-point ( = | U10m - U_ice | ) (masked) 
    453453#if defined key_vectopt_loop 
     
    605605      !!   9.0  !  05-08  (L. Brodeau) Rewriting and optimization 
    606606      !!---------------------------------------------------------------------- 
    607       !! * Arguments 
    608  
    609607      REAL(wp), INTENT(in) :: zu                 ! altitude of wind measurement       [m] 
    610608      REAL(wp), INTENT(in), DIMENSION(jpi,jpj) ::  & 
     
    646644         grav   = 9.8,          &  ! gravity                        
    647645         kappa  = 0.4              ! von Karman s constant 
    648  
     646      !!---------------------------------------------------------------------- 
    649647      !! * Start 
    650648      !! Air/sea differences 
     
    770768         grav   = 9.8,      &  ! gravity                        
    771769         kappa  = 0.4          ! von Karman's constant 
    772  
     770      !!---------------------------------------------------------------------- 
    773771      !!  * Start 
    774772 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbccpl.F90

    r2287 r2370  
    44   !! Surface Boundary Condition :  momentum, heat and freshwater fluxes in coupled mode 
    55   !!====================================================================== 
    6    !! History :  2.0  !  06-2007  (R. Redler, N. Keenlyside, W. Park) Original code split into flxmod & taumod 
    7    !!            3.0  !  02-2008  (G. Madec, C Talandier)  surface module 
    8    !!            3.1  !  02-2009  (G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface 
     6   !! History :  2.0  ! 2007-06  (R. Redler, N. Keenlyside, W. Park) Original code split into flxmod & taumod 
     7   !!            3.0  ! 2008-02  (G. Madec, C Talandier)  surface module 
     8   !!            3.1  ! 2009_02  (G. Madec, S. Masson, E. Maisonave, A. Caubel) generic coupled interface 
    99   !!---------------------------------------------------------------------- 
    1010#if defined key_oasis3 || defined key_oasis4 
     
    788788      !! ** Method  :   transform the received stress from the atmosphere into 
    789789      !!             an atmosphere-ice stress in the (i,j) ocean referencial 
    790       !!             and at the velocity point of the sea-ice model (cigr_type): 
     790      !!             and at the velocity point of the sea-ice model (cp_ice_msh): 
    791791      !!                'C'-grid : i- (j-) components given at U- (V-) point  
    792       !!                'B'-grid : both components given at I-point  
     792      !!                'I'-grid : B-grid lower-left corner: both components given at I-point  
    793793      !! 
    794794      !!                The received stress are : 
     
    803803      !!                 first  as  2 components on the sphere  
    804804      !!                 second as  2 components oriented along the local grid 
    805       !!                 third  as  2 components on the cigr_type point  
     805      !!                 third  as  2 components on the cp_ice_msh point  
    806806      !! 
    807807      !!                In 'oce and ice' case, only one vector stress field  
     
    809809      !!             so that it is now defined as (i,j) components given at U- 
    810810      !!             and V-points, respectively. Therefore, here only the third 
    811       !!             transformation is done and only if the ice-grid is a 'B'-grid.  
    812       !! 
    813       !! ** Action  :   return ptau_i, ptau_j, the stress over the ice at cigr_type point 
     811      !!             transformation is done and only if the ice-grid is a 'I'-grid.  
     812      !! 
     813      !! ** Action  :   return ptau_i, ptau_j, the stress over the ice at cp_ice_msh point 
    814814      !!---------------------------------------------------------------------- 
    815815      REAL(wp), INTENT(out), DIMENSION(jpi,jpj) ::   p_taui   ! i- & j-components of atmos-ice stress [N/m2] 
     
    872872         !     
    873873         !                                                  j+1   j     -----V---F 
    874          ! ice stress on ice velocity point (cigr_type)                  !       | 
     874         ! ice stress on ice velocity point (cp_ice_msh)                 !       | 
    875875         ! (C-grid ==>(U,V) or B-grid ==> I or F)                 j      |   T   U 
    876876         !                                                               |       | 
     
    879879         !                                                              i-1  i   i 
    880880         !                                                               i      i+1 (for I) 
    881          SELECT CASE ( cigr_type ) 
     881         SELECT CASE ( cp_ice_msh ) 
    882882            ! 
    883883         CASE( 'I' )                                         ! B-grid ==> I 
     
    12581258            END DO 
    12591259         CASE( 'weighted oce and ice' )    
    1260             SELECT CASE ( cigr_type ) 
     1260            SELECT CASE ( cp_ice_msh ) 
    12611261            CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T 
    12621262               DO jj = 2, jpjm1 
     
    12931293            CALL lbc_lnk( zitx1, 'T', -1. )   ;   CALL lbc_lnk( zity1, 'T', -1. ) 
    12941294         CASE( 'mixed oce-ice'        ) 
    1295             SELECT CASE ( cigr_type ) 
     1295            SELECT CASE ( cp_ice_msh ) 
    12961296            CASE( 'C' )                      ! Ocean and Ice on C-grid ==> T 
    12971297               DO jj = 2, jpjm1 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim.F90

    r2287 r2370  
    55   !!       &           covered area using LIM sea-ice model 
    66   !! Sea-Ice model  :  LIM 3.0 Sea ice model time-stepping 
    7    !!====================================================================== 
    8    !! History :  2.0   !  2006-12  (M. Vancoppenolle) Original code 
    9    !!            3.0   !  2008-02  (C. Talandier)  Surface module from icestp.F90 
    10    !!            9.0   !  2008-04  (G. Madec)  sltyle and lim_ctl routine 
     7   !!===================================================================== 
     8   !! History :  2.0  ! 2006-12  (M. Vancoppenolle) Original code 
     9   !!            3.0  ! 2008-02  (C. Talandier)  Surface module from icestp.F90 
     10   !!             -   ! 2008-04  (G. Madec)  sltyle and lim_ctl routine 
     11   !!            3.3  ! 2010-11  (G. Madec) ice-ocean stress always computed at each ocean time-step 
    1112   !!---------------------------------------------------------------------- 
    1213#if defined key_lim3 
     
    2021   USE oce             ! ocean dynamics and tracers 
    2122   USE dom_oce         ! ocean space and time domain 
    22    USE lib_mpp 
     23   USE lib_mpp         ! MPP library 
    2324   USE par_ice         ! sea-ice parameters 
    24    USE ice 
    25    USE iceini 
    26    USE dom_ice 
     25   USE ice             ! LIM-3: ice variables 
     26   USE iceini          ! LIM-3: ice initialisation 
     27   USE dom_ice         ! LIM-3: ice domain 
    2728 
    2829   USE sbc_oce         ! Surface boundary condition: ocean fields 
     
    3031   USE sbcblk_core     ! Surface boundary condition: CORE bulk 
    3132   USE sbcblk_clio     ! Surface boundary condition: CLIO bulk 
    32    USE albedo 
     33   USE albedo          ! ocean & ice albedo 
    3334 
    3435   USE phycst          ! Define parameters for the routines 
     
    4647   USE limvar          ! Ice variables switch 
    4748 
    48    USE lbclnk 
     49   USE lbclnk          ! lateral boundary condition - MPP link 
    4950   USE iom             ! I/O manager library 
    5051   USE in_out_manager  ! I/O manager 
     
    5657   PUBLIC sbc_ice_lim  ! routine called by sbcmod.F90 
    5758    
    58    CHARACTER(len=1) ::   cl_grid = 'C'     ! type of grid used in ice dynamics 
    59  
    6059   !! * Substitutions 
    6160#  include "domzgr_substitute.h90" 
     
    6463   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    6564   !! $Id$ 
    66    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    67    !!---------------------------------------------------------------------- 
    68  
     65   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     66   !!---------------------------------------------------------------------- 
    6967CONTAINS 
    7068 
    71    SUBROUTINE sbc_ice_lim( kt, kblk, kico ) 
     69   SUBROUTINE sbc_ice_lim( kt, kblk ) 
    7270      !!--------------------------------------------------------------------- 
    7371      !!                  ***  ROUTINE sbc_ice_lim  *** 
     
    9189      INTEGER, INTENT(in) ::   kt      ! ocean time step 
    9290      INTEGER, INTENT(in) ::   kblk    ! type of bulk (=3 CLIO, =4 CORE) 
    93       INTEGER, INTENT(in) ::   kico    ! ice-ocean stress treatment 
    9491      !! 
    9592      INTEGER  ::   jl                 ! loop index 
     
    142139               &                      qla_ice   , dqns_ice  , dqla_ice  ,               & 
    143140               &                      tprecip   , sprecip   ,                           & 
    144                &                      fr1_i0    , fr2_i0    , cl_grid, jpl  ) 
     141               &                      fr1_i0    , fr2_i0    , cp_ice_msh, jpl  ) 
    145142            !          
    146143         CASE( 4 )                                       ! CORE bulk formulation 
     
    149146               &                      qla_ice   , dqns_ice  , dqla_ice  ,               & 
    150147               &                      tprecip   , sprecip   ,                           & 
    151                &                      fr1_i0    , fr2_i0    , cl_grid, jpl  ) 
     148               &                      fr1_i0    , fr2_i0    , cp_ice_msh, jpl  ) 
    152149         END SELECT 
    153150 
     
    160157         !                                           ! Store previous ice values 
    161158!!gm : remark   old_...   should becomes ...b  as tn versus tb   
    162          old_a_i(:,:,:)   = a_i(:,:,:)     ! ice area 
    163          old_e_i(:,:,:,:) = e_i(:,:,:,:)   ! ice thermal energy 
    164          old_v_i(:,:,:)   = v_i(:,:,:)     ! ice volume 
    165          old_v_s(:,:,:)   = v_s(:,:,:)     ! snow volume  
    166          old_e_s(:,:,:,:) = e_s(:,:,:,:)   ! snow thermal energy 
    167          old_smv_i(:,:,:) = smv_i(:,:,:)   ! salt content 
    168          old_oa_i(:,:,:)  = oa_i(:,:,:)    ! areal age content 
     159         old_a_i  (:,:,:)   = a_i  (:,:,:)     ! ice area 
     160         old_e_i  (:,:,:,:) = e_i  (:,:,:,:)   ! ice thermal energy 
     161         old_v_i  (:,:,:)   = v_i  (:,:,:)     ! ice volume 
     162         old_v_s  (:,:,:)   = v_s  (:,:,:)     ! snow volume  
     163         old_e_s  (:,:,:,:) = e_s  (:,:,:,:)   ! snow thermal energy 
     164         old_smv_i(:,:,:)   = smv_i(:,:,:)     ! salt content 
     165         old_oa_i (:,:,:)   = oa_i (:,:,:)     ! areal age content 
    169166 
    170167         !                                           ! intialisation to zero    !!gm is it truly necessary ??? 
    171          d_a_i_thd(:,:,:)   = 0.e0 ; d_a_i_trp(:,:,:)   = 0.e0 
    172          d_v_i_thd(:,:,:)   = 0.e0 ; d_v_i_trp(:,:,:)   = 0.e0 
    173          d_e_i_thd(:,:,:,:) = 0.e0 ; d_e_i_trp(:,:,:,:) = 0.e0 
    174          d_v_s_thd(:,:,:)   = 0.e0 ; d_v_s_trp(:,:,:)   = 0.e0 
    175          d_e_s_thd(:,:,:,:) = 0.e0 ; d_e_s_trp(:,:,:,:) = 0.e0 
    176          d_smv_i_thd(:,:,:) = 0.e0 ; d_smv_i_trp(:,:,:) = 0.e0 
    177          d_oa_i_thd(:,:,:)  = 0.e0 ; d_oa_i_trp(:,:,:)  = 0.e0 
    178          ! 
    179          fseqv(:,:)    = 0.e0 
    180          fsbri(:,:)     = 0.e0     ; fsalt_res(:,:) = 0.e0 
     168         d_a_i_thd  (:,:,:)   = 0.e0   ;   d_a_i_trp  (:,:,:)   = 0.e0 
     169         d_v_i_thd  (:,:,:)   = 0.e0   ;   d_v_i_trp  (:,:,:)   = 0.e0 
     170         d_e_i_thd  (:,:,:,:) = 0.e0   ;   d_e_i_trp  (:,:,:,:) = 0.e0 
     171         d_v_s_thd  (:,:,:)   = 0.e0   ;   d_v_s_trp  (:,:,:)   = 0.e0 
     172         d_e_s_thd  (:,:,:,:) = 0.e0   ;   d_e_s_trp  (:,:,:,:) = 0.e0 
     173         d_smv_i_thd(:,:,:)   = 0.e0   ;   d_smv_i_trp(:,:,:)  = 0.e0 
     174         d_oa_i_thd (:,:,:)   = 0.e0   ;   d_oa_i_trp (:,:,:)   = 0.e0 
     175         ! 
     176         fseqv    (:,:) = 0.e0 
     177         fsbri    (:,:) = 0.e0     ;  fsalt_res(:,:) = 0.e0 
    181178         fsalt_rpo(:,:) = 0.e0 
    182          fhmec(:,:)     = 0.e0     ; fhbri(:,:)    = 0.e0 
    183          fmmec(:,:)     = 0.e0     ; fheat_res(:,:) = 0.e0 
    184          fheat_rpo(:,:) = 0.e0     ; focea2D(:,:)  = 0.e0 
    185          fsup2D(:,:)    = 0.e0 
     179         fhmec    (:,:) = 0.e0     ;   fhbri    (:,:) = 0.e0 
     180         fmmec    (:,:) = 0.e0     ;  fheat_res(:,:) = 0.e0 
     181         fheat_rpo(:,:) = 0.e0     ;   focea2D  (:,:) = 0.e0 
     182         fsup2D   (:,:) = 0.e0 
    186183         !  
    187          diag_sni_gr(:,:) = 0.e0   ; diag_lat_gr(:,:) = 0.e0 
    188          diag_bot_gr(:,:) = 0.e0   ; diag_dyn_gr(:,:) = 0.e0 
    189          diag_bot_me(:,:) = 0.e0   ; diag_sur_me(:,:) = 0.e0 
     184         diag_sni_gr(:,:) = 0.e0   ;   diag_lat_gr(:,:) = 0.e0 
     185         diag_bot_gr(:,:) = 0.e0   ;   diag_dyn_gr(:,:) = 0.e0 
     186         diag_bot_me(:,:) = 0.e0   ;   diag_sur_me(:,:) = 0.e0 
    190187         ! dynamical invariants 
    191          delta_i(:,:) = 0.e0       ; divu_i (:,:) = 0.e0       ;    shear_i(:,:) = 0.e0 
     188         delta_i(:,:) = 0.e0       ;   divu_i(:,:) = 0.e0       ;   shear_i(:,:) = 0.e0 
    192189 
    193190                          CALL lim_rst_opn( kt )     ! Open Ice restart file 
     
    196193         ! 
    197194#if ! defined key_c1d 
    198          ! Ice dynamics & transport (not in 1D case) 
     195                                                     ! Ice dynamics & transport (not in 1D case) 
    199196                          CALL lim_dyn( kt )              ! Ice dynamics    ( rheology/dynamics ) 
    200197                          CALL lim_trp( kt )              ! Ice transport   ( Advection/diffusion ) 
     
    204201                          CALL lim_itd_me                 ! Mechanical redistribution ! (ridging/rafting) 
    205202#endif 
    206          ! 
    207203         !                                           ! Ice thermodynamics  
    208204                          CALL lim_var_glo2eqv            ! equivalent variables 
     
    216212                          CALL lim_itd_th( kt )           !  Remap ice categories, lateral accretion  ! 
    217213         ! 
    218          !                                           ! Global variables update | 
     214         !                                           ! Global variables update 
    219215                          CALL lim_var_agg( 1 )           ! requested by limupdate 
    220216                          CALL lim_update                 ! Global variables update 
     
    223219         IF( ln_nicep )   CALL lim_prt_state( jiindx, jjindx, 2, ' - Final state - ' )   ! control print 
    224220         ! 
    225          !                                           ! Fluxes of mass and heat to the ocean | 
    226                          CALL lim_sbc_flx( kt )           ! Ice/Ocean heat freshwater/salt fluxes 
    227          IF( ln_limdyn .AND. kico == 0 )   &              ! Ice/Ocean stresses (only in ice-dynamic case) 
    228             &            CALL lim_sbc_tau( kt, kico )     ! otherwise the atm.-ocean stresses are used everywhere 
     221                          CALL lim_sbc_flx( kt )     ! Update surface ocean mass, heat and salt fluxes 
    229222         ! 
    230223         IF( ln_nicep )   CALL lim_prt_state( jiindx, jjindx, 3, ' - Final state lim_sbc - ' )   ! control print 
     
    239232         IF( ln_nicep )   CALL lim_ctl               ! alerts in case of model crash 
    240233         ! 
    241       ENDIF                                          ! End sea-ice time step only 
    242  
    243       !                                              !--------------------------! 
    244       ! Ice/Ocean stresses (nn_ico_cpl=1 or 2 cases) !  at all ocean time step  ! 
    245       !                                              !--------------------------! 
    246       IF( ln_limdyn .AND. kico /= 0 )   & 
    247          &                CALL lim_sbc_tau( kt, kico )  
    248 !!gm   remark, in this case the ocean-ice stress is not saved in diag call above .....  find a solution!!! 
     234      ENDIF                                    ! End sea-ice time step only 
     235 
     236      !                                        !--------------------------! 
     237      !                                        !  at all ocean time step  ! 
     238      !                                        !--------------------------! 
     239      !                                                
     240      !                                              ! Update surface ocean stresses (only in ice-dynamic case) 
     241      !                                                   ! otherwise the atm.-ocean stresses are used everywhere 
     242      IF( ln_limdyn )     CALL lim_sbc_tau( kt, ub(:,:,1), vb(:,:,1) )  ! using before instantaneous surf. currents 
     243       
     244!!gm   remark, the ocean-ice stress is not saved in ice diag call above .....  find a solution!!! 
    249245      ! 
    250246   END SUBROUTINE sbc_ice_lim 
     
    664660   !!---------------------------------------------------------------------- 
    665661CONTAINS 
    666    SUBROUTINE sbc_ice_lim ( kt, kblk, kico )     ! Dummy routine 
    667       WRITE(*,*) 'sbc_ice_lim: You should not have seen this print! error?', kt, kblk, kico 
     662   SUBROUTINE sbc_ice_lim ( kt, kblk )     ! Dummy routine 
     663      WRITE(*,*) 'sbc_ice_lim: You should not have seen this print! error?', kt, kblk 
    668664   END SUBROUTINE sbc_ice_lim 
    669665#endif 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcice_lim_2.F90

    r2319 r2370  
    22   !!====================================================================== 
    33   !!                       ***  MODULE  sbcice_lim_2  *** 
    4    !! Surface module :  update surface ocean boundary condition over ice 
    5    !!                   covered area using LIM sea-ice model 
    6    !! Sea-Ice model  :  LIM 2.0 Sea ice model time-stepping 
     4   !! Surface module :  update surface ocean boundary condition over ice covered area using LIM sea-ice model 
     5   !! Sea-Ice model  :  LIM-2 Sea ice model time-stepping 
    76   !!====================================================================== 
    87   !! History :  1.0   !  06-2006  (G. Madec)  from icestp_2.F90 
     
    1211#if defined key_lim2 
    1312   !!---------------------------------------------------------------------- 
    14    !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    15    !!---------------------------------------------------------------------- 
    16    !!   sbc_ice_lim_2  : sea-ice model time-stepping and 
    17    !!                    update ocean sbc over ice-covered area 
    18    !!---------------------------------------------------------------------- 
    19    USE oce             ! ocean dynamics and tracers 
    20    USE dom_oce         ! ocean space and time domain 
    21    USE lib_mpp 
     13   !!   'key_lim2' :                                    LIM-2 sea-ice model 
     14   !!---------------------------------------------------------------------- 
     15   !!   sbc_ice_lim_2   : sea-ice model time-stepping and update ocean sbc over ice-covered area 
     16   !!---------------------------------------------------------------------- 
     17   USE oce              ! ocean dynamics and tracers 
     18   USE dom_oce          ! ocean space and time domain 
    2219   USE ice_2 
    2320   USE par_ice_2 
     
    2522   USE dom_ice_2 
    2623 
    27    USE sbc_oce         ! Surface boundary condition: ocean fields 
    28    USE sbc_ice         ! Surface boundary condition: ice   fields 
    29    USE sbcblk_core     ! Surface boundary condition: CORE bulk 
    30    USE sbcblk_clio     ! Surface boundary condition: CLIO bulk 
    31    USE sbccpl          ! Surface boundary condition: coupled interface 
     24   USE sbc_oce          ! Surface boundary condition: ocean fields 
     25   USE sbc_ice          ! Surface boundary condition: ice   fields 
     26   USE sbcblk_core      ! Surface boundary condition: CORE bulk 
     27   USE sbcblk_clio      ! Surface boundary condition: CLIO bulk 
     28   USE sbccpl           ! Surface boundary condition: coupled interface 
    3229   USE albedo 
    3330 
    34    USE phycst          ! Define parameters for the routines 
    35    USE eosbn2          ! equation of state 
     31   USE phycst           ! Define parameters for the routines 
     32   USE eosbn2           ! equation of state 
    3633   USE limdyn_2 
    3734   USE limtrp_2 
    3835   USE limdmp_2 
    3936   USE limthd_2 
    40    USE limsbc_2        ! sea surface boundary condition 
     37   USE limsbc_2         ! sea surface boundary condition 
    4138   USE limdia_2 
    4239   USE limwri_2 
    4340   USE limrst_2 
    4441 
    45    USE lbclnk 
    46    USE iom             ! I/O manager library 
    47    USE in_out_manager  ! I/O manager 
    48    USE prtctl          ! Print control 
     42   USE lbclnk           ! lateral boundary condition - MPP link 
     43   USE lib_mpp          ! MPP library 
     44   USE iom              ! I/O manager library 
     45   USE in_out_manager   ! I/O manager 
     46   USE prtctl           ! Print control 
    4947 
    5048   IMPLICIT NONE 
     
    5957   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    6058   !! $Id$ 
    61    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    62    !!---------------------------------------------------------------------- 
    63  
     59   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     60   !!---------------------------------------------------------------------- 
    6461CONTAINS 
    6562 
     
    9794         IF(lwp) WRITE(numout,*) 'sbc_ice_lim_2 : update ocean surface boudary condition'  
    9895         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~   via Louvain la Neuve Ice Model (LIM) time stepping' 
    99  
     96         ! 
    10097         CALL ice_init_2 
    101  
    10298      ENDIF 
    10399 
    104       IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN 
    105          ! 
     100      !                                        !----------------------! 
     101      IF( MOD( kt-1, nn_fsbc ) == 0 ) THEN     !  Ice time-step only  ! 
     102         !                                     !----------------------! 
     103         !  Bulk Formulea ! 
     104         !----------------! 
    106105         ! ... mean surface ocean current at ice dynamics point 
    107          !     B-grid dynamics :  I-point  
    108          DO jj = 2, jpj 
    109             DO ji = 2, jpi   ! B grid : no vector opt. 
    110                u_oce(ji,jj) = 0.5 * ( ssu_m(ji-1,jj  ) + ssu_m(ji-1,jj-1) ) * tmu(ji,jj) 
    111                v_oce(ji,jj) = 0.5 * ( ssv_m(ji  ,jj-1) + ssv_m(ji-1,jj-1) ) * tmu(ji,jj) 
     106         SELECT CASE( cp_ice_msh ) 
     107         CASE( 'I' )                  !== B-grid ice dynamics :   I-point (i.e. F-point with sea-ice indexation) 
     108            DO jj = 2, jpj 
     109               DO ji = 2, jpi   ! NO vector opt. possible 
     110                  u_oce(ji,jj) = 0.5_wp * ( ssu_m(ji-1,jj  ) + ssu_m(ji-1,jj-1) ) * tmu(ji,jj) 
     111                  v_oce(ji,jj) = 0.5_wp * ( ssv_m(ji  ,jj-1) + ssv_m(ji-1,jj-1) ) * tmu(ji,jj) 
     112               END DO 
    112113            END DO 
    113          END DO 
    114          CALL lbc_lnk( u_oce, 'I', -1. )   ! I-point (i.e. F-point with ice indices) 
    115          CALL lbc_lnk( v_oce, 'I', -1. )   ! I-point (i.e. F-point with ice indices) 
     114            CALL lbc_lnk( u_oce, 'I', -1. )   ! I-point (i.e. F-point with ice indices) 
     115            CALL lbc_lnk( v_oce, 'I', -1. )   ! I-point (i.e. F-point with ice indices) 
     116            ! 
     117         CASE( 'C' )                  !== C-grid ice dynamics :   U & V-points (same as ocean) 
     118            u_oce(:,:) = ssu_m(:,:)                     ! mean surface ocean current at ice velocity point 
     119            v_oce(:,:) = ssv_m(:,:) 
     120            ! 
     121         END SELECT 
    116122 
    117123         ! ... masked sea surface freezing temperature [Kelvin] (set to rt0 over land) 
     
    142148               &                      qla_ice    , dqns_ice   , dqla_ice   ,            & 
    143149               &                      tprecip    , sprecip    ,                         & 
    144                &                      fr1_i0     , fr2_i0     , cl_grid    , jpl  ) 
     150               &                      fr1_i0     , fr2_i0     , cp_ice_msh , jpl  ) 
    145151 
    146152         CASE( 4 )           ! CORE bulk formulation 
     
    149155               &                      qla_ice    , dqns_ice   , dqla_ice   ,            & 
    150156               &                      tprecip    , sprecip    ,                         & 
    151                &                      fr1_i0     , fr2_i0     , cl_grid    , jpl  ) 
     157               &                      fr1_i0     , fr2_i0     , cp_ice_msh , jpl  ) 
    152158         CASE( 5 )           ! Coupled formulation : atmosphere-ice stress only (fluxes provided after ice dynamics) 
    153159            CALL sbc_cpl_ice_tau( utau_ice , vtau_ice ) 
     
    170176         !  Ice model step  ! 
    171177         ! ---------------- ! 
    172          numit = numit + nn_fsbc                                             !Ice model time step 
    173  
    174                                         CALL lim_rst_opn_2  ( kt )           ! Open Ice restart file 
     178         numit = numit + nn_fsbc                           ! Ice model time step 
     179 
     180                           CALL lim_rst_opn_2  ( kt )  ! Open Ice restart file 
    175181#if ! defined key_c1d 
    176             ! Ice dynamics & transport (not in 1D case) 
    177                                         CALL lim_dyn_2      ( kt )           ! Ice dynamics    ( rheology/dynamics ) 
    178                                         CALL lim_trp_2      ( kt )           ! Ice transport   ( Advection/diffusion ) 
    179             IF( ln_limdmp )             CALL lim_dmp_2      ( kt )           ! Ice damping  
     182                                                       ! Ice dynamics & transport (except in 1D case) 
     183                           CALL lim_dyn_2      ( kt )      ! Ice dynamics    ( rheology/dynamics ) 
     184                           CALL lim_trp_2      ( kt )      ! Ice transport   ( Advection/diffusion ) 
     185         IF( ln_limdmp )   CALL lim_dmp_2      ( kt )      ! Ice damping  
    180186#endif 
    181187#if defined key_coupled 
    182          IF( ksbc == 5    )             CALL sbc_cpl_ice_flx( reshape( frld, (/jpi,jpj,1/) ),        & 
    183       &                                                       qns_tot, qns_ice, qsr_tot , qsr_ice,   & 
    184       &                                                       emp_tot, emp_ice, dqns_ice, sprecip,   & 
     188         !                                             ! Ice surface fluxes in coupled mode  
     189         IF( ksbc == 5 )   CALL sbc_cpl_ice_flx( reshape( frld, (/jpi,jpj,1/) ),                 & 
     190      &                                                   qns_tot, qns_ice, qsr_tot , qsr_ice,   & 
     191      &                                                   emp_tot, emp_ice, dqns_ice, sprecip,   & 
    185192      !                                      optional arguments, used only in 'mixed oce-ice' case 
    186       &                                                       palbi = zalb_ice_cs, psst = sst_m, pist = sist ) 
     193      &                                                   palbi = zalb_ice_cs, psst = sst_m, pist = sist ) 
    187194#endif 
    188                                         CALL lim_thd_2      ( kt )      ! Ice thermodynamics  
    189                                         CALL lim_sbc_2      ( kt )      ! Ice/Ocean Mass & Heat fluxes  
     195                           CALL lim_thd_2      ( kt )      ! Ice thermodynamics  
     196                           CALL lim_sbc_flx_2  ( kt )      ! update surface ocean mass, heat & salt fluxes  
    190197 
    191198         IF( ( MOD( kt+nn_fsbc-1, ninfo ) == 0 .OR. ntmoy == 1 ) .AND. .NOT. lk_mpp )   & 
    192             &                           CALL lim_dia_2      ( kt )      ! Ice Diagnostics 
     199            &              CALL lim_dia_2      ( kt )      ! Ice Diagnostics 
    193200# if ! defined key_iomput 
    194                                         CALL lim_wri_2      ( kt )      ! Ice outputs 
     201                           CALL lim_wri_2      ( kt )      ! Ice outputs 
    195202# endif 
    196          IF( lrst_ice )                 CALL lim_rst_write_2( kt )      ! Ice restart file 
     203         IF( lrst_ice  )   CALL lim_rst_write_2( kt )      ! Ice restart file 
    197204         ! 
    198       ENDIF 
     205      ENDIF                                    ! End sea-ice time step only 
     206      ! 
     207      !                                        !--------------------------! 
     208      !                                        !  at all ocean time step  ! 
     209      !                                        !--------------------------! 
     210      !                                                
     211      !                                              ! Update surface ocean stresses (only in ice-dynamic case) 
     212      !                                                   ! otherwise the atm.-ocean stresses are used everywhere 
     213      IF( ln_limdyn    )   CALL lim_sbc_tau_2( kt, ub(:,:,1), vb(:,:,1) )  ! using before instantaneous surf. currents 
    199214      ! 
    200215   END SUBROUTINE sbc_ice_lim_2 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcmod.F90

    r2305 r2370  
    44   !! Surface module :  provide to the ocean its surface boundary condition 
    55   !!====================================================================== 
    6    !! History :  3.0  !  2006-07  (G. Madec)  Original code 
    7    !!            3.1  !  2008-08  (S. Masson, A. Caubel, E. Maisonnave, G. Madec) coupled interface 
    8    !!            3.3  !  2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps 
    9    !!            3.3  !  2010-10  (S. Masson)  add diurnal cycle 
    10    !!            3.3   !  09-2010  (D. Storkey) add ice boundary conditions (BDY) 
     6   !! History :  3.0  ! 2006-07  (G. Madec)  Original code 
     7   !!            3.1  ! 2008-08  (S. Masson, A. Caubel, E. Maisonnave, G. Madec) coupled interface 
     8   !!            3.3  ! 2010-04  (M. Leclair, G. Madec)  Forcing averaged over 2 time steps 
     9   !!            3.3  ! 2010-10  (S. Masson)  add diurnal cycle 
     10   !!            3.3  ! 2010-09  (D. Storkey) add ice boundary conditions (BDY) 
     11   !!             -   ! 2010-11  (G. Madec) ice-ocean stress always computed at each ocean time-step 
    1112   !!---------------------------------------------------------------------- 
    1213 
     
    7475      !! 
    7576      NAMELIST/namsbc/ nn_fsbc, ln_ana  , ln_flx, ln_blk_clio, ln_blk_core, ln_cpl    ,   & 
    76          &             nn_ice , ln_dm2dc, ln_rnf, ln_ssr     , nn_fwb     , nn_ico_cpl 
     77         &             nn_ice , ln_dm2dc, ln_rnf, ln_ssr     , nn_fwb 
    7778      !!---------------------------------------------------------------------- 
    7879 
     
    107108         WRITE(numout,*) '           Misc. options of sbc : ' 
    108109         WRITE(numout,*) '              ice management in the sbc (=0/1/2/3)       nn_ice      = ', nn_ice  
    109          WRITE(numout,*) '              ice-ocean stress computation (=0/1/2)      nn_ico_cpl  = ', nn_ico_cpl 
    110110         WRITE(numout,*) '              daily mean to diurnal cycle qsr            ln_dm2dc    = ', ln_dm2dc  
    111111         WRITE(numout,*) '              runoff / runoff mouths                     ln_rnf      = ', ln_rnf 
     
    238238       
    239239      SELECT CASE( nn_ice )                                     ! Update heat and freshwater fluxes over sea-ice areas 
    240       CASE(  1 )   ;       CALL sbc_ice_if   ( kt )                   ! Ice-cover climatology ("Ice-if" model) 
     240      CASE(  1 )   ;       CALL sbc_ice_if   ( kt )                  ! Ice-cover climatology ("Ice-if" model) 
    241241         !                                                       
    242       CASE(  2 )   ;       CALL sbc_ice_lim_2( kt, nsbc )             ! LIM 2.0 ice model 
    243                            IF( lk_bdy ) CALL bdy_ice_frs( kt ) 
     242      CASE(  2 )   ;       CALL sbc_ice_lim_2( kt, nsbc )            ! LIM-2 ice model 
     243         IF( lk_bdy )      CALL bdy_ice_frs  ( kt )                  ! BDY boundary condition 
    244244         !                                                      
    245       CASE(  3 )   ;       CALL sbc_ice_lim  ( kt, nsbc, nn_ico_cpl)  ! LIM 3.0 ice model 
     245      CASE(  3 )   ;       CALL sbc_ice_lim  ( kt, nsbc )            ! LIM-3 ice model 
    246246      END SELECT                                               
    247247 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/SBC/sbcrnf.F90

    r2350 r2370  
    44   !! Ocean forcing:  river runoff 
    55   !!===================================================================== 
    6    !! History :  OPA  !  2000-11  (R. Hordoir, E. Durand)  NetCDF FORMAT 
    7    !!   NEMO     1.0  !  2002-09  (G. Madec)  F90: Free form and module 
    8    !!            3.0  !  2006-07  (G. Madec)  Surface module  
    9    !!            3.2  !  2009-04  (B. Lemaire)  Introduce iom_put 
     6   !! History :  OPA  ! 2000-11  (R. Hordoir, E. Durand)  NetCDF FORMAT 
     7   !!   NEMO     1.0  ! 2002-09  (G. Madec)  F90: Free form and module 
     8   !!            3.0  ! 2006-07  (G. Madec)  Surface module  
     9   !!            3.2  ! 2009-04  (B. Lemaire)  Introduce iom_put 
     10   !!            3.3  ! 2010-10  (R. Furner, G. Madec) runoff distributed over ocean levels 
    1011   !!---------------------------------------------------------------------- 
    1112 
     
    4950   REAL(wp), PUBLIC, DIMENSION(jpk)     ::   rnfmsk_z     !: river mouth mask (vert.) 
    5051 
    51    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_rnf       !: structure of input river runoff (file information, fields read) 
    52  
    53    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_s_rnf     !: structure of input river runoff salinity (file information, fields read)   
    54    TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_t_rnf     !: structure of input river runoff temperature (file information, fields read)   
     52   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_rnf       !: structure: river runoff (file information, fields read) 
     53 
     54   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_s_rnf     !: structure: river runoff salinity (file information, fields read)   
     55   TYPE(FLD), ALLOCATABLE, DIMENSION(:) ::   sf_t_rnf     !: structure: river runoff temperature (file information, fields read)   
    5556  
    5657   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   h_rnf        !: depth of runoff in m 
    5758   INTEGER,  PUBLIC, DIMENSION(jpi,jpj) ::   nk_rnf       !: depth of runoff in model levels 
    5859 
    59    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpts) :: rnf_tsc_b, rnf_tsc  !: before and now  
    60    !                                                                 !: temp. & sal. content of river runoffs   [K.m/s & PSU.m/s] 
     60   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpts) :: rnf_tsc_b, rnf_tsc  !: before and now T & S contents of runoffs  [K.m/s & PSU.m/s] 
    6161 
    6262   !! * Substitutions   
     
    6565   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    6666   !! $Id$ 
    67    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     67   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    6868   !!---------------------------------------------------------------------- 
    69  
    7069CONTAINS 
    7170 
     
    171170         CALL iom_rstput( kt, nitrst, numrow, 'rnf_sc_b', rnf_tsc(:,:,jp_sal) ) 
    172171      ENDIF 
    173  
    174172      ! 
    175173   END SUBROUTINE sbc_rnf 
     174 
    176175 
    177176   SUBROUTINE sbc_rnf_div( phdivn ) 
  • branches/nemo_v3_3_beta/NEMOGCM/NEMO/OPA_SRC/step.F90

    r2348 r2370  
    2626 
    2727   !!---------------------------------------------------------------------- 
    28    !!   stp            : OPA system time-stepping 
     28   !!   stp             : OPA system time-stepping 
    2929   !!---------------------------------------------------------------------- 
    3030   USE step_oce         ! time stepping definition modules  
    3131#if defined key_top 
    32    USE trcstp          ! passive tracer time-stepping      (trc_stp routine) 
    33 #endif 
    34    USE asminc          ! assimilation increments    (tra_asm_inc, dyn_asm_inc routines) 
    35    USE stpctl          ! time stepping control            (stp_ctl routine) 
    36    USE restart         ! ocean restart                    (rst_wri routine) 
    37    USE prtctl          ! Print control                    (prt_ctl routine) 
     32   USE trcstp           ! passive tracer time-stepping      (trc_stp routine) 
     33#endif 
     34   USE asminc           ! assimilation increments    (tra_asm_inc, dyn_asm_inc routines) 
     35   USE stpctl           ! time stepping control            (stp_ctl routine) 
     36   USE restart          ! ocean restart                    (rst_wri routine) 
     37   USE prtctl           ! Print control                    (prt_ctl routine) 
    3838 
    3939#if defined key_agrif 
     
    4444   PRIVATE 
    4545 
    46    PUBLIC   stp        ! called by opa.F90 
     46   PUBLIC   stp   ! called by opa.F90 
    4747 
    4848   !! * Substitutions 
     
    5252   !! NEMO/OPA 3.3 , NEMO Consortium (2010) 
    5353   !! $Id$ 
    54    !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    55    !!---------------------------------------------------------------------- 
    56  
     54   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     55   !!---------------------------------------------------------------------- 
    5756CONTAINS 
    5857 
     
    149148         IF( ln_zps )    CALL zps_hde( kstp, jpts, tsb, gtsu, gtsv,  &    ! Partial steps: before horizontal gradient 
    150149            &                                      rhd, gru , grv  )      ! of t, s, rd at the last ocean level 
    151          IF ( .NOT. ln_traldf_grif ) THEN 
    152                          CALL ldf_slp( kstp, rhd, rn2b )        ! before slope of the lateral mixing 
     150         IF( ln_traldf_grif ) THEN                           ! before slope for Griffies operator 
     151                         CALL ldf_slp_grif( kstp ) 
    153152         ELSE 
    154                          CALL ldf_slp_grif( kstp ) 
    155                          IF ( ln_dynldf_bilap .OR. ln_dynldf_iso ) CALL ldf_slp( kstp, rhd, rn2b ) 
    156       ENDIF 
     153                         CALL ldf_slp( kstp, rhd, rn2b )     ! before slope for Madec operator 
     154         ENDIF 
    157155      ENDIF 
    158156#if defined key_traldf_c2d 
  • branches/nemo_v3_3_beta/NEMOGCM/TOOLS/COMPILE/cfg.txt

    r2281 r2370  
    44ORCA2_OFF_PISCES OPA_SRC OFF_SRC TOP_SRC 
    55POMME OPA_SRC NST_SRC 
     6ORCA2_LIM3 OPA_SRC LIM_SRC_3 
    67ORCA2_LIM OPA_SRC LIM_SRC_2 NST_SRC C1D_SRC 
Note: See TracChangeset for help on using the changeset viewer.