Wiki

Clone wiki

templates / Modifying the Coupling

This page explains the steps to either fully or partially manipulate the land-atmosphere coupling in WRF-LIS-CABLE. As code versions can change, this is a guide on which parts of the model code you will need to modify.

Full vs. Partial Uncoupled Simulations

What's the difference and why would it matter? This will depend on what research question you are hoping to answer.

For questions that are about general land-atmosphere coupling for a location one of the most common experimental designs is the Global Land-Atmosphere Coupling Experiment (GLACE) described by Koster et al. 2006. Without going into too much detail, essentially this involves running two sets of experiments:

  1. with full coupling involving two-way land-atmosphere interactions between the land surface model and the atmospheric model
  2. with the soil (moisture) states prescribed to effectively decouple the land surface model from the variability in the atmospheric model. Often mean soil moisture or anomalously dry or wet soil moisture conditions are used.

Once these sets of simulations are complete, one can compare the outputs, generally of the atmospheric variables, to evaluate the impact of the two-way interactions.

There are several metrics however for a GLACE-type experiment evaluating the coupling over Australia using WRF-LIS-CABLE see Hirsch et al. 2014.

For partial de-coupling, the concept is similar to the full de-coupling however instead of prescribing soil moisture everywhere, this is done selectively over a region and/or time of interest. For an example of how this has been done in WRF-LIS-CABLE please see Hirsch and King 2020.

Necessary Model Code Modifications for either Full or Partial Decoupling

Note that to see the implementation of these modifications please see the EFforcing branch of the nuwrf bitbucket.



TIP: when modifying code for the first time, it helps to see how its been done for other things. So if adding a new namelist option in the lis.config file, see how its been done for one of the CABLE switches that is already there.



Because codes do change between versions, this guide tries to limit use of line numbers and instead refers to the relevant subroutine and/or module.

Conceptually, what one does is add a namelist option so that you can toggle between running an uncoupled experiment prescribing all or part of the land surface state to running fully coupled. This requires one to add a new input file that contains the data that you want to prescribe.

Add a new namelist switch

To begin, let's add the new namelist option to the lis.config file. In this example, we'll create a new namelist switch called PRESCRIBE_SWITCH (this is very generic, but feel free to use something more creative!)

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_lsmMod.F90 - this file is where you can define new namelist options. In the section starting with type, public :: cable_type_dec add: logical :: PRESCRIBE_SWITCH

This will create a namelist option that is either .TRUE. or .FALSE.

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_readcrd.F90 - this file is where we tell LIS how to read the option from the lis.config file. Scrolling through this file you will see some code snippets of a similar format for each namelist option already available. Add the following within the cable_readcrd subroutine:

  call ESMF_ConfigFindLabel(LIS_config,                            &
       "CABLE PRESCRIBE_SWITCH:",rc=rc)
  do n=1,LIS_rc%nnest
     call ESMF_ConfigGetAttribute(LIS_config,                      &
          cable_struc(n)%prescribe_switch,rc=rc)
     call LIS_verify(rc,'CABLE PRESCRIBE_SWITCH: not defined')
  enddo

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_setup.F90 - this file allows you to transfer the namelist option from LIS to CABLE.

In the subroutine called cable_setup there is a section to set cable_user add the following: cable_user%PRESCRIBE_SWITCH = cable_struc(1)%PRESCRIBE_SWITCH

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/core/biogeophys/cable_common.F90 - where you add switches for CABLE via cable_user. Where other logical kbl_user_switches are defined add PRESCRIBE_SWITCH = .FALSE.

Add new input file

Now, we need to add a new input file to be read. This is a little more involved modification:

<path_to_model_code>/nuwrf/WRFV3/lis/core/LIS_readConfig.F90 - here one needs to add code in a few different sections. Line numbers are indicative of the general area. If you search for useslopemap you will get a general sense of the location.

Near line 240 add: 

allocate(LIS_rc%useprescribemap(LIS_rc%nnest))
allocate(LIS_rc%prescribefile(LIS_rc%nnest))

Near line 275 add:

LIS_rc%useprescribemap = "none"

The section starting near line 300 and ending near 470 add the following:

  call ESMF_ConfigFindLabel(LIS_config,"Prescribe land data file:",&
       rc=rc)
  do i=1,LIS_rc%nnest
     call ESMF_ConfigGetAttribute(LIS_config,LIS_rc%prescribefile(i),rc=rc)
     call LIS_verify(rc,"Prescribe land data file: not defined")
  enddo
 
  call ESMF_ConfigFindLabel(LIS_config,"Prescribe land data source:",rc=rc)
  do i=1,LIS_rc%nnest
     call ESMF_ConfigGetAttribute(LIS_config,LIS_rc%useprescribemap(i),rc=rc)
     call LIS_verify(rc,"Prescribe land data source: not defined")
  enddo

Towards the end near line 1000 where one set's the time, add:

LIS_rc%prescribetime = 0.0

<path_to_model_code>/nuwrf/WRFV3/lis/core/LIS_PRIV_rcMod.F90 - this is a module file for specifying independent variables in LIS. Again, search for useslopemap for an idea on what content needs to be added. In our example search for lisrcdec and in this section add:

character*50, allocatable  :: useprescribemap(:)
integer                    :: prescribeflag
character*100, allocatable :: prescribefile(:)
real*8                     :: prescribetime

<path_to_model_code>/nuwrf/WRFV3/lis/core/LIS_prescribeMod.F90

This is a new code file to read in the new input that you will need to create. I did this by copying the following file <path_to_model_code>/nuwrf/WRFV3/lis/core/LIS_vegDataMod.F90 and adapt what is done for LAI for the new input. To see an implementation see this code snippet noting that the word efforce is used instead of prescribe.

<path_to_model_code>/nuwrf/WRFV3/lis/core/LIS_histDataMod.F90 - where you can define new output variables. Here search for SWNET to understand placement of the following code additions:

public :: LIS_MOC_PRESCRIBE
integer :: LIS_MOC_PRESCRIBE   = -9999
    call ESMF_ConfigFindLabel(modelSpecConfig,"PRESCRIBE:",rc=rc)
    call get_moc_attributes(modelSpecConfig, LIS_histData(n)%head_lsm_list, &
         "PRESCRIBE",&
         "Quantity_used_to_prescribe",&
         "Quantity used to prescribe",rc)
    if ( rc == 1 ) then
       call register_dataEntry(LIS_MOC_LSM_COUNT,LIS_MOC_PRESCRIBE,&
            LIS_histData(n)%head_lsm_list,&
            n,1,ntiles,(/"-"/),1,(/"-"/),1,1,1,&
            model_patch=.true.)
    endif

<path_to_model_code>/nuwrf/WRFV3/lis/core/LIS_paramsMod.F90 - this is an interface to manage different sources of parameter datasets, it's where we can get LIS to read in some data that will be used to prescribe land surface state variables.

Before the implicit none after the module LIS_paramsMod statement we need to tell it to use the new module file that we have created (called LIS_prescribeMod.F90) by adding:

use LIS_prescribeMod

Then within the subroutine LIS_param_init() add:

call LIS_prescribe_setup

Then within the subroutine LIS_setDynparams(n) add:

call LIS_read_prescribe(n)

Then within the subroutine LIS_param_finalize() add:

call LIS_prescribe_finalize

Then within the subroutine LIS_param_reset() add:

call LIS_prescribe_reset()

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_driv_mod.F90 - here we need CABLE to understand what to do for water points, following the similar syntax in the subroutine called cable_watervals add:

cable_struc(n)%cable(t)%prescribe = LIS_rc%udef

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_dynsetup.F90 - to handle time-dependent variables like the quantity you want to prescribe!

Where modules that are used are specified add:

use LIS_prescribeMod

Then after the point where the LAI is updated add:

  if (LIS_rc%useprescribemap(n).ne."none") then
     do t = 1,LIS_rc%npatch(n,LIS_rc%lsm_index)
        tid = LIS_surface(n,LIS_rc%lsm_index)%tile(t)%tile_id
        cable_struc(n)%cable(t)%prescribe = LIS_prescribe(n)%prescribe(tid)
     enddo
  endif

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_setveg.F90 - worth checking but the above can also be added to where the vegetation parameters are read in a similar format

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_module.F90 - this is where new state variables are defined. Pending what you plan to prescribe you will need to add:

  • for a 1D variable
    real :: prescribe
    
  • or for a 2D variable (e.g. varies with soil layers)
    real :: prescribe(ms)
    

For the next steps the method depends on what you are trying to achieve.

How to do Fully Uncoupled / Prescribed Soil Moisture Everywhere

Run your fully coupled simulations first as per your domain and period of interest BUT before you run off to do that first make sure that you appropriately set the output interval of the land surface state variables of interest. Here one can take advantage of the WRF-LIS-CABLE model that the outputs and restarts of the atmospheric model (i.e. WRF) can be at a different frequency to that of the land surface model (i.e. CABLE).

When I have run GLACE-type experiments in the past, I write the restart files every timestep as the restart files generally save the state variables that would be of interest in this type of experiment. This does increase the i/o and hence the cost of the run so if in the end you only want to prescribe a climatological value you could instead save the restarts daily and use these to calculate your climatological estimate to be read in.

Now once you have these, modify the code according to the previous section. There may be different quantities that you may want to prescribe. Here, we will focus on prescribing soil moisture for all soil layers except the top layer, however, feel free to adapt for other variables should you so wish.

For prescribing sub-surface soil moisture at every time step you will need to have a file containing the soil moisture values that you want to prescribe. This guide assumes that you have soil moisture values for every time step of your simulation, with the time variable in the file coincident with your period of interest.

NOTE: It is advantageous to make sure that the format of your input file that is read by LIS_prescribeMod.F90 looks like the lis_input files created by LDT.

If you have this, then the way to prescribe the soil moisture is to allow CABLE to run through its calculations for the timestep, and then overwrite the soil moisture value with the value you want to prescribe. Then on the next timestep-integration, that soil moisture value will be used to determine the surface energy fluxes.

The place to do that is after the call to cbm in: <path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_driver.F90 in a subroutine called cable_2cable.

Here, <path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_driv_mod.F90

In the subroutine cable_2cable search for wb which is the variable name for soil moisture in CABLE. You will see a line:

cable_struc(n)%cable(t)%wb(:) = ssnow%wb(1,:)

After this line add the following:

if (LIS_rc%useprescribemap(n).ne."none")) then
  cable_struc(n)%cable(t)%wb(2:ms) = cable_struc(n)%cable(t)%prescribe(2:ms)
endif

How to do Partially Uncoupled / Prescribed Soil Moisture Over a region and/or portion of time

For this kind of set-up, the main modifications covered in the Necessary Model Code Modifications section previously are done with further modifications that will be detailed shortly!

Here, I'll describe how things were done for Hirsch and King 2020 noting that to produce the input file can be an involved process. Briefly:

Step 1: Run fully coupled simulations with WRF-LIS-CABLE for your location/period of interest

Step 2: Run the LAGRANTO back-trajectory codes on the WRF outputs to track air flow over the days preceding the event of interest. To do this, Dr. Malcolm King from the Monash CLEX node has experience in how to do this for WRF outputs.

Step 3: I wrote some code on how to take the back trajectories, along with their height and information on the WRF PBL height to create a forcing file. These are available at the zenodo record linked to the Hirsch and King 2020 paper. You will want to look at the python notebook called create_ef_forcing.ipynb. There are options to either create a file where the perturbation is only at one instance, or that when the parcel enters the boundary layer is starts and is sustained, as previous but not within the seed region (e.g. the area where the event of interest took place) to understand the impact of remote conditions.

Step 4: With this input file, as well as updating the lis.config files with the new setup run the simulations with the partial uncoupled set-up

In the Hirsch and King 2020 paper we aimed to prescribe the energy partitioning between the sensible and latent heat fluxes at the land surface for only parts of Australia in the days leading up to a heatwave event according to when an air mass entered the boundary layer (hence the need for the back-trajectory calculations). As most land-atmosphere coupled/uncoupled experiments aim to perturb the soil moisture to effectively change the energy partitioning this is a more direct way of doing so. To do this, one uses a quantity called the Evaporative Fraction (EF) which is calculated as:

EF = latent_heat_flux / (latent_heat_flux + sensible_heat_flux)

Since the partitioning is done within CABLE, further code modifications to CABLE itself are necessary:

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/core/biogeophys/cable_define_types.F90 we need to define a new variable that is the quantity we want to prescribe. In this implementation, it has been done within the CABLE canopy types so within this file search for TYPE canopy_type you can pick a variable like, for example, fev to understand the placement of the following code additions:

prescribe ! The data that you want to prescribe
ALLOCATE( var% prescribe(mp) )
DEALLOCATE( var% prescribe )

NOTE I have used the word PRESCRIBE for many of these modifications, this is intentional to avoid using a word that may already be in use in the model code.

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable/cable_driv_mod.F90 in the subroutine called cable_2canopy add the following in the style similar to other canopy% variables:

canopy%prescribe = cable_struc(n)%cable(t)%prescribe

Finally, to modify the energy partitioning this is applied within CABLE such that once this is applied, all the variables within the land surface model (particularly the soil moisture) are consistent with this change in the energy partitioning.

<path_to_model_code>/nuwrf/WRFV3/lis/surfacemodels/land/cable.2.3.4/core/biogeophys/cable_cbm.F90

where local variables are defined add:

    REAL, DIMENSION(mp) ::                                                      &
      ffess,ffes,ffhs,ffevc,ffhvc,ffevw,ffhvw

and later after the line CALL define_canopy(bal,rad,rough,air,met,dels,ssnow,soil,veg, canopy,climate) and before the call to the soil hydrology routines add:

   if (cable_user%prescribe_switch) then
      where (canopy%prescribe.gt.0.1)
         ffes = canopy%prescribe * (canopy%fes + canopy%fhs)
         ffhs = (1.0 - canopy%prescribe) * (canopy%fes + canopy%fhs)
         ! NB: fes = fess + fesp
         ! I assume that delta fes = delta fess and that there is no change in
         ! fesp
         ffess = canopy%prescribe * (canopy%fes + canopy%fhs)
         ffevc = canopy%prescribe * (canopy%fevc + canopy%fhv - canopy%fhvw)
         ffhvc = (1.0 - canopy%prescribe) * (canopy%fevc + canopy%fhv - canopy%fhvw)
         ffevw = canopy%prescribe * (canopy%fev + canopy%fhv)
         ffhvw = (1.0 - canopy%prescribe) * (canopy%fev + canopy%fhv)
         canopy%fess = ffess
         canopy%fes = ffes
         canopy%fhs = ffhs
         canopy%fevc = ffevc
         canopy%fevw = ffevw
         canopy%fhvw = ffhvw
         ! Calculate total latent heat flux:
         canopy%fev  = canopy%fevc + canopy%fevw
         canopy%fe = canopy%fev + canopy%fes
         ! Calculate total sensible heat flux:
         canopy%fhv  = ffhvc + canopy%fhvw
         canopy%fh = canopy%fhv + canopy%fhs
      endwhere
   endif

In this way the soil moisture will be updated to account for any changes in extraction due to the imposed changes in the latent heat flux.

Modifications to the templates and python setup codes!

It's worth updating the templates and codes for handling the new namelist option and. In particular, you probably want the prescribing added to the template_lis.config in the section where parameters are defined following:

Prescribe land data file:         ${PRESCRIBEFILE}
Prescribe land data source:       ${PRESCRIBESOURCE}

And later in the CABLE section of the namelist:

CABLE PRESCRIBE_SWITCH: ${PRESCRIBESWITCH} # .TRUE., .FALSE.

In the template_runwrf_raijin.deck where the input files are either linked or transferred to the run directory add the following:

    rm -f my_prescribe_data_file.nc
    ln -sf ${BDYDATA_PATH}/my_prescribe_data_file.nc .

For the python codes, you want to add the following in codes/prepare_LIS.py, search for RUNMODE and add near there the following:

    changes['PRESCRIBEFILE'] = 'none'
    changes['PRESCRIBESOURCE'] = 'none'
    changes['PRESCRIBESWITCH'] = '.FALSE.'

and also add in codes/prepare_WRF.py, search for RUNMODE and add near there the following:

    changes['PRESCRIBEFILE'] = 'my_prescribe_data_file.nc'
    changes['PRESCRIBESOURCE'] = 'EF'
    changes['PRESCRIBESWITCH'] = '.TRUE.'

Updated