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

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

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

Update NEMOGCM from branch nemo_v3_3_beta

Location:
trunk/NEMOGCM/NEMO/LIM_SRC_3
Files:
29 edited

Legend:

Unmodified
Added
Removed

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/dom_ice.F90

    r1608 r2528  
    2626 
    2727   !!---------------------------------------------------------------------- 
    28    !! NEMO/LIM 3.2, UCL-ASTR-LOCEAN-IPSL (2009) 
     28   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    2929   !! $Id$ 
    30    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     30   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    3131   !!====================================================================== 
    3232END MODULE dom_ice 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/ice.F90

    r1471 r2528  
    11MODULE ice 
     2   !!====================================================================== 
     3   !!                        ***  MODULE ice  *** 
     4   !! LIM-3 Sea Ice physics:  diagnostics variables of ice defined in memory 
     5   !!===================================================================== 
     6   !! History :  3.0  ! 2008-03  (M. Vancoppenolle) : original code LIM-3 
     7   !!---------------------------------------------------------------------- 
    28#if defined key_lim3 
    39   !!---------------------------------------------------------------------- 
    410   !!   'key_lim3' :                                   LIM3 sea-ice model 
    511   !!---------------------------------------------------------------------- 
    6    !! History : 
    7    !!   2.0  !  03-08  (C. Ethe)  F90: Free form and module 
    8    !!   3.0  !  08-03  (M. Vancoppenolle) : LIM3 ! 
    9    !!---------------------------------------------------------------------- 
    10    !!  LIM 3.0, UCL-LOCEAN-IPSL (2005) 
    11    !! $Id$ 
    12    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
    13    !!---------------------------------------------------------------------- 
    14    !! * Modules used 
    1512   USE par_ice          ! LIM sea-ice parameters 
    1613 
    1714   IMPLICIT NONE 
    1815   PRIVATE 
    19    !! 
     16 
    2017   !!====================================================================== 
    21    !!                        ***  MODULE ice  *** 
    22    !! 
    23    !!                           ************** 
    24    !!                           * L I M  3.0 * 
    25    !!                           ************** 
    26    !! 
    27    !!                         ''in ice we trust'' 
    28    !! 
    29    !!                   This module contains the sea ice  
    30    !!                 diagnostics variables of ice defined  
    31    !!                             in memory 
    32    !! 
    33    !!====================================================================== 
    34    !! 
    3518   !! LIM3 by the use of sweat, agile fingers and sometimes brain juice,  
    3619   !!  was developed in Louvain-la-Neuve by :  
    … …  
    4831   !! 
    4932   !!    * Gurvan Madec, Claude Talandier, Christian Ethe  
    50    !!      and Rachid Benshila (LOCEAN-IPSL, France) 
     33   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    5134   !!    * Xavier Fettweis (UCL-ASTR), Ralph Timmermann (AWI, Germany) 
    5235   !!    * Bill Lipscomb (LANL), Cecilia Bitz (UWa)  
    … …  
    6447   !!    * Bouillon et al., in prep for 2008. 
    6548   !! 
    66    !!    Or the reference manual, that should be available by 2009 
    67    !! 
     49   !!    Or the reference manual, that should be available by 2010 
    6850   !!====================================================================== 
    6951   !!                                                                     | 
    70    !!            *****************************************                | 
    71    !!            *                                       *                | 
    72    !! ************ I C E   S T A T E   V A R I A B L E S **************** | 
    73    !!            *                                       *                | 
    74    !!            *****************************************                | 
     52   !!              I C E   S T A T E   V A R I A B L E S                  | 
    7553   !!                                                                     | 
    7654   !! Introduction :                                                      | 
    7755   !! --------------                                                      | 
    78    !!                                                                     | 
    7956   !! Every ice-covered grid cell is characterized by a series of state   | 
    8057   !! variables. To account for unresolved spatial variability in ice     | 
    … …  
    130107   !!                                                                     | 
    131108   !! ** Global variables                                                 | 
    132    !!                                                                     | 
    133109   !!-------------|-------------|---------------------------------|-------| 
    134110   !! a_i         | a_i_b       |    Ice concentration            |       | 
    … …  
    145121   !!                                                                     | 
    146122   !! ** Equivalent variables                                             | 
    147    !!                                                                     | 
    148123   !!-------------|-------------|---------------------------------|-------| 
    149124   !!                                                                     | 
    … …  
    179154   !! et_s        !      -      !    Total snow enthalpy          | 10^9 J|  
    180155   !! bv_i        !      -      !    Mean relative brine volume   | ???   |  
    181    !!                                                                     | 
    182    !!                                                                     | 
    183156   !!===================================================================== 
    184157 
    185    LOGICAL, PUBLIC ::    & 
    186       con_i = .false.           ! switch for conservation test 
     158   LOGICAL, PUBLIC ::   con_i = .false.   ! switch for conservation test 
    187159 
    188160   !!-------------------------------------------------------------------------- 
    189161   !! * Share Module variables 
    190162   !!-------------------------------------------------------------------------- 
    191    REAL(wp), PUBLIC ::   rdt_ice      !: ice time step 
     163   INTEGER , PUBLIC ::   nstart    !: iteration number of the begining of the run  
     164   INTEGER , PUBLIC ::   nlast     !: iteration number of the end of the run  
     165   INTEGER , PUBLIC ::   nitrun    !: number of iteration 
     166   INTEGER , PUBLIC ::   numit     !: iteration number 
     167   REAL(wp), PUBLIC ::   tpstot    !: time of the run in seconds 
     168   REAL(wp), PUBLIC ::   rdt_ice   !: ice time step 
    192169 
    193170   INTEGER , PUBLIC ::   &     !!: ** ice-dynamic namelist (namicedyn) ** 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/iceini.F90

    r2477 r2528  
    44   !!   Sea-ice model : LIM Sea ice model Initialization 
    55   !!====================================================================== 
     6   !! History :  3.0  ! 2008-03  (M. Vancoppenolle) LIM-3 original code 
     7   !!            3.3  ! 2010-12  (G. Madec) add call to lim_thd_init and lim_thd_sal_init 
     8   !!---------------------------------------------------------------------- 
    69#if defined key_lim3 
    710   !!---------------------------------------------------------------------- 
    … …  
    1013   !!   ice_init       : sea-ice model initialization 
    1114   !!---------------------------------------------------------------------- 
    12    USE dom_oce 
    13    USE in_out_manager 
    14    USE sbc_oce         ! Surface boundary condition: ocean fields 
    15    USE sbc_ice         ! Surface boundary condition: ice fields 
    16    USE phycst          ! Define parameters for the routines 
    17    USE ice 
    18    USE limmsh 
    19    USE limistate 
    20    USE limthd         ! LIM: ice thermodynamics  
    21    USE limthd_sal     ! LIM: ice thermodynamics: salinity  
    22    USE limrst 
    23    USE par_ice 
    24    USE limvar 
    25    USE lib_mpp 
     15   USE phycst         ! physical constants 
     16   USE dom_oce        ! ocean domain 
     17   USE sbc_oce        ! Surface boundary condition: ocean fields 
     18   USE sbc_ice        ! Surface boundary condition: ice fields 
     19   USE ice            ! LIM: sea-ice variables 
     20   USE limmsh         ! LIM: mesh 
     21   USE limistate      ! LIM: initial state 
     22   USE limrst         ! LIM: restart 
     23   USE limthd         ! LIM: ice thermodynamics 
     24   USE limthd_sal     ! LIM: ice thermodynamics: salinity 
     25   USE par_ice        ! LIM: sea-ice parameters 
     26   USE limvar         ! LIM: variables 
     27   USE in_out_manager ! I/O manager 
     28   USE lib_mpp        ! MPP library 
    2629 
    2730   IMPLICIT NONE 
    2831   PRIVATE 
    2932 
    30    PUBLIC ice_init                 ! called by opa.F90 
    31    !!---------------------------------------------------------------------- 
    32    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
     33   PUBLIC   ice_init   ! called by opa.F90 
     34 
     35   !!---------------------------------------------------------------------- 
     36   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3337   !! $Id$ 
    34    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
    35    !!---------------------------------------------------------------------- 
    36  
     38   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     39   !!---------------------------------------------------------------------- 
    3740CONTAINS 
    3841 
    … …  
    4245      !! 
    4346      !! ** purpose :    
    44       !! 
    45       !! History : 
    46       !!   2.0  !  02-08  (G. Madec)  F90: Free form and modules 
    47       !!   3.0  !  08-03  (M. Vancop) ITD, salinity, EVP-C 
    4847      !!---------------------------------------------------------------------- 
    49  
    50       ! Open the namelist file  
     48      ! 
     49      !                                ! Open the namelist file  
    5150      CALL ctl_opn( numnam_ice, 'namelist_ice', 'OLD', 'FORMATTED', 'SEQUENTIAL', -1, numout, lwp ) 
    52  
    53       CALL ice_run                    !  read in namelist some run parameters 
     51      ! 
     52      CALL ice_run                     ! namelist read some ice run parameters 
    5453      ! 
    5554      CALL lim_thd_init                ! namelist read ice thermodynics parameters 
    … …  
    6362      CALL lim_itd_ini                 ! initialize the ice thickness distribution 
    6463 
    65      !                                ! Initial sea-ice state 
     64      !                                ! Initial sea-ice state 
    6665      IF( .NOT.ln_rstart ) THEN              ! start from rest 
    6766         numit = 0 
    … …  
    7978      fr_i(:,:) = at_i(:,:)           ! initialisation of sea-ice fraction 
    8079      ! 
    81       nstart = numit  + nn_fsbc 
    82       nitrun = nitend - nit000 + 1 
    83       nlast  = numit  + nitrun 
     80      nstart = numit  + nn_fsbc       
     81      nitrun = nitend - nit000 + 1  
     82      nlast  = numit  + nitrun  
    8483      ! 
    8584      IF( nstock == 0  )   nstock = nlast + 1 
    86  
     85      ! 
    8786   END SUBROUTINE ice_init 
     87 
    8888 
    8989   SUBROUTINE ice_run 
    … …  
    9797      !! 
    9898      !! ** input   :   Namelist namicerun 
    99       !! 
    100       !! history : 
    101       !!   2.0  !  03-08 (C. Ethe)  Original code 
    102       !!   3.0  !  08-03 (M. Vancop) LIM3 
    10399      !!------------------------------------------------------------------- 
    104100      NAMELIST/namicerun/ cn_icerst_in, cn_icerst_out, ln_limdyn, acrit, hsndif, hicdif, cai, cao, ln_nicep 
    105101      !!------------------------------------------------------------------- 
    106  
    107       !                                           ! Read Namelist namicerun  
    108       REWIND ( numnam_ice ) 
    109       READ   ( numnam_ice , namicerun ) 
    110  
     102      !                     
     103      REWIND( numnam_ice )                ! Read Namelist namicerun  
     104      READ  ( numnam_ice , namicerun ) 
     105      ! 
    111106      IF( lk_mpp .AND. ln_nicep ) THEN 
    112107         ln_nicep = .FALSE. 
    113108         CALL ctl_warn( 'ice_run : specific control print for LIM3 desactivated with MPI' ) 
    114109      ENDIF        
    115  
    116       IF(lwp) THEN 
     110      ! 
     111      IF(lwp) THEN                        ! control print 
    117112         WRITE(numout,*) 
    118113         WRITE(numout,*) 'ice_run : ice share parameters for dynamics/advection/thermo of sea-ice' 
    … …  
    124119         WRITE(numout,*) '   atmospheric drag over sea ice                           = ', cai 
    125120         WRITE(numout,*) '   atmospheric drag over ocean                             = ', cao 
    126          WRITE(numout,*) '   Several ice points in the ice or not in ocean.output = ', ln_nicep 
    127       ENDIF 
    128  
     121         WRITE(numout,*) '   Several ice points in the ice or not in ocean.output    = ', ln_nicep 
     122      ENDIF 
     123      ! 
    129124   END SUBROUTINE ice_run 
     125 
    130126 
    131127   SUBROUTINE lim_itd_ini 
    132128      !!------------------------------------------------------------------ 
    133129      !!                ***  ROUTINE lim_itd_ini *** 
    134       !! ** Purpose : 
    135       !!            Initializes the ice thickness distribution 
    136       !! ** Method  : 
    137       !!            Very simple. Currently there are no ice types in the 
    138       !!            model... 
    139       !! 
    140       !! ** Arguments : 
    141       !!           kideb , kiut : Starting and ending points on which the 
    142       !!                         the computation is applied 
    143       !! 
    144       !! ** Inputs / Ouputs : (global commons) 
    145       !! 
    146       !! ** External : 
    147       !! 
    148       !! ** References : 
    149       !! 
    150       !! ** History : 
    151       !!           (12-2005) Martin Vancoppenolle 
    152       !! 
     130      !! 
     131      !! ** Purpose :   Initializes the ice thickness distribution 
     132      !! ** Method  :   ... 
    153133      !!------------------------------------------------------------------ 
    154       !! * Arguments 
    155  
    156       !! * Local variables 
    157       INTEGER ::   jl,       &   ! ice category dummy loop index 
    158          jm            ! ice types    dummy loop index 
    159  
    160       REAL(wp)  ::           &  ! constant values 
    161          zeps      =  1.0e-10,   & ! 
    162          zc1                 ,   & ! 
    163          zc2                 ,   & ! 
    164          zc3                 ,   & ! 
    165          zx1 
    166  
     134      INTEGER  ::   jl, jm               ! dummy loop index 
     135      REAL(wp) ::   zc1, zc2, zc3, zx1   ! local scalars 
     136      !!------------------------------------------------------------------ 
     137 
     138      IF(lwp) WRITE(numout,*) 
    167139      IF(lwp) WRITE(numout,*) 'lim_itd_ini : Initialization of ice thickness distribution ' 
    168140      IF(lwp) WRITE(numout,*) '~~~~~~~~~~~~' 
    169  
    170       !!-- End of declarations 
    171       !!------------------------------------------------------------------------------ 
    172141 
    173142      !------------------------------------------------------------------------------! 
    … …  
    207176      !- Thickness categories boundaries  
    208177      !---------------------------------- 
    209       hi_max(:) = 0.0 
    210       hi_max_typ(:,:) = 0.0 
     178      hi_max(:) = 0._wp 
     179      hi_max_typ(:,:) = 0._wp 
    211180 
    212181      !- Type 1 - undeformed ice 
    213       zc1 = 3./REAL(ice_cat_bounds(1,2)-ice_cat_bounds(1,1)+1) 
    214       zc2 = 10.0*zc1 
    215       zc3 = 3.0 
     182      zc1 =  3._wp / REAL( ice_cat_bounds(1,2) - ice_cat_bounds(1,1) + 1 , wp ) 
     183      zc2 = 10._wp * zc1 
     184      zc3 =  3._wp 
    216185 
    217186      DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
    218          zx1 = REAL(jl-1) / REAL(ice_cat_bounds(1,2)-ice_cat_bounds(1,1)+1) 
    219          hi_max(jl) = hi_max(jl-1) + zc1 + zc2 * (1.0 + TANH ( zc3 * (zx1 - 1.0 ) ) ) 
     187         zx1 = REAL( jl-1 , wp ) / REAL( ice_cat_bounds(1,2) - ice_cat_bounds(1,1) + 1 , wp ) 
     188         hi_max(jl) = hi_max(jl-1) + zc1 + zc2 * (1._wp + TANH( zc3 * (zx1 - 1._wp ) ) ) 
    220189      END DO 
    221190 
    222191      !- Fill in the hi_max_typ vector, useful in other circumstances 
    223       ! Tricky trick 
    224       ! hi_max_typ is actually not used in the code and will be removed in a 
    225       ! next flyspray at this time, the tricky trick will also be removed 
    226       ! Martin, march 08 
     192      ! Tricky trick: hi_max_typ is actually not used in the code and will be removed in a 
     193      ! next flyspray at this time, the tricky trick will also be removed (Martin, march 08) 
    227194      DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
    228195         hi_max_typ(jl,1) = hi_max(jl) 
    … …  
    239206         END DO 
    240207      ENDIF 
    241  
     208      ! 
    242209      DO jl = 1, jpl 
    243          hi_mean(jl) = ( hi_max(jl) + hi_max(jl-1) ) / 2.0 
    244       END DO 
    245  
     210         hi_mean(jl) = ( hi_max(jl) + hi_max(jl-1) ) * 0.5_wp 
     211      END DO 
     212      ! 
    246213      tn_ice(:,:,:) = t_su(:,:,:) 
    247  
     214      ! 
    248215   END SUBROUTINE lim_itd_ini 
    249216 
    … …  
    255222   SUBROUTINE ice_init        ! Empty routine 
    256223   END SUBROUTINE ice_init 
    257  
    258    SUBROUTINE lim_itd_ini 
    259    END SUBROUTINE lim_itd_ini 
    260224#endif 
    261225 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limadv.F90

    r1530 r2528  
    3434#  include "vectopt_loop_substitute.h90" 
    3535   !!---------------------------------------------------------------------- 
    36    !! NEMO/LIM 3.2,  UCL-LOCEAN-IPSL (2009)  
     36   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3737   !! $Id$ 
    38    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     38   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    3939   !!---------------------------------------------------------------------- 
    4040 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limcons.F90

    r1465 r2528  
    3636   !! * Module variables 
    3737   !!---------------------------------------------------------------------- 
    38    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2005) 
     38   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3939   !! $Id$ 
    40    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     40   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    4141   !!---------------------------------------------------------------------- 
    4242 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limdia.F90

    r2477 r2528  
    1414   !!   lim_dia_init   : initialization and namelist read 
    1515   !!---------------------------------------------------------------------- 
    16    !! * Modules used 
    17    USE phycst 
    18    USE in_out_manager 
    19    USE par_ice         ! ice parameters 
    20    USE sbc_ice         ! ice variables 
    21    USE daymod 
    22    USE dom_ice 
    23    USE ice 
    24    USE dom_oce 
    25    USE sbc_oce         ! Surface boundary condition: ocean fields 
    26    USE dom_oce 
    27    USE lib_mpp 
    28    USE in_out_manager 
    29   
     16   USE ice             ! LIM-3: sea-ice variable 
     17   USE par_ice         ! LIM-3: ice parameters 
     18   USE dom_ice         ! LIM-3: sea-ice domain 
     19   USE dom_oce         ! ocean domain 
     20   USE sbc_oce         ! surface boundary condition: ocean fields 
     21   USE daymod          ! model calendar 
     22   USE phycst          ! physical constant 
     23   USE in_out_manager  ! I/O manager 
     24   USE lib_mpp         ! MPP library 
     25    
    3026   IMPLICIT NONE 
    3127   PRIVATE 
    … …  
    4642   INTEGER  ::   nbvt                ! number of time variables 
    4743   INTEGER  ::   naveg               ! number of step for accumulation before averaging 
    48    REAL(wp) ::   epsi06 = 1.e-6   ! small number 
     44   REAL(wp) ::   epsi06 = 1.e-6_wp   ! small number 
    4945 
    5046   CHARACTER(len= 8) ::   fmtinf = '1PE13.5 '   ! format of the output values   
    … …  
    8480      ! 0) date from the minimum of ice extent 
    8581      !--------------------------------------- 
    86       zday_min = 273.        ! zday_min = date of minimum extent, here September 30th 
    87       zday = REAL(numit-nit000) * rdt_ice / ( 86400. * REAL(nn_fsbc) ) 
     82      zday_min = 273._wp        ! zday_min = date of minimum extent, here September 30th 
     83      zday = REAL(numit-nit000,wp) * rdt_ice / ( 86400._wp * REAL(nn_fsbc,wp) ) 
    8884      ! 
    8985      IF( zday > zday_min ) THEN   ;   zshift_date  =  zday - zday_min 
    … …  
    9793 
    9894      DO jv = nbvt + 1, nvinfo      ! put everything to zero 
    99          vinfor(jv) = 0. 
     95         vinfor(jv) = 0._wp 
    10096      END DO 
    10197 
    … …  
    108104         DO ji = fs_2, fs_jpim1   ! vector opt. 
    109105            IF( tms(ji,jj) == 1 ) THEN 
    110                vinfor(3)  = vinfor(3)  + at_i(ji,jj)*aire(ji,jj) * 1.e-12 !ice area 
    111                IF (at_i(ji,jj).GT.0.15) vinfor(5) = vinfor(5) + aire(ji,jj) * 1.e-12 !ice extent 
    112                vinfor(7)  = vinfor(7)  + vt_i(ji,jj)*aire(ji,jj) * 1.e-12 !ice volume 
    113                vinfor(9)  = vinfor(9)  + vt_s(ji,jj)*aire(ji,jj) * 1.e-12 !snow volume 
    114                vinfor(15) = vinfor(15) + ot_i(ji,jj) *vt_i(ji,jj)*aire(ji,jj) * 1.e-12 !mean age 
    115                vinfor(29) = vinfor(29) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12 !mean salinity 
     106               vinfor(3)  = vinfor(3)  + at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !ice area 
     107               IF (at_i(ji,jj).GT.0.15) vinfor(5) = vinfor(5) + aire(ji,jj) * 1.e-12_wp !ice extent 
     108               vinfor(7)  = vinfor(7)  + vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !ice volume 
     109               vinfor(9)  = vinfor(9)  + vt_s(ji,jj)*aire(ji,jj) * 1.e-12_wp !snow volume 
     110               vinfor(15) = vinfor(15) + ot_i(ji,jj) *vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !mean age 
     111               vinfor(29) = vinfor(29) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !mean salinity 
    116112               ! the computation of this diagnostic is not reliable 
    117113               vinfor(31) = vinfor(31) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
    118114                  v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12  
    119                vinfor(53) = vinfor(53) + emps(ji,jj)*aire(ji,jj) * 1.e-12 !salt flux 
    120                vinfor(55) = vinfor(55) + fsbri(ji,jj)*aire(ji,jj) * 1.e-12 !brine drainage flux 
    121                vinfor(57) = vinfor(57) + fseqv(ji,jj)*aire(ji,jj) * 1.e-12 !equivalent salt flux 
    122                vinfor(59) = vinfor(59) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) * 1.e-12  !SST 
    123                vinfor(61) = vinfor(61) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) * 1.e-12  !SSS 
    124                vinfor(65) = vinfor(65) + et_s(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12  ! snow temperature 
    125                vinfor(67) = vinfor(67) + et_i(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12       ! ice heat content 
    126                vinfor(69) = vinfor(69) + v_i(ji,jj,1)*aire(ji,jj) * 1.e-12 !ice volume 
    127                vinfor(71) = vinfor(71) + v_i(ji,jj,2)*aire(ji,jj) * 1.e-12 !ice volume 
    128                vinfor(73) = vinfor(73) + v_i(ji,jj,3)*aire(ji,jj) * 1.e-12 !ice volume 
    129                vinfor(75) = vinfor(75) + v_i(ji,jj,4)*aire(ji,jj) * 1.e-12 !ice volume 
    130                vinfor(77) = vinfor(77) + v_i(ji,jj,5)*aire(ji,jj) * 1.e-12 !ice volume 
     115               vinfor(53) = vinfor(53) + emps(ji,jj)*aire(ji,jj) * 1.e-12_wp !salt flux 
     116               vinfor(55) = vinfor(55) + fsbri(ji,jj)*aire(ji,jj) * 1.e-12_wp !brine drainage flux 
     117               vinfor(57) = vinfor(57) + fseqv(ji,jj)*aire(ji,jj) * 1.e-12_wp !equivalent salt flux 
     118               vinfor(59) = vinfor(59) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp  !SST 
     119               vinfor(61) = vinfor(61) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp  !SSS 
     120               vinfor(65) = vinfor(65) + et_s(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12_wp  ! snow temperature 
     121               vinfor(67) = vinfor(67) + et_i(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12_wp       ! ice heat content 
     122               vinfor(69) = vinfor(69) + v_i(ji,jj,1)*aire(ji,jj) * 1.e-12_wp !ice volume 
     123               vinfor(71) = vinfor(71) + v_i(ji,jj,2)*aire(ji,jj) * 1.e-12_wp !ice volume 
     124               vinfor(73) = vinfor(73) + v_i(ji,jj,3)*aire(ji,jj) * 1.e-12_wp !ice volume 
     125               vinfor(75) = vinfor(75) + v_i(ji,jj,4)*aire(ji,jj) * 1.e-12_wp !ice volume 
     126               vinfor(77) = vinfor(77) + v_i(ji,jj,5)*aire(ji,jj) * 1.e-12_wp !ice volume 
    131127               vinfor(79) = 0.0 
    132                vinfor(81) = vinfor(81) + emp(ji,jj)*aire(ji,jj) * 1.e-12 ! mass flux 
     128               vinfor(81) = vinfor(81) + emp(ji,jj)*aire(ji,jj) * 1.e-12_wp ! mass flux 
    133129            ENDIF 
    134130         END DO 
    … …  
    139135            DO ji = fs_2, fs_jpim1   ! vector opt. 
    140136               IF( tms(ji,jj) == 1 ) THEN 
    141                   vinfor(11) = vinfor(11) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 !undef def ice volume 
     137                  vinfor(11) = vinfor(11) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !undef def ice volume 
    142138               ENDIF 
    143139            END DO 
    … …  
    145141      END DO 
    146142 
    147       vinfor(13) = 0. 
     143      vinfor(13) = 0._wp 
    148144 
    149145      vinfor(15) = vinfor(15) / MAX(vinfor(7),epsi06) ! these have to be divided by total ice volume to have the 
    … …  
    169165            DO ji = fs_2, fs_jpim1   ! vector opt. 
    170166               IF( tms(ji,jj) == 1 ) THEN 
    171                   vinfor(33) = vinfor(33) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) * 1.e-12 !ice volume 
    172                   vinfor(35) = vinfor(35) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) * 1.e-12 !ice volume 
     167                  vinfor(33) = vinfor(33) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !ice volume 
     168                  vinfor(35) = vinfor(35) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !ice volume 
    173169               ENDIF 
    174170            END DO 
    … …  
    179175         DO ji = fs_2, fs_jpim1   ! vector opt. 
    180176            IF( tms(ji,jj) == 1 ) THEN 
    181                vinfor(37) = vinfor(37) + diag_sni_gr(ji,jj)*aire(ji,jj) * 1.e-12 !th growth rates 
    182                vinfor(39) = vinfor(39) + diag_lat_gr(ji,jj)*aire(ji,jj) * 1.e-12  
    183                vinfor(41) = vinfor(41) + diag_bot_gr(ji,jj)*aire(ji,jj) * 1.e-12 
    184                vinfor(43) = vinfor(43) + diag_dyn_gr(ji,jj)*aire(ji,jj) * 1.e-12  
    185                vinfor(45) = vinfor(45) + dv_dt_thd(ji,jj,5)*aire(ji,jj) * 1.e-12 
    186                vinfor(47) = vinfor(47) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12 / rdt_ice ! volume acc in OW 
     177               vinfor(37) = vinfor(37) + diag_sni_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp !th growth rates 
     178               vinfor(39) = vinfor(39) + diag_lat_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp  
     179               vinfor(41) = vinfor(41) + diag_bot_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp 
     180               vinfor(43) = vinfor(43) + diag_dyn_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp  
     181               vinfor(45) = vinfor(45) + dv_dt_thd(ji,jj,5)*aire(ji,jj) * 1.e-12_wp 
     182               vinfor(47) = vinfor(47) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12_wp / rdt_ice ! volume acc in OW 
    187183            ENDIF 
    188184         END DO 
    … …  
    193189            DO ji = fs_2, fs_jpim1   ! vector opt. 
    194190               IF( tms(ji,jj) == 1 ) THEN 
    195                   vinfor(63) = vinfor(63) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 
     191                  vinfor(63) = vinfor(63) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp 
    196192               ENDIF 
    197193            END DO 
    … …  
    209205               DO jl = 1, jpl 
    210206                  IF ((o_i(ji,jj,jl) - zshift_date).LT.0.0) THEN 
    211                      vinfor(17) = vinfor(17) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 ! FY ice area 
    212                      vinfor(25) = vinfor(25) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 ! FY ice volume 
    213                      vinfor(49) = vinfor(49) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 !FY ice salinity 
     207                     vinfor(17) = vinfor(17) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp ! FY ice area 
     208                     vinfor(25) = vinfor(25) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp ! FY ice volume 
     209                     vinfor(49) = vinfor(49) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !FY ice salinity 
    214210                     zafy = zafy + a_i(ji,jj,jl) 
    215211                  ENDIF 
    216212                  IF ((o_i(ji,jj,jl) - zshift_date).GT.0.0) THEN 
    217                      vinfor(19) = vinfor(19) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12    ! MY ice area 
    218                      vinfor(27) = vinfor(27) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 ! MY ice volume 
    219                      vinfor(51) = vinfor(51) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 !MY ice salinity 
     213                     vinfor(19) = vinfor(19) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp    ! MY ice area 
     214                     vinfor(27) = vinfor(27) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp ! MY ice volume 
     215                     vinfor(51) = vinfor(51) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !MY ice salinity 
    220216                     zamy = zamy + a_i(ji,jj,jl) 
    221217                  ENDIF 
    222218               END DO 
    223219               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN 
    224                   vinfor(21) = vinfor(21) + aire(ji,jj) * 1.e-12 ! Seasonal ice extent 
     220                  vinfor(21) = vinfor(21) + aire(ji,jj) * 1.e-12_wp ! Seasonal ice extent 
    225221               ENDIF 
    226222               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN 
    227                   vinfor(23) = vinfor(23) + aire(ji,jj) * 1.e-12 ! Perennial ice extent 
     223                  vinfor(23) = vinfor(23) + aire(ji,jj) * 1.e-12_wp ! Perennial ice extent 
    228224               ENDIF 
    229225            ENDIF 
    … …  
    245241      DO ji = 134, 138 
    246242         vinfor(83) = vinfor(83) - v_ice(ji,jj) * &  
    247             e1t(ji,jj)*at_i(ji,jj)*rdt_ice * 1.e-12 
     243            e1t(ji,jj)*at_i(ji,jj)*rdt_ice * 1.e-12_wp 
    248244         vinfor(84) = vinfor(84) - v_ice(ji,jj) * &  
    249             e1t(ji,jj)*vt_i(ji,jj)*rdt_ice * 1.e-12 
     245            e1t(ji,jj)*vt_i(ji,jj)*rdt_ice * 1.e-12_wp 
    250246      END DO 
    251247 
    … …  
    258254         DO ji = fs_2, fs_jpim1   ! vector opt. 
    259255            IF( tms(ji,jj) == 1 ) THEN 
    260                vinfor(4)  = vinfor(4)  + at_i(ji,jj)*aire(ji,jj) * 1.e-12 !ice area 
    261                IF (at_i(ji,jj).GT.0.15) vinfor(6) = vinfor(6) + aire(ji,jj) * 1.e-12 !ice extent 
    262                vinfor(8)  = vinfor(8)  + vt_i(ji,jj)*aire(ji,jj) * 1.e-12 !ice volume 
    263                vinfor(10) = vinfor(10) + vt_s(ji,jj)*aire(ji,jj) * 1.e-12 !snow volume 
    264                vinfor(16) = vinfor(16) + ot_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12 !mean age 
    265                vinfor(30) = vinfor(30) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12 !mean salinity 
     256               vinfor(4)  = vinfor(4)  + at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !ice area 
     257               IF (at_i(ji,jj).GT.0.15) vinfor(6) = vinfor(6) + aire(ji,jj) * 1.e-12_wp !ice extent 
     258               vinfor(8)  = vinfor(8)  + vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !ice volume 
     259               vinfor(10) = vinfor(10) + vt_s(ji,jj)*aire(ji,jj) * 1.e-12_wp !snow volume 
     260               vinfor(16) = vinfor(16) + ot_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !mean age 
     261               vinfor(30) = vinfor(30) + smt_i(ji,jj)*vt_i(ji,jj)*aire(ji,jj) * 1.e-12_wp !mean salinity 
    266262               ! this diagnostic is not well computed (weighted by vol instead 
    267263               ! of area) 
    268264               vinfor(32) = vinfor(32) + vt_i(ji,jj)*( u_ice(ji,jj)*u_ice(ji,jj) + &  
    269265                  v_ice(ji,jj)*v_ice(ji,jj) )*aire(ji,jj)/1.0e12 !ice vel 
    270                vinfor(54) = vinfor(54) + at_i(ji,jj)*emps(ji,jj)*aire(ji,jj) * 1.e-12 ! Total salt flux 
    271                vinfor(56) = vinfor(56) + at_i(ji,jj)*fsbri(ji,jj)*aire(ji,jj) * 1.e-12 ! Brine drainage salt flux 
    272                vinfor(58) = vinfor(58) + at_i(ji,jj)*fseqv(ji,jj)*aire(ji,jj) * 1.e-12 ! Equivalent salt flux 
    273                vinfor(60) = vinfor(60) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) * 1.e-12  !SST 
    274                vinfor(62) = vinfor(62) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) * 1.e-12  !SSS 
    275                vinfor(66) = vinfor(66) + et_s(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12 ! snow temperature 
    276                vinfor(68) = vinfor(68) + et_i(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12 ! ice enthalpy 
    277                vinfor(70) = vinfor(70) + v_i(ji,jj,1)*aire(ji,jj) * 1.e-12 !ice volume 
    278                vinfor(72) = vinfor(72) + v_i(ji,jj,2)*aire(ji,jj) * 1.e-12 !ice volume 
    279                vinfor(74) = vinfor(74) + v_i(ji,jj,3)*aire(ji,jj) * 1.e-12 !ice volume 
    280                vinfor(76) = vinfor(76) + v_i(ji,jj,4)*aire(ji,jj) * 1.e-12 !ice volume 
    281                vinfor(78) = vinfor(78) + v_i(ji,jj,5)*aire(ji,jj) * 1.e-12 !ice volume 
     266               vinfor(54) = vinfor(54) + at_i(ji,jj)*emps(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Total salt flux 
     267               vinfor(56) = vinfor(56) + at_i(ji,jj)*fsbri(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Brine drainage salt flux 
     268               vinfor(58) = vinfor(58) + at_i(ji,jj)*fseqv(ji,jj)*aire(ji,jj) * 1.e-12_wp ! Equivalent salt flux 
     269               vinfor(60) = vinfor(60) +(sst_m(ji,jj)+rt0)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp  !SST 
     270               vinfor(62) = vinfor(62) + sss_m(ji,jj)*at_i(ji,jj)*aire(ji,jj) * 1.e-12_wp  !SSS 
     271               vinfor(66) = vinfor(66) + et_s(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12_wp ! snow temperature 
     272               vinfor(68) = vinfor(68) + et_i(ji,jj)/1.0e9*aire(ji,jj) * 1.e-12_wp ! ice enthalpy 
     273               vinfor(70) = vinfor(70) + v_i(ji,jj,1)*aire(ji,jj) * 1.e-12_wp !ice volume 
     274               vinfor(72) = vinfor(72) + v_i(ji,jj,2)*aire(ji,jj) * 1.e-12_wp !ice volume 
     275               vinfor(74) = vinfor(74) + v_i(ji,jj,3)*aire(ji,jj) * 1.e-12_wp !ice volume 
     276               vinfor(76) = vinfor(76) + v_i(ji,jj,4)*aire(ji,jj) * 1.e-12_wp !ice volume 
     277               vinfor(78) = vinfor(78) + v_i(ji,jj,5)*aire(ji,jj) * 1.e-12_wp !ice volume 
    282278               vinfor(80) = 0.0 
    283                vinfor(82) = vinfor(82) + emp(ji,jj)*aire(ji,jj) * 1.e-12 ! mass flux 
     279               vinfor(82) = vinfor(82) + emp(ji,jj)*aire(ji,jj) * 1.e-12_wp ! mass flux 
    284280            ENDIF 
    285281         END DO 
    … …  
    289285         DO jj = 2, njeqm1 
    290286            DO ji = fs_2, fs_jpim1   ! vector opt. 
    291                vinfor(12) = vinfor(12) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 !undef def ice volume 
     287               vinfor(12) = vinfor(12) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !undef def ice volume 
    292288            END DO 
    293289         END DO 
    … …  
    320316            DO ji = fs_2, fs_jpim1   ! vector opt. 
    321317               IF( tms(ji,jj) == 1 ) THEN 
    322                   vinfor(34) = vinfor(34) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) * 1.e-12 !ice volume 
    323                   vinfor(36) = vinfor(36) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) * 1.e-12 !ice volume 
     318                  vinfor(34) = vinfor(34) + d_v_i_trp(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !ice volume 
     319                  vinfor(36) = vinfor(36) + d_v_i_thd(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !ice volume 
    324320               ENDIF 
    325321            END DO 
    … …  
    330326         DO ji = fs_2, fs_jpim1   ! vector opt. 
    331327            IF( tms(ji,jj) == 1 ) THEN 
    332                vinfor(38) = vinfor(38) + diag_sni_gr(ji,jj)*aire(ji,jj) * 1.e-12 !th growth rates 
    333                vinfor(40) = vinfor(40) + diag_lat_gr(ji,jj)*aire(ji,jj) * 1.e-12  
    334                vinfor(42) = vinfor(42) + diag_bot_gr(ji,jj)*aire(ji,jj) * 1.e-12 
    335                vinfor(44) = vinfor(44) + diag_dyn_gr(ji,jj)*aire(ji,jj) * 1.e-12  
    336                vinfor(46) = vinfor(46) + dv_dt_thd(ji,jj,5)*aire(ji,jj) * 1.e-12 
    337                vinfor(48) = vinfor(48) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12 / rdt_ice ! volume acc in OW 
     328               vinfor(38) = vinfor(38) + diag_sni_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp !th growth rates 
     329               vinfor(40) = vinfor(40) + diag_lat_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp  
     330               vinfor(42) = vinfor(42) + diag_bot_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp 
     331               vinfor(44) = vinfor(44) + diag_dyn_gr(ji,jj)*aire(ji,jj) * 1.e-12_wp  
     332               vinfor(46) = vinfor(46) + dv_dt_thd(ji,jj,5)*aire(ji,jj) * 1.e-12_wp 
     333               vinfor(48) = vinfor(48) + v_newice(ji,jj) *aire(ji,jj) * 1.e-12_wp / rdt_ice ! volume acc in OW 
    338334            ENDIF 
    339335         END DO 
    … …  
    344340            DO ji = fs_2, fs_jpim1   ! vector opt. 
    345341               IF( tms(ji,jj) == 1 ) THEN 
    346                   vinfor(64) = vinfor(64) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 
     342                  vinfor(64) = vinfor(64) + t_su(ji,jj,jl)*a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp 
    347343               ENDIF 
    348344            END DO 
    … …  
    356352         DO ji = fs_2, fs_jpim1   ! vector opt. 
    357353            IF( tms(ji,jj) == 1 ) THEN 
    358                zafy = 0. 
    359                zamy = 0. 
     354               zafy = 0._wp 
     355               zamy = 0._wp 
    360356               DO jl = 1, jpl 
    361                   IF( (o_i(ji,jj,jl) - zshift_date) < 0. ) THEN 
    362                      vinfor(18) = vinfor(18) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 ! FY ice area 
    363                      vinfor(26) = vinfor(26) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 ! FY ice volume 
     357                  IF( (o_i(ji,jj,jl) - zshift_date) < 0._wp ) THEN 
     358                     vinfor(18) = vinfor(18) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp ! FY ice area 
     359                     vinfor(26) = vinfor(26) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp ! FY ice volume 
    364360                     zafy = zafy + a_i(ji,jj,jl) 
    365                      vinfor(50) = vinfor(50) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 !FY ice salinity 
     361                     vinfor(50) = vinfor(50) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !FY ice salinity 
    366362                  ENDIF 
    367                   IF( (o_i(ji,jj,jl) - zshift_date) > 0. ) THEN 
    368                      vinfor(20) = vinfor(20) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12    ! MY ice area 
    369                      vinfor(28) = vinfor(28) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 
    370                      vinfor(52) = vinfor(52) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12 !FY ice salinity 
     363                  IF( (o_i(ji,jj,jl) - zshift_date) > 0._wp ) THEN 
     364                     vinfor(20) = vinfor(20) + a_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp    ! MY ice area 
     365                     vinfor(28) = vinfor(28) + v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp 
     366                     vinfor(52) = vinfor(52) + sm_i(ji,jj,jl)*v_i(ji,jj,jl)*aire(ji,jj) * 1.e-12_wp !FY ice salinity 
    371367                     zamy = zamy + a_i(ji,jj,jl) 
    372368                  ENDIF 
    373369               END DO ! jl 
    374370               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.GT.zamy)) THEN 
    375                   vinfor(22) = vinfor(22) + aire(ji,jj) * 1.e-12 ! Seasonal ice extent 
     371                  vinfor(22) = vinfor(22) + aire(ji,jj) * 1.e-12_wp ! Seasonal ice extent 
    376372               ENDIF 
    377373               IF ((at_i(ji,jj).GT.0.15).AND.(zafy.LE.zamy)) THEN 
    378                   vinfor(24) = vinfor(24) + aire(ji,jj) * 1.e-12 ! Perennial ice extent 
     374                  vinfor(24) = vinfor(24) + aire(ji,jj) * 1.e-12_wp ! Perennial ice extent 
    379375               ENDIF 
    380376            ENDIF ! tms 
    … …  
    397393      naveg = 0 
    398394      DO jv = 1, nvinfo 
    399          vinfom(jv) = 0. 
     395         vinfom(jv) = 0._wp 
    400396      END DO 
    401397      !MV      ENDIF 
    … …  
    568564      !--Initialisation of the arrays for the accumulation 
    569565      DO jv = 1, nvinfo 
    570          vinfom(jv) = 0. 
     566         vinfom(jv) = 0._wp 
    571567      END DO 
    572568      naveg = 0 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limdyn.F90

    r2477 r2528  
    44   !!   Sea-Ice dynamics :   
    55   !!====================================================================== 
     6   !! history :  1.0  ! 2002-08 (C. Ethe, G. Madec)  original VP code  
     7   !!            3.0  ! 2007-03 (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle)  LIM3: EVP-Cgrid 
     8   !!---------------------------------------------------------------------- 
    69#if defined key_lim3 
    710   !!---------------------------------------------------------------------- 
    … …  
    1114   !!    lim_dyn_init : initialization and namelist read 
    1215   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    1416   USE phycst 
    1517   USE in_out_manager  ! I/O manager 
    … …  
    2830   PRIVATE 
    2931 
    30    !! * Accessibility 
    31    PUBLIC lim_dyn  ! routine called by ice_step 
     32   PUBLIC   lim_dyn   ! routine called by ice_step 
    3233 
    3334   !! * Substitutions 
    3435#  include "vectopt_loop_substitute.h90" 
    35  
    36    !! * Module variables 
    37    REAL(wp)  ::  rone    = 1.e0   ! constant value 
    38  
    39    !!---------------------------------------------------------------------- 
    40    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
     36   !!---------------------------------------------------------------------- 
     37   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4138   !! $Id$ 
    42    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    43    !!---------------------------------------------------------------------- 
    44  
     39   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     40   !!---------------------------------------------------------------------- 
    4541CONTAINS 
    4642 
    … …  
    5753      !!              - computation of the stress at the ocean surface          
    5854      !!              - treatment of the case if no ice dynamic 
    59       !! History : 
    60       !!   1.0  !  01-04   (LIM)  Original code 
    61       !!   2.0  !  02-08   (C. Ethe, G. Madec)  F90, mpp 
    62       !!   3.0  !  2007-03 (M.A. Morales Maqueda, S. Bouillon, M. Vancoppenolle)  
    63       !!                   LIM3, EVP, C-grid 
    6455      !!------------------------------------------------------------------------------------ 
    6556      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    66       !! * Local variables 
     57      !! 
    6758      INTEGER  ::   ji, jj, jl, ja    ! dummy loop indices 
    6859      INTEGER  ::   i_j1, i_jpj       ! Starting/ending j-indices for rheology 
    69       REAL(wp) ::   zcoef             ! temporary scalar 
     60      REAL(wp) ::   zcoef             ! local scalar 
    7061      REAL(wp), DIMENSION(jpj)     ::   zind           ! i-averaged indicator of sea-ice 
    7162      REAL(wp), DIMENSION(jpj)     ::   zmsk           ! i-averaged of tmask 
    … …  
    9182         ! --------------------------------------------------- 
    9283 
    93          IF( lk_mpp .OR. nbit_cmp == 1 ) THEN                    ! mpp: compute over the whole domain 
     84         IF( lk_mpp .OR. lk_mpp_rep ) THEN                    ! mpp: compute over the whole domain 
    9485            i_j1 = 1 
    9586            i_jpj = jpj 
    … …  
    156147         zv_io(:,:) = v_ice(:,:) - ssv_m(:,:) 
    157148         ! frictional velocity at T-point 
     149         zcoef = 0.5 * cw 
    158150         DO jj = 2, jpjm1  
    159151            DO ji = fs_2, fs_jpim1   ! vector opt. 
    160                ust2s(ji,jj) = 0.5 * cw                                                          & 
    161                   &         * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
    162                   &            + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
     152               ust2s(ji,jj) = zcoef * (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
     153                  &                    + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1)   ) * tms(ji,jj) 
    163154            END DO 
    164155         END DO 
    … …  
    169160         DO jj = 2, jpjm1 
    170161            DO ji = fs_2, fs_jpim1   ! vector opt. 
    171                ust2s(ji,jj) = zcoef * tms(ji,jj) * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
    172                   &                                     + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) 
     162               ust2s(ji,jj) = zcoef * SQRT(  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
     163                  &                        + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1)   ) * tms(ji,jj) 
    173164            END DO 
    174165         END DO 
    175166         ! 
    176167      ENDIF 
    177  
    178168      CALL lbc_lnk( ust2s, 'T',  1. )   ! T-point 
    179169 
    … …  
    216206         END DO 
    217207      ENDIF 
    218  
     208      ! 
    219209   END SUBROUTINE lim_dyn 
     210 
    220211 
    221212   SUBROUTINE lim_dyn_init 
    … …  
    230221      !! 
    231222      !! ** input   :   Namelist namicedyn 
    232       !! 
    233       !! history : 
    234       !!  8.5  ! 03-08 (C. Ethe) original code 
    235       !!  9.0  ! 07-03 (MA Morales Maqueda, S. Bouillon, M. Vancoppenolle) 
    236       !!               EVP-Cgrid-LIM3 
    237223      !!------------------------------------------------------------------- 
    238224      NAMELIST/namicedyn/ epsd, alpha,     & 
    … …  
    242228      !!------------------------------------------------------------------- 
    243229 
    244       ! Define the initial parameters 
    245       ! ------------------------- 
    246  
    247       ! Read Namelist namicedyn 
    248       REWIND ( numnam_ice ) 
    249       READ   ( numnam_ice  , namicedyn ) 
    250       IF(lwp) THEN 
     230      REWIND( numnam_ice )                ! Read Namelist namicedyn 
     231      READ  ( numnam_ice  , namicedyn ) 
     232       
     233      IF(lwp) THEN                        ! control print 
    251234         WRITE(numout,*) 
    252235         WRITE(numout,*) 'lim_dyn_init : ice parameters for ice dynamics ' 
    … …  
    270253         WRITE(numout,*) '   timescale for elastic waves                      telast = ', telast 
    271254         WRITE(numout,*) '   coefficient for the solution of int. stresses  alphaevp = ', alphaevp 
    272  
    273       ENDIF 
    274  
    275       usecc2 = 1.0 / ( ecc * ecc ) 
    276       rhoco  = rau0 * cw 
     255      ENDIF 
     256      ! 
     257      IF( angvg /= 0._wp ) THEN 
     258         CALL ctl_warn( 'lim_dyn_init: turning angle for oceanic stress not properly coded for EVP ',   & 
     259            &           '(see limsbc module). We force  angvg = 0._wp'  ) 
     260         angvg = 0._wp 
     261      ENDIF 
     262       
     263      usecc2 = 1._wp / ( ecc * ecc ) 
     264      rhoco  = rau0  * cw 
    277265      angvg  = angvg * rad 
    278266      sangvg = SIN( angvg ) 
    279267      cangvg = COS( angvg ) 
    280       pstarh = pstar / 2.0 
     268      pstarh = pstar * 0.5_wp 
    281269 
    282270      !  Diffusion coefficients. 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limhdf.F90

    r1465 r2528  
    3232#  include "vectopt_loop_substitute.h90" 
    3333   !!---------------------------------------------------------------------- 
    34    !!   LIM 2.0,  UCL-LOCEAN-IPSL (2005)  
     34   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3535   !! $Id$ 
    36    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt  
     36   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    3737   !!---------------------------------------------------------------------- 
    3838 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90

    r1471 r2528  
    44   !!              Initialisation of diagnostics ice variables 
    55   !!====================================================================== 
     6   !! History :  2.0  ! 2004-01 (C. Ethe, G. Madec)  Original code 
     7   !!---------------------------------------------------------------------- 
    68#if defined key_lim3 
    79   !!---------------------------------------------------------------------- 
    … …  
    1113   !!   lim_istate_init :  initialization of ice state and namelist read 
    1214   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    14    USE phycst 
    15    USE oce             ! dynamics and tracers variables 
    16    USE dom_oce 
    17    USE sbc_oce         ! Surface boundary condition: ocean fields 
    18    USE par_ice         ! ice parameters 
    19    USE eosbn2          ! equation of state 
    20    USE in_out_manager 
    21    USE dom_ice 
    22    USE ice 
    23    USE lbclnk 
     15   USE phycst           ! physical constant 
     16   USE oce              ! dynamics and tracers variables 
     17   USE dom_oce          ! ocean domain 
     18   USE sbc_oce          ! Surface boundary condition: ocean fields 
     19   USE eosbn2           ! equation of state 
     20   USE ice              ! sea-ice variables 
     21   USE par_ice          ! ice parameters 
     22   USE dom_ice          ! sea-ice domain 
     23   USE in_out_manager   ! I/O manager 
     24   USE lbclnk           ! lateral boundary condition - MPP exchanges 
    2425 
    2526   IMPLICIT NONE 
    2627   PRIVATE 
    2728 
    28    !! * Accessibility 
    29    PUBLIC lim_istate      ! routine called by lim_init.F90 
    30  
    31    !! * Module variables 
    32    REAL(wp) ::             & !!! ** init namelist (namiceini) ** 
    33       ttest    = 2.0  ,    &  ! threshold water temperature for initial sea ice 
    34       hninn    = 0.5  ,    &  ! initial snow thickness in the north 
    35       hginn_u  = 2.5  ,    &  ! initial ice thickness in the north 
    36       aginn_u  = 0.7  ,    &  ! initial leads area in the north 
    37       hginn_d  = 5.0  ,    &  ! initial ice thickness in the north 
    38       aginn_d  = 0.25 ,    &  ! initial leads area in the north 
    39       hnins    = 0.1  ,    &  ! initial snow thickness in the south 
    40       hgins_u  = 1.0  ,    &  ! initial ice thickness in the south 
    41       agins_u  = 0.7  ,    &  ! initial leads area in the south 
    42       hgins_d  = 2.0  ,    &  ! initial ice thickness in the south 
    43       agins_d  = 0.2  ,    &  ! initial leads area in the south 
    44       sinn     = 6.301 ,   &  ! initial salinity  
    45       sins     = 6.301 
    46  
    47    REAL(wp)  ::            &  ! constant values 
    48       zzero   = 0.0     ,  & 
    49       zone    = 1.0 
    50  
    51    !!---------------------------------------------------------------------- 
    52    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008) 
     29   PUBLIC   lim_istate      ! routine called by lim_init.F90 
     30 
     31   !                                  !!** init namelist (namiceini) ** 
     32   REAL(wp) ::   ttest    = 2.0_wp     ! threshold water temperature for initial sea ice 
     33   REAL(wp) ::   hninn    = 0.5_wp     ! initial snow thickness in the north 
     34   REAL(wp) ::   hginn_u  = 2.5_wp     ! initial ice thickness in the north 
     35   REAL(wp) ::   aginn_u  = 0.7_wp     ! initial leads area in the north 
     36   REAL(wp) ::   hginn_d  = 5.0_wp     ! initial ice thickness in the north 
     37   REAL(wp) ::   aginn_d  = 0.25_wp    ! initial leads area in the north 
     38   REAL(wp) ::   hnins    = 0.1_wp     ! initial snow thickness in the south 
     39   REAL(wp) ::   hgins_u  = 1.0_wp     ! initial ice thickness in the south 
     40   REAL(wp) ::   agins_u  = 0.7_wp     ! initial leads area in the south 
     41   REAL(wp) ::   hgins_d  = 2.0_wp     ! initial ice thickness in the south 
     42   REAL(wp) ::   agins_d  = 0.2_wp     ! initial leads area in the south 
     43   REAL(wp) ::   sinn     = 6.301_wp   ! initial salinity  
     44   REAL(wp) ::   sins     = 6.301_wp   ! 
     45 
     46   !!---------------------------------------------------------------------- 
     47   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    5348   !! $Id$ 
    54    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    55    !!---------------------------------------------------------------------- 
    56  
     49   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     50   !!---------------------------------------------------------------------- 
    5751CONTAINS 
    5852 
    … …  
    6559      !! ** Method  :   restart from a state defined in a binary file 
    6660      !!                or from arbitrary sea-ice conditions 
    67       !! 
    68       !! History : 
    69       !!   2.0  !  01-04  (C. Ethe, G. Madec)  Original code 
    70       !!-------------------------------------------------------------------- 
    71  
    72       !! * Local variables 
     61      !!------------------------------------------------------------------- 
    7362      INTEGER  ::   ji, jj, jk, jl             ! dummy loop indices 
    74  
    75       REAL(wp) ::       &  ! temporary scalar 
    76          zeps6, zeps, ztmelts, & 
    77          epsi06 
    78       REAL(wp), DIMENSION(jpm) ::   & 
    79          zgfactorn, zhin, & 
    80          zgfactors, zhis 
    81       REAL(wp) ::  & 
    82          zvol, zare, zh, zh1, zh2, zh3, zan, zbn, zas, zbs  
    83       REAL(wp), DIMENSION(jpi,jpj) ::   zidto    ! ice indicator 
     63      REAL(wp) ::   zeps6, zeps, ztmelts, epsi06   ! local scalars 
     64      REAL(wp) ::  zvol, zare, zh, zh1, zh2, zh3, zan, zbn, zas, zbs  
     65      REAL(wp), DIMENSION(jpm)     ::   zgfactorn, zhin  
     66      REAL(wp), DIMENSION(jpm)     ::   zgfactors, zhis 
     67      REAL(wp), DIMENSION(jpi,jpj) ::   zidto      ! ice indicator 
    8468      !-------------------------------------------------------------------- 
    8569 
    … …  
    8771      ! 1) Preliminary things  
    8872      !-------------------------------------------------------------------- 
    89       epsi06 = 1.0e-6 
     73      epsi06 = 1.e-6_wp 
    9074 
    9175      CALL lim_istate_init     !  reading the initials parameters of the ice 
    … …  
    116100 
    117101      ! constants for heat contents 
    118       zeps   = 1.0d-20 
    119       zeps6  = 1.0d-06 
     102      zeps   = 1.e-20_wp 
     103      zeps6  = 1.e-06_wp 
    120104 
    121105      ! zgfactor for initial ice distribution 
    122       zgfactorn(:) = 0.0 
    123       zgfactors(:) = 0.0 
     106      zgfactorn(:) = 0._wp 
     107      zgfactors(:) = 0._wp 
    124108 
    125109      ! first ice type 
    126110      DO jl = ice_cat_bounds(1,1), ice_cat_bounds(1,2) 
    127          zhin (1)     = ( hi_max(jl-1) + hi_max(jl) ) / 2.0 
    128          zgfactorn(1) = zgfactorn(1) + exp(-(zhin(1)-hginn_u)*(zhin(1)-hginn_u)/2.0) 
    129          zhis (1)     = ( hi_max(jl-1) + hi_max(jl) ) / 2.0 
    130          zgfactors(1) = zgfactors(1) + exp(-(zhis(1)-hgins_u)*(zhis(1)-hgins_u)/2.0) 
     111         zhin (1)     = ( hi_max(jl-1) + hi_max(jl) ) * 0.5_wp 
     112         zgfactorn(1) = zgfactorn(1) + exp(-(zhin(1)-hginn_u)*(zhin(1)-hginn_u) * 0.5_wp ) 
     113         zhis (1)     = ( hi_max(jl-1) + hi_max(jl) ) * 0.5_wp 
     114         zgfactors(1) = zgfactors(1) + exp(-(zhis(1)-hgins_u)*(zhis(1)-hgins_u) * 0.5_wp ) 
    131115      END DO ! jl 
    132116      zgfactorn(1) = aginn_u / zgfactorn(1) 
    … …  
    135119      ! ------------- 
    136120      ! new distribution, polynom of second order, conserving area and volume 
    137       zh1 = 0.0 
    138       zh2 = 0.0 
    139       zh3 = 0.0 
     121      zh1 = 0._wp 
     122      zh2 = 0._wp 
     123      zh3 = 0._wp 
    140124      DO jl = 1, jpl 
    141          zh = ( hi_max(jl-1) + hi_max(jl) ) / 2.0 
     125         zh = ( hi_max(jl-1) + hi_max(jl) ) * 0.5_wp 
    142126         zh1 = zh1 + zh 
    143          zh2 = zh2 + zh*zh 
    144          zh3 = zh3 + zh*zh*zh 
     127         zh2 = zh2 + zh * zh 
     128         zh3 = zh3 + zh * zh * zh 
    145129      END DO 
    146130      IF(lwp) WRITE(numout,*) ' zh1 : ', zh1 
    … …  
    148132      IF(lwp) WRITE(numout,*) ' zh3 : ', zh3 
    149133 
    150       zvol = aginn_u*hginn_u 
     134      zvol = aginn_u * hginn_u 
    151135      zare = aginn_u 
    152       IF ( jpl .GE. 2 ) THEN 
     136      IF( jpl >= 2 ) THEN 
    153137         zbn = ( zvol*zh2 - zare*zh3 ) / ( zh2*zh2 - zh1*zh3) 
    154138         zan = ( zare - zbn*zh1 ) / zh2 
    … …  
    160144      IF(lwp) WRITE(numout,*) ' zan : ', zan  
    161145 
    162       zvol = agins_u*hgins_u 
     146      zvol = agins_u * hgins_u 
    163147      zare = agins_u 
    164       IF ( jpl .GE. 2 ) THEN 
     148      IF( jpl >= 2 ) THEN 
    165149         zbs = ( zvol*zh2 - zare*zh3 ) / ( zh2*zh2 - zh1*zh3) 
    166150         zas = ( zare - zbs*zh1 ) / zh2 
    … …  
    205189            !--- Northern hemisphere 
    206190            !---------------------------------------------------------------- 
    207             IF( fcor(ji,jj) >= 0.e0 ) THEN     
     191            IF( fcor(ji,jj) >= 0._wp ) THEN     
    208192 
    209193               !----------------------- 
    … …  
    453437            ENDIF ! on fcor 
    454438 
    455          ENDDO 
    456       ENDDO 
     439         END DO 
     440      END DO 
    457441 
    458442      !-------------------------------------------------------------------- 
    … …  
    494478 
    495479      DO jl = 1, jpl 
    496  
    497480         CALL lbc_lnk( a_i(:,:,jl)  , 'T', 1. ) 
    498481         CALL lbc_lnk( v_i(:,:,jl)  , 'T', 1. ) 
    … …  
    500483         CALL lbc_lnk( smv_i(:,:,jl), 'T', 1. ) 
    501484         CALL lbc_lnk( oa_i(:,:,jl) , 'T', 1. ) 
    502  
     485         ! 
    503486         CALL lbc_lnk( ht_i(:,:,jl) , 'T', 1. ) 
    504487         CALL lbc_lnk( ht_s(:,:,jl) , 'T', 1. ) 
    … …  
    514497            CALL lbc_lnk(e_i(:,:,jk,jl), 'T', 1. ) 
    515498         END DO 
    516  
    517          a_i (:,:,jl) = tms(:,:) * a_i(:,:,jl) 
    518  
     499         ! 
     500         a_i(:,:,jl) = tms(:,:) * a_i(:,:,jl) 
    519501      END DO 
    520502 
    521503      CALL lbc_lnk( at_i , 'T', 1. ) 
    522504      at_i(:,:) = tms(:,:) * at_i(:,:)                       ! put 0 over land 
    523  
     505      ! 
    524506      CALL lbc_lnk( fsbbq  , 'T', 1. ) 
    525  
     507      ! 
    526508   END SUBROUTINE lim_istate 
     509 
    527510 
    528511   SUBROUTINE lim_istate_init 
    … …  
    532515      !! ** Purpose : Definition of initial state of the ice  
    533516      !! 
    534       !! ** Method : Read the namiceini namelist and check the parameter  
    535       !!       values called at the first timestep (nit000) 
     517      !! ** Method :   Read the namiceini namelist and check the parameter  
     518      !!             values called at the first timestep (nit000) 
    536519      !! 
    537       !! ** input :  
    538       !!        Namelist namiceini 
    539       !! 
    540       !! history : 
    541       !!  8.5  ! 03-08 (C. Ethe) original code 
     520      !! ** input  :   namelist namiceini 
    542521      !!----------------------------------------------------------------------------- 
    543       NAMELIST/namiceini/ ttest, hninn, hginn_u, aginn_u, hginn_d, aginn_d, hnins, & 
    544          hgins_u, agins_u, hgins_d, agins_d, sinn, sins 
     522      NAMELIST/namiceini/ ttest, hninn, hginn_u, aginn_u, hginn_d, aginn_d, hnins,   & 
     523         &                hgins_u, agins_u, hgins_d, agins_d, sinn, sins 
    545524      !!----------------------------------------------------------------------------- 
    546  
    547       ! Define the initial parameters 
    548       ! ------------------------- 
    549  
    550       ! Read Namelist namiceini  
    551       REWIND ( numnam_ice ) 
     525      ! 
     526      REWIND ( numnam_ice )               ! Read Namelist namiceini  
    552527      READ   ( numnam_ice , namiceini ) 
    553       IF(lwp) THEN 
     528      ! 
     529      IF(lwp) THEN                        ! control print 
    554530         WRITE(numout,*) 
    555531         WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation ' 
    … …  
    569545         WRITE(numout,*) '   initial  ice salinity       in the south     sins       = ', sins 
    570546      ENDIF 
    571  
     547      ! 
    572548   END SUBROUTINE lim_istate_init 
    573549 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limitd_me.F90

    r2477 r2528  
    8686 
    8787   !!---------------------------------------------------------------------- 
    88    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008) 
     88   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    8989   !! $Id$ 
    90    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     90   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    9191   !!---------------------------------------------------------------------- 
    9292 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limitd_th.F90

    r2477 r2528  
    11MODULE limitd_th 
    2 #if defined key_lim3 
    3    !!---------------------------------------------------------------------- 
    4    !!   'key_lim3' :                                   LIM3 sea-ice model 
    5    !!---------------------------------------------------------------------- 
    62   !!====================================================================== 
    73   !!                       ***  MODULE limitd_th *** 
    … …  
    95   !!                   computation of changes in g(h)       
    106   !!====================================================================== 
    11  
     7#if defined key_lim3 
    128   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
     9   !!   'key_lim3' :                                   LIM3 sea-ice model 
     10   !!---------------------------------------------------------------------- 
     11   !!---------------------------------------------------------------------- 
    1412   USE dom_ice 
    1513   USE par_oce          ! ocean parameters 
    … …  
    1715   USE phycst           ! physical constants (ocean directory)  
    1816   USE thd_ice 
    19    USE in_out_manager 
    2017   USE ice 
    2118   USE par_ice 
    … …  
    2421   USE limcons 
    2522   USE prtctl           ! Print control 
     23   USE in_out_manager 
    2624   USE lib_mpp  
    2725 
    … …  
    2927   PRIVATE 
    3028 
    31    !! * Routine accessibility 
    3229   PUBLIC lim_itd_th        ! called by ice_stp 
    3330   PUBLIC lim_itd_th_rem 
    … …  
    3633   PUBLIC lim_itd_shiftice 
    3734 
    38    !! * Module variables 
    3935   REAL(wp)  ::           &  ! constant values 
    4036      epsi20 = 1e-20   ,  & 
    … …  
    4440 
    4541   !!---------------------------------------------------------------------- 
    46    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008) 
     42   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4743   !! $Id$ 
    48    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     44   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    4945   !!---------------------------------------------------------------------- 
    5046 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limmsh.F90

    r1923 r2528  
    2424 
    2525   !!---------------------------------------------------------------------- 
    26    !! NEMO/LIM 3.2,  UCL-ASTR-LOCEAN-IPSL (2009)  
     26   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    2727   !! $Id$ 
    28    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     28   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    2929   !!---------------------------------------------------------------------- 
    3030 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limrhg.F90

    r2477 r2528  
    77   !!            3.0  !  2008-03  (M. Vancoppenolle) LIM3 
    88   !!             -   !  2008-11  (M. Vancoppenolle, S. Bouillon, Y. Aksenov) add surface tilt in ice rheolohy  
     9   !!            3.3  !  2009-05  (G.Garric) addition of the lim2_evp cas 
    910   !!---------------------------------------------------------------------- 
    10 #if defined key_lim3 
     11#if defined key_lim3 || (  defined key_lim2 && ! defined key_lim2_vp ) 
    1112   !!---------------------------------------------------------------------- 
    12    !!   'key_lim3'                                      LIM3 sea-ice model 
     13   !!   'key_lim3'               OR                     LIM-3 sea-ice model 
     14   !!   'key_lim2' AND NOT 'key_lim2_vp'             VP LIM-2 sea-ice model 
    1315   !!---------------------------------------------------------------------- 
    1416   !!   lim_rhg   : computes ice velocities 
    1517   !!---------------------------------------------------------------------- 
    16    !! * Modules used 
    17    USE phycst 
    18    USE par_oce 
    19    USE dom_oce 
    20    USE dom_ice 
    21    USE sbc_oce         ! Surface boundary condition: ocean fields 
    22    USE sbc_ice         ! Surface boundary condition: ice fields 
    23    USE ice 
    24    USE lbclnk 
    25    USE lib_mpp 
    26    USE in_out_manager  ! I/O manager 
    27    USE limitd_me 
    28    USE prtctl          ! Print control 
    29  
     18   USE phycst           ! Physical constant 
     19   USE par_oce          ! Ocean parameters 
     20   USE dom_oce          ! Ocean domain 
     21   USE sbc_oce          ! Surface boundary condition: ocean fields 
     22   USE sbc_ice          ! Surface boundary condition: ice fields 
     23   USE lbclnk           ! Lateral Boundary Condition / MPP link 
     24   USE lib_mpp          ! MPP library 
     25   USE in_out_manager   ! I/O manager 
     26   USE prtctl           ! Print control 
     27#if defined key_lim3 
     28   USE ice              ! LIM-3: ice variables 
     29   USE dom_ice          ! LIM-3: ice domain 
     30   USE limitd_me        ! LIM-3:  
     31#else 
     32   USE ice_2            ! LIM2: ice variables 
     33   USE dom_ice_2        ! LIM2: ice domain 
     34#endif 
    3035 
    3136   IMPLICIT NONE 
    3237   PRIVATE 
    3338 
    34    !! * Routine accessibility 
    35    PUBLIC lim_rhg  ! routine called by lim_dyn 
    36  
     39   PUBLIC   lim_rhg   ! routine called by lim_dyn (or lim_dyn_2) 
     40 
     41   REAL(wp) ::   rzero   = 0._wp   ! constant values 
     42   REAL(wp) ::   rone    = 1._wp   ! constant values 
     43       
    3744   !! * Substitutions 
    3845#  include "vectopt_loop_substitute.h90" 
    39  
    40    !! * Module variables 
    41    REAL(wp)  ::           &  ! constant values 
    42       rzero   = 0.e0   ,  & 
    43       rone    = 1.e0 
    4446   !!---------------------------------------------------------------------- 
    45    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008)  
     47   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4648   !! $Id$ 
    47    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     49   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    4850   !!---------------------------------------------------------------------- 
    49  
    5051CONTAINS 
    5152 
    5253   SUBROUTINE lim_rhg( k_j1, k_jpj ) 
    53  
    5454      !!------------------------------------------------------------------- 
    5555      !!                 ***  SUBROUTINE lim_rhg  *** 
    … …  
    9898      !!                 e.g. in the Canadian Archipelago 
    9999      !! 
    100       !! ** References : Hunke and Dukowicz, JPO97 
    101       !!                 Bouillon et al., 08, in prep (update this when 
    102       !!                 published) 
    103       !!                 Vancoppenolle et al., OM08 
     100      !! References : Hunke and Dukowicz, JPO97 
     101      !!              Bouillon et al., Ocean Modelling 2009 
     102      !!              Vancoppenolle et al., Ocean Modelling 2008 
     103      !!------------------------------------------------------------------- 
     104      INTEGER, INTENT(in) ::   k_j1    ! southern j-index for ice computation 
     105      INTEGER, INTENT(in) ::   k_jpj   ! northern j-index for ice computation 
    104106      !! 
    105       !!------------------------------------------------------------------- 
    106       ! * Arguments 
    107       ! 
    108       INTEGER, INTENT(in) :: & 
    109          k_j1 ,                      & !: southern j-index for ice computation 
    110          k_jpj                         !: northern j-index for ice computation 
    111  
    112       ! * Local variables 
    113       INTEGER ::   ji, jj              !: dummy loop indices 
    114  
    115       INTEGER  :: & 
    116          jter                          !: temporary integers 
    117  
     107      INTEGER ::   ji, jj   ! dummy loop indices 
     108      INTEGER ::   jter     ! local integers 
    118109      CHARACTER (len=50) ::   charout 
    119  
    120       REAL(wp) :: & 
    121          zt11, zt12, zt21, zt22,     & !: temporary scalars 
    122          ztagnx, ztagny,             & !: wind stress on U/V points                        
    123          delta                         ! 
    124  
    125       REAL(wp) :: & 
    126          za,                         & !: 
    127          zstms,                      & !: temporary scalar for ice strength 
    128          zsang,                      & !: temporary scalar for coriolis term 
    129          zmask                         !: mask for the computation of ice mass 
     110      REAL(wp) ::   zt11, zt12, zt21, zt22, ztagnx, ztagny, delta                         ! 
     111      REAL(wp) ::   za, zstms, zsang, zmask   ! local scalars 
    130112 
    131113      REAL(wp),DIMENSION(jpi,jpj) :: & 
    … …  
    145127 
    146128      REAL(wp) :: & 
    147          dtevp,                      & !: time step for subcycling 
    148          dtotel,                     & !: 
    149          ecc2,                       & !: square of yield ellipse eccenticity 
    150          z0,                         & !: temporary scalar 
    151          zr,                         & !: temporary scalar 
    152          zcca, zccb,                 & !: temporary scalars 
    153          zu_ice2,                    & !:  
    154          zv_ice1,                    & !: 
    155          zddc, zdtc,                 & !: temporary array for delta on corners 
    156          zdst,                       & !: temporary array for delta on centre 
    157          zdsshx, zdsshy,             & !: term for the gradient of ocean surface 
    158          sigma1, sigma2                !: internal ice stress 
     129         dtevp,                      & ! time step for subcycling 
     130         dtotel,                     & ! 
     131         ecc2,                       & ! square of yield ellipse eccenticity 
     132         z0,                         & ! temporary scalar 
     133         zr,                         & ! temporary scalar 
     134         zcca, zccb,                 & ! temporary scalars 
     135         zu_ice2,                    & !  
     136         zv_ice1,                    & ! 
     137         zddc, zdtc,                 & ! temporary array for delta on corners 
     138         zdst,                       & ! temporary array for delta on centre 
     139         zdsshx, zdsshy,             & ! term for the gradient of ocean surface 
     140         sigma1, sigma2                ! internal ice stress 
     141 
     142      REAL(wp),DIMENSION(jpi,jpj) ::   zf1, zf2   ! arrays for internal stresses 
    159143 
    160144      REAL(wp),DIMENSION(jpi,jpj) :: & 
    161          zf1, zf2                      !: arrays for internal stresses 
    162  
    163       REAL(wp),DIMENSION(jpi,jpj) :: & 
    164          zdd, zdt,                   & !: Divergence and tension at centre of grid cells 
    165          zds,                        & !: Shear on northeast corner of grid cells 
    166          deltat,                     & !: Delta at centre of grid cells 
    167          deltac,                     & !: Delta on corners 
    168          zs1, zs2,                   & !: Diagonal stress tensor components zs1 and zs2  
    169          zs12                          !: Non-diagonal stress tensor component zs12 
     145         zdd, zdt,                   & ! Divergence and tension at centre of grid cells 
     146         zds,                        & ! Shear on northeast corner of grid cells 
     147         deltat,                     & ! Delta at centre of grid cells 
     148         deltac,                     & ! Delta on corners 
     149         zs1, zs2,                   & ! Diagonal stress tensor components zs1 and zs2  
     150         zs12                          ! Non-diagonal stress tensor component zs12 
    170151 
    171152      REAL(wp) :: & 
    172          zresm            ,          & !: Maximal error on ice velocity 
    173          zindb            ,          & !: ice (1) or not (0)       
    174          zdummy                        !: dummy argument 
    175  
    176       REAL(wp),DIMENSION(jpi,jpj) :: & 
    177          zu_ice           ,          & !: Ice velocity on previous time step 
    178          zv_ice           ,          & 
    179          zresr                         !: Local error on velocity 
    180  
     153         zresm            ,          & ! Maximal error on ice velocity 
     154         zindb            ,          & ! ice (1) or not (0)       
     155         zdummy                        ! dummy argument 
     156 
     157      REAL(wp),DIMENSION(jpi,jpj) ::   zu_ice, zv_ice, zresr   ! Local error on velocity 
     158      !!------------------------------------------------------------------- 
     159#if  defined key_lim2 && ! defined key_lim2_vp 
     160# if defined key_agrif 
     161     USE ice_2, vt_s => hsnm 
     162     USE ice_2, vt_i => hicm 
     163# else 
     164     vt_s => hsnm 
     165     vt_i => hicm 
     166# endif 
     167     at_i(:,:) = 1. - frld(:,:) 
     168#endif 
    181169      ! 
    182170      !------------------------------------------------------------------------------! 
    … …  
    185173      ! 
    186174      ! Put every vector to 0 
    187       zpresh(:,:) = 0.0 ; zc1(:,:) = 0.0 
    188       zpreshc(:,:) = 0.0 
    189       u_ice2(:,:)  = 0.0 ; v_ice1(:,:)  = 0.0 
    190       zdd(:,:)     = 0.0 ; zdt(:,:)     = 0.0 ; zds(:,:)     = 0.0 
    191  
    192       ! Ice strength on T-points 
    193       CALL lim_itd_me_icestrength(ridge_scheme_swi) 
    194  
    195       ! Ice mass and temp variables 
    196 !CDIR NOVERRCHK 
    197       DO jj = k_j1 , k_jpj 
     175      zpresh (:,:) = 0._wp   ;   zc1   (:,:) = 0._wp 
     176      zpreshc(:,:) = 0._wp 
     177      u_ice2 (:,:) = 0._wp   ;   v_ice1(:,:) = 0._wp 
     178      zdd    (:,:) = 0._wp   ;   zdt   (:,:) = 0._wp   ;   zds(:,:) = 0._wp 
     179 
     180#if defined key_lim3 
     181      CALL lim_itd_me_icestrength( ridge_scheme_swi )      ! LIM-3: Ice strength on T-points 
     182#endif 
     183 
     184!CDIR NOVERRCHK 
     185      DO jj = k_j1 , k_jpj       ! Ice mass and temp variables 
    198186!CDIR NOVERRCHK 
    199187         DO ji = 1 , jpi 
    200188            zc1(ji,jj)    = tms(ji,jj) * ( rhosn * vt_s(ji,jj) + rhoic * vt_i(ji,jj) ) 
    201             zpresh(ji,jj) = tms(ji,jj) *  strength(ji,jj) / 2. 
     189#if defined key_lim3 
     190            zpresh(ji,jj) = tms(ji,jj) *  strength(ji,jj) * 0.5_wp 
     191#endif 
    202192            ! tmi = 1 where there is ice or on land 
    203             tmi(ji,jj)    = 1.0 - ( 1.0 - MAX( 0.0 , SIGN ( 1.0 , vt_i(ji,jj) - & 
    204                epsd ) ) ) * tms(ji,jj) 
     193            tmi(ji,jj)    = 1._wp - ( 1._wp - MAX( 0._wp , SIGN ( 1._wp , vt_i(ji,jj) - epsd ) ) ) * tms(ji,jj) 
    205194         END DO 
    206195      END DO 
    … …  
    247236         DO ji = fs_2, fs_jpim1 
    248237 
    249             zt11 = tms(ji,jj)*e1t(ji,jj) 
    250             zt12 = tms(ji+1,jj)*e1t(ji+1,jj) 
    251             zt21 = tms(ji,jj)*e2t(ji,jj) 
    252             zt22 = tms(ji,jj+1)*e2t(ji,jj+1) 
     238            zt11 = tms(ji  ,jj) * e1t(ji  ,jj) 
     239            zt12 = tms(ji+1,jj) * e1t(ji+1,jj) 
     240            zt21 = tms(ji,jj  ) * e2t(ji,jj  ) 
     241            zt22 = tms(ji,jj+1) * e2t(ji,jj+1) 
    253242 
    254243            ! Leads area. 
    255             zfrld1(ji,jj) = ( zt12 * ( 1.0 - at_i(ji,jj) ) + & 
    256                &                        zt11 * ( 1.0 - at_i(ji+1,jj) ) ) / ( zt11 + zt12 + epsd ) 
    257             zfrld2(ji,jj) = ( zt22 * ( 1.0 - at_i(ji,jj) ) + & 
    258                &                        zt21 * ( 1.0 - at_i(ji,jj+1) ) ) / ( zt21 + zt22 + epsd ) 
     244            zfrld1(ji,jj) = ( zt12 * ( 1.0 - at_i(ji,jj) ) + zt11 * ( 1.0 - at_i(ji+1,jj) ) ) / ( zt11 + zt12 + epsd ) 
     245            zfrld2(ji,jj) = ( zt22 * ( 1.0 - at_i(ji,jj) ) + zt21 * ( 1.0 - at_i(ji,jj+1) ) ) / ( zt21 + zt22 + epsd ) 
    259246 
    260247            ! Mass, coriolis coeff. and currents 
    261248            zmass1(ji,jj) = ( zt12*zc1(ji,jj) + zt11*zc1(ji+1,jj) ) / (zt11+zt12+epsd) 
    262249            zmass2(ji,jj) = ( zt22*zc1(ji,jj) + zt21*zc1(ji,jj+1) ) / (zt21+zt22+epsd) 
    263             zcorl1(ji,jj) = zmass1(ji,jj) * ( e1t(ji+1,jj)*fcor(ji,jj) + & 
    264                e1t(ji,jj)*fcor(ji+1,jj) ) & 
    265                / (e1t(ji,jj) + e1t(ji+1,jj) + epsd ) 
    266             zcorl2(ji,jj) = zmass2(ji,jj) * ( e2t(ji,jj+1)*fcor(ji,jj) + & 
    267                e2t(ji,jj)*fcor(ji,jj+1) ) & 
    268                / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 
     250            zcorl1(ji,jj) = zmass1(ji,jj) * ( e1t(ji+1,jj)*fcor(ji,jj) + e1t(ji,jj)*fcor(ji+1,jj) )   & 
     251               &                          / ( e1t(ji,jj) + e1t(ji+1,jj) + epsd ) 
     252            zcorl2(ji,jj) = zmass2(ji,jj) * ( e2t(ji,jj+1)*fcor(ji,jj) + e2t(ji,jj)*fcor(ji,jj+1) )   & 
     253               &                          / ( e2t(ji,jj+1) + e2t(ji,jj) + epsd ) 
    269254            ! 
    270255            u_oce1(ji,jj)  = u_oce(ji,jj) 
    … …  
    290275            ! include it later 
    291276 
    292             zdsshx =  (ssh_m(ji+1,jj) - ssh_m(ji,jj))/e1u(ji,jj) 
    293             zdsshy =  (ssh_m(ji,jj+1) - ssh_m(ji,jj))/e2v(ji,jj) 
     277            zdsshx =  ( ssh_m(ji+1,jj) - ssh_m(ji,jj) ) / e1u(ji,jj) 
     278            zdsshy =  ( ssh_m(ji,jj+1) - ssh_m(ji,jj) ) / e2v(ji,jj) 
    294279 
    295280            za1ct(ji,jj) = ztagnx - zmass1(ji,jj) * grav * zdsshx 
    … …  
    306291      ! Time step for subcycling 
    307292      dtevp  = rdt_ice / nevp 
    308       dtotel = dtevp / ( 2.0 * telast ) 
     293      dtotel = dtevp / ( 2._wp * telast ) 
    309294 
    310295      !-ecc2: square of yield ellipse eccenticrity (reminder: must become a namelist parameter) 
    311       ecc2 = ecc*ecc 
     296      ecc2 = ecc * ecc 
    312297 
    313298      !-Initialise stress tensor  
    314       zs1(:,:)  = stress1_i(:,:)  
    315       zs2(:,:)  = stress2_i(:,:) 
     299      zs1 (:,:) = stress1_i (:,:)  
     300      zs2 (:,:) = stress2_i (:,:) 
    316301      zs12(:,:) = stress12_i(:,:) 
    317302 
    318       !----------------------! 
     303      !                                               !----------------------! 
    319304      DO jter = 1 , nevp                              !    loop over jter    ! 
    320          !----------------------!         
     305         !                                            !----------------------!         
    321306         DO jj = k_j1, k_jpj-1 
    322             zu_ice(:,jj) = u_ice(:,jj) ! velocity at previous time step 
     307            zu_ice(:,jj) = u_ice(:,jj)    ! velocity at previous time step 
    323308            zv_ice(:,jj) = v_ice(:,jj) 
    324309         END DO 
    … …  
    385370            END DO 
    386371         END DO 
    387          CALL lbc_lnk( v_ice1(:,:), 'U', -1. ) 
    388          CALL lbc_lnk( u_ice2(:,:), 'V', -1. ) 
     372         CALL lbc_lnk( v_ice1, 'U', -1. )   ;   CALL lbc_lnk( u_ice2, 'V', -1. )      ! lateral boundary cond. 
    389373 
    390374!CDIR NOVERRCHK 
    … …  
    394378 
    395379               !- Calculate Delta at centre of grid cells 
    396                zdst       = (  e2u( ji  , jj   ) * v_ice1(ji,jj)            & 
    397                   &          - e2u( ji-1, jj   ) * v_ice1(ji-1,jj)          & 
    398                   &          + e1v( ji  , jj   ) * u_ice2(ji,jj)            & 
    399                   &          - e1v( ji  , jj-1 ) * u_ice2(ji,jj-1)          & 
    400                   &          )                                              & 
     380               zdst       = (  e2u(ji  , jj) * v_ice1(ji  ,jj)          & 
     381                  &          - e2u(ji-1, jj) * v_ice1(ji-1,jj)          & 
     382                  &          + e1v(ji, jj  ) * u_ice2(ji,jj  )          & 
     383                  &          - e1v(ji, jj-1) * u_ice2(ji,jj-1)          & 
     384                  &          )                                          & 
    401385                  &         / area(ji,jj) 
    402386 
    403                delta = SQRT( zdd(ji,jj)*zdd(ji,jj) +                        &  
    404                   &                       ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 )   
    405                deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 +                    & 
    406                   (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 
     387               delta = SQRT( zdd(ji,jj)*zdd(ji,jj) + ( zdt(ji,jj)*zdt(ji,jj) + zdst*zdst ) * usecc2 )   
     388               deltat(ji,jj) = MAX( SQRT(zdd(ji,jj)**2 + (zdt(ji,jj)**2 + zdst**2)*usecc2), creepl ) 
    407389 
    408390               !-Calculate stress tensor components zs1 and zs2  

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limrst.F90

    r2477 r2528  
    44   !! Ice restart :  write the ice restart file 
    55   !!====================================================================== 
     6   !! History:   -   ! 2005-04 (M. Vancoppenolle) Original code 
     7   !!           3.0  ! 2008-03 (C. Ethe) restart files in using IOM interface 
     8   !!---------------------------------------------------------------------- 
    69#if defined key_lim3 
    710   !!---------------------------------------------------------------------- 
    … …  
    1215   !!   lim_rst_read    : read  the restart file  
    1316   !!---------------------------------------------------------------------- 
    14    !! * Modules used 
    15    USE ice 
    16    USE par_ice 
    17    USE in_out_manager 
    18    USE dom_oce 
    19    USE sbc_oce         ! Surface boundary condition: ocean fields 
    20    USE sbc_ice         ! Surface boundary condition: ice fields 
    21    USE iom 
     17   USE ice              ! sea-ice variables 
     18   USE par_ice          ! sea-ice parameters 
     19   USE dom_oce          ! ocean domain 
     20   USE sbc_oce          ! Surface boundary condition: ocean fields 
     21   USE sbc_ice          ! Surface boundary condition: ice fields 
     22   USE in_out_manager   ! I/O manager 
     23   USE iom              ! I/O library 
    2224 
    2325   IMPLICIT NONE 
    2426   PRIVATE 
    2527 
    26    !! * Accessibility 
    27    PUBLIC lim_rst_opn    ! routine called by icestep.F90 
    28    PUBLIC lim_rst_write  ! routine called by icestep.F90 
    29    PUBLIC lim_rst_read   ! routine called by iceinit.F90 
     28   PUBLIC   lim_rst_opn    ! routine called by icestep.F90 
     29   PUBLIC   lim_rst_write  ! routine called by icestep.F90 
     30   PUBLIC   lim_rst_read   ! routine called by iceini.F90 
    3031 
    3132   LOGICAL, PUBLIC ::   lrst_ice         !: logical to control the ice restart write  
    … …  
    3334 
    3435   !!---------------------------------------------------------------------- 
    35    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008) 
     36   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3637   !! $Id$ 
    37    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    38    !!---------------------------------------------------------------------- 
    39  
     38   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     39   !!---------------------------------------------------------------------- 
    4040CONTAINS 
    4141 
    … …  
    4848      INTEGER, INTENT(in) ::   kt       ! number of iteration 
    4949      ! 
    50       CHARACTER(LEN=20)   ::   clkt     ! ocean time-step deine as a character 
     50      CHARACTER(LEN=20)   ::   clkt     ! ocean time-step define as a character 
    5151      CHARACTER(LEN=50)   ::   clname   ! ice output restart file name 
    5252      !!---------------------------------------------------------------------- 
    … …  
    5454      IF( kt == nit000 )   lrst_ice = .FALSE.   ! default definition 
    5555 
    56       ! to get better performances with NetCDF format: 
    57       ! we open and define the ice restart file one ice time step before writing the data (-> at nitrst - 2*nn_fsbc + 1) 
    58       ! except if we write ice restart files every ice time step or if an ice restart file was writen at nitend - 2*nn_fsbc + 1 
    59       IF( kt == nitrst - 2*nn_fsbc + 1 .OR. nstock == nn_fsbc .OR. ( kt == nitend - nn_fsbc + 1 .AND. .NOT. lrst_ice ) ) THEN 
     56      ! in order to get better performances with NetCDF format, we open and define the ice restart file  
     57      ! one ice time step before writing the data (-> at nitrst - 2*nn_fsbc + 1), except if we write ice  
     58      ! restart files every ice time step or if an ice restart file was writen at nitend - 2*nn_fsbc + 1 
     59      IF( kt == nitrst - 2*nn_fsbc + 1 .OR. nstock == nn_fsbc    & 
     60         &                             .OR. ( kt == nitend - nn_fsbc + 1 .AND. .NOT. lrst_ice ) ) THEN 
    6061         ! beware of the format used to write kt (default is i8.8, that should be large enough...) 
    6162         IF( nitrst > 99999999 ) THEN   ;   WRITE(clkt, *       ) nitrst 
    … …  
    7576            ENDIF 
    7677         ENDIF 
    77  
     78         ! 
    7879         CALL iom_open( clname, numriw, ldwrt = .TRUE., kiolib = jprstlib ) 
    7980         lrst_ice = .TRUE. 
    … …  
    8283   END SUBROUTINE lim_rst_opn 
    8384 
     85 
    8486   SUBROUTINE lim_rst_write( kt ) 
    8587      !!---------------------------------------------------------------------- 
    … …  
    8789      !! 
    8890      !! ** purpose  :   output of sea-ice variable in a netcdf file 
     91      !!---------------------------------------------------------------------- 
     92      INTEGER, INTENT(in) ::   kt     ! number of iteration 
    8993      !! 
    90       !!---------------------------------------------------------------------- 
    91       ! Arguments : 
    92       INTEGER, INTENT(in) ::   kt     ! number of iteration 
    93  
    94       ! Local variables : 
     94      INTEGER ::   ji, jj, jk ,jl   ! dummy loop indices 
     95      INTEGER ::   iter 
     96      CHARACTER(len=15) ::   znam 
     97      CHARACTER(len=1)  ::   zchar, zchar1 
    9598      REAL(wp), DIMENSION(jpi,jpj) :: z2d 
    96       INTEGER :: ji, jj, jk ,jl 
    97       INTEGER :: iter 
    98       CHARACTER(len=15) :: znam 
    99       CHARACTER(len=1)  :: zchar, zchar1 
    10099      !!---------------------------------------------------------------------- 
    101100 
    … …  
    111110      ! ------------------  
    112111      !                                                                        ! calendar control 
    113       CALL iom_rstput( iter, nitrst, numriw, 'nn_fsbc', REAL( nn_fsbc, wp) )      ! time-step  
    114       CALL iom_rstput( iter, nitrst, numriw, 'kt_ice' , REAL( iter , wp) )        ! date 
     112      CALL iom_rstput( iter, nitrst, numriw, 'nn_fsbc', REAL( nn_fsbc, wp ) )      ! time-step  
     113      CALL iom_rstput( iter, nitrst, numriw, 'kt_ice' , REAL( iter   , wp ) )      ! date 
    115114 
    116115      ! Prognostic variables  
    … …  
    288287      ENDIF 
    289288      ! 
    290  
    291289   END SUBROUTINE lim_rst_write 
     290 
    292291 
    293292   SUBROUTINE lim_rst_read 
    … …  
    297296      !! ** purpose  :   read of sea-ice variable restart in a netcdf file 
    298297      !!---------------------------------------------------------------------- 
    299       ! Local variables 
    300298      INTEGER :: ji, jj, jk, jl, indx 
    301299      REAL(wp) ::   zfice, ziter 
    302       REAL(wp) :: & !parameters for the salinity profile 
    303          zs_inf, z_slope_s, zsmax, zsmin, zalpha, zindb 
    304       REAL(wp), DIMENSION(nlay_i) :: zs_zero  
    305       REAL(wp), DIMENSION(jpi,jpj) :: z2d 
    306       CHARACTER(len=15) :: znam 
    307       CHARACTER(len=1)  :: zchar, zchar1 
    308       INTEGER           :: jlibalt = jprstlib 
    309       LOGICAL           :: llok 
     300      REAL(wp) ::   zs_inf, z_slope_s, zsmax, zsmin, zalpha, zindb   ! local scalars used for the salinity profile 
     301      REAL(wp), DIMENSION(nlay_i)  ::   zs_zero  
     302      REAL(wp), DIMENSION(jpi,jpj) ::   z2d 
     303      CHARACTER(len=15) ::   znam 
     304      CHARACTER(len=1)  ::   zchar, zchar1 
     305      INTEGER           ::   jlibalt = jprstlib 
     306      LOGICAL           ::   llok 
    310307      !!---------------------------------------------------------------------- 
    311308 
    … …  
    384381      END DO 
    385382 
    386       ! Salinity profile 
    387       !----------------- 
    388       IF( num_sal == 2 ) THEN 
    389          !     CALL lim_var_salprof 
     383      IF( num_sal == 2 ) THEN      ! Salinity profile 
    390384         DO jl = 1, jpl 
    391385            DO jk = 1, nlay_i 
    … …  
    393387                  DO ji = 1, jpi 
    394388                     zs_inf        = sm_i(ji,jj,jl) 
    395                      z_slope_s     = 2.0*sm_i(ji,jj,jl)/MAX(0.01,ht_i(ji,jj,jl)) 
     389                     z_slope_s     = 2._wp * sm_i(ji,jj,jl) / MAX( 0.01_wp , ht_i(ji,jj,jl) ) 
    396390                     !- slope of the salinity profile 
    397                      zs_zero(jk)   = z_slope_s * ( FLOAT(jk)-1.0/2.0 ) * & 
    398                         ht_i(ji,jj,jl) / FLOAT(nlay_i) 
    399                      zsmax = 4.5 
    400                      zsmin = 3.5 
     391                     zs_zero(jk)   = z_slope_s * ( REAL(jk,wp) - 0.5_wp ) * ht_i(ji,jj,jl) / REAL(nlay_i,wp) 
     392                     zsmax = 4.5_wp 
     393                     zsmin = 3.5_wp 
    401394                     IF( sm_i(ji,jj,jl) .LT. zsmin ) THEN 
    402                         zalpha = 1.0 
     395                        zalpha = 1._wp 
    403396                     ELSEIF( sm_i(ji,jj,jl) .LT.zsmax ) THEN 
    404                         zalpha = sm_i(ji,jj,jl) / (zsmin-zsmax) + zsmax / (zsmax-zsmin) 
     397                        zalpha = sm_i(ji,jj,jl) / ( zsmin - zsmax ) + zsmax / ( zsmax - zsmin ) 
    405398                     ELSE 
    406                         zalpha = 0.0 
     399                        zalpha = 0._wp 
    407400                     ENDIF 
    408                      s_i(ji,jj,jk,jl) = zalpha*zs_zero(jk) + ( 1.0 - zalpha )*zs_inf 
     401                     s_i(ji,jj,jk,jl) = zalpha * zs_zero(jk) + ( 1._wp - zalpha ) * zs_inf 
    409402                  END DO 
    410403               END DO 
    … …  
    558551         END DO 
    559552      END DO 
    560  
     553      ! 
    561554      CALL iom_close( numrir ) 
    562  
     555      ! 
    563556   END SUBROUTINE lim_rst_read 
    564  
    565557 
    566558#else 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limsbc.F90

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

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limtab.F90

    r1156 r2528  
    2222 
    2323   !!---------------------------------------------------------------------- 
    24    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
     24   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    2525   !! $Id$ 
    26    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     26   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    2727   !!---------------------------------------------------------------------- 
    2828CONTAINS 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90

    r2477 r2528  
    77   !!            2.0  ! 2002-07 (C. Ethe, G. Madec)  LIM-2 (F90 rewriting) 
    88   !!            3.0  ! 2005-11 (M. Vancoppenolle)  LIM-3 : Multi-layer thermodynamics + salinity variations 
    9    !!             -   ! 2007-04 (M. Vancoppenolle) add lim_thd_glohec and lim_thd_con_dif 
     9   !!             -   ! 2007-04 (M. Vancoppenolle) add lim_thd_glohec, lim_thd_con_dh and lim_thd_con_dif 
    1010   !!            3.2  ! 2009-07 (M. Vancoppenolle, Y. Aksenov, G. Madec) bug correction in rdmsnif 
     11   !!            3.3  ! 2010-11 (G. Madec) corrected snow melting heat (due to factor betas) 
    1112   !!---------------------------------------------------------------------- 
    1213#if defined key_lim3 
    … …  
    1415   !!   'key_lim3'                                      LIM3 sea-ice model 
    1516   !!---------------------------------------------------------------------- 
    16    !!   lim_thd      : thermodynamic of sea ice 
     17   !!   lim_thd        : thermodynamic of sea ice 
     18   !!   lim_thd_init   : initialisation of sea-ice thermodynamic 
    1719   !!---------------------------------------------------------------------- 
    1820   USE phycst          ! physical constants 
    1921   USE dom_oce         ! ocean space and time domain variables 
    20    USE ice             ! LIM sea-ice variables 
    21    USE par_ice 
     22   USE ice             ! LIM: sea-ice variables 
     23   USE par_ice         ! LIM: sea-ice parameters 
    2224   USE sbc_oce         ! Surface boundary condition: ocean fields 
    2325   USE sbc_ice         ! Surface boundary condition: ice fields 
    2426   USE thd_ice         ! LIM thermodynamic sea-ice variables 
    2527   USE dom_ice         ! LIM sea-ice domain 
    26    USE domvvl 
    27    USE limthd_dif 
    28    USE limthd_dh 
    29    USE limthd_sal 
    30    USE limthd_ent 
    31    USE limtab 
    32    USE limvar 
     28   USE domvvl          ! domain: variable volume level 
     29   USE limthd_dif      ! LIM: thermodynamics, vertical diffusion 
     30   USE limthd_dh       ! LIM: thermodynamics, ice and snow thickness variation 
     31   USE limthd_sal      ! LIM: thermodynamics, ice salinity 
     32   USE limthd_ent      ! LIM: thermodynamics, ice enthalpy redistribution 
     33   USE limtab          ! LIM: 1D <==> 2D transformation 
     34   USE limvar          ! LIM: sea-ice variables 
     35   USE lbclnk          ! lateral boundary condition - MPP links 
     36   USE lib_mpp         ! MPP library 
    3337   USE in_out_manager  ! I/O manager 
    3438   USE prtctl          ! Print control 
    35    USE lbclnk 
    36    USE lib_mpp 
    3739 
    3840   IMPLICIT NONE 
    3941   PRIVATE 
    4042 
    41    PUBLIC   lim_thd    ! called by lim_step 
    42  
    43    REAL(wp) ::   epsi20 = 1e-20   ! constant values 
    44    REAL(wp) ::   epsi16 = 1e-16   ! 
    45    REAL(wp) ::   epsi06 = 1e-06   ! 
    46    REAL(wp) ::   epsi04 = 1e-04   ! 
    47    REAL(wp) ::   zzero  = 0.e0    ! 
    48    REAL(wp) ::   zone   = 1.e0    ! 
     43   PUBLIC   lim_thd        ! called by limstp module 
     44   PUBLIC   lim_thd_init   ! called by iceini module 
     45 
     46   REAL(wp) ::   epsi20 = 1e-20_wp   ! constant values 
     47   REAL(wp) ::   epsi16 = 1e-16_wp   ! 
     48   REAL(wp) ::   epsi06 = 1e-06_wp   ! 
     49   REAL(wp) ::   epsi04 = 1e-04_wp   ! 
     50   REAL(wp) ::   zzero  = 0._wp      ! 
     51   REAL(wp) ::   zone   = 1._wp      ! 
    4952 
    5053   !! * Substitutions 
    … …  
    5255#  include "vectopt_loop_substitute.h90" 
    5356   !!---------------------------------------------------------------------- 
    54    !! NEMO/LIM  3.2 ,  UCL-LOCEAN-IPSL (2009)  
     57   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    5558   !! $Id$ 
    56    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     59   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    5760   !!---------------------------------------------------------------------- 
    58  
    5961CONTAINS 
    6062 
    … …  
    8486      REAL(wp) ::   zfric_umax = 2e-02    ! upper bound for the friction velocity 
    8587      REAL(wp) ::   zinda, zindb, zthsnice, zfric_u    ! temporary scalar 
    86       REAL(wp) ::   zfntlat, zpareff                   !    -         - 
     88      REAL(wp) ::   zfntlat, zpareff   !    -         - 
    8789      REAL(wp) ::   zeps, zareamin, zcoef 
    8890      REAL(wp), DIMENSION(jpi,jpj) ::   zqlbsbq   ! link with lead energy budget qldif 
    8991      !!------------------------------------------------------------------- 
    90  
    9192       
    9293      !------------------------------------------------------------------------------! 
    … …  
    99100      !-------------------- 
    100101      ! Change the units of heat content; from global units to J.m3 
    101  
    102102      DO jl = 1, jpl 
    103103         DO jk = 1, nlay_i 
    … …  
    166166            pfrld(ji,jj)   = 1.0 - at_i(ji,jj) 
    167167            zinda          = 1.0 - MAX( zzero , SIGN( zone , - ( 1.0 - pfrld(ji,jj) ) ) ) 
    168  
     168            ! 
    169169            !           !  solar irradiance transmission at the mixed layer bottom and used in the lead heat budget 
    170170            !           !  practically no "direct lateral ablation" 
    … …  
    181181            qdtcn(ji,jj)  = zindb * fdtcn(ji,jj) * (1.0 - at_i(ji,jj)) * rdt_ice 
    182182            !                        
    183  
    184             ! still need to be updated : fdtcn !!!! 
    185             !           !-- Lead heat budget (part 1, next one is in limthd_dh 
    186             !           !-- qldif -- (or qldif_1d in 1d routines) 
    187             qldif(ji,jj) =  tms(ji,jj) * rdt_ice  * (                             &  
    188                &   pfrld(ji,jj)        * (  qsr(ji,jj)                            &   ! solar heat  
    189                &                            + qns(ji,jj)                          &   ! non solar heat  
    190                &                            + fdtcn(ji,jj)                        &   ! turbulent ice-ocean heat  
    191                &                            + fsbbq(ji,jj) * ( 1.0 - zindb )  )   &   ! residual heat from previous step  
     183            !           !-- Lead heat budget, qldif (part 1, next one is in limthd_dh)  
     184            !           !   caution: exponent betas used as more snow can fallinto leads 
     185            qldif(ji,jj) =  tms(ji,jj) * rdt_ice  * (                             & 
     186               &   pfrld(ji,jj)        * (  qsr(ji,jj)                            &   ! solar heat 
     187               &                            + qns(ji,jj)                          &   ! non solar heat 
     188               &                            + fdtcn(ji,jj)                        &   ! turbulent ice-ocean heat 
     189               &                            + fsbbq(ji,jj) * ( 1.0 - zindb )  )   &   ! residual heat from previous step 
    192190               & - pfrld(ji,jj)**betas * sprecip(ji,jj) * lfus                    )   ! latent heat of sprecip melting 
    193              
     191            ! 
    194192            ! Positive heat budget is used for bottom ablation 
    195193            zfntlat        = 1.0 - MAX( zzero , SIGN( zone ,  - qldif(ji,jj) ) ) 
    … …  
    198196            != 0 if ice and positive heat budget and 1 if one of those two is false 
    199197            zqlbsbq(ji,jj) = qldif(ji,jj) * ( 1.0 - zpareff ) / MAX( at_i(ji,jj) * rdt_ice , epsi16 ) 
    200  
     198            ! 
    201199            ! Heat budget of the lead, energy transferred from ice to ocean 
    202200            qldif  (ji,jj) = zpareff * qldif(ji,jj) 
    203201            qdtcn  (ji,jj) = zpareff * qdtcn(ji,jj) 
    204  
     202            ! 
    205203            ! Energy needed to bring ocean surface layer until its freezing (qcmif, limflx) 
    206204            qcmif  (ji,jj) =  rau0 * rcp * fse3t(ji,jj,1) * ( t_bo(ji,jj) - (sst_m(ji,jj) + rt0) ) * ( 1. - zinda ) 
    207  
     205            ! 
    208206            ! oceanic heat flux (limthd_dh) 
    209207            fbif   (ji,jj) = zindb * (  fsbbq(ji,jj) / MAX( at_i(ji,jj) , epsi20 ) + fdtcn(ji,jj) ) 
    … …  
    223221         ENDIF 
    224222 
    225          zareamin = 1.0e-10 
     223         zareamin = 1.e-10 
    226224         nbpb = 0 
    227225         DO jj = 1, jpj 
    … …  
    359357            CALL tab_1d_2d( nbpb, rdvonif, npb, dvnbq_1d  (1:nbpb), jpi, jpj )  
    360358            CALL tab_1d_2d( nbpb, fseqv  , npb, fseqv_1d  (1:nbpb), jpi, jpj ) 
    361  
     359            ! 
    362360            IF( num_sal == 2 ) THEN 
    363361               CALL tab_1d_2d( nbpb, fsbri, npb, fsbri_1d(1:nbpb), jpi, jpj ) 
    364362               CALL tab_1d_2d( nbpb, fhbri, npb, fhbri_1d(1:nbpb), jpi, jpj ) 
    365363            ENDIF 
    366  
     364            ! 
    367365            !+++++ 
    368366            !temporary stuff for a dummy version 
    … …  
    375373            CALL tab_1d_2d( nbpb, qns_ice(:,:,jl), npb, qnsr_ice_1d(1:nbpb), jpi, jpj) 
    376374            !+++++ 
    377  
     375            ! 
    378376            IF( lk_mpp )   CALL mpp_comm_free( ncomm_ice ) !RB necessary ?? 
    379377         ENDIF 
    … …  
    388386      ! 5.1) Ice heat content               
    389387      !------------------------ 
    390  
    391388      ! Enthalpies are global variables we have to readjust the units 
    392389      zcoef = 1.e0 / ( unit_fac * REAL(nlay_i) ) 
    … …  
    401398      ! 5.2) Snow heat content               
    402399      !------------------------ 
    403  
    404400      ! Enthalpies are global variables we have to readjust the units 
    405401      zcoef = 1.e0 / ( unit_fac * REAL(nlay_s) ) 
    … …  
    424420      IF( con_i )   fbif(:,:) = fbif(:,:) + zqlbsbq(:,:) 
    425421 
    426       IF(ln_ctl) THEN   ! Control print 
     422      IF(ln_ctl) THEN            ! Control print 
    427423         CALL prt_ctl_info(' ') 
    428424         CALL prt_ctl_info(' - Cell values : ') 
    … …  
    454450            END DO 
    455451         END DO 
    456  
    457452      ENDIF 
    458  
     453      ! 
    459454   END SUBROUTINE lim_thd 
    460455 
    … …  
    524519      END DO 
    525520      ! layer by layer 
    526       dq_i_layer(:,:)    = q_i_layer_fin(:,:) - q_i_layer_in(:,:) 
     521      dq_i_layer(:,:) = q_i_layer_fin(:,:) - q_i_layer_in(:,:) 
    527522 
    528523      !---------------------------------------- 
    529524      ! Atmospheric heat flux, ice heat budget 
    530525      !---------------------------------------- 
    531  
    532       DO ji = kideb, kiut 
    533          zji = MOD( npb(ji) - 1, jpi ) + 1 
    534          zjj = ( npb(ji) - 1 ) / jpi + 1 
    535  
    536          fatm(ji,jl) = qnsr_ice_1d(ji) + (1.0-i0(ji))*qsr_ice_1d(ji) 
    537  
    538          sum_fluxq(ji,jl) = fc_su(ji) - fc_bo_i(ji) + qsr_ice_1d(ji)*i0(ji) - fstroc(zji,zjj,jl) 
     526      DO ji = kideb, kiut 
     527         zji = MOD( npb(ji) - 1 , jpi ) + 1 
     528         zjj =    ( npb(ji) - 1 ) / jpi + 1 
     529         fatm     (ji,jl) = qnsr_ice_1d(ji) + ( 1._wp - i0(ji) ) * qsr_ice_1d(ji) 
     530         sum_fluxq(ji,jl) = fc_su(ji) - fc_bo_i(ji) + qsr_ice_1d(ji) * i0(ji) - fstroc(zji,zjj,jl) 
    539531      END DO 
    540532 
    … …  
    542534      ! Conservation error 
    543535      !-------------------- 
    544  
    545536      DO ji = kideb, kiut 
    546537         cons_error(ji,jl) = ABS( dq_i(ji,jl) / rdt_ice + sum_fluxq(ji,jl) ) 
    547538      END DO 
    548539 
    549       numce = 0 
     540      numce  = 0 
    550541      meance = 0.0 
    551542      DO ji = kideb, kiut 
    … …  
    554545            meance = meance + cons_error(ji,jl) 
    555546         ENDIF 
    556       ENDDO 
    557       IF (numce .GT. 0 ) meance = meance / numce 
     547      END DO 
     548      IF( numce .GT. 0 ) meance = meance / numce 
    558549 
    559550      WRITE(numout,*) ' Maximum tolerated conservation error : ', max_cons_err 
    … …  
    565556      ! Surface error due to imbalance between Fatm and Fcsu 
    566557      !------------------------------------------------------- 
    567       numce  = 0.0 
     558      numce  = 0 
    568559      meance = 0.0 
    569560 
    570561      DO ji = kideb, kiut 
    571562         surf_error(ji,jl) = ABS ( fatm(ji,jl) - fc_su(ji) ) 
    572          IF ( ( t_su_b(ji) .LT. rtt ) .AND. ( surf_error(ji,jl) .GT. & 
    573             max_surf_err ) ) THEN 
     563         IF( ( t_su_b(ji) .LT. rtt ) .AND. ( surf_error(ji,jl) .GT. max_surf_err ) ) THEN 
    574564            numce = numce + 1  
    575565            meance = meance + surf_error(ji,jl) 
    576566         ENDIF 
    577567      ENDDO 
    578       IF (numce .GT. 0 ) meance = meance / numce 
     568      IF( numce .GT. 0 ) meance = meance / numce 
    579569 
    580570      WRITE(numout,*) ' Maximum tolerated surface error : ', max_surf_err 
    … …  
    590580      ! Write ice state in case of big errors 
    591581      !--------------------------------------- 
    592  
    593582      DO ji = kideb, kiut 
    594583         IF ( ( ( t_su_b(ji) .LT. rtt ) .AND. ( surf_error(ji,jl) .GT. max_surf_err ) ) .OR. & 
    … …  
    596585            zji                 = MOD( npb(ji) - 1, jpi ) + 1 
    597586            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    598  
     587            ! 
    599588            WRITE(numout,*) ' alerte 1     ' 
    600589            WRITE(numout,*) ' Untolerated conservation / surface error after ' 
    … …  
    612601            !        WRITE(numout,*) ' qt_i_fin  : ', qt_i_fin(ji,jl) 
    613602            !        WRITE(numout,*) ' qt_s_fin  : ', qt_s_fin(ji,jl) 
    614             !        WRITE(numout,*) ' qt        : ', qt_i_fin(ji,jl) + & 
    615             !                                         qt_s_fin(ji,jl) 
     603            !        WRITE(numout,*) ' qt        : ', qt_i_fin(ji,jl) + qt_s_fin(ji,jl) 
    616604            WRITE(numout,*) ' ht_i       : ', ht_i_b(ji) 
    617605            WRITE(numout,*) ' ht_s       : ', ht_s_b(ji) 
    … …  
    640628            WRITE(numout,*) 
    641629            WRITE(numout,*) ' Layer by layer ... ' 
    642             WRITE(numout,*) ' dq_snow : ', ( qt_s_fin(ji,jl) - & 
    643                qt_s_in(ji,jl) )  &  
    644                / rdt_ice 
    645             WRITE(numout,*) ' dfc_snow  : ', fc_s(ji,1) -      & 
    646                fc_s(ji,0) 
     630            WRITE(numout,*) ' dq_snow : ', ( qt_s_fin(ji,jl) - qt_s_in(ji,jl) ) / rdt_ice 
     631            WRITE(numout,*) ' dfc_snow  : ', fc_s(ji,1) - fc_s(ji,0) 
    647632            DO jk = 1, nlay_i 
    648633               WRITE(numout,*) ' layer  : ', jk 
    649634               WRITE(numout,*) ' dq_ice : ', dq_i_layer(ji,jk) / rdt_ice   
    650635               WRITE(numout,*) ' radab  : ', radab(ji,jk) 
    651                WRITE(numout,*) ' dfc_i  : ', fc_i(ji,jk) -               & 
    652                   fc_i(ji,jk-1) 
    653                WRITE(numout,*) ' tot f  : ', fc_i(ji,jk) -               & 
    654                   fc_i(ji,jk-1) - radab(ji,jk) 
     636               WRITE(numout,*) ' dfc_i  : ', fc_i(ji,jk) - fc_i(ji,jk-1) 
     637               WRITE(numout,*) ' tot f  : ', fc_i(ji,jk) - fc_i(ji,jk-1) - radab(ji,jk) 
    655638            END DO 
    656639 
    … …  
    662645 
    663646 
    664    SUBROUTINE lim_thd_con_dh(kideb,kiut,jl) 
     647   SUBROUTINE lim_thd_con_dh( kideb, kiut, jl ) 
    665648      !!----------------------------------------------------------------------- 
    666649      !!                   ***  ROUTINE lim_thd_con_dh  ***  
    667650      !!                  
    668651      !! ** Purpose :   Test energy conservation after enthalpy redistr. 
    669       !! 
    670       !! history : 
    671       !!  9.9  ! 07-04 (M.Vancoppenolle) original code 
    672652      !!----------------------------------------------------------------------- 
    673653      INTEGER, INTENT(in) ::        & 
    … …  
    745725      ! Write ice state in case of big errors 
    746726      !--------------------------------------- 
    747  
    748727      DO ji = kideb, kiut 
    749728         IF ( cons_error(ji,jl) .GT. max_cons_err  ) THEN 
    750729            zji                 = MOD( npb(ji) - 1, jpi ) + 1 
    751730            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    752  
     731            ! 
    753732            WRITE(numout,*) ' alerte 1 - category : ', jl 
    754733            WRITE(numout,*) ' Untolerated conservation error after limthd_ent ' 
    … …  
    771750            WRITE(numout,*) ' qt_s_in    : ', qt_s_in(ji,jl) / rdt_ice 
    772751            WRITE(numout,*) ' qt_i_in    : ', qt_i_in(ji,jl) / rdt_ice 
    773             WRITE(numout,*) ' qt_in      : ', ( qt_i_in(ji,jl) + & 
    774                qt_s_in(ji,jl) ) / rdt_ice 
     752            WRITE(numout,*) ' qt_in      : ', ( qt_i_in(ji,jl) + qt_s_in(ji,jl) ) / rdt_ice 
    775753            WRITE(numout,*) ' qt_s_fin   : ', qt_s_fin(ji,jl) / rdt_ice 
    776754            WRITE(numout,*) ' qt_i_fin   : ', qt_i_fin(ji,jl) / rdt_ice 
    777             WRITE(numout,*) ' qt_fin     : ', ( qt_i_fin(ji,jl) + & 
    778                qt_s_fin(ji,jl) ) / rdt_ice 
     755            WRITE(numout,*) ' qt_fin     : ', ( qt_i_fin(ji,jl) + qt_s_fin(ji,jl) ) / rdt_ice 
    779756            WRITE(numout,*) ' * ' 
    780757            WRITE(numout,*) ' Ice variables --- : ' 
    … …  
    785762            WRITE(numout,*) ' dh_i_surf : ', dh_i_surf(ji) 
    786763            WRITE(numout,*) ' dh_i_bott : ', dh_i_bott(ji) 
    787  
    788764         ENDIF 
    789765         ! 
    … …  
    800776      !! 
    801777      !! ** Method  :   Formula (Bitz and Lipscomb, 1999) 
    802       !! 
    803778      !!------------------------------------------------------------------- 
    804779      INTEGER, INTENT(in) ::   kideb, kiut   ! bounds for the spatial loop 
    … …  
    827802 
    828803   SUBROUTINE lim_thd_init 
    829  
    830804      !!----------------------------------------------------------------------- 
    831805      !!                   ***  ROUTINE lim_thd_init ***  
    … …  
    838812      !! 
    839813      !! ** input   :   Namelist namicether 
    840        !!------------------------------------------------------------------- 
     814      !!------------------------------------------------------------------- 
    841815      NAMELIST/namicethd/ hmelt , hiccrit, fraz_swi, maxfrazb, vfrazb, Cfrazb,   & 
    842816         &                hicmin, hiclim, amax  ,                                & 
    … …  
    845819         &                kappa_i, nconv_i_thd, maxer_i_thd, thcon_i_swi 
    846820      !!------------------------------------------------------------------- 
    847  
     821      ! 
    848822      IF(lwp) THEN 
    849823         WRITE(numout,*) 
    … …  
    851825         WRITE(numout,*) '~~~~~~~' 
    852826      ENDIF 
    853  
     827      ! 
    854828      REWIND( numnam_ice )                  ! read Namelist numnam_ice 
    855829      READ  ( numnam_ice , namicethd ) 
    856        
     830      ! 
    857831      IF(lwp) THEN                          ! control print 
    858832         WRITE(numout,*) 
    … …  
    892866#else 
    893867   !!---------------------------------------------------------------------- 
    894    !!   Default option                               NO  LIM3 sea-ice model 
     868   !!   Default option         Dummy module          NO  LIM3 sea-ice model 
    895869   !!---------------------------------------------------------------------- 
    896 CONTAINS 
    897    SUBROUTINE lim_thd         ! Empty routine 
    898    END SUBROUTINE lim_thd 
    899    SUBROUTINE lim_thd_con_dif 
    900    END SUBROUTINE lim_thd_con_dif 
    901870#endif 
    902871 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_dh.F90

    r2477 r2528  
    1717   USE phycst           ! physical constants (OCE directory)  
    1818   USE sbc_oce          ! Surface boundary condition: ocean fields 
     19   USE ice 
     20   USE par_ice 
    1921   USE thd_ice 
    2022   USE in_out_manager 
    21    USE ice 
    22    USE par_ice 
    2323   USE lib_mpp 
    2424 
    … …  
    3535 
    3636   !!---------------------------------------------------------------------- 
    37    !! NEMO/LIM 3.2,  UCL-LOCEAN-IPSL (2009) 
     37   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3838   !! $Id$ 
    39    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     39   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    4040   !!---------------------------------------------------------------------- 
    4141 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_dif.F90

    r2477 r2528  
    11MODULE limthd_dif 
    2 #if defined key_lim3 
    3    !!---------------------------------------------------------------------- 
    4    !!   'key_lim3'                                      LIM3 sea-ice model 
    5    !!---------------------------------------------------------------------- 
    62   !!====================================================================== 
    73   !!                       ***  MODULE limthd_dif *** 
    … …  
    95   !!                   computation of surface and inner T   
    106   !!====================================================================== 
    11  
    127   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
     8#if defined key_lim3 
     9   !!---------------------------------------------------------------------- 
     10   !!   'key_lim3'                                      LIM3 sea-ice model 
     11   !!---------------------------------------------------------------------- 
    1412   USE par_oce          ! ocean parameters 
    1513   USE phycst           ! physical constants (ocean directory)  
    … …  
    2321   PRIVATE 
    2422 
    25    !! * Routine accessibility 
    26    PUBLIC lim_thd_dif        ! called by lim_thd 
    27  
    28    !! * Module variables 
     23   PUBLIC   lim_thd_dif   ! called by lim_thd 
     24 
    2925   REAL(wp)  ::           &  ! constant values 
    3026      epsi20 = 1e-20   ,  & 
    … …  
    3430 
    3531   !!---------------------------------------------------------------------- 
    36    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2005) 
     32   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3733   !! $Id$ 
    38    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     34   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    3935   !!---------------------------------------------------------------------- 
    40  
    4136CONTAINS 
    4237 
    … …  
    870865   END SUBROUTINE lim_thd_dif 
    871866#endif 
     867   !!====================================================================== 
    872868END MODULE limthd_dif 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_ent.F90

    r2477 r2528  
    11MODULE limthd_ent 
    2 #if defined key_lim3 
    3    !!---------------------------------------------------------------------- 
    4    !!   'key_lim3'                                      LIM3 sea-ice model 
    5    !!---------------------------------------------------------------------- 
    62   !!====================================================================== 
    73   !!                       ***  MODULE limthd_ent   *** 
    … …  
    106   !!                       after vertical growth/decay 
    117   !!====================================================================== 
     8#if defined key_lim3 
     9   !!---------------------------------------------------------------------- 
     10   !!   'key_lim3'                                      LIM3 sea-ice model 
     11   !!---------------------------------------------------------------------- 
    1212   !!   lim_thd_ent : ice redistribution of enthalpy 
    13  
    14    !! * Modules used 
     13   !!---------------------------------------------------------------------- 
    1514   USE par_oce          ! ocean parameters 
    1615   USE dom_oce 
    … …  
    2726   PRIVATE 
    2827 
    29    !! * Routine accessibility 
    30    PUBLIC  lim_thd_ent     ! called by lim_thd 
    31  
    32    !! * Module variables 
     28   PUBLIC   lim_thd_ent     ! called by lim_thd 
     29 
    3330   REAL(wp)  ::           &  ! constant values 
    3431      epsi20 = 1.e-20  ,  & 
    … …  
    3734      zone   = 1.e0    ,  & 
    3835      epsi10 = 1.0e-10 
    39  
    4036   !!---------------------------------------------------------------------- 
    41    !!   LIM 3.0,  UCL-LOCEAN-IPSL (2005) 
     37   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4238   !! $Id$ 
    43    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     39   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
    4440   !!---------------------------------------------------------------------- 
    45  
    4641CONTAINS 
    4742 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_lac.F90

    r2477 r2528  
    11MODULE limthd_lac 
    2    !!---------------------------------------------------------------------- 
    3    !!   'key_lim3'                                      LIM3 sea-ice model 
    4    !!---------------------------------------------------------------------- 
    52   !!====================================================================== 
    63   !!                       ***  MODULE limthd_lac   *** 
    … …  
    96#if defined key_lim3 
    107   !!---------------------------------------------------------------------- 
     8   !!   'key_lim3'                                      LIM3 sea-ice model 
     9   !!---------------------------------------------------------------------- 
    1110   !!   lim_lat_acr    : lateral accretion of ice 
    12    !! * Modules used 
     11   !!---------------------------------------------------------------------- 
    1312   USE par_oce          ! ocean parameters 
    1413   USE dom_oce 
    … …  
    4241 
    4342   !!---------------------------------------------------------------------- 
    44    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008)  
     43   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    4544   !! $Id$ 
    46    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     45   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    4746   !!---------------------------------------------------------------------- 
    4847 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limthd_sal.F90

    r2477 r2528  
    11MODULE limthd_sal 
    2    !!---------------------------------------------------------------------- 
    3    !!   'key_lim3'                                      LIM3 sea-ice model 
    4    !!---------------------------------------------------------------------- 
    52   !!====================================================================== 
    63   !!                       ***  MODULE limthd_sal *** 
    7    !!                   computation salinity variations in 
    8    !!                               the ice 
     4   !! LIM-3 sea-ice :  computation of salinity variations in the ice 
    95   !!====================================================================== 
     6   !! History :   -   ! 2003-05 (M. Vancoppenolle) UCL-ASTR first coding for LIM3-1D 
     7   !!            3.0  ! 2005-12 (M. Vancoppenolle) adapted to the 3-D version 
     8   !!--------------------------------------------------------------------- 
    109#if defined key_lim3 
    1110   !!---------------------------------------------------------------------- 
     11   !!   'key_lim3'                                      LIM-3 sea-ice model 
     12   !!---------------------------------------------------------------------- 
    1213   !!   lim_thd_sal : salinity variations in the ice 
    13    !! * Modules used 
     14   !!---------------------------------------------------------------------- 
    1415   USE par_oce          ! ocean parameters 
    1516   USE phycst           ! physical constants (ocean directory) 
    1617   USE sbc_oce          ! Surface boundary condition: ocean fields 
    17    USE thd_ice 
    18    USE ice 
    19    USE in_out_manager 
    20    USE limvar 
    21    USE par_ice 
     18   USE ice              ! LIM: sea-ice variables 
     19   USE par_ice          ! LIM: sea-ice parameters 
     20   USE thd_ice          ! LIM: sea-ice thermodynamics 
     21   USE limvar           ! LIM: sea-ice variables 
     22   USE in_out_manager   ! I/O manager 
    2223 
    2324   IMPLICIT NONE 
    2425   PRIVATE 
    2526 
    26    !! * Routine accessibility 
    27    PUBLIC lim_thd_sal        ! called by lim_thd 
    28    PUBLIC lim_thd_sal_init   ! called by lim_thd 
    29  
    30    !! * Module variables 
    31    REAL(wp)  ::           &  ! constant values 
    32       epsi20 = 1e-20   ,  & 
    33       epsi13 = 1e-13   ,  & 
    34       zzero  = 0.e0    ,  & 
    35       zone   = 1.e0 
    36  
    37    !!---------------------------------------------------------------------- 
    38    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
     27   PUBLIC   lim_thd_sal        ! called by limthd module 
     28   PUBLIC   lim_thd_sal_init   ! called by iceini module 
     29 
     30   !!---------------------------------------------------------------------- 
     31   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    3932   !! $Id$ 
    40    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    41    !!---------------------------------------------------------------------- 
    42  
     33   !! Software governed by the CeCILL licence     (NEMOGCM/NEMO_CeCILL.txt) 
     34   !!---------------------------------------------------------------------- 
    4335CONTAINS 
    4436 
    45    SUBROUTINE lim_thd_sal(kideb,kiut) 
     37   SUBROUTINE lim_thd_sal( kideb, kiut ) 
    4638      !!------------------------------------------------------------------- 
    47       !!                ***  ROUTINE lim_thd_sal  ***        
    48       !! ** Purpose : 
    49       !!        This routine computes new salinities in the ice 
     39      !!                ***  ROUTINE lim_thd_sal  ***     
     40      !!    
     41      !! ** Purpose :  computes new salinities in the ice 
    5042      !! 
    5143      !! ** Method  :  4 possibilities 
    … …  
    5446      !!               -> num_sal = 3 -> S = S(z)   [multiyear ice] 
    5547      !!               -> num_sal = 4 -> S = S(h)   [Cox and Weeks 74] 
    56       !!            
    57       !! ** Steps 
    58       !! 
    59       !! ** Arguments 
    60       !! 
    61       !! ** Inputs / Outputs 
    62       !! 
    63       !! ** External 
    64       !! 
    65       !! ** References 
    66       !! 
    67       !! ** History  :  
    68       !! 
    69       !! "Je ne suis pour l'instant qu'a 80% de ma condition, mais c'est  
    70       !!  les 30% qui restent qui seront les plus difficiles" 
    71       !!                                           E. Mpenza 
    72       !! 
    73       !!------------------------------------------------------------------- 
    74       !! History : 
    75       !!   ori  !  03-05 M. Vancoppenolle UCL-ASTR first coding for LIM-1D 
    76       !!   3.0  !  05-12 Routine rewritten for the 3-D version 
    7748      !!--------------------------------------------------------------------- 
    78       !! 
    79       !! * Local variables 
    80       INTEGER, INTENT(in) :: & 
    81          kideb, kiut             !: thickness category index 
    82  
    83       INTEGER ::             & 
    84          ji, jk ,            &   !: geographic and layer index  
    85          zji, zjj 
    86  
    87       REAL(wp) ::            & 
    88          zsold,              &   !: old salinity 
    89          zeps=1.0e-06   ,    &   !: very small 
    90          iflush         ,    &   !: flushing (1) or not (0) 
    91          iaccrbo        ,    &   !: bottom accretion (1) or not (0) 
    92          igravdr        ,    &   !: gravity drainage or not 
    93          isnowic        ,    &   !: gravity drainage or not 
    94          i_ice_switch   ,    &   !: ice thickness above a certain treshold or not 
    95          ztmelts        ,    &   !: freezing point of sea ice 
    96          zaaa           ,    &   !: dummy factor 
    97          zbbb           ,    &   !: dummy factor 
    98          zccc           ,    &   !: dummy factor 
    99          zdiscrim                !: dummy factor 
    100  
    101       REAL(wp), DIMENSION(jpij) ::          & 
    102          ze_init        ,    &   !initial total enthalpy 
    103          zhiold         ,    & 
    104          zsiold 
    105  
     49      INTEGER, INTENT(in) ::  kideb, kiut   ! thickness category index 
     50      ! 
     51      INTEGER  ::   ji, jk     ! dummy loop indices  
     52      INTEGER  ::   zji, zjj   ! local integers 
     53      REAL(wp) ::   zsold, zeps, iflush, iaccrbo, igravdr, isnowic, i_ice_switch,  ztmelts   ! local scalars 
     54      REAL(wp) ::   zaaa, zbbb, zccc, zdiscrim   ! local scalars 
     55      REAL(wp), DIMENSION(jpij) ::   ze_init, zhiold, zsiold   ! 1D workspace 
    10656      !!--------------------------------------------------------------------- 
    10757 
     58      zeps=1.0e-06_wp 
     59 
    10860      !------------------------------------------------------------------------------| 
    10961      ! 1) Constant salinity, constant in time                                       | 
    11062      !------------------------------------------------------------------------------| 
    11163 
    112       IF (num_sal.eq.1) THEN 
    113  
    114          !         WRITE(numout,*) 
    115          !         WRITE(numout,*) 'lim_thd_sal : Ice salinity computation module ', & 
    116          !         num_sal 
    117          !         WRITE(numout,*) '~~~~~~~~~~~~' 
    118  
     64      IF( num_sal == 1 ) THEN 
    11965         DO jk = 1, nlay_i 
    12066            DO ji = kideb, kiut 
    … …  
    12268            END DO ! ji 
    12369         END DO ! jk 
    124  
     70         ! 
    12571         DO ji = kideb, kiut 
    12672            sm_i_b(ji)      =  bulk_sal  
    12773         END DO ! ji 
    128  
    129       ENDIF ! num_sal .EQ. 1 
     74         ! 
     75      ENDIF 
    13076 
    13177      !------------------------------------------------------------------------------| 
    … …  
    174120            ! isnowic : 1 if snow ice formation 
    175121            i_ice_switch = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - ht_i_b(ji) + 1.0e-2 ) ) 
    176             isnowic      = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - dh_snowice(ji) ) ) * & 
    177                i_ice_switch 
     122            isnowic      = 1.0 - MAX ( 0.0 , SIGN ( 1.0 , - dh_snowice(ji) ) ) * i_ice_switch 
    178123 
    179124            !--------------------- 
    … …  
    182127 
    183128            ! drainage by gravity drainage 
    184             dsm_i_gd_1d(ji) = - igravdr *                                     &  
    185                MAX( sm_i_b(ji) - sal_G , 0.0 ) /             & 
    186                time_G * rdt_ice  
     129            dsm_i_gd_1d(ji) = - igravdr * MAX( sm_i_b(ji) - sal_G , 0._wp ) / time_G * rdt_ice  
    187130 
    188131            ! drainage by flushing   
    189             dsm_i_fl_1d(ji)  = - iflush *                                     & 
    190                MAX( sm_i_b(ji) - sal_F , 0.0 ) /             &  
    191                time_F * rdt_ice 
     132            dsm_i_fl_1d(ji)  = - iflush * MAX( sm_i_b(ji) - sal_F , 0._wp ) / time_F * rdt_ice 
    192133 
    193134            !----------------- 
    … …  
    199140            ! to conserve energy 
    200141            zsiold(ji) = sm_i_b(ji) 
    201             sm_i_b(ji) = sm_i_b(ji)                                           & 
    202                + dsm_i_fl_1d(ji) + dsm_i_gd_1d(ji)  !                 & 
     142            sm_i_b(ji) = sm_i_b(ji) + dsm_i_fl_1d(ji) + dsm_i_gd_1d(ji) 
    203143 
    204144            ! if no ice, salinity eq 0.1 
    205             i_ice_switch  = 1.0 - MAX ( 0.0, SIGN (1.0 , - ht_i_b(ji) ) ) 
    206             sm_i_b(ji)    = i_ice_switch*sm_i_b(ji) + s_i_min * ( 1.0 -       & 
    207                i_ice_switch ) 
     145            i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) ) 
     146            sm_i_b(ji)   = i_ice_switch*sm_i_b(ji) + s_i_min * ( 1._wp - i_ice_switch ) 
    208147         END DO ! ji 
    209148 
    210149         ! Salinity profile 
    211          CALL lim_var_salprof1d(kideb,kiut) 
     150         CALL lim_var_salprof1d( kideb, kiut ) 
    212151 
    213152         !---------------------------- 
    … …  
    217156         DO ji = kideb, kiut 
    218157            ! iflush  : 1 if summer  
    219             iflush       =  MAX( 0.0 , SIGN ( 1.0 , t_su_b(ji) - rtt ) )  
     158            iflush  =  MAX( 0._wp , SIGN ( 1._wp , t_su_b(ji) - rtt ) )  
    220159            ! igravdr : 1 if t_su lt t_bo 
    221             igravdr      =  MAX( 0.0 , SIGN ( 1.0 , t_bo_b(ji) - t_su_b(ji) ) )  
     160            igravdr =  MAX( 0._wp , SIGN ( 1._wp , t_bo_b(ji) - t_su_b(ji) ) )  
    222161            ! iaccrbo : 1 if bottom accretion 
    223             iaccrbo      =  MAX( 0.0 , SIGN ( 1.0 , dh_i_bott(ji) ) ) 
    224  
    225             fhbri_1d(ji) = 0.0 
     162            iaccrbo =  MAX( 0._wp , SIGN ( 1._wp , dh_i_bott(ji) ) ) 
     163            ! 
     164            fhbri_1d(ji) = 0._wp 
    226165         END DO ! ji 
    227166 
    … …  
    230169         !---------------------------- 
    231170         DO ji = kideb, kiut 
    232             i_ice_switch  = 1.0 - MAX ( 0.0, SIGN (1.0 , - ht_i_b(ji) ) ) 
    233             fsbri_1d(ji) = fsbri_1d(ji) -  & 
    234                i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji) *  & 
    235                ( MAX(dsm_i_gd_1d(ji) + dsm_i_fl_1d(ji), & 
    236                sm_i_b(ji) - zsiold(ji) ) ) / rdt_ice 
    237             IF ( num_sal .EQ. 4 ) fsbri_1d(ji) = 0.0 
    238  
     171            i_ice_switch = 1._wp - MAX ( 0._wp, SIGN( 1._wp , - ht_i_b(ji) ) ) 
     172            fsbri_1d(ji) = fsbri_1d(ji) - i_ice_switch * rhoic * a_i_b(ji) * ht_i_b(ji)         & 
     173               &         * ( MAX(dsm_i_gd_1d(ji) + dsm_i_fl_1d(ji), sm_i_b(ji) - zsiold(ji) ) ) / rdt_ice 
     174            IF( num_sal == 4 ) fsbri_1d(ji) = 0._wp 
    239175         END DO ! ji 
    240176 
    … …  
    244180         !-------------------- 
    245181         DO jk = 1, nlay_i 
    246  
    247182            DO ji = kideb, kiut 
    248  
    249183               ztmelts    =  -tmut*s_i_b(ji,jk) + rtt 
    250184               !Conversion q(S,T) -> T (second order equation) 
    251185               zaaa         =  cpic 
    252                zbbb         =  ( rcp - cpic ) * ( ztmelts - rtt ) + & 
    253                   q_i_b(ji,jk) / rhoic - lfus 
     186               zbbb         =  ( rcp - cpic ) * ( ztmelts - rtt ) + q_i_b(ji,jk) / rhoic - lfus 
    254187               zccc         =  lfus * ( ztmelts - rtt ) 
    255188               zdiscrim     =  SQRT( MAX(zbbb*zbbb - 4.0*zaaa*zccc,0.0) ) 
    256                t_i_b(ji,jk) =  rtt - ( zbbb + zdiscrim ) / &  
    257                   ( 2.0 *zaaa ) 
     189               t_i_b(ji,jk) =  rtt - ( zbbb + zdiscrim ) / ( 2.0 *zaaa ) 
    258190            END DO !ji 
    259  
    260191         END DO !jk 
    261  
     192         ! 
    262193      ENDIF ! num_sal .EQ. 2 
    263194 
    … …  
    266197      !------------------------------------------------------------------------------| 
    267198 
    268       IF ( num_sal .EQ. 3 ) THEN 
     199      IF( num_sal .EQ. 3 ) THEN 
    269200 
    270201         WRITE(numout,*) 
    … …  
    320251            zji                 = MOD( npb(ji) - 1, jpi ) + 1 
    321252            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    322             fseqv_1d(ji) = fseqv_1d(ji)              + &  
    323                ( sss_m(zji,zjj) - bulk_sal    ) * &  
    324                rhoic * a_i_b(ji) * & 
    325                MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 
     253            fseqv_1d(ji) = fseqv_1d(ji) + ( sss_m(zji,zjj) - bulk_sal    )               & 
     254               &                        * rhoic * a_i_b(ji) * MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 
    326255         END DO 
    327256      ELSE 
    … …  
    329258            zji                 = MOD( npb(ji) - 1, jpi ) + 1 
    330259            zjj                 = ( npb(ji) - 1 ) / jpi + 1 
    331             fseqv_1d(ji) = fseqv_1d(ji)              + &  
    332                ( sss_m(zji,zjj) - s_i_new(ji) ) * &  
    333                rhoic * a_i_b(ji) * & 
    334                MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 
     260            fseqv_1d(ji) = fseqv_1d(ji) + ( sss_m(zji,zjj) - s_i_new(ji) )               & 
     261               &                        * rhoic * a_i_b(ji) * MAX( dh_i_bott(ji) , 0.0 ) / rdt_ice 
    335262         END DO ! ji 
    336263      ENDIF 
    337  
    338       !-- End of salinity computations 
     264      ! 
    339265   END SUBROUTINE lim_thd_sal 
    340    !============================================================================== 
     266 
    341267 
    342268   SUBROUTINE lim_thd_sal_init 
    … …  
    346272      !! ** Purpose :   initialization of ice salinity parameters 
    347273      !! 
    348       !! ** Method  : Read the namicesal namelist and check the parameter 
    349       !!       values called at the first timestep (nit000) 
     274      !! ** Method  :   Read the namicesal namelist and check the parameter 
     275      !!              values called at the first timestep (nit000) 
    350276      !! 
    351277      !! ** input   :   Namelist namicesal 
    352       !! 
    353       !! history : 
    354       !!   3.0  !  July 2005 M. Vancoppenolle  Original code 
    355278      !!------------------------------------------------------------------- 
    356       NAMELIST/namicesal/ num_sal, bulk_sal, sal_G, time_G, sal_F, time_F, & 
    357          s_i_max, s_i_min, s_i_0, s_i_1 
     279      NAMELIST/namicesal/ num_sal, bulk_sal, sal_G, time_G, sal_F, time_F,   & 
     280         &                s_i_max, s_i_min, s_i_0, s_i_1 
    358281      !!------------------------------------------------------------------- 
    359  
    360       ! Read Namelist namicesal 
    361       REWIND ( numnam_ice ) 
    362       READ   ( numnam_ice  , namicesal ) 
    363       IF(lwp) THEN 
     282      ! 
     283      REWIND( numnam_ice )                   ! Read Namelist namicesal 
     284      READ  ( numnam_ice  , namicesal ) 
     285      ! 
     286      IF(lwp) THEN                           ! control print 
    364287         WRITE(numout,*) 
    365288         WRITE(numout,*) 'lim_thd_sal_init : Ice parameters for salinity ' 
    … …  
    376299         WRITE(numout,*) ' 2nd salinity for salinity profile  : ', s_i_1 
    377300      ENDIF 
    378  
     301      ! 
    379302   END SUBROUTINE lim_thd_sal_init 
    380303 
    381304#else 
    382305   !!---------------------------------------------------------------------- 
    383    !!   Default option         Empty Module                No sea-ice model 
    384    !!---------------------------------------------------------------------- 
    385 CONTAINS 
    386    SUBROUTINE lim_thd_sal        ! Empty routine 
    387    END SUBROUTINE lim_thd_sal 
     306   !!   Default option         Dummy Module          No LIM-3 sea-ice model 
     307   !!---------------------------------------------------------------------- 
    388308#endif 
    389309   !!====================================================================== 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limtrp.F90

    r2477 r2528  
    1111   !!   lim_trp_init : initialization and namelist read 
    1212   !!---------------------------------------------------------------------- 
    13    !! * Modules used 
    1413   USE phycst 
    1514   USE dom_oce 
    … …  
    4746#  include "vectopt_loop_substitute.h90" 
    4847   !!---------------------------------------------------------------------- 
    49    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
     48   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    5049   !! $Id$ 
    51    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     50   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    5251   !!---------------------------------------------------------------------- 
    5352 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limupdate.F90

    r2477 r2528  
    1313   !!   'key_lim3'                                      LIM3 sea-ice model 
    1414   !!---------------------------------------------------------------------- 
    15    !!    lim_update   : computes update of sea-ice global variables  
    16    !!                   from trend terms 
     15   !!    lim_update   : computes update of sea-ice global variables from trend terms 
    1716   !!---------------------------------------------------------------------- 
    18    !! * Modules used 
    19    USE limistate 
    2017   USE limrhg          ! ice rheology 
    2118   USE lbclnk 
    … …  
    5451 
    5552   !!---------------------------------------------------------------------- 
    56    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
     53   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    5754   !! $Id$ 
    58    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     55   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    5956   !!---------------------------------------------------------------------- 
    6057 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limvar.F90

    r1465 r2528  
    6363 
    6464   !!---------------------------------------------------------------------- 
    65    !!   LIM 3.0,  UCL-ASTR-LOCEAN-IPSL (2008) 
     65   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    6666   !! $Id$ 
    67    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     67   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    6868   !!---------------------------------------------------------------------- 
    6969 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limwri.F90

    r2477 r2528  
    11MODULE limwri 
    2    !!---------------------------------------------------------------------- 
    3    !!   'key_lim3'                                      LIM3 sea-ice model 
    4    !!---------------------------------------------------------------------- 
    52   !!====================================================================== 
    63   !!                     ***  MODULE  limwri  *** 
    … …  
    96#if defined key_lim3 
    107   !!---------------------------------------------------------------------- 
    11    !!  LIM 3.0, UCL-LOCEAN-IPSL (2008) 
    12    !! $Id$ 
    13    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     8   !!   'key_lim3'                                      LIM3 sea-ice model 
    149   !!---------------------------------------------------------------------- 
    1510   !!   lim_wri      : write of the diagnostics variables in ouput file  
    1611   !!   lim_wri_init : initialization and namelist read 
    1712   !!---------------------------------------------------------------------- 
    18    !! * Modules used 
    1913   USE ioipsl 
    2014   USE dianam          ! build name of file (routine) 
    … …  
    6357      zzero  = 0.e0     ,  & 
    6458      zone   = 1.e0 
    65  
     59       
     60   !!---------------------------------------------------------------------- 
     61   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
     62   !! $Id$ 
     63   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
     64   !!---------------------------------------------------------------------- 
    6665CONTAINS 
     66 
    6767#if defined key_dimgout 
    68  
    6968# include "limwri_dimg.h90" 
    70  
    7169#else 
    7270 
    … …  
    145143         CALL dia_nam ( clhstnam, nwrite, 'icemod' ) 
    146144         CALL histbeg ( clhstnam, jpi, glamt, jpj, gphit, 1, jpi, 1, jpj, niter, zjulian, rdt_ice,   & 
    147             &           nhorid, nice, domain_id=nidom ) 
     145            &           nhorid, nice, domain_id=nidom, snc4chunks=snc4set ) 
    148146         CALL histvert( nice, "deptht", "Vertical T levels", "m", 1, zdept, ndepid, "down") 
    149147         CALL wheneq  ( jpij , tmask(:,:,1), 1, 1., ndex51, ndim) 
    … …  
    159157         END DO 
    160158 
    161          CALL histend(nice) 
     159         CALL histend(nice, snc4set) 
    162160 
    163161         !----------------- 
    … …  
    174172            niter, zjulian, rdt_ice,   & ! time 
    175173            nhorida,                   & ! ? linked with horizontal ... 
    176             nicea , domain_id=nidom)                  ! file  
     174            nicea , domain_id=nidom, snc4chunks=snc4set)                  ! file  
    177175         CALL histvert( nicea, "icethi", "L levels",               & 
    178176            "m", ipl , hi_mean , nz ) 
    … …  
    194192         CALL histdef( nicea, "iice_etd", "Brine volume distr. "               , "%"    ,   &   
    195193            jpi, jpj, nhorida, jpl, 1, jpl, nz, 15, clop, zsto, zout ) 
    196          CALL histend(nicea) 
     194         CALL histend(nicea, snc4set) 
    197195      ENDIF 
    198196 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/limwri_dimg.h90

    r1694 r2528  
    11SUBROUTINE lim_wri 
    22   !!---------------------------------------------------------------------- 
    3    !!  LIM 3.0, UCL-LOCEAN-IPSL (2008) 
     3   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    44   !! $Id$ 
    5    !! This software is governed by the CeCILL licence see modipsl/doc/NEMO_CeCILL.txt 
     5   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    66   !!---------------------------------------------------------------------- 
    77   !!------------------------------------------------------------------- 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/par_ice.F90

    r1608 r2528  
    2020 
    2121   !!---------------------------------------------------------------------- 
    22    !! NEMO/LIM 3.2, UCL-ASTR-LOCEAN-IPSL (2009) 
     22   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    2323   !! $Id$ 
    24    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     24   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    2525   !!====================================================================== 
    2626END MODULE par_ice 

  • trunk/NEMOGCM/NEMO/LIM_SRC_3/thd_ice.F90

    r1465 r2528  
    77   !!   2.0  !  02-11  (C. Ethe)  F90: Free form and module 
    88   !!---------------------------------------------------------------------- 
    9    !!   LIM 3.0, UCL-LOCEAN-IPSL (2008) 
     9   !! NEMO/LIM3 3.3 , UCL - NEMO Consortium (2010) 
    1010   !! $Id$ 
    11    !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
     11   !! Software governed by the CeCILL licence (NEMOGCM/NEMO_CeCILL.txt) 
    1212   !!---------------------------------------------------------------------- 
    1313   !! * Modules used 
Note: See TracChangeset for help on using the changeset viewer.