Changeset 1026 for LMDZ4/trunk


Ignore:
Timestamp:
Oct 24, 2008, 5:09:39 PM (16 years ago)
Author:
lmdzadmin
Message:

Modifs thermiques
FH/IM

Location:
LMDZ4/trunk/libf/phylmd
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/trunk/libf/phylmd/calltherm.F90

    r987 r1026  
    88     &      ,fm_therm,entr_therm,detr_therm,zqasc,clwcon0,lmax,ratqscth,  &
    99     &       ratqsdiff,zqsatth,Ale_bl,Alp_bl,lalim_conv,wght_th, &
    10      &       zmax0,f0)
     10     &       zmax0,f0,zw2,fraca)
    1111
    1212      USE dimphy
    … …  
    3939      REAL wght_th(klon,klev)
    4040      INTEGER lalim_conv(klon)
     41      REAL zw2(klon,klev+1),fraca(klon,klev+1)
    4142
    4243!FH Update Thermiques
    … …  
    5051      real fmc_therm(klon,klev+1),zqasc(klon,klev)
    5152      real zqla(klon,klev)
     53      real zqta(klon,klev)
    5254      real wmax_sec(klon)
    5355      real zmax_sec(klon)
    … …  
    187189     &      ,r_aspect_thermals,l_mix_thermals,w2di_thermals  &
    188190     &      ,tau_thermals)
    189 !            CALL calcul_sec_entr(klon,klev,zdt
    190 !     s      ,pplay,paprs,pphi,zlev,debut
    191 !     s      ,u_seri,v_seri,t_seri,q_seri               
    192 !     s      ,zmax_sec,wmax_sec,zw_sec,lmix_sec
    193 !     s      ,r_aspect_thermals,l_mix_thermals,w2di_thermals
    194 !     s      ,tau_thermals,3)
    195 !           CALL thermcell_pluie_detr(klon,klev,zdt  &
    196 !    &      ,pplay,paprs,pphi,zlev,debut  &
    197 !    &      ,u_seri,v_seri,t_seri,q_seri  &
    198 !    &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
    199 !    &      ,zfm_therm,zentr_therm,zqla,lmax  &
    200 !    &      ,zmax_sec,wmax_sec,zw_sec,lmix_sec  &
    201 !    &      ,ratqscth,ratqsdiff,zqsatth  &
    202 !    &      ,r_aspect_thermals,l_mix_thermals,w2di_thermals  &
    203 !    &      ,tau_thermals)
    204191          else if (iflag_thermals.ge.13) then
    205192            CALL thermcell_main(itap,klon,klev,zdt  &
    … …  
    207194     &      ,u_seri,v_seri,t_seri,q_seri  &
    208195     &      ,d_u_the,d_v_the,d_t_the,d_q_the  &
    209      &      ,zfm_therm,zentr_therm,zdetr_therm,zqla,lmax  &
     196     &      ,zfm_therm,zentr_therm,zdetr_therm,zqasc,zqla,lmax  &
    210197     &      ,ratqscth,ratqsdiff,zqsatth  &
    211198     &      ,r_aspect_thermals,l_mix_thermals  &
    212199     &      ,tau_thermals,Ale,Alp,lalim_conv,wght_th &
    213      &      ,zmax0,f0)
     200     &      ,zmax0,f0,zw2,fraca)
    214201         endif
    215202
    … …  
    231218            entr_therm(:,k)=entr_therm(:,k)  &
    232219     &       +zentr_therm(:,k)*fact(:)
     220            detr_therm(:,k)=detr_therm(:,k)  &
     221     &       +zdetr_therm(:,k)*fact(:)
    233222      ENDDO
    234223       fm_therm(:,klev+1)=0.

  • LMDZ4/trunk/libf/phylmd/thermcell.h

    r987 r1026  
    33      integer w2di_thermals,isplit
    44      integer iflag_coupl,iflag_clos,iflag_wake
     5      integer iflag_thermals_ed,iflag_thermals_optflux
    56
    67      common/ctherm1/iflag_thermals,nsplit_thermals
    … …  
    89      common/ctherm3/w2di_thermals
    910      common/ctherm4/iflag_coupl,iflag_clos,iflag_wake
     11      common/ctherm5/iflag_thermals_ed,iflag_thermals_optflux
     12
    1013!$OMP THREADPRIVATE(/ctherm1/,/ctherm2/,/ctherm3/,/ctherm4/)

  • LMDZ4/trunk/libf/phylmd/thermcell_closure.F90

    r972 r1026  
    11      SUBROUTINE thermcell_closure(ngrid,nlay,r_aspect,ptimestep,rho,  &
    2      &   zlev,lalim,alim_star,zmax_sec,wmax_sec,zmax,wmax,f,lev_out)
     2     &   zlev,lalim,alim_star,alim_star_tot,zmax_sec,wmax_sec,zmax,wmax,f,lev_out)
    33
    44!-------------------------------------------------------------------------
    … …  
    88
    99#include "iniprint.h"
     10#include "thermcell.h"
    1011      INTEGER ngrid,nlay
    1112      INTEGER ig,k       
    … …  
    1516      INTEGER lalim(ngrid)
    1617      REAL alim_star(ngrid,nlay)
     18      REAL alim_star_tot(ngrid)
    1719      REAL rho(ngrid,nlay)
    1820      REAL zlev(ngrid,nlay)
    … …  
    4850                stop
    4951             endif
    50              if ((zmax_sec(ig).gt.1.e-10).and.(1.eq.1)) then
    51              f(ig)=wmax_sec(ig)/(max(500.,zmax_sec(ig))*r_aspect  &
     52             if ((zmax_sec(ig).gt.1.e-10).and.(iflag_thermals_ed.eq.0)) then
     53             f(ig)=wmax_sec(ig)*alim_star_tot(ig)/(max(500.,zmax_sec(ig))*r_aspect  &
    5254     &             *alim_star2(ig))
    5355!            f(ig)=f(ig)+(f0(ig)-f(ig))*exp((-ptimestep/  &
    5456!    &                     zmax_sec(ig))*wmax_sec(ig))
     57             print*,'closure dry',f(ig),wmax_sec(ig),alim_star_tot(ig),zmax_sec(ig)
    5558             else
    56              f(ig)=wmax(ig)/zdenom
     59             f(ig)=wmax(ig)*alim_star_tot(ig)/zdenom
    5760!            f(ig)=f(ig)+(f0(ig)-f(ig))*exp((-ptimestep/  &
    5861!     &                     zmax(ig))*wmax(ig))
     62             print*,'closure moist',f(ig),wmax(ig),alim_star_tot(ig),zmax(ig)
    5963             endif
    6064         endif

  • LMDZ4/trunk/libf/phylmd/thermcell_flux2.F90

    r987 r1026  
    1212      IMPLICIT NONE
    1313#include "iniprint.h"
     14#include "thermcell.h"
    1415     
    1516      INTEGER ig,l
    … …  
    4546      REAL fomass_max,alphamax
    4647      save fomass_max,alphamax
    47 !$OMP THREADPRIVATE(fomass_max,alphamax)
    4848
    4949      fomass_max=0.5
    … …  
    202202!Test sur fraca croissant
    203203!-------------------------------------------------------------------------
    204       if (1.eq.1) then
     204      if (iflag_thermals_optflux==0) then
    205205!     do l=1,klev
    206206         do ig=1,ngrid
    … …  
    228228!test sur flux de masse croissant
    229229!-------------------------------------------------------------------------
    230       if (1.eq.1) then
     230      if (iflag_thermals_optflux==0) then
    231231!     do l=1,klev
    232232         do ig=1,ngrid

  • LMDZ4/trunk/libf/phylmd/thermcell_height.F90

    r938 r1026  
    77      IMPLICIT NONE
    88#include "iniprint.h"
     9#include "thermcell.h"
    910
    1011      INTEGER ig,l
    … …  
    2324      REAL zmax0(ngrid)
    2425      REAL zmix(ngrid)
     26      REAL num(ngrid)
     27      REAL denom(ngrid)
    2528
    2629      REAL zlevinter(ngrid)
    … …  
    7073         zlevinter(ig)=zlev(ig,1)
    7174      enddo
     75
     76      if (iflag_thermals_ed.ge.1) then
     77
     78         num(:)=0.
     79         denom(:)=0.
     80         do ig=1,ngrid
     81          do l=1,nlay
     82             num(ig)=num(ig)+zw2(ig,l)*zlev(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     83             denom(ig)=denom(ig)+zw2(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     84          enddo
     85       enddo
     86       do ig=1,ngrid
     87       if (denom(ig).gt.1.e-10) then
     88          zmax(ig)=2.*num(ig)/denom(ig)
     89          zmax0(ig)=zmax(ig)
     90       endif
     91       enddo
     92
     93       else
     94
    7295      do  ig=1,ngrid
    7396! calcul de zlevinter
    … …  
    80103       zmax0(ig)=zmax(ig)
    81104      enddo
     105
     106
     107      endif
     108!endif iflag_thermals_ed
    82109!
    83110! def de  zmix continu (profil parabolique des vitesses)

  • LMDZ4/trunk/libf/phylmd/thermcell_init.F90

    r972 r1026  
    77      IMPLICIT NONE
    88#include "iniprint.h"
     9#include "thermcell.h"
    910
    1011      INTEGER l,ig
    … …  
    3435         lalim(ig)=1
    3536      enddo
     37
     38      if (iflag_thermals_ed.ge.1) then
     39!si la première couche est instable, on declenche un thermique
     40         do ig=1,ngrid
     41            if (ztv(ig,1).gt.ztv(ig,2)) then
     42               lmin(ig)=1
     43               lalim(ig)=2
     44               alim_star(ig,1)=1.
     45               alim_star_tot(ig)=alim_star(ig,1)
     46               print*,'init',alim_star(ig,1),alim_star_tot(ig)
     47            else
     48                lmin(ig)=1
     49                lalim(ig)=1
     50                alim_star(ig,1)=0.
     51                alim_star_tot(ig)=0.
     52            endif
     53         enddo
     54 
     55         else
     56!else iflag_thermals_ed=0 ancienne def de l alim
    3657
    3758!on ne considere que les premieres couches instables
    … …  
    131152      enddo
    132153       
     154!on remet alim_star_tot a 1
     155      do ig=1,ngrid
     156         alim_star_tot(ig)=1.
     157      enddo
     158
     159      endif
     160!endif iflag_thermals_ed
    133161      return
    134162      end 

  • LMDZ4/trunk/libf/phylmd/thermcell_main.F90

    r987 r1026  
    66     &                  ,pu,pv,pt,po  &
    77     &                  ,pduadj,pdvadj,pdtadj,pdoadj  &
    8      &                  ,fm0,entr0,detr0,zqla,lmax  &
     8     &                  ,fm0,entr0,detr0,zqta,zqla,lmax  &
    99     &                  ,ratqscth,ratqsdiff,zqsatth  &
    1010     &                  ,r_aspect,l_mix,tau_thermals &
    1111     &                  ,Ale_bl,Alp_bl,lalim_conv,wght_th &
    12      &                  ,zmax0, f0)
     12     &                  ,zmax0, f0,zw2,fraca)
    1313
    1414      USE dimphy
     15      USE comgeomphy , ONLY:rlond,rlatd
    1516      IMPLICIT NONE
    1617
    … …  
    4142
    4243#include "dimensions.h"
    43 !#include "dimphy.h"
    4444#include "YOMCST.h"
    4545#include "YOETHF.h"
    … …  
    8282      INTEGER lmax(klon),lmin(klon),lalim(klon)
    8383      INTEGER lmix(klon)
     84      INTEGER lmix_bis(klon)
    8485      real linter(klon)
    8586      real zmix(klon)
    86       real zmax(klon),zw2(klon,klev+1),ztva(klon,klev)
     87      real zmax(klon),zw2(klon,klev+1),ztva(klon,klev),zw_est(klon,klev+1)
     88!      real fraca(klon,klev)
     89
    8790      real zmax_sec(klon)
    88       real w_est(klon,klev+1)
    8991!on garde le zmax du pas de temps precedent
    9092      real zmax0(klon)
    … …  
    118120      real vardiff
    119121      real ratqsdiff(klon,klev)
    120       integer isplit,nsplit
    121       parameter (nsplit=10)
    122       data isplit/0/
    123       save isplit
    124 !$OMP THREADPRIVATE(isplit)
    125122
    126123      logical sorties
    … …  
    183180         entr0=0.
    184181         detr0=0.
     182
     183
     184! #define wrgrads_thermcell
     185#undef wrgrads_thermcell
     186#ifdef wrgrads_thermcell
     187! Initialisation des sorties grads pour les thermiques.
     188! Pour l'instant en 1D sur le point igout.
     189! Utilise par thermcell_out3d.h
     190         str10='therm'
     191         call inigrads(1,1,rlond(igout),1.,-180.,180.,jjm, &
     192     &   rlatd(igout),-90.,90.,1.,llm,pplay(igout,:),1.,   &
     193     &   ptimestep,str10,'therm ')
     194#endif
     195
     196
     197
    185198      endif
    186199
    … …  
    409422!IM 140508   CALL thermcell_plume(ngrid,nlay,ptimestep,ztv,zthl,po,zl,rhobarz,  &
    410423      CALL thermcell_plume(itap,ngrid,nlay,ptimestep,ztv,zthl,po,zl,rhobarz,  &
    411      &           zlev,pplev,pphi,zpspsk,l_mix,r_aspect,alim_star,  &
     424     &           zlev,pplev,pphi,zpspsk,l_mix,r_aspect,alim_star,alim_star_tot,  &
    412425     &           lalim,zmax_sec,f0,detr_star,entr_star,f_star,ztva,  &
    413      &           ztla,zqla,zqta,zha,zw2,zqsatth,lmix,linter,lev_out)
     426     &           ztla,zqla,zqta,zha,zw2,zw_est,zqsatth,lmix,lmix_bis,linter &
     427     &            ,lev_out,lunout1,igout)
    414428      if (prt_level.ge.1) print*,'apres thermcell_plume ',lev_out
    415429
    … …  
    445459! Fermeture,determination de f
    446460!-------------------------------------------------------------------------------
     461!
     462!avant closure: on redéfinit lalim, alim_star_tot et alim_star
     463!       do ig=1,klon
     464!       do l=2,lalim(ig)
     465!       alim_star(ig,l)=entr_star(ig,l)
     466!       entr_star(ig,l)=0.
     467!       enddo
     468!       enddo
    447469
    448470      CALL thermcell_closure(ngrid,nlay,r_aspect,ptimestep,rho,  &
    449      &   zlev,lalim,alim_star,zmax_sec,wmax_sec,zmax,wmax,f,lev_out)
     471     &   zlev,lalim,alim_star,alim_star_tot,zmax_sec,wmax_sec,zmax,wmax,f,lev_out)
    450472
    451473      if(prt_level.ge.1)print*,'thermcell_closure apres thermcell_closure'
    … …  
    757779!     print*,'15 OK convect8'
    758780
    759       isplit=isplit+1       
    760 
    761 
    762 #ifdef und
    763       if (prt_level.ge.1) print*,'thermcell_main sorties 1D'
    764 #include "thermcell_out1d.h"
    765 #endif
    766 
    767 
    768 #define troisD
    769 #undef troisD
    770781      if (prt_level.ge.1) print*,'thermcell_main sorties 3D'
    771 #ifdef troisD
     782#ifdef wrgrads_thermcell
    772783#include "thermcell_out3d.h"
    773784#endif

  • LMDZ4/trunk/libf/phylmd/thermcell_plume.F90

    r972 r1026  
    11      SUBROUTINE thermcell_plume(itap,ngrid,klev,ptimestep,ztv,zthl,po,zl,rhobarz,  &
    2      &           zlev,pplev,pphi,zpspsk,l_mix,r_aspect,alim_star,  &
     2     &           zlev,pplev,pphi,zpspsk,l_mix,r_aspect,alim_star,alim_star_tot,  &
    33     &           lalim,zmax_sec,f0,detr_star,entr_star,f_star,ztva,  &
    4      &           ztla,zqla,zqta,zha,zw2,zqsatth,lmix,linter,lev_out)
     4     &           ztla,zqla,zqta,zha,zw2,w_est,zqsatth,lmix,lmix_bis,linter &
     5     &           ,lev_out,lunout1,igout)
    56
    67!--------------------------------------------------------------------------
    … …  
    1415#include "FCTTRE.h"
    1516#include "iniprint.h"
     17#include "thermcell.h"
    1618
    1719      INTEGER itap
    18 
     20      INTEGER lunout1,igout
    1921      INTEGER ngrid,klev
    2022      REAL ptimestep
    … …  
    3638      integer lev_out                           ! niveau pour les print
    3739      real zcon2(ngrid)
     40   
     41      real alim_star_tot(ngrid)
    3842
    3943      REAL ztva(ngrid,klev)
    … …  
    6771      REAL zqla_est(ngrid,klev)
    6872      REAL zqsatth(ngrid,klev)
     73      REAL zta_est(ngrid,klev)
    6974
    7075      REAL linter(ngrid)
    … …  
    8186      PARAMETER (DDT0=.01)
    8287      logical Zsat
     88      REAL fact_gamma,fact_epsilon
     89      REAL c2(ngrid,klev)
    8390
    8491      Zsat=.false.
    … …  
    8693      RLvCp = RLVTT/RCPD
    8794     
     95      if (iflag_thermals_ed==0) then
     96         fact_gamma=1.
     97         fact_epsilon=1.
     98      else if (iflag_thermals_ed==1)  then
     99         fact_gamma=1.
     100         fact_epsilon=1.
     101      else if (iflag_thermals_ed==2)  then
     102         fact_gamma=1.
     103         fact_epsilon=2.
     104      endif
     105
    88106      do l=1,klev
    89107         do ig=1,ngrid
    … …  
    94112      enddo
    95113
     114!CR: attention test couche alim
     115!     do l=2,klev
     116!     do ig=1,ngrid
     117!        alim_star(ig,l)=0.
     118!     enddo
     119!     enddo
    96120!AM:initialisations du thermique
    97121      do k=1,klev
    … …  
    141165         linter(ig)=1.
    142166         lmix(ig)=1
     167         lmix_bis(ig)=2
    143168         wmaxa(ig)=0.
    144169      enddo
    … …  
    212237     &      /(f_star(ig,2)+alim_star(ig,2))**2+ &
    213238     &      2.*RG*(ztva(ig,2)-ztv(ig,2))/ztv(ig,2) &
     239!     &      *1./3. &
    214240     &      *(zlev(ig,3)-zlev(ig,2))
    215241       endif
    … …  
    253279        ztva_est(ig,l) = ztla(ig,l-1)*zpspsk(ig,l)+RLvCp*zqla_est(ig,l)
    254280        ztva_est(ig,l) = ztva_est(ig,l)/zpspsk(ig,l)
     281        zta_est(ig,l)=ztva_est(ig,l)
    255282        ztva_est(ig,l) = ztva_est(ig,l)*(1.+RETV*(zqta(ig,l-1)  &
    256283     &      -zqla_est(ig,l))-zqla_est(ig,l))
    … …  
    260287     &                   /(f_star(ig,l)+alim_star(ig,l))**2+  &
    261288     &                   2.*RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l)  &
     289!     &                   *1./3. &
    262290     &                   *(zlev(ig,l+1)-zlev(ig,l))
    263291             if (w_est(ig,l+1).lt.0.) then
    … …  
    268296!=======================
    269297
     298!CR:on vire les modifs
     299         if (iflag_thermals_ed==0) then
    270300
    271301! Modifications du calcul du detrainement.
    … …  
    523553        if (prt_level.ge.20) print*,'coucou calcul detr 450: ig, l', ig, l
    524554
     555       endif  ! iflag_thermals_ed==0
     556
     557!CR:nvlle def de entr_star et detr_star
     558      if (iflag_thermals_ed>=1) then
     559!      if (l.lt.lalim(ig)) then
     560!      if (l.lt.2) then
     561!        entr_star(ig,l)=0.
     562!        detr_star(ig,l)=0.
     563!      else
     564!      if (0.001.gt.(RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l))/(2.*w_est(ig,l+1)))) then
     565!         entr_star(ig,l)=0.001*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     566!      else
     567!         entr_star(ig,l)=  &
     568!     &                f_star(ig,l)/(2.*w_est(ig,l+1))        &
     569!     &                *RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l))   &
     570!     &                *(zlev(ig,l+1)-zlev(ig,l))
     571
     572 
     573         entr_star(ig,l)=MAX(0.*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l)),  &         
     574     &                f_star(ig,l)/(2.*w_est(ig,l+1))        &
     575     &                *RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l)   &
     576     &                *(zlev(ig,l+1)-zlev(ig,l))) &
     577     &                +0.0001*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     578
     579        if (((ztva_est(ig,l)-ztv(ig,l)).gt.1.e-10).and.(l.le.lmix_bis(ig))) then
     580            alim_star_tot(ig)=alim_star_tot(ig)+entr_star(ig,l)
     581            lalim(ig)=lmix_bis(ig)
     582            print*,'alim_star_tot',alim_star_tot(ig),entr_star(ig,l)
     583        endif
     584
     585        if (((ztva_est(ig,l)-ztv(ig,l)).gt.1.e-10).and.(l.le.lmix_bis(ig))) then
     586!        c2(ig,l)=2500000.*zqla_est(ig,l)/(1004.*zta_est(ig,l))
     587         c2(ig,l)=0.001
     588         detr_star(ig,l)=MAX(0.*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l)),  &
     589     &                c2(ig,l)*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l)) &
     590     &                -f_star(ig,l)/(2.*w_est(ig,l+1))       &
     591     &                *RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l)       &
     592     &                *(zlev(ig,l+1)-zlev(ig,l)))                    &
     593     &                +0.0001*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     594
     595       else
     596!         c2(ig,l)=2500000.*zqla_est(ig,l)/(1004.*zta_est(ig,l))
     597          c2(ig,l)=0.003
     598
     599         detr_star(ig,l)=MAX(0.*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l)),  &
     600     &                c2(ig,l)*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l)) &
     601     &                -f_star(ig,l)/(2.*w_est(ig,l+1))       &
     602     &                *RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l)       &
     603     &                *(zlev(ig,l+1)-zlev(ig,l))) &
     604     &                +0.0002*f_star(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     605       endif
     606         
     607           
     608!        detr_star(ig,l)=detr_star(ig,l)*3.
     609!        if (l.lt.lalim(ig)) then
     610!          entr_star(ig,l)=0.
     611!        endif
     612!        if (l.lt.2) then
     613!          entr_star(ig,l)=0.
     614!          detr_star(ig,l)=0.
     615!        endif
     616
     617
     618!      endif
     619!      else if ((ztva_est(ig,l)-ztv(ig,l)).gt.1.e-10) then
     620!      entr_star(ig,l)=MAX(0.,0.8*f_star(ig,l)/(2.*w_est(ig,l+1))        &
     621!     &                *RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l))   &
     622!     &                *(zlev(ig,l+1)-zlev(ig,l))
     623!      detr_star(ig,l)=0.002*f_star(ig,l)                         &
     624!     &                *(zlev(ig,l+1)-zlev(ig,l))
     625!      else
     626!      entr_star(ig,l)=0.001*f_star(ig,l)                         &
     627!     &                *(zlev(ig,l+1)-zlev(ig,l))
     628!      detr_star(ig,l)=MAX(0.,-0.2*f_star(ig,l)/(2.*w_est(ig,l+1))       &
     629!     &                *RG*(ztva_est(ig,l)-ztv(ig,l))/ztv(ig,l))       &
     630!     &                *(zlev(ig,l+1)-zlev(ig,l))                      &
     631!     &                +0.002*f_star(ig,l)                             &
     632!     &                *(zlev(ig,l+1)-zlev(ig,l))
     633!      endif
     634
     635      endif   ! iflag_thermals_ed==1
     636
    525637!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    526638! FH inutile si on conserve comme on l'a fait plus haut entr=detr
    … …  
    533645!
    534646!prise en compte du detrainement et de l entrainement dans le calcul du flux
     647
    535648      f_star(ig,l+1)=f_star(ig,l)+alim_star(ig,l)+entr_star(ig,l)  &
    536649     &              -detr_star(ig,l)
    … …  
    597710           zw2(ig,l+1)=zw2(ig,l)*  &
    598711     &                 ((f_star(ig,l))**2)  &
    599      &                 /(f_star(ig,l+1)+detr_star(ig,l))**2+             &
     712!  Tests de Catherine
     713!     &                 /(f_star(ig,l+1)+detr_star(ig,l))**2+             &
     714     &      /(f_star(ig,l+1)+detr_star(ig,l)-entr_star(ig,l)*(1.-fact_epsilon))**2+ &
    600715     &                 2.*RG*(ztva(ig,l)-ztv(ig,l))/ztv(ig,l)  &
     716     &                 *fact_gamma &
    601717     &                 *(zlev(ig,l+1)-zlev(ig,l))
    602                
     718!prise en compte des forces de pression que qd flottabilité<0
     719!              zw2(ig,l+1)=zw2(ig,l)*  &
     720!     &            1./(1.+2.*entr_star(ig,l)/f_star(ig,l)) + &       
     721!     &                 (f_star(ig,l))**2 &
     722!     &                 /(f_star(ig,l)+entr_star(ig,l))**2+ &
     723!     &                 (f_star(ig,l)-2.*entr_star(ig,l))**2/(f_star(ig,l)+2.*entr_star(ig,l))**2+  &       
     724!     &                 /(f_star(ig,l+1)+detr_star(ig,l)-entr_star(ig,l)*(1.-2.))**2+ &
     725!     &                 /(f_star(ig,l)**2+2.*2.*detr_star(ig,l)*f_star(ig,l)+2.*entr_star(ig,l)*f_star(ig,l))+ &
     726!     &                 2.*RG*(ztva(ig,l)-ztv(ig,l))/ztv(ig,l)  &
     727!     &                 *1./3. &
     728!     &                 *(zlev(ig,l+1)-zlev(ig,l))
     729         
     730!        write(30,*),l+1,zw2(ig,l+1)-zw2(ig,l), &
     731!     &              -2.*entr_star(ig,l)/f_star(ig,l)*zw2(ig,l), &
     732!     &               2.*RG*(ztva(ig,l)-ztv(ig,l))/ztv(ig,l)*(zlev(ig,l+1)-zlev(ig,l))
     733
     734 
     735!             zw2(ig,l+1)=zw2(ig,l)*  &
     736!     &                 (2.-2.*entr_star(ig,l)/f_star(ig,l)) & 
     737!     &                 -zw2(ig,l-1)+  &       
     738!     &                 2.*RG*(ztva(ig,l)-ztv(ig,l))/ztv(ig,l)  &
     739!     &                 *1./3. &
     740!     &                 *(zlev(ig,l+1)-zlev(ig,l))             
     741
    603742            endif
    604743        endif
    … …  
    614753            endif
    615754
    616 
     755!        if ((zw2(ig,l).gt.0.).and. (zw2(ig,l+1).le.0.)) then
    617756        if (zw2(ig,l+1).lt.0.) then
    618757           linter(ig)=(l*(zw2(ig,l+1)-zw2(ig,l))  &
    … …  
    625764        if (wa_moy(ig,l+1).gt.wmaxa(ig)) then
    626765!   lmix est le niveau de la couche ou w (wa_moy) est maximum
     766!on rajoute le calcul de lmix_bis
     767            if (zqla(ig,l).lt.1.e-10) then
     768               lmix_bis(ig)=l+1
     769            endif
    627770            lmix(ig)=l+1
    628771            wmaxa(ig)=wa_moy(ig,l+1)
    … …  
    631774      enddo
    632775
     776!on remplace a* par e* ds premiere couche
     777!      if (iflag_thermals_ed.ge.1) then
     778!       do ig=1,ngrid
     779!       do l=2,klev
     780!          if (l.lt.lalim(ig)) then
     781!             alim_star(ig,l)=entr_star(ig,l)
     782!          endif
     783!       enddo
     784!       enddo
     785!       do ig=1,ngrid
     786!          lalim(ig)=lmix_bis(ig)
     787!       enddo
     788!      endif
     789       if (iflag_thermals_ed.ge.1) then
     790          do ig=1,ngrid
     791             do l=2,lalim(ig)
     792                alim_star(ig,l)=entr_star(ig,l)
     793                entr_star(ig,l)=0.
     794             enddo
     795           enddo
     796       endif
    633797        if (prt_level.ge.20) print*,'coucou calcul detr 470: ig, l', ig, l
    634 #ifdef troisD
    635       call wrgradsfi(1,klev,zqla_est,'ql_es     ','ql_es     ')
    636       call wrgradsfi(1,klev,entr_star1,'e_st1     ','e_st1     ')
    637       call wrgradsfi(1,klev,entr_star2,'e_st2     ','e_st2     ')
    638       call wrgradsfi(1,klev,detr_stara,'d_sta     ','d_sta     ')
    639       call wrgradsfi(1,klev,detr_starb,'d_stb     ','d_stb     ')
    640       call wrgradsfi(1,klev,detr_starc,'d_stc     ','d_stc     ')
    641       call wrgradsfi(1,klev,zqla0     ,'zqla0     ','zqla0     ')
    642       call wrgradsfi(1,klev,detr_star0,'d_st0     ','d_st0     ')
    643       call wrgradsfi(1,klev,detr_star1,'d_st1     ','d_st1     ')
    644       call wrgradsfi(1,klev,detr_star2,'d_st2     ','d_st2     ')
    645 #endif
    646798
    647799!     print*,'thermcell_plume OK'
Note: See TracChangeset for help on using the changeset viewer.