Changeset 2535 for trunk/LMDZ.VENUS


Ignore:
Timestamp:
Jun 10, 2021, 11:42:13 AM (3 years ago)
Author:
emillour
Message:

Venus GCM:
Update yamada4 to use q2 computed from previous time step (or read from startphy) so that it inow conforms to 1+1=2. Some code tidying in clmain and yamada4 along the way.
VB + EM

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

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.VENUS/libf/phyvenus/clmain.F

    r2048 r2535  
    1212     .                  solsw, sollw, sollwdown, fder,
    1313     .                  rlon, rlat, dx, dy,
     14     .                  q2,
    1415     .                  debut, lafin,
    1516     .                  d_t,d_u,d_v,d_ts,
     
    3233cAA il faudra sortir ces memes champs en leur ajoutant une dimension,
    3334cAA c'est a dire nbsrf (nbre de subsurface).
    34       USE ioipsl
    35       USE interface_surf
    36       use dimphy
     35!      USE ioipsl
     36!      USE interface_surf
     37      use dimphy, only: klon, klev
    3738      use mod_grid_phy_lmdz, only: nbp_lev
    3839      use cpdet_phy_mod, only: t2tpot
    3940      use turb_mod, only :yustar
     41     
     42#ifdef CPP_XIOS     
     43      use xios_output_mod, only: send_xios_field
     44#endif     
     45     
    4046      IMPLICIT none
    4147c======================================================================
     
    5662c dy-----input-R- resolution des mailles en y (m)
    5763c
     64c q2-----inoutput-R- $q^2$ TKE at inter-layers
     65c
    5866c d_t------output-R- le changement pour "t"
    5967c d_u------output-R- le changement pour "u"
     
    6977c======================================================================
    7078c$$$ PB ajout pour soil
    71 #include "dimsoil.h"
    72 #include "iniprint.h"
    73 #include "clesphys.h"
    74 #include "compbl.h"
    75 c
    76       REAL dtime
    77       integer itap
    78       REAL t(klon,klev)
    79       REAL u(klon,klev), v(klon,klev)
    80       REAL paprs(klon,klev+1), pplay(klon,klev), radsol(klon)
    81 ! ADAPTATION GCM POUR CP(T)
    82       real ppk(klon,klev)
    83       REAL rlon(klon), rlat(klon), dx(klon), dy(klon)
    84       REAL d_t(klon, klev)
    85       REAL d_u(klon, klev), d_v(klon, klev)
    86       REAL flux_t(klon,klev)
    87       REAL dflux_t(klon)
    88 
    89       REAL flux_u(klon,klev), flux_v(klon,klev)
    90       REAL cdragh(klon), cdragm(klon)
    91       real rmu0(klon)         ! cosinus de l'angle solaire zenithal
    92       LOGICAL debut, lafin
    93 c
    94       REAL ts(klon)
    95       REAL d_ts(klon)
    96       REAL albe(klon)
    97 C
    98       REAL fder(klon)
    99       REAL sollw(klon), solsw(klon), sollwdown(klon)
     79      include "dimsoil.h"
     80      include "iniprint.h"
     81      include "clesphys.h"
     82      include "compbl.h"
     83c
     84c Parametres d'entree
     85      REAL, INTENT(IN) :: dtime
     86      INTEGER, INTENT(IN) :: itap
     87      REAL, INTENT(IN) :: t(klon,klev)
     88      REAL, INTENT(IN) :: u(klon,klev), v(klon,klev)
     89      REAL, INTENT(IN) :: paprs(klon,klev+1), pplay(klon,klev)
     90! ADAPTATION GCM POUR CP(T)
     91      REAL, INTENT(IN) :: ppk(klon,klev)
     92      REAL, INTENT(IN) :: rlon(klon), rlat(klon), dx(klon), dy(klon)
     93      REAL, INTENT(IN) :: rmu0(klon)         ! cosine of solar zenith angle
     94      LOGICAL, INTENT(IN) :: debut, lafin
     95      REAL, INTENT(IN) :: ts(klon)
     96      REAL, INTENT(IN) :: sollw(klon), solsw(klon), sollwdown(klon)
     97     
     98c Paramètres IN/OUT
     99      REAL, INTENT(INOUT) :: fder(klon)
     100      REAL, INTENT(INOUT) :: flux_u(klon,klev), flux_v(klon,klev)
     101      REAL, INTENT(INOUT) :: radsol(klon)
     102      REAL, INTENT(INOUT) :: q2(klon,klev+1)
     103     
     104c Parametres de sorties
     105      REAL, INTENT(OUT) :: d_t(klon, klev)
     106      REAL, INTENT(OUT) :: d_u(klon, klev), d_v(klon, klev)
     107      REAL, INTENT(OUT) :: flux_t(klon,klev)
     108      REAL, INTENT(OUT) :: dflux_t(klon)
     109      REAL, INTENT(OUT) :: cdragh(klon), cdragm(klon)
     110      REAL, INTENT(OUT) :: d_ts(klon)
     111      REAL, INTENT(OUT) :: albe(klon)
    100112cAA
    101       REAL zcoefh(klon,klev)
    102       REAL zu1(klon)
    103       REAL zv1(klon)
     113      REAL, INTENT(OUT) :: zcoefh(klon,klev)
     114      REAL, INTENT(OUT) :: zu1(klon)
     115      REAL, INTENT(OUT) :: zv1(klon)
    104116cAA
    105117c$$$ PB ajout pour soil
    106       REAL ftsoil(klon,nsoilmx)
     118      REAL, INTENT(INOUT) :: ftsoil(klon,nsoilmx)
     119
     120c======================================================================
     121      EXTERNAL clqh, clvent, coefkz
     122c======================================================================
     123c Parametre locaux
    107124      REAL ytsoil(klon,nsoilmx)
    108 c======================================================================
    109       EXTERNAL clqh, clvent, coefkz
    110 c======================================================================
    111125      REAL yts(klon)
    112126      REAL yalb(klon)
     
    135149      real y_cd_m(klon),y_cd_h(klon)
    136150c
    137 #include "YOMCST.h"
     151      include "YOMCST.h"
    138152      REAL u1lay(klon), v1lay(klon)
    139153      REAL delp(klon,klev)
    140154      INTEGER i, k
    141       INTEGER ni(klon), knon, j
    142      
     155      INTEGER ni(klon), knon, j     
    143156c======================================================================
    144157      REAL zx_alf1, zx_alf2 !valeur ambiante par extrapola.
     
    200213      y_flux_u = 0.0
    201214      y_flux_v = 0.0
     215      ycoefh=0.0
     216      ycoefm=0.0
    202217C$$ PB
    203218      y_dflux_t = 0.0
     
    360375            call yamada4(knon,dtime,rg,rd,ypaprs,yt
    361376     s      ,yzlev,yzlay,yu,yv,yteta
    362      s      ,y_cd_m,ykmm,ykmn,ykmq,yustar,
     377     s      ,y_cd_m,q2,ykmm,ykmn,ykmq,yustar,
    363378     s      iflag_pbl)
    364379         endif
     
    383398      CALL clvent(knon,dtime,yu1,yv1,ycoefm,yt,yv,ypaprs,ypplay,ydelp,
    384399     s            y_d_v,y_flux_v)
    385 
     400 
    386401c pour le couplage
    387402      ytaux = y_flux_u(:,1)
     
    447462      ENDDO
    448463
    449 c --------------------
    450 c TEST!!!!! PAS DE MELANGE PAR TURBULENCE !!!
    451 c       d_u = 0.
    452 c       d_v = 0.
    453 c       flux_u = 0.
    454 c       flux_v = 0.
    455 c --------------------
    456 
    457 c     print*,"y_d_t apres clqh=",y_d_t(klon/2,:)
    458 
    459       RETURN
    460464      END
    461465
     
    478482     s                flux_t, dflux_s)
    479483
    480       USE interface_surf
    481       use dimphy
     484      use interface_surf, only: interfsurf_hq
     485      use dimphy, only: klon, klev
    482486      use mod_grid_phy_lmdz, only: nbp_lon, nbp_lat, nbp_lev
    483487      use cpdet_phy_mod, only: t2tpot,tpot2t,cpdet
     
    488492c Objet: diffusion verticale de "h"
    489493c======================================================================
    490 #include "YOMCST.h"
    491 #include "dimsoil.h"
    492 #include "iniprint.h"
     494      include "YOMCST.h"
     495      include "dimsoil.h"
     496      include "iniprint.h"
    493497
    494498c Arguments:
    495       INTEGER knon
    496       REAL dtime              ! intervalle du temps (s)
    497       REAL u1lay(klon)        ! vitesse u de la 1ere couche (m/s)
    498       REAL v1lay(klon)        ! vitesse v de la 1ere couche (m/s)
    499       REAL coef(klon,klev)    ! le coefficient d'echange (m**2/s)
     499c Parametres d'entree
     500      INTEGER, INTENT(IN) :: knon
     501      REAL, INTENT(IN) :: dtime              ! intervalle du temps (s)
     502      REAL, INTENT(IN) :: u1lay(klon)        ! vitesse u de la 1ere couche (m/s)
     503      REAL, INTENT(IN) :: v1lay(klon)        ! vitesse v de la 1ere couche (m/s)
     504      REAL, INTENT(IN) :: coef(klon,klev)    ! le coefficient d'echange (m**2/s)
    500505c                               multiplie par le cisaillement du
    501506c                               vent (dV/dz); la premiere valeur
    502507c                               indique la valeur de Cdrag (sans unite)
    503       REAL t(klon,klev)       ! temperature (K)
    504       REAL ts(klon)           ! temperature du sol (K)
    505       REAL paprs(klon,klev+1) ! pression a inter-couche (Pa)
    506       REAL pplay(klon,klev)   ! pression au milieu de couche (Pa)
    507 ! ADAPTATION GCM POUR CP(T)
    508       REAL ppk(klon,klev)     ! fonction d'Exner milieu de couche
    509       REAL delp(klon,klev)    ! epaisseur de couche en pression (Pa)
    510       REAL radsol(klon)       ! ray. net au sol (Solaire+IR) W/m2
    511       REAL albedo(klon)       ! albedo de la surface
    512       real rmu0(klon)         ! cosinus de l'angle solaire zenithal
    513       real rlon(klon), rlat(klon), dx(klon), dy(klon)
    514 c
    515       REAL d_t(klon,klev)     ! incrementation de "t"
    516       REAL d_ts(klon)         ! incrementation de "ts"
    517       REAL flux_t(klon,klev)  ! (diagnostic) flux de la chaleur
     508      REAL, INTENT(IN) :: t(klon,klev)       ! temperature (K)
     509      REAL, INTENT(IN) :: ts(klon)           ! temperature du sol (K)
     510      REAL, INTENT(IN) :: paprs(klon,klev+1) ! pression a inter-couche (Pa)
     511      REAL, INTENT(IN) :: pplay(klon,klev)   ! pression au milieu de couche (Pa)
     512! ADAPTATION GCM POUR CP(T)
     513      REAL, INTENT(IN) :: ppk(klon,klev)     ! fonction d'Exner milieu de couche
     514      REAL, INTENT(IN) :: delp(klon,klev)    ! epaisseur de couche en pression (Pa)
     515      REAL, INTENT(INOUT) :: radsol(klon)       ! ray. net au sol (Solaire+IR) W/m2
     516      REAL, INTENT(IN) :: rmu0(klon)         ! cosinus de l'angle solaire zenithal
     517      REAL, INTENT(IN) :: rlon(klon), rlat(klon), dx(klon), dy(klon)
     518      INTEGER, INTENT(IN) :: itime
     519      LOGICAL, INTENT(IN) :: debut, lafin
     520      REAL, INTENT(INOUT) :: fder(klon)
     521      REAL, INTENT(IN) :: taux(klon), tauy(klon)
     522      REAL, INTENT(IN) :: sollw(klon), sollwdown(klon)
     523      REAL, INTENT(IN) :: swnet(klon)
     524     
     525     
     526c$$$C PB ajout pour soil
     527      LOGICAL, INTENT(IN) :: soil_model
     528      REAL, INTENT(IN) :: tsoil(klon, nsoilmx)
     529c
     530c Parametre de sorties
     531      REAL, INTENT(OUT) :: albedo(klon)       ! albedo de la surface
     532      REAL, INTENT(OUT) :: d_t(klon,klev)     ! incrementation de "t"
     533      REAL, INTENT(OUT) :: d_ts(klon)         ! incrementation de "ts"
     534      REAL, INTENT(OUT) :: flux_t(klon,klev)  ! (diagnostic) flux de la chaleur
    518535c                               sensible, flux de Cp*T, positif vers
    519536c                               le bas: j/(m**2 s) c.a.d.: W/m2
    520       REAL dflux_s(klon) ! derivee du flux sensible dF/dTs
    521 c======================================================================
     537      REAL, INTENT(OUT) :: dflux_s(klon) ! derivee du flux sensible dF/dTs
     538c======================================================================
     539c Variables locales
    522540      INTEGER i, k
    523541      REAL zx_ch(klon,klev)
     
    527545      REAL local_h(klon,klev) ! enthalpie potentielle
    528546      REAL local_ts(klon)
     547      REAL swdown(klon)
    529548c======================================================================
    530549c contre-gradient pour la chaleur sensible: Kelvin/metre
     
    534553      REAL zdelz
    535554c======================================================================
    536 #include "compbl.h"
     555      include "compbl.h"
    537556c======================================================================
    538557c Rajout pour l'interface
    539       integer itime
    540       logical debut, lafin
    541558      real zlev1(klon)
    542       real fder(klon), taux(klon), tauy(klon)
    543559      real temp_air(klon)
    544560      real epot_air(klon)
    545561      real tq_cdrag(klon), petAcoef(klon)
    546562      real petBcoef(klon)
    547       real sollw(klon), sollwdown(klon), swnet(klon), swdown(klon)
    548       real p1lay(klon),pkh1(klon)
    549 c$$$C PB ajout pour soil
    550       LOGICAL soil_model
    551       REAL tsoil(klon, nsoilmx)
     563      real p1lay(klon), pkh1(klon)
    552564
    553565! Parametres de sortie
     
    777789c
    778790
    779       RETURN
    780791      END
    781792     
     
    791802     s                  d_ven,flux_v)
    792803
    793       use dimphy
     804      use dimphy, only: klon, klev
    794805      IMPLICIT none
    795806c======================================================================
     
    814825c flux_v---output-R- (diagnostic) flux du vent: (kg m/s)/(m**2 s)
    815826c======================================================================
    816 #include "iniprint.h"
    817       INTEGER knon
    818       REAL dtime
    819       REAL u1lay(klon), v1lay(klon)
    820       REAL coef(klon,klev)
    821       REAL t(klon,klev), ven(klon,klev)
    822       REAL paprs(klon,klev+1), pplay(klon,klev), delp(klon,klev)
    823       REAL d_ven(klon,klev)
    824       REAL flux_v(klon,klev)
    825 c======================================================================
    826 #include "YOMCST.h"
    827 c======================================================================
     827      include "iniprint.h"
     828
     829c Parametres d'entree
     830      INTEGER, INTENT(IN) :: knon
     831      REAL, INTENT(IN) :: dtime
     832      REAL, INTENT(IN) :: u1lay(klon) , v1lay(klon)
     833      REAL, INTENT(IN) :: coef(klon, klev)
     834      REAL, INTENT(IN) :: t(klon, klev), ven(klon, klev)
     835      REAL, INTENT(IN) :: paprs(klon, klev+1), pplay(klon, klev)
     836      REAL, INTENT(IN) :: delp(klon, klev)
     837     
     838c Parametres de sorties
     839      REAL, INTENT(OUT) :: d_ven(klon, klev)
     840      REAL, INTENT(OUT) :: flux_v(klon, klev)
     841c======================================================================
     842      include "YOMCST.h"
     843c======================================================================
     844c Parametres locaux
    828845      INTEGER i, k
    829846      REAL zx_cv(klon,2:klev)
     
    905922      ENDDO
    906923c
    907       RETURN
    908924      END
    909925     
     
    919935     .                  pcfm, pcfh)
    920936
    921       use dimphy
     937      use dimphy, only: klon, klev
    922938      use cpdet_phy_mod, only: cpdet,t2tpot
    923939#ifdef CPP_XIOS     
     
    944960c pcfh-----output-R- coefficients a calculer (chaleur et humidite)
    945961c======================================================================
    946 #include "YOMCST.h"
    947 #include "iniprint.h"
    948 #include "compbl.h"
    949 #include "clesphys.h"
     962      include "YOMCST.h"
     963      include "iniprint.h"
     964      include "compbl.h"
     965      include "clesphys.h"
    950966c
    951967c Arguments:
    952968c
    953       INTEGER knon
    954       REAL ts(klon)
    955       REAL paprs(klon,klev+1), pplay(klon,klev)
    956 ! ADAPTATION GCM POUR CP(T)
    957       real ppk(klon,klev)
    958       REAL u(klon,klev), v(klon,klev), t(klon,klev)
    959 c
    960       REAL pcfm(klon,klev), pcfh(klon,klev)
    961       INTEGER itop(klon)
     969c Parametres d'entree
     970      INTEGER, INTENT(IN) :: knon
     971      REAL, INTENT(IN) :: ts(klon)
     972      REAL, INTENT(IN) :: pplay(klon, klev)
     973      REAL, INTENT(IN) ::paprs(klon, klev+1)
     974! ADAPTATION GCM POUR CP(T)
     975      REAL, INTENT(IN) :: ppk(klon, klev)
     976      REAL, INTENT(IN) :: u(klon, klev), v(klon, klev), t(klon, klev)
     977     
     978c Parametre de sorties
     979      REAL, INTENT(OUT) :: pcfm(klon, klev), pcfh(klon, klev)
    962980c
    963981c Quelques constantes et options:
    964982c
    965       REAL cepdu2, ckap, cb, cc, cd, clam
     983      REAL, PARAMETER :: cepdu2=((1.e-5)**2)
     984c     REAL, PARAMETER :: cepdu2 =(0.1)**2
     985      REAL, PARAMETER :: ckap=0.4
     986      REAL, PARAMETER :: cb=5.0
     987      REAL, PARAMETER :: cc=5.0
     988      REAL, PARAMETER :: cd=5.0
     989      REAL, PARAMETER :: clam=160.0
     990      REAL, PARAMETER :: ric=0.4 ! nombre de Richardson critique
     991      REAL, PARAMETER :: prandtl=0.4
     992      INTEGER isommet ! le sommet de la couche limite
     993     
     994      LOGICAL, PARAMETER :: tvirtu=.TRUE. ! calculer Ri d'une maniere plus performante
     995      LOGICAL, PARAMETER :: opt_ec=.FALSE. ! formule du Centre Europeen dans l'atmosphere
     996     
     997c      REAL cepdu2, ckap, cb, cc, cd, clam
    966998c TEST VENUS
    967999c     PARAMETER (cepdu2 =(0.1)**2)
    968       PARAMETER (cepdu2 =(1.e-5)**2)
    969 
    970       PARAMETER (CKAP=0.4)
    971       PARAMETER (cb=5.0)
    972       PARAMETER (cc=5.0)
    973       PARAMETER (cd=5.0)
    974       PARAMETER (clam=160.0)
    975       REAL ric ! nombre de Richardson critique
    976       PARAMETER(ric=0.4)
    977       REAL prandtl
    978       PARAMETER (prandtl=0.4)
    979       INTEGER isommet ! le sommet de la couche limite
    980 
    981       LOGICAL tvirtu ! calculer Ri d'une maniere plus performante
    982       PARAMETER (tvirtu=.TRUE.)
    983       LOGICAL opt_ec ! formule du Centre Europeen dans l'atmosphere
    984       PARAMETER (opt_ec=.FALSE.)
     1000c      PARAMETER (cepdu2 =(1.e-5)**2)
     1001c      PARAMETER (CKAP=0.4)
     1002c      PARAMETER (cb=5.0)
     1003c      PARAMETER (cc=5.0)
     1004c      PARAMETER (cd=5.0)
     1005c      PARAMETER (clam=160.0)
     1006c      REAL ric ! nombre de Richardson critique
     1007c      PARAMETER(ric=0.4)
     1008c      REAL prandtl
     1009c      PARAMETER (prandtl=0.4)
     1010c      INTEGER isommet ! le sommet de la couche limite
     1011
     1012c      LOGICAL tvirtu ! calculer Ri d'une maniere plus performante
     1013c      PARAMETER (tvirtu=.TRUE.)
     1014c      LOGICAL opt_ec ! formule du Centre Europeen dans l'atmosphere
     1015c      PARAMETER (opt_ec=.FALSE.)
    9851016
    9861017c
    9871018c Variables locales:
    9881019c
     1020      INTEGER itop(klon)
    9891021      INTEGER i, k
    9901022      REAL zgeop(klon,klev)
     
    10011033      REAL z2geomf, zalh2, zalm2, zscfh, zscfm
    10021034cIM
    1003       LOGICAL check
    1004       PARAMETER (check=.false.)
     1035      LOGICAL, PARAMETER :: check=.false.
    10051036c
    10061037c contre-gradient pour la chaleur sensible: Kelvin/metre
     
    10081039      REAL gamh(2:klev)
    10091040c
    1010       LOGICAL appel1er
    1011       SAVE appel1er
     1041      LOGICAL, SAVE :: appel1er = .TRUE.
    10121042c
    10131043c Fonctions thermodynamiques et fonctions d'instabilite
    10141044      REAL fsta, fins, x
    1015       LOGICAL zxli ! utiliser un jeu de fonctions simples
    1016       PARAMETER (zxli=.FALSE.)
     1045      LOGICAL, PARAMETER :: zxli=.FALSE. ! utiliser un jeu de fonctions simples PARAMETER (zxli=.FALSE.)
    10171046c
    10181047      fsta(x) = 1.0 / (1.0+10.0*x*(1+8.0*x))
    10191048      fins(x) = SQRT(1.0-18.0*x)
    1020 c
    1021       DATA appel1er /.TRUE./
    10221049c
    10231050c----------------------
     
    12281255c#endif
    12291256
    1230       RETURN
    12311257      END
    12321258
     
    12391265     .                  pcfm, pcfh)
    12401266
    1241       use dimphy
     1267      use dimphy, only: klon, klev
    12421268      use mod_grid_phy_lmdz, only: nbp_lev
    12431269      use cpdet_phy_mod, only: cpdet
     
    12571283c pcfh-----output-R- coefficients a calculer (chaleur et humidite)
    12581284c======================================================================
    1259 #include "YOMCST.h"
    1260 #include "iniprint.h"
     1285      include "YOMCST.h"
     1286      include "iniprint.h"
    12611287c
    12621288c Arguments:
    12631289c
    1264       INTEGER knon
    1265       REAL paprs(klon,klev+1), pplay(klon,klev)
    1266       REAL t(klon,klev)
    1267 c
    1268       REAL pcfm(klon,klev), pcfh(klon,klev)
     1290      INTEGER,INTENT(IN) :: knon
     1291      REAL,INTENT(IN) :: paprs(klon,klev+1)
     1292      REAL,INTENT(IN) :: pplay(klon,klev)
     1293      REAL,INTENT(IN) :: t(klon,klev)
     1294c
     1295      REAL,INTENT(OUT) :: pcfm(klon,klev), pcfh(klon,klev)
    12691296c
    12701297c Variables locales:
     
    13041331      ENDDO
    13051332c
    1306       RETURN
    13071333      END
    13081334
  • trunk/LMDZ.VENUS/libf/phyvenus/physiq_mod.F

    r2523 r2535  
    12771277     $            paprs,pplay,ppk,radsol,falbe,
    12781278     e            solsw, sollw, sollwdown, fder,
    1279      e            longitude_deg, latitude_deg, dx, dy,   
     1279     e            longitude_deg, latitude_deg, dx, dy,
     1280     &            q2,
    12801281     e            debut, lafin,
    12811282     s            d_t_vdf,d_u_vdf,d_v_vdf,d_ts,
  • trunk/LMDZ.VENUS/libf/phyvenus/yamada4.F

    r2534 r2535  
    33!
    44      SUBROUTINE yamada4(ngrid,dt,g,rconst,plev,temp
    5 c     s   ,zlev,zlay,u,v,teta,cd,q2,km,kn,kq,ustar
    6      s   ,zlev,zlay,u,v,teta,cd,km,kn,kq,ustar
     5     s   ,zlev,zlay,u,v,teta,cd,q2,km,kn,kq,ustar
    76     s   ,iflag_pbl)
    87c.......................................................................
    9       use dimphy
    10       use turb_mod, only: q2,l0
     8      use dimphy, only: klon, klev
     9      use turb_mod, only: l0
     10         
     11#ifdef CPP_XIOS     
     12      use xios_output_mod, only: send_xios_field
     13#endif     
     14   
    1115      IMPLICIT NONE
    1216c.......................................................................
     
    4044
    4145c.......................................................................
    42       REAL dt,g,rconst
    43       real plev(klon,klev+1),temp(klon,klev)
    44       real ustar(klon)
     46      REAL, INTENT(IN) :: dt,g,rconst
     47      REAL, INTENT(IN) :: plev(klon,klev+1),temp(klon,klev)
     48      REAL, INTENT(IN) :: ustar(klon)
     49      REAL, INTENT(INOUT) :: zlev(klon,klev+1)
     50      REAL, INTENT(IN) :: zlay(klon,klev)
     51      REAL, INTENT(IN) :: u(klon,klev)
     52      REAL, INTENT(IN) :: v(klon,klev)
     53      REAL, INTENT(IN) :: teta(klon,klev)
     54      REAL, INTENT(IN) :: cd(klon)
     55      INTEGER, INTENT(IN) :: iflag_pbl,ngrid
     56     
     57      REAL, INTENT(OUT) :: km(klon,klev+1)
     58      REAL, INTENT(OUT) :: kn(klon,klev+1)
     59      REAL, INTENT(OUT) :: kq(klon,klev+1)
     60     
     61      REAL, INTENT(INOUT) :: q2(klon,klev+1)
     62c Variables locales:
    4563      real kmin,qmin,pblhmin(klon),coriol(klon)
    46       REAL zlev(klon,klev+1)
    47       REAL zlay(klon,klev)
    48       REAL u(klon,klev)
    49       REAL v(klon,klev)
    50       REAL teta(klon,klev)
    51       REAL cd(klon)
    5264      REAL qpre
    5365      REAL unsdz(klon,klev)
    5466      REAL unsdzdec(klon,klev+1)
    5567
    56       REAL km(klon,klev+1)
    5768      REAL kmpre(klon,klev+1),tmp2
    5869      REAL mpre(klon,klev+1)
    59       REAL kn(klon,klev+1)
    60       REAL kq(klon,klev+1)
     70     
    6171      real ff(klon,klev+1),delta(klon,klev+1)
    6272      real aa(klon,klev+1),aa0,aa1
    63       integer iflag_pbl,ngrid
    64 
    65 
     73     
    6674      integer nlay,nlev
    6775
    68       logical first
    69       integer ipas
    70       save first,ipas
    71       data first,ipas/.true.,0/
     76      logical,save :: first=.false. ! not neede any more
     77      integer,save :: ipas=0
    7278
    7379
     
    8288      real m2cstat,mcstat,kmcstat
    8389      real l(klon,klev+1)
    84       !real,save,allocatable :: l0(:)
    85 c  ATTENTION! mis ici car j'ai enlevé q2 des arguments...
    86 c   sinon, c'est au-dessus que ça se passe...
    87       !REAL,save,allocatable :: q2(:,:)
    8890
    8991      real sq(klon),sqz(klon),zz(klon,klev+1)
     
    9799      real frif,falpha,fsm
    98100      real fl,zzz,zl0,zq2,zn2
    99 
     101     
    100102c     real rino(klon,klev+1),smyam(klon,klev),styam(klon,klev)
    101103c    s  ,lyam(klon,klev),knyam(klon,klev)
     
    109111     s     ,0.5*sqrt(q2(ig,k))/sqrt(max(n2(ig,k),1.e-10))) ,1.)
    110112
     113   
    111114      if (.not.(iflag_pbl.ge.6.and.iflag_pbl.le.9)) then
    112115           stop'probleme de coherence dans appel a MY'
     
    118121      nlev=klev+1
    119122
    120       if (first) then
    121 !        IF (.not.ALLOCATED(l0)) allocate(l0(klon))
    122 !        IF (.not.ALLOCATED(q2)) allocate(q2(klon,klevp1))
    123 
    124 c (surtout pour k=1, à cause diagnostiques...)
    125         q2 = 0.
    126       endif
     123c      if (first) then
     124c        IF (.not.ALLOCATED(l0)) allocate(l0(klon))
     125c        IF (.not.ALLOCATED(q2)) allocate(q2(klon,klev+1))
     126c      endif
    127127c===================================
    128      
     128       
    129129      ipas=ipas+1
    130       if (0.eq.1.and.first) then
    131       do ig=1,1000
    132          ri=(ig-800.)/500.
    133          if (ri.lt.ric) then
    134             zrif=frif(ri)
    135          else
    136             zrif=rifc
    137          endif
    138          if(zrif.lt.0.16) then
    139             zalpha=falpha(zrif)
    140             zsm=fsm(zrif)
    141          else
    142             zalpha=1.12
    143             zsm=0.085
    144          endif
     130!      if (0.eq.1.and.first) then
     131!      do ig=1,1000
     132!         ri=(ig-800.)/500.
     133!         if (ri.lt.ric) then
     134!            zrif=frif(ri)
     135!         else
     136!            zrif=rifc
     137!         endif
     138!         if(zrif.lt.0.16) then
     139!            zalpha=falpha(zrif)
     140!            zsm=fsm(zrif)
     141!         else
     142!            zalpha=1.12
     143!            zsm=0.085
     144!         endif
    145145c     print*,ri,rif,zalpha,zsm
    146       enddo
    147       endif
    148 
     146c      enddo
     147c      endif
     148     
    149149c.......................................................................
    150150c  les increments verticaux
    151151c.......................................................................
    152152c
    153 c!!!!! allerte !!!!!c
    154 c!!!!! zlev n'est pas declare a nlev !!!!!c
    155 c!!!!! ---->
    156                                                       DO ig=1,ngrid
    157             zlev(ig,nlev)=zlay(ig,nlay)
    158      &             +( zlay(ig,nlay) - zlev(ig,nlev-1) )
    159                                                       ENDDO
    160 c!!!!! <----
    161 c!!!!! allerte !!!!!c
    162 c
    163153      DO k=1,nlay
    164                                                       DO ig=1,ngrid
     154       DO ig=1,ngrid
    165155        unsdz(ig,k)=1.E+0/(zlev(ig,k+1)-zlev(ig,k))
    166                                                       ENDDO
     156       ENDDO
    167157      ENDDO
    168                                                       DO ig=1,ngrid
     158      DO ig=1,ngrid
    169159      unsdzdec(ig,1)=1.E+0/(zlay(ig,1)-zlev(ig,1))
    170                                                       ENDDO
     160      ENDDO
    171161      DO k=2,nlay
    172                                                       DO ig=1,ngrid
     162       DO ig=1,ngrid
    173163        unsdzdec(ig,k)=1.E+0/(zlay(ig,k)-zlay(ig,k-1))
    174                                                      ENDDO
     164       ENDDO
    175165      ENDDO
    176                                                       DO ig=1,ngrid
     166      DO ig=1,ngrid
    177167      unsdzdec(ig,nlay+1)=1.E+0/(zlev(ig,nlay+1)-zlay(ig,nlay))
    178                                                      ENDDO
    179 c
    180 c.......................................................................
    181 
     168      ENDDO
     169c
     170c.......................................................................
     171 
    182172c===================================
    183173c INITIALISATIONS (surtout pour k=1, à cause diagnostiques...)
    184         dz = 0.
    185         m2 = 0.
    186         dtetadz = 0.
    187         n2 = 0.
    188         rif = 0.
    189         alpha = 0.
    190         sm = 0.
    191         zz = 0.
    192         l = 0.
     174        dz(:,:) = 0.
     175        m2(:,:) = 0.
     176        dtetadz(:,:) = 0.
     177        n2(:,:) = 0.
     178        rif(:,:) = 0.
     179        alpha(:,:) = 0.
     180        sm(:,:) = 0.
     181        zz(:,:) = 0.
     182        l(:,:) = 0.
     183        km(:,:) = 0.
     184        kn(:,:) = 0.
     185        kmpre(:,:) = 0.
     186        mpre(:,:) = 0.
     187
    193188c===================================
    194189      do k=2,klev
    195                                                           do ig=1,ngrid
     190       do ig=1,ngrid
    196191         dz(ig,k)=zlay(ig,k)-zlay(ig,k-1)
    197192         m2(ig,k)=((u(ig,k)-u(ig,k-1))**2+(v(ig,k)-v(ig,k-1))**2)
     
    216211c     print*,'RIF L=',k,rif(ig,k),ri*alpha(ig,k)
    217212
    218 
    219                                                           enddo
    220       enddo
    221 
     213       enddo
     214      enddo
    222215
    223216c====================================================================
     
    259252
    260253c     print*,'Fin de l initialisation de q2 et l0'
    261 
    262       endif ! first
    263 
     254      endif ! of if (first.or.iflag_pbl.eq.6)
     255   
    264256c====================================================================
    265257c  Calcul de la longueur de melange.
     
    292284                                                          enddo
    293285      enddo
    294 
    295286c====================================================================
    296287c   Yamada 2.0
     
    402393
    403394      endif ! Fin du cas 8
    404 
    405 c     print*,'OK8'
    406 
    407395c====================================================================
    408396c   Calcul des coefficients de mélange
     
    466454      endif
    467455
    468 c     print*,'YAMADA4 1'
    469 
    470456c   Estimations de w'2 et T'2 d'apres Abdela et McFarlane
    471457
     
    481467c     print*,'OKFIN'
    482468      first=.false.
    483       return
     469     
    484470      end
Note: See TracChangeset for help on using the changeset viewer.