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 8379 for branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90 – NEMO

Ignore:
Timestamp:
2017-07-26T17:35:49+02:00 (7 years ago)
Author:
clem
Message:

more cleaning

File:
1 edited

Legend:

Unmodified
Added
Removed
  • branches/2017/dev_r8183_ICEMODEL/NEMOGCM/NEMO/LIM_SRC_3/limthd.F90

    r8378 r8379  
    8383      INTEGER  :: ji, jj, jk, jl   ! dummy loop indices 
    8484      REAL(wp) :: zfric_u, zqld, zqfr, zqfr_neg 
    85       REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b  
     85      REAL(wp) :: zvi_b, zsmv_b, zei_b, zfs_b, zfw_b, zft_b ! conservation check 
    8686      REAL(wp), PARAMETER :: zfric_umin = 0._wp           ! lower bound for the friction velocity (cice value=5.e-04) 
    8787      REAL(wp), PARAMETER :: zch        = 0.0057_wp       ! heat transfer coefficient 
     
    218218         END DO 
    219219 
    220          ! debug point to follow 
    221          jiindex_1d = 0 
    222          IF( ln_limctl ) THEN 
    223             DO ji = mi0(iiceprt), mi1(iiceprt) 
    224                DO jj = mj0(jiceprt), mj1(jiceprt) 
    225                   jiindex_1d = (jj - 1) * jpi + ji 
    226                   WRITE(numout,*) ' lim_thd : Category no : ', jl  
    227                END DO 
    228             END DO 
    229          ENDIF 
    230  
    231220         IF( lk_mpp )         CALL mpp_ini_ice( nidx , numout ) 
    232221 
     
    234223            !                                                                 
    235224                              CALL lim_thd_1d2d( jl, 1 )            ! --- Move to 1D arrays --- ! 
    236             ! 
    237             DO jk = 1, nlay_i                                       ! --- Change units from J/m2 to J/m3 --- ! 
    238                WHERE( ht_i_1d(1:nidx)>0._wp ) e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) / (ht_i_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_i 
    239             ENDDO 
    240             DO jk = 1, nlay_s 
    241                WHERE( ht_s_1d(1:nidx)>0._wp ) e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) / (ht_s_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_s 
    242             ENDDO 
     225            !                                                       ! --- & Change units of e_i, e_s from J/m2 to J/m3 --- ! 
    243226            ! 
    244227            s_i_new   (1:nidx) = 0._wp ; dh_s_tot (1:nidx) = 0._wp  ! --- some init --- !  (important to have them here)  
     
    264247            IF( ln_limdA )    CALL lim_thd_da                       ! --- lateral melting --- ! 
    265248            ! 
    266             DO jk = 1, nlay_i                                       ! --- Change units from J/m3 to J/m2 --- ! 
    267                e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) * ht_i_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_i 
    268             ENDDO 
    269             DO jk = 1, nlay_s 
    270                e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) * ht_s_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_s 
    271             ENDDO 
    272             ! 
    273             ! Change thickness to volume 
    274             v_i_1d(1:nidx)   = ht_i_1d(1:nidx) * a_i_1d(1:nidx) 
    275             v_s_1d(1:nidx)   = ht_s_1d(1:nidx) * a_i_1d(1:nidx) 
    276             smv_i_1d(1:nidx) = sm_i_1d(1:nidx) * v_i_1d(1:nidx) 
    277             ! 
    278                               CALL lim_thd_1d2d( jl, 2 )            ! --- Move to 2D arrays --- ! 
     249                              CALL lim_thd_1d2d( jl, 2 )            ! --- Change units of e_i, e_s from J/m3 to J/m2 --- ! 
     250            !                                                       ! --- & Move to 2D arrays --- ! 
    279251            ! 
    280252            IF( lk_mpp )      CALL mpp_comm_free( ncomm_ice ) !RB necessary ?? 
     
    282254         ! 
    283255      END DO 
    284       at_i(:,:) = SUM( a_i(:,:,:), dim=3 ) 
    285  
    286256      ! update ice age (in case a_i changed, i.e. becomes 0 or lateral melting) 
    287       DO jl  = 1, jpl 
    288          DO jj = 1, jpj 
    289             DO ji = 1, jpi 
    290                rswitch = MAX( 0._wp , SIGN( 1._wp, a_i_b(ji,jj,jl) - epsi10 ) ) 
    291                oa_i(ji,jj,jl) = rswitch * oa_i(ji,jj,jl) * a_i(ji,jj,jl) / MAX( a_i_b(ji,jj,jl), epsi10 ) 
    292             END DO 
    293          END DO 
    294       END DO 
     257      oa_i(:,:,:) = o_i(:,:,:) * a_i(:,:,:) 
    295258 
    296259      IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'limthd', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    297  
    298       CALL lim_var_zapsmall 
    299  
     260      ! 
     261                         CALL lim_var_zapsmall           ! --- remove very small ice concentration (<1e-10) --- ! 
     262      !                                                  !     & make sure at_i=SUM(a_i) & ato_i=1 where at_i=0 
     263      !                    
     264      IF( jpl > 1 )      CALL lim_itd_th_rem( kt )       ! --- Transport ice between thickness categories --- ! 
     265      ! 
     266      IF( ln_limdO )     CALL lim_thd_lac                ! --- frazil ice growing in leads --- ! 
     267      ! 
    300268      IF( ln_limctl )    CALL lim_prt( kt, iiceprt, jiceprt, 1, ' - ice thermodyn. - ' )   ! control print 
    301       ! 
    302       !------------------------------------------------! 
    303       !  Transport ice between thickness categories 
    304       !------------------------------------------------! 
    305       IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    306  
    307       IF( jpl > 1 )      CALL lim_itd_th_rem( kt ) 
    308  
    309       IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'limitd_th_rem', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    310  
    311       !------------------------------------------------! 
    312       !  Add frazil ice growing in leads 
    313       !------------------------------------------------! 
    314       IF( ln_limdiachk ) CALL lim_cons_hsm(0, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    315  
    316       IF( ln_limdO )     CALL lim_thd_lac 
    317        
    318       IF( ln_limdiachk ) CALL lim_cons_hsm(1, 'limthd_lac', zvi_b, zsmv_b, zei_b, zfw_b, zfs_b, zft_b) 
    319  
    320       IF( ln_ctl )       CALL lim_prt3D( 'limthd' )  ! Control print 
     269      IF( ln_ctl )       CALL lim_prt3D( 'limthd' )      ! Control print 
    321270      ! 
    322271      IF( nn_timing == 1 )  CALL timing_stop('limthd') 
     
    480429         CALL tab_2d_1d( nidx, idxice(1:nidx), sst_1d(1:nidx), sst_m ) 
    481430         CALL tab_2d_1d( nidx, idxice(1:nidx), sss_1d(1:nidx), sss_m ) 
     431 
     432         ! --- Change units of e_i, e_s from J/m2 to J/m3 --- ! 
     433         DO jk = 1, nlay_i 
     434            WHERE( ht_i_1d(1:nidx)>0._wp ) e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) / (ht_i_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_i 
     435         ENDDO 
     436         DO jk = 1, nlay_s 
     437            WHERE( ht_s_1d(1:nidx)>0._wp ) e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) / (ht_s_1d(1:nidx) * a_i_1d(1:nidx)) * nlay_s 
     438         ENDDO 
    482439         ! 
    483440      CASE( 2 )            ! from 1D to 2D 
    484441         ! 
     442         ! --- Change units of e_i, e_s from J/m3 to J/m2 --- ! 
     443         DO jk = 1, nlay_i 
     444            e_i_1d(1:nidx,jk) = e_i_1d(1:nidx,jk) * ht_i_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_i 
     445         ENDDO 
     446         DO jk = 1, nlay_s 
     447            e_s_1d(1:nidx,jk) = e_s_1d(1:nidx,jk) * ht_s_1d(1:nidx) * a_i_1d(1:nidx) * r1_nlay_s 
     448         ENDDO 
     449         ! 
     450         ! Change thickness to volume 
     451         v_i_1d(1:nidx)   = ht_i_1d(1:nidx) * a_i_1d(1:nidx) 
     452         v_s_1d(1:nidx)   = ht_s_1d(1:nidx) * a_i_1d(1:nidx) 
     453         smv_i_1d(1:nidx) = sm_i_1d(1:nidx) * v_i_1d(1:nidx) 
     454          
    485455         CALL tab_1d_2d( nidx, idxice(1:nidx), at_i_1d(1:nidx), at_i             ) 
    486456         CALL tab_1d_2d( nidx, idxice(1:nidx), a_i_1d (1:nidx), a_i(:,:,jl)      ) 
Note: See TracChangeset for help on using the changeset viewer.