Ignore:
Timestamp:
May 3, 2023, 6:21:08 PM (19 months ago)
Author:
evignon
Message:

merge de la branche blowing snow vers la trunk
premiere tentative
Etienne

Location:
LMDZ6/trunk
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • LMDZ6/trunk

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

    r4522 r4523  
    7979    USE wxios, ONLY: g_ctx, wxios_set_context
    8080#endif
    81     USE atke_turbulence_ini_mod, ONLY : atke_ini
    82     USE lscp_ini_mod, ONLY : lscp_ini
    8381    USE lscp_mod, ONLY : lscp
    8482    USE wake_ini_mod, ONLY : wake_ini
    8583    USE yamada_ini_mod, ONLY : yamada_ini
     84    USE atke_turbulence_ini_mod, ONLY : atke_ini
    8685    USE thermcell_ini_mod, ONLY : thermcell_ini
     86    USE blowing_snow_ini_mod, ONLY : blowing_snow_ini , qbst_bs
     87    USE lscp_ini_mod, ONLY : lscp_ini
    8788
    8889    !USE cmp_seri_mod
    … …  
    183184       ! [Variables internes non sauvegardees de la physique]
    184185       ! Variables locales pour effectuer les appels en serie
    185        t_seri,q_seri,ql_seri,qs_seri,u_seri,v_seri,tr_seri,rneb_seri, &
     186       t_seri,q_seri,ql_seri,qs_seri,qbs_seri,u_seri,v_seri,tr_seri,rneb_seri, &
    186187       ! Dynamic tendencies (diagnostics)
    187        d_t_dyn,d_q_dyn,d_ql_dyn,d_qs_dyn,d_u_dyn,d_v_dyn,d_tr_dyn,d_rneb_dyn, &
    188        d_q_dyn2d,d_ql_dyn2d,d_qs_dyn2d, &
     188       d_t_dyn,d_q_dyn,d_ql_dyn,d_qs_dyn,d_qbs_dyn,d_u_dyn,d_v_dyn,d_tr_dyn,d_rneb_dyn, &
     189       d_q_dyn2d,d_ql_dyn2d,d_qs_dyn2d,d_qbs_dyn2d, &
    189190       ! Physic tendencies
    190191       d_t_con,d_q_con,d_u_con,d_v_con, &
    … …  
    203204       plul_st,plul_th, &
    204205       !
    205        d_t_vdf,d_q_vdf,d_u_vdf,d_v_vdf,d_t_diss, &
     206       d_t_vdf,d_q_vdf, d_qbs_vdf, d_u_vdf,d_v_vdf,d_t_diss, &
    206207       d_t_vdf_x, d_t_vdf_w, &
    207208       d_q_vdf_x, d_q_vdf_w, &
    208209       d_ts, &
     210       !
     211       d_t_bs,d_q_bs,d_qbs_bs, &
    209212       !
    210213!       d_t_oli,d_u_oli,d_v_oli, &
    … …  
    254257       cldh, cldl,cldm, cldq, cldt,      &
    255258       JrNt,                             &
    256        dthmin, evap, fder, plcl, plfc,   &
    257        prw, prlw, prsw,                  &
     259       dthmin, evap, snowerosion,fder, plcl, plfc,   &
     260       prw, prlw, prsw, prbsw,                  &
    258261       s_lcl, s_pblh, s_pblt, s_therm,   &
    259262       cdragm, cdragh,                   &
    … …  
    343346       fsolsw, wfbils, wfbilo,  &
    344347       wfevap, wfrain, wfsnow,  & 
    345        prfl, psfl, fraca, Vprecip,  &
     348       prfl, psfl,bsfl, fraca, Vprecip,  &
    346349       zw2,  &
    347350       !
    … …  
    526529    !======================================================================
    527530    !
    528     ! indices de traceurs eau vapeur, liquide, glace, fraction nuageuse LS (optional)
    529     INTEGER,SAVE :: ivap, iliq, isol, irneb
    530 !$OMP THREADPRIVATE(ivap, iliq, isol, irneb)
     531    ! indices de traceurs eau vapeur, liquide, glace, fraction nuageuse LS (optional), blowing snow (optional)
     532    INTEGER,SAVE :: ivap, iliq, isol, irneb, ibs
     533!$OMP THREADPRIVATE(ivap, iliq, isol, irneb, ibs)
    531534    !
    532535    !
    … …  
    905908    REAL dialiq(klon,klev)  ! eau liquide nuageuse
    906909    REAL diafra(klon,klev)  ! fraction nuageuse
    907     REAL cldliq(klon,klev)  ! eau liquide nuageuse
     910    REAL radocond(klon,klev)  ! eau condensee nuageuse
    908911    !
    909912    !XXX PB
    910913    REAL fluxq(klon,klev, nbsrf)   ! flux turbulent d'humidite
     914    REAL fluxqbs(klon,klev, nbsrf)   ! flux turbulent de neige soufflee
    911915    !
    912916    REAL zxfluxt(klon, klev)
    913917    REAL zxfluxq(klon, klev)
     918    REAL zxfluxqbs(klon,klev)
    914919    REAL zxfluxu(klon, klev)
    915920    REAL zxfluxv(klon, klev)
    … …  
    10091014    !
    10101015    ! tendance nulles
    1011     REAL, dimension(klon,klev):: du0, dv0, dt0, dq0, dql0, dqi0
     1016    REAL, dimension(klon,klev):: du0, dv0, dt0, dq0, dql0, dqi0, dqbs0
    10121017#ifdef ISO
    10131018    REAL, dimension(ntraciso,klon,klev):: dxt0, dxtl0, dxti0
    … …  
    11561161    REAL ztsol(klon)
    11571162    REAL q2m(klon,nbsrf)  ! humidite a 2m
     1163    REAL fsnowerosion(klon,nbsrf) ! blowing snow flux at surface
     1164    REAL qbsfra  ! blowing snow fraction
    11581165#ifdef ISO
    11591166    REAL d_xtw(ntraciso),d_xtl(ntraciso), d_xts(ntraciso)
    … …  
    13721379       isol = strIdx(tracers(:)%name, addPhase('H2O', 's'))
    13731380       irneb= strIdx(tracers(:)%name, addPhase('H2O', 'r'))
     1381       ibs  = strIdx(tracers(:)%name, addPhase('H2O', 'b'))
    13741382       CALL init_etat0_limit_unstruct
    13751383       IF (.NOT. create_etat0_limit) CALL init_limit_read(days_elapsed)
    … …  
    14211429       ENDIF
    14221430
    1423        IF (ok_ice_sursat.AND.(nqo.NE.4)) THEN
     1431       IF (ok_ice_sursat.AND.(nqo.LT.4)) THEN
    14241432          WRITE (lunout, *) ' ok_ice_sursat=y requires 4 H2O tracers ', &
    14251433               '(H2O_g, H2O_l, H2O_s, H2O_r) but nqo=', nqo, '. Might as well stop here.'
    … …  
    14401448       ENDIF
    14411449
     1450        IF (ok_bs) THEN
     1451          abort_message='blowing snow cannot be activated with water isotopes yet'
     1452          CALL abort_physic(modname,abort_message, 1)
     1453         IF ((ok_ice_sursat.AND.nqo .LT.5).OR.(.NOT.ok_ice_sursat.AND.nqo.LT.4)) THEN
     1454             WRITE (lunout, *) 'activation of blowing snow needs a specific H2O tracer', &
     1455                               'but nqo=', nqo
     1456             abort_message='see above'
     1457             CALL abort_physic(modname,abort_message, 1)
     1458         ENDIF
     1459        ENDIF
    14421460       Ncvpaseq1 = 0
    14431461       dnwd0=0.0
    … …  
    18931911           CALL lscp_ini(pdtphys,ok_ice_sursat)
    18941912       endif
     1913       CALL blowing_snow_ini(prt_level,lunout, &
     1914                             RCPD, RLSTT, RLVTT, RLMLT, &
     1915                             RVTMP2, RTT,RD,RG)
    18951916
    18961917!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    … …  
    19191940       CALL phys_output_write(itap, pdtphys, paprs, pphis,                    &
    19201941                              pplay, lmax_th, aerosol_couple,                 &
    1921                               ok_ade, ok_aie, ok_volcan, ivap, iliq, isol, ok_sync,&
     1942                              ok_ade, ok_aie, ok_volcan, ivap, iliq, isol, ibs,  ok_sync,&
    19221943                              ptconv, read_climoz, clevSTD,                   &
    19231944                              ptconvth, d_u, d_t, qx, d_qx, zmasse,           &
    … …  
    24062427    dql0(:,:)=0.
    24072428    dqi0(:,:)=0.
     2429    dqbs0(:,:)=0.
    24082430#ifdef ISO
    24092431      dxt0(:,:,:)=0.
    … …  
    24632485          q_seri(i,k)  = qx(i,k,ivap)
    24642486          ql_seri(i,k) = qx(i,k,iliq)
     2487          qbs_seri(i,k) = 0.
    24652488          !CR: ATTENTION, on rajoute la variable glace
    24662489          IF (nqo.EQ.2) THEN             !--vapour and liquid only
    … …  
    24702493             qs_seri(i,k) = qx(i,k,isol)
    24712494             rneb_seri(i,k) = 0.
    2472           ELSE IF (nqo.EQ.4) THEN        !--vapour, liquid, ice and rneb
     2495          ELSE IF (nqo.GE.4) THEN        !--vapour, liquid, ice and rneb and blowing snow
    24732496             qs_seri(i,k) = qx(i,k,isol)
     2497             IF (ok_ice_sursat) THEN
    24742498             rneb_seri(i,k) = qx(i,k,irneb)
     2499             ENDIF
     2500             IF (ok_bs) THEN
     2501             qbs_seri(i,k)= qx(i,k,ibs)
     2502             ENDIF
     2503
    24752504          ENDIF
     2505
     2506
    24762507       ENDDO
    24772508    ENDDO
    … …  
    25152546    qql1(:)=0.0
    25162547    DO k = 1, klev
    2517       qql1(:)=qql1(:)+(q_seri(:,k)+ql_seri(:,k)+qs_seri(:,k))*zmasse(:,k)
     2548      qql1(:)=qql1(:)+(q_seri(:,k)+ql_seri(:,k)+qs_seri(:,k)+qbs_seri(:,k))*zmasse(:,k)
    25182549    ENDDO
    25192550#ifdef ISO
    … …  
    26312662       d_ql_dyn(:,:) = (ql_seri(:,:)-ql_ancien(:,:))/phys_tstep
    26322663       d_qs_dyn(:,:) = (qs_seri(:,:)-qs_ancien(:,:))/phys_tstep
     2664       d_qbs_dyn(:,:) = (qbs_seri(:,:)-qbs_ancien(:,:))/phys_tstep
    26332665       CALL water_int(klon,klev,q_seri,zmasse,zx_tmp_fi2d)
    26342666       d_q_dyn2d(:)=(zx_tmp_fi2d(:)-prw_ancien(:))/phys_tstep
    … …  
    26372669       CALL water_int(klon,klev,qs_seri,zmasse,zx_tmp_fi2d)
    26382670       d_qs_dyn2d(:)=(zx_tmp_fi2d(:)-prsw_ancien(:))/phys_tstep
     2671       CALL water_int(klon,klev,qbs_seri,zmasse,zx_tmp_fi2d)
     2672       d_qbs_dyn2d(:)=(zx_tmp_fi2d(:)-prbsw_ancien(:))/phys_tstep
    26392673       ! !! RomP >>>   td dyn traceur
    26402674       IF (nqtot > nqo) d_tr_dyn(:,:,:)=(tr_seri(:,:,:)-tr_ancien(:,:,:))/phys_tstep
    … …  
    27232757       d_ql_dyn2d(:) = 0.0
    27242758       d_qs_dyn2d(:) = 0.0
     2759       d_qbs_dyn2d(:)= 0.0
    27252760
    27262761#ifdef ISO
    … …  
    27462781       ! !! RomP <<<
    27472782       d_rneb_dyn(:,:)=0.0
     2783       d_qbs_dyn(:,:)=0.0
    27482784       ancien_ok = .TRUE.
    27492785    ENDIF
    … …  
    28782914
    28792915     CALL add_phys_tend &
    2880             (du0,dv0,d_t_eva,d_q_eva,d_ql_eva,d_qi_eva,paprs,&
     2916            (du0,dv0,d_t_eva,d_q_eva,d_ql_eva,d_qi_eva,dqbs0,paprs,&
    28812917               'eva',abortphy,flag_inhib_tend,itap,0 &
    28822918#ifdef ISO
    … …  
    30473083            longitude_deg, latitude_deg, rugoro,  zrmu0,      &
    30483084            sollwdown,    cldt,      &
    3049             rain_fall, snow_fall, solsw,   solswfdiff, sollw,     &
     3085            rain_fall, snow_fall, bs_fall, solsw,   solswfdiff, sollw,     &
    30503086            gustiness,                                &
    3051             t_seri,    q_seri,    u_seri,  v_seri,    &
     3087            t_seri,    q_seri,   qbs_seri, u_seri,  v_seri,    &
    30523088                                !nrlmd+jyg<
    30533089            wake_deltat, wake_deltaq, wake_cstar, wake_s, &
    … …  
    30603096                                !albedo SB >>>
    30613097                                ! albsol1,   albsol2,   sens,    evap,      &
    3062             albsol_dir,   albsol_dif,   sens,    evap,   & 
     3098            albsol_dir,   albsol_dif,   sens,    evap, snowerosion, & 
    30633099                                !albedo SB <<<
    30643100            albsol3_lic,runoff,   snowhgt,   qsnow, to_ice, sissnow, &
    30653101            zxtsol,    zxfluxlat, zt2m,    qsat2m,  zn2mout, &
    3066             d_t_vdf,   d_q_vdf,   d_u_vdf, d_v_vdf, d_t_diss, &
     3102            d_t_vdf,   d_q_vdf, d_qbs_vdf,  d_u_vdf, d_v_vdf, d_t_diss, &
    30673103                                !nrlmd<
    30683104                                !jyg<
    … …  
    30903126            fluxt,   fluxu,  fluxv, &
    30913127            dsens,     devap,     zxsnow, &
    3092             zxfluxt,   zxfluxq,   q2m,     fluxq, pbl_tke, &
     3128            zxfluxt,   zxfluxq,  zxfluxqbs,  q2m, fluxq, fluxqbs, pbl_tke, &
    30933129                                !nrlmd+jyg<
    30943130            wake_delta_pbl_TKE, &
    … …  
    32023238       IF (klon_glo==1) THEN
    32033239          CALL add_pbl_tend &
    3204                (d_u_vdf,d_v_vdf,d_t_vdf+d_t_diss,d_q_vdf,dql0,dqi0,paprs,&
     3240               (d_u_vdf,d_v_vdf,d_t_vdf+d_t_diss,d_q_vdf,dql0,dqi0,d_qbs_vdf,paprs,&
    32053241               'vdf',abortphy,flag_inhib_tend,itap &
    32063242#ifdef ISO
    … …  
    32103246       ELSE
    32113247          CALL add_phys_tend &
    3212                (d_u_vdf,d_v_vdf,d_t_vdf+d_t_diss,d_q_vdf,dql0,dqi0,paprs,&
     3248               (d_u_vdf,d_v_vdf,d_t_vdf+d_t_diss,d_q_vdf,dql0,dqi0,d_qbs_vdf,paprs,&
    32133249               'vdf',abortphy,flag_inhib_tend,itap,0 &
    32143250#ifdef ISO
    … …  
    32723308
    32733309    ENDIF
     3310
     3311    ! ==================================================================
     3312    ! Blowing snow sublimation and sedimentation
     3313
     3314    d_t_bs(:,:)=0.
     3315    d_q_bs(:,:)=0.
     3316    d_qbs_bs(:,:)=0.
     3317    bsfl(:,:)=0.
     3318    bs_fall(:)=0.
     3319    IF (ok_bs) THEN
     3320
     3321     CALL call_blowing_snow_sublim_sedim(klon,klev,phys_tstep,t_seri,q_seri,qbs_seri,pplay,paprs, &
     3322                                        d_t_bs,d_q_bs,d_qbs_bs,bsfl,bs_fall)
     3323
     3324     CALL add_phys_tend &
     3325               (du0,dv0,d_t_bs,d_q_bs,dql0,dqi0,d_qbs_bs,paprs,&
     3326               'bs',abortphy,flag_inhib_tend,itap,0  &
     3327#ifdef ISO                                       
     3328       &    ,dxt0,dxtl0,dxti0 &                   
     3329#endif                                       
     3330       &   )
     3331
     3332    ENDIF
     3333
    32743334    ! =================================================================== c
    32753335    !   Calcul de Qsat
    … …  
    39594019!!
    39604020!!
    3961     CALL add_phys_tend(d_u_con, d_v_con, d_t_con, d_q_con, dql0, dqi0, paprs, &
     4021    CALL add_phys_tend(d_u_con, d_v_con, d_t_con, d_q_con, dql0, dqi0, dqbs0, paprs, &
    39624022         'convection',abortphy,flag_inhib_tend,itap,0 &
    39634023#ifdef ISO
    … …  
    41914251       !-----------------------------------------------------------------------
    41924252       ! ajout des tendances des poches froides
    4193        CALL add_phys_tend(du0,dv0,d_t_wake,d_q_wake,dql0,dqi0,paprs,'wake', &
     4253       CALL add_phys_tend(du0,dv0,d_t_wake,d_q_wake,dql0,dqi0,dqbs0,paprs,'wake', &
    41944254            abortphy,flag_inhib_tend,itap,0 &
    41954255#ifdef ISO
    … …  
    44704530          !
    44714531          CALL add_phys_tend(d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs,  &
    4472                              dql0,dqi0,paprs,'thermals', abortphy,flag_inhib_tend,itap,0 &
     4532                             dql0,dqi0,dqbs0,paprs,'thermals', abortphy,flag_inhib_tend,itap,0 &
    44734533#ifdef ISO
    44744534     &    ,d_xt_ajs,dxtl0,dxti0 &
    … …  
    45864646          !--------------------------------------------------------------------
    45874647          ! ajout des tendances de l'ajustement sec ou des thermiques
    4588           CALL add_phys_tend(du0,dv0,d_t_ajsb,d_q_ajsb,dql0,dqi0,paprs, &
     4648          CALL add_phys_tend(du0,dv0,d_t_ajsb,d_q_ajsb,dql0,dqi0,dqbs0,paprs, &
    45894649               'ajsb',abortphy,flag_inhib_tend,itap,0 &
    45904650#ifdef ISO
    … …  
    47304790    CALL lscp(klon,klev,phys_tstep,missing_val,paprs,pplay, &
    47314791         t_seri, q_seri,ptconv,ratqs, &
    4732          d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, rneb_seri, &
    4733          cldliq, picefra, rain_lsc, snow_lsc, &
     4792         d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, rneblsvol, rneb_seri, & 
     4793         radocond, picefra, rain_lsc, snow_lsc, &
    47344794         frac_impa, frac_nucl, beta_prec_fisrt, &
    47354795         prfl, psfl, rhcl,  &
    … …  
    47434803    CALL fisrtilp(phys_tstep,paprs,pplay, &
    47444804         t_seri, q_seri,ptconv,ratqs, &
    4745          d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, cldliq, &
     4805         d_t_lsc, d_q_lsc, d_ql_lsc, d_qi_lsc, rneb, radocond, &
    47464806         rain_lsc, snow_lsc, &
    47474807         pfrac_impa, pfrac_nucl, pfrac_1nucl, &
    … …  
    47914851!    write(*,9000) "rcpv","rcw",rcpv,rcw,rcs,t_seri(1,1)
    47924852!-JLD
    4793     CALL add_phys_tend(du0,dv0,d_t_lsc,d_q_lsc,d_ql_lsc,d_qi_lsc,paprs, &
     4853    CALL add_phys_tend(du0,dv0,d_t_lsc,d_q_lsc,d_ql_lsc,d_qi_lsc,dqbs0,paprs, &
    47944854         'lsc',abortphy,flag_inhib_tend,itap,0 &
    47954855#ifdef ISO
    … …  
    48284888    ENDIF
    48294889
    4830     !---------------------------------------------------------------------------
     4890
     4891!---------------------------------------------------------------------------
    48314892    DO k = 1, klev
    48324893       DO i = 1, klon
    48334894          cldfra(i,k) = rneb(i,k)
    48344895          !CR: a quoi ca sert? Faut-il ajouter qs_seri?
    4835           IF (.NOT.new_oliq) cldliq(i,k) = ql_seri(i,k)
     4896          !EV: en effet etrange, j'ajouterais aussi qs_seri
     4897          !    plus largement, je nettoierais (enleverrais) ces lignes
     4898          IF (.NOT.new_oliq) radocond(i,k) = ql_seri(i,k)
    48364899       ENDDO
    48374900    ENDDO
    48384901
    48394902
    4840 
     4903    ! Option to activate the radiative effect of blowing snow (ok_rad_bs)
     4904    ! makes sense only if the new large scale condensation scheme is active
     4905    ! with the ok_icefra_lscp flag active as well
     4906
     4907    IF (ok_bs .AND. ok_rad_bs) THEN
     4908       IF (ok_new_lscp .AND. ok_icefra_lscp) THEN
     4909           DO k=1,klev
     4910             DO i=1,klon
     4911                radocond(i,k)=radocond(i,k)+qbs_seri(i,k)
     4912                picefra(i,k)=(radocond(i,k)*picefra(i,k)+qbs_seri(i,k))/(radocond(i,k))
     4913                qbsfra=min(qbs_seri(i,k)/qbst_bs,1.0)
     4914                cldfra(i,k)=max(cldfra(i,k),qbsfra)
     4915             ENDDO
     4916           ENDDO
     4917       ELSE
     4918          WRITE(lunout,*)"PAY ATTENTION, you try to activate the radiative effect of blowing snow"
     4919          WRITE(lunout,*)"with ok_new_lscp=false and/or ok_icefra_lscp=false"
     4920          abort_message='inconsistency in cloud phase for blowing snow'
     4921          CALL abort_physic(modname,abort_message,1)
     4922       ENDIF
     4923
     4924    ENDIF
    48414925#ifdef ISO     
    48424926!#ifdef ISOVERIF
    … …  
    50015085          DO i = 1, klon
    50025086             IF (diafra(i,k).GT.cldfra(i,k)) THEN
    5003                 cldliq(i,k) = dialiq(i,k)
     5087                radocond(i,k) = dialiq(i,k)
    50045088                cldfra(i,k) = diafra(i,k)
    50055089             ENDIF
    … …  
    50385122                DO i=1,klon
    50395123                   IF (ptconv(i,k).AND.ptconvth(i,k)) THEN
    5040                       cldliq(i,k)=cldliq(i,k)+rnebcon(i,k)*clwcon(i,k)
     5124                      radocond(i,k)=radocond(i,k)+rnebcon(i,k)*clwcon(i,k)
    50415125                      cldfra(i,k)=min(cldfra(i,k)+rnebcon(i,k),1.)
    50425126                   ELSE IF (ptconv(i,k)) THEN
    50435127                      cldfra(i,k)=rnebcon(i,k)
    5044                       cldliq(i,k)=rnebcon(i,k)*clwcon(i,k)
     5128                      radocond(i,k)=rnebcon(i,k)*clwcon(i,k)
    50455129                   ENDIF
    50465130                ENDDO
    … …  
    50515135                DO i=1,klon
    50525136                   cldfra(i,k)=min(cldfra(i,k)+rnebcon(i,k),1.)
    5053                    cldliq(i,k)=cldliq(i,k)+rnebcon(i,k)*clwcon(i,k)
     5137                   radocond(i,k)=radocond(i,k)+rnebcon(i,k)*clwcon(i,k)
    50545138                ENDDO
    50555139             ENDDO
    … …  
    50695153                   IF (ptconv(i,k).AND. .NOT.ptconvth(i,k)) THEN
    50705154                      cldfra(i,k)=rnebcon(i,k)
    5071                       cldliq(i,k)=rnebcon(i,k)*clwcon(i,k)
     5155                      radocond(i,k)=rnebcon(i,k)*clwcon(i,k)
    50725156                   ENDIF
    50735157                ENDDO
    … …  
    50805164          ! Ancienne version
    50815165          cldfra(:,:)=min(max(cldfra(:,:),rnebcon(:,:)),1.)
    5082           cldliq(:,:)=cldliq(:,:)+rnebcon(:,:)*clwcon(:,:)
     5166          radocond(:,:)=radocond(:,:)+rnebcon(:,:)*clwcon(:,:)
    50835167       ENDIF
    50845168
    … …  
    51005184          DO i = 1, klon
    51015185             IF (diafra(i,k).GT.cldfra(i,k)) THEN
    5102                 cldliq(i,k) = dialiq(i,k)
     5186                radocond(i,k) = dialiq(i,k)
    51035187                cldfra(i,k) = diafra(i,k)
    51045188             ENDIF
    … …  
    54755559          ENDIF
    54765560          CALL newmicro (flag_aerosol, ok_cdnc, bl95_b0, bl95_b1, &
    5477                paprs, pplay, t_seri, cldliq, picefra, cldfra, &
     5561               paprs, pplay, t_seri, radocond, picefra, cldfra, &
    54785562               cldtau, cldemi, cldh, cldl, cldm, cldt, cldq, &
    54795563               flwp, fiwp, flwc, fiwc, &
    … …  
    54835567       ELSE
    54845568          CALL nuage (paprs, pplay, &
    5485                t_seri, cldliq, picefra, cldfra, cldtau, cldemi, &
     5569               t_seri, radocond, picefra, cldfra, cldtau, cldemi, &
    54865570               cldh, cldl, cldm, cldt, cldq, &
    54875571               ok_aie, &
    … …  
    57925876    ENDDO
    57935877
    5794     CALL add_phys_tend(du0,dv0,d_t_swr,dq0,dql0,dqi0,paprs,'SW',abortphy,flag_inhib_tend,itap,0 &
     5878    CALL add_phys_tend(du0,dv0,d_t_swr,dq0,dql0,dqi0,dqbs0,paprs,'SW',abortphy,flag_inhib_tend,itap,0 &
    57955879#ifdef ISO
    57965880     &    ,dxt0,dxtl0,dxti0 &
    … …  
    57985882     &  )
    57995883    CALL prt_enerbil('SW',itap)
    5800     CALL add_phys_tend(du0,dv0,d_t_lwr,dq0,dql0,dqi0,paprs,'LW',abortphy,flag_inhib_tend,itap,0 &
     5884    CALL add_phys_tend(du0,dv0,d_t_lwr,dq0,dql0,dqi0,dqbs0,paprs,'LW',abortphy,flag_inhib_tend,itap,0 &
    58015885#ifdef ISO
    58025886     &    ,dxt0,dxtl0,dxti0 &
    … …  
    58465930          ! -> condition on zrel_oro can deactivate the drag on tilted planar terrains
    58475931          !    such as ice sheets (work by V. Wiener)
    5848           ! zpmm_orodr_t and zstd_orodr_t are activation thresholds set by F. Lott to 
     5932          ! zpmm_orodr_t and zstd_orodr_t are activation thresholds set by F. Lott to
    58495933          ! earn computation time but they are not physical.
    58505934          IF (((zpic(i)-zmea(i)).GT.zpmm_orodr_t).AND.(zstd(i).GT.zstd_orodr_t).AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
    … …  
    58775961       !-----------------------------------------------------------------------
    58785962       ! ajout des tendances de la trainee de l'orographie
    5879        CALL add_phys_tend(d_u_oro,d_v_oro,d_t_oro,dq0,dql0,dqi0,paprs,'oro', &
     5963       CALL add_phys_tend(d_u_oro,d_v_oro,d_t_oro,dq0,dql0,dqi0,dqbs0,paprs,'oro', &
    58805964            abortphy,flag_inhib_tend,itap,0 &
    58815965#ifdef ISO
    … …  
    59326016
    59336017       ! ajout des tendances de la portance de l'orographie
    5934        CALL add_phys_tend(d_u_lif, d_v_lif, d_t_lif, dq0, dql0, dqi0, paprs, &
     6018       CALL add_phys_tend(d_u_lif, d_v_lif, d_t_lif, dq0, dql0, dqi0, dqbs0,paprs, &
    59356019            'lif', abortphy,flag_inhib_tend,itap,0 &
    59366020#ifdef ISO
    … …  
    59616045       d_t_hin(:, :)=0.
    59626046       CALL add_phys_tend(du_gwd_hines, dv_gwd_hines, d_t_hin, dq0, dql0, &
    5963             dqi0, paprs, 'hin', abortphy,flag_inhib_tend,itap,0 &
     6047            dqi0, dqbs0,paprs, 'hin', abortphy,flag_inhib_tend,itap,0 &
    59646048#ifdef ISO
    59656049     &    ,dxt0,dxtl0,dxti0 &
    … …  
    59836067       ENDDO
    59846068
    5985        CALL add_phys_tend(du_gwd_front, dv_gwd_front, dt0, dq0, dql0, dqi0, &
     6069       CALL add_phys_tend(du_gwd_front, dv_gwd_front, dt0, dq0, dql0, dqi0, dqbs0, &
    59866070            paprs, 'front_gwd_rando', abortphy,flag_inhib_tend,itap,0 &
    59876071#ifdef ISO
    … …  
    59966080            rain_fall + snow_fall, zustr_gwd_rando, zvstr_gwd_rando, &
    59976081            du_gwd_rando, dv_gwd_rando, east_gwstress, west_gwstress)
    5998        CALL add_phys_tend(du_gwd_rando, dv_gwd_rando, dt0, dq0, dql0, dqi0, &
     6082       CALL add_phys_tend(du_gwd_rando, dv_gwd_rando, dt0, dq0, dql0, dqi0, dqbs0, &
    59996083            paprs, 'flott_gwd_rando', abortphy,flag_inhib_tend,itap,0 &
    60006084#ifdef ISO
    … …  
    60556139       d_xt_ch4_dtime(:,:,:) = d_xt_ch4(:,:,:)*phys_tstep
    60566140#endif
    6057        CALL add_phys_tend(du0, dv0, dt0, d_q_ch4_dtime, dql0, dqi0, paprs, &
     6141       CALL add_phys_tend(du0, dv0, dt0, d_q_ch4_dtime, dql0, dqi0, dqbs0, paprs, &
    60586142            'q_ch4', abortphy,flag_inhib_tend,itap,0 &
    60596143#ifdef ISO
    … …  
    61216205! car on peut s'attendre a ce que les petites echelles produisent aussi de la TKE
    61226206! Mais attention, cela ne va pas dans le sens de la conservation de l'energie!
    6123           IF ((zstd(i).GT.1.0).AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
     6207          IF ((zstd(i).GT.1.0) .AND.(zrel_oro(i).LE.zrel_oro_t)) THEN
    61246208             itest(i)=1
    61256209             igwd=igwd+1
    … …  
    63816465    !
    63826466
    6383     IF (type_trac=='repr') THEN
     6467    IF (type_trac == 'repr') THEN
    63846468!MM pas d'impact, car on recupere q_seri,tr_seri,t_seri via phys_local_var_mod
    63856469!MM                               dans Reprobus
    … …  
    64316515         presnivs, pphis,     pphi,     albsol1, &
    64326516         sh_in,   ch_in,    rhcl,      cldfra,   rneb, &
    6433          diafra,   cldliq,    itop_con, ibas_con, &
     6517         diafra,   radocond,    itop_con, ibas_con, &
    64346518         pmflxr,   pmflxs,    prfl,     psfl, &
    64356519         da,       phi,       mp,       upwd, &
    … …  
    64486532
    64496533     CALL add_phys_tend &
    6450             (du0,dv0,dt0,d_q_rep,d_ql_rep,d_qi_rep,paprs,&
     6534            (du0,dv0,dt0,d_q_rep,d_ql_rep,d_qi_rep,dqbs0,paprs,&
    64516535             'rep',abortphy,flag_inhib_tend,itap,0)
    64526536        IF (abortphy==1) Print*,'ERROR ABORT REP'
    … …  
    65256609    !   prlw = colonne eau liquide
    65266610    !   prlw = colonne eau solide
     6611    !   prbsw = colonne neige soufflee
    65276612    prw(:) = 0.
    65286613    prlw(:) = 0.
    65296614    prsw(:) = 0.
     6615    prbsw(:) = 0.
    65306616    DO k = 1, klev
    65316617       prw(:)  = prw(:)  + q_seri(:,k)*zmasse(:,k)
    65326618       prlw(:) = prlw(:) + ql_seri(:,k)*zmasse(:,k)
    65336619       prsw(:) = prsw(:) + qs_seri(:,k)*zmasse(:,k)
     6620       prbsw(:)= prbsw(:) + qbs_seri(:,k)*zmasse(:,k)
    65346621    ENDDO
    65356622
    … …  
    66026689          ENDIF
    66036690          !--ice_sursat: nqo=4, on ajoute rneb
    6604           IF (nqo == 4) THEN
     6691          IF (nqo.ge.4 .and. ok_ice_sursat) THEN
    66056692             d_qx(i,k,irneb) = ( rneb_seri(i,k) - qx(i,k,irneb) ) / phys_tstep
    66066693          ENDIF
     6694
     6695           IF (nqo.ge.4 .and. ok_bs) THEN
     6696             d_qx(i,k,ibs) = ( qbs_seri(i,k) - qx(i,k,ibs) ) / phys_tstep
     6697          ENDIF
     6698
     6699
    66076700       ENDDO
    66086701    ENDDO
    … …  
    66846777    ql_ancien(:,:) = ql_seri(:,:)
    66856778    qs_ancien(:,:) = qs_seri(:,:)
     6779    qbs_ancien(:,:) = qbs_seri(:,:)
    66866780    rneb_ancien(:,:) = rneb_seri(:,:)
    66876781#ifdef ISO
    … …  
    66936787    CALL water_int(klon,klev,ql_ancien,zmasse,prlw_ancien)
    66946788    CALL water_int(klon,klev,qs_ancien,zmasse,prsw_ancien)
     6789    CALL water_int(klon,klev,qbs_ancien,zmasse,prbsw_ancien)
    66956790    ! !! RomP >>>
    66966791    IF (nqtot > nqo) tr_ancien(:,:,:) = tr_seri(:,:,:)
    … …  
    68216916    CALL phys_output_write(itap, pdtphys, paprs, pphis,  &
    68226917         pplay, lmax_th, aerosol_couple,                 &
    6823          ok_ade, ok_aie, ok_volcan, ivap, iliq, isol,    &
     6918         ok_ade, ok_aie, ok_volcan, ivap, iliq, isol, ibs,   &
    68246919         ok_sync, ptconv, read_climoz, clevSTD,          &
    68256920         ptconvth, d_u, d_t, qx, d_qx, zmasse,           &
Note: See TracChangeset for help on using the changeset viewer.