#include <misc.h>
#include <params.h>
subroutine dyn(irow ,grlps1 ,grt1 ,grz1 ,grd1 , & 1,9
grfu1 ,grfv1 ,grut1 ,grvt1 ,grrh1 , &
grlps2 ,grt2 ,grz2 ,grd2 ,grfu2 , &
grfv2 ,grut2 ,grvt2 ,grrh2 )
!-----------------------------------------------------------------------
!
! Combine undifferentiated and longitudinally differentiated Fourier
! coefficient terms for later use in the Gaussian quadrature
!
! Computational note: Index "2*m-1" refers to the real part of the
! complex coefficient, and "2*m" to the imaginary.
!
! The naming convention is as follows:
! - t, q, d, z refer to temperature, specific humidity, divergence
! and vorticity
! - "1" suffix to an array => symmetric component of current latitude pair
! - "2" suffix to an array => antisymmetric component
!
!---------------------------Code history--------------------------------
!
! Original version: J. Rosinski
! Standardized: J. Rosinski, June 1992
! Reviewed: D. Williamson, B. Boville, J. Hack, August 1992
! Reviewed: D. Williamson, March 1996
! Modified: P. Worley, September 2002
!
!-----------------------------------------------------------------------
!
! $Id: dyn.F90,v 1.4.2.5 2004/04/28 23:43:54 eaton Exp $
! $Author: eaton $
!
!-----------------------------------------------------------------------
use shr_kind_mod
, only: r8 => shr_kind_r8
use pmgrid
use pspect
use rgrid
use comspe
use commap
use dynconst
, only: rearth
use time_manager
, only: get_step_size, is_first_step
implicit none
!
! Input arguments
!
integer irow ! latitude pair index
!
! Input/output arguments
!
real(r8) grlps1(2*maxm) ! sym. surface pressure equation term
real(r8) grt1(2*maxm,plev) ! sym. undifferentiated term in t eqn.
real(r8) grz1(2*maxm,plev) ! sym. undifferentiated term in z eqn.
real(r8) grd1(2*maxm,plev) ! sym. undifferentiated term in d eqn.
real(r8) grfu1(2*maxm,plev) ! sym. nonlinear terms in u eqn.
real(r8) grfv1(2*maxm,plev) ! sym. nonlinear terms in v eqn.
real(r8) grut1(2*maxm,plev) ! sym. lambda derivative term in t eqn.
real(r8) grvt1(2*maxm,plev) ! sym. mu derivative term in t eqn.
real(r8) grrh1(2*maxm,plev) ! sym. RHS of divergence eqn (del^2 term)
real(r8) grlps2(2*maxm) ! antisym. surface pressure equation term
real(r8) grt2(2*maxm,plev) ! antisym. undifferentiated term in t eqn.
real(r8) grz2(2*maxm,plev) ! antisym. undifferentiated term in z eqn.
real(r8) grd2(2*maxm,plev) ! antisym. undifferentiated term in d eqn.
real(r8) grfu2(2*maxm,plev) ! antisym. nonlinear terms in u eqn.
real(r8) grfv2(2*maxm,plev) ! antisym. nonlinear terms in v eqn.
real(r8) grut2(2*maxm,plev) ! antisym. lambda derivative term in t eqn.
real(r8) grvt2(2*maxm,plev) ! antisym. mu derivative term in t eqn.
real(r8) grrh2(2*maxm,plev) ! antisym. RHS of divergence eqn (del^2 term)
!
!---------------------------Local workspace-----------------------------
!
real(r8) tmp1,tmp2 ! temporaries
real(r8) zxm(pmmax) ! m*2dt/(a*cos(lat)**2)
real(r8) zrcsj ! 1./(a*cos(lat)**2)
real(r8) dtime ! timestep size [seconds]
real(r8) ztdtrc ! 2dt/(a*cos(lat)**2) 1dt/..... at nstep=0
integer lm, mlength ! local Fourier wavenumber index
! and number of local indices
integer k ! level index
!DIR$ NOSTREAM
!
! Set constants
!
mlength = numm(iam)
dtime = get_step_size
()
zrcsj = 1./(cs(irow)*rearth)
if (is_first_step()) then
ztdtrc = dtime*zrcsj
else
ztdtrc = 2.0*dtime*zrcsj
end if
!
! Combine constants with Fourier wavenumber m
!
!cdir altcode=(loopcnt)
do lm=1,mlength
zxm(lm) = ztdtrc*xm(locm(lm,iam))
end do
!
! Combine undifferentiated and longitudinal derivative terms for
! later use in Gaussian quadrature
!
do k=1,plev
!cdir loopchg
do lm=1,mlength
grt1(2*lm-1,k) = grt1(2*lm-1,k) + zxm(lm)*grut1(2*lm,k)
grt1(2*lm,k) = grt1(2*lm,k) - zxm(lm)*grut1(2*lm-1,k)
grd1(2*lm-1,k) = grd1(2*lm-1,k) - zxm(lm)*grfu1(2*lm,k)
grd1(2*lm,k) = grd1(2*lm,k) + zxm(lm)*grfu1(2*lm-1,k)
grz1(2*lm-1,k) = grz1(2*lm-1,k) - zxm(lm)*grfv1(2*lm,k)
grz1(2*lm,k) = grz1(2*lm,k) + zxm(lm)*grfv1(2*lm-1,k)
!
grt2(2*lm-1,k) = grt2(2*lm-1,k) + zxm(lm)*grut2(2*lm,k)
grt2(2*lm,k) = grt2(2*lm,k) - zxm(lm)*grut2(2*lm-1,k)
grd2(2*lm-1,k) = grd2(2*lm-1,k) - zxm(lm)*grfu2(2*lm,k)
grd2(2*lm,k) = grd2(2*lm,k) + zxm(lm)*grfu2(2*lm-1,k)
grz2(2*lm-1,k) = grz2(2*lm-1,k) - zxm(lm)*grfv2(2*lm,k)
grz2(2*lm,k) = grz2(2*lm,k) + zxm(lm)*grfv2(2*lm-1,k)
end do
end do
return
end subroutine dyn