DART PROGRAM model_mod_check

DART project logo

Jump to DART Documentation Main Index
version information for this file:
$Id: model_mod_check.html 6380 2013-08-05 23:47:11Z nancy $

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 a DART ics file and writes out restart files, netCDF 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 
   input_file           = "dart_ics",        
   output_file          = "check_me",
   advance_time_present = .FALSE.,
   x_ind                = -1,
   loc_of_interest      = -1.0, -1.0, -1.0,
   kind_of_interest     = 'ANY'
   verbose              = .FALSE.,
   /


Item Type Description
input_file character(len=129) Name of a file containing DART initial conditions for the model.
output_file character(len=129) base portion of the name of the test netCDF file that will exercise the DART routines that create the True_State.nc, Prior_Diag.nc, and Posterior_Diag.nc files. The proper file extension will be added.
advance_time_present logical Flag to indicate if the DART restart file has the advance time present in the file.
x_ind integer 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. Tests the routine to find the closest gridpoint.
kind_of_interest character(len=32) Since there are usually many state variables on the same grid, it may be useful to restrict the search for a location of interest to include a particular kind of state variable.
verbose logical Print extra info about the model_mod_check run.


[top]

OTHER MODULES USED

assim_model_mod
types_mod
location_mod
model_mod
null_mpi_utilities_mod
obs_def_mod
obs_kind_mod
random_nr_mod
random_seq_mod
time_manager_mod
model_mod_check
utilities_mod
[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:

assim_model/assim_model_mod.f90
common/types_mod.f90
location/threed_sphere/location_mod.f90
models/your_model/model_mod.f90
mpi_utilities/null_mpi_utilities_mod.f90
obs_def/obs_def_mod.f90
obs_kind/obs_kind_mod.f90
random_nr/random_nr_mod.f90
random_seq/random_seq_mod.f90
time_manager/time_manager_mod.f90
utilities/model_mod_check.f90
utilities/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. write a trivial restart file,
  3. read either a restart file,
  4. check the netCDF routines used to create the diagnostic output files,
  5. check the metadata, and
  6. [optionally] run a test of the model interpolation routine.

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.

call initialize_utilities(progname='model_mod_check', output_flag=verbose)

call find_namelist_in_file("input.nml", "model_mod_check_nml", iunit)
read(iunit, nml = model_mod_check_nml, iostat = io)
call check_namelist_read(iunit, io, "model_mod_check_nml")

call static_init_model() ! Exercise the initialization process

call get_gridsize(numlons, numlats, numlevs)
write(*,'(''nlons, nlats, nlevs'',3(1x,i10))') numlons,numlats,numlevs

x_size = get_model_size()
write(*,'(''state vector has length'',i10)') x_size
allocate(statevector(x_size))

Writing a trivial restart file:

It's almost inconceivable that this cannot work, but sometimes it is just a nice sanity check to see how big a DART restart file for one ensemble member should be.

statevector = 1.0_r8;
model_time  = set_time(21600, 149446)   ! 06Z 4 March 2010 in the Gregorian calendar

iunit = open_restart_write('allones.ics')
call awrite_state_restart(model_time, statevector, iunit)
call close_restart(iunit)

Reading a restart file:

I generally take the allones.ics file into Matlab and stuff something interesting into the appropriate slots and write it out again so I can read it in for a more rigorous test of model_mod:aread_state_restart(). Generally, I'll try to fill up what I *intend* to be the locations of the first variable with some constant value or some nice field from Matlab. The second variable gets its own field, etc. That way, when it's read it in, I know what I'm supposed to find in the output. aread_state_restart() is used to read restart files as well as the intermediate files created by DART right before the model advance. Those files (created internally by DART) have an additional data record specifying the "advance_to_time" of the model state. Common practice has the dart_to_model program reading this intermediate file and communicating the advance_to_time to the model control mechanism and then parsing the DART state vector into a form compatible with the model. TIP: if you write the restart files in ASCII, you can trivially add another DART time record in the first line to test the case in which the "advance_to_time" is present.

iunit = open_restart_read(input_file)
if ( advance_time_present ) then
   call aread_state_restart(model_time, statevector, iunit, adv_to_time)
else
   call aread_state_restart(model_time, statevector, iunit)
endif

call close_restart(iunit)
call print_date( model_time,'model_mod_check:model date')
call print_time( model_time,'model_mod_check:model time')

Checking the diagnostic output netCDF routines:

This block happens after a call to aread_state_restart(), so, depending on what was in the restart file (presumably, once you get model_to_dart working, you have converted a real model state to a DART restart and are using that), you can fine-tune what gets put into the DART True_State.nc, Prior_Diag.nc, and Posterior_Diag.nc diagnostic files. Only one ensemble member is needed to test the routines (hence the hardcoded 1 in the test block).

state_meta(1) = 'restart test'
ncFileID = init_diag_output(trim(output_file),'just testing a restart', 1, state_meta)
call aoutput_diagnostics(ncFileID, model_time, statevector, 1)
call nc_check( finalize_diag_output(ncFileID), 'model_mod_check:main', 'finalize')

Check the metadata, and ...

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. The find_closest_gridpoint() routine is designed to ensure that your variable layout is as you expect. "closest" in this context is close in the horizontal only - all vertical levels will be reported.

call check_meta_data( x_ind )
call find_closest_gridpoint( loc_of_interest )

... [optionally] run a test of the model interpolation routine.

If you like, and you have a test_interpolate subroutine, this would be the place to run it.

call test_interpolate(statevector, test_pressure = 500.0_r8, &
                                   start_lon = 142.5_r8)
[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 2004 - 2013 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: 6380 $
Source: $URL: https://svn-dares-dart.cgd.ucar.edu/DART/releases/Lanai/utilities/model_mod_check.html $
Change Date: $Date: 2013-08-05 17:47:11 -0600 (Mon, 05 Aug 2013) $
Change history:  try "svn log" or "svn diff"