c c a test program. The comment statements contained therein should be c pretty self-explanatory. The WVPARAM subroutine is a temporary makeshift c one valid for December solstice. It provides eastward acceleration in c m/s**2 that should be utilized as a momentum source term in the TGCM c to simulate the effects of wave dissipation (instead of your current c Rayleigh friction). This approach should be viewed for the moment as c exploratory. Let me know what happens. Jeff c c SUBROUTINE WVPARAM (Z, XLAT, VB, FUB) c This version of WVPARAM is a temporary one, valid only for the month of c December, and intended only for interim numerical experimentation. It is c based on table look-up in the meridional wind model of Yu.I. Portnyagin c (YP), which is based on meteor and SAD radar data (see Miyahara et al., c JGR, 96, p. 1226-1227,1991, for a description). The YP table extends from 70 c to 110 km. An exponential fix is introduced down to 60 km and up to 120 km c to avoid discontinuities. This subroutine will be replaced by an analytic c one based on Legendre Polynomial decomposition of the full YP tables for c 12 months of the year. c Input: z = altitude in km c xlat = latitude in degrees c Output: c vb = prevailing northward wind (m/s) c fub = eastward acceleration (m/s2) c = -f*vb c where f = coriolis parameter c Caution: The assumption here is that the YP model represents true zonal c mean winds. The possibility exists that some correction needs c to be introduced to correct the YP winds for geostrophic c meridional winds, say due to wavenumber 1 or 2 stationary c waves. This is under examination. c c J. Forbes, November, 1991 subroutine wvparam(z,xlat,vb,fub) dimension iv(37,21) c Tabulated prevailing northward winds from Portnyagin model c c lat = - 90 ht = 70,72,................110 c - 85 c . c . c . c + 90 data ((iv(l,iz),iz=1,21),l=1,10)/ + 1,1,1,1,1,2,2,2,2,2,1,1,0,0,0,0,0,1,1,1,1, + 1,1,1,1,1,2,2,2,2,2,1,1,0,0,0,0,0,1,1,1,1, + 3,3,3,3,3,3,3,3,2,2,2,2,2,2,2,2,2,2,1,1,1, + 4,4,4,4,4,4,4,3,3,3,3,3,3,3,3,3,3,2,2,2,1, + 3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,3,2,2,2,2, + 1,1,1,1,0,0,1,1,1,2,2,2,2,2,2,2,2,2,2,2,2, + -2,-2,-2,-2,-3,-3,-3,-2,-1,-1,0,0,0,0,1,1,1,2,2,2,3, + -4,-5,-5,-5,-5,-5,-5,-4,-4,-3,-2,-1,-1,-1,0,1,1,2,2,2,3, + -6,-6,-6,-6,-6,-6,-6,-5,-3,-4,-3,-3,-2,-1,0,1,1,2,2,2,3, + -6,-6,-6,-6,-5,-5,-5,-4,-4,-3,-3,-2,-1,0,1,1,2,2,2,2,2/ data ((iv(l,iz),iz=1,21),l=11,20)/ + -5,-4,-4,-3,-3,-2,-2,-1,-1,-1,-1,0,0,1,1,2,2,2,2,2,2, + -2,-2,-1,0,1,1,2,2,2,3,3,2,2,2,2,2,2,2,2,2,2, + 1,1,2,3,3,4,5,5,6,6,6,5,4,3,2,2,1,1,1,1,1, + 2,3,3,4,5,5,6,7,7,8,8,7,5,4,2,1,1,0,1,1,1, + 3,3,4,4,4,5,5,6,7,7,7,7,6,4,2,1,0,0,0,0,1, + 3,3,3,3,3,3,3,4,5,5,6,6,5,4,2,0,-1,-1,-1,0,1, + 2,2,2,2,1,1,1,1,2,3,3,4,3,2,1,-1,-2,-2,-2,0,1, + 3,2,2,2,1,1,1,0,1,1,1,2,2,1,-1,-2,-3,-3,-2,0,1, + 4,4,3,3,3,3,2,2,1,1,0,0,0,-2,-3,-4,-5,-4,-3,0,1, + 5,6,6,6,6,6,6,5,4,3,2,0,-1,-4,-5,-6,-6,-5,-2,0,2/ data ((iv(l,iz),iz=1,21),l=21,30)/ + 7,8,9,10,10,10,10,9,8,6,4,2,-2,-5,-7,-8,-7,-5,-2,1,2, + 9,10,11,12,13,13,13,12,11,9,7,4,-1,-5,-8,-9,-7,-4,-1,2,2, + 8,9,11,12,13,14,14,14,13,11,9,5,0,-5,-7,-8,-6,-3,0,2,3, + 5,6,8,9,10,11,12,12,12,11,10,6,1,-3,-5,-6,-4,-1,2,3,3, + 1,1,2,3,5,6,7,8,8,8,8,6,3,0,-1,-2,-1,1,3,4,4, + -5,-5,-4,-3,-2,0,1,2,4,5,6,5,4,4,3,3,3,4,4,4,5, + -9,-10,-9,-9,-8,-6,-4,-3,-1,1,3,4,6,6,7,6,6,6,5,5,5, + -13,-13,-13,-12,-11,-9,-8,-6,-3,1,1,4,6,8,9,9,8,7,6,5,5, + -14,-14,-14,-13,-12,-10,-8,-6,-4,-1,1,4,7,9,9,9,8,7,6,5,5, + -13,-13,-12,-12,-10,-8,-6,-4,-2,0,2,4,6,8,8,8,7,6,5,5,5/ data ((iv(l,iz),iz=1,21),l=31,37)/ + -11,-11,-10,-9,-7,-5,-4,-2,1,2,4,5,6,6,6,6,5,5,5,4,4, + -9,-8,-6,-5,-4,-2,0,1,3,4,5,5,5,4,4,4,3,3,3,3,3, + -6,-5,-4,-3,-1,0,2,3,4,5,5,5,4,3,3,2,2,2,2,2,2, + -5,-4,-2,-1,0,2,3,4,5,5,5,4,3,3,2,2,2,1,1,1,1, + -4,-3,-1,0,1,2,3,4,4,4,4,4,3,3,2,2,2,1,1,1,1, + -4,-2,1,1,2,3,4,4,5,5,4,4,4,4,3,3,2,2,2,1,1, + -4,-2,1,1,2,3,4,4,5,5,4,4,4,4,3,3,2,2,2,1,1/ vb = 0.0 fub = 0.0 if(z.lt.60.or.z.gt.120.) go to 99 fact = 1.0 if(z.le.70.) fact = exp(-((z-70.)/5.)**2) if(z.ge.110.) fact = exp(-((z-110.)/5.)**2) f = 2.*.72722e-05*sin(xlat*.0174533) c Since TIGCM is set up on 5 degree grid, forget about interpolating c in latitude: l = ((xlat+90.001)/5.) + 1 c Interpolate in height: jz = ((z-70.)/2. + .00001) dz = (z-70.-float(jz)*2.)/2. if(z.le.70.) jz=0 if(z.ge.108.) jz=19 vb = float(iv(l,jz+1))*(1.-dz) + dz*float(iv(l,jz+2)) vb = vb*fact fub = -f*vb 99 return end