Changeset 2047


Ignore:
Timestamp:
Nov 29, 2018, 4:47:07 PM (6 years ago)
Author:
slebonnois
Message:

SL: VENUS, ajout des modifs apportees par Thomas Navarro pour la parametrisation des ondes de gravite orographiques

Location:
trunk/LMDZ.VENUS/libf/phyvenus
Files:
8 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.VENUS/libf/phyvenus/YOEGWD.h

    r780 r2047  
    66      real    :: GFRCRIT,GKWAKE,GRCRIT,GVCRIT,GKDRAG,GKLIFT
    77      real    :: GHMAX,GRAHILO,GSIGCR,GSSEC,GTSEC,GVSEC
     8      real    :: TAUBS
     9      integer :: LEVBS
    810      COMMON/YOEGWD/ GFRCRIT,GKWAKE,GRCRIT,GVCRIT,GKDRAG,GKLIFT         &
    9      &   ,GHMAX,GRAHILO,GSIGCR,NKTOPG,NTOP,GSSEC,GTSEC,GVSEC
     11     &   ,GHMAX,GRAHILO,GSIGCR,NKTOPG,NTOP,GSSEC,GTSEC,GVSEC            &
     12     &   ,TAUBS,LEVBS
    1013

  • trunk/LMDZ.VENUS/libf/phyvenus/drag_noro.F

    r1530 r2047  
    77     e                   t, u, v,
    88     s                   pulow, pvlow, pustr, pvstr,
    9      s                   d_t, d_u, d_v)
     9     s                   d_t, d_u, d_v,
     10     s                   blustr,blvstr,pnlow,zeff,zbl,
     11     s                   ptau,tau0,knu2,kbreak)
    1012c
    1113      use dimphy
    … …  
    8688      REAL t(nlon,nlev), u(nlon,nlev), v(nlon,nlev)
    8789      REAL d_t(nlon,nlev), d_u(nlon,nlev), d_v(nlon,nlev)
     90      REAL blustr(nlon),blvstr(nlon),pnlow(nlon),zeff(nlon),zbl(nlon)
     91      REAL knu2(nlon),kbreak(nlon)
     92      REAL ztau(klon,klev+1), ptau(klon,klev), tau0(klon)
    8893c
    8994      INTEGER i, k, kgwd,  kdx(nlon), ktest(nlon)
     95      INTEGER ikenvh(nlon)
     96      INTEGER iknu2(nlon)
     97      INTEGER ikbreak(nlon)
    9098c
    9199c LOCAL VARIABLES:
    … …  
    135143      DO i = 1, klon
    136144         zgeom(i,k) = pgeop(i,klev-k+1)/RG
    137         zn2(i,k)   = pn2(i,klev-k+1)
     145        zn2(i,k)   = pn2(i,klev-k+1)
    138146      ENDDO
    139147      ENDDO
    … …  
    147155     .            pmea, pstd, psig, pgam, pthe, ppic,pval,
    148156     .            pulow,pvlow,
    149      .            pdudt,pdvdt,pdtdt)
     157     .            pdudt,pdvdt,pdtdt,
     158     .            blustr,blvstr,pnlow,zeff,ikenvh,
     159     .            ztau,iknu2,ikbreak)
     160
     161      zbl(:) = real(klev-ikenvh(:))
     162      knu2(:) = real(klev-iknu2(:))
     163      kbreak(:) = real(klev-ikbreak(:))
     164      tau0 = ztau(:,klev+1)
     165
    150166C
    151167C COMPUTE INCREMENTS AND STRESS FROM TENDENCIES
    … …  
    153169      DO k = 1, klev
    154170      DO i = 1, klon
     171         ptau(i,klev+1-k) = ztau(i,k)
    155172         d_u(i,klev+1-k) = dtime*pdudt(i,k)
    156173         d_v(i,klev+1-k) = dtime*pdvdt(i,k)

  • trunk/LMDZ.VENUS/libf/phyvenus/gwprofil.F

    r1530 r2047  
    33     *         , kgwd ,kdx  , ktest
    44     *         , kkcrit, kkcrith, kcrit ,  kkenvh, kknu,kknu2
     5     *         , kkbreak
    56     *         , paphm1, prho   , pstab , ptfr , pvph , pri , ptau
    67     *         , pdmod   , pnu   , psig ,pgamma, pstd, ppic,pval)
    … …  
    5051      integer kkcrit(nlon),kkcrith(nlon),kcrit(nlon)
    5152     *       ,kdx(nlon),ktest(nlon)
    52      *       ,kkenvh(nlon),kknu(nlon),kknu2(nlon)
     53     *       ,kkenvh(nlon),kknu(nlon),kknu2(nlon),kkbreak(nlon)
    5354C
    5455      real paphm1(nlon,nlev+1), pstab(nlon,nlev+1),
    … …  
    128129C                 level
    129130C
    130                          
     131
     132
    131133      do 440 jk=klev,1,-1
    132134
    … …  
    212214 531  continue       
    213215
     216c Yo, this is Venus.
     217      do jl=kidia,kfdia
     218        do jk=klev,1,-1
     219          if(ktest(jl).eq.1)then
     220            if(jk.lt.kkbreak(jl)) ptau(jl,jk)=0.0
     221          endif
     222        enddo
     223      enddo
     224
     225                         
     226! Venus: resolve waves
     227      do jk=klev,1,-1
     228      do jl=kidia,kfdia
     229      if(ktest(jl).eq.1)then
     230      ! if surface stress greater than threshold
     231      if (ztau(jl,klev+1) .ge. taubs) then
     232           ! then enforce same stress in the atmosphere up to the predefined level
     233           if  ((jk.gt.levbs)) then
     234             ptau(jl,jk) = ztau(jl,klev+1)
     235           ! and zero above
     236           elseif (jk.le.levbs) then
     237             ptau(jl,jk) = 0.
     238           endif
     239!      else
     240          !if (jk.le.klev-1) ptau(jl,jk) = 0.
     241!          ptau(jl,jk) = 0.
     242      endif
     243      endif
     244      enddo
     245      enddo
     246
     247
    214248
    215249 123   format(i4,1x,20(f6.3,1x))

  • trunk/LMDZ.VENUS/libf/phyvenus/gwstress.F

    r1530 r2047  
    55     *         , prho  , pstab , pvph  , pstd, psig
    66     *         , pmea , ppic , pval  , ptfr  , ptau 
    7      *         , pgeom1 , pgamma , pd1  , pd2   , pdmod , pnu )
     7     *         , pgeom1 , pgamma , pd1  , pd2   , pdmod , pnu
     8     *         , zeff )
    89c
    910c**** *gwstress*
    … …  
    7475      real pmea(nlon),ppic(nlon),pval(nlon)
    7576      real pdmod(nlon)
     77      real zeff(nlon) ! effective height seen by the flow when there is blocking
    7678c
    7779c-----------------------------------------------------------------------
    … …  
    7981c*       0.2   local arrays
    8082c              ------------
    81 c  zeff--real: effective height seen by the flow when there is blocking
    8283
    8384      integer jl
    84       real zeff 
    8585c
    8686c-----------------------------------------------------------------------
    … …  
    101101c
    102102c
     103      zeff = 0.
    103104      do 301 jl=kidia,kfdia
    104105      if(ktest(jl).eq.1) then
    … …  
    106107c  effective mountain height above the blocked flow
    107108 
    108          zeff=ppic(jl)-pval(jl)
     109         zeff(jl)=ppic(jl)-pval(jl)
    109110         if(kkenvh(jl).lt.klev)then
    110          zeff=amin1(GFRCRIT*pvph(jl,klev+1)/sqrt(pstab(jl,klev+1))
    111      c              ,zeff)
     111         zeff(jl)=amin1(GFRCRIT*pvph(jl,klev+1)/sqrt(pstab(jl,klev+1))
     112     c              ,zeff(jl))
    112113         endif
    113114
    … …  
    116117     *     *psig(jl)*pdmod(jl)/4./pstd(jl)
    117118     *     *pvph(jl,klev+1)*sqrt(pstab(jl,klev+1))
    118      *     *zeff**2
     119     *     *zeff(jl)**2
    119120
    120121

  • trunk/LMDZ.VENUS/libf/phyvenus/orodrag.F

    r1530 r2047  
    66c outputs
    77     r                 , pulow,pvlow
    8      r                 , pvom,pvol,pte )
     8     r                 , pvom,pvol,pte
     9     r                 , blustr,blvstr,pnlow,zeff,ikenvh
     10c 3D and temporary outputs
     11     r                 , ztau,iknu2,ikbreak)
    912     
    1013      use dimphy
    … …  
    110113     *      pgeom1(nlon,nlev),pn2m1(nlon,nlev),
    111114     *      papm1(nlon,nlev),
    112      *      paphm1(nlon,nlev+1)
     115     *      paphm1(nlon,nlev+1),
     116     *      pnlow(nlon), ! low level stability
     117     *      blustr(nlon),blvstr(nlon), ! blocked stress
     118     *      zeff(nlon) ! effective height
     119
    113120c
    114121      integer  kdx(nlon),ktest(nlon)
    … …  
    123130     *         iknu(klon),
    124131     *         iknu2(klon),
     132     *         ikbreak(klon),
    125133     *         ikcrit(klon),
    126134     *         ikhlim(klon)
    … …  
    199207     *     ( nlon, nlev , ktest
    200208     *     , ikcrit, ikcrith, icrit, isect, ikhlim, ikenvh,iknu,iknu2
     209     *     , ikbreak
    201210     *     , paphm1, papm1 , pum1   , pvm1 , ptm1 , pgeom1, zstab, pstd
    202211     *     , zrho  , zri   , ztau , zvph , zpsi, zzdep
    203212     *     , pulow, pvlow
    204213     *     , pthe,pgam,pmea,ppic,pval,znu  ,zd1,  zd2,  zdmod )
     214
     215
     216      pnlow(:) = sqrt(zstab(:,klev+1))
    205217
    206218c
    … …  
    221233     *    , zrho  , zstab, zvph  , pstd,  psig, pmea, ppic, pval
    222234     *    , ztfr   , ztau
    223      *    , pgeom1,pgam,zd1,zd2,zdmod,znu)
     235     *    , pgeom1,pgam,zd1,zd2,zdmod,znu,zeff)
    224236
    225237c
    … …  
    236248     *       (  nlon , nlev
    237249     *       , kgwd   , kdx  , ktest
    238      *       , ikcrit, ikcrith, icrit  , ikenvh, iknu
    239      *       ,iknu2 , paphm1, zrho   , zstab , ztfr  , zvph
     250     *       , ikcrit, ikcrith, icrit , ikenvh, iknu
     251     *       ,iknu2 , ikbreak, paphm1, zrho , zstab , ztfr , zvph
    240252     *       , zri   , ztau
    241253 
    … …  
    259271      zdvdt(jl)=0.0
    260272      zdtdt(jl)=0.0
     273      blustr(jl)=0.0
     274      blvstr(jl)=0.0
    261275  510 continue
    262276c
    … …  
    294308      rover=0.25
    295309      zforc=sqrt(zdudt(ji)**2+zdvdt(ji)**2)       
    296       ztend=sqrt(pum1(ji,jk)**2+pvm1(ji,jk)**2)/ztmst                     
     310      ztend=sqrt(pum1(ji,jk)**2+pvm1(ji,jk)**2)/ztmst
    297311
    298312      if(zforc.ge.rover*ztend)then
    … …  
    330344         zdudt(ji)=zdudt(ji)*(zbet/(1.+zbet))
    331345         zdvdt(ji)=zdvdt(ji)*(zbet/(1.+zbet))
     346
     347c output blocked flow stress
     348         blustr(ji)  = blustr(ji)
     349     .              +zdudt(ji)*(paphm1(ji,jk+1)-paphm1(ji,jk))/rg
     350         blvstr(ji)  = blvstr(ji)
     351     .              +zdvdt(ji)*(paphm1(ji,jk+1)-paphm1(ji,jk))/rg
     352
     353
    332354      end if
    333355      pvom(ji,jk)=zdudt(ji)

  • trunk/LMDZ.VENUS/libf/phyvenus/orosetup.F

    r1530 r2047  
    22     *         ( nlon   , nlev  , ktest
    33     *         , kkcrit, kkcrith, kcrit, ksect , kkhlim
    4      *         , kkenvh, kknu  , kknu2
     4     *         , kkenvh, kknu  , kknu2, kkbreak
    55     *         , paphm1, papm1 , pum1 , pvm1, ptm1, pgeom1, pstab, pstd
    66     *         , prho  , pri   , ptau, pvph, ppsi, pzdep
    … …  
    109109      integer kkcrit(nlon),kkcrith(nlon),kcrit(nlon),ksect(nlon),
    110110     *        kkhlim(nlon),ktest(nlon),kkenvh(nlon)
    111 
     111      integer kkbreak(nlon)
    112112c
    113113      real paphm1(nlon,klev+1),papm1(nlon,klev),pum1(nlon,klev),
    … …  
    179179      kknub(jl)   =klev
    180180      kknul(jl)   =klev
     181      kkbreak(jl) =klev + 1
    181182      pgam(jl) =max(pgam(jl),gtsec)
    182183      ll1(jl,klev+1)=.false.
    … …  
    190191      ENDDO
    191192      ENDDO
     193
     194c      VENUS: define break for subcritical stress
     195c      ----------------------------
     196      do jk=klev,ilevh,-1
     197      do jl=kidia,kfdia
     198      if(ktest(jl).eq.1) then
     199      !zhgeo=pgeom1(jl,jk)/rg
     200      !!zhcrit(jl,jk)=ppic(jl)
     201      !zhcrit(jl,jk)=amin1(ppic(jl)-pmea(jl),pmea(jl)-pval(jl))
     202      !ll1(jl,jk)=(zhgeo.gt.zhcrit(jl,jk))
     203      !if(ll1(jl,jk) .neqv. ll1(jl,jk+1)) then
     204      !  kkbreak(jl)=jk
     205      !endif
     206
     207      !if (paphm1(jl,jk) .ge. 7.e6) kkbreak(jl)=jk
     208      !kkbreak(jl) = klev - 2 ! gwd1103
     209      !kkbreak(jl) = klev - 4 ! gwd1104
     210      !kkbreak(jl) = klev - 3 ! gwd1105
     211
     212      endif     
     213      enddo
     214      enddo
     215
    192216c
    193217c*      define top of low level flow
    … …  
    214238      do 2005 jl=kidia,kfdia
    215239      if(ktest(jl).eq.1) then
    216       zhcrit(jl,jk)=ppic(jl)-pmea(jl)
     240!      zhcrit(jl,jk)=ppic(jl)-pmea(jl)
    217241      zhgeo=pgeom1(jl,jk)/rg
    218       ll1(jl,jk)=(zhgeo.gt.zhcrit(jl,jk))
     242!      ll1(jl,jk)=(zhgeo.gt.zhcrit(jl,jk))
     243      ll1(jl,jk)=(zhgeo.gt.0.5*pstd(jl)) ! TN : do not consider outlier peaks
     244!      ll1(jl,jk)=(zhgeo.gt.pstd(jl)) ! TN : do not consider outlier peaks
     245!      ll1(jl,jk)=(zhgeo.gt.2*pstd(jl)) ! TN : do not consider outlier peaks
    219246      if(ll1(jl,jk) .neqv. ll1(jl,jk+1)) then
    220247        kknu2(jl)=jk
    … …  
    224251 2005 continue
    225252 2004 continue
     253
     254!      do 2104 jk=klev,ilevh,-1
     255!      do 2105 jl=kidia,kfdia
     256!      if(ktest(jl).eq.1) then
     257!      zhgeo=pgeom1(jl,jk)/rg
     258!      ll1(jl,jk)=(zhgeo.gt.pstd(jl)) ! TN : do not consider outlier peaks
     259!      if(ll1(jl,jk) .neqv. ll1(jl,jk+1)) then
     260!        kknul(jl)=jk
     261!      endif
     262!      if(.not.ll1(jl,ilevh))kknul(jl)=ilevh
     263!      endif
     264! 2105 continue
     265! 2104 continue
     266
     267
    226268      do 2006 jk=klev,ilevh,-1
    227269      do 2007 jl=kidia,kfdia
    … …  
    244286      kknu2(jl)=min(kknu2(jl),nktopg)
    245287      kknub(jl)=min(kknub(jl),nktopg)
    246       kknul(jl)=klev
     288!      kknul(jl)=klev
    247289      endif
    248290 2010 continue     

  • trunk/LMDZ.VENUS/libf/phyvenus/physiq_mod.F

    r1726 r2047  
    189189      REAL zustrli(klon), zvstrli(klon)
    190190      REAL zustrhi(klon), zvstrhi(klon)
     191      REAL zublstrdr(klon), zvblstrdr(klon)
     192      REAL znlow(klon), zeff(klon)
     193      REAL zbl(klon), knu2(klon),kbreak(nlon)
     194      REAL tau0(klon), ztau(klon,klev)
    191195
    192196c Pour calcul GW drag oro et nonoro: CALCUL de N2:
    … …  
    854858      END DO
    855859
     860c======
     861c GEOP CORRECTION
    856862c
    857863c Ajouter le geopotentiel du sol:
    … …  
    862868      ENDDO
    863869      ENDDO
     870
     871c............................
     872c CETTE CORRECTION VA DE PAIR AVEC DES MODIFS DE LEAPFROG(_p)
     873c ELLE MARCHE A 50 NIVEAUX (si mmean constante...)
     874c MAIS PAS A 78 NIVEAUX (quand mmean varie...)
     875c A ANALYSER PLUS EN DETAIL AVANT D'UTILISER
     876c............................
     877c zphi is recomputed (pphi is not ok if mean molecular mass varies)
     878c with     dphi = RT/mmean d(ln p) [evaluated at interface]
     879
     880c     DO i = 1, klon
     881c       zphi(i,1) = pphis(i) + R*t_seri(i,1)/mmean(i,1)*1000.
     882c    *                *( log(paprs(i,1)) - log(pplay(i,1)) )     
     883c       DO k = 2, klev
     884c        zphi(i,k) = zphi(i,k-1)
     885c    *      + R*500.*(t_seri(i,k)/mmean(i,k)+t_seri(i,k-1)/mmean(i,k-1))
     886c    *          * (log(pplay(i,k-1)) - log(pplay(i,k)))
     887c       ENDDO
     888c     ENDDO
     889c............................
     890c=====
    864891
    865892c   calcul du geopotentiel aux niveaux intercouches
    … …  
    13861413     e             paprs, pplay,ftsol, t_seri)
    13871414
     1415c albedo variations: test for Yeon Joo Lee
     1416c +12% in 4 Vd / increment to increase it for 20 Vd => +80%
     1417c       heat(:,:)=heat(:,:)*(1.+0.12*(rjourvrai+gmtime)/4.)*1.12**4
    13881418
    13891419c       CO2 near infrared absorption
    … …  
    15311561c====================================================================
    15321562c
    1533       if (ok_orodr.or.ok_gw_nonoro) then
     1563c     if (ok_orodr.or.ok_gw_nonoro) then
    15341564
    15351565c  CALCUL DE N2   
    … …  
    15551585       enddo
    15561586
    1557       endif
     1587c     endif
    15581588     
    15591589c ----------------------------ORODRAG
    … …  
    15731603c        igwdim=MAX(1,igwd)
    15741604c
    1575 c A ADAPTER POUR VENUS!!!
     1605c A ADAPTER POUR VENUS!!!  [ TN: c'est fait ! ]
    15761606        CALL drag_noro(klon,klev,dtime,paprs,pplay,pphi,zn2,
    15771607     e                   zmea,zstd, zsig, zgam, zthe,zpic,zval,
    … …  
    15791609     e                   t_seri, u_seri, v_seri,
    15801610     s                   zulow, zvlow, zustrdr, zvstrdr,
    1581      s                   d_t_oro, d_u_oro, d_v_oro)
     1611     s                   d_t_oro, d_u_oro, d_v_oro,
     1612     s                   zublstrdr,zvblstrdr,znlow,zeff,zbl,
     1613     s                   ztau,tau0,knu2,kbreak)
    15821614
    15831615c       print*,"d_u_oro=",d_u_oro(klon/2,:)
    … …  
    15961628         zustrdr = 0.
    15971629         zvstrdr = 0.
     1630         zublstrdr = 0.
     1631         zvblstrdr = 0.
     1632         znlow = 0.
     1633         zeff = 0.
     1634         zbl = 0
     1635         knu2 = 0
     1636         kbreak = 0
     1637         ztau = 0
     1638         tau0 = 0.
    15981639c
    15991640      ENDIF ! fin de test sur ok_orodr
    … …  
    18821923      CALL send_xios_field("mmean",mmean)
    18831924      CALL send_xios_field("rho",rho)
     1925      CALL send_xios_field("BV2",zn2)
    18841926
    18851927      CALL send_xios_field("dudyn",d_u_dyn)

  • trunk/LMDZ.VENUS/libf/phyvenus/sugwd.F

    r1530 r2047  
    140140
    141141c valeurs dans les routines Mars
    142       GKDRAG=0.1
    143       GRAHILO=1.0   
    144       GRCRIT=0.25
    145       GFRCRIT=1.00
    146       GKWAKE=1.0
    147 C
     142c      GKDRAG=0.1
     143c      GRAHILO=1.0   
     144c      GRCRIT=0.25
     145c      GFRCRIT=1.00
     146c      GKWAKE=1.0
     147C
     148C VENUS
     149      GKDRAG=0.5      ! G
     150      GRAHILO=1.0     ! beta - useless
     151      GRCRIT=0.25     ! Ric  - useless
     152      GFRCRIT=1.0     ! Hnc
     153      GKWAKE=1.0      ! Cd
     154      TAUBS=2.0       ! VENUS: stress threshold is 2 Pa
     155      !TAUBS=1.0       ! VENUS: stress threshold is 1 Pa
     156      !TAUBS=0.5       ! VENUS: stress threshold is 0.5 Pa
     157      LEVBS=nlev-9   ! VENUS: level release is 9
     158      !LEVBS=nlev-19   ! VENUS: level release is 19
     159      !LEVBS=nlev-13   ! VENUS: level release is 13
     160
    148161      GKLIFT=0.25
    149162      GVCRIT =0.0
    … …  
    151164      WRITE(UNIT=6,FMT='('' *** SSO essential constants ***'')')
    152165      WRITE(UNIT=6,FMT='('' *** SPECIFIED IN SUGWD ***'')')
    153       WRITE(UNIT=6,FMT='('' Gravity wave ct '',E13.7,'' '')')GKDRAG
    154       WRITE(UNIT=6,FMT='('' Trapped/total wave dag '',E13.7,'' '')')
    155      S      GRAHILO
    156       WRITE(UNIT=6,FMT='('' Critical Richardson   = '',E13.7,'' '')')
     166      WRITE(UNIT=6,FMT='('' Gravity wave ct '',E14.7,'' '')')GKDRAG
     167      WRITE(UNIT=6,FMT='('' Trapped/total wave dag '',E14.7,'' '')')
     168     S                  GRAHILO
     169      WRITE(UNIT=6,FMT='('' Critical Richardson   = '',E14.7,'' '')')
    157170     S                  GRCRIT
    158       WRITE(UNIT=6,FMT='('' Critical Froude'',e13.7)') GFRCRIT
    159       WRITE(UNIT=6,FMT='('' Low level Wake bluff cte'',e13.7)') GKWAKE
    160       WRITE(UNIT=6,FMT='('' Low level lift  cte'',e13.7)') GKLIFT
    161 
     171      WRITE(UNIT=6,FMT='('' Critical Froude'',e14.7)') GFRCRIT
     172      WRITE(UNIT=6,FMT='('' Low level Wake bluff cte'',e14.7)') GKWAKE
     173      WRITE(UNIT=6,FMT='('' Low level lift  cte'',e14.7)') GKLIFT
     174
     175      WRITE(UNIT=6,FMT='('' VENUS: Mountain stress threshold'',E14.7)')
     176     S                  TAUBS
     177      WRITE(UNIT=6,FMT='('' VENUS: Level release'',I5)') nlev - LEVBS
    162178C
    163179C
Note: See TracChangeset for help on using the changeset viewer.