Ignore:
Timestamp:
Jul 15, 2011, 2:55:17 PM (14 years ago)
Author:
aslmd
Message:

MESOSCALE/LMDZ.MARS.new
--> modified to impact last changes

MESOSCALE/LMD_MM_MARS/makemeso
MESOSCALE/LMD_MM_MARS/SRC/WRFV2/call_meso_physiq?.inc
MESOSCALE/LMD_MM_MARS/SRC/WRFV2/call_meso_inifis?.inc
MESOSCALE/LMD_MM_MARS/SRC/WRFV2/phys/module_lmd_driver.F
--> modified to get rid of ecri_phys

and make changes related to meso_physiq and meso_inifis

LMDZ.MARS/libf/phymars
--> see LMDZ.MARS/README

15/07/2011 == AS

  • Modified the mesoscale part so that the previous change by EM does not imply an error in the mesoscale case. More development is needed though to get the "varying z0" capability in the mesoscale model.
  • Worked on versions of meso_physiq and meso_inifis as close as possible to physiq and inifis for more continuity in the process of impacting changes (and even possibly to reach a common version of physiq and inifis).

    The main point is to make the mesoscale significant specific parts

    coded into include files in meso_inc so that meso_physiq and meso_inifis looks very close to physiq and inifis.

    This is completely transparent for GCM users who does not need the

    contents of meso_inc.
  • Slight cosmetic changes to physiq.f and inifis.F --- some of them e.g. to prepare convergence between meso_physiq and physiq
Location:
trunk/MESOSCALE/LMDZ.MARS.new/in_lmdz_mars_newphys
Files:
2 edited

Legend:

Unmodified
Added
Removed

  • trunk/MESOSCALE/LMDZ.MARS.new/in_lmdz_mars_newphys/compile

    r117 r226  
    2929nz=25
    3030##################
    31 #tracers=2
    32 #nx=64
    33 #ny=48
    34 #nz=25
     31tracers=2
     32nx=64
     33ny=48
     34nz=25
    3535##################
    3636

  • trunk/MESOSCALE/LMDZ.MARS.new/in_lmdz_mars_newphys/physiq.F

    r56 r226  
    1       SUBROUTINE physiq(ngrid,nlayer,nq,
    2      $            firstcall,lastcall,
    3      $            pday,ptime,ptimestep,
    4      $            pplev,pplay,pphi,
    5      $            pu,pv,pt,pq,
    6      $            pw,
    7      $            pdu,pdv,pdt,pdq,pdpsrf,tracerdyn)
    8 
     1      SUBROUTINE physiq(
     2     $            ngrid,nlayer,nq
     3     $            ,firstcall,lastcall
     4     $            ,pday,ptime,ptimestep
     5     $            ,pplev,pplay,pphi
     6     $            ,pu,pv,pt,pq
     7     $            ,pw
     8     $            ,pdu,pdv,pdt,pdq,pdpsrf,tracerdyn
     9     $            )
    910
    1011      IMPLICIT NONE
    … …  
    187188                                     !   of the size distribution
    188189c     Albedo of deposited surface ice
    189       REAL, PARAMETER :: alb_surfice = 0.4 ! 0.45
     190      !!REAL, PARAMETER :: alb_surfice = 0.4 ! 0.45
     191      REAL, PARAMETER :: alb_surfice = 0.45 !!TESTS_JB
    190192
    191193      SAVE day_ini, icount
    … …  
    278280      REAL ccn(ngridmx,nlayermx)   ! Cloud condensation nuclei
    279281                                   !   (particules kg-1)
     282      SAVE ccn  !! in case iradia != 1
    280283      real rdust(ngridmx,nlayermx) ! dust geometric mean radius (m)
    281284      real qtot1,qtot2 ! total aerosol mass
    … …  
    298301
    299302      REAL time_phys
     303
     304c Variables from thermal
     305
     306      REAL lmax_th_out(ngridmx),zmax_th(ngridmx)
     307      REAL wmax_th(ngridmx)
     308      REAL ,SAVE :: hfmax_th(ngridmx)
     309      REAL pdu_th(ngridmx,nlayermx),pdv_th(ngridmx,nlayermx)
     310      REAL pdt_th(ngridmx,nlayermx),pdq_th(ngridmx,nlayermx,nqmx)
     311      INTEGER lmax_th(ngridmx)
     312      REAL dtke_th(ngridmx,nlayermx+1)
     313      REAL dummycol(ngridmx)
    300314
    301315c=======================================================================
    … …  
    513527c          Radiative transfer
    514528c          ------------------
     529
    515530           CALL callradite(icount,ngrid,nlayer,nq,zday,zls,pq,albedo,
    516531     $     emis,mu0,pplev,pplay,pt,tsurf,fract,dist_sol,igout,
    … …  
    549564           ENDIF
    550565
    551 
    552 
    553566        ENDIF ! of if(mod(icount-1,iradia).eq.0)
    554567
    … …  
    565578               fluxrad(ig)=fluxrad_sky(ig)-zplanck(ig)
    566579           ENDDO
    567 
    568580
    569581         DO l=1,nlayer
    … …  
    597609c    4. Vertical diffusion (turbulent mixing):
    598610c    -----------------------------------------
    599 c
     611
    600612      IF (calldifv) THEN
    601 
    602613
    603614         DO ig=1,ngrid
    … …  
    612623            enddo
    613624         enddo
    614          
     625
    615626c        Calling vdif (Martian version WITH CO2 condensation)
    616627         CALL vdifc(ngrid,nlayer,nq,co2ice,zpopsk,
    … …  
    633644         ENDDO
    634645
    635          DO ig=1,ngrid
    636             zdtsurf(ig)=zdtsurf(ig)+zdtsdif(ig)
    637          ENDDO
     646          DO ig=1,ngrid
     647             zdtsurf(ig)=zdtsurf(ig)+zdtsdif(ig)
     648          ENDDO
    638649
    639650         if (tracer) then
    … …  
    659670      ENDIF ! of IF (calldifv)
    660671
     672c-----------------------------------------------------------------------
     673c   TEST. Thermals :
     674c HIGHLY EXPERIMENTAL, BEWARE !!
     675c   -----------------------------
     676 
     677      if(calltherm) then
     678 
     679        call calltherm_interface(firstcall,
     680     $ long,lati,zzlev,zzlay,
     681     $ ptimestep,pu,pv,pt,pq,pdu,pdv,pdt,pdq,q2,
     682     $ pplay,pplev,pphi,zpopsk,
     683     $ pdu_th,pdv_th,pdt_th,pdq_th,lmax_th,zmax_th,
     684     $ dtke_th,hfmax_th,wmax_th)
     685 
     686         DO l=1,nlayer
     687           DO ig=1,ngrid
     688              pdu(ig,l)=pdu(ig,l)+pdu_th(ig,l)
     689              pdv(ig,l)=pdv(ig,l)+pdv_th(ig,l)
     690              pdt(ig,l)=pdt(ig,l)+pdt_th(ig,l)
     691              q2(ig,l)=q2(ig,l)+dtke_th(ig,l)*ptimestep
     692           ENDDO
     693        ENDDO
     694 
     695        DO ig=1,ngrid
     696          q2(ig,nlayer+1)=q2(ig,nlayer+1)+dtke_th(ig,nlayer+1)*ptimestep
     697        ENDDO     
     698   
     699        if (tracer) then
     700        DO iq=1,nq
     701         DO l=1,nlayer
     702           DO ig=1,ngrid
     703             pdq(ig,l,iq)=pdq(ig,l,iq)+pdq_th(ig,l,iq)
     704           ENDDO
     705         ENDDO
     706        ENDDO
     707        endif
     708
     709        else   !of if calltherm
     710        lmax_th(:)=0
     711        end if
    661712
    662713c-----------------------------------------------------------------------
    … …  
    676727
    677728         CALL convadj(ngrid,nlayer,nq,ptimestep,
    678      $                pplay,pplev,zpopsk,
     729     $                pplay,pplev,zpopsk,lmax_th,
    679730     $                pu,pv,zh,pq,
    680731     $                pdu,pdv,zdh,pdq,
    681732     $                zduadj,zdvadj,zdhadj,
    682733     $                zdqadj)
     734
    683735
    684736         DO l=1,nlayer
    … …  
    713765     $              co2ice,albedo,emis,
    714766     $              zdtc,zdtsurfc,pdpsrf,zduc,zdvc,zdqc,
    715      $              fluxsurf_sw,zls)
     767     $              fluxsurf_sw,zls)
    716768
    717769         DO l=1,nlayer
    … …  
    10211073      ENDDO
    10221074
    1023 c    Compute surface stress : (NB: z0 is a common in planete.h)
     1075c    Compute surface stress : (NB: z0 is a common in surfdat.h)
    10241076c     DO ig=1,ngrid
    1025 c        cd = (0.4/log(zzlay(ig,1)/z0))**2
     1077c        cd = (0.4/log(zzlay(ig,1)/z0(ig)))**2
    10261078c        zstress(ig) = rho(ig,1)*cd*(zu(ig,1)**2 + zv(ig,1)**2)
    10271079c     ENDDO
    … …  
    10291081c     Sum of fluxes in solar spectral bands (for output only)
    10301082      DO ig=1,ngrid
    1031              fluxtop_sw_tot(ig)=fluxtop_sw(ig,1) + fluxtop_sw(ig,2)
    1032              fluxsurf_sw_tot(ig)=fluxsurf_sw(ig,1) + fluxsurf_sw(ig,2)
     1083             fluxtop_sw_tot(ig)=fluxtop_sw(ig,1) + fluxtop_sw(ig,2)
     1084             fluxsurf_sw_tot(ig)=fluxsurf_sw(ig,1) + fluxsurf_sw(ig,2)
    10331085      ENDDO
    10341086c ******* TEST ******************************************************
    … …  
    10701122
    10711123      IF (ngrid.NE.1) THEN
    1072          print*,'Ls =',zls*180./pi,
    1073      &   ' tauref(700 Pa,lat=0) =',tauref(ngrid/2),
    1074      &   ' tau(Viking1) =',tau(ig_vl1,1)
     1124         print*,'Ls =',zls*180./pi
     1125     &   ,' tauref(700 Pa,lat=0) =',tauref(ngrid/2)
     1126     &   ,' tau(Viking1) =',tau(ig_vl1,1)
    10751127
    10761128
    … …  
    11421194c        which can later be used to make the statistic files of the run:
    11431195c        "stats")          only possible in 3D runs !
    1144 
    11451196         
    11461197         IF (callstats) THEN
    … …  
    12571308c        WRITEDIAGFI can ALSO be called from any other subroutines
    12581309c        for any variables !!
    1259         call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2,
     1310         call WRITEDIAGFI(ngrid,"emis","Surface emissivity","w.m-1",2,
    12601311     &                  emis)
     1312!         call WRITEDIAGFI(ngrid,"pplay","Pressure","Pa",3,zplay)
     1313!         call WRITEDIAGFI(ngrid,"pplev","Pressure","Pa",3,zplev)
    12611314         call WRITEDIAGFI(ngrid,"tsurf","Surface temperature","K",2,
    12621315     &                  tsurf)
    … …  
    12641317        call WRITEDIAGFI(ngrid,"co2ice","co2 ice thickness","kg.m-2",2,
    12651318     &                  co2ice)
     1319c         call WRITEDIAGFI(ngrid,"temp7","temperature in layer 7",
     1320c     &                  "K",2,zt(1,7))
    12661321c         call WRITEDIAGFI(ngrid,"fluxsurf_lw","fluxsurf_lw","W.m-2",2,
    12671322c     &                  fluxsurf_lw)
    … …  
    12731328c     &                  fluxtop_sw_tot)
    12741329         call WRITEDIAGFI(ngrid,"temp","temperature","K",3,zt)
    1275 c        call WRITEDIAGFI(ngrid,"tau","tau"," ",2,tau)
    12761330        call WRITEDIAGFI(ngrid,"u","Zonal wind","m.s-1",3,zu)
    12771331        call WRITEDIAGFI(ngrid,"v","Meridional wind","m.s-1",3,zv)
    … …  
    12791333c         call WRITEDIAGFI(ngrid,"rho","density","none",3,rho)
    12801334c        call WRITEDIAGFI(ngrid,"q2","q2","kg.m-3",3,q2)
    1281 c        call WRITEDIAGFI(ngrid,'Teta','T potentielle','K',3,zh)
     1335!        call WRITEDIAGFI(ngrid,'Teta','T potentielle','K',3,zh)
    12821336c        call WRITEDIAGFI(ngrid,"pressure","Pressure","Pa",3,pplay)
    12831337c        call WRITEDIAGFI(ngrid,"ssurf","Surface stress","N.m-2",2,
    … …  
    12871341c        call WRITEDIAGFI(ngridmx,'lw_htrt','lw heat. rate',
    12881342c    &                   'w.m-2',3,zdtlw)
     1343c        CALL WRITEDIAGFI(ngridmx,'tauTESap',
     1344c     &                         'tau abs 825 cm-1',
     1345c     &                         '',2,tauTES)
    12891346
    12901347!!!!!!!!!!!!!!!!!!!!!!!!SOIL SOIL SOIL
    … …  
    13241381
    13251382            !!!! waterice = q01, voir readmeteo.F90
    1326             call WRITEDIAGFI(ngridmx,'q01',noms(iq),
     1383            call WRITEDIAGFI(ngridmx,'q01',noms(igcm_h2o_ice),
    13271384     &                      'kg/kg',3,
    13281385     &                       zq(1:ngridmx,1:nlayermx,igcm_h2o_ice))
    13291386            !!!! watervapor = q02, voir readmeteo.F90
    1330             call WRITEDIAGFI(ngridmx,'q02',noms(iq),
     1387            call WRITEDIAGFI(ngridmx,'q02',noms(igcm_h2o_vap),
    13311388     &                      'kg/kg',3,
    13321389     &                       zq(1:ngridmx,1:nlayermx,igcm_h2o_vap))
    1333 
    1334 c            call WRITEDIAGFI(ngridmx,'qsurf'//str2,noms(iq),
    1335 c     &                     'kg.m-2',2,qsurf(1,iq))
    1336 
    1337 
    1338 c             CALL WRITEDIAGFI(ngridmx,'mtot',
    1339 c     &                       'total mass of water vapor',
    1340 c     &                       'kg/m2',2,mtot)
    1341 c             CALL WRITEDIAGFI(ngridmx,'icetot',
    1342 c     &                       'total mass of water ice',
    1343 c     &                       'kg/m2',2,icetot)
    1344 cc            vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice)
    1345 cc    &                *mugaz/mmol(igcm_h2o_ice)
    1346 cc            call WRITEDIAGFI(ngridmx,'vmr_h2oice','h2o ice vmr',
    1347 cc    &                       'mol/mol',3,vmr)
    1348 c             CALL WRITEDIAGFI(ngridmx,'reffice',
    1349 c     &                       'Mean reff',
    1350 c     &                       'm',2,rave)
    1351 cc            call WRITEDIAGFI(ngridmx,'rice','Ice particle size',
    1352 cc    &                       'm',3,rice)
    1353 cc            If activice is true, tauTES is computed in aeropacity.F;
    1354 c             if (.not.activice) then
    1355 c               CALL WRITEDIAGFI(ngridmx,'tauTESap',
    1356 c     &                         'tau abs 825 cm-1',
    1357 c     &                         '',2,tauTES)
    1358 c             endif
    1359 c             call WRITEDIAGFI(ngridmx,'h2o_ice_s',
    1360 c     &                       'surface h2o_ice',
    1361 c     &                       'kg.m-2',2,qsurf(1,igcm_h2o_ice))
     1390            !!!! surface waterice qsurf02 (voir readmeteo)
     1391            call WRITEDIAGFI(ngridmx,'qsurf02','surface tracer',
     1392     &                      'kg.m-2',2,
     1393     &                       qsurf(1:ngridmx,igcm_h2o_ice))
     1394
     1395             CALL WRITEDIAGFI(ngridmx,'mtot',
     1396     &                       'total mass of water vapor',
     1397     &                       'kg/m2',2,mtot)
     1398             CALL WRITEDIAGFI(ngridmx,'icetot',
     1399     &                       'total mass of water ice',
     1400     &                       'kg/m2',2,icetot)
     1401c            vmr=zq(1:ngridmx,1:nlayermx,igcm_h2o_ice)
     1402c    &                *mugaz/mmol(igcm_h2o_ice)
     1403c            call WRITEDIAGFI(ngridmx,'vmr_h2oice','h2o ice vmr',
     1404c    &                       'mol/mol',3,vmr)
     1405             CALL WRITEDIAGFI(ngridmx,'reffice',
     1406     &                       'Mean reff',
     1407     &                       'm',2,rave)
     1408c            call WRITEDIAGFI(ngridmx,'rice','Ice particle size',
     1409c    &                       'm',3,rice)
     1410c            If activice is true, tauTES is computed in aeropacity.F;
     1411             if (.not.activice) then
     1412               CALL WRITEDIAGFI(ngridmx,'tauTESap',
     1413     &                         'tau abs 825 cm-1',
     1414     &                         '',2,tauTES)
     1415             endif
     1416             call WRITEDIAGFI(ngridmx,'h2o_ice_s',
     1417     &                       'surface h2o_ice',
     1418     &                       'kg.m-2',2,qsurf(1,igcm_h2o_ice))
    13621419           endif !(water)
    13631420
    … …  
    14011458             call WRITEDIAGFI(ngridmx,'dustq','Dust mass mr',
    14021459     &                        'kg/kg',3,pq(1,1,igcm_dust_mass))
    1403             call WRITEDIAGFI(ngridmx,'dustN','Dust number',
     1460             call WRITEDIAGFI(ngridmx,'dustN','Dust number',
    14041461     &                        'part/kg',3,pq(1,1,igcm_dust_number))
    14051462           else
    … …  
    14191476
    14201477c        ----------------------------------------------------------
     1478c        Outputs of thermals
     1479c        ----------------------------------------------------------
     1480         if (calltherm) then
     1481
     1482!        call WRITEDIAGFI(ngrid,'dtke',
     1483!     &              'tendance tke thermiques','m**2/s**2',
     1484!     &                         3,dtke_th)
     1485!        call WRITEDIAGFI(ngrid,'d_u_ajs',
     1486!     &              'tendance u thermiques','m/s',
     1487!     &                         3,pdu_th*ptimestep)
     1488!        call WRITEDIAGFI(ngrid,'d_v_ajs',
     1489!     &              'tendance v thermiques','m/s',
     1490!     &                         3,pdv_th*ptimestep)
     1491!        if (tracer) then
     1492!        if (nq .eq. 2) then
     1493!        call WRITEDIAGFI(ngrid,'deltaq_th',
     1494!     &              'delta q thermiques','kg/kg',
     1495!     &                         3,ptimestep*pdq_th(:,:,2))
     1496!        endif
     1497!        endif
     1498
     1499        lmax_th_out(:)=real(lmax_th(:))
     1500
     1501        call WRITEDIAGFI(ngridmx,'lmax_th',
     1502     &              'hauteur du thermique','K',
     1503     &                         2,lmax_th_out)
     1504        call WRITEDIAGFI(ngridmx,'hfmax_th',
     1505     &              'maximum TH heat flux','K.m/s',
     1506     &                         2,hfmax_th)
     1507        call WRITEDIAGFI(ngridmx,'wmax_th',
     1508     &              'maximum TH vertical velocity','m/s',
     1509     &                         2,wmax_th)
     1510
     1511         endif
     1512
     1513c        ----------------------------------------------------------
    14211514c        Output in netcdf file "diagsoil.nc" for subterranean
    14221515c          variables (output every "ecritphy", as for writediagfi)
    … …  
    14521545c        CALL writeg1d(ngrid,nlayer,pw,'w','m.s-1')
    14531546
     1547! THERMALS STUFF 1D
     1548
     1549         if(calltherm) then
     1550
     1551        lmax_th_out(:)=real(lmax_th(:))
     1552
     1553        call WRITEDIAGFI(ngridmx,'lmax_th',
     1554     &              'hauteur du thermique','point',
     1555     &                         0,lmax_th_out)
     1556        call WRITEDIAGFI(ngridmx,'hfmax_th',
     1557     &              'maximum TH heat flux','K.m/s',
     1558     &                         0,hfmax_th)
     1559        call WRITEDIAGFI(ngridmx,'wmax_th',
     1560     &              'maximum TH vertical velocity','m/s',
     1561     &                         0,wmax_th)
     1562
     1563
     1564         co2col(:)=0.
     1565         dummycol(:)=0.
     1566         if (tracer) then
     1567         do l=1,nlayermx
     1568           do ig=1,ngrid
     1569             co2col(ig)=co2col(ig)+
     1570     &                  zq(ig,l,1)*(pplev(ig,l)-pplev(ig,l+1))/g
     1571         if (nqmx .gt. 1) then
     1572             iq=2 ! to avoid out-of-bounds spotting by picky compilers
     1573                  ! when gcm is compiled with only one tracer
     1574             dummycol(ig)=dummycol(ig)+
     1575     &                  zq(ig,l,iq)*(pplev(ig,l)-pplev(ig,l+1))/g
     1576         endif
     1577         enddo
     1578         enddo
     1579
     1580         end if
     1581         call WRITEDIAGFI(ngrid,'co2col','integrated co2 mass'          &
     1582     &                                      ,'kg/m-2',0,co2col)
     1583         call WRITEDIAGFI(ngrid,'dummycol','integrated dummy mass'      &
     1584     &                                      ,'kg/m-2',0,dummycol)
     1585         endif
     1586         call WRITEDIAGFI(ngrid,'w','vertical velocity'                 &
     1587     &                              ,'m/s',1,pw)
     1588         call WRITEDIAGFI(ngrid,"q2","q2","kg.m-3",1,q2)
     1589         call WRITEDIAGFI(ngrid,"tsurf","Surface temperature","K",0,
     1590     &                  tsurf)
     1591
     1592         call WRITEDIAGFI(ngrid,"pplay","Pressure","Pa",1,zplay)
     1593         call WRITEDIAGFI(ngrid,"pplev","Pressure","Pa",1,zplev)
    14541594! or output in diagfi.nc (for testphys1d)
    14551595         call WRITEDIAGFI(ngridmx,'ps','Surface pressure','Pa',0,ps)
Note: See TracChangeset for help on using the changeset viewer.