Day 2 Adding a New Scheme to the CCPP
Day 2 Adding a New Scheme to the CCPP carson Mon, 03/04/2019 - 15:14Adding a new scheme to the CCPP
In this exercise, we will implement a new scheme in CCPP that calculates the Modified Red Flag Threat Index (RFTIM) following https://www.weather.gov/media/ama/rftimod_primer.pdf
We will also be adding the RFTIM to the diagnostic output, run a short forecast and display the results.
The exercise is designed to touch all the relevant bits and pieces for writing a new scheme, adding interstitial variables and a diagnostic variable, and outputting the diagnostic variable to disk.
Prepare for the GitHub workflow
Prepare for the GitHub workflow carson Mon, 03/04/2019 - 14:51Set up your development forks
In order to submit a pull request (PR) for the code changes at a later stage, it is necessary to create a GitHub user fork of the authoritative repositories for which changes are made. For this exercise, changes are made to the umbrella repository NEMSfv3gfs and the two submodules ccpp-physics and FV3. Creating a fork is a one-off task.
- Navigate to http://github.com/NCAR/NEMSfv3gfs, locate and click on the the "fork" button towards the top right of the page.
- Navigate to http://github.com/NCAR/FV3 locate and click on the the "fork" button
- Navigate to http://github.com/NCAR/ccpp-physics locate and click on the the "fork" button
If you already forked this repository in the past, you will see a message "You've already forked NEMSfv3gfs" and a link that you can click on to get to your fork (http://github.com/YOUR_GITHUB_USERNAME/NEMSfv3gfs). Otherwise, follow the instructions to create your user fork. Repeat for ccpp-physics and FV3.
Set up the directory structure for this hands-on session and checkout the source code.
mkdir ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training
cd ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training
git clone --branch=emc_training_march_2019 --recursive https://github.com/NCAR/NEMSfv3gfs develop
Note. Potential errors while checking out the submodule NCEPLIBS-pyprodutil (e.g., "gerrit: Name or service not known") can be ignored.
At this point, you should have a directory tree named "develop" which contains sub-directories for FV3, FMS, NEMS, ccpp/physics, ccpp/framework. Each of these contains the tag "emc_training_march_2019", as a starting point for these exercises.
Configure remotes and create branches
Next, configure the remote destinations (remotes) so that your development is pushed to your fork and not (inadvertently) to the GMTB repository at NCAR. Also, create branches in the three repositories where changes will be made:
cd develop
# Submodule ccpp-physics
cd ccpp/physics
git remote rename origin upstream
git remote add origin https://github.com/YOUR_GITHUB_USERNAME/ccpp-physics
git remote -v show
git remote update
git checkout -b my_rftim_implementation
cd ../../
# Submodule FV3
cd FV3
git remote rename origin upstream
git remote add origin https://github.com/YOUR_GITHUB_USERNAME/FV3
git remote -v show
git remote update
git checkout -b my_rftim_implementation
cd ../
# Umbrella repository NEMSfv3gfs
git remote rename origin upstream
git remote add origin https://github.com/YOUR_GITHUB_USERNAME/NEMSfv3gfs
git remote -v show
git remote update
git checkout -b my_rftim_implementation
Preparing to Calculate the Modified Red Flag Threat Index
Preparing to Calculate the Modified Red Flag Threat Index carson Mon, 03/04/2019 - 14:53The calculation of the RFTIM is based on the 2-m relative humidity in percent (rh2m) and on the 6-m wind speed in miles per hour (ws6m). It uses a geographically-dependent climatology of thresholds for the two quantities to return an integer value between 0 and 5 for each of the two variables. In our simplified version, we use the climatological values derived for KAMA (Rick Husband, Amarillo International Airport, Amarillo, TX), independent of the geographical location. This calculation could be made more sophisticated by using a spatially-varying climatology and by considering additional variables, such as snow cover.
Add the new variables
This calculation is done best in the existing scheme sfc_diag (ccpp/physics/physics/sfc_diag.f), because all the necessary quantities to derive relative humidity and wind speed are computed in the scheme. Before making changes to sfc_diag.f, note that the required input variables, rh2m and ws6m, are not provided by the model. The first task therefore is to add these two variables to FV3/gfsphysics/GFS_layer/GFS_typedefs.F90.
First let’s work with rh2m and ws6m. Both rh2m and ws6m are pure interstitial variables, i.e. they do not need to be persistent in memory from one time step to another. It is thus recommended to add them to the GFS_interstitial data type.
Variable rh2m:
standard_name : relative_humidity_at_2m
units : percent
rank : 1
type : real
kind : kind_phys
Variable ws6m:
standard_name : wind_speed_at_6m_in_miles_per_hour
units : mi h-1
rank : 1
type : real
kind : kind_phys
FYI: When variables have non-standard units, such as miles per hour, the variables standard name are decorated with information about the units, This is the case with variable ws6m above.
Each variable that is added to FV3/gfsphysics/GFS_layer/GFS_typedefs.F90 requires a minimum of 4 additions to the file:
- a metadata table entry
- a variable declaration inside one of the derived type definitions
- an allocation statement within a *_create method for the derived type
- an initialization statement.
For expediency, you can duplicate the entries pertaining to the variable zt1d for rh2m and ws6m, making the necessary adjustments to use the metadata above, that is, use the proper standard_name and units for each variable added. You will note that zt1d appears a 5th time, which is a print statement.
FYI: For some variables added to the GFS_interstitial_type or GFS_diag_type types, a few more additions are necessary for things like resetting their values in a *_reset method and adding a line to a *_print method.
Hint 1: Refer to presentation day2_5_Host-Side_Coding, section “Adding a variable to the host+CCPP”, “Case 1: the new variable is a member of a derived data type that already exists on the host model side”.
Hint 2: Using the standard_name and units exactly as listed above is critical for this development to work, These are the two aspects of metadata used in communication between the CCPP and the host model. Regarding the longname, you can use any description of your preference.
Hint 3: Use optional=F for all variable declarations in GFS_typedefs.F90
Next let’s work with rftim. Variable rftim also needs to be added to FV3/gfsphysics/GFS_layer/GFS_typedefs.F90 in 4 locations: metadata table, variable declaration, variable allocation, and variable initialization.
Variable rftim:
standard_name : modified_red_flag_threat_index
units : index
rank : 1
type : real
kind : kind_phys
Since variable rftim will be output by the model, it needs to be added to the GFS_diag_type data type. You need to modify the metadata table where the various GFS_diag_type variables are listed and add a line for rftim. For example, you could add a line at the end of the list, right after the line for atmosphere_momentum_diffusivity_for_mynnpbl. Note: if adding to the end of the metadata table list, be sure to leave a blank line starting with !! at the bottom of the table -- this is needed for parsing by the CCPP framework.
Next, you need to declare variable rftim. Add the following to the list of variables in section “Output - only in physics”. You can place this after the declaration for “MP quantities for 3D”
!--- Modified Read Flag Thread Index
real (kind=kind_phys), pointer :: rftim(:) => null() !<
You need to allocate the variable rftim in subroutine diag_create and then initialize it. You can enter the text below right before the line call Diag%rad_zero (Model)
! Modified Red Flag Threat Index calculation
allocate (Diag%rftim (IM))
Diag%rftim = clear_val
Calculate the New Variables
To calculate these variables in sfc_diag, you can use the following equations:
!> - For Modified Red Flag Threat Index calculation,
!! relative humidity 2m above ground in percent
rh2m(i) = 100.0*q2m(i)/qss
!> - Wind speed 6m above ground in miles per hour,
!! conversion factor from 10m to 6m is 0.886
!! (Lindley et al., 2011), and m/s to mi/h is 2.23694
ws6m(i) = sqrt(u10m(i)**2 + v10m(i)**2) * 0.886 * 2.23694
Edit file ccpp/physics/physics/sfc_diag.f. You will need to do the following modifications:
- Add the new variables (rh2m and ws6m) to the metadata table.
- Add the new variables as arguments to the subroutine.
- Declare the variables as real(kind=kind_phys), dimension(im), intent(out)
- Compute the variables using the equations provided above
- Print out the variable (optional)
To make these changes, you can follow what is done to variable u10m, since it also goes through the 5 steps above.
Hint 1: Don’t forget to add the new variables to the metadata table. You can copy/paste the entries from GFS_typdefs.F90 and remember to change the intent to OUT.
Hint 2: Don’t forget to pass your variables in and out of the subroutine. Arguments in the subroutine call must be in the same order as in the table of metadata. Also don’t forgot to declare them appropriately in the declaration section of the subroutine using intent(out) and dimension (im).
Calculate the Index and Add it to the Output
Calculate the Index and Add it to the Output carson Mon, 03/04/2019 - 14:57Add diagnostic variables for output
In order to add RFTIM to the output files, we need to add a new diagnostic variable to NEMSfv3gfs that stores the RFTIM. Because the NEMSfv3gfs diagnostic output doesn't work with integers, we need to declare rftim as a real array constituent of GFS_diag_type:
Variable rftim:
standard_name : modified_red_flag_threat_index
units : index
rank : 1
type : real
kind : kind_phys
Let's register this variable in the GFS diagnostics module. Open FV3/gfsphysics/GFS_layer/GFS_diagnostics.F90 and add the following blob to the end of routine GFS_externaldiag_populate (be sure to enter this after the last #endif preprocessor directive):
idx = idx + 1
ExtDiag(idx)%axes = 2
ExtDiag(idx)%name = 'rftim'
ExtDiag(idx)%desc = 'modified red flag threat index'
ExtDiag(idx)%unit = 'n/a'
ExtDiag(idx)%mod_name = 'gfs_sfc'
ExtDiag(idx)%intpl_method = 'bilinear'
allocate (ExtDiag(idx)%data(nblks))
do nb = 1,nblks
ExtDiag(idx)%data(nb)%var2 => IntDiag(nb)%rftim(:)
enddo
Calculate the index
Now we can begin the task of calculating the RFTIM (variable name rftim in the following). A skeleton CCPP scheme called rtfim, as well as a module containing the simplified climatology to map rh2m and ws6m to red flag threat indices, are provided for you to copy into your development codebase. Both files need to be copied into ccpp/physics/physics:
/scratch4/BMC/gmtb/emc_training_march_2019/handson_day2/templates/module_rftim_clima.F90
/scratch4/BMC/gmtb/emc_training_march_2019/handson_day2/templates/rftim.F90
Please take a look at module_rftim_clima.F90, but note that you do not need to make any modification in this file. You will see two functions in this file. The functions rfti_rh_clima and rfti_ws_clima can be called for arrays or scalars, and return integers. The interfaces are for relative humidity and wind speed are:
elemental function rfti_rh_clima(lat, lon, rh) result(rfti_rh)
real(kind_phys), intent(in) :: lat
real(kind_phys), intent(in) :: lon
real(kind_phys), intent(in) :: rh
integer :: rfti_rh
end function rfti_rh_clima
elemental function rfti_ws_clima(lat, lon, ws) result(rfti_ws)
real(kind_phys), intent(in) :: lat
real(kind_phys), intent(in) :: lon
real(kind_phys), intent(in) :: ws
integer :: rfti_ws
end function rfti_ws_clima
Next you need to edit rftim.F90 to compute RFTIM.
Note that in file rftim.F90 there are 3 subroutines wrapped in a module: rftim_init, rftim_finalize, and rftim_run. This is typical of CCPP schemes. No changes are needed in rftim_init and rftim_finalize, but you need to edit rftim_run and do the following:
-
Update the metadata table in the module
-
Declare local variable i
-
Calculate RFTIM
The argument list for subroutine rftim_run is:
subroutine rftim_run(im, lat, lon, islmsk, rh2, ws6, rftim, errmsg, errflg)
The metadata table needs to be populated with the same variables and in the same order as the subroutine call. The local variable names, standard names, intent, and optional columns are already pre-populated. Your job is to add the remaining metadata. Note that units are critical because they are used during the build system to compare the variables requested by the physics against those supplied by the host model. Therefore, the units needs to be exactly what the host model provides.
Variable im:
standard_name : horizontal_loop_extent
units : count
rank : 0
type : integer
kind :
Variable lat:
standard_name : latitude
units : radians
rank : 1
type : real
kind : kind_phys
Variable lon:
standard_name : longitude
units : radians
rank : 1
type : real
kind : kind_phys
Variable islmsk:
standard_name : sea_land_ice_mask
units : flag
rank : 1
type : integer
kind :
Variable rh2m:
standard_name : relative_humidity_at_2m
units : percent
rank : 1
type : real
kind : kind_phys
Variable ws6m:
standard_name : wind_speed_at_6m_in_miles_per_hour
units : mi h-1
rank : 1
type : real
kind : kind_phys
Variable rftim:
standard_name : modified_red_flag_threat_index
units : index
rank : 1
type : real
kind : kind_phys
Variable errmsg:
standard_name : ccpp_error_message
units : none
rank : 0
type : character
kind : len=*
Variable errflg:
standard_name : ccpp_error_flag
units : flag
rank : 0
type : integer
kind :
Hint: For unit “mi h-1” make sure you have the white space in the middle. If not, there will be a mismatch between the variable provided by the host model and what scheme rtfim requests.
Now that you are done with the metadata table, add to the subroutine
integer :: i
do i=1,im
rftim(i) = real(rfti_rh_clima(lat(i), lon(i), rh2(i)) + rfti_ws_clima(lat(i), lon(i), ws6(i)), kind_phys)
enddo
Scientific Documentation
Scientific Documentation carson Wed, 03/06/2019 - 12:58Prepare the Code
In this part of the training you will see how to add documentation to the new scheme, rtfim. The documentation uses Doxygen markup, so it can later be converted for display as an html or PDF document. The doxygen markup have been added to the template. At the end of this section, students are expected to
- finish the metadata table
- insert the calculation of RFTIM
- add the citation in library.bib file.
Subroutines rftim_init() and rftim_finalize() are just empty placeholders, we do not need to document them. Since rftim_run() is the entry point of the RFTIM scheme, Doxygen documentation should include:
1. In the first line of the Fortran file, add a brief one-sentence overview of the file purpose following “\file”, e.g.,
!>\file rftim.F90
!! This file contains the modified Red Flag Threat Index calculation.
The Doxygen code block begins with “!>”, and subsequent lines begin with “!!”.
2. Define a parent module using “\defgroup”, e.g., this module implements a particular scheme functionality. Also add a one-line description of the scheme with reference using “\cite” if applicable. Add the following lines after “end subroutine rftim_finalize”
!>\defgroup rftim_cal The Modified Red Flag Threat Index calculation
!!\brief This module calculates RFTIM using the 2-m relative humidity in percent
!! and 6m wind speed in miles per hour Lindley et al. (2011) \cite lindley_et_al_2011.
Doxygen can handle in-line paper citations and link to an automatically created bibliography page. A ccpp/physics/physics/docs/library.bib file in BibTex format is provided in the repository. Insert the following item in the file (order is not important):
@article {lindley_et_al_2011,
Author = {T.T. Lindley and J.D. Vitale and W.S. Burgett and M.-J. Beierle},
Title = {Proximity meteorological observations for wind-driven grassland wildfire starts on the southern High Plains},
Journal = {Electronic J. Severe Storms Meteor.},
Year = {2011},
Volume= {6},
Number = {1},
Pages = {1-27}}
3. In rftim_run(), which is the entry point to the scheme, further documentation will include:
An argument table with:
- Local variable name “local_name”
- CF-standard compliant “standard_name”: this is used as a key to communicate variables between the CCPP and the host model
- Short description “long_name”
- Units (format follows “unit exponent”, i.e., m2 s-2 for m2/s2)
- Rank (0 for scalar, 1 for 1-D array, 2 for 2-D array, etc)
Part of the argument table is intentionally left blank. Student are expected to fill in the remaining metadata.
The scheme general algorithm:
- Start with “\section sectionID XYZ General Algorithm”
- List in-line calculation steps by using “-#” markers
4. For those child subroutines that are not an entry point to the scheme, we can associate a subroutine or function with its parent module by using “\ingroup” to enable the auto-generation of a calling hierarchy within the scheme.
!>\ingroup rftim_cal
!!This function calculates the climatological relative humidity to integer.
!!...(detailed description could go here)...
!!\param lon real, longitude in radiance
!!\param rh real, relative humidity in percentage
!!\param rfti_rh integer, the climatological relative humidity
elemental function rfti_rh_clima(lat,lon,rh) result(rfti_rh)
Following a brief description of this function or subroutine, use “\param” to define each parameter with local name, provide a short description and list the unit (optional but highly recommended).
Run Doxygen
- Get set up to use Doxygen v1.8.10 by adding the following line into .cshrc file under your home directory:
alias doxygen /scratch4/BMC/gmtb/doxygen-1.8.10/bin/doxygen
- Source your .cshrc file.
- Edit file ccpplatex_dox and add the new files to the INPUT variable:
../calpreciptype.f90 \
### EMC_CCPP_TRAINING rftim
../rftim.F90 \
.. module_rftim_clima.F90 \
Here RFTIM scheme is inserted after calpreciptype.f90
Under ccpp/physics/physics/docs, run Doxygen by typing:
bash:
doxygen ccpplatex_dox 2>&1 | tee ccppdoc.log
tcsh:
doxygen ccpplatex_dox | & tee ccppdoc.log
- The generated html will be located under ccpp-physics/physics/doc/html/.
- In order to visualize the output, later you can transfer the html/ directory to your local machine. It can be viewed by pointing a browser to the index.html file.
Add the scheme to CCPP
Add the scheme to CCPP carson Mon, 03/04/2019 - 15:00Add the scheme to CCPP
Once the CCPP scheme rftim is complete, it is time to add the scheme rftim.F90 and its prerequisite module_rftim_clima.F90 to the CCPP prebuild configuration in ccpp/config/ccpp_prebuild_config.py.
Hint: Follow the procedure outlined in talk day2_3_prebuild.ppt section “Modifying CCPP prebuild config”. Note that scheme rftim will be added to the “slow physics” group. rftim.F90 goes in SCHEME_FILES and module_rftim_clima.F90 goes in SCHEME_FILES_DEPENDENCIES.
The next step is to create a new Suite Definition File and edit it so that it runs the scheme rftim just after the sfc_diag scheme:
cd ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/develop
cp ccpp/suites/suite_FV3_GFS_2017_updated_gfdlmp.xml ccpp/suites/suite_FV3_GFS_2017_rftim.xml
Edit ccpp/suites/suite_FV3_GFS_2017_rftim.xml and
- change the suite name to 'FV3_GFS_2017_rftim'
- add the scheme "rftim" after sfc_diag.
Copy the development directory to a test directory and change into that directory:
rsync -av ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/develop/ ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/compile/
cd ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/compile/tests
Compile the code (this example is for a dynamic build since the STATIC variable is not specified):
- for bash shell:
./compile.sh $PWD/../FV3 theia.intel “CCPP=Y HYBRID=N” '' NO NO 2>&1 | tee compile.log
- for t/csh shell:
./compile.sh $PWD/../FV3 theia.intel “CCPP=Y HYBRID=N” '' NO NO |& tee compile.log
Hint 1: After HYBRID=NO, there is a double quote plus two single quotes.
Hint 2: In file compile.log, look for the following statement: INFO: CCPP prebuild step completed successfully. If not found, some problem occurred in the CCPP prebuild stage. If found, prebuild was successful and you can look at other messages in compile.log pertaining to prebuild.
If the entire prebuild and build process was successful, this creates an executable 'fv3.exe' in the same directory.
- The compile.sh script takes up to 6 arguments.
PATHTR BUILD_TARGET [ MAKE_OPT [ BUILD_NR ] [ clean_before ] [ clean_after ]
- so, in the example command-line given above,
- PATHTR = $PWD/../FV3
- BUILD_TARGET = theia.intel
- MAKE_OPT = "CCPP=Y HYBRID=N"
- BUILD_NR (name of the fv3.exe executable) = null (so no additional string, just fv3.exe)
- clean_before = NO
- clean_after = NO
If successful, this creates an executable 'fv3.exe' in the same directory.
Run the model
Run the model carson Mon, 03/04/2019 - 15:00Run the model
Now copy the run directory template to your directory:
rsync -av /scratch4/BMC/gmtb/emc_training_march_2019/handson_day2/fv3_rftim/ ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/run/
cd ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/run/
Edit job_card to:
- Adjust the account name in line "# PBS -A ...".
- Adjust CODEDIR accordingly - this will also adjust the LD_LIBRARY_PATH later in the job_card script.
Take a quick look at the additional steps required for CCPP:
-
Link the Suite Definition File you previously created (develop/ccpp/suites/suite_FV3_GFS_2017_rftim.xml) to the run directory
ln -s ${YOUR_WORK_DIR}/develop/ccpp/suites/suite_FV3_GFS_2017_rftim.xml ccpp_suite.xml
-
Edit file input.nml and set it to use your Suite Definition File
ccpp_suite = 'ccpp_suite.xml'
Note that in this example the CCPP library will be loaded at runtime since the compilation command was for the dynamic CCPP build.
Before we run the case, we need to add the new diagnostic variable to the list in file diag_table. Edit diag_table and add the following line just below the 'qrain' line:
"gfs_sfc", "rftim", "rftim", "fv3_history2d", "all", .false., "none", 2
Submit the job script
qsub job_card && touch out && tail -f out
When finished, plot the new diagnostic variable rftim:
module load intel/15.1.133 ncl/6.5.0
ncl plot_rftim.ncl
Contribute your changes to GitHub
Contribute your changes to GitHub carson Mon, 03/04/2019 - 15:02Assuming you are happy with your development, commit the changes and push them to your GitHub fork.
cd ${YOUR_WORK_DIR}/NEMSfv3gfs_CCPP_training/develop
# ccpp-physics
cd ccpp/physics
git add .
git commit -m "Add scheme rftim and dependencies"
git push origin my_rftim_implementation
cd ../..
# FV3
cd FV3
git add .
git commit -m "Add variables and configure diagnostic output for RFTIM calculation"
git push origin my_rftim_implementation
cd ..
# NEMSfv3gfs
git add .
git commit -m "Add scheme rftim and dependencies to CCPP prebuild config, add SDF"
git push origin my_rftim_implementation
At this point, your changes have been committed to your forks of the repositories.
To propose these changes for inclusion into the GMTB authoritative repository on NCAR's GitHub organization, you need to create three pull requests (one for each of the repositories). For the example of ccpp-physics:
- Navigate to http://github.com/NCAR/ccpp-physics and click on 'new pull request'
- Select 'compare across forks', set base repository to 'NCAR/ccpp-physics', base to 'master', head repository to 'YOUR_GITHUB_USERNAME/ccpp-physics', compare to 'my_rftim_implementation'
- Add a descriptive title and short description in the text boxes
- This concludes the exercise
- If this were a real development for inclusion in master, you would now click on 'Create pull request' . PLEASE DO NOT DO THAT FOR THIS EXERCISE!