Welcome to the MPAS tutorial practice guide

Throughout this guide, shell commands are written assuming the 'bash' shell; if you're familiar with another shell and can translate the shell command from bash to another shell, please feel free to switch shells.

In order to compile Parallel-NetCDF, PIO, and MPAS, we'll need to ensure that there are working Fortran and C compilers in our path, as well as the MPI compiler wrappers 'mpif90' and 'mpicc'. On the server used for this tutorial, the GNU C and Fortran compilers should already be in your path, as should the MPI compiler wrappers; you may verify this with:

$ which gfortran
$ which gcc

At this point, the GNU Fortran and C compilers should be in our path, and the 'which' commands should return paths to the gfortran and gcc commands. If so, we can confirm that the MPI compiler wrappers are also available:

$ which mpif90
$ which mpicc

The 'which' commands should return paths to the mpif90 and mpicc commands.

Later in the tutorial, we will be using the ncdump, ncview, and ncl to examine NetCDF files created by the model, its pre-processors, and post-processing programs. You may verify that these commands are all in your path with:

$ which ncdump
$ which ncview
$ which ncl

If all of the above commands were successful, your shell path should now contain the compilers that will be used to build Parallel-NetCDF, PIO, and MPAS, as well as some useful visualization commands!

Now that our shell path contains the Fortran and C compilers, as well as the MPI compiler wrappers, we can begin to install the first of two I/O libraries that are needed by MPAS.

In this tutorial, we will install the Parallel-NetCDF and PIO libraries in our own directory. To begin, we'll unpack the Parallel-NetCDF source code:

$ cd ${HOME}/MPAS
$ tar xjvf /home/wrf/mpas/source/parallel-netcdf-1.6.0.tar.bz2
$ cd parallel-netcdf-1.6.0

Before configuring and compiling the code, we must set a few environment variables to tell the configuration script which compilers should be used to build the library.

$ export FC=gfortran
$ export F77=gfortran
$ export CC=gcc
$ export MPIF77=mpif90
$ export MPIF90=mpif90
$ export MPICC=mpicc

Next, we can run the configuration script; when we do this, we add two options: the first is used to disable the C++ interface to the library (since we don't need it for PIO or MPAS) and the second is used to specify where the library should be installed. Recall that we will be installing the library in our own directory; to keep our directory tidy, though, we'll actually install the library in a sub-directory named 'pnetcdf'.

$ ./configure --disable-cxx --prefix=${HOME}/MPAS/pnetcdf

This configuration step generally completes in less than a minute. Once it finishes, we can compile and install the library:

$ make
$ make install

Compilation of the library itself can take several minutes (though probably not more than five or so). After installing the Parallel-NetCDF library, we need to set an environment variable, $PNETCDF, so that MPAS will be able to locate the library. Setting this environment variable is very important!

$ export PNETCDF=${HOME}/MPAS/pnetcdf

Now, we can move on to the installation of the PIO library.

The PIO library is a high-level parallel I/O library originally developed for the Community Earth System Model (CESM). The purpose of the PIO library is to provide a higher-level I/O interface for models, and to handle the problem of shuffling data from the decomposition used by the model to a decomposition that is optimal for performing file input and output.

The PIO library does not perform the file I/O itself, but rather, it passes the shuffled field onto other I/O libraries, such as serial NetCDF, Parallel-NetCDF, MPI-IO, etc. In MPAS, we generally always read and write fields through the Parallel-NetCDF library via PIO, so before configuring PIO, we need to set an environment variable to the installation directory for Parallel-NetCDF:

$ export PNETCDF_PATH=${HOME}/MPAS/pnetcdf

In other words, for the purposes of installing PIO, the environment variable $PNETCDF_PATH must be set to the same value as the environment variable $PNETCDF. However, after PIO has been configured and compiled, there is no longer a need for this environment variable to be set.

Before running the configure script for PIO, we need to unpack the source code and set one additional environment variable (in addition to the FC, F77, CC, MPIF77, MPIF90, and MPICC environment variables that were already set when installing Parallel-NetCDF):

$ cd ${HOME}/MPAS
$ tar xzvf /home/wrf/mpas/source/pio1_7_1.tgz
$ cd pio1_7_1/pio
$ export MPIFC=mpif90

Now, we can run the PIO configuration script with one option, --prefix, to tell the script where PIO should be installed. Additionally, we'll add an option to disable linking to the serial NetCDF library in PIO, since we won't be using the serial output option in MPAS. As with Parallel-NetCDF, we'll install PIO in a subdirectory of our home directory:

$ ./configure --prefix=${HOME}/MPAS/pio --disable-netcdf

Configuration should take just a minute, and when it finishes successfully, we can compile and install the library:

$ make
$ make install

After installing the PIO library, we need to set an environment variable, PIO, so that MPAS will be able to locate the library. Setting this environment variable is very important!

$ export PIO=${HOME}/MPAS/pio

Now, we can move on to compiling the MPAS init_atmosphere and atmosphere cores.

Before beginning to compile any part of MPAS, it's a good idea to verify that the environment variables $PNETCDF and $PIO are properly set:

$ echo $PNETCDF
$ echo $PIO

Does printing these environment variables correctly show the paths where the Parallel-NetCDF and PIO libraries were installed? As a second check, we can verify that these paths actually contain the Parallel-NetCDF and PIO libraries:

$ ls $PNETCDF/lib/libpnetcdf.a
$ ls $PIO/lib/libpio.a

If the libraries themselves are present, we're ready to compile MPAS code.

The first thing we need to do is to unpack the MPAS source code. For this tutorial, a copy of the most recent MPAS release, v4.0, is already available on the server; in general, you can obtain the latest MPAS source code through the "MPAS-Atmosphere download" link on the MPAS homepage: http://mpas-dev.github.io/. We can unpack the MPAS source code in our home directory:

$ cd ${HOME}/MPAS
$ tar xzvf /home/wrf/mpas/source/MPAS-Release-4.0.tar.gz
$ cd MPAS-Release-4.0

After changing directories to MPAS-Release-4.0, it's worth looking around the directory structure and noting the following:

• The 'Makefile' file is used to build any MPAS core, including the MPAS ocean and land-ice cores. There is no code in the Makefile that is specific to MPAS-Atmosphere.
• Inside the 'Makefile' are sets of rules (targets) for compiling MPAS with different sets of compilers; the target that we will use in this tutorial is the 'gfortran' target. Generally, MPAS should be compiled with the same compilers as were used to install the I/O libraries (as well as the MPI implementation, e.g., MPICH).
• Inside the src/ directory are several directories named core_*, each of which corresponds to code that is particular to one MPAS core (e.g., the atmosphere core or the ocean core); in this tutorial, we'll only be concerned with the init_atmosphere and atmosphere cores.
• The src/framework directory contains model infrastructure code that is shared by different MPAS cores.

In order to compile the MPAS init_atmosphere core with the GNU Fortran and C compilers, we will run the following command:

$ make gfortran CORE=init_atmosphere

Compilation should take several minutes. If compilation completes with no errors, the following files should now exist in the MPAS-Release-4.0 directory:

• init_atmosphere_model — the actual executable program for the init_atmosphere core.
• namelist.init_atmosphere — the Fortran namelist for the init_atmosphere core.
• streams.init_atmosphere — the XML I/O configuration file used to specify input and output files for the init_atmosphere core.

Additionally, there will be a directory named default_inputs; the files in this directory just contain the default settings for the init_atmosphere core, and can be used as a reference or to get back to a default state after making changes to the namelist.init_atmosphere or streams.init_atmosphere files.

If the files listed above were created in the MPAS-Release-4.0 directory, then compilation was successful.

After compiling the init_atmosphere core, we can then compile the atmosphere core, which will create an executable for the MPAS-Atmosphere model itself.

It might seem reasonable to run the following make command

$ make gfortran CORE=atmosphere

However, doing so would generate the following error message:

*******************************************************************************
 The MPAS infrastructure is currently built for the init_atmosphere_model core.
 Before building the atmosphere core, please do one of the following.


 To remove the init_atmosphere_model_model executable and clean the MPAS infrastructure, run:
      make clean CORE=init_atmosphere_model

 To preserve all executables except atmosphere_model and clean the MPAS infrastructure, run:
      make clean CORE=atmosphere

 Alternatively, AUTOCLEAN=true can be appended to the make command to force a clean,
 build a new atmosphere_model executable, and preserve all other executables.

*******************************************************************************

We cannot immediately compile the atmosphere core, because the infrastructure code in the src/framework directory was previously compiled for the init_atmosphere core. Before compiling the atmosphere core, we must first clean any compiled infrastructure that will be needed by the atmosphere core with the command

$ make clean CORE=atmosphere

Note that we are specifying the CORE as 'atmosphere' and not 'init_atmosphere'! This is because we don't want to remove the init_atmosphere executable, but instead, we want to simply clean any infrastructure files that are shared by both the init_atmosphere and atmosphere cores.

After running the 'make clean' command, above, the init_atmosphere_model program should still exist. Then, we can compile the atmosphere core with the command

$ make gfortran CORE=atmosphere

Compilation may take around 10 minutes, and if the compilation finishes without errors, the following new files should be created in the MPAS-Release-4.0 directory:

• atmosphere_model — the actual executable program for the MPAS-Atmosphere model core.
• namelist.atmosphere — the Fortran namelist for the MPAS-Atmosphere model core.
• streams.atmosphere — the XML I/O configuration file used to specify input and output files for MPAS-Atmosphere model core.
• stream_list.atmosphere.output — a list of primarily 3-d fields that will be written to the model output (sometimes called "history") files.
• stream_list.atmosphere.diagnostics — a list of primarily 2-d diagnostic fields that will be written to a separate set of files.
• stream_list.atmosphere.surface — a list of fields that may be periodically read to update the SST and sea-ice fields as the model is integrating.

Observe that several symbolic links have also been created in the MPAS-Release-4.0 directory, e.g., VEGPARM.TBL, LANDUSE.TBL, RRTMG_LW_DATA_DBL, RRTMG_SW_DATA_DBL, etc. These links point to look-up tables that are used by the WRF radiation and land-surface physics schemes that are used in MPAS-Atmosphere.

If the files mentioned above were created in the MPAS-Release-4.0 directory, then compilation was successful.

The process of generating a new variable-resolution mesh for MPAS involves a long series of iterations that can take quite some time to converge; accordingly, generating a new variable-resolution mesh is an expensive process, and whenever possible, it is preferable to "re-use" an existing variable-resolution mesh by simply shifting or rotating the refinement region from its previous location to the desired location.

Rotating an MPAS mesh is accomplished using the grid_rotate utility program that may be downloaded from the same page as the MPAS meshes. For this tutorial, the grid_rotate source code is available on the server in /home/wrf/mpas/source/grid_rotate.tar.gz . We'll begin by compiling this utility program from source in a directory parallel to our MPAS-Release-4.0 directory. First, we need to change directories and unpack the source code:

$ cd ${HOME}/MPAS
$ tar xzvf /home/wrf/mpas/source/grid_rotate.tar.gz
$ cd grid_rotate

The grid_rotate program relies on the standard (serial) NetCDF library, and it looks for this library in the path indicated by the $NETCDF environment variable; we'll next this variable and compile the code:

$ export NETCDF=/home/wrf/mpas/local
$ make

Compiling the program should take just a few seconds. If compilation was successful, the grid_rotate program should now be present. In general, it may be necessary to edit the Makefile to select the proper compiler, compiler flags, or library paths before compiling.

As described in Section 4.2 of the MPAS-Atmosphere Users' Guide, all that is needed to run the grid_rotate utility is to edit the namelist.input file to describe how the refinement region should be rotated, then to run the utility itself.

On the server, a variable-resolution mesh that refines from about 60-km resolution over the globe down to about 15-km resolution over an area useful for hurricane or typhoon forecasting in the Western Pacific is available in /home/wrf/mpas/meshes/x4.535554.grid.nc. Before proceeding, we'll create a symbolic link to this file (there's no need to copy the original) in our grid_rotate directory:

$ ln -s /home/wrf/mpas/meshes/x4.535554.grid.nc .

This SCVT grid file has the refinement region centered over 0.0° latitude, 0.0° longitude, as shown in the figure below.

As an aside, MPAS-Atmosphere meshes that are available through the MPAS-Atmosphere download page are typically named according to the refinement factor of the mesh ("x1" implies no refinement, while "x4" implies refinement by a factor of four) and the number of grid cells in the horizontal mesh. So, the "x4.535554" mesh refines by a factor of four and has a total of 535,554 grid cells.

For the purposes of this tutorial, you may rotate the refinement region wherever you would like. For example, setting up the grid_rotate program's namelist.input file like this

&input
   config_original_latitude_degrees = 0
   config_original_longitude_degrees = 0

   config_new_latitude_degrees = 25
   config_new_longitude_degrees = 150

   config_birdseye_rotation_counter_clockwise_degrees = 0
/

and running the grid_rotate program like this

$ ./grid_rotate x4.535554.grid.nc x4.535554.rotated.grid.nc

should yield the following rotated grid:

In order to visualize the location of the refined region in the rotated mesh, you can make use of the mesh_resolution.ncl NCL script provided in the grid_rotate utility directory. To run this script, you will need to first set the environment variable FNAME to the name of the rotated grid file that you created with the grid_rotate program, e.g.,

$ export FNAME=x4.535554.rotated.grid.nc
$ ncl mesh_resolution.ncl

where "x4.535554.rotated.grid.nc" is replaced with the name of your rotated grid file. To display the PDF file that is created by the script, we can use the command

$ display -rotate -90 mesh_resolution.pdf

Rotating the mesh doesn't typically require much time, so feel free to try moving the refinement region in a few different ways. In subsequent steps, we'll assume that the rotated mesh file is named x4.535554.rotated.grid.nc; you can either rename your rotated mesh using this filename, or make the appropriate substitutions in the namelist, I/O control, and other files.

Regardless of whether a variable-resolution mesh or a uniform-resolution mesh will be used with MPAS-Atmosphere, the first step in preparing real-data initial conditions is to interpolate static (i.e., not changing in time), geographical fields onto the mesh. The geographical datasets used in MPAS-Atmosphere are the same as those used in the Weather Research and Forecasting (WRF) model's preprocessing system (WPS). These datasets can be downloaded from the WRF model download page, or through the link in Section 7.2.1 in the MPAS-Atmosphere Users' Guide. However, for this tutorial, the geographical datasets are already available on the server in the directory /home/wrf/geog/.

The process of interpolating geographical data can in principle be done in any working directory provided the necessary files are linked or copied there; however, for simplicity, we'll be running in our MPAS-Release-4.0 directory, where the init_atmosphere_model program was compiled. Before proceeding, we'll change to this directory:

$ cd ${HOME}/MPAS/MPAS-Release-4.0

Important note: The 60—15-km mesh that we rotated in the previous section would generate unnecessarily large initial condition files for the purposes of this tutorial, and it would take longer to process than necessary; so, we'll assume that at this point, we understand how to prepare a variable-resolution mesh, and for the rest of the process of creating initial conditions, we'll switch to a coarser, uniform-resolution mesh. The uniform mesh that we will use has 10242 horizontal grid cells, with a resolution of approximately 240 km. Before proceeding, we should link this grid file into our working directory:

$ ln -s /home/wrf/mpas/meshes/x1.10242.grid.nc .

Following Section 7.2.1 of the Users' Guide, we begin the process of interpolating the geographical data onto the MPAS mesh by editing the namelist.init_atmosphere file in the MPAS-Release-4.0 directory. The key namelist options that must be set are shown below; other namelist options can be ignored.

&nhyd_model
    config_init_case = 7
/

&dimensions
    config_nvertlevels = 1
    config_nsoillevels = 1
/

&data_sources
    config_geog_data_path = '/home/wrf/geog/'
    config_landuse_data = 'USGS'
/

&preproc_stages
    config_static_interp = true
    config_vertical_grid = false
    config_met_interp = false
    config_input_sst = false
    config_frac_seaice = false
/

After editing the namelist.init_atmosphere file, it will also be necessary to tell the init_atmosphere_model program the name of our input grid file, as well as the name of the "static" output file that we would like to create. We do this by editing the streams.init_atmosphere file, where we first set the name of the grid file for the "input" stream:

<immutable_stream name="input"
                  type="input"
                  filename_template="x1.10242.grid.nc"
                  input_interval="initial_only"/>

and then set the name of the static file to be created for the "output" stream:

<immutable_stream name="output"
                  type="output"
                  filename_template="x1.10242.static.nc"
                  packages="initial_conds"
                  output_interval="initial_only" />

When interpolating geographical fields to a mesh, the "surface" stream can be ignored.

After editing the namelist.init_atmosphere and streams.init_atmosphere files, ensure that the x1.10242.grid.nc file is moved or linked into the MPAS-Release-4.0 directory, where we will run the init_atmosphere_model program to interpolate geographical data.

$ ./init_atmosphere_model

If the init_atmosphere_model program successfully ran, the log.0000.out file should contain timing information:

$ cat log.0000.out 
   timer_name                total    calls         min           max           avg    percent  efficiency
 0 total time             231.75340        1     231.75340     231.75340     231.75340
 1  initialize             77.30472        1      77.30472      77.30472      77.30472   0.33   0.00

and the log.0000.err file should contain no error messages. Additionally, the file x1.10242.static.nc should have been created. To check that the interpolation of the static, geographical fields worked as expected, you can use the NCL script from /home/wrf/mpas/ncl_scripts/plot_terrain.ncl to plot the terrain field after setting the FNAME environment variable to "x1.10242.static.nc":

$ cp /home/wrf/mpas/ncl_scripts/plot_terrain.ncl .
$ export FNAME=x1.10242.static.nc
$ ncl plot_terrain.ncl

If the init_atmosphere_model program ran successfully and produced reasonable geographical fields, you're ready to generate a vertical grid and to interpolate meteorological initial conditions to the grid.

It should be noted that, because the geographical fields (e.g., terrain, land use, etc) do not change in time, the interpolation of these fields only needs to be done once for a particular mesh. The same "static" file can be re-used for many different MPAS-Atmosphere simulations

The "static" file (e.g., x1.10242.static.nc) created by the init_atmosphere model has the terrain, land use category, soil category, and other geographic iterpolated to the horizontal MPAS mesh, but the file does not contain any information about the distribution of vertical levels in the atmosphere, nor does it contain any atmospheric fields. The next step is to generate a vertical grid, and to then interpolate atmospheric initial conditions to full 3-d MPAS grid.

Following Section 6.2.3 of the Users' Guide, we must edit the namelist.init_atmosphere file to specify the particular time of the atmospheric initial condition fields, the number of vertical levels to be generated, and the height of the top of the model. The key namelist options that must be set are shown below; other namelist options can be ignored.

&nhyd_model
    config_init_case = 7
    config_start_time = '2014-09-10_00:00:00'
    config_theta_adv_order = 3
/

&dimensions
    config_nvertlevels = 41
    config_nsoillevels = 4
    config_nfglevels = 38
    config_nfgsoillevels = 4
/

&data_sources
    config_met_prefix = 'GFS'
/

&vertical_grid
    config_ztop = 30000.0
    config_nsmterrain = 1
    config_smooth_surfaces = true
/

&preproc_stages
    config_static_interp = false
    config_vertical_grid = true
    config_met_interp = true
    config_input_sst = false
    config_frac_seaice = false
/

Besides the namelist.init_atmosphere file, we must also edit the XML I/O configuration file, streams.init_atmosphere, to specify the name of the static file that will serve as input to the init_atmosphere_model program, as well as the name of the MPAS-Atmosphere initial condition file to be created. Specifically, we must set the filename_template for the "input" stream to the name of our static file:

<immutable_stream name="input"
                  type="input"
                  filename_template="x1.10242.static.nc"
                  input_interval="initial_only"/>

and we must set the name of the initial condition file to be created in the "output" stream:

<immutable_stream name="output"
                  type="output"
                  filename_template="x1.10242.init.nc"
                  packages="initial_conds"
                  output_interval="initial_only" />

Before running the init_atmosphere_model program, we must first link a WPS intermediate format file into our working directory; this intermediate file contains the atmospheric and land-surface initial conditions for the model run, and its name must match the parameters, config_met_prefix and config_start_time specified in the namelist.init_atmosphere file.

Ordinarily, an intermediate file must be prepared using the WPS ungrib.exe program described in the WRF ARW Users' Guide; however, for this tutorial, intermediate data have already been processed, and can be found in the directory /home/wrf/mpas/met_data.

$ ln -s /home/wrf/mpas/met_data/GFS:2014-09-10_00 .

After editing the namelist.init_atmosphere and streams.atmosphere files, and linking the intermediate file to our run directory, we can run the init_atmosphere_model program.

$ ./init_atmosphere_model

Because a larger volume of data is being processed than in the geographical data interpolation step, this execution of the init_atmosphere core can take a few minutes. You can check on the progress of the program by 'tailing' the log.0000.err file from another terminal window, e.g.,

$ tail -f log.0000.err

Once the program finishes, if the log.0000.out file contains some timing statistics and the log.0000.err file contains no error messages, then the init_atmosphere_model program probably ran correctly. The program should have created the output file specified in the streams.init_atmosphere file, x1.10242.init.nc. To check that reasonable atmospheric initial conditions exist in the x1.10242.init.nc file, you can use the NCL script from /home/wrf/mpas/ncl_scripts/plot_qv.ncl to plot the water vapor mixing ratio in the lowest model level after setting the FNAME environment variable to the name of this initial condition file.

$ cp /home/wrf/mpas/ncl_scripts/plot_qv.ncl .
$ export FNAME=x1.10242.init.nc
$ ncl plot_qv.ncl

The x1.10242.init.nc file is the MPAS-Atmosphere initial condition file. However, because performing an MPAS-Atmosphere simulation on this mesh could take quite some time, we will actually use a coarser mesh for which initial conditions have already been prepared when it comes time to run the model itself.

Because MPAS-Atmosphere — at least, when run as a stand-alone model — does not contain prognostic equations for the SST and sea-ice fraction, these fields would remain constant if not updated from an external source; this is, of course, not realistic, and it will generally impact the quality of longer model simulations. Consequently, for MPAS-Atmosphere simulations longer than roughly a week, it is typically necessary to periodically update the sea-surface temperature (SST) field in the model. For real-time simulations, this is generally not an option, but it is feasible for retrospective simulations, where we have observed SST analyses available to us.

In order to create a "surface update" file for MPAS-Atmosphere, we will make use of a sequence of intermediate files containing SST and sea-ice analyses that are available on the server in /home/wrf/mpas/met_data. Before proceeding, we'll link all of these intermediate files into our working directory:

$ ln -sf /home/wrf/mpas/met_data/GFS* .

Following Section 7.2.2 of the MPAS-Atmosphere Users' Guide, we must edit the namelist.init_atmosphere file to specify the range of dates for which we have SST data, as well as the frequency at which the SST intermediate files are available. The key namelist options that must be set are shown below; other options can be ignored.

Since the NCEP GFS analysis files contain an SST and a sea-ice field, we can actually use the GFS intermediate files as a source of SST data for the purposes of creating a surface update file.

&nhyd_model
    config_init_case = 8
    config_start_time = '2014-09-10_00:00:00'
    config_stop_time = '2014-09-20_00:00:00'
/

&data_sources
    config_sfc_prefix = 'GFS'
    config_fg_interval = 86400
/

&preproc_stages
    config_static_interp = false
    config_vertical_grid = false
    config_met_interp = false
    config_input_sst = true
    config_frac_seaice = true
/

Note in particular that we have set the config_init_case variable to 8! This is the initialization case used to create surface update files, instead of real-data initial condition files.

As before, we also need to edit the streams.init_atmosphere file, this time setting the name of the SST update file to be created by the "surface" stream, as well as the frequency at which this update file should contain records:

<immutable_stream name="surface"
                  type="output"
                  filename_template="x1.10242.sfc_update.nc"
                  filename_interval="none"
                  packages="sfc_update"
                  output_interval="86400"/>

Note that we have set both the filename_template and output_interval attributes of the "surface" stream. The output_interval should match the interval specified in the namelist for the config_fg_interval variable. The other streams ("input" and "output") can remain unchanged — the input file should still be set to the name of the static file.

After setting up the namelist.init_atmosphere and streams.atmosphere file, we're ready to run the init_atmosphere_model program:

$ ./init_atmosphere_model

After the job completes, as before, the log.0000.out file should contain timing statistics, and the log.0000.err file should contain no error messages. Additionally, the file "x1.10242.sfc_udpate.nc" file should have been created.

We can confirm that the x1.10242.sfc_udpate.nc contains the requested valid times for the SST and sea-ice fields by printing the "xtime" variable from the file; in MPAS, the "xtime" variable in a NetCDF file is always used to store the times at which records in the file are valid.

$ ncdump -v xtime x1.10242.sfc_update.nc

This ncdump command should print the following times:

 xtime =
  "2014-09-10_00:00:00                                             ",
  "2014-09-11_00:00:00                                             ",
  "2014-09-12_00:00:00                                             ",
  "2014-09-13_00:00:00                                             ",
  "2014-09-14_00:00:00                                             ",
  "2014-09-15_00:00:00                                             ",
  "2014-09-16_00:00:00                                             ",
  "2014-09-17_00:00:00                                             ",
  "2014-09-18_00:00:00                                             ",
  "2014-09-19_00:00:00                                             ",
  "2014-09-20_00:00:00                                             " ;

If the times shown above were printed, then we have successfully created an SST and sea-ice update file for MPAS-Atmosphere. As a final check, we can use the NCL script from /home/wrf/mpas/ncl_scripts/plot_delta_sst.ncl to plot the difference in the SST field between the last time in the surface update file and the first time in the file. Because the surface update file contains no information about the latitude and longitude of grid cells, we need to specify two environment variables when running this script: one to give the name of the surface update file, and the second to give the name of the grid file, which contains cell latitude and longitude fields:

$ cp /home/wrf/mpas/ncl_scripts/plot_delta_sst.ncl .
$ export FNAME=x1.10242.sfc_update.nc
$ export GNAME=x1.10242.static.nc
$ ncl plot_delta_sst.ncl

In this plot, we have masked out the SST differences over land, since the values of the field over land are not representative of actual SST differences, but may represent differences in, e.g., skin temperature or 2-m air temperature, depending on the source of the SST analyses.

Assuming that an initial condition file, and, optionally, a surface update file have been created, we're ready to run the MPAS-Atmosphere model itself!

Important note: In order for us to be able to simulate several days in a reasonable amount of time, rather than using the 240-km initial conditions and SST update file that we created in previous sections, we will actually be running the model using an even coarser, uniform grid for which initial conditions and an SST (surface) update file have already been prepared. This mesh has 2562 horizontal grid cells and an approximate resolution of 480 km. The input files are located on the server in the directory /home/wrf/mpas/input.

At this point, we can either delete or move to another directory the 240-km files, that were created in previous sections, since we will be actually running the model on a coarser mesh. Doing this will create a cleaner working directory for us to more easily see the files created by the model:

$ mkdir ../mpas_initialization
$ mv x1.10242* ../mpas_initialization

To begin running the model, we'll need to link (copying is unnecessary) the initial condition file, the surface update file, and the mesh decomposition files from /home/wrf/mpas/input to our MPAS-Release-4.0 directory:

$ ln -s /home/wrf/mpas/input/x1.2562.init.nc .
$ ln -s /home/wrf/mpas/input/x1.2562.sfc_update.nc .
$ ln -s /home/wrf/mpas/meshes/x1.2562.graph.info.part.4 .

Next, we will need to edit the namelist.atmosphere file, which is the namelist for the MPAS-Atmosphere model itself. The default namelist.atmosphere is set up for a 120-km mesh; in order to run the model on a different mesh, there are several key parameters that depend on the model resolution:

• config_dt — the model timestep (delta-t) in seconds.
• config_len_disp — the length-scale for explicit horizontal diffusion, in meters.

To tell the model the date and time at which integration will begin (important, e.g., for computing solar radiation parameters), and to specify the length of the integration to be run, there are two other parameters that must be set for each simulation:

• config_start_time — the start time of the integration.
• config_run_duration — the length of the integration.

Lastly, when running the MPAS-Atmosphere model in parallel, we must tell the model the prefix of the filenames that contain mesh partitioning information for different MPI task counts:

• config_block_decomp_file_prefix — the file prefix (to be suffixed with processor count) of the file containging mesh partitioning information.

Accordingly, we must edit at least these parameters in the namelist.atmosphere file; other parameters are described in Appendix B of the MPAS-Atmosphere Users' Guide, and do not necessarily need to be changed:

&nhyd_model
    config_dt = 1800.0
    config_start_time = '2014-09-10_00:00:00'
    config_run_duration = '2_00:00:00'
    config_len_disp = 480000.0
/

&decomposition
    config_block_decomp_file_prefix = 'x1.2562.graph.info.part.'
/

After changing the parameters shown above in the namelist.atmosphere file, we must also set the name of the initial condition file, the name of the surface update file, and the interval at which we would like to read the surface update file in the streams.atmosphere file:

<immutable_stream name="input"
                  type="input"
                  filename_template="x1.2562.init.nc"
                  input_interval="initial_only"/>

<stream name="surface"
        type="input"
        filename_template="x1.2562.sfc_update.nc"
        filename_interval="none"
        input_interval="86400">

    <file name="stream_list.atmosphere.surface"/>

</stream>

Having set up the namelist.atmosphere and streams.atmosphere files, we're ready to run the model in parallel on the cluster. As a general rule, a mesh with N grid columns should be run on at most N/160 processors in MPAS-Atmosphere to make efficient use of the processors. So, for the mesh in this exercise, which has 2562 grid columns, we should not use more than about 16 MPI tasks in order to efficiently make use of computational resources. For this exercise, using just four processors should be sufficient:

$ mpiexec -n 4 ./atmosphere_model

Once the model begins to run, it will write information about its progress to the log.0000.err file. You can "tail" this file with the "-f" option to follow the updates to the file as the model is running:

$ tail -f log.0000.err

If the model has started up successfully and begun to run, we should see messages like this

 Begin timestep 2014-09-10_00:00:00
 --- time to update background surface albedo, greeness fraction.
 --- time to run the LW radiation scheme L_RADLW = T
 --- time to run the SW radiation scheme L_RADSW = T
 --- time to run the convection scheme L_CONV    = T
 --- time to apply limit to accumulated rainc and rainnc L_ACRAIN   = F
 --- time to apply limit to accumulated radiation diags. L_ACRADT   = F
 --- time to calculate additional physics_diagnostics               = F
 
 global min, max w  -0.128822076860735       0.213503433301904     
 global min, max u   -114.854987176182        115.563325406686   

written for each timestep that the model takes. If not, the log.0000.err file may have some messages to indicate what might have gone wrong.

When the model finishes simulating for the two days that we requested in the namelist.atmosphere file, we should have several different sets of output files:

• diag.2014.09.*.nc — These files contain mostly 2-d diagnostic fields that were listed in the stream_list.atmosphere.diagnostics file.
• history.2014.09.*.nc — These files contain mostly 3-d prognostic and diagnostic files from the model.
• restart.2014.09.*.nc — These are model restart files that are essentially checkpoints of the model state; a simulation can be re-started from any of these checkpoint/restart files.

If you see these files, then you've successfully performed a global simulation using MPAS-Atmosphere!

When we ran the MPAS-Atmosphere model, the model produced restart files (named restart.2014-09-*.nc) once per day. The frequency with which the model creates these restart files is controlled in the XML I/O configuration file streams.atmosphere in the definition of the "restart" stream:

<immutable_stream name="restart"
                  type="input;output"
                  filename_template="restart.$Y-$M-$D_$h.$m.$s.nc"
                  input_interval="initial_only"
                  output_interval="1_00:00:00"/>

Notice in particular that the "output_interval" option for the "restart" stream was set to one day, which is the interval at which we now have restart files.

If we would like to continue the simulation beyond the original two days in our simulation, we could restart the model from the restart.2014-09-12_00.00.00.nc file. There are many reasons why we might want or need to restart the model — perhaps the cluster that we are running on only allows jobs to run for a certain number of hours, and to complete a long MPAS-Atmosphere simulation, we need to complete the simulation in a series of shorter runs using restart files, for example.

To try out the restart capability in MPAS-Atmosphere, we will restart the model from 2014-09-12_00:00:00. Since we have a restart file valid at this time, all that we need to do is to edit three variables in the namelist.atmosphere file:

&nhyd_model
    config_start_time = '2014-09-12_00:00:00'
    config_run_duration = '1_00:00:00'
/
&restart
    config_do_restart = true
/

In the namelist.atmosphere file, we are specifying the time from which we would like to restart the model, the duration of the integration, and, most imporantly, that we would like the simulation to restart from an existing restart file.

After making these three changes in the namelist.atmosphere file, we can run the model, again using four processors:

$ mpiexec -n 4 ./atmosphere_model

Once the model begins to run, note the valid time of the timesteps in the log.0000.err file; the model has started simulating from 2014-09-12_00:00:00.

After the model finishes running, we now have model simulation output for a total of three days — from 2014-09-10_00:00:00 through 2014-09-13_00:00:00. A key property of model restarts is that they produce bitwise-identical results; that is, running the model continuously for three days, or restarting after two days and simulating for an additional day, should both give identical results.

A new feature in MPAS v3.0 is the ability to write additional model output files without needing to edit the model's Registry.xml file and re-compile the model. To see how this capability works, we will make another restart run from our last restart file, restart.2014-09-13_00.00.00.nc, but before doing so, we will define a new stream in the streams.atmosphere file that will cause the 10-m diagnosed winds to be written every 60 simulated minutes to a new set of output files.

A more complete description of input and output streams is given in Chapter 5 of the MPAS-Atmosphere Users' Guide; for this exercise, we will simply describe the changes to be made to the streams.atmosphere file.

Notice that the first and last lines in the streams.atmosphere file are opening and closing "streams" tags:

<streams>

      [... some stream definitions ...]

</streams>

All streams in MPAS must be defined between these opening and closing tags!

We can define a new stream, which we will call the "sfc_winds" stream, by adding the following to the end of the streams.atmosphere file just above the closing </streams> tag:

<stream name="sfc_winds"
        type="output"
        filename_template="sfc_winds.$Y-$M-$D_$h.$m.$s.nc"
        output_interval="1:00:00">

    <var name="u10"/>
    <var name="v10"/>

</stream>

Observe that the only two variables that will be written to this stream are the "u10" and "v10" fields, which represent the diagnosed 10-m wind zonal and meridional components. Before running the model, we will also set the start time in the namelist to "2014-09-13_00:00:00":

&nhyd_model
    config_start_time = '2014-09-13_00:00:00'
/

The config_do_restart variable should already be set to "true" from when we restart the model previously. After updating the namelist and adding the stream definition to the streams.atmosphere file, we can run the model again using four processors:

$ mpiexec -n 4 ./atmosphere_model

After the model has successfully completed running, in addition to the usual history.*.nc, diag.*.nc, and restart.*.nc files, we should also have a new set of files named sfc_winds.2014-09-??_??.00.00.nc. As requested in the stream definition, these files are written every simulated hour (i.e., with an interval of "1:00:00").

One of the simplest types of plots to create from MPAS-Atmosphere output using NCL is a horizontal contour plot. In this exercise, we'll make a horizontal contour plot using the example script in /home/wrf/mpas/ncl_scripts/plot_contours.ncl. Since we will need to edit this script, we'll begin by copying it to our working directory:

$ cp /home/wrf/mpas/ncl_scripts/plot_contours.ncl .

In order to work with higher-resolution model output than the 480-km simulation that we performed earlier, we'll make use of some pre-computed model output that is available on the server; for simplicity we will symbolically link this output to our working directory:

$ ln -s /home/wrf/mpas/output/output.2014*.nc .

The example "plot_contours.ncl" script uses the environment variable $FNAME to specify the name of the MPAS-Atmosphere output file that will be plotted, and the environment variable $T to specify the time record to be plotted; in NCL, indexing of variables begins at 0, so setting T to 0 would specify the first time in the file. In the script, the default field to plotted is the "qv" (water vapor mixing ratio) field at the lowest model level. To plot the lowest model level "qv" field, then, we just need to:

$ export FNAME=output.2014-01-11_00.00.00.nc
$ export T=0
$ ncl plot_contours.ncl

The output of the script should appear in a file named "contours.pdf" and should look something like this:

In addition to horizontal contour plots of model fields, it can also be useful when analyzing a simulation to produce vertical cross-sections of fields. Using the example script /home/wrf/mpas/ncl_scripts/plot_xsec.ncl, we'll see how we can create vertical cross sections with endpoints of our choosing.

To begin, we'll make a copy of the plot_xsec.ncl script into our working directory. By default, this script will plot a cross-section of the "qv" (water vapor mixing ratio) field, with the cross-section starting at the point (40.00° N, 140° W) and ending at the point (40.00° N, 80° W).

We can change the starting and ending points of the cross-section so that it crosses over Taiwan by editing the following variables in the script:

  start_lat = 23.58
  start_lon = 90.0
  end_lat = 23.58
  end_lon = 150.0

After making the changes to these variables, we can run the script to produce a cross-section of the "qv" field that should look something like the following:

The unstructured horizontal mesh used by MPAS can often be problematic for some visualization tools to work with, and makes it more difficult for us to quickly check model input or output fields when trying to track down problems in the model. For these reasons, it can be helpful to interpolate MPAS fields from their native, unstructured mesh onto a regular, rectangular mesh such as a regular latitude-longitude grid.

We can use the NCL script found in /home/wrf/mpas/ncl_scripts/mpas_to_latlon.ncl to remap (i.e., interpolate) any cell-based field from an MPAS mesh onto a latitude-longitude grid. To begin, let's first copy this script to our MPAS-Release-3.1 directory, since we will be making changes to the script:

$ cp /home/wrf/mpas/ncl_scripts/mpas_to_latlon.ncl .

The "mpas_to_latlon.ncl" script is generally run in two stages:

1. In the first step, we specify the resolution of the latitude-longitude grid to which we would like to interpolate, and we run the script to generate the remapping weights (i.e., interpolation weights) that will be used in the second step to actually interpolate from the MPAS mesh to the latitude-longitude grid.
2. In the second step, after having generated remapping weights for a particular MPAS mesh and a particular resolution of latitude-longitude grid, we can run the script to remap any cell-centered field of our choosing from the MPAS mesh to the latitude-longitude grid.

For the first step, we'll edit the mpas_to_latlon.ncl script to specify a 1° x 1° latitude-longitude grid by changing this line

   opt@DstGridType = "0.5x0.5" ; destination grid

to

   opt@DstGridType = "1.0x1.0" ; destination grid

To ensure that the script will attempt to generate remapping weights, rather than to actually remap a field, also verify that the variable "gen_weights" near the top of the script is set to "True":

gen_weights = True

By default, the mpas_to_latlon.ncl script will attempt to generate remapping weights from MPAS output found in files named "history*ncl"; for this exercise, we'll need to change the names of the output files to "x1.40962.output.*.nc" in the script; we'll locate this line in the script:

   files = systemfunc("ls history.*.nc")

and change it to:

   files = systemfunc("ls output.2014*.nc")

After making these changes to the script, we can make symbolic links to pre-computed MPAS output on an approximately 120-km mesh to our working directory, then run the script to generate remapping weights:

$ ln -sf /home/wrf/mpas/output/output.2014*.nc .
$ ncl mpas_to_latlon.ncl

Generating the remapping weights may take a minute. After the script finishes running, if the generation of the remapping weights was successful, we should now have the following files in our directory: src_grid.nc, dst_grid.nc, and weights.nc. If we have these files, we're ready to remap an MPAS field in the second step of the process.

To remap a field, we must first edit the mpas_to_latlon.ncl script and set the "gen_weights" variable to "False"; doing so will cause the script to remap a field, rather than to generate remapping weights.

By default, the script is set up to remap the "qv" (water vapor mixing ratio) field. We'll try to change this variable later; for now, after setting "gen_weights = False", we can run the script a second time to create a NetCDF file containing the "qv" field on a 1° x 1° latitude-longitude grid:

$ ncl mpas_to_latlon.ncl 

The remapping should take very little time to finish; if it is successful, a new file, named "latlon.nc", should be created. We can now use the "ncview" command to quickly scan the "qv" field at different levels.

$ ncview latlon.nc

As a final part of this exercise, we'll try to change the field to be remapped from "qv" to the zonal velocity field, "uReconstructZonal". Edit the mpas_to_latlon.ncl script, locate the occurrences of the variable "qv", and change them to "uReconstructZonal". After doing this, remove the "latlon.nc" output file and re-run the script. Note that, because the remapping weights depend only on the MPAS mesh and our chosen latitude-longitude grid, changing variables does not require us to re-compute the remapping weights. So, we can just run the mpas_to_latlon.ncl script with "gen_weights = False" to create a new "latlon.nc" file with the remapped "uReconstructZonal" field.

Because MPAS-Atmosphere may perform simulations on variable-resolution meshes, care must be taken when computing global averages of fields: if values aren't weighted in proportion to the areas of the grid cells that they represent, it's possible to arrive at very biased values for the mean of the field.

In this exercise, we will compute the global mean precipitation amount during a six-hour period in a variable-resolution simulation in two different ways: the first, incorrect method will give equal weight to every grid cell, while the second, correct method will weight the precipitation in an MPAS grid cell in proportion to that cell's area.

On the server, the files /home/wrf/mpas/output/x4.153762.output.2012-12-01_*.00.00.nc contain output from a variable-resolution MPAS-Atmosphere simulation. We will begin by linking this output files to our working directory:

$ ln -s /home/wrf/mpas/output/x4.153762.output.2012-12-01_* .

The total accumulated precipitation from an MPAS-Atmosphere simulation is given by the sum of the explicit and the parameterized precipitation, i.e., the sum of the variables "rainnc" and "rainc". To get the precipitation over some period, we can simply subtract the precipitation at the beginning of the period from the precipitation at the end of the period. We can use the NCL script from /home/wrf/mpas/ncl_scripts/plot_6hr_precip.ncl to plot the 6-hour accumulated precipitation between 2012-12-01 1200 UTC and 2012-12-01 1800 UTC:

$ cp /home/wrf/mpas/ncl_scripts/plot_6hr_precip.ncl .
$ ncl plot_6hr_precip.ncl

The precipitation field should look something like the following:

Now, to compute the global mean 6-hourly precipitation, we can launch an interactive NCL session and run the following commands:

$ ncl
 Copyright (C) 1995-2014 - All Rights Reserved
 University Corporation for Atmospheric Research
 NCAR Command Language Version 6.2.0
 The use of this software is governed by a License Agreement.
 See http://www.ncl.ucar.edu/ for more details.
ncl 0> load "$NCARG_ROOT/lib/ncarg/nclscripts/csm/contributed.ncl"
ncl 1> f1 = addfile("x4.153762.output.2012-12-01_12.00.00.nc","r")
ncl 2> f2 = addfile("x4.153762.output.2012-12-01_18.00.00.nc","r")
ncl 3> fld = (f2->rainc(0,:) + f2->rainnc(0,:)) - (f1->rainc(0,:) + f1->rainnc(0,:))
ncl 4> dims = dimsizes(fld)
ncl 5> nCells = dims(0)
ncl 6> print(sum(fld)/int2flt(nCells))
(0)	1.127321613497147

So, without weighting the precipitation amounts in each grid cell by the relative area of the cell, we find a global 6-hour precipitation amount of 1.1273 mm.

With the NCL session still open, we can now compute the global mean in the correct way, by multiplying the precipitation value from each cell by the fraction of the total area of the mesh represented by that cell:

ncl 7> print(sum(fld * f1->areaCell(:)) / sum(f1->areaCell(:)))
(0)	0.8885738885384865
ncl 8> quit

Weighting the precipitation by the relative area of each grid cell yields a (correct) value of 0.88857 mm, which is significantly less than the value computed earlier! The key point of this exercise is to demonstrate that, when working on a variable-resolution MPAS mesh, it is very important to compute averages of fields by taking into account the relative areas of the cells in the mesh.

Often, it is necessary to efficiently find the cell in an MPAS mesh that is nearest to some specified latitude and longitude location. A simple but inefficient method to find this nearest cell would be to scan all the cells in the mesh, computing the distance from each cell to the specified location, then choosing the cell that gave a minimum distance. When there are only a few locations for which we need the nearest cell, this method is acceptable, but when there are many such locations, a more efficient method is often needed.

Recall from the lecture earlier today that we can efficiently find a cell nearest to some specified location using a simple search; the algorithm was described as:

We can see an illustration of the search path taken by this algorithm in the figure below, assuming that the solid blue cell represents our starting cell, and the small green circle represents the latitude and longitude location that we are searching for:

In order to try out this algorithm ourselves, we can copy the script from /home/wrf/mpas/ncl_scripts/nearest_cell.ncl to our working directory

$ cp /home/wrf/mpas/ncl_scripts/nearest_cell.ncl .

and edit the location to search for in these two lines:

  search_lat = 40.0
  search_lon = -140.0

For this exercise, feel free to set the location to any that you'd like, for example, to the latitude and longitude of Taipei.

After editing the script to set the target for the search, we'll need to specify the name of an MPAS input or output file that contains the mesh in which we would like to search using the environment variable $FNAME; for example, if we still had the variable-resolution mesh output from the previous section in our working directory, we could run:

$ export FNAME=x4.153762.output.2012-12-01_18.00.00.nc
$ ncl nearest_cell.ncl

After running the script, does it find a cell whose latitude and longitude is closest among those in the mesh to your specified target location? Remember: the latitude and longitude of the closest cell should be within one grid distance of your target location!