Changeset 290 for trunk/LMDZ.MARS/libf


Ignore:
Timestamp:
Sep 9, 2011, 4:24:40 PM (13 years ago)
Author:
acolaitis
Message:

--- AC 07/09/2011 ---

This is a major update for thermals and richardson layer parametrization. This update stabilizes thermals (further
studies might show that we can reduce the value of nsplit in gcm. To be tested.), and prevent downdrafts from descending into
the surface layer, which was acting as a cold feedback on the thermals. The Richardson surface layer now features
different gustiness coefficients for Richardson, heat and momentum so that u* and t* are correctly represented.

Upcoming updates will change surflayer_interpol.F90 to implement those changes in the interpolation scheme as well.

*
IMPORTANT : several parts of the vdifc code might want to use these new definitions for gustiness, u* and t*. exemple : dust devil routines
that recompute u* ? lifting routines ? TODO !

M 289 libf/phymars/thermcell_main_mars.F90
-----------------> Added iterations to the velocity / buoyancy / entrainment / detrainment

computation to ensure correct convergence. Iteration number is for now set to
4, which is probably too much. This will be changed once several cases are tested.
The minimum is probably 2 iterations.

M 289 libf/phymars/vdifc.F
-----------------> Subgrid gustiness parametrization now uses different gustiness beta coefficients

for heat and momentum. Comparisons with LES at Exomars landing site, matching u*
and teta* suggests values of beta=0.7 for momentum and beta=1.2 for heat, where
the formula for large scale horizontal winds in the first layer is :

U02 = pu2 + pv2 + (beta*wmax_th)2

M 289 libf/phymars/vdif_cd.F
-----------------> LES data suggests that Richardson number distribution during daytime has a very large

standart deviation due to highly negative Richardsons on several gridpoints. As a consequence,
the mean Richardson in daytime in the LES is much more negative than in LES. In the gcm,
parametrized subgrid turbulence prevents such cases, which can be dangerous in nearly unstable conditions.
Hence, we use beta=0.5 instead of one, so that we keep the anti-crazy-hfx function of beta and we increase the
norm of negative Richardsons in general for daytime conditions. This is set in conjonction with
beta settings for heat and momentum in vdifc.

M 289 libf/phymars/meso_inc/meso_inc_les.F
-----------------> HFX and USTM computations now uses the different betas for heat and momentum.

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
4 edited

Legend:

Unmodified
Added
Removed

  • trunk/LMDZ.MARS/libf/phymars/meso_inc/meso_inc_les.F

    r284 r290  
    2020            hfx(ig) = (pplay(ig,1)/(r*pt(ig,1)))*cpp
    2121     &        *sqrt(pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1)
    22      &        + (1.0*wmax_th(ig))**2)
     22     &        + (1.2*wmax_th(ig))**2)
    2323     &        *zcdh(ig)*(tsurf(ig)-zh(ig,1))
    24 
    2524
    2625! New SL parametrization, ust is more accurately computed in vdif_cd :
    2726            ust(ig) = sqrt(zcdv(ig)*
    28      &   (pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) + (1.0*wmax_th(ig))**2)
     27     &   (pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) + (0.7*wmax_th(ig))**2)
    2928     &                     )
    3029

  • trunk/LMDZ.MARS/libf/phymars/thermcell_main_mars.F90

    r190 r290  
    114114      LOGICAL active(ngridmx),activetmp(ngridmx)
    115115      REAL a1,b1,ae,be,ad,bd
     116      INTEGER tic
    116117
    117118! ==========================================
    … …  
    334335               alim_star(ig,1)=MAX((ztv(ig,1)-ztv(ig,2)),0.)  &
    335336     &                       *sqrt(zlev(ig,2))
     337!     &                       /sqrt(zlev(ig,2))
    336338!      &                       *zlev(ig,2)
    337339               lalim(ig)=2
    … …  
    381383! dans un panache conservatif
    382384          f_star(ig,1)=0.
     385         
    383386          f_star(ig,2)=alim_star(ig,1)
     387
    384388          zw2(ig,2)=2.*g*(ztv(ig,1)-ztv(ig,2))/ztv(ig,2)  &
    385389      &                     *(zlev(ig,2)-zlev(ig,1))  &
    386       &                     *0.4*pphi(ig,1)/(pphi(ig,2)-pphi(ig,1))
     390      &                    *0.4*pphi(ig,1)/(pphi(ig,2)-pphi(ig,1))
    387391          w_est(ig,2)=zw2(ig,2)
    388392
    … …  
    420424
    421425!                if(l .lt. lalim(ig)) then
    422 !                ztva_est(ig,l)=(f_star(ig,l)*ztla(ig,l-1)+  &
     426!               ztva_est(ig,l)=(f_star(ig,l)*ztla(ig,l-1)+  &
    423427!     &            alim_star(ig,l)*ztv(ig,l))  &
    424428!     &            /(f_star(ig,l)+alim_star(ig,l))
    … …  
    522526!---------------------------------------------------------------------------
    523527
     528      DO tic=0,3
     529
    524530      activetmp(:)=active(:) .and. f_star(:,l+1)>1.e-10
    525531      do ig=1,ngridmx
    … …  
    547553      enddo
    548554
     555! ================ RECOMPUTE ENTR, DETR, and F FROM NEW W2 ===================
     556
     557
     558      do ig=1,ngridmx
     559        if (activetmp(ig)) then
     560
     561          zw2m=zw2(ig,l+1)
     562          if(zw2m .gt. 0) then
     563          if((a1*(zbuoy(ig,l)/zw2m)-b1) .gt. 0.) then
     564          entr_star(ig,l)=f_star(ig,l)*zdz*  &
     565        &   MAX(0.,ae*(a1*(zbuoy(ig,l)/zw2m)-b1)**be)
     566          else
     567          entr_star(ig,l)=0.
     568          endif
     569
     570          if(zbuoy(ig,l) .gt. 0.) then
     571             if(l .lt. lalim(ig)) then
     572                detr_star(ig,l)=0.
     573             else
     574                 detr_star(ig,l) = f_star(ig,l)*zdz*              &
     575            &  ad
     576             endif
     577          else
     578          detr_star(ig,l)=f_star(ig,l)*zdz*                        &
     579                &     bd*zbuoy(ig,l)/zw2m
     580          endif
     581          else
     582          entr_star(ig,l)=0.
     583          detr_star(ig,l)=0.
     584          endif
     585
     586! En dessous de lalim, on prend le max de alim_star et entr_star pour
     587! alim_star et 0 sinon
     588
     589        if (l.lt.lalim(ig)) then
     590          alim_star(ig,l)=max(alim_star(ig,l),entr_star(ig,l))
     591          entr_star(ig,l)=0.
     592        endif
     593
     594! Calcul du flux montant normalise
     595
     596      f_star(ig,l+1)=f_star(ig,l)+alim_star(ig,l)+entr_star(ig,l)  &
     597     &              -detr_star(ig,l)
     598
     599      endif
     600      enddo
     601
     602      ENDDO
    549603!---------------------------------------------------------------------------
    550604!initialisations pour le calcul de la hauteur du thermique, de l'inversion et de la vitesse verticale max
    … …  
    11211175      do ig=1,ngridmx
    11221176         if(lmax(ig) .gt. 1) then
    1123          do k=1,lmax(ig)
     1177! No downdraft in the very-near surface layer, we begin at k=3
     1178         do k=3,lmax(ig)
    11241179            zdthladj_down(ig,k)=(1./masse(ig,k))*(fm_down(ig,k+1)*ztv(ig,k+1)- &
    11251180     & fm_down(ig,k)*ztv(ig,k)+fm_down(ig,k)*ztvd(ig,k)-fm_down(ig,k+1)*ztvd(ig,k+1))

  • trunk/LMDZ.MARS/libf/phymars/vdif_cd.F

    r284 r290  
    149149!         zu2=pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1)
    150150!         zu2=MAX(pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1),wmax(ig)**2)
    151        zu2(ig)= pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) + (1.*wmax(ig))**2
     151       zu2(ig)=pu(ig,1)*pu(ig,1) + pv(ig,1)*pv(ig,1) + (0.5*wmax(ig))**2
    152152
    153153               ! we add the wmax to simulate

  • trunk/LMDZ.MARS/libf/phymars/vdifc.F

    r284 r290  
    234234     &             ,zcdv_true,zcdh_true)
    235235
    236 
    237         IF (callrichsl .eq. .true.) then
    238            zu2(:)=pu(:,1)*pu(:,1)+pv(:,1)*pv(:,1)
    239      &             +(1.*wmax(:))**2        ! subgrid gustiness added by quadratic interpolation with a coeff beta, here beta=1. (tuned from
    240                                             ! LES comparisons. This parameter is linked to the thermals settings)
     236        zu2(:)=pu(:,1)*pu(:,1)+pv(:,1)*pv(:,1)
     237
     238        IF (callrichsl .eq. .true.) THEN
     239        zcdv(:)=zcdv_true(:)*sqrt(zu2(:)+(0.7*wmax(:))**2)     ! subgrid gustiness added by quadratic interpolation with a coeff beta
     240        zcdh(:)=zcdh_true(:)*sqrt(zu2(:)+(1.2*wmax(:))**2)     ! LES comparisons. This parameter is linked to the thermals settings)
    241241        ELSE
    242            zu2(:)=pu(:,1)*pu(:,1)+pv(:,1)*pv(:,1)
    243         ENDIF       
    244 
    245242        zcdv(:)=zcdv_true(:)*sqrt(zu2(:))     ! 1 / bulk aerodynamic momentum conductance
    246243        zcdh(:)=zcdh_true(:)*sqrt(zu2(:))     ! 1 / bulk aerodynamic heat conductance
    247                                              ! these are the quantities to be looked at when comparing surface layers of different models
     244        ENDIF
     245
    248246
    249247! Some usefull diagnostics for the new surface layer parametrization :
    … …  
    257255!        call WRITEDIAGFI(ngridmx,'u*',
    258256!     &              'friction velocity','m/s',
    259 !     &                         2,sqrt(zcdv_true(:))*sqrt(zu2(:)))
     257!     &             2,sqrt(zcdv_true(:))*sqrt(zu2(:)+(0.7*wmax(:))**2))
    260258!        call WRITEDIAGFI(ngridmx,'T*',
    261259!     &              'friction temperature','K',
    262 !     &          2,-zcdh_true(:)*(ptsrf(:)-ph(:,1))/sqrt(zcdv_true(:)))
     260!     &   2,zcdh_true(:)*sqrt(zu2(:)+(1.2*wmax(:))**2)*
     261!     &   -(ptsrf(:)-ph(:,1))
     262!     & /(sqrt(zcdv_true(:))*sqrt(zu2(:)+(0.7*wmax(:))**2)))
    263263!        call WRITEDIAGFI(ngridmx,'rm-1',
    264264!     &              'aerodyn momentum conductance','m/s',
Note: See TracChangeset for help on using the changeset viewer.