Ignore:
Timestamp:
Jan 14, 2010, 3:35:30 PM (14 years ago)
Author:
Laurent Fairhead
Message:

Modifications pour la nouvelle version des thermiques (2009/2010) CR et FH

File:
1 edited

Legend:

Unmodified
Added
Removed

  • LMDZ4/branches/LMDZ4V5.0-dev/libf/phylmd/thermcell_main.F90

    r1146 r1294  
    2222!   de "thermiques" explicitement representes avec processus nuageux
    2323!
    24 !   Réécriture à partir d'un listing papier à Habas, le 14/02/00
    25 !
    26 !   le thermique est supposé homogène et dissipé par mélange avec
    27 !   son environnement. la longueur l_mix contrôle l'efficacité du
    28 !   mélange
    29 !
    30 !   Le calcul du transport des différentes espèces se fait en prenant
     24!   Reecriture a partir d'un listing papier a Habas, le 14/02/00
     25!
     26!   le thermique est suppose homogene et dissipe par melange avec
     27!   son environnement. la longueur l_mix controle l'efficacite du
     28!   melange
     29!
     30!   Le calcul du transport des differentes especes se fait en prenant
    3131!   en compte:
    3232!     1. un flux de masse montant
    … …  
    8585      real linter(klon)
    8686      real zmix(klon)
    87       real zmax(klon),zw2(klon,klev+1),ztva(klon,klev),zw_est(klon,klev+1)
     87      real zmax(klon),zw2(klon,klev+1),ztva(klon,klev),zw_est(klon,klev+1),ztva_est(klon,klev)
    8888!      real fraca(klon,klev)
    8989
    … …  
    115115! FH probleme de dimensionnement avec l'allocation dynamique
    116116!     common/comtherm/thetath2,wth2
     117      real wq(klon,klev)
     118      real wthl(klon,klev)
     119      real wthv(klon,klev)
    117120   
    118121      real ratqscth(klon,klev)
    … …  
    142145      real f_star(klon,klev+1),entr_star(klon,klev)
    143146      real detr_star(klon,klev)
    144       real alim_star_tot(klon),alim_star2(klon)
     147      real alim_star_tot(klon)
    145148      real alim_star(klon,klev)
     149      real alim_star_clos(klon,klev)
    146150      real f(klon), f0(klon)
    147151!FH/IM     save f0
    … …  
    149153      logical debut
    150154       real seuil
     155      real csc(klon,klev)
    151156
    152157!
    … …  
    220225      ENDIF
    221226!
    222 !Initialisation
    223 !
    224 !    IF (1.eq.0) THEN
    225 !     do ig=1,klon     
    226 !FH/IM 130308     if ((debut).or.((.not.debut).and.(f0(ig).lt.1.e-10))) then
    227 !     if ((.not.debut).and.(f0(ig).lt.1.e-10)) then
    228 !           f0(ig)=1.e-5
    229 !           zmax0(ig)=40.
    230 !v1d        therm=.false.
    231 !     endif
    232 !     enddo
    233 !    ENDIF !(1.eq.0) THEN
    234      if (prt_level.ge.10)write(lunout,*)                                &
    235     &     'WARNING thermcell_main f0=max(f0,1.e-2)'
     227     write(lunout,*)'WARNING thermcell_main f0=max(f0,1.e-2)'
    236228     do ig=1,klon
    237229      if (prt_level.ge.20) then
    … …  
    239231      endif
    240232         f0(ig)=max(f0(ig),1.e-2)
     233         zmax0(ig)=max(zmax0(ig),40.)
    241234!IMmarche pas ?!       if (f0(ig)<1.e-2) f0(ig)=1.e-2
    242235     enddo
    … …  
    364357
    365358!------------------------------------------------------------------
    366 !  1. alim_star est le profil vertical de l'alimentation à la base du
    367 !     panache thermique, calculé à partir de la flotabilité de l'air sec
     359!  1. alim_star est le profil vertical de l'alimentation a la base du
     360!     panache thermique, calcule a partir de la flotabilite de l'air sec
    368361!  2. lmin et lalim sont les indices inferieurs et superieurs de alim_star
    369362!------------------------------------------------------------------
    370363!
    371364      entr_star=0. ; detr_star=0. ; alim_star=0. ; alim_star_tot=0.
    372       CALL thermcell_init(ngrid,nlay,ztv,zlay,zlev,  &
    373      &                    lalim,lmin,alim_star,alim_star_tot,lev_out)
    374 
    375 call test_ltherm(ngrid,nlay,pplev,pplay,lmin,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_init lmin  ')
    376 call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_init lalim ')
    377 
    378 
    379       if (prt_level.ge.1) print*,'thermcell_main apres thermcell_init'
    380       if (prt_level.ge.10) then
    381          write(lunout1,*) 'Dans thermcell_main 1'
    382          write(lunout1,*) 'lmin ',lmin(igout)
    383          write(lunout1,*) 'lalim ',lalim(igout)
    384          write(lunout1,*) ' ig l alim_star thetav'
    385          write(lunout1,'(i6,i4,2e15.5)') (igout,l,alim_star(igout,l) &
    386      &   ,ztv(igout,l),l=1,lalim(igout)+4)
    387       endif
    388 
    389 !v1d      do ig=1,klon
    390 !v1d     if (alim_star(ig,1).gt.1.e-10) then
    391 !v1d     therm=.true.
    392 !v1d     endif
    393 !v1d      enddo
     365      lmin=1
     366
    394367!-----------------------------------------------------------------------------
    395368!  3. wmax_sec et zmax_sec sont les vitesses et altitudes maximum d'un
    396369!     panache sec conservatif (e=d=0) alimente selon alim_star
    397370!     Il s'agit d'un calcul de type CAPE
    398 !     zmax_sec est utilisé pour déterminer la géométrie du thermique.
     371!     zmax_sec est utilise pour determiner la geometrie du thermique.
    399372!------------------------------------------------------------------------------
    400 !
     373!---------------------------------------------------------------------------------
     374!calcul du melange et des variables dans le thermique
     375!--------------------------------------------------------------------------------
     376!
     377      if (prt_level.ge.1) print*,'avant thermcell_plume ',lev_out
     378!IM 140508   CALL thermcell_plume(ngrid,nlay,ptimestep,ztv,zthl,po,zl,rhobarz,  &
     379      CALL thermcell_plume(itap,ngrid,nlay,ptimestep,ztv,zthl,po,zl,rhobarz,  &
     380     &           zlev,pplev,pphi,zpspsk,alim_star,alim_star_tot,  &
     381     &           lalim,f0,detr_star,entr_star,f_star,csc,ztva,  &
     382     &           ztla,zqla,zqta,zha,zw2,zw_est,ztva_est,zqsatth,lmix,lmix_bis,linter &
     383     &            ,lev_out,lunout1,igout)
     384      if (prt_level.ge.1) print*,'apres thermcell_plume ',lev_out
     385
     386      call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_plum lalim ')
     387      call test_ltherm(ngrid,nlay,pplev,pplay,lmix ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_plum lmix  ')
     388
     389      if (prt_level.ge.1) print*,'thermcell_main apres thermcell_plume'
     390      if (prt_level.ge.10) then
     391         write(lunout1,*) 'Dans thermcell_main 2'
     392         write(lunout1,*) 'lmin ',lmin(igout)
     393         write(lunout1,*) 'lalim ',lalim(igout)
     394         write(lunout1,*) ' ig l alim_star entr_star detr_star f_star '
     395         write(lunout1,'(i6,i4,4e15.5)') (igout,l,alim_star(igout,l),entr_star(igout,l),detr_star(igout,l) &
     396     &    ,f_star(igout,l+1),l=1,nint(linter(igout))+5)
     397      endif
     398
     399!-------------------------------------------------------------------------------
     400! Calcul des caracteristiques du thermique:zmax,zmix,wmax
     401!-------------------------------------------------------------------------------
     402!
     403      CALL thermcell_height(ngrid,nlay,lalim,lmin,linter,lmix,zw2,  &
     404     &           zlev,lmax,zmax,zmax0,zmix,wmax,lev_out)
     405
     406
     407
     408      call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lalim ')
     409      call test_ltherm(ngrid,nlay,pplev,pplay,lmin ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmin  ')
     410      call test_ltherm(ngrid,nlay,pplev,pplay,lmix ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmix  ')
     411      call test_ltherm(ngrid,nlay,pplev,pplay,lmax ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmax  ')
     412
     413      if (prt_level.ge.1) print*,'thermcell_main apres thermcell_height'
     414
     415!-------------------------------------------------------------------------------
     416! Fermeture,determination de f
     417!-------------------------------------------------------------------------------
     418!
     419!
     420      write(lunout,*)'THERM NOUVEAU RIO2009 31B'
    401421      CALL thermcell_dry(ngrid,nlay,zlev,pphi,ztv,alim_star,  &
    402      &                      lalim,lmin,zmax_sec,wmax_sec,lev_out)
     422    &                      lalim,lmin,zmax_sec,wmax_sec,lev_out)
    403423
    404424call test_ltherm(ngrid,nlay,pplev,pplay,lmin,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_dry  lmin  ')
    … …  
    417437
    418438
    419 !---------------------------------------------------------------------------------
    420 !calcul du melange et des variables dans le thermique
    421 !--------------------------------------------------------------------------------
    422 !
    423       if (prt_level.ge.1) print*,'avant thermcell_plume ',lev_out
    424 !IM 140508   CALL thermcell_plume(ngrid,nlay,ptimestep,ztv,zthl,po,zl,rhobarz,  &
    425       CALL thermcell_plume(itap,ngrid,nlay,ptimestep,ztv,zthl,po,zl,rhobarz,  &
    426      &           zlev,pplev,pphi,zpspsk,l_mix,r_aspect,alim_star,alim_star_tot,  &
    427      &           lalim,zmax_sec,f0,detr_star,entr_star,f_star,ztva,  &
    428      &           ztla,zqla,zqta,zha,zw2,zw_est,zqsatth,lmix,lmix_bis,linter &
    429      &            ,lev_out,lunout1,igout)
    430       if (prt_level.ge.1) print*,'apres thermcell_plume ',lev_out
    431 
    432       call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_plum lalim ')
    433       call test_ltherm(ngrid,nlay,pplev,pplay,lmix ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_plum lmix  ')
    434 
    435       if (prt_level.ge.1) print*,'thermcell_main apres thermcell_plume'
    436       if (prt_level.ge.10) then
    437          write(lunout1,*) 'Dans thermcell_main 2'
    438          write(lunout1,*) 'lmin ',lmin(igout)
    439          write(lunout1,*) 'lalim ',lalim(igout)
    440          write(lunout1,*) ' ig l alim_star entr_star detr_star f_star '
    441          write(lunout1,'(i6,i4,4e15.5)') (igout,l,alim_star(igout,l),entr_star(igout,l),detr_star(igout,l) &
    442      &    ,f_star(igout,l+1),l=1,nint(linter(igout))+5)
    443       endif
    444 
    445 !-------------------------------------------------------------------------------
    446 ! Calcul des caracteristiques du thermique:zmax,zmix,wmax
    447 !-------------------------------------------------------------------------------
    448 !
    449       CALL thermcell_height(ngrid,nlay,lalim,lmin,linter,lmix,zw2,  &
    450      &           zlev,lmax,zmax,zmax0,zmix,wmax,lev_out)
    451 
    452 
    453       call test_ltherm(ngrid,nlay,pplev,pplay,lalim,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lalim ')
    454       call test_ltherm(ngrid,nlay,pplev,pplay,lmin ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmin  ')
    455       call test_ltherm(ngrid,nlay,pplev,pplay,lmix ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmix  ')
    456       call test_ltherm(ngrid,nlay,pplev,pplay,lmax ,seuil,ztv,po,ztva,zqla,f_star,zw2,'thermcell_heig lmax  ')
    457 
    458       if (prt_level.ge.1) print*,'thermcell_main apres thermcell_height'
    459 
    460 !-------------------------------------------------------------------------------
    461 ! Fermeture,determination de f
    462 !-------------------------------------------------------------------------------
    463 !
    464 !avant closure: on redéfinit lalim, alim_star_tot et alim_star
    465 !       do ig=1,klon
    466 !       do l=2,lalim(ig)
    467 !       alim_star(ig,l)=entr_star(ig,l)
    468 !       entr_star(ig,l)=0.
    469 !       enddo
    470 !       enddo
    471 
     439print*,'THERM 26JJJ'
     440
     441! Choix de la fonction d'alimentation utilisee pour la fermeture.
     442! Apparemment sans importance
     443      alim_star_clos(:,:)=alim_star(:,:)
     444      alim_star_clos(:,:)=entr_star(:,:)+alim_star(:,:)
     445
     446! Appel avec la version seche
    472447      CALL thermcell_closure(ngrid,nlay,r_aspect,ptimestep,rho,  &
    473      &   zlev,lalim,alim_star,alim_star_tot,zmax_sec,wmax_sec,zmax,wmax,f,lev_out)
     448     &   zlev,lalim,alim_star_clos,f_star,zmax_sec,wmax_sec,f,lev_out)
     449
     450!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
     451! Appel avec les zmax et wmax tenant compte de la condensation
     452! Semble moins bien marcher
     453!     CALL thermcell_closure(ngrid,nlay,r_aspect,ptimestep,rho,  &
     454!    &   zlev,lalim,alim_star,f_star,zmax,wmax,f,lev_out)
     455!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    474456
    475457      if(prt_level.ge.1)print*,'thermcell_closure apres thermcell_closure'
    … …  
    484466! Test valable seulement en 1D mais pas genant
    485467      if (.not. (f0(1).ge.0.) ) then
    486            stop 'Dans thermcell_main'
     468           stop'Dans thermcell_main'
    487469      endif
    488470
    … …  
    511493         fm0=(1.-lambda)*fm+lambda*fm0
    512494         entr0=(1.-lambda)*entr+lambda*entr0
    513 !        detr0=(1.-lambda)*detr+lambda*detr0
     495         detr0=(1.-lambda)*detr+lambda*detr0
    514496      else
    515497         fm0=fm
    … …  
    560542     &    ,fraca,zmax  &
    561543     &    ,zu,zv,pduadj,pdvadj,zua,zva,lev_out)
    562 !IM 050508    &    ,zu,zv,pduadj,pdvadj,zua,zva,igout,lev_out)
     544
    563545      else
    564546
    … …  
    636618!
    637619      if (prt_level.ge.10) print*,'14f OK convect8 ig,l,zha zh zpspsk ',ig,l,zha(ig,l),zh(ig,l),zpspsk(ig,l)
    638             thetath2(ig,l)=zf2*(zha(ig,l)-zh(ig,l)/zpspsk(ig,l))**2
     620            thetath2(ig,l)=zf2*(ztla(ig,l)-zthl(ig,l))**2
    639621            if(zw2(ig,l).gt.1.e-10) then
    640622             wth2(ig,l)=zf2*(zw2(ig,l))**2
    … …  
    651633         enddo
    652634      enddo
     635!calcul des flux: q, thetal et thetav
     636      do l=1,nlay
     637         do ig=1,ngrid
     638      wq(ig,l)=fraca(ig,l)*zw2(ig,l)*(zqta(ig,l)*1000.-po(ig,l)*1000.)
     639      wthl(ig,l)=fraca(ig,l)*zw2(ig,l)*(ztla(ig,l)-zthl(ig,l))
     640      wthv(ig,l)=fraca(ig,l)*zw2(ig,l)*(ztva(ig,l)-ztv(ig,l))
     641         enddo
     642      enddo
    653643!calcul de ale_bl et alp_bl
    654 !pour le calcul d'une valeur intégrée entre la surface et lmax
     644!pour le calcul d'une valeur integree entre la surface et lmax
    655645      do ig=1,ngrid
    656646      alp_int(ig)=0.
    … …  
    782772!     print*,'15 OK convect8'
    783773
     774#ifdef wrgrads_thermcell
    784775      if (prt_level.ge.1) print*,'thermcell_main sorties 3D'
    785 #ifdef wrgrads_thermcell
    786776#include "thermcell_out3d.h"
    787777#endif
Note: See TracChangeset for help on using the changeset viewer.