Ignore:
Timestamp:
Nov 2, 2011, 5:32:28 PM (13 years ago)
Author:
acolaitis
Message:

10/10/2011 == AC

*
This commit aims at increasing the thermals speed, especially for large tracer number configurations. The idea behind this commit is to advect non-active conserved variables outside of the sub-timestep of the thermals. Because these variables are not used in thermals computation, we can decouple them:

momentum: can be decoupled because we assume a constant ratio between horizontal velocity in alimentation layer and maximum vertical velocity in the thermals.

tracer: can be decoupled because we do not take condensation of any tracer into account and hence do not liberate latent heat nor form clouds in the thermals.

temperature: cannot be decoupled (of course)
*

D 336 libf/phymars/thermcell_dqupdown.F90
---------------- Deleted and replaced by a simpler version. Notes about downdraft advection are still available from revision 336 of SVN in thermcell_dqupdown.

A 0 libf/phymars/thermcell_dqup.F90
---------------- New upward advection for tracers and momentum in thermals. Several changes are done compared to the old approach:

  • Updraft quantities are not longer computed by making hypothesis on the amount of advected air.
    • In general, the formalism for updraft computation is much simpler and clearer.
  • Tracer tendancies are no longer computed using the conservation equation. Instead, we use the divergence of an approximated turbulent flux of the concerned quantity, where downdraft are also neglected.

M 336 libf/phymars/thermcell_main_mars.F90
---------------- The Main does not call anymore thermcell_dqupdown, which it was doing 2+tracer_number times per subtimestep (140 times per physical step for a 2 tracer config)

M 336 libf/phymars/calltherm_mars.F90
---------------- Entrainment, detrainment and mass-fluxes are recomputed in the sub-timestep loop. Their final value after iterations is used by the new advection routine to compute tracer and momentum fluxes.

* Results

  • Conservation of tracers has been assessed over 1 yr in 1D and found to be comparable to that obtained with the simple convective adjustment. (it actually seems to be better by a factor of 10%!)
  • GCM speed-up is of about 20% compared to the old thermal configuration, for a 2 tracer case.
  • Advection of sharp tracer profiles has been successfully observed, similar to the old method.
Location:
trunk/LMDZ.MARS/libf/phymars
Files:
1 added
2 edited

Legend:

Unmodified
Added
Removed

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

    r336 r337  
    22! $Id: calltherm.F90 1428 2010-09-13 08:43:37Z fairhead $
    33!
    4       subroutine calltherm_mars(dtime,zzlev,zzlay  &
    5      &      ,pplay,paprs,pphi  &
     4      subroutine calltherm_mars(ptimestep,zzlev,zzlay  &
     5     &      ,pplay,pplev,pphi  &
    66     &      ,u_seri,v_seri,t_seri,pq_therm,q2_therm  &
    77     &      ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs,dq2_therm  &
    … …  
    1515#include "dimensions.h"
    1616#include "dimphys.h"
    17 
    18       REAL dtime
     17#include "comcstfi.h"
     18
     19      REAL ptimestep
    1920      LOGICAL logexpr0, logexpr2(ngridmx,nlayermx), logexpr1(ngridmx)
    2021      REAL fact
    21       INTEGER nbptspb
     22      INTEGER nbptspb,iq,l
    2223
    2324      REAL, INTENT(IN) :: zzlay(ngridmx,nlayermx)
    … …  
    2728      REAL t_seri(ngridmx,nlayermx),pq_therm(ngridmx,nlayermx,nqmx)
    2829      REAL q2_therm(ngridmx,nlayermx)
    29       REAL paprs(ngridmx,nlayermx+1)
     30      REAL pplev(ngridmx,nlayermx+1)
    3031      REAL pplay(ngridmx,nlayermx)
    3132      REAL pphi(ngridmx,nlayermx)
    … …  
    4142      real fm_therm(ngridmx,nlayermx+1)
    4243      real entr_therm(ngridmx,nlayermx),detr_therm(ngridmx,nlayermx)
     44      REAL masse(ngridmx,nlayermx)
    4345
    4446!********************************************************
    … …  
    5153      real lmax_real(ngridmx)
    5254      real zmax(ngridmx),zmaxth(ngridmx)
     55      REAL zdz(ngridmx,nlayermx)
     56 
    5357
    5458!nouvelles variables pour la convection
    … …  
    7680      real zbuoyancyEst(ngridmx,nlayermx)
    7781
    78       character (len=20) :: modname='calltherm'
     82      character (len=20) :: modname
    7983      character (len=80) :: abort_message
    8084
    … …  
    107111         call getin("r_aspect_thermals",r_aspect_thermals)
    108112
    109 !         fm_therm(:,:)=0.
    110 !         detr_therm(:,:)=0.
    111 !         entr_therm(:,:)=0.
     113         fm_therm(:,:)=0.
     114         detr_therm(:,:)=0.
     115         entr_therm(:,:)=0.
    112116
    113117         heatFlux(:,:)=0.
    … …  
    120124         lmax_real(:)=0.
    121125
    122          zdt=dtime/REAL(nsplit_thermals)
     126         zdt=ptimestep/REAL(nsplit_thermals)
    123127
    124128         do isplit=1,nsplit_thermals
    … …  
    130134! cas de splitting
    131135
    132 !         zfm_therm(:,:)=0.
    133 !         zentr_therm(:,:)=0.
    134 !         zdetr_therm(:,:)=0.
     136         zfm_therm(:,:)=0.
     137         zentr_therm(:,:)=0.
     138         zdetr_therm(:,:)=0.
    135139!
    136140         zheatFlux(:,:)=0.
    … …  
    153157             CALL thermcell_main_mars(zdt  &
    154158!             CALL thermcell_main_mars_coupled_v2(zdt  &
    155      &      ,pplay,paprs,pphi,zzlev,zzlay  &
     159     &      ,pplay,pplev,pphi,zzlev,zzlay  &
    156160     &      ,u_seri,v_seri,t_seri,pq_therm,q2_therm  &
    157161     &      ,d_u_the,d_v_the,d_t_the,d_q_the,dq2_the  &
    … …  
    165169!  transformation de la derivee en tendance
    166170
    167             d_t_the(:,:)=d_t_the(:,:)*dtime*fact
    168             d_u_the(:,:)=d_u_the(:,:)*fact
    169             d_v_the(:,:)=d_v_the(:,:)*fact
     171            d_t_the(:,:)=d_t_the(:,:)*ptimestep*fact
     172!            d_u_the(:,:)=d_u_the(:,:)*fact
     173!            d_v_the(:,:)=d_v_the(:,:)*fact
    170174!            dq2_the(:,:)=dq2_the(:,:)*fact           
    171175
    172             if (nqmx .ne. 0) then
    173                d_q_the(:,:,:)=d_q_the(:,:,:)*fact
    174             endif
     176!            if (nqmx .ne. 0) then
     177!               d_q_the(:,:,:)=d_q_the(:,:,:)*fact
     178!            endif
    175179
    176180             zmaxth(:)=zmaxth(:)+zmax(:)*fact
    177181             lmax_real(:)=lmax_real(:)+float(lmax(:))*fact
    178 !            fm_therm(:,:)=fm_therm(:,:)  &
    179 !     &      +zfm_therm(:,:)*fact
    180 !            entr_therm(:,:)=entr_therm(:,:)  &
    181 !     &       +zentr_therm(:,:)*fact
    182 !            detr_therm(:,:)=detr_therm(:,:)  &
    183 !     &       +zdetr_therm(:,:)*fact
     182            fm_therm(:,:)=fm_therm(:,:)  &
     183     &      +zfm_therm(:,:)*fact
     184            entr_therm(:,:)=entr_therm(:,:)  &
     185     &       +zentr_therm(:,:)*fact
     186            detr_therm(:,:)=detr_therm(:,:)  &
     187     &       +zdetr_therm(:,:)*fact
    184188
    185189            heatFlux(:,:)=heatFlux(:,:) &
    … …  
    197201     
    198202            d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
    199             d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
    200             d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
    201             d_q_ajs(:,:,:)=d_q_ajs(:,:,:)+d_q_the(:,:,:)
     203!           d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
     204!           d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
     205!            d_q_ajs(:,:,:)=d_q_ajs(:,:,:)+d_q_the(:,:,:)
    202206!            dq2_therm(:,:)=dq2_therm(:,:)+dq2_the(:,:)
    203207!  incrementation des variables meteo
    204208     
    205209            t_seri(:,:) = t_seri(:,:) + d_t_the(:,:)
    206             u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
    207             v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
    208             pq_therm(:,:,:) = pq_therm(:,:,:) + d_q_the(:,:,:)
     210!            u_seri(:,:) = u_seri(:,:) + d_u_the(:,:)
     211!            v_seri(:,:) = v_seri(:,:) + d_v_the(:,:)
     212!            pq_therm(:,:,:) = pq_therm(:,:,:) + d_q_the(:,:,:)
    209213!            q2_therm(:,:) = q2_therm(:,:) + dq2_therm(:,:)
    210214
    … …  
    218222!****************************************************************
    219223
    220 !          do i=1,ngridmx
    221 !             do k=1,nlayermx
    222 !                if (ztla(i,k) .lt. 1.e-10) fraca(i,k) =0.
    223 !               print*,'youpi je sers a quelque chose !'
    224 !             enddo
    225 !          enddo
    226        
    227           DO i=1,ngridmx
    228             hfmax(i)=MAXVAL(heatFlux(i,:)+heatFlux_down(i,:))
    229             wmax(i)=MAXVAL(zw2(i,:))
    230           ENDDO
    231  
    232          lmax(:)=nint(lmax_real(:))
     224! Now that we have computed total entrainment and detrainment, we can
     225! advect u, v, and q in thermals. (theta already advected). We can do
     226! that separatly because u,v,and q are not used in thermcell_main for
     227! any thermals-related computation : they are purely passive.
     228
     229!calcul de la masse
     230      do l=1,nlayermx
     231         masse(:,l)=(pplev(:,l)-pplev(:,l+1))/g
     232      enddo
     233
     234!calcul de l'epaisseur des couches
     235      do l=1,nlayermx
     236         zdz(:,l)=zzlev(:,l+1)-zzlev(:,l)
     237      enddo
     238
     239
     240      modname='momentum'
     241      call thermcell_dqup(ngridmx,nlayermx,ptimestep                &
     242     &      ,fm_therm,entr_therm,detr_therm,  &
     243     &     masse,u_seri,d_u_ajs,modname,zdz)
     244
     245      call thermcell_dqup(ngridmx,nlayermx,ptimestep    &
     246     &       ,fm_therm,entr_therm,detr_therm,  &
     247     &     masse,v_seri,d_v_ajs,modname,zdz)
     248
     249      if (nqmx .ne. 0.) then
     250      modname='tracer'
     251      DO iq=1,nqmx
     252      call thermcell_dqup(ngridmx,nlayermx,ptimestep     &
     253     &     ,fm_therm,entr_therm,detr_therm,  &
     254     &    masse,pq_therm(:,:,iq),d_q_ajs(:,:,iq),modname,zdz)
     255
     256      ENDDO
     257      endif
     258
     259      DO i=1,ngridmx
     260         hfmax(i)=MAXVAL(heatFlux(i,:)+heatFlux_down(i,:))
     261         wmax(i)=MAXVAL(zw2(i,:))
     262      ENDDO
     263
     264      lmax(:)=nint(lmax_real(:))
    233265         
    234266      return

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

    r336 r337  
    3838
    3939      REAL, INTENT(OUT) :: pdtadj(ngridmx,nlayermx)
    40       REAL, INTENT(OUT) :: pduadj(ngridmx,nlayermx)
    41       REAL, INTENT(OUT) :: pdvadj(ngridmx,nlayermx)
    42       REAL, INTENT(OUT) :: pdqadj(ngridmx,nlayermx,nqmx)
     40      REAL :: pduadj(ngridmx,nlayermx)
     41      REAL :: pdvadj(ngridmx,nlayermx)
     42      REAL :: pdqadj(ngridmx,nlayermx,nqmx)
    4343!      REAL, INTENT(OUT) :: pdq2adj(ngridmx,nlayermx)
    4444      REAL :: pdq2adj(ngridmx,nlayermx)
    … …  
    185185      detr(:,:)=0.
    186186      fm(:,:)=0.
    187       zu(:,:)=pu(:,:)
    188       zv(:,:)=pv(:,:)
     187!      zu(:,:)=pu(:,:)
     188!      zv(:,:)=pv(:,:)
    189189      ztv(:,:)=pt(:,:)/zpopsk(:,:)
    190190
    … …  
    13081308!                 gamma(ig,k)=gamma0(ig,k)
    13091309!   On choisit une relaxation quadratique.
    1310                   gamma(ig,k)=gamma0(ig,k)*sqrt(dua(ig,k)**2+dva(ig,k)**2)
     1310                gamma(ig,k)=gamma0(ig,k)*sqrt(dua(ig,k)**2+dva(ig,k)**2)
    13111311                  zua(ig,k)=(fm(ig,k)*zua(ig,k-1)  &
    13121312     &               +(zf2*entr(ig,k)+gamma(ig,k))*zu(ig,k))  &
    … …  
    13771377      else
    13781378
    1379       modname='momentum'
    1380       call thermcell_dqupdown(ngridmx,nlayermx,ptimestep,fm,entr,detr,  &
    1381      &     masse,zu,pduadj,ztvd,fm_down,ztv,modname,lmax)
    1382 
    1383       call thermcell_dqupdown(ngridmx,nlayermx,ptimestep,fm,entr,detr,  &
    1384      &     masse,zv,pdvadj,ztvd,fm_down,ztv,modname,lmax)
     1379!      modname='momentum'
     1380!      call thermcell_dqupdown(ngridmx,nlayermx,ptimestep,fm,entr,detr,  &
     1381!     &     masse,zu,pduadj,ztvd,fm_down,ztv,modname,lmax)
     1382!
     1383!      call thermcell_dqupdown(ngridmx,nlayermx,ptimestep,fm,entr,detr,  &
     1384!     &     masse,zv,pdvadj,ztvd,fm_down,ztv,modname,lmax)
    13851385
    13861386      endif
    … …  
    14001400!------------------------------------------------------------------
    14011401
    1402       if (nqmx .ne. 0.) then
    1403       modname='tracer'
    1404       DO iq=1,nqmx
    1405           call thermcell_dqupdown(ngridmx,nlayermx,ptimestep,fm,entr,detr,  &
    1406           &      masse,pq(:,:,iq),pdqadj(:,:,iq),ztvd,fm_down,ztv,modname,lmax)
    1407 
    1408       ENDDO
    1409       endif
     1402!      if (nqmx .ne. 0.) then
     1403!      modname='tracer'
     1404!      DO iq=1,nqmx
     1405!          call thermcell_dqupdown(ngridmx,nlayermx,ptimestep,fm,entr,detr,  &
     1406!          &      masse,pq(:,:,iq),pdqadj(:,:,iq),ztvd,fm_down,ztv,modname,lmax)
     1407!
     1408!      ENDDO
     1409!      endif
    14101410
    14111411!------------------------------------------------------------------
Note: See TracChangeset for help on using the changeset viewer.