11 



3t 



1 3 /UV 



y „ 3x \ H 



1 3 hi 



Py 9y V H 





-1-7^(^1 + -7^1^)+ H(vu.)„=o 1+ I^U .51111+ J_(,^„.,,J 



3p ^ 3H 

 -^ da + g -— 

 3y ^ 3y 



Py 3y P( 



sy'^by' 



(/ 



pdo+op 



Hdc + (H.D.), 



:i4) 



where w^^Edx/da and pyidy/dy are the stretching coefficients, a)=da/dt is 

 the vertical velocity in the stretched coordinate, and (H.D.)^ and 

 (H.D.) are the horizontal diffusion terms. Notice that the bottom 

 stresses (t^j^, j^^) are determined from the latest three-dimensional 

 velocity profiles available from the internal mode computation, and 

 hence are more accurate than the traditional vertically-integrated 

 models which assume the bottom stress is proportional to the local 

 vertically-integrated velocity or its square. 



Treating implicitly all the terms in the continuity equation, 

 while only the time derivatives and the surface slopes in the momentum 

 equations, one can obtain the following finite-difference equations: 



n+1 



(15) 



where 



= [r+(r-4.)Xj^+(r-4,)x ]w + Ato 



X ,, AY X' 



p^AX 





D = D 



16) 



where (^x.Ay) are the spatial grids, ^t is time step, D^ and Dy are 

 terms in Eqs. (13) and (14) excluding the time derivatives ana the 

 surface slopes, superscripts n+1 and n indicate present and previous 

 time step of integration, 6^, and 6 are central difference spatial 

 operators, and 4) is a weighting factor, 0<({]<1. If <P=0, Eq. (15) 

 reduces to a two-step explicit scheme. If 4)>0 the resulting schemes 

 are implicit, with (})=l/2 corresponding to the Crank-Nicholson scheme 

 and 0=1 corresponding to the fully implicit scheme. Eq. (15) can be 

 factorized such that solution can be obtained by consecutive 

 tridiagonal matrix inversions in the x-direction ,and y-.direction. 

 Further, we employ a method that solves only two variables during each 

 sweep. This method allows very large time step to be used and has been 

 found to be more stable than the traditional ADI method. Courant 

 number (^%sed on the maximum propagation speed of surface gravity wave, 

 (gH^gj^) • At/Ax, may now be as large as 100, compared to the limit of 

 1 for the explicit method. The maximum step is now governed by the CFL 

 condition based on vertically-averaged advection speed in the system. 



In the full three-dimensional model, the external mode computation 



272 



Sheng 



Am 



