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 7355 for branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90 – NEMO

Ignore:
Timestamp:
2016-11-28T18:21:42+01:00 (8 years ago)
Author:
flavoni
Message:

merge branch dev_CNRS_2016 & dev_CNRS_AGRIF_LIM3

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2016/dev_CNRS_2016/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90

    r7280 r7355  
    2222   USE phycst         ! physical constants 
    2323   USE dom_oce        ! ocean space and time domain variables 
    24    USE ice            ! LIM: sea-ice variables 
     24   USE ice            ! sea-ice variables 
    2525   USE sbc_oce        ! Surface boundary condition: ocean fields 
    2626   USE sbc_ice        ! Surface boundary condition: ice fields 
    27    USE dom_ice        ! LIM: sea-ice domain 
    28    USE thd_ice        ! LIM: thermodynamic sea-ice variables 
    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 limthd_lac     ! LIM: lateral accretion 
    34    USE limitd_th      ! LIM: remapping thickness distribution 
    35    USE limtab         ! LIM: 1D <==> 2D transformation 
    36    USE limvar         ! LIM: sea-ice variables 
    37    USE limcons        ! LIM: conservation tests 
    38    USE limctl         ! LIM: control print 
     27   USE thd_ice        ! thermodynamic sea-ice variables 
     28   USE limthd_dif     ! vertical diffusion 
     29   USE limthd_dh      ! ice-snow growth and melt 
     30   USE limthd_da      ! lateral melting 
     31   USE limthd_sal     ! ice salinity 
     32   USE limthd_ent     ! ice enthalpy redistribution 
     33   USE limthd_lac     ! lateral accretion 
     34   USE limitd_th      ! remapping thickness distribution 
     35   USE limtab         ! 1D <==> 2D transformation 
     36   USE limvar         ! 
     37   USE limcons        ! conservation tests 
     38   USE limctl         ! control print 
    3939   ! 
    4040   USE in_out_manager ! I/O manager 
    41    USE prtctl         ! Print control 
    4241   USE lbclnk         ! lateral boundary condition - MPP links 
    4342   USE lib_mpp        ! MPP library 
     
    8887      REAL(wp), PARAMETER :: zfric_umin = 0._wp           ! lower bound for the friction velocity (cice value=5.e-04) 
    8988      REAL(wp), PARAMETER :: zch        = 0.0057_wp       ! heat transfer coefficient 
     89      REAL(wp), POINTER, DIMENSION(:,:) ::   zu_io, zv_io, zfric   ! ice-ocean velocity (m/s) and frictional velocity (m2/s2) 
     90      ! 
    9091      !!------------------------------------------------------------------- 
    9192 
    9293      IF( nn_timing == 1 )   CALL timing_start('limthd') 
    9394 
     95      CALL wrk_alloc( jpi,jpj, zu_io, zv_io, zfric ) 
     96 
     97      IF( kt == nit000 .AND. lwp ) THEN 
     98         WRITE(numout,*)''  
     99         WRITE(numout,*)' lim_thd ' 
     100         WRITE(numout,*)' ~~~~~~~~' 
     101      ENDIF 
     102       
    94103      ! conservation test 
    95       IF( ln_limdiahsb )   CALL lim_cons_hsm( 0, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b ) 
     104      IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    96105 
    97106      CALL lim_var_glo2eqv 
    98       !------------------------------------------------------------------------! 
    99       ! 1) Initialization of some variables                                    ! 
    100       !------------------------------------------------------------------------! 
     107 
     108      !---------------------------------------------! 
     109      ! computation of friction velocity at T points 
     110      !---------------------------------------------! 
     111      IF( ln_limdyn ) THEN 
     112         zu_io(:,:) = u_ice(:,:) - ssu_m(:,:) 
     113         zv_io(:,:) = v_ice(:,:) - ssv_m(:,:) 
     114         DO jj = 2, jpjm1  
     115            DO ji = fs_2, fs_jpim1 
     116               zfric(ji,jj) = rn_cio * ( 0.5_wp *  & 
     117                  &                    (  zu_io(ji,jj) * zu_io(ji,jj) + zu_io(ji-1,jj) * zu_io(ji-1,jj)   & 
     118                  &                     + zv_io(ji,jj) * zv_io(ji,jj) + zv_io(ji,jj-1) * zv_io(ji,jj-1) ) ) * tmask(ji,jj,1) 
     119            END DO 
     120         END DO 
     121      ELSE      !  if no ice dynamics => transmit directly the atmospheric stress to the ocean 
     122         DO jj = 2, jpjm1 
     123            DO ji = fs_2, fs_jpim1 
     124               zfric(ji,jj) = r1_rau0 * SQRT( 0.5_wp *  & 
     125                  &                         (  utau(ji,jj) * utau(ji,jj) + utau(ji-1,jj) * utau(ji-1,jj)   & 
     126                  &                          + vtau(ji,jj) * vtau(ji,jj) + vtau(ji,jj-1) * vtau(ji,jj-1) ) ) * tmask(ji,jj,1) 
     127            END DO 
     128         END DO 
     129      ENDIF 
     130      CALL lbc_lnk( zfric, 'T',  1. ) 
     131      ! 
     132      !----------------------------------! 
     133      ! Initialization and units change 
     134      !----------------------------------! 
    101135      ftr_ice(:,:,:) = 0._wp  ! part of solar radiation transmitted through the ice 
    102136 
    103       !-------------------- 
    104       ! 1.2) Heat content     
    105       !-------------------- 
    106137      ! Change the units of heat content; from J/m2 to J/m3 
    107138      DO jl = 1, jpl 
     
    109140            DO jj = 1, jpj 
    110141               DO ji = 1, jpi 
    111                   !0 if no ice and 1 if yes 
    112142                  rswitch = MAX(  0._wp , SIGN( 1._wp , v_i(ji,jj,jl) - epsi20 )  ) 
    113143                  !Energy of melting q(S,T) [J.m-3] 
     
    119149            DO jj = 1, jpj 
    120150               DO ji = 1, jpi 
    121                   !0 if no ice and 1 if yes 
    122151                  rswitch = MAX(  0._wp , SIGN( 1._wp , v_s(ji,jj,jl) - epsi20 )  ) 
    123152                  !Energy of melting q(S,T) [J.m-3] 
     
    128157      END DO 
    129158 
    130       ! 2) Partial computation of forcing for the thermodynamic sea ice model.      ! 
    131       !-----------------------------------------------------------------------------! 
     159      !--------------------------------------------------------------------! 
     160      ! Partial computation of forcing for the thermodynamic sea ice model 
     161      !--------------------------------------------------------------------! 
    132162      DO jj = 1, jpj 
    133163         DO ji = 1, jpi 
     
    148178 
    149179            ! --- Energy from the turbulent oceanic heat flux (W/m2) --- ! 
    150             zfric_u      = MAX( SQRT( ust2s(ji,jj) ), zfric_umin )  
     180            zfric_u      = MAX( SQRT( zfric(ji,jj) ), zfric_umin )  
    151181            fhtur(ji,jj) = MAX( 0._wp, rswitch * rau0 * rcp * zch  * zfric_u * ( ( sst_m(ji,jj) + rt0 ) - t_bo(ji,jj) ) ) ! W.m-2 
    152182            fhtur(ji,jj) = rswitch * MIN( fhtur(ji,jj), - zqfr * r1_rdtice / MAX( at_i(ji,jj), epsi10 ) ) 
     
    166196            ENDIF 
    167197            ! 
    168             ! ----------------------------------------- 
    169             ! Net heat flux on top of ice-ocean [W.m-2] 
    170             ! ----------------------------------------- 
     198            ! Net heat flux on top of the ice-ocean [W.m-2] 
     199            ! --------------------------------------------- 
    171200            hfx_in(ji,jj) = qns_tot(ji,jj) + qsr_tot(ji,jj)  
    172  
    173             ! ----------------------------------------------------------------------------- 
    174             ! Net heat flux on top of the ocean after ice thermo (1st step) [W.m-2] 
    175             ! ----------------------------------------------------------------------------- 
    176             !     First  step here              :  non solar + precip - qlead - qturb 
    177             !     Second step in limthd_dh      :  heat remaining if total melt (zq_rema)  
    178             !     Third  step in limsbc         :  heat from ice-ocean mass exchange (zf_mass) + solar 
     201         END DO 
     202      END DO 
     203       
     204      ! In case we bypass open-water ice formation 
     205      IF( .NOT. ln_limdO )  qlead(:,:) = 0._wp 
     206      ! In case we bypass growing/melting from top and bottom: we suppose ice is impermeable => ocean is isolated from atmosphere 
     207      IF( .NOT. ln_limdH )  hfx_in(:,:) = pfrld(:,:) * ( qns_oce(:,:) + qsr_oce(:,:) ) + qemp_oce(:,:) 
     208      IF( .NOT. ln_limdH )  fhtur (:,:) = 0._wp  ;  fhld  (:,:) = 0._wp 
     209 
     210      ! --------------------------------------------------------------------- 
     211      ! Net heat flux on top of the ocean after ice thermo (1st step) [W.m-2] 
     212      ! --------------------------------------------------------------------- 
     213      !     First  step here              :  non solar + precip - qlead - qturb 
     214      !     Second step in limthd_dh      :  heat remaining if total melt (zq_rema)  
     215      !     Third  step in limsbc         :  heat from ice-ocean mass exchange (zf_mass) + solar 
     216      DO jj = 1, jpj 
     217         DO ji = 1, jpi 
    179218            hfx_out(ji,jj) =   pfrld(ji,jj) * qns_oce(ji,jj) + qemp_oce(ji,jj)  &  ! Non solar heat flux received by the ocean                
    180219               &             - qlead(ji,jj) * r1_rdtice                         &  ! heat flux taken from the ocean where there is open water ice formation 
     
    186225 
    187226      !------------------------------------------------------------------------------! 
    188       ! 3) Select icy points and fulfill arrays for the vectorial grid.             
     227      ! Thermodynamic computation (only on grid points covered by ice) 
    189228      !------------------------------------------------------------------------------! 
    190229 
    191230      DO jl = 1, jpl      !loop over ice categories 
    192231 
    193          IF( kt == nit000 .AND. lwp ) THEN 
    194             WRITE(numout,*) ' lim_thd : transfer to 1D vectors. Category no : ', jl  
    195             WRITE(numout,*) ' ~~~~~~~~' 
    196          ENDIF 
    197  
     232         ! select ice covered grid points 
    198233         nbpb = 0 
    199234         DO jj = 1, jpj 
     
    208243         ! debug point to follow 
    209244         jiindex_1d = 0 
    210          IF( ln_icectl ) THEN 
     245         IF( ln_limctl ) THEN 
    211246            DO ji = mi0(iiceprt), mi1(iiceprt) 
    212247               DO jj = mj0(jiceprt), mj1(jiceprt) 
     
    217252         ENDIF 
    218253 
    219          !------------------------------------------------------------------------------! 
    220          ! 4) Thermodynamic computation 
    221          !------------------------------------------------------------------------------! 
    222  
    223          IF( lk_mpp )   CALL mpp_ini_ice( nbpb , numout ) 
     254         IF( lk_mpp )         CALL mpp_ini_ice( nbpb , numout ) 
    224255 
    225256         IF( nbpb > 0 ) THEN  ! If there is no ice, do nothing. 
    226             ! 
    227             CALL lim_thd_1d2d( nbpb, jl, 1 )                ! --- Move to 1D arrays ---! 
    228             ! 
    229             CALL lim_thd_dif ( 1, nbpb )                    ! --- Ice/Snow Temperature profile --- ! 
    230             ! 
    231             CALL lim_thd_dh  ( 1, nbpb )                    ! --- Ice/Snow thickness ---! 
    232             ! 
    233             CALL lim_thd_ent ( 1, nbpb, q_i_1d(1:nbpb,:) )  ! --- Ice enthalpy remapping --- ! 
    234             ! 
    235             CALL lim_thd_sal ( 1, nbpb )                    ! --- Ice salinity ---            ! 
    236             ! 
    237             CALL lim_thd_temp( 1, nbpb )                    ! --- temperature update ---      ! 
    238             ! 
    239             !                                               ! --- lateral melting if monocat --- ! 
    240             IF ( ( nn_monocat == 1 .OR. nn_monocat == 4 ) .AND. jpl == 1 ) THEN 
    241                CALL lim_thd_lam( 1, nbpb ) 
     257            !                                                                 
     258            s_i_new   (:) = 0._wp ; dh_s_tot (:) = 0._wp                     ! --- some init --- ! 
     259            dh_i_surf (:) = 0._wp ; dh_i_bott(:) = 0._wp 
     260            dh_snowice(:) = 0._wp ; dh_i_sub (:) = 0._wp 
     261 
     262                              CALL lim_thd_1d2d( nbpb, jl, 1 )               ! --- Move to 1D arrays --- ! 
     263            ! 
     264            IF( ln_limdH )    CALL lim_thd_dif( 1, nbpb )                    ! --- Ice/Snow Temperature profile --- ! 
     265            ! 
     266            IF( ln_limdH )    CALL lim_thd_dh( 1, nbpb )                     ! --- Ice/Snow thickness --- !     
     267            ! 
     268            IF( ln_limdH )    CALL lim_thd_ent( 1, nbpb, q_i_1d(1:nbpb,:) )  ! --- Ice enthalpy remapping --- ! 
     269            ! 
     270                              CALL lim_thd_sal( 1, nbpb )                    ! --- Ice salinity --- !     
     271            ! 
     272                              CALL lim_thd_temp( 1, nbpb )                   ! --- temperature update --- ! 
     273            ! 
     274            IF( ln_limdH ) THEN 
     275               IF ( ( nn_monocat == 1 .OR. nn_monocat == 4 ) .AND. jpl == 1 ) THEN 
     276                              CALL lim_thd_lam( 1, nbpb )                    ! --- extra lateral melting if monocat --- ! 
     277               END IF 
    242278            END IF 
    243279            ! 
    244             CALL lim_thd_1d2d( nbpb, jl, 2 )                ! --- Move to 2D arrays --- 
    245             ! 
    246             IF( lk_mpp )   CALL mpp_comm_free( ncomm_ice ) !RB necessary ?? 
     280                              CALL lim_thd_1d2d( nbpb, jl, 2 )               ! --- Move to 2D arrays --- ! 
     281            ! 
     282            IF( lk_mpp )      CALL mpp_comm_free( ncomm_ice ) !RB necessary ?? 
    247283         ENDIF 
    248284         ! 
    249285      END DO !jl 
    250286 
    251       !------------------------------------------------------------------------------! 
    252       ! 5) Global variables, diagnostics 
    253       !------------------------------------------------------------------------------! 
    254  
    255       !------------------------ 
    256       ! Ice heat content               
    257       !------------------------ 
     287      IF( ln_limdA)           CALL lim_thd_da                                ! --- lateral melting --- ! 
     288 
    258289      ! Enthalpies are global variables we have to readjust the units (heat content in J/m2) 
    259290      DO jl = 1, jpl 
     
    261292            e_i(:,:,jk,jl) = e_i(:,:,jk,jl) * a_i(:,:,jl) * ht_i(:,:,jl) * r1_nlay_i 
    262293         END DO 
    263       END DO 
    264  
    265       !------------------------ 
    266       ! Snow heat content               
    267       !------------------------ 
    268       ! Enthalpies are global variables we have to readjust the units (heat content in J/m2) 
    269       DO jl = 1, jpl 
    270294         DO jk = 1, nlay_s 
    271295            e_s(:,:,jk,jl) = e_s(:,:,jk,jl) * a_i(:,:,jl) * ht_s(:,:,jl) * r1_nlay_s 
     
    273297      END DO 
    274298  
    275       !---------------------------------- 
    276299      ! Change thickness to volume 
    277       !---------------------------------- 
    278300      v_i(:,:,:)   = ht_i(:,:,:) * a_i(:,:,:) 
    279301      v_s(:,:,:)   = ht_s(:,:,:) * a_i(:,:,:) 
     
    292314      CALL lim_var_zapsmall 
    293315 
    294       !-------------------------------------------- 
    295       ! Diagnostic thermodynamic growth rates 
    296       !-------------------------------------------- 
    297       IF( ln_icectl )   CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermodyn. - ' )   ! control print 
    298  
    299       IF(ln_ctl) THEN            ! Control print 
    300          CALL prt_ctl_info(' ') 
    301          CALL prt_ctl_info(' - Cell values : ') 
    302          CALL prt_ctl_info('   ~~~~~~~~~~~~~ ') 
    303          CALL prt_ctl(tab2d_1=e1e2t, clinfo1=' lim_thd  : cell area :') 
    304          CALL prt_ctl(tab2d_1=at_i , clinfo1=' lim_thd  : at_i      :') 
    305          CALL prt_ctl(tab2d_1=vt_i , clinfo1=' lim_thd  : vt_i      :') 
    306          CALL prt_ctl(tab2d_1=vt_s , clinfo1=' lim_thd  : vt_s      :') 
    307          DO jl = 1, jpl 
    308             CALL prt_ctl_info(' ') 
    309             CALL prt_ctl_info(' - Category : ', ivar1=jl) 
    310             CALL prt_ctl_info('   ~~~~~~~~~~') 
    311             CALL prt_ctl(tab2d_1=a_i   (:,:,jl)   , clinfo1= ' lim_thd  : a_i      : ') 
    312             CALL prt_ctl(tab2d_1=ht_i  (:,:,jl)   , clinfo1= ' lim_thd  : ht_i     : ') 
    313             CALL prt_ctl(tab2d_1=ht_s  (:,:,jl)   , clinfo1= ' lim_thd  : ht_s     : ') 
    314             CALL prt_ctl(tab2d_1=v_i   (:,:,jl)   , clinfo1= ' lim_thd  : v_i      : ') 
    315             CALL prt_ctl(tab2d_1=v_s   (:,:,jl)   , clinfo1= ' lim_thd  : v_s      : ') 
    316             CALL prt_ctl(tab2d_1=e_s   (:,:,1,jl) , clinfo1= ' lim_thd  : e_s      : ') 
    317             CALL prt_ctl(tab2d_1=t_su  (:,:,jl)   , clinfo1= ' lim_thd  : t_su     : ') 
    318             CALL prt_ctl(tab2d_1=t_s   (:,:,1,jl) , clinfo1= ' lim_thd  : t_snow   : ') 
    319             CALL prt_ctl(tab2d_1=sm_i  (:,:,jl)   , clinfo1= ' lim_thd  : sm_i     : ') 
    320             CALL prt_ctl(tab2d_1=smv_i (:,:,jl)   , clinfo1= ' lim_thd  : smv_i    : ') 
    321             DO jk = 1, nlay_i 
    322                CALL prt_ctl_info(' ') 
    323                CALL prt_ctl_info(' - Layer : ', ivar1=jk) 
    324                CALL prt_ctl_info('   ~~~~~~~') 
    325                CALL prt_ctl(tab2d_1=t_i(:,:,jk,jl) , clinfo1= ' lim_thd  : t_i      : ') 
    326                CALL prt_ctl(tab2d_1=e_i(:,:,jk,jl) , clinfo1= ' lim_thd  : e_i      : ') 
    327             END DO 
    328          END DO 
    329       ENDIF 
    330       ! 
    331       ! 
    332       IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    333  
    334       !------------------------------------------------------------------------------| 
    335       !  6) Transport of ice between thickness categories.                           | 
    336       !------------------------------------------------------------------------------| 
     316      ! control checks 
     317      IF( ln_limctl )    CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermodyn. - ' )   ! control print 
     318      ! 
     319      IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     320 
     321      !------------------------------------------------! 
     322      !  Transport ice between thickness categories 
     323      !------------------------------------------------! 
    337324      ! Given thermodynamic growth rates, transport ice between thickness categories. 
    338       IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     325      IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    339326 
    340327      IF( jpl > 1 )      CALL lim_itd_th_rem( 1, jpl, kt ) 
    341328 
    342       IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    343  
    344       !------------------------------------------------------------------------------| 
    345       !  7) Add frazil ice growing in leads. 
    346       !------------------------------------------------------------------------------| 
    347       IF( ln_limdiahsb ) CALL lim_cons_hsm(0, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    348  
    349       CALL lim_thd_lac 
     329      IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     330 
     331      !------------------------------------------------! 
     332      !  Add frazil ice growing in leads 
     333      !------------------------------------------------! 
     334      IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     335 
     336      IF( ln_limdO )     CALL lim_thd_lac 
    350337       
    351       IF( ln_limdiahsb ) CALL lim_cons_hsm(1, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
     338      ! conservation test 
     339      IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    352340 
    353341      ! Control print 
    354       IF(ln_ctl) THEN 
    355          CALL lim_var_glo2eqv 
    356  
    357          CALL prt_ctl_info(' ') 
    358          CALL prt_ctl_info(' - Cell values : ') 
    359          CALL prt_ctl_info('   ~~~~~~~~~~~~~ ') 
    360          CALL prt_ctl(tab2d_1=e1e2t, clinfo1=' lim_itd_th  : cell area :') 
    361          CALL prt_ctl(tab2d_1=at_i , clinfo1=' lim_itd_th  : at_i      :') 
    362          CALL prt_ctl(tab2d_1=vt_i , clinfo1=' lim_itd_th  : vt_i      :') 
    363          CALL prt_ctl(tab2d_1=vt_s , clinfo1=' lim_itd_th  : vt_s      :') 
    364          DO jl = 1, jpl 
    365             CALL prt_ctl_info(' ') 
    366             CALL prt_ctl_info(' - Category : ', ivar1=jl) 
    367             CALL prt_ctl_info('   ~~~~~~~~~~') 
    368             CALL prt_ctl(tab2d_1=a_i   (:,:,jl)   , clinfo1= ' lim_itd_th  : a_i      : ') 
    369             CALL prt_ctl(tab2d_1=ht_i  (:,:,jl)   , clinfo1= ' lim_itd_th  : ht_i     : ') 
    370             CALL prt_ctl(tab2d_1=ht_s  (:,:,jl)   , clinfo1= ' lim_itd_th  : ht_s     : ') 
    371             CALL prt_ctl(tab2d_1=v_i   (:,:,jl)   , clinfo1= ' lim_itd_th  : v_i      : ') 
    372             CALL prt_ctl(tab2d_1=v_s   (:,:,jl)   , clinfo1= ' lim_itd_th  : v_s      : ') 
    373             CALL prt_ctl(tab2d_1=e_s   (:,:,1,jl) , clinfo1= ' lim_itd_th  : e_s      : ') 
    374             CALL prt_ctl(tab2d_1=t_su  (:,:,jl)   , clinfo1= ' lim_itd_th  : t_su     : ') 
    375             CALL prt_ctl(tab2d_1=t_s   (:,:,1,jl) , clinfo1= ' lim_itd_th  : t_snow   : ') 
    376             CALL prt_ctl(tab2d_1=sm_i  (:,:,jl)   , clinfo1= ' lim_itd_th  : sm_i     : ') 
    377             CALL prt_ctl(tab2d_1=smv_i (:,:,jl)   , clinfo1= ' lim_itd_th  : smv_i    : ') 
    378             DO jk = 1, nlay_i 
    379                CALL prt_ctl_info(' ') 
    380                CALL prt_ctl_info(' - Layer : ', ivar1=jk) 
    381                CALL prt_ctl_info('   ~~~~~~~') 
    382                CALL prt_ctl(tab2d_1=t_i(:,:,jk,jl) , clinfo1= ' lim_itd_th  : t_i      : ') 
    383                CALL prt_ctl(tab2d_1=e_i(:,:,jk,jl) , clinfo1= ' lim_itd_th  : e_i      : ') 
    384             END DO 
    385          END DO 
    386       ENDIF 
    387       ! 
    388       IF( nn_timing == 1 )   CALL timing_stop('limthd') 
    389       ! 
     342      IF( ln_ctl )       CALL lim_prt3D( 'limthd' ) 
     343      ! 
     344      CALL wrk_dealloc( jpi,jpj, zu_io, zv_io, zfric ) 
     345      ! 
     346      IF( nn_timing == 1 )  CALL timing_stop('limthd') 
     347 
    390348   END SUBROUTINE lim_thd  
    391349 
     
    449407            zda_mel     = rswitch * a_i_1d(ji) * zdh_mel / ( 2._wp * MAX( zhi_bef, epsi20 ) ) 
    450408            a_i_1d(ji)  = MAX( epsi20, a_i_1d(ji) + zda_mel )  
    451              ! adjust thickness 
     409            ! adjust thickness 
    452410            ht_i_1d(ji) = zvi / a_i_1d(ji)             
    453411            ht_s_1d(ji) = zvs / a_i_1d(ji)             
     
    613571      !!------------------------------------------------------------------- 
    614572      INTEGER  ::   ios                 ! Local integer output status for namelist read 
    615       !! 
    616       NAMELIST/namicethd/ rn_hnewice, ln_frazil, rn_maxfrazb, rn_vfrazb, rn_Cfrazb,                & 
    617          &                rn_himin, rn_betas, rn_kappa_i, nn_conv_dif, rn_terr_dif, nn_ice_thcon,  & 
    618          &                nn_monocat, ln_it_qnsice 
     573      NAMELIST/namicethd/ rn_kappa_i, nn_conv_dif, rn_terr_dif, nn_ice_thcon,ln_it_qnsice,nn_monocat,  & 
     574         &                ln_limdH, rn_betas,                                                          & 
     575         &                ln_limdA, rn_beta, rn_dmin,                                                  & 
     576         &                ln_limdO, rn_hnewice, ln_frazil, rn_maxfrazb, rn_vfrazb, rn_Cfrazb, rn_himin 
    619577      !!------------------------------------------------------------------- 
    620       ! 
    621       IF(lwp) THEN 
    622          WRITE(numout,*) 
    623          WRITE(numout,*) 'lim_thd_init : Ice Thermodynamics initialization' 
    624          WRITE(numout,*) '~~~~~~~~~~~~' 
    625       ENDIF 
    626578      ! 
    627579      REWIND( numnam_ice_ref )              ! Namelist namicethd in reference namelist : Ice thermodynamics 
     
    635587      ! 
    636588      IF(lwp) THEN                          ! control print 
    637          WRITE(numout,*)'   Namelist of ice parameters for ice thermodynamic computation ' 
     589         WRITE(numout,*) 'lim_thd_init : Ice Thermodynamics' 
     590         WRITE(numout,*) '~~~~~~~~~~~~~' 
     591         WRITE(numout,*)'   -- limthd_dif --' 
     592         WRITE(numout,*)'      extinction radiation parameter in sea ice               rn_kappa_i   = ', rn_kappa_i 
     593         WRITE(numout,*)'      maximal n. of iter. for heat diffusion computation      nn_conv_dif  = ', nn_conv_dif 
     594         WRITE(numout,*)'      maximal err. on T for heat diffusion computation        rn_terr_dif  = ', rn_terr_dif 
     595         WRITE(numout,*)'      switch for comp. of thermal conductivity in the ice     nn_ice_thcon = ', nn_ice_thcon 
     596         WRITE(numout,*)'      iterate the surface non-solar flux (T) or not (F)       ln_it_qnsice = ', ln_it_qnsice 
     597         WRITE(numout,*)'      virtual ITD mono-category parameterizations (1) or not  nn_monocat   = ', nn_monocat 
     598         WRITE(numout,*)'   -- limthd_dh --' 
     599         WRITE(numout,*)'      activate ice thick change from top/bot (T) or not (F)   ln_limdH     = ', ln_limdH 
     600         WRITE(numout,*)'      coefficient for ice-lead partition of snowfall          rn_betas     = ', rn_betas 
     601         WRITE(numout,*)'   -- limthd_da --' 
     602         WRITE(numout,*)'      activate lateral melting (T) or not (F)                 ln_limdA     = ', ln_limdA 
     603         WRITE(numout,*)'      Coef. beta for lateral melting param.                   rn_beta      = ', rn_beta 
     604         WRITE(numout,*)'      Minimum floe diameter for lateral melting param.        rn_dmin      = ', rn_dmin 
     605         WRITE(numout,*)'   -- limthd_lac --' 
     606         WRITE(numout,*)'      activate ice growth in open-water (T) or not (F)        ln_limdO     = ', ln_limdO 
    638607         WRITE(numout,*)'      ice thick. for lateral accretion                        rn_hnewice   = ', rn_hnewice 
    639608         WRITE(numout,*)'      Frazil ice thickness as a function of wind or not       ln_frazil    = ', ln_frazil 
     
    641610         WRITE(numout,*)'      Thresold relative drift speed for collection of frazil  rn_vfrazb    = ', rn_vfrazb 
    642611         WRITE(numout,*)'      Squeezing coefficient for collection of frazil          rn_Cfrazb    = ', rn_Cfrazb 
     612         WRITE(numout,*)'   -- limitd_th --' 
    643613         WRITE(numout,*)'      minimum ice thickness                                   rn_himin     = ', rn_himin  
    644          WRITE(numout,*)'      numerical carac. of the scheme for diffusion in ice ' 
    645          WRITE(numout,*)'      coefficient for ice-lead partition of snowfall          rn_betas     = ', rn_betas 
    646          WRITE(numout,*)'      extinction radiation parameter in sea ice               rn_kappa_i   = ', rn_kappa_i 
    647          WRITE(numout,*)'      maximal n. of iter. for heat diffusion computation      nn_conv_dif  = ', nn_conv_dif 
    648          WRITE(numout,*)'      maximal err. on T for heat diffusion computation        rn_terr_dif  = ', rn_terr_dif 
    649          WRITE(numout,*)'      switch for comp. of thermal conductivity in the ice     nn_ice_thcon = ', nn_ice_thcon 
    650614         WRITE(numout,*)'      check heat conservation in the ice/snow                 con_i        = ', con_i 
    651          WRITE(numout,*)'      virtual ITD mono-category parameterizations (1) or not  nn_monocat   = ', nn_monocat 
    652          WRITE(numout,*)'      iterate the surface non-solar flux (T) or not (F)       ln_it_qnsice = ', ln_it_qnsice 
    653615      ENDIF 
    654616      IF( jpl > 1 .AND. nn_monocat == 1 ) THEN  
Note: See TracChangeset for help on using the changeset viewer.