subroutine scandyn (ztodt, etadot, etamid, grlps1, grt1, & grz1, grd1, grfu1, grfv1, grut1, & grvt1, grrh1, grlps2, grt2, grz2, & grd2, grfu2, grfv2, grut2, grvt2, & grrh2, vcour, vmax2d, vmax2dt, detam, & cwava, flx_net, t2, fu, fv) !----------------------------------------------------------------------- ! ! Purpose: ! ! Method: ! "After coupling" gaussian latitude scan for which some of the physics ! and nonlinear dynamics calculations are completed. The main loop over ! latitude in this routine is multitasked. ! ! Note: the "ifdef" constructs in this routine are associated with the ! message-passing version of CAM. Messages are sent which ! have no relevance to the shared-memory case. ! ! Author: ! Original version: CCM3 !----------------------------------------------------------------------- ! ! $Id$ ! $Author$ ! !----------------------------------------------------------------------- use shr_kind_mod, only: r8 => shr_kind_r8 use pmgrid, only: plon, plat, plev, beglat, endlat, plevp use prognostics, only: u3, v3, q3, t3, div, vort, phis, omga, dpsl, & dpsm, ps, n3m1, n3, n3m2, qminus, pdeld use constituents, only: pcnst use scanslt, only: hw1lat use rgrid, only: nlon use comspe, only: maxm use linemsdyn, only: linemsdyn_bft, linemsdyn_fft, linemsdyn_aft, & plondfft use commap, only: w use qmassa, only: qmassarun use perf_mod !----------------------------------------------------------------------- implicit none !----------------------------------------------------------------------- ! ! Input arguments ! real(r8), intent(in) :: ztodt ! two delta t unless nstep =0 real(r8), intent(inout) :: etadot(plon,plevp,beglat:endlat) ! vertical motion (slt) real(r8), intent(in) :: etamid(plev) ! hybrd coord value at levels real(r8), intent(in) :: detam(plev) ! ! Fourier coefficient arrays which have a latitude index on them for ! multitasking. These arrays are defined in LINEMSDYN and and used in QUAD ! to compute spectral coefficients. They contain a latitude index so ! that the sums over latitude can be performed in a specified order. ! real(r8), intent(in) :: cwava(plat) ! weight applied to global integrals real(r8), intent(in) :: flx_net(plon,beglat:endlat) ! net flx from physics real(r8), intent(inout) :: t2(plon,plev,beglat:endlat) ! tot dT/dt to to physics real(r8), intent(inout) :: fu(plon,plev,beglat:endlat) ! u wind tend real(r8), intent(inout) :: fv(plon,plev,beglat:endlat) ! v wind tend ! ! Output arguments ! real(r8), intent(out) :: grlps1(2*maxm,(plat+1)/2) ! sym. undiff. term in lnps eqn. real(r8), intent(out) :: grlps2(2*maxm,(plat+1)/2) ! antisym undiff. term in lnps eqn. real(r8), intent(out) :: grt1(2*maxm,plev,(plat+1)/2) ! sym. undiff. term in t eqn. real(r8), intent(out) :: grt2(2*maxm,plev,(plat+1)/2) ! antisym. undiff. term in t eqn. real(r8), intent(out) :: grz1(2*maxm,plev,(plat+1)/2) ! sym. undiff. term in z eqn. real(r8), intent(out) :: grz2(2*maxm,plev,(plat+1)/2) ! antisym. undiff. term in z eqn. real(r8), intent(out) :: grd1(2*maxm,plev,(plat+1)/2) ! sym. undiff. term in d eqn. real(r8), intent(out) :: grd2(2*maxm,plev,(plat+1)/2) ! antisym. undiff. term in d eqn. real(r8), intent(out) :: grfu1(2*maxm,plev,(plat+1)/2) ! sym. nonlinear terms in u eqn. real(r8), intent(out) :: grfu2(2*maxm,plev,(plat+1)/2) ! antisym. nonlinear terms in u eqn. real(r8), intent(out) :: grfv1(2*maxm,plev,(plat+1)/2) ! sym. nonlinear terms in v eqn. real(r8), intent(out) :: grfv2(2*maxm,plev,(plat+1)/2) ! antisym. nonlinear terms in v eqn. real(r8), intent(out) :: grut1(2*maxm,plev,(plat+1)/2) ! sym. lambda deriv. term in t eqn. real(r8), intent(out) :: grut2(2*maxm,plev,(plat+1)/2) ! antisym. lambda deriv. term in t eqn. real(r8), intent(out) :: grvt1(2*maxm,plev,(plat+1)/2) ! sym. mu derivative term in t eqn. real(r8), intent(out) :: grvt2(2*maxm,plev,(plat+1)/2) ! antisym. mu deriv. term in t eqn. real(r8), intent(out) :: grrh1(2*maxm,plev,(plat+1)/2) ! sym. del**2 term in d eqn. real(r8), intent(out) :: grrh2(2*maxm,plev,(plat+1)/2) ! antisym. del**2 term in d eqn. real(r8), intent(out) :: vcour(plev,plat) ! maximum Courant number in vert. real(r8), intent(out) :: vmax2d(plev,plat) ! max. wind at each level, latitude real(r8), intent(out) :: vmax2dt(plev,plat) ! max. truncated wind at each lvl,lat ! Local variables integer irow ! latitude pair index integer lat,latn,lats ! latitude indices integer nlon_fft_in ! FFT work array inner dimension integer nlon_fft_out ! FFT work array inner dimension real(r8) pmid(plon,plev) ! pressure at model levels real(r8) pint(plon,plevp) ! pressure at interfaces real(r8) pdel(plon,plev) ! pressure difference between integer :: m ! constituent index ! ! FFT buffers ! real(r8), allocatable:: fftbuf_in(:,:,:,:) ! fftbuf_in(nlon_fft_in,9,plev,beglat:endlat) real(r8), allocatable:: fftbuf_out(:,:,:,:) ! fftbuf_out(nlon_fft_out,9,plev,plat) ! call t_startf ('scandyn_alloc') nlon_fft_in = plondfft allocate(fftbuf_in(nlon_fft_in,9,plev,beglat:endlat)) #if ( defined SPMD ) #ifdef NEC_SX nlon_fft_out = 2*maxm + 1 #else nlon_fft_out = 2*maxm #endif allocate(fftbuf_out(nlon_fft_out,9,plev,plat)) #else nlon_fft_out = 1 allocate(fftbuf_out(1,1,1,1)) #endif call t_stopf ('scandyn_alloc') ! call t_startf ('linemsdyn_bft') #ifdef OUTER_OMP !$OMP PARALLEL DO PRIVATE (LAT) #endif do lat=beglat,endlat call linemsdyn_bft (lat, nlon(lat), nlon_fft_in, & ps(1,lat,n3m1), ps(1,lat,n3m2), u3(1,1,lat,n3m1), & u3(1,1,lat,n3m2), v3(1,1,lat,n3m1), v3(1,1,lat,n3m2), t3(1,1,lat,n3m1), t3(1,1,lat,n3m2), & q3(1,1,1,lat,n3m1), etadot(1,1,lat), etamid, & ztodt, vcour(1,lat), vmax2d(1,lat), vmax2dt(1,lat), & detam, t2(1,1,lat), fu(1,1,lat), fv(1,1,lat), & div(1,1,lat,n3m1), vort(1,1,lat,n3m2), div(1,1,lat,n3m2), vort(1,1,lat,n3m1), & phis(1,lat), dpsl(1,lat), dpsm(1,lat), omga(1,1,lat), & cwava(lat), flx_net(1,lat), fftbuf_in(1,1,1,lat) ) end do call t_stopf ('linemsdyn_bft') call t_startf ('linemsdyn_fft') call linemsdyn_fft (nlon_fft_in,nlon_fft_out,fftbuf_in,fftbuf_out) call t_stopf ('linemsdyn_fft') call t_startf ('linemsdyn_aft') !$OMP PARALLEL DO PRIVATE (IROW, LATN, LATS) do irow=1,plat/2 lats = irow latn = plat - irow + 1 #if ( defined SPMD ) call linemsdyn_aft (irow, nlon_fft_out, fftbuf_out(1,1,1,lats), fftbuf_out(1,1,1,latn), & grlps1(1,irow), grt1(1,1,irow), grz1(1,1,irow), grd1(1,1,irow), & grfu1(1,1,irow), grfv1(1,1,irow), & grut1(1,1,irow), grvt1(1,1,irow), grrh1(1,1,irow), grlps2(1,irow),grt2(1,1,irow), & grz2(1,1,irow), grd2(1,1,irow), grfu2(1,1,irow), grfv2(1,1,irow), grut2(1,1,irow), & grvt2(1,1,irow), grrh2(1,1,irow) ) #else call linemsdyn_aft (irow, nlon_fft_in, fftbuf_in(1,1,1,lats), fftbuf_in(1,1,1,latn), & grlps1(1,irow), grt1(1,1,irow), grz1(1,1,irow), grd1(1,1,irow), & grfu1(1,1,irow), grfv1(1,1,irow), & grut1(1,1,irow), grvt1(1,1,irow), grrh1(1,1,irow), grlps2(1,irow),grt2(1,1,irow), & grz2(1,1,irow), grd2(1,1,irow), grfu2(1,1,irow), grfv2(1,1,irow), grut2(1,1,irow), & grvt2(1,1,irow), grrh2(1,1,irow) ) #endif end do call t_stopf ('linemsdyn_aft') ! call t_startf ('scandyn_dealloc') deallocate(fftbuf_in) deallocate(fftbuf_out) call t_stopf ('scandyn_dealloc') ! call t_startf ('moisture_mass') ! ! Initialize moisture mass integrals. ! hw1lat = 0.0_r8 ! ! Calculate total mass of moisture in fields advected ! #ifdef OUTER_OMP !$OMP PARALLEL DO PRIVATE (LAT, IROW) #endif do lat=beglat,endlat if(lat.le.plat/2) then irow = lat else irow = plat + 1 - lat end if ! ! Only pdel is needed pint and pmid are not. ! call plevs0 (nlon(lat),plon,plev,ps(1,lat,n3m2), pint, pmid, pdel) ! ! Calculate mass of moisture in field being advected ! ! q3 is plon,plev,pcnst,beglat:endlat,ptimelevs ! qminus is plon,plev,pcnst,beglat:endlat call qmassarun (cwava(lat),w(irow) ,qminus(1,1,1,lat),pdel , & hw1lat(1,lat),nlon(lat), q3(1,1,1,lat,n3m2), lat, & pdeld(:,:,lat,n3m2 )) end do call t_stopf ('moisture_mass') return end subroutine scandyn