model code for linear one-and-half layer ocean

[vinu]

A model code for a simple 1-1/2-layer hydrodynamic
ocean is given (attached file- rgmodel.f)

This code may be useful for experimental runs for understanding basic
Ocean dynamics, waves in the equatorial betaplane, Sverdrup
circulation, planetary wave propagation etc. This code may not be
suitable for any particular research purpose.
(written in Fortran-77)

Model

The linear 1-1/2-layer reduced gravity model (http://www.oceanographers.net/forums/showthread.php?p=237) with rigid-lid approximation,

dU/dt - f.V = -g'.dh/dx + tx/(rho.H) + Ah.d^2U/dx^2 + Ah.d^2U/dy^2
dV/dt + f.U = -g'.dh/dy + ty/(rho.H) + Ah.d^2V/dx^2 + Ah.d^2V/dy^2
dh/dt = -H(dU/dx + dV/dy) + Ah.d^2h/dx^2 + Ah.d^2h/dy^2

where U,V are velocities of the active layer.
H is the un-disturbed (initial) layer thickness of the active layer
and 'h' is the time-dependent layer thickness. g' is the reduced gravity.
tx and ty are wind stress vectors. The term containing Ah is
the viscosity and diffusion for momentum and layer thickness
respectively.

A simple A-grid (http://www.oceanographers.net/forums/showthread.php?p=237) with beta-plane implimentation is adopted.
Initial condition (wind or layer thickness) can be specified inside the
model. A leapfrog time step (centered in time) with centered in space
is used. A box-shape boundary with no-slip boundary condition
is implimented.

to compile and run the model (under linux/unix)

f77 rgmodel.f
./a.out

(Also compatible with windows fortran compiler).

The output is written in binary format (out.dat)
A grads control file is also attached (rgmodel.ctl) will be useful to
make plots.



-Thanks
------------------------------------------------------------------------------------------

[Malte]

Dear Vinu,

first at all thanks for the code and the description in the other thread. I understand most of it but have some problems understanding the discretization of the differential equations. To specify: I don´t understand clearly how the equations:

dU/dt - f.V = -g'.dh/dx + tx/(rho.H) + Ah.d^2U/dx^2 + Ah.d^2U/dy^2
dV/dt + f.U = -g'.dh/dy + ty/(rho.H) + Ah.d^2V/dx^2 + Ah.d^2V/dy^2
dh/dt = -H(dU/dx + dV/dy) + Ah.d^2h/dx^2 + Ah.d^2h/dy^2

translate into the fortran code:

u(i,j,lp) = u(i,j,lm) + (
& + f(j)*v(i,j,l)
& - rg*(h(ip,j,l) - h(im,j,l))/d2x
& + tau(i,j)/hinitial/rho
& + ah*( (u(ip,j,l)+u(im,j,l)-2.0*u(i,j,l))/dx2 +
& (u(i,jp,l)+u(i,jm,l)-2.0*u(i,j,l))/dy2)
& )*delt
v(i,j,lp) = v(i,j,lm) + (
& - f(j)*u(i,j,l)
& - rg*(h(i,jp,l) - h(i,jm,l))/d2y
& + ah*( (v(ip,j,l)+v(im,j,l)-2.0*v(i,j,l))/dx2 +
& (v(i,jp,l)+v(i,jm,l)-2.0*v(i,j,l))/dy2)
& )*delt

h(i,j,lp) = h(i,j,lm) + (
& -hinitial*((u(ip,j,l)-u(im,j,l))/d2x
& +(v(i,jp,l)-v(i,jm,l))/d2y)
& + ah*( (h(ip,j,l)+h(im,j,l)-2.0*h(i,j,l))/dx2 +
& (h(i,jp,l)+h(i,jm,l)-2.0*h(i,j,l))/dy2)
& )*delt

I have only used box models so far and the grid aspect is quite new to me.
As an exercise I tried to modify your code for the use of the Stommel equations but have problems with the discretization.
Do you have some suggestions (links or pdf sources) which explain discretization and especially the translation into fortran code?

Thank you very much in advance

[vinu]

Hello Malte,

An easy way of capturing the method of descritization
would be, to draw a grid box in a paper and assign the
grid edges with values i,ip and im, such that
ip = i + 1
im = i - 1.

Now the same way, assign for the jp, j and jm.

In the given example, an A-grid is assumed. That means
the location of U,V and eta is at the same location.
So you have to draw four grid squar boxes to locate
ip, i, im, j, jp, jm in the cartesian coordinate.
In this four-grid squar boxes, the location
of i and j is at the center.

Then a first order
differential (i.e. the presure gradient term) at 'i' is descritized
as, (ip - im)/2dx in the zonal direction and (jp-jm)/2dy in
the meridional direction. Remember that the ip, i and im
refers to the position only. The actual values at the position
'i' is h(i,j) and ip is h(ip,j) in the 2-D fortan array.

The finite difference form of
diffusion (second order partial differential is)
d2/dx2 = (ip + im -2.i)/dx2.

Referring to the stommels equation, (second derivatives invloved)
the boundary values (east and west) are important. I think we
are free to put psi (stream function) = 0 and these boundaries
because, in a gyre, water flows northward will return to south and no-normal
flow at east and west boundary.


I will post if I find proper finite difference form of Stommel's equation.


Thanks
-Vinu

[vinu]

Hello Malte,

I could find the Stommel's problem and its method of solution
in Book of Numerical models of ocean and oceanic processes by
Kantha and Clayson (2000), International Geophysical series-66.
Section 2.12.4

Since it involves the laplacian equation the numerical solution
is found by matrix form. The same book section 2.4 describes
the methods in details. I havent tried myself solving it
and I cant be much help.

-Vinu

[Malte]

Hello Vinu,

first of all thank you very much for the answer! I will try to do the discretization again within the next days with your tips.
Unfortunately our universtiy library hasn´t the book you recommended in stock.

Cheers,
Malte

PS: Can you recommend the book (Numerical models of ocean and oceanic processes by
Kantha and Clayson (2000)) as an introduction to the field for someone who doesn`t have a strong background in physics and phsical oceanography, or are there better introductional books for the topic?

[vinu]

Hi Malte,

For the modelling, yes. I recommend this book for some one
"who wants to learn/do some
fundamental numerical modeling in physcal oceanography" although
havent learned much about this subject a prior. This book is unique in its
presentation by starting with the fundamentals of ocean dynamics (in brief)
and then proceeding to the modeling (in much detail). A sister book of this
by the same authors, "Small scale processes in Geophyscial fluid dynamics"
published an year after as Vol-67 encapsulates the small scale dynamic
modeling such as mixing, turbulance etc.. in detail. Thus combining these
two volumes, it becomes a wealth of information for an ocean modeller.
Essentially these two books referes to those who want to learn modelling.

For those who want to learn Physical Oceanography as a beginner, the
above books may not be appropriate, but should follow books such as,

Descritptive Physical Oceanography by Pickard and Emery
Regional Oceanography by Tomczak and Godfrey
Introductory Dynamical Oceanography by Pond and Pickard
Atmopshere-Ocean Dynamics by A. E. Gill
and so on.


-Vinu

[cd264]

I supply another version writen with matlab for people use matlab more often than fortran.The version is found from the website http://ic.ucsc.edu
It uses C-grid and can produce a realistic tropical Pacific Ocean circulation.more details see the Attached Files.
In addition,i would like to say thank you to Vinu,he helps me a lot!

[cd264]

The NLOM model is a good application of the rgmodel.It contains the muti-layers of the rgmodel.
But i don't know whether the model is still downloadable,can anyone tell me that?
A users guide of the model is available,named< The Navy Layered Ocean Model Users Guide.>.