Changeset 566


Ignore:
Timestamp:
Mar 7, 2012, 3:48:18 PM (13 years ago)
Author:
acolaitis
Message:

Something wierd happened in the last commit... So back to where it was and retry... Update to pbl_parameter so that several values of z_out can be sought in one pass. Pbl_parameter can also now be called when calltherm=false, as long as callrichsl=true.

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

Legend:

Unmodified
Added
Removed

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

    r565 r566  
    249249     $              "1-> new correction",
    250250     $              "(matters only if callnirco2=T)"
    251          nircorr=0      !default value
     251         nircorr=1      !default value
    252252         call getin("nircorr",nircorr)
    253253         write(*,*) " nircorr = ",nircorr
    … …  
    290290         endif
    291291
    292 #ifndef MESOSCALE
    293292         if (calladj .and. callrichsl .and. (.not. calltherm)) then
    294293          print*,'You should not be calling the convective adjustment
    … …  
    301300          stop
    302301         endif
    303 #endif
     302
    304303         write(*,*) "call CO2 condensation ?"
    305304         callcond=.true. ! default value

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

    r529 r566  
    11      SUBROUTINE pbl_parameters(ngrid,nlay,ps,pplay,pz0,
    2      & pg,zzlay,zzlev,pu,pv,wstar_in,hfmax,zmax,pts,ph,z_out,
     2     & pg,zzlay,zzlev,pu,pv,wstar_in,hfmax,zmax,pts,ph,z_out,n_out,
    33     & Teta_out,u_out,ustar,tstar,L_mo,vhf,vvv)
    44      IMPLICIT NONE
    … …  
    3030!     pts(ngrid)       surface temperature
    3131!     ph(ngrid,nlay)   potential temperature T*(p/ps)^kappa
    32 !     z_out            height of interpolation
     32!     z_out(n_out)     heights of interpolation
     33!     n_out            number of points for interpolation
    3334!
    3435!   outputs:
    3536!   ------
    3637!
    37 !     Teta_out(ngrid)  interpolated teta
    38 !     u_out(ngrid)     interpolated u
     38!     Teta_out(ngrid,n_out)  interpolated teta
     39!     u_out(ngrid,n_out)     interpolated u
    3940!     ustar(ngrid)     friction velocity
    4041!     tstar(ngrid)     friction temperature
    … …  
    5556!   ----------
    5657
    57       INTEGER, INTENT(IN) :: ngrid,nlay
     58      INTEGER, INTENT(IN) :: ngrid,nlay,n_out
    5859      REAL, INTENT(IN) :: pz0(ngrid),ps(ngrid),pplay(ngrid,nlay)
    5960      REAL, INTENT(IN) :: pg,zzlay(ngrid,nlay),zzlev(ngrid,nlay)
    … …  
    6162      REAL, INTENT(IN) :: wstar_in(ngrid),hfmax(ngrid),zmax(ngrid)
    6263      REAL, INTENT(IN) :: pts(ngrid),ph(ngrid,nlay)
    63       REAL, INTENT(IN) :: z_out
     64      REAL, INTENT(IN) :: z_out(n_out)
    6465
    6566!    Outputs:
    6667!    --------
    6768
    68       REAL, INTENT(OUT) :: Teta_out(ngrid),u_out(ngrid)
    69       REAL T_out(ngrid)
     69      REAL, INTENT(OUT) :: Teta_out(ngrid,n_out),u_out(ngrid,n_out)
     70      REAL T_out(ngrid,n_out)
    7071      REAL, INTENT(OUT) :: ustar(ngrid), tstar(ngrid)
    7172      REAL wstar(ngrid)
    … …  
    7576!   ------
    7677
    77       INTEGER ig,k
     78      INTEGER ig,k,n
    7879      REAL karman,nu
    7980      DATA karman,nu/.41,0.001/
    … …  
    118119!------------------------------------------------------------------------
    119120
     121      DO n=1,n_out
     122
    120123c Initialisation :
    121124
    … …  
    123126      ustar(:)=0.
    124127      tstar(:)=0.
    125       zout=z_out
     128      zout=z_out(n)
    126129      reynolds(:)=0.
    127130      pz0t = 0.
    … …  
    247250      DO ig=1,ngrid
    248251        IF(zout .lt. pz0tcomp(ig)) THEN
    249            u_out(ig)=0.
    250            Teta_out(ig)=pts(ig)
     252           u_out(ig,n)=0.
     253           Teta_out(ig,n)=pts(ig)
    251254        ELSEIF (L_mo(ig) .gt. 0.) THEN
    252            u_out(ig)=(ustar(ig)/karman)*log(zout/pz0(ig)) +
     255           u_out(ig,n)=(ustar(ig)/karman)*log(zout/pz0(ig)) +
    253256     &            5.*(ustar(ig)/(karman*L_mo(ig)))*(zout-pz0(ig))
    254            Teta_out(ig)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))
     257           Teta_out(ig,n)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))
    255258     &                            *log(zout/pz0tcomp(ig)) +
    256259     &            5.*(tstar(ig)/(prandtl(ig)*karman*L_mo(ig)))
    … …  
    260263          IF(L_mo(ig) .gt. -1000.) THEN
    261264           
    262            u_out(ig)=(ustar(ig)/karman)*(
     265           u_out(ig,n)=(ustar(ig)/karman)*(
    263266     &  2.*atan((1.-16.*zout/L_mo(ig))**0.25)
    264267     & -2.*atan((1.-16.*pz0(ig)/L_mo(ig))**0.25)
    … …  
    270273     &                                  )
    271274
    272            Teta_out(ig)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))*(
     275           Teta_out(ig,n)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))*(
    273276     &  2.*log(1.+sqrt(1.-16.*pz0tcomp(ig)/L_mo(ig)))
    274277     & -2.*log(1.+sqrt(1.-16.*zout/L_mo(ig)))
    … …  
    285288! (we do that to avoid using r*4 precision, otherwise, we get -inf values)         
    286289
    287            u_out(ig)=(ustar(ig)/karman)*(
     290           u_out(ig,n)=(ustar(ig)/karman)*(
    288291     &     (4./L_mo(ig))*(zout-pz0(ig))
    289292     &   + (20./(L_mo(ig))**2)*(zout**2-pz0(ig)**2)
    … …  
    292295     &                                   )
    293296
    294            Teta_out(ig)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))*(
     297           Teta_out(ig,n)=pts(ig)+(tstar(ig)/(prandtl(ig)*karman))*(
    295298     &     (8./L_mo(ig))*(zout-pz0tcomp(ig))
    296299     &   + (48./(L_mo(ig))**2)*(zout**2-pz0tcomp(ig)**2)
    … …  
    301304          ENDIF
    302305        ELSE
    303            u_out(ig)=0.
    304            Teta_out(ig)=pts(ig)
     306           u_out(ig,n)=0.
     307           Teta_out(ig,n)=pts(ig)
    305308        ENDIF
    306309        IF(zout .lt. pz0(ig)) THEN
    307            u_out(ig)=0.
     310           u_out(ig,n)=0.
    308311        ENDIF
    309312      ENDDO
    … …  
    314317
    315318      IF ((.not.calltherm).and.(calladj)) THEN
    316            Teta_out(:)=ph(:,1)
     319           Teta_out(:,n)=ph(:,1)
    317320      ENDIF
    318321
    319               T_out(:) = Teta_out(:)*(exp(
     322              T_out(:,n) = Teta_out(:,n)*(exp(
    320323     &   (zout/zzlay(:,1))*(log(pplay(:,1)/ps))
    321324     &                             )
    322325     &                         )**rcp
    323326
     327      ENDDO   !of n=1,n_out
    324328
    325329!------------------------------------------------------------------------
    … …  
    333337
    334338! Nearest index for the pbl height
     339
     340      IF (calltherm) THEN
    335341
    336342      pbl_height_index(:)=1
    … …  
    361367      ENDDO
    362368
     369! Recompute wstar
    363370! We follow Spiga et. al 2010 (QJRMS)
    364371! ------------
    … …  
    410417      vvv(:) = dvvv(:)*(wstar(:))**2
    411418
    412 !------------------------------------------------------------------------
    413 !------------------------------------------------------------------------
    414 ! OUTPUTS
    415 !------------------------------------------------------------------------
    416 !------------------------------------------------------------------------
    417 
    418       IF (ngrid .eq. 1) THEN
    419          dimout=0
    420       ELSE
    421          dimout=2
    422       ENDIF
    423 
    424          call WRITEDIAGFI(ngrid,'Teta_out',
    425      &              'potential temperature at z_out','K',
    426      &                         dimout,Teta_out)
    427          call WRITEDIAGFI(ngrid,'u_out',
    428      &              'horizontal velocity norm at z_out','m/s',
    429      &                         dimout,u_out)
    430          call WRITEDIAGFI(ngrid,'u_star',
    431      &              'friction velocity','m/s',
    432      &                         dimout,ustar)
    433          call WRITEDIAGFI(ngrid,'teta_star',
    434      &              'friction potential temperature','K',
    435      &                         dimout,tstar)
    436          call WRITEDIAGFI(ngrid,'L',
    437      &              'Monin Obukhov length','m',
    438      &                         dimout,L_mo)
    439 !         call WRITEDIAGFI(ngrid,'w_star',
    440 !     &              'Free convection velocity','m',
    441 !     &                         dimout,wstar)
    442 !         call WRITEDIAGFI(ngrid,'z_i',
    443 !     &              'PBL height','m',
    444 !     &                         dimout,zmax)
    445 !         call WRITEDIAGFI(ngrid,'hf_max',
    446 !     &              'Maximum vertical heat flux','m',
    447 !     &                         dimout,hfmax)
    448          call WRITEDIAGFI(ngrid,'vvv',
    449      &              'Vertical velocity variance at zout','m',
    450      &                         dimout,vvv)
    451          call WRITEDIAGFI(ngrid,'vhf',
    452      &              'Vertical heat flux at zout','m',
    453      &                         dimout,vhf)
     419      ENDIF ! of if calltherm
    454420
    455421      RETURN

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

    r556 r566  
    345345      REAL pdu_th(ngridmx,nlayermx),pdv_th(ngridmx,nlayermx)
    346346      REAL pdt_th(ngridmx,nlayermx),pdq_th(ngridmx,nlayermx,nqmx)
    347       INTEGER lmax_th(ngridmx)
     347      INTEGER lmax_th(ngridmx),dimout,n_out,n
     348      CHARACTER(50) zstring
    348349      REAL dtke_th(ngridmx,nlayermx+1)
    349350      REAL zcdv(ngridmx), zcdh(ngridmx)
    350       REAL Teta_out(ngridmx),u_out(ngridmx)  ! Interpolated teta and u at z_out
    351       REAL z_out                          ! height of interpolation between z0 and z1 [meters]
     351      REAL, ALLOCATABLE, DIMENSION(:,:) :: Teta_out
     352      REAL, ALLOCATABLE, DIMENSION(:,:) :: u_out ! Interpolated teta and u at z_out
     353      REAL, ALLOCATABLE, DIMENSION(:) :: z_out     ! height of interpolation between z0 and z1 [meters]
    352354      REAL ustar(ngridmx),tstar(ngridmx)  ! friction velocity and friction potential temp
    353355      REAL L_mo(ngridmx),wstarpbl(ngridmx),vhf(ngridmx),vvv(ngridmx)
    … …  
    838840
    839841c-----------------------------------------------------------------------
    840 c   TEST. Thermals :
    841 c HIGHLY EXPERIMENTAL, BEWARE !!
     842c   5. Thermals :
    842843c   -----------------------------
    843  
     844
    844845      if(calltherm) then
    845846 
    … …  
    13771378      ENDIF ! of IF (lwrite)
    13781379
     1380c        ----------------------------------------------------------
     1381c        ----------------------------------------------------------
     1382c        INTERPOLATIONS IN THE SURFACE-LAYER
     1383c        ----------------------------------------------------------
     1384c        ----------------------------------------------------------
     1385
     1386         IF (1 .eq. 0.) THEN
     1387         IF (callrichsl) THEN
     1388            n_out=5
     1389
     1390            ALLOCATE(z_out(n_out))
     1391            ALLOCATE(Teta_out(ngrid,n_out))
     1392            ALLOCATE(u_out(ngrid,n_out))
     1393
     1394            z_out(:)=[0.001,0.05,0.1,0.5,1.]
     1395            u_out(:,:)=0.
     1396            Teta_out(:,:)=0.
     1397
     1398            call pbl_parameters(ngrid,nlayer,ps,zplay,z0,
     1399     & g,zzlay,zzlev,zu,zv,wstar,hfmax_th,zmax_th,tsurf,zh,z_out,n_out,
     1400     & Teta_out,u_out,ustar,tstar,wstarpbl,L_mo,vhf,vvv)
     1401
     1402#ifndef MESOSCALE
     1403            IF (ngrid .eq. 1) THEN
     1404               dimout=0
     1405            ELSE
     1406               dimout=2
     1407            ENDIF
     1408            DO n=1,n_out
     1409               write(zstring, '(F9.6)') z_out(n)
     1410               call WRITEDIAGFI(ngrid,'Teta_out_'//trim(zstring),
     1411     &   'potential temperature at z_out','K',dimout,Teta_out(:,n))
     1412               call WRITEDIAGFI(ngrid,'u_out_'//trim(zstring),
     1413     &   'horizontal velocity norm at z_out','m/s',dimout,u_out(:,n))
     1414            ENDDO
     1415            call WRITEDIAGFI(ngrid,'u_star',
     1416     &   'friction velocity','m/s',dimout,ustar)
     1417            call WRITEDIAGFI(ngrid,'teta_star',
     1418     &   'friction potential temperature','K',dimout,tstar)
     1419            call WRITEDIAGFI(ngrid,'L',
     1420     &   'Monin Obukhov length','m',dimout,L_mo)
     1421            call WRITEDIAGFI(ngrid,'vvv',
     1422     &   'Vertical velocity variance at zout','m',dimout,vvv)
     1423            call WRITEDIAGFI(ngrid,'vhf',
     1424     &   'Vertical heat flux at zout','m',dimout,vhf)
     1425#endif
     1426
     1427         ENDIF
     1428         ENDIF   ! of pbl interpolation outputs
     1429
     1430c        ----------------------------------------------------------
     1431c        ----------------------------------------------------------
     1432c        END OF SURFACE LAYER INTERPOLATIONS
     1433c        ----------------------------------------------------------
     1434c        ----------------------------------------------------------
     1435
    13791436      IF (ngrid.NE.1) THEN
    13801437
    … …  
    18631920c        ----------------------------------------------------------
    18641921
    1865 
    1866 c        ----------------------------------------------------------
    1867 c        Outputs of surface layer
    1868 c        ----------------------------------------------------------
    1869 
    1870 
    1871          z_out=0.
    1872          if (calltherm .and. (z_out .gt. 0.)) then
    1873 
    1874          call pbl_parameters(ngrid,nlayer,ps,zplay,z0,
    1875      & g,zzlay,zzlev,zu,zv,wstar,hfmax_th,zmax_th,tsurf,zh,z_out,
    1876      & Teta_out,u_out,ustar,tstar,wstarpbl,L_mo,vhf,vvv)
    1877 
    1878          else
    1879            if((.not. calltherm).and.(z_out .gt. 0.)) then
    1880             print*, 'WARNING : no interpolation in surface-layer :'
    1881             print*, 'Outputing surface-layer quantities without thermals
    1882      & does not make sense'
    1883            endif
    1884          endif
    1885 
    18861922c        ----------------------------------------------------------
    18871923c        Outputs of thermals
    … …  
    19611997
    19621998! THERMALS STUFF 1D
    1963 
    1964          z_out=1.
    1965          if (calltherm .and. (z_out .gt. 0.)) then
    1966 
    1967          call pbl_parameters(ngrid,nlayer,ps,zplay,z0,
    1968      & g,zzlay,zzlev,zu,zv,wstar,hfmax_th,zmax_th,tsurf,zh,z_out,
    1969      & Teta_out,u_out,ustar,tstar,wstarpbl,L_mo,vhf,vvv)
    1970 
    1971          else
    1972            if((.not. calltherm).and.(z_out .gt. 0.)) then
    1973             print*, 'WARNING : no interpolation in surface-layer :'
    1974             print*, 'Outputing surface-layer quantities without thermals
    1975      & does not make sense'
    1976            endif
    1977          endif
    1978 
    19791999         if(calltherm) then
    19802000
Note: See TracChangeset for help on using the changeset viewer.