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

Ignore:
Timestamp:
2017-09-05T19:53:41+02:00 (7 years ago)
Author:
clem
Message:

changes in style - part2 -

File:
1 edited

Legend:

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

    r8496 r8498  
    157157      !!------------------------------------------------------------------ 
    158158      INTEGER  ::   ji, jj, jk, jl   ! dummy loop indices 
    159       REAL(wp) ::   ze_i, z1_CpR, zdiscrim, zaaa, z1_2aaa             ! local scalars 
    160       REAL(wp) ::   ze_s, zL_Cp , ztmelts , zbbb, zccc                !   -      - 
    161       REAL(wp) ::   zhmax, z1_zhmax, zsm_i, zcpMcpic, zt_i   !   -      - 
    162       REAL(wp) ::   zlay_i, zlay_s   !   -      - 
     159      REAL(wp) ::   ze_i, z1_cp, z1_2cp             ! local scalars 
     160      REAL(wp) ::   ze_s, ztmelts, zbbb, zccc       !   -      - 
     161      REAL(wp) ::   zhmax, z1_zhmax, zsm_i          !   -      - 
     162      REAL(wp) ::   zlay_i, zlay_s                  !   -      - 
    163163      REAL(wp), DIMENSION(jpi,jpj,jpl) ::   z1_a_i, z1_v_i 
    164164      !!------------------------------------------------------------------ 
     
    181181      ht_i(:,:,:) = v_i (:,:,:) * z1_a_i(:,:,:)    !--- ice thickness 
    182182 
    183      zhmax    =          hi_max(jpl) 
     183      zhmax    =          hi_max(jpl) 
    184184      z1_zhmax =  1._wp / hi_max(jpl)                
    185185      WHERE( ht_i(:,:,jpl) > zhmax )               !--- bound ht_i by hi_max (i.e. 99 m) with associated update of ice area 
     
    205205      !------------------- 
    206206      zlay_i   = REAL( nlay_i , wp )    ! number of layers 
    207       zaaa     = cpic                   ! Conversion q(S,T) -> T (second order equation) 
    208       z1_2aaa  = 1._wp / ( 2._wp * zaaa ) 
    209       zcpMcpic = rcp - cpic 
     207      z1_2cp  = 1._wp / ( 2._wp * cpic ) 
    210208      DO jl = 1, jpl 
    211209         DO jk = 1, nlay_i 
     
    214212                  IF ( v_i(ji,jj,jl) > epsi20 ) THEN     !--- icy area  
    215213                     ! 
    216                      ze_i    =   e_i(ji,jj,jk,jl) / v_i(ji,jj,jl) * zlay_i               ! Energy of melting e(S,T) [J.m-3] 
    217                      ztmelts = - s_i(ji,jj,jk,jl) * tmut                                 ! Ice layer melt temperature [C] 
    218                      ! 
    219                      zbbb     =  zcpMcpic * ztmelts + ze_i * r1_rhoic - lfus 
    220                      zccc     =  lfus * ztmelts 
    221                      zdiscrim =  SQRT(  MAX( zbbb*zbbb - 4._wp*zaaa*zccc , 0._wp)  ) 
    222                      zt_i     = - ( zbbb + zdiscrim ) * z1_2aaa                          ! [C] 
    223                      t_i(ji,jj,jk,jl) = MAX( -100._wp , MIN( zt_i , ztmelts )  ) + rt0   ! [K] with bounds: -100 < zt_i < ztmelts 
     214                     ze_i             =   e_i(ji,jj,jk,jl) / v_i(ji,jj,jl) * zlay_i               ! Energy of melting e(S,T) [J.m-3] 
     215                     ztmelts          = - s_i(ji,jj,jk,jl) * tmut                                 ! Ice layer melt temperature [C] 
     216                     ! Conversion q(S,T) -> T (second order equation) 
     217                     zbbb             = ( rcp - cpic ) * ztmelts + ze_i * r1_rhoic - lfus 
     218                     zccc             = SQRT( MAX( zbbb * zbbb - 4._wp * cpic * lfus * ztmelts , 0._wp) ) 
     219                     t_i(ji,jj,jk,jl) = MAX( -100._wp , MIN( -( zbbb + zccc ) * z1_2cp , ztmelts ) ) + rt0   ! [K] with bounds: -100 < t_i < ztmelts 
    224220                     ! 
    225221                  ELSE                                !--- no ice 
    226                      t_i(ji,jj,jk,jl) =  rt0 - 100._wp                                   ! impose 173.15 K (i.e. -100 C) 
     222                     t_i(ji,jj,jk,jl) = rt0 - 100._wp                                   ! impose 173.15 K (i.e. -100 C) 
    227223!!clem: I think we should impose rt0 instead 
    228224                  ENDIF 
     
    235231      ! Snow temperature   [K]   (with a minimum value (rt0 - 100.) imposed everywhere) 
    236232      !-------------------- 
    237       z1_CpR = 1._wp / ( cpic * rhosn ) 
    238       zL_Cp  = lfus  /   cpic 
    239233      zlay_s = REAL( nlay_s , wp ) 
     234      z1_cp  = 1._wp / cpic 
    240235      DO jk = 1, nlay_s 
    241236         WHERE( v_s(:,:,:) > epsi20 )        !--- icy area 
    242             t_s(:,:,jk,:) = MAX(  -100._wp , MIN( - z1_CpR * ( e_s(:,:,jk,:)/v_s(:,:,:)*zlay_s ) + zL_Cp , 0._wp )  ) + rt0        
     237            t_s(:,:,jk,:) = MAX( -100._wp , MIN( z1_cp * ( -r1_rhosn * (e_s(:,:,jk,:)/v_s(:,:,:)*zlay_s) + lfus ) , 0._wp ) ) + rt0 
    243238         ELSEWHERE                           !--- no ice 
    244239            t_s(:,:,jk,:) =  rt0 - 100._wp         ! impose 173.15 K (i.e. -100 C) 
    245240         END WHERE 
    246241      END DO 
    247 !!gm perhaps better like this (?) :  (if this compile, yes! one test instead of nlay_s tests) 
    248 !      WHERE( v_s(:,:,:) > epsi20 )        !--- icy area 
    249 !         DO jk = 1, nlay_s 
    250 !            t_s(:,:,jk,:) = MAX(  -100._wp , MIN( - z1_CpR * ( e_s(:,:,jk,:)/v_s(:,:,:)*zlay_s ) + zL_Cp , 0._wp )  ) + rt0        
    251 !         END DO 
    252 !      ELSEWHERE                           !--- no ice 
    253 !         DO jk = 1, nlay_s 
    254 !            t_s(:,:,jk,:) =  rt0 - 100._wp      ! impose 173.15 K (i.e. -100 C) 
    255 !         END DO 
    256 !      END WHERE 
    257 !!gm end better ? 
    258242 
    259243      ! integrated values  
Note: See TracChangeset for help on using the changeset viewer.