!############################# Change Log ################################## ! 5.0.0 ! !########################################################################### ! Copyright (C) 1990, 1995, 1999, 2000, 2003 - All Rights Reserved ! Regional Atmospheric Modeling System - RAMS !########################################################################### module ModTimestep contains subroutine timestep(OneGrid) use ModMessageSet, only: & PostRecvSendMsgs, & WaitRecvMsgs use ModAcoust, only: acoustic_new use ModGrid, only: & Grid use mem_basic, only: & basic_g ! INTENT(INOUT) use node_mod, only: & mzp, mxp, myp, & ! INTENT(IN) ia, iz, ja, jz, & ! INTENT(IN) i0, j0, & ! INTENT(IN) izu, jzv, & ! INTENT(IN) mynum, & ! INTENT(IN) ibcon, & ! INTENT(IN) nmachs ! INTENT(IN) use mem_cuparm, only: & NNQPARM, & ! INTENT(IN) IF_CUINV ! INTENT(IN) use mem_varinit, only: & NUD_TYPE ! INTENT(IN) use mem_turb, only: & IF_URBAN_CANOPY, & ! INTENT(IN) ihorgrad ! INTENT(IN) use mem_oda, only: & if_oda ! INTENT(IN) use micphys, only: & mcphys_type, &! INTENT(IN) level ! INTENT(IN) use mem_grid, only: & ngrids, & ! INTENT(IN) ngrid, & ! INTENT(IN) npatch, & ! INTENT(IN) time, & ! INTENT(IN) dtlong, & ! INTENT(IN) dtlongn, & ! INTENT(IN) iyear1, & ! INTENT(IN) imonth1, & ! INTENT(IN) idate1, & ! INTENT(IN) grid_g, & ! INTENT(INOUT) nxtnest, & ! INTENT(IN) if_adap, & ! INTENT(IN) dtlt, & ! INTENT(IN) istp, & ! INTENT(IN) jdim, & ! INTENT(IN) nzp, & ! INTENT(IN) f_thermo_e, & ! INTENT(IN) f_thermo_w, & ! INTENT(IN) f_thermo_s, & ! INTENT(IN) f_thermo_n, & ! INTENT(IN) zt, & zm, & dzt, & nzpmax, & itime1, & vveldamp ! INTENT(IN) use shcu_vars_const, only: & ! For Shallow Cumulus Paramet. NNSHCU ! INTENT(IN) use mem_scalar, only: & ! For SiB scalar_g ! INTENT(IN) use mem_leaf, only: & ! For SiB ISFCL ! INTENT(IN) ! TEB_SPM use teb_spm_start, only: & TEB_SPM ! INTENT(IN) use mem_emiss, only: & ichemi, & ! INTENT(IN) isource ! INTENT(IN) !ALF- Necessary in new advection scheme use advect_kit, only : & calc_advec ! Subroutine ! For specific optimization depending the type of machine use machine_arq, only: & machine ! INTENT(IN) use rconstants, only: & g, & ! (IN) cp, & ! (IN) cpor, & ! (IN) p00, & ! (IN) rgas, & ! (IN) pi180 ! (IN) use ccatt_start, only: & ccatt ! (IN) use mem_chem1, only: & nvert_src=>chem1_src_z_dim_g, & ! (IN) chem1_g, & ! (INOUT) nsrc, & ! (IN) chem1_src_g, & ! %sc_src(INOUT) chemistry, & ! (IN) split_method, & ! (IN) n_dyn_chem, & ntimes_src use mem_aer1, only: & aerosol, &! (IN) aer1_g, &! %sc_src(INOUT) aer2_g, &! %sc_src(INOUT) aer_nvert_src=>aer1_src_z_dim_g ! (IN) use mem_plume_chem1, only: plume_mean_g ! %flam_frac(IN), %fire_size(IN) use mem_stilt, only: & iexev, & ! (IN) imassflx, & ! (IN) stilt_g ! %dnp (IN) use mem_radiate, only: radiate_g use chem_sources, only : & alloc_emiss_cycle, & ! Subroutine init_actual_time_index, & ! Subroutine emiss_cycle, & ! (INOUT) emiss_cycle_alloc, & srcmapfn ! (IN) use ChemSourcesDriver, only: sources_driver ! Subroutine use ChemDryDepDriver , only: drydep_driver ! Subroutine use module_chemistry_driver, only: chemistry_driver ! Subroutine use radiation, only: radiate ! Subroutine use ModTimeStamp, only: SynchronizedTimeStamp, TimeStamp use cuparm_grell3, only: cuparm_grell3_catt ! subroutine use digitalFilter, only: & applyDigitalFilter, & ! subroutine fileNameDF, & ! intent(inout) - file control dfVars, & applyDF USE monotonic_adv, only: & advmnt_driver, & ! subroutine advmnt USE DriverMatrix, ONLY: MatrixDriver !Matrix Aerosol Model USE rams_microphysics_2M, only: micro_2M_rams60,negadj1_2M_rams60 implicit none type(Grid), pointer :: OneGrid ! execution time instrumentation include "tsNames.h" INTEGER, PARAMETER :: acoshdp = 0 ! +-------------------------------------------------------------+ ! | Timestep driver for the hybrid non-hydrostatic time-split | ! | model. | ! +-------------------------------------------------------------+ ! Zero out all tendency arrays. !-------------------------------- call TEND0() ! Thermodynamic diagnosis !-------------------------------- if (mcphys_type <= 1 .and. level/=3) then call THERMO(mzp,mxp,myp,ia,iz,ja,jz,'SUPSAT') endif ! Radiation parameterization !-------------------------------- !call RADIATE(mzp,mxp,myp,ia,iz,ja,jz,mynum) !---------------------------------------- if (CCATT==1 .and. chemistry >= 10000000) then!tmp !- emission !- allocation for diurnal cycle of emission arrays and !- actual_time_index on nodes if( (.not. emiss_cycle_alloc) .and. (chemistry >= 0) .and. & trim(srcmapfn) .ne. 'NONE' .and. trim(srcmapfn) .ne. 'none') then call alloc_emiss_cycle(mxp,myp,ngrids,nsrc) call init_actual_time_index(nsrc,ntimes_src) endif ! !-LFR Sea salt Aerossol inline source call SeaSaltDriver(ia,iz,ja,jz,ngrid,mxp,myp) !plume_mean_g(:,:) instead of plume_mean_g(:,ngrid) to avoid memory errors. !emiss_cycle(:,:) instead of emiss_cycle(:,ngrid) to avoid memory errors. !the same for the others var call sources_driver(ngrid, mzp,mxp,myp,ia,iz,ja,jz, & g,cp,cpor,p00,rgas,pi180, & radiate_g(ngrid)%cosz,basic_g(ngrid)%theta, & basic_g(ngrid)%pp,basic_g(ngrid)%pi0,basic_g(ngrid)%rv, & basic_g(ngrid)%dn0,basic_g(ngrid)%up,basic_g(ngrid)%vp, & time,iyear1,imonth1,idate1,itime1,dtlt, & grid_g(ngrid)%rtgt,grid_g(ngrid)%lpw,grid_g(ngrid)%glat, & grid_g(ngrid)%glon,zt,zm,dzt,nzpmax, & nvert_src (:,:), & chem1_g (:,:), & chem1_src_g (:,:,:,:), & aer1_g (:,:,:), & aer_nvert_src(:,:), & plume_mean_g (:,:), & stilt_g(ngrid)%dnp, & emiss_cycle (:,:), & aer2_g (:,:) ) !- call dry deposition and sedimentation routines !call drydep_driver(mzp,mxp,myp,ia,iz,ja,jz) !- call Matrix Aerosol Model (for parallel spliting) !---------------------------------------- !if(AEROSOL==2) then ! CALL MatrixDriver(ia,iz,ja,jz,mzp,mxp,myp) !endif endif !------------ transport ------------------------------------ ! ! Send boundaries to adjoining nodes !-------------------------------------------------- if (nmachs > 1) then call PostRecvSendMsgs(OneGrid%SelectedGhostZoneSend, OneGrid%SelectedGhostZoneRecv) endif ! Analysis nudging and boundary condition !--------------------------------------------------- if (NUD_TYPE>0) call DATASSIM() ! Get the overlap region between parallel nodes !--------------------------------------------------- if (nmachs > 1) then call WaitRecvMsgs(OneGrid%SelectedGhostZoneSend, OneGrid%SelectedGhostZoneRecv) endif ! Sub-grid diffusion terms !---------------------------------------- ! call diffuse() !---------------------------------------- IF(advmnt == 1) THEN CALL AdvMnt_driver('T',mzp,mxp,myp,ia,iz,ja,jz,izu,jzv,mynum) ELSE ! CALL ADVECTc('T',mzp,mxp,myp,ia,iz,ja,jz,izu,jzv,mynum) !----- for R-K time integration --- ! Scalar advection ! output: scalart (=>scalar_tab%var_t) CALL advectc_rk('T',mzp,mxp,myp,ia,iz,ja,jz,izu,jzv,mynum) !----- for R-Ks time integration --- ENDIF !------------ end of transport ------------------------------- ! Update scalars !---------------------------------------- call PREDTR() ! Apply scalar b.c.'s !---------------------------------------- !call TRSETS() !--------> tmp srf !- Chemistry !---------------------------------------- if (ccatt==1) THEN !- call Matrix Aerosol Model !---------------------------------------- !- using symmetric/sequential spliting operator if(AEROSOL==2) CALL MatrixDriver(ia,iz,ja,jz,mzp,mxp,myp) endif if (mcphys_type == 2 .or. mcphys_type == 3 ) then ! G. Thompson microphysics !call micro_thompson() !--------> tmp srf endif end subroutine timestep end module ModTimestep !************************************************************************* subroutine mass_flux(n1,n2,n3,m1,m2,m3,up,vp,wp & ,dn0,rtgu,rtgv,dyu,dxv,pp,pi0) use mem_grid use rconstants implicit none integer :: n1,n2,n3,m1,m2,m3 real :: up(m1,m2,m3),vp(m1,m2,m3),wp(m1,m2,m3) & ,dn0(n1,n2,n3),rtgu(n2,n3),dyu(n2,n3),dxv(n2,n3) & ,rtgv(n2,n3),pp(m1,m2,m3),pi0(n1,n2,n3) real, save :: aintmass=0. integer :: i,j,k real :: wmass,emass,smass,nmass,prtot,tmass,ppp,area !cc if (mod(time,300.).gt..1) return ! west/east bound wmass=0. emass=0. do j=2,nyp-1 do k=2,nzp-1 i=1 wmass=wmass + & up(k,i,j)*rtgu(i,j)/(dyu(i,j)*dzt(k)) & *(dn0(k,i,j)+dn0(k,i+1,j))*.5 i=nxp-1 emass=emass - & up(k,i,j)*rtgu(i,j)/(dyu(i,j)*dzt(k)) & *(dn0(k,i,j)+dn0(k,i+1,j))*.5 enddo enddo ! north/south bound smass=0. nmass=0. do i=2,nxp-1 do k=2,nzp-1 j=1 smass=smass + & vp(k,i,j)*rtgv(i,j)/(dxv(i,j)*dzt(k)) & *(dn0(k,i,j)+dn0(k,i,j+1))*.5 j=nyp-1 nmass=nmass - & vp(k,i,j)*rtgv(i,j)/(dxv(i,j)*dzt(k)) & *(dn0(k,i,j)+dn0(k,i,j+1))*.5 enddo enddo k=2 prtot=0. do j=2,nyp-1 do i=2,nxp-1 ppp= ( (pp(k,i,j)+pi0(k,i,j))/cp )**cpor*p00 prtot=prtot+ppp/(dyu(i,j)*dxv(i,j)) enddo enddo tmass=wmass+emass+smass+nmass aintmass=aintmass+tmass*dtlong area=(nxp-2)*deltax*(nyp-2)*deltay print*,'===============================' print*,' Mass flux - W, E, S, N' print*, wmass,emass,smass,nmass print*, 'total (kg/(m2 s):',tmass/area print*, 'total (kg/m2):',aintmass/area print*, 'total pr change (pa):',aintmass/area*9.8 print*, 'computed mean press:',prtot/area print*,'===============================' return end subroutine mass_flux ! ***************************************************************** subroutine w_damping(mzp,mxp,myp,ia,iz,ja,jz,mynum) use mem_basic, only: & basic_g ! intent(in) use mem_grid, only: & ngrid, & ! intent(in) grid_g ! intent(in) use mem_tend, only: tend implicit none integer, intent(in) :: mzp integer, intent(in) :: mxp integer, intent(in) :: myp integer, intent(in) :: ia integer, intent(in) :: iz integer, intent(in) :: ja integer, intent(in) :: jz integer, intent(in) :: mynum call apply_wdamp(mzp,mxp,myp,ia,iz,ja,jz,mynum & ,basic_g(ngrid)%up ,basic_g(ngrid)%vp & ,basic_g(ngrid)%wp ,grid_g(ngrid)%rtgt & ,grid_g(ngrid)%f13t ,grid_g(ngrid)%f23t & ,grid_g(ngrid)%dxt ,grid_g(ngrid)%dyt & ,tend%ut ,tend%vt & ,tend%wt ) end subroutine w_damping ! ********************************************************************** subroutine apply_wdamp(m1,m2,m3,ia,iz,ja,jz,mynum,up,vp,wp,rtgt,f13t,f23t,dxt,dyt,ut,vt,wt) use mem_grid, only: & jdim, & ! intent(in) ngrid, & ! intent(in) dtlt, & ! intent(in) ht, & ! intent(in) dzt ! intent(in) implicit none integer, intent(in) :: m1 integer, intent(in) :: m2 integer, intent(in) :: m3 integer, intent(in) :: ia integer, intent(in) :: iz integer, intent(in) :: ja integer, intent(in) :: jz integer, intent(in) :: mynum real, intent(in) :: up(m1,m2,m3) real, intent(in) :: vp(m1,m2,m3) real, intent(in) :: wp(m1,m2,m3) real, intent(in) :: rtgt(m2,m3) real, intent(in) :: f13t(m2,m3) real, intent(in) :: f23t(m2,m3) real, intent(in) :: dxt(m2,m3) real, intent(in) :: dyt(m2,m3) real, intent(inout) :: ut(m1,m2,m3) real, intent(inout) :: vt(m1,m2,m3) real, intent(inout) :: wt(m1,m2,m3) integer :: i,j,k,ifm,icm,innest real :: c1x,c1y,c1z,cflnumh,cflnumv,cflz real :: vctr1(m1) real :: vctr2(m1) real :: vctr3(m1) real , parameter ::gama_w=0.3 !m/sē ! This routine damps the vertical velocity when CFLZ is exceed ! (Actually check on 80% of CFL) cflnumh = .80 cflnumv = .50 cflz=0.0 ! !c1x=0.0 !vctr3(:)=0.0 do j = ja,jz do i = ia,iz do k = 2,m1-1 !cflx !vctr1(k) = .5*(up(k,i,j)+up(k,i-1,j))*dtlt*dxt(i,j) !cfly !vctr2(k) = .5*(vp(k,i,j)+vp(k,i,j-jdim))*dtlt*dyt(i,j) !cflz vctr3(k) = ((wp(k,i,j)+wp(k-1,i,j)) & +(up(k,i,j)+up(k,i-1,j))*f13t(i,j)*ht(k)*rtgt(i,j) & +(vp(k,i,j)+vp(k,i,j-jdim))*f23t(i,j)*ht(k)*rtgt(i,j) & )*.5*dtlt*dzt(k) c1z=abs(vctr3(k)) if(c1z > cflnumv) then wt(k,i,j) = wt(k,i,j) -gama_w*sign(1.,wp(k,i,j))*(c1z-cflnumv) print*,'wdamp applied at=',k,i,j,mynum,c1z,wp(k,i,j) call flush(6) endif enddo !do k = 2,m1-1 ! c1x = abs(vctr1(k)) ! c1y = abs(vctr2(k)) ! c1z = abs(vctr3(k)) ! ! if (c1x .gt. cflxy(ngrid)) cflxy(ngrid) = c1x ! if (c1y .gt. cflxy(ngrid)) cflxy(ngrid) = c1y ! if (c1z .gt. cflz) cflz = c1z !enddo enddo enddo !print*,'at wdamp2-max cflz',mynum,cflz !call flush(6) end subroutine apply_wdamp ! ***************************************************************************