4. MITgcm Example Experiments¶
The full MITgcm distribution comes with a set of pre-configured
numerical experiments. Some of these example experiments are tests of
individual parts of the model code, but many are fully fledged
numerical simulations. Full tutorials exist for a few of the examples,
and are documented in sections Section 4.2 -
sec_eg_tank
. The other examples follow the same general
structure as the tutorial examples. However, they only include brief
instructions in a text file called {it README}. The examples are
located in subdirectories under the directory texttt{verification}.
Each example is briefly described below.
4.1. Full list of model examples¶
Once you have chosen the example you want to run, you are ready to compile the code.
4.2. Barotropic Gyre MITgcm Example¶
(in directory: verification/tutorial_barotropic_gyre/)
This example experiment demonstrates using the MITgcm to simulate a Barotropic, wind-forced, ocean gyre circulation. The files for this experiment can be found in the verification directory verification/tutorial_barotropic_gyre. The experiment is a numerical rendition of the gyre circulation problem similar to the problems described analytically by Stommel in 1966 [Sto48] and numerically in Holland et. al [HL5a].
In this experiment the model is configured to represent a rectangular enclosed box of fluid, \(1200 \times 1200\) km in lateral extent. The fluid is 5 km deep and is forced by a constant in time zonal wind stress, \(\tau_x\), that varies sinusoidally in the ‘north-south’ direction. Topologically the grid is Cartesian and the coriolis parameter \(f\) is defined according to a mid-latitude beta-plane equation
where \(y\) is the distance along the ‘north-south’ axis of the simulated domain. For this experiment \(f_{0}\) is set to \(10^{-4}s^{-1}\) in (1) and \(\beta = 10^{-11}s^{-1}m^{-1}\).
The sinusoidal wind-stress variations are defined according to
where \(L_{y}\) is the lateral domain extent (1200~km) and \(\tau_0\) is set to \(0.1N m^{-2}\).
Figure 4.1 summarizes the configuration simulated.
4.2.1. Equations Solved¶
The model is configured in hydrostatic form. The implicit free surface form of the pressure equation described in [MHPA97] is employed. A horizontal Laplacian operator \(\nabla_{h}^2\) provides viscous dissipation. The wind-stress momentum input is added to the momentum equation for the ‘zonal flow’, \(u\). Other terms in the model are explicitly switched off for this experiment configuration (see section Section 4.2.3 ), yielding an active set of equations solved in this configuration as follows
where \(u\) and \(v\) and the \(x\) and \(y\) components of the flow vector \(\vec{u}\).
4.2.2. Discrete Numerical Configuration¶
The domain is discretised with a uniform grid spacing in the horizontal set to \(\Delta x=\Delta y=20\) km, so that there are sixty grid cells in the \(x\) and \(y\) directions. Vertically the model is configured with a single layer with depth, \(\Delta z\), of \(5000\) m.
4.2.2.1. Numerical Stability Criteria¶
The Laplacian dissipation coefficient, \(A_{h}\), is set to \(400 m s^{-1}\). This value is chosen to yield a Munk layer width [Adc95],
of \(\approx\) 100km. This is greater than the model resolution \(\Delta x\), ensuring that the frictional boundary layer is well resolved.
The model is stepped forward with a time step \(\delta t=1200\) secs. With this time step the stability parameter to the horizontal Laplacian friction [Adc95]
evaluates to 0.012, which is well below the 0.3 upper limit for stability.
The numerical stability for inertial oscillations [Adc95]
evaluates to \(0.0144\) , which is well below the 0.5 upper limit for stability.
The advective CFL [Adc95] for an extreme maximum horizontal flow speed of \({|\vec{u}|} = 2 ms^{-1}\)
evaluates to 0.12. This is approaching the stability limit of 0.5 and limits \(\delta t\) to 1200 s.
4.2.3. Code Configuration¶
The model configuration for this experiment resides under the directory verification/tutorial_barotropic_gyre/
.
The experiment files
- input/data
- input/data.pkg
- input/eedata
- input/windx.sin_y
- input/topog.box
- code/CPP_EEOPTIONS.h
- code/CPP_OPTIONS.h
- code/SIZE.h
contain the code customizations and parameter settings for this experiments. Below we describe the customizations to these files associated with this experiment.
4.2.3.1. File input/data¶
This file, reproduced completely below, specifies the main parameters for the experiment. The parameters that are significant for this configuration are
Line 7
- viscAh=4.E2,
- this line sets the Laplacian friction coefficient to \(400 m^2s^{-1}\)
Line 10
- beta=1.E-11,
- this line sets \(\beta\) (the gradient of the coriolis parameter, \(f\)) to \(10^{-11} s^{-1}m^{-1}\)
Lines 15 and 16
- rigidLid=.FALSE.,
- implicitFreeSurface=.TRUE.,
- these lines suppress the rigid lid formulation of the surface pressure inverter and activate the implicit free surface form of the pressure inverter.
Line 27
- startTime=0,
- this line indicates that the experiment should start from \(t=0\) and implicitly suppresses searching for checkpoint files associated with restarting an numerical integration from a previously saved state.
Line 29
- endTime=12000,
- this line indicates that the experiment should start finish at \(t=12000s\). A restart file will be written at this time that will enable the simulation to be continued from this point.
Line 30
- deltaTmom=1200,
- This line sets the momentum equation timestep to \(1200s\).
Line 39
- usingCartesianGrid=.TRUE.,
- This line requests that the simulation be performed in a Cartesian coordinate system.
Line 41
- delX=60*20E3,
- This line sets the horizontal grid spacing between each x-coordinate line in the discrete grid. The syntax indicates that the discrete grid should be comprise of $60$ grid lines each separated by \(20 \times 10^{3}m\) (20 km).
Line 42
- delY=60*20E3,
- This line sets the horizontal grid spacing between each y-coordinate line in the discrete grid to \(20 \times 10^{3}m\) (20 km).
Line 43
- delZ=5000,
- This line sets the vertical grid spacing between each z-coordinate line in the discrete grid to 5000m (5 km).
Line 46
- bathyFile=’topog.box’
- This line specifies the name of the file from which the domain bathymetry is read. This file is a two-dimensional (\(x,y\)) map of depths. This file is assumed to contain 64-bit binary numbers giving the depth of the model at each grid cell, ordered with the x coordinate varying fastest. The points are ordered from low coordinate to high coordinate for both axes. The units and orientation of the depths in this file are the same as used in the MITgcm code. In this experiment, a depth of 0 m indicates a solid wall and a depth of -5000 m indicates open ocean. The matlab program input/gendata.m shows an example of how to generate a bathymetry file.
Line 49
- zonalWindFile=’windx.sin_y’
- This line specifies the name of the file from which the x-direction surface wind stress is read. This file is also a two-dimensional (\(x,y\)) map and is enumerated and formatted in the same manner as the bathymetry file. The matlab program input/gendata.m includes example code to generate a valid zonalWindFile file.
other lines in the file input/data are standard values that are described in the MITgcm Getting Started and MITgcm Parameters notes.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 | # Model parameters
# Continuous equation parameters
&PARM01
tRef=20.,
sRef=10.,
viscAz=1.E-2,
viscAh=4.E2,
diffKhT=4.E2,
diffKzT=1.E-2,
beta=1.E-11,
tAlpha=2.E-4,
sBeta =0.,
gravity=9.81,
gBaro=9.81,
rigidLid=.FALSE.,
implicitFreeSurface=.TRUE.,
eosType='LINEAR',
readBinaryPrec=64,
&
# Elliptic solver parameters
&PARM02
cg2dMaxIters=1000,
cg2dTargetResidual=1.E-7,
&
# Time stepping parameters
&PARM03
startTime=0,
#endTime=311040000,
endTime=12000.0,
deltaTmom=1200.0,
deltaTtracer=1200.0,
abEps=0.1,
pChkptFreq=2592000.0,
chkptFreq=120000.0,
dumpFreq=2592000.0,
monitorSelect=2,
monitorFreq=1.,
&
# Gridding parameters
&PARM04
usingCartesianGrid=.TRUE.,
usingSphericalPolarGrid=.FALSE.,
delX=60*20E3,
delY=60*20E3,
delZ=5000.,
&
# Input datasets
&PARM05
bathyFile='topog.box',
hydrogThetaFile=,
hydrogSaltFile=,
zonalWindFile='windx.sin_y',
meridWindFile=,
&
|
4.2.3.2. File input/data.pkg¶
This file uses standard default values and does not contain customizations for this experiment.
4.2.3.3. File input/eedata¶
This file uses standard default values and does not contain customizations for this experiment.
4.2.3.4. File input/windx.sin_y¶
The input/windx.sin_y file specifies a two-dimensional (\(x,y\)) map of wind stress, \(\tau_{x}\), values. The units used are \(Nm^{-2}\). Although \(\tau_{x}\) is only a function of \(y\) in this experiment this file must still define a complete two-dimensional map in order to be compatible with the standard code for loading forcing fields in MITgcm. The included matlab program input/gendata.m gives a complete code for creating the input/windx.sin_y file.
4.2.3.5. File input/topog.box¶
The input/topog.box file specifies a two-dimensional (\(x,y\)) map of depth values. For this experiment values are either 0 m or \(-delZ\) m, corresponding respectively to a wall or to deep ocean. The file contains a raw binary stream of data that is enumerated in the same way as standard MITgcm two-dimensional, horizontal arrays. The included matlab program input/gendata.m gives a completecode for creating the input/topog.box file.
4.2.3.6. File code/SIZE.h¶
Two lines are customized in this file for the current experiment
- Line 39
- sNx=60,
- this line sets the lateral domain extent in grid points for the axis aligned with the x-coordinate.
- Line 40
- sNy=60,
- this line sets the lateral domain extent in grid points for the axis aligned with the y-coordinate.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 | C $Header$
C $Name$
C
C /==========================================================\
C | SIZE.h Declare size of underlying computational grid. |
C |==========================================================|
C | The design here support a three-dimensional model grid |
C | with indices I,J and K. The three-dimensional domain |
C | is comprised of nPx*nSx blocks of size sNx along one axis|
C | nPy*nSy blocks of size sNy along another axis and one |
C | block of size Nz along the final axis. |
C | Blocks have overlap regions of size OLx and OLy along the|
C | dimensions that are subdivided. |
C \==========================================================/
C Voodoo numbers controlling data layout.
C sNx - No. X points in sub-grid.
C sNy - No. Y points in sub-grid.
C OLx - Overlap extent in X.
C OLy - Overlat extent in Y.
C nSx - No. sub-grids in X.
C nSy - No. sub-grids in Y.
C nPx - No. of processes to use in X.
C nPy - No. of processes to use in Y.
C Nx - No. points in X for the total domain.
C Ny - No. points in Y for the total domain.
C Nr - No. points in R for full process domain.
INTEGER sNx
INTEGER sNy
INTEGER OLx
INTEGER OLy
INTEGER nSx
INTEGER nSy
INTEGER nPx
INTEGER nPy
INTEGER Nx
INTEGER Ny
INTEGER Nr
PARAMETER (
& sNx = 30,
& sNy = 30,
& OLx = 2,
& OLy = 2,
& nSx = 2,
& nSy = 2,
& nPx = 1,
& nPy = 1,
& Nx = sNx*nSx*nPx,
& Ny = sNy*nSy*nPy,
& Nr = 1)
C MAX_OLX - Set to the maximum overlap region size of any array
C MAX_OLY that will be exchanged. Controls the sizing of exch
C routine buufers.
INTEGER MAX_OLX
INTEGER MAX_OLY
PARAMETER ( MAX_OLX = OLx,
& MAX_OLY = OLy )
|
4.2.3.7. File code/CPP_OPTIONS.h¶
This file uses standard default values and does not contain customizations for this experiment.
4.2.3.8. File code/CPP_EEOPTIONS.h¶
This file uses standard default values and does not contain customizations for this experiment.
4.3. A Rotating Tank in Cylindrical Coordinates¶
(in directory:verification/rotating_tank/
)
4.3.1. Overview¶
This example configuration demonstrates using the MITgcm to simulate a
laboratory demonstration using a differentially heated rotating
annulus of water. The simulation is configured for a laboratory scale
on a \(3^{\circ}\times1\mathrm{cm}\) cyclindrical grid with twenty-nine
vertical levels of 0.5cm each. This is a typical laboratory setup for
illustration principles of GFD, as well as for a laboratory data
assimilation project. The files for this experiment can be found in
the verification directory under rotating_tank
.
example illustration from GFD lab here
4.3.2. Equations Solved¶
4.3.3. Discrete Numerical Configuration¶
The domain is discretised with a uniform cylindrical grid spacing in the horizontal set to \(\Delta a=1`~cm and :math:\)Delta phi=3^{circ}`, so that there are 120 grid cells in the azimuthal direction and thirty-one grid cells in the radial, representing a tank 62cm in diameter. The bathymetry file sets the depth=0 in the nine lowest radial rows to represent the central of the annulus. Vertically the model is configured with twenty-nine layers of uniform 0.5cm thickness.
something about heat flux
4.3.4. Code Configuration¶
The model configuration for this experiment resides under the
directory verification/rotatingi_tank/
. The experiment files
input/data
input/data.pkg
input/eedata
input/bathyPol.bin
input/thetaPol.bin
code/CPP\_EEOPTIONS.h
code/CPP\_OPTIONS.h
code/SIZE.h
contain the code customizations and parameter settings for this experiments. Below we describe the customizations to these files associated with this experiment.
4.3.4.1. File input/data¶
This file, reproduced completely below, specifies the main parameters for the experiment. The parameters that are significant for this configuration are
- Lines 9-10,
- viscAh=5.0E-6,
- viscAz=5.0E-6,
These lines set the Laplacian friction coefficient in the horizontal and vertical, respectively. Note that they are several orders of magnitude smaller than the other examples due to the small scale of this example.
- Lines 13-16,
- diffKhT=2.5E-6,
- diffKzT=2.5E-6,
- diffKhS=1.0E-6,
- diffKzS=1.0E-6,
These lines set horizontal and vertical diffusion coefficients for temperature and salinity. Similarly to the friction coefficients, the values are a couple of orders of magnitude less than most
configurations.
- Line 17, f0=0.5, this line sets the
coriolis term, and represents a tank spinning at about 2.4 rpm.
- Lines 23 and 24
- rigidLid=.TRUE.,
- implicitFreeSurface=.FALSE.,
These lines activate the rigid lid formulation of the surface pressure inverter and suppress the implicit free surface form of the pressure inverter.
- Line 40,
- nIter=0,
This line indicates that the experiment should start from $t=0$ and implicitly suppresses searching for checkpoint files associated with restarting an numerical integration from a previously saved state. Instead, the file thetaPol.bin will be loaded to initialized the temperature fields as indicated below, and other variables will be initialized to their defaults.
- Line 43,
- deltaT=0.1,
This line sets the integration timestep to $0.1s$. This is an unsually small value among the examples due to the small physical scale of the experiment. Using the ensemble Kalman filter to produce input fields can necessitate even shorter timesteps.
- Line 56,
- usingCylindricalGrid=.TRUE.,
This line requests that the simulation be performed in a cylindrical coordinate system.
- Line 57,
- dXspacing=3,
This line sets the azimuthal grid spacing between each $x$-coordinate line in the discrete grid. The syntax indicates that the discrete grid should be comprised of $120$ grid lines each separated by $3^{circ}$.
- Line 58,
- dYspacing=0.01,
This line sets the radial cylindrical grid spacing between each \(a\)-coordinate line in the discrete grid to \(1cm\).
- Line 59,
- delZ=29*0.005,
This line sets the vertical grid spacing between each of 29 z-coordinate lines in the discrete grid to $0.005m$ ($5$~mm).
- Line 64,
- bathyFile=’bathyPol.bin’,
This line specifies the name of the file from which the domain ‘bathymetry’ (tank depth) is read. This file is a two-dimensional (\(a,\phi\)) map of depths. This file is assumed to contain 64-bit binary numbers giving the depth of the model at each grid cell, ordered with the $phi$ coordinate varying fastest. The points are ordered from low coordinate to high coordinate for both axes. The units and orientation of the depths in this file are the same as used in the MITgcm code. In this experiment, a depth of $0m$ indicates an area outside of the tank and a depth f \(-0.145m\) indicates the tank itself.
- Line 65,
- hydrogThetaFile=’thetaPol.bin’,
This line specifies the name of the file from which the initial values of temperature are read. This file is a three-dimensional (\(x,y,z\)) map and is enumerated and formatted in the same manner as the bathymetry file.
- Lines 66 and 67
- tCylIn = 0
- tCylOut = 20
These line specify the temperatures in degrees Celsius of the interior and exterior walls of the tank – typically taken to be icewater on the inside and room temperature on the outside.
Other lines in the file input/data are standard values that are described in the MITgcm Getting Started and MITgcm Parameters notes.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60 61 62 63 64 65 66 67 | # ====================
# | Model parameters |
# ====================
#
# Continuous equation parameters
&PARM01
tRef=29*20.0,
sRef=29*35.0,
viscAh=5.0E-6,
viscAz=5.0E-6,
no_slip_sides=.FALSE.,
no_slip_bottom=.FALSE.,
diffKhT=2.5E-6,
diffKzT=2.5E-6,
diffKhS=1.0E-6,
diffKzS=1.0E-6,
f0=0.5,
eosType='LINEAR',
sBeta =0.,
gravity=9.81,
rhoConst=1000.0,
rhoNil=1000.0,
#heatCapacity_Cp=3900.0,
rigidLid=.TRUE.,
implicitFreeSurface=.FALSE.,
nonHydrostatic=.TRUE.,
readBinaryPrec=32,
&
# Elliptic solver parameters
&PARM02
cg2dMaxIters=1000,
cg2dTargetResidual=1.E-7,
cg3dMaxIters=10,
cg3dTargetResidual=1.E-9,
&
# Time stepping parameters
&PARM03
nIter0=0,
nTimeSteps=20,
#nTimeSteps=36000000,
deltaT=0.1,
abEps=0.1,
pChkptFreq=2.0,
#chkptFreq=2.0,
dumpFreq=2.0,
monitorSelect=2,
monitorFreq=0.1,
&
# Gridding parameters
&PARM04
usingCylindricalGrid=.TRUE.,
dXspacing=3.,
dYspacing=0.01,
delZ=29*0.005,
ygOrigin=0.07,
&
# Input datasets
&PARM05
hydrogThetaFile='thetaPolR.bin',
bathyFile='bathyPolR.bin',
tCylIn = 0.,
tCylOut = 20.,
&
|
4.3.4.2. File input/data.pkg¶
This file uses standard default values and does not contain customizations for this experiment.
4.3.4.3. File input/eedata¶
This file uses standard default values and does not contain customizations for this experiment.
4.3.4.4. File input/thetaPol.bin¶
The {it input/thetaPol.bin} file specifies a three-dimensional ($x,y,z$) map of initial values of $theta$ in degrees Celsius. This particular experiment is set to random values x around 20C to provide initial perturbations.
4.3.4.5. File input/bathyPol.bin¶
The {it input/bathyPol.bin} file specifies a two-dimensional ($x,y$) map of depth values. For this experiment values are either $0m$ or {bf -delZ}m, corresponding respectively to outside or inside of the tank. The file contains a raw binary stream of data that is enumerated in the same way as standard MITgcm two-dimensional, horizontal arrays.
4.3.4.6. File code/SIZE.h¶
Two lines are customized in this file for the current experiment
- Line 39, - sNx=120,
this line sets the lateral domain extent in grid points for the axis aligned with the x-coordinate.
- Line 40, - sNy=31,
this line sets the lateral domain extent in grid points for the axis aligned with the y-coordinate.
1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 | C $Header$
C $Name$
C
C /==========================================================\
C | SIZE.h Declare size of underlying computational grid. |
C |==========================================================|
C | The design here support a three-dimensional model grid |
C | with indices I,J and K. The three-dimensional domain |
C | is comprised of nPx*nSx blocks of size sNx along one axis|
C | nPy*nSy blocks of size sNy along another axis and one |
C | block of size Nz along the final axis. |
C | Blocks have overlap regions of size OLx and OLy along the|
C | dimensions that are subdivided. |
C \==========================================================/
C Voodoo numbers controlling data layout.
C sNx - No. X points in sub-grid.
C sNy - No. Y points in sub-grid.
C OLx - Overlap extent in X.
C OLy - Overlat extent in Y.
C nSx - No. sub-grids in X.
C nSy - No. sub-grids in Y.
C nPx - No. of processes to use in X.
C nPy - No. of processes to use in Y.
C Nx - No. points in X for the total domain.
C Ny - No. points in Y for the total domain.
C Nr - No. points in Z for full process domain.
INTEGER sNx
INTEGER sNy
INTEGER OLx
INTEGER OLy
INTEGER nSx
INTEGER nSy
INTEGER nPx
INTEGER nPy
INTEGER Nx
INTEGER Ny
INTEGER Nr
PARAMETER (
& sNx = 30,
& sNy = 23,
& OLx = 3,
& OLy = 3,
& nSx = 4,
& nSy = 1,
& nPx = 1,
& nPy = 1,
& Nx = sNx*nSx*nPx,
& Ny = sNy*nSy*nPy,
& Nr = 29)
C MAX_OLX - Set to the maximum overlap region size of any array
C MAX_OLY that will be exchanged. Controls the sizing of exch
C routine buufers.
INTEGER MAX_OLX
INTEGER MAX_OLY
PARAMETER ( MAX_OLX = OLx,
& MAX_OLY = OLy )
|
4.3.4.7. File code/CPP_OPTIONS.h¶
This file uses standard default values and does not contain customizations for this experiment.
4.3.4.8. File code/CPP_EEOPTIONS.h¶
This file uses standard default values and does not contain customizations for this experiment.