There is no emphasis for description of the particular main 3D time

stepping procedure used in ROMS v. 1.8....2.2 in the Shchepetkin

McWilliams,2005 ROMS paper because this algorithm was always considered

as a provisional variant to be replaced with a more refined one.

The main 3D procedure of v. 1.8....2.2 can be classified as an AB3-TR

Generalized forward-backward scheme [AB3 = Adams--Bashforth 3rd order

(open parabolic rule); TR = trapezoidal rule) and is summarized as

follows:

**Code:**

u^{n+1} = u^n +dt*[(23/12)*RHSU^n - (4/3)*RHSU^{n-1} + (5/12)*RHSU^{n-2}]

T^{n+1} = T^n -dt*div[ (1/2)*(U^{n+1}+U^n)*T^{n+1/2} ]

where T^{n+1/2} = (5/12)*T^{n+1,*} + (2/3)*T^n - (1/12)*T^{n-1} (the

interpolation coefficients are recognized as the ones for Adams--Moulton

3rd-order closed parabolic integration rule). In its turn T^{n+1,*} comes

from LF predictor sub-step,

**Code:**

T^{n+1,*}=T^{n-1} -2*dt*div[ U^n * T^n ]

which is implemented in practice as a pseudo-compressible finite-volume

step with artificial continuity equation, and, The LF-step is also

combined with AM3 interpolation (this is done in pre_step.F), so that

you actually do not see a naked LF explicitly in the code, but rather

you see

**Code:**

T^{n+1/2} = (5/12)*[ T^{n-1} -2*dt*div(U^n * T^n)]

+ (2/3)*T^n -(1/12)*T^{n-1}

or

**Code:**

T^{n+1/2} = (1/3)*T^{n-1} + (2/3)*T^n -(5/6)*dt*div(U^n * T^n)

where (1/3) appear in "pre_step.F" as "1/2-Gamma"; (2/3) as "1/2+Gamma";

(5/6) as "1-Gamma"; and Gamma is set to 1/6.

Above for simplicity I have omitted multiplication and division by Hz

and metric factors. Lowercase "u" means velocity; uppercase "U" fluxes

(which map onto Huon and Hvom in the code.

Setting of fluxes for corrector sub-step of tracer equations via

Trapezoidal Rule (TR) occurs in "step3d_uv.F", see

**Code:**

Huon(i,j,k)=0.5_r8*(Huon(i,j,k)+u(i,j,k,nnew)*DC(i,k))

where u(i,j,k,nnew) is u^{n+1}, DC(i,k) is new-time step cross-section

(Hz^{n+1}/pn or pm for v-component), and Huon(i,j,k) is previously

computed flux "U^n" (in the actual code Huon(i,j,k) is initially computed

within "set_massflux.F" (in ROMS 2.2; I believe that older codes may have

different file name), and then modified in "step3d_uv.F".

Assuming that the stability limit is governed primarily by the phase speed

of the fastest internal mode (basically the first baroclinic mode), the

details of LF-AM3 stepping for tracer become unimportant because the rate

of change of T is dominated by vertical gradient of T multiplied by

vertical velocity. In this case the whole algorithm maps onto Eq. (2.49)

for Shchepetkin McWilliams, 2005 with beta=5/12, delta=1/2,

gamma=epsilon=0. Its stability limit alpha_max=1.1441551 and its

placement of characteristic roots relatively to unit circle looks similar

to one of Fig. 13, upper left panel (that panel shows AB3-AM3, rather

than AB3-TR version).

The exact figure showing AB3-TR characteristic roots can be seen on page 22 in

http://marine.rutgers.edu/po/Workshops/ ... petkin.pdf
which is presentation on Venice, 2004 workshop.