Changeset 342 for trunk/LMDZ.MARS/libf


Ignore:
Timestamp:
Nov 3, 2011, 3:16:48 PM (13 years ago)
Author:
acolaitis
Message:

M 341 libf/phymars/calltherm_interface.F90
D 341 libf/phymars/calltherm_mars.F90
---------------- Merged calltherm_mars with calltherm_interface for simplicity.

Cleaned up variables and code

M 341 libf/phymars/thermcell_dqup.F90
---------------- dqup now output derivatives and not tendancies

Location:
trunk/LMDZ.MARS/libf/phymars
Files:
1 deleted
2 edited

Legend:

Unmodified
Added
Removed

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

    r319 r342  
    66     & ptimestep,pu,pv,pt,pq,pdu,pdv,pdt,pdq,q2, &
    77     & pplay,pplev,pphi,zpopsk, &
    8      & pdu_th,pdv_th,pdt_th,pdq_th,lmax_th,zmax_th,pbl_dtke,hfmax,wmax)
    9 
    10        !USE ioipsl_getincom
     8     & pdu_th,pdv_th,pdt_th,pdq_th,lmax,zmaxth,pbl_dtke,hfmax,wmax)
     9
     10       USE ioipsl_getincom
    1111
    1212      implicit none
    … …  
    1414#include "dimensions.h"
    1515#include "dimphys.h"
    16 
    17 !--------------------------------------------------------
    18 ! Variables d'entree
     16#include "comcstfi.h"
     17
     18!--------------------------------------------------------
     19! Input Variables
    1920!--------------------------------------------------------
    2021
    … …  
    3435
    3536!--------------------------------------------------------
    36 ! Variables de sortie (ou entree/sortie)
    37 !--------------------------------------------------------
    38 
    39       REAL pdu_th(ngridmx,nlayermx),pdv_th(ngridmx,nlayermx)
    40       REAL pdt_th(ngridmx,nlayermx),pdq_th(ngridmx,nlayermx,nqmx)
    41       INTEGER lmax_th(ngridmx)
    42       REAL zmax_th(ngridmx)
    43       REAL pbl_dtke(ngridmx,nlayermx+1)
    44 
    45 !--------------------------------------------------------
    46 ! Variables du thermique
    47 !--------------------------------------------------------
    48       REAL u_seri(ngridmx,nlayermx), v_seri(ngridmx,nlayermx)
    49       REAL t_seri(ngridmx,nlayermx)
     37! Output Variables
     38!--------------------------------------------------------
     39
     40      REAL, INTENT(OUT) :: pdu_th(ngridmx,nlayermx)
     41      REAL, INTENT(OUT) :: pdv_th(ngridmx,nlayermx)
     42      REAL, INTENT(OUT) :: pdt_th(ngridmx,nlayermx)
     43      REAL, INTENT(OUT) :: pdq_th(ngridmx,nlayermx,nqmx)
     44      INTEGER, INTENT(OUT) :: lmax(ngridmx)
     45      REAL, INTENT(OUT) :: zmaxth(ngridmx)
     46      REAL, INTENT(OUT) :: pbl_dtke(ngridmx,nlayermx+1)
     47
     48!--------------------------------------------------------
     49! Thermals local variables
     50!--------------------------------------------------------
     51      REAL zu(ngridmx,nlayermx), zv(ngridmx,nlayermx)
     52      REAL zt(ngridmx,nlayermx)
    5053      REAL d_t_ajs(ngridmx,nlayermx)
    5154      REAL d_u_ajs(ngridmx,nlayermx), d_q_ajs(ngridmx,nlayermx,nqmx)
    … …  
    5962      REAL dq_therm(ngridmx,nlayermx), dq_thermdown(ngridmx,nlayermx)
    6063      REAL q2_therm(ngridmx,nlayermx), dq2_therm(ngridmx,nlayermx)
    61 
     64      REAL lmax_real(ngridmx)
     65      REAL masse(ngridmx,nlayermx)
     66      REAL zdz(ngridmx,nlayermx)
    6267      LOGICAL qtransport_thermals,dtke_thermals
    63 
    6468      INTEGER l,ig,iq
    65 
    66 ! Variable de diagnostique : flux de chaleur vertical
     69      CHARACTER (LEN=20) :: modname
     70
     71!--------------------------------------------------------
     72! Local variables for sub-timestep
     73!--------------------------------------------------------
     74
     75      REAL d_t_the(ngridmx,nlayermx), d_q_the(ngridmx,nlayermx,nqmx)
     76      REAL d_u_the(ngridmx,nlayermx),d_v_the(ngridmx,nlayermx)
     77      REAL dq2_the(ngridmx,nlayermx)
     78      INTEGER isplit,nsplit_thermals
     79      REAL r_aspect_thermals
     80      REAL fact
     81      REAL zfm_therm(ngridmx,nlayermx+1),zdt
     82      REAL zentr_therm(ngridmx,nlayermx),zdetr_therm(ngridmx,nlayermx)
     83      REAL zheatFlux(ngridmx,nlayermx)
     84      REAL zheatFlux_down(ngridmx,nlayermx)
     85      REAL zbuoyancyOut(ngridmx,nlayermx)
     86      REAL zbuoyancyEst(ngridmx,nlayermx)
     87      REAL zzw2(ngridmx,nlayermx+1)
     88      REAL zmax(ngridmx)
     89
     90!--------------------------------------------------------
     91! Diagnostics
     92!--------------------------------------------------------
    6793
    6894      REAL heatFlux(ngridmx,nlayermx)
    … …  
    7399
    74100!---------------------------------------------------------
    75 !---------------------------------------------------------
    76 ! **********************************************************************
    77 ! Thermique
    78 ! **********************************************************************
    79 
    80 ! Initialisation des sorties
    81 
    82       lmax_th(:)=1
     101
     102
     103! **********************************************************************
     104! Initialization
     105! **********************************************************************
     106
     107      lmax(:)=0.
    83108      pdu_th(:,:)=0.
    84109      pdv_th(:,:)=0.
    … …  
    98123         pdq_th(:,:,:)=0.
    99124      end if
    100 
    101 ! Dans le model terrestres, les seri sont des q+dq tendances déja cumulées. Il n'y a donc pas de
    102 ! cumulage à l'intérieur de la routine comme dans le model martien. On le fait ici :
    103 
    104             u_seri(:,:)=pu(:,:)+pdu(:,:)*ptimestep
    105             v_seri(:,:)=pv(:,:)+pdv(:,:)*ptimestep
    106             t_seri(:,:)=pt(:,:)+pdt(:,:)*ptimestep
    107 
    108             pq_therm(:,:,:)=0.
    109             qtransport_thermals=.true. !! default setting
    110             !call getin("qtransport_thermals",qtransport_thermals)
    111             if(qtransport_thermals) then
    112             if(tracer) then
    113             pq_therm(:,:,:)=pq(:,:,:)+pdq(:,:,:)*ptimestep
    114             endif
    115             endif
    116 
    117             d_t_ajs(:,:)=0.
    118             d_u_ajs(:,:)=0.
    119             d_v_ajs(:,:)=0.
    120             d_q_ajs(:,:,:)=0.
    121             heatFlux(:,:)=0.
    122             heatFlux_down(:,:)=0.
    123             buoyancyOut(:,:)=0.
    124             buoyancyEst(:,:)=0.
    125 
    126        dtke_thermals=.false. !! default setting
    127        !call getin("dtke_thermals",dtke_thermals)
    128          if(dtke_thermals) then
    129 
    130          DO l=1,nlayermx
    131               q2_therm(:,l)=0.5*(q2(:,l)+q2(:,l+1))
    132          ENDDO
     125      d_t_ajs(:,:)=0.
     126      d_u_ajs(:,:)=0.
     127      d_v_ajs(:,:)=0.
     128      d_q_ajs(:,:,:)=0.
     129      heatFlux(:,:)=0.
     130      heatFlux_down(:,:)=0.
     131      buoyancyOut(:,:)=0.
     132      buoyancyEst(:,:)=0.
     133      zmaxth(:)=0.
     134      lmax_real(:)=0.
     135
     136
     137! **********************************************************************
     138! Preparing inputs for the thermals
     139! **********************************************************************
     140
     141       zu(:,:)=pu(:,:)+pdu(:,:)*ptimestep
     142       zv(:,:)=pv(:,:)+pdv(:,:)*ptimestep
     143       zt(:,:)=pt(:,:)+pdt(:,:)*ptimestep
     144
     145       pq_therm(:,:,:)=0.
     146       qtransport_thermals=.true. !! default setting
     147       !call getin("qtransport_thermals",qtransport_thermals)
     148
     149       if(qtransport_thermals) then
     150          if(tracer) then
     151                pq_therm(:,:,:)=pq(:,:,:)+pdq(:,:,:)*ptimestep
     152          endif
     153       endif
     154
     155!       dtke_thermals=.false. !! default setting
     156!       !call getin("dtke_thermals",dtke_thermals)
     157!
     158!       IF(dtke_thermals) THEN
     159!          DO l=1,nlayermx
     160!              q2_therm(:,l)=0.5*(q2(:,l)+q2(:,l+1))
     161!          ENDDO
     162!       ENDIF
     163
     164! **********************************************************************
     165! **********************************************************************
     166! **********************************************************************
     167! CALLTHERM
     168! **********************************************************************
     169! **********************************************************************
     170! **********************************************************************
     171
     172!         r_aspect_thermals     ! ultimately conrols the amount of mass going through the thermals
     173                                ! decreasing it increases the thermals effect. Tests at gcm resolution
     174                                ! shows that too low values destabilize the model
     175                                ! when changing this value, one should check that the surface layer model
     176                                ! outputs the correct Cd*u and Ch*u through changing the gustiness coefficient beta
     177
     178#ifdef MESOSCALE
     179         !! valid for timesteps < 200s
     180         nsplit_thermals=2
     181         r_aspect_thermals=0.7
     182#else
     183         nsplit_thermals=35
     184         r_aspect_thermals=1.5
     185#endif
     186         call getin("nsplit_thermals",nsplit_thermals)
     187         call getin("r_aspect_thermals",r_aspect_thermals)
     188
     189! **********************************************************************
     190! SUB-TIMESTEP LOOP
     191! **********************************************************************
     192
     193         zdt=ptimestep/REAL(nsplit_thermals)
     194
     195         DO isplit=1,nsplit_thermals
     196
     197! Initialization of intermediary variables
     198
     199         zfm_therm(:,:)=0.
     200         zentr_therm(:,:)=0.
     201         zdetr_therm(:,:)=0.
     202         zheatFlux(:,:)=0.
     203         zheatFlux_down(:,:)=0.
     204!         zbuoyancyOut(:,:)=0.
     205!         zbuoyancyEst(:,:)=0.
     206         zzw2(:,:)=0.
     207         zmax(:)=0.
     208         lmax(:)=0.
     209         d_t_the(:,:)=0.
     210         d_u_the(:,:)=0.
     211         d_v_the(:,:)=0.
     212         dq2_the(:,:)=0.
     213         if (nqmx .ne. 0) then
     214            d_q_the(:,:,:)=0.
    133215         endif
    134216
    135          CALL calltherm_mars(ptimestep,zzlev,zzlay &
    136      &      ,pplay,pplev,pphi &
    137      &      ,u_seri,v_seri,t_seri,pq_therm, q2_therm &
    138      &      ,d_u_ajs,d_v_ajs,d_t_ajs,d_q_ajs, dq2_therm &
    139      &      ,fm_therm,entr_therm,detr_therm &
    140      &      ,lmax_th,zmax_th &
    141      &      ,zw2,fraca &
    142      &      ,zpopsk,ztla,heatFlux,heatFlux_down &
    143      &      ,buoyancyOut,buoyancyEst,hfmax,wmax)
    144 
    145 ! Accumulation des  tendances. On n'accumule pas les quantités de traceurs car celle ci n'a pas du changer
    146 ! étant donné qu'on ne prends en compte que q_seri de la vap d'eau = 0
    147 
    148 ! INCREMENTATION : les d_u_ sont des tendances alors que les pdu sont des dérivees, attention !
    149 
    150            pdu_th(:,:)=d_u_ajs(:,:)/ptimestep
    151            pdv_th(:,:)=d_v_ajs(:,:)/ptimestep
     217             CALL thermcell_main_mars(zdt  &
     218     &      ,pplay,pplev,pphi,zzlev,zzlay  &
     219     &      ,zu,zv,zt,pq_therm,q2_therm  &
     220     &      ,d_u_the,d_v_the,d_t_the,d_q_the,dq2_the  &
     221     &      ,zfm_therm,zentr_therm,zdetr_therm,lmax,zmax  &
     222     &      ,r_aspect_thermals &
     223     &      ,zzw2,fraca,zpopsk &
     224     &      ,ztla,zheatFlux,zheatFlux_down &
     225     &      ,zbuoyancyOut,zbuoyancyEst)
     226
     227      fact=1./REAL(nsplit_thermals)
     228
     229            d_t_the(:,:)=d_t_the(:,:)*ptimestep*fact
     230!            d_u_the(:,:)=d_u_the(:,:)*fact
     231!            d_v_the(:,:)=d_v_the(:,:)*fact
     232!            dq2_the(:,:)=dq2_the(:,:)*fact
     233!            if (nqmx .ne. 0) then
     234!               d_q_the(:,:,:)=d_q_the(:,:,:)*fact
     235!            endif
     236
     237             zmaxth(:)=zmaxth(:)+zmax(:)*fact
     238             lmax_real(:)=lmax_real(:)+float(lmax(:))*fact
     239            fm_therm(:,:)=fm_therm(:,:)  &
     240     &      +zfm_therm(:,:)*fact
     241            entr_therm(:,:)=entr_therm(:,:)  &
     242     &       +zentr_therm(:,:)*fact
     243            detr_therm(:,:)=detr_therm(:,:)  &
     244     &       +zdetr_therm(:,:)*fact
     245
     246            heatFlux(:,:)=heatFlux(:,:) &
     247     &       +zheatFlux(:,:)*fact
     248            heatFlux_down(:,:)=heatFlux_down(:,:) &
     249     &       +zheatFlux_down(:,:)*fact
     250!            buoyancyOut(:,:)=buoyancyOut(:,:) &
     251!     &       +zbuoyancyOut(:,:)*fact
     252!            buoyancyEst(:,:)=buoyancyEst(:,:) &
     253!     &       +zbuoyancyEst(:,:)*fact
     254
     255            zw2(:,:)=zw2(:,:) + zzw2(:,:)*fact
     256
     257!  accumulation de la tendance
     258
     259            d_t_ajs(:,:)=d_t_ajs(:,:)+d_t_the(:,:)
     260!           d_u_ajs(:,:)=d_u_ajs(:,:)+d_u_the(:,:)
     261!           d_v_ajs(:,:)=d_v_ajs(:,:)+d_v_the(:,:)
     262!            d_q_ajs(:,:,:)=d_q_ajs(:,:,:)+d_q_the(:,:,:)
     263!            dq2_therm(:,:)=dq2_therm(:,:)+dq2_the(:,:)
     264!  incrementation des variables meteo
     265
     266            zt(:,:) = zt(:,:) + d_t_the(:,:)
     267!            zu(:,:) = zu(:,:) + d_u_the(:,:)
     268!            zv(:,:) = zv(:,:) + d_v_the(:,:)
     269!            pq_therm(:,:,:) = pq_therm(:,:,:) + d_q_the(:,:,:)
     270!            q2_therm(:,:) = q2_therm(:,:) + dq2_therm(:,:)
     271
     272
     273         ENDDO ! isplit
     274!****************************************************************
     275
     276! Now that we have computed total entrainment and detrainment, we can
     277! advect u, v, and q in thermals. (theta already advected). We can do
     278! that separatly because u,v,and q are not used in thermcell_main for
     279! any thermals-related computation : they are purely passive.
     280
     281! mass of cells
     282      do l=1,nlayermx
     283         masse(:,l)=(pplev(:,l)-pplev(:,l+1))/g
     284      enddo
     285
     286! thickness of layers
     287      do l=1,nlayermx
     288         zdz(:,l)=zzlev(:,l+1)-zzlev(:,l)
     289      enddo
     290
     291      modname='momentum'
     292      call thermcell_dqup(ngridmx,nlayermx,ptimestep                &
     293     &      ,fm_therm,entr_therm,detr_therm,  &
     294     &     masse,zu,d_u_ajs,modname,zdz)
     295
     296      call thermcell_dqup(ngridmx,nlayermx,ptimestep    &
     297     &       ,fm_therm,entr_therm,detr_therm,  &
     298     &     masse,zv,d_v_ajs,modname,zdz)
     299
     300      if (nqmx .ne. 0.) then
     301      modname='tracer'
     302      DO iq=1,nqmx
     303      call thermcell_dqup(ngridmx,nlayermx,ptimestep     &
     304     &     ,fm_therm,entr_therm,detr_therm,  &
     305     &    masse,pq_therm(:,:,iq),d_q_ajs(:,:,iq),modname,zdz)
     306
     307      ENDDO
     308      endif
     309
     310      DO ig=1,ngridmx
     311         hfmax(ig)=MAXVAL(heatFlux(ig,:)+heatFlux_down(ig,:))
     312         wmax(ig)=MAXVAL(zw2(ig,:))
     313      ENDDO
     314
     315      lmax(:)=nint(lmax_real(:))
     316
     317! **********************************************************************
     318! **********************************************************************
     319! **********************************************************************
     320! CALLTHERM END
     321! **********************************************************************
     322! **********************************************************************
     323! **********************************************************************
     324
     325
     326! **********************************************************************
     327! Preparing outputs
     328! **********************************************************************
     329
     330! Winds and tracers PDU, PDV, and PDQ are DERIVATIVES
     331
     332           pdu_th(:,:)=d_u_ajs(:,:)
     333           pdv_th(:,:)=d_v_ajs(:,:)
     334
     335           if(qtransport_thermals) then
     336              if(tracer) then
     337                  pdq_th(:,:,:)=d_q_ajs(:,:,:)
     338              endif
     339           endif
     340
     341!           IF(dtke_thermals) THEN
     342!              DO l=2,nlayermx
     343!                 pbl_dtke(:,l)=0.5*(dq2_therm(:,l-1)+dq2_therm(:,l))
     344!              ENDDO
     345! 
     346!              pbl_dtke(:,1)=0.5*dq2_therm(:,1)
     347!              pbl_dtke(:,nlayermx+1)=0.
     348!           ENDIF
     349
     350
     351! Temperature PDT is a TENDANCY
    152352           pdt_th(:,:)=d_t_ajs(:,:)/ptimestep
    153            if(qtransport_thermals) then
    154            if(tracer) then
    155            pdq_th(:,:,:)=d_q_ajs(:,:,:)/ptimestep
    156            endif
    157            endif
    158 
    159 
    160          DO l=2,nlayermx
    161               pbl_dtke(:,l)=0.5*(dq2_therm(:,l-1)+dq2_therm(:,l))/ptimestep
    162          ENDDO
    163 
    164          pbl_dtke(:,1)=0.5*dq2_therm(:,1)/ptimestep
    165          pbl_dtke(:,nlayermx+1)=0.
    166 !! DIAGNOSTICS
     353
     354! **********************************************************************
     355! Diagnostics
     356! **********************************************************************
    167357       
    168358        if(outptherm) then
    … …  
    189379     &         'tendance temp TH','K',1,d_t_ajs)
    190380        call WRITEDIAGFI(ngridmx,'zmax',  &
    191      &         'pbl height','m',0,zmax_th)
     381     &         'pbl height','m',0,zmaxth)
    192382      else
    193383

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

    r337 r342  
    2929! ============================ OUTPUTS ===========================
    3030
    31       REAL, INTENT(OUT) :: dq_therm(ngridmx,nlayermx)
     31      REAL, INTENT(OUT) :: dq_therm(ngridmx,nlayermx)  ! dq/dt -> derivative
    3232
    3333! ============================ LOCAL =============================
    … …  
    9494        do ig=1, ngridmx
    9595           do k=1,nlayermx-1
    96               dq_therm(ig,k)=-(ptimestep/masse0(ig,k))*(  &
     96              dq_therm(ig,k)=-(1./masse0(ig,k))*(  &
    9797     &           fm0(ig,k+1)*(qa(ig,k+1)-q(ig,k+1)) -   &
    9898     &           fm0(ig,k)*(qa(ig,k)-q(ig,k))          ) &
Note: See TracChangeset for help on using the changeset viewer.