Changeset 4818


Ignore:
Timestamp:
Feb 14, 2024, 4:44:16 PM (3 months ago)
Author:
evignon
Message:

modifications suite à l'atelier nuages du jour
le nouveau schema de precip est full operationnel!
(+ petite commentaire d'un print dans wake)

Location:
LMDZ6/trunk/libf
Files:
7 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk/libf/phylmd/lmdz_lscp.F90

    r4803 r4818  
    946946                    zqn(i) = zq(i)
    947947                    rneb(i,k) = 1.0
    948                     zcond(i) = MAX(0.0,zqn(i)-gammasat(i)*zqs(i))
     948                    zcond(i) = MAX(0.0,zqn(i)-zqs(i))
    949949                    rhcl(i,k)=1.0
    950950                ELSE
    951                     zcond(i) = MAX(0.0,zqn(i)-gammasat(i)*zqs(i))*rneb(i,k)
     951                    zcond(i) = MAX(0.0,zqn(i)-zqs(i))*rneb(i,k)
    952952                    ! following line is very strange and probably wrong:
    953953                    rhcl(i,k)=(zqs(i)+zq(i))/2./zqs(i)

  • LMDZ6/trunk/libf/phylmd/lmdz_lscp_ini.F90

    r4809 r4818  
    66  !--------------------
    77 
    8   REAL RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG
    9   !$OMP THREADPRIVATE(RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG)
     8  REAL RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG, RPI
     9  !$OMP THREADPRIVATE(RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG, RPI)
    1010
    1111  REAL, SAVE, PROTECTED :: seuil_neb=0.001      ! cloud fraction threshold: a cloud really exists when exceeded
    … …  
    164164  !$OMP THREADPRIVATE(r_snow)
    165165
    166   REAL, SAVE, PROTECTED :: Eff_rain_liq=1.0                  ! A COMMENTER TODO [-]
    167   !$OMP THREADPRIVATE(Eff_rain_liq)
    168 
    169   REAL, SAVE, PROTECTED :: Eff_snow_ice=0.5                ! A COMMENTER TODO [-]
    170   !$OMP THREADPRIVATE(Eff_snow_ice)
    171 
    172   REAL, SAVE, PROTECTED :: Eff_snow_liq=1.0              ! A COMMENTER TODO [-]
    173   !$OMP THREADPRIVATE(Eff_snow_liq)
    174 
    175166  REAL, SAVE, PROTECTED :: tau_auto_snow_min=1800.          ! A COMMENTER TODO [s]
    176167  !$OMP THREADPRIVATE(tau_auto_snow_min)
    … …  
    181172  REAL, SAVE, PROTECTED :: eps=1.E-10                       ! A COMMENTER TODO [-]
    182173  !$OMP THREADPRIVATE(eps)
     174
     175  REAL, SAVE, PROTECTED :: air_thermal_conduct=2.4e-2      ! A COMMENTER TODO [-]
     176  !$OMP THREADPRIVATE(air_thermal_conduct)
     177
     178  REAL, SAVE, PROTECTED :: coef_ventil=1.                   ! A COMMENTER TODO [-]
     179  !$OMP THREADPRIVATE(coef_ventil)
     180
     181  REAL, SAVE, PROTECTED :: alpha_freez=4.                 ! A COMMENTER TODO [-]
     182  !$OMP THREADPRIVATE(alpha_freez)
     183
     184  REAL, SAVE, PROTECTED :: gamma_freez=0.1               ! A COMMENTER TODO [m-3.s-1]
     185  !$OMP THREADPRIVATE(gamma_freez)
    183186  !--End of the parameters for poprecip
    184187
    … …  
    190193
    191194SUBROUTINE lscp_ini(dtime,lunout_in,prt_level_in,ok_ice_sursat, iflag_ratqs, fl_cor_ebil_in, RCPD_in, RLSTT_in, RLVTT_in, RLMLT_in, &
    192                     RVTMP2_in, RTT_in,RD_in,RG_in)
     195                    RVTMP2_in, RTT_in,RD_in,RG_in,RPI_in)
    193196
    194197
    … …  
    202205
    203206   REAL, INTENT(IN)      :: RCPD_in, RLSTT_in, RLVTT_in, RLMLT_in
    204    REAL, INTENT(IN)      ::  RVTMP2_in, RTT_in, RD_in, RG_in
     207   REAL, INTENT(IN)      ::  RVTMP2_in, RTT_in, RD_in, RG_in, RPI_in
    205208   character (len=20) :: modname='lscp_ini_mod'
    206209   character (len=80) :: abort_message
    … …  
    220223    RG=RG_in
    221224    RVTMP2=RVTMP2_in
     225    RPI=RPI_in
    222226
    223227
    … …  
    262266    CALL getin_p('gamma_agg',gamma_agg)
    263267    CALL getin_p('gamma_col',gamma_col)
     268    CALL getin_p('gamma_rim',gamma_rim)
    264269
    265270
    … …  
    303308    WRITE(lunout,*) 'lscp_ini, gamma_agg:', gamma_agg
    304309    WRITE(lunout,*) 'lscp_ini, gamma_col:', gamma_col
     310    WRITE(lunout,*) 'lscp_ini, gamma_rim:', gamma_rim
    305311
    306312

  • LMDZ6/trunk/libf/phylmd/lmdz_lscp_poprecip.F90

    r4809 r4818  
    232232           precipfracclr, precipfraccld, &
    233233           rain, rainclr, raincld, snow, snowclr, snowcld, &
    234            dqrauto,dqrcol,dqrmelt,dqrfreez,dqsauto,dqsagg,dqsrim,dqsmelt,dqsfreez)
     234           dqrauto, dqrcol, dqrmelt, dqrfreez, dqsauto, dqsagg, &
     235           dqsrim, dqsmelt, dqsfreez)
    235236
    236237USE lmdz_lscp_ini, ONLY : prt_level, lunout
    237 USE lmdz_lscp_ini, ONLY : RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG
     238USE lmdz_lscp_ini, ONLY : RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT, RD, RG, RPI
    238239USE lmdz_lscp_tools, ONLY : calc_qsat_ecmwf
    239240
    … …  
    241242                          cld_lc_lsc, cld_tau_lsc, cld_expo_lsc, rain_int_min, &
    242243                          thresh_precip_frac, gamma_col, gamma_agg, gamma_rim, &
    243                           rho_rain, rho_snow, r_rain, r_snow, Eff_rain_liq,    &
    244                           Eff_snow_ice, Eff_snow_liq, tau_auto_snow_min,       &
    245                           tau_auto_snow_max, thresh_precip_frac, eps,          &
     244                          rho_rain, rho_snow, r_rain, r_snow,                  &
     245                          tau_auto_snow_min, tau_auto_snow_max,                &
     246                          thresh_precip_frac, eps, air_thermal_conduct,        &
     247                          coef_ventil, alpha_freez, gamma_freez, temp_nowater, &
    246248                          iflag_cloudth_vert, iflag_rain_incloud_vol
     249
    247250
    248251IMPLICIT NONE
    … …  
    292295
    293296INTEGER :: i
    294 
    295 REAL :: hum_to_flux
     297REAL, DIMENSION(klon) :: hum_to_flux
     298REAL, DIMENSION(klon) :: qtot !--includes vap, liq, ice and precip
     299
     300! partition of the fluxes
    296301REAL :: dcldfra
    297302REAL :: precipfractot
    … …  
    299304REAL :: drainclr, dsnowclr
    300305REAL :: draincld, dsnowcld
     306
     307! collection, aggregation and riming
    301308REAL :: eff_cldfra
    302 REAL :: coef_col, coef_agg, coef_tmp, qrain
     309REAL :: coef_col, coef_agg, coef_rim, coef_tmp, qrain
     310REAL :: Eff_rain_liq, Eff_snow_ice, Eff_snow_liq
     311REAL :: dqlcol           ! loss of liquid cloud content due to collection by rain [kg/kg/s]
     312REAL :: dqiagg           ! loss of ice cloud content due to collection by aggregation [kg/kg/s]
     313REAL :: dqlrim           ! loss of liquid cloud content due to riming on snow[kg/kg/s]
     314
     315! autoconversion
    303316REAL :: qthresh_auto_rain, tau_auto_rain, expo_auto_rain
    304317REAL :: qthresh_auto_snow, tau_auto_snow, expo_auto_snow
    305 REAL ::  dqlcol           ! loss of liquid cloud content due to collection by rain [kg/kg/s]
    306 REAL ::  dqiagg           ! loss of ice cloud content due to collection by aggregation [kg/kg/s]
    307 REAL ::  dqlauto          ! loss of liquid cloud content due to autoconversion to rain [kg/kg/s]
    308 REAL ::  dqiauto          ! loss of ice cloud content due to autoconversion to snow [kg/kg/s]
    309 REAL ::  dqlrim           ! loss of liquid cloud content due to riming on snow[kg/kg/s]
     318REAL :: dqlauto          ! loss of liquid cloud content due to autoconversion to rain [kg/kg/s]
     319REAL :: dqiauto          ! loss of ice cloud content due to autoconversion to snow [kg/kg/s]
     320
     321! melting
     322REAL :: dqsmelt_max
     323REAL :: fallice_clr, fallice_cld
     324REAL :: nb_snowflake_clr, nb_snowflake_cld
     325REAL :: capa_snowflake, temp_wetbulb
     326REAL :: dqsclrmelt, dqscldmelt, dqstotmelt
     327REAL, DIMENSION(klon) :: qzero, qsat, dqsat
     328
     329! freezing
     330REAL :: dqrfreez_max
     331REAL :: tau_freez
     332REAL :: dqrclrfreez, dqrcldfreez, dqrtotfreez
    310333
    311334
    312335!--Initialisation of variables
    313336
    314 dqrcol(:) = 0.
    315 dqsagg(:) = 0.
    316 dqsauto(:) = 0.
    317 dqrauto(:) = 0.
    318 dqsrim(:) = 0.
    319 dqrmelt(:) = 0.
    320 dqsmelt(:) = 0.
     337qzero(:)    = 0.
     338
     339dqrcol(:)   = 0.
     340dqsagg(:)   = 0.
     341dqsauto(:)  = 0.
     342dqrauto(:)  = 0.
     343dqsrim(:)   = 0.
     344dqrmelt(:)  = 0.
     345dqsmelt(:)  = 0.
    321346dqrfreez(:) = 0.
    322347dqsfreez(:) = 0.
    … …  
    325350DO i = 1, klon
    326351
    327 
    328352  ! variables initialisation
    329   dqlrim = 0.0
    330   dqlcol = 0.0
    331   dqiagg = 0.0
    332   dqiauto = 0.0
    333   dqlauto = 0.0
     353  dqlrim  = 0.
     354  dqlcol  = 0.
     355  dqiagg  = 0.
     356  dqiauto = 0.
     357  dqlauto = 0.
     358
     359  !--hum_to_flux: coef to convert a specific quantity to a flux
     360  !-- hum_to_flux = rho * dz/dt = 1 / g * dP/dt
     361  hum_to_flux(i) = ( paprsdn(i) - paprsup(i) ) / RG / dtime
     362  qtot(i) = qvap(i) + qliq(i) + qice(i) &
     363          + ( raincld(i) + rainclr(i) + snowcld(i) + snowclr(i) ) / hum_to_flux(i)
    334364
    335365  !------------------------------------------------------------
    … …  
    344374
    345375  !--Initialisation
    346   !--hum_to_flux: coef to convert a specific quantity to a flux
    347   !-- hum_to_flux = rho * dz/dt = 1 / g * dP/dt
    348   hum_to_flux = ( paprsdn(i) - paprsup(i) ) / RG / dtime
    349376  precipfractot = precipfracclr(i) + precipfraccld(i)
    350377
    … …  
    446473    !--get in-cloud mean quantities. The two divisions by eff_cldfra are
    447474    !--then simplified.
    448 
     475    Eff_rain_liq = 1.
    449476    coef_col = gamma_col * 3. / 4. / rho_rain / r_rain * Eff_rain_liq
    450477    IF ((raincld(i) .GT. 0.) .AND. (coef_col .GT. 0.)) THEN
    451        !--Explicit version
    452        !dqlcol = - coef_col * qliq(i) * raincld(i) / precipfraccld(i) *dtime
    453        !--Semi-implicit version
    454        !dqlcol = qliq(i) * ( 1. / ( 1. + coef_col * raincld(i) / precipfraccld(i)*dtime ) - 1. )
    455        !--Implicit version
    456        !qrain = raincld(i) / hum_to_flux
    457        !coef_tmp = coef_col * dtime *( qrain / precipfraccld(i) + qliq(i) / eff_cldfra )
    458        !dqlcol = qliq(i) * ( 1. / ( 1. + 0.5 * ( coef_tmp - 1. + SQRT( &
    459        !    ( 1. - coef_tmp )**2. + 4. * coef_col * dtime *qrain / precipfraccld(i) ) &
    460        !                   ) ) - 1. )
    461        !dqlcol=max(dqlcol,-qliq(i))
    462        ! Exact version
    463        dqlcol=qliq(i)*(exp(-dtime * coef_col * raincld(i) / precipfraccld(i))-1.)
     478      !--Explicit version
     479      !dqlcol = - coef_col * qliq(i) * raincld(i) / precipfraccld(i) *dtime
     480      !--Semi-implicit version
     481      !dqlcol = qliq(i) * ( 1. / ( 1. + coef_col * raincld(i) / precipfraccld(i)*dtime ) - 1. )
     482      !--Implicit version
     483      !qrain = raincld(i) / hum_to_flux(i)
     484      !coef_tmp = coef_col * dtime * ( qrain / precipfraccld(i) + qliq(i) / eff_cldfra )
     485      !dqlcol = qliq(i) * ( 1. / ( 1. + 0.5 * ( coef_tmp - 1. + SQRT( &
     486      !    ( 1. - coef_tmp )**2. + 4. * coef_col * dtime * qrain / precipfraccld(i) ) &
     487      !                   ) ) - 1. )
     488      !--Barriers so that the processes do not consume more liquid/ice than
     489      !--available.
     490      !dqlcol = MAX( - qliq(i), dqlcol )
     491      !--Exact version
     492      dqlcol = qliq(i) * ( EXP( - dtime * coef_col * raincld(i) / precipfraccld(i) ) - 1. )
     493
     494      !--Add tendencies
     495      qliq(i) = qliq(i) + dqlcol
     496      raincld(i) = raincld(i) - dqlcol * hum_to_flux(i)
     497
     498      !--Outputs
     499      dqrcol(i)  = - dqlcol  / dtime
    464500    ENDIF
    465501
    466502    !--Same as for aggregation
    467     coef_agg=gamma_agg * 3. / 4. / rho_snow / r_snow * Eff_snow_ice
     503    !--Following Milbrandt and Yau 2005, it s a product of a collection
     504    !--efficiency and a sticking efficiency
     505    Eff_snow_ice = 0.05 * EXP( 0.1 * ( temp(i) - RTT ) )
     506    coef_agg = gamma_agg * 3. / 4. / rho_snow / r_snow * Eff_snow_ice
    468507    IF ((snowcld(i) .GT. 0.) .AND. (coef_agg .GT. 0.)) THEN
    469         !--Explicit version     
    470         !dqiagg = - coef_agg &
    471         !      * qice(i) * snowcld(i) / precipfraccld(i) * dtime
    472         ! Exact version
    473         dqiagg=qice(i)*(exp(-dtime * coef_agg * snowcld(i) / precipfraccld(i))-1.)
    474     ENDIF
    475     !--Barriers so that the processes do not consume more liquid/ice than
    476     !--available.
    477     dqlcol = MAX( - qliq(i), dqlcol )
    478     dqiagg = MAX( - qice(i), dqiagg )
    479 
    480     !--Add tendencies
    481     qliq(i) = qliq(i) + dqlcol
    482     qice(i) = qice(i) + dqiagg
    483     raincld(i) = raincld(i) - dqlcol * hum_to_flux
    484     snowcld(i) = snowcld(i) - dqiagg * hum_to_flux
     508      !--Explicit version
     509      !dqiagg = - coef_agg &
     510      !      * qice(i) * snowcld(i) / precipfraccld(i) * dtime
     511      !--Barriers so that the processes do not consume more liquid/ice than
     512      !--available.
     513      !dqiagg = MAX( - qice(i), dqiagg )
     514      !--Exact version
     515      dqiagg = qice(i) * ( EXP( - dtime * coef_agg * snowcld(i) / precipfraccld(i) ) - 1. )
     516
     517      !--Add tendencies
     518      qice(i) = qice(i) + dqiagg
     519      snowcld(i) = snowcld(i) - dqiagg * hum_to_flux(i)
     520
     521      !--Outputs
     522      dqsagg(i)  = - dqiagg  / dtime
     523    ENDIF
    485524
    486525
    … …  
    532571    qice(i) = qice(i) + dqiauto
    533572
    534     raincld(i) = raincld(i) - dqlauto * hum_to_flux
    535     snowcld(i) = snowcld(i) - dqiauto * hum_to_flux
     573    raincld(i) = raincld(i) - dqlauto * hum_to_flux(i)
     574    snowcld(i) = snowcld(i) - dqiauto * hum_to_flux(i)
     575
     576    !--Outputs
     577    dqsauto(i) = - dqiauto / dtime
     578    dqrauto(i) = - dqlauto / dtime
    536579
    537580
    … …  
    543586    !---------------------------------------------------------
    544587
    545     dqlrim=0.0
    546 
    547     ! remplacer la premiere ligne par "coef_rim" ?
    548     IF (snowcld(i) .GT. 0.) THEN
    549        dqlrim = - gamma_rim * 3. / 4. / rho_snow / r_snow * Eff_snow_liq &
    550               * qliq(i) * snowcld(i) /  precipfraccld(i) * dtime
    551     ENDIF
    552     dqlrim = MAX( - qliq(i), dqlrim )
    553 
    554     qliq(i) = qliq(i) + dqlrim
    555 
    556     snowcld(i) = snowcld(i) - dqlrim * hum_to_flux
    557 
     588    !--Following Seifert and Beheng 2006, assuming a cloud droplet diameter
     589    !--of 20 microns.
     590    Eff_snow_liq = 0.2
     591    coef_rim = gamma_rim * 3. / 4. / rho_snow / r_snow * Eff_snow_liq
     592    IF ((snowcld(i) .GT. 0.) .AND. (coef_rim .GT. 0.)) THEN
     593      !--Explicit version
     594      !dqlrim = - gamma_rim * 3. / 4. / rho_snow / r_snow * Eff_snow_liq &
     595      !       * qliq(i) * snowcld(i) /  precipfraccld(i) * dtime
     596      !--Barriers so that the processes do not consume more liquid than
     597      !--available.
     598      !dqlrim = MAX( - qliq(i), dqlrim )
     599      !--Exact version
     600      dqlrim = qliq(i) * ( EXP( - dtime * coef_col * snowcld(i) / precipfraccld(i) ) - 1. )
     601
     602      qliq(i) = qliq(i) + dqlrim
     603      snowcld(i) = snowcld(i) - dqlrim * hum_to_flux(i)
     604
     605      ! Latent heat of melting with precipitation thermalization
     606      temp(i) = temp(i) - dqlrim * RLMLT / RCPD &
     607                        / ( 1. + RVTMP2 * qtot(i) )
     608
     609      !--Outputs
     610      dqsrim(i)  = - dqlrim  / dtime
     611    ENDIF
    558612
    559613  ENDIF ! rneb .GE. seuil_neb
    560614
     615ENDDO ! iteration on klon
     616
     617
     618! Calculation of saturation specific humidity
     619! depending on temperature:
     620CALL calc_qsat_ecmwf(klon, temp, qzero, pplay, RTT, 0, .FALSE., qsat, dqsat)
     621
     622DO i = 1, klon
     623
     624  !---------------------------------------------------------
     625  !--                  MELTING
     626  !---------------------------------------------------------
     627
     628  IF ( ( snowclr(i) + snowcld(i) ) .GT. 0. ) THEN
     629    dqsmelt_max = MIN(0., ( RTT - temp(i) ) / RLMLT * RCPD &
     630                        * ( 1. + RVTMP2 * qtot(i) ))
     631   
     632    dqsclrmelt = 0.
     633    dqscldmelt = 0.
     634
     635    capa_snowflake = r_snow
     636    ! ATTENTION POUR L'INSTANT C'EST UN WBULB SELON FORMULE ECMWF
     637    ! ATTENTION EST-CE QVAP ?????
     638    temp_wetbulb = temp(i) - ( qsat(i) - qvap(i) ) &
     639                 * ( 1329.31 + 0.0074615 * ( pplay(i) - 0.85e5 ) &
     640                 - 40.637 * ( temp(i) - 275. ) )
     641
     642    IF ( snowclr(i) .GT. 0. ) THEN
     643      ! ATTENTION ATTENTION ATTENTION
     644      fallice_clr = 1.
     645      nb_snowflake_clr = snowclr(i) / precipfracclr(i) / fallice_clr &
     646                       / ( 4. / 3. * RPI * r_snow**3. * rho_snow )
     647      dqsclrmelt = - nb_snowflake_clr * 4. * RPI * air_thermal_conduct &
     648                   * capa_snowflake / RLMLT * coef_ventil &
     649                   * MAX(0., temp_wetbulb - RTT) * dtime
     650    ENDIF
     651
     652    IF ( snowcld(i) .GT. 0. ) THEN
     653      ! ATTENTION ATTENTION ATTENTION
     654      fallice_cld = 1.
     655      nb_snowflake_cld = snowcld(i) / precipfraccld(i) / fallice_cld &
     656                       / ( 4. / 3. * RPI * r_snow**3. * rho_snow )
     657      dqscldmelt = - nb_snowflake_cld * 4. * RPI * air_thermal_conduct &
     658                   * capa_snowflake / RLMLT * coef_ventil &
     659                   * MAX(0., temp_wetbulb - RTT) * dtime
     660    ENDIF
     661   
     662    ! barrier
     663    ! lower than bec. negative values
     664    dqstotmelt = dqsclrmelt + dqscldmelt
     665    IF ( dqstotmelt .LT. dqsmelt_max ) THEN
     666      dqsclrmelt = dqsmelt_max * dqsclrmelt / dqstotmelt
     667      dqscldmelt = dqsmelt_max * dqscldmelt / dqstotmelt
     668      dqstotmelt = dqsmelt_max
     669    ENDIF
     670
     671    ! update of rainfall and snowfall due to melting
     672    rainclr(i) = rainclr(i) - dqsclrmelt * hum_to_flux(i)
     673    raincld(i) = raincld(i) - dqscldmelt * hum_to_flux(i)
     674    snowclr(i) = snowclr(i) + dqsclrmelt * hum_to_flux(i)
     675    snowcld(i) = snowcld(i) + dqscldmelt * hum_to_flux(i)
     676
     677    ! Latent heat of melting with precipitation thermalization
     678    temp(i) = temp(i) + dqstotmelt * RLMLT / RCPD &
     679                      / ( 1. + RVTMP2 * qtot(i) )
     680
     681    dqrmelt(i) = - dqstotmelt / dtime
     682    dqsmelt(i) = dqstotmelt / dtime
     683
     684  ENDIF
    561685
    562686
    … …  
    565689  !---------------------------------------------------------
    566690
    567   !dqrfree_max = MINOUMAX(0., ( RTT - temp(i) ) / RLMLT * RCPD / ( 1. + RVTMP2 * ( qtot(i) + qprecip(i) ) ))
    568 
    569   !---------------------------------------------------------
    570   !--                  MELTING
    571   !---------------------------------------------------------
    572 
    573   !flux = velocity * N_0 * 4. / 3. * PI * r_snow**3. * rho_snow
    574 
    575   !IF ( ( snowclr(i) + snowcld(i) ) .GT. 0. ) THEN
    576   !  dqsmelt_max = MIN(0., ( RTT - temp(i) ) / RLMLT * RCPD / ( 1. + RVTMP2 * ( qtot(i) + qprecip(i) ) ))
    577   !  dsnowtotmelt_max = dqsmelt_max * hum_to_flux(i)
    578   ! 
    579   !  dsnowtotmelt = - nb_snowflake * 4. * PI * mol_diff_vap * snowflake_capa / RLMLT * coef_ventil &
    580   !              * MAX(0., ticebulb - RTT) &
    581   !              * ( paprsdn(i) - paprsup(i) ) / RG
    582   !  ! max bec. negative values
    583   !  dsnowtotmelt = MAX(dsnowtotmelt, dsnowmelt_max)
    584   !  dsnowclrmelt = dsnowtotmelt * snowclr(i) / ( snowclr(i) + snowcld(i) )
    585   !  dsnowcldmelt = dsnowtotmelt - dsnowclrmelt
    586 
    587 
    588   !  ! update of rainfall and snowfall due to melting
    589   !  rainclr(i) = rainclr(i) - dsnowclrmelt(i)
    590   !  raincld(i) = raincld(i) - dsnowcldmelt(i)
    591   !  snowclr(i) = snowclr(i) + dsnowclrmelt(i)
    592   !  snowcld(i) = snowcld(i) + dsnowcldmelt(i)
    593   !ENDIF
    594   !
    595   !! Latent heat of melting with precipitation thermalization
    596   !zt(i)=zt(i)-zifl(i)*zmelt*(RG*dtime)/(paprs(i,k)-paprs(i,k+1)) &
    597   !*RLMLT/RCPD/(1.0+RVTMP2*(zq(i)+zmqc(i)))
     691  IF ( ( rainclr(i) + raincld(i) ) .GT. 0. ) THEN
     692
     693    dqrfreez_max = MIN(0., ( temp(i) - RTT ) / RLMLT * RCPD &
     694                         * ( 1. + RVTMP2 * qtot(i) ))
     695   
     696    tau_freez = 1. / ( gamma_freez &
     697              * EXP( - alpha_freez * ( temp(i) - temp_nowater ) / ( RTT - temp_nowater ) ) )
     698    dqrtotfreez = raincld(i) / hum_to_flux(i) * ( EXP( - dtime / tau_freez ) - 1. )
     699
     700    ! barrier
     701    ! max bec. those are negative values
     702    dqrtotfreez = MAX(dqrtotfreez, dqrfreez_max)
     703
     704    ! partition between clear and cloudy air
     705    ! proportionnal to the rain fluxes in clear / cloud
     706    dqrclrfreez = dqrtotfreez * rainclr(i) / ( rainclr(i) + raincld(i) )
     707    dqrcldfreez = dqrtotfreez - dqrclrfreez
     708
     709    ! update of rainfall and snowfall due to melting
     710    rainclr(i) = rainclr(i) + dqrclrfreez * hum_to_flux(i)
     711    raincld(i) = raincld(i) + dqrcldfreez * hum_to_flux(i)
     712    snowclr(i) = snowclr(i) - dqrclrfreez * hum_to_flux(i)
     713    snowcld(i) = snowcld(i) - dqrcldfreez * hum_to_flux(i)
     714
     715    ! Latent heat of melting with precipitation thermalization
     716    temp(i) = temp(i) - dqrtotfreez * RLMLT / RCPD &
     717                      / ( 1. + RVTMP2 * qtot(i) )
     718
     719    dqrfreez(i) = dqrtotfreez / dtime
     720    dqsfreez(i) = - dqrtotfreez / dtime
     721
     722  ENDIF
    598723
    599724
    … …  
    609734! write output tendencies for rain and snow
    610735
    611   dqsrim(i)  = -dqlrim/dtime
    612   dqrcol(i)  = -dqlcol/dtime
    613   dqsagg(i)  = -dqiagg/dtime
    614   dqsauto(i) = -dqiauto/dtime
    615   dqrauto(i) = -dqlauto/dtime
    616 
    617 
    618 
    619736ENDDO
    620737

  • LMDZ6/trunk/libf/phylmd/lmdz_lscp_tools.F90

    r4665 r4818  
    359359            ! 'Cirrus regime'
    360360               fcirrus=acirrus-temp(i)/bcirrus
    361                IF (fcirrus .LT. qsl(i)/qsi(i)) THEN
     361               IF (fcirrus .GT. qsl(i)/qsi(i)) THEN
    362362                  gammasat(i)=qsl(i)/qsi(i)
    363363                  dgammasatdt(i)=(dqsl(i)*qsi(i)-dqsi(i)*qsl(i))/qsi(i)/qsi(i)

  • LMDZ6/trunk/libf/phylmd/lmdz_wake.F90

    r4816 r4818  
    24852485
    24862486    ! -     1/ Pressure of the level where dth changes sign.
    2487     print*,'WAKE LJYF'
     2487    !print*,'WAKE LJYF'
    24882488
    24892489    DO i = 1, klon

  • LMDZ6/trunk/libf/phylmd/physiq_mod.F90

    r4806 r4818  
    18571857   &    RG,RD,RCPD,RKAPPA,RLVTT,RETV)
    18581858       CALL ratqs_ini(klon,klev,iflag_thermals,lunout,nbsrf,is_lic,is_ter,RG,RV,RD,RCPD,RLSTT,RLVTT,RTT)
    1859        CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_sursat,iflag_ratqs,fl_cor_ebil,RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT,RD,RG)
     1859       CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_sursat,iflag_ratqs,fl_cor_ebil,RCPD,RLSTT,RLVTT,RLMLT,RVTMP2,RTT,RD,RG,RPI)
    18601860       CALL blowing_snow_ini(prt_level,lunout, &
    18611861                             RCPD, RLSTT, RLVTT, RLMLT, &

  • LMDZ6/trunk/libf/phylmdiso/physiq_mod.F90

    r4803 r4818  
    19491949   &    RG,RD,RCPD,RKAPPA,RLVTT,RETV)
    19501950       CALL ratqs_ini(klon,klev,iflag_thermals,lunout,nbsrf,is_lic,is_ter,RG,RV,RD,RCPD,RLSTT,RLVTT,RTT)
    1951        CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_sursat,iflag_ratqs,fl_cor_ebil,RCPD, RLSTT, RLVTT, RLMLT, RVTMP2, RTT,RD,RG)
     1951       CALL lscp_ini(pdtphys,lunout,prt_level,ok_ice_sursat,iflag_ratqs,fl_cor_ebil,RCPD,RLSTT,RLVTT,RLMLT,RVTMP2,RTT,RD,RG,RPI)
    19521952       CALL blowing_snow_ini(prt_level,lunout, &
    19531953                             RCPD, RLSTT, RLVTT, RLMLT, &
Note: See TracChangeset for help on using the changeset viewer.