subroutine statsdump use modfields, only : um,up,vm,wm,svm,qtm,thlm,pres0,ncstaty,ncstatxy,ncstatyt,ncstattke,ncstatmint,& ncstatkslice,ncstatislice,ncstatjslice,t_t,t_v,t_p,t_sgs,d_sgs,p_b,p_t,adv,& IIc,IIu,IIv,IIw,IIuw,IIvw,IIct,IIwt,IIut,IIvt,IIuwt,IIuv,& IIcs,IIws,IIus,IIvs,IIuws,IIvws,IIuvs,& vyt,uyt,wyt,thlyt,qtyt,& sca1yt,sca2yt,sca3yt,thlsgsyt,qtsgsyt,sv1sgsyt,sv2sgsyt,sv3sgsyt,usgsyt,wsgsyt,& usgsxyt,thlsgsxyt,vsgsxyt,uwtik,vwtjk,uvtij,utik,wtik,wtjk,vtjk,utij,vtij,& wthltk,wqttk,thlthlt,qtqtt,sv1sv1t,sv2sv2t,sv3sv3t,sv4sv4t,wmt,thltk,qttk,thlt,uxyt,vxyt,wxyt,thlxyt,& ncstatxyt,qtxyt,pxyt,ncstatt,uutc,vvtc,wwtc,utc,vtc,wtc,& umt,vmt,sv1t,sv2t,sv3t,sv4t,sv1tk,sv2tk,sv3tk,sv4tk,wsv1tk,wsv2tk,wsv3tk,wsv4tk,& sv1sgst,sv2sgst,sv3sgst,sv4sgst,qtt,pt,PSSt,& !,sv1max,sv2max,sv3max,sv4max ncstattr,tr_u,tr_ut,tr_v,tr_vt,tr_w,tr_wt,tr_thl,tr_thlt,tr_qt,tr_qtR,& tr_qtA,tr_qtt,tr_qtRt,tr_qtAt,tr_sv,tr_sv1t, PSSt, tr_sv2t,tr_omega,tr_omegat use modglobal, only : ib,ie,ih,ihc,xf,xh,jb,je,jhc,jgb,jge,dy,dyi,jh,ke,kb,kh,khc,rk3step,& timee,cexpnr,tsample,tstatsdump,jtot,imax,jmax,dzf,& ltempeq,zh,dxf,dzf,dzh2i,lprofforc,lscasrcl,& lkslicedump,lislicedump,ljslicedump,lchem,dzhi,dzfi,dzhiq,dxhi,lmoist,nsv,& k1,JNO2,lchem,kslice,islice,jslice,isliceloc,jsliceloc,islicerank,jslicerank,& ltreedump ! use modsubgriddata, only : ekm,sbshr use modstat_nc, only : writestat_nc,writestat_1D_nc use modmpi, only : myid,cmyid,my_real,mpi_sum,avey_ibm,mpierr,& comm3d,avexy_ibm,nprocs use modsurfdata, only : thls use modsubgrid, only : ekh,ekm use modstatistics, only : genstats,tkestats implicit none !> Create fields to be used in statistics ! interpolated fields ! real, dimension(ib:ie,jb:je,kb:ke) :: umc ! real, dimension(ib:ie,jb:je,kb:ke) :: vmc ! real, dimension(ib:ie,jb:je,kb:ke) :: wmc ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: thlk ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: qtk ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: uik ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wik ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: vjk ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wjk ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: uij ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: vij ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: uc ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: vc ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wc real, allocatable :: thlk(:,:,:) real, allocatable :: qtk(:,:,:) real, allocatable :: uik(:,:,:) real, allocatable :: wik(:,:,:) real, allocatable :: vjk(:,:,:) real, allocatable :: wjk(:,:,:) real, allocatable :: uij(:,:,:) real, allocatable :: vij(:,:,:) real, allocatable :: uc(:,:,:) real, allocatable :: vc(:,:,:) real, allocatable :: wc(:,:,:) ! SGS fluxes ! real, dimension(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh) :: thlsgs ! real, dimension(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh) :: qtsgs ! real, dimension(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh) :: usgs ! real, dimension(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh) :: vsgs ! real, dimension(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh) :: wsgs real, allocatable :: thlsgs(:,:,:) real, allocatable :: qtsgs(:,:,:) real, allocatable :: usgs(:,:,:) real, allocatable :: vsgs(:,:,:) real, allocatable :: wsgs(:,:,:) ! t-averaged fields ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv1k ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv2k ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv3k ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv4k ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wpsv1p ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wpsv2p ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wpsv3p ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wpsv4p ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv1sgs ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv2sgs ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv3sgs ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: sv4sgs real, allocatable :: sv1k(:,:,:) real, allocatable :: sv2k(:,:,:) real, allocatable :: sv3k(:,:,:) real, allocatable :: sv4k(:,:,:) real, allocatable :: sv1sgs(:,:,:) real, allocatable :: sv2sgs(:,:,:) real, allocatable :: sv3sgs(:,:,:) real, allocatable :: sv4sgs(:,:,:) real, allocatable :: wpsv1p(:,:,:) real, allocatable :: wpsv2p(:,:,:) real, allocatable :: wpsv3p(:,:,:) real, allocatable :: wpsv4p(:,:,:) real, allocatable :: sv1psv1pt(:,:,:) real, allocatable :: sv2psv2pt(:,:,:) real, allocatable :: sv3psv3pt(:,:,:) real, allocatable :: sv4psv4pt(:,:,:) real, allocatable :: PSS(:,:,:) ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: upwptik ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: vpwptjk ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: upvptij ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wpthlptk ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: thlpthlpt ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: upuptc ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: vpvptc ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: wpwptc ! real, dimension(ib:ie,jb:je,kb:ke+kh) :: tketc real, allocatable :: upwptik(:,:,:) real, allocatable :: vpwptjk(:,:,:) real, allocatable :: upvptij(:,:,:) real, allocatable :: wpthlptk(:,:,:) real, allocatable :: thlpthlpt(:,:,:) real, allocatable :: upuptc(:,:,:) real, allocatable :: vpvptc(:,:,:) real, allocatable :: wpwptc(:,:,:) real, allocatable :: tketc(:,:,:) ! y-averaged fields real, dimension(ib:ie,kb:ke) :: uy real, dimension(ib:ie,kb:ke) :: vy real, dimension(ib:ie,kb:ke) :: wy real, dimension(ib:ie,kb:ke) :: thly real, dimension(ib:ie,kb:ke) :: qty real, dimension(ib:ie,kb:ke) :: sca1y real, dimension(ib:ie,kb:ke) :: sca2y real, dimension(ib:ie,kb:ke) :: sca3y real, dimension(ib:ie,kb:ke) :: usgsy real, dimension(ib:ie,kb:ke) :: wsgsy real, dimension(ib:ie,kb:ke) :: thlsgsy real, dimension(ib:ie,kb:ke) :: qtsgsy real, dimension(ib:ie,kb:ke) :: sv1sgsy real, dimension(ib:ie,kb:ke) :: sv2sgsy real, dimension(ib:ie,kb:ke) :: sv3sgsy real, dimension(ib:ie,kb:ke) :: uwyik real, dimension(ib:ie,kb:ke) :: wthlyk real, dimension(ib:ie,kb:ke) :: wqtyk ! real, dimension(ib:ie,kb:ke) :: wsv1yk ! real, dimension(ib:ie,kb:ke) :: wsv2yk real, dimension(ib:ie,kb:ke) :: wyik real, dimension(ib:ie,kb:ke) :: uyik real, dimension(ib:ie,kb:ke) :: thlyk real, dimension(ib:ie,kb:ke) :: upwpyik real, dimension(ib:ie,kb:ke) :: wpthlpyk ! ty-averaged fluxes real, dimension(ib:ie,kb:ke) :: upwptyik real, dimension(ib:ie,kb:ke) :: wpthlptyk real, dimension(ib:ie,kb:ke) :: wpqtptyk real, dimension(ib:ie,kb:ke) :: wpsv1ptyk real, dimension(ib:ie,kb:ke) :: wpsv2ptyk real, dimension(ib:ie,kb:ke) :: wpsv3ptyk real, dimension(ib:ie,kb:ke) :: upuptyc real, dimension(ib:ie,kb:ke) :: vpvptyc real, dimension(ib:ie,kb:ke) :: wpwptyc real, dimension(ib:ie,kb:ke) :: qtpqtpty real, dimension(ib:ie,kb:ke) :: thlpthlpty real, dimension(ib:ie,kb:ke) :: sv1psv1pty real, dimension(ib:ie,kb:ke) :: sv2psv2pty real, dimension(ib:ie,kb:ke) :: sv3psv3pty real, dimension(ib:ie,kb:ke) :: uwtyik real, dimension(ib:ie,kb:ke) :: wthltyk real, dimension(ib:ie,kb:ke) :: wqttyk real, dimension(ib:ie,kb:ke) :: wsv1tyk real, dimension(ib:ie,kb:ke) :: wsv2tyk real, dimension(ib:ie,kb:ke) :: wsv3tyk ! xy-averaged fields real, dimension(kb:ke+kh) :: uxy real, dimension(kb:ke+kh) :: vxy real, dimension(kb:ke+kh) :: wxy real, dimension(kb:ke+kh) :: thlxy real, dimension(kb:ke+kh) :: qtxy real, dimension(kb:ke+kh) :: pxy real, dimension(kb:ke+kh) :: usgsxy real, dimension(kb:ke+kh) :: thlsgsxy real, dimension(kb:ke+kh) :: vsgsxy real, dimension(kb:ke+kh) :: sca1xy real, dimension(kb:ke+kh) :: uwxyik real, dimension(kb:ke+kh) :: vwxyjk real, dimension(kb:ke+kh) :: wthlxyk real, dimension(kb:ke+kh) :: thlxyk real, dimension(kb:ke+kh) :: wxyik real, dimension(kb:ke+kh) :: uxyik real, dimension(kb:ke+kh) :: vxyjk real, dimension(kb:ke+kh) :: wxyjk real, dimension(kb:ke+kh) :: upwpxyik real, dimension(kb:ke+kh) :: wpthlpxyk real, dimension(kb:ke+kh) :: vpwpxyjk ! txy-averaged fields real, dimension(kb:ke+kh) :: upwptxyik real, dimension(kb:ke+kh) :: wpthlptxyk real, dimension(kb:ke+kh) :: thlpthlptxy real, dimension(kb:ke+kh) :: upuptxyc real, dimension(kb:ke+kh) :: vpvptxyc real, dimension(kb:ke+kh) :: wpwptxyc real, dimension(kb:ke+kh) :: tketxyc real, dimension(kb:ke+kh) :: vpwptxyjk real, dimension(kb:ke+kh) :: upvptxyij real, dimension(kb:ke+kh) :: uwtxyik real, dimension(kb:ke+kh) :: wthltxyk real, dimension(kb:ke+kh) :: vwtxyjk real, dimension(kb:ke+kh) :: wwtxyk real, dimension(kb:ke+kh) :: uvtxyij real, allocatable :: field(:,:), varsy(:,:,:),varsyt(:,:,:),varstke(:,:),varsxy(:,:),& varkslice(:,:,:),varislice(:,:,:),varjslice(:,:,:),varsxyt(:,:),varst(:,:,:,:),varstr(:,:,:,:),varsmint(:,:,:,:) real :: tstatsdumppi,emom integer :: i,j,k,ip,im,jp,jm,kp,km integer :: writecounter = 1 integer :: reclength if (.not.(lytdump .or. lydump .or. lxydump .or. lxytdump .or. ltdump .or. lmintdump & .or. lkslicedump.or. lislicedump.or. ljslicedump)) return allocate(thlk(ib:ie,jb:je,kb:ke+kh)) allocate(qtk(ib:ie,jb:je,kb:ke+kh)) allocate(uik(ib:ie,jb:je,kb:ke+kh)) allocate(wik(ib:ie,jb:je,kb:ke+kh)) allocate(vjk(ib:ie,jb:je,kb:ke+kh)) allocate(wjk(ib:ie,jb:je,kb:ke+kh)) allocate(uij(ib:ie,jb:je,kb:ke+kh)) allocate(vij(ib:ie,jb:je,kb:ke+kh)) allocate(uc(ib:ie,jb:je,kb:ke+kh)) allocate(vc(ib:ie,jb:je,kb:ke+kh)) allocate(wc(ib:ie,jb:je,kb:ke+kh)) allocate(thlsgs(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh)) allocate(qtsgs(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh)) allocate(usgs(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh)) allocate(vsgs(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh)) allocate(wsgs(ib-ih:ie+ih,jb-jh:je+jh,kb:ke+kh)) allocate(sv1k(ib:ie,jb:je,kb:ke+kh)) allocate(sv2k(ib:ie,jb:je,kb:ke+kh)) allocate(sv3k(ib:ie,jb:je,kb:ke+kh)) allocate(sv4k(ib:ie,jb:je,kb:ke+kh)) allocate(sv1sgs(ib:ie,jb:je,kb:ke+kh)) allocate(sv2sgs(ib:ie,jb:je,kb:ke+kh)) allocate(sv3sgs(ib:ie,jb:je,kb:ke+kh)) allocate(sv4sgs(ib:ie,jb:je,kb:ke+kh)) allocate(wpsv1p(ib:ie,jb:je,kb:ke+kh)) allocate(wpsv2p(ib:ie,jb:je,kb:ke+kh)) allocate(wpsv3p(ib:ie,jb:je,kb:ke+kh)) allocate(wpsv4p(ib:ie,jb:je,kb:ke+kh)) allocate(sv1psv1pt(ib:ie,jb:je,kb:ke+kh)) allocate(sv2psv2pt(ib:ie,jb:je,kb:ke+kh)) allocate(sv3psv3pt(ib:ie,jb:je,kb:ke+kh)) allocate(sv4psv4pt(ib:ie,jb:je,kb:ke+kh)) allocate(PSS(ib:ie,jb:je,kb:ke+kh)) allocate(upwptik(ib:ie,jb:je,kb:ke+kh)) allocate(vpwptjk(ib:ie,jb:je,kb:ke+kh)) allocate(upvptij(ib:ie,jb:je,kb:ke+kh)) allocate(wpthlptk(ib:ie,jb:je,kb:ke+kh)) allocate(thlpthlpt(ib:ie,jb:je,kb:ke+kh)) allocate(upuptc(ib:ie,jb:je,kb:ke+kh)) allocate(vpvptc(ib:ie,jb:je,kb:ke+kh)) allocate(wpwptc(ib:ie,jb:je,kb:ke+kh)) allocate(tketc(ib:ie,jb:je,kb:ke+kh)) ! thlk=0.;qtk=0.;uik=0.;wik=0.;vjk=0.;wjk=0.;uij=0.;vij=0.;uc=0.;vc=0.;wc=0.;thlsgs=0.;qtsgs=0.;usgs=0.;vsgs=0.;wsgs=0.;sv1k=0.;sv2k=0.;sv3k=0.;sv4k=0.;sv1sgs=0.;sv2sgs=0.;sv3sgs=0.;sv4sgs=0.;wpsv1p=0.;wpsv2p=0. ! wpsv3p=0.;wpsv4p=0.;sv1psv1pt=0.;sv2psv2pt=0.;sv3psv3pt=0.;sv4psv4pt=0.;PSS=0.;upwptik=0.;vpwptjk=0.;upvptij=0.;wpthlptk=0.;thlpthlpt=0.;upuptc=0.;vpvptc=0.;wpwptc=0.;tketc=0. ! upwptyik=0.;wpthlptyk=0.;wpqtptyk=0.;wpsv1ptyk=0.;wpsv2ptyk=0.;wpsv3ptyk=0.;uwtyik=0.;wthltyk=0.;wqttyk=0.;wsv1tyk=0.;wsv2tyk=0.;wsv3tyk=0.;upuptyc=0.;wpwptyc=0.;thlpthlpty=0. ! qtpqtpty=0.;sv1psv1pty=0.;sv2psv2pty=0.;sv3psv3pty=0. if (.not. rk3step==3) return if (tsamplep > tsample) then if (lytdump .or. lydump .or. lxydump .or. lxytdump .or. ltdump .or. lmintdump) then ! wpthlptyk=0.;wpqtptyk=0.;wpsv1ptyk=0.;wpsv2ptyk=0. tstatsdumppi = 1./tstatsdumpp !> Perform required interpolations for flux calculations ! tg3315 for non-equidistant x and z-grids this needs to change do k=kb,ke+kh do j=jb,je do i=ib,ie uik(i,j,k) = 0.5*dzhi(k)*(um(i,j,k)*dzf(k-1) + um(i,j,k-1)*dzf(k)) wik(i,j,k) = 0.5*dxhi(i)*(wm(i,j,k)*dxf(i-1) + wm(i-1,j,k)*dxf(i)) vjk(i,j,k) = 0.5*dzhi(k)*(vm(i,j,k)*dzf(k-1) + vm(i,j,k-1)*dzf(k)) wjk(i,j,k) = 0.5* (wm(i,j,k) + wm(i,j-1,k)) uij(i,j,k) = 0.5* (um(i,j,k) + um(i,j-1,k)) vij(i,j,k) = 0.5*dxhi(i)*(vm(i,j,k)*dxf(i-1) + vm(i-1,j,k)*dxf(i)) uc (i,j,k) = 0.5*dxhi(i)*(um(i,j,k)*dxf(i-1) + um(i-1,j,k)*dxf(i)) vc (i,j,k) = 0.5* (vm(i,j,k) + vm(i,j-1,k)) wc (i,j,k) = 0.5*dzhi(k)*(wm(i,j,k)*dzf(k-1) + wm(i,j,k-1)*dzf(k)) ! SGS fluxes ! interps ekm to cell corner (uw) emom = ( dzf(k-1) * ( ekm(i,j,k)*dxf(i-1) + ekm(i-1,j,k)*dxf(i) ) + & dzf(k) * ( ekm(i,j,k-1)*dxf(i-1) + ekm(i-1,j,k-1)*dxf(i) ) )*dxhi(i) * dzhiq(k) usgs(i,j,k) = emom * ( (um(i,j,k)-um(i,j,k-1)) *dzhi(k) & +(wm(i,j,k)-wm(i-1,j,k)) *dxhi(i)) ! interps ekm to cell corner (vw) emom = ( dzf(k-1) * ( ekm(i,j,k) + ekm(i,j-1,k) ) + & dzf(k) * ( ekm(i,j,k-1) + ekm(i,j-1,k-1) ) ) * dzhiq(k) vsgs(i,j,k) = emom * ( (vm(i,j,k)-vm(i,j,k-1)) *dzhi(k) & +(wm(i,j,k)-wm(i,j-1,k)) *dyi) end do end do end do do k=kb,ke do j=jb,je do i=ib,ie wsgs(i,j,k) = ( ekm(i,j,k) * (wm(i,j,k+1)-wm(i,j,k)) *dzfi(k) & -ekm(i,j,k-1)* (wm(i,j,k)-wm(i,j,k-1)) *dzfi(k-1) ) * 2. & * dzhi(k) ! tg3315 check this end do end do end do if (ltempeq) then do k=kb,ke+kh do j=jb,je do i=ib,ie thlk(i,j,k) = 0.5*dzhi(k)*(thlm(i,j,k)*dzf(k-1) + thlm(i,j,k-1)*dzf(k)) end do end do end do do k=kb,ke !> SGS fluxes thlsgs(:,:,k) = 0.5 * (dzf(k-1)*ekh(:,:,k) + dzf(k)*ekh(:,:,k-1)) & * (thlm(:,:,k)-thlm(:,:,k-1)) * dzh2i(k) end do end if if (lmoist) then do k=kb,ke+kh do j=jb,je do i=ib,ie qtk(i,j,k) = 0.5*dzhi(k)*(qtm(i,j,k)*dzf(k-1) + qtm(i,j,k-1)*dzf(k)) end do end do end do do k=kb,ke !> SGS fluxes qtsgs(:,:,k) = 0.5 * (dzf(k-1)*ekh(:,:,k) + dzf(k)*ekh(:,:,k-1)) & * (qtm(:,:,k)-qtm(:,:,k-1)) * dzh2i(k) end do end if if (nsv>0) then do k=kb,ke do j=jb,je do i=ib,ie sv1k(i,j,k) = 0.5*dzhi(k)*(svm(i,j,k,1)*dzf(k-1) + svm(i,j,k-1,1)*dzf(k)) end do end do end do do k=kb,ke sv1sgs(ib:ie,jb:je,k) = 0.5 * (dzf(k-1)*ekh(ib:ie,jb:je,k) + dzf(k)*ekh(ib:ie,jb:je,k-1)) & * (svm(ib:ie,jb:je,k,1)-svm(ib:ie,jb:je,k-1,1)) * dzh2i(k) end do end if if (nsv>1) then do k=kb,ke+kh do j=jb,je do i=ib,ie sv2k(i,j,k) = 0.5*dzhi(k)*(svm(i,j,k,2)*dzf(k-1) + svm(i,j,k-1,2)*dzf(k)) end do end do end do do k=kb,ke sv2sgs(ib:ie,jb:je,k) = 0.5 * (dzf(k-1)*ekh(ib:ie,jb:je,k) + dzf(k)*ekh(ib:ie,jb:je,k-1)) & * (svm(ib:ie,jb:je,k,2)-svm(ib:ie,jb:je,k-1,2)) * dzh2i(k) end do end if if (nsv>2) then do k=kb,ke+kh do j=jb,je do i=ib,ie sv3k(i,j,k) = 0.5*dzhi(k)*(svm(i,j,k,3)*dzf(k-1) + svm(i,j,k-1,3)*dzf(k)) end do end do end do do k=kb,ke sv3sgs(ib:ie,jb:je,k) = 0.5 * (dzf(k-1)*ekh(ib:ie,jb:je,k) + dzf(k)*ekh(ib:ie,jb:je,k-1)) & * (svm(ib:ie,jb:je,k,3)-svm(ib:ie,jb:je,k-1,3)) * dzh2i(k) end do end if if (nsv>3) then do k=kb,ke+kh do j=jb,je do i=ib,ie sv4k(i,j,k) = 0.5*dzhi(k)*(svm(i,j,k,4)*dzf(k-1) + svm(i,j,k-1,4)*dzf(k)) end do end do end do do k=kb,ke sv4sgs(ib:ie,jb:je,k) = 0.5 * (dzf(k-1)*ekh(ib:ie,jb:je,k) + dzf(k)*ekh(ib:ie,jb:je,k-1)) & * (svm(ib:ie,jb:je,k,4)-svm(ib:ie,jb:je,k-1,4)) * dzh2i(k) end do end if if ((nsv>2) .and. (lchem .eqv. .true.)) then do k=kb,ke do j=jb,je do i=ib,ie if ((ABS(svm(i,j,k,2)) .gt. 1.e-40) .and. (IIc(i,j,k)==1)) then PSS(i,j,k) = ( ( (k1*(svm(i,j,k,1)/30.)*(svm(i,j,k,3)/48.))/(JNO2*(svm(i,j,k,2)/46.)) ) - 1 ) * 100 end if end do end do end do end if !!>> CALCS FOR INST. STATS !> Note: More computationally efficient to spatially average mean quantities first & ! for time dependant stats, hence the .or.s. Assuming homogeneity in y. !> Average in y-direction if (lydump .or. lytdump) then uy=0.;vy=0.;wy=0.;uwyik=0.;usgsy=0.;wsgsy=0.;thly=0.;wthlyk=0.;thlsgsy=0. qty=0.;wqtyk=0.;qtsgsy=0.;sca1y=0.;sv1sgsy=0.;sv2sgsy=0.;sca2y=0.;sca3y=0.;sv3sgsy=0. call avey_ibm(uy,um(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIu(ib:ie,jb:je,kb:ke),IIut(ib:ie,kb:ke)) call avey_ibm(vy,vm(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIv(ib:ie,jb:je,kb:ke),IIvt(ib:ie,kb:ke)) call avey_ibm(wy,wm(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(uwyik,uik(ib:ie,jb:je,kb:ke)*wik(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) call avey_ibm(usgsy,usgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) call avey_ibm(wsgsy,wsgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) if (ltempeq) then call avey_ibm(thly,thlm(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) call avey_ibm(wthlyk,wm(ib:ie,jb:je,kb:ke)*thlk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(thlsgsy,thlsgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) end if if (lmoist) then call avey_ibm(qty,qtm(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) call avey_ibm(wqtyk,wm(ib:ie,jb:je,kb:ke)*qtk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(qtsgsy,qtsgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) end if if(nsv>0) then call avey_ibm(sca1y,svm(ib:ie,jb:je,kb:ke,1),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) ! call avey_ibm(wsv1yk,wm(ib:ie,jb:je,kb:ke)*sv1k(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(sv1sgsy,sv1sgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) end if if (nsv>1) then call avey_ibm(sca2y,svm(ib:ie,jb:je,kb:ke,2),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) ! call avey_ibm(wsv2yk,wm(ib:ie,jb:je,kb:ke)*sv2k(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(sv2sgsy,sv2sgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) end if if (nsv>2) then call avey_ibm(sca3y,svm(ib:ie,jb:je,kb:ke,3),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) call avey_ibm(sv3sgsy,sv3sgs(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) end if end if ! lydump .or. lytdump if (lydump) then uwyik=0.;wthlyk=0.;uyik=0.;wyik=0.;thlyk=0.;wpthlpyk=0. call avey_ibm(uwyik,uik(ib:ie,jb:je,kb:ke)*wik(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) call avey_ibm(uyik,uik(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) call avey_ibm(wyik,wik(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) if (ltempeq) then call avey_ibm(wthlyk,wm(ib:ie,jb:je,kb:ke)*thlk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(thlyk,thlk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) end if upwpyik = uwyik - uyik*wyik if (ltempeq) then wpthlpyk = wthlyk - wy*thlyk where (IIwt==0) wpthlpyk = -999.0 endwhere end if where (IIuwt==0) upwpyik = -999.0 endwhere end if ! lydump !> tg3315 10.07.18 - in any case where averaging spatially can assume homogeneity and therefore average ! spatially first? Perhaps not due to UCL...? Would save space but goes against triple decomposition ! definition !> Average in x and y-direction if (lxydump .or. lxytdump) then uxy=0.;vxy=0.;wxy=0.;thlxy=0.;qtxy=0.;pxy=0.;usgsxy=0.;thlsgsxy=0.;sca1xy=0.;vsgsxy=0. !> Spatial averages of mean quantities call avexy_ibm(uxy(kb:ke+kh),um(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIu(ib:ie,jb:je,kb:ke+kh),IIus(kb:ke+kh),.false.) call avexy_ibm(vxy(kb:ke+kh),vm(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIv(ib:ie,jb:je,kb:ke+kh),IIvs(kb:ke+kh),.false.) call avexy_ibm(wxy(kb:ke+kh),wm(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.false.) if (ltempeq) then call avexy_ibm(thlxy(kb:ke+kh),thlm(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) call avexy_ibm(thlsgsxy(kb:ke+kh),thlsgs(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.false.) end if if (lmoist) then call avexy_ibm(qtxy(kb:ke+kh),qtm(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) end if call avexy_ibm(pxy(kb:ke+kh),pres0(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) call avexy_ibm(usgsxy(kb:ke+kh),usgs(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuw(ib:ie,jb:je,kb:ke+kh),IIuws(kb:ke+kh),.false.) call avexy_ibm(vsgsxy(kb:ke+kh),vsgs(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIvw(ib:ie,jb:je,kb:ke+kh),IIvws(kb:ke+kh),.false.) end if ! lxydump .or. lxytdump if (lxydump) then uwxyik=0.;vwxyjk=0.;uxyik=0.;wxyik=0.;vxyjk=0.;wxyjk=0.;wthlxyk=0.;thlxyk=0.;wpthlpxyk=0. call avexy_ibm(uwxyik(kb:ke+kh),uik(ib:ie,jb:je,kb:ke+kh)*wik(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuw(ib:ie,jb:je,kb:ke+kh),IIuws(kb:ke+kh),.true.) call avexy_ibm(vwxyjk(kb:ke+kh),vjk(ib:ie,jb:je,kb:ke+kh)*wjk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIvw(ib:ie,jb:je,kb:ke+kh),IIvws(kb:ke+kh),.true.) call avexy_ibm(uxyik(kb:ke+kh),uik(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuw(ib:ie,jb:je,kb:ke+kh),IIuws(kb:ke+kh),.true.) call avexy_ibm(wxyik(kb:ke+kh),wik(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuw(ib:ie,jb:je,kb:ke+kh),IIuws(kb:ke+kh),.true.) call avexy_ibm(wxyjk(kb:ke+kh),wjk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIvw(ib:ie,jb:je,kb:ke+kh),IIvws(kb:ke+kh),.true.) call avexy_ibm(vxyjk(kb:ke+kh),vjk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIvw(ib:ie,jb:je,kb:ke+kh),IIvws(kb:ke+kh),.true.) if (ltempeq) then call avexy_ibm(wthlxyk(kb:ke+kh),wm(ib:ie,jb:je,kb:ke+kh)*thlk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.true.) call avexy_ibm(thlxyk(kb:ke+kh),thlk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.true.) end if upwpxyik = uwxyik - uxyik*wxyik vpwpxyjk = vwxyjk - vxyjk*wxyjk if (ltempeq) then wpthlpxyk = wthlxyk - wxy*thlxyk end if end if ! lxydump !!>> CALCS FOR TIME DEPENDANT (AVERAGED) STATS !> Average 1-D fields in time if (lxytdump) then uxyt(kb:ke+kh) = (uxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + uxy(kb:ke+kh)*tsamplep)*tstatsdumppi vxyt(kb:ke+kh) = (vxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + vxy(kb:ke+kh)*tsamplep)*tstatsdumppi wxyt(kb:ke+kh) = (wxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + wxy(kb:ke+kh)*tsamplep)*tstatsdumppi qtxyt(kb:ke+kh) = (qtxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + qtxy(kb:ke+kh)*tsamplep)*tstatsdumppi pxyt(kb:ke+kh) = (pxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + pxy(kb:ke+kh)*tsamplep)*tstatsdumppi usgsxyt(kb:ke+kh) = (usgsxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + usgsxy(kb:ke+kh)*tsamplep)*tstatsdumppi vsgsxyt(kb:ke+kh) = (vsgsxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + vsgsxy(kb:ke+kh)*tsamplep)*tstatsdumppi if (ltempeq) then thlsgsxyt(kb:ke+kh) = (thlsgsxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + thlsgsxy(kb:ke+kh)*tsamplep)*tstatsdumppi thlxyt(kb:ke+kh) = (thlxyt(kb:ke+kh)*(tstatsdumpp-tsamplep) + thlxy(kb:ke+kh)*tsamplep)*tstatsdumppi end if end if ! lxytdump !> Average 2-D fields in time if (lytdump) then if (myid==0) then vyt(ib:ie,kb:ke) = (vyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + vy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi uyt(ib:ie,kb:ke) = (uyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + uy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi wyt(ib:ie,kb:ke) = (wyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + wy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi usgsyt(ib:ie,kb:ke) = (usgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + usgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi wsgsyt(ib:ie,kb:ke) = (wsgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + wsgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi if (ltempeq) then thlyt(ib:ie,kb:ke) = (thlyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + thly(ib:ie,kb:ke)*tsamplep)*tstatsdumppi thlsgsyt(ib:ie,kb:ke) = (thlsgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + thlsgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi end if if (lmoist) then qtyt(ib:ie,kb:ke) = (qtyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + qty(ib:ie,kb:ke)*tsamplep)*tstatsdumppi qtsgsyt(ib:ie,kb:ke) = (qtsgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + qtsgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi end if if (nsv>0) then sca1yt(ib:ie,kb:ke) = (sca1yt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + sca1y(ib:ie,kb:ke)*tsamplep)*tstatsdumppi sv1sgsyt(ib:ie,kb:ke) = (sv1sgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + sv1sgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi end if if (nsv>1) then sca2yt(ib:ie,kb:ke) = (sca2yt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + sca2y(ib:ie,kb:ke)*tsamplep)*tstatsdumppi sv2sgsyt(ib:ie,kb:ke) = (sv2sgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + sv2sgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi end if if (nsv>2) then sca3yt(ib:ie,kb:ke) = (sca3yt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + sca3y(ib:ie,kb:ke)*tsamplep)*tstatsdumppi sv3sgsyt(ib:ie,kb:ke) = (sv3sgsyt(ib:ie,kb:ke)*(tstatsdumpp-tsamplep) + sv3sgsy(ib:ie,kb:ke)*tsamplep)*tstatsdumppi end if end if ! myid end if !lytdump ! Average 3-D fields in time ! tg3315 may be possible to do less calculations by splitting up ! some calcs not necessary for xyt or yt... if (lxytdump .or. lytdump .or. ltdump .or. lmintdump) then uwtik(:,:,kb:ke+kh) = (uwtik(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + wik(:,:,kb:ke+kh)*uik(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi vwtjk(:,:,kb:ke+kh) = (vwtjk(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + wjk(:,:,kb:ke+kh)*vjk(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi uvtij(:,:,kb:ke+kh) = (uvtij(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + uij(:,:,kb:ke+kh)*vij(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi uutc(ib:ie,jb:je,kb:ke+kh) = (uutc(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + uc(ib:ie,jb:je,kb:ke+kh)*uc(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi vvtc(ib:ie,jb:je,kb:ke+kh) = (vvtc(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + vc(ib:ie,jb:je,kb:ke+kh)*vc(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wwtc(ib:ie,jb:je,kb:ke+kh) = (wwtc(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wc(ib:ie,jb:je,kb:ke+kh)*wc(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi utik(:,:,kb:ke+kh) = (utik(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + uik(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi wtik(:,:,kb:ke+kh) = (wtik(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + wik(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi vtjk(:,:,kb:ke+kh) = (vtjk(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + vjk(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi wtjk(:,:,kb:ke+kh) = (wtjk(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + wjk(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi utij(:,:,kb:ke+kh) = (utij(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + uij(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi vtij(:,:,kb:ke+kh) = (vtij(:,:,kb:ke+kh)*(tstatsdumpp-tsamplep) + vij(:,:,kb:ke+kh)*tsamplep)*tstatsdumppi umt(ib:ie,jb:je,kb:ke+kh) = (umt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + um(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi vmt(ib:ie,jb:je,kb:ke+kh) = (vmt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + vm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wmt(ib:ie,jb:je,kb:ke+kh) = (wmt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi utc(ib:ie,jb:je,kb:ke+kh) = (utc(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + uc(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi vtc(ib:ie,jb:je,kb:ke+kh) = (vtc(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + vc(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wtc(ib:ie,jb:je,kb:ke+kh) = (wtc(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wc(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi pt(ib:ie,jb:je,kb:ke+kh) = (pt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + pres0(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi if (ltempeq) then wthltk(ib:ie,jb:je,kb:ke+kh) = (wthltk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + thlk(ib:ie,jb:je,kb:ke+kh)*wm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi thlthlt(ib:ie,jb:je,kb:ke+kh) = (thlthlt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + thlm(ib:ie,jb:je,kb:ke+kh)*thlm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi thltk(ib:ie,jb:je,kb:ke+kh) = (thltk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + thlk(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi thlt(ib:ie,jb:je,kb:ke+kh) = (thlt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + thlm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi end if if (lmoist) then wqttk(ib:ie,jb:je,kb:ke+kh) = (wqttk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + qtk(ib:ie,jb:je,kb:ke+kh)*wm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi qtqtt(ib:ie,jb:je,kb:ke+kh) = (qtqtt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + qtm(ib:ie,jb:je,kb:ke+kh)*qtm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi qttk(ib:ie,jb:je,kb:ke+kh) = (qttk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + qtk(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi qtt(ib:ie,jb:je,kb:ke+kh) = (qtt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + qtm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi end if if (nsv>0) then sv1t(ib:ie,jb:je,kb:ke+kh) = (sv1t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,1)*tsamplep)*tstatsdumppi sv1tk(ib:ie,jb:je,kb:ke+kh) = (sv1tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv1k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wsv1tk(ib:ie,jb:je,kb:ke+kh) = (wsv1tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv1k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv1sgst(ib:ie,jb:je,kb:ke+kh) = (sv1sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv1sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv1sv1t(ib:ie,jb:je,kb:ke+kh) = (sv1sv1t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,1)*svm(ib:ie,jb:je,kb:ke+kh,1)*tsamplep)*tstatsdumppi ! sv1max(ib:ie,jb:je,kb:ke) = max(sv1max(ib:ie,jb:je,kb:ke),svm(ib:ie,jb:je,kb:ke,1)) end if if ((lchem .eqv. .true.) .and. (nsv>2)) then PSSt(ib:ie,jb:je,kb:ke+kh) = (PSSt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + PSS(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi end if if (nsv>1) then sv2t(ib:ie,jb:je,kb:ke+kh) = (sv2t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,2)*tsamplep)*tstatsdumppi sv2tk(ib:ie,jb:je,kb:ke+kh) = (sv2tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv2k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wsv2tk(ib:ie,jb:je,kb:ke+kh) = (wsv2tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv2k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv2sgst(ib:ie,jb:je,kb:ke+kh) = (sv2sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv2sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv2sv2t(ib:ie,jb:je,kb:ke+kh) = (sv2sv2t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,2)*svm(ib:ie,jb:je,kb:ke+kh,2)*tsamplep)*tstatsdumppi ! sv2max(ib:ie,jb:je,kb:ke) = max(sv2max(ib:ie,jb:je,kb:ke),svm(ib:ie,jb:je,kb:ke,2)) end if if (nsv>2) then sv3t(ib:ie,jb:je,kb:ke+kh) = (sv3t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,3)*tsamplep)*tstatsdumppi sv3tk(ib:ie,jb:je,kb:ke+kh) = (sv3tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv3k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wsv3tk(ib:ie,jb:je,kb:ke+kh) = (wsv3tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv3k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv3sgst(ib:ie,jb:je,kb:ke+kh) = (sv3sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv3sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv3sv3t(ib:ie,jb:je,kb:ke+kh) = (sv3sv3t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,3)*svm(ib:ie,jb:je,kb:ke+kh,3)*tsamplep)*tstatsdumppi ! sv3max(ib:ie,jb:je,kb:ke) = max(sv3max(ib:ie,jb:je,kb:ke),svm(ib:ie,jb:je,kb:ke,3)) end if if (nsv>3) then sv4t(ib:ie,jb:je,kb:ke+kh) = (sv4t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,4)*tsamplep)*tstatsdumppi sv4tk(ib:ie,jb:je,kb:ke+kh) = (sv4tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv4k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi wsv4tk(ib:ie,jb:je,kb:ke+kh) = (wsv4tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv4k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv4sgst(ib:ie,jb:je,kb:ke+kh) = (sv4sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv4sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi sv4sv4t(ib:ie,jb:je,kb:ke+kh) = (sv4sv4t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,4)*svm(ib:ie,jb:je,kb:ke+kh,4)*tsamplep)*tstatsdumppi ! sv4max(ib:ie,jb:je,kb:ke) = max(sv4max(ib:ie,jb:je,kb:ke),svm(ib:ie,jb:je,kb:ke,4)) end if end if !lxytdump .or. lytdump .or. ltdump ! Other 3-D fields specifically for tdump !if (ltdump) then ! bss116 already calculated above ! wmt(ib:ie,jb:je,kb:ke+kh) = (wmt(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv1t(ib:ie,jb:je,kb:ke+kh) = (sv1t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,1)*tsamplep)*tstatsdumppi ! sv2t(ib:ie,jb:je,kb:ke+kh) = (sv2t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,2)*tsamplep)*tstatsdumppi ! sv3t(ib:ie,jb:je,kb:ke+kh) = (sv3t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,3)*tsamplep)*tstatsdumppi ! ! sv4t(ib:ie,jb:je,kb:ke+kh) = (sv4t(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + svm(ib:ie,jb:je,kb:ke+kh,4)*tsamplep)*tstatsdumppi ! sv1tk(ib:ie,jb:je,kb:ke+kh) = (sv1tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv1k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv2tk(ib:ie,jb:je,kb:ke+kh) = (sv2tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv2k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv3tk(ib:ie,jb:je,kb:ke+kh) = (sv3tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv3k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv4tk(ib:ie,jb:je,kb:ke+kh) = (sv4tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv4k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! wsv1tk(ib:ie,jb:je,kb:ke+kh) = (wsv1tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv1k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! wsv2tk(ib:ie,jb:je,kb:ke+kh) = (wsv2tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv2k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! wsv3tk(ib:ie,jb:je,kb:ke+kh) = (wsv3tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv3k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! wsv4tk(ib:ie,jb:je,kb:ke+kh) = (wsv4tk(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + wm(ib:ie,jb:je,kb:ke+kh)*sv4k(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv1sgst(ib:ie,jb:je,kb:ke+kh) = (sv1sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv1sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv2sgst(ib:ie,jb:je,kb:ke+kh) = (sv2sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv2sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv3sgst(ib:ie,jb:je,kb:ke+kh) = (sv3sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv3sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi ! sv4sgst(ib:ie,jb:je,kb:ke+kh) = (sv4sgst(ib:ie,jb:je,kb:ke+kh)*(tstatsdumpp-tsamplep) + sv4sgs(ib:ie,jb:je,kb:ke+kh)*tsamplep)*tstatsdumppi !end if ! ltdump if (ltreedump) then tr_ut(ib:ie,jb:je,kb:ke) = (tr_ut(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_u(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi tr_vt(ib:ie,jb:je,kb:ke) = (tr_vt(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_v(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi tr_wt(ib:ie,jb:je,kb:ke) = (tr_wt(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_w(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi if (ltempeq) then tr_thlt(ib:ie,jb:je,kb:ke) = (tr_thlt(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_thl(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi end if if (lmoist) then tr_qtt(ib:ie,jb:je,kb:ke) = (tr_qtt(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_qt(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi tr_qtRt(ib:ie,jb:je,kb:ke) = (tr_qtRt(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_qtR(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi tr_qtAt(ib:ie,jb:je,kb:ke) = (tr_qtAt(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_qtA(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi tr_omegat(ib:ie,jb:je,kb:ke) = (tr_omegat(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_omega(ib:ie,jb:je,kb:ke)*tsamplep)*tstatsdumppi end if if (nsv>0) then tr_sv1t(ib:ie,jb:je,kb:ke) = (tr_sv1t(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_sv(ib:ie,jb:je,kb:ke,1)*tsamplep)*tstatsdumppi tr_sv2t(ib:ie,jb:je,kb:ke) = (tr_sv2t(ib:ie,jb:je,kb:ke)*(tstatsdumpp-tsamplep) + tr_sv(ib:ie,jb:je,kb:ke,2)*tsamplep)*tstatsdumppi end if end if ! where (IIuwt==0) ! upwpyik = -999 ! upwpytik = -999 ! endwhere end if ! lytdump .or. lydump .or. lxydump .or. lxytdump .or. ltdump .or. lmintdump ! Write y-averaged statistics every tsample if (lydump) then if (myid == 0) then allocate(field(ib:ie,kb:ke)) allocate(varsy(imax,khigh-klow+1,nstaty)) varsy = 0. varsy(:,:,1) = uy(ib:ie,kb:ke) varsy(:,:,2) = vy(ib:ie,kb:ke) varsy(:,:,3) = wy(ib:ie,kb:ke) varsy(:,:,4) = thly(ib:ie,kb:ke) varsy(:,:,5) = qty(ib:ie,kb:ke) varsy(:,:,6) = sca1y(ib:ie,kb:ke) varsy(:,:,7) = sca2y(ib:ie,kb:ke) varsy(:,:,8) = sca3y(ib:ie,kb:ke) varsy(:,:,9) = upwpyik(ib:ie,kb:ke) varsy(:,:,10) = wpthlpyk(ib:ie,kb:ke) varsy(:,:,11) = usgsy(ib:ie,kb:ke) varsy(:,:,12) = thlsgsy(ib:ie,kb:ke) varsy(:,:,13) = uwyik(ib:ie,kb:ke) varsy(:,:,14) = wthlyk(ib:ie,kb:ke) call writestat_nc(ncidy,1,tncstaty,(/timee/),nrecy,.true.) call writestat_nc(ncidy,nstaty,ncstaty,varsy,nrecy,imax,khigh-klow+1) deallocate(field,varsy) endif !myid endif !lydump ! Write xy-averaged statistics every tsample if (lxydump) then if (myid == 0) then call writestat_nc(ncidxy,1,tncstatxy,(/timee/),nrecxy,.true.) allocate(varsxy(khigh-klow+1,nstatxy)) varsxy(:,1) = uxy(kb:ke) varsxy(:,2) = vxy(kb:ke) varsxy(:,3) = wxy(kb:ke) varsxy(:,4) = thlxy(kb:ke) varsxy(:,5) = qtxy(kb:ke) varsxy(:,6) = pxy(kb:ke) varsxy(:,7) = upwpxyik(kb:ke) varsxy(:,8) = wpthlpxyk(kb:ke) varsxy(:,9) = vpwpxyjk(kb:ke) varsxy(:,10) = usgsxy(kb:ke) varsxy(:,11) = thlsgsxy(kb:ke) !wdthldtc(kb:ke) varsxy(:,12) = vsgsxy(kb:ke) varsxy(:,13) = uwxyik(kb:ke) varsxy(:,14) = wthlxyk(kb:ke) varsxy(:,15) = vwxyjk(kb:ke) call writestat_1D_nc(ncidxy,nstatxy,ncstatxy,varsxy,nrecxy,khigh-klow+1) end if !myid end if !lxydump if (lkslicedump) then allocate(varkslice(imax,jmax,nstatkslice)) call writestat_nc(ncidkslice,1,tncstatkslice,(/timee/),nreckslice,.true.) varkslice(:,:,1) = um(ib:ie,jb:je,kslice) varkslice(:,:,2) = vm(ib:ie,jb:je,kslice) varkslice(:,:,3) = 0.5*(wm(ib:ie,jb:je,kslice)+wm(ib:ie,jb:je,kslice+1)) ! assumes equidistant varkslice(:,:,4) = thlm(ib:ie,jb:je,kslice) varkslice(:,:,5) = qtm(ib:ie,jb:je,kslice) call writestat_nc(ncidkslice,nstatkslice,ncstatkslice,varkslice,nreckslice,imax,jmax) endif if (lislicedump) then if (islicerank) then allocate(varislice(jmax,khigh-klow+1,nstatislice)) call writestat_nc(ncidislice,1,tncstatislice,(/timee/),nrecislice,.true.) varislice(:,:,1) = 0.5*(um(isliceloc,jb:je,kb:ke)+um(isliceloc+1,jb:je,kb:ke)) varislice(:,:,2) = vm(isliceloc,jb:je,kb:ke) varislice(:,:,3) = wm(isliceloc,jb:je,kb:ke) varislice(:,:,4) = thlm(isliceloc,jb:je,kb:ke) varislice(:,:,5) = qtm(isliceloc,jb:je,kb:ke) call writestat_nc(ncidislice,nstatislice,ncstatislice,varislice,nrecislice,jmax,khigh-klow+1) endif endif if (ljslicedump) then if (jslicerank) then allocate(varjslice(imax,khigh-klow+1,nstatjslice)) call writestat_nc(ncidjslice,1,tncstatjslice,(/timee/),nrecjslice,.true.) varjslice(:,:,1) = um(ib:ie,jsliceloc,kb:ke) varjslice(:,:,2) = 0.5*(vm(ib:ie,jsliceloc,kb:ke)+vm(ib:ie,jsliceloc+1,kb:ke)) varjslice(:,:,3) = wm(ib:ie,jsliceloc,kb:ke) varjslice(:,:,4) = thlm(ib:ie,jsliceloc,kb:ke) varjslice(:,:,5) = qtm(ib:ie,jsliceloc,kb:ke) call writestat_nc(ncidjslice,nstatjslice,ncstatjslice,varjslice,nrecjslice,imax,khigh-klow+1) endif endif if (ltkedump) then !call genstats(tsamplep,tstatsdumpp,umc,vmc,wmc) endif tsamplep = dt else !timestatsdumpp < tsample tsamplep = tsamplep + dt endif if (tstatsdumpp > tstatsdump) then ! Final calculations and write xyt-averaged statistics every tsample if (lxytdump) then wthltxyk=0.;uwtxyik=0.;vwtxyjk=0.;wwtxyk=0.;uvtxyij=0.;wpthlptxyk=0.;upwptxyik=0.;vpwptxyjk=0.;upvptxyij=0.;thlpthlptxy=0.;upuptxyc=0.;vpvptxyc=0.;wpwptxyc=0.;tketxyc=0. !> Advective flux if (ltempeq) then call avexy_ibm(wthltxyk(kb:ke+kh),wmt(ib:ie,jb:je,kb:ke+kh)*thltk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.false.) end if call avexy_ibm(uwtxyik(kb:ke+kh),utik(ib:ie,jb:je,kb:ke+kh)*wtik(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuw(ib:ie,jb:je,kb:ke+kh),IIuws(kb:ke+kh),.false.) call avexy_ibm(vwtxyjk(kb:ke+kh),vtjk(ib:ie,jb:je,kb:ke+kh)*wtjk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIvw(ib:ie,jb:je,kb:ke+kh),IIvws(kb:ke+kh),.false.) call avexy_ibm(wwtxyk(kb:ke+kh),wmt(ib:ie,jb:je,kb:ke+kh)*wmt(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.false.) call avexy_ibm(uvtxyij(kb:ke+kh),utij(ib:ie,jb:je,kb:ke+kh)*vtij(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuv(ib:ie,jb:je,kb:ke+kh),IIuvs(kb:ke+kh),.false.) !> Turbulent fluxes if (ltempeq) then call avexy_ibm(wpthlptxyk(kb:ke+kh),wthltk(ib:ie,jb:je,kb:ke+kh)-wmt(ib:ie,jb:je,kb:ke+kh)*thltk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIw(ib:ie,jb:je,kb:ke+kh),IIws(kb:ke+kh),.false.) end if call avexy_ibm(upwptxyik(kb:ke+kh),uwtik(ib:ie,jb:je,kb:ke+kh)-utik(ib:ie,jb:je,kb:ke+kh)*wtik(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuw(ib:ie,jb:je,kb:ke+kh),IIuws(kb:ke+kh),.false.) call avexy_ibm(vpwptxyjk(kb:ke+kh),vwtjk(ib:ie,jb:je,kb:ke+kh)-vtjk(ib:ie,jb:je,kb:ke+kh)*wtjk(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIvw(ib:ie,jb:je,kb:ke+kh),IIvws(kb:ke+kh),.false.) call avexy_ibm(upvptxyij(kb:ke+kh),uvtij(ib:ie,jb:je,kb:ke+kh)-utij(ib:ie,jb:je,kb:ke+kh)*vtij(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIuv(ib:ie,jb:je,kb:ke+kh),IIuvs(kb:ke+kh),.false.) !> Variances and TKE if (ltempeq) then call avexy_ibm(thlpthlptxy(kb:ke+kh),thlthlt(ib:ie,jb:je,kb:ke+kh)-thlt(ib:ie,jb:je,kb:ke+kh)*thlt(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) end if call avexy_ibm(upuptxyc(kb:ke+kh),uutc(ib:ie,jb:je,kb:ke+kh)-utc(ib:ie,jb:je,kb:ke+kh)*utc(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) call avexy_ibm(vpvptxyc(kb:ke+kh),vvtc(ib:ie,jb:je,kb:ke+kh)-vtc(ib:ie,jb:je,kb:ke+kh)*vtc(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) call avexy_ibm(wpwptxyc(kb:ke+kh),wwtc(ib:ie,jb:je,kb:ke+kh)-wtc(ib:ie,jb:je,kb:ke+kh)*wtc(ib:ie,jb:je,kb:ke+kh),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) call avexy_ibm(tketxyc(kb:ke+kh),0.5*((wwtc(ib:ie,jb:je,kb:ke+kh)-wtc(ib:ie,jb:je,kb:ke+kh)*wtc(ib:ie,jb:je,kb:ke+kh))+(vvtc(ib:ie,jb:je,kb:ke+kh)-vtc(ib:ie,jb:je,kb:ke+kh)*vtc(ib:ie,jb:je,kb:ke+kh))+(uutc(ib:ie,jb:je,kb:ke+kh)-utc(ib:ie,jb:je,kb:ke+kh)*utc(ib:ie,jb:je,kb:ke+kh))),ib,ie,jb,je,kb,ke,ih,jh,kh,IIc(ib:ie,jb:je,kb:ke+kh),IIcs(kb:ke+kh),.false.) if (myid == 0) then call writestat_nc(ncidxyt,1,tncstatxyt,(/timee/),nrecxyt,.true.) allocate(varsxyt(khigh-klow+1,nstatxyt)) varsxyt(:,1) = uxyt(kb:ke) varsxyt(:,2) = vxyt(kb:ke) varsxyt(:,3) = wxyt(kb:ke) varsxyt(:,4) = thlxyt(kb:ke) varsxyt(:,5) = qtxyt(kb:ke) varsxyt(:,6) = pxyt(kb:ke) varsxyt(:,7) = upwptxyik(kb:ke) varsxyt(:,8) = wpthlptxyk(kb:ke) varsxyt(:,9) = vpwptxyjk(kb:ke) varsxyt(:,10) = upvptxyij(kb:ke) varsxyt(:,11) = uwtxyik(kb:ke) varsxyt(:,12) = wthltxyk(kb:ke) !wdthldtc(kb:ke) varsxyt(:,13) = uvtxyij(kb:ke) varsxyt(:,14) = vwtxyjk(kb:ke) varsxyt(:,15) = wwtxyk(kb:ke) varsxyt(:,16) = usgsxyt(kb:ke) !wdthldtw(kb:ke) varsxyt(:,17) = thlsgsxyt(kb:ke) varsxyt(:,18) = vsgsxyt(kb:ke) varsxyt(:,19) = thlpthlptxy(kb:ke) varsxyt(:,20) = upuptxyc(kb:ke) varsxyt(:,21) = vpvptxyc(kb:ke) varsxyt(:,22) = wpwptxyc(kb:ke) varsxyt(:,23) = tketxyc(kb:ke) call writestat_1D_nc(ncidxyt,nstatxyt,ncstatxyt,varsxyt,nrecxyt,khigh-klow+1) end if !myid end if !lxytdump ! Final calculations and write yt-averaged statistics every tsample if (lytdump) then ! call MPI_BCAST(sca1yt ,(ke+kh-(kb-kh))*(ie+ih-(ib-ih)),MY_REAL ,7,comm3d,mpierr) ! Turbulent flux call avey_ibm(upwptyik,uwtik(ib:ie,jb:je,kb:ke)-utik(ib:ie,jb:je,kb:ke)*wtik(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) call avey_ibm(uwtyik,utik(ib:ie,jb:je,kb:ke)*wtik(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIuw(ib:ie,jb:je,kb:ke),IIuwt(ib:ie,kb:ke)) call avey_ibm(upuptyc,uutc(ib:ie,jb:je,kb:ke)-utc(ib:ie,jb:je,kb:ke)*utc(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) call avey_ibm(vpvptyc,vvtc(ib:ie,jb:je,kb:ke)-vtc(ib:ie,jb:je,kb:ke)*vtc(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) call avey_ibm(wpwptyc,wwtc(ib:ie,jb:je,kb:ke)-wtc(ib:ie,jb:je,kb:ke)*wtc(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) if (ltempeq) then call avey_ibm(wpthlptyk,wthltk(ib:ie,jb:je,kb:ke)-wmt(ib:ie,jb:je,kb:ke)*thltk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(wthltyk,wmt(ib:ie,jb:je,kb:ke)*thltk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(thlpthlpty,thlthlt(ib:ie,jb:je,kb:ke)-thlt(ib:ie,jb:je,kb:ke)*thlt(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) end if if (lmoist) then call avey_ibm(wpqtptyk,wqttk(ib:ie,jb:je,kb:ke)-wmt(ib:ie,jb:je,kb:ke)*qttk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(wqttyk,wmt(ib:ie,jb:je,kb:ke)*qttk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(qtpqtpty,qtqtt(ib:ie,jb:je,kb:ke)-qtt(ib:ie,jb:je,kb:ke)*qtt(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) end if if (nsv>0) then call avey_ibm(wpsv1ptyk,wsv1tk(ib:ie,jb:je,kb:ke)-wmt(ib:ie,jb:je,kb:ke)*sv1tk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(wsv1tyk,wmt(ib:ie,jb:je,kb:ke)*sv1tk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(sv1psv1pty,sv1sv1t(ib:ie,jb:je,kb:ke)-sv1t(ib:ie,jb:je,kb:ke)*sv1t(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) end if if (nsv>1) then call avey_ibm(wpsv2ptyk,wsv2tk(ib:ie,jb:je,kb:ke)-wmt(ib:ie,jb:je,kb:ke)*sv2tk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(wsv2tyk,wmt(ib:ie,jb:je,kb:ke)*sv2tk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(sv2psv2pty,sv2sv2t(ib:ie,jb:je,kb:ke)-sv2t(ib:ie,jb:je,kb:ke)*sv2t(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) end if if (nsv>2) then call avey_ibm(wpsv3ptyk,wsv3tk(ib:ie,jb:je,kb:ke)-wmt(ib:ie,jb:je,kb:ke)*sv3tk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(wsv3tyk,wmt(ib:ie,jb:je,kb:ke)*sv3tk(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIw(ib:ie,jb:je,kb:ke),IIwt(ib:ie,kb:ke)) call avey_ibm(sv3psv3pty,sv3sv3t(ib:ie,jb:je,kb:ke)-sv3t(ib:ie,jb:je,kb:ke)*sv3t(ib:ie,jb:je,kb:ke),ib,ie,jb,je,kb,ke,IIc(ib:ie,jb:je,kb:ke),IIct(ib:ie,kb:ke)) end if if (myid == 0) then allocate(varsyt(imax,khigh-klow+1,nstatyt)) varsyt = 0. call writestat_nc(ncidyt,1,tncstatyt,(/timee/),nrecyt,.true.) varsyt(:,:,1) = uyt(ib:ie,kb:ke) varsyt(:,:,2) = vyt(ib:ie,kb:ke) varsyt(:,:,3) = wyt(ib:ie,kb:ke) varsyt(:,:,4) = thlyt(ib:ie,kb:ke) varsyt(:,:,5) = qtyt(ib:ie,kb:ke) varsyt(:,:,6) = sca1yt(ib:ie,kb:ke) varsyt(:,:,7) = sca2yt(ib:ie,kb:ke) varsyt(:,:,8) = sca3yt(ib:ie,kb:ke) varsyt(:,:,9) = upwptyik(ib:ie,kb:ke) varsyt(:,:,10) = wpthlptyk(ib:ie,kb:ke) varsyt(:,:,11) = wpqtptyk(ib:ie,kb:ke) varsyt(:,:,12) = wpsv1ptyk(ib:ie,kb:ke) varsyt(:,:,13) = wpsv2ptyk(ib:ie,kb:ke) varsyt(:,:,14) = wpsv3ptyk(ib:ie,kb:ke) varsyt(:,:,15) = uwtyik(ib:ie,kb:ke) varsyt(:,:,16) = wthltyk(ib:ie,kb:ke) varsyt(:,:,17) = wqttyk(ib:ie,kb:ke) varsyt(:,:,18) = wsv1tyk(ib:ie,kb:ke) varsyt(:,:,19) = wsv2tyk(ib:ie,kb:ke) varsyt(:,:,20) = wsv3tyk(ib:ie,kb:ke) varsyt(:,:,21) = upuptyc(ib:ie,kb:ke) varsyt(:,:,22) = wpwptyc(ib:ie,kb:ke) varsyt(:,:,23) = thlpthlpty(ib:ie,kb:ke) varsyt(:,:,24) = qtpqtpty(ib:ie,kb:ke) varsyt(:,:,25) = sv1psv1pty(ib:ie,kb:ke) varsyt(:,:,26) = sv2psv2pty(ib:ie,kb:ke) varsyt(:,:,27) = sv3psv3pty(ib:ie,kb:ke) varsyt(:,:,28) = usgsyt(ib:ie,kb:ke) varsyt(:,:,29) = wsgsyt(ib:ie,kb:ke) varsyt(:,:,30) = thlsgsyt(ib:ie,kb:ke) varsyt(:,:,31) = qtsgsyt(ib:ie,kb:ke) varsyt(:,:,32) = sv1sgsyt(ib:ie,kb:ke) varsyt(:,:,33) = sv2sgsyt(ib:ie,kb:ke) varsyt(:,:,34) = sv3sgsyt(ib:ie,kb:ke) call writestat_nc(ncidyt,nstatyt,ncstatyt,varsyt,nrecyt,imax,khigh-klow+1) end if !myid end if !lytdump ! Final calculations and write t-averaged statistics every tsample if (ltdump) then ! wpsv1p = wsv1tk - wmt*sv1tk ! wpsv2p = wsv2tk - wmt*sv2tk ! wpsv3p = wsv3tk - wmt*sv3tk ! wpsv4p = wsv4tk - wmt*sv4tk wpthlptk=0.;thlpthlpt=0. !> Turbulent fluxes upwptik = uwtik - utik*wtik vpwptjk = vwtjk - vtjk*wtjk upvptij = uvtij - utij*vtij if (ltempeq) then wpthlptk = wthltk - wmt*thlk end if if (nsv>0) then wpsv1p = wsv1tk - wmt*sv1tk end if if (nsv>1) then wpsv2p = wsv2tk - wmt*sv2tk end if if (nsv>2) then wpsv3p = wsv3tk - wmt*sv3tk end if if (nsv>3) then wpsv4p = wsv4tk - wmt*sv4tk end if !> Variances and TKE if (ltempeq) then thlpthlpt = thlthlt - thlt*thlt end if upuptc = uutc - utc*utc vpvptc = vvtc - vtc*vtc wpwptc = wwtc - wtc*wtc tketc = 0.5*(upuptc + vpvptc + wpwptc) if (nsv>0) then sv1psv1pt = sv1sv1t - sv1t*sv1t end if if (nsv>1) then sv2psv2pt = sv2sv2t - sv2t*sv2t end if if (nsv>2) then sv3psv3pt = sv3sv3t - sv3t*sv3t end if if (nsv>3) then sv4psv4pt = sv4sv4t - sv4t*sv4t end if ! if (myid == 0) then allocate(varst(imax,jmax,khigh-klow+1,nstatt)) call writestat_nc(ncidt,1,tncstatt,(/timee/),nrect,.true.) varst(:,:,:,1) = umt(ib:ie,jb:je,kb:ke) varst(:,:,:,2) = vmt(ib:ie,jb:je,kb:ke) varst(:,:,:,3) = wmt(ib:ie,jb:je,kb:ke) varst(:,:,:,4) = thlt(ib:ie,jb:je,kb:ke) varst(:,:,:,5) = qtt(ib:ie,jb:je,kb:ke) varst(:,:,:,6) = pt(ib:ie,jb:je,kb:ke) varst(:,:,:,7) = sv1t(ib:ie,jb:je,kb:ke) varst(:,:,:,8) = sv2t(ib:ie,jb:je,kb:ke) varst(:,:,:,9) = sv3t(ib:ie,jb:je,kb:ke) varst(:,:,:,10) = sv4t(ib:ie,jb:je,kb:ke) varst(:,:,:,11) = PSSt(ib:ie,jb:je,kb:ke) varst(:,:,:,12) = upwptik(ib:ie,jb:je,kb:ke) varst(:,:,:,13) = vpwptjk(ib:ie,jb:je,kb:ke) varst(:,:,:,14) = upvptij(ib:ie,jb:je,kb:ke) varst(:,:,:,15) = wpthlptk(ib:ie,jb:je,kb:ke) varst(:,:,:,16) = wpsv1p(ib:ie,jb:je,kb:ke) varst(:,:,:,17) = wpsv2p(ib:ie,jb:je,kb:ke) varst(:,:,:,18) = wpsv3p(ib:ie,jb:je,kb:ke) varst(:,:,:,19) = wpsv4p(ib:ie,jb:je,kb:ke) varst(:,:,:,20) = thlpthlpt(ib:ie,jb:je,kb:ke) varst(:,:,:,21) = upuptc(ib:ie,jb:je,kb:ke) varst(:,:,:,22) = vpvptc(ib:ie,jb:je,kb:ke) varst(:,:,:,23) = wpwptc(ib:ie,jb:je,kb:ke) varst(:,:,:,24) = tketc(ib:ie,jb:je,kb:ke) varst(:,:,:,25) = sv1psv1pt(ib:ie,jb:je,kb:ke) varst(:,:,:,26) = sv2psv2pt(ib:ie,jb:je,kb:ke) varst(:,:,:,27) = sv3psv3pt(ib:ie,jb:je,kb:ke) varst(:,:,:,28) = sv4psv4pt(ib:ie,jb:je,kb:ke) varst(:,:,:,29) = sv1sgst(ib:ie,jb:je,kb:ke) varst(:,:,:,30) = sv2sgst(ib:ie,jb:je,kb:ke) varst(:,:,:,31) = sv3sgst(ib:ie,jb:je,kb:ke) varst(:,:,:,32) = sv4sgst(ib:ie,jb:je,kb:ke) ! varst(:,:,:,33) = sv1max(ib:ie,jb:je,kb:ke) ! varst(:,:,:,34) = sv2max(ib:ie,jb:je,kb:ke) ! varst(:,:,:,35) = sv3max(ib:ie,jb:je,kb:ke) ! varst(:,:,:,36) = sv4max(ib:ie,jb:je,kb:ke) call writestat_nc(ncidt,nstatt,ncstatt,varst,nrect,imax,jmax,khigh-klow+1) ! end if !myid deallocate(varst) end if !ltdump if (lmintdump) then ! if (myid == 0) then allocate(varsmint(imax,jmax,khigh-klow+1,nstatmint)) call writestat_nc(ncidmint,1,tncstatmint,(/timee/),nrecmint,.true.) varsmint(:,:,:,1) = umt(ib:ie,jb:je,kb:ke) varsmint(:,:,:,2) = vmt(ib:ie,jb:je,kb:ke) varsmint(:,:,:,3) = wmt(ib:ie,jb:je,kb:ke) varsmint(:,:,:,4) = thlt(ib:ie,jb:je,kb:ke) varsmint(:,:,:,5) = qtt(ib:ie,jb:je,kb:ke) varsmint(:,:,:,6) = pt(ib:ie,jb:je,kb:ke) call writestat_nc(ncidmint,nstatmint,ncstatmint,varsmint,nrecmint,imax,jmax,khigh-klow+1) ! end if !myid deallocate(varsmint) end if !lmintdump ! Final calculations and write t-averaged statistics for the trees if (ltreedump) then ! if (myid == 0) then allocate(varstr(imax,jmax,khigh-klow+1,nstattr)) call writestat_nc(ncidtr,1,tncstattr,(/timee/),nrectr,.true.) varstr(:,:,:,1) = tr_ut(ib:ie,jb:je,kb:ke) varstr(:,:,:,2) = tr_vt(ib:ie,jb:je,kb:ke) varstr(:,:,:,3) = tr_wt(ib:ie,jb:je,kb:ke) varstr(:,:,:,4) = tr_thlt(ib:ie,jb:je,kb:ke) varstr(:,:,:,5) = tr_qtt(ib:ie,jb:je,kb:ke) varstr(:,:,:,6) = tr_qtRt(ib:ie,jb:je,kb:ke) varstr(:,:,:,7) = tr_qtAt(ib:ie,jb:je,kb:ke) varstr(:,:,:,8) = tr_sv1t(ib:ie,jb:je,kb:ke) varstr(:,:,:,9) = tr_sv2t(ib:ie,jb:je,kb:ke) varstr(:,:,:,10) = tr_omegat(ib:ie,jb:je,kb:ke) call writestat_nc(ncidtr,nstattr,ncstattr,varstr,nrectr,imax,jmax,khigh-klow+1) ! end if !myid deallocate(varstr) end if !ltdump if (ltkedump) then call tkestatsdump if (myid == 0) then call writestat_nc(ncidtke,1,tncstattke,(/timee/),nrectke,.true.) allocate(varstke(khigh-klow+1,nstattke)) varstke(:,1) = p_b(kb:ke+kh) varstke(:,2) = t_p(kb:ke+kh) varstke(:,3) = adv(kb:ke+kh) varstke(:,4) = t_t(kb:ke+kh) varstke(:,5) = t_sgs(kb:ke+kh) varstke(:,6) = p_t(kb:ke+kh) varstke(:,7) = t_v(kb:ke+kh) varstke(:,8) = d_sgs(kb:ke+kh) call writestat_1D_nc(ncidtke,nstattke,ncstattke,varstke,nrectke,khigh-klow+1) end if !myid endif !ltkedump tstatsdumpp = dt else !tstatsdumpp < tstatsdump tstatsdumpp = tstatsdumpp + dt endif deallocate(thlk,qtk,uik,wik,vjk,wjk,uij,vij,uc,vc,wc) deallocate(thlsgs,qtsgs,usgs,vsgs,wsgs) deallocate(sv1k,sv2k,sv3k,sv4k,sv1sgs,sv2sgs,sv3sgs,sv4sgs,PSS,wpsv1p,wpsv2p,wpsv3p,wpsv4p,sv1psv1pt,sv2psv2pt,sv3psv3pt,sv4psv4pt) deallocate(upwptik,vpwptjk,upvptij,wpthlptk,thlpthlpt,upuptc,vpvptc,wpwptc,tketc) end subroutine statsdump