New URL for NEMO forge!   http://forge.nemo-ocean.eu

Since March 2022 along with NEMO 4.2 release, the code development moved to a self-hosted GitLab.
This present forge is now archived and remained online for history.
Changeset 1855 – NEMO

Changeset 1855


Ignore:
Timestamp:
2010-04-30T17:49:04+02:00 (14 years ago)
Author:
gm
Message:

ticket:#665 style change only, with the suppression of thd_ice_2 (merged in ice_2)

Location:
branches/DEV_r1837_mass_heat_salt_fluxes/NEMO
Files:
1 deleted
21 edited

Legend:

Unmodified
Added
Removed

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/dom_ice_2.F90

    r1228 r1855  
    22   !!====================================================================== 
    33   !!                   ***  MODULE  dom_ice  *** 
    4    !! LIM 2.0 Sea Ice :   Domain  variables 
     4   !! LIM-2 Sea Ice :   Domain  variables 
    55   !!====================================================================== 
    6    !! History :   2.0  !  03-08  (C. Ethe)  Free form and module 
     6   !! History :   2.0  !  2003-08  (C. Ethe)  Free form and module 
    77   !!---------------------------------------------------------------------- 
    88#if defined key_lim2 
    99   !!---------------------------------------------------------------------- 
    10    !!   LIM 2.0, UCL-LOCEAN-IPSL (2005) 
    11    !! $Id$ 
    12    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     10   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1311   !!---------------------------------------------------------------------- 
    1412   USE par_ice_2 
    … …  
    1816 
    1917   LOGICAL, PUBLIC ::   l_jeq     = .TRUE.     !: Equator inside the domain flag 
     18   INTEGER, PUBLIC ::   njeq , njeqm1          !: j-index of the equator if it is inside the domain 
    2019 
    21    INTEGER, PUBLIC ::   njeq , njeqm1          !: j-index of the equator if it is inside the domain 
    22       !                                        !  (otherwise = jpj+10 (SH) or -10 (SH) ) 
     20   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)         ::   fs2cor, fcor     !: coriolis factor and coeficient 
     21   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)         ::   covrai           !: sine of geographic latitude 
     22   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)         ::   area             !: surface of grid cell  
     23   REAL(wp), PUBLIC, DIMENSION(jpi,jpj)         ::   tms   , tmu      !: temperature and velocity points masks 
     24   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2)     ::   wght             !: weight of the 4 neighbours to compute averages 
     25   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2)     ::   akappa, bkappa   !: first and third group of metric coefficients 
     26   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2,2,2) ::   alambd           !: second group of metric coefficients 
    2327 
    24    REAL(wp), PUBLIC, DIMENSION(jpi,jpj)         ::   fs2cor , fcor,   &  !: coriolis factor and coeficient 
    25       &                                              covrai ,         &  !: sine of geographic latitude 
    26       &                                              area   ,         &  !: surface of grid cell  
    27       &                                              tms    , tmu        !: temperature and velocity points masks 
    28    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2)     ::   wght   ,         &  !: weight of the 4 neighbours to compute averages 
    29       &                                              akappa , bkappa     !: first and third group of metric coefficients 
    30    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,2,2,2,2) ::   alambd   !: second group of metric coefficients 
     28#else 
     29   !!---------------------------------------------------------------------- 
     30   !!   Default option          Dummy module       NO LIM 2.0 sea-ice model 
     31   !!---------------------------------------------------------------------- 
     32#endif 
    3133 
     34   !!-------- ------------------------------------------------------------- 
     35   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
     36   !! $Id$ 
     37   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3238   !!====================================================================== 
    33 #endif 
    3439END MODULE dom_ice_2 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/ice_2.F90

    r1756 r1855  
    44   !! Sea Ice physics:  diagnostics variables of ice defined in memory 
    55   !!===================================================================== 
    6    !! History :  2.0  !  03-08  (C. Ethe)  F90: Free form and module 
     6   !! History :  2.0  !  2003-08 (C. Ethe)  F90: Free form and module 
     7   !!             -   !  2010-05 (G. Madec) suppression of thd_ice_2 
    78   !!---------------------------------------------------------------------- 
    89#if defined key_lim2 
    … …  
    1011   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1112   !!---------------------------------------------------------------------- 
    12    !!  LIM 2.0, UCL-LOCEAN-IPSL (2006) 
    13    !! $Id$ 
    14    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    15    !!---------------------------------------------------------------------- 
    16    !! * Modules used 
    1713   USE par_ice_2          ! LIM sea-ice parameters 
    1814 
    … …  
    2016   PRIVATE 
    2117 
    22    !!* Share parameters namelist (namicerun read in iceini) * 
    23    CHARACTER(len=32)     , PUBLIC ::   cn_icerst_in  = "restart_ice_in"   !: suffix of ice restart name (input) 
    24    CHARACTER(len=32)     , PUBLIC ::   cn_icerst_out = "restart_ice"      !: suffix of ice restart name (output) 
    25    LOGICAL               , PUBLIC ::   ln_limdyn     = .TRUE.             !: flag for ice dynamics (T) or not (F) 
    26    LOGICAL               , PUBLIC ::   ln_limdmp     = .FALSE.            !: Ice damping 
    27    REAL(wp)              , PUBLIC ::   hsndif        = 0.e0               !: computation of temp. in snow (0) or not (9999) 
    28    REAL(wp)              , PUBLIC ::   hicdif        = 0.e0               !: computation of temp. in ice (0) or not (9999) 
    29    REAL(wp), DIMENSION(2), PUBLIC ::   acrit = (/ 1.e-06 , 1.e-06 /)      !: minimum fraction for leads in  
    30    !                                                                      !: north and south hemisphere 
    31    !!* ice-dynamic namelist (namicedyn) * 
    32    INTEGER , PUBLIC ::   nbiter = 1         !: number of sub-time steps for relaxation 
    33    INTEGER , PUBLIC ::   nbitdr = 250       !: maximum number of iterations for relaxation 
    34    REAL(wp), PUBLIC ::   rdt_ice            !: ice time step 
    35    REAL(wp), PUBLIC ::   epsd   = 1.0e-20   !: tolerance parameter for dynamic 
    36    REAL(wp), PUBLIC ::   alpha  = 0.5       !: coefficient for semi-implicit coriolis 
    37    REAL(wp), PUBLIC ::   dm     = 0.6e+03   !: diffusion constant for dynamics 
    38    REAL(wp), PUBLIC ::   om     = 0.5       !: relaxation constant 
    39    REAL(wp), PUBLIC ::   resl   = 5.0e-05   !: maximum value for the residual of relaxation 
    40    REAL(wp), PUBLIC ::   cw     = 5.0e-03   !: drag coefficient for oceanic stress 
    41    REAL(wp), PUBLIC ::   angvg  = 0.e0      !: turning angle for oceanic stress 
    42    REAL(wp), PUBLIC ::   pstar  = 1.0e+04   !: first bulk-rheology parameter 
    43    REAL(wp), PUBLIC ::   c_rhg  = 20.e0     !: second bulk-rhelogy parameter 
    44    REAL(wp), PUBLIC ::   etamn  = 0.e+07    !: minimun value for viscosity 
    45    REAL(wp), PUBLIC ::   creepl = 2.e-08    !: creep limit 
    46    REAL(wp), PUBLIC ::   ecc    = 2.e0      !: eccentricity of the elliptical yield curve 
    47    REAL(wp), PUBLIC ::   ahi0   = 350.e0    !: sea-ice hor. eddy diffusivity coeff. (m2/s) 
     18   !---------------------------------------------------! 
     19   !                 Sea-Ice namelist                  ! 
     20   !---------------------------------------------------! 
     21   !                                                                !!! **  Namelist namicerun  ** 
     22   CHARACTER(len=32), PUBLIC ::   cn_icerst_in  = "restart_ice_in"   !: suffix of ice restart name (input) 
     23   CHARACTER(len=32), PUBLIC ::   cn_icerst_out = "restart_ice"      !: suffix of ice restart name (output) 
     24   LOGICAL          , PUBLIC ::   ln_limdyn     = .TRUE.             !: flag for ice dynamics (T) or not (F) 
     25   LOGICAL          , PUBLIC ::   ln_limdmp     = .FALSE.            !: Ice damping 
     26   REAL(wp)         , PUBLIC ::   hsndif        = 0.e0               !: computation of temp. in snow (0) or not (9999) 
     27   REAL(wp)         , PUBLIC ::   hicdif        = 0.e0               !: computation of temp. in ice (0) or not (9999) 
     28   REAL(wp)         , PUBLIC ::   acrit(2) = (/ 1.e-06 , 1.e-06 /)   !: min lead fraction in north and south hemispheres 
     29    
     30   !                                         !!! **  Namelist  namicedyn  ** 
     31   INTEGER , PUBLIC ::   nbiter =   1         !: number of sub-time steps for relaxation 
     32   INTEGER , PUBLIC ::   nbitdr = 250         !: maximum number of iterations for relaxation 
     33   REAL(wp), PUBLIC ::   rdt_ice              !: ice time step 
     34   REAL(wp), PUBLIC ::   epsd   =   1.0e-20   !: tolerance parameter for dynamic 
     35   REAL(wp), PUBLIC ::   alpha  =   0.5       !: coefficient for semi-implicit coriolis 
     36   REAL(wp), PUBLIC ::   dm     =   0.6e+03   !: diffusion constant for dynamics 
     37   REAL(wp), PUBLIC ::   om     =   0.5       !: relaxation constant 
     38   REAL(wp), PUBLIC ::   resl   =   5.0e-05   !: maximum value for the residual of relaxation 
     39   REAL(wp), PUBLIC ::   cw     =   5.0e-03   !: drag coefficient for oceanic stress 
     40   REAL(wp), PUBLIC ::   angvg  =   0.e0      !: turning angle for oceanic stress 
     41   REAL(wp), PUBLIC ::   pstar  =   1.0e+04   !: first bulk-rheology parameter 
     42   REAL(wp), PUBLIC ::   c_rhg  =  20.e0      !: second bulk-rhelogy parameter 
     43   REAL(wp), PUBLIC ::   etamn  =   0.e+07    !: minimun value for viscosity 
     44   REAL(wp), PUBLIC ::   creepl =   2.e-08    !: creep limit 
     45   REAL(wp), PUBLIC ::   ecc    =   2.e0      !: eccentricity of the elliptical yield curve 
     46   REAL(wp), PUBLIC ::   ahi0   = 350.e0      !: sea-ice hor. eddy diffusivity coeff. (m2/s) 
    4847 
    49    REAL(wp), PUBLIC ::   usecc2             !:  = 1.0 / ( ecc * ecc ) 
    50    REAL(wp), PUBLIC ::   rhoco              !: = rau0 * cw 
    51    REAL(wp), PUBLIC ::   sangvg, cangvg     !: sin and cos of the turning angle for ocean stress 
    52    REAL(wp), PUBLIC ::   pstarh             !: pstar / 2.0 
     48   !                                              !!! ** Namelist  namicethd  ** 
     49   REAL(wp), PUBLIC ::   hmelt      = -0.15        !: maximum melting at the bottom 
     50   REAL(wp), PUBLIC ::   hicmin     =  0.2         !: ice th. corr. to max. ener. in brine pocket 
     51   REAL(wp), PUBLIC ::   hiclim     =  0.05        !: minimum ice thickness 
     52   REAL(wp), PUBLIC ::   amax       =  0.999       !: maximum lead fraction 
     53   REAL(wp), PUBLIC ::   swiqst     =  1.0         !: energy stored in brine pocket (1) or not (0) 
     54   REAL(wp), PUBLIC ::   sbeta      =  1.0         !: numerical scheme for diffusion in ice  
     55   REAL(wp), PUBLIC ::   parlat     =  0.0         !: percent. of energy used for lateral ablation 
     56   REAL(wp), PUBLIC ::   hakspl     =  0.5         !: slope of distr. for Hakkinen-Mellro's lat. melt 
     57   REAL(wp), PUBLIC ::   hibspl     =  0.5         !: slope of distribution for Hibler's lat. melt 
     58   REAL(wp), PUBLIC ::   exld       =  2.0         !: exponent for leads-closure rate 
     59   REAL(wp), PUBLIC ::   hakdif     =  1.0         !: coefficient for diffusions of ice and snow 
     60   REAL(wp), PUBLIC ::   thth       =  0.2         !: thick. for comp. of eq. thermal conduct 
     61   REAL(wp), PUBLIC ::   hnzst      =  0.1         !: thick. of the surf. layer in temp. comp. 
     62   REAL(wp), PUBLIC ::   parsub     =  1.0         !: switch for snow sublimation or not 
     63   REAL(wp), PUBLIC ::   alphs      =  1.0         !: coef. for snow density when snow-ice formation 
     64   REAL(wp), PUBLIC ::   hiccrit(2) = (/0.3,0.3/)  !: ice th. for lateral accretion in the NH (SH) (m) 
     65 
     66   !---------------------------------------------------! 
     67   !                Sea-Ice variables                  ! 
     68   !---------------------------------------------------! 
     69 
     70   REAL(wp), PUBLIC ::   usecc2           !:  = 1.0 / ( ecc * ecc ) 
     71   REAL(wp), PUBLIC ::   rhoco            !: = rau0 * cw 
     72   REAL(wp), PUBLIC ::   sangvg, cangvg   !: sin and cos of the turning angle for ocean stress 
     73   REAL(wp), PUBLIC ::   pstarh           !: pstar / 2.0 
     74   REAL(wp), PUBLIC ::   uscomi           !: inverse of minimum lead fraction      !!gm to be suppress 
     75   REAL(wp), PUBLIC ::   cnscg            !: ratio  rcpsn/rcpic                    !!gm to be suppress 
    5376 
    5477   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ahiu , ahiv   !: hor. diffusivity coeff. at ocean U- and V-points (m2/s) 
    5578   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   pahu , pahv   !: ice hor. eddy diffusivity coef. at ocean U- and V-points 
    5679   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hsnm , hicm   !: mean snow and ice thicknesses 
    57    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ust2s                 !: friction velocity 
     80   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ust2s         !: friction velocity 
    5881 
    59    !!* diagnostic quantities 
    60    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sist          !: Sea-Ice Surface Temperature (Kelvin) 
    61    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   tfu           !: Freezing/Melting point temperature of sea water at SSS 
    62    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hicif         !: Ice thickness 
    63    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hsnif         !: Snow thickness 
    64    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   hicifp        !: Ice production/melting 
    65    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   frld          !: Leads fraction = 1-a/totalarea 
    66    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   phicif        !: ice thickness  at previous time  
    67    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   pfrld         !: Leads fraction at previous time   
    68    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qstoif        !: Energy stored in the brine pockets 
    69    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fbif          !: Heat flux at the ice base 
    70    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdmsnif       !: Variation of snow mass 
    71    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   rdmicif       !: Variation of ice mass 
    72    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qldif         !: heat balance of the lead (or of the open ocean) 
    73    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qcmif         !: Energy needed to bring the ocean surface layer until its freezing  
    74    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fdtcn         !: net downward heat flux from the ice to the ocean 
    75    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qdtcn         !: energy from the ice to the ocean point (at a factor 2) 
    76    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   thcm          !: part of the solar energy used in the lead heat budget 
    77    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fstric        !: Solar flux transmitted trough the ice 
    78    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   ffltbif       !: Array linked with the max heat contained in brine pockets (?) 
    79    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fscmbq        !: Linked with the solar flux below the ice (?) 
    80    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fsbbq         !: Also linked with the solar flux below the ice (?) 
    81    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   qfvbq         !: Array used to store energy in case of toral lateral ablation (?) 
    82    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   dmgwi         !: Variation of the mass of snow ice 
     82   INTEGER , PUBLIC, DIMENSION(jpij) ::   npb         !: number of points where computations has to be done 
     83   INTEGER , PUBLIC, DIMENSION(jpij) ::   npac        !: correspondance between the points 
    8384 
    84    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   u_ice, v_ice   !: two components of the ice   velocity at I-point (m/s) 
    85    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   u_oce, v_oce   !: two components of the ocean velocity at I-point (m/s) 
     85   !                                                         !!* sea-ice variables * 
     86   REAL(wp), PUBLIC ::   sist   (jpi,jpj), sist_1d   (jpij)   !: Sea-Ice Surface Temperature (Kelvin) 
     87   REAL(wp), PUBLIC ::   tfu    (jpi,jpj), tfu_1d    (jpij)   !: Freezing/Melting point temperature of sea water at SSS 
     88   REAL(wp), PUBLIC ::   hicif  (jpi,jpj), h_ice_1d  (jpij)   !: ice thickness 
     89   REAL(wp), PUBLIC ::   phicif (jpi,jpj)                     !: ice thickness  at previous time  
     90   REAL(wp), PUBLIC ::   hsnif  (jpi,jpj), h_snow_1d (jpij)   !: Snow thickness 
     91   REAL(wp), PUBLIC ::   hicifp (jpi,jpj)                     !: Ice production/melting 
     92   REAL(wp), PUBLIC ::   frld   (jpi,jpj), frld_1d   (jpij)   !: Leads fraction = 1-a/totalarea 
     93   REAL(wp), PUBLIC ::   pfrld  (jpi,jpj)                     !: Leads fraction at previous time   
     94   REAL(wp), PUBLIC ::   qstoif (jpi,jpj), qstbif_1d (jpij)   !: Energy stored in the brine pockets         !!gm  
     95   REAL(wp), PUBLIC ::   fbif   (jpi,jpj), fbif_1d   (jpij)   !: Heat flux at the ice base 
     96   REAL(wp), PUBLIC ::   rdmsnif(jpi,jpj), rdmsnif_1d(jpij)   !: Variation of snow mass 
     97   REAL(wp), PUBLIC ::   rdmicif(jpi,jpj), rdmicif_1d(jpij)   !: Variation of ice mass 
     98   REAL(wp), PUBLIC ::   qldif  (jpi,jpj), qldif_1d  (jpij)   !: heat balance of the lead (or of the open ocean) 
     99   REAL(wp), PUBLIC ::   qcmif  (jpi,jpj), qcmif_1d  (jpij)   !: Energy needed to freeze the ocean surface layer  
     100   REAL(wp), PUBLIC ::   fdtcn  (jpi,jpj)                     !: net downward heat flux from the ice to the ocean 
     101   REAL(wp), PUBLIC ::   qdtcn  (jpi,jpj)                     !: energy from the ice to the ocean point (at a factor 2) 
     102   REAL(wp), PUBLIC ::   thcm   (jpi,jpj), thcm_1d   (jpij)   !: part of the solar energy used in the lead heat budget 
     103   REAL(wp), PUBLIC ::   fstric (jpi,jpj), fstric_1d (jpij)   !: Solar flux transmitted trough the ice 
     104   REAL(wp), PUBLIC ::   ffltbif(jpi,jpj), fltbif_1d (jpij)   !: Array linked with the max heat contained in brine pockets (?) 
     105   REAL(wp), PUBLIC ::   fscmbq (jpi,jpj), fscbq_1d  (jpij)   !: Linked with the solar flux below the ice (?) 
     106   REAL(wp), PUBLIC ::   fsbbq  (jpi,jpj)                     !: Also linked with the solar flux below the ice (?) 
     107   REAL(wp), PUBLIC ::   qfvbq  (jpi,jpj), qfvbq_1d  (jpij)   !: Array used to store energy in case of total lateral ablation (?) 
     108   REAL(wp), PUBLIC ::   dmgwi(jpi,jpj), dmgwi_1d(jpij)       !: Variation of the mass of snow ice 
    86109 
    87    REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jplayersp1) ::   tbif  !: Temperature inside the ice/snow layer 
     110   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   u_ice, v_ice     !: two components of the ice   velocity at I-point (m/s) 
     111   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   u_oce, v_oce     !: two components of the ocean velocity at I-point (m/s) 
    88112 
    89    !!* moment used in the advection scheme 
     113   REAL(wp), PUBLIC ::   tbif(jpi,jpj,jplayersp1), tbif_1d(jpij,jplayersp1)  !: Temperature inside the ice+snow layer 
     114 
     115   ! Surface forcing transford into  1D field  (2D field defined in /OPA_SRC/SBC/sbc_ice.F90) 
     116   REAL(wp), PUBLIC, DIMENSION(jpij) ::   qns_ice_1d   !:   qns_ice    
     117   REAL(wp), PUBLIC, DIMENSION(jpij) ::   qsr_ice_1d   !:   qsr_ice 
     118   REAL(wp), PUBLIC, DIMENSION(jpij) ::   qla_ice_1d   !:   qla_ice 
     119   REAL(wp), PUBLIC, DIMENSION(jpij) ::   dqla_ice_1d  !:   dqla_ice 
     120   REAL(wp), PUBLIC, DIMENSION(jpij) ::   dqns_ice_1d  !:   dqns_ice 
     121   REAL(wp), PUBLIC, DIMENSION(jpij) ::   fr1_i0_1d    !:   fr1_i0 
     122   REAL(wp), PUBLIC, DIMENSION(jpij) ::   fr2_i0_1d    !:   fr2_i0 
     123 
     124   ! Surface forcing transformed into 1D field  (2D field defined in /OPA_SRC/SBC/sbc_oce.F90) 
     125   REAL(wp), PUBLIC, DIMENSION(jpij) ::   sprecip_1d   !:   sprecip 
     126    
     127   ! local variable  transformed into 1D field  (2D field defined in /LIM_SRC_2/limthd_2.F90) 
     128   REAL(wp), PUBLIC, DIMENSION(jpij) ::   qlbbq_1d     !:   zqlbsbq 
     129   REAL(wp), PUBLIC, DIMENSION(jpij) ::   dvsbq_1d     !:   zdvosif 
     130   REAL(wp), PUBLIC, DIMENSION(jpij) ::   dvbbq_1d     !:   zdvobif 
     131   REAL(wp), PUBLIC, DIMENSION(jpij) ::   dvlbq_1d     !:   zdvolif 
     132   REAL(wp), PUBLIC, DIMENSION(jpij) ::   dvnbq_1d     !:   zdvonif 
     133    
     134    
     135   REAL(wp), PUBLIC, DIMENSION(jpij) ::   rdvomif_1d    !:   rdvomif       <<<== !gm should be suppressed ! 
     136 
     137   !---------------------------------------------------! 
     138   !   Prather' advection scheme : 1st & 2nd moments   ! 
     139   !---------------------------------------------------! 
    90140   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxice, syice, sxxice, syyice, sxyice   !: for ice  volume 
    91141   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   sxsn,  sysn,  sxxsn,  syysn,  sxysn    !: for snow volume 
    … …  
    102152#endif 
    103153 
     154   !!-------- ------------------------------------------------------------- 
     155   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
     156   !! $Id$ 
     157   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    104158   !!====================================================================== 
    105159END MODULE ice_2 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/iceini_2.F90

    r1581 r1855  
    44   !!   Sea-ice model : LIM 2.0 Sea ice model Initialization 
    55   !!====================================================================== 
    6    !! History :   1.0  !  02-08  (G. Madec)  F90: Free form and modules 
    7    !!             2.0  !  03-08  (C. Ethe)  add ice_run 
     6   !! History :  1.0  !  2002-08  (G. Madec)  F90: Free form and modules 
     7   !!            2.0  !  2003-08  (C. Ethe)  add ice_run 
    88   !!---------------------------------------------------------------------- 
    99#if defined key_lim2 
    … …  
    1111   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1212   !!---------------------------------------------------------------------- 
     13   !!   ice_init_2      : sea-ice model initialization 
    1314   !!---------------------------------------------------------------------- 
    14    !!   ice_init_2       : sea-ice model initialization 
    15    !!   ice_run_2        : Definition some run parameter for ice model 
    16    !!---------------------------------------------------------------------- 
    17    USE dom_oce 
    18    USE dom_ice_2 
    19    USE sbc_oce         ! surface boundary condition: ocean 
    20    USE sbc_ice         ! surface boundary condition: ice 
    21    USE phycst          ! Define parameters for the routines 
    22    USE ice_2 
    23    USE limmsh_2 
    24    USE limistate_2 
    25    USE limrst_2    
    26    USE in_out_manager 
     15   USE dom_oce          ! ocean domain 
     16   USE dom_ice_2        ! LIM-2 ice domain 
     17   USE sbc_oce          ! surface boundary condition: ocean 
     18   USE sbc_ice          ! surface boundary condition: sea ice 
     19   USE phycst           ! physical constant 
     20   USE ice_2            ! LIM-2 ice varibles 
     21   USE limmsh_2         ! LIM-2 mesh 
     22   USE limistate_2      ! LIM-2 initial state 
     23   USE limrst_2         ! LIM-2 restart 
     24   USE in_out_manager   ! I/O manager 
    2725       
    2826   IMPLICIT NONE 
    2927   PRIVATE 
    3028 
    31    PUBLIC   ice_init_2               ! called by sbcice_lim_2.F90 
     29   PUBLIC   ice_init_2   ! called by sbcice_lim_2.F90 
    3230 
    3331   !!---------------------------------------------------------------------- 
    34    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     32   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    3533   !! $Id$  
    36    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     34   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3735   !!---------------------------------------------------------------------- 
    3836 
    … …  
    4543      !! ** purpose :    
    4644      !!---------------------------------------------------------------------- 
    47       ! 
    48       ! Open the namelist file  
    49       CALL ctl_opn( numnam_ice, 'namelist_ice', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )       
    50       CALL ice_run_2                    !  read in namelist some run parameters 
    51                   
    52       ! Louvain la Neuve Ice model 
    53       rdt_ice = nn_fsbc * rdttra(1) 
    54  
    55       CALL lim_msh_2                  ! ice mesh initialization 
    56       
    57       ! Initial sea-ice state 
    58       IF( .NOT.ln_rstart ) THEN 
    59          CALL lim_istate_2            ! start from rest: sea-ice deduced from sst 
    60       ELSE 
    61          CALL lim_rst_read_2          ! start from a restart file 
    62       ENDIF 
    63        
    64       tn_ice(:,:,1) = sist(:,:)         ! initialisation of ice temperature    
    65       fr_i  (:,:) = 1.0 - frld(:,:)   ! initialisation of sea-ice fraction     
    66       ! 
    67    END SUBROUTINE ice_init_2 
    68  
    69  
    70    SUBROUTINE ice_run_2 
    71       !!------------------------------------------------------------------- 
    72       !!                  ***  ROUTINE ice_run_2 *** 
    73       !!                  
    74       !! ** Purpose :   Definition some run parameter for ice model 
    75       !! 
    76       !! ** Method  :   Read the namicerun namelist and check the parameter  
    77       !!       values called at the first timestep (nit000) 
    78       !! 
    79       !! ** input   :   Namelist namicerun 
    80       !!------------------------------------------------------------------- 
    8145      NAMELIST/namicerun/ cn_icerst_in, cn_icerst_out, ln_limdyn, ln_limdmp, acrit, hsndif, hicdif 
    8246      !!------------------------------------------------------------------- 
    83       !                     
    84       REWIND ( numnam_ice )                       ! Read Namelist namicerun  
     47      ! 
     48      !                                      ! Open the ice namelist file  
     49      CALL ctl_opn( numnam_ice, 'namelist_ice', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp )     
     50         
     51      REWIND ( numnam_ice )                  ! Read namicerun namelist 
    8552      READ   ( numnam_ice , namicerun ) 
    86  
    87       IF(lwp) THEN 
     53      ! 
     54      IF(lwp) THEN                           ! control print 
    8855         WRITE(numout,*) 
    8956         WRITE(numout,*) 'ice_run : ice share parameters for dynamics/advection/thermo of sea-ice' 
    … …  
    9562         WRITE(numout,*) '   computation of temp. in ice  (=0) or not (=9999) hicdif = ', hicdif 
    9663      ENDIF 
     64                        
     65      rdt_ice = nn_fsbc * rdttra(1)          ! ice time step  
     66 
     67      CALL lim_msh_2                         ! ice mesh initialization 
     68      
     69      !                                      ! Initial sea-ice state 
     70      IF( .NOT.ln_rstart ) THEN   ;   CALL lim_istate_2        ! start from rest: sea-ice deduced from sst 
     71      ELSE                        ;   CALL lim_rst_read_2      ! restart from a file 
     72      ENDIF 
     73 
     74      !                                      ! ice variables see by the ocean       
     75      tn_ice(:,:,1) = sist(:,:)                 ! surface ice temperature    
     76      fr_i  (:,:)   = 1.0 - frld(:,:)           ! sea-ice fraction     
    9777      ! 
    98    END SUBROUTINE ice_run_2 
     78   END SUBROUTINE ice_init_2 
    9979 
    10080#else 
    … …  
    10282   !!   Default option :        Empty module       NO LIM 2.0 sea-ice model 
    10383   !!---------------------------------------------------------------------- 
    104 CONTAINS 
    105    SUBROUTINE ice_init_2      ! Empty routine 
    106    END SUBROUTINE ice_init_2 
    10784#endif 
    10885 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limadv_2.F90

    r1530 r1855  
    22   !!====================================================================== 
    33   !!                       ***  MODULE limadv_2   *** 
    4    !! LIM 2.0 sea-ice model : sea-ice advection 
     4   !! LIM-2 sea-ice model : sea-ice advection 
    55   !!====================================================================== 
    66   !! History :  OPA  ! 2000-01 (LIM)  Original code 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limdia_2.F90

    r1715 r1855  
    1212   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1313   !!---------------------------------------------------------------------- 
    14    !!---------------------------------------------------------------------- 
    1514   !!   lim_dia_2      : computation of the time evolution of keys var. 
    1615   !!   lim_dia_init_2 : initialization and namelist read 
    … …  
    2827   PRIVATE 
    2928 
    30    PUBLIC               lim_dia_2          ! called by sbc_ice_lim_2 
    31    INTEGER, PUBLIC ::   ntmoy   = 1 ,   &  !: instantaneous values of ice evolution or averaging ntmoy 
    32       &                 ninfo   = 1        !: frequency of ouputs on file ice_evolu in case of averaging 
    33  
    34    INTEGER, PARAMETER ::   &  ! Parameters for outputs to files "evolu" 
    35       jpinfmx = 100         ,    &  ! maximum number of key variables 
    36       jpchinf = 5           ,    &  ! ??? 
    37       jpchsep = jpchinf + 2         ! ??? 
    38  
    39    INTEGER ::   & 
    40       nfrinf  = 4 ,     &  ! number of variables written in one line  
    41       nferme ,          &  ! last time step at which the var. are written on file 
    42       nvinfo ,          &  ! number of total variables  
    43       nbvt   ,          &  ! number of time variables 
    44       naveg                ! number of step for accumulation before averaging 
    45  
    46    CHARACTER(len= 8) ::   fmtinf  = '1PE13.5 ' ! format of the output values   
    47    CHARACTER(len=30) ::   fmtw  ,           &  ! formats 
    48       &                   fmtr  ,           &  ! ??? 
    49       &                   fmtitr               ! ??? 
    50    CHARACTER(len=jpchsep), DIMENSION(jpinfmx) ::   titvar               ! title of key variables 
     29   PUBLIC   lim_dia_2   ! called by sbc_ice_lim_2 
     30    
     31   INTEGER, PUBLIC ::   ntmoy   = 1   !: instantaneous values of ice evolution or averaging ntmoy 
     32   INTEGER, PUBLIC ::   ninfo   = 1   !: frequency of ouputs on file ice_evolu in case of averaging 
     33 
     34   !                                              !!! Parameters for outputs to files "evolu" 
     35   INTEGER, PARAMETER ::   jpinfmx = 100           ! maximum number of key variables 
     36   INTEGER, PARAMETER ::   jpchinf = 5             ! ??? 
     37   INTEGER, PARAMETER ::   jpchsep = jpchinf + 2   ! ??? 
     38 
     39   INTEGER ::   nfrinf  = 4   ! number of variables written in one line  
     40   INTEGER ::   nferme        ! last time step at which the var. are written on file 
     41   INTEGER ::   nvinfo        ! number of total variables  
     42   INTEGER ::   nbvt          ! number of time variables 
     43   INTEGER ::   naveg         ! number of step for accumulation before averaging 
     44 
     45   CHARACTER(len= 8) ::   fmtinf  = '1PE13.5 '   ! format of the output values   
     46   CHARACTER(len=30) ::   fmtw, fmtr, fmtitr     ! various formats 
     47   CHARACTER(len=jpchsep), DIMENSION(jpinfmx) ::   titvar   ! title of key variables 
    5148  
    5249   REAL(wp)                     ::   epsi06 = 1.e-06      ! ??? 
    … …  
    5754#  include "vectopt_loop_substitute.h90" 
    5855   !!---------------------------------------------------------------------- 
    59    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     56   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    6057   !! $Id$ 
    6158   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    7168      !!      the temporal evolution of some key variables 
    7269      !!------------------------------------------------------------------- 
    73       INTEGER, INTENT(in) ::   kt     ! number of iteration 
     70      INTEGER, INTENT(in) ::   kt   ! number of iteration 
    7471      !! 
    7572      INTEGER  ::   jv,ji, jj   ! dummy loop indices 
    7673      INTEGER  ::   nv          ! indice of variable  
    77       REAL(wp) ::   zarea    , zldarea  ,    &  ! sea-ice and leads area 
    78          &          zextent15, zextent85,    &  ! sea-ice extent (15% and 85%) 
    79          &          zicevol  , zsnwvol  ,    &  ! sea-ice and snow volume volume 
    80          &          zicespd                     ! sea-ice velocity 
     74      REAL(wp) ::   zarea    , zldarea     ! sea-ice and leads area 
     75      REAL(wp) ::   zextent15, zextent85   ! sea-ice extent (15% and 85%) 
     76      REAL(wp) ::   zicevol  , zsnwvol     ! sea-ice and snow volume volume 
     77      REAL(wp) ::   zicespd                     ! sea-ice velocity 
    8178      REAL(wp), DIMENSION(jpinfmx) ::   vinfor  ! temporary working space  
    8279      !!------------------------------------------------------------------- 
    … …  
    120117      vinfor(nv+13) = SQRT( vinfor(nv+13) / MAX( vinfor(nv+9) , epsi06 ) ) 
    121118 
    122  
    123119      ! variables in southern Hemis 
    124120       nv = nv + 1 
    … …  
    177173       INTEGER  ::   nv            ! indice of variable  
    178174       REAL(wp) ::   zxx0, zxx1    ! temporary scalars 
    179  
     175       !! 
    180176       NAMELIST/namicedia/fmtinf, nfrinf, ninfo, ntmoy 
    181177       !!------------------------------------------------------------------- 
    182178 
    183        ! Read Namelist namicedia 
    184        REWIND ( numnam_ice ) 
     179       REWIND ( numnam_ice )                 ! Read Namelist namicedia 
    185180       READ   ( numnam_ice  , namicedia ) 
    186181        
    187        IF(lwp) THEN 
     182       IF(lwp) THEN                          ! control print 
    188183          WRITE(numout,*) 
    189184          WRITE(numout,*) 'lim_dia_init_2 : ice parameters for ice diagnostics ' 
    … …  
    274269   !!   Default option :                           NO LIM 2.0 sea-ice model 
    275270   !!---------------------------------------------------------------------- 
    276 CONTAINS 
    277    SUBROUTINE lim_dia_2         ! Empty routine 
    278    END SUBROUTINE lim_dia_2 
    279271#endif 
    280272 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limdmp_2.F90

    r1715 r1855  
    44   !!  Ice model : restoring Ice thickness and Fraction leads 
    55   !!====================================================================== 
    6    !! History :   2.0  !  04-04 (S. Theetten) Original code 
     6   !! History :   2.0  ! 2004-04 (S. Theetten) Original code 
    77   !!---------------------------------------------------------------------- 
    88#if defined key_lim2   &&   defined key_tradmp 
    … …  
    1111   !!   'key_tradmp'                                                Damping 
    1212   !!---------------------------------------------------------------------- 
     13   !!   lim_dmp_2       : ice model damping 
    1314   !!---------------------------------------------------------------------- 
    14    !!   lim_dmp_2      : ice model damping 
    15    !!---------------------------------------------------------------------- 
    16    USE in_out_manager  ! I/O manager 
    17    USE phycst          ! physical constants 
    18    USE ice_2 
    19    USE tradmp 
    20    USE dom_oce 
    21    USE oce 
    22    USE iom 
     15   USE oce              ! ocean variables 
     16   USE dom_oce          ! ocean domain 
     17   USE phycst           ! physical constants 
     18   USE ice_2            ! LIM-2 variables 
     19   USE tradmp           ! traceur damping 
     20   USE in_out_manager   ! I/O manager 
     21   USE iom              ! 
    2322    
    2423   IMPLICIT NONE 
    … …  
    2726   PUBLIC   lim_dmp_2     ! called by ice_step_2 
    2827    
    29    INTEGER                        ::   nice1, nice2,  &  ! first and second record used 
    30       &                                inumice_dmp       ! logical unit for ice variables (damping) 
    31    REAL(wp), DIMENSION(jpi,jpj)   ::   hicif_dta  ,   &  ! ice thickness at a given time 
    32       &                                frld_dta          ! fraction lead at a given time 
    33    REAL(wp), DIMENSION(jpi,jpj,2) ::   hicif_data ,   &  ! ice thickness data at two consecutive times 
    34       &                                frld_data         ! fraction lead data at two consecutive times 
     28   INTEGER                        ::   nice1, nice2   ! first and second record used 
     29   INTEGER                        ::   inumice_dmp    ! logical unit for ice variables (damping) 
     30   REAL(wp), DIMENSION(jpi,jpj)   ::   hicif_dta      ! ice thickness at a given time 
     31   REAL(wp), DIMENSION(jpi,jpj)   ::   frld_dta       ! fraction lead at a given time 
     32   REAL(wp), DIMENSION(jpi,jpj,2) ::   hicif_data     ! ice thickness data at two consecutive times 
     33   REAL(wp), DIMENSION(jpi,jpj,2) ::   frld_data      ! fraction lead data at two consecutive times 
    3534 
    3635   !! * Substitution 
    3736#  include "vectopt_loop_substitute.h90" 
    3837   !!---------------------------------------------------------------------- 
    39    !!   LIM 2.0 , UCL-LOCEAN-IPSL  (2006) 
     38   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    4039   !! $Id$ 
    4140   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    4443CONTAINS 
    4544 
    46    SUBROUTINE lim_dmp_2(kt) 
     45   SUBROUTINE lim_dmp_2( kt ) 
    4746      !!------------------------------------------------------------------- 
    4847      !!                   ***  ROUTINE lim_dmp_2 *** 
    … …  
    5352      !! ** method  : the key_tradmp must be used to compute resto(:,:) coef. 
    5453      !!--------------------------------------------------------------------- 
    55       INTEGER, INTENT(in) ::   kt     ! ocean time-step 
    56       ! 
    57       INTEGER             ::   ji, jj         ! dummy loop indices 
     54      INTEGER, INTENT(in) ::   kt   ! ocean time-step 
     55      !! 
     56      INTEGER ::   ji, jj   ! dummy loop indices 
    5857      !!--------------------------------------------------------------------- 
    5958      ! 
    6059      CALL dta_lim_2( kt ) 
    61  
     60      ! 
    6261      DO jj = 2, jpjm1 
    6362         DO ji = fs_2, fs_jpim1   ! vector opt. 
    6463            hicif(ji,jj) = hicif(ji,jj) - rdt_ice * resto(ji,jj,1) * ( hicif(ji,jj) - hicif_dta(ji,jj) ) 
    65             frld(ji,jj)  = frld (ji,jj) - rdt_ice * resto(ji,jj,1) * ( frld(ji,jj) - frld_dta (ji,jj) )   
     64            frld(ji,jj)  = frld (ji,jj) - rdt_ice * resto(ji,jj,1) * ( frld (ji,jj) - frld_dta (ji,jj) )   
    6665         END DO 
    6766      END DO 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limdyn_2.F90

    r1694 r1855  
    1717   !!---------------------------------------------------------------------- 
    1818   USE dom_oce        ! ocean space and time domain 
    19    USE sbc_oce        ! 
    20    USE phycst         ! 
    21    USE ice_2          ! 
    22    USE dom_ice_2      ! 
    23    USE limistate_2    ! 
    24    USE limrhg_2       ! ice rheology 
     19   USE sbc_oce        ! surface boundary condition : ocean 
     20   USE phycst         ! physical constants 
     21   USE ice_2          ! LIM-2 sea-ice variables 
     22   USE dom_ice_2      ! LIM-2 domain 
     23   USE limistate_2    ! LIM-2 initial state 
     24   USE limrhg_2       ! LIM-2 rheology 
    2525 
    2626   USE lbclnk         ! 
    … …  
    3434   PUBLIC   lim_dyn_2 ! routine called by sbc_ice_lim 
    3535 
    36    !! * Module variables 
    3736   REAL(wp)  ::  rone    = 1.e0   ! constant value 
    3837 
    3938#  include "vectopt_loop_substitute.h90" 
    4039   !!---------------------------------------------------------------------- 
    41    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     40   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    4241   !! $Id$ 
    4342   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    5958      !!              - treatment of the case if no ice dynamic 
    6059      !!--------------------------------------------------------------------- 
    61       INTEGER, INTENT(in) ::   kt     ! number of iteration 
    62       !! 
    63       INTEGER  ::   ji, jj             ! dummy loop indices 
    64       INTEGER  ::   i_j1, i_jpj        ! Starting/ending j-indices for rheology 
    65       REAL(wp) ::   zcoef              ! temporary scalar 
     60      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) ::   zcoef         ! temporary scalar 
    6665      REAL(wp), DIMENSION(jpj)     ::   zind           ! i-averaged indicator of sea-ice 
    6766      REAL(wp), DIMENSION(jpj)     ::   zmsk           ! i-averaged of tmask 
    … …  
    144143         ENDIF 
    145144 
    146          IF(ln_ctl)   CALL prt_ctl(tab2d_1=u_ice , clinfo1=' lim_dyn  : u_ice :', tab2d_2=v_ice , clinfo2=' v_ice :') 
     145         IF(ln_ctl)   CALL prt_ctl( tab2d_1=u_ice, clinfo1=' lim_dyn  : u_ice :', tab2d_2=v_ice, clinfo2=' v_ice :' ) 
    147146          
    148147         ! computation of friction velocity 
    … …  
    179178      CALL lbc_lnk( ust2s, 'T',  1. )   ! T-point 
    180179      ! 
    181       IF(ln_ctl)   CALL prt_ctl(tab2d_1=ust2s , clinfo1=' lim_dyn  : ust2s :') 
     180      IF(ln_ctl)   CALL prt_ctl( tab2d_1=ust2s, clinfo1=' lim_dyn  : ust2s :' ) 
    182181 
    183182   END SUBROUTINE lim_dyn_2 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limhdf_2.F90

    r1465 r1855  
    44   !! LIM 2.0 ice model : horizontal diffusion of sea-ice quantities 
    55   !!====================================================================== 
     6   !! History :  LIM  ! 2000-01 (UCL)  Original code 
     7   !!            1.0  ! 2001-05 (G. Madec, R. Hordoir) opa norm 
     8   !!            2.0  ! 2002-08 (C. Ethe)  F90, free form 
     9   !!---------------------------------------------------------------------- 
    610#if defined key_lim2 
    711   !!---------------------------------------------------------------------- 
    … …  
    1014   !!   lim_hdf_2  : diffusion trend on sea-ice variable 
    1115   !!---------------------------------------------------------------------- 
    12    !! * Modules used 
    13    USE dom_oce 
    14    USE in_out_manager 
    15    USE ice_2 
    16    USE lbclnk 
    17    USE lib_mpp 
    18    USE prtctl          ! Print control 
     16   USE dom_oce          ! ocean domain 
     17   USE ice_2            ! LIM-2 variables 
     18   USE prtctl           ! Print control 
     19   USE lbclnk           ! 
     20   USE lib_mpp          ! 
     21   USE in_out_manager   ! I/O manager 
    1922 
    2023   IMPLICIT NONE 
    2124   PRIVATE 
    2225 
    23    !! * Routine accessibility 
    24    PUBLIC lim_hdf_2    ! called by lim_tra_2 
     26   PUBLIC   lim_hdf_2   ! called by lim_tra_2 
    2527 
    26    !! * Module variables 
    27    LOGICAL  ::   linit = .TRUE.              ! ??? 
     28   LOGICAL  ::   linit = .TRUE.              ! indictor of initialisation 
    2829   REAL(wp) ::   epsi04 = 1e-04              ! constant 
    29    REAL(wp), DIMENSION(jpi,jpj) ::   zfact   ! ??? 
     30   REAL(wp), DIMENSION(jpi,jpj) ::   zfact   ! metric term 
    3031 
    3132   !! * Substitution  
    3233#  include "vectopt_loop_substitute.h90" 
    3334   !!---------------------------------------------------------------------- 
    34    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     35   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    3536   !! $Id$ 
    36    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     37   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3738   !!---------------------------------------------------------------------- 
    3839 
    … …  
    4344      !!                  ***  ROUTINE lim_hdf_2  *** 
    4445      !! 
    45       !! ** purpose :   Compute and add the diffusive trend on sea-ice 
    46       !!      variables 
     46      !! ** purpose :   Compute and add the diffusive trend on sea-ice variables 
    4747      !! 
    4848      !! ** method  :   Second order diffusive operator evaluated using a 
    49       !!      Cranck-Nicholson time Scheme. 
     49      !!              Cranck-Nicholson time Scheme. 
    5050      !! 
    51       !! ** Action  :    update ptab with the diffusive contribution 
     51      !! ** Action  :   update ptab with the diffusive contribution 
     52      !!------------------------------------------------------------------- 
     53      REAL(wp), INTENT(inout), DIMENSION(jpi,jpj) ::   ptab   ! Field on which the diffusion is applied   
    5254      !! 
    53       !! History : 
    54       !!        !  00-01 (LIM) Original code 
    55       !!        !  01-05 (G. Madec, R. Hordoir) opa norm 
    56       !!        !  02-08 (C. Ethe)  F90, free form 
    57       !!------------------------------------------------------------------- 
    58       ! * Arguments 
    59       REAL(wp), DIMENSION(jpi,jpj), INTENT( inout ) ::   & 
    60          ptab                 ! Field on which the diffusion is applied   
    61       REAL(wp), DIMENSION(jpi,jpj) ::   & 
    62          ptab0                ! ??? 
    63  
    64       ! * Local variables 
    6555      INTEGER ::  ji, jj      ! dummy loop indices 
    66       INTEGER ::  & 
    67          its, iter            ! temporary integers 
     56      INTEGER ::  its, iter   ! temporary integers 
    6857      CHARACTER (len=55) :: charout 
    69       REAL(wp) ::  & 
    70          zalfa, zrlxint, zconv, zeps   ! temporary scalars 
    71       REAL(wp), DIMENSION(jpi,jpj) ::  &  
    72          zrlx, zflu, zflv, &  ! temporary workspaces 
    73          zdiv0, zdiv          !    "           " 
     58      REAL(wp) ::   zalfa, zrlxint, zconv, zeps   ! temporary scalars 
     59      REAL(wp), DIMENSION(jpi,jpj) ::   zrlx, zflu, zflv, zdiv0, zdiv  ! temporary workspaces 
     60      REAL(wp), DIMENSION(jpi,jpj) ::   ztab0                ! ??? 
    7461      !!------------------------------------------------------------------- 
    7562 
    … …  
    8269 
    8370      ! Arrays initialization 
    84       ptab0 (:, : ) = ptab(:,:) 
    85 !bug  zflu (:,jpj) = 0.e0 
    86 !bug  zflv (:,jpj) = 0.e0 
     71      ztab0(:, : ) = ptab(:,:) 
    8772      zdiv0(:, 1 ) = 0.e0 
    8873      zdiv0(:,jpj) = 0.e0 
    … …  
    9883         DO jj = 2, jpjm1   
    9984            DO ji = fs_2 , fs_jpim1   ! vector opt. 
    100                zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj  ) + e1v(ji,jj) + e1v(ji,jj-1) ) & 
    101                   &          / ( e1t(ji,jj) * e2t(ji,jj) ) 
     85               zfact(ji,jj) = ( e2u(ji,jj) + e2u(ji-1,jj) + e1v(ji,jj) + e1v(ji,jj-1) ) / ( e1t(ji,jj) * e2t(ji,jj) ) 
    10286            END DO 
    10387         END DO 
    … …  
    11094      iter  = 0 
    11195 
    112       !                                                   !=================== 
    113       DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) )   ! Sub-time step loop 
    114          !                                                !=================== 
    115          ! incrementation of the sub-time step number 
    116          iter = iter + 1 
     96      !                                                   !======================! 
     97      DO WHILE ( ( zconv > zeps ) .AND. (iter <= its) )   !  Sub-time step loop  ! 
     98         !                                                !======================! 
     99         iter = iter + 1               ! incrementation of the sub-time step number 
    117100 
    118          ! diffusive fluxes in U- and V- direction 
    119          DO jj = 1, jpjm1 
     101         DO jj = 1, jpjm1              ! diffusive fluxes in U- and V- direction 
    120102            DO ji = 1 , fs_jpim1   ! vector opt. 
    121103               zflu(ji,jj) = pahu(ji,jj) * e2u(ji,jj) / e1u(ji,jj) * ( ptab(ji+1,jj) - ptab(ji,jj) ) 
    … …  
    123105            END DO 
    124106         END DO 
    125  
    126          ! diffusive trend : divergence of the fluxes 
    127          DO jj= 2, jpjm1 
     107         DO jj= 2, jpjm1               ! diffusive trend : divergence of the fluxes 
    128108            DO ji = fs_2 , fs_jpim1   ! vector opt.  
    129109               zdiv (ji,jj) = (  zflu(ji,jj) - zflu(ji-1,jj  )   & 
    … …  
    131111            END DO 
    132112         END DO 
    133  
    134          ! save the first evaluation of the diffusive trend in zdiv0 
     113         !                             ! save the first evaluation of the diffusive trend in zdiv0 
    135114         IF( iter == 1 )   zdiv0(:,:) = zdiv(:,:)        
    136115 
    137          ! XXXX iterative evaluation????? 
    138          DO jj = 2, jpjm1 
     116         DO jj = 2, jpjm1              ! XXXX iterative evaluation????? 
    139117            DO ji = fs_2 , fs_jpim1   ! vector opt. 
    140                zrlxint = (   ptab0(ji,jj)    & 
     118               zrlxint = (   ztab0(ji,jj)    & 
    141119                  &       +  rdt_ice * (           zalfa   * ( zdiv(ji,jj) + zfact(ji,jj) * ptab(ji,jj) )   & 
    142120                  &                      + ( 1.0 - zalfa ) *   zdiv0(ji,jj) )  )                             &  
    … …  
    145123            END DO 
    146124         END DO 
    147  
    148          ! lateral boundary condition on ptab 
    149125         CALL lbc_lnk( zrlx, 'T', 1. ) 
    150126 
    151          ! convergence test 
    152          zconv = 0.e0 
     127         zconv = 0.e0                  ! convergence test 
    153128         DO jj = 2, jpjm1 
    154129            DO ji = 2, jpim1 
    … …  
    157132         END DO 
    158133         IF( lk_mpp )   CALL mpp_max( zconv )   ! max over the global domain 
    159  
    160          ptab(:,:) = zrlx(:,:) 
    161  
    162          !                                         !========================== 
    163       END DO                                       ! end of sub-time step loop 
    164       !                                            !========================== 
     134         ! 
     135         ptab(:,:) = zrlx(:,:)         ! update value 
     136         !                                         !=============================! 
     137      END DO                                       !  end of sub-time step loop  ! 
     138      !                                            !=============================! 
    165139 
    166140      IF(ln_ctl)   THEN 
    167          zrlx(:,:) = ptab(:,:) - ptab0(:,:) 
     141         zrlx(:,:) = ptab(:,:) - ztab0(:,:) 
    168142         WRITE(charout,FMT="(' lim_hdf  : zconv =',D23.16, ' iter =',I4,2X)") zconv, iter 
    169          CALL prt_ctl(tab2d_1=zrlx, clinfo1=charout) 
     143         CALL prt_ctl( tab2d_1=zrlx, clinfo1=charout ) 
    170144      ENDIF 
    171  
     145      ! 
    172146   END SUBROUTINE lim_hdf_2 
    173147 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limistate_2.F90

    r1471 r1855  
    44   !!              Initialisation of diagnostics ice variables 
    55   !!====================================================================== 
    6    !! History :   1.0  !  01-04  (C. Ethe, G. Madec)  Original code 
    7    !!             2.0  !  03-08  (G. Madec)  add lim_istate_init 
    8    !!                  !  04-04  (S. Theetten) initialization from a file 
    9    !!                  !  06-07  (S. Masson)  IOM to read the restart 
    10    !!                  !  07-10  (G. Madec)  surface module 
     6   !! History :   1.0  !  2001-04  (C. Ethe, G. Madec)  Original code 
     7   !!             2.0  !  2003-08  (G. Madec)  add lim_istate_init 
     8   !!              -   !  2004-04  (S. Theetten) initialization from a file 
     9   !!              -   !  2006-07  (S. Masson)  IOM to read the restart 
     10   !!              -   !  2007-10  (G. Madec)  surface module 
    1111   !!-------------------------------------------------------------------- 
    1212#if defined key_lim2 
    … …  
    1414   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1515   !!---------------------------------------------------------------------- 
    16    !!---------------------------------------------------------------------- 
    1716   !!   lim_istate_2      :  Initialisation of diagnostics ice variables 
    1817   !!   lim_istate_init_2 :  initialization of ice state and namelist read 
    1918   !!---------------------------------------------------------------------- 
    20    USE phycst 
    21    USE par_ice_2       ! ice parameters 
    22    USE dom_ice_2 
    23    USE eosbn2          ! equation of state 
    24    USE lbclnk 
    25    USE oce 
    26    USE ice_2 
    27    USE iom 
    28    USE in_out_manager 
     19   USE oce              ! ocean variables 
     20   USE ice_2            ! LIM-2 variables 
     21   USE par_ice_2        ! LIM-2 ice parameters 
     22   USE dom_ice_2        ! LIM-2 domain 
     23   USE phycst           ! physical constants 
     24   USE eosbn2           ! equation of state 
     25   USE lbclnk           !  
     26   USE iom              ! 
     27   USE in_out_manager   ! 
    2928 
    3029   IMPLICIT NONE 
    3130   PRIVATE 
    3231 
    33    PUBLIC lim_istate_2      ! routine called by lim_init_2.F90 
    34  
    35 !!! ** init namelist (namiceini) ** 
    36    LOGICAL  ::   ln_limini = .FALSE.  !: Ice initialization state 
     32   PUBLIC   lim_istate_2   ! routine called by lim_init_2.F90 
     33 
     34   !                                 !!! ** init namelist (namiceini) ** 
     35   LOGICAL  ::   ln_limini = .FALSE.  ! Ice initialization state 
    3736   REAL(wp) ::   ttest     = 2.0      ! threshold water temperature for initial sea ice 
    3837   REAL(wp) ::   hninn     = 0.5      ! initial snow thickness in the north 
    … …  
    4544   REAL(wp) ::   zero      = 0.e0     ! constant value = 0 
    4645   REAL(wp) ::   zone      = 1.e0     ! constant value = 1 
    47    !!---------------------------------------------------------------------- 
    48    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     46    
     47   !!---------------------------------------------------------------------- 
     48   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    4949   !! $Id$ 
    5050   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    6262      !!                or from arbitrary sea-ice conditions 
    6363      !!-------------------------------------------------------------------- 
    64       INTEGER  ::   ji, jj, jk                ! dummy loop indices 
    65       REAL(wp) ::   zidto                     ! temporary scalar 
     64      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     65      REAL(wp) ::   zidto        ! temporary scalar 
    6666      !-------------------------------------------------------------------- 
    6767  
    … …  
    6969 
    7070      IF( .NOT. ln_limini ) THEN   
    71           
     71         ! 
    7272         tfu(:,:) = tfreez( sn(:,:,1) ) * tmask(:,:,1)       ! freezing/melting point of sea water [Celcius] 
    73  
     73         ! 
    7474         DO jj = 1, jpj 
    7575            DO ji = 1, jpi 
    … …  
    126126 
    127127      !-- lateral boundary conditions 
    128       CALL lbc_lnk( hicif, 'T', 1. ) 
    129       CALL lbc_lnk( frld , 'T', 1. ) 
     128      CALL lbc_lnk( hicif, 'T', 1. )   ;   CALL lbc_lnk( frld , 'T', 1. ) 
    130129 
    131130      ! C A U T I O N  frld = 1 over land and lbc_lnk put zero along  
    … …  
    140139      CALL lbc_lnk( fsbbq  , 'T', 1. ) 
    141140      CALL lbc_lnk( qstoif , 'T', 1. ) 
    142  
     141      ! 
    143142   END SUBROUTINE lim_istate_2 
    144143 
    … …  
    151150      !! 
    152151      !! ** Method  :   Read the namiceini namelist and check the parameter  
    153       !!       values called at the first timestep (nit000) 
     152      !!              values called at the first timestep (nit000) 
    154153      !! 
    155154      !! ** input   :   Namelist namiceini 
    156155      !!------------------------------------------------------------------- 
    157       INTEGER :: inum_ice 
    158       INTEGER :: ji,jj 
    159  
     156      INTEGER ::   ji,jj      ! dummy loop indices 
     157      INTEGER ::   inum_ice   ! temporary integer 
     158      !! 
    160159      NAMELIST/namiceini/ ln_limini, ttest, hninn, hginn, alinn, & 
    161          &                hnins, hgins, alins 
     160         &                                  hnins, hgins, alins 
    162161      !!------------------------------------------------------------------- 
    163162      ! 
    … …  
    165164      READ   ( numnam_ice , namiceini ) 
    166165      ! 
    167       IF(lwp) THEN 
     166      IF(lwp) THEN                        ! control print 
    168167         WRITE(numout,*) 
    169168         WRITE(numout,*) 'lim_istate_init_2 : ice parameters inititialisation ' 
    … …  
    179178      ENDIF 
    180179 
    181       IF( ln_limini ) THEN                      ! Ice initialization using input file 
     180      IF( ln_limini ) THEN                ! Ice initialization using input file 
    182181         ! 
    183182         CALL iom_open( 'Ice_initialization.nc', inum_ice ) 
    … …  
    186185            IF(lwp) WRITE(numout,*) 
    187186            IF(lwp) WRITE(numout,*) '                  ice state initialization with : Ice_initialization.nc' 
    188              
     187            ! 
    189188            CALL iom_get( inum_ice, jpdom_data, 'hicif', hicif )       
    190189            CALL iom_get( inum_ice, jpdom_data, 'hsnif', hsnif )       
    … …  
    192191            CALL iom_get( inum_ice, jpdom_data, 'ts'   , sist  ) 
    193192            CALL iom_get( inum_ice, jpdom_unknown, 'tbif', tbif(1:nlci,1:nlcj,:),   & 
    194                  &        kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,jplayersp1 /) ) 
     193               &          kstart = (/ mig(1),mjg(1),1 /), kcount = (/ nlci,nlcj,jplayersp1 /) ) 
    195194            ! put some values in the extra-halo... 
    196195            DO jj = nlcj+1, jpj   ;   tbif(1:nlci,jj,:) = tbif(1:nlci,nlej,:)   ;   END DO 
    197196            DO ji = nlci+1, jpi   ;   tbif(ji    ,: ,:) = tbif(nlei  ,:   ,:)   ;   END DO 
    198  
     197            ! 
    199198            CALL iom_close( inum_ice) 
    200199            ! 
    … …  
    208207   !!   Default option :         Empty module      NO LIM 2.0 sea-ice model 
    209208   !!---------------------------------------------------------------------- 
    210 CONTAINS 
    211    SUBROUTINE lim_istate_2        ! Empty routine 
    212    END SUBROUTINE lim_istate_2 
    213209#endif 
    214210 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limmsh_2.F90

    r1694 r1855  
    44   !! LIM 2.0 ice model :   definition of the ice mesh parameters 
    55   !!====================================================================== 
     6   !! ** History :  LIM  ! 2001-04 (UCL)  original code 
     7   !!               2.0  ! 2002-08 (C. Ethe, G. Madec) F90, module 
     8   !!---------------------------------------------------------------------- 
    69#if defined key_lim2 
    710   !!---------------------------------------------------------------------- 
    … …  
    1013   !!   lim_msh_2   : definition of the ice mesh 
    1114   !!---------------------------------------------------------------------- 
    12    !! * Modules used 
    1315   USE phycst 
    1416   USE dom_oce 
    … …  
    2022   PRIVATE 
    2123 
    22    !! * Accessibility 
    23    PUBLIC lim_msh_2      ! routine called by ice_ini_2.F90 
    24  
    25    !!---------------------------------------------------------------------- 
    26    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     24   PUBLIC   lim_msh_2   ! routine called by ice_ini_2.F90 
     25 
     26   !!---------------------------------------------------------------------- 
     27   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    2728   !! $Id$ 
    28    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     29   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    2930   !!---------------------------------------------------------------------- 
    3031 
    … …  
    4243      !!              - Initialization of the ice masks (tmsk, umsk) 
    4344      !!  
    44       !! ** Refer.  : Deleersnijder et al. Ocean Modelling 100, 7-10  
    45       !! 
    46       !! ** History : 
    47       !!         original    : 01-04 (LIM) 
    48       !!         addition    : 02-08 (C. Ethe, G. Madec) 
     45      !! Reference  : Deleersnijder et al. Ocean Modelling 100, 7-10  
    4946      !!---------------------------------------------------------------------  
    50       !! * Local variables 
    51       INTEGER :: ji, jj      ! dummy loop indices 
    52  
    53       REAL(wp), DIMENSION(jpi,jpj) ::  & 
    54          zd2d1 , zd1d2       ! Derivative of zh2 (resp. zh1) in the x direction 
    55          !                   ! (resp. y direction) (defined at the center) 
    56       REAL(wp) ::         & 
    57          zh1p  , zh2p   , &  ! Idem zh1, zh2 for the bottom left corner of the grid 
    58          zd2d1p, zd1d2p , &  ! Idem zd2d1, zd1d2 for the bottom left corner of the grid 
    59          zusden, zusden2     ! temporary scalars 
     47      INTEGER  ::   ji, jj   ! dummy loop indices 
     48      REAL(wp) ::   zh1p, zd1d2p, zusden    ! temporary scalars 
     49      REAL(wp) ::   zh2p, zd2d1p, zusden2   !    -         - 
     50      REAL(wp), DIMENSION(jpi,jpj) ::   zd2d1 , zd1d2   ! 2D workspace 
    6051      !!--------------------------------------------------------------------- 
    6152 
    … …  
    112103      !------------------- 
    113104!!ibug ??? 
    114       akappa(:,:,:,:) = 0.e0 
    115       wght(:,:,:,:) = 0.e0 
     105      akappa(:,:,:,:)     = 0.e0 
     106      wght  (:,:,:,:)    = 0.e0 
    116107      alambd(:,:,:,:,:,:) = 0.e0 
    117       tmu(:,:) = 0.e0 
     108      tmu   (:,:)        = 0.e0 
    118109!!i 
    119110       
    … …  
    238229         END DO 
    239230      END DO 
    240        
    241       !--lateral boundary conditions     
    242       CALL lbc_lnk( tmu(:,:), 'I', 1. ) 
    243        
    244       ! unmasked and masked area of T-grid cell 
    245       area(:,:) = e1t(:,:) * e2t(:,:) 
     231      CALL lbc_lnk( tmu(:,:), 'I', 1. )      !--lateral boundary conditions     
     232       
     233      area(:,:) = e1t(:,:) * e2t(:,:)        ! unmasked and masked area of T-grid cell 
    246234       
    247235   END SUBROUTINE lim_msh_2 
    … …  
    251239   !!   Default option            Dummy Module         NO LIM sea-ice model 
    252240   !!---------------------------------------------------------------------- 
    253 CONTAINS 
    254    SUBROUTINE lim_msh_2           ! Dummy routine 
    255    END SUBROUTINE lim_msh_2 
    256241#endif 
    257242 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limrhg_2.F90

    r1774 r1855  
    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 
     6   !! History :  0.0  !  1993-12  (M.A. Morales Maqueda.)  Original code 
     7   !!            1.0  !  1994-12  (H. Goosse)  
     8   !!            2.0  !  1903-12  (C. Ethe, G. Madec)  F90, mpp 
     9   !!             -   !  2006-08  (G. Madec)  surface module, ice-stress at I-point 
     10   !!             -   !  2009-09  (G. Madec)  Huge verctor optimisation 
    1111   !!---------------------------------------------------------------------- 
    1212#if defined key_lim2 
    … …  
    1414   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    1515   !!---------------------------------------------------------------------- 
    16    !!---------------------------------------------------------------------- 
    17    !!   lim_rhg_2   : computes ice velocities 
     16   !!   lim_rhg_2     : computes ice velocities 
    1817   !!---------------------------------------------------------------------- 
    1918   USE par_oce        ! ocean parameter 
    … …  
    3433   PUBLIC   lim_rhg_2 ! routine called by lim_dyn 
    3534 
    36    REAL(wp) ::   rzero   = 0.e0   ! constant value: zero 
    37    REAL(wp) ::   rone    = 1.e0   !            and  one 
     35   REAL(wp) ::   rzero = 0.e0   ! constant value: zero 
     36   REAL(wp) ::   rone  = 1.e0   !            and  one 
    3837 
    3938   !! * Substitutions 
    4039#  include "vectopt_loop_substitute.h90" 
    4140   !!---------------------------------------------------------------------- 
    42    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     41   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    4342   !! $Id$ 
    4443   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    5655      !!  viscous-plastic law including shear strength and a bulk rheology. 
    5756      !! 
    58       !! ** Action  : - compute u_ice, v_ice the sea-ice velocity defined 
    59       !!              at I-point 
     57      !! ** Action  : - compute u_ice, v_ice the sea-ice velocity defined at I-point 
    6058      !!------------------------------------------------------------------- 
    6159      INTEGER, INTENT(in) ::   k_j1    ! southern j-index for ice computation 
    … …  
    582580   !!   Default option          Dummy module       NO 2.0 LIM sea-ice model 
    583581   !!---------------------------------------------------------------------- 
    584 CONTAINS 
    585    SUBROUTINE lim_rhg_2( k1 , k2 )       ! Dummy routine 
    586       WRITE(*,*) 'lim_rhg_2: You should not have seen this print! error?', k1, k2 
    587    END SUBROUTINE lim_rhg_2 
    588582#endif 
    589583 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limrst_2.F90

    r1715 r1855  
    22   !!====================================================================== 
    33   !!                     ***  MODULE  limrst_2  *** 
    4    !! Ice restart :  write the ice restart file 
     4   !! LIM-2 ice model : open/read/write the restart file 
    55   !!====================================================================== 
    6    !! History :  2.0  !  01-04  (C. Ethe, G. Madec)  Original code 
    7    !!                 !  06-07  (S. Masson)  use IOM for restart read/write 
     6   !! History :  2.0  ! 2001-04  (C. Ethe, G. Madec)  Original code 
     7   !!                 ! 2006-07  (S. Masson)  use IOM for restart read/write 
    88   !!---------------------------------------------------------------------- 
    99#if defined key_lim2 
    1010   !!---------------------------------------------------------------------- 
    1111   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    12    !!---------------------------------------------------------------------- 
    1312   !!---------------------------------------------------------------------- 
    1413   !!   lim_rst_opn_2   : open ice restart file 
    … …  
    2019   USE sbc_oce 
    2120   USE sbc_ice 
    22  
    2321   USE in_out_manager 
    2422   USE iom 
    … …  
    3533 
    3634   !!---------------------------------------------------------------------- 
    37    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     35   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    3836   !! $Id$ 
    3937   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    8482   END SUBROUTINE lim_rst_opn_2 
    8583 
     84 
    8685   SUBROUTINE lim_rst_write_2( kt ) 
    8786      !!---------------------------------------------------------------------- 
    … …  
    9190      !!---------------------------------------------------------------------- 
    9291      INTEGER, INTENT(in) ::   kt   ! number of iteration 
    93       ! 
     92      !! 
    9493      INTEGER ::   iter   ! kt + nn_fsbc -1 
    9594      !!---------------------------------------------------------------------- 
    … …  
    181180      ENDIF 
    182181 
    183       IF ( jprstlib == jprstdimg ) THEN 
     182      IF( jprstlib == jprstdimg ) THEN 
    184183        ! eventually read netcdf file (monobloc)  for restarting on different number of processors 
    185184        ! if {cn_icerst_in}.nc exists, then set jlibalt to jpnf90 
    … …  
    188187      ENDIF 
    189188 
    190       CALL iom_open ( cn_icerst_in, numrir, kiolib = jlibalt ) 
     189      CALL iom_open( cn_icerst_in, numrir, kiolib = jlibalt ) 
    191190 
    192191      CALL iom_get( numrir, 'kt_ice' , ziter ) 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limsbc_2.F90

    r1756 r1855  
    22   !!====================================================================== 
    33   !!                       ***  MODULE limsbc_2   *** 
    4    !!           computation of the flux at the sea ice/ocean interface 
     4   !!   LIM 2 ice model :   heat, salt, mass and momentum fluxes at the sea ice/ocean interface 
    55   !!====================================================================== 
    6    !! History : 00-01 (H. Goosse) Original code 
    7    !!           02-07 (C. Ethe, G. Madec) re-writing F90 
    8    !!           06-07 (G. Madec) surface module 
     6   !! History :  LIM  ! 2000-01 (H. Goosse) Original code 
     7   !!            2.0  ! 2002-07 (C. Ethe, G. Madec) re-writing F90 
     8   !!             -   ! 2006-07 (G. Madec) surface module 
    99   !!---------------------------------------------------------------------- 
    1010#if defined key_lim2 
    1111   !!---------------------------------------------------------------------- 
    1212   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    13    !!---------------------------------------------------------------------- 
    1413   !!---------------------------------------------------------------------- 
    1514   !!   lim_sbc_2  : flux at the ice / ocean interface 
    … …  
    2019   USE sbc_oce          ! surface boundary condition 
    2120   USE phycst           ! physical constants 
    22    USE ice_2            ! LIM sea-ice variables 
     21   USE ice_2            ! LIM 2 sea-ice variables 
    2322 
    2423   USE lbclnk           ! ocean lateral boundary condition 
    … …  
    3332   PRIVATE 
    3433 
    35    PUBLIC lim_sbc_2     ! called by sbc_ice_lim_2 
     34   PUBLIC   lim_sbc_2   ! called by sbc_ice_lim_2 
    3635 
    3736   REAL(wp)  ::   epsi16 = 1.e-16  ! constant values 
    … …  
    4443#  include "vectopt_loop_substitute.h90" 
    4544   !!---------------------------------------------------------------------- 
    46    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2006)  
     45   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    4746   !! $Id$ 
    4847   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    109108         sice_r(:,:) = sice 
    110109         ! 
    111          IF( cp_cfg == "orca"  .AND. jp_cfg == 2 ) THEN 
    112             !                                        ! ======================= 
    113             !                                        !  ORCA_R2 configuration 
    114             !                                        ! ======================= 
     110         IF( cp_cfg == "orca"  .AND. jp_cfg == 2 ) THEN             !  ORCA_R2 configuration 
    115111            ii0 = 145   ;   ii1 = 180        ! Baltic Sea 
    116112            ij0 = 113   ;   ij1 = 130   ;   soce_r(mi0(ii0):mi1(ii1) , mj0(ij0):mj1(ij1) ) = 4.e0 
    … …  
    194190      END DO 
    195191 
    196       CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:) )       
    197       CALL iom_put( 'qns_io_cea', qns(:,:) - zqnsoce(:,:) * pfrld(:,:) )       
    198       CALL iom_put( 'qsr_io_cea', fstric(:,:) * (1. - pfrld(:,:)) ) 
     192      CALL iom_put( 'hflx_ice_cea', - fdtcn(:,:)                         )       
     193      CALL iom_put( 'qns_io_cea'  , qns(:,:) - zqnsoce(:,:) * pfrld(:,:) )       
     194      CALL iom_put( 'qsr_io_cea'  , fstric(:,:) * (1. - pfrld(:,:))      ) 
    199195 
    200196      !------------------------------------------! 
    … …  
    302298      !-----------------------------------------------! 
    303299 
    304       IF ( lk_cpl ) THEN            
    305          ! Ice surface temperature  
     300      IF( lk_cpl ) THEN            
    306301         tn_ice(:,:,1) = sist(:,:)          ! sea-ice surface temperature        
    307          ! Computation of snow/ice and ocean albedo 
     302         !                                  ! snow/ice and ocean albedos 
    308303         CALL albedo_ice( tn_ice, reshape( hicif, (/jpi,jpj,1/) ), reshape( hsnif, (/jpi,jpj,1/) ), zalbp, zalb ) 
    309304         alb_ice(:,:,1) =  0.5 * ( zalbp(:,:,1) + zalb (:,:,1) )   ! Ice albedo (mean clear and overcast skys) 
    … …  
    312307 
    313308      IF(ln_ctl) THEN 
    314          CALL prt_ctl(tab2d_1=qsr   , clinfo1=' lim_sbc: qsr    : ', tab2d_2=qns   , clinfo2=' qns     : ') 
    315          CALL prt_ctl(tab2d_1=emp   , clinfo1=' lim_sbc: emp    : ', tab2d_2=emps  , clinfo2=' emps    : ') 
    316          CALL prt_ctl(tab2d_1=utau  , clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
    317             &         tab2d_2=vtau  , clinfo2=' vtau    : '        , mask2=vmask ) 
    318          CALL prt_ctl(tab2d_1=fr_i  , clinfo1=' lim_sbc: fr_i   : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice  : ') 
     309         CALL prt_ctl( tab2d_1=qsr , clinfo1=' lim_sbc: qsr    : ', tab2d_2=qns  , clinfo2=' qns     : ' ) 
     310         CALL prt_ctl( tab2d_1=emp , clinfo1=' lim_sbc: emp    : ', tab2d_2=emps , clinfo2=' emps    : ' ) 
     311         CALL prt_ctl( tab2d_1=utau, clinfo1=' lim_sbc: utau   : ', mask1=umask,   & 
     312            &          tab2d_2=vtau, clinfo2=' vtau    : '        , mask2=vmask ) 
     313         CALL prt_ctl( tab2d_1=fr_i, clinfo1=' lim_sbc: fr_i   : ', tab2d_2=tn_ice(:,:,1), clinfo2=' tn_ice  : ' ) 
    319314      ENDIF  
    320     
    321     END SUBROUTINE lim_sbc_2 
     315      ! 
     316   END SUBROUTINE lim_sbc_2 
    322317 
    323318#else 
    … …  
    325320   !!   Default option :        Dummy module       NO LIM 2.0 sea-ice model 
    326321   !!---------------------------------------------------------------------- 
    327 CONTAINS 
    328    SUBROUTINE lim_sbc_2         ! Dummy routine 
    329    END SUBROUTINE lim_sbc_2 
    330322#endif  
    331323 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limtab_2.F90

    r1156 r1855  
    22   !!====================================================================== 
    33   !!                       ***  MODULE limtab_2   *** 
    4    !!              transform 1D (2D) array to a 2D (1D) table 
     4   !!   LIM 2 ice model : transform 1D (2D) array to a 2D (1D) array 
    55   !!====================================================================== 
    66#if defined key_lim2 
    77   !!---------------------------------------------------------------------- 
    8    !!   tab_2d_1d  : 2-D to 1-D 
    9    !!   tab_1d_2d  : 1-D to 2-D 
     8   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    109   !!---------------------------------------------------------------------- 
    11    !! * Modules used 
     10   !!   tab_2d_1d_2  : 2-D to 1-D 
     11   !!   tab_1d_2d_2  : 1-D to 2-D 
     12   !!---------------------------------------------------------------------- 
    1213   USE par_kind 
    1314 
    … …  
    1516   PRIVATE 
    1617 
    17    !! * Routine accessibility 
    18    PUBLIC tab_2d_1d_2  ! called by lim_ther 
    19    PUBLIC tab_1d_2d_2  ! called by lim_ther 
     18   PUBLIC   tab_2d_1d_2   ! called by lim_thd_2 
     19   PUBLIC   tab_1d_2d_2   ! called by lim_thd_2 
    2020 
    2121   !!---------------------------------------------------------------------- 
    22    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     22   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    2323   !! $Id$ 
    2424   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    2626CONTAINS 
    2727 
    28    SUBROUTINE tab_2d_1d_2 ( ndim1d, tab1d, tab2d, ndim2d_x, ndim2d_y, tab_ind ) 
    29  
    30       INTEGER, INTENT(in) :: & 
    31          ndim1d, ndim2d_x, ndim2d_y 
    32  
    33       REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT(in) ::  & 
    34          tab2d 
    35  
    36       INTEGER, DIMENSION ( ndim1d), INTENT ( in) :: & 
    37          tab_ind 
    38  
    39       REAL(wp), DIMENSION(ndim1d), INTENT ( out) ::  &  
    40          tab1d 
    41  
    42       INTEGER ::  & 
    43          jn , jid, jjd 
    44          
    45       DO jn = 1, ndim1d 
    46          jid        = MOD( tab_ind(jn) - 1, ndim2d_x ) + 1 
    47          jjd        = ( tab_ind(jn) - 1 ) / ndim2d_x + 1 
    48          tab1d( jn) = tab2d( jid, jjd) 
     28   SUBROUTINE tab_2d_1d_2( kdim1d, ptab1d, ptab2d, kdim2d_i, kdim2d_j, ptab_ind ) 
     29      !!------------------------------------------------------------------- 
     30      INTEGER , INTENT(in   )                                ::   kdim1d 
     31      INTEGER , INTENT(in   )                                ::   kdim2d_i, kdim2d_j 
     32      REAL(wp), INTENT(in   ), DIMENSION(kdim2d_i, kdim2d_j) ::   ptab2d 
     33      INTEGER , INTENT(in   ), DIMENSION(kdim1d)             ::   ptab_ind 
     34      REAL(wp), INTENT(  out), DIMENSION(kdim1d)             ::   ptab1d 
     35      !! 
     36      INTEGER ::   jn , jid, jjd   ! dummy loop indices 
     37      !!------------------------------------------------------------------- 
     38      !  
     39      DO jn = 1, kdim1d 
     40         jid = MOD( ptab_ind(jn) - 1, kdim2d_i ) + 1 
     41         jjd = ( ptab_ind(jn) - 1 ) / kdim2d_i + 1 
     42         ptab1d(jn) = ptab2d(jid,jjd) 
    4943      END DO  
    50  
     44      ! 
    5145   END SUBROUTINE tab_2d_1d_2 
    5246 
    5347 
    54    SUBROUTINE tab_1d_2d_2 ( ndim1d, tab2d, tab_ind, tab1d, ndim2d_x, ndim2d_y ) 
    55  
    56       INTEGER, INTENT ( in) :: & 
    57          ndim1d, ndim2d_x, ndim2d_y 
    58  
    59       INTEGER, DIMENSION (ndim1d) , INTENT (in) :: & 
    60          tab_ind 
    61  
    62       REAL(wp), DIMENSION(ndim1d), INTENT (in) ::  & 
    63          tab1d   
    64  
    65       REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT ( out) :: & 
    66          tab2d 
    67  
    68       INTEGER :: & 
    69          jn, jid, jjd 
    70  
    71       DO jn = 1, ndim1d 
    72          jid             = MOD( tab_ind(jn) - 1, ndim2d_x) + 1 
    73          jjd             =    ( tab_ind(jn) - 1 ) / ndim2d_x  + 1 
    74          tab2d(jid, jjd) = tab1d( jn) 
     48   SUBROUTINE tab_1d_2d_2( kdim1d, ptab2d, ptab_ind, ptab1d, kdim2d_i, kdim2d_j ) 
     49      !!------------------------------------------------------------------- 
     50      INTEGER , INTENT(in   )                                ::   kdim1d 
     51      INTEGER , INTENT(in   )                                ::   kdim2d_i, kdim2d_j 
     52      INTEGER , INTENT(in   ), DIMENSION(kdim1d)             ::   ptab_ind 
     53      REAL(wp), INTENT(in   ), DIMENSION(kdim1d)             ::   ptab1d 
     54      REAL(wp), INTENT(  out), DIMENSION(kdim2d_i, kdim2d_j) ::   ptab2d 
     55      !! 
     56      INTEGER ::   jn, jid, jjd   ! dummy loop indices 
     57      !!------------------------------------------------------------------- 
     58      ! 
     59      DO jn = 1, kdim1d 
     60         jid = MOD( ptab_ind(jn) - 1, kdim2d_i) + 1 
     61         jjd =    ( ptab_ind(jn) - 1 ) / kdim2d_i  + 1 
     62         ptab2d(jid,jjd) = ptab1d(jn) 
    7563      END DO 
    76  
     64      ! 
    7765   END SUBROUTINE tab_1d_2d_2 
    7866 
     67#else 
     68   !!---------------------------------------------------------------------- 
     69   !!   Default option          Dummy module       NO LIM 2.0 sea-ice model 
     70   !!---------------------------------------------------------------------- 
    7971#endif 
     72 
     73   !!====================================================================== 
    8074END MODULE limtab_2 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limthd_2.F90

    r1758 r1855  
    22   !!====================================================================== 
    33   !!                  ***  MODULE limthd_2   *** 
    4    !!              LIM thermo ice model : ice thermodynamic 
     4   !!   LIM 2 ice model : thermodynamics 
    55   !!====================================================================== 
    6    !! History :  1.0  ! 2000-01 (LIM) 
     6   !! History :  LIM  ! 2000-01 (UCL) Original code 
    77   !!            2.0  ! 2002-07 (C. Ethe, G. Madec) F90 
    8    !!            2.0  ! 2003-08 (C. Ethe)  add lim_thd_init 
    9    !!             -   ! 2008-2008  (A. Caubel, G. Madec, E. Maisonnave, S. Masson ) generic coupled interface 
     8   !!             -   ! 2003-08 (C. Ethe)  add lim_thd_init 
     9   !!             -   ! 2008-08  (A. Caubel, G. Madec, E. Maisonnave, S. Masson ) generic coupled interface 
    1010   !!--------------------------------------------------------------------- 
    1111#if defined key_lim2 
    … …  
    1313   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
    1414   !!---------------------------------------------------------------------- 
    15    !!   lim_thd_2      : thermodynamic of sea ice 
     15   !!   lim_thd_2      : thermodynamics of sea ice 
    1616   !!   lim_thd_init_2 : initialisation of sea-ice thermodynamic 
    1717   !!---------------------------------------------------------------------- 
    … …  
    2323   USE lib_mpp 
    2424   USE iom             ! IOM library 
    25    USE ice_2           ! LIM sea-ice variables 
     25   USE ice_2           ! LIM 2 sea-ice variables 
    2626   USE sbc_oce         !  
    2727   USE sbc_ice         !  
    28    USE thd_ice_2       ! LIM thermodynamic sea-ice variables 
    29    USE dom_ice_2       ! LIM sea-ice domain 
     28   USE dom_ice_2       ! LIM 2 sea-ice domain 
    3029   USE limthd_zdf_2 
    3130   USE limthd_lac_2 
    … …  
    4948#  include "domzgr_substitute.h90" 
    5049#  include "vectopt_loop_substitute.h90" 
    51    !!-------- ------------------------------------------------------------- 
    52    !! NEMO/LIM 3.2,  UCL-LOCEAN-IPSL (2009)  
     50   !!---------------------------------------------------------------------- 
     51   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    5352   !! $Id$ 
    5453   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    566565   !!   Default option          Dummy module       NO LIM 2.0 sea-ice model 
    567566   !!---------------------------------------------------------------------- 
    568 CONTAINS 
    569    SUBROUTINE lim_thd_2         ! Dummy routine 
    570    END SUBROUTINE lim_thd_2 
    571567#endif 
    572568 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limthd_lac_2.F90

    r1697 r1855  
    11MODULE limthd_lac_2 
    2 #if defined key_lim2 
    32   !!====================================================================== 
    43   !!                       ***  MODULE limthd_lac_2   *** 
    5    !!                lateral thermodynamic growth of the ice  
     4   !!   LIM 2 ice model :   thermodynamics -- lateral accretion of the ice  
    65   !!====================================================================== 
    7  
    8    !!---------------------------------------------------------------------- 
    9    !!   lim_lat_acr_2    : lateral accretion of ice 
    10    !! * Modules used 
     6   !! History :  LIM  !  2001-04 (UCL)  original code 
     7   !!            2.0  !  2002-08 (C. Ethe, G. Madec)  F90, mpp 
     8   !!---------------------------------------------------------------------- 
     9#if defined key_lim2 
     10   !!---------------------------------------------------------------------- 
     11   !!   'key_lim2' :                                  LIM 2.0 sea-ice model 
     12   !!---------------------------------------------------------------------- 
     13   !!   lim_lat_acr_2   : lateral accretion of ice 
     14   !!---------------------------------------------------------------------- 
    1115   USE par_oce          ! ocean parameters 
    12    USE phycst 
    13    USE thd_ice_2 
    14    USE ice_2 
    15    USE limistate_2  
     16   USE phycst           ! physical constants 
     17   USE ice_2            ! LIM 2 sea-ice variables 
     18   USE limistate_2      ! LIM 2 initial state 
    1619      
    1720   IMPLICIT NONE 
    1821   PRIVATE 
    1922 
    20    !! * Routine accessibility 
    21    PUBLIC lim_thd_lac_2   ! called by lim_thd_2 
    22  
    23    !! * Module variables 
    24    REAL(wp)  ::           &  ! constant values 
    25       epsi20 = 1.e-20  ,  & 
    26       epsi13 = 1.e-13  ,  & 
    27       zzero  = 0.e0    ,  & 
    28       zone   = 1.e0 
    29    !!---------------------------------------------------------------------- 
    30    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     23   PUBLIC   lim_thd_lac_2   ! called by lim_thd_2 
     24 
     25   REAL(wp) ::   epsi20 = 1.e-20   ! constant values 
     26   REAL(wp) ::   epsi13 = 1.e-13   ! 
     27   REAL(wp) ::   rzero  = 0.e0     ! 
     28   REAL(wp) ::   rone   = 1.e0     ! 
     29 
     30   !!---------------------------------------------------------------------- 
     31   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    3132   !! $Id$ 
    32    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     33   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    3334   !!---------------------------------------------------------------------- 
    3435CONTAINS 
    … …  
    6061      !!               update h_snow_1d, h_ice_1d and tbif_1d(:,:)       
    6162      !!  
    62       !! ** References : 
    63       !!      M. Maqueda, 1995, PhD Thesis, Univesidad Complutense Madrid 
    64       !!      Fichefet T. and M. Maqueda 1997, J. Geo. Res., 102(C6),  
    65       !!                                                12609 -12646    
    66       !! History : 
    67       !!   1.0  !  01-04 (LIM)  original code 
    68       !!   2.0  !  02-08 (C. Ethe, G. Madec)  F90, mpp 
     63      !! References :   M. Maqueda, 1995, PhD Thesis, Univesidad Complutense Madrid 
     64      !!                Fichefet T. and M. Maqueda 1997, J. Geo. Res., 102(C6), 12609 -12646    
    6965      !!------------------------------------------------------------------- 
    70       !! * Arguments 
    71       INTEGER , INTENT(IN)::  & 
    72          kideb          ,   &  ! start point on which the the computation is applied 
    73          kiut                  ! end point on which the the computation is applied 
    74  
    75       !! * Local variables 
    76       INTEGER ::            & 
    77          ji             ,   &  !  dummy loop indices 
    78          iicefr         ,   &  !  1 = existing ice ; 0 = no ice 
    79          iiceform       ,   &  !  1 = ice formed   ; 0 = no ice formed 
    80          ihemis                !  dummy indice 
    81       REAL(wp), DIMENSION(jpij) :: & 
    82          zqbgow           ,  &  !  heat budget of the open water (negative) 
    83          zfrl_old         ,  &  !  previous sea/ice fraction 
    84          zhice_old        ,  &  !  previous ice thickness 
    85          zhice0           ,  &  !  thickness of newly formed ice in leads 
    86          zfrlmin          ,  &  !  minimum fraction for leads 
    87          zdhicbot               !  part of thickness of newly formed ice in leads which  
    88                                 !  has been already used in transport for example 
    89       REAL(wp)  ::  & 
    90          zhemis           ,  &  !  hemisphere (0 = North, 1 = South) 
    91          zhicenew         ,  &  !  new ice thickness 
    92          zholds2          ,  &  !  ratio of previous ice thickness and 2  
    93          zhnews2          ,  &  !  ratio of new ice thickness and 2  
    94          zfrlnew          ,  &  !  new sea/ice fraction 
    95          zfrld            ,  &  !  ratio of sea/ice fraction and minimum fraction for leads 
    96          zfrrate          ,  &  !  leads-closure rate 
    97          zdfrl                  !  sea-ice fraction increment 
    98       REAL(wp)  ::  & 
    99          zdh1 , zdh2 , zdh3 , zdh4, zdh5   , &   ! tempory scalars 
    100          ztint , zta1 , zta2 , zta3 , zta4 , & 
    101          zah, zalpha , zbeta 
     66      INTEGER, INTENT(in) ::   kideb   ! start point on which the the computation is applied 
     67      INTEGER, INTENT(in) ::   kiut    ! end point on which the the computation is applied 
     68      !! 
     69      INTEGER ::   ji         !  dummy loop indices 
     70      INTEGER ::   iicefr     !  1 = existing ice ; 0 = no ice 
     71      INTEGER ::   iiceform   !  1 = ice formed   ; 0 = no ice formed 
     72      INTEGER ::   ihemis     !  hemisphere indicator  
     73      REAL(wp), DIMENSION(jpij) ::   zqbgow      !  heat budget of the open water (negative) 
     74      REAL(wp), DIMENSION(jpij) ::   zfrl_old    !  previous sea/ice fraction 
     75      REAL(wp), DIMENSION(jpij) ::   zhice_old   !  previous ice thickness 
     76      REAL(wp), DIMENSION(jpij) ::   zhice0      !  thickness of newly formed ice in leads 
     77      REAL(wp), DIMENSION(jpij) ::   zfrlmin     !  minimum fraction for leads 
     78      REAL(wp), DIMENSION(jpij) ::   zdhicbot    !  part of thickness of newly formed ice in leads which  
     79      !                                          !  has been already used in transport for example 
     80      REAL(wp) ::   zhemis     !  hemisphere (0 = North, 1 = South) 
     81      REAL(wp) ::   zhicenew   !  new ice thickness 
     82      REAL(wp) ::   zholds2    !  ratio of previous ice thickness and 2  
     83      REAL(wp) ::   zhnews2    !  ratio of new ice thickness and 2  
     84      REAL(wp) ::   zfrlnew    !  new sea/ice fraction 
     85      REAL(wp) ::   zfrld      !  ratio of sea/ice fraction and minimum fraction for leads 
     86      REAL(wp) ::   zfrrate    !  leads-closure rate 
     87      REAL(wp) ::   zdfrl      !  sea-ice fraction increment 
     88      REAL(wp) ::   zdh1 , zdh2, zdh3, zdh4, zdh5   ! tempory scalars 
     89      REAL(wp) ::   ztint, zta1, zta2, zta3, zta4   !    -       - 
     90      REAL(wp) ::   zah, zalpha, zbeta              !    -       - 
    10291      !!---------------------------------------------------------------------       
    10392                    
    10493      !-------------------------------------------------------------- 
    10594      !   Computation of the heat budget of the open water (negative) 
    106       !-------------------------------------------------------------- 
    107        
     95      !--------------------------------------------------------------  
    10896      DO ji = kideb , kiut       
    10997         zqbgow(ji) = qldif_1d(ji) - qcmif_1d(ji) 
    … …  
    131119      !     Morales Maqueda, 1995 - Fichefet and Morales Maqueda, 1997 
    132120      !--------------------------------------------------------------------- 
    133        
    134121!CDIR NOVERRCHK 
    135122      DO ji = kideb , kiut 
    … …  
    145132         !--computation of the remaining part of ice thickness which has been already used 
    146133         zdhicbot(ji) =  ( frld_1d(ji) - zfrlnew ) * zhice0(ji) / ( 1.0 - zfrlmin(ji) ) &  
    147                       -  (  ( 1.0 - zfrrate ) / ( 1.0 - frld_1d(ji) ) )  * ( zqbgow(ji) / xlic )  
     134            &         -  (  ( 1.0 - zfrrate ) / ( 1.0 - frld_1d(ji) ) )  * ( zqbgow(ji) / xlic )  
    148135      END DO 
    149136  
    … …  
    155142      !         (1-Anew) * hsnew = (1-Aold) * hsold             
    156143      !-------------------------------------------------------------------------------------------- 
    157        
    158144      DO ji = kideb , kiut 
    159145         iicefr       = 1 - MAX( 0, INT( SIGN( 1.5 * zone , zfrl_old(ji) - 1.0 + epsi13 ) ) ) 
    … …  
    169155      !   Ajustement of ice internal temperatures 
    170156      !------------------------------------------------------- 
    171        
    172157      DO ji = kideb , kiut 
    173158         iicefr      = 1 - MAX( 0, INT( SIGN( 1.5 * zone , zfrl_old(ji) - 1.0 + epsi13 ) ) ) 
    … …  
    213198      !           dV = Vnew - Vold 
    214199      !------------------------------------------------------------- 
    215        
    216200      DO ji = kideb , kiut 
    217201         dvlbq_1d  (ji) = ( 1. - frld_1d(ji) ) * h_ice_1d(ji) - ( 1. - zfrl_old(ji) ) * zhice_old(ji) 
    218202         rdmicif_1d(ji) = rdmicif_1d(ji) + rhoic * dvlbq_1d(ji) 
    219203      END DO 
    220        
     204      ! 
    221205   END SUBROUTINE lim_thd_lac_2 
     206    
    222207#else 
     208   !!---------------------------------------------------------------------- 
     209   !!   Default option          Dummy module       NO LIM 2.0 sea-ice model 
     210   !!---------------------------------------------------------------------- 
     211#endif 
     212 
    223213   !!====================================================================== 
    224    !!                       ***  MODULE limthd_lac_2   *** 
    225    !!                           no sea ice model 
    226    !!====================================================================== 
    227 CONTAINS 
    228    SUBROUTINE lim_thd_lac_2           ! Empty routine 
    229    END SUBROUTINE lim_thd_lac_2 
    230 #endif 
    231214END MODULE limthd_lac_2 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limthd_zdf_2.F90

    r1756 r1855  
    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 
     6   !! History :  1.0  !  2001-04 (LIM) Original code 
     7   !!            2.0  !  2002-08 (C. Ethe, G. Madec) F90 
    88   !!---------------------------------------------------------------------- 
    99#if defined key_lim2 
    … …  
    1111   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    1212   !!---------------------------------------------------------------------- 
    13    !!---------------------------------------------------------------------- 
    1413   !!   lim_thd_zdf_2 : vertical accr./abl. and lateral ablation of sea ice 
    1514   !!---------------------------------------------------------------------- 
    16    !! * Modules used 
    1715   USE par_oce          ! ocean parameters 
    18    USE phycst           ! ??? 
    19    USE thd_ice_2 
    20    USE ice_2 
    21    USE limistate_2 
    22    USE in_out_manager 
     16   USE phycst           ! physical constants 
     17   USE ice_2            ! LIM 2 sea-ice variables 
     18   USE limistate_2      ! LIM 2 initial state 
     19   USE in_out_manager   ! I/O manager 
     20   !! 
    2321   USE cpl_oasis3, ONLY : lk_cpl 
    2422       
    … …  
    2624   PRIVATE 
    2725 
    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 
     26   PUBLIC   lim_thd_zdf_2   ! called by lim_thd_2 
     27 
     28   REAL(wp) ::   epsi20 = 1.e-20   ! constant values 
     29   REAL(wp) ::   epsi13 = 1.e-13   ! 
     30   REAL(wp) ::   rzero  = 0.e0     ! 
     31   REAL(wp) ::   rone   = 1.e0     ! 
     32 
    3433   !!---------------------------------------------------------------------- 
    35    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     34   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    3635   !! $Id$ 
    3736   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    7372      !! 
    7473      INTEGER ::   ji       ! dummy loop indices 
    75       REAL(wp), DIMENSION(jpij,2) ::   zqcmlt        ! energy due to surface( /1 ) and bottom melting( /2 ) 
    76       REAL(wp), DIMENSION(jpij) ::  & 
    77          ztsmlt      &    ! snow/ice surface melting temperature 
    78          ,ztbif      &    ! int. temp. at the mid-point of the 1st layer of the snow/ice sys.  
    79          ,zksn       &    ! effective conductivity of snow 
    80          ,zkic       &    ! effective conductivity of ice 
    81          ,zksndh     &    ! thermal cond. at the mid-point of the 1st layer of the snow/ice sys.  
    82          , zfcsu     &    ! conductive heat flux at the surface of the snow/ice system  
    83          , zfcsudt   &    ! = zfcsu * dt 
    84          , zi0       &    ! frac. of the net SW rad. which is not absorbed at the surface 
    85          , z1mi0     &    ! fraction of the net SW radiation absorbed at the surface 
    86          , zqmax     &    ! maximum energy stored in brine pockets 
    87          , zrcpdt    &    ! h_su*rho_su*cp_su/dt(h_su being the thick. of surf. layer) 
    88          , zts_old   &    ! previous surface temperature 
    89          , zidsn , z1midsn , zidsnic ! tempory variables 
    90       REAL(wp), DIMENSION(jpij) ::   & 
    91           zfnet       &  ! net heat flux at the top surface( incl. conductive heat flux) 
    92           , zsprecip  &    ! snow accumulation 
    93           , zhsnw_old &    ! previous snow thickness 
    94           , zdhictop  &    ! change in ice thickness due to top surf ablation/accretion 
    95           , zdhicbot  &    ! change in ice thickness due to bottom surf abl/acc 
    96           , zqsup     &    ! energy transmitted to ocean (coming from top surf abl/acc) 
    97           , zqocea    &    ! energy transmitted to ocean (coming from bottom sur abl/acc) 
    98           , zfrl_old  &    ! previous sea/ice fraction 
    99           , zfrld_1d    &    ! new sea/ice fraction 
    100           , zep            ! internal temperature of the 2nd layer of the snow/ice system 
    101        REAL(wp), DIMENSION(3) :: &  
    102           zplediag  &    ! principle diagonal, subdiag. and supdiag. of the  
    103           , zsubdiag  &    ! tri-diagonal matrix coming from the computation 
    104           , zsupdiag  &    ! of the temperatures inside the snow-ice system 
    105           , zsmbr          ! second member 
    106        REAL(wp) :: &  
    107           zhsu     &     ! thickness of surface layer 
    108           , zhe      &     ! effective thickness for compu. of equ. thermal conductivity 
    109           , zheshth  &     ! = zhe / thth 
    110           , zghe     &     ! correction factor of the thermal conductivity 
    111           , zumsb    &     ! parameter for numerical method to solve heat-diffusion eq. 
    112           , zkhsn    &     ! conductivity at the snow layer 
    113           , zkhic    &     ! conductivity at the ice layers 
    114           , zkint    &     ! equivalent conductivity at the snow-ice interface 
    115           , zkhsnint &     ! = zkint*dt / (hsn*rhosn*cpsn)   
    116           , zkhicint &     ! = 2*zkint*dt / (hic*rhoic*cpic) 
    117           , zpiv1 , zpiv2  &       ! tempory scalars used to solve the tri-diagonal system 
    118           , zb2 , zd2 , zb3 , zd3 & 
    119           , ztint          ! equivalent temperature at the snow-ice interface 
    120        REAL(wp) :: &  
    121           zexp      &     ! exponential function of the ice thickness 
    122           , zfsab     &     ! part of solar radiation stored in brine pockets 
    123           , zfts      &     ! value of energy balance function when the temp. equal surf. temp. 
    124           , zdfts     &     ! value of derivative of ztfs when the temp. equal surf. temp. 
    125           , zdts      &     ! surface temperature increment 
    126           , zqsnw_mlt &     ! energy needed to melt snow 
    127           , zdhsmlt   &     ! change in snow thickness due to melt 
    128           , zhsn      &     ! snow thickness (previous+accumulation-melt) 
    129           , zqsn_mlt_rem &  ! remaining heat coming from snow melting 
    130           , zqice_top_mlt & ! energy used to melt ice at top surface 
    131           , zdhssub      &  ! change in snow thick. due to sublimation or evaporation 
    132           , zdhisub      &  ! change in ice thick. due to sublimation or evaporation     
    133           , zdhsn        &  ! snow ice thickness increment 
    134           , zdtsn        &  ! snow internal temp. increment 
    135           , zdtic        &  ! ice internal temp. increment 
    136           , zqnes          ! conductive energy due to ice melting in the first ice layer 
    137        REAL(wp) :: &  
    138           ztbot     &      ! temperature at the bottom surface 
    139           , zfcbot    &      ! conductive heat flux at bottom surface 
    140           , zqice_bot &      ! energy used for bottom melting/growing 
    141           , zqice_bot_mlt &  ! energy used for bottom melting 
    142           , zqstbif_bot  &  ! part of energy stored in brine pockets used for bottom melting 
    143           , zqstbif_old  &  ! tempory var. for zqstbif_bot 
    144           , zdhicmlt      &  ! change in ice thickness due to bottom melting 
    145           , zdhicm        &  ! change in ice thickness var.  
    146           , zdhsnm        &  ! change in snow thickness var.  
    147           , zhsnfi        &  ! snow thickness var.  
    148           , zc1, zpc1, zc2, zpc2, zp1, zp2 & ! tempory variables 
    149           , ztb2, ztb3 
    150        REAL(wp) :: &  
    151           zdrmh         &   ! change in snow/ice thick. after snow-ice formation 
    152           , zhicnew       &   ! new ice thickness 
    153           , zhsnnew       &   ! new snow thickness 
    154           , zquot , ztneq &   ! tempory temp. variables 
    155           , zqice, zqicetot & ! total heat inside the snow/ice system 
    156           , zdfrl         &   ! change in ice concentration 
    157           , zdvsnvol      &   ! change in snow volume 
    158           , zdrfrl1, zdrfrl2 &  ! tempory scalars 
    159           , zihsn, zidhb, zihic, zihe, zihq, ziexp, ziqf, zihnf, zibmlt, ziqr, zihgnew, zind 
    160        !!---------------------------------------------------------------------- 
    161  
    162        !----------------------------------------------------------------------- 
    163        !  1. Boundaries conditions for snow/ice system internal temperature 
    164        !       - If tbif_1d(ji,1) > rt0_snow, tbif_1d(ji,1) = rt0_snow  
    165        !       - If tbif_1d(ji,2/3) > rt0_ice, tbif_1d(ji,2/3) = rt0_ice  
    166        !     Computation of energies due to surface and bottom melting  
    167        !----------------------------------------------------------------------- 
    168         
    169        DO ji = kideb , kiut 
    170           zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 
    171           zihic = MAX( zzero , SIGN( zone , hicdif - h_ice_1d(ji) ) ) 
    172           !--computation of energy due to surface melting 
    173           zqcmlt(ji,1) = ( MAX ( zzero ,  & 
    174              &                   rcpsn * h_snow_1d(ji) * ( tbif_1d(ji,1) - rt0_snow ) ) ) * ( 1.0 - zihsn ) 
    175           !--computation of energy due to bottom melting 
    176           zqcmlt(ji,2) = ( MAX( zzero , & 
    177              &                  rcpic * ( tbif_1d(ji,2) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 
    178              &           + MAX( zzero , & 
    179              &                  rcpic * ( tbif_1d(ji,3) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 
    180              &           ) * ( 1.0 - zihic  ) 
    181           !--limitation of  snow/ice system internal temperature 
    182           tbif_1d(ji,1)   = MIN( rt0_snow, tbif_1d(ji,1) ) 
    183           tbif_1d(ji,2)   = MIN( rt0_ice , tbif_1d(ji,2) ) 
    184           tbif_1d(ji,3)   = MIN( rt0_ice , tbif_1d(ji,3) ) 
    185        END DO 
    186  
    187        !------------------------------------------- 
    188        !  2. Calculate some intermediate variables.   
    189        !------------------------------------------- 
    190         
    191        ! initialisation of the thickness of surface layer 
    192        zhsu = hnzst   
    193  
    194        DO ji = kideb , kiut 
    195           zind   = MAX( zzero , SIGN( zone , zhsu - h_snow_1d(ji) ) ) 
    196           zihsn  = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 
    197           zihsn  = MAX( zihsn , zind ) 
    198           zihic  = MAX( zzero , sign( zone , hicdif - h_ice_1d(ji) ) ) 
    199           !     2.1. Computation of surface melting temperature 
    200           !---------------------------------------------------- 
    201           zind  = MAX( zzero , SIGN( zone , -h_snow_1d(ji) ) ) 
    202           ztsmlt(ji) = ( 1.0 - zind ) * rt0_snow + zind * rt0_ice 
    203           ! 
    204           !     2.2. Effective conductivity of snow and ice 
    205           !----------------------------------------------- 
    206  
    207           !---computation of the correction factor on the thermal conductivity 
    208           !-- (Morales Maqueda, 1995 ; Fichefet and Morales Maqueda, 1997) 
    209           zhe      =  ( rcdsn / ( rcdsn + rcdic ) ) * h_ice_1d(ji)   & 
    210              &     + ( rcdic / ( rcdsn + rcdic ) ) * h_snow_1d(ji)  
    211           zihe     = MAX( zzero , SIGN( zone , 2.0 * zhe - thth ) ) 
    212           zheshth  = zhe / thth 
    213           zghe     = ( 1.0 - zihe ) * zheshth * ( 2.0 - zheshth )   & 
    214              &     +         zihe   * 0.5 * ( 1.5 + LOG( 2.0 * zheshth ) ) 
    215  
    216           !---effective conductivities  
    217           zksn(ji)  = zghe * rcdsn   
    218           zkic(ji)  = zghe * rcdic 
    219  
    220           ! 
    221           !     2.3. Computation of the conductive heat flux from the snow/ice  
    222           !          system interior toward the top surface 
    223           !------------------------------------------------------------------ 
    224  
    225           !---Thermal conductivity at the mid-point of the first snow/ice system layer 
    226           zksndh(ji) =   ( ( 1.0 - zihsn ) * 2.0 * zksn(ji) + zihsn * 4.0 * zkic(ji) )   & 
    227              &         / ( ( 1.0 - zihsn ) *  h_snow_1d(ji)                              & 
    228              &           +        zihsn   *  ( ( 1.0 + 3.0 * zihic ) * h_ice_1d(ji)      & 
    229              &           + 4.0 * zkic(ji)/zksn(ji) * h_snow_1d(ji) ) ) 
    230  
    231           !---internal temperature at the mid-point of the first snow/ice system layer 
    232           ztbif(ji)  = ( 1.0 - zihsn ) * tbif_1d(ji,1)                       & 
    233              &       +         zihsn   * ( ( 1.0 - zihic ) * tbif_1d(ji,2)   & 
    234              &       +         zihic   * tfu_1d(ji)   ) 
    235           !---conductive heat flux  
    236           zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 
    237  
    238        END DO 
    239  
    240        !-------------------------------------------------------------------- 
    241        !  3. Calculate :  
    242        !     - fstbif_1d, part of solar radiation absorbing inside the ice 
    243        !       assuming an exponential absorption (Grenfell and Maykut, 1977) 
    244        !     - zqmax,  maximum energy stored in brine pockets 
    245        !     - qstbif_1d, total energy stored in brine pockets (updating) 
    246        !------------------------------------------------------------------- 
    247  
    248        DO ji = kideb , kiut 
    249           zihsn  = MAX( zzero , SIGN (zone , -h_snow_1d(ji) ) ) 
    250           zihic  = MAX( zzero , 1.0 - ( h_ice_1d(ji) / zhsu ) )      
    251           zind   = MAX( zzero , SIGN (zone , hicdif - h_ice_1d(ji) ) ) 
    252           !--Computation of the fraction of the net shortwave radiation which 
    253           !--penetrates inside the ice cover ( See Forcat) 
    254           zi0(ji)  = zihsn * ( fr1_i0_1d(ji) + zihic * fr2_i0_1d(ji) ) 
    255           zexp     = MIN( zone , EXP( -1.5 * ( h_ice_1d(ji) - zhsu ) ) ) 
    256           fstbif_1d(ji) = zi0(ji) * qsr_ice_1d(ji) * zexp 
    257           !--Computation of maximum energy stored in brine pockets zqmax and update 
    258           !--the total energy stored in brine pockets, if less than zqmax 
    259           zqmax(ji) = MAX( zzero , 0.5 * xlic * ( h_ice_1d(ji) - hicmin ) ) 
    260           zfsab   = zi0(ji) * qsr_ice_1d(ji) * ( 1.0 - zexp ) 
    261           zihq    = ( 1.0 - zind ) * MAX(zzero, SIGN( zone , qstbif_1d(ji) - zqmax(ji) ) ) & 
    262              &    +         zind   * zone 
    263           qstbif_1d(ji) = ( qstbif_1d(ji) + ( 1.0 - zihq ) * zfsab * rdt_ice ) * swiqst 
    264           !--fraction of shortwave radiation absorbed at surface 
    265           ziexp = zihq * zexp + ( 1.0 - zihq ) * ( swiqst + ( 1.0 - swiqst ) * zexp ) 
    266           z1mi0(ji) = 1.0 - zi0(ji) * ziexp 
    267        END DO 
    268  
    269        !-------------------------------------------------------------------------------- 
    270        !  4. Computation of the surface temperature : determined by considering the  
    271        !     budget of a thin layer of thick. zhsu at the top surface (H. Grenier, 1995)  
    272        !     and based on a surface energy balance :  
    273        !     hsu * rcp * dT/dt = Fsr + Fnsr(T) + Fcs(T), 
    274        !     where - Fsr is the net absorbed solar radiation,  
    275        !           - Fnsr is the total non solar radiation (incoming and outgoing long-wave, 
    276        !             sensible and latent heat fluxes) 
    277        !           - Fcs the conductive heat flux at the top of surface 
    278        !------------------------------------------------------------------------------ 
    279  
    280        !     4.1. Computation of intermediate values 
    281        !--------------------------------------------- 
    282        DO ji = kideb, kiut 
    283           zrcpdt(ji) = ( rcpsn * MIN( h_snow_1d(ji) , zhsu )    & 
    284              &       + rcpic * MAX( zhsu - h_snow_1d(ji) , zzero ) ) / rdt_ice 
    285           zts_old(ji) =  sist_1d(ji) 
    286        END DO 
    287  
    288        !     4.2. Computation of surface temperature by expanding the eq. of energy balance 
    289        !          with Ts = Tp + DT. One obtain , F(Tp) + DT * DF(Tp) = 0 
    290        !          where  - F(Tp) = Fsr + Fnsr(Tp) + Fcs(Tp)  
    291        !                 - DF(Tp)= (dFnsr(Tp)/dT) + (dFcs(Tp)/dT) - hsu*rcp/dt 
    292        !--------------------------------------------------------------------------------- 
    293  
    294        DO ji = kideb, kiut 
    295           !---computation of the derivative of energy balance function  
    296           zdfts    =  zksndh(ji)   & ! contribution of the conductive heat flux 
    297              &      + zrcpdt(ji)   & ! contribution of hsu * rcp / dt 
    298              &      - dqns_ice_1d (ji)     ! contribution of the total non solar radiation  
    299           !---computation of the energy balance function  
    300           zfts    = - z1mi0 (ji) * qsr_ice_1d(ji)   & ! net absorbed solar radiation 
    301              &      - qns_ice_1d(ji)                & ! total non solar radiation 
    302              &      - zfcsu (ji)                      ! conductive heat flux from the surface 
    303           !---computation of surface temperature increment   
    304           zdts    = -zfts / zdfts 
    305           !---computation of the new surface temperature  
    306           sist_1d(ji) = sist_1d(ji) + zdts 
    307        END DO 
    308  
    309        !---------------------------------------------------------------------------- 
    310        !  5. Boundary condition at the top surface 
    311        !--    IF Tsb < Tmelt, Fnet = Fcs (the net heat flux equal the conductive heat flux) 
    312        !      Otherwise Tsb = Tmelt and Qnet(Tmelt) > 0  
    313        !      Fnet(Tmelt) is therefore the net surface flux needed for melting 
    314        !---------------------------------------------------------------------------- 
     74      REAL(wp), DIMENSION(jpij,2) ::   zqcmlt   ! energy due to surface( /1 ) and bottom melting( /2 ) 
     75      REAL(wp), DIMENSION(jpij) ::   ztsmlt     ! snow/ice surface melting temperature 
     76      REAL(wp), DIMENSION(jpij) ::   ztbif      ! int. temp. at the mid-point of the 1st layer of the snow/ice sys.  
     77      REAL(wp), DIMENSION(jpij) ::   zksn       ! effective conductivity of snow 
     78      REAL(wp), DIMENSION(jpij) ::   zkic       ! effective conductivity of ice 
     79      REAL(wp), DIMENSION(jpij) ::   zksndh     ! thermal cond. at the mid-point of the 1st layer of the snow/ice sys.  
     80      REAL(wp), DIMENSION(jpij) ::   zfcsu      ! conductive heat flux at the surface of the snow/ice system  
     81      REAL(wp), DIMENSION(jpij) ::   zfcsudt    ! = zfcsu * dt 
     82      REAL(wp), DIMENSION(jpij) ::   zi0        ! frac. of the net SW rad. which is not absorbed at the surface 
     83      REAL(wp), DIMENSION(jpij) ::   z1mi0      ! fraction of the net SW radiation absorbed at the surface 
     84      REAL(wp), DIMENSION(jpij) ::   zqmax      ! maximum energy stored in brine pockets 
     85      REAL(wp), DIMENSION(jpij) ::   zrcpdt     ! h_su*rho_su*cp_su/dt(h_su being the thick. of surf. layer) 
     86      REAL(wp), DIMENSION(jpij) ::   zts_old    ! previous surface temperature 
     87      REAL(wp), DIMENSION(jpij) ::   zidsn , z1midsn , zidsnic ! tempory variables 
     88      !! 
     89      REAL(wp), DIMENSION(jpij) ::   zfnet       ! net heat flux at the top surface( incl. conductive heat flux) 
     90      REAL(wp), DIMENSION(jpij) ::   zsprecip    ! snow accumulation 
     91      REAL(wp), DIMENSION(jpij) ::   zhsnw_old   ! previous snow thickness 
     92      REAL(wp), DIMENSION(jpij) ::   zdhictop    ! change in ice thickness due to top surf ablation/accretion 
     93      REAL(wp), DIMENSION(jpij) ::   zdhicbot    ! change in ice thickness due to bottom surf abl/acc 
     94      REAL(wp), DIMENSION(jpij) ::   zqsup       ! energy transmitted to ocean (coming from top surf abl/acc) 
     95      REAL(wp), DIMENSION(jpij) ::   zqocea      ! energy transmitted to ocean (coming from bottom sur abl/acc) 
     96      REAL(wp), DIMENSION(jpij) ::   zfrl_old    ! previous sea/ice fraction 
     97      REAL(wp), DIMENSION(jpij) ::   zfrld_1d    ! new sea/ice fraction 
     98      REAL(wp), DIMENSION(jpij) ::   zep         ! internal temperature of the 2nd layer of the snow/ice system 
     99      !! 
     100      REAL(wp), DIMENSION(3) ::   zplediag       ! principle diagonal, subdiag. and supdiag. of the  
     101      REAL(wp), DIMENSION(3) ::   zsubdiag       ! tri-diagonal matrix coming from the computation 
     102      REAL(wp), DIMENSION(3) ::   zsupdiag       ! of the temperatures inside the snow-ice system 
     103      REAL(wp), DIMENSION(3) ::   zsmbr          ! second member 
     104      REAL(wp) ::   zhsu        ! thickness of surface layer 
     105      REAL(wp) ::   zhe         ! effective thickness for compu. of equ. thermal conductivity 
     106      REAL(wp) ::   zheshth     ! = zhe / thth 
     107      REAL(wp) ::   zghe        ! correction factor of the thermal conductivity 
     108      REAL(wp) ::   zumsb       ! parameter for numerical method to solve heat-diffusion eq. 
     109      REAL(wp) ::   zkhsn       ! conductivity at the snow layer 
     110      REAL(wp) ::   zkhic       ! conductivity at the ice layers 
     111      REAL(wp) ::   zkint       ! equivalent conductivity at the snow-ice interface 
     112      REAL(wp) ::   zkhsnint    ! = zkint*dt / (hsn*rhosn*cpsn)   
     113      REAL(wp) ::   zkhicint    ! = 2*zkint*dt / (hic*rhoic*cpic) 
     114      REAL(wp) ::   zpiv1 , zpiv2        ! tempory scalars used to solve the tri-diagonal system 
     115      REAL(wp) ::   zb2, zd2, zb3, zd3 
     116      REAL(wp) ::   ztint          ! equivalent temperature at the snow-ice interface 
     117      REAL(wp) ::   zexp         ! exponential function of the ice thickness 
     118      REAL(wp) ::   zfsab        ! part of solar radiation stored in brine pockets 
     119      REAL(wp) ::   zfts         ! value of energy balance function when the temp. equal surf. temp. 
     120      REAL(wp) ::   zdfts        ! value of derivative of ztfs when the temp. equal surf. temp. 
     121      REAL(wp) ::   zdts         ! surface temperature increment 
     122      REAL(wp) ::   zqsnw_mlt    ! energy needed to melt snow 
     123      REAL(wp) ::   zdhsmlt      ! change in snow thickness due to melt 
     124      REAL(wp) ::   zhsn         ! snow thickness (previous+accumulation-melt) 
     125      REAL(wp) ::   zqsn_mlt_rem   ! remaining heat coming from snow melting 
     126      REAL(wp) ::   zqice_top_mlt   ! energy used to melt ice at top surface 
     127      REAL(wp) ::   zdhssub        ! change in snow thick. due to sublimation or evaporation 
     128      REAL(wp) ::   zdhisub        ! change in ice thick. due to sublimation or evaporation     
     129      REAL(wp) ::   zdhsn          ! snow ice thickness increment 
     130      REAL(wp) ::   zdtsn          ! snow internal temp. increment 
     131      REAL(wp) ::   zdtic          ! ice internal temp. increment 
     132      REAL(wp) ::   zqnes          ! conductive energy due to ice melting in the first ice layer 
     133      !! 
     134      REAL(wp) ::   ztbot           ! temperature at the bottom surface 
     135      REAL(wp) ::   zfcbot          ! conductive heat flux at bottom surface 
     136      REAL(wp) ::   zqice_bot       ! energy used for bottom melting/growing 
     137      REAL(wp) ::   zqice_bot_mlt   ! energy used for bottom melting 
     138      REAL(wp) ::   zqstbif_bot    ! part of energy stored in brine pockets used for bottom melting 
     139      REAL(wp) ::   zqstbif_old    ! tempory var. for zqstbif_bot 
     140      REAL(wp) ::   zdhicmlt        ! change in ice thickness due to bottom melting 
     141      REAL(wp) ::   zdhicm          ! change in ice thickness var.  
     142      REAL(wp) ::   zdhsnm          ! change in snow thickness var.  
     143      REAL(wp) ::   zhsnfi          ! snow thickness var.  
     144      REAL(wp) ::   zc1, zpc1, zc2, zpc2, zp1, zp2   ! temporary scalars 
     145      REAL(wp) ::   ztb2, ztb3 
     146      !! 
     147      REAL(wp) ::   zdrmh           ! change in snow/ice thick. after snow-ice formation 
     148      REAL(wp) ::   zhicnew         ! new ice thickness 
     149      REAL(wp) ::   zhsnnew         ! new snow thickness 
     150      REAL(wp) ::   zquot , ztneq   ! tempory temp. variables 
     151      REAL(wp) ::   zqice, zqicetot & ! total heat inside the snow/ice system 
     152      REAL(wp) ::   zdfrl           ! change in ice concentration 
     153      REAL(wp) ::   zdvsnvol        ! change in snow volume 
     154      REAL(wp) ::   zdrfrl1, zdrfrl2   ! tempory scalars 
     155      REAL(wp) ::   zihsn, zidhb, zihic, zihe, zihq, ziexp, ziqf, zihnf, zibmlt, ziqr, zihgnew, zind 
     156      !!---------------------------------------------------------------------- 
     157 
     158      !----------------------------------------------------------------------- 
     159      !  1. Boundaries conditions for snow/ice system internal temperature 
     160      !       - If tbif_1d(ji,1) > rt0_snow, tbif_1d(ji,1) = rt0_snow  
     161      !       - If tbif_1d(ji,2/3) > rt0_ice, tbif_1d(ji,2/3) = rt0_ice  
     162      !     Computation of energies due to surface and bottom melting  
     163      !----------------------------------------------------------------------- 
     164      DO ji = kideb , kiut 
     165         zihsn = MAX( zzero , SIGN( zone , hsndif - h_snow_1d(ji) ) ) 
     166         zihic = MAX( zzero , SIGN( zone , hicdif - h_ice_1d(ji) ) ) 
     167         !--computation of energy due to surface melting 
     168         zqcmlt(ji,1) = ( MAX ( zzero ,  & 
     169            &                   rcpsn * h_snow_1d(ji) * ( tbif_1d(ji,1) - rt0_snow ) ) ) * ( 1.0 - zihsn ) 
     170         !--computation of energy due to bottom melting 
     171         zqcmlt(ji,2) = ( MAX( zzero , & 
     172            &                  rcpic * ( tbif_1d(ji,2) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 
     173            &           + MAX( zzero , & 
     174            &                  rcpic * ( tbif_1d(ji,3) - rt0_ice ) * ( h_ice_1d(ji) / 2. ) ) & 
     175            &           ) * ( 1.0 - zihic  ) 
     176         !--limitation of  snow/ice system internal temperature 
     177         tbif_1d(ji,1)   = MIN( rt0_snow, tbif_1d(ji,1) ) 
     178         tbif_1d(ji,2)   = MIN( rt0_ice , tbif_1d(ji,2) ) 
     179         tbif_1d(ji,3)   = MIN( rt0_ice , tbif_1d(ji,3) ) 
     180      END DO 
     181 
     182      !------------------------------------------- 
     183      !  2. Calculate some intermediate variables.   
     184      !------------------------------------------- 
     185      zhsu = hnzst         ! initialisation of the thickness of surface layer 
     186      ! 
     187      DO ji = kideb , kiut 
     188         zind   = MAX(  zzero , SIGN( zone , zhsu   - h_snow_1d(ji) )  ) 
     189         zihsn  = MAX(  zzero , SIGN( zone , hsndif - h_snow_1d(ji) )  ) 
     190         zihsn  = MAX(  zihsn , zind                                   ) 
     191         zihic  = MAX(  zzero , sign( zone , hicdif - h_ice_1d (ji) )  ) 
     192         ! 
     193         !     2.1. Computation of surface melting temperature 
     194         !---------------------------------------------------- 
     195         zind  = MAX( zzero , SIGN( zone , -h_snow_1d(ji) ) ) 
     196         ztsmlt(ji) = ( 1.0 - zind ) * rt0_snow + zind * rt0_ice 
     197         ! 
     198         !     2.2. Effective conductivity of snow and ice 
     199         !----------------------------------------------- 
     200 
     201         !---computation of the correction factor on the thermal conductivity 
     202         !-- (Morales Maqueda, 1995 ; Fichefet and Morales Maqueda, 1997) 
     203         zhe      =  ( rcdsn / ( rcdsn + rcdic ) ) * h_ice_1d(ji)   & 
     204            &     + ( rcdic / ( rcdsn + rcdic ) ) * h_snow_1d(ji)  
     205         zihe     = MAX( zzero , SIGN( zone , 2.0 * zhe - thth ) ) 
     206         zheshth  = zhe / thth 
     207         zghe     = ( 1.0 - zihe ) * zheshth * ( 2.0 - zheshth )   & 
     208            &     +         zihe   * 0.5 * ( 1.5 + LOG( 2.0 * zheshth ) ) 
     209 
     210         !---effective conductivities  
     211         zksn(ji)  = zghe * rcdsn   
     212         zkic(ji)  = zghe * rcdic 
     213         ! 
     214         !     2.3. Computation of the conductive heat flux from the snow/ice  
     215         !          system interior toward the top surface 
     216         !------------------------------------------------------------------ 
     217 
     218         !---Thermal conductivity at the mid-point of the first snow/ice system layer 
     219         zksndh(ji) =   ( ( 1.0 - zihsn ) * 2.0 * zksn(ji) + zihsn * 4.0 * zkic(ji) )   & 
     220            &         / ( ( 1.0 - zihsn ) *  h_snow_1d(ji)                              & 
     221            &           +        zihsn   *  ( ( 1.0 + 3.0 * zihic ) * h_ice_1d(ji)      & 
     222            &           + 4.0 * zkic(ji)/zksn(ji) * h_snow_1d(ji) ) ) 
     223 
     224         !---internal temperature at the mid-point of the first snow/ice system layer 
     225         ztbif(ji)  = ( 1.0 - zihsn ) * tbif_1d(ji,1)                       & 
     226            &       +         zihsn   * ( ( 1.0 - zihic ) * tbif_1d(ji,2)   & 
     227            &       +         zihic   * tfu_1d(ji)   ) 
     228         !---conductive heat flux  
     229         zfcsu(ji) = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 
     230 
     231      END DO 
     232 
     233      !-------------------------------------------------------------------- 
     234      !  3. Calculate :  
     235      !     - fstbif_1d, part of solar radiation absorbing inside the ice 
     236      !       assuming an exponential absorption (Grenfell and Maykut, 1977) 
     237      !     - zqmax,  maximum energy stored in brine pockets 
     238      !     - qstbif_1d, total energy stored in brine pockets (updating) 
     239      !------------------------------------------------------------------- 
     240 
     241      DO ji = kideb , kiut 
     242         zihsn  = MAX( zzero , SIGN (zone , -h_snow_1d(ji) ) ) 
     243         zihic  = MAX( zzero , 1.0 - ( h_ice_1d(ji) / zhsu ) )      
     244         zind   = MAX( zzero , SIGN (zone , hicdif - h_ice_1d(ji) ) ) 
     245         !--Computation of the fraction of the net shortwave radiation which 
     246         !--penetrates inside the ice cover ( See Forcat) 
     247         zi0(ji)  = zihsn * ( fr1_i0_1d(ji) + zihic * fr2_i0_1d(ji) ) 
     248         zexp     = MIN( zone , EXP( -1.5 * ( h_ice_1d(ji) - zhsu ) ) ) 
     249         fstbif_1d(ji) = zi0(ji) * qsr_ice_1d(ji) * zexp 
     250         !--Computation of maximum energy stored in brine pockets zqmax and update 
     251         !--the total energy stored in brine pockets, if less than zqmax 
     252         zqmax(ji) = MAX( zzero , 0.5 * xlic * ( h_ice_1d(ji) - hicmin ) ) 
     253         zfsab   = zi0(ji) * qsr_ice_1d(ji) * ( 1.0 - zexp ) 
     254         zihq    = ( 1.0 - zind ) * MAX(zzero, SIGN( zone , qstbif_1d(ji) - zqmax(ji) ) ) & 
     255            &    +         zind   * zone 
     256         qstbif_1d(ji) = ( qstbif_1d(ji) + ( 1.0 - zihq ) * zfsab * rdt_ice ) * swiqst 
     257         !--fraction of shortwave radiation absorbed at surface 
     258         ziexp = zihq * zexp + ( 1.0 - zihq ) * ( swiqst + ( 1.0 - swiqst ) * zexp ) 
     259         z1mi0(ji) = 1.0 - zi0(ji) * ziexp 
     260      END DO 
     261 
     262      !-------------------------------------------------------------------------------- 
     263      !  4. Computation of the surface temperature : determined by considering the  
     264      !     budget of a thin layer of thick. zhsu at the top surface (H. Grenier, 1995)  
     265      !     and based on a surface energy balance :  
     266      !     hsu * rcp * dT/dt = Fsr + Fnsr(T) + Fcs(T), 
     267      !     where - Fsr is the net absorbed solar radiation,  
     268      !           - Fnsr is the total non solar radiation (incoming and outgoing long-wave, 
     269      !             sensible and latent heat fluxes) 
     270      !           - Fcs the conductive heat flux at the top of surface 
     271      !------------------------------------------------------------------------------ 
     272 
     273      !     4.1. Computation of intermediate values 
     274      !--------------------------------------------- 
     275      DO ji = kideb, kiut 
     276         zrcpdt(ji) = ( rcpsn * MIN( h_snow_1d(ji) , zhsu )    & 
     277            &       + rcpic * MAX( zhsu - h_snow_1d(ji) , zzero ) ) / rdt_ice 
     278         zts_old(ji) =  sist_1d(ji) 
     279      END DO 
     280 
     281      !     4.2. Computation of surface temperature by expanding the eq. of energy balance 
     282      !          with Ts = Tp + DT. One obtain , F(Tp) + DT * DF(Tp) = 0 
     283      !          where  - F(Tp) = Fsr + Fnsr(Tp) + Fcs(Tp)  
     284      !                 - DF(Tp)= (dFnsr(Tp)/dT) + (dFcs(Tp)/dT) - hsu*rcp/dt 
     285      !--------------------------------------------------------------------------------- 
     286 
     287      DO ji = kideb, kiut 
     288         !---computation of the derivative of energy balance function  
     289         zdfts    =  zksndh(ji)   & ! contribution of the conductive heat flux 
     290            &      + zrcpdt(ji)   & ! contribution of hsu * rcp / dt 
     291            &      - dqns_ice_1d (ji)     ! contribution of the total non solar radiation  
     292         !---computation of the energy balance function  
     293         zfts    = - z1mi0 (ji) * qsr_ice_1d(ji)   & ! net absorbed solar radiation 
     294            &      - qns_ice_1d(ji)                & ! total non solar radiation 
     295            &      - zfcsu (ji)                      ! conductive heat flux from the surface 
     296         !---computation of surface temperature increment   
     297         zdts    = -zfts / zdfts 
     298         !---computation of the new surface temperature  
     299         sist_1d(ji) = sist_1d(ji) + zdts 
     300      END DO 
     301 
     302      !---------------------------------------------------------------------------- 
     303      !  5. Boundary condition at the top surface 
     304      !--    IF Tsb < Tmelt, Fnet = Fcs (the net heat flux equal the conductive heat flux) 
     305      !      Otherwise Tsb = Tmelt and Qnet(Tmelt) > 0  
     306      !      Fnet(Tmelt) is therefore the net surface flux needed for melting 
     307      !---------------------------------------------------------------------------- 
    315308        
    316309        
    317        !     5.1.  Limitation of surface temperature and update total non solar fluxes, 
    318        !          latent heat flux and conductive flux at the top surface  
    319        !----------------------------------------------------------------------   
     310      !     5.1.  Limitation of surface temperature and update total non solar fluxes, 
     311      !          latent heat flux and conductive flux at the top surface  
     312      !----------------------------------------------------------------------   
    320313                      
    321        IF ( .NOT. lk_cpl ) THEN   ! duplicate the loop for performances issues 
    322           DO ji = kideb, kiut 
    323              sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 
    324              qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
    325              qla_ice_1d(ji) = qla_ice_1d(ji) + dqla_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
    326              zfcsu(ji)  = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 
    327           END DO 
    328        ELSE 
    329           DO ji = kideb, kiut 
    330              sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 
    331              zfcsu(ji)  = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 
    332           END DO 
    333        ENDIF 
    334  
    335        !     5.2. Calculate available heat for surface ablation.  
    336        !--------------------------------------------------------------------- 
    337  
    338        DO ji = kideb, kiut 
    339           zfnet(ji) = qns_ice_1d(ji) + z1mi0(ji) * qsr_ice_1d(ji) + zfcsu(ji)           
    340           zfnet(ji) = MAX( zzero , zfnet(ji) ) 
    341           zfnet(ji) = zfnet(ji) * MAX( zzero , SIGN( zone , sist_1d(ji) - ztsmlt(ji) ) ) 
    342        END DO 
    343  
    344        !------------------------------------------------------------------------- 
    345        !  6. Calculate changes in internal temperature due to vertical diffusion    
    346        !     processes. The evolution of this temperature is governed by the one- 
    347        !     dimensionnal heat-diffusion equation.  
    348        !     Given the temperature tbif(1/2/3), at time m we solve a set 
    349        !     of finite difference equations to obtain new tempe. Each tempe is coupled 
    350        !     to the temp. immediatly above and below by heat conduction terms. Thus  
    351        !     we have a set of equations of the form A * T = B, where A is a tridiagonal 
    352        !     matrix, T a vector whose components are the unknown new temp. 
    353        !------------------------------------------------------------------------- 
     314      IF ( .NOT. lk_cpl ) THEN   ! duplicate the loop for performances issues 
     315         DO ji = kideb, kiut 
     316            sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 
     317            qns_ice_1d(ji) = qns_ice_1d(ji) + dqns_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
     318            qla_ice_1d(ji) = qla_ice_1d(ji) + dqla_ice_1d(ji) * ( sist_1d(ji) - zts_old(ji) ) 
     319            zfcsu(ji)  = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 
     320         END DO 
     321      ELSE 
     322         DO ji = kideb, kiut 
     323            sist_1d(ji) = MIN( ztsmlt(ji) , sist_1d(ji) ) 
     324            zfcsu(ji)  = zksndh(ji) * ( ztbif(ji) - sist_1d(ji) ) 
     325         END DO 
     326      ENDIF 
     327 
     328      !     5.2. Calculate available heat for surface ablation.  
     329      !--------------------------------------------------------------------- 
     330 
     331      DO ji = kideb, kiut 
     332         zfnet(ji) = qns_ice_1d(ji) + z1mi0(ji) * qsr_ice_1d(ji) + zfcsu(ji)           
     333         zfnet(ji) = MAX( zzero , zfnet(ji) ) 
     334         zfnet(ji) = zfnet(ji) * MAX( zzero , SIGN( zone , sist_1d(ji) - ztsmlt(ji) ) ) 
     335      END DO 
     336 
     337      !------------------------------------------------------------------------- 
     338      !  6. Calculate changes in internal temperature due to vertical diffusion    
     339      !     processes. The evolution of this temperature is governed by the one- 
     340      !     dimensionnal heat-diffusion equation.  
     341      !     Given the temperature tbif(1/2/3), at time m we solve a set 
     342      !     of finite difference equations to obtain new tempe. Each tempe is coupled 
     343      !     to the temp. immediatly above and below by heat conduction terms. Thus  
     344      !     we have a set of equations of the form A * T = B, where A is a tridiagonal 
     345      !     matrix, T a vector whose components are the unknown new temp. 
     346      !------------------------------------------------------------------------- 
    354347        
    355        !--parameter for the numerical methode use to solve the heat-diffusion equation 
    356        !- implicit, explicit or Crank-Nicholson 
    357        zumsb = 1.0 - sbeta   
    358        DO ji = kideb, kiut 
    359           zidsn(ji)   = MAX ( zzero, SIGN( zone, hsndif - h_snow_1d(ji) ) )  
    360           z1midsn(ji) = 1.0 - zidsn(ji) 
    361           zihic       = MAX ( zzero, SIGN( zone, hicdif - h_ice_1d(ji) ) )  
    362           zidsnic(ji) = zidsn(ji) *  zihic  
    363           zfcsudt(ji) = zfcsu(ji) * rdt_ice  
    364        END DO  
     348      !--parameter for the numerical methode use to solve the heat-diffusion equation 
     349      !- implicit, explicit or Crank-Nicholson 
     350      zumsb = 1.0 - sbeta   
     351      DO ji = kideb, kiut 
     352         zidsn(ji)   = MAX ( zzero, SIGN( zone, hsndif - h_snow_1d(ji) ) )  
     353         z1midsn(ji) = 1.0 - zidsn(ji) 
     354         zihic       = MAX ( zzero, SIGN( zone, hicdif - h_ice_1d(ji) ) )  
     355         zidsnic(ji) = zidsn(ji) *  zihic  
     356         zfcsudt(ji) = zfcsu(ji) * rdt_ice  
     357      END DO  
    365358    
    366        DO ji = kideb, kiut 
    367  
    368           !     6.1 Calculate intermediate variables. 
    369           !---------------------------------------- 
    370  
    371           !--conductivity at the snow surface 
    372           zkhsn = 2.0 * zksn(ji) * rdt_ice / rcpsn 
    373           !--conductivity at the ice surface 
    374           zkhic = 4.0 * zkic(ji) * rdt_ice / MAX( h_ice_1d(ji) * h_ice_1d(ji) * rcpic , epsi20 ) 
    375           !--conductivity at the snow/ice interface  
    376           zkint = 4.0 * zksn(ji) * zkic(ji)  & 
    377              &        / ( zksn(ji) * h_ice_1d(ji) + 2.0 * zkic(ji) * h_snow_1d(ji) * z1midsn(ji))  
    378           zkhsnint = zkint * rdt_ice / rcpsn 
    379           zkhicint = zkint * 2.0 * rdt_ice / MAX( h_ice_1d(ji) * rcpic , epsi20 ) 
    380            
    381           !     6.2. Fulfill the linear system matrix. 
    382           !----------------------------------------- 
     359      DO ji = kideb, kiut 
     360 
     361         !     6.1 Calculate intermediate variables. 
     362         !---------------------------------------- 
     363 
     364         !--conductivity at the snow surface 
     365         zkhsn = 2.0 * zksn(ji) * rdt_ice / rcpsn 
     366         !--conductivity at the ice surface 
     367         zkhic = 4.0 * zkic(ji) * rdt_ice / MAX( h_ice_1d(ji) * h_ice_1d(ji) * rcpic , epsi20 ) 
     368         !--conductivity at the snow/ice interface  
     369         zkint = 4.0 * zksn(ji) * zkic(ji)  & 
     370            &        / ( zksn(ji) * h_ice_1d(ji) + 2.0 * zkic(ji) * h_snow_1d(ji) * z1midsn(ji))  
     371         zkhsnint = zkint * rdt_ice / rcpsn 
     372         zkhicint = zkint * 2.0 * rdt_ice / MAX( h_ice_1d(ji) * rcpic , epsi20 ) 
     373           
     374         !     6.2. Fulfill the linear system matrix. 
     375         !----------------------------------------- 
    383376!$$$          zplediag(1) = 1 + sbeta * z1midsn(ji) * ( zkhsn + zkhsnint )        
    384           zplediag(1) =   zidsn(ji) + z1midsn(ji) * h_snow_1d(ji)   & 
    385              &          + sbeta * z1midsn(ji) * zkhsnint  
    386           zplediag(2) = 1 + sbeta * ( z1midsn(ji) * zkhicint + zkhic )  
    387           zplediag(3) = 1 + 3.0 * sbeta * zkhic    
    388  
    389           zsubdiag(1) =  0.e0               
    390           zsubdiag(2) = -1.e0 * z1midsn(ji) * sbeta * zkhicint 
    391           zsubdiag(3) = -1.e0 * sbeta * zkhic  
    392  
    393           zsupdiag(1) = -1.e0 * z1midsn(ji) * sbeta * zkhsnint  
    394           zsupdiag(2) = zsubdiag(3) 
    395           zsupdiag(3) =  0.e0 
    396            
    397           !     6.3. Fulfill the idependent term vector. 
    398           !------------------------------------------- 
     377         zplediag(1) =   zidsn(ji) + z1midsn(ji) * h_snow_1d(ji)   & 
     378            &          + sbeta * z1midsn(ji) * zkhsnint  
     379         zplediag(2) = 1 + sbeta * ( z1midsn(ji) * zkhicint + zkhic )  
     380         zplediag(3) = 1 + 3.0 * sbeta * zkhic    
     381 
     382         zsubdiag(1) =  0.e0               
     383         zsubdiag(2) = -1.e0 * z1midsn(ji) * sbeta * zkhicint 
     384         zsubdiag(3) = -1.e0 * sbeta * zkhic  
     385 
     386         zsupdiag(1) = -1.e0 * z1midsn(ji) * sbeta * zkhsnint  
     387         zsupdiag(2) = zsubdiag(3) 
     388         zsupdiag(3) =  0.e0 
     389           
     390         !     6.3. Fulfill the idependent term vector. 
     391         !------------------------------------------- 
    399392           
    400393!$$$          zsmbr(1) = zidsn(ji) * sist_1d(ji) + z1midsn(ji) *   & 
    … …  
    402395!$$$             &         - zumsb * ( zkhsn * tbif_1d(ji,1) 
    403396!$$$             &                   + zkhsnint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) ) ) 
    404           zsmbr(1) = zidsn(ji) * sist_1d(ji) + z1midsn(ji) *    & 
    405              &       ( h_snow_1d(ji) * tbif_1d(ji,1) - ( zfcsudt(ji) / rcpsn )  & 
    406              &       - zumsb * zkhsnint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) ) 
    407  
    408           zsmbr(2) =  tbif_1d(ji,2)  & 
    409              &      - zidsn(ji) * ( 1.0 - zidsnic(ji) ) & 
    410              &        * ( zfcsudt(ji) / MAX( h_ice_1d(ji) * rcpic , epsi20 ) ) & 
    411              &      + zumsb * ( zkhicint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) & 
    412              &                   - zkhic * ( tbif_1d(ji,2) - tbif_1d(ji,3) )  ) 
    413  
    414           zsmbr(3) =  tbif_1d(ji,3)  & 
    415              &      + zkhic * ( 2.0 * tfu_1d(ji) & 
    416              &                + zumsb * ( tbif_1d(ji,2) - 3.0 * tbif_1d(ji,3) ) ) 
    417            
    418           !     6.4. Solve the system (Gauss elimination method). 
    419           !---------------------------------------------------- 
    420            
    421           zpiv1 = zsubdiag(2) / zplediag(1)  
    422           zb2   = zplediag(2) - zpiv1 * zsupdiag(1) 
    423           zd2   = zsmbr(2) - zpiv1 * zsmbr(1) 
    424  
    425           zpiv2 = zsubdiag(3) / zb2 
    426           zb3   = zplediag(3) - zpiv2 * zsupdiag(2) 
    427           zd3   = zsmbr(3) - zpiv2 * zd2 
    428  
    429           tbif_1d(ji,3) = zd3 / zb3 
    430           tbif_1d(ji,2) = ( zd2 - zsupdiag(2) * tbif_1d(ji,3) ) / zb2 
    431           tbif_1d(ji,1) = ( zsmbr(1) - zsupdiag(1) * tbif_1d(ji,2) ) / zplediag(1)             
    432  
    433           !--- taking into account the particular case of  zidsnic(ji) = 1 
    434           ztint =  (  zkic(ji) * h_snow_1d(ji) * tfu_1d (ji)    & 
    435              &      + zksn(ji) * h_ice_1d(ji) * sist_1d(ji) )   & 
    436              &   / ( zkic(ji) * h_snow_1d(ji) + zksn(ji) * h_ice_1d(ji) )  
    437  
    438           tbif_1d(ji,1) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,1)   & 
    439              &                + zidsnic(ji)   * ( ztint + sist_1d(ji) ) / 2.0 
    440           tbif_1d(ji,2) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,2)   & 
    441              &                + zidsnic(ji)   * ( 3.0 * ztint + tfu_1d(ji) ) / 4.0 
    442           tbif_1d(ji,3) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,3)   & 
    443              &                + zidsnic(ji)   * ( ztint + 3.0 * tfu_1d(ji) ) / 4.0      
    444        END DO 
     397         zsmbr(1) = zidsn(ji) * sist_1d(ji) + z1midsn(ji) *    & 
     398            &       ( h_snow_1d(ji) * tbif_1d(ji,1) - ( zfcsudt(ji) / rcpsn )  & 
     399            &       - zumsb * zkhsnint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) ) 
     400 
     401         zsmbr(2) =  tbif_1d(ji,2)  & 
     402            &      - zidsn(ji) * ( 1.0 - zidsnic(ji) ) & 
     403            &        * ( zfcsudt(ji) / MAX( h_ice_1d(ji) * rcpic , epsi20 ) ) & 
     404            &      + zumsb * ( zkhicint * ( tbif_1d(ji,1) - tbif_1d(ji,2) ) & 
     405            &                   - zkhic * ( tbif_1d(ji,2) - tbif_1d(ji,3) )  ) 
     406 
     407         zsmbr(3) =  tbif_1d(ji,3)  & 
     408            &      + zkhic * ( 2.0 * tfu_1d(ji) & 
     409            &                + zumsb * ( tbif_1d(ji,2) - 3.0 * tbif_1d(ji,3) ) ) 
     410           
     411         !     6.4. Solve the system (Gauss elimination method). 
     412         !---------------------------------------------------- 
     413           
     414         zpiv1 = zsubdiag(2) / zplediag(1)  
     415         zb2   = zplediag(2) - zpiv1 * zsupdiag(1) 
     416         zd2   = zsmbr(2) - zpiv1 * zsmbr(1) 
     417 
     418         zpiv2 = zsubdiag(3) / zb2 
     419         zb3   = zplediag(3) - zpiv2 * zsupdiag(2) 
     420         zd3   = zsmbr(3) - zpiv2 * zd2 
     421 
     422         tbif_1d(ji,3) = zd3 / zb3 
     423         tbif_1d(ji,2) = ( zd2 - zsupdiag(2) * tbif_1d(ji,3) ) / zb2 
     424         tbif_1d(ji,1) = ( zsmbr(1) - zsupdiag(1) * tbif_1d(ji,2) ) / zplediag(1)             
     425 
     426         !--- taking into account the particular case of  zidsnic(ji) = 1 
     427         ztint =  (  zkic(ji) * h_snow_1d(ji) * tfu_1d (ji)    & 
     428            &      + zksn(ji) * h_ice_1d(ji) * sist_1d(ji) )   & 
     429            &   / ( zkic(ji) * h_snow_1d(ji) + zksn(ji) * h_ice_1d(ji) )  
     430 
     431         tbif_1d(ji,1) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,1)   & 
     432            &                + zidsnic(ji)   * ( ztint + sist_1d(ji) ) / 2.0 
     433         tbif_1d(ji,2) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,2)   & 
     434            &                + zidsnic(ji)   * ( 3.0 * ztint + tfu_1d(ji) ) / 4.0 
     435         tbif_1d(ji,3) = ( 1.0 - zidsnic(ji) ) * tbif_1d(ji,3)   & 
     436            &                + zidsnic(ji)   * ( ztint + 3.0 * tfu_1d(ji) ) / 4.0      
     437      END DO 
    445438  
    446        !---------------------------------------------------------------------- 
    447        !  9. Take into account surface ablation and bottom accretion-ablation.| 
    448        !---------------------------------------------------------------------- 
     439      !---------------------------------------------------------------------- 
     440      !  9. Take into account surface ablation and bottom accretion-ablation.| 
     441      !---------------------------------------------------------------------- 
    449442        
    450        !---Snow accumulation in one thermodynamic time step 
    451        zsprecip(kideb:kiut) = sprecip_1d(kideb:kiut) * rdt_ice / rhosn 
    452  
    453  
    454        DO ji = kideb, kiut 
    455            
    456           !      9.1. Surface ablation and update of snow thickness and qstbif_1d 
    457           !-------------------------------------------------------------------- 
     443      !---Snow accumulation in one thermodynamic time step 
     444      zsprecip(kideb:kiut) = sprecip_1d(kideb:kiut) * rdt_ice / rhosn 
     445 
     446 
     447      DO ji = kideb, kiut 
     448           
     449         !      9.1. Surface ablation and update of snow thickness and qstbif_1d 
     450         !-------------------------------------------------------------------- 
    458451           
    459452          !-------------------------------------------------------------------------- 
    … …  
    759752          qstbif_1d(ji) = zdrfrl2 * qstbif_1d(ji) 
    760753          frld_1d(ji)    = zfrld_1d(ji) 
    761           ! 
    762        END DO 
    763        !  
    764     END SUBROUTINE lim_thd_zdf_2 
     754         ! 
     755      END DO 
     756      !  
     757   END SUBROUTINE lim_thd_zdf_2 
    765758 
    766759#else 
    … …  
    768761   !!   Default Option                                     NO sea-ice model 
    769762   !!---------------------------------------------------------------------- 
    770 CONTAINS 
    771    SUBROUTINE lim_thd_zdf_2          ! Empty routine 
    772    END SUBROUTINE lim_thd_zdf_2 
    773763#endif 
    774764 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limtrp_2.F90

    r1715 r1855  
    44   !! LIM 2.0 transport ice model : sea-ice advection/diffusion 
    55   !!====================================================================== 
     6   !! History :  LIM  ! 2000-01 (LIM)  Original code 
     7   !!            1.0  ! 2001-05 (G. Madec, R. Hordoir) opa norm 
     8   !!            2.0  ! 2004-01 (G. Madec, C. Ethe)  F90, mpp 
     9   !!---------------------------------------------------------------------- 
    610#if defined key_lim2 
    711   !!---------------------------------------------------------------------- 
    … …  
    1115   !!   lim_trp_init_2 : initialization and namelist read 
    1216   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    1417   USE phycst 
    1518   USE dom_oce 
    … …  
    2629   PRIVATE 
    2730 
    28    !! * Routine accessibility 
    29    PUBLIC lim_trp_2     ! called by sbc_ice_lim_2 
    30  
    31    !! * Shared module variables 
    32    REAL(wp), PUBLIC  ::   &  !: 
    33       bound  = 0.e0          !: boundary condit. (0.0 no-slip, 1.0 free-slip) 
    34  
    35    !! * Module variables 
    36    REAL(wp)  ::           &  ! constant values 
    37       epsi06 = 1.e-06  ,  & 
    38       epsi03 = 1.e-03  ,  & 
    39       epsi16 = 1.e-16  ,  & 
    40       rzero  = 0.e0    ,  & 
    41       rone   = 1.e0 
     31   PUBLIC   lim_trp_2   ! called by sbc_ice_lim_2 
     32 
     33   REAL(wp), PUBLIC  ::   bound  = 0.e0   !: boundary condit. (0.0 no-slip, 1.0 free-slip) 
     34 
     35   REAL(wp) ::   epsi06 = 1.e-06   ! constant values 
     36   REAL(wp) ::   epsi03 = 1.e-03   ! 
     37   REAL(wp) ::   epsi16 = 1.e-16   ! 
     38   REAL(wp) ::   rzero  = 0.e0     ! 
     39   REAL(wp) ::   rone   = 1.e0     ! 
    4240 
    4341   !! * Substitution 
    4442#  include "vectopt_loop_substitute.h90" 
    4543   !!---------------------------------------------------------------------- 
    46    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     44   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    4745   !! $Id$ 
    48    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     46   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    4947   !!---------------------------------------------------------------------- 
    5048 
    … …  
    6260      !! 
    6361      !! ** action : 
    64       !! 
    65       !! History : 
    66       !!   1.0  !  00-01 (LIM)  Original code 
    67       !!        !  01-05 (G. Madec, R. Hordoir) opa norm 
    68       !!   2.0  !  04-01 (G. Madec, C. Ethe)  F90, mpp 
    6962      !!--------------------------------------------------------------------- 
    7063      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    71  
    72       INTEGER  ::   ji, jj, jk,   &  ! dummy loop indices 
    73          &          initad           ! number of sub-timestep for the advection 
    74  
    75       REAL(wp) ::  &                               
    76          zindb  ,  & 
    77          zacrith, & 
    78          zindsn , & 
    79          zindic , & 
    80          zusvosn, & 
    81          zusvoic, & 
    82          zignm  , & 
    83          zindhe , & 
    84          zvbord , & 
    85          zcfl   , & 
    86          zusnit , & 
    87          zrtt, ztsn, ztic1, ztic2 
    88  
    89       REAL(wp), DIMENSION(jpi,jpj)  ::   &  ! temporary workspace 
    90          zui_u , zvi_v , zsm   ,         & 
    91          zs0ice, zs0sn , zs0a  ,         & 
    92          zs0c0 , zs0c1 , zs0c2 ,         & 
    93          zs0st 
     64      !! 
     65      INTEGER  ::   ji, jj, jk   ! dummy loop indices 
     66      INTEGER  ::   initad       ! number of sub-timestep for the advection 
     67      REAL(wp) ::   zindb , zacrith, zindsn , zignm , zvbord, zrtt, ztic1, zusnit 
     68      REAL(wp) ::   zindic, zusvosn, zusvoic, zindhe, zcfl  , ztsn, ztic2 
     69      REAL(wp), DIMENSION(jpi,jpj) ::   zui_u , zs0sn, zs0c0, zs0a , zsm      ! 2D workspace 
     70      REAL(wp), DIMENSION(jpi,jpj) ::   zvi_v , zs0st, zs0c1, zs0c2, zs0ice   !  -      - 
    9471      !--------------------------------------------------------------------- 
    9572 
    … …  
    9875      zsm(:,:) = area(:,:) 
    9976       
    100       IF( ln_limdyn ) THEN 
    101          !-------------------------------------! 
    102          !   Advection of sea ice properties   ! 
    103          !-------------------------------------! 
     77      !                             !-------------------------------------! 
     78      IF( ln_limdyn ) THEN          !   Advection of sea ice properties   ! 
     79         !                          !-------------------------------------! 
    10480 
    10581         ! ice velocities at ocean U- and V-points (zui_u,zvi_v) 
    … …  
    11389            END DO 
    11490         END DO 
    115          ! Lateral boundary conditions on zui_u, zvi_v 
    116          CALL lbc_lnk( zui_u, 'U', -1. ) 
    117          CALL lbc_lnk( zvi_v, 'V', -1. ) 
     91         CALL lbc_lnk( zui_u, 'U', -1. )   ;   CALL lbc_lnk( zvi_v, 'V', -1. )         ! Lateral boundary conditions 
    11892 
    11993         ! CFL test for stability 
    … …  
    12397         zcfl  = MAX( zcfl, MAXVAL( ABS( zvi_v( :     ,1:jpjm1) ) * rdt_ice / e2v( :     ,1:jpjm1) ) ) 
    12498 
    125          IF (lk_mpp ) CALL mpp_max(zcfl) 
    126  
    127          IF ( zcfl > 0.5 .AND. lwp )   WRITE(numout,*) 'lim_trp_2 : violation of cfl criterion the ',nday,'th day, cfl = ',zcfl 
     99         IF( lk_mpp ) CALL mpp_max( zcfl ) 
     100 
     101         IF( zcfl > 0.5 .AND. lwp )   WRITE(numout,*) 'lim_trp_2 : CFL violation at the ',nday,'th day, cfl = ',zcfl 
    128102 
    129103         ! content of properties 
    … …  
    144118         zusnit = 1.0 / REAL( initad )  
    145119          
    146          IF ( MOD( nday , 2 ) == 0) THEN 
     120         IF( MOD( nday , 2 ) == 0) THEN 
    147121            DO jk = 1,initad 
    148122               CALL lim_adv_x_2( zusnit, zui_u, rone , zsm, zs0ice, sxice, sxxice, syice, syyice, sxyice ) 
    … …  
    228202         !   Up-dating and limitation of sea ice properties after transport   ! 
    229203         ! -------------------------------------------------------------------! 
    230  
    231          ! Up-dating and limitation of sea ice properties after transport. 
    232204         DO jj = 1, jpj 
    233 !!!iii      zindhe = REAL( MAX( 0, isign(1, jj - njeq ) ) )              !ibug mpp  !!bugmpp  njeq! 
    234205            zindhe = MAX( 0.e0, SIGN( 1.e0, fcor(1,jj) ) )              ! = 0 for SH, =1 for NH 
    235206            DO ji = 1, jpi 
    236  
    237                ! Recover mean values over the grid squares. 
     207               !                                                        ! Recover mean values over the grid squares. 
    238208               zs0sn (ji,jj) = MAX( rzero, zs0sn (ji,jj)/area(ji,jj) ) 
    239209               zs0ice(ji,jj) = MAX( rzero, zs0ice(ji,jj)/area(ji,jj) ) 
    … …  
    243213               zs0c2 (ji,jj) = MAX( rzero, zs0c2 (ji,jj)/area(ji,jj) ) 
    244214               zs0st (ji,jj) = MAX( rzero, zs0st (ji,jj)/area(ji,jj) ) 
    245  
    246                ! Recover in situ values. 
     215               !                                                        ! Recover in situ values. 
    247216               zindb         = MAX( rzero, SIGN( rone, zs0a(ji,jj) - epsi06 ) ) 
    248217               zacrith       = 1.0 - ( zindhe * acrit(1) + ( 1.0 - zindhe ) * acrit(2) ) 
    … …  
    262231               zrtt          = 173.15 * rone  
    263232               ztsn          =          zignm   * tbif(ji,jj,1)  & 
    264                               + ( 1.0 - zignm ) * MIN( MAX( zrtt, rt0_snow * zusvosn * zs0c0(ji,jj)) , tfu(ji,jj) )  
     233                  &           + ( 1.0 - zignm ) * MIN( MAX( zrtt, rt0_snow * zusvosn * zs0c0(ji,jj)) , tfu(ji,jj) )  
    265234               ztic1          = MIN( MAX( zrtt, rt0_ice * zusvoic * zs0c1(ji,jj) ) , tfu(ji,jj) ) 
    266235               ztic2          = MIN( MAX( zrtt, rt0_ice * zusvoic * zs0c2(ji,jj) ) , tfu(ji,jj) ) 
    267   
     236               ! 
    268237               tbif(ji,jj,1) = zindsn * ztsn  + ( 1.0 - zindsn ) * tfu(ji,jj)                
    269238               tbif(ji,jj,2) = zindic * ztic1 + ( 1.0 - zindic ) * tfu(ji,jj) 
    … …  
    272241            END DO 
    273242         END DO 
    274           
     243         ! 
    275244      ENDIF 
    276        
     245      ! 
    277246   END SUBROUTINE lim_trp_2 
    278247 
    … …  
    288257      !! 
    289258      !! ** input   :   Namelist namicetrp 
    290       !! 
    291       !! history : 
    292       !!   2.0  !  03-08 (C. Ethe)  Original code 
    293259      !!------------------------------------------------------------------- 
    294260      NAMELIST/namicetrp/ bound 
    295261      !!------------------------------------------------------------------- 
    296  
    297       ! Read Namelist namicetrp 
    298       REWIND ( numnam_ice ) 
     262      ! 
     263      REWIND ( numnam_ice )                  ! Read Namelist namicetrp 
    299264      READ   ( numnam_ice  , namicetrp ) 
    300       IF(lwp) THEN 
     265      IF(lwp) THEN                           ! control print 
    301266         WRITE(numout,*) 
    302267         WRITE(numout,*) 'lim_trp_init_2 : Ice parameters for advection ' 
    … …  
    304269         WRITE(numout,*) '   boundary conditions (0. no-slip, 1. free-slip) bound  = ', bound 
    305270      ENDIF 
    306              
     271      !  
    307272   END SUBROUTINE lim_trp_init_2 
    308273 
    … …  
    311276   !!   Default option         Empty Module                No sea-ice model 
    312277   !!---------------------------------------------------------------------- 
    313 CONTAINS 
    314    SUBROUTINE lim_trp_2        ! Empty routine 
    315    END SUBROUTINE lim_trp_2 
    316278#endif 
    317279 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_2/limwri_2.F90

    r1818 r1855  
    1212   !!   'key_lim2'                                    LIM 2.0 sea-ice model 
    1313   !!---------------------------------------------------------------------- 
    14    !!---------------------------------------------------------------------- 
    15    !!   lim_wri_2      : write of the diagnostics variables in ouput file  
    16    !!   lim_wri_init_2 : initialization and namelist read 
    17    !!   lim_wri_state_2 : write for initial state or/and abandon: 
    18    !!                     > output.init.nc (if ninist = 1 in namelist) 
    19    !!                     > output.abort.nc 
     14   !!   lim_wri_2       : write of the diagnostics variables in ouput file  
     15   !!   lim_wri_init_2  : initialization and namelist read 
     16   !!   lim_wri_state_2 : write for initial state (output.init.nc if ninist=1) or/and in the abort file 
    2017   !!---------------------------------------------------------------------- 
    2118   USE phycst 
    … …  
    4239   INTEGER, PARAMETER                       ::   jpnoumax = 40   ! maximum number of variable for ice output 
    4340   INTEGER                                  ::   noumef          ! number of fields 
    44    REAL(wp)           , DIMENSION(jpnoumax) ::   cmulti ,     &  ! multiplicative constant 
    45       &                                          cadd            ! additive constant 
     41   REAL(wp)           , DIMENSION(jpnoumax) ::   cmulti         ! multiplicative constant 
     42   REAL(wp)           , DIMENSION(jpnoumax) ::   cadd            ! additive constant 
    4643   CHARACTER(len = 35), DIMENSION(jpnoumax) ::   titn            ! title of the field 
    4744   CHARACTER(len = 8 ), DIMENSION(jpnoumax) ::   nam             ! name of the field 
    … …  
    5249   INTEGER , DIMENSION( jpij ) ::   ndex51              ! ???? 
    5350 
    54    REAL(wp)  ::            &  ! constant values 
    55       epsi16 = 1.e-16   ,  & 
    56       zzero  = 0.e0     ,  & 
    57       zone   = 1.e0 
     51   REAL(wp) ::   epsi16 = 1.e-16   ! constant values 
     52   REAL(wp) ::   rzero  = 0.e0     ! 
     53   REAL(wp) ::   rone   = 1.e0     ! 
    5854 
    5955   !! * Substitutions 
    6056#   include "vectopt_loop_substitute.h90" 
    6157   !!---------------------------------------------------------------------- 
    62    !!  LIM 2.0, UCL-LOCEAN-IPSL (2006) 
     58   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    6359   !! $Id$ 
    6460   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    8177      !! 
    8278      !! ** Method  :   computes the average of some variables and write 
    83       !!      it in the NetCDF ouput files 
    84       !!      CAUTION: the sea-ice time-step must be an integer fraction 
    85       !!      of a day 
     79      !!              it in the NetCDF ouput files 
     80      !! CAUTION: the sea-ice time-step must be an integer fraction of a day 
    8681      !!------------------------------------------------------------------- 
    8782      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    8883      !! 
    89       INTEGER  ::   ji, jj, jf                      ! dummy loop indices 
     84      INTEGER  ::   ji, jj, jf   ! dummy loop indices 
    9085      CHARACTER(len = 40)  ::   clhstnam, clop 
    91       REAL(wp) ::   zsto, zjulian, zout,   &  ! temporary scalars 
    92          &          zindh, zinda, zindb, ztmu 
     86      REAL(wp) ::   zsto, zjulian, zout          ! temporary scalars 
     87      REAL(wp) ::   zindh, zinda, zindb, ztmu 
    9388      REAL(wp), DIMENSION(1)                ::   zdept 
    9489      REAL(wp), DIMENSION(jpi,jpj)          ::   zfield 
    … …  
    132127      DO jj = 2 , jpjm1 
    133128         DO ji = 1 , jpim1   ! NO vector opt. 
    134             zindh  = MAX( zzero , SIGN( zone , hicif(ji,jj) * (1.0 - frld(ji,jj) ) - 0.10 ) ) 
    135             zinda  = MAX( zzero , SIGN( zone , ( 1.0 - frld(ji,jj) ) - 0.10 ) ) 
     129            zindh  = MAX( rzero , SIGN( rone , hicif(ji,jj) * (1.0 - frld(ji,jj) ) - 0.10 ) ) 
     130            zinda  = MAX( rzero , SIGN( rone , ( 1.0 - frld(ji,jj) ) - 0.10 ) ) 
    136131            zindb  = zindh * zinda 
    137             ztmu   = MAX( 0.5 * zone , ( tmu(ji,jj) + tmu(ji+1,jj) + tmu(ji,jj+1) + tmu(ji+1,jj+1) ) )  
     132            ztmu   = MAX( 0.5 * rone , ( tmu(ji,jj) + tmu(ji+1,jj) + tmu(ji,jj+1) + tmu(ji+1,jj+1) ) )  
    138133            zcmo(ji,jj,1)  = hsnif (ji,jj) 
    139134            zcmo(ji,jj,2)  = hicif (ji,jj) 
    … …  
    143138            zcmo(ji,jj,6)  = fbif  (ji,jj) 
    144139            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) ) & 
    146                                   / ztmu  
     140               &                      + u_ice(ji,jj+1) * tmu(ji,jj+1) + u_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     141               &                  / ztmu  
    147142 
    148143            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) ) & 
    150                                   / ztmu 
     144               &                      + v_ice(ji,jj+1) * tmu(ji,jj+1) + v_ice(ji+1,jj+1) * tmu(ji+1,jj+1) ) & 
     145               &                  / ztmu 
    151146            zcmo(ji,jj,9)  = sst_m(ji,jj) 
    152147            zcmo(ji,jj,10) = sss_m(ji,jj) 
    … …  
    200195      !! ** input   :   Namelist namicewri 
    201196      !!------------------------------------------------------------------- 
    202       INTEGER ::   nf      ! ??? 
     197      INTEGER ::   jf      ! dummy loop indices 
    203198      TYPE FIELD  
    204199         CHARACTER(len = 35) :: ztitle  
    … …  
    209204         REAL                :: zcadd         
    210205      END TYPE FIELD 
    211       TYPE(FIELD) ::  & 
    212          field_1 , field_2 , field_3 , field_4 , field_5 , field_6 ,   & 
    213          field_7 , field_8 , field_9 , field_10, field_11, field_12,   & 
    214          field_13, field_14, field_15, field_16, field_17, field_18,   & 
    215          field_19 
     206      TYPE(FIELD) ::   field_1 , field_2 , field_3 , field_4 , field_5 , field_6   
     207      TYPE(FIELD) ::   field_7 , field_8 , field_9 , field_10, field_11, field_12 
     208      TYPE(FIELD) ::   field_13, field_14, field_15, field_16, field_17, field_18, field_19 
    216209      TYPE(FIELD) , DIMENSION(jpnoumax) :: zfield 
    217  
    218       NAMELIST/namiceout/ noumef, & 
    219          field_1 , field_2 , field_3 , field_4 , field_5 , field_6 ,   & 
    220          field_7 , field_8 , field_9 , field_10, field_11, field_12,   & 
    221          field_13, field_14, field_15, field_16, field_17, field_18,   & 
    222          field_19 
     210      !! 
     211      NAMELIST/namiceout/ noumef,                                          & 
     212         field_1 , field_2 , field_3 , field_4 , field_5 , field_6 ,       & 
     213         field_7 , field_8 , field_9 , field_10, field_11, field_12,       & 
     214         field_13, field_14, field_15, field_16, field_17, field_18, field_19 
    223215      !!------------------------------------------------------------------- 
    224216 
    … …  
    246238      zfield(19) = field_19 
    247239       
    248       DO nf = 1, noumef 
    249          titn  (nf) = zfield(nf)%ztitle 
    250          nam   (nf) = zfield(nf)%zname 
    251          uni   (nf) = zfield(nf)%zunit 
    252          nc    (nf) = zfield(nf)%znc 
    253          cmulti(nf) = zfield(nf)%zcmulti 
    254          cadd  (nf) = zfield(nf)%zcadd 
     240      DO jf = 1, noumef 
     241         titn  (jf) = zfield(jf)%ztitle 
     242         nam   (jf) = zfield(jf)%zname 
     243         uni   (jf) = zfield(jf)%zunit 
     244         nc    (jf) = zfield(jf)%znc 
     245         cmulti(jf) = zfield(jf)%zcmulti 
     246         cadd  (jf) = zfield(jf)%zcadd 
    255247      END DO 
    256248 
    … …  
    262254         WRITE(numout,*) '           title                            name     unit      Saving (1/0) ',   & 
    263255            &            '    multiplicative constant       additive constant ' 
    264          DO nf = 1 , noumef          
    265             WRITE(numout,*) '   ', titn(nf), '   ', nam(nf),'      ', uni(nf),'  ', nc(nf),'        ', cmulti(nf),   & 
    266                &       '        ', cadd(nf) 
     256         DO jf = 1 , noumef          
     257            WRITE(numout,*) '   ', titn(jf), '   ', nam(jf),'      ', uni(jf),'  ', nc(jf),'        ', cmulti(jf),   & 
     258               &       '        ', cadd(jf) 
    267259         END DO 
    268260      ENDIF 
    … …  
    281273      !!        Used to find errors in the initial state or save the last 
    282274      !!      ocean state in case of abnormal end of a simulation 
    283       !! 
    284       !! History : 
    285       !!   2.0  !  2009-06  (B. Lemaire) 
    286275      !!---------------------------------------------------------------------- 
    287276      INTEGER, INTENT( in ) ::   kt               ! ocean time-step index) 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/LIM_SRC_3/limtab.F90

    r1156 r1855  
    22   !!====================================================================== 
    33   !!                       ***  MODULE limtab   *** 
    4    !!              transform 1D (2D) array to a 2D (1D) table 
     4   !!   LIM-3 ice model : transform 1D (2D) array to a 2D (1D) array 
    55   !!====================================================================== 
    66#if defined key_lim3 
    77   !!---------------------------------------------------------------------- 
    8    !!   'key_lim3'                                      LIM3 sea-ice model 
     8   !!   'key_lim3' :                                  LIM 2.0 sea-ice model 
    99   !!---------------------------------------------------------------------- 
    1010   !!   tab_2d_1d  : 2-D to 1-D 
    1111   !!   tab_1d_2d  : 1-D to 2-D 
    1212   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    1413   USE par_kind 
    1514 
    … …  
    1716   PRIVATE 
    1817 
    19    !! * Routine accessibility 
    20    PUBLIC tab_2d_1d  ! called by lim_ther 
    21    PUBLIC tab_1d_2d  ! called by lim_ther 
     18   PUBLIC   tab_2d_1d   ! called by lim_thd 
     19   PUBLIC   tab_1d_2d   ! called by lim_thd 
    2220 
    2321   !!---------------------------------------------------------------------- 
    24    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
     22   !! NEMO/LIM 3.3,  UCL-LOCEAN-IPSL (2010)  
    2523   !! $Id$ 
    2624   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    … …  
    2826CONTAINS 
    2927 
    30    SUBROUTINE tab_2d_1d ( ndim1d, tab1d, tab2d, ndim2d_x, ndim2d_y, tab_ind ) 
    31  
    32       INTEGER, INTENT(in) :: & 
    33          ndim1d, ndim2d_x, ndim2d_y 
    34  
    35       REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT(in) ::  & 
    36          tab2d 
    37  
    38       INTEGER, DIMENSION ( ndim1d), INTENT ( in) :: & 
    39          tab_ind 
    40  
    41       REAL(wp), DIMENSION(ndim1d), INTENT ( out) ::  &  
    42          tab1d 
    43  
    44       INTEGER ::  & 
    45          jn , jid, jjd 
    46  
    47       DO jn = 1, ndim1d 
    48          jid        = MOD( tab_ind(jn) - 1, ndim2d_x ) + 1 
    49          jjd        = ( tab_ind(jn) - 1 ) / ndim2d_x + 1 
    50          tab1d( jn) = tab2d( jid, jjd) 
    51       END DO 
    52  
     28   SUBROUTINE tab_2d_1d( kdim1d, ptab1d, ptab2d, kdim2d_i, kdim2d_j, ptab_ind ) 
     29      !!------------------------------------------------------------------- 
     30      INTEGER , INTENT(in   )                                ::   kdim1d 
     31      INTEGER , INTENT(in   )                                ::   kdim2d_i, kdim2d_j 
     32      REAL(wp), INTENT(in   ), DIMENSION(kdim2d_i, kdim2d_j) ::   ptab2d 
     33      INTEGER , INTENT(in   ), DIMENSION(kdim1d)             ::   ptab_ind 
     34      REAL(wp), INTENT(  out), DIMENSION(kdim1d)             ::   ptab1d 
     35      !! 
     36      INTEGER ::   jn , jid, jjd   ! dummy loop indices 
     37      !!------------------------------------------------------------------- 
     38      !  
     39      DO jn = 1, kdim1d 
     40         jid = MOD( ptab_ind(jn) - 1, kdim2d_i ) + 1 
     41         jjd = ( ptab_ind(jn) - 1 ) / kdim2d_i + 1 
     42         ptab1d(jn) = ptab2d(jid,jjd) 
     43      END DO  
     44      ! 
    5345   END SUBROUTINE tab_2d_1d 
    5446 
    5547 
    56    SUBROUTINE tab_1d_2d ( ndim1d, tab2d, tab_ind, tab1d, ndim2d_x, ndim2d_y ) 
    57  
    58       INTEGER, INTENT ( in) :: & 
    59          ndim1d, ndim2d_x, ndim2d_y 
    60  
    61       INTEGER, DIMENSION (ndim1d) , INTENT (in) :: & 
    62          tab_ind 
    63  
    64       REAL(wp), DIMENSION(ndim1d), INTENT (in) ::  & 
    65          tab1d   
    66  
    67       REAL(wp), DIMENSION (ndim2d_x, ndim2d_y), INTENT ( out) :: & 
    68          tab2d 
    69  
    70       INTEGER :: & 
    71          jn, jid, jjd 
    72  
    73       DO jn = 1, ndim1d 
    74          jid             = MOD( tab_ind(jn) - 1, ndim2d_x) + 1 
    75          jjd             =    ( tab_ind(jn) - 1 ) / ndim2d_x  + 1 
    76          tab2d(jid, jjd) = tab1d( jn) 
     48   SUBROUTINE tab_1d_2d( kdim1d, ptab2d, ptab_ind, ptab1d, kdim2d_i, kdim2d_j ) 
     49      !!------------------------------------------------------------------- 
     50      INTEGER , INTENT(in   )                                ::   kdim1d 
     51      INTEGER , INTENT(in   )                                ::   kdim2d_i, kdim2d_j 
     52      INTEGER , INTENT(in   ), DIMENSION(kdim1d)             ::   ptab_ind 
     53      REAL(wp), INTENT(in   ), DIMENSION(kdim1d)             ::   ptab1d 
     54      REAL(wp), INTENT(  out), DIMENSION(kdim2d_i, kdim2d_j) ::   ptab2d 
     55      !! 
     56      INTEGER ::   jn, jid, jjd   ! dummy loop indices 
     57      !!------------------------------------------------------------------- 
     58      ! 
     59      DO jn = 1, kdim1d 
     60         jid = MOD( ptab_ind(jn) - 1, kdim2d_i) + 1 
     61         jjd =    ( ptab_ind(jn) - 1 ) / kdim2d_i  + 1 
     62         ptab2d(jid,jjd) = ptab1d(jn) 
    7763      END DO 
    78  
     64      ! 
    7965   END SUBROUTINE tab_1d_2d 
    8066 
     67#else 
     68   !!---------------------------------------------------------------------- 
     69   !!   Default option          Dummy module       NO LIM 2.0 sea-ice model 
     70   !!---------------------------------------------------------------------- 
    8171#endif 
     72 
     73   !!====================================================================== 
    8274END MODULE limtab 

  • branches/DEV_r1837_mass_heat_salt_fluxes/NEMO/OPA_SRC/SBC/sbc_ice.F90

    r1482 r1855  
    4141   REAL(wp), PUBLIC, DIMENSION(jpi,jpj,jpl) ::   alb_ice   !: albedo of ice 
    4242 
    43    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   utau_ice    !: u-stress over ice (I-point for LIM2 or U,V-point for LIM3)   [N/m2] 
    44    REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   vtau_ice    !: v-stress over ice (I-point for LIM2 or U,V-point for LIM3)   [N/m2] 
    4543   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fr1_i0      !: 1st fraction of sol. rad. which penetrate inside the ice cover 
    4644   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   fr2_i0      !: 2nd fraction of sol. rad. which penetrate inside the ice cover 
    4745   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   emp_ice     !: solid freshwater budget over ice: sublivation - snow 
     46   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   utau_ice    !: u-stress over ice (I-point for LIM2 or U,V-point for LIM3)   [N/m2] 
     47   REAL(wp), PUBLIC, DIMENSION(jpi,jpj) ::   vtau_ice    !: v-stress over ice (I-point for LIM2 or U,V-point for LIM3)   [N/m2] 
    4848 
    4949# if defined key_lim3 
Note: See TracChangeset for help on using the changeset viewer.