2. Basic Features

The core of gcmfaces lies in (1) its representation of connected arrays, or faces, as objects of the class defined in the @gcmfaces/ directory (see Section 2.1) and (2) its handling of C-Grid specifications via the mygrid global structure (Section 2.2). Other basic features include functions that exchange data between faces (Section 2.3), or overload common operations (Section 2.4), as well as I/O routines (Section 2.5). These and other gcmfaces functions are generally documented through help sections that are readily accessible via the Matlab or Octave command window.

2.1. The gcmfaces Class

Figure 1 illustrates four types of grids that are have been used in ocean general circulation models. Despite evident design differences, these grids can all be represented as sets of connected arrays, or faces, as shown in Figure 2.1 in the case of the LLC90 grid. gcmfaces simply takes advantage of this fact by defining a class for these objects, within @gcmfaces/, that represents gridded earth variables generically as sets of connected arrays.

Grid specifics, such as the number of faces and the size of each face, are embedded within the gcmfaces objects (see Gcmfaces object structure). Objects of the gcmfaces class can thus be manipulated simply through compact and generic expressions such as a+b that are robust to changes in grid design (Figure 1). The gcmfaces class inherits many of its basic operations from the double class as illustrated for the Overloaded + function (see Section 2.4 for details).

Ocean topography on the LLC90 grid displayed face by face

Figure 2.1 Ocean topography on the LLC90 grid (Figure 1, bottom right) displayed face by face (going from 1 to 5). This plot generated using example_display(1) illustrates how gcmfaces organizes data in memory (see Gcmfaces object structure). Within each face, grid point indices increase from left to right and bottom to top.

Gcmfaces object structure

An Earth variable on the LLC90 grid (Figure 1, bottom right) stored as a gcmfaces object called fld has the data structure depicted below. In this example, fld is a two dimensional field, and the five face arrays plotted in Figure 2.1 are denoted as f1 to f5.

  nFaces: 5
  f1: [90x270 double]
  f2: [90x270 double]
  f3: [90x90 double]
  f4: [270x90 double]
  f5: [270x90 double]

2.2. Handling C-Grids

In practice gcmfaces gets activated by adding, to the least, the @gcmfaces/ directory to the Matlab path and then loading a grid to memory (Section 1.3). The default grid is LLC90, which can be loaded to memory by calling grid_load.m without any argument. Section 2.5 and help grid_load; provide additional information regarding, respectively, and supported file formats and grid_load.m arguments. As an alternative to grid_load.m, MITgcm input grid files can be read grid_load_native.m as shown here (see README and demo_grids.m).

Both grid_load.m and grid_load_native.m store all C-grid variables at once in a global variable named mygrid (Table 2.1). gcmfaces functions then rely on mygrid that they get access to by calling gcmfaces_global.m which also returns system information via myenv. If these global variables get deleted at some point, for example by a call to clear all;, user may need to rerun grid_load.m or grid_load_native.m. In such situtations, any call to gcmfaces_global.m will generate a warning that mygrid has not yet been loaded to memory.

Table 2.1 List of grid variables available via the mygrid global variable. The naming convention is directly inherited from the MITgcm naming convention [1].
XC : [1x1 gcmfaces] longitude (tracer)
YC : [1x1 gcmfaces] latitude (tracer)
RC : [50x1 double] depth (tracer)
XG : [1x1 gcmfaces] longitude (vorticity)
YG : [1x1 gcmfaces] latitude (vorticity)
RF : [51x1 double] depth (velocity along 3rd dim)
DXC : [1x1 gcmfaces] grid spacing (tracer, 1st dim)
DYC : [1x1 gcmfaces] grid spacing (tracer, 2nd dim)
DRC : [50x1 double] grid spacing (tracer, 3nd dim)
RAC : [1x1 gcmfaces] grid cell area (tracer)
DXG : [1x1 gcmfaces] grid spacing (vorticity, 1st dim)
DYG : [1x1 gcmfaces] grid spacing (vorticity, 2nd dim)
DRF : [50x1 double] grid spacing (velocity, 3nd dim)
RAZ : [1x1 gcmfaces] grid cell area (vorticity)
AngleCS : [1x1 gcmfaces] grid orientation (tracer, cosine)
AngleSN : [1x1 gcmfaces] grid orientation (tracer, cosine)
Depth : [1x1 gcmfaces] ocean bottom depth (tracer)
hFacC : [1x1 gcmfaces] partial cell factor (tracer)
hFacS : [1x1 gcmfaces] partial cell factor (velocity, 2nd dim)
hFacW : [1x1 gcmfaces] partial cell factor (velocity, 1rst dim)

The C-grid variable names listed in Table 2.1 derive from the MITgcm naming convention [1]. In brief, XC, YC, and RC denote longitude, latitude, and vertical position of tracer variable locations. DXC, DYC, DRC and RAC are the corresponding grid spacings, in m, and grid cell areas, in m\(^2\). A different set of such variables (XG, YG, RF, DXG, DYG, DRF, RAZ) corresponds to velocity and vorticity variables that are staggered in the C-grid approach [1].

Indexing and vector orientation conventions also derive from MITgcm conventions [1]. The indexing convention is illustrated in Figure 2.1. For vector fields, the first component (U) is directed toward the right of the page and the second component (V) toward the top of the page. As compared with tracers, velocity variable locations are shifted by half a grid point to the left of the page (U components) or the bottom of the page (V components) following the C-grid approach [1].

2.3. Exchange Functions

Many computations of interest (e.g., gradients and flow convergences) involve values from nearby grid points on neighboring faces. In practice rows and columns need to be appended at each face edge that are exchanged between neighboring faces – e.g., rows and columns from faces #2, #3, and #5 need to be appended at the face #1 edges in Figure 2.1. Exchanges are operated by exch_T_N.m for tracer-type variables and by exch_UV_N.m for velocity-type variables. These are notably used to compute gradients (calc_T_grad.m) and flow convergences (calc_UV_conv.m).

2.4. Overloaded Functions

As in the case of the Overloaded + function, common operations and functions are overloaded as part of the gcmfaces class definition within the @gcmfaces/ directory:

  1. Logical operators: and, eq, ge, gt, isnan, le, lt, ne, not, or.
  2. Numerical operators: abs, angle, cat, cos, cumsum, diff, exp, imag, log2, max, mean, median, min, minus, mrdivide, mtimes, nanmax, nanmean, nanmedian, nanmin, nanstd, nansum, plus, power, rdivide, real, sin, sqrt, std, sum, tan, times, uminus, uplus.
  3. Indexing operators: subsasgn, subsref, find, get, set, squeeze, repmat.

It may be worth highlighting @gcmfaces/subsasgn.m (subscripted assignment) and @gcmfaces/subsref.m (subscripted reference) since they overload some of the most commonly used Matlab functions. For example, if fld is of the double class then tmp2=fld(1); and fld(1)=1; call subsref.m and subsasgn.m, respectively. If fld instead is of the gcmfaces class then @gcmfaces/subsref.m behaves as follows:

fld{n}     returns the n^{th} face data (i.e., an array).
fld(:,:,n) returns the n^{th} vertical level (i.e., a gcmfaces object).

and @gcmfaces/subsasgn.m behaves similarly but for assignments.

Overloaded + function

function r = plus(p,q)
%overloaded gcmfaces `+' function :
%  simply calls double `+' function for each face data
%  if any of the two arguments is a gcmfaces object
if isa(p,'gcmfaces'); r=p; else; r=q; end;
for iFace=1:r.nFaces;
   if isa(p,'gcmfaces')&isa(q,'gcmfaces');
       eval(['r.f' iF '=p.f' iF '+q.f' iF ';']);
   elseif isa(p,'gcmfaces')&isa(q,'double');
       eval(['r.f' iF '=p.f' iF '+q;']);
   elseif isa(p,'double')&isa(q,'gcmfaces');
       eval(['r.f' iF '=p+q.f' iF ';']);
      error('gcmfaces plus: types are incompatible')

2.5. Input / Output Files

Objects of the gcmfaces class can readily be saved to file using Matlab’s proprietary I/O format (.mat files). Reloading them in a later Matlab session works seamlessly as long as the gcmfaces class has been defined by including @gcmfaces/ to the Matlab path.

Alternatively, gcmfaces variables can be written to files in the nctiles format [FCH+15]. Illustrations in this user guide rely upon ECCO version 4 fields which are distributed in this format (see Section 1.2; Demo Directories Organization and Demo Directory Downloads). The associated I/O functions provided in gcmfaces (write2nctiles.m and read_nctiles.m) reformat data on the fly.

Finally, gcmfaces can read MITgcm binary output in the mds format [2]. The provided I/O functions (rdmds2gcmfaces.m and read_bin.m) rely on convert2gcmfaces.m to convert mds output to gcmfaces objects on the fly. The reverse conversion occurs when convert2gcmfaces.m is called with a gcmfaces input argument. This approach provides a unified framework to analyze MITgcm output or prepare MITgcm input for all known grids (see README and demo_grids.m).

[1](1, 2, 3, 4, 5) For details, see sections 2.11 and 6.2.4 in http://mitgcm.org/public/r2_manual/latest/online_documents/manual.pdf
[2]For details, see section 7.3 in http://mitgcm.org/public/r2_manual/latest/online_documents/manual.pdf