Ignore:
Timestamp:
Jan 11, 2013, 10:19:19 AM (12 years ago)
Author:
Laurent Fairhead
Message:

Version testing basée sur la r1706


Testing release based on r1706

Location:
LMDZ5/branches/testing
Files:
2 edited

Legend:

Unmodified
Added
Removed
  • LMDZ5/branches/testing

  • LMDZ5/branches/testing/libf/phydev/physiq.F90

    r1665 r1707  
    1111     &            , PVteta)
    1212
    13       USE dimphy
    14       USE infotrac
    15       USE comgeomphy
     13      USE dimphy, only : klon,klev
     14      USE infotrac, only : nqtot
     15      USE comgeomphy, only : rlatd
     16      USE comcstphy, only : rg
     17      USE iophy, only : histbeg_phy,histwrite_phy
     18      USE ioipsl, only : getin,histvert,histdef,histend,ymds2ju
     19      USE mod_phys_lmdz_para, only : jj_nb
     20      USE phys_state_var_mod, only : phys_state_var_init
    1621
    1722      IMPLICIT none
    18 !======================================================================
    19 ! Objet: Moniteur general de la physique du modele
    20 !======================================================================
     23#include "dimensions.h"
     24
     25      integer,parameter :: jjmp1=jjm+1-1/jjm
     26      integer,parameter :: iip1=iim+1
    2127!
    22 !  Arguments:
     28! Routine argument:
    2329!
    24 ! nlon----input-I-nombre de points horizontaux
    25 ! nlev----input-I-nombre de couches verticales, doit etre egale a klev
    26 ! debut---input-L-variable logique indiquant le premier passage
    27 ! lafin---input-L-variable logique indiquant le dernier passage
    28 ! jD_cur       -R-jour courant a l'appel de la physique (jour julien)
    29 ! jH_cur       -R-heure courante a l'appel de la physique (jour julien)
    30 ! pdtphys-input-R-pas d'integration pour la physique (seconde)
    31 ! paprs---input-R-pression pour chaque inter-couche (en Pa)
    32 ! pplay---input-R-pression pour le mileu de chaque couche (en Pa)
    33 ! pphi----input-R-geopotentiel de chaque couche (g z) (reference sol)
    34 ! pphis---input-R-geopotentiel du sol
    35 ! presnivs-input_R_pressions approximat. des milieux couches ( en PA)
    36 ! u-------input-R-vitesse dans la direction X (de O a E) en m/s
    37 ! v-------input-R-vitesse Y (de S a N) en m/s
    38 ! t-------input-R-temperature (K)
    39 ! qx------input-R-humidite specifique (kg/kg) et d'autres traceurs
    40 ! d_t_dyn-input-R-tendance dynamique pour "t" (K/s)
    41 ! d_q_dyn-input-R-tendance dynamique pour "q" (kg/kg/s)
    42 ! flxmass_w -input-R- flux de masse verticale
    43 ! d_u-----output-R-tendance physique de "u" (m/s/s)
    44 ! d_v-----output-R-tendance physique de "v" (m/s/s)
    45 ! d_t-----output-R-tendance physique de "t" (K/s)
    46 ! d_qx----output-R-tendance physique de "qx" (kg/kg/s)
    47 ! d_ps----output-R-tendance physique de la pression au sol
    48 !IM
    49 ! PVteta--output-R-vorticite potentielle a des thetas constantes
    50 !======================================================================
    51 #include "dimensions.h"
    52 #include "comcstphy.h"
     30      integer,intent(in) :: nlon ! number of atmospheric colums
     31      integer,intent(in) :: nlev ! number of vertical levels (should be =klev)
     32      real,intent(in) :: jD_cur ! current day number (Julian day)
     33      real,intent(in) :: jH_cur ! current time of day (as fraction of day)
     34      logical,intent(in) :: debut ! signals first call to physics
     35      logical,intent(in) :: lafin ! signals last call to physics
     36      real,intent(in) :: pdtphys ! physics time step (s)
     37      real,intent(in) :: paprs(klon,klev+1) ! interlayer pressure (Pa)
     38      real,intent(in) :: pplay(klon,klev) ! mid-layer pressure (Pa)
     39      real,intent(in) :: pphi(klon,klev) ! geopotential at mid-layer
     40      real,intent(in) :: pphis(klon) ! surface geopotential
     41      real,intent(in) :: presnivs(klev) ! pseudo-pressure (Pa) of mid-layers
     42      integer,parameter :: longcles=20
     43      real,intent(in) :: clesphy0(longcles) ! Not used
     44      real,intent(in) :: u(klon,klev) ! eastward zonal wind (m/s)
     45      real,intent(in) :: v(klon,klev) ! northward meridional wind (m/s)
     46      real,intent(in) :: t(klon,klev) ! temperature (K)
     47      real,intent(in) :: qx(klon,klev,nqtot) ! tracers (.../kg_air)
     48      real,intent(in) :: flxmass_w(klon,klev) ! vertical mass flux
     49      real,intent(out) :: d_u(klon,klev) ! physics tendency on u (m/s/s)
     50      real,intent(out) :: d_v(klon,klev) ! physics tendency on v (m/s/s)
     51      real,intent(out) :: d_t(klon,klev) ! physics tendency on t (K/s)
     52      real,intent(out) :: d_qx(klon,klev,nqtot) ! physics tendency on tracers
     53      real,intent(out) :: d_ps(klon) ! physics tendency on surface pressure
     54      real,intent(in) :: dudyn(iim+1,jjmp1,klev) ! Not used
     55!FH! REAL PVteta(klon,nbteta)
     56!      REAL PVteta(klon,1)
     57      real,intent(in) :: PVteta(klon,3) ! Not used ; should match definition
     58                                        ! in calfis.F
    5359
    54       integer jjmp1
    55       parameter (jjmp1=jjm+1-1/jjm)
    56       integer iip1
    57       parameter (iip1=iim+1)
     60integer,save :: itau=0 ! counter to count number of calls to physics
     61!$OMP THREADPRIVATE(itau)
     62real :: temp_newton(klon,klev)
     63integer :: k
     64logical, save :: first=.true.
     65!$OMP THREADPRIVATE(first)
    5866
    59       INTEGER ivap          ! indice de traceurs pour vapeur d'eau
    60       PARAMETER (ivap=1)
    61       INTEGER iliq          ! indice de traceurs pour eau liquide
    62       PARAMETER (iliq=2)
     67! For I/Os
     68integer :: itau0
     69real :: zjulian
     70real :: dtime
     71integer :: nhori ! horizontal coordinate ID
     72integer,save :: nid_hist ! output file ID
     73!$OMP THREADPRIVATE(nid_hist)
     74integer :: zvertid ! vertical coordinate ID
     75integer,save :: iwrite_phys ! output every iwrite_phys physics step
     76!$OMP THREADPRIVATE(iwrite_phys)
     77real :: t_ops ! frequency of the IOIPSL operations (eg average over...)
     78real :: t_wrt ! frequency of the IOIPSL outputs
    6379
    64 !
    65 !
    66 ! Variables argument:
    67 !
    68       INTEGER nlon
    69       INTEGER nlev
    70       REAL, intent(in):: jD_cur, jH_cur
     80! initializations
     81if (debut) then ! Things to do only for the first call to physics
     82! load initial conditions for physics (including the grid)
     83  call phys_state_var_init() ! some initializations, required before calling phyetat0
     84  call phyetat0("startphy.nc")
    7185
    72       REAL pdtphys
    73       LOGICAL debut, lafin
    74       REAL paprs(klon,klev+1)
    75       REAL pplay(klon,klev)
    76       REAL pphi(klon,klev)
    77       REAL pphis(klon)
    78       REAL presnivs(klev)
    79       REAL znivsig(klev)
    80       real pir
     86! Initialize outputs:
     87  itau0=0
     88  iwrite_phys=1 !default: output every physics timestep
     89  call getin("iwrite_phys",iwrite_phys)
     90  t_ops=pdtphys*iwrite_phys ! frequency of the IOIPSL operation
     91  t_wrt=pdtphys*iwrite_phys ! frequency of the outputs in the file
     92  ! compute zjulian for annee0=1979 and month=1 dayref=1 and hour=0.0
     93  !CALL ymds2ju(annee0, month, dayref, hour, zjulian)
     94  call ymds2ju(1979, 1, 1, 0.0, zjulian)
     95  dtime=pdtphys
     96  call histbeg_phy("histins.nc",itau0,zjulian,dtime,nhori,nid_hist)
     97!$OMP MASTER
     98  ! define vertical coordinate
     99  call histvert(nid_hist,"presnivs","Vertical levels","Pa",klev, &
     100                presnivs,zvertid,'down')
     101  ! define variables which will be written in "histins.nc" file
     102  call histdef(nid_hist,'temperature','Atmospheric temperature','K', &
     103               iim,jj_nb,nhori,klev,1,klev,zvertid,32, &
     104               'inst(X)',t_ops,t_wrt)
     105  call histdef(nid_hist,'u','Eastward Zonal Wind','m/s', &
     106               iim,jj_nb,nhori,klev,1,klev,zvertid,32, &
     107               'inst(X)',t_ops,t_wrt)
     108  call histdef(nid_hist,'v','Northward Meridional Wind','m/s', &
     109               iim,jj_nb,nhori,klev,1,klev,zvertid,32, &
     110               'inst(X)',t_ops,t_wrt)
     111  call histdef(nid_hist,'ps','Surface Pressure','Pa', &
     112               iim,jj_nb,nhori,1,1,1,zvertid,32, &
     113               'inst(X)',t_ops,t_wrt)
     114  ! end definition sequence
     115  call histend(nid_hist)
     116!$OMP END MASTER
     117endif ! of if (debut)
    81118
    82       REAL u(klon,klev)
    83       REAL v(klon,klev)
    84       REAL t(klon,klev),theta(klon,klev)
    85       REAL qx(klon,klev,nqtot)
    86       REAL flxmass_w(klon,klev)
    87       REAL omega(klon,klev) ! vitesse verticale en Pa/s
    88       REAL d_u(klon,klev)
    89       REAL d_v(klon,klev)
    90       REAL d_t(klon,klev)
    91       REAL d_qx(klon,klev,nqtot)
    92       REAL d_ps(klon)
    93       real da(klon,klev),phi(klon,klev,klev),mp(klon,klev)
    94 !IM definition dynamique o_trac dans phys_output_open
    95 !      type(ctrl_out) :: o_trac(nqtot)
    96 !FH! REAL PVteta(klon,nbteta)
    97       REAL PVteta(klon,1)
    98       REAL dudyn(iim+1,jjmp1,klev)
     119! increment counter itau
     120itau=itau+1
    99121
    100     INTEGER        longcles
    101     PARAMETER    ( longcles = 20 )
     122! set all tendencies to zero
     123d_u(1:klon,1:klev)=0.
     124d_v(1:klon,1:klev)=0.
     125d_t(1:klon,1:klev)=0.
     126d_qx(1:klon,1:klev,1:nqtot)=0.
     127d_ps(1:klon)=0.
    102128
    103 real temp_newton(klon,klev)
    104 integer k
    105 logical, save :: first=.true.
    106 
    107       REAL clesphy0( longcles      )
    108 
    109 d_u=0.
    110 d_v=0.
    111 d_t=0.
    112 d_qx=0.
    113 d_ps=0.
    114 
    115      d_u(:,1)=-u(:,1)/86400.
    116      do k=1,klev
    117         temp_newton(:,k)=280.+cos(rlatd(:))*40.-pphi(:,k)/rg*6.e-3
    118         d_t(:,k)=(temp_newton(:,k)-t(:,k))/1.e5
    119      enddo
     129! compute tendencies to return to the dynamics:
     130! "friction" on the first layer
     131d_u(1:klon,1)=-u(1:klon,1)/86400.
     132d_v(1:klon,1)=-v(1:klon,1)/86400.
     133! newtonian relaxation towards temp_newton()
     134do k=1,klev
     135  temp_newton(1:klon,k)=280.+cos(rlatd(1:klon))*40.-pphi(1:klon,k)/rg*6.e-3
     136  d_t(1:klon,k)=(temp_newton(1:klon,k)-t(1:klon,k))/1.e5
     137enddo
    120138
    121139
    122       print*,'COUCOU PHYDEV'
    123       return
    124       end
     140!print*,'PHYDEV: itau=',itau
     141
     142! write some outputs:
     143if (modulo(itau,iwrite_phys)==0) then
     144  call histwrite_phy(nid_hist,.false.,"temperature",itau,t)
     145  call histwrite_phy(nid_hist,.false.,"u",itau,u)
     146  call histwrite_phy(nid_hist,.false.,"v",itau,v)
     147  call histwrite_phy(nid_hist,.false.,"ps",itau,paprs(:,1))
     148endif
     149
     150! if lastcall, then it is time to write "restartphy.nc" file
     151if (lafin) then
     152  call phyredem("restartphy.nc")
     153endif
     154
     155end
Note: See TracChangeset for help on using the changeset viewer.