PDAFlocal - the localization module of PDAF

Overview

PDAFlocal is optional functionality and was introduced with PDAF V2.3.

PDAFlocal provides an improved handling of the localization of state vectors as used in the loca filters LESTKF, LETKF, LNETF, LSEIK, and LKNETF. A leading motivation to implement PDAFlocal was that it leads to much better performance when the Python interface pyPDAF is used. However, it is also useful for native implementations in e.g. Fortran because it also provides additional functionality for e.g. vertical localization or weakly coupled assimilation.

Compared to not using PDAFlocal removes the need to implement the routines g2l_state_pdaf and 'l2g_state_pdaf' (see e.g. ImplementAnalysisLocal). For this, PDAFlocal provides three initialization routines and a set of alternative routines that are called for the analysis step.

Please note that if you have an existing implementation there is no direct need to switch to using PDAFlocal. If might be useful if you plan to use the additional functionality of PDAFlocal_set_increment_weights (see below).

PDAFlocal provides the initialization routines

• PDAFlocal_set indices : This routine stores the index array that describes the indices of the element of a local state vector in the global state vector. It is used for mapping between both the local and global vectors.
• PDAFlocal_set_increment_weights : This routine allows to provide PDAFlocal with a vector of local increment weights. When the local state vector is written back to the global state vector after a local analysis update, this allows to prescribe increment weights for each element of the local state vector. With this, one can e.g. implement a vertical localization if the local state vector is a column of the model grid.
• PDAFlocal_clear_increment_weights : This routine clears and deallocates increment weights which have been set before with PDAFlocal_set_increment_weights. After executing this routine, the weights are reset to one for all elements of the local state vector (i.e. unweighted increment)

These routines are called in init_dim_l_pdaf.

Initialization routines

Here we provide an overview of the interfaces of the 3 initialization routines.

PDAFlocal_set_indices

This routine is called in init_dim_l_pdaf when the dimension of the local state vector has to be determined. In addition to this dimension, the implementation guide recommends to also initialize the coordinates of the local analysis domain and the indices of the elements ofthe local state vector in the non-localized (global or domain-decomposed) state vector state_p. Thise indices are then used in the mapping from the global to the local state vector before the analysis and back after the analysis update in the local analysis loop. With PDAFlocal, one uses PDAFlocal_set_indices to pass the vector of indices to PDAFlocal, which later performs the mapping for all ensemble states.

The interface is:

SUBROUTINE PDAFlocal_set_indices(dim_l, map)

  INTEGER, INTENT(in) :: dim_l          ! Dimension of local state vector
  INTEGER, INTENT(in) :: map(dim_l)     ! Index array for mapping

Hints:

• For complex cases it might feel more intuitive to perform the mapping between the global and local state vectors in an explicit loop. However, this might lead to repeated computation of the indices and can hence be less efficient than computing the indices before.
• The initialization of hthe index vector map is analogous to a loop that directly performs the initialization of a local state vector. However, here only the indices are stored.

PDAFlocal_set_increment_weights

This routine provides the optional functionality to prescribe the assimilation increment of each element of the local state vector a different weight. This can be used, e.g. to implement a vertical localization in the case that the local state vector is a full vertical column of the model grid. In this case, one can make the increment weight depending on the height (or depth) of a grid point. Another application is to implement weakly-coupled assimilation in which the local state vector contains all variables, but only a subset of them is updated. This is achieved by givening those element that should not be updated the weight 0.

The routine can be called in init_dim_l_pdaf when the dimension of the local state vector has to be determined. In addition to this dimension, the implementation guide recommends to also initialize the coordinates of the local analysis domain and the indices of the elements ofthe local state vector in the non-localized (global or domain-decomposed) state vector state_p. In addition one can initialize a vector of increment weights and provide it to PDAFlocal by calling this routine.

The interface is:

SUBROUTINE PDAFlocal_set_increment_weights(dim_l, weights)

  INTEGER, INTENT(in) :: dim_l          ! Dimension of local state vector
  REAL, INTENT(in) :: weights(dim_l)    ! Weights array

Hints:

• The loop to initialize weights can be analogous to the loop that initializes the index array for mapping between the global and local state vector.
• It is possible to clear the weights array by calling PDAFlocal_clear_increment_weights (see below). Afterwards, unit weights are used.
• One can call this routine again to change the weights. One does not need to call PDAFlocal_clear_increment_weights before. In any case, one has to be careful to ensure that the size of weights remains consistent when the size of the local state vector changes.

PDAFlocal_clear_increment_weights

If PDAFlocal_set_increment_weights was used to prescribe the assimilation increment of each element of the local state vector a different weight this routine can be used to clear and deallocate the weights array. Afterwards a weight of one is used for the increment of all elements of the local state vector when initializing the global state vector after the analysis update in the local analysis loop.

The interface is without arguments:

SUBROUTINE PDAFlocal_clear_increment_weights()

Hint:

• After calling this routine, the local increment weights are one until one calls PDAFlocal_set_increment_weights with prescribed weights again.

Routines to perform the analysis step with OMI

Generally the routines to be called for the analysis step are analogous to those when not used PDAFlocal (e.g. in implementations before the release of PDAF V2.3). The difference is in the name and in the fact that the routines g2l_state_pdaf and l2g_state_pdaf are not present in the interface of the PDAFlocal routines.

Ensemble filters

for diagonal R matrix

Depending on the fully parallel are flexible implementation choose the routine

fully parallel	flexible
PDAFlocalomi_assimilate	PDAFlocalomi_put_state

For non-diagonal R matrix

See OMI_nondiagonal_observation_error_covariance_matrices for information on using non-diagonal R-matrices with OMI. The routines are only partly generic depending on the needed observation-specific routine:

filter	fully parallel	flexible
LESTKF LETKF LSEIK	PDAFlocalomi_assimilate_nondiagR	PDAFlocalomi_put_state_nondiagR
LNETF	PDAFlocalomi_assimilate_lnetf_nondiagR	PDAFlocalomi_put_state_lnetf_nondiagR
LLNETF	PDAFlocalomi_assimilate_lknetf_nondiagR	PDAFlocalomi_put_state_lknetf_nondiagR

3D-Var

Only En3DVar and hybrid 3D-Var use a local filter. The routines are the following.

for diagonal R matrix

Method	fully parallel	flexible
En3DVar	PDAFlocalomi_assimilate_en3dvar_lestkf	PDAFlocalomi_put_state_en3dvar_lestkf
Hyb3DVar	PDAFlocalomi_assimilate_hyb3dvar_lestkf	PDAFlocalomi_put_state_hyb3dvar_lestkf

for non-diagonal R matrix

See OMI_nondiagonal_observation_error_covariance_matrices for information on using non-diagonal R-matrices with OMI.

Method	fully parallel	flexible
En3DVar	PDAFlocalomi_assimilate_en3dvar_lestkf_nondiagR	PDAFlocalomi_put_state_en3dvar_lestkf_nondiagR
Hyb3DVar	PDAFlocalomi_assimilate_hyb3dvar_lestkf_nondiagR	PDAFlocalomi_put_state_hyb3dvar_lestkf_nondiagR

Routines to perform the analysis step not using OMI

Thesse are the routines for the older observation handling without OMI. There is no distinction with regard to diagonal or non-diagonal R because without PDAF-OMI one always needs to implement all observation routines on the user side.

Generally the routines to be called for the analysis step are analogous to those when not used PDAFlocal (e.g. in implementations before the release of PDAF V2.3). The difference is in the name and in the fact that the routines g2l_state_pdaf and l2g_state_pdaf are not present in the interface of the PDAFlocal routines.

Ensemble filters

Depending on the fully parallel are flexible implementation choose the routine

Filter	fully parallel	flexible
LESTKF	PDAFlocal_assimilate_lestkf	PDAFlocal_put_state_lestkf
LETKF	PDAFlocal_assimilate_letkf	PDAFlocal_put_state_letkf
LSEIK	PDAFlocal_assimilate_lseik	PDAFlocal_put_state_lseik
LNETF	PDAFlocal_assimilate_lnetf	PDAFlocal_put_state_lnetf
LKNETF	PDAFlocal_assimilate_lknetf	PDAFlocal_put_state_lknetf

3D-Var

Only En3DVar and hybrid 3D-Var use a local filter. The routines are the following.

for diagonal R matrix

Method	fully parallel	flexible
En3DVar	PDAFlocal_assimilate_en3dvar_lestkf	PDAFlocal_put_state_en3dvar_lestkf
Hyb3DVar	PDAFlocal_assimilate_hyb3dvar_lestkf	PDAFlocal_put_state_hyb3dvar_lestkf

Porting to PDAFlocal

If you like to port you existing code to using PDAFlocal, the steps are the following

1. Change the call to the analysis routine from PDAFomi_ to PDAFlocalomi_ (or PDAF_ to PDAFlocal_ for non-OMI implementations)
2. Remove g2l_state_pdaf and l2g_state_pdaf from the interface of the PDAFlocalomi_ routine
3. In init_dim_l_pdaf ensure that the index array for the mapping between the global and local state vectors is implemented (this should exist if the implementation followed the tutorial and templates where it is the array id_lstate_in_pstate).
4. Insert the call to PDAFlocal_set_indices.

After these steps it should be possible to recompile and run. If this was successful, one can safely remove the files holding g2l_state_pdaf and l2g_state_pdaf.