Changeset 888


Ignore:
Timestamp:
2008-04-11T19:05:03+02:00 (13 years ago)
Author:
ctlod
Message:

merge dev_001_SBC branche with the trunk to include the New Surface Module package, see ticket: #113

Location:
trunk/NEMO
Files:
21 added
28 deleted
91 edited

Legend:

Unmodified
Added
Removed
  • trunk/NEMO/C1D_SRC/diawri1d.F90

    r833 r888  
    1313   USE dom_oce         ! ocean space and time domain 
    1414   USE zdf_oce         ! ocean vertical physics 
     15   USE sbc_oce         ! surface boundary condition: ocean 
     16   USE sbc_ice         ! surface boundary condition: ice 
    1517   USE zdftke          ! TKE vertical mixing 
    1618   USE zdfkpp          ! KPP vertical mixing 
     
    1921   USE phycst          ! physical constants 
    2022   USE ocfzpt          ! ??? 
    21    USE ocesbc          ! surface thermohaline fluxes 
    22    USE taumod          ! surface stress 
    23    USE flxrnf          ! ??? 
    2423   USE zdfmxl          ! mixed layer 
    2524   USE daymod          ! calendar 
     
    4948   !!---------------------------------------------------------------------- 
    5049   !!   OPA 9.0 , LOCEAN-IPSL  (2005) 
    51    !! $Header$  
     50   !! $Id$ 
    5251   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    5352   !!---------------------------------------------------------------------- 
     
    194193            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    195194 
    196 #if ! defined key_dynspg_rl && defined key_lim3 
    197          ! sowaflup = sowaflep + sorunoff + sowafldp + a term associated to 
    198          !    internal damping to Levitus that can be diagnosed from others 
    199          ! sowaflcd = sowaflep + sorunoff + sowafldp + iowaflup 
    200          CALL histdef( nid_T, "iowaflup", "Ice=>ocean net freshwater"          , "kg/m2/s",   &  ! fsalt 
    201             &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    202          CALL histdef( nid_T, "sowaflep", "atmos=>ocean net freshwater"        , "kg/m2/s",   &  ! fmass 
    203             &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    204 #endif 
     195!!$#if ! defined key_dynspg_rl && ( defined key_lim2 || defined key_lim2 ) 
     196!!$         ! sowaflup = sowaflep + sorunoff + sowafldp + a term associated to 
     197!!$         !    internal damping to Levitus that can be diagnosed from others 
     198!!$         ! sowaflcd = sowaflep + sorunoff + sowafldp + iowaflup 
     199!!$         CALL histdef( nid_T, "iowaflup", "Ice=>ocean net freshwater"          , "kg/m2/s",   &  ! fsalt 
     200!!$            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     201!!$         CALL histdef( nid_T, "sowaflep", "atmos=>ocean net freshwater"        , "kg/m2/s",   &  ! fmass 
     202!!$            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     203!!$#endif 
    205204         CALL histdef( nid_T, "sowaflup", "Net Upward Water Flux"              , "Kg/m2/s",   &  ! emp 
    206205            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     
    211210         CALL histdef( nid_T, "sosalflx", "Surface Salt Flux"                  , "Kg/m2/s",   &  ! emps * sn 
    212211            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    213          CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux"             , "W/m2"   ,   &  ! qt 
     212         CALL histdef( nid_T, "sohefldo", "Net Downward Heat Flux"             , "W/m2"   ,   &  ! qsr + qns 
    214213            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
    215214         CALL histdef( nid_T, "soshfldo", "Shortwave Radiation"                , "W/m2"   ,   &  ! qsr 
     
    238237#endif 
    239238 
    240 #if ( defined key_coupled && ! defined key_lim3 )  
     239#if ( defined key_coupled && ! ( defined key_lim3 || defined key_lim2 ) )  
    241240         CALL histdef( nid_T, "sohefldp", "Surface Heat Flux: Damping"         , "W/m2"   ,   &  ! qrp 
    242241            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     
    260259#endif 
    261260 
    262 #if defined key_lim3 && defined key_coupled 
     261#if ( defined key_lim3 || defined key_lim2 ) && defined key_coupled 
    263262         CALL histdef( nid_T,"soicetem" , "Ice Surface Temperature"            , "K"      ,   &  ! tn_ice 
    264263            &          jpi, jpj, nh_T, 1  , 1, 1  , -99 , 32, clop, zsto, zout ) 
     
    275274#endif 
    276275         !                                                                                      !!! nid_U : 2D 
    277          CALL histdef( nid_T, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! taux 
     276         CALL histdef( nid_T, "sozotaux", "Wind Stress along i-axis"           , "N/m2"   ,   &  ! utau 
    278277            &          jpi, jpj, nh_T, 1  , 1, 1  , - 99, 32, clop, zsto, zout ) 
    279278 
     
    286285#endif 
    287286         !                                                                                      !!! nid_V : 2D 
    288          CALL histdef( nid_T, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! tauy 
     287         CALL histdef( nid_T, "sometauy", "Wind Stress along j-axis"           , "N/m2"   ,   &  ! vtau 
    289288            &          jpi, jpj, nh_T, 1  , 1, 1  , - 99, 32, clop, zsto, zout ) 
    290289#if defined key_zdftke 
     
    365364      CALL histwrite( nid_T, "sosstsst", it, tn(:,:,1)     , ndim_hT, ndex_hT )   ! sea surface temperature 
    366365      CALL histwrite( nid_T, "sosaline", it, sn(:,:,1)     , ndim_hT, ndex_hT )   ! sea surface salinity 
    367 #if ! defined key_dynspg_rl && defined key_lim3 
     366#if ! defined key_dynspg_rl && ( defined key_lim3 || defined key_lim2 ) 
    368367      CALL histwrite( nid_T, "iowaflup", it, fsalt(:,:)    , ndim_hT, ndex_hT )   ! ice=>ocean water flux 
    369368      CALL histwrite( nid_T, "sowaflep", it, fmass(:,:)    , ndim_hT, ndex_hT )   ! atmos=>ocean water flux 
     
    374373      zw2d(:,:) = emps(:,:) * sn(:,:,1) * tmask(:,:,1) 
    375374      CALL histwrite( nid_T, "sosalflx", it, zw2d          , ndim_hT, ndex_hT )   ! c/d salt flux 
    376       CALL histwrite( nid_T, "sohefldo", it, qt            , ndim_hT, ndex_hT )   ! total heat flux 
     375      CALL histwrite( nid_T, "sohefldo", it, qsr + qns     , ndim_hT, ndex_hT )   ! total heat flux 
    377376      CALL histwrite( nid_T, "soshfldo", it, qsr           , ndim_hT, ndex_hT )   ! solar heat flux 
    378377      CALL histwrite( nid_T, "somxl010", it, hmlp          , ndim_hT, ndex_hT )   ! mixed layer depth 
     
    397396      CALL histwrite( nid_T, "sosafldp", it, zw2d          , ndim_hT, ndex_hT )   ! salt flux damping 
    398397#endif 
    399 #if ( defined key_coupled && ! defined key_lim3 )  
     398#if ( defined key_coupled && ! ( defined key_lim3 || defined key_lim2 ) )  
    400399      CALL histwrite( nid_T, "sohefldp", it, qrp           , ndim_hT, ndex_hT )   ! heat flux damping 
    401400      CALL histwrite( nid_T, "sowafldp", it, erp           , ndim_hT, ndex_hT )   ! freshwater flux damping 
     
    412411      CALL histwrite( nid_T, "sohtc300", it, htc3          , ndim_hT, ndex_hT )   ! first 300m heaat content 
    413412#endif 
    414 #if defined key_lim3 &&  defined key_coupled  
     413#if ( defined key_lim3 || defined key_lim2 ) &&  defined key_coupled  
    415414      CALL histwrite( nid_T, "soicetem", it, tn_ice        , ndim_hT, ndex_hT )   ! surf. ice temperature 
    416415      CALL histwrite( nid_T, "soicealb", it, alb_ice       , ndim_hT, ndex_hT )   ! ice albedo 
     
    418417 
    419418      CALL histwrite( nid_T, "vozocrtx", it, un            , ndim_T , ndex_T )    ! i-current 
    420       CALL histwrite( nid_T, "sozotaux", it, taux          , ndim_hT, ndex_hT )   ! i-wind stress 
     419      CALL histwrite( nid_T, "sozotaux", it, utau          , ndim_hT, ndex_hT )   ! i-wind stress 
    421420      CALL histwrite( nid_T, "vomecrty", it, vn            , ndim_T , ndex_T  )   ! j-current 
    422       CALL histwrite( nid_T, "sometauy", it, tauy          , ndim_hT, ndex_hT )   ! j-wind stress 
     421      CALL histwrite( nid_T, "sometauy", it, vtau          , ndim_hT, ndex_hT )   ! j-wind stress 
    423422#if defined key_zdftke 
    424423      CALL histwrite( nid_T, "votlsdis", it, e_dis         , ndim_T , ndex_T )    ! Diss. Turb. lenght scale 
  • trunk/NEMO/C1D_SRC/icestp1d.F90

    r833 r888  
    66   !! History :   9.0  !  04-10  (C. Ethe)  from icestp, 1D configuration 
    77   !!---------------------------------------------------------------------- 
    8 #if defined key_cfg_1d && defined key_lim3 
     8#if defined key_cfg_1d && ( defined key_lim3 || defined key_lim2 ) 
    99   !!---------------------------------------------------------------------- 
    1010   !!   'key_cfg_1d'  .AND.                                1D Configuration 
    11    !!   'key_lim3'                                     Lim sea-ice model 
     11   !!   'key_lim2' OR 'key_lim3' :             LIM 2.0 or 3.0 sea-ice model 
    1212   !!---------------------------------------------------------------------- 
    1313   !!---------------------------------------------------------------------- 
     
    1818   USE in_out_manager  ! I/O manager 
    1919   USE ice_oce         ! ice variables 
    20    USE flx_oce         ! forcings variables 
    21    USE dom_ice         ! LIM sea-ice domain 
    22    USE cpl_oce         ! coupled ocean-atmosphere variables 
    23    USE blk_oce         ! bulk variables 
     20   USE dom_ice_2       ! LIM sea-ice domain 
     21   USE sbc_oce         ! surface boundary condition: ocean 
     22   USE sbc_ice         ! surface boundary condition: ice 
    2423   USE daymod          ! calendar 
    2524   USE phycst          ! Define parameters for the routines 
    26    USE taumod          ! surface stress forcing 
    27    USE ice             ! ice variables 
     25   USE ice_2           ! ice variables 
    2826   USE lbclnk          ! ocean lateral boundary conditions (or mpp link) 
    29    USE limthd 
    30    USE limflx 
    31    USE limwri 
    32    USE limrst 
    33  
    34    USE ocesbc          ! thermohaline fluxes 
    35    USE flxmod          ! thermohaline forcing 
    36    USE flxrnf          ! runoffs forcing 
     27   USE limthd_2 
     28   USE limwri_2 
     29   USE limrst_2 
     30 
    3731   USE tradmp          ! damping salinity trend 
    3832   USE dtatem          ! ocean temperature data 
     
    5246   !!---------------------------------------------------------------------- 
    5347   !!   LIM 2.0 , UCL-LOCEAN-IPSL (2006)  
    54    !! $Header$  
     48   !! $Id$ 
    5549   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5650   !!---------------------------------------------------------------------- 
     
    109103         u_io  (:,:) = u_io  (:,:) / FLOAT( nfice ) 
    110104         v_io  (:,:) = v_io  (:,:) / FLOAT( nfice ) 
    111          gtaux (:,:) = taux  (:,:) 
    112          gtauy (:,:) = tauy  (:,:) 
     105         gtaux (:,:) = utau  (:,:) 
     106         gtauy (:,:) = vtau  (:,:) 
    113107 
    114108         zsss_io (:,:) = SQRT( sss_io(:,:) )  
     
    220214      IF( kt == nit000 ) THEN      
    221215         qsr    (:,:) = 0.e0 
    222          qt     (:,:) = 0.e0 
     216         qns    (:,:) = 0.e0 
    223217         qrp    (:,:) = 0.e0 
    224218         emp    (:,:) = 0.e0 
     
    238232      ! ----------------- 
    239233       
    240       qt  (:,:) = fnsolar(:,:) + fsolar(:,:)     ! non solar heat flux + solar flux 
     234      qns (:,:) = fnsolar(:,:)                    ! non solar heat flux 
    241235      qsr (:,:) = fsolar(:,:)                     ! solar flux 
    242236       
     
    261255         DO ji = 1, fs_jpim1   ! vertor opt. 
    262256            ztxy        = freezn(ji,jj)             ! ice/ocean indicator at T-points 
    263             taux(ji,jj) = (1.-ztxy) * taux(ji,jj) + ztxy * ftaux(ji,jj)    ! stress at the ocean surface 
    264             tauy(ji,jj) = (1.-ztxy) * tauy(ji,jj) + ztxy * ftauy(ji,jj) 
    265          END DO 
    266       END DO 
    267        
    268       ! boundary condition on the stress (taux,tauy) 
    269       CALL lbc_lnk( taux, 'U', -1. ) 
    270       CALL lbc_lnk( tauy, 'V', -1. ) 
     257            utau(ji,jj) = (1.-ztxy) * utau(ji,jj) + ztxy * ftaux(ji,jj)    ! stress at the ocean surface 
     258            vtau(ji,jj) = (1.-ztxy) * vtau(ji,jj) + ztxy * ftauy(ji,jj) 
     259         END DO 
     260      END DO 
     261       
     262      ! boundary condition on the stress (utau,vtau) 
     263      CALL lbc_lnk( utau, 'U', -1. ) 
     264      CALL lbc_lnk( vtau, 'V', -1. ) 
    271265       
    272266      ! Re-initialization of fluxes 
  • trunk/NEMO/C1D_SRC/step1d.F90

    r719 r888  
    1515   USE dom_oce         ! ocean space and time domain variables  
    1616   USE zdf_oce         ! ocean vertical physics variables 
     17   USE sbc_oce         ! surface boundary condition: ocean 
    1718   USE ldftra_oce 
    1819   USE ldfdyn_oce 
     
    2425   USE dtatem          ! ocean temperature data           (dta_tem routine) 
    2526   USE dtasal          ! ocean salinity    data           (dta_sal routine) 
    26    USE dtasst          ! ocean sea surface temerature     (dta_sst routine) 
    27    USE taumod          ! surface stress                   (tau     routine) 
    28    USE flxmod          ! thermohaline fluxes              (flx     routine) 
    29    USE ocesbc          ! thermohaline fluxes              (oce_sbc routine) 
    30    USE flxrnf          ! runoffs                          (flx_rnf routine) 
    31    USE flxfwb          ! freshwater budget correction     (flx_fwb routine) 
    3227   USE ocfzpt          ! surface ocean freezing point    (oc_fz_pt routine) 
    3328 
     
    7570   !!---------------------------------------------------------------------- 
    7671   !!   OPA 9.0 , LOCEAN-IPSL (2005)  
    77    !! $Header$  
     72   !! $Id$ 
    7873   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    7974   !!---------------------------------------------------------------------- 
     
    157152         CALL prt_ctl(tab2d_1=emp    , clinfo1=' emp  -   : ', mask1=tmask, ovlap=1) 
    158153         CALL prt_ctl(tab2d_1=emps   , clinfo1=' emps -   : ', mask1=tmask, ovlap=1) 
    159          CALL prt_ctl(tab2d_1=qt     , clinfo1=' qt   -   : ', mask1=tmask, ovlap=1) 
     154         CALL prt_ctl(tab2d_1=qns    , clinfo1=' qns  -   : ', mask1=tmask, ovlap=1) 
    160155         CALL prt_ctl(tab2d_1=qsr    , clinfo1=' qsr  -   : ', mask1=tmask, ovlap=1) 
    161          CALL prt_ctl(tab2d_1=runoff , clinfo1=' runoff   : ', mask1=tmask, ovlap=1) 
    162156         CALL prt_ctl(tab3d_1=tmask  , clinfo1=' tmask    : ', mask1=tmask, ovlap=1, kdim=jpk) 
    163157         CALL prt_ctl(tab3d_1=tn     , clinfo1=' sst  -   : ', mask1=tmask, ovlap=1, kdim=1) 
    164158         CALL prt_ctl(tab3d_1=sn     , clinfo1=' sss  -   : ', mask1=tmask, ovlap=1, kdim=1) 
    165          CALL prt_ctl(tab2d_1=taux   , clinfo1=' tau  - x : ', mask1=umask, & 
    166             &         tab2d_2=tauy   , clinfo2='      - y : ', mask2=vmask, ovlap=1) 
     159         CALL prt_ctl(tab2d_1=utau   , clinfo1=' tau  - u : ', mask1=umask, & 
     160            &         tab2d_2=vtau   , clinfo2='      - v : ', mask2=vmask, ovlap=1) 
    167161      ENDIF 
    168162 
  • trunk/NEMO/LIM_SRC_2/dom_ice_2.F90

    r823 r888  
    11MODULE dom_ice_2 
    2 #if defined key_lim2 
    32   !!====================================================================== 
    43   !!                   ***  MODULE  dom_ice  *** 
     
    76   !! History :   2.0  !  03-08  (C. Ethe)  Free form and module 
    87   !!---------------------------------------------------------------------- 
    9  
     8#if defined key_lim2 
    109   !!---------------------------------------------------------------------- 
    1110   !!   LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    12    !! $Header$ 
     11   !! $ Id: $ 
    1312   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    1413   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/ice_2.F90

    r823 r888  
    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 
     7   !!---------------------------------------------------------------------- 
    68#if defined key_lim2 
    79   !!---------------------------------------------------------------------- 
    810   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    911   !!---------------------------------------------------------------------- 
    10    !! History : 
    11    !!   2.0  !  03-08  (C. Ethe)  F90: Free form and module 
    12    !!---------------------------------------------------------------------- 
    13    !!  LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    14    !! $Header$ 
    15    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     12   !!  LIM 2.0, UCL-LOCEAN-IPSL (2006) 
     13   !! $ Id: $ 
     14   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    1615   !!---------------------------------------------------------------------- 
    1716   !! * Modules used 
     
    2120   PRIVATE 
    2221 
    23    !! * Share Module variables 
    24    INTEGER , PUBLIC ::   & !!: ** ice-dynamic namelist (namicedyn) ** 
    25       nbiter = 1      ,  &  !: number of sub-time steps for relaxation 
    26       nbitdr = 250          !: maximum number of iterations for relaxation 
     22   !!* ice-dynamic namelist (namicedyn) * 
     23   INTEGER , PUBLIC ::   nbiter = 1         !: number of sub-time steps for relaxation 
     24   INTEGER , PUBLIC ::   nbitdr = 250       !: maximum number of iterations for relaxation 
     25   REAL(wp), PUBLIC ::   epsd   = 1.0e-20   !: tolerance parameter for dynamic 
     26   REAL(wp), PUBLIC ::   alpha  = 0.5       !: coefficient for semi-implicit coriolis 
     27   REAL(wp), PUBLIC ::   dm     = 0.6e+03   !: diffusion constant for dynamics 
     28   REAL(wp), PUBLIC ::   om     = 0.5       !: relaxation constant 
     29   REAL(wp), PUBLIC ::   resl   = 5.0e-05   !: maximum value for the residual of relaxation 
     30   REAL(wp), PUBLIC ::   cw     = 5.0e-03   !: drag coefficient for oceanic stress 
     31   REAL(wp), PUBLIC ::   angvg  = 0.e0      !: turning angle for oceanic stress 
     32   REAL(wp), PUBLIC ::   pstar  = 1.0e+04   !: first bulk-rheology parameter 
     33   REAL(wp), PUBLIC ::   c_rhg  = 20.e0     !: second bulk-rhelogy parameter 
     34   REAL(wp), PUBLIC ::   etamn  = 0.e+07    !: minimun value for viscosity 
     35   REAL(wp), PUBLIC ::   creepl = 2.e-08    !: creep limit 
     36   REAL(wp), PUBLIC ::   ecc    = 2.e0      !: eccentricity of the elliptical yield curve 
     37   REAL(wp), PUBLIC ::   ahi0   = 350.e0    !: sea-ice hor. eddy diffusivity coeff. (m2/s) 
    2738 
    28    REAL(wp), PUBLIC ::   & !!: ** ice-dynamic namelist (namicedyn) ** 
    29       epsd   = 1.0e-20,  &  !: tolerance parameter for dynamic 
    30       alpha  = 0.5    ,  &  !: coefficient for semi-implicit coriolis 
    31       dm     = 0.6e+03,  &  !: diffusion constant for dynamics 
    32       om     = 0.5    ,  &  !: relaxation constant 
    33       resl   = 5.0e-05,  &  !: maximum value for the residual of relaxation 
    34       cw     = 5.0e-03,  &  !: drag coefficient for oceanic stress 
    35       angvg  = 0.e0   ,  &  !: turning angle for oceanic stress 
    36       pstar  = 1.0e+04,  &  !: first bulk-rheology parameter 
    37       c_rhg  = 20.e0  ,  &  !: second bulk-rhelogy parameter 
    38       etamn  = 0.e+07,   &  !: minimun value for viscosity 
    39       creepl = 2.e-08,   &  !: creep limit 
    40       ecc    = 2.e0   ,  &  !: eccentricity of the elliptical yield curve 
    41       ahi0   = 350.e0       !: sea-ice hor. eddy diffusivity coeff. (m2/s) 
     39   REAL(wp), PUBLIC ::   usecc2             !:  = 1.0 / ( ecc * ecc ) 
     40   REAL(wp), PUBLIC ::   rhoco              !: = rau0 * cw 
     41   REAL(wp), PUBLIC ::   sangvg, cangvg     !: sin and cos of the turning angle for ocean stress 
     42   REAL(wp), PUBLIC ::   pstarh             !: pstar / 2.0 
    4243 
    43    REAL(wp), PUBLIC ::   &  !: 
    44       usecc2          ,  &  !:  = 1.0 / ( ecc * ecc ) 
    45       rhoco           ,  &  !: = rau0 * cw 
    46       sangvg, cangvg  ,  &  !: sin and cos of the turning angle for ocean stress 
    47       pstarh                !: pstar / 2.0 
     44   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ahiu , ahiv   !: hor. diffusivity coeff. at ocean U- and V-points (m2/s) 
     45   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   pahu , pahv   !: ice hor. eddy diffusivity coef. at ocean U- and V-points 
     46   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hsnm , hicm   !: mean snow and ice thicknesses 
     47   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ust2s                 !: friction velocity 
    4848 
    49    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::  &  !: 
    50       u_oce, v_oce,      &  !: surface ocean velocity used in ice dynamics 
    51       ahiu , ahiv ,      &  !: hor. diffusivity coeff. at ocean U- and V-points (m2/s) 
    52       pahu , pahv ,      &  !: ice hor. eddy diffusivity coef. at ocean U- and V-points 
    53       hsnm , hicm ,      &  !: mean snow and ice thicknesses 
    54       ust2s                 !: friction velocity 
     49   !!* diagnostic quantities 
     50!! REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   firic         !: IR flux over the ice (only used for outputs) 
     51!! REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fcsic         !: Sensible heat flux over the ice (only used for outputs) 
     52!! REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fleic         !: Latent heat flux over the ice (only used for outputs) 
     53!! REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qlatic        !: latent flux (only used for outputs) 
     54   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdvosif       !: Variation of volume at surface (only used for outputs) 
     55   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdvobif       !: Variation of ice volume at the bottom ice (only used for outputs) 
     56   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fdvolif       !: Total variation of ice volume (only used for outputs) 
     57   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdvonif       !: Lateral Variation of ice volume (only used for outputs) 
    5558 
    56    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::  &  !: 
    57         sst_ini,         &  !: sst read from a file for ice model initialization  
    58         sss_ini             !: sss read from a file for ice model initialization  
     59   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sist          !: Sea-Ice Surface Temperature (Kelvin ??? degree ??? I don't know) 
     60   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tfu           !: Freezing/Melting point temperature of sea water at SSS 
     61   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hicif         !: Ice thickness 
     62   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hsnif         !: Snow thickness 
     63   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hicifp        !: Ice production/melting 
     64   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   frld          !: Leads fraction = 1-a/totalarea 
     65   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   phicif        !: ice thickness  at previous time  
     66   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   pfrld         !: Leads fraction at previous time   
     67   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qstoif        !: Energy stored in the brine pockets 
     68   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fbif          !: Heat flux at the ice base 
     69   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdmsnif       !: Variation of snow mass 
     70   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdmicif       !: Variation of ice mass 
     71   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qldif         !: heat balance of the lead (or of the open ocean) 
     72   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qcmif         !: Energy needed to bring the ocean surface layer until its freezing  
     73   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fdtcn         !: net downward heat flux from the ice to the ocean 
     74   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qdtcn         !: energy from the ice to the ocean point (at a factor 2) 
     75   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   thcm          !: part of the solar energy used in the lead heat budget 
     76   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fstric        !: Solar flux transmitted trough the ice 
     77   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ffltbif       !: Array linked with the max heat contained in brine pockets (?) 
     78   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fscmbq        !: Linked with the solar flux below the ice (?) 
     79   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fsbbq         !: Also linked with the solar flux below the ice (?) 
     80   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qfvbq         !: Array used to store energy in case of toral lateral ablation (?) 
     81   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   dmgwi         !: Variation of the mass of snow ice 
    5982 
    60    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &  !: 
    61       firic  ,   &  !: IR flux over the ice (only used for outputs) 
    62       fcsic  ,   &  !: Sensible heat flux over the ice (only used for outputs) 
    63       fleic  ,   &  !: Latent heat flux over the ice (only used for outputs) 
    64       qlatic ,   &  !: latent flux 
    65       rdvosif,   &  !: Variation of volume at surface (only used for outputs) 
    66       rdvobif,   &  !: Variation of ice volume at the bottom ice (only used for outputs) 
    67       fdvolif,   &  !: Total variation of ice volume (only used for outputs) 
    68       rdvonif,   &  !: Lateral Variation of ice volume (only used for outputs) 
    69       sist   ,   &  !: Sea-Ice Surface Temperature (Kelvin ??? degree ??? I don't know) 
    70       tfu    ,   &  !: Melting point temperature of sea water 
    71       hsnif  ,   &  !: Snow thickness 
    72       hicif  ,   &  !: Ice thickness 
    73       hicifp ,   &  !: Ice production/melting 
    74       frld   ,   &  !: Leads fraction = 1-a/totalarea 
    75       phicif ,   &  !: ice thickness  at previous time  
    76       pfrld  ,   &  !: Leads fraction at previous time   
    77       qstoif ,   &  !: Energy stored in the brine pockets 
    78       fbif   ,   &  !: Heat flux at the ice base 
    79       rdmsnif,   &  !: Variation of snow mass 
    80       rdmicif,   &  !: Variation of ice mass 
    81       qldif  ,   &  !: heat balance of the lead (or of the open ocean) 
    82       qcmif  ,   &  !: Energy needed to bring the ocean surface layer until its freezing  
    83       fdtcn  ,   &  !: net downward heat flux from the ice to the ocean 
    84       qdtcn  ,   &  !: energy from the ice to the ocean 
    85       !             !  point (at a factor 2) 
    86       thcm   ,   &  !: part of the solar energy used in the lead heat budget 
    87       fstric ,   &  !: Solar flux transmitted trough the ice 
    88       ffltbif,   &  !: Array linked with the max heat contained in brine pockets (?) 
    89       fscmbq ,   &  !: Linked with the solar flux below the ice (?) 
    90       fsbbq  ,   &  !: Also linked with the solar flux below the ice (?) 
    91       qfvbq  ,   &  !: Array used to store energy in case of toral lateral ablation (?) 
    92       dmgwi         !: Variation of the mass of snow ice 
     83   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   albege        !: Albedo of the snow or ice (only for outputs) 
     84   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   albecn        !: Albedo of the ocean (only for outputs) 
     85   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tauc          !: Cloud optical depth 
    9386 
    94    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &  !: 
    95       albege ,   &  !: Albedo of the snow or ice (only for outputs) 
    96       albecn ,   &  !: Albedo of the ocean (only for outputs) 
    97       tauc   ,   &  !: Cloud optical depth 
    98       sdvt          !: u*^2/(Stress/density) 
     87   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ui_ice, vi_ice   !: two components of the ice   velocity at I-point (m/s) 
     88   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ui_oce, vi_oce   !: two components of the ocean velocity at I-point (m/s) 
    9989 
     90   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpsmax)     ::   scal0   !: ??? 
     91   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jplayersp1) ::   tbif  !: Temperature inside the ice/snow layer 
    10092 
    101    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &  !: 
    102       u_ice, v_ice,   &  !: two components of the ice velocity (m/s) 
    103       tio_u, tio_v       !: two components of the ice-ocean stress (N/m2) 
     93!! REAL(wp), DIMENSION(jpi,jpj,0:jpkmax+1) ::   reslum        !: Relative absorption of solar radiation in each ocean level 
    10494 
    105    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpsmax) ::   &  !: 
    106       scal0              !: ??? 
    107  
    108    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jplayersp1) ::   &  !: 
    109       tbif          !: Temperature inside the ice/snow layer 
    110  
    111    REAL(wp), DIMENSION(jpi,jpj,0:jpkmax+1) ::    &  !: 
    112       reslum        !: Relative absorption of solar radiation in each ocean level 
    113  
    114    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   &  !: 
    115          sxice, syice, sxxice, syyice, sxyice,      &  !: moments for advection 
    116          sxsn,  sysn,  sxxsn,  syysn,  sxysn,       &  !: 
    117          sxa,   sya,   sxxa,   syya,   sxya,        &  !: 
    118          sxc0,  syc0,  sxxc0,  syyc0,  sxyc0,       &  !: 
    119          sxc1,  syc1,  sxxc1,  syyc1,  sxyc1,       &  !: 
    120          sxc2,  syc2,  sxxc2,  syyc2,  sxyc2,       &  !: 
    121          sxst,  syst,  sxxst,  syyst,  sxyst           !: 
     95   !!* moment used in the advection scheme 
     96   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxice, syice, sxxice, syyice, sxyice   !: for ice  volume 
     97   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxsn,  sysn,  sxxsn,  syysn,  sxysn    !: for snow volume 
     98   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxa,   sya,   sxxa,   syya,   sxya     !: for ice cover area 
     99   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxc0,  syc0,  sxxc0,  syyc0,  sxyc0    !: for heat content of snow 
     100   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxc1,  syc1,  sxxc1,  syyc1,  sxyc1    !: for heat content of 1st ice layer 
     101   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxc2,  syc2,  sxxc2,  syyc2,  sxyc2    !: for heat content of 2nd ice layer 
     102   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxst,  syst,  sxxst,  syyst,  sxyst    !: for heat content of brine pockets 
    122103 
    123104#else 
  • trunk/NEMO/LIM_SRC_2/iceini_2.F90

    r823 r888  
    1717   USE dom_oce 
    1818   USE dom_ice_2 
    19    USE in_out_manager 
    2019   USE ice_oce         ! ice variables 
    21    USE flx_oce 
     20   USE sbc_oce         ! surface boundary condition: ocean 
     21   USE sbc_ice         ! surface boundary condition: ice 
    2222   USE phycst          ! Define parameters for the routines 
    2323   USE ocfzpt 
     
    2727   USE limrst_2    
    2828   USE ini1d           ! initialization of the 1D configuration 
     29   USE in_out_manager 
    2930       
    3031   IMPLICIT NONE 
     
    4041   !!---------------------------------------------------------------------- 
    4142   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    42    !! $Header$  
     43   !! $ Id: $  
    4344   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    4445   !!---------------------------------------------------------------------- 
     
    6263                  
    6364      ! Louvain la Neuve Ice model 
    64       IF( nacc == 1 ) THEN 
    65           dtsd2   = nfice * rdtmin * 0.5 
    66           rdt_ice = nfice * rdtmin 
    67       ELSE 
    68           dtsd2   = nfice * rdt * 0.5 
    69           rdt_ice = nfice * rdt 
    70       ENDIF 
     65      dtsd2   = nn_fsbc * rdttra(1) * 0.5 
     66      rdt_ice = nn_fsbc * rdttra(1) 
    7167 
    7268      CALL lim_msh_2                  ! ice mesh initialization 
  • trunk/NEMO/LIM_SRC_2/limadv_2.F90

    r823 r888  
    3333   !!---------------------------------------------------------------------- 
    3434   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    35    !! $Header$  
     35   !! $ Id: $  
    3636   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    3737   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/limdia_2.F90

    r823 r888  
    1919   USE par_ice_2       ! ice parameters 
    2020   USE ice_oce         ! ice variables 
     21   USE sbc_oce         ! surface boundary condition variables 
    2122   USE daymod          ! 
    2223   USE dom_ice_2       ! 
     
    2829   PRIVATE 
    2930 
    30    PUBLIC               lim_dia_2          ! called by ice_step 
     31   PUBLIC               lim_dia_2          ! called by sbc_ice_lim_2 
    3132   INTEGER, PUBLIC ::   ntmoy   = 1 ,   &  !: instantaneous values of ice evolution or averaging ntmoy 
    3233      &                 ninfo   = 1        !: frequency of ouputs on file ice_evolu in case of averaging 
     
    5859   !!---------------------------------------------------------------------- 
    5960   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    60    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limdia.F90,v 1.9 2007/06/29 17:03:12 opalod Exp $  
     61   !! $ Id: $ 
    6162   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    6263   !!---------------------------------------------------------------------- 
     
    8788 
    8889      nv = 1  
    89       vinfor(nv) = REAL( kt + nfice - 1 ) 
     90      vinfor(nv) = REAL( kt + nn_fsbc - 1 ) 
    9091      nv = nv + 1 
    9192      vinfor(nv) = nyear 
     
    107108               zicevol = zarea   * hicif(ji,jj) 
    108109               zsnwvol = zarea   * hsnif(ji,jj) 
    109                zicespd = zicevol * ( u_ice(ji,jj) * u_ice(ji,jj)   & 
    110                   &                + v_ice(ji,jj) * v_ice(ji,jj) ) 
     110               zicespd = zicevol * ( ui_ice(ji,jj) * ui_ice(ji,jj)   & 
     111                  &                + vi_ice(ji,jj) * vi_ice(ji,jj) ) 
    111112               vinfor(nv+ 1) = vinfor(nv+ 1) + zarea 
    112113               vinfor(nv+ 3) = vinfor(nv+ 3) + zextent15 
     
    133134                zicevol = zarea   * hicif(ji,jj) 
    134135                zsnwvol = zarea   * hsnif(ji,jj) 
    135                 zicespd = zicevol * ( u_ice(ji,jj) * u_ice(ji,jj)   & 
    136                    &                + v_ice(ji,jj) * v_ice(ji,jj) ) 
     136                zicespd = zicevol * ( ui_ice(ji,jj) * ui_ice(ji,jj)   & 
     137                   &                + vi_ice(ji,jj) * vi_ice(ji,jj) ) 
    137138                vinfor(nv+ 1) = vinfor(nv+ 1) + zarea 
    138139                vinfor(nv+ 3) = vinfor(nv+ 3) + zextent15 
     
    154155     
    155156       ! oututs on file ice_evolu     
    156        IF( MOD( kt + nfice - 1, ninfo ) == 0 ) THEN 
     157       IF( MOD( kt + nn_fsbc - 1, ninfo ) == 0 ) THEN 
    157158          WRITE(numevo_ice,fmtw) ( titvar(jv), vinfom(jv)/naveg, jv = 1, nvinfo ) 
    158159          naveg = 0 
     
    227228 
    228229       ! Definition et Ecriture de l'entete : nombre d'enregistrements  
    229        ndeb   = ( nit000 - 1 + nfice - 1 ) / ninfo 
    230        IF( nit000 - 1 + nfice == 1 ) ndeb = -1 
    231  
    232        nferme = ( nitend + nfice - 1 ) / ninfo ! nit000 - 1 + nfice - 1 + nitend - nit000 + 1 
     230       ndeb   = ( nit000 - 1 + nn_fsbc - 1 ) / ninfo 
     231       IF( nit000 - 1 + nn_fsbc == 1 ) ndeb = -1 
     232 
     233       nferme = ( nitend + nn_fsbc - 1 ) / ninfo ! nit000 - 1 + nn_fsbc - 1 + nitend - nit000 + 1 
    233234       ntot   = nferme - ndeb 
    234235       ndeb   = ninfo * ( 1 + ndeb ) 
  • trunk/NEMO/LIM_SRC_2/limdmp_2.F90

    r823 r888  
    4040   !!---------------------------------------------------------------------- 
    4141   !!   LIM 2.0 , UCL-LOCEAN-IPSL  (2006) 
    42    !! $Header$ 
     42   !! $ Id: $ 
    4343   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4444   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/limdyn_2.F90

    r823 r888  
    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 
     10   !!--------------------------------------------------------------------- 
    611#if defined key_lim2 
    712   !!---------------------------------------------------------------------- 
     
    1116   !!    lim_dyn_init_2 : initialization and namelist read 
    1217   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    14    USE phycst 
    15    USE in_out_manager  ! I/O manager 
    16    USE dom_ice_2 
    17    USE dom_oce         ! ocean space and time domain 
    18    USE ice_2 
    19    USE ice_oce 
    20    USE iceini_2 
    21    USE limistate_2 
    22    USE limrhg_2        ! ice rheology 
    23    USE lbclnk 
    24    USE lib_mpp 
    25    USE prtctl          ! Print control 
     18   USE dom_oce        ! ocean space and time domain 
     19   USE sbc_oce        ! 
     20   USE phycst         ! 
     21   USE ice_2          ! 
     22   USE ice_oce        ! 
     23   USE dom_ice_2      ! 
     24   USE iceini_2       ! 
     25   USE limistate_2    ! 
     26   USE limrhg_2       ! ice rheology 
     27 
     28   USE lbclnk         ! 
     29   USE lib_mpp        ! 
     30   USE in_out_manager ! I/O manager 
     31   USE prtctl         ! Print control 
    2632 
    2733   IMPLICIT NONE 
    2834   PRIVATE 
    2935 
    30    !! * Accessibility 
    31    PUBLIC lim_dyn_2  ! routine called by ice_step 
     36   PUBLIC   lim_dyn_2 ! routine called by sbc_ice_lim 
    3237 
    3338   !! * Module variables 
    3439   REAL(wp)  ::  rone    = 1.e0   ! constant value 
    3540 
    36    !!---------------------------------------------------------------------- 
    37    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    38    !! $Header$  
    39    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     41#  include "vectopt_loop_substitute.h90" 
     42   !!---------------------------------------------------------------------- 
     43   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     44   !! $ Id: $ 
     45   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4046   !!---------------------------------------------------------------------- 
    4147 
     
    4652      !!               ***  ROUTINE lim_dyn_2  *** 
    4753      !!                
    48       !! ** Purpose :   compute ice velocity and ocean-ice stress 
     54      !! ** Purpose :   compute ice velocity and ocean-ice friction velocity 
    4955      !!                 
    5056      !! ** Method  :  
     
    5258      !! ** Action  : - Initialisation 
    5359      !!              - Call of the dynamic routine for each hemisphere 
    54       !!              - computation of the stress at the ocean surface          
     60      !!              - computation of the friction velocity at the sea-ice base 
    5561      !!              - treatment of the case if no ice dynamic 
    56       !! History : 
    57       !!   1.0  !  01-04  (LIM)  Original code 
    58       !!   2.0  !  02-08  (C. Ethe, G. Madec)  F90, mpp 
    5962      !!--------------------------------------------------------------------- 
    6063      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    61  
    62       INTEGER ::   ji, jj             ! dummy loop indices 
    63       INTEGER ::   i_j1, i_jpj        ! Starting/ending j-indices for rheology 
    64       REAL(wp) ::   & 
    65          ztairx, ztairy,           &  ! tempory scalars 
    66          zsang , zmod,             & 
    67          ztglx , ztgly ,           & 
    68          zt11, zt12, zt21, zt22 ,  & 
    69          zustm, zsfrld, zsfrldm4,  & 
    70          zu_ice, zv_ice, ztair2 
    71       REAL(wp),DIMENSION(jpj) ::   & 
    72          zind,                     &  ! i-averaged indicator of sea-ice 
    73          zmsk                         ! i-averaged of tmask 
     64      !! 
     65      INTEGER  ::   ji, jj             ! dummy loop indices 
     66      INTEGER  ::   i_j1, i_jpj        ! Starting/ending j-indices for rheology 
     67      REAL(wp) ::   zcoef              ! temporary scalar 
     68      REAL(wp), DIMENSION(jpj)     ::   zind           ! i-averaged indicator of sea-ice 
     69      REAL(wp), DIMENSION(jpj)     ::   zmsk           ! i-averaged of tmask 
     70      REAL(wp), DIMENSION(jpi,jpj) ::   zu_io, zv_io   ! ice-ocean velocity 
    7471      !!--------------------------------------------------------------------- 
    7572 
    76       IF( kt == nit000  )   CALL lim_dyn_init_2   ! Initialization (first time-step only) 
     73      IF( kt == nit000 )   CALL lim_dyn_init_2   ! Initialization (first time-step only) 
    7774       
    78       IF ( ln_limdyn ) THEN 
    79  
     75      IF( ln_limdyn ) THEN 
     76         ! 
    8077         ! Mean ice and snow thicknesses.           
    8178         hsnm(:,:)  = ( 1.0 - frld(:,:) ) * hsnif(:,:) 
    8279         hicm(:,:)  = ( 1.0 - frld(:,:) ) * hicif(:,:) 
    83  
    84          u_oce(:,:)  = u_io(:,:) * tmu(:,:) 
    85          v_oce(:,:)  = v_io(:,:) * tmu(:,:) 
    86         
    87          !                                         ! Rheology (ice dynamics) 
    88          !                                         ! ======== 
     80         ! 
     81         !                                     ! Rheology (ice dynamics) 
     82         !                                     ! ======== 
    8983          
    9084         !  Define the j-limits where ice rheology is computed 
     
    9488            i_j1 = 1    
    9589            i_jpj = jpj 
    96             IF(ln_ctl)    THEN 
    97                CALL prt_ctl_info('lim_dyn  :    i_j1 = ', ivar1=i_j1, clinfo2=' ij_jpj = ', ivar2=i_jpj) 
    98             ENDIF 
     90            IF(ln_ctl)   CALL prt_ctl_info( 'lim_dyn  :    i_j1 = ', ivar1=i_j1, clinfo2=' ij_jpj = ', ivar2=i_jpj ) 
    9991            CALL lim_rhg_2( i_j1, i_jpj ) 
    100  
     92            ! 
    10193         ELSE                                 ! optimization of the computational area 
    102  
     94            ! 
    10395            DO jj = 1, jpj 
    10496               zind(jj) = SUM( frld (:,jj  ) )   ! = FLOAT(jpj) if ocean everywhere on a j-line 
    10597               zmsk(jj) = SUM( tmask(:,jj,1) )   ! = 0          if land  everywhere on a j-line 
    106    !!i         write(numout,*) narea, 'limdyn' , jj, zind(jj), zmsk(jj) 
    107             END DO 
    108  
     98            END DO 
     99            ! 
    109100            IF( l_jeq ) THEN                     ! local domain include both hemisphere 
    110101               !                                 ! Rheology is computed in each hemisphere 
     
    118109               i_j1 = MAX( 1, i_j1-1 ) 
    119110               IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : NH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 
    120      
     111               !  
    121112               CALL lim_rhg_2( i_j1, i_jpj ) 
    122      
     113               !  
    123114               ! Southern hemisphere 
    124115               i_j1  =  1  
     
    129120               i_jpj = MIN( jpj, i_jpj+2 ) 
    130121               IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : SH i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 
    131      
     122               !  
    132123               CALL lim_rhg_2( i_j1, i_jpj ) 
    133      
     124               !  
    134125            ELSE                                 ! local domain extends over one hemisphere only 
    135126               !                                 ! Rheology is computed only over the ice cover 
     
    148139     
    149140               IF(ln_ctl)   WRITE(numout,*) 'lim_dyn : one hemisphere: i_j1 = ', i_j1, ' ij_jpj = ', i_jpj 
    150      
     141               !  
    151142               CALL lim_rhg_2( i_j1, i_jpj ) 
    152  
     143               ! 
    153144            ENDIF 
    154  
     145            ! 
    155146         ENDIF 
    156147 
    157          IF(ln_ctl)   THEN  
    158             CALL prt_ctl(tab2d_1=u_oce , clinfo1=' lim_dyn  : u_oce :', tab2d_2=v_oce , clinfo2=' v_oce :') 
    159             CALL prt_ctl(tab2d_1=u_ice , clinfo1=' lim_dyn  : u_ice :', tab2d_2=v_ice , clinfo2=' v_ice :') 
    160          ENDIF 
    161           
    162          !                                         ! Ice-Ocean stress 
    163          !                                         ! ================ 
     148         IF(ln_ctl)   CALL prt_ctl(tab2d_1=ui_ice , clinfo1=' lim_dyn  : ui_ice :', tab2d_2=vi_ice , clinfo2=' vi_ice :') 
     149          
     150         ! computation of friction velocity 
     151         ! -------------------------------- 
     152         ! ice-ocean velocity at U & V-points (ui_ice vi_ice at I-point ; ssu_m, ssv_m at U- & V-points) 
     153          
     154         DO jj = 1, jpjm1 
     155            DO ji = 1, fs_jpim1   ! vector opt. 
     156               zu_io(ji,jj) = 0.5 * ( ui_ice(ji+1,jj+1) + ui_ice(ji+1,jj  ) ) - ssu_m(ji,jj) 
     157               zv_io(ji,jj) = 0.5 * ( vi_ice(ji+1,jj+1) + vi_ice(ji  ,jj+1) ) - ssv_m(ji,jj) 
     158            END DO 
     159         END DO 
     160         ! frictional velocity at T-point 
    164161         DO jj = 2, jpjm1 
    165             zsang  = SIGN(1.e0, gphif(1,jj-1) ) * sangvg 
    166             DO ji = 2, jpim1 
    167                ! computation of wind stress over ocean in X and Y direction 
    168 #if defined key_coupled && defined key_lim_cp1 
    169                ztairx =  frld(ji-1,jj  ) * gtaux(ji-1,jj  ) + frld(ji,jj  ) * gtaux(ji,jj  )      & 
    170                   &    + frld(ji-1,jj-1) * gtaux(ji-1,jj-1) + frld(ji,jj-1) * gtaux(ji,jj-1) 
    171  
    172                ztairy =  frld(ji-1,jj  ) * gtauy(ji-1,jj  ) + frld(ji,jj  ) * gtauy(ji,jj  )      & 
    173                   &    + frld(ji-1,jj-1) * gtauy(ji-1,jj-1) + frld(ji,jj-1) * gtauy(ji,jj-1) 
    174 #else 
    175                zsfrld  = frld(ji,jj) + frld(ji-1,jj) + frld(ji-1,jj-1) + frld(ji,jj-1) 
    176                ztairx  = zsfrld * gtaux(ji,jj) 
    177                ztairy  = zsfrld * gtauy(ji,jj) 
    178 #endif 
    179                zsfrldm4 = 4 - frld(ji,jj) - frld(ji-1,jj) - frld(ji-1,jj-1) - frld(ji,jj-1) 
    180                zu_ice   = u_ice(ji,jj) - u_oce(ji,jj) 
    181                zv_ice   = v_ice(ji,jj) - v_oce(ji,jj) 
    182                zmod     = SQRT( zu_ice * zu_ice + zv_ice * zv_ice )  
    183                ztglx   = zsfrldm4 * rhoco * zmod * ( cangvg * zu_ice - zsang * zv_ice )  
    184                ztgly   = zsfrldm4 * rhoco * zmod * ( cangvg * zv_ice + zsang * zu_ice )  
    185  
    186                tio_u(ji,jj) = - ( ztairx + 1.0 * ztglx ) / ( 4 * rau0 ) 
    187                tio_v(ji,jj) = - ( ztairy + 1.0 * ztgly ) / ( 4 * rau0 ) 
     162            DO ji = fs_2, fs_jpim1   ! vector opt. 
     163               ust2s(ji,jj) = 0.5 * cw                                                          & 
     164                  &         * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
     165                  &            + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
    188166            END DO 
    189167         END DO 
    190           
    191          ! computation of friction velocity 
     168         ! 
     169      ELSE      ! no ice dynamics : transmit directly the atmospheric stress to the ocean 
     170         ! 
     171         zcoef = SQRT( 0.5 ) / rau0 
    192172         DO jj = 2, jpjm1 
    193             DO ji = 2, jpim1 
    194  
    195                zu_ice   = u_ice(ji-1,jj-1) - u_oce(ji-1,jj-1) 
    196                zv_ice   = v_ice(ji-1,jj-1) - v_oce(ji-1,jj-1) 
    197                zt11  = rhoco * ( zu_ice * zu_ice + zv_ice * zv_ice ) 
    198  
    199                zu_ice   = u_ice(ji-1,jj) - u_oce(ji-1,jj) 
    200                zv_ice   = v_ice(ji-1,jj) - v_oce(ji-1,jj) 
    201                zt12  = rhoco * ( zu_ice * zu_ice + zv_ice * zv_ice )  
    202  
    203                zu_ice   = u_ice(ji,jj-1) - u_oce(ji,jj-1) 
    204                zv_ice   = v_ice(ji,jj-1) - v_oce(ji,jj-1) 
    205                zt21  = rhoco * ( zu_ice * zu_ice + zv_ice * zv_ice )  
    206  
    207                zu_ice   = u_ice(ji,jj) - u_oce(ji,jj) 
    208                zv_ice   = v_ice(ji,jj) - v_oce(ji,jj) 
    209                zt22  = rhoco * ( zu_ice * zu_ice + zv_ice * zv_ice )  
    210  
    211                ztair2 = gtaux(ji,jj) * gtaux(ji,jj) + gtauy(ji,jj) * gtauy(ji,jj) 
    212  
    213                zustm =  ( 1 - frld(ji,jj) ) * 0.25 * ( zt11 + zt12 + zt21 + zt22 )        & 
    214                   &  +        frld(ji,jj)   * SQRT( ztair2 ) 
    215  
    216                ust2s(ji,jj) = ( zustm / rau0 ) * ( rone + sdvt(ji,jj) ) * tms(ji,jj) 
     173            DO ji = fs_2, fs_jpim1   ! vector opt. 
     174               ust2s(ji,jj) = zcoef * tms(ji,jj) * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
     175                  &                                     + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) 
    217176            END DO 
    218177         END DO 
    219  
    220        ELSE      ! no ice dynamics : transmit directly the atmospheric stress to the ocean 
    221                      
    222           DO jj = 2, jpjm1 
    223              DO ji = 2, jpim1 
    224 #if defined key_coupled && defined key_lim_cp1 
    225                 tio_u(ji,jj) = - (  gtaux(ji  ,jj  ) + gtaux(ji-1,jj  )       & 
    226                    &              + gtaux(ji-1,jj-1) + gtaux(ji  ,jj-1) ) / ( 4 * rau0 ) 
    227  
    228                 tio_v(ji,jj) = - (  gtauy(ji  ,jj )  + gtauy(ji-1,jj  )       & 
    229                    &              + gtauy(ji-1,jj-1) + gtauy(ji  ,jj-1) ) / ( 4 * rau0 ) 
    230 #else 
    231                 tio_u(ji,jj) = - gtaux(ji,jj) / rau0 
    232                 tio_v(ji,jj) = - gtauy(ji,jj) / rau0  
    233 #endif 
    234                 ztair2       = gtaux(ji,jj) * gtaux(ji,jj) + gtauy(ji,jj) * gtauy(ji,jj) 
    235                 zustm        = SQRT( ztair2  ) 
    236  
    237                 ust2s(ji,jj) = ( zustm / rau0 ) * ( rone + sdvt(ji,jj) ) * tms(ji,jj) 
    238             END DO 
    239          END DO 
    240  
     178         ! 
    241179      ENDIF 
    242  
     180      ! 
    243181      CALL lbc_lnk( ust2s, 'T',  1. )   ! T-point 
    244       CALL lbc_lnk( tio_u, 'I', -1. )   ! I-point (i.e. ice U-V point) 
    245       CALL lbc_lnk( tio_v, 'I', -1. )   ! I-point (i.e. ice U-V point) 
    246  
    247       IF(ln_ctl) THEN  
    248             CALL prt_ctl(tab2d_1=tio_u , clinfo1=' lim_dyn  : tio_u :', tab2d_2=tio_v , clinfo2=' tio_v :') 
    249             CALL prt_ctl(tab2d_1=ust2s , clinfo1=' lim_dyn  : ust2s :') 
    250       ENDIF 
     182      ! 
     183      IF(ln_ctl)   CALL prt_ctl(tab2d_1=ust2s , clinfo1=' lim_dyn  : ust2s :') 
    251184 
    252185   END SUBROUTINE lim_dyn_2 
     
    257190      !!                  ***  ROUTINE lim_dyn_init_2  *** 
    258191      !! 
    259       !! ** Purpose : Physical constants and parameters linked to the ice 
    260       !!      dynamics 
    261       !! 
    262       !! ** Method  :  Read the namicedyn namelist and check the ice-dynamic 
    263       !!       parameter values called at the first timestep (nit000) 
     192      !! ** Purpose :   Physical constants and parameters linked to the ice 
     193      !!              dynamics 
     194      !! 
     195      !! ** Method  :   Read the namicedyn namelist and check the ice-dynamic 
     196      !!              parameter values 
    264197      !! 
    265198      !! ** input   :   Namelist namicedyn 
    266       !! 
    267       !! history : 
    268       !!  8.5  ! 03-08 (C. Ethe) original code 
    269199      !!------------------------------------------------------------------- 
    270200      NAMELIST/namicedyn/ epsd, alpha,     & 
     
    273203      !!------------------------------------------------------------------- 
    274204 
    275       ! Define the initial parameters 
    276       ! ------------------------- 
    277  
    278       ! Read Namelist namicedyn 
    279       REWIND ( numnam_ice ) 
     205      REWIND ( numnam_ice )                       ! Read Namelist namicedyn 
    280206      READ   ( numnam_ice  , namicedyn ) 
    281       IF(lwp) THEN 
     207 
     208      IF(lwp) THEN                                ! Control print 
    282209         WRITE(numout,*) 
    283210         WRITE(numout,*) 'lim_dyn_init_2: ice parameters for ice dynamics ' 
     
    291218         WRITE(numout,*) '       maximum value for the residual of relaxation     resl   = ', resl 
    292219         WRITE(numout,*) '       drag coefficient for oceanic stress              cw     = ', cw 
    293          WRITE(numout,*) '       turning angle for oceanic stress                 angvg  = ', angvg 
     220         WRITE(numout,*) '       turning angle for oceanic stress                 angvg  = ', angvg, ' degrees' 
    294221         WRITE(numout,*) '       first bulk-rheology parameter                    pstar  = ', pstar 
    295222         WRITE(numout,*) '       second bulk-rhelogy parameter                    c_rhg  = ', c_rhg 
     
    303230      usecc2 = 1.0 / ( ecc * ecc ) 
    304231      rhoco  = rau0 * cw 
    305       angvg  = angvg * rad 
     232      angvg  = angvg * rad      ! convert angvg from degree to radian 
    306233      sangvg = SIN( angvg ) 
    307234      cangvg = COS( angvg ) 
    308235      pstarh = pstar / 2.0 
    309       sdvt(:,:) = 0.e0 
    310  
    311       !  Diffusion coefficients. 
    312       ahiu(:,:) = ahi0 * umask(:,:,1) 
     236      ! 
     237      ahiu(:,:) = ahi0 * umask(:,:,1)            ! Ice eddy Diffusivity coefficients. 
    313238      ahiv(:,:) = ahi0 * vmask(:,:,1) 
    314  
     239      ! 
    315240   END SUBROUTINE lim_dyn_init_2 
    316241 
  • trunk/NEMO/LIM_SRC_2/limhdf_2.F90

    r823 r888  
    66#if defined key_lim2 
    77   !!---------------------------------------------------------------------- 
    8    !!   'key_lim2'  i                                 LIM 2.0 sea-ice model 
     8   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    99   !!---------------------------------------------------------------------- 
    1010   !!   lim_hdf_2  : diffusion trend on sea-ice variable 
     
    3434   !!---------------------------------------------------------------------- 
    3535   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    36    !! $Header$  
     36   !! $ Id: $ 
    3737   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    3838   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/limistate_2.F90

    r823 r888  
    44   !!              Initialisation of diagnostics ice variables 
    55   !!====================================================================== 
    6    !! History :   2.0  !  01-04  (C. Ethe, G. Madec)  Original code 
     6   !! History :   1.0  !  01-04  (C. Ethe, G. Madec)  Original code 
     7   !!             2.0  !  03-08  (G. Madec)  add lim_istate_init 
    78   !!                  !  04-04  (S. Theetten) initialization from a file 
    89   !!                  !  06-07  (S. Masson)  IOM to read the restart 
     10   !!                  !  07-10  (G. Madec)  surface module 
    911   !!-------------------------------------------------------------------- 
    1012#if defined key_lim2 
     
    1820   USE phycst 
    1921   USE ocfzpt 
    20    USE oce             ! dynamics and tracers variables      !!gm used??? 
    21    USE dom_oce                                                     !!gm used??? 
    2222   USE par_ice_2       ! ice parameters 
    2323   USE ice_oce         ! ice variables 
    2424   USE dom_ice_2 
    2525   USE lbclnk 
     26   USE oce 
    2627   USE ice_2 
    2728   USE iom 
     
    4748   !!---------------------------------------------------------------------- 
    4849   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
    49    !! $Header$  
     50   !! $ Id: $ 
    5051   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5152   !!---------------------------------------------------------------------- 
     
    6667      REAL(wp), DIMENSION(jpi,jpj) ::   ztn   ! workspace 
    6768      !-------------------------------------------------------------------- 
    68  
    69        CALL lim_istate_init_2   !  reading the initials parameters of the ice 
    70  
    71       !-- Initialisation of sst,sss,u,v do i=1,jpi 
    72       u_io(:,:)  = 0.e0       ! ice velocity in x direction 
    73       v_io(:,:)  = 0.e0       ! ice velocity in y direction 
    74  
    75       IF( ln_limini ) THEN    !  
    76          
    77          ! Initialisation at tn if no ice or sst_ini if ice 
    78          ! Idem for salinity 
    79  
    80       !--- Criterion for presence (zidto=1.) or absence (zidto=0.) of ice 
    81          DO jj = 1 , jpj 
    82             DO ji = 1 , jpi 
    83                 
    84                zidto = MAX(zzero, - SIGN(1.,frld(ji,jj) - 1.)) 
    85                 
    86                sst_io(ji,jj) = ( nfice - 1 ) * (zidto * sst_ini(ji,jj)  + &   ! use the ocean initial values 
    87                     &          (1.0 - zidto ) * ( tn(ji,jj,1) + rt0 ))        ! tricky trick *(nfice-1) ! 
    88                sss_io(ji,jj) = ( nfice - 1 ) * (zidto * sss_ini(ji,jj) + & 
    89                     &          (1.0 - zidto ) *  sn(ji,jj,1) ) 
    90  
    91                ! to avoid the the melting of ice, several layers (mixed layer) should be 
    92                ! set to sst_ini (sss_ini) if there is ice 
    93                ! example for one layer  
    94                ! tn(ji,jj,1) = zidto * ( sst_ini(ji,jj) - rt0 )  + (1.0 - zidto ) *  tn(ji,jj,1) 
    95                ! sn(ji,jj,1) = zidto * sss_ini(ji,jj)  + (1.0 - zidto ) *  sn(ji,jj,1) 
    96                ! tb(ji,jj,1) = tn(ji,jj,1) 
    97                ! sb(ji,jj,1) = sn(ji,jj,1) 
    98             END DO 
    99          END DO 
    100           
    101           
    102          !  tfu: Melting point of sea water 
    103          tfu(:,:)  = ztf    
    104           
    105          tfu(:,:)  = ABS ( rt0 - 0.0575       * sss_ini(:,:)                               & 
    106               &                    + 1.710523e-03 * sss_ini(:,:) * SQRT( sss_ini(:,:) )    & 
    107               &                    - 2.154996e-04 * sss_ini(:,:) * sss_ini(:,:) ) 
    108       ELSE                     ! 
    109  
     69  
     70      CALL lim_istate_init_2     !  reading the initials parameters of the ice 
     71 
     72      IF( .NOT. ln_limini ) THEN   
    11073          
    11174         ! Initialisation at tn or -2 if ice 
     
    11679            END DO 
    11780         END DO 
    118           
    119          u_io  (:,:) = 0.e0 
    120          v_io  (:,:) = 0.e0 
    121          sst_io(:,:) = ( nfice - 1 ) * ( tn(:,:,1) + rt0 )   ! use the ocean initial values 
    122          sss_io(:,:) = ( nfice - 1 ) *   sn(:,:,1)           ! tricky trick *(nfice-1) ! 
    123           
    124          ! reference salinity 34psu 
     81                   
     82         !  tfu: Melting point of sea water [Kelvin] 
    12583         zs0 = 34.e0 
    126          ztf = ABS ( rt0 - 0.0575       * zs0                           & 
    127               &                    + 1.710523e-03 * zs0 * SQRT( zs0 )   & 
    128               &                    - 2.154996e-04 * zs0 *zs0          ) 
    129           
    130          !  tfu: Melting point of sea water 
    131          tfu(:,:)  = ztf    
     84         ztf = rt0 + ( - 0.0575 + 1.710523e-3 * SQRT( zs0 ) - 2.154996e-4 * zs0 ) * zs0 
     85         tfu(:,:) = ztf 
    13286          
    13387         DO jj = 1, jpj 
     
    152106         tbif  (:,:,2) = tfu(:,:) 
    153107         tbif  (:,:,3) = tfu(:,:) 
    154        
     108 
    155109      ENDIF 
     110      
    156111      fsbbq (:,:)   = 0.e0 
    157112      qstoif(:,:)   = 0.e0 
    158       u_ice (:,:)   = 0.e0 
    159       v_ice (:,:)   = 0.e0 
     113      ui_ice(:,:)   = 0.e0 
     114      vi_ice(:,:)   = 0.e0 
    160115# if defined key_coupled 
    161116      albege(:,:)   = 0.8 * tms(:,:) 
     
    191146 
    192147      CALL lbc_lnk( hsnif, 'T', 1. ) 
    193       CALL lbc_lnk( sist , 'T', 1. ) 
     148      CALL lbc_lnk( sist , 'T', 1. , pval = rt0 )      ! set rt0 on closed boundary (required by bulk formulation) 
    194149      DO jk = 1, jplayersp1 
    195150         CALL lbc_lnk(tbif(:,:,jk), 'T', 1. ) 
     
    197152      CALL lbc_lnk( fsbbq  , 'T', 1. ) 
    198153      CALL lbc_lnk( qstoif , 'T', 1. ) 
    199       CALL lbc_lnk( sss_io , 'T', 1. ) 
    200       ! 
     154 
    201155   END SUBROUTINE lim_istate_2 
    202156 
     
    209163      !! 
    210164      !! ** Method  :   Read the namiceini namelist and check the parameter  
    211       !!                values called at the first timestep (nit000) 
    212       !!                or 
    213       !!                Read 7 variables from a previous restart file 
    214       !!                sst, sst, hicif, hsnif, frld, ts & tbif 
     165      !!       values called at the first timestep (nit000) 
    215166      !! 
    216167      !! ** input   :   Namelist namiceini 
     
    222173         &                hnins, hgins, alins 
    223174      !!------------------------------------------------------------------- 
    224        
    225       ! Read Namelist namiceini  
    226       REWIND ( numnam_ice ) 
     175      ! 
     176      REWIND ( numnam_ice )               ! Read Namelist namiceini  
    227177      READ   ( numnam_ice , namiceini ) 
    228        
    229       IF(.NOT. ln_limini) THEN  
    230          IF(lwp) THEN 
    231             WRITE(numout,*) 
    232             WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation ' 
    233             WRITE(numout,*) '~~~~~~~~~~~~~~~' 
    234             WRITE(numout,*) '         threshold water temp. for initial sea-ice    ttest      = ', ttest 
    235             WRITE(numout,*) '         initial snow thickness in the north          hninn      = ', hninn 
    236             WRITE(numout,*) '         initial ice thickness in the north           hginn      = ', hginn  
    237             WRITE(numout,*) '         initial leads area in the north              alinn      = ', alinn             
    238             WRITE(numout,*) '         initial snow thickness in the south          hnins      = ', hnins  
    239             WRITE(numout,*) '         initial ice thickness in the south           hgins      = ', hgins 
    240             WRITE(numout,*) '         initial leads area in the south              alins      = ', alins 
    241          ENDIF 
     178      ! 
     179      IF(lwp) THEN 
     180         WRITE(numout,*) 
     181         WRITE(numout,*) 'lim_istate_init_2 : ice parameters inititialisation ' 
     182         WRITE(numout,*) '~~~~~~~~~~~~~~~~~' 
     183         WRITE(numout,*) '         threshold water temp. for initial sea-ice    ttest      = ', ttest 
     184         WRITE(numout,*) '         initial snow thickness in the north          hninn      = ', hninn 
     185         WRITE(numout,*) '         initial ice thickness in the north           hginn      = ', hginn  
     186         WRITE(numout,*) '         initial leads area in the north              alinn      = ', alinn             
     187         WRITE(numout,*) '         initial snow thickness in the south          hnins      = ', hnins  
     188         WRITE(numout,*) '         initial ice thickness in the south           hgins      = ', hgins 
     189         WRITE(numout,*) '         initial leads area in the south              alins      = ', alins 
     190         WRITE(numout,*) '         Ice state initialization using input file    ln_limini  = ', ln_limini 
    242191      ENDIF 
    243192 
    244193      IF( ln_limini ) THEN                      ! Ice initialization using input file 
    245  
     194         ! 
    246195         CALL iom_open( 'Ice_initialization.nc', inum_ice ) 
    247  
     196         ! 
    248197         IF( inum_ice > 0 ) THEN 
    249             IF(lwp) THEN 
    250                WRITE(numout,*) ' ' 
    251                WRITE(numout,*) 'lim_istate_init : ice state initialization with : Ice_initialization.nc' 
    252                WRITE(numout,*) '~~~~~~~~~~~~~~~' 
    253                WRITE(numout,*) '         Ice state initialization using input file    ln_limini  = ', ln_limini 
    254                WRITE(numout,*) ' ' 
    255             ENDIF 
     198            IF(lwp) WRITE(numout,*) 
     199            IF(lwp) WRITE(numout,*) '                  ice state initialization with : Ice_initialization.nc' 
    256200             
    257             CALL iom_get( inum_ice, jpdom_data, 'sst'  , sst_ini(:,:) )         
    258             CALL iom_get( inum_ice, jpdom_data, 'sss'  , sss_ini(:,:) )        
    259             CALL iom_get( inum_ice, jpdom_data, 'hicif', hicif  (:,:) )       
    260             CALL iom_get( inum_ice, jpdom_data, 'hsnif', hsnif  (:,:) )       
    261             CALL iom_get( inum_ice, jpdom_data, 'frld' , frld   (:,:) )      
    262             CALL iom_get( inum_ice, jpdom_data, 'ts'   , sist   (:,:) ) 
     201            CALL iom_get( inum_ice, jpdom_data, 'hicif', hicif )       
     202            CALL iom_get( inum_ice, jpdom_data, 'hsnif', hsnif )       
     203            CALL iom_get( inum_ice, jpdom_data, 'frld' , frld  )      
     204            CALL iom_get( inum_ice, jpdom_data, 'ts'   , sist  ) 
    263205            CALL iom_get( inum_ice, jpdom_unknown, 'tbif', tbif(1:nlci,1:nlcj,:),   & 
    264206                 &        kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,jplayersp1 /) ) 
     
    268210 
    269211            CALL iom_close( inum_ice) 
    270              
     212            ! 
    271213         ENDIF 
    272214      ENDIF 
    273       ! 
     215      !      
    274216   END SUBROUTINE lim_istate_init_2 
    275217 
  • trunk/NEMO/LIM_SRC_2/limmsh_2.F90

    r823 r888  
    2525   !!---------------------------------------------------------------------- 
    2626   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    27    !! $Header$  
     27   !! $ Id: $ 
    2828   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    2929   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/limrhg_2.F90

    r823 r888  
    44   !!   Ice rheology :  performs sea ice rheology 
    55   !!====================================================================== 
     6   !! History :  0.0  !  93-12  (M.A. Morales Maqueda.)  Original code 
     7   !!            1.0  !  94-12  (H. Goosse)  
     8   !!            2.0  !  03-12  (C. Ethe, G. Madec)  F90, mpp 
     9   !!            " "  !  06-08  (G. Madec)  surface module, ice-stress at I-point 
     10   !!            " "  !  09-09  (G. Madec)  Huge verctor optimisation 
     11   !!---------------------------------------------------------------------- 
    612#if defined key_lim2 
    713   !!---------------------------------------------------------------------- 
    814   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    915   !!---------------------------------------------------------------------- 
     16   !!---------------------------------------------------------------------- 
    1017   !!   lim_rhg_2   : computes ice velocities 
    1118   !!---------------------------------------------------------------------- 
    12    !! * Modules used 
    13    USE phycst 
    14    USE par_oce 
    15    USE ice_oce         ! ice variables 
    16    USE dom_ice_2 
    17    USE ice_2 
    18    USE lbclnk 
    19    USE lib_mpp 
    20    USE in_out_manager  ! I/O manager 
    21    USE prtctl          ! Print control 
     19   USE par_oce        ! ocean parameter 
     20   USE ice_oce        ! ice variables 
     21   USE sbc_ice        ! surface boundary condition: ice variables 
     22   USE dom_ice_2      ! domaine: ice variables 
     23   USE phycst         ! physical constant 
     24   USE ice_2          ! ice variables 
     25   USE lbclnk         ! lateral boundary condition 
     26   USE lib_mpp        ! MPP library 
     27   USE in_out_manager ! I/O manager 
     28   USE prtctl         ! Print control 
    2229 
    2330   IMPLICIT NONE 
    2431   PRIVATE 
    2532 
    26    !! * Routine accessibility 
    27    PUBLIC lim_rhg_2  ! routine called by lim_dyn_2 
    28  
    29    !! * Module variables 
    30    REAL(wp)  ::           &  ! constant values 
    31       rzero   = 0.e0   ,  & 
    32       rone    = 1.e0 
    33    !!---------------------------------------------------------------------- 
    34    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    35    !! $Header$  
    36    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     33   PUBLIC   lim_rhg_2 ! routine called by lim_dyn 
     34 
     35   REAL(wp) ::   rzero   = 0.e0   ! constant value: zero 
     36   REAL(wp) ::   rone    = 1.e0   !            and  one 
     37 
     38   !! * Substitutions 
     39#  include "vectopt_loop_substitute.h90" 
     40   !!---------------------------------------------------------------------- 
     41   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     42   !! $ Id: $ 
     43   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3744   !!---------------------------------------------------------------------- 
    3845 
     
    4855      !!  viscous-plastic law including shear strength and a bulk rheology. 
    4956      !! 
    50       !! ** Action  : - compute u_ice, v_ice the sea-ice velocity 
     57      !! ** Action  : - compute ui_ice, vi_ice the sea-ice velocity defined 
     58      !!              at I-point 
     59      !!------------------------------------------------------------------- 
     60      INTEGER, INTENT(in) ::   k_j1    ! southern j-index for ice computation 
     61      INTEGER, INTENT(in) ::   k_jpj   ! northern j-index for ice computation 
    5162      !! 
    52       !! History : 
    53       !!   0.0  !  93-12  (M.A. Morales Maqueda.)  Original code 
    54       !!   1.0  !  94-12  (H. Goosse)  
    55       !!   2.0  !  03-12  (C. Ethe, G. Madec)  F90, mpp 
     63      INTEGER ::   ji, jj              ! dummy loop indices 
     64      INTEGER ::   iter, jter          ! temporary integers 
     65      CHARACTER (len=50) ::   charout 
     66      REAL(wp) ::   ze11  , ze12  , ze22  , ze21       ! temporary scalars 
     67      REAL(wp) ::   zt11  , zt12  , zt21  , zt22       !    "         " 
     68      REAL(wp) ::   zvis11, zvis21, zvis12, zvis22     !    "         " 
     69      REAL(wp) ::   zgphsx, ztagnx, zunw, zur, zusw    !    "         " 
     70      REAL(wp) ::   zgphsy, ztagny, zvnw, zvr          !    "         " 
     71      REAL(wp) ::   zresm,  za, zac, zmod 
     72      REAL(wp) ::   zmpzas, zstms, zindu, zusdtp, zmassdt, zcorlal 
     73      REAL(wp) ::   ztrace2, zdeter, zdelta, zmask, zdgp, zdgi, zdiag 
     74      REAL(wp) ::   za1, zb1, zc1, zd1 
     75      REAL(wp) ::   za2, zb2, zc2, zd2, zden 
     76      REAL(wp) ::   zs11_11, zs11_12, zs11_21, zs11_22 
     77      REAL(wp) ::   zs12_11, zs12_12, zs12_21, zs12_22 
     78      REAL(wp) ::   zs21_11, zs21_12, zs21_21, zs21_22 
     79      REAL(wp) ::   zs22_11, zs22_12, zs22_21, zs22_22 
     80      REAL(wp), DIMENSION(jpi,  jpj  ) ::   zfrld, zmass, zcorl 
     81      REAL(wp), DIMENSION(jpi,  jpj  ) ::   za1ct, za2ct, zresr 
     82      REAL(wp), DIMENSION(jpi,  jpj  ) ::   zc1u, zc1v, zc2u, zc2v 
     83      REAL(wp), DIMENSION(jpi,  jpj  ) ::   zsang 
     84      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zu0, zv0 
     85      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zu_n, zv_n 
     86      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zu_a, zv_a 
     87      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zviszeta, zviseta 
     88      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zzfrld, zztms 
     89      REAL(wp), DIMENSION(jpi,0:jpj+1) ::   zi1, zi2, zmasst, zpresh 
     90 
    5691      !!------------------------------------------------------------------- 
    57       ! * Arguments 
    58       INTEGER, INTENT(in) ::   k_j1 ,  &  ! southern j-index for ice computation 
    59          &                     k_jpj      ! northern j-index for ice computation 
    60  
    61       ! * Local variables 
    62       INTEGER ::   ji, jj              ! dummy loop indices 
    63  
    64       INTEGER  :: & 
    65          iim1, ijm1, iip1 , ijp1   , & ! temporary integers 
    66          iter, jter                    !    "          " 
    67  
    68       CHARACTER (len=50) :: charout 
    69  
    70       REAL(wp) :: & 
    71          ze11  , ze12  , ze22  , ze21  ,   &  ! temporary scalars 
    72          zt11  , zt12  , zt21  , zt22  ,   &  !    "         " 
    73          zvis11, zvis21, zvis12, zvis22,   &  !    "         " 
    74          zgphsx, ztagnx, zusw  ,           &  !    "         " 
    75          zgphsy, ztagny                       !    "         " 
    76       REAL(wp) :: & 
    77          zresm, zunw, zvnw, zur, zvr, zmod, za, zac, & 
    78          zmpzas, zstms, zindu, zindu1, zusdtp, zmassdt, zcorlal,  & 
    79          ztrace2, zdeter, zdelta, zsang, zmask, zdgp, zdgi, zdiag 
    80       REAL(wp),DIMENSION(jpi,jpj) ::   & 
    81          zpresh, zfrld, zmass, zcorl,     & 
    82          zu0, zv0, zviszeta, zviseta,     & 
    83          zc1u, zc1v, zc2u, zc2v, za1ct, za2ct, za1, za2, zb1, zb2,  & 
    84          zc1, zc2, zd1, zd2, zden, zu_ice, zv_ice, zresr 
    85       REAL(wp),DIMENSION(jpi,jpj,2,2) :: & 
    86          zs11, zs12, zs22, zs21 
    87       !!------------------------------------------------------------------- 
     92 
     93!!bug 
     94!!    ui_oce(:,:) = 0.e0 
     95!!    vi_oce(:,:) = 0.e0 
     96!!    write(*,*) 'rhg min, max u & v', maxval(ui_oce), minval(ui_oce), maxval(vi_oce), minval(vi_oce) 
     97!!bug 
    8898       
    8999      !  Store initial velocities 
    90       !  ------------------------ 
    91       zu0(:,:) = u_ice(:,:) 
    92       zv0(:,:) = v_ice(:,:) 
     100      !  ---------------- 
     101      zztms(:,0    ) = 0.e0       ;    zzfrld(:,0    ) = 0.e0 
     102      zztms(:,jpj+1) = 0.e0       ;    zzfrld(:,jpj+1) = 0.e0 
     103      zu0(:,0    ) = 0.e0         ;    zv0(:,0    ) = 0.e0 
     104      zu0(:,jpj+1) = 0.e0         ;    zv0(:,jpj+1) = 0.e0 
     105      zztms(:,1:jpj) = tms(:,:)   ;    zzfrld(:,1:jpj) = frld(:,:) 
     106      zu0(:,1:jpj) = ui_ice(:,:)   ;    zv0(:,1:jpj) = vi_ice(:,:) 
     107 
     108      zu_a(:,:)    = zu0(:,:)     ;   zv_a(:,:) = zv0(:,:) 
     109      zu_n(:,:)    = zu0(:,:)     ;   zv_n(:,:) = zv0(:,:) 
     110 
     111!i 
     112      zi1   (:,:) = 0.e0 
     113      zi2   (:,:) = 0.e0 
     114      zpresh(:,:) = 0.e0 
     115      zmasst(:,:) = 0.e0 
     116!i 
     117!!gm violant 
     118      zfrld(:,:) =0.e0 
     119      zcorl(:,:) =0.e0 
     120      zmass(:,:) =0.e0 
     121      za1ct(:,:) =0.e0 
     122      za2ct(:,:) =0.e0 
     123!!gm end 
     124 
     125      zviszeta(:,:) = 0.e0 
     126      zviseta (:,:) = 0.e0 
     127 
     128!i    zviszeta(:,0    ) = 0.e0    ;    zviseta(:,0    ) = 0.e0 
     129!i    zviszeta(:,jpj  ) = 0.e0    ;    zviseta(:,jpj  ) = 0.e0 
     130!i    zviszeta(:,jpj+1) = 0.e0    ;    zviseta(:,jpj+1) = 0.e0 
     131 
    93132 
    94133      ! Ice mass, ice strength, and wind stress at the center            | 
     
    96135      !------------------------------------------------------------------- 
    97136 
     137!CDIR NOVERRCHK 
    98138      DO jj = k_j1 , k_jpj-1 
     139!CDIR NOVERRCHK 
    99140         DO ji = 1 , jpi 
    100             za1(ji,jj)    = tms(ji,jj) * ( rhosn * hsnm(ji,jj) + rhoic * hicm(ji,jj) ) 
     141            ! only the sinus changes its sign with the hemisphere 
     142            zsang(ji,jj)  = SIGN( 1.e0, fcor(ji,jj) ) * sangvg   ! only the sinus changes its sign with the hemisphere 
     143            ! 
     144            zmasst(ji,jj) = tms(ji,jj) * ( rhosn * hsnm(ji,jj) + rhoic * hicm(ji,jj) ) 
    101145            zpresh(ji,jj) = tms(ji,jj) *  pstarh * hicm(ji,jj) * EXP( -c_rhg * frld(ji,jj) ) 
    102 #if defined key_lim_cp1 && defined key_coupled 
    103             zb1(ji,jj)    = tms(ji,jj) * gtaux(ji,jj) * ( 1.0 - frld(ji,jj) ) 
    104             zb2(ji,jj)    = tms(ji,jj) * gtauy(ji,jj) * ( 1.0 - frld(ji,jj) ) 
    105 #else 
    106             zb1(ji,jj)    = tms(ji,jj) * ( 1.0 - frld(ji,jj) ) 
    107             zb2(ji,jj)    = tms(ji,jj) * ( 1.0 - frld(ji,jj) ) 
    108 #endif 
     146!!gm  :: stress given at I-point (F-point for the ocean) only compute the ponderation with the ice fraction (1-frld) 
     147            zi1(ji,jj)    = tms(ji,jj) * ( 1.0 - frld(ji,jj) ) 
     148            zi2(ji,jj)    = tms(ji,jj) * ( 1.0 - frld(ji,jj) ) 
    109149         END DO 
    110150      END DO 
     
    117157          
    118158      DO jj = k_j1+1, k_jpj-1 
    119          DO ji = 2, jpi 
    120             zstms = tms(ji,jj  ) * wght(ji,jj,2,2) + tms(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    121                &  + tms(ji,jj-1) * wght(ji,jj,2,1) + tms(ji-1,jj-1) * wght(ji,jj,1,1) 
     159         DO ji = fs_2, jpi 
     160            zstms = zztms(ji,jj  ) * wght(ji,jj,2,2) + zztms(ji-1,jj  ) * wght(ji,jj,1,2)   & 
     161               &  + zztms(ji,jj-1) * wght(ji,jj,2,1) + zztms(ji-1,jj-1) * wght(ji,jj,1,1) 
    122162            zusw  = 1.0 / MAX( zstms, epsd ) 
    123163 
    124             zt11 = tms(ji  ,jj  ) * frld(ji  ,jj  )  
    125             zt12 = tms(ji-1,jj  ) * frld(ji-1,jj  )  
    126             zt21 = tms(ji  ,jj-1) * frld(ji  ,jj-1)  
    127             zt22 = tms(ji-1,jj-1) * frld(ji-1,jj-1) 
     164            zt11 = zztms(ji  ,jj  ) * zzfrld(ji  ,jj  )  
     165            zt12 = zztms(ji-1,jj  ) * zzfrld(ji-1,jj  )  
     166            zt21 = zztms(ji  ,jj-1) * zzfrld(ji  ,jj-1)  
     167            zt22 = zztms(ji-1,jj-1) * zzfrld(ji-1,jj-1) 
    128168 
    129169            ! Leads area. 
     
    131171               &             + zt21 * wght(ji,jj,2,1) + zt22 * wght(ji,jj,1,1) ) * zusw 
    132172 
    133             ! Mass and coriolis coeff. 
    134             zmass(ji,jj) = ( za1(ji,jj  ) * wght(ji,jj,2,2) + za1(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    135                &           + za1(ji,jj-1) * wght(ji,jj,2,1) + za1(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw 
     173            ! Mass and coriolis coeff. at I-point 
     174            zmass(ji,jj) = ( zmasst(ji,jj  ) * wght(ji,jj,2,2) + zmasst(ji-1,jj  ) * wght(ji,jj,1,2)   & 
     175               &           + zmasst(ji,jj-1) * wght(ji,jj,2,1) + zmasst(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw 
    136176            zcorl(ji,jj) = zmass(ji,jj) * fcor(ji,jj) 
    137177 
    138178            ! Wind stress. 
    139 #if defined key_lim_cp1 && defined key_coupled 
    140             ztagnx = ( zb1(ji,jj  ) * wght(ji,jj,2,2) + zb1(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    141                &     + zb1(ji,jj-1) * wght(ji,jj,2,1) + zb1(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw 
    142             ztagny = ( zb2(ji,jj  ) * wght(ji,jj,2,2) + zb2(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    143                &     + zb2(ji,jj-1) * wght(ji,jj,2,1) + zb2(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw 
    144 #else 
    145             ztagnx = ( zb1(ji,jj  ) * wght(ji,jj,2,2) + zb1(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    146                &     + zb1(ji,jj-1) * wght(ji,jj,2,1) + zb1(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw * gtaux(ji,jj) 
    147             ztagny = ( zb2(ji,jj  ) * wght(ji,jj,2,2) + zb2(ji-1,jj  ) * wght(ji,jj,1,2)   & 
    148                &     + zb2(ji,jj-1) * wght(ji,jj,2,1) + zb2(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw * gtauy(ji,jj) 
    149 #endif 
     179            ! always provide stress at I-point (ocean F-point) 
     180            ztagnx = ( zi1(ji,jj  ) * wght(ji,jj,2,2) + zi1(ji-1,jj  ) * wght(ji,jj,1,2)   & 
     181               &     + zi1(ji,jj-1) * wght(ji,jj,2,1) + zi1(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw * utaui_ice(ji,jj) 
     182            ztagny = ( zi2(ji,jj  ) * wght(ji,jj,2,2) + zi2(ji-1,jj  ) * wght(ji,jj,1,2)   & 
     183               &     + zi2(ji,jj-1) * wght(ji,jj,2,1) + zi2(ji-1,jj-1) * wght(ji,jj,1,1) ) * zusw * vtaui_ice(ji,jj) 
    150184 
    151185            ! Gradient of ice strength 
     
    161195 
    162196            ! Computation of the velocity field taking into account the ice-ice interaction.                                  
    163             ! Terms that are independent of the velocity field. 
    164             za1ct(ji,jj) = ztagnx - zcorl(ji,jj) * v_oce(ji,jj) - zgphsx 
    165             za2ct(ji,jj) = ztagny + zcorl(ji,jj) * u_oce(ji,jj) - zgphsy 
     197            ! Terms that are independent of the ice velocity field. 
     198            za1ct(ji,jj) = ztagnx - zcorl(ji,jj) * vi_oce(ji,jj) - zgphsx 
     199            za2ct(ji,jj) = ztagny + zcorl(ji,jj) * ui_oce(ji,jj) - zgphsy 
    166200         END DO 
    167201      END DO 
    168  
    169 !! inutile!! 
    170 !!??    CALL lbc_lnk( za1ct, 'I', -1. ) 
    171 !!??    CALL lbc_lnk( za2ct, 'I', -1. ) 
    172202 
    173203 
     
    182212         ! Computation of free drift field for free slip boundary conditions. 
    183213 
    184            DO jj = k_j1, k_jpj-1 
    185               DO ji = 1, jpim1 
    186                  !- Rate of strain tensor. 
    187                  zt11 =   akappa(ji,jj,1,1) * ( u_ice(ji+1,jj) + u_ice(ji+1,jj+1) - u_ice(ji,jj  ) - u_ice(ji  ,jj+1) )  & 
    188                     &   + akappa(ji,jj,1,2) * ( v_ice(ji+1,jj) + v_ice(ji+1,jj+1) + v_ice(ji,jj  ) + v_ice(ji  ,jj+1) ) 
    189                  zt12 = - akappa(ji,jj,2,2) * ( u_ice(ji  ,jj) + u_ice(ji+1,jj  ) - u_ice(ji,jj+1) - u_ice(ji+1,jj+1) )  & 
    190                     &   - akappa(ji,jj,2,1) * ( v_ice(ji  ,jj) + v_ice(ji+1,jj  ) + v_ice(ji,jj+1) + v_ice(ji+1,jj+1) ) 
    191                  zt22 = - akappa(ji,jj,2,2) * ( v_ice(ji  ,jj) + v_ice(ji+1,jj  ) - v_ice(ji,jj+1) - v_ice(ji+1,jj+1) )  & 
    192                     &   + akappa(ji,jj,2,1) * ( u_ice(ji  ,jj) + u_ice(ji+1,jj  ) + u_ice(ji,jj+1) + u_ice(ji+1,jj+1) ) 
    193                  zt21 =   akappa(ji,jj,1,1) * ( v_ice(ji+1,jj) + v_ice(ji+1,jj+1) - v_ice(ji,jj  ) - v_ice(ji  ,jj+1) )  & 
    194                     &   - akappa(ji,jj,1,2) * ( u_ice(ji+1,jj) + u_ice(ji+1,jj+1) + u_ice(ji,jj  ) + u_ice(ji  ,jj+1) ) 
    195  
    196                  !- Rate of strain tensor.  
    197                  zdgp = zt11 + zt22 
    198                  zdgi = zt12 + zt21 
    199                  ztrace2 = zdgp * zdgp  
    200                  zdeter  = zt11 * zt22 - 0.25 * zdgi * zdgi 
    201  
    202                  !  Creep limit depends on the size of the grid. 
    203                  zdelta = MAX( SQRT( ztrace2 + ( ztrace2 - 4.0 * zdeter ) * usecc2),  creepl) 
    204  
    205                  !-  Computation of viscosities. 
    206                  zviszeta(ji,jj) = MAX( zpresh(ji,jj) / zdelta, etamn ) 
    207                  zviseta (ji,jj) = zviszeta(ji,jj) * usecc2 
    208               END DO 
    209            END DO 
    210 !!??       CALL lbc_lnk( zviszeta, 'I', -1. )  ! or T point???   semble reellement inutile 
    211 !!??       CALL lbc_lnk( zviseta , 'I', -1. ) 
    212  
    213  
    214            !-  Determination of zc1u, zc2u, zc1v and zc2v. 
    215            DO jj = k_j1+1, k_jpj-1 
    216               DO ji = 2, jpim1 
    217                  ze11   =  akappa(ji-1,jj-1,1,1) 
    218                  ze12   = +akappa(ji-1,jj-1,2,2) 
    219                  ze22   =  akappa(ji-1,jj-1,2,1) 
    220                  ze21   = -akappa(ji-1,jj-1,1,2) 
    221                  zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
    222                  zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
    223                  zvis12 =       zviseta (ji-1,jj-1) + dm 
    224                  zvis21 =       zviseta (ji-1,jj-1) 
    225  
    226                  zdiag = zvis22 * ( ze11 + ze22 ) 
    227                  zs11(ji,jj,1,1) =  zvis11 * ze11 + zdiag 
    228                  zs12(ji,jj,1,1) =  zvis12 * ze12 + zvis21 * ze21 
    229                  zs22(ji,jj,1,1) =  zvis11 * ze22 + zdiag 
    230                  zs21(ji,jj,1,1) =  zvis12 * ze21 + zvis21 * ze12 
    231  
    232                  ze11   = -akappa(ji,jj-1,1,1) 
    233                  ze12   = +akappa(ji,jj-1,2,2) 
    234                  ze22   =  akappa(ji,jj-1,2,1) 
    235                  ze21   = -akappa(ji,jj-1,1,2) 
    236                  zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
    237                  zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
    238                  zvis12 =       zviseta (ji,jj-1) + dm 
    239                  zvis21 =       zviseta (ji,jj-1) 
    240  
    241                  zdiag = zvis22 * ( ze11 + ze22 ) 
    242                  zs11(ji,jj,2,1) =  zvis11 * ze11 + zdiag 
    243                  zs12(ji,jj,2,1) =  zvis12 * ze12 + zvis21 * ze21 
    244                  zs22(ji,jj,2,1) =  zvis11 * ze22 + zdiag 
    245                  zs21(ji,jj,2,1) =  zvis12 * ze21 + zvis21 * ze12 
    246  
    247                  ze11   =  akappa(ji-1,jj,1,1) 
    248                  ze12   = -akappa(ji-1,jj,2,2) 
    249                  ze22   =  akappa(ji-1,jj,2,1) 
    250                  ze21   = -akappa(ji-1,jj,1,2) 
    251                  zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
    252                  zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
    253                  zvis12 =       zviseta (ji-1,jj) + dm 
    254                  zvis21 =       zviseta (ji-1,jj) 
    255  
    256                  zdiag = zvis22 * ( ze11 + ze22 )  
    257                  zs11(ji,jj,1,2) =  zvis11 * ze11 + zdiag 
    258                  zs12(ji,jj,1,2) =  zvis12 * ze12 + zvis21 * ze21 
    259                  zs22(ji,jj,1,2) =  zvis11 * ze22 + zdiag 
    260                  zs21(ji,jj,1,2) =  zvis12 * ze21 + zvis21 * ze12 
    261  
    262                  ze11   = -akappa(ji,jj,1,1) 
    263                  ze12   = -akappa(ji,jj,2,2) 
    264                  ze22   =  akappa(ji,jj,2,1) 
    265                  ze21   = -akappa(ji,jj,1,2) 
    266                  zvis11 = 2.0 * zviseta (ji,jj) + dm 
    267                  zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
    268                  zvis12 =       zviseta (ji,jj) + dm 
    269                  zvis21 =       zviseta (ji,jj) 
    270  
    271                  zdiag = zvis22 * ( ze11 + ze22 ) 
    272                  zs11(ji,jj,2,2) =  zvis11 * ze11 + zdiag 
    273                  zs12(ji,jj,2,2) =  zvis12 * ze12 + zvis21 * ze21 
    274                  zs22(ji,jj,2,2) =  zvis11 * ze22 + zdiag  
    275                  zs21(ji,jj,2,2) =  zvis12 * ze21 + zvis21 * ze12 
    276               END DO 
    277            END DO 
    278  
    279            DO jj = k_j1+1, k_jpj-1 
    280               DO ji = 2, jpim1 
    281                  zc1u(ji,jj) =   & 
    282                     + alambd(ji,jj,2,2,2,1) * zs11(ji,jj,2,1) + alambd(ji,jj,2,2,2,2) * zs11(ji,jj,2,2)   & 
    283                     - alambd(ji,jj,2,2,1,1) * zs11(ji,jj,1,1) - alambd(ji,jj,2,2,1,2) * zs11(ji,jj,1,2)   & 
    284                     - alambd(ji,jj,1,1,2,1) * zs12(ji,jj,2,1) - alambd(ji,jj,1,1,1,1) * zs12(ji,jj,1,1)   & 
    285                     + alambd(ji,jj,1,1,2,2) * zs12(ji,jj,2,2) + alambd(ji,jj,1,1,1,2) * zs12(ji,jj,1,2)   & 
    286                     + alambd(ji,jj,1,2,1,1) * zs21(ji,jj,1,1) + alambd(ji,jj,1,2,2,1) * zs21(ji,jj,2,1)   & 
    287                     + alambd(ji,jj,1,2,1,2) * zs21(ji,jj,1,2) + alambd(ji,jj,1,2,2,2) * zs21(ji,jj,2,2)   & 
    288                     - alambd(ji,jj,2,1,1,1) * zs22(ji,jj,1,1) - alambd(ji,jj,2,1,2,1) * zs22(ji,jj,2,1)   & 
    289                     - alambd(ji,jj,2,1,1,2) * zs22(ji,jj,1,2) - alambd(ji,jj,2,1,2,2) * zs22(ji,jj,2,2) 
    290                   
    291                  zc2u(ji,jj) =   & 
    292                     + alambd(ji,jj,2,2,2,1) * zs21(ji,jj,2,1) + alambd(ji,jj,2,2,2,2) * zs21(ji,jj,2,2)   & 
    293                     - alambd(ji,jj,2,2,1,1) * zs21(ji,jj,1,1) - alambd(ji,jj,2,2,1,2) * zs21(ji,jj,1,2)   & 
    294                     - alambd(ji,jj,1,1,2,1) * zs22(ji,jj,2,1) - alambd(ji,jj,1,1,1,1) * zs22(ji,jj,1,1)   & 
    295                     + alambd(ji,jj,1,1,2,2) * zs22(ji,jj,2,2) + alambd(ji,jj,1,1,1,2) * zs22(ji,jj,1,2)   & 
    296                     - alambd(ji,jj,1,2,1,1) * zs11(ji,jj,1,1) - alambd(ji,jj,1,2,2,1) * zs11(ji,jj,2,1)   & 
    297                     - alambd(ji,jj,1,2,1,2) * zs11(ji,jj,1,2) - alambd(ji,jj,1,2,2,2) * zs11(ji,jj,2,2)   & 
    298                     + alambd(ji,jj,2,1,1,1) * zs12(ji,jj,1,1) + alambd(ji,jj,2,1,2,1) * zs12(ji,jj,2,1)   & 
    299                     + alambd(ji,jj,2,1,1,2) * zs12(ji,jj,1,2) + alambd(ji,jj,2,1,2,2) * zs12(ji,jj,2,2) 
    300              END DO 
    301            END DO 
    302  
    303            DO jj = k_j1+1, k_jpj-1 
    304               DO ji = 2, jpim1 
    305                  !  zc1v , zc2v. 
    306                  ze11   =  akappa(ji-1,jj-1,1,2) 
    307                  ze12   = -akappa(ji-1,jj-1,2,1) 
    308                  ze22   = +akappa(ji-1,jj-1,2,2) 
    309                  ze21   =  akappa(ji-1,jj-1,1,1) 
    310                  zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
    311                  zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
    312                  zvis12 =       zviseta (ji-1,jj-1) + dm 
    313                  zvis21 =       zviseta (ji-1,jj-1) 
    314  
    315                  zdiag = zvis22 * ( ze11 + ze22 ) 
    316                  zs11(ji,jj,1,1) =  zvis11 * ze11 + zdiag 
    317                  zs12(ji,jj,1,1) =  zvis12 * ze12 + zvis21 * ze21 
    318                  zs22(ji,jj,1,1) =  zvis11 * ze22 + zdiag 
    319                  zs21(ji,jj,1,1) =  zvis12 * ze21 + zvis21 * ze12 
     214!CDIR NOVERRCHK 
     215         DO jj = k_j1, k_jpj-1 
     216!CDIR NOVERRCHK 
     217            DO ji = 1, fs_jpim1 
     218               !- Rate of strain tensor. 
     219               zt11 =   akappa(ji,jj,1,1) * ( zu_a(ji+1,jj) + zu_a(ji+1,jj+1) - zu_a(ji,jj  ) - zu_a(ji  ,jj+1) )  & 
     220                  &   + akappa(ji,jj,1,2) * ( zv_a(ji+1,jj) + zv_a(ji+1,jj+1) + zv_a(ji,jj  ) + zv_a(ji  ,jj+1) ) 
     221               zt12 = - akappa(ji,jj,2,2) * ( zu_a(ji  ,jj) + zu_a(ji+1,jj  ) - zu_a(ji,jj+1) - zu_a(ji+1,jj+1) )  & 
     222                  &   - akappa(ji,jj,2,1) * ( zv_a(ji  ,jj) + zv_a(ji+1,jj  ) + zv_a(ji,jj+1) + zv_a(ji+1,jj+1) ) 
     223               zt22 = - akappa(ji,jj,2,2) * ( zv_a(ji  ,jj) + zv_a(ji+1,jj  ) - zv_a(ji,jj+1) - zv_a(ji+1,jj+1) )  & 
     224                  &   + akappa(ji,jj,2,1) * ( zu_a(ji  ,jj) + zu_a(ji+1,jj  ) + zu_a(ji,jj+1) + zu_a(ji+1,jj+1) ) 
     225               zt21 =   akappa(ji,jj,1,1) * ( zv_a(ji+1,jj) + zv_a(ji+1,jj+1) - zv_a(ji,jj  ) - zv_a(ji  ,jj+1) )  & 
     226                  &   - akappa(ji,jj,1,2) * ( zu_a(ji+1,jj) + zu_a(ji+1,jj+1) + zu_a(ji,jj  ) + zu_a(ji  ,jj+1) ) 
     227 
     228               !- Rate of strain tensor.  
     229               zdgp = zt11 + zt22 
     230               zdgi = zt12 + zt21 
     231               ztrace2 = zdgp * zdgp  
     232               zdeter  = zt11 * zt22 - 0.25 * zdgi * zdgi 
     233 
     234               !  Creep limit depends on the size of the grid. 
     235               zdelta = MAX( SQRT( ztrace2 + ( ztrace2 - 4.0 * zdeter ) * usecc2 ),  creepl) 
     236 
     237               !-  Computation of viscosities. 
     238               zviszeta(ji,jj) = MAX( zpresh(ji,jj) / zdelta, etamn ) 
     239               zviseta (ji,jj) = zviszeta(ji,jj) * usecc2 
     240            END DO 
     241         END DO 
     242 
     243         !-  Determination of zc1u, zc2u, zc1v and zc2v. 
     244         DO jj = k_j1+1, k_jpj-1 
     245            DO ji = fs_2, fs_jpim1 
     246               !* zc1u , zc2v 
     247               zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
     248               zvis12 =       zviseta (ji-1,jj-1) + dm 
     249               zvis21 =       zviseta (ji-1,jj-1) 
     250               zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
     251               zdiag  = zvis22 * ( akappa(ji-1,jj-1,1,1) + akappa(ji-1,jj-1,2,1) ) 
     252               zs11_11 =  zvis11 * akappa(ji-1,jj-1,1,1) + zdiag 
     253               zs12_11 =  zvis12 * akappa(ji-1,jj-1,2,2) - zvis21 * akappa(ji-1,jj-1,1,2) 
     254               zs21_11 = -zvis12 * akappa(ji-1,jj-1,1,2) + zvis21 * akappa(ji-1,jj-1,2,2) 
     255               zs22_11 =  zvis11 * akappa(ji-1,jj-1,2,1) + zdiag 
     256 
     257               zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
     258               zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     259               zvis12 =       zviseta (ji,jj-1) + dm 
     260               zvis21 =       zviseta (ji,jj-1) 
     261               zdiag = zvis22 * ( -akappa(ji,jj-1,1,1) + akappa(ji,jj-1,2,1) ) 
     262               zs11_21 = -zvis11 * akappa(ji,jj-1,1,1) + zdiag 
     263               zs12_21 =  zvis12 * akappa(ji,jj-1,2,2) - zvis21 * akappa(ji,jj-1,1,2) 
     264               zs22_21 =  zvis11 * akappa(ji,jj-1,2,1) + zdiag 
     265               zs21_21 = -zvis12 * akappa(ji,jj-1,1,2) + zvis21 * akappa(ji,jj-1,2,2) 
     266 
     267               zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
     268               zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     269               zvis12 =       zviseta (ji-1,jj) + dm 
     270               zvis21 =       zviseta (ji-1,jj) 
     271               zdiag = zvis22 * ( akappa(ji-1,jj,1,1) + akappa(ji-1,jj,2,1) ) 
     272               zs11_12 =  zvis11 * akappa(ji-1,jj,1,1) + zdiag 
     273               zs12_12 = -zvis12 * akappa(ji-1,jj,2,2) - zvis21 * akappa(ji-1,jj,1,2) 
     274               zs22_12 =  zvis11 * akappa(ji-1,jj,2,1) + zdiag 
     275               zs21_12 = -zvis12 * akappa(ji-1,jj,1,2) - zvis21 * akappa(ji-1,jj,2,2) 
     276 
     277               zvis11 = 2.0 * zviseta (ji,jj) + dm 
     278               zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
     279               zvis12 =       zviseta (ji,jj) + dm 
     280               zvis21 =       zviseta (ji,jj) 
     281               zdiag = zvis22 * ( -akappa(ji,jj,1,1) + akappa(ji,jj,2,1) ) 
     282               zs11_22 = -zvis11 * akappa(ji,jj,1,1) + zdiag 
     283               zs12_22 = -zvis12 * akappa(ji,jj,2,2) - zvis21 * akappa(ji,jj,1,2) 
     284               zs22_22 =  zvis11 * akappa(ji,jj,2,1) + zdiag 
     285               zs21_22 = -zvis12 * akappa(ji,jj,1,2) - zvis21 * akappa(ji,jj,2,2) 
     286 
     287               zc1u(ji,jj) = + alambd(ji,jj,2,2,2,1) * zs11_21 + alambd(ji,jj,2,2,2,2) * zs11_22   & 
     288                  &          - alambd(ji,jj,2,2,1,1) * zs11_11 - alambd(ji,jj,2,2,1,2) * zs11_12   & 
     289                  &          - alambd(ji,jj,1,1,2,1) * zs12_21 - alambd(ji,jj,1,1,1,1) * zs12_11   & 
     290                  &          + alambd(ji,jj,1,1,2,2) * zs12_22 + alambd(ji,jj,1,1,1,2) * zs12_12   & 
     291                  &          + alambd(ji,jj,1,2,1,1) * zs21_11 + alambd(ji,jj,1,2,2,1) * zs21_21   & 
     292                  &          + alambd(ji,jj,1,2,1,2) * zs21_12 + alambd(ji,jj,1,2,2,2) * zs21_22   & 
     293                  &          - alambd(ji,jj,2,1,1,1) * zs22_11 - alambd(ji,jj,2,1,2,1) * zs22_21   & 
     294                  &          - alambd(ji,jj,2,1,1,2) * zs22_12 - alambd(ji,jj,2,1,2,2) * zs22_22 
     295 
     296               zc2u(ji,jj) = + alambd(ji,jj,2,2,2,1) * zs21_21 + alambd(ji,jj,2,2,2,2) * zs21_22   & 
     297                  &          - alambd(ji,jj,2,2,1,1) * zs21_11 - alambd(ji,jj,2,2,1,2) * zs21_12   & 
     298                  &          - alambd(ji,jj,1,1,2,1) * zs22_21 - alambd(ji,jj,1,1,1,1) * zs22_11   & 
     299                  &          + alambd(ji,jj,1,1,2,2) * zs22_22 + alambd(ji,jj,1,1,1,2) * zs22_12   & 
     300                  &          - alambd(ji,jj,1,2,1,1) * zs11_11 - alambd(ji,jj,1,2,2,1) * zs11_21   & 
     301                  &          - alambd(ji,jj,1,2,1,2) * zs11_12 - alambd(ji,jj,1,2,2,2) * zs11_22   & 
     302                  &          + alambd(ji,jj,2,1,1,1) * zs12_11 + alambd(ji,jj,2,1,2,1) * zs12_21   & 
     303                  &          + alambd(ji,jj,2,1,1,2) * zs12_12 + alambd(ji,jj,2,1,2,2) * zs12_22 
     304 
     305               !* zc1v , zc2v. 
     306               zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
     307               zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
     308               zvis12 =       zviseta (ji-1,jj-1) + dm 
     309               zvis21 =       zviseta (ji-1,jj-1) 
     310               zdiag = zvis22 * ( akappa(ji-1,jj-1,1,2) + akappa(ji-1,jj-1,2,2) ) 
     311               zs11_11 =  zvis11 * akappa(ji-1,jj-1,1,2) + zdiag 
     312               zs12_11 = -zvis12 * akappa(ji-1,jj-1,2,1) + zvis21 * akappa(ji-1,jj-1,1,1) 
     313               zs22_11 =  zvis11 * akappa(ji-1,jj-1,2,2) + zdiag 
     314               zs21_11 =  zvis12 * akappa(ji-1,jj-1,1,1) - zvis21 * akappa(ji-1,jj-1,2,1) 
    320315  
    321                  ze11   =  akappa(ji,jj-1,1,2) 
    322                  ze12   = -akappa(ji,jj-1,2,1) 
    323                  ze22   = +akappa(ji,jj-1,2,2) 
    324                  ze21   = -akappa(ji,jj-1,1,1) 
    325                  zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
    326                  zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
    327                  zvis12 =       zviseta (ji,jj-1) + dm 
    328                  zvis21 =       zviseta (ji,jj-1) 
    329  
    330                  zdiag = zvis22 * ( ze11 + ze22 ) 
    331                  zs11(ji,jj,2,1) =  zvis11 * ze11 + zdiag 
    332                  zs12(ji,jj,2,1) =  zvis12 * ze12 + zvis21 * ze21 
    333                  zs22(ji,jj,2,1) =  zvis11 * ze22 + zdiag 
    334                  zs21(ji,jj,2,1) =  zvis12 * ze21 + zvis21 * ze12 
    335  
    336                  ze11   =  akappa(ji-1,jj,1,2) 
    337                  ze12   = -akappa(ji-1,jj,2,1) 
    338                  ze22   = -akappa(ji-1,jj,2,2) 
    339                  ze21   =  akappa(ji-1,jj,1,1) 
    340                  zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
    341                  zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
    342                  zvis12 =       zviseta (ji-1,jj) + dm 
    343                  zvis21 =       zviseta (ji-1,jj) 
    344  
    345                  zdiag = zvis22 * ( ze11 + ze22 ) 
    346                  zs11(ji,jj,1,2) =  zvis11 * ze11 + zdiag 
    347                  zs12(ji,jj,1,2) =  zvis12 * ze12 + zvis21 * ze21 
    348                  zs22(ji,jj,1,2) =  zvis11 * ze22 + zdiag 
    349                  zs21(ji,jj,1,2) =  zvis12 * ze21 + zvis21 * ze12 
    350  
    351                  ze11   =  akappa(ji,jj,1,2) 
    352                  ze12   = -akappa(ji,jj,2,1) 
    353                  ze22   = -akappa(ji,jj,2,2) 
    354                  ze21   = -akappa(ji,jj,1,1) 
    355                  zvis11 = 2.0 * zviseta (ji,jj) + dm 
    356                  zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
    357                  zvis12 =       zviseta (ji,jj) + dm 
    358                  zvis21 =       zviseta (ji,jj) 
    359  
    360                  zdiag = zvis22 * ( ze11 + ze22 ) 
    361                  zs11(ji,jj,2,2) =  zvis11 * ze11 + zdiag 
    362                  zs12(ji,jj,2,2) =  zvis12 * ze12 + zvis21 * ze21 
    363                  zs22(ji,jj,2,2) =  zvis11 * ze22 + zdiag 
    364                  zs21(ji,jj,2,2) =  zvis12 * ze21 + zvis21 * ze12 
    365  
    366               END DO 
    367            END DO 
    368  
    369            DO jj = k_j1+1, k_jpj-1 
    370               DO ji = 2, jpim1 
    371                  zc1v(ji,jj) =   & 
    372                     + alambd(ji,jj,2,2,2,1) * zs11(ji,jj,2,1) + alambd(ji,jj,2,2,2,2) * zs11(ji,jj,2,2)   & 
    373                     - alambd(ji,jj,2,2,1,1) * zs11(ji,jj,1,1) - alambd(ji,jj,2,2,1,2) * zs11(ji,jj,1,2)   & 
    374                     - alambd(ji,jj,1,1,2,1) * zs12(ji,jj,2,1) - alambd(ji,jj,1,1,1,1) * zs12(ji,jj,1,1)   & 
    375                     + alambd(ji,jj,1,1,2,2) * zs12(ji,jj,2,2) + alambd(ji,jj,1,1,1,2) * zs12(ji,jj,1,2)   & 
    376                     + alambd(ji,jj,1,2,1,1) * zs21(ji,jj,1,1) + alambd(ji,jj,1,2,2,1) * zs21(ji,jj,2,1)   & 
    377                     + alambd(ji,jj,1,2,1,2) * zs21(ji,jj,1,2) + alambd(ji,jj,1,2,2,2) * zs21(ji,jj,2,2)   & 
    378                     - alambd(ji,jj,2,1,1,1) * zs22(ji,jj,1,1) - alambd(ji,jj,2,1,2,1) * zs22(ji,jj,2,1)   & 
    379                     - alambd(ji,jj,2,1,1,2) * zs22(ji,jj,1,2) - alambd(ji,jj,2,1,2,2) * zs22(ji,jj,2,2) 
    380                  zc2v(ji,jj) =   & 
    381                     + alambd(ji,jj,2,2,2,1) * zs21(ji,jj,2,1) + alambd(ji,jj,2,2,2,2) * zs21(ji,jj,2,2)   & 
    382                     - alambd(ji,jj,2,2,1,1) * zs21(ji,jj,1,1) - alambd(ji,jj,2,2,1,2) * zs21(ji,jj,1,2)   & 
    383                     - alambd(ji,jj,1,1,2,1) * zs22(ji,jj,2,1) - alambd(ji,jj,1,1,1,1) * zs22(ji,jj,1,1)   & 
    384                     + alambd(ji,jj,1,1,2,2) * zs22(ji,jj,2,2) + alambd(ji,jj,1,1,1,2) * zs22(ji,jj,1,2)   & 
    385                     - alambd(ji,jj,1,2,1,1) * zs11(ji,jj,1,1) - alambd(ji,jj,1,2,2,1) * zs11(ji,jj,2,1)   & 
    386                     - alambd(ji,jj,1,2,1,2) * zs11(ji,jj,1,2) - alambd(ji,jj,1,2,2,2) * zs11(ji,jj,2,2)   & 
    387                     + alambd(ji,jj,2,1,1,1) * zs12(ji,jj,1,1) + alambd(ji,jj,2,1,2,1) * zs12(ji,jj,2,1)   & 
    388                     + alambd(ji,jj,2,1,1,2) * zs12(ji,jj,1,2) + alambd(ji,jj,2,1,2,2) * zs12(ji,jj,2,2) 
    389               END DO 
    390            END DO 
    391  
    392          ! Relaxation. 
    393             
    394 iflag:   DO jter = 1 , nbitdr 
    395  
    396             !  Store previous drift field.    
    397             DO jj = k_j1, k_jpj-1 
    398                zu_ice(:,jj) = u_ice(:,jj) 
    399                zv_ice(:,jj) = v_ice(:,jj) 
     316               zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
     317               zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     318               zvis12 =       zviseta (ji,jj-1) + dm 
     319               zvis21 =       zviseta (ji,jj-1) 
     320               zdiag = zvis22 * ( akappa(ji,jj-1,1,2) + akappa(ji,jj-1,2,2) ) 
     321               zs11_21 =  zvis11 * akappa(ji,jj-1,1,2) + zdiag 
     322               zs12_21 = -zvis12 * akappa(ji,jj-1,2,1) - zvis21 * akappa(ji,jj-1,1,1) 
     323               zs22_21 =  zvis11 * akappa(ji,jj-1,2,2) + zdiag 
     324               zs21_21 = -zvis12 * akappa(ji,jj-1,1,1) - zvis21 * akappa(ji,jj-1,2,1) 
     325 
     326               zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
     327               zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     328               zvis12 =       zviseta (ji-1,jj) + dm 
     329               zvis21 =       zviseta (ji-1,jj) 
     330               zdiag = zvis22 * ( akappa(ji-1,jj,1,2) - akappa(ji-1,jj,2,2) ) 
     331               zs11_12 =  zvis11 * akappa(ji-1,jj,1,2) + zdiag 
     332               zs12_12 = -zvis12 * akappa(ji-1,jj,2,1) + zvis21 * akappa(ji-1,jj,1,1) 
     333               zs22_12 = -zvis11 * akappa(ji-1,jj,2,2) + zdiag 
     334               zs21_12 =  zvis12 * akappa(ji-1,jj,1,1) - zvis21 * akappa(ji-1,jj,2,1) 
     335 
     336               zvis11 = 2.0 * zviseta (ji,jj) + dm 
     337               zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
     338               zvis12 =       zviseta (ji,jj) + dm 
     339               zvis21 =       zviseta (ji,jj) 
     340               zdiag = zvis22 * ( akappa(ji,jj,1,2) - akappa(ji,jj,2,2) ) 
     341               zs11_22 =  zvis11 * akappa(ji,jj,1,2) + zdiag 
     342               zs12_22 = -zvis12 * akappa(ji,jj,2,1) - zvis21 * akappa(ji,jj,1,1) 
     343               zs22_22 = -zvis11 * akappa(ji,jj,2,2) + zdiag 
     344               zs21_22 = -zvis12 * akappa(ji,jj,1,1) - zvis21 * akappa(ji,jj,2,1) 
     345 
     346               zc1v(ji,jj) = + alambd(ji,jj,2,2,2,1) * zs11_21 + alambd(ji,jj,2,2,2,2) * zs11_22   & 
     347                  &          - alambd(ji,jj,2,2,1,1) * zs11_11 - alambd(ji,jj,2,2,1,2) * zs11_12   & 
     348                  &          - alambd(ji,jj,1,1,2,1) * zs12_21 - alambd(ji,jj,1,1,1,1) * zs12_11   & 
     349                  &          + alambd(ji,jj,1,1,2,2) * zs12_22 + alambd(ji,jj,1,1,1,2) * zs12_12   & 
     350                  &          + alambd(ji,jj,1,2,1,1) * zs21_11 + alambd(ji,jj,1,2,2,1) * zs21_21   & 
     351                  &          + alambd(ji,jj,1,2,1,2) * zs21_12 + alambd(ji,jj,1,2,2,2) * zs21_22   & 
     352                  &          - alambd(ji,jj,2,1,1,1) * zs22_11 - alambd(ji,jj,2,1,2,1) * zs22_21   & 
     353                  &          - alambd(ji,jj,2,1,1,2) * zs22_12 - alambd(ji,jj,2,1,2,2) * zs22_22 
     354 
     355               zc2v(ji,jj) = + alambd(ji,jj,2,2,2,1) * zs21_21 + alambd(ji,jj,2,2,2,2) * zs21_22   & 
     356                  &          - alambd(ji,jj,2,2,1,1) * zs21_11 - alambd(ji,jj,2,2,1,2) * zs21_12   & 
     357                  &          - alambd(ji,jj,1,1,2,1) * zs22_21 - alambd(ji,jj,1,1,1,1) * zs22_11   & 
     358                  &          + alambd(ji,jj,1,1,2,2) * zs22_22 + alambd(ji,jj,1,1,1,2) * zs22_12   & 
     359                  &          - alambd(ji,jj,1,2,1,1) * zs11_11 - alambd(ji,jj,1,2,2,1) * zs11_21   & 
     360                  &          - alambd(ji,jj,1,2,1,2) * zs11_12 - alambd(ji,jj,1,2,2,2) * zs11_22   & 
     361                  &          + alambd(ji,jj,2,1,1,1) * zs12_11 + alambd(ji,jj,2,1,2,1) * zs12_21   & 
     362                  &          + alambd(ji,jj,2,1,1,2) * zs12_12 + alambd(ji,jj,2,1,2,2) * zs12_22 
    400363            END DO 
    401  
     364         END DO 
     365 
     366         ! GAUSS-SEIDEL method 
     367         !                                                      ! ================ ! 
     368iflag:   DO jter = 1 , nbitdr                                   !    Relaxation    ! 
     369            !                                                   ! ================ ! 
     370!CDIR NOVERRCHK 
    402371            DO jj = k_j1+1, k_jpj-1 
    403                zsang  = SIGN( 1.e0, fcor(1,jj) ) * sangvg   ! only the sinus changes its sign with the hemisphere 
    404                DO ji = 2, jpim1 
    405                  zur     = u_ice(ji,jj) - u_oce(ji,jj) 
    406                  zvr     = v_ice(ji,jj) - v_oce(ji,jj) 
    407                  zmod    = SQRT( zur * zur + zvr * zvr) * ( 1.0 - zfrld(ji,jj) ) 
    408                  za      = rhoco * zmod 
    409                  zac     = za * cangvg 
    410                   zmpzas  = alpha * zcorl(ji,jj) + za * zsang 
     372!CDIR NOVERRCHK 
     373               DO ji = fs_2, fs_jpim1 
     374                  ! 
     375                  ze11 =   akappa(ji,jj-1,1,1) * zu_a(ji+1,jj) + akappa(ji,jj-1,1,2) * zv_a(ji+1,jj) 
     376                  ze12 = + akappa(ji,jj-1,2,2) * zu_a(ji+1,jj) - akappa(ji,jj-1,2,1) * zv_a(ji+1,jj) 
     377                  ze22 = + akappa(ji,jj-1,2,2) * zv_a(ji+1,jj) + akappa(ji,jj-1,2,1) * zu_a(ji+1,jj) 
     378                  ze21 =   akappa(ji,jj-1,1,1) * zv_a(ji+1,jj) - akappa(ji,jj-1,1,2) * zu_a(ji+1,jj) 
     379                  zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
     380                  zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     381                  zvis12 =       zviseta (ji,jj-1) + dm 
     382                  zvis21 =       zviseta (ji,jj-1) 
     383                  zdiag = zvis22 * ( ze11 + ze22 ) 
     384                  zs11_21 =  zvis11 * ze11 + zdiag 
     385                  zs12_21 =  zvis12 * ze12 + zvis21 * ze21 
     386                  zs22_21 =  zvis11 * ze22 + zdiag 
     387                  zs21_21 =  zvis12 * ze21 + zvis21 * ze12 
     388 
     389                  ze11 =   akappa(ji-1,jj,1,1) * ( zu_a(ji  ,jj+1) - zu_a(ji-1,jj+1) )   & 
     390                     &   + akappa(ji-1,jj,1,2) * ( zv_a(ji  ,jj+1) + zv_a(ji-1,jj+1) ) 
     391                  ze12 = + akappa(ji-1,jj,2,2) * ( zu_a(ji-1,jj+1) + zu_a(ji  ,jj+1) )   & 
     392                     &   - akappa(ji-1,jj,2,1) * ( zv_a(ji-1,jj+1) + zv_a(ji  ,jj+1) ) 
     393                  ze22 = + akappa(ji-1,jj,2,2) * ( zv_a(ji-1,jj+1) + zv_a(ji  ,jj+1) )   & 
     394                     &   + akappa(ji-1,jj,2,1) * ( zu_a(ji-1,jj+1) + zu_a(ji  ,jj+1) ) 
     395                  ze21 =   akappa(ji-1,jj,1,1) * ( zv_a(ji  ,jj+1) - zv_a(ji-1,jj+1) )   & 
     396                     &   - akappa(ji-1,jj,1,2) * ( zu_a(ji  ,jj+1) + zu_a(ji-1,jj+1) ) 
     397                  zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
     398                  zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     399                  zvis12 =       zviseta (ji-1,jj) + dm 
     400                  zvis21 =       zviseta (ji-1,jj) 
     401                  zdiag = zvis22 * ( ze11 + ze22 ) 
     402                  zs11_12 =  zvis11 * ze11 + zdiag 
     403                  zs12_12 =  zvis12 * ze12 + zvis21 * ze21 
     404                  zs22_12 =  zvis11 * ze22 + zdiag 
     405                  zs21_12 =  zvis12 * ze21 + zvis21 * ze12 
     406 
     407                  ze11 =   akappa(ji,jj,1,1) * ( zu_a(ji+1,jj) + zu_a(ji+1,jj+1) - zu_a(ji  ,jj+1) )   & 
     408                     &   + akappa(ji,jj,1,2) * ( zv_a(ji+1,jj) + zv_a(ji+1,jj+1) + zv_a(ji  ,jj+1) ) 
     409                  ze12 = - akappa(ji,jj,2,2) * ( zu_a(ji+1,jj) - zu_a(ji  ,jj+1) - zu_a(ji+1,jj+1) )   & 
     410                     &   - akappa(ji,jj,2,1) * ( zv_a(ji+1,jj) + zv_a(ji  ,jj+1) + zv_a(ji+1,jj+1) ) 
     411                  ze22 = - akappa(ji,jj,2,2) * ( zv_a(ji+1,jj) - zv_a(ji  ,jj+1) - zv_a(ji+1,jj+1) )   & 
     412                     &   + akappa(ji,jj,2,1) * ( zu_a(ji+1,jj) + zu_a(ji  ,jj+1) + zu_a(ji+1,jj+1) ) 
     413                  ze21 =   akappa(ji,jj,1,1) * ( zv_a(ji+1,jj) + zv_a(ji+1,jj+1) - zv_a(ji  ,jj+1) )   & 
     414                     &   - akappa(ji,jj,1,2) * ( zu_a(ji+1,jj) + zu_a(ji+1,jj+1) + zu_a(ji  ,jj+1) ) 
     415                  zvis11 = 2.0 * zviseta (ji,jj) + dm 
     416                  zvis22 =       zviszeta(ji,jj) - zviseta(ji,jj) 
     417                  zvis12 =       zviseta (ji,jj) + dm 
     418                  zvis21 =       zviseta (ji,jj) 
     419                  zdiag = zvis22 * ( ze11 + ze22 ) 
     420                  zs11_22 =  zvis11 * ze11 + zdiag 
     421                  zs12_22 =  zvis12 * ze12 + zvis21 * ze21 
     422                  zs22_22 =  zvis11 * ze22 + zdiag 
     423                  zs21_22 =  zvis12 * ze21 + zvis21 * ze12 
     424 
     425            ! 2nd part 
     426                  ze11 =   akappa(ji-1,jj-1,1,1) * ( zu_a(ji  ,jj-1) - zu_a(ji-1,jj-1) - zu_a(ji-1,jj) )   & 
     427                     &   + akappa(ji-1,jj-1,1,2) * ( zv_a(ji  ,jj-1) + zv_a(ji-1,jj-1) + zv_a(ji-1,jj) ) 
     428                  ze12 = - akappa(ji-1,jj-1,2,2) * ( zu_a(ji-1,jj-1) + zu_a(ji  ,jj-1) - zu_a(ji-1,jj) )   & 
     429                     &   - akappa(ji-1,jj-1,2,1) * ( zv_a(ji-1,jj-1) + zv_a(ji  ,jj-1) + zv_a(ji-1,jj) ) 
     430                  ze22 = - akappa(ji-1,jj-1,2,2) * ( zv_a(ji-1,jj-1) + zv_a(ji  ,jj-1) - zv_a(ji-1,jj) )   & 
     431                     &   + akappa(ji-1,jj-1,2,1) * ( zu_a(ji-1,jj-1) + zu_a(ji  ,jj-1) + zu_a(ji-1,jj) ) 
     432                  ze21 =   akappa(ji-1,jj-1,1,1) * ( zv_a(ji  ,jj-1) - zv_a(ji-1,jj-1) - zv_a(ji-1,jj) )   & 
     433                     &  -  akappa(ji-1,jj-1,1,2) * ( zu_a(ji  ,jj-1) + zu_a(ji-1,jj-1) + zu_a(ji-1,jj) ) 
     434                  zvis11 = 2.0 * zviseta (ji-1,jj-1) + dm 
     435                  zvis22 =       zviszeta(ji-1,jj-1) - zviseta(ji-1,jj-1) 
     436                  zvis12 =       zviseta (ji-1,jj-1) + dm 
     437                  zvis21 =       zviseta (ji-1,jj-1) 
     438                  zdiag = zvis22 * ( ze11 + ze22 ) 
     439                  zs11_11 =  zvis11 * ze11 + zdiag 
     440                  zs12_11 =  zvis12 * ze12 + zvis21 * ze21 
     441                  zs22_11 =  zvis11 * ze22 + zdiag 
     442                  zs21_11 =  zvis12 * ze21 + zvis21 * ze12 
     443 
     444                  ze11 =   akappa(ji,jj-1,1,1) * ( zu_a(ji+1,jj-1) - zu_a(ji  ,jj-1) )   & 
     445                     &   + akappa(ji,jj-1,1,2) * ( zv_a(ji+1,jj-1) + zv_a(ji  ,jj-1) ) 
     446                  ze12 = - akappa(ji,jj-1,2,2) * ( zu_a(ji  ,jj-1) + zu_a(ji+1,jj-1) )   & 
     447                     &   - akappa(ji,jj-1,2,1) * ( zv_a(ji  ,jj-1) + zv_a(ji+1,jj-1) ) 
     448                  ze22 = - akappa(ji,jj-1,2,2) * ( zv_a(ji  ,jj-1) + zv_a(ji+1,jj-1) )   & 
     449                     &   + akappa(ji,jj-1,2,1) * ( zu_a(ji  ,jj-1) + zu_a(ji+1,jj-1) ) 
     450                  ze21 =   akappa(ji,jj-1,1,1) * ( zv_a(ji+1,jj-1) - zv_a(ji  ,jj-1) )   & 
     451                     &   - akappa(ji,jj-1,1,2) * ( zu_a(ji+1,jj-1) + zu_a(ji  ,jj-1) ) 
     452                  zvis11 = 2.0 * zviseta (ji,jj-1) + dm 
     453                  zvis22 =       zviszeta(ji,jj-1) - zviseta(ji,jj-1) 
     454                  zvis12 =       zviseta (ji,jj-1) + dm 
     455                  zvis21 =       zviseta (ji,jj-1) 
     456                  zdiag = zvis22 * ( ze11 + ze22 ) 
     457                  zs11_21 =  zs11_21 + zvis11 * ze11 + zdiag 
     458                  zs12_21 =  zs12_21 + zvis12 * ze12 + zvis21 * ze21 
     459                  zs22_21 =  zs22_21 + zvis11 * ze22 + zdiag 
     460                  zs21_21 =  zs21_21 + zvis12 * ze21 + zvis21 * ze12 
     461 
     462                  ze11 = - akappa(ji-1,jj,1,1) * zu_a(ji-1,jj) + akappa(ji-1,jj,1,2) * zv_a(ji-1,jj) 
     463                  ze12 = - akappa(ji-1,jj,2,2) * zu_a(ji-1,jj) - akappa(ji-1,jj,2,1) * zv_a(ji-1,jj) 
     464                  ze22 = - akappa(ji-1,jj,2,2) * zv_a(ji-1,jj) + akappa(ji-1,jj,2,1) * zu_a(ji-1,jj) 
     465                  ze21 = - akappa(ji-1,jj,1,1) * zv_a(ji-1,jj) - akappa(ji-1,jj,1,2) * zu_a(ji-1,jj) 
     466                  zvis11 = 2.0 * zviseta (ji-1,jj) + dm 
     467                  zvis22 =       zviszeta(ji-1,jj) - zviseta(ji-1,jj) 
     468                  zvis12 =       zviseta (ji-1,jj) + dm 
     469                  zvis21 =       zviseta (ji-1,jj) 
     470                  zdiag = zvis22 * ( ze11 + ze22 ) 
     471                  zs11_12 =  zs11_12 + zvis11 * ze11 + zdiag 
     472                  zs12_12 =  zs12_12 + zvis12 * ze12 + zvis21 * ze21 
     473                  zs22_12 =  zs22_12 + zvis11 * ze22 + zdiag 
     474                  zs21_12 =  zs21_12 + zvis12 * ze21 + zvis21 * ze12 
     475 
     476                  zd1 = + alambd(ji,jj,2,2,2,1) * zs11_21 + alambd(ji,jj,2,2,2,2) * zs11_22  & 
     477                     &  - alambd(ji,jj,2,2,1,1) * zs11_11 - alambd(ji,jj,2,2,1,2) * zs11_12  & 
     478                     &  - alambd(ji,jj,1,1,2,1) * zs12_21 - alambd(ji,jj,1,1,1,1) * zs12_11  & 
     479                     &  + alambd(ji,jj,1,1,2,2) * zs12_22 + alambd(ji,jj,1,1,1,2) * zs12_12  & 
     480                     &  + alambd(ji,jj,1,2,1,1) * zs21_11 + alambd(ji,jj,1,2,2,1) * zs21_21  & 
     481                     &  + alambd(ji,jj,1,2,1,2) * zs21_12 + alambd(ji,jj,1,2,2,2) * zs21_22  & 
     482                     &  - alambd(ji,jj,2,1,1,1) * zs22_11 - alambd(ji,jj,2,1,2,1) * zs22_21  & 
     483                     &  - alambd(ji,jj,2,1,1,2) * zs22_12 - alambd(ji,jj,2,1,2,2) * zs22_22 
     484 
     485                  zd2 = + alambd(ji,jj,2,2,2,1) * zs21_21 + alambd(ji,jj,2,2,2,2) * zs21_22  & 
     486                     &  - alambd(ji,jj,2,2,1,1) * zs21_11 - alambd(ji,jj,2,2,1,2) * zs21_12  & 
     487                     &  - alambd(ji,jj,1,1,2,1) * zs22_21 - alambd(ji,jj,1,1,1,1) * zs22_11  & 
     488                     &  + alambd(ji,jj,1,1,2,2) * zs22_22 + alambd(ji,jj,1,1,1,2) * zs22_12  & 
     489                     &  - alambd(ji,jj,1,2,1,1) * zs11_11 - alambd(ji,jj,1,2,2,1) * zs11_21  & 
     490                     &  - alambd(ji,jj,1,2,1,2) * zs11_12 - alambd(ji,jj,1,2,2,2) * zs11_22  & 
     491                     &  + alambd(ji,jj,2,1,1,1) * zs12_11 + alambd(ji,jj,2,1,2,1) * zs12_21  & 
     492                     &  + alambd(ji,jj,2,1,1,2) * zs12_12 + alambd(ji,jj,2,1,2,2) * zs12_22 
     493 
     494                  zur     = zu_a(ji,jj) - ui_oce(ji,jj) 
     495                  zvr     = zv_a(ji,jj) - vi_oce(ji,jj) 
     496!!!! 
     497                  zmod    = SQRT( zur*zur + zvr*zvr ) * ( 1.0 - zfrld(ji,jj) ) 
     498                  za      = rhoco * zmod 
     499!!!! 
     500!!gm chg resul    za      = rhoco * SQRT( zur*zur + zvr*zvr ) * ( 1.0 - zfrld(ji,jj) ) 
     501                  zac     = za * cangvg 
     502                  zmpzas  = alpha * zcorl(ji,jj) + za * zsang(ji,jj) 
    411503                  zmassdt = zusdtp * zmass(ji,jj) 
    412504                  zcorlal = ( 1.0 - alpha ) * zcorl(ji,jj) 
    413505 
    414                   za1(ji,jj) =  zmassdt * zu0(ji,jj) + zcorlal * zv0(ji,jj) + za1ct(ji,jj)   & 
    415                      &        + za * ( cangvg * u_oce(ji,jj) - zsang * v_oce(ji,jj) ) 
    416  
    417                   za2(ji,jj) =  zmassdt * zv0(ji,jj) - zcorlal * zu0(ji,jj) + za2ct(ji,jj)   & 
    418                      &        + za * ( cangvg * v_oce(ji,jj) + zsang * u_oce(ji,jj) ) 
    419  
    420                   zb1(ji,jj)  = zmassdt + zac - zc1u(ji,jj) 
    421                   zb2(ji,jj)  = zmpzas  - zc2u(ji,jj) 
    422                   zc1(ji,jj)  = zmpzas  + zc1v(ji,jj) 
    423                   zc2(ji,jj)  = zmassdt + zac  - zc2v(ji,jj)  
    424                   zdeter      = zc1(ji,jj) * zb2(ji,jj) + zc2(ji,jj) * zb1(ji,jj) 
    425                   zden(ji,jj) = SIGN( rone, zdeter) / MAX( epsd , ABS( zdeter ) ) 
     506                  za1 =  zmassdt * zu0(ji,jj) + zcorlal * zv0(ji,jj) + za1ct(ji,jj)   & 
     507                     &        + za * ( cangvg * ui_oce(ji,jj) - zsang(ji,jj) * vi_oce(ji,jj) ) 
     508                  za2 =  zmassdt * zv0(ji,jj) - zcorlal * zu0(ji,jj) + za2ct(ji,jj)   & 
     509                     &        + za * ( cangvg * vi_oce(ji,jj) + zsang(ji,jj) * ui_oce(ji,jj) ) 
     510                  zb1    = zmassdt + zac - zc1u(ji,jj) 
     511                  zb2    = zmpzas        - zc2u(ji,jj) 
     512                  zc1    = zmpzas        + zc1v(ji,jj) 
     513                  zc2    = zmassdt + zac - zc2v(ji,jj) 
     514                  zdeter = zc1 * zb2 + zc2 * zb1 
     515                  zden   = SIGN( rone, zdeter) / MAX( epsd , ABS( zdeter ) ) 
     516                  zunw   = (  ( za1 + zd1 ) * zc2 + ( za2 + zd2 ) * zc1 ) * zden 
     517                  zvnw   = (  ( za2 + zd2 ) * zb1 - ( za1 + zd1 ) * zb2 ) * zden 
     518                  zmask  = ( 1.0 - MAX( rzero, SIGN( rone , 1.0 - zmass(ji,jj) ) ) ) * tmu(ji,jj) 
     519 
     520                  zu_n(ji,jj) = ( zu_a(ji,jj) + om * ( zunw - zu_a(ji,jj) ) * tmu(ji,jj) ) * zmask 
     521                  zv_n(ji,jj) = ( zv_a(ji,jj) + om * ( zvnw - zv_a(ji,jj) ) * tmu(ji,jj) ) * zmask 
    426522               END DO 
    427523            END DO 
    428524 
    429             ! The computation of ice interaction term is splitted into two parts 
    430             !------------------------------------------------------------------------- 
    431  
    432             ! Terms that do not involve already up-dated velocities. 
    433           
    434             DO jj = k_j1+1, k_jpj-1 
    435                DO ji = 2, jpim1 
    436                   iim1 = ji 
    437                   ijm1 = jj - 1 
    438                   iip1 = ji + 1 
    439                   ijp1 = jj 
    440                   ze11 =   akappa(iim1,ijm1,1,1) * u_ice(iip1,ijp1) + akappa(iim1,ijm1,1,2) * v_ice(iip1,ijp1) 
    441                   ze12 = + akappa(iim1,ijm1,2,2) * u_ice(iip1,ijp1) - akappa(iim1,ijm1,2,1) * v_ice(iip1,ijp1) 
    442                   ze22 = + akappa(iim1,ijm1,2,2) * v_ice(iip1,ijp1) + akappa(iim1,ijm1,2,1) * u_ice(iip1,ijp1) 
    443                   ze21 =   akappa(iim1,ijm1,1,1) * v_ice(iip1,ijp1) - akappa(iim1,ijm1,1,2) * u_ice(iip1,ijp1) 
    444                   zvis11 = 2.0 * zviseta (iim1,ijm1) + dm 
    445                   zvis22 =       zviszeta(iim1,ijm1) - zviseta(iim1,ijm1) 
    446                   zvis12 =       zviseta (iim1,ijm1) + dm 
    447                   zvis21 =       zviseta (iim1,ijm1) 
    448                   zdiag = zvis22 * ( ze11 + ze22 ) 
    449                   zs11(ji,jj,2,1) =  zvis11 * ze11 + zdiag 
    450                   zs12(ji,jj,2,1) =  zvis12 * ze12 + zvis21 * ze21 
    451                   zs22(ji,jj,2,1) =  zvis11 * ze22 + zdiag 
    452                   zs21(ji,jj,2,1) =  zvis12 * ze21 + zvis21 * ze12 
    453  
    454  
    455                   iim1 = ji - 1 
    456                   ijm1 = jj 
    457                   iip1 = ji 
    458                   ijp1 = jj + 1                    
    459                   ze11 =   akappa(iim1,ijm1,1,1) * ( u_ice(iip1,ijp1) - u_ice(iim1,ijp1) )   & 
    460                      &   + akappa(iim1,ijm1,1,2) * ( v_ice(iip1,ijp1) + v_ice(iim1,ijp1) ) 
    461                   ze12 = + akappa(iim1,ijm1,2,2) * ( u_ice(iim1,ijp1) + u_ice(iip1,ijp1) )   & 
    462                      &   - akappa(iim1,ijm1,2,1) * ( v_ice(iim1,ijp1) + v_ice(iip1,ijp1) ) 
    463                   ze22 = + akappa(iim1,ijm1,2,2) * ( v_ice(iim1,ijp1) + v_ice(iip1,ijp1) )   & 
    464                      &   + akappa(iim1,ijm1,2,1) * ( u_ice(iim1,ijp1) + u_ice(iip1,ijp1) ) 
    465                   ze21 =   akappa(iim1,ijm1,1,1) * ( v_ice(iip1,ijp1) - v_ice(iim1,ijp1) )   & 
    466                      &   - akappa(iim1,ijm1,1,2) * ( u_ice(iip1,ijp1) + u_ice(iim1,ijp1) ) 
    467                   zvis11 = 2.0 * zviseta (iim1,ijm1) + dm 
    468                   zvis22 =       zviszeta(iim1,ijm1) - zviseta(iim1,ijm1) 
    469                   zvis12 =       zviseta (iim1,ijm1) + dm 
    470                   zvis21 =       zviseta (iim1,ijm1) 
    471                   zdiag = zvis22 * ( ze11 + ze22 ) 
    472                   zs11(ji,jj,1,2) =  zvis11 * ze11 + zdiag 
    473                   zs12(ji,jj,1,2) =  zvis12 * ze12 + zvis21 * ze21 
    474                   zs22(ji,jj,1,2) =  zvis11 * ze22 + zdiag 
    475                   zs21(ji,jj,1,2) =  zvis12 * ze21 + zvis21 * ze12 
    476  
    477                   iim1 = ji 
    478                   ijm1 = jj 
    479                   iip1 = ji + 1 
    480                   ijp1 = jj + 1 
    481                   ze11 =   akappa(iim1,ijm1,1,1) * ( u_ice(iip1,ijm1) + u_ice(iip1,ijp1) - u_ice(iim1,ijp1) )   & 
    482                      &   + akappa(iim1,ijm1,1,2) * ( v_ice(iip1,ijm1) + v_ice(iip1,ijp1) + v_ice(iim1,ijp1) ) 
    483                   ze12 = - akappa(iim1,ijm1,2,2) * ( u_ice(iip1,ijm1) - u_ice(iim1,ijp1) - u_ice(iip1,ijp1) )   & 
    484                      &   - akappa(iim1,ijm1,2,1) * ( v_ice(iip1,ijm1) + v_ice(iim1,ijp1) + v_ice(iip1,ijp1) ) 
    485                   ze22 = - akappa(iim1,ijm1,2,2) * ( v_ice(iip1,ijm1) - v_ice(iim1,ijp1) - v_ice(iip1,ijp1) )   & 
    486                      &   + akappa(iim1,ijm1,2,1) * ( u_ice(iip1,ijm1) + u_ice(iim1,ijp1) + u_ice(iip1,ijp1) ) 
    487                   ze21 =   akappa(iim1,ijm1,1,1) * ( v_ice(iip1,ijm1) + v_ice(iip1,ijp1) - v_ice(iim1,ijp1) )   & 
    488                      &   - akappa(iim1,ijm1,1,2) * ( u_ice(iip1,ijm1) + u_ice(iip1,ijp1) + u_ice(iim1,ijp1) )  
    489                   zvis11 = 2.0 * zviseta (iim1,ijm1) + dm 
    490                   zvis22 =       zviszeta(iim1,ijm1) - zviseta(iim1,ijm1) 
    491                   zvis12 =       zviseta (iim1,ijm1) + dm 
    492                   zvis21 =       zviseta (iim1,ijm1) 
    493  
    494                   zdiag = zvis22 * ( ze11 + ze22 ) 
    495                   zs11(ji,jj,2,2) =  zvis11 * ze11 + zdiag 
    496                   zs12(ji,jj,2,2) =  zvis12 * ze12 + zvis21 * ze21 
    497                   zs22(ji,jj,2,2) =  zvis11 * ze22 + zdiag 
    498                   zs21(ji,jj,2,2) =  zvis12 * ze21 + zvis21 * ze12 
    499  
    500                END DO 
     525            CALL lbc_lnk( zu_n(:,1:jpj), 'I', -1. ) 
     526            CALL lbc_lnk( zv_n(:,1:jpj), 'I', -1. ) 
     527 
     528            ! Test of Convergence 
     529            DO jj = k_j1+1 , k_jpj-1 
     530               zresr(:,jj) = MAX( ABS( zu_a(:,jj) - zu_n(:,jj) ) , ABS( zv_a(:,jj) - zv_n(:,jj) ) ) 
    501531            END DO 
    502  
    503             ! Terms involving already up-dated velocities. 
    504             !-Using the arrays zu_ice and zv_ice in the computation of the terms ze leads to JACOBI's method;  
    505             ! Using arrays u and v in the computation of the terms ze leads to GAUSS-SEIDEL method. 
    506               
    507             DO jj = k_j1+1, k_jpj-1 
    508                DO ji = 2, jpim1 
    509                   iim1 = ji - 1 
    510                   ijm1 = jj - 1 
    511                   iip1 = ji 
    512                   ijp1 = jj 
    513                   ze11 =   akappa(iim1,ijm1,1,1) * ( zu_ice(iip1,ijm1) - zu_ice(iim1,ijm1) - zu_ice(iim1,ijp1) )   & 
    514                      &   + akappa(iim1,ijm1,1,2) * ( zv_ice(iip1,ijm1) + zv_ice(iim1,ijm1) + zv_ice(iim1,ijp1) ) 
    515                   ze12 = - akappa(iim1,ijm1,2,2) * ( zu_ice(iim1,ijm1) + zu_ice(iip1,ijm1) - zu_ice(iim1,ijp1) )   & 
    516                      &   - akappa(iim1,ijm1,2,1) * ( zv_ice(iim1,ijm1) + zv_ice(iip1,ijm1) + zv_ice(iim1,ijp1) ) 
    517                   ze22 = - akappa(iim1,ijm1,2,2) * ( zv_ice(iim1,ijm1) + zv_ice(iip1,ijm1) - zv_ice(iim1,ijp1) )   & 
    518                      &   + akappa(iim1,ijm1,2,1) * ( zu_ice(iim1,ijm1) + zu_ice(iip1,ijm1) + zu_ice(iim1,ijp1) ) 
    519                   ze21 =   akappa(iim1,ijm1,1,1) * ( zv_ice(iip1,ijm1) - zv_ice(iim1,ijm1) - zv_ice(iim1,ijp1) )   & 
    520                      &  -  akappa(iim1,ijm1,1,2) * ( zu_ice(iip1,ijm1) + zu_ice(iim1,ijm1) + zu_ice(iim1,ijp1) ) 
    521                   zvis11 = 2.0 * zviseta (iim1,ijm1) + dm 
    522                   zvis22 =       zviszeta(iim1,ijm1) - zviseta(iim1,ijm1) 
    523                   zvis12 =       zviseta (iim1,ijm1) + dm 
    524                   zvis21 =       zviseta (iim1,ijm1) 
    525  
    526                   zdiag = zvis22 * ( ze11 + ze22 ) 
    527                   zs11(ji,jj,1,1) =  zvis11 * ze11 + zdiag 
    528                   zs12(ji,jj,1,1) =  zvis12 * ze12 + zvis21 * ze21 
    529                   zs22(ji,jj,1,1) =  zvis11 * ze22 + zdiag 
    530                   zs21(ji,jj,1,1) =  zvis12 * ze21 + zvis21 * ze12 
    531  
    532 #if defined key_agrif 
    533              END DO 
    534           END DO 
    535  
    536           DO jj = k_j1+1, k_jpj-1 
    537              DO ji = 2, jpim1 
    538 #endif 
    539  
    540                   iim1 = ji 
    541                   ijm1 = jj - 1 
    542                   iip1 = ji + 1 
    543                   ze11 =   akappa(iim1,ijm1,1,1) * ( zu_ice(iip1,ijm1) - zu_ice(iim1,ijm1) )   & 
    544                      &   + akappa(iim1,ijm1,1,2) * ( zv_ice(iip1,ijm1) + zv_ice(iim1,ijm1) ) 
    545                   ze12 = - akappa(iim1,ijm1,2,2) * ( zu_ice(iim1,ijm1) + zu_ice(iip1,ijm1) )   & 
    546                      &   - akappa(iim1,ijm1,2,1) * ( zv_ice(iim1,ijm1) + zv_ice(iip1,ijm1) ) 
    547                   ze22 = - akappa(iim1,ijm1,2,2) * ( zv_ice(iim1,ijm1) + zv_ice(iip1,ijm1) )   & 
    548                      &   + akappa(iim1,ijm1,2,1) * ( zu_ice(iim1,ijm1) + zu_ice(iip1,ijm1) ) 
    549                   ze21 =   akappa(iim1,ijm1,1,1) * ( zv_ice(iip1,ijm1) - zv_ice(iim1,ijm1) )   & 
    550                      &   - akappa(iim1,ijm1,1,2) * ( zu_ice(iip1,ijm1) + zu_ice(iim1,ijm1) ) 
    551                   zvis11 = 2.0 * zviseta (iim1,ijm1) + dm 
    552                   zvis22 =       zviszeta(iim1,ijm1) - zviseta(iim1,ijm1) 
    553                   zvis12 =       zviseta (iim1,ijm1) + dm 
    554                   zvis21 =       zviseta (iim1,ijm1) 
    555  
    556                   zdiag = zvis22 * ( ze11 + ze22 ) 
    557                   zs11(ji,jj,2,1) =  zs11(ji,jj,2,1) + zvis11 * ze11 + zdiag 
    558                   zs12(ji,jj,2,1) =  zs12(ji,jj,2,1) + zvis12 * ze12 + zvis21 * ze21 
    559                   zs22(ji,jj,2,1) =  zs22(ji,jj,2,1) + zvis11 * ze22 + zdiag 
    560                   zs21(ji,jj,2,1) =  zs21(ji,jj,2,1) + zvis12 * ze21 + zvis21 * ze12 
    561  
    562  
    563                   iim1 = ji - 1 
    564                   ijm1 = jj  
    565                   ze11 = - akappa(iim1,ijm1,1,1) * zu_ice(iim1,ijm1) + akappa(iim1,ijm1,1,2) * zv_ice(iim1,ijm1) 
    566                   ze12 = - akappa(iim1,ijm1,2,2) * zu_ice(iim1,ijm1) - akappa(iim1,ijm1,2,1) * zv_ice(iim1,ijm1) 
    567                   ze22 = - akappa(iim1,ijm1,2,2) * zv_ice(iim1,ijm1) + akappa(iim1,ijm1,2,1) * zu_ice(iim1,ijm1) 
    568                   ze21 = - akappa(iim1,ijm1,1,1) * zv_ice(iim1,ijm1) - akappa(iim1,ijm1,1,2) * zu_ice(iim1,ijm1) 
    569                   zvis11 = 2.0 * zviseta (iim1,ijm1) + dm 
    570                   zvis22 =       zviszeta(iim1,ijm1) - zviseta(iim1,ijm1) 
    571                   zvis12 =       zviseta (iim1,ijm1) + dm 
    572                   zvis21 =       zviseta (iim1,ijm1) 
    573  
    574                   zdiag = zvis22 * ( ze11 + ze22 ) 
    575                   zs11(ji,jj,1,2) =  zs11(ji,jj,1,2) + zvis11 * ze11 + zdiag  
    576                   zs12(ji,jj,1,2) =  zs12(ji,jj,1,2) + zvis12 * ze12 + zvis21 * ze21 
    577                   zs22(ji,jj,1,2) =  zs22(ji,jj,1,2) + zvis11 * ze22 + zdiag 
    578                   zs21(ji,jj,1,2) =  zs21(ji,jj,1,2) + zvis12 * ze21 + zvis21 * ze12 
    579  
    580 #if defined key_agrif 
    581              END DO 
    582           END DO 
    583  
    584           DO jj = k_j1+1, k_jpj-1 
    585              DO ji = 2, jpim1 
    586 #endif 
    587                   zd1(ji,jj) =   & 
    588                      + alambd(ji,jj,2,2,2,1) * zs11(ji,jj,2,1) + alambd(ji,jj,2,2,2,2) * zs11(ji,jj,2,2)  & 
    589                      - alambd(ji,jj,2,2,1,1) * zs11(ji,jj,1,1) - alambd(ji,jj,2,2,1,2) * zs11(ji,jj,1,2)  & 
    590                      - alambd(ji,jj,1,1,2,1) * zs12(ji,jj,2,1) - alambd(ji,jj,1,1,1,1) * zs12(ji,jj,1,1)  & 
    591                      + alambd(ji,jj,1,1,2,2) * zs12(ji,jj,2,2) + alambd(ji,jj,1,1,1,2) * zs12(ji,jj,1,2)  & 
    592                      + alambd(ji,jj,1,2,1,1) * zs21(ji,jj,1,1) + alambd(ji,jj,1,2,2,1) * zs21(ji,jj,2,1)  & 
    593                      + alambd(ji,jj,1,2,1,2) * zs21(ji,jj,1,2) + alambd(ji,jj,1,2,2,2) * zs21(ji,jj,2,2)  & 
    594                      - alambd(ji,jj,2,1,1,1) * zs22(ji,jj,1,1) - alambd(ji,jj,2,1,2,1) * zs22(ji,jj,2,1)  & 
    595                      - alambd(ji,jj,2,1,1,2) * zs22(ji,jj,1,2) - alambd(ji,jj,2,1,2,2) * zs22(ji,jj,2,2) 
    596                   zd2(ji,jj) =   & 
    597                      + alambd(ji,jj,2,2,2,1) * zs21(ji,jj,2,1) + alambd(ji,jj,2,2,2,2) * zs21(ji,jj,2,2)  & 
    598                      - alambd(ji,jj,2,2,1,1) * zs21(ji,jj,1,1) - alambd(ji,jj,2,2,1,2) * zs21(ji,jj,1,2)  & 
    599                      - alambd(ji,jj,1,1,2,1) * zs22(ji,jj,2,1) - alambd(ji,jj,1,1,1,1) * zs22(ji,jj,1,1)  & 
    600                      + alambd(ji,jj,1,1,2,2) * zs22(ji,jj,2,2) + alambd(ji,jj,1,1,1,2) * zs22(ji,jj,1,2)  & 
    601                      - alambd(ji,jj,1,2,1,1) * zs11(ji,jj,1,1) - alambd(ji,jj,1,2,2,1) * zs11(ji,jj,2,1)  & 
    602                      - alambd(ji,jj,1,2,1,2) * zs11(ji,jj,1,2) - alambd(ji,jj,1,2,2,2) * zs11(ji,jj,2,2)  & 
    603                      + alambd(ji,jj,2,1,1,1) * zs12(ji,jj,1,1) + alambd(ji,jj,2,1,2,1) * zs12(ji,jj,2,1)  & 
    604                      + alambd(ji,jj,2,1,1,2) * zs12(ji,jj,1,2) + alambd(ji,jj,2,1,2,2) * zs12(ji,jj,2,2) 
    605                END DO 
     532            zresm = MAXVAL( zresr(1:jpi,k_j1+1:k_jpj-1) ) 
     533!!!! this should be faster on scalar processor 
     534!           zresm = MAXVAL(  MAX( ABS( zu_a(1:jpi,k_j1+1:k_jpj-1) - zu_n(1:jpi,k_j1+1:k_jpj-1) ),   & 
     535!              &                  ABS( zv_a(1:jpi,k_j1+1:k_jpj-1) - zv_n(1:jpi,k_j1+1:k_jpj-1) ) )  ) 
     536!!!! 
     537            IF( lk_mpp )   CALL mpp_max( zresm )   ! max over the global domain 
     538 
     539            DO jj = k_j1, k_jpj 
     540               zu_a(:,jj) = zu_n(:,jj) 
     541               zv_a(:,jj) = zv_n(:,jj) 
    606542            END DO 
    607543 
    608             DO jj = k_j1+1, k_jpj-1 
    609                DO ji = 2, jpim1 
    610                   zunw = (  ( za1(ji,jj) + zd1(ji,jj) ) * zc2(ji,jj)        & 
    611                      &    + ( za2(ji,jj) + zd2(ji,jj) ) * zc1(ji,jj) ) * zden(ji,jj) 
    612  
    613                   zvnw = (  ( za2(ji,jj) + zd2(ji,jj) ) * zb1(ji,jj)        & 
    614                      &    - ( za1(ji,jj) + zd1(ji,jj) ) * zb2(ji,jj) ) * zden(ji,jj) 
    615  
    616                   zmask = ( 1.0 - MAX( rzero, SIGN( rone , 1.0 - zmass(ji,jj) ) ) ) * tmu(ji,jj) 
    617  
    618                   u_ice(ji,jj) = ( u_ice(ji,jj) + om * ( zunw - u_ice(ji,jj) ) * tmu(ji,jj) ) * zmask 
    619                   v_ice(ji,jj) = ( v_ice(ji,jj) + om * ( zvnw - v_ice(ji,jj) ) * tmu(ji,jj) ) * zmask 
    620                END DO 
     544            IF( zresm <= resl )   EXIT   iflag 
     545 
     546            !                                                   ! ================ ! 
     547         END DO    iflag                                        !  end Relaxation  ! 
     548         !                                                      ! ================ ! 
     549 
     550         IF( zindu == 0 ) THEN      ! even iteration 
     551            DO jj = k_j1 , k_jpj-1 
     552               zu0(:,jj) = zu_a(:,jj) 
     553               zv0(:,jj) = zv_a(:,jj) 
    621554            END DO 
    622  
    623             CALL lbc_lnk( u_ice, 'I', -1. ) 
    624             CALL lbc_lnk( v_ice, 'I', -1. ) 
    625  
    626             !---  5.2.5.4. Convergence test. 
    627             DO jj = k_j1+1 , k_jpj-1 
    628                zresr(:,jj) = MAX( ABS( u_ice(:,jj) - zu_ice(:,jj) ) , ABS( v_ice(:,jj) - zv_ice(:,jj) ) ) 
    629             END DO 
    630             zresm = MAXVAL( zresr( 1:jpi , k_j1+1:k_jpj-1 ) ) 
    631             IF( lk_mpp )   CALL mpp_max( zresm )   ! max over the global domain 
    632  
    633             IF ( zresm <= resl) EXIT iflag 
    634  
    635          END DO iflag 
    636  
    637          zindu1 = 1.0 - zindu 
    638          DO jj = k_j1 , k_jpj-1 
    639             zu0(:,jj) = zindu * zu0(:,jj) + zindu1 * u_ice(:,jj) 
    640             zv0(:,jj) = zindu * zv0(:,jj) + zindu1 * v_ice(:,jj) 
    641          END DO 
    642       !                                                   ! ==================== ! 
     555         ENDIF 
     556         !                                                ! ==================== ! 
    643557      END DO                                              !  end loop over iter  ! 
    644558      !                                                   ! ==================== ! 
     559 
     560      ui_ice(:,:) = zu_a(:,1:jpj) 
     561      vi_ice(:,:) = zv_a(:,1:jpj) 
    645562 
    646563      IF(ln_ctl) THEN 
    647564         WRITE(charout,FMT="('lim_rhg  : res =',D23.16, ' iter =',I4)") zresm, jter 
    648565         CALL prt_ctl_info(charout) 
    649          CALL prt_ctl(tab2d_1=u_ice, clinfo1=' lim_rhg  : u_ice :', tab2d_2=v_ice, clinfo2=' v_ice :') 
     566         CALL prt_ctl(tab2d_1=ui_ice, clinfo1=' lim_rhg  : ui_ice :', tab2d_2=vi_ice, clinfo2=' vi_ice :') 
    650567      ENDIF 
    651568 
  • trunk/NEMO/LIM_SRC_2/limrst_2.F90

    r823 r888  
    1717   !!---------------------------------------------------------------------- 
    1818   USE ice_2 
    19    USE dom_oce 
    20    USE ice_oce         ! ice variables 
     19   USE sbc_oce 
     20   USE sbc_ice 
    2121   USE daymod 
    2222 
     
    2727   PRIVATE 
    2828 
    29    PUBLIC   lim_rst_opn_2     ! routine called by ??? module 
    30    PUBLIC   lim_rst_write_2   ! routine called by ??? module 
    31    PUBLIC   lim_rst_read_2    ! routine called by ??? module 
     29   PUBLIC   lim_rst_opn_2     ! routine called by sbcice_lim_2.F90 
     30   PUBLIC   lim_rst_write_2   ! routine called by sbcice_lim_2.F90 
     31   PUBLIC   lim_rst_read_2    ! routine called by iceini_2.F90 
    3232 
    3333   LOGICAL, PUBLIC ::   lrst_ice         !: logical to control the ice restart write  
     
    3636   !!---------------------------------------------------------------------- 
    3737   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
    38    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limrst.F90,v 1.15 2007/06/29 14:54:06 opalod Exp $  
     38   !! $ Id: $ 
    3939   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4040   !!---------------------------------------------------------------------- 
     
    5757       
    5858      ! to get better performances with NetCDF format: 
    59       ! we open and define the ice restart file one ice time step before writing the data (-> at nitrst - 2*nfice + 1) 
    60       ! except if we write ice restart files every ice time step or if an ice restart file was writen at nitend - 2*nfice + 1 
    61       IF( kt == nitrst - 2*nfice + 1 .OR. nstock == nfice .OR. ( kt == nitend - nfice + 1 .AND. .NOT. lrst_ice ) ) THEN 
     59      ! we open and define the ice restart file one ice time step before writing the data (-> at nitrst - 2*nn_fsbc + 1) 
     60      ! except if we write ice restart files every ice time step or if an ice restart file was writen at nitend - 2*nn_fsbc + 1 
     61      IF( kt == nitrst - 2*nn_fsbc + 1 .OR. nstock == nn_fsbc .OR. ( kt == nitend - nn_fsbc + 1 .AND. .NOT. lrst_ice ) ) THEN 
    6262         ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
    6363         IF( nitrst > 99999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     
    7272            CASE DEFAULT         ;   WRITE(numout,*) '             open ice restart NetCDF file: '//clname 
    7373            END SELECT 
    74             IF( kt == nitrst - 2*nfice + 1 ) THEN    
    75                WRITE(numout,*)         '             kt = nitrst - 2*nfice + 1 = ', kt,' date= ', ndastp 
    76             ELSE   ;   WRITE(numout,*) '             kt = '                       , kt,' date= ', ndastp 
     74            IF( kt == nitrst - 2*nn_fsbc + 1 ) THEN    
     75               WRITE(numout,*)         '             kt = nitrst - 2*nn_fsbc + 1 = ', kt,' date= ', ndastp 
     76            ELSE   ;   WRITE(numout,*) '             kt = '                         , kt,' date= ', ndastp 
    7777            ENDIF 
    7878         ENDIF 
     
    9090      !! ** purpose  :   output of sea-ice variable in a netcdf file 
    9191      !!---------------------------------------------------------------------- 
    92       INTEGER, INTENT(in) ::   kt     ! number of iteration 
    93       !! 
    94       INTEGER ::   iter     ! kt + nfice -1 
    95       !!---------------------------------------------------------------------- 
    96  
    97       iter = kt + nfice - 1   ! ice restarts are written at kt == nitrst - nfice + 1 
     92      INTEGER, INTENT(in) ::   kt   ! number of iteration 
     93      ! 
     94      INTEGER ::   iter   ! kt + nn_fsbc -1 
     95      !!---------------------------------------------------------------------- 
     96 
     97      iter = kt + nn_fsbc - 1   ! ice restarts are written at kt == nitrst - nn_fsbc + 1 
    9898 
    9999      IF( iter == nitrst ) THEN 
    100100         IF(lwp) WRITE(numout,*) 
    101101         IF(lwp) WRITE(numout,*) 'lim_rst_write_2 : write ice restart file  kt =', kt 
    102          IF(lwp) WRITE(numout,*) '~~~~~~~'          
     102         IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~~~~' 
    103103      ENDIF 
    104104 
     
    106106      ! ------------------  
    107107      !                                                                     ! calendar control 
    108       CALL iom_rstput( iter, nitrst, numriw, 'nfice' , REAL( nfice, wp) )      ! time-step  
    109       CALL iom_rstput( iter, nitrst, numriw, 'kt_ice', REAL( iter , wp) )      ! date 
     108      CALL iom_rstput( iter, nitrst, numriw, 'kt_ice', REAL( iter, wp) )  
    110109       
    111110      CALL iom_rstput( iter, nitrst, numriw, 'hicif' , hicif (:,:)   )      ! prognostic variables  
     
    119118      CALL iom_rstput( iter, nitrst, numriw, 'tbif2' , tbif  (:,:,2) ) 
    120119      CALL iom_rstput( iter, nitrst, numriw, 'tbif3' , tbif  (:,:,3) ) 
    121       CALL iom_rstput( iter, nitrst, numriw, 'u_ice' , u_ice (:,:)   ) 
    122       CALL iom_rstput( iter, nitrst, numriw, 'v_ice' , v_ice (:,:)   ) 
    123       CALL iom_rstput( iter, nitrst, numriw, 'gtaux' , gtaux (:,:)   ) 
    124       CALL iom_rstput( iter, nitrst, numriw, 'gtauy' , gtauy (:,:)   ) 
     120      CALL iom_rstput( iter, nitrst, numriw, 'ui_ice', ui_ice(:,:)   ) 
     121      CALL iom_rstput( iter, nitrst, numriw, 'vi_ice', vi_ice(:,:)   ) 
    125122      CALL iom_rstput( iter, nitrst, numriw, 'qstoif', qstoif(:,:)   ) 
    126123      CALL iom_rstput( iter, nitrst, numriw, 'fsbbq' , fsbbq (:,:)   ) 
     
    175172      !! ** purpose  :   read of sea-ice variable restart in a netcdf file 
    176173      !!---------------------------------------------------------------------- 
    177       ! 
    178       REAL(wp) ::   zfice, ziter 
     174      REAL(wp) ::   ziter 
    179175      !!---------------------------------------------------------------------- 
    180176 
     
    182178         WRITE(numout,*) 
    183179         WRITE(numout,*) 'lim_rst_read_2 : read ice NetCDF restart file' 
    184          WRITE(numout,*) '~~~~~~~~' 
     180         WRITE(numout,*) '~~~~~~~~~~~~~~' 
    185181      ENDIF 
    186182 
    187183      CALL iom_open ( 'restart_ice_in', numrir, kiolib = jprstlib ) 
    188184 
    189       CALL iom_get( numrir, 'nfice' , zfice ) 
    190       CALL iom_get( numrir, 'kt_ice', ziter )     
    191       IF(lwp) WRITE(numout,*) '   read ice restart file at time step    : ', ziter 
     185      CALL iom_get( numrir, 'kt_ice' , ziter ) 
     186      IF(lwp) WRITE(numout,*) '   read ice restart file at time step    : ', INT( ziter ) 
    192187      IF(lwp) WRITE(numout,*) '   in any case we force it to nit000 - 1 : ', nit000 - 1 
    193188 
     
    196191      IF( ( nit000 - INT(ziter) ) /= 1 .AND. ABS( nrstdt ) == 1 )   & 
    197192         &     CALL ctl_stop( 'lim_rst_read ===>>>> : problem with nit000 in ice restart',  & 
    198          &                   '   verify the file or rerun with the value 0 for the',        & 
    199          &                   '   control of time parameter  nrstdt' ) 
    200       IF( INT(zfice) /= nfice          .AND. ABS( nrstdt ) == 1 )   & 
    201          &     CALL ctl_stop( 'lim_rst_read ===>>>> : problem with nfice in ice restart',  & 
    202193         &                   '   verify the file or rerun with the value 0 for the',        & 
    203194         &                   '   control of time parameter  nrstdt' ) 
     
    213204      CALL iom_get( numrir, jpdom_autoglo, 'tbif2' , tbif(:,:,2) )     
    214205      CALL iom_get( numrir, jpdom_autoglo, 'tbif3' , tbif(:,:,3) )     
    215       CALL iom_get( numrir, jpdom_autoglo, 'u_ice' , u_ice  )     
    216       CALL iom_get( numrir, jpdom_autoglo, 'v_ice' , v_ice  )     
    217       CALL iom_get( numrir, jpdom_autoglo, 'gtaux' , gtaux  )     
    218       CALL iom_get( numrir, jpdom_autoglo, 'gtauy' , gtauy  )     
     206      CALL iom_get( numrir, jpdom_autoglo, 'ui_ice', ui_ice )     
     207      CALL iom_get( numrir, jpdom_autoglo, 'vi_ice', vi_ice )     
    219208      CALL iom_get( numrir, jpdom_autoglo, 'qstoif', qstoif )     
    220209      CALL iom_get( numrir, jpdom_autoglo, 'fsbbq' , fsbbq  )     
  • trunk/NEMO/LIM_SRC_2/limtab_2.F90

    r823 r888  
    2121   !!---------------------------------------------------------------------- 
    2222   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    23    !! $Header$  
    24    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     23   !! $ Id: $ 
     24   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    2525   !!---------------------------------------------------------------------- 
    2626CONTAINS 
  • trunk/NEMO/LIM_SRC_2/limthd_2.F90

    r823 r888  
    44   !!              LIM thermo ice model : ice thermodynamic 
    55   !!====================================================================== 
     6   !! History :  1.0  !  00-01 (LIM) 
     7   !!            2.0  !  02-07 (C. Ethe, G. Madec) F90 
     8   !!            2.0  !  03-08 (C. Ethe)  add lim_thd_init 
     9   !!--------------------------------------------------------------------- 
    610#if defined key_lim2 
    711   !!---------------------------------------------------------------------- 
     
    1822   USE ice_2           ! LIM sea-ice variables 
    1923   USE ice_oce         ! sea-ice/ocean variables 
    20    USE flx_oce         ! sea-ice/ocean forcings variables  
     24   USE sbc_oce         !  
     25   USE sbc_ice         !  
    2126   USE thd_ice_2       ! LIM thermodynamic sea-ice variables 
    2227   USE dom_ice_2       ! LIM sea-ice domain 
     
    3035   PRIVATE 
    3136 
    32    !! * Routine accessibility 
    33    PUBLIC lim_thd_2       ! called by lim_step_2 
    34  
    35    !! * Module variables 
    36    REAL(wp)  ::            &  ! constant values 
    37       epsi20 = 1.e-20   ,  & 
    38       epsi16 = 1.e-16   ,  & 
    39       epsi04 = 1.e-04   ,  & 
    40       zzero  = 0.e0     ,  & 
    41       zone   = 1.e0 
     37   PUBLIC   lim_thd_2  ! called by lim_step 
     38 
     39   REAL(wp)  ::   epsi20 = 1.e-20   ,  &  ! constant values 
     40      &           epsi16 = 1.e-16   ,  & 
     41      &           epsi04 = 1.e-04   ,  & 
     42      &           zzero  = 0.e0     ,  & 
     43      &           zone   = 1.e0 
    4244 
    4345   !! * Substitutions 
     
    4648   !!-------- ------------------------------------------------------------- 
    4749   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    48    !! $Header$  
    49    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     50   !! $ Id: $ 
     51   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5052   !!---------------------------------------------------------------------- 
    5153 
     
    6870      !!             - back to the geographic grid 
    6971      !! 
    70       !! ** References : 
    71       !!       H. Goosse et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90 
    72       !! 
    73       !! History : 
    74       !!   1.0  !  00-01 (LIM) 
    75       !!   2.0  !  02-07 (C. Ethe, G. Madec) F90 
     72      !! References :   Goosse et al. 1996, Bul. Soc. Roy. Sc. Liege, 65, 87-90 
    7673      !!--------------------------------------------------------------------- 
    7774      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    78  
     75      !! 
    7976      INTEGER  ::   ji, jj,    &   ! dummy loop indices 
    8077         nbpb  ,               &   ! nb of icy pts for thermo. cal. 
     
    9289         zfontn             ,  &   ! heat flux from snow thickness 
    9390         zfntlat, zpareff          ! test. the val. of lead heat budget 
    94       REAL(wp), DIMENSION(jpi,jpj) :: & 
    95          zhicifp            ,  &   ! ice thickness for outputs 
    96          zqlbsbq                   ! link with lead energy budget qldif 
    97       REAL(wp), DIMENSION(jpi,jpj,jpk) :: & 
    98          zmsk                      ! working array 
     91      REAL(wp), DIMENSION(jpi,jpj) ::   zhicifp,   &  ! ice thickness for outputs 
     92         &                              zqlbsbq       ! link with lead energy budget qldif 
     93      REAL(wp), DIMENSION(jpi,jpj,jpk) ::   zmsk      ! working array 
    9994      !!------------------------------------------------------------------- 
    10095 
    101       IF( kt == nit000  )   CALL lim_thd_init_2  ! Initialization (first time-step only) 
     96      IF( kt == nit000 )   CALL lim_thd_init_2  ! Initialization (first time-step only) 
    10297    
    10398      !-------------------------------------------! 
     
    173168      !-------------------------------------------------------------------------- 
    174169 
     170      sst_m(:,:) = sst_m(:,:) + rt0 
     171 
    175172      !CDIR NOVERRCHK 
    176173      DO jj = 1, jpj 
     
    188185            !  temperature and turbulent mixing (McPhee, 1992) 
    189186            zfric_u        = MAX ( MIN( SQRT( ust2s(ji,jj) ) , zfric_umax ) , zfric_umin )  ! friction velocity 
    190             fdtcn(ji,jj)  = zindb * rau0 * rcp * 0.006  * zfric_u * ( sst_io(ji,jj) - tfu(ji,jj) )  
     187            fdtcn(ji,jj)  = zindb * rau0 * rcp * 0.006  * zfric_u * ( sst_m(ji,jj) - tfu(ji,jj) )  
    191188            qdtcn(ji,jj)  = zindb * fdtcn(ji,jj) * frld(ji,jj) * rdt_ice 
    192189                         
    193190            !  partial computation of the lead energy budget (qldif) 
    194191            zfontn         = ( sprecip(ji,jj) / rhosn ) * xlsn  !   energy for melting 
    195             zfnsol         = qnsr_oce(ji,jj)  !  total non solar flux 
    196             qldif(ji,jj)   = tms(ji,jj) * ( qsr_oce(ji,jj) * ( 1.0 - thcm(ji,jj) )   & 
     192            zfnsol         = qns(ji,jj)                         !  total non solar flux over the ocean 
     193            qldif(ji,jj)   = tms(ji,jj) * ( qsr(ji,jj) * ( 1.0 - thcm(ji,jj) )   & 
    197194               &                               + zfnsol + fdtcn(ji,jj) - zfontn     & 
    198195               &                               + ( 1.0 - zindb ) * fsbbq(ji,jj) )   & 
     
    206203             
    207204            !  energy needed to bring ocean surface layer until its freezing 
    208             qcmif  (ji,jj) =  rau0 * rcp * fse3t(ji,jj,1) * ( tfu(ji,jj) - sst_io(ji,jj) ) * ( 1 - zinda ) 
     205            qcmif  (ji,jj) =  rau0 * rcp * fse3t(ji,jj,1) * ( tfu(ji,jj) - sst_m(ji,jj) ) * ( 1 - zinda ) 
    209206             
    210207            !  calculate oceanic heat flux. 
     
    216213      END DO 
    217214       
     215      sst_m(:,:) = sst_m(:,:) - rt0 
    218216       
    219217      !         Select icy points and fulfill arrays for the vectorial grid. 
     
    258256         CALL tab_2d_1d_2( nbpb, fr1_i0_1d  (1:nbpb)     , fr1_i0     , jpi, jpj, npb(1:nbpb) ) 
    259257         CALL tab_2d_1d_2( nbpb, fr2_i0_1d  (1:nbpb)     , fr2_i0     , jpi, jpj, npb(1:nbpb) ) 
    260          CALL tab_2d_1d_2( nbpb, qnsr_ice_1d(1:nbpb)     , qnsr_ice   , jpi, jpj, npb(1:nbpb) ) 
     258         CALL tab_2d_1d_2( nbpb, qns_ice_1d (1:nbpb)     , qns_ice    , jpi, jpj, npb(1:nbpb) ) 
    261259#if ! defined key_coupled 
    262260         CALL tab_2d_1d_2( nbpb, qla_ice_1d (1:nbpb)     , qla_ice    , jpi, jpj, npb(1:nbpb) ) 
     
    404402         CALL prt_ctl(tab2d_1=qstoif, clinfo1=' lim_thd: qstoif  : ', tab2d_2=fsbbq , clinfo2=' fsbbq  : ') 
    405403      ENDIF 
    406  
     404       ! 
    407405    END SUBROUTINE lim_thd_2 
    408406 
     
    419417      !! 
    420418      !! ** input   :   Namelist namicether 
    421       !! 
    422       !! history : 
    423       !!  8.5  ! 03-08 (C. Ethe) original code 
    424419      !!------------------------------------------------------------------- 
    425420      NAMELIST/namicethd/ hmelt , hiccrit, hicmin, hiclim, amax  ,        & 
  • trunk/NEMO/LIM_SRC_2/limthd_lac_2.F90

    r823 r888  
    3030   !!---------------------------------------------------------------------- 
    3131   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    32    !! $Header$  
     32   !! $ Id: $ 
    3333   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    3434   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/limthd_zdf_2.F90

    r823 r888  
    44   !!                thermodynamic growth and decay of the ice  
    55   !!====================================================================== 
     6   !! History :  1.0  !  01-04 (LIM) Original code 
     7   !!            2.0  !  02-08 (C. Ethe, G. Madec) F90 
     8   !!---------------------------------------------------------------------- 
    69#if defined key_lim2 
    710   !!---------------------------------------------------------------------- 
    811   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
     12   !!---------------------------------------------------------------------- 
    913   !!---------------------------------------------------------------------- 
    1014   !!   lim_thd_zdf_2 : vertical accr./abl. and lateral ablation of sea ice 
     
    2226   PRIVATE 
    2327 
    24    !! * Routine accessibility 
    25    PUBLIC lim_thd_zdf_2      ! called by lim_thd_2 
    26  
    27    !! * Module variables 
    28    REAL(wp)  ::           &  ! constant values 
    29       epsi20 = 1.e-20  ,  & 
    30       epsi13 = 1.e-13  ,  & 
    31       zzero  = 0.e0    ,  & 
    32       zone   = 1.e0 
     28   PUBLIC   lim_thd_zdf_2        ! called by lim_thd_2 
     29 
     30   REAL(wp) ::   epsi20 = 1.e-20  ,  &  ! constant values 
     31      &          epsi13 = 1.e-13  ,  & 
     32      &          zzero  = 0.e0    ,  & 
     33      &          zone   = 1.e0 
    3334   !!---------------------------------------------------------------------- 
    3435   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    35    !! $Header$  
    36    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    37    !!---------------------------------------------------------------------- 
     36   !! $ Id: $ 
     37   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     38   !!---------------------------------------------------------------------- 
     39 
    3840CONTAINS 
    3941 
     
    6466      !!              - Performs lateral ablation 
    6567      !! 
    66       !! References : 
    67       !!   Fichefet T. and M. Maqueda 1997, J. Geophys. Res., 102(C6), 12609-12646    
    68       !!   Fichefet T. and M. Maqueda 1999, Clim. Dyn, 15(4), 251-268   
     68      !! References : Fichefet T. and M. Maqueda 1997, J. Geophys. Res., 102(C6), 12609-12646    
     69      !!              Fichefet T. and M. Maqueda 1999, Clim. Dyn, 15(4), 251-268   
     70      !!------------------------------------------------------------------ 
     71      INTEGER, INTENT(in) ::   kideb    ! Start point on which the  the computation is applied 
     72      INTEGER, INTENT(in) ::   kiut     ! End point on which the  the computation is applied 
    6973      !! 
    70       !! History : 
    71       !!   original    : 01-04 (LIM) 
    72       !!   addition    : 02-08 (C. Ethe, G. Madec) 
    73       !!------------------------------------------------------------------ 
    74       !! * Arguments 
    75       INTEGER , INTENT (in) ::  & 
    76          kideb ,  &  ! Start point on which the  the computation is applied 
    77          kiut        ! End point on which the  the computation is applied 
    78  
    79       !! * Local variables 
    8074      INTEGER ::   ji       ! dummy loop indices 
    81  
    82       REAL(wp) , DIMENSION(jpij,2) ::  & 
    83          zqcmlt        ! energy due to surface( /1 ) and bottom melting( /2 ) 
    84  
     75      REAL(wp), DIMENSION(jpij,2) ::   zqcmlt        ! energy due to surface( /1 ) and bottom melting( /2 ) 
    8576      REAL(wp), DIMENSION(jpij) ::  & 
    8677         ztsmlt      &    ! snow/ice surface melting temperature 
     
    9788         , zts_old   &    ! previous surface temperature 
    9889         , zidsn , z1midsn , zidsnic ! tempory variables 
    99  
    100       REAL(wp), DIMENSION(jpij) :: & 
     90      REAL(wp), DIMENSION(jpij) ::   & 
    10191          zfnet       &  ! net heat flux at the top surface( incl. conductive heat flux) 
    10292          , zsprecip  &    ! snow accumulation 
     
    10999          , zfrld_1d    &    ! new sea/ice fraction 
    110100          , zep            ! internal temperature of the 2nd layer of the snow/ice system 
    111  
    112101       REAL(wp), DIMENSION(3) :: &  
    113102          zplediag  &    ! principle diagonal, subdiag. and supdiag. of the  
     
    115104          , zsupdiag  &    ! of the temperatures inside the snow-ice system 
    116105          , zsmbr          ! second member 
    117  
    118106       REAL(wp) :: &  
    119107          zhsu     &     ! thickness of surface layer 
     
    130118          , zb2 , zd2 , zb3 , zd3 & 
    131119          , ztint          ! equivalent temperature at the snow-ice interface 
    132  
    133120       REAL(wp) :: &  
    134121          zexp      &     ! exponential function of the ice thickness 
     
    148135          , zdtic        &  ! ice internal temp. increment 
    149136          , zqnes          ! conductive energy due to ice melting in the first ice layer 
    150  
    151137       REAL(wp) :: &  
    152138          ztbot     &      ! temperature at the bottom surface 
     
    162148          , zc1, zpc1, zc2, zpc2, zp1, zp2 & ! tempory variables 
    163149          , ztb2, ztb3 
    164  
    165150       REAL(wp) :: &  
    166151          zdrmh         &   ! change in snow/ice thick. after snow-ice formation 
     
    181166       !     Computation of energies due to surface and bottom melting  
    182167       !----------------------------------------------------------------------- 
    183         
    184168        
    185169       DO ji = kideb , kiut 
     
    201185       END DO 
    202186 
    203  
    204187       !------------------------------------------- 
    205188       !  2. Calculate some intermediate variables.   
     
    265248       !     - qstbif_1d, total energy stored in brine pockets (updating) 
    266249       !------------------------------------------------------------------- 
    267  
    268250 
    269251       DO ji = kideb , kiut 
     
    288270       END DO 
    289271 
    290  
    291272       !-------------------------------------------------------------------------------- 
    292273       !  4. Computation of the surface temperature : determined by considering the  
     
    333314          !---computation of the energy balance function  
    334315          zfts    = - z1mi0 (ji) * qsr_ice_1d(ji)   & ! net absorbed solar radiation 
    335              &      - qnsr_ice_1d(ji)                & ! total non solar radiation 
    336              &      - zfcsu (ji)                  ! conductive heat flux from the surface 
     316             &      - qns_ice_1d(ji)                & ! total non solar radiation 
     317             &      - zfcsu (ji)                      ! conductive heat flux from the surface 
    337318          !---computation of surface temperature increment   
    338319          zdts    = -zfts / zdfts 
     
    360341          sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 
    361342#if ! defined key_coupled 
    362           qnsr_ice_1d(ji) = qnsr_ice_1d(ji) + dqns_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
     343          qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
    363344          qla_ice_1d (ji) = qla_ice_1d (ji) + dqla_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
    364345#endif 
     
    366347       END DO 
    367348 
    368       
    369  
    370349       !     5.2. Calculate available heat for surface ablation.  
    371350       !--------------------------------------------------------------------- 
    372351 
    373352       DO ji = kideb, kiut 
    374           zfnet(ji) = qnsr_ice_1d(ji) + z1mi0(ji) * qsr_ice_1d(ji) + zfcsu(ji)           
     353          zfnet(ji) = qns_ice_1d(ji) + z1mi0(ji) * qsr_ice_1d(ji) + zfcsu(ji)           
    375354          zfnet(ji) = MAX( zzero , zfnet(ji) ) 
    376355          zfnet(ji) = zfnet(ji) * MAX( zzero , SIGN( zone , sist_1d(ji) - ztsmlt(ji) ) ) 
     
    730709          dvnbq_1d(ji) = ( 1.0 - frld_1d(ji) ) * ( zhicnew - h_ice_1d(ji) ) 
    731710          dmgwi_1d(ji) = dmgwi_1d(ji) + ( 1.0 -frld_1d(ji) ) * ( h_snow_1d(ji) - zhsnnew ) * rhosn 
    732           !  case of natural freshwater flux 
    733 #if defined key_lim_fdd   
    734           rdmicif_1d(ji) = rdmicif_1d(ji) + ( 1.0 - frld_1d(ji) )   * ( zhicnew - h_ice_1d(ji) )  * rhoic 
     711          !---  volume change of ice and snow (used for ocean-ice freshwater flux computation) 
     712          rdmicif_1d(ji) = rdmicif_1d(ji) + ( 1.0 - frld_1d(ji) )   * ( zhicnew - h_ice_1d (ji) ) * rhoic 
    735713          rdmsnif_1d(ji) = rdmsnif_1d(ji) + ( 1.0 - frld_1d(ji) )   * ( zhsnnew - h_snow_1d(ji) ) * rhosn 
    736714 
    737 #else 
    738           rdmicif_1d(ji) = rdmicif_1d(ji) + ( 1.0 - frld_1d(ji) ) * (  ( zhicnew - h_ice_1d(ji) ) * rhoic   & 
    739              &                                                    + ( zhsnnew - h_snow_1d(ji) ) * rhosn ) 
    740 #endif 
    741  
    742715          !---  Actualize new snow and ice thickness. 
    743  
    744716          h_snow_1d(ji)  = zhsnnew 
    745           h_ice_1d(ji)  = zhicnew 
     717          h_ice_1d (ji)  = zhicnew 
    746718 
    747719       END DO 
     
    799771          qstbif_1d(ji) = zdrfrl2 * qstbif_1d(ji) 
    800772          frld_1d(ji)    = zfrld_1d(ji) 
    801  
    802        END DO 
    803         
     773          ! 
     774       END DO 
     775       !  
    804776    END SUBROUTINE lim_thd_zdf_2 
     777 
    805778#else 
    806    !!====================================================================== 
    807    !!                       ***  MODULE limthd_zdf_2   *** 
    808    !!                              no sea ice model   
    809    !!====================================================================== 
     779   !!---------------------------------------------------------------------- 
     780   !!   Default Option                                     NO sea-ice model 
     781   !!---------------------------------------------------------------------- 
    810782CONTAINS 
    811783   SUBROUTINE lim_thd_zdf_2          ! Empty routine 
    812784   END SUBROUTINE lim_thd_zdf_2 
    813785#endif 
    814  END MODULE limthd_zdf_2 
     786 
     787   !!====================================================================== 
     788END MODULE limthd_zdf_2 
  • trunk/NEMO/LIM_SRC_2/limtrp_2.F90

    r823 r888  
    3030 
    3131   !! * Routine accessibility 
    32    PUBLIC lim_trp_2     ! called by ice_step_2 
     32   PUBLIC lim_trp_2     ! called by sbc_ice_lim_2 
    3333 
    3434   !! * Shared module variables 
     
    4848   !!---------------------------------------------------------------------- 
    4949   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    50    !! $Header$  
     50   !! $ Id: $ 
    5151   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    5252   !!---------------------------------------------------------------------- 
     
    112112         DO jj = 1, jpjm1 
    113113            DO ji = 1, jpim1 
    114                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 ) ) 
    115                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 ) ) 
     114               zui_u(ji,jj) = ( ui_ice(ji+1,jj  ) + ui_ice(ji+1,jj+1) ) / ( MAX( tmu(ji+1,jj  ) + tmu(ji+1,jj+1), zvbord ) ) 
     115               zvi_v(ji,jj) = ( vi_ice(ji  ,jj+1) + vi_ice(ji+1,jj+1) ) / ( MAX( tmu(ji  ,jj+1) + tmu(ji+1,jj+1), zvbord ) ) 
    116116            END DO 
    117117         END DO 
     
    128128         IF (lk_mpp ) CALL mpp_max(zcfl) 
    129129 
    130          IF ( zcfl > 0.5 .AND. lwp )   WRITE(numout,*) 'lim_trp : violation of cfl criterion the ',nday,'th day, cfl = ',zcfl 
     130         IF ( zcfl > 0.5 .AND. lwp )   WRITE(numout,*) 'lim_trp_2 : violation of cfl criterion the ',nday,'th day, cfl = ',zcfl 
    131131 
    132132         ! content of properties 
  • trunk/NEMO/LIM_SRC_2/limwri_2.F90

    r823 r888  
    99#if defined key_lim2 
    1010   !!---------------------------------------------------------------------- 
    11    !!   'key_lim2'  i                                 LIM 2.0 sea-ice model 
     11   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    1212   !!---------------------------------------------------------------------- 
    1313   !!---------------------------------------------------------------------- 
     
    1515   !!   lim_wri_init_2 : initialization and namelist read 
    1616   !!---------------------------------------------------------------------- 
    17    USE ioipsl 
    18    USE dianam    ! build name of file (routine) 
    1917   USE phycst 
    2018   USE dom_oce 
    2119   USE daymod 
    22    USE in_out_manager 
    2320   USE ice_oce         ! ice variables 
    24    USE flx_oce 
     21   USE sbc_oce 
     22   USE sbc_ice 
    2523   USE dom_ice_2 
    2624   USE ice_2 
     25 
    2726   USE lbclnk 
     27   USE dianam          ! build name of file (routine) 
     28   USE in_out_manager 
     29   USE ioipsl 
    2830 
    2931   IMPLICIT NONE 
    3032   PRIVATE 
    3133 
    32    PUBLIC   lim_wri_2        ! routine called by lim_step_2.F90 
     34   PUBLIC   lim_wri_2      ! routine called by sbc_ice_lim_2 
    3335 
    3436   INTEGER, PARAMETER                       ::   jpnoumax = 40   ! maximum number of variable for ice output 
     
    4951      zone   = 1.e0 
    5052 
    51    !!---------------------------------------------------------------------- 
    52    !!  LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    53    !! $Header$ 
     53   !! * Substitutions 
     54#   include "vectopt_loop_substitute.h90" 
     55   !!---------------------------------------------------------------------- 
     56   !!  LIM 2.0, UCL-LOCEAN-IPSL (2006) 
     57   !! $ Id: $ 
    5458   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5559   !!---------------------------------------------------------------------- 
     
    7983      !!------------------------------------------------------------------- 
    8084      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    81  
     85      !! 
    8286      INTEGER  ::   ji, jj, jf                      ! dummy loop indices 
    8387      CHARACTER(len = 40)  ::   clhstnam, clop 
     
    9094 
    9195      !                                          !--------------------! 
    92       IF ( kt == nit000 ) THEN                !   Initialisation   ! 
     96      IF( kt == nit000 ) THEN                    !   Initialisation   ! 
    9397         !                                       !--------------------! 
    9498         CALL lim_wri_init_2  
     
    97101!!Chris         clop     = "ave(only(x))"      !ibug  namelist parameter a ajouter 
    98102         clop     = "ave(x)" 
    99          zout     = nwrite * rdt_ice / nfice 
     103         zout     = nwrite * rdt_ice / nn_fsbc 
    100104         zsec     = 0. 
    101105         niter    = 0 
     
    110114          
    111115         DO jf = 1, noumef 
    112             IF ( nc(jf) == 1 )   CALL histdef( nice, nam(jf), titn(jf), uni(jf), jpi, jpj   & 
    113                   &                                , nhorid, 1, 1, 1, -99, 32, clop, zsto, zout ) 
     116            IF( nc(jf) == 1 )   CALL histdef( nice, nam(jf), titn(jf), uni(jf), jpi, jpj   & 
     117                                               , nhorid, 1, 1, 1, -99, 32, clop, zsto, zout ) 
    114118         END DO 
    115119         CALL histend( nice ) 
    116           
     120         ! 
    117121      ENDIF 
    118122      !                                          !--------------------! 
     
    120124      !                                          !--------------------! 
    121125 
    122 !!gm  change the print below to have it only at output time step or when nitend =< 100 
    123       IF(lwp) THEN 
    124          WRITE(numout,*) 
    125          WRITE(numout,*) 'lim_wri_2 : write ice outputs in NetCDF files at time : ', nyear, nmonth, nday, kt + nfice - 1 
    126          WRITE(numout,*) '~~~~~~~~~ ' 
    127       ENDIF 
    128  
    129       !-- calculs des valeurs instantanees 
     126      !-- Store instantaneous values in zcmo 
    130127       
    131128      zcmo(:,:, 1:jpnoumax ) = 0.e0  
    132129      DO jj = 2 , jpjm1 
    133          DO ji = 2 , jpim1 
     130         DO ji = fs_2 , fs_jpim1 
    134131            zindh  = MAX( zzero , SIGN( zone , hicif(ji,jj) * (1.0 - frld(ji,jj) ) - 0.10 ) ) 
    135132            zinda  = MAX( zzero , SIGN( zone , ( 1.0 - frld(ji,jj) ) - 0.10 ) ) 
     
    142139            zcmo(ji,jj,5)  = sist  (ji,jj) 
    143140            zcmo(ji,jj,6)  = fbif  (ji,jj) 
    144             zcmo(ji,jj,7)  = zindb * (  u_ice(ji,jj  ) * tmu(ji,jj  ) + u_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
    145                                       + u_ice(ji,jj+1) * tmu(ji,jj+1) + u_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     141            zcmo(ji,jj,7)  = zindb * (  ui_ice(ji,jj  ) * tmu(ji,jj  ) + ui_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
     142                                      + ui_ice(ji,jj+1) * tmu(ji,jj+1) + ui_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
    146143                                  / ztmu  
    147144 
    148             zcmo(ji,jj,8)  = zindb * (  v_ice(ji,jj  ) * tmu(ji,jj  ) + v_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
    149                                       + v_ice(ji,jj+1) * tmu(ji,jj+1) + v_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     145            zcmo(ji,jj,8)  = zindb * (  vi_ice(ji,jj  ) * tmu(ji,jj  ) + vi_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
     146                                      + vi_ice(ji,jj+1) * tmu(ji,jj+1) + vi_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
    150147                                  / ztmu 
    151             zcmo(ji,jj,9)  = sst_io(ji,jj) 
    152             zcmo(ji,jj,10) = sss_io(ji,jj) 
    153  
    154             zcmo(ji,jj,11) = fnsolar(ji,jj) + fsolar(ji,jj) 
    155             zcmo(ji,jj,12) = fsolar (ji,jj) 
    156             zcmo(ji,jj,13) = fnsolar(ji,jj) 
     148            zcmo(ji,jj,9)  = sst_m(ji,jj) 
     149            zcmo(ji,jj,10) = sss_m(ji,jj) 
     150            zcmo(ji,jj,11) = qns(ji,jj) + qsr(ji,jj) 
     151            zcmo(ji,jj,12) = qsr(ji,jj) 
     152            zcmo(ji,jj,13) = qns(ji,jj) 
    157153            ! See thersf for the coefficient 
    158             zcmo(ji,jj,14) = - fsalt(ji,jj) * rday * ( sss_io(ji,jj) + epsi16 ) / soce 
    159             zcmo(ji,jj,15) = gtaux(ji,jj) 
    160             zcmo(ji,jj,16) = gtauy(ji,jj) 
    161             zcmo(ji,jj,17) = ( 1.0 - frld(ji,jj) ) * qsr_ice (ji,jj) + frld(ji,jj) * qsr_oce (ji,jj) 
    162             zcmo(ji,jj,18) = ( 1.0 - frld(ji,jj) ) * qnsr_ice(ji,jj) + frld(ji,jj) * qnsr_oce(ji,jj) 
     154            zcmo(ji,jj,14) = - emps(ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce    !!gm ??? 
     155            zcmo(ji,jj,15) = utaui_ice(ji,jj) 
     156            zcmo(ji,jj,16) = vtaui_ice(ji,jj) 
     157            zcmo(ji,jj,17) = qsr_ice(ji,jj) 
     158            zcmo(ji,jj,18) = qns_ice(ji,jj) 
    163159            zcmo(ji,jj,19) = sprecip(ji,jj) 
    164160         END DO 
     
    175171         END DO 
    176172          
    177          IF ( jf == 7  .OR. jf == 8  .OR. jf == 15 .OR. jf == 16 ) THEN  
     173         IF( jf == 7  .OR. jf == 8  .OR. jf == 15 .OR. jf == 16 ) THEN 
    178174            CALL lbc_lnk( zfield, 'T', -1. ) 
    179175         ELSE  
     
    181177         ENDIF 
    182178          
    183          IF ( nc(jf) == 1 ) CALL histwrite( nice, nam(jf), niter, zfield, ndim, ndex51 ) 
     179         IF( nc(jf) == 1 )  CALL histwrite( nice, nam(jf), niter, zfield, ndim, ndex51 ) 
    184180          
    185181      END DO 
    186182       
    187       IF ( ( nfice * niter + nit000 - 1 ) >= nitend ) THEN 
    188          CALL histclo( nice )  
    189       ENDIF 
     183      IF( ( nn_fsbc * niter + nit000 - 1 ) >= nitend )   CALL histclo( nice )  
    190184      ! 
    191185   END SUBROUTINE lim_wri_2 
     
    225219         field_13, field_14, field_15, field_16, field_17, field_18,   & 
    226220         field_19 
    227 !!gm      NAMELIST/namiceout/ noumef, & 
    228 !!           zfield( 1), zfield( 2), zfield( 3), zfield( 4), zfield( 5),   & 
    229 !!           zfield( 6), zfield( 7), zfield( 8), zfield( 9), zfield(10),   & 
    230 !!           zfield(11), zfield(12), zfield(13), zfield(14), zfield(15),   & 
    231 !!gm         zfield(16), zfield(17), zfield(18), zfield(19) 
    232       !!------------------------------------------------------------------- 
    233  
    234       ! Read Namelist namicewri 
    235       REWIND ( numnam_ice ) 
     221      !!------------------------------------------------------------------- 
     222 
     223      REWIND ( numnam_ice )                ! Read Namelist namicewri 
    236224      READ   ( numnam_ice  , namiceout ) 
    237225       
    238       zfield(1) = field_1 
    239       zfield(2) = field_2 
    240       zfield(3) = field_3 
    241       zfield(4) = field_4 
    242       zfield(5) = field_5 
    243       zfield(6) = field_6 
    244       zfield(7) = field_7 
    245       zfield(8) = field_8 
    246       zfield(9) = field_9 
     226      zfield( 1) = field_1 
     227      zfield( 2) = field_2 
     228      zfield( 3) = field_3 
     229      zfield( 4) = field_4 
     230      zfield( 5) = field_5 
     231      zfield( 6) = field_6 
     232      zfield( 7) = field_7 
     233      zfield( 8) = field_8 
     234      zfield( 9) = field_9 
    247235      zfield(10) = field_10 
    248236      zfield(11) = field_11 
     
    274262         DO nf = 1 , noumef          
    275263            WRITE(numout,*) '   ', titn(nf), '   ', nam(nf),'      ', uni(nf),'  ', nc(nf),'        ', cmulti(nf),   & 
    276                '        ', cadd(nf) 
     264               &       '        ', cadd(nf) 
    277265         END DO 
    278266      ENDIF 
  • trunk/NEMO/LIM_SRC_2/limwri_dimg_2.h90

    r823 r888  
    22   !!---------------------------------------------------------------------- 
    33   !!  LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    4    !! $Header$ 
     4   !! $ Id: $ 
    55   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    66   !!---------------------------------------------------------------------- 
     
    8282 
    8383       zsto     = rdt_ice 
    84        zout     = nwrite * rdt_ice / nfice 
     84       zout     = nwrite * rdt_ice / nn_fsbc 
    8585       zsec     = 0. 
    8686       niter    = 0 
     
    106106          zcmo(ji,jj,5)  = sist  (ji,jj) 
    107107          zcmo(ji,jj,6)  = fbif  (ji,jj) 
    108           zcmo(ji,jj,7)  = zindb * (  u_ice(ji,jj  ) * tmu(ji,jj  ) + u_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
    109                + u_ice(ji,jj+1) * tmu(ji,jj+1) + u_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     108          zcmo(ji,jj,7)  = zindb * (  ui_ice(ji,jj  ) * tmu(ji,jj  ) + ui_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
     109             &                      + ui_ice(ji,jj+1) * tmu(ji,jj+1) + ui_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
    110110               / ztmu  
    111111 
    112           zcmo(ji,jj,8)  = zindb * (  v_ice(ji,jj  ) * tmu(ji,jj  ) + v_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
    113                + v_ice(ji,jj+1) * tmu(ji,jj+1) + v_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     112          zcmo(ji,jj,8)  = zindb * (  vi_ice(ji,jj  ) * tmu(ji,jj  ) + vi_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
     113             &                      + vi_ice(ji,jj+1) * tmu(ji,jj+1) + vi_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
    114114               / ztmu 
    115           zcmo(ji,jj,9)  = sst_io(ji,jj) 
    116           zcmo(ji,jj,10) = sss_io(ji,jj) 
    117  
    118           zcmo(ji,jj,11) = fnsolar(ji,jj) + fsolar(ji,jj) 
    119           zcmo(ji,jj,12) = fsolar (ji,jj) 
    120           zcmo(ji,jj,13) = fnsolar(ji,jj) 
     115          zcmo(ji,jj,9)  = sst_m(ji,jj) 
     116          zcmo(ji,jj,10) = sss_m(ji,jj) 
     117 
     118          zcmo(ji,jj,11) = qns(ji,jj) + qsr(ji,jj) 
     119          zcmo(ji,jj,12) = qsr(ji,jj) 
     120          zcmo(ji,jj,13) = qns(ji,jj) 
    121121          ! See thersf for the coefficient 
    122           zcmo(ji,jj,14) = - fsalt(ji,jj) * rday * ( sss_io(ji,jj) + epsi16 ) / soce 
    123           zcmo(ji,jj,15) = gtaux(ji,jj) 
    124           zcmo(ji,jj,16) = gtauy(ji,jj) 
    125           zcmo(ji,jj,17) = ( 1.0 - frld(ji,jj) ) * qsr_ice (ji,jj) + frld(ji,jj) * qsr_oce (ji,jj) 
    126           zcmo(ji,jj,18) = ( 1.0 - frld(ji,jj) ) * qnsr_ice(ji,jj) + frld(ji,jj) * qnsr_oce(ji,jj) 
     122          zcmo(ji,jj,14) = - emps(ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce 
     123          zcmo(ji,jj,15) = utaui_ice(ji,jj) 
     124          zcmo(ji,jj,16) = vtaui_ice(ji,jj) 
     125          zcmo(ji,jj,17) = qsr_ice(ji,jj) 
     126          zcmo(ji,jj,18) = qns_ice(ji,jj) 
    127127          zcmo(ji,jj,19) = sprecip(ji,jj) 
    128128       END DO 
     
    132132    nmoyice = nmoyice + 1  
    133133    ! compute mean value if it is time to write on file 
    134     IF ( MOD(kt+nfice-1-nit000+1,nwrite) == 0 ) THEN 
     134    IF ( MOD(kt+nn_fsbc-1-nit000+1,nwrite) == 0 ) THEN 
    135135       rcmoy(:,:,:) = rcmoy(:,:,:) / FLOAT(nmoyice) 
    136136#else   
    137        IF ( MOD(kt-nfice-1-nit000+1,nwrite) == 0 ) THEN  
     137       IF ( MOD(kt-nn_fsbc-1-nit000+1,nwrite) == 0 ) THEN  
    138138          !  case of instantaneaous output rcmoy(:,:, 1:jpnoumax ) = 0.e0 
    139139          DO jj = 2 , jpjm1 
     
    149149                rcmoy(ji,jj,5)  = sist  (ji,jj) 
    150150                rcmoy(ji,jj,6)  = fbif  (ji,jj) 
    151                 rcmoy(ji,jj,7)  = zindb * (  u_ice(ji,jj  ) * tmu(ji,jj  ) + u_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
    152                      + u_ice(ji,jj+1) * tmu(ji,jj+1) + u_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     151                rcmoy(ji,jj,7)  = zindb * (  ui_ice(ji,jj  ) * tmu(ji,jj  ) + ui_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
     152                   &                       + ui_ice(ji,jj+1) * tmu(ji,jj+1) + ui_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
    153153                     / ztmu 
    154154 
    155                 rcmoy(ji,jj,8)  = zindb * (  v_ice(ji,jj  ) * tmu(ji,jj  ) + v_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
    156                      + v_ice(ji,jj+1) * tmu(ji,jj+1) + v_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     155                rcmoy(ji,jj,8)  = zindb * (  vi_ice(ji,jj  ) * tmu(ji,jj  ) + vi_ice(ji+1,jj  ) * tmu(ji+1,jj  )   & 
     156                   &                       + vi_ice(ji,jj+1) * tmu(ji,jj+1) + vi_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
    157157                     / ztmu 
    158                 rcmoy(ji,jj,9)  = sst_io(ji,jj) 
    159                 rcmoy(ji,jj,10) = sss_io(ji,jj) 
    160  
    161                 rcmoy(ji,jj,11) = fnsolar(ji,jj) + fsolar(ji,jj) 
    162                 rcmoy(ji,jj,12) = fsolar (ji,jj) 
    163                 rcmoy(ji,jj,13) = fnsolar(ji,jj) 
     158                rcmoy(ji,jj,9)  = sst_m(ji,jj) 
     159                rcmoy(ji,jj,10) = sss_m(ji,jj) 
     160 
     161                rcmoy(ji,jj,11) = qns(ji,jj) + qsr(ji,jj) 
     162                rcmoy(ji,jj,12) = qsr(ji,jj) 
     163                rcmoy(ji,jj,13) = qns(ji,jj) 
    164164                ! See thersf for the coefficient 
    165                 rcmoy(ji,jj,14) = - fsalt(ji,jj) * rday * ( sss_io(ji,jj) + epsi16 ) / soce 
    166                 rcmoy(ji,jj,15) = gtaux(ji,jj) 
    167                 rcmoy(ji,jj,16) = gtauy(ji,jj) 
    168                 rcmoy(ji,jj,17) = ( 1.0 - frld(ji,jj) ) * qsr_ice (ji,jj) + frld(ji,jj) * qsr_oce (ji,jj) 
    169                 rcmoy(ji,jj,18) = ( 1.0 - frld(ji,jj) ) * qnsr_ice(ji,jj) + frld(ji,jj) * qnsr_oce(ji,jj) 
     165                rcmoy(ji,jj,14) = - emps(ji,jj) * rday * ( sss_m(ji,jj) + epsi16 ) / soce 
     166                rcmoy(ji,jj,15) = utaui_ice(ji,jj) 
     167                rcmoy(ji,jj,16) = vtaui_ice(ji,jj) 
     168                rcmoy(ji,jj,17) = qsr_ice(ji,jj) 
     169                rcmoy(ji,jj,18) = qns_ice(ji,jj) 
    170170                rcmoy(ji,jj,19) = sprecip(ji,jj) 
    171171             END DO 
     
    201201          rcmoy(:,:,:) = 0.0 
    202202          nmoyice = 0  
    203        END IF     !  MOD(kt+nfice-1-nit000+1, nwrite == 0 ) ! 
     203       END IF     !  MOD(kt+nn_fsbc-1-nit000+1, nwrite == 0 ) ! 
    204204 
    205205     END SUBROUTINE lim_wri_2 
  • trunk/NEMO/LIM_SRC_2/par_ice_2.F90

    r823 r888  
    77   !!---------------------------------------------------------------------- 
    88   !!  LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    9    !! $Header$ 
     9   !! $ Id: $ 
    1010   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    1111   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_2/thd_ice_2.F90

    r823 r888  
    99   !!---------------------------------------------------------------------- 
    1010   !!   LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    11    !! $Header$ 
     11   !! $ Id: $ 
    1212   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    1313   !!---------------------------------------------------------------------- 
     
    5757      fr1_i0_1d   ,     &  !:    "                  "      fr1_i0 
    5858      fr2_i0_1d   ,     &  !:    "                  "      fr2_i0 
    59       qnsr_ice_1d ,     &  !:    "                  "      qns_ice 
     59      qns_ice_1d ,     &  !:    "                  "      qns_ice 
    6060      qfvbq_1d    ,     &  !:    "                  "      qfvbq 
    6161      sist_1d     ,     &  !:    "                  "      sist 
  • trunk/NEMO/LIM_SRC_3/ice.F90

    r834 r888  
    99   !!---------------------------------------------------------------------- 
    1010   !!  LIM 3.0, UCL-LOCEAN-IPSL (2005) 
    11    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/ice.F90,v 1.4 2005/03/27 18:34:41 opalod Exp $ 
     11   !! $ Id: $ 
    1212   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    1313   !!---------------------------------------------------------------------- 
     
    493493      diag_bot_me,                           & ! vertical bottom melt  
    494494      diag_sur_me                              ! vertical surface melt 
    495    INTEGER , PUBLIC ::   &                      !: indexes of the debugging 
    496       jiindex,           &                      !  point 
    497       jjindex 
     495   INTEGER , PUBLIC ::   jiindx, jjindx        !: indexes of the debugging point 
    498496 
    499497#else 
  • trunk/NEMO/LIM_SRC_3/iceini.F90

    r862 r888  
    1313   USE in_out_manager 
    1414   USE ice_oce         ! ice variables 
    15    USE flx_oce 
     15   USE sbc_oce         ! Surface boundary condition: ocean fields 
     16   USE sbc_ice         ! Surface boundary condition: ice fields 
    1617   USE phycst          ! Define parameters for the routines 
    1718   USE ocfzpt 
     
    5051   !!---------------------------------------------------------------------- 
    5152   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
    52    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/iceini.F90,v 1.4 2005/03/27 18:34:41 opalod Exp $  
     53   !! $ Id: $ 
    5354   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    5455   !!---------------------------------------------------------------------- 
     
    7475      ! Louvain la Neuve Ice model 
    7576      IF( nacc == 1 ) THEN 
    76           dtsd2   = nfice * rdtmin * 0.5 
    77           rdt_ice = nfice * rdtmin 
     77          dtsd2   = nn_fsbc * rdtmin * 0.5 
     78          rdt_ice = nn_fsbc * rdtmin 
    7879      ELSE 
    79           dtsd2   = nfice * rdt * 0.5 
    80           rdt_ice = nfice * rdt 
     80          dtsd2   = nn_fsbc * rdt * 0.5 
     81          rdt_ice = nn_fsbc * rdt 
    8182      ENDIF 
    8283 
     
    104105      freeze(:,:) = at_i(:,:)   ! initialisation of sea/ice cover     
    105106# if defined key_coupled 
    106       alb_ice(:,:) = albege(:,:)      ! sea-ice albedo 
     107      Must be adpated to LIM3  
     108      alb_ice(:,:,:) = albege(:,:)      ! sea-ice albedo 
    107109# endif 
    108110       
    109       nstart = numit  + nfice       
     111      nstart = numit  + nn_fsbc       
    110112      nitrun = nitend - nit000 + 1  
    111113      nlast  = numit  + nitrun  
     
    188190 
    189191       WRITE(numout,*) 'lim_itd_ini : Initialization of ice thickness distribution ' 
    190        WRITE(numout,*) '~~~~~~~~~~~~~~~' 
     192       WRITE(numout,*) '~~~~~~~~~~~~' 
    191193 
    192194!!-- End of declarations 
  • trunk/NEMO/LIM_SRC_3/limadv.F90

    r868 r888  
    3636   !!---------------------------------------------------------------------- 
    3737   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    38    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limadv.F90,v 1.4 2005/03/27 18:34:41 opalod Exp $  
     38   !! $ Id: $ 
    3939   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    4040   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_3/limdia.F90

    r869 r888  
    77   !! 1) in lim_dia : add its definition for both hemispheres if wished 
    88   !! 2) add the new titles in lim_dia_init 
    9    !! 
     9   !!---------------------------------------------------------------------- 
    1010#if defined key_lim3 
    1111   !!---------------------------------------------------------------------- 
     
    2626   USE limistate 
    2727   USE dom_oce 
     28   USE sbc_oce         ! Surface boundary condition: ocean fields 
    2829 
    2930   IMPLICIT NONE 
     
    7374   !!---------------------------------------------------------------------- 
    7475   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005) 
    75    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limdia.F90,v 1.5 2005/03/27 18:34:41 opalod Exp $ 
     76   !! $ Id: $ 
    7677   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    7778   !!---------------------------------------------------------------------- 
     
    107108       !--------------------------------------- 
    108109       zday_min = 273.0        ! zday_min = date of minimum extent, here September 30th 
    109        zday = FLOAT(numit-nit000) * rdt_ice / ( 86400.0 * FLOAT(nfice) ) 
     110       zday = FLOAT(numit-nit000) * rdt_ice / ( 86400.0 * FLOAT(nn_fsbc) ) 
    110111       IF (zday.GT.zday_min) THEN  
    111112          zshift_date  =  zday - zday_min 
     
    142143                vinfor(31) = vinfor(31) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
    143144                                                        v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12  
    144                 vinfor(53) = vinfor(53) + fsalt(ji,jj)*aire(ji,jj) / 1.0e12 !salt flux 
     145                vinfor(53) = vinfor(53) + emps(ji,jj)*aire(ji,jj) / 1.0e12 !salt flux 
    145146                vinfor(55) = vinfor(55) + fsbri(ji,jj)*aire(ji,jj) / 1.0e12 !brine drainage flux 
    146147                vinfor(57) = vinfor(57) + fseqv(ji,jj)*aire(ji,jj) / 1.0e12 !equivalent salt flux 
    147                 vinfor(59) = vinfor(59) + sst_io(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
    148                 vinfor(61) = vinfor(61) + sss_io(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
     148                vinfor(59) = vinfor(59) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
     149                vinfor(61) = vinfor(61) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
    149150                vinfor(65) = vinfor(65) + et_s(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12  ! snow temperature 
    150151                vinfor(67) = vinfor(67) + et_i(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12       ! ice heat content 
     
    155156                vinfor(77) = vinfor(77) + v_i(ji,jj,5)*aire(ji,jj) / 1.0e12 !ice volume 
    156157                vinfor(79) = 0.0 
    157                 vinfor(81) = vinfor(81) + fmass(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
     158                vinfor(81) = vinfor(81) + emp(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
    158159             ENDIF 
    159160          END DO 
     
    293294                vinfor(32) = vinfor(32) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
    294295                                                        v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12 !ice vel 
    295                 vinfor(54) = vinfor(54) + at_i(ji,jj)*fsalt(ji,jj)*aire(ji,jj) / 1.0e12 ! Total salt flux 
     296                vinfor(54) = vinfor(54) + at_i(ji,jj)*emps(ji,jj)*aire(ji,jj) / 1.0e12 ! Total salt flux 
    296297                vinfor(56) = vinfor(56) + at_i(ji,jj)*fsbri(ji,jj)*aire(ji,jj) / 1.0e12 ! Brine drainage salt flux 
    297298                vinfor(58) = vinfor(58) + at_i(ji,jj)*fseqv(ji,jj)*aire(ji,jj) / 1.0e12 ! Equivalent salt flux 
    298                 vinfor(60) = vinfor(60) + sst_io(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
    299                 vinfor(62) = vinfor(62) + sss_io(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
     299                vinfor(60) = vinfor(60) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SST 
     300                vinfor(62) = vinfor(62) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) / 1.0e12  !SSS 
    300301                vinfor(66) = vinfor(66) + et_s(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12 ! snow temperature 
    301302                vinfor(68) = vinfor(68) + et_i(ji,jj)/1.0e9*aire(ji,jj) / 1.0e12 ! ice enthalpy 
     
    306307                vinfor(78) = vinfor(78) + v_i(ji,jj,5)*aire(ji,jj) / 1.0e12 !ice volume 
    307308                vinfor(80) = 0.0 
    308                 vinfor(82) = vinfor(82) + fmass(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
     309                vinfor(82) = vinfor(82) + emp(ji,jj)*aire(ji,jj) / 1.0e12 ! mass flux 
    309310             ENDIF 
    310311          END DO 
  • trunk/NEMO/LIM_SRC_3/limdyn.F90

    r869 r888  
    1616   USE dom_ice 
    1717   USE dom_oce         ! ocean space and time domain 
    18    USE taumod 
    1918   USE ice 
    2019   USE par_ice 
     20   USE sbc_ice         ! Surface boundary condition: ice fields 
    2121   USE ice_oce 
    2222   USE iceini 
     
    4141   !!---------------------------------------------------------------------- 
    4242   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
    43    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limdyn.F90,v 1.5 2005/03/27 18:34:41 opalod Exp $ 
    44    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     43   !! $ Id: $ 
     44   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4545   !!---------------------------------------------------------------------- 
    4646 
     
    9090      IF ( ln_limdyn ) THEN 
    9191 
    92          ! ocean velocity 
    93          u_oce(:,:)  = u_io(:,:) * tmu(:,:) 
    94          v_oce(:,:)  = v_io(:,:) * tmv(:,:) 
    95           
    9692         old_u_ice(:,:) = u_ice(:,:) * tmu(:,:) 
    9793         old_v_ice(:,:) = v_ice(:,:) * tmv(:,:) 
     
    162158         ENDIF 
    163159 
    164          ! Ice-Ocean stress 
    165          ! ================ 
    166          DO jj = 2, jpjm1 
    167             zsang  = SIGN(1.e0, gphif(1,jj-1) ) * sangvg 
    168  
    169             DO ji = fs_2, fs_jpim1 
    170                ! computation of wind stress over ocean in X and Y direction 
    171 #if defined key_coupled && defined key_lim_cp1 
    172 !              ztairx =  ( 1.0 - at_i(ji-1,jj)   ) * gtaux(ji-1,jj)   + & 
    173 !                        ( 1.0 - at_i(ji,jj)     ) * gtaux(ji,jj  )   + & 
    174 !                        ( 1.0 - at_i(ji-1,jj-1) ) * gtaux(ji-1,jj-1) + &  
    175 !                        ( 1.0 - at_i(ji,jj-1)   ) * gtaux(ji,jj-1) 
    176  
    177 !              ztairy =  ( 1.0 - at_i(ji-1,jj)   ) * gtauy(ji-1,jj  ) + & 
    178 !                        ( 1.0 - at_i(ji,jj  )   ) * gtauy(ji,jj    ) + & 
    179 !                        ( 1.0 - at_i(ji-1,jj-1) ) * gtauy(ji-1,jj-1) + &  
    180 !                        ( 1.0 - at_i(ji,jj-1)   ) * gtauy(ji,jj-1) 
    181 #else 
    182                ztairx =  ( 2.0 - at_i(ji,jj) - at_i(ji+1,jj) ) * gtaux(ji,jj) / cai * cao 
    183                ztairy =  ( 2.0 - at_i(ji,jj) - at_i(ji,jj+1) ) * gtauy(ji,jj) / cai * cao 
    184  
    185                zsfrldmx2 = at_i(ji,jj) + at_i(ji+1,jj) 
    186                zsfrldmy2 = at_i(ji,jj) + at_i(ji,jj+1) 
    187  
    188 #endif 
    189                zu_ice   = u_ice(ji,jj) - u_oce(ji,jj) 
    190                zv_ice   = v_ice(ji,jj) - v_oce(ji,jj) 
    191                zmod     = SQRT( zu_ice * zu_ice + zv_ice * zv_ice )  
    192  
    193                ! quadratic drag formulation 
    194                ztglx   = zsfrldmx2 * rhoco * zmod * ( cangvg * zu_ice - zsang * zv_ice )  
    195                ztgly   = zsfrldmy2 * rhoco * zmod * ( cangvg * zv_ice + zsang * zu_ice )  
    196 ! 
    197 !              ! IMPORTANT 
    198 !              ! these lignes are bound to prevent numerical oscillations 
    199 !              ! in the ice-ocean stress 
    200 !              ! They are physically ill-based. There is a cleaner solution 
    201 !              ! to try (remember discussion in Paris Gurvan) 
    202 ! 
    203                ztglx   = ztglx * exp( - zmod / 0.5 ) 
    204                ztgly   = ztglx * exp( - zmod / 0.5 ) 
    205  
    206                tio_u(ji,jj) = - ( ztairx + 1.0 * ztglx ) / ( 2. * rau0 ) 
    207                tio_v(ji,jj) = - ( ztairy + 1.0 * ztgly ) / ( 2. * rau0 ) 
    208             END DO 
    209          END DO 
    210           
    211160         ! computation of friction velocity 
    212161         DO jj = 2, jpjm1 
    213162            DO ji = fs_2, fs_jpim1 
    214163 
    215                zu_ice   = u_ice(ji,jj) - u_io(ji,jj) 
     164               zu_ice   = u_ice(ji,jj) - u_oce(ji,jj) 
    216165               zt11  = rhoco * zu_ice * zu_ice 
    217166 
    218                zu_ice   = u_ice(ji-1,jj) - u_io(ji-1,jj) 
     167               zu_ice   = u_ice(ji-1,jj) - u_oce(ji-1,jj) 
    219168               zt12  = rhoco * zu_ice * zu_ice 
    220169 
    221                zv_ice   = v_ice(ji,jj) - v_io(ji,jj) 
     170               zv_ice   = v_ice(ji,jj) - v_oce(ji,jj) 
    222171               zt21  = rhoco * zv_ice * zv_ice 
    223172 
    224                zv_ice   = v_ice(ji,jj-1) - v_io(ji,jj-1) 
     173               zv_ice   = v_ice(ji,jj-1) - v_oce(ji,jj-1) 
    225174               zt22  = rhoco * zv_ice * zv_ice 
    226                ztair2 = ( ( gtaux(ji,jj) + gtaux(ji-1,jj) ) / 2. )**2 + & 
    227                         ( ( gtauy(ji,jj) + gtauy(ji,jj-1) ) / 2. )**2 
     175               ztair2 = ( ( utaui_ice(ji,jj) + utaui_ice(ji-1,jj) ) / 2. )**2 + & 
     176                        ( ( vtaui_ice(ji,jj) + vtaui_ice(ji,jj-1) ) / 2. )**2 
    228177 
    229178               ! should not be weighted 
     
    241190          DO jj = 2, jpjm1 
    242191             DO ji = fs_2, fs_jpim1 
    243 #if defined key_coupled && defined key_lim_cp1 
    244                 tio_u(ji,jj) = - (  gtaux(ji  ,jj  ) + gtaux(ji-1,jj  )       & 
    245                    &              + gtaux(ji-1,jj-1) + gtaux(ji  ,jj-1) ) / ( 4 * rau0 ) 
    246  
    247                 tio_v(ji,jj) = - (  gtauy(ji  ,jj )  + gtauy(ji-1,jj  )       & 
    248                    &              + gtauy(ji-1,jj-1) + gtauy(ji  ,jj-1) ) / ( 4 * rau0 ) 
    249 #else 
    250                 tio_u(ji,jj) = - gtaux(ji,jj) / cai * cao / rau0 
    251                 tio_v(ji,jj) = - gtauy(ji,jj) / cai * cao / rau0  
    252 #endif 
    253                 ztair2 = ( ( gtaux(ji,jj) + gtaux(ji-1,jj) ) / 2. )**2 + & 
    254                          ( ( gtauy(ji,jj) + gtauy(ji,jj-1) ) / 2. )**2 
     192                ztair2 = ( ( utaui_ice(ji,jj) + utaui_ice(ji-1,jj) ) / 2. )**2 + & 
     193                         ( ( vtaui_ice(ji,jj) + vtaui_ice(ji,jj-1) ) / 2. )**2 
    255194                zustm        = SQRT( ztair2  ) 
    256195 
     
    262201 
    263202      CALL lbc_lnk( ust2s, 'T',  1. )   ! T-point 
    264       CALL lbc_lnk( tio_u, 'U', -1. )   ! I-point (i.e. ice U-V point) 
    265       CALL lbc_lnk( tio_v, 'V', -1. )   ! I-point (i.e. ice U-V point) 
    266203 
    267204      IF(ln_ctl) THEN   ! Control print 
     
    269206         CALL prt_ctl_info(' - Cell values : ') 
    270207         CALL prt_ctl_info('   ~~~~~~~~~~~~~ ') 
    271          CALL prt_ctl(tab2d_1=tio_u     , clinfo1=' lim_dyn  : tio_u     :', tab2d_2=tio_v , clinfo2=' tio_v :') 
    272208         CALL prt_ctl(tab2d_1=ust2s     , clinfo1=' lim_dyn  : ust2s     :') 
    273209         CALL prt_ctl(tab2d_1=divu_i    , clinfo1=' lim_dyn  : divu_i    :') 
  • trunk/NEMO/LIM_SRC_3/limhdf.F90

    r868 r888  
    3434   !!---------------------------------------------------------------------- 
    3535   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
    36    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limhdf.F90,v 1.5 2005/03/27 18:34:41 opalod Exp $  
     36   !! $ Id: $ 
    3737   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    3838   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_3/limistate.F90

    r869 r888  
    1616   USE oce             ! dynamics and tracers variables 
    1717   USE dom_oce 
     18   USE sbc_oce         ! Surface boundary condition: ocean fields 
    1819   USE par_ice         ! ice parameters 
    1920   USE ice_oce         ! ice variables 
     
    5152   !!---------------------------------------------------------------------- 
    5253   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005) 
    53    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limistate.F90,v 1.3 2005/03/27 18:34:41 opalod Exp $ 
    54    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     54   !! $ Id: $ 
     55   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5556   !!---------------------------------------------------------------------- 
    5657 
     
    9394      CALL lim_istate_init     !  reading the initials parameters of the ice 
    9495 
    95       !-- Initialisation of sst,sss,u,v do i=1,jpi 
    96       u_io(:,:)  = 0.e0       ! ice velocity in x direction 
    97       v_io(:,:)  = 0.e0       ! ice velocity in y direction 
    98  
    9996      ! Initialisation at tn or -2 if ice 
    10097      DO jj = 1, jpj 
     
    104101         END DO 
    105102      END DO 
    106  
    107       u_io  (:,:) = 0. 
    108       v_io  (:,:) = 0. 
    109       sst_io(:,:) = ( nfice - 1 ) * ( tn(:,:,1) + rt0 )   ! use the ocean initial values 
    110       sss_io(:,:) = ( nfice - 1 ) *   sn(:,:,1)           ! tricky trick *(nfice-1) ! 
    111103 
    112104      !-------------------------------------------------------------------- 
     
    280272                  !--------------- 
    281273                  sm_i(ji,jj,jl)   = zidto * sinn  + ( 1.0 - zidto ) * 0.1 
    282                   smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_io(ji,jj) ) * v_i(ji,jj,jl) 
     274                  smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) 
    283275 
    284276                  !---------- 
     
    405397 
    406398                  sm_i(ji,jj,jl)   = zidto * sins  + ( 1.0 - zidto ) * 0.1 
    407                   smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_io(ji,jj) ) * v_i(ji,jj,jl) 
     399                  smv_i(ji,jj,jl)  = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) 
    408400 
    409401                  !---------- 
     
    538530 
    539531      CALL lbc_lnk( fsbbq  , 'T', 1. ) 
    540       CALL lbc_lnk( sss_io , 'T', 1. ) 
    541532 
    542533   END SUBROUTINE lim_istate 
  • trunk/NEMO/LIM_SRC_3/limitd_me.F90

    r869 r888  
    2020   USE phycst           ! physical constants (ocean directory)  
    2121   USE ice_oce          ! ice variables 
     22   USE sbc_oce          ! Surface boundary condition: ocean fields 
    2223   USE thd_ice 
    2324   USE limistate 
     
    743744      ! Temporal smoothing 
    744745      !-------------------- 
    745       IF ( numit .EQ. nit000 + nfice - 1 ) THEN 
     746      IF ( numit .EQ. nit000 + nn_fsbc - 1 ) THEN 
    746747         strp1(:,:) = 0.0             
    747748         strp2(:,:) = 0.0             
     
    11941195      IF ( con_i ) THEN 
    11951196         CALL lim_column_sum (jpl,   v_i, vice_init ) 
    1196          WRITE(numout,*) ' vice_init  : ', vice_init(jiindex,jjindex) 
     1197         WRITE(numout,*) ' vice_init  : ', vice_init(jiindx,jjindx) 
    11971198         CALL lim_column_sum_energy (jpl, nlay_i,  e_i, eice_init ) 
    1198          WRITE(numout,*) ' eice_init  : ', eice_init(jiindex,jjindex) 
     1199         WRITE(numout,*) ' eice_init  : ', eice_init(jiindx,jjindx) 
    11991200      ENDIF 
    12001201 
     
    13631364            ! Salinity 
    13641365            !------------- 
    1365             smsw(ji,jj)       = sss_io(ji,jj) * vsw(ji,jj) * ridge_por  
     1366            smsw(ji,jj)       = sss_m(ji,jj) * vsw(ji,jj) * ridge_por  
    13661367 
    13671368            ! salinity of new ridge 
     
    14471448                                        - eirft(ji,jj,jk) 
    14481449            ! sea water heat content 
    1449             ztmelts          = - tmut * sss_io(ji,jj) + rtt 
     1450            ztmelts          = - tmut * sss_m(ji,jj) + rtt 
    14501451            ! heat content per unit volume 
    1451             zdummy0          = - rcp * ( sst_io(ji,jj) - rtt ) * vsw(ji,jj) 
     1452            zdummy0          = - rcp * ( sst_m(ji,jj) + rt0 - rtt ) * vsw(ji,jj) 
    14521453 
    14531454            ! corrected sea water salinity 
     
    16161617         fieldid = ' v_i : limitd_me ' 
    16171618         CALL lim_cons_check (vice_init, vice_final, 1.0e-6, fieldid)  
    1618          WRITE(numout,*) ' vice_init  : ', vice_init(jiindex,jjindex) 
    1619          WRITE(numout,*) ' vice_final : ', vice_final(jiindex,jjindex) 
     1619         WRITE(numout,*) ' vice_init  : ', vice_init(jiindx,jjindx) 
     1620         WRITE(numout,*) ' vice_final : ', vice_final(jiindx,jjindx) 
    16201621 
    16211622         CALL lim_column_sum_energy (jpl, nlay_i,  e_i, eice_final ) 
    16221623         fieldid = ' e_i : limitd_me ' 
    16231624         CALL lim_cons_check (eice_init, eice_final, 1.0e-2, fieldid)  
    1624          WRITE(numout,*) ' eice_init  : ', eice_init(jiindex,jjindex) 
    1625          WRITE(numout,*) ' eice_final : ', eice_final(jiindex,jjindex) 
     1625         WRITE(numout,*) ' eice_init  : ', eice_init(jiindx,jjindx) 
     1626         WRITE(numout,*) ' eice_final : ', eice_final(jiindx,jjindx) 
    16261627      ENDIF 
    16271628 
     
    18391840!           fresh_hist(i,j) = fresh_hist(i,j) + xtmp 
    18401841 
    1841 !           fsalt_res(ji,jj)  = fsalt_res(ji,jj) + ( sss_io(ji,jj)                  ) * &  
     1842!           fsalt_res(ji,jj)  = fsalt_res(ji,jj) + ( sss_m(ji,jj)                  ) * &  
    18421843!                               rhosn * v_s(ji,jj,jl) / rdt_ice 
    18431844 
    1844 !           fsalt_res(ji,jj)  = fsalt_res(ji,jj) + ( sss_io(ji,jj) - sm_i(ji,jj,jl) ) * &  
     1845!           fsalt_res(ji,jj)  = fsalt_res(ji,jj) + ( sss_m(ji,jj) - sm_i(ji,jj,jl) ) * &  
    18451846!                               rhoic * v_i(ji,jj,jl) / rdt_ice 
    18461847 
    1847 !           fsalt(i,j)      = fsalt(i,j)      + xtmp 
     1848!           emps(i,j)      = emps(i,j)      + xtmp 
    18481849!           fsalt_hist(i,j) = fsalt_hist(i,j) + xtmp 
    18491850 
  • trunk/NEMO/LIM_SRC_3/limmsh.F90

    r869 r888  
    2525   !!---------------------------------------------------------------------- 
    2626   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
    27    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limmsh.F90,v 1.5 2005/03/27 18:34:42 opalod Exp $  
     27   !! $ Id: $ 
    2828   !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    2929   !!---------------------------------------------------------------------- 
  • trunk/NEMO/LIM_SRC_3/limrhg.F90

    r869 r888  
    1616   USE dom_oce 
    1717   USE dom_ice 
     18   USE sbc_ice         ! Surface boundary condition: ice fields 
    1819   USE ice 
    1920   USE iceini 
     
    4041   !!---------------------------------------------------------------------- 
    4142   !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008)  
    42    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limrhg.F90,v 1.5 2005/03/27 18:34:42 opalod Exp $  
    43    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     43   !! $ Id: $ 
     44   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4445   !!---------------------------------------------------------------------- 
    4546 
     
    268269                                          / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 
    269270            ! 
    270             u_oce1(ji,jj)  = u_io(ji,jj) 
    271             v_oce2(ji,jj)  = v_io(ji,jj) 
     271            u_oce1(ji,jj)  = u_oce(ji,jj) 
     272            v_oce2(ji,jj)  = v_oce(ji,jj) 
    272273 
    273274            ! Ocean has no slip boundary condition 
    274             v_oce1(ji,jj)  = 0.5*( (v_io(ji,jj)+v_io(ji,jj-1))*e1t(ji,jj)    & 
    275                 &                 +(v_io(ji+1,jj)+v_io(ji+1,jj-1))*e1t(ji+1,jj)) & 
     275            v_oce1(ji,jj)  = 0.5*( (v_oce(ji,jj)+v_oce(ji,jj-1))*e1t(ji,jj)    & 
     276                &                 +(v_oce(ji+1,jj)+v_oce(ji+1,jj-1))*e1t(ji+1,jj)) & 
    276277                &               /(e1t(ji+1,jj)+e1t(ji,jj)) * tmu(ji,jj)   
    277278 
    278             u_oce2(ji,jj)  = 0.5*((u_io(ji,jj)+u_io(ji-1,jj))*e2t(ji,jj)     & 
    279                 &                 +(u_io(ji,jj+1)+u_io(ji-1,jj+1))*e2t(ji,jj+1)) & 
     279            u_oce2(ji,jj)  = 0.5*((u_oce(ji,jj)+u_oce(ji-1,jj))*e2t(ji,jj)     & 
     280                &                 +(u_oce(ji,jj+1)+u_oce(ji-1,jj+1))*e2t(ji,jj+1)) & 
    280281                &                / (e2t(ji,jj+1)+e2t(ji,jj)) * tmv(ji,jj) 
    281282 
    282283            ! Wind stress. 
    283             ztagnx = ( 1. - zfrld1(ji,jj) ) * gtaux(ji,jj) 
    284             ztagny = ( 1. - zfrld2(ji,jj) ) * gtauy(ji,jj) 
     284            ztagnx = ( 1. - zfrld1(ji,jj) ) * utaui_ice(ji,jj) 
     285            ztagny = ( 1. - zfrld2(ji,jj) ) * vtaui_ice(ji,jj) 
    285286 
    286287            ! Computation of the velocity field taking into account the ice internal interaction. 
     
    621622            zdummy = zindb * vt_i(ji,jj) / MAX(at_i(ji,jj) , 1.0e-06 ) 
    622623            IF ( zdummy .LE. 5.0e-2 ) THEN 
    623                u_ice(ji,jj) = u_io(ji,jj) 
    624                v_ice(ji,jj) = v_io(ji,jj) 
     624               u_ice(ji,jj) = u_oce(ji,jj) 
     625               v_ice(ji,jj) = v_oce(ji,jj) 
    625626            ENDIF ! zdummy 
    626627         END DO 
  • trunk/NEMO/LIM_SRC_3/limrst.F90

    r838 r888  
    1818   USE dom_oce 
    1919   USE ice_oce         ! ice variables 
     20   USE sbc_oce         ! Surface boundary condition: ocean fields 
     21   USE sbc_ice         ! Surface boundary condition: ice fields 
    2022   USE daymod 
    2123   USE iom 
     
    3436   !!---------------------------------------------------------------------- 
    3537   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005) 
    36    !! $Id:$ 
     38   !! $ Id: $ 
    3739   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3840   !!---------------------------------------------------------------------- 
     
    5557       
    5658      ! to get better performances with NetCDF format: 
    57       ! we open and define the ice restart file one ice time step before writing the data (-> at nitrst - 2*nfice + 1) 
    58       ! except if we write ice restart files every ice time step or if an ice restart file was writen at nitend - 2*nfice + 1 
    59       IF( kt == nitrst - 2*nfice + 1 .OR. nstock == nfice .OR. ( kt == nitend - nfice + 1 .AND. .NOT. lrst_ice ) ) THEN 
     59      ! we open and define the ice restart file one ice time step before writing the data (-> at nitrst - 2*nn_fsbc + 1) 
     60      ! except if we write ice restart files every ice time step or if an ice restart file was writen at nitend - 2*nn_fsbc + 1 
     61      IF( kt == nitrst - 2*nn_fsbc + 1 .OR. nstock == nn_fsbc .OR. ( kt == nitend - nn_fsbc + 1 .AND. .NOT. lrst_ice ) ) THEN 
    6062         ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
    6163         IF( nitrst > 99999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
     
    7072            CASE DEFAULT         ;   WRITE(numout,*) '             open ice restart NetCDF file: '//clname 
    7173            END SELECT 
    72             IF( kt == nitrst - 2*nfice + 1 ) THEN    
    73                WRITE(numout,*)         '             kt = nitrst - 2*nfice + 1 = ', kt,' date= ', ndastp 
    74             ELSE   ;   WRITE(numout,*) '             kt = '                       , kt,' date= ', ndastp 
     74            IF( kt == nitrst - 2*nn_fsbc + 1 ) THEN    
     75               WRITE(numout,*)         '             kt = nitrst - 2*nn_fsbc + 1 = ', kt,' date= ', ndastp 
     76            ELSE   ;   WRITE(numout,*) '             kt = '                         , kt,' date= ', ndastp 
    7577            ENDIF 
    7678         ENDIF 
     
    100102      !!---------------------------------------------------------------------- 
    101103    
    102       iter = kt + nfice - 1   ! ice restarts are written at kt == nitrst - nfice + 1 
     104      iter = kt + nn_fsbc - 1   ! ice restarts are written at kt == nitrst - nn_fsbc + 1 
    103105 
    104106      IF( iter == nitrst ) THEN 
     
    111113      ! ------------------  
    112114      !                                                                        ! calendar control 
    113       CALL iom_rstput( iter, nitrst, numriw, 'nfice' , REAL( nfice, wp) )      ! time-step  
    114       CALL iom_rstput( iter, nitrst, numriw, 'kt_ice', REAL( iter , wp) )      ! date 
     115      CALL iom_rstput( iter, nitrst, numriw, 'nn_fsbc', REAL( nn_fsbc, wp) )      ! time-step  
     116      CALL iom_rstput( iter, nitrst, numriw, 'kt_ice' , REAL( iter , wp) )        ! date 
    115117 
    116118      ! Prognostic variables  
     
    158160      CALL iom_rstput( iter, nitrst, numriw, 'u_ice'     , u_ice      ) 
    159161      CALL iom_rstput( iter, nitrst, numriw, 'v_ice'     , v_ice      ) 
    160       CALL iom_rstput( iter, nitrst, numriw, 'gtaux'     , gtaux      ) 
    161       CALL iom_rstput( iter, nitrst, numriw, 'gtauy'     , gtauy      ) 
     162      CALL iom_rstput( iter, nitrst, numriw, 'utaui_ice' , utaui_ice  ) 
     163      CALL iom_rstput( iter, nitrst, numriw, 'vtaui_ice' , vtaui_ice  ) 
    162164      CALL iom_rstput( iter, nitrst, numriw, 'fsbbq'     , fsbbq      ) 
    163165      CALL iom_rstput( iter, nitrst, numriw, 'stress1_i' , stress1_i  ) 
     
    299301               WRITE(numout,*) ' ~~~ Arctic' 
    300302    
    301                ji = jiindex 
    302                jj = jjindex 
     303               ji = jiindx 
     304               jj = jjindx 
    303305    
    304306               WRITE(numout,*) ' ji, jj ', ji, jj 
     
    387389      !!---------------------------------------------------------------------- 
    388390      ! Local variables 
    389       INTEGER :: ji, jj, jk, jl, index 
     391      INTEGER :: ji, jj, jk, jl, indx 
    390392      REAL(wp) ::   zfice, ziter 
    391393      REAL(wp) :: & !parameters for the salinity profile 
     
    405407      CALL iom_open ( 'restart_ice_in', numrir, kiolib = jprstlib ) 
    406408 
    407       CALL iom_get( numrir, 'nfice' , zfice ) 
    408       CALL iom_get( numrir, 'kt_ice', ziter )     
     409      CALL iom_get( numrir, 'nn_fsbc', zfice ) 
     410      CALL iom_get( numrir, 'kt_ice' , ziter )     
    409411      IF(lwp) WRITE(numout,*) '   read ice restart file at time step    : ', ziter 
    410412      IF(lwp) WRITE(numout,*) '   in any case we force it to nit000 - 1 : ', nit000 - 1 
     
    416418         &                   '   verify the file or rerun with the value 0 for the',        & 
    417419         &                   '   control of time parameter  nrstdt' ) 
    418       IF( INT(zfice) /= nfice          .AND. ABS( nrstdt ) == 1 )   & 
    419          &     CALL ctl_stop( 'lim_rst_read ===>>>> : problem with nfice in ice restart',  & 
    420          &                   '   verify the file or rerun with the value 0 for the',        & 
     420      IF( INT(zfice) /= nn_fsbc          .AND. ABS( nrstdt ) == 1 )   & 
     421         &     CALL ctl_stop( 'lim_rst_read ===>>>> : problem with nn_fsbc in ice restart',  & 
     422         &                   '   verify the file or rerun with the value 0 for the',         & 
    421423         &                   '   control of time parameter  nrstdt' ) 
    422424 
     
    512514      CALL iom_get( numrir, jpdom_autoglo, 'u_ice'     , u_ice      ) 
    513515      CALL iom_get( numrir, jpdom_autoglo, 'v_ice'     , v_ice      ) 
    514       CALL iom_get( numrir, jpdom_autoglo, 'gtaux'     , gtaux      ) 
    515       CALL iom_get( numrir, jpdom_autoglo, 'gtauy'     , gtauy      ) 
     516      CALL iom_get( numrir, jpdom_autoglo, 'utaui_ice' , utaui_ice  ) 
     517      CALL iom_get( numrir, jpdom_autoglo, 'vtaui_ice' , vtaui_ice  ) 
    516518      CALL iom_get( numrir, jpdom_autoglo, 'fsbbq'     , fsbbq      ) 
    517519      CALL iom_get( numrir, jpdom_autoglo, 'stress1_i' , stress1_i  ) 
     
    650652               WRITE(numout,*) ' ~~~ Arctic' 
    651653    
    652                index = 1 
     654               indx = 1 
    653655               ji = 24 
    654656               jj = 24 
  • trunk/NEMO/LIM_SRC_3/limtab.F90

    r834 r888  
    2323   !!---------------------------------------------------------------------- 
    2424   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
    25    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limtab.F90,v 1.2 2005/03/27 18:34:42 opalod Exp $  
    26    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     25   !! $ Id: $ 
     26   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    2727   !!---------------------------------------------------------------------- 
    2828CONTAINS 
  • trunk/NEMO/LIM_SRC_3/limthd.F90

    r869 r888  
    1818   USE ice             ! LIM sea-ice variables 
    1919   USE ice_oce         ! sea-ice/ocean variables 
    20    USE flx_oce         ! sea-ice/ocean forcings variables  
     20   USE sbc_oce         ! Surface boundary condition: ocean fields 
     21   USE sbc_ice         ! Surface boundary condition: ice fields 
    2122   USE thd_ice         ! LIM thermodynamic sea-ice variables 
    2223   USE dom_ice         ! LIM sea-ice domain 
     
    5253   !!---------------------------------------------------------------------- 
    5354   !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005) 
    54    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limthd.F90,v 1.6 2005/03/27 18:34:42 opalod Exp $ 
    55    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     55   !! $ Id: $ 
     56   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5657   !!---------------------------------------------------------------------- 
    5758 
     
    8485      !!--------------------------------------------------------------------- 
    8586      !! * Local variables 
    86       INTEGER  ::  ji, jj, jk, jl,  & 
    87                    zji  , zjj,      &   ! dummy loop indices 
    88                    nbpb ,           &   ! nb of icy pts for thermo. cal. 
    89                    index 
     87      INTEGER  ::  ji, jj, jk, jl, nbpb   ! nb of icy pts for thermo. cal. 
    9088 
    9189      REAL(wp) ::  & 
     
    211209 
    212210            ! here the drag will depend on ice thickness and type (0.006) 
    213             fdtcn(ji,jj)  = zindb * rau0 * rcp * 0.006  * zfric_u * ( sst_io(ji,jj) - t_bo(ji,jj) )  
     211            fdtcn(ji,jj)  = zindb * rau0 * rcp * 0.006  * zfric_u * ( (sst_m(ji,jj) + rt0) - t_bo(ji,jj) )  
    214212            ! also category dependent 
    215213!           !-- Energy from the turbulent oceanic heat flux heat flux coming in the lead  
     
    220218!           !-- Lead heat budget (part 1, next one is in limthd_dh 
    221219!           !-- qldif -- (or qldif_1d in 1d routines) 
    222             zfontn         = sprecip(ji,jj) * lfus              !   energy of melting 
    223             zfnsol         = qnsr_oce(ji,jj)  ! total non solar flux 
    224             qldif(ji,jj)   = tms(ji,jj) * ( qsr_oce(ji,jj)                          & 
     220            zfontn         = sprecip(ji,jj) * lfus              ! energy of melting 
     221            zfnsol         = qns(ji,jj)                         ! total non solar flux 
     222            qldif(ji,jj)   = tms(ji,jj) * ( qsr(ji,jj)                              & 
    225223               &                               + zfnsol + fdtcn(ji,jj) - zfontn     & 
    226224               &                               + ( 1.0 - zindb ) * fsbbq(ji,jj) )   & 
     
    242240            ! Energy needed to bring ocean surface layer until its freezing 
    243241            ! qcmif, limflx 
    244             qcmif  (ji,jj) =  rau0 * rcp * fse3t(ji,jj,1) * ( t_bo(ji,jj) - sst_io(ji,jj) ) * ( 1. - zinda ) 
     242            qcmif  (ji,jj) =  rau0 * rcp * fse3t(ji,jj,1) * ( t_bo(ji,jj) - (sst_m(ji,jj) + rt0) ) * ( 1. - zinda ) 
    245243 
    246244            !  calculate oceanic heat flux (limthd_dh) 
     
    271269               ENDIF 
    272270               ! debug point to follow 
    273                IF ( (ji.eq.jiindex).AND.(jj.eq.jjindex) ) THEN 
     271               IF ( (ji.eq.jiindx).AND.(jj.eq.jjindx) ) THEN 
    274272                   jiindex_1d = nbpb 
    275273               ENDIF 
     
    310308            CALL tab_2d_1d( nbpb, fr1_i0_1d  (1:nbpb)     , fr1_i0     , jpi, jpj, npb(1:nbpb) ) 
    311309            CALL tab_2d_1d( nbpb, fr2_i0_1d  (1:nbpb)     , fr2_i0     , jpi, jpj, npb(1:nbpb) ) 
    312             CALL tab_2d_1d( nbpb, qnsr_ice_1d(1:nbpb)     , qnsr_ice(:,:,jl)  , jpi, jpj, npb(1:nbpb) ) 
     310            CALL tab_2d_1d( nbpb, qnsr_ice_1d(1:nbpb)     , qns_ice(:,:,jl) , jpi, jpj, npb(1:nbpb) ) 
    313311 
    314312#if ! defined key_coupled 
     
    360358 
    361359                                          !---------------------------------! 
    362             CALL lim_thd_sal(1,nbpb,jl)   ! Ice salinity computation        ! 
     360            CALL lim_thd_sal(1,nbpb)      ! Ice salinity computation        ! 
    363361                                          !---------------------------------! 
    364362 
     
    415413            CALL tab_1d_2d( nbpb, s_i_newice , npb, s_i_new  (1:nbpb)      , jpi, jpj ) 
    416414            CALL tab_1d_2d( nbpb, izero(:,:,jl) , npb, i0    (1:nbpb)      , jpi, jpj ) 
    417             CALL tab_1d_2d( nbpb, qnsr_ice(:,:,jl), npb, qnsr_ice_1d(1:nbpb), jpi, jpj) 
     415            CALL tab_1d_2d( nbpb, qns_ice(:,:,jl), npb, qnsr_ice_1d(1:nbpb), jpi, jpj) 
    418416            !+++++ 
    419417 
     
    543541 
    544542      INTEGER  :: & 
    545          ji,jj,jk            ! loop indices 
     543         ji,jk               ! loop indices 
    546544 
    547545      !!----------------------------------------------------------------------- 
     
    598596                                       !  is violated 
    599597      INTEGER  :: & 
    600          ji,jj,jk,                  &  !: loop indices 
     598         ji,jk,                     &  !: loop indices 
    601599         zji, zjj 
    602600      !!--------------------------------------------------------------------- 
     
    726724         WRITE(numout,*) ' foc        : ', fbif_1d(ji) 
    727725         WRITE(numout,*) ' fstroc     : ', fstroc   (zji,zjj,jl) 
    728          WRITE(numout,*) ' i0        : ', i0(ji) 
    729          WRITE(numout,*) ' fsolar     : ', (1.0-i0(ji))*qsr_ice_1d(ji) 
    730          WRITE(numout,*) ' fnsolar    : ', qnsr_ice_1d(ji) 
     726         WRITE(numout,*) ' i0         : ', i0(ji) 
     727         WRITE(numout,*) ' qsr_ice    : ', (1.0-i0(ji))*qsr_ice_1d(ji) 
     728         WRITE(numout,*) ' qns_ice    : ', qnsr_ice_1d(ji) 
    731729         WRITE(numout,*) ' Conduction fluxes : ' 
    732730         WRITE(numout,*) ' fc_s      : ', fc_s(ji,0:nlay_s) 
     
    778776         numce                         !: number of points for which conservation 
    779777                                       !  is violated 
    780       INTEGER  :: & 
    781          ji,jj,jk,                  &  !: loop indices 
    782          zji, zjj 
    783  
     778      INTEGER  ::  ji, zji, zjj        ! loop indices 
    784779      !!--------------------------------------------------------------------- 
    785780 
  • trunk/NEMO/LIM_SRC_3/limthd_dh.F90

    r869 r888  
    1616   USE phycst           ! physical constants (OCE directory)  
    1717   USE ice_oce          ! ice variables 
     18   USE sbc_oce          ! Surface boundary condition: ocean fields 
    1819   USE thd_ice 
    1920   USE iceini 
     
    338339            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    339340            zfsalt_melt(ji)     = zfsalt_melt(ji) +                           & 
    340                                   ( sss_io(zji,zjj) - sm_i_b(ji)   ) *        & 
     341                                  ( sss_m(zji,zjj) - sm_i_b(ji)   ) *         & 
    341342                                  a_i_b(ji) *                                 & 
    342343                                  MIN( zdeltah(ji,jk) , 0.0 ) * rhoic / rdt_ice  
     
    368369            WRITE(numout,*) ' qlbbq_1d: ', qlbbq_1d(ji) 
    369370            WRITE(numout,*) ' s_i_new : ', s_i_new(ji) 
    370             WRITE(numout,*) ' sss_io  : ', sss_io(zji,zjj) 
     371            WRITE(numout,*) ' sss_m   : ', sss_m(zji,zjj) 
    371372         ENDIF 
    372373 
     
    494495                           zswi2  * 0.26 /  & 
    495496                           ( 0.26 + 0.74 * EXP ( - 724300.0 * zgrr ) )  
    496                   zds         = zfracs*sss_io(zji,zjj) - s_i_new(ji) 
    497                   s_i_new(ji) = zfracs * sss_io(zji,zjj) 
     497                  zds         = zfracs*sss_m(zji,zjj) - s_i_new(ji) 
     498                  s_i_new(ji) = zfracs * sss_m(zji,zjj) 
    498499               ENDIF ! fc_bo_i 
    499500            END DO ! ji 
     
    567568                  zjj             = ( npb(ji) - 1 ) / jpi + 1 
    568569                  zfsalt_melt(ji) = zfsalt_melt(ji) +                         & 
    569                                    ( sss_io(zji,zjj) - sm_i_b(ji)   ) *       & 
     570                                   ( sss_m(zji,zjj) - sm_i_b(ji)   ) *        & 
    570571                                   a_i_b(ji) * & 
    571572                                   MIN( zdeltah(ji,jk) , 0.0 ) * rhoic / rdt_ice  
     
    596597                WRITE(numout,*) ' qlbbq_1d: ', qlbbq_1d(ji) 
    597598                WRITE(numout,*) ' s_i_new : ', s_i_new(ji) 
    598                 WRITE(numout,*) ' sss_io  : ', sss_io(zji,zjj) 
     599                WRITE(numout,*) ' sss_m   : ', sss_m(zji,zjj) 
    599600                WRITE(numout,*) ' dh_i_bott : ', dh_i_bott(ji) 
    600601                WRITE(numout,*) ' innermelt : ', innermelt(ji) 
     
    701702         fseqv_1d(ji)  = fseqv_1d(ji) + zihgnew * zfsalt_melt(ji) +           & 
    702703                          (1.0 - zihgnew) * rdmicif_1d(ji) *                  & 
    703                           ( sss_io(zji,zjj) - sm_i_b(ji) ) / rdt_ice 
     704                          ( sss_m(zji,zjj) - sm_i_b(ji) ) / rdt_ice 
    704705         ! new lines 
    705706         IF ( num_sal .EQ. 4 ) & 
    706707         fseqv_1d(ji)  = fseqv_1d(ji) + zihgnew * zfsalt_melt(ji) +           & 
    707708                          (1.0 - zihgnew) * rdmicif_1d(ji) *                  & 
    708                           ( sss_io(zji,zjj) - bulk_sal ) / rdt_ice 
     709                          ( sss_m(zji,zjj) - bulk_sal ) / rdt_ice 
    709710         ! Heat flux 
    710711         ! excessive bottom ablation energy (fsup) - 0 except if jpl = 1 
     
    762763                             *(ht_s_b(ji)-zhnnew)*rhosn 
    763764 
    764 #if defined key_lim_fdd  
    765 !(presently Activated)  
    766765         rdmicif_1d(ji) = rdmicif_1d(ji) + a_i_b(ji) &  
    767766                                         * ( zhgnew(ji) - ht_i_b(ji) )*rhoic  
     
    775774 
    776775         zsm_snowice  = ( rhoic - rhosn ) / rhoic *            & 
    777                         sss_io(zji,zjj)  
     776                        sss_m(zji,zjj)  
    778777 
    779778         IF ( num_sal .NE. 2 ) zsm_snowice = sm_i_b(ji) 
     
    781780         IF ( num_sal .NE. 4 ) & 
    782781         fseqv_1d(ji)   = fseqv_1d(ji)   + & 
    783                           ( sss_io(zji,zjj) - zsm_snowice ) * & 
     782                          ( sss_m(zji,zjj) - zsm_snowice ) * & 
    784783                            a_i_b(ji)   * & 
    785784                          ( zhgnew(ji) - ht_i_b(ji) ) * rhoic / rdt_ice 
     
    787786         IF ( num_sal .EQ. 4 ) & 
    788787         fseqv_1d(ji)   = fseqv_1d(ji)   + & 
    789                           ( sss_io(zji,zjj) - bulk_sal    ) * & 
     788                          ( sss_m(zji,zjj) - bulk_sal    ) * & 
    790789                            a_i_b(ji)   * & 
    791790                          ( zhgnew(ji) - ht_i_b(ji) ) * rhoic / rdt_ice 
     
    801800                            - sm_i_b(ji) ) * isnowic      
    802801 
    803 #else 
    804          rdmicif_1d(ji) = rdmicif_1d(ji) + a_i_b(ji) & 
    805                                          * ( zhgnew(ji) - ht_i_b(ji) ) * rhoic & 
    806                                          + ( zhnnew - ht_s_b(ji) ) * rhosn ) 
    807 #endif 
    808802!  Actualize new snow and ice thickness. 
    809803         ht_s_b(ji)  = zhnnew 
  • trunk/NEMO/LIM_SRC_3/limthd_lac.F90

    r865 r888  
    11MODULE limthd_lac 
    2 #if defined key_lim3 
    32   !!---------------------------------------------------------------------- 
    43   !!   'key_lim3'                                      LIM3 sea-ice model 
     
    87   !!                lateral thermodynamic growth of the ice  
    98   !!====================================================================== 
    10  
     9#if defined key_lim3 
    1110   !!---------------------------------------------------------------------- 
    1211   !!   lim_lat_acr    : lateral accretion of ice 
     
    1716   USE phycst 
    1817   USE ice_oce         ! ice variables 
     18   USE sbc_oce         ! Surface boundary condition: ocean fields 
     19   USE sbc_ice         ! Surface boundary condition: ice fields 
    1920   USE thd_ice 
    2021   USE dom_ice 
     
    2324   USE iceini 
    2425   USE limtab 
    25    USE taumod 
    26    USE blk_oce 
    2726   USE limcons 
    2827      
     
    4645   !!---------------------------------------------------------------------- 
    4746   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
    48    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limthd_lac.F90,v 1.5 2005/03/27 18:34:42 opalod Exp $  
    49    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     47   !! $ Id: $ 
     48   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5049   !!---------------------------------------------------------------------- 
    5150 
     
    181180         vt_s_init, vt_s_final,   &  !  snow volume summed over categories 
    182181         et_i_init, et_i_final,   &  !  ice energy summed over categories 
    183          et_s_init, et_s_final       !  snow energy summed over categories 
     182         et_s_init                   !  snow energy summed over categories 
    184183 
    185184      REAL(wp) ::            & 
     
    267266            !------------- 
    268267            ! C-grid wind stress components 
    269             ztaux         = ( gtaux(ji-1,jj  ) * tmu(ji-1,jj  ) & 
    270                           +   gtaux(ji  ,jj  ) * tmu(ji  ,jj  ) ) / 2.0 
    271             ztauy         = ( gtauy(ji  ,jj-1) * tmv(ji  ,jj-1) & 
    272                           +   gtauy(ji  ,jj  ) * tmv(ji  ,jj  ) ) / 2.0 
     268            ztaux         = ( utaui_ice(ji-1,jj  ) * tmu(ji-1,jj  ) & 
     269                          +   utaui_ice(ji  ,jj  ) * tmu(ji  ,jj  ) ) / 2.0 
     270            ztauy         = ( vtaui_ice(ji  ,jj-1) * tmv(ji  ,jj-1) & 
     271                          +   vtaui_ice(ji  ,jj  ) * tmv(ji  ,jj  ) ) / 2.0 
    273272            ! Square root of wind stress 
    274273            ztenagm       =  SQRT( SQRT( ztaux * ztaux + ztauy * ztauy ) ) 
     
    343342               nbpac = nbpac + 1 
    344343               npac( nbpac ) = (jj - 1) * jpi + ji 
    345                IF ( (ji.eq.jiindex).AND.(jj.eq.jjindex) ) THEN 
     344               IF ( (ji.eq.jiindx).AND.(jj.eq.jjindx) ) THEN 
    346345                  jiindex_1d = nbpac 
    347346               ENDIF 
     
    418417              zji            =   MOD( npac(ji) - 1, jpi ) + 1 
    419418              zjj            =   ( npac(ji) - 1 ) / jpi + 1 
    420               zs_newice(ji)  =   MIN( 0.5*sss_io(zji,zjj) , zs_newice(ji) ) 
     419              zs_newice(ji)  =   MIN( 0.5*sss_m(zji,zjj) , zs_newice(ji) ) 
    421420           END DO ! jl 
    422421 
     
    476475              zjj            = ( npac(ji) - 1 ) / jpi + 1 
    477476              fseqv_1d(ji)   = fseqv_1d(ji) +                                     & 
    478                                ( sss_io(zji,zjj) - bulk_sal      ) * rhoic *      & 
     477                               ( sss_m(zji,zjj) - bulk_sal      ) * rhoic *       & 
    479478                               zv_newice(ji) / rdt_ice 
    480479           END DO 
     
    484483              zjj            = ( npac(ji) - 1 ) / jpi + 1 
    485484              fseqv_1d(ji)   = fseqv_1d(ji) +                                     & 
    486                                ( sss_io(zji,zjj) - zs_newice(ji) ) * rhoic *      & 
     485                               ( sss_m(zji,zjj) - zs_newice(ji) ) * rhoic *       & 
    487486                               zv_newice(ji) / rdt_ice 
    488487           END DO ! ji 
     
    617616        END DO 
    618617 
    619         WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindex, 1:jpl) 
     618        WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindx, 1:jpl) 
    620619        DO jl = ice_cat_bounds(jm,1), ice_cat_bounds(jm,2) 
    621620           DO ji = 1, nbpac 
     
    626625           END DO ! ji 
    627626        END DO ! jl 
    628         WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindex, 1:jpl) 
     627        WRITE(numout,*) ' zv_i_ac : ', zv_i_ac(jiindx, 1:jpl) 
    629628 
    630629        !--------------------------------- 
     
    796795!     CALL lim_cons_check (et_s_init, et_s_final, 1.0e-3, fieldid)  
    797796 
    798       WRITE(numout,*) ' vt_i_init : ', vt_i_init(jiindex,jjindex) 
    799       WRITE(numout,*) ' vt_i_final: ', vt_i_final(jiindex,jjindex) 
    800       WRITE(numout,*) ' et_i_init : ', et_i_init(jiindex,jjindex) 
    801       WRITE(numout,*) ' et_i_final: ', et_i_final(jiindex,jjindex) 
     797      WRITE(numout,*) ' vt_i_init : ', vt_i_init(jiindx,jjindx) 
     798      WRITE(numout,*) ' vt_i_final: ', vt_i_final(jiindx,jjindx) 
     799      WRITE(numout,*) ' et_i_init : ', et_i_init(jiindx,jjindx) 
     800      WRITE(numout,*) ' et_i_final: ', et_i_final(jiindx,jjindx) 
    802801 
    803802      ENDIF 
  • trunk/NEMO/LIM_SRC_3/limthd_sal.F90

    r869 r888  
    11MODULE limthd_sal 
    2 #if defined key_lim3 
    32   !!---------------------------------------------------------------------- 
    43   !!   'key_lim3'                                      LIM3 sea-ice model 
     
    98   !!                               the ice 
    109   !!====================================================================== 
    11  
     10#if defined key_lim3 
    1211   !!---------------------------------------------------------------------- 
    1312   !!   lim_thd_sal : salinity variations in the ice 
     
    1615   USE phycst           ! physical constants (ocean directory) 
    1716   USE ice_oce          ! ice variables 
     17   USE sbc_oce          ! Surface boundary condition: ocean fields 
    1818   USE thd_ice 
    1919   USE iceini 
     
    4040   CONTAINS 
    4141 
    42    SUBROUTINE lim_thd_sal(kideb,kiut,jl) 
     42   SUBROUTINE lim_thd_sal(kideb,kiut) 
    4343      !!------------------------------------------------------------------- 
    4444      !!                ***  ROUTINE lim_thd_sal  ***        
     
    7676      !! * Local variables 
    7777      INTEGER, INTENT(in) :: & 
    78          kideb, kiut, jl         !: thickness category index 
     78         kideb, kiut             !: thickness category index 
    7979 
    8080      INTEGER ::             & 
     
    318318            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    319319            fseqv_1d(ji) = fseqv_1d(ji)              + &  
    320                            ( sss_io(zji,zjj) - bulk_sal    ) * &  
     320                           ( sss_m(zji,zjj) - bulk_sal    ) * &  
    321321                           rhoic * a_i_b(ji) * & 
    322322                           MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 
     
    327327            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    328328            fseqv_1d(ji) = fseqv_1d(ji)              + &  
    329                            ( sss_io(zji,zjj) - s_i_new(ji) ) * &  
     329                           ( sss_m(zji,zjj) - s_i_new(ji) ) * &  
    330330                             rhoic * a_i_b(ji) * & 
    331331                             MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 
  • trunk/NEMO/LIM_SRC_3/limtrp.F90

    r868 r888  
    1717   USE in_out_manager  ! I/O manager 
    1818   USE ice_oce         ! ice variables 
     19   USE sbc_oce         ! Surface boundary condition: ocean fields 
    1920   USE dom_ice 
    2021   USE ice 
     
    5152   !!---------------------------------------------------------------------- 
    5253   !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
    53    !! $Header: /home/opalod/NEMOCVSROOT/NEMO/LIM_SRC/limtrp.F90,v 1.5 2005/03/27 18:34:42 opalod Exp $ 
    54    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     54   !! $ Id: $ 
     55   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5556   !!---------------------------------------------------------------------- 
    5657CONTAINS 
     
    519520 
    520521                  ! Ice salinity and age 
    521                   zsal            = MAX( MIN( (rhoic-rhosn)/rhoic*sss_io(ji,jj)  , & 
     522                  zsal            = MAX( MIN( (rhoic-rhosn)/rhoic*sss_m(ji,jj)  , & 
    522523                                            zusvoic * zs0sm(ji,jj,jl) ), s_i_min ) * & 
    523524                                            v_i(ji,jj,jl) 
  • trunk/NEMO/LIM_SRC_3/limupdate.F90

    r869 r888  
    2525   USE in_out_manager 
    2626   USE ice_oce         ! ice variables 
    27    USE flx_oce         ! forcings variables 
     27   USE sbc_oce         ! Surface boundary condition: ocean fields 
     28   USE sbc_ice         ! Surface boundary condition: ice fields 
    2829   USE dom_ice 
    2930   USE daymod 
    3031   USE phycst          ! Define parameters for the routines 
    31    USE taumod 
    3232   USE ice 
    3333   USE iceini 
    34    USE ocesbc 
    3534   USE lbclnk 
    3635   USE limdyn 
    3736   USE limtrp 
    3837   USE limthd 
    39    USE limflx 
     38   USE limsbc 
    4039   USE limdia 
    4140