CCPP and SCM Online Tutorial | Exercise 3: Add a new CCPP-compliant scheme

Adding a New CCPP-Compliant Scheme

In this exercise you will learn how  to implement a new scheme in CCPP. This will be a simple case in which all variables that need to be passed to and from the physics are already defined in the host model and available for the physics. You will create a new CCPP-compliant scheme to modify the sensible and latent heat fluxes by adding a tunable amount to each of them. You will then implement the new scheme, which we will call fix_sys_bias_sfc, in the SCM_RRFS_v1beta suite. Note that this procedure is provided only as an exercise to examine the impact of the surface fluxes on the prediction. While correcting the surface fluxes can help fix systematic biases,  there is no scientific basis for arbitrarily modifying the fluxes in this case.

This exercise has four steps:

  • Prepare a CCPP-compliant scheme 
  • Add the scheme to the Cmake list of files
  • Add the scheme to the Suite Definition File
  • Run the model and inspect the output

Step 1: Prepare a CCPP-compliant Scheme

The first step is to prepare a new CCPP-compliant scheme, which in this case will be called fix_sys_bias_sfc and will be responsible for “tuning” the surface fluxes. A skeleton CCPP-compliant scheme is provided for you to download into your codebase and continue developing. Two skeleton files (fix_sys_bias_sfc.F90 and fix_sys_bias_sfc.meta), which is the scheme itself and the metadata file with the description of the variables passed to and from the scheme, need to be copied from tutorial_files/ directory and into the CCPP physics source code directory:

cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/fix_sys_bias_sfc_skeleton.F90 $SCM_WORK/ccpp-scm/ccpp/physics/physics/fix_sys_bias_sfc.F90
cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/fix_sys_bias_sfc_skeleton.meta $SCM_WORK/ccpp-scm/ccpp/physics/physics/fix_sys_bias_sfc.meta

In the file fix_sys_bias_sfc.F90, you will need to add the following modifications to the subroutine  fix_sys_bias_sfc_run:

  • Add hflx_r and qflx_r which correspond to the sensible and latent heat flux, respectively,  to the argument list in subroutine fix_sys_bias_sfc_run with 

subroutine fix_sys_bias_sfc_run (im, con_cp, con_rd, con_hvap, p1, t1, hflx_r, qflx_r, errmsg, errflg)

  • Declare and allocate variables hflx_r and qflx_r

real(kind=kind_phys), intent(inout) :: hflx_r(:), qflx_r(:)
...
  • Add code to modify the surface fluxes by adding a tunable amount of each of them or, in case of negative fluxes, setting the fluxes to zero. Note that the modification factors are converted to kinematic units since the variables hflx_r and qflx_r are in kinematic units, not W m2.
do i=1, im
   rho = p1(i)/(con_rd*t1(i))
   !convert mod_factor to kinematic units and add
   hflx_r(i) = MAX(sens_mod_factor/(rho*con_cp) + hflx_r(i), 0.0)
   qflx_r(i) = MAX(lat_mod_factor/(rho*con_hvap) + qflx_r(i), 0.0)

end do

Where con_rd (287.05 J kg-1 K-1; ideal gas constant for dry air), con_cp (1004.6 J kg-1 K-1; specific heat of dry air at constant pressure), and con_hvap (2.5106 J kg-1; latent heat of evaporation and sublimation) are constants defined in the host model ($SCM_WORK/ccpp-scm/scm/src/scm_physical_constants.F90), and  sens_mod_factor and lat_mod_factor are arbitrary constants defined locally:

   real(kind=kind_phys), parameter :: sens_mod_factor = 0 !W m-2              
   real(kind=kind_phys), parameter :: lat_mod_factor = 40 !W m-2              

Given the value of the constants, this results in an increase of 40 W m-2 to the surface latent heat flux, while the sensible heat flux is not altered.                                                                                                            

After saving your changes, edit the corresponding metadata file: fix_sys_bias_sfc.meta to add variables hflx_r and qflx_r after the entry for variable t1. The order of the variables on the table should match the order of the arguments in the subroutine call.

The metadata file needs to be populated with the same variables as the subroutine fix_sys_bias_sfc_run call:

[t1]
  standard_name = air_temperature_at_surface_adjacent_layer
  long_name = mean temperature at lowest model layer
  units = K
  dimensions = (horizontal_loop_extent)
  type = real
  kind = kind_phys
  intent = in
[qflx_r]
  standard_name = surface_upward_specific_humidity_flux
  long_name = kinematic surface upward latent heat flux
  units = kg kg-1 m s-1
  dimensions = (horizontal_loop_extent)
  type = real
  kind = kind_phys
  intent = inout
[hflx_r]
  standard_name = surface_upward_temperature_flux
  long_name = kinematic surface upward sensible heat flux
  units = K m s-1
  dimensions = (horizontal_loop_extent)
  type = real
  kind = kind_phys
  intent = inout

[errmsg]

Note: If you don’t want to type, you can copy over completed versions of these files:

cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/fix_sys_bias_sfc.F90 $SCM_WORK/ccpp-scm/ccpp/physics/physics
cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/fix_sys_bias_sfc.meta $SCM_WORK/ccpp-scm/ccpp/physics/physics

Step 2: Add the Scheme To the Cmake List of Files

Now that the CCPP-compliant scheme fix_sys_bias_sfc.F90 and the accompanying fix_sys_bias_sfc.meta are complete, it is the time to add the scheme to the CCPP prebuild configuration. 

Edit file $SCM_WORK/ccpp-scm/ccpp/config/ccpp_prebuild_config.py and include scheme fix_sys_bias_sfc.F90 in the SCHEME_FILES section:

SCHEME_FILES = [
...
    'ccpp/physics/physics/mynnsfc_wrapper.F90'                     ,
    'ccpp/physics/physics/fix_sys_bias_sfc.F90'                    ,
...  

Note: If you don’t want to type, you can copy over completed version:

cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/ccpp_prebuild_config.py $SCM_WORK/ccpp-scm/ccpp/config

Step 3: Add the Scheme To the Suite Definition File

The next step is to create a new Suite Definition File (SDF) and edit it so that it lists the scheme fix_sys_bias_sfc just after the GFS_surface_generic_post interstitial scheme:

First, create a new SDF from an existing SDF.

cd $SCM_WORK/ccpp-scm/ccpp/suites/
cp suite_SCM_RRFS_v1beta_sas.xml suite_SCM_RRFS_v1beta_sas_sfcmod.xml

Next, the new SDF, suite_SCM_RRFS_v1beta_sas_sfcmod.xml, needs the following modifications:

  • Change the suite name to “SCM_RRFS_v1beta_sas_sfcmod
  • Add the scheme “fix_sys_bias_sfc” after GFS_surface_generic_post

Note: If you don’t want to type, you can copy the completed SDF and modified supported  suite_list:

cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/suite_SCM_RRFS_v1beta_sas_sfcmod.xml $SCM_WORK/ccpp-scm/ccpp/suites
cp $SCM_WORK/ccpp-scm/tutorial_files/create_new_suite/suite_info.py $SCM_WORK/ccpp-scm/scm/src/suite_info.py

Step 4: Run the Model and Inspect the Output

After creating the new SDF, follow these commands to rebuild the SCM executable:

cd $SCM_WORK/ccpp-scm/scm/bin
cmake ../src -DCCPP_SUITES=ALL
make -j4

The following two commands will run the SCM for the TWP-ICE case using the modified suite containing the fix_sys_bias_sfc scheme, respectively.

Note: if you are using a Docker container, add -d at the end of the command

./run_scm.py -c twpice -s SCM_RRFS_v1beta_sas_sfcmod

Upon completion of the run, the following subdirectory will be created containing a log file and results in NetCDF format in file output.nc

$SCM_WORK/ccpp-scm/scm/run/output_twpice_SCM_RRFS_v1beta_sas_sfcmod

Now you will create plots to visualize the results. Edit $SCM_WORK/ccpp-scm/scm/etc/scripts/plot_configs/twpice_scm_tutorial.ini to define the results that will be displayed: the RRFS_v1beta_sas results (which use the saSAS scheme), and the new results obtained with the surface fluxes modification:

scm_datasets = output_twpice_SCM_RRFS_v1beta_sas/output.nc,
output_twpice_SCM_RRFS_v1beta_sas_sfcmod/output.nc
scm_datasets_labels = RRFSv1beta_sas,RRFSv1beta_sas_sfcmod
plot_dir = plots_twpice_sas_sfcmod/

Note: If you don’t want to type, you can copy over a pre-generated plot configuration file:

cp $SCM_WORK/ccpp-scm/tutorial_files/add_new_scheme/twpice_scm_tutorial.ini $SCM_WORK/ccpp-scm/scm/etc/scripts/plot_configs

Change into the following directory:

cd $SCM_WORK/ccpp-scm/scm/run/

Create the plots for the new suite with the command below.

Note: if you are using a Docker container, add -d at the end of the command

./scm_analysis.py twpice_scm_tutorial.ini

The new plots are under: 

$SCM_WORK/ccpp-scm/scm/run/plots_twpice_R_sas_sfcmod/comp/active/

Inspect the SCM results

Q: Examine the time series of latent and sensible heat fluxes in files time_series_lhf.pdf and time_series_shf.pdf, respectively. Are the modified latent and sensible heat fluxes as expected?

 

A: Compared to RRFSv1beta_sas, the surface latent heat flux in RRFSv1beta_sas_sfcmod is increased as expected, and the surface sensible heat flux is adjusted to a lower value.

 

Q: Examine the vertical profile of convective tendencies in file profiles_mean_dT_dt_conv.pdf. How  does the increased latent heat flux affect the atmospheric convection?

The evaporation process associated with the increased surface latent heat flux makes more moisture available in the atmosphere, which triggers stronger convection.

A: This exercise shows us that if some properties of the surface energy budget system change in such a way that the surface latent heat flux goes up, the system will rebalance to compensate for this change. In this case, we note that the sensible heat flux was reduced after the second day of forecast, even though we did not modify it directly. The sensible heat flux depends on the difference in temperature between the sea surface temperature and the air, the wind speed, and the surface heat exchange coefficient. Given that sea surface temperature changes very slowly, the changes in sensible heat flux can be attributed to lowering of the air temperature or increase in the wind speed. Those are non-linear effects deriving from the increase in latent heat flux.