program model_mod_check

DART project logo

Jump to DART Documentation Main Index
version information for this file:
$Id: model_mod_check.html 11997 2017-10-17 22:26:56Z nancy@ucar.edu $

NAMELIST / MODULES / FILES / USAGE / REFERENCES / ERRORS / PLANS / TERMS OF USE

Overview

model_mod_check tests some of the more fundamental routines in any model_mod. This is intended to be used when adding a new model to DART - test the pieces as they are written. As such, this program is meant to be hacked up and customized to your own purpose. Right now, it reads in model netCDF file(s) - one per domain/nest/whatever - and writes out files, queries the metdata, etc. It also exercises static_init_model(), which is the first routine to get right ...

[top]

NAMELIST

This namelist is read from the file input.nml. Namelists start with an ampersand '&' and terminate with a slash '/'. Character strings that contain a '/' must be enclosed in quotes to prevent them from prematurely terminating the namelist.

&model_mod_check 
   num_ens               = 1
   single_file           = .FALSE.
   input_state_files     = 'null'
   output_state_files    = 'null'
   all_metadata_file     = 'metadata.txt'

   x_ind                 = -1
   loc_of_interest       = -1.0, -1.0, -1.0
   kind_of_interest      = 'ANY'

   interp_test_dlon      = 10.0
   interp_test_dlat      = 10.0
   interp_test_dvert     = 10.0

   interp_test_lonrange  = 0.0, 120.0
   interp_test_latrange  = 0.0, 120.0
   interp_test_vertrange = 0.0, 100.0

   interp_test_dx        = -888888.0
   interp_test_dy        = -888888.0
   interp_test_dz        = -888888.0

   interp_test_xrange    = -888888.0, -888888.0
   interp_test_yrange    = -888888.0, -888888.0
   interp_test_zrange    = -888888.0, -888888.0

   interp_test_vertcoord = 'VERTISHEIGHT'

   test1thru             = 15
   run_tests             = -1
   verbose               = .FALSE.
   /


Item Type Description
num_ens integer Provided for future use. Must be 1. Ultimately, The number of ensemble members you would like to read in for testing.
single_file logical If .TRUE. all members are stored in a single restart file.
input_state_files(:) character(len=256) The name(s) of the NetCDF file(s) containing the model states, one per domain.
output_state_files(:) character(len=256) The name(s) of the output NetCDF file(s) for testing IO, one per domain.
all_metadata_file character(len=256) Test 6 produces an exhaustive list of metadata for EVERY element in the DART state vector. The metadata get written to this file name.
x_ind integer(i8) An integer index into the DART state vector. This will be used to test the metadata routines. Answers questions about location, what variable type is stored there, etc.
loc_of_interest real(r8), dimension(3) The lat/lon/level for a particular location. Used in Test 4, the single-point interpolation test. Indirectly tests the routine to find the closest gridpoint.
kind_of_interest character(len=32) Specifies the QUANTITY of the model state to use in Test 4.
interp_test_dlon real(r8) The distance (measured in degrees) on the longitude interpolation grid. Ignored if interpolating with cartesian coordinates. Used in Test 5.
interp_test_dlat real(r8) The distance (measured in degrees) on the latitude interpolation grid. Ignored if interpolating with cartesian coordinates. Used in Test 5.
interp_test_dvert real(r8) The distance (measured in interp_vertcoord) on the vertical interpolation grid. Ignored if interpolating with cartesian coordinates. Used in Test 5.
interp_test_lonrange real(r8) The range of y to be tested with model_interpolate, with spacing interp_test_dlon. Ignored if interpolating with cartesian coordinates. Used in Test 5.
interp_test_latrange real(r8) The range of y to be tested with model_interpolate, with spacing interp_test_dlat. Ignored if interpolating with cartesian coordinates. Used in Test 5.
interp_test_vertrange real(r8) The range in the vertical direction to be tested with model_interpolate, with spacing interp_test_dvert. Ignored if interpolating with cartesian coordinates. Used in Test 5.
interp_test_dx real(r8) The interval on the x axis of the interpolation grid. This is used in Test 5 for models with threed_cartesian coordinates.
interp_test_dy real(r8) The interval on the y axis of the interpolation grid. This is used in Test 5 for models with threed_cartesian coordinates.
interp_test_dz real(r8) The interval on the z axis of the interpolation grid. This is used in Test 5 for models with threed_cartesian coordinates.
interp_test_xrange real(r8) The range of x to be tested with model_interpolate in Test 5, with spacing interp_test_dx.
interp_test_yrange real(r8) The range of y to be tested with model_interpolate in Test 5, with spacing interp_test_dy.
interp_test_zrange real(r8) The range in the vertical direction to be tested with model_interpolate in Test 5, with spacing interp_test_dz.
interp_test_vertcoord character(len=32) Specifies the vertical coordinate system to use during the interpolation tests. Valid values are: 'VERTISHEIGHT','VERTISPRESSURE','VERTISLEVEL', and 'VERTISSCALEHEIGHT'.
test1thru integer If test1thru > 0, specifies the last test to be performed. All tests get performed sequentially. If test1thru < 0, run_tests is used to specify the tests to perform.
  1. Test static_init_model() by calling get_model_size().
  2. Reads and writes a restart file.
  3. Tests get_state_meta_data() for a single index into the DART state. Helps determine if the state vector is constructed correctly.
  4. Tests model_interpolate for a single point.
  5. Tests model_interpolate for a range of interpolation points.
  6. Long, expensive test to return the metadata for every element of the state vector. May be useful to decide on known locations for subsequent testing.
run_tests(:) integer Specifies a list of tests to be performed. Same test numbers as described in test1thru. There are some dependencies. Tests 4 and 5 require a valid model state - which is read by Test 2. If a required test is not specified, the required test is enabled and run. A value of -1 means that test1thru will be used.
verbose logical Print extra info about the model_mod_check run.

A more typical namelist for a single ensemble member for a model with an outer grid and a single nested grid is shown below.

&model_mod_check_nml
   input_state_files     = 'dart_vector1.nc','dart_vector2.nc'
   output_state_files    = 'check_me1.nc', 'check_me2.nc'
   all_metadata_file     = 'metadata.txt'
   verbose               = .TRUE.
   test1thru             = 5
   run_tests             = -1
   loc_of_interest       = 243.72386169, 52.78578186, 10.0
   x_ind                 = 12666739
   kind_of_interest      = 'QTY_POTENTIAL_TEMPERATURE'
   interp_test_lonrange  = 144.0, 326.0
   interp_test_dlon      = 1.0
   interp_test_latrange  = -5.0, 80.0
   interp_test_dlat      = 1.0
   interp_test_vertrange = 100.0, 11000.0
   interp_test_dvert     = 200.0
   interp_test_vertcoord = 'VERTISHEIGHT'
  /

[top]

OTHER MODULES USED

assimilation_code/location/threed_sphere/location_mod.f90
assimilation_code/location/utilities/default_location_mod.f90
assimilation_code/location/utilities/location_io_mod.f90
assimilation_code/modules/assimilation/adaptive_inflate_mod.f90
assimilation_code/modules/assimilation/assim_model_mod.f90
assimilation_code/modules/assimilation/assim_tools_mod.f90
assimilation_code/modules/assimilation/cov_cutoff_mod.f90
assimilation_code/modules/assimilation/filter_mod.f90
assimilation_code/modules/assimilation/obs_model_mod.f90
assimilation_code/modules/assimilation/quality_control_mod.f90
assimilation_code/modules/assimilation/reg_factor_mod.f90
assimilation_code/modules/assimilation/sampling_error_correction_mod.f90
assimilation_code/modules/assimilation/smoother_mod.f90
assimilation_code/modules/io/dart_time_io_mod.f90
assimilation_code/modules/io/direct_netcdf_mod.f90
assimilation_code/modules/io/io_filenames_mod.f90
assimilation_code/modules/io/state_structure_mod.f90
assimilation_code/modules/io/state_vector_io_mod.f90
assimilation_code/modules/observations/forward_operator_mod.f90
assimilation_code/modules/observations/obs_kind_mod.f90
assimilation_code/modules/observations/obs_sequence_mod.f90
assimilation_code/modules/utilities/distributed_state_mod.f90
assimilation_code/modules/utilities/ensemble_manager_mod.f90
assimilation_code/modules/utilities/netcdf_utilities_mod.f90
assimilation_code/modules/utilities/null_mpi_utilities_mod.f90
assimilation_code/modules/utilities/null_win_mod.f90
assimilation_code/modules/utilities/obs_impact_mod.f90
assimilation_code/modules/utilities/options_mod.f90
assimilation_code/modules/utilities/parse_args_mod.f90
assimilation_code/modules/utilities/random_seq_mod.f90
assimilation_code/modules/utilities/sort_mod.f90
assimilation_code/modules/utilities/time_manager_mod.f90
assimilation_code/modules/utilities/types_mod.f90
assimilation_code/modules/utilities/utilities_mod.f90
assimilation_code/programs/model_mod_check/model_mod_check.f90
models/your_model_here/model_mod.f90
models/model_mod_tools/test_interpolate_threed_sphere.f90
models/utilities/default_model_mod.f90
observations/forward_operators/obs_def_mod.f90
observations/forward_operators/obs_def_utilities_mod.f90

Items highlighted may change based on which model is being tested.

[top]

FILES

[top]

USAGE

Normal circumstances indicate that you are trying to put a new model into DART, so to be able to build and run model_mod_check, you will need to create a path_names_model_mod_check file with the following contents:

assimilation_code/location/threed_sphere/location_mod.f90
assimilation_code/location/utilities/default_location_mod.f90
assimilation_code/location/utilities/location_io_mod.f90
assimilation_code/modules/assimilation/adaptive_inflate_mod.f90
assimilation_code/modules/assimilation/assim_model_mod.f90
assimilation_code/modules/assimilation/assim_tools_mod.f90
assimilation_code/modules/assimilation/cov_cutoff_mod.f90
assimilation_code/modules/assimilation/filter_mod.f90
assimilation_code/modules/assimilation/obs_model_mod.f90
assimilation_code/modules/assimilation/quality_control_mod.f90
assimilation_code/modules/assimilation/reg_factor_mod.f90
assimilation_code/modules/assimilation/sampling_error_correction_mod.f90
assimilation_code/modules/assimilation/smoother_mod.f90
assimilation_code/modules/io/dart_time_io_mod.f90
assimilation_code/modules/io/direct_netcdf_mod.f90
assimilation_code/modules/io/io_filenames_mod.f90
assimilation_code/modules/io/state_structure_mod.f90
assimilation_code/modules/io/state_vector_io_mod.f90
assimilation_code/modules/observations/forward_operator_mod.f90
assimilation_code/modules/observations/obs_kind_mod.f90
assimilation_code/modules/observations/obs_sequence_mod.f90
assimilation_code/modules/utilities/distributed_state_mod.f90
assimilation_code/modules/utilities/ensemble_manager_mod.f90
assimilation_code/modules/utilities/netcdf_utilities_mod.f90
assimilation_code/modules/utilities/null_mpi_utilities_mod.f90
assimilation_code/modules/utilities/null_win_mod.f90
assimilation_code/modules/utilities/obs_impact_mod.f90
assimilation_code/modules/utilities/options_mod.f90
assimilation_code/modules/utilities/parse_args_mod.f90
assimilation_code/modules/utilities/random_seq_mod.f90
assimilation_code/modules/utilities/sort_mod.f90
assimilation_code/modules/utilities/time_manager_mod.f90
assimilation_code/modules/utilities/types_mod.f90
assimilation_code/modules/utilities/utilities_mod.f90
assimilation_code/programs/model_mod_check/model_mod_check.f90
models/your_model_here/model_mod.f90
models/model_mod_tools/test_interpolate_threed_sphere.f90
models/utilities/default_model_mod.f90
observations/forward_operators/obs_def_mod.f90
observations/forward_operators/obs_def_utilities_mod.f90
as well as a mkmf_model_mod_check script. You should be able to look at any other mkmf_xxxx script and figure out what to change. Once they exist:

[~/DART/models/yourmodel/work] % csh mkmf_model_mod_check
[~/DART/models/yourmodel/work] % make
[~/DART/models/yourmodel/work] % ./model_mod_check

Unlike other DART components, you are expected to modify model_mod_check.f90 to suit your needs as you develop your model_mod. The code is roughly divided into the following categories:

  1. Check the geometry information,
  2. Read/write a restart file,
  3. Check the construction of the state vector ... i.e. the metadata,
  4. Interpolate at a single point,
  5. Interpolate for a range of points.

1. Checking the Geometry Information:

The first block of code in model_mod_check is intended to test some of the most basic routines, especially static_init_model - which generally sets the geometry of the grid, the number of state variables and their shape, etc. Virtually everything requires knowledge of the grid and state vector, so this block should never be skipped.

2. Read/writing a restart file:

This directly reads and write state variables from the model netCDF file. This is a nice sanity check to make sure that the DART state vector is being read in properly.

3. Check the construction of the state vector:

It is critical to return the correct metadata for any given index into the DART state vector. This code block tests the two most common features of the metadata. As a bonus, this routine is also quite useful to determine EXACTLY where to place your first test observation. If you test precisely at a grid location, you should be able to really get a handle on debugging your model_interpolate() routine.

4. Test interpolation on a single point.

This tests your model's interpolation routine on a single point and returns the interpolated value. This requires that Test 2 works - it needs a valid model state with data.

5. Test interpolation on a range of values.

This tests your model's interpolation routine on a range of values returns the interpolated grid in check_me_interptest.nc and check_me_interptest.m which can be read in Matlab and used to visualize the result.

6. Exhaustively test the construction of the state vector.

This can be a long test, depending on the size of your state vector. This returns the same data as in Test 3 - for every element in the state vector. The metadata are written to a file specified by all_metadata_file and check_me_interptest.m which can be read in Matlab and used to visualize the result.

[top]

REFERENCES

[top]

ERROR CODES and CONDITIONS

There are no error conditions to check. This program is intended to demonstrate simple checks that will allow you to proceed with improving and testing the model_mod. There will be plenty of run-time errors, I suggest compiling your code with "bounds checking" turned on - at a minimum.

KNOWN BUGS

none at this time

[top]

FUTURE PLANS

Expanded instructions on how to add a model, and how to methodically test piece-by-piece.

[top]

Terms of Use

DART software - Copyright UCAR. This open source software is provided by UCAR, "as is", without charge, subject to all terms of use at http://www.image.ucar.edu/DAReS/DART/DART_download

Contact: Tim Hoar
Revision: $Revision: 11997 $
Source: $URL: https://svn-dares-dart.cgd.ucar.edu/DART/releases/Manhattan/assimilation_code/programs/model_mod_check/model_mod_check.html $
Change Date: $Date: 2017-10-17 16:26:56 -0600 (Tue, 17 Oct 2017) $
Change history:  try "svn log" or "svn diff"