Changeset 1236 for trunk/LMDZ.MARS/libf


Ignore:
Timestamp:
May 5, 2014, 11:38:51 AM (11 years ago)
Author:
aslmd
Message:

MESOSCALE. A necessary complement commit to r1234 where a upgraded interface making use of modules was proposed. Completed the new formulation for module_lmd_driver for newphys with improved interface with both ini/bdy conditions and physical parameterizations. Changed the Registry accordingly. Finished changes about I/O with the LMD physics (see LMDZ.MARS/README). Made all those changes compatible for old interface, and LES runs (checked with test cases), as well as old input files. Changed makemeso to account for full flexibility on changin nx ny ntracers nproc with newphys. Cleaned the now obsolete bits of code used in LMD physics shared with the GCM. ----- Everything is now ready to properly code both restart runs and nesting for mesoscale runs with new physics.

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
1 added
6 deleted
5 edited

Legend:

Unmodified
Added
Removed
  • trunk/LMDZ.MARS/libf/phymars/meso_inc/meso_inc_les.F

    r790 r1236  
    22       
    33             DO ig=1,ngrid
    4 !! sensible heat flux in W/m2
    5 
    6              sensheat(ig) = zflubid(ig)-capcal(ig)*zdtsdif(ig)
     4!!! sensible heat flux in W/m2
     5!
     6!             sensheat(ig) = zflubid(ig)-capcal(ig)*zdtsdif(ig)
    77
    88!! u star in similarity theory in m/s
     
    1616            DO ig=1,ngrid
    1717
    18 ! New SL parametrization, correct formulation for sensheat :
    19 
    20             sensheat(ig) = (pplay(ig,1)/(r*pt(ig,1)))*cpp
    21      &        *sqrt(pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1)
    22      &        + (log(1.+0.7*wstar(ig) + 2.3*wstar(ig)**2))**2)
    23      &        *zcdh(ig)*(tsurf(ig)-zh(ig,1))
     18!! New SL parametrization, correct formulation for sensheat :
     19!
     20!            sensheat(ig) = (pplay(ig,1)/(r*pt(ig,1)))*cpp
     21!     &        *sqrt(pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1)
     22!     &        + (log(1.+0.7*wstar(ig) + 2.3*wstar(ig)**2))**2)
     23!     &        *zcdh(ig)*(tsurf(ig)-zh(ig,1))
    2424
    2525! New SL parametrization, ustar is more accurately computed in vdif_cd :
     
    4040!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!!
    4141!!! LES LES
    42        IF (flag_LES) THEN       
     42       IF (turb_resolved) THEN       
    4343
    4444         write (*,*) '************************************************'
     
    4949         write (*,*) '************************************************'
    5050
    51          DO ig=1,ngrid
    52           wstar(ig)=0.  !! no additional gustiness needed in surface layer (see vdifc.F)
    53           DO l=1,nlayer
    54             zdvdif(ig,l) = 0.
    55             zdudif(ig,l) = 0.
    56             zdhdif(ig,l) = 0.
    57           ENDDO
    58          ENDDO
    5951         IF (lifting .and. doubleq) THEN
    6052         !! lifted dust is injected in the first layer.
  • trunk/LMDZ.MARS/libf/phymars/physiq.F

    r1233 r1236  
    66     $            ,pu,pv,pt,pq
    77     $            ,pw
    8      $            ,pdu,pdv,pdt,pdq,pdpsrf,tracerdyn
    9 #ifdef MESOSCALE
    10      $            ,output_tab2d, output_tab3d, flag_LES
    11 #endif
    12      $            )
     8     $            ,pdu,pdv,pdt,pdq,pdpsrf,tracerdyn)
    139
    1410      use tracer_mod, only: noms, mmol, igcm_co2, igcm_n2,
     
    3228      use dimradmars_mod, only: tauscaling, aerosol,
    3329     &                          dtrad, fluxrad_sky, fluxrad, albedo
    34       use turb_mod, only: q2, wstar, hfmax_th
     30      use turb_mod, only: q2, wstar, ustar, sensibFlux,
     31     &                    zmax_th, hfmax_th, turb_resolved
    3532      use planete_h, only: aphelie, periheli, year_day, peri_day,
    3633     &                     obliquit
     
    3936      use comsoil_h, only: mlayer,layer
    4037      use surfdat_h, only: z0_default
     38      use comm_wrf
    4139#else
    4240      use phyredem, only: physdem0, physdem1
     
    170168!#include "slope.h"
    171169
    172 #ifdef MESOSCALE
    173 #include "wrf_output_2d.h"
    174 #include "wrf_output_3d.h"
    175 !#include "advtrac.h"   !!! this is necessary for tracers (in dyn3d)
    176 #include "meso_inc/meso_inc_var.F"
    177 #endif
    178 
    179170c Arguments :
    180171c -----------
     
    359350c Variables for PBL
    360351      REAL zz1(ngrid)
    361       REAL lmax_th_out(ngrid),zmax_th(ngrid)
     352      REAL lmax_th_out(ngrid)
    362353      REAL pdu_th(ngrid,nlayer),pdv_th(ngrid,nlayer)
    363354      REAL pdt_th(ngrid,nlayer),pdq_th(ngrid,nlayer,nq)
     
    371362      REAL T_out1(ngrid)
    372363      REAL, ALLOCATABLE, DIMENSION(:) :: z_out     ! height of interpolation between z0 and z1 [meters]
    373       REAL ustar(ngrid),tstar(ngrid)  ! friction velocity and friction potential temp
     364      REAL tstar(ngrid)  ! friction velocity and friction potential temp
    374365      REAL L_mo(ngrid),vhf(ngrid),vvv(ngrid)
    375366!      REAL zu2(ngrid)
    376       REAL sensibFlux(ngrid)
    377367
    378368c=======================================================================
     
    812802         enddo
    813803
    814 
    815 #ifdef MESOSCALE
    816       IF (.not.flag_LES) THEN
    817 #endif
    818804c ----------------------
    819805c Treatment of a special case : using new surface layer (Richardson based)
     
    822808c a unit subgrid gustiness. Remember that thermals should be used we using the
    823809c Richardson based surface layer model.
    824         IF ( .not.calltherm .and. callrichsl ) THEN
     810        IF ( .not.calltherm
     811     .       .and. callrichsl
     812     .       .and. .not.turb_resolved) THEN
    825813          DO ig=1, ngrid
    826814             IF (zh(ig,1) .lt. tsurf(ig)) THEN
     
    834822        ENDIF
    835823c ----------------------
    836 #ifdef MESOSCALE
    837       ENDIF
    838 #endif
    839824
    840825         IF (tke_heat_flux .ne. 0.) THEN
     
    851836     $        zdum1,zdum2,zdh,pdq,zflubid,
    852837     $        zdudif,zdvdif,zdhdif,zdtsdif,q2,
    853      &        zdqdif,zdqsdif,wstar,zcdv,zcdh,hfmax_th,sensibFlux
    854 #ifdef MESOSCALE
    855      &        ,flag_LES
    856 #endif
    857      &        )
    858 
     838     &        zdqdif,zdqsdif,wstar,zcdv,zcdh,hfmax_th,sensibFlux)
    859839
    860840#ifdef MESOSCALE
    861841#include "meso_inc/meso_inc_les.F"
    862842#endif
    863          DO l=1,nlayer
     843
     844          DO ig=1,ngrid
     845             zdtsurf(ig)=zdtsurf(ig)+zdtsdif(ig)
     846          ENDDO
     847
     848         IF (.not.turb_resolved) THEN
     849          DO l=1,nlayer
    864850            DO ig=1,ngrid
    865851               pdv(ig,l)=pdv(ig,l)+zdvdif(ig,l)
     
    869855               zdtdif(ig,l)=zdhdif(ig,l)*zpopsk(ig,l) ! for diagnostic only
    870856            ENDDO
    871          ENDDO
    872 
    873           DO ig=1,ngrid
    874              zdtsurf(ig)=zdtsurf(ig)+zdtsdif(ig)
    875857          ENDDO
    876858
    877          if (tracer) then
    878 #ifdef MESOSCALE
    879          IF (.not.flag_LES) THEN
    880 #endif
     859          if (tracer) then
    881860           DO iq=1, nq
    882861            DO l=1,nlayer
     
    891870              ENDDO
    892871           ENDDO
    893 #ifdef MESOSCALE
     872          end if ! of if (tracer)
    894873         ENDIF
    895 #endif
    896          end if ! of if (tracer)
    897874
    898875      ELSE   
     
    901878     s      (fluxrad(ig)+fluxgrd(ig))/capcal(ig)
    902879         ENDDO
    903 #ifdef MESOSCALE
    904          IF (flag_LES) THEN
    905             write(*,*) 'LES mode !'
     880         IF (turb_resolved) THEN
     881            write(*,*) 'Turbulent-resolving mode !'
    906882            write(*,*) 'Please set calldifv to T in callphys.def'
    907883            STOP
    908884         ENDIF
    909 #endif
    910885      ENDIF ! of IF (calldifv)
    911886
     
    914889c   -----------------------------
    915890
    916       if(calltherm) then
     891      if(calltherm .and. .not.turb_resolved) then
    917892 
    918893        call calltherm_interface(ngrid,nlayer,nq,
     
    949924        lmax_th_out(:)=real(lmax_th(:))
    950925
    951         else   !of if calltherm
     926      else   !of if calltherm
    952927        lmax_th(:)=0
    953928        wstar(:)=0.
    954929        hfmax_th(:)=0.
    955930        lmax_th_out(:)=0.
    956         end if
     931      end if
    957932
    958933c-----------------------------------------------------------------------
     
    11431118         END IF   ! of IF (sedimentation)
    11441119       
    1145 c Add lifted dust to tendancies after sedimentation in the LES
    1146 #ifdef MESOSCALE
    1147 #include "meso_inc/meso_inc_lift_les.F"
    1148 #endif
     1120c Add lifted dust to tendancies after sedimentation in the LES (AC)
     1121      IF (turb_resolved) THEN
     1122           DO iq=1, nq
     1123            DO l=1,nlayer
     1124              DO ig=1,ngrid
     1125                 pdq(ig,l,iq)=pdq(ig,l,iq)+ zdqdif(ig,l,iq)
     1126              ENDDO
     1127            ENDDO
     1128           ENDDO
     1129           DO iq=1, nq
     1130              DO ig=1,ngrid
     1131                 dqsurf(ig,iq)=dqsurf(ig,iq) + zdqsdif(ig,iq)
     1132              ENDDO
     1133           ENDDO
     1134      ENDIF
    11491135
    11501136c
     
    18851871!endif of ifndef MESOSCALE
    18861872
    1887 
    1888 
    18891873#ifdef MESOSCALE
    1890       !!!
    1891       !!! OUTPUT FIELDS
    1892       !!!
    1893       !wtsurf(1:ngrid) = tsurf(1:ngrid)    !! surface temperature
    1894       !wco2ice(1:ngrid) = co2ice(1:ngrid)  !! co2 ice
    1895       TAU_lay(:)=tau(:,1)!!true opacity (not a reference like tauref)
    1896       IF (tracer) THEN
    1897       mtot(1:ngrid) = mtot(1:ngrid) * 1.e6 / rho_ice
    1898       icetot(1:ngrid) = icetot(1:ngrid) * 1.e6 / rho_ice
    1899       !! JF
     1874     
     1875      !! see comm_wrf.
     1876      !! not needed when an array is already in a shared module.
     1877      !! --> example : hfmax_th, zmax_th
     1878
     1879      !state  real  HR_SW     ikj   misc  1  -  h  "HR_SW"     "HEATING RATE SW"                 "K/s"
     1880      comm_HR_SW(1:ngrid,1:nlayer) = zdtsw(1:ngrid,1:nlayer)
     1881      !state  real  HR_LW     ikj   misc  1  -  h  "HR_LW"     "HEATING RATE LW"                 "K/s"
     1882      comm_HR_LW(1:ngrid,1:nlayer) = zdtlw(1:ngrid,1:nlayer)
     1883      !state  real  SWDOWNZ    ij   misc  1  -  h  "SWDOWNZ"   "DOWNWARD SW FLUX AT SURFACE"     "W m-2"
     1884      comm_SWDOWNZ(1:ngrid) = fluxsurf_sw_tot(1:ngrid)
     1885      !state  real  TAU_DUST   ij   misc  1  -  h  "TAU_DUST"  "REFERENCE VISIBLE DUST OPACITY"  ""
     1886      comm_TAU_DUST(1:ngrid) = tauref(1:ngrid)
     1887      !state  real  RDUST     ikj   misc  1  -  h  "RDUST"     "DUST RADIUS"                     "m"
     1888      comm_RDUST(1:ngrid,1:nlayer) = rdust(1:ngrid,1:nlayer)
     1889      !state  real  QSURFDUST  ij   misc  1  -  h  "QSURFDUST" "DUST MASS AT SURFACE"            "kg m-2"
    19001890      IF (igcm_dust_mass .ne. 0) THEN
    1901         qsurfdust(1:ngrid) = qsurf(1:ngrid,igcm_dust_mass)
     1891        comm_QSURFDUST(1:ngrid) = qsurf(1:ngrid,igcm_dust_mass)
     1892      ELSE
     1893        comm_QSURFDUST(1:ngrid) = 0.
    19021894      ENDIF
    1903       IF (igcm_h2o_ice .ne. 0) THEN     
    1904         qsurfice(1:ngrid) = qsurf(1:ngrid,igcm_h2o_ice)
    1905         vmr=1.e6 * zq(1:ngrid,1:nlayer,igcm_h2o_ice)
    1906      .           *mmean(1:ngrid,1:nlayer) / mmol(igcm_h2o_ice)
     1895      !state  real  MTOT       ij   misc  1  -  h  "MTOT"      "TOTAL MASS WATER VAPOR in pmic"  "pmic"
     1896      comm_MTOT(1:ngrid) = mtot(1:ngrid) * 1.e6 / rho_ice
     1897      !state  real  ICETOT     ij   misc  1  -  h  "ICETOT"    "TOTAL MASS WATER ICE"            "kg m-2"
     1898      comm_ICETOT(1:ngrid) = icetot(1:ngrid) * 1.e6 / rho_ice
     1899      !state  real  VMR_ICE   ikj   misc  1  -  h  "VMR_ICE"   "VOL. MIXING RATIO ICE"           "ppm"
     1900      IF (igcm_h2o_ice .ne. 0) THEN
     1901        comm_VMR_ICE(1:ngrid,1:nlayer) = 1.e6
     1902     .    * zq(1:ngrid,1:nlayer,igcm_h2o_ice)
     1903     .    * mmean(1:ngrid,1:nlayer) / mmol(igcm_h2o_ice)
     1904      ELSE
     1905        comm_VMR_ICE(1:ngrid,1:nlayer) = 0.
    19071906      ENDIF
    1908       !! Dust quantity integration along the vertical axe
    1909       dustot(:)=0
    1910       IF (igcm_dust_mass .ne. 0) THEN
    1911       do ig=1,ngrid
    1912        do l=1,nlayer
    1913         dustot(ig) = dustot(ig) +
    1914      &               zq(ig,l,igcm_dust_mass)
    1915      &               *  (zplev(ig,l) - zplev(ig,l+1)) / g
    1916        enddo
    1917       enddo
    1918       ENDIF
    1919       ENDIF
    1920       !! TAU water ice as seen by TES
    1921       if (activice) tauTES = taucloudtes
    1922 c AUTOMATICALLY GENERATED FROM REGISTRY
    1923 #include "fill_save.inc"
     1907      !state  real  TAU_ICE    ij   misc  1  -  h  "TAU_ICE"   "CLOUD OD at 825 cm-1 TES"        ""
     1908      comm_TAU_ICE(1:ngrid) = taucloudtes(1:ngrid)
     1909      !state  real  RICE      ikj   misc  1  -  h  "RICE"      "ICE RADIUS"                      "m"
     1910      comm_RICE(1:ngrid,1:nlayer) = rice(1:ngrid,1:nlayer)
    19241911#else
    19251912#ifndef MESOINI
  • trunk/LMDZ.MARS/libf/phymars/turb_mod.F90

    r1224 r1236  
    44  REAL,SAVE,ALLOCATABLE :: q2(:,:)    ! Turbulent Kinetic Energy
    55  REAL,allocatable,SAVE :: l0(:)
     6  REAL,SAVE,ALLOCATABLE :: ustar(:)
    67  REAL,SAVE,ALLOCATABLE :: wstar(:)
    78  REAL,SAVE,ALLOCATABLE :: hfmax_th(:)
     9  REAL,SAVE,ALLOCATABLE :: zmax_th(:)
     10  REAL,SAVE,ALLOCATABLE :: sensibFlux(:)
     11  LOGICAL :: turb_resolved = .false.
     12      ! this is a flag to say 'turbulence is resolved'
     13      ! mostly for LES use. default is FALSE (for GCM and mesoscale)
    814
    915contains
     
    1824    allocate(l0(ngrid))
    1925    allocate(wstar(ngrid))
     26    allocate(ustar(ngrid))
    2027    allocate(hfmax_th(ngrid))
    21  
     28    allocate(zmax_th(ngrid))
     29    allocate(sensibFlux(ngrid))
     30
    2231  end subroutine ini_turb_mod
    2332
  • trunk/LMDZ.MARS/libf/phymars/vdif_cd.F

    r1226 r1236  
    11      SUBROUTINE vdif_cd(ngrid,nlay,pz0,
    2      &           pg,pz,pu,pv,wstar,pts,ph,pcdv,pcdh
    3 #ifdef MESOSCALE
    4      &                ,flag_LES
    5 #endif
    6      &                  )
     2     &           pg,pz,pu,pv,wstar,pts,ph,pcdv,pcdh)
    73      USE comcstfi_h
     4      use turb_mod, only: turb_resolved
    85      IMPLICIT NONE
    96c=======================================================================
     
    5249      REAL, INTENT(IN) :: wstar(ngrid)
    5350      REAL, INTENT(OUT) :: pcdv(ngrid),pcdh(ngrid) ! momentum and heat drag coefficient
    54 
    55 #ifdef MESOSCALE
    56       LOGICAL, INTENT(IN) :: flag_LES     !! pour LES avec isfflx!=0
    57 #endif
    5851
    5952c   Local:
     
    160153          zu2(ig)=pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1)
    161154     &     + (log(1.+0.7*wstar(ig) + 2.3*wstar(ig)**2))**2
    162 #ifdef MESOSCALE
    163           if(flag_LES) then
     155          if(turb_resolved) then
    164156             zu2(ig)=MAX(zu2(ig),1.)
    165157          endif
    166 #endif
    167158!       zu2(ig)=pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) + (0.5*wstar(ig))**2
    168159
  • trunk/LMDZ.MARS/libf/phymars/vdifc.F

    r1226 r1236  
    66     $                pdudif,pdvdif,pdhdif,pdtsrf,pq2,
    77     $                pdqdif,pdqsdif,wstar,zcdv_true,zcdh_true,
    8      $                hfmax,sensibFlux
    9 #ifdef MESOSCALE
    10      &                ,flag_LES
    11 #endif
    12      &                )
     8     $                hfmax,sensibFlux)
     9
    1310      use tracer_mod, only: noms, igcm_dust_mass, igcm_dust_number,
    1411     &                      igcm_dust_submicron, igcm_h2o_vap,
     
    1613      use surfdat_h, only: watercaptag, frost_albedo_threshold, dryness
    1714      USE comcstfi_h
     15      use turb_mod, only: turb_resolved
    1816      IMPLICIT NONE
    1917
     
    156154      REAL,INTENT(OUT) :: sensibFlux(ngrid)
    157155
    158 #ifdef MESOSCALE
    159       LOGICAL flag_LES     !! pour LES avec isfflx!=0
    160 #endif
    161156c    ** un petit test de coherence
    162157c       --------------------------
     
    330325
    331326      CALL vdif_cd(ngrid,nlay,pz0,g,pzlay,pu,pv,wstar,ptsrf,ph
    332      &             ,zcdv_true,zcdh_true
    333 #ifdef MESOSCALE
    334      &                ,flag_LES
    335 #endif
    336      &           )
     327     &             ,zcdv_true,zcdh_true)
    337328
    338329        zu2(:)=pu(:,1)*pu(:,1)+pv(:,1)*pv(:,1)
     
    554545
    555546      llnt(:)=1
    556 #ifdef MESOSCALE
    557       IF (.not.flag_LES) THEN
    558 #endif
     547      IF (.not.turb_resolved) THEN
    559548      IF (callcond) THEN
    560549       DO ig=1,ngrid
     
    570559      ENDIF
    571560
    572 #ifdef MESOSCALE
    573561      ENDIF
    574 #endif
    575562
    576563      DO ig=1,ngrid
     
    717704             end do
    718705           else
     706#endif
    719707            call dustlift(ngrid,nlay,nq,rho,zcdh_true,zcdh,co2ice,
    720708     &                    pdqsdif)
     709#ifndef MESOSCALE
    721710           endif !doubleq.AND.submicron
    722 #else
    723             call dustlift(ngrid,nlay,nq,rho,zcdh_true,zcdh,co2ice,
    724      &                    pdqsdif)
    725711#endif
    726712        else
Note: See TracChangeset for help on using the changeset viewer.