source: trunk/LMDZ.MARS/libf/phymars/calltherm_interface.F90 @ 278

!
! AC 2011-01-05
!
      SUBROUTINE calltherm_interface (firstcall, &
     & long,lati,zzlev,zzlay, &
     & ptimestep,pu,pv,pt,pq,pdu,pdv,pdt,pdq,q2, &
     & pplay,pplev,pphi,zpopsk, &
     & pdu_th,pdv_th,pdt_th,pdq_th,lmax_th,zmax_th,pbl_dtke,hfmax,wmax)

       USE ioipsl_getincom

      implicit none
#include "callkeys.h"
#include "dimensions.h"
#include "dimphys.h"

!--------------------------------------------------------
! Variables d'entree
!--------------------------------------------------------

      REAL, INTENT(IN) :: ptimestep
      REAL, INTENT(IN) :: pplev(ngridmx,nlayermx+1),pplay(ngridmx,nlayermx)
      REAL, INTENT(IN) :: pphi(ngridmx,nlayermx)
      REAL, INTENT(IN) :: pu(ngridmx,nlayermx),pv(ngridmx,nlayermx)
      REAL, INTENT(IN) :: pt(ngridmx,nlayermx),pq(ngridmx,nlayermx,nqmx)
      REAL, INTENT(IN) :: zzlay(ngridmx,nlayermx)
      REAL, INTENT(IN) :: zzlev(ngridmx,nlayermx+1) 
      LOGICAL, INTENT(IN) :: firstcall
      REAL, INTENT(IN) :: pdu(ngridmx,nlayermx),pdv(ngridmx,nlayermx)
      REAL, INTENT(IN) :: pdq(ngridmx,nlayermx,nqmx),pdt(ngridmx,nlayermx)
      REAL, INTENT(IN) :: q2(ngridmx,nlayermx+1)
      REAL, INTENT(IN) :: long(ngridmx),lati(ngridmx)
      REAL, INTENT(IN) :: zpopsk(ngridmx,nlayermx)

!--------------------------------------------------------
! Variables de sortie (ou entree/sortie)
!--------------------------------------------------------

      REAL pdu_th(ngridmx,nlayermx),pdv_th(ngridmx,nlayermx)
      REAL pdt_th(ngridmx,nlayermx),pdq_th(ngridmx,nlayermx,nqmx)
      INTEGER lmax_th(ngridmx)
      REAL zmax_th(ngridmx)
      REAL pbl_dtke(ngridmx,nlayermx+1)

!--------------------------------------------------------
! Variables du thermique
!--------------------------------------------------------
      REAL u_seri(ngridmx,nlayermx), v_seri(ngridmx,nlayermx)
      REAL t_seri(ngridmx,nlayermx)
      REAL d_t_ajs(ngridmx,nlayermx)
      REAL d_u_ajs(ngridmx,nlayermx), d_q_ajs(ngridmx,nlayermx,nqmx)
      REAL d_v_ajs(ngridmx,nlayermx) 
      REAL fm_therm(ngridmx,nlayermx+1), entr_therm(ngridmx,nlayermx)
      REAL detr_therm(ngridmx,nlayermx)
      REAL zw2(ngridmx,nlayermx+1)
      REAL fraca(ngridmx,nlayermx+1)
      REAL ztla(ngridmx,nlayermx)
      REAL q_therm(ngridmx,nlayermx), pq_therm(ngridmx,nlayermx,nqmx)
      REAL dq_therm(ngridmx,nlayermx), dq_thermdown(ngridmx,nlayermx)
      REAL q2_therm(ngridmx,nlayermx), dq2_therm(ngridmx,nlayermx)

      LOGICAL qtransport_thermals,dtke_thermals

      INTEGER l,ig,iq

! Variable de diagnostique : flux de chaleur vertical

      REAL heatFlux(ngridmx,nlayermx)
      REAL heatFlux_down(ngridmx,nlayermx)
      REAL buoyancyOut(ngridmx,nlayermx)
      REAL buoyancyEst(ngridmx,nlayermx)
      REAL hfmax(ngridmx),wmax(ngridmx)

      REAL tstart,tstop

!---------------------------------------------------------
!---------------------------------------------------------
! **********************************************************************
! Thermique
! **********************************************************************

! Initialisation des sorties

      lmax_th(:)=1
      pdu_th(:,:)=0.
      pdv_th(:,:)=0.
      pdt_th(:,:)=0.
      entr_therm(:,:)=0.
      detr_therm(:,:)=0.
      q2_therm(:,:)=0.
      dq2_therm(:,:)=0.
      dq_therm(:,:)=0.
      dq_thermdown(:,:)=0.
      ztla(:,:)=0.
      pbl_dtke(:,:)=0.
      fm_therm(:,:)=0.
      zw2(:,:)=0.
      fraca(:,:)=0.
      if (tracer) then
         pdq_th(:,:,:)=0.
      end if

! Dans le model terrestres, les seri sont des q+dq tendances déja cumulées. Il n'y a donc pas de
! cumulage à l'intérieur de la routine comme dans le model martien. On le fait ici :

            u_seri(:,:)=pu(:,:)+pdu(:,:)*ptimestep
            v_seri(:,:)=pv(:,:)+pdv(:,:)*ptimestep
            t_seri(:,:)=pt(:,:)+pdt(:,:)*ptimestep

            pq_therm(:,:,:)=0.
            qtransport_thermals=.true.
            call getin("qtransport_thermals",qtransport_thermals)
            if(qtransport_thermals) then
            if(tracer) then
            pq_therm(:,:,:)=pq(:,:,:)+pdq(:,:,:)*ptimestep
            endif
            endif

            d_t_ajs(:,:)=0.
            d_u_ajs(:,:)=0.
            d_v_ajs(:,:)=0.
            d_q_ajs(:,:,:)=0.
            heatFlux(:,:)=0.
            heatFlux_down(:,:)=0.
            buoyancyOut(:,:)=0.
            buoyancyEst(:,:)=0.

       dtke_thermals=.false.
       call getin("dtke_thermals",dtke_thermals)
         if(dtke_thermals) then

         DO l=1,nlayermx
              q2_therm(:,l)=0.5*(q2(:,l)+q2(:,l+1))
         ENDDO
         endif

         call cpu_time(tstart)

         CALL calltherm_mars(ptimestep,zzlev,zzlay &
     &      ,pplay,pplev,pphi &
     &      ,u_seri,v_seri,t_seri,pq_therm, q2_therm &
     &      ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs, dq2_therm &
     &      ,fm_therm,entr_therm,detr_therm &
     &      ,lmax_th,zmax_th &
     &      ,zw2,fraca &
     &      ,zpopsk,ztla,heatFlux,heatFlux_down &
     &      ,buoyancyOut,buoyancyEst,hfmax,wmax)
          call cpu_time(tstop)
         print*,'TOTAL elapsed time in thermals : ',tstop-tstart


! Accumulation des  tendances. On n'accumule pas les quantités de traceurs car celle ci n'a pas du changer
! étant donné qu'on ne prends en compte que q_seri de la vap d'eau = 0

! INCREMENTATION : les d_u_ sont des tendances alors que les pdu sont des dérivees, attention !

           pdu_th(:,:)=d_u_ajs(:,:)/ptimestep
           pdv_th(:,:)=d_v_ajs(:,:)/ptimestep
           pdt_th(:,:)=d_t_ajs(:,:)/ptimestep
           if(qtransport_thermals) then
           if(tracer) then
           pdq_th(:,:,:)=d_q_ajs(:,:,:)/ptimestep
           endif
           endif


         DO l=2,nlayermx
              pbl_dtke(:,l)=0.5*(dq2_therm(:,l-1)+dq2_therm(:,l))/ptimestep
         ENDDO

         pbl_dtke(:,1)=0.5*dq2_therm(:,1)/ptimestep
         pbl_dtke(:,nlayermx+1)=0.
!! DIAGNOSTICS
        
        if(outptherm) then
        if (ngridmx .eq. 1) then
        call WRITEDIAGFI(ngridmx,'entr_therm','entrainement thermique',&
     &                       'kg/m-2',1,entr_therm)
        call WRITEDIAGFI(ngridmx,'detr_therm','detrainement thermique',&
     &                       'kg/m-2',1,detr_therm)
        call WRITEDIAGFI(ngridmx,'fm_therm','flux masse thermique',&
     &                       'kg/m-2',1,fm_therm)
        call WRITEDIAGFI(ngridmx,'zw2','vitesse verticale thermique',&
     &                       'm/s',1,zw2)
        call WRITEDIAGFI(ngridmx,'heatFlux_up','heatFlux_updraft',&
     &                       'SI',1,heatFlux)
       call WRITEDIAGFI(ngridmx,'heatFlux_down','heatFlux_downdraft',&
     &                       'SI',1,heatFlux_down)
        call WRITEDIAGFI(ngridmx,'fraca','fraction coverage',&
     &                       'percent',1,fraca)
        call WRITEDIAGFI(ngridmx,'buoyancyOut','buoyancyOut',&
     &                       'm.s-2',1,buoyancyOut)
        call WRITEDIAGFI(ngridmx,'buoyancyEst','buoyancyEst',&
     &                       'm.s-2',1,buoyancyEst)
        call WRITEDIAGFI(ngridmx,'d_t_th',  &
     &         'tendance temp TH','K',1,d_t_ajs)
        call WRITEDIAGFI(ngridmx,'zmax',  &
     &         'pbl height','m',0,zmax_th)
      else

        call WRITEDIAGFI(ngridmx,'entr_therm','entrainement thermique',&
     &                       'kg/m-2',3,entr_therm)
        call WRITEDIAGFI(ngridmx,'detr_therm','detrainement thermique',&
     &                       'kg/m-2',3,detr_therm)
        call WRITEDIAGFI(ngridmx,'fm_therm','flux masse thermique',&
     &                       'kg/m-2',3,fm_therm)
        call WRITEDIAGFI(ngridmx,'zw2','vitesse verticale thermique',&
     &                       'm/s',3,zw2)
        call WRITEDIAGFI(ngridmx,'heatFlux','heatFlux',&
     &                       'SI',3,heatFlux)
        call WRITEDIAGFI(ngridmx,'buoyancyOut','buoyancyOut',&
     &                       'SI',3,buoyancyOut)
        call WRITEDIAGFI(ngridmx,'d_t_th',  &
     &         'tendance temp TH','K',3,d_t_ajs)

      endif
      endif

       END