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 5581 for branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90 – NEMO

Ignore:
Timestamp:
2015-07-10T13:28:53+02:00 (9 years ago)
Author:
timgraham
Message:

Merged head of trunk into branch

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2014/dev_r4765_CNRS_agrif/NEMOGCM/NEMO/LIM_SRC_3/limistate.F90

    r4793 r5581  
    66   !! History :  2.0  ! 2004-01 (C. Ethe, G. Madec)  Original code 
    77   !!            4.0  ! 2011-02 (G. Madec) dynamical allocation 
    8    !!             -   ! 2012    (C. Rousset) add par_oce (for jp_sal)...bug? 
     8   !!             -   ! 2014    (C. Rousset) add N/S initializations 
    99   !!---------------------------------------------------------------------- 
    1010#if defined key_lim3 
     
    2222   USE eosbn2           ! equation of state 
    2323   USE ice              ! sea-ice variables 
    24    USE par_ice          ! ice parameters 
    2524   USE par_oce          ! ocean parameters 
    2625   USE dom_ice          ! sea-ice domain 
     
    2928   USE lib_fortran      ! Fortran utilities (allows no signed zero when 'key_nosignedzero' defined)   
    3029   USE wrk_nemo         ! work arrays 
    31    USE cpl_oasis3, ONLY : lk_cpl 
    3230 
    3331   IMPLICIT NONE 
     
    3634   PUBLIC   lim_istate      ! routine called by lim_init.F90 
    3735 
    38    !! * Module variables 
    3936   !                          !!** init namelist (namiceini) ** 
    40    REAL(wp) ::   thres_sst   ! threshold water temperature for initial sea ice 
    41    REAL(wp) ::   hts_ini_n   ! initial snow thickness in the north 
    42    REAL(wp) ::   hts_ini_s   ! initial snow thickness in the south 
    43    REAL(wp) ::   hti_ini_n   ! initial ice thickness in the north 
    44    REAL(wp) ::   hti_ini_s   ! initial ice thickness in the south 
    45    REAL(wp) ::   ati_ini_n   ! initial leads area in the north 
    46    REAL(wp) ::   ati_ini_s   ! initial leads area in the south 
    47    REAL(wp) ::   smi_ini_n   ! initial salinity  
    48    REAL(wp) ::   smi_ini_s   ! initial salinity 
    49    REAL(wp) ::   tmi_ini_n   ! initial temperature 
    50    REAL(wp) ::   tmi_ini_s   ! initial temperature 
    51  
    52    LOGICAL  ::  ln_limini    ! initialization or not 
     37   REAL(wp) ::   rn_thres_sst   ! threshold water temperature for initial sea ice 
     38   REAL(wp) ::   rn_hts_ini_n   ! initial snow thickness in the north 
     39   REAL(wp) ::   rn_hts_ini_s   ! initial snow thickness in the south 
     40   REAL(wp) ::   rn_hti_ini_n   ! initial ice thickness in the north 
     41   REAL(wp) ::   rn_hti_ini_s   ! initial ice thickness in the south 
     42   REAL(wp) ::   rn_ati_ini_n   ! initial leads area in the north 
     43   REAL(wp) ::   rn_ati_ini_s   ! initial leads area in the south 
     44   REAL(wp) ::   rn_smi_ini_n   ! initial salinity  
     45   REAL(wp) ::   rn_smi_ini_s   ! initial salinity 
     46   REAL(wp) ::   rn_tmi_ini_n   ! initial temperature 
     47   REAL(wp) ::   rn_tmi_ini_s   ! initial temperature 
     48 
     49   LOGICAL  ::  ln_iceini    ! initialization or not 
    5350   !!---------------------------------------------------------------------- 
    5451   !!   LIM 3.0,  UCL-LOCEAN-IPSL (2008) 
     
    5653   !! Software governed by the CeCILL licence (modipsl/doc/NEMO_CeCILL.txt) 
    5754   !!---------------------------------------------------------------------- 
    58  
    5955CONTAINS 
    6056 
     
    8076      !! 
    8177      !! ** Notes   : o_i, t_su, t_s, t_i, s_i must be filled everywhere, even 
    82       !!              where there is no ice (clem: I do not know why but it is mandatory)  
     78      !!              where there is no ice (clem: I do not know why, is it mandatory?)  
    8379      !! 
    8480      !! History : 
     
    9086      !! * Local variables 
    9187      INTEGER    :: ji, jj, jk, jl             ! dummy loop indices 
    92       REAL(wp)   :: epsi20, ztmelts, zdh 
     88      REAL(wp)   :: ztmelts, zdh 
    9389      INTEGER    :: i_hemis, i_fill, jl0   
    9490      REAL(wp)   :: ztest_1, ztest_2, ztest_3, ztest_4, ztests, zsigma, zarg, zA, zV, zA_cons, zV_cons, zconv 
     
    104100      CALL wrk_alloc(   2,      zht_i_ini, zat_i_ini, zvt_i_ini, zht_s_ini, zsm_i_ini, ztm_i_ini ) 
    105101 
    106       epsi20   = 1.e-20_wp 
    107  
    108102      IF(lwp) WRITE(numout,*) 
    109103      IF(lwp) WRITE(numout,*) 'lim_istate : Ice initialization ' 
     
    116110      CALL lim_istate_init     !  reading the initials parameters of the ice 
    117111 
    118 # if defined key_coupled 
    119       albege(:,:)   = 0.8 * tms(:,:) 
    120 # endif 
    121  
    122112      ! surface temperature 
    123113      DO jl = 1, jpl ! loop over categories 
    124          t_su  (:,:,jl) = rtt * tms(:,:) 
    125          tn_ice(:,:,jl) = rtt * tms(:,:) 
    126       END DO 
    127       ! Basal temperature is set to the freezing point of seawater in Kelvin 
    128       t_bo(:,:) = ( tfreez( tsn(:,:,1,jp_sal) ) + rt0 ) * tms(:,:)  
    129  
    130       IF( ln_limini ) THEN 
     114         t_su  (:,:,jl) = rt0 * tmask(:,:,1) 
     115         tn_ice(:,:,jl) = rt0 * tmask(:,:,1) 
     116      END DO 
     117 
     118      ! basal temperature (considered at freezing point) 
     119      CALL eos_fzp( sss_m(:,:), t_bo(:,:) ) 
     120      t_bo(:,:) = ( t_bo(:,:) + rt0 ) * tmask(:,:,1)  
     121 
     122      IF( ln_iceini ) THEN 
    131123 
    132124      !-------------------------------------------------------------------- 
    133125      ! 2) Basal temperature, ice mask and hemispheric index 
    134126      !-------------------------------------------------------------------- 
    135       ! ice if sst <= t-freez + thres_sst 
    136       DO jj = 1, jpj                                        
     127 
     128      DO jj = 1, jpj                                       ! ice if sst <= t-freez + ttest 
    137129         DO ji = 1, jpi 
    138             IF( ( tsn(ji,jj,1,jp_tem)  - ( t_bo(ji,jj) - rt0 ) ) * tms(ji,jj) >= thres_sst ) THEN  
    139                zswitch(ji,jj) = 0._wp * tms(ji,jj)    ! no ice 
     130            IF( ( sst_m(ji,jj)  - ( t_bo(ji,jj) - rt0 ) ) * tmask(ji,jj,1) >= rn_thres_sst ) THEN  
     131               zswitch(ji,jj) = 0._wp * tmask(ji,jj,1)    ! no ice 
    140132            ELSE                                                                                    
    141                zswitch(ji,jj) = 1._wp * tms(ji,jj)    !    ice 
     133               zswitch(ji,jj) = 1._wp * tmask(ji,jj,1)    !    ice 
    142134            ENDIF 
    143135         END DO 
     
    146138 
    147139      ! Hemispheric index 
    148       ! MV 2011 new initialization 
    149140      DO jj = 1, jpj 
    150141         DO ji = 1, jpi 
     
    156147         END DO 
    157148      END DO 
    158       ! END MV 2011 new initialization 
    159149 
    160150      !-------------------------------------------------------------------- 
     
    166156      !----------------------------- 
    167157      ! assign initial thickness, concentration, snow depth and salinity to an hemisphere-dependent array 
    168       zht_i_ini(1) = hti_ini_n ; zht_i_ini(2) = hti_ini_s  ! ice thickness 
    169       zht_s_ini(1) = hts_ini_n ; zht_s_ini(2) = hts_ini_s  ! snow depth 
    170       zat_i_ini(1) = ati_ini_n ; zat_i_ini(2) = ati_ini_s  ! ice concentration 
    171       zsm_i_ini(1) = smi_ini_n ; zsm_i_ini(2) = smi_ini_s  ! bulk ice salinity 
    172       ztm_i_ini(1) = tmi_ini_n ; ztm_i_ini(2) = tmi_ini_s  ! temperature (ice and snow) 
     158      zht_i_ini(1) = rn_hti_ini_n ; zht_i_ini(2) = rn_hti_ini_s  ! ice thickness 
     159      zht_s_ini(1) = rn_hts_ini_n ; zht_s_ini(2) = rn_hts_ini_s  ! snow depth 
     160      zat_i_ini(1) = rn_ati_ini_n ; zat_i_ini(2) = rn_ati_ini_s  ! ice concentration 
     161      zsm_i_ini(1) = rn_smi_ini_n ; zsm_i_ini(2) = rn_smi_ini_s  ! bulk ice salinity 
     162      ztm_i_ini(1) = rn_tmi_ini_n ; ztm_i_ini(2) = rn_tmi_ini_s  ! temperature (ice and snow) 
    173163 
    174164      zvt_i_ini(:) = zht_i_ini(:) * zat_i_ini(:)   ! ice volume 
     
    205195               !--- Ice thicknesses in the i_fill - 1 first categories 
    206196               DO jl = 1, i_fill - 1 
    207                   zh_i_ini(jl,i_hemis)    = 0.5 * ( hi_max(jl) + hi_max(jl-1) ) 
     197                  zh_i_ini(jl,i_hemis) = hi_mean(jl) 
    208198               END DO 
    209199                
    210200               !--- jl0: most likely index where cc will be maximum 
    211201               DO jl = 1, jpl 
    212                   IF ( ( zht_i_ini(i_hemis) .GT. hi_max(jl-1) ) .AND. & 
    213                      ( zht_i_ini(i_hemis) .LE. hi_max(jl)   ) ) THEN 
     202                  IF ( ( zht_i_ini(i_hemis) > hi_max(jl-1) ) .AND. & 
     203                     & ( zht_i_ini(i_hemis) <= hi_max(jl)   ) ) THEN 
    214204                     jl0 = jl 
    215205                  ENDIF 
     
    275265 
    276266            ! Test 3: thickness of the last category is in-bounds ? 
    277             IF ( zh_i_ini(i_fill, i_hemis) .GT. hi_max(i_fill-1) ) THEN 
     267            IF ( zh_i_ini(i_fill, i_hemis) > hi_max(i_fill-1) ) THEN 
    278268               ztest_3 = 1 
    279269            ELSE 
     
    325315            DO ji = 1, jpi 
    326316               a_i(ji,jj,jl)   = zswitch(ji,jj) * za_i_ini (jl,zhemis(ji,jj))  ! concentration 
    327                ht_i(ji,jj,jl)  = zswitch(ji,jj) * zh_i_ini(jl,zhemis(ji,jj))  ! ice thickness 
     317               ht_i(ji,jj,jl)  = zswitch(ji,jj) * zh_i_ini(jl,zhemis(ji,jj))   ! ice thickness 
    328318               ht_s(ji,jj,jl)  = ht_i(ji,jj,jl) * ( zht_s_ini( zhemis(ji,jj) ) / zht_i_ini( zhemis(ji,jj) ) )  ! snow depth 
    329                sm_i(ji,jj,jl)  = zswitch(ji,jj) * zsm_i_ini(zhemis(ji,jj)) !+ ( 1._wp - zswitch(ji,jj) ) * s_i_min ! salinity 
    330                o_i(ji,jj,jl)   = zswitch(ji,jj) * 1._wp + ( 1._wp - zswitch(ji,jj) ) ! age 
    331                t_su(ji,jj,jl)  = zswitch(ji,jj) * ztm_i_ini(zhemis(ji,jj)) + ( 1._wp - zswitch(ji,jj) ) * rtt ! surf temp 
     319               sm_i(ji,jj,jl)  = zswitch(ji,jj) * zsm_i_ini(zhemis(ji,jj))     ! salinity 
     320               o_i(ji,jj,jl)   = zswitch(ji,jj) * 1._wp                        ! age (1 day) 
     321               t_su(ji,jj,jl)  = zswitch(ji,jj) * ztm_i_ini(zhemis(ji,jj)) + ( 1._wp - zswitch(ji,jj) ) * rt0 ! surf temp 
    332322 
    333323               ! This case below should not be used if (ht_s/ht_i) is ok in namelist 
     
    337327               ! recompute ht_i, ht_s avoiding out of bounds values 
    338328               ht_i(ji,jj,jl) = MIN( hi_max(jl), ht_i(ji,jj,jl) + zdh ) 
    339                ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic / rhosn ) 
     329               ht_s(ji,jj,jl) = MAX( 0._wp, ht_s(ji,jj,jl) - zdh * rhoic * r1_rhosn ) 
    340330 
    341331               ! ice volume, salt content, age content 
     
    344334               smv_i(ji,jj,jl) = MIN( sm_i(ji,jj,jl) , sss_m(ji,jj) ) * v_i(ji,jj,jl) ! salt content 
    345335               oa_i(ji,jj,jl)  = o_i(ji,jj,jl) * a_i(ji,jj,jl)               ! age content 
    346             END DO ! ji 
    347          END DO ! jj 
    348       END DO ! jl 
     336            END DO 
     337         END DO 
     338      END DO 
    349339 
    350340      ! Snow temperature and heat content 
     
    353343            DO jj = 1, jpj 
    354344               DO ji = 1, jpi 
    355                    t_s(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(zhemis(ji,jj)) + ( 1._wp - zswitch(ji,jj) ) * rtt 
     345                   t_s(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(zhemis(ji,jj)) + ( 1._wp - zswitch(ji,jj) ) * rt0 
    356346                   ! Snow energy of melting 
    357                    e_s(ji,jj,jk,jl) = zswitch(ji,jj) * rhosn * ( cpic * ( rtt - t_s(ji,jj,jk,jl) ) + lfus ) 
    358                    ! Change dimensions 
    359                    e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) / unit_fac 
    360                    ! Multiply by volume, so that heat content in Joules 
    361                    e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * area(ji,jj) * v_s(ji,jj,jl) / nlay_s 
    362                END DO ! ji 
    363             END DO ! jj 
    364          END DO ! jl 
    365       END DO ! jk 
     347                   e_s(ji,jj,jk,jl) = zswitch(ji,jj) * rhosn * ( cpic * ( rt0 - t_s(ji,jj,jk,jl) ) + lfus ) 
     348 
     349                   ! Mutliply by volume, and divide by number of layers to get heat content in J/m2 
     350                   e_s(ji,jj,jk,jl) = e_s(ji,jj,jk,jl) * v_s(ji,jj,jl) * r1_nlay_s 
     351               END DO 
     352            END DO 
     353         END DO 
     354      END DO 
    366355 
    367356      ! Ice salinity, temperature and heat content 
     
    370359            DO jj = 1, jpj 
    371360               DO ji = 1, jpi 
    372                    t_i(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(zhemis(ji,jj)) + ( 1._wp - zswitch(ji,jj) ) * rtt  
    373                    s_i(ji,jj,jk,jl) = zswitch(ji,jj) * zsm_i_ini(zhemis(ji,jj)) !+ ( 1._wp - zswitch(ji,jj) ) * s_i_min 
    374                    ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rtt !Melting temperature in K 
     361                   t_i(ji,jj,jk,jl) = zswitch(ji,jj) * ztm_i_ini(zhemis(ji,jj)) + ( 1._wp - zswitch(ji,jj) ) * rt0  
     362                   s_i(ji,jj,jk,jl) = zswitch(ji,jj) * zsm_i_ini(zhemis(ji,jj)) !+ ( 1._wp - zswitch(ji,jj) ) * rn_simin 
     363                   ztmelts          = - tmut * s_i(ji,jj,jk,jl) + rt0 !Melting temperature in K 
    375364 
    376365                   ! heat content per unit volume 
    377366                   e_i(ji,jj,jk,jl) = zswitch(ji,jj) * rhoic * (   cpic    * ( ztmelts - t_i(ji,jj,jk,jl) ) & 
    378                       +   lfus    * ( 1._wp - (ztmelts-rtt) / MIN((t_i(ji,jj,jk,jl)-rtt),-epsi20) ) & 
    379                       -   rcp     * ( ztmelts - rtt ) ) 
    380  
    381                    ! Correct dimensions to avoid big values 
    382                    e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) / unit_fac  
    383  
    384                    ! Mutliply by ice volume, and divide by number of layers to get heat content in J 
    385                    e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * area(ji,jj) * v_i(ji,jj,jl) / nlay_i 
    386                END DO ! ji 
    387             END DO ! jj 
    388          END DO ! jl 
    389       END DO ! jk 
     367                      +   lfus    * ( 1._wp - (ztmelts-rt0) / MIN((t_i(ji,jj,jk,jl)-rt0),-epsi20) ) & 
     368                      -   rcp     * ( ztmelts - rt0 ) ) 
     369 
     370                   ! Mutliply by ice volume, and divide by number of layers to get heat content in J/m2 
     371                   e_i(ji,jj,jk,jl) = e_i(ji,jj,jk,jl) * v_i(ji,jj,jl) * r1_nlay_i 
     372               END DO 
     373            END DO 
     374         END DO 
     375      END DO 
    390376 
    391377      tn_ice (:,:,:) = t_su (:,:,:) 
    392378 
    393379      ELSE  
    394          ! if ln_limini=false 
     380         ! if ln_iceini=false 
    395381         a_i  (:,:,:) = 0._wp 
    396382         v_i  (:,:,:) = 0._wp 
     
    408394         DO jl = 1, jpl 
    409395            DO jk = 1, nlay_i 
    410                t_i(:,:,jk,jl) = rtt * tms(:,:) 
     396               t_i(:,:,jk,jl) = rt0 * tmask(:,:,1) 
    411397            END DO 
    412398            DO jk = 1, nlay_s 
    413                t_s(:,:,jk,jl) = rtt * tms(:,:) 
     399               t_s(:,:,jk,jl) = rt0 * tmask(:,:,1) 
    414400            END DO 
    415401         END DO 
    416402       
    417       ENDIF ! ln_limini 
     403      ENDIF ! ln_iceini 
    418404       
    419405      at_i (:,:) = 0.0_wp 
     
    489475      !!  8.5  ! 07-11 (M. Vancoppenolle) rewritten initialization 
    490476      !!----------------------------------------------------------------------------- 
    491       NAMELIST/namiceini/ ln_limini, thres_sst, hts_ini_n, hts_ini_s, hti_ini_n, hti_ini_s,  & 
    492          &                                      ati_ini_n, ati_ini_s, smi_ini_n, smi_ini_s, tmi_ini_n, tmi_ini_s 
     477      NAMELIST/namiceini/ ln_iceini, rn_thres_sst, rn_hts_ini_n, rn_hts_ini_s, rn_hti_ini_n, rn_hti_ini_s,  & 
     478         &                                      rn_ati_ini_n, rn_ati_ini_s, rn_smi_ini_n, rn_smi_ini_s, rn_tmi_ini_n, rn_tmi_ini_s 
    493479      INTEGER :: ios                 ! Local integer output status for namelist read 
    494480      !!----------------------------------------------------------------------------- 
     
    510496         WRITE(numout,*) 'lim_istate_init : ice parameters inititialisation ' 
    511497         WRITE(numout,*) '~~~~~~~~~~~~~~~' 
    512          WRITE(numout,*) '   initialization with ice (T) or not (F)       ln_limini   = ', ln_limini 
    513          WRITE(numout,*) '   threshold water temp. for initial sea-ice    thres_sst  = ', thres_sst 
    514          WRITE(numout,*) '   initial snow thickness in the north          hts_ini_n  = ', hts_ini_n 
    515          WRITE(numout,*) '   initial snow thickness in the south          hts_ini_s  = ', hts_ini_s  
    516          WRITE(numout,*) '   initial ice thickness  in the north          hti_ini_n  = ', hti_ini_n 
    517          WRITE(numout,*) '   initial ice thickness  in the south          hti_ini_s  = ', hti_ini_s 
    518          WRITE(numout,*) '   initial ice concentr.  in the north          ati_ini_n  = ', ati_ini_n 
    519          WRITE(numout,*) '   initial ice concentr.  in the north          ati_ini_s  = ', ati_ini_s 
    520          WRITE(numout,*) '   initial  ice salinity  in the north          smi_ini_n  = ', smi_ini_n 
    521          WRITE(numout,*) '   initial  ice salinity  in the south          smi_ini_s  = ', smi_ini_s 
    522          WRITE(numout,*) '   initial  ice/snw temp  in the north          tmi_ini_n  = ', tmi_ini_n 
    523          WRITE(numout,*) '   initial  ice/snw temp  in the south          tmi_ini_s  = ', tmi_ini_s 
     498         WRITE(numout,*) '   initialization with ice (T) or not (F)       ln_iceini     = ', ln_iceini 
     499         WRITE(numout,*) '   threshold water temp. for initial sea-ice    rn_thres_sst  = ', rn_thres_sst 
     500         WRITE(numout,*) '   initial snow thickness in the north          rn_hts_ini_n  = ', rn_hts_ini_n 
     501         WRITE(numout,*) '   initial snow thickness in the south          rn_hts_ini_s  = ', rn_hts_ini_s  
     502         WRITE(numout,*) '   initial ice thickness  in the north          rn_hti_ini_n  = ', rn_hti_ini_n 
     503         WRITE(numout,*) '   initial ice thickness  in the south          rn_hti_ini_s  = ', rn_hti_ini_s 
     504         WRITE(numout,*) '   initial ice concentr.  in the north          rn_ati_ini_n  = ', rn_ati_ini_n 
     505         WRITE(numout,*) '   initial ice concentr.  in the north          rn_ati_ini_s  = ', rn_ati_ini_s 
     506         WRITE(numout,*) '   initial  ice salinity  in the north          rn_smi_ini_n  = ', rn_smi_ini_n 
     507         WRITE(numout,*) '   initial  ice salinity  in the south          rn_smi_ini_s  = ', rn_smi_ini_s 
     508         WRITE(numout,*) '   initial  ice/snw temp  in the north          rn_tmi_ini_n  = ', rn_tmi_ini_n 
     509         WRITE(numout,*) '   initial  ice/snw temp  in the south          rn_tmi_ini_s  = ', rn_tmi_ini_s 
    524510      ENDIF 
    525511 
Note: See TracChangeset for help on using the changeset viewer.