| 57 | | In our study, Nerger et al. (2015), we discussed that applying localization can lead to stability issues of the ENSRF. The filter performs a loop over all single observations. With localization the assimilation result depends on the order in which these observations are assimilated. This actually leads to the effect that the assimilation result at some grid point does not only depend on the observations within the localization radius **r**, but also on observations further away. This happens if these observation at larger distance are assimilated before the observations within the radius **r** and if they influence the state close to distance **r**. In this case the innovation for the state close to distance **r** is already changed when assimilating the observation within the radius and hence leads to a different result. This effect has implications on the parallelization since keeping the observation order constant over the full model domain leads to a partial serialization of the algorithm. In the implementation of the serial observation processing in PDAF we use the parallelization approach that does not guarantee the same order of the observations. Usually, the differences when changing the observation order are small, but the benefit is a better scaling because the serialization is avoided. Nonetheless, we generally recommend using LESTKF or LETKF, or their global variants ESTKF or ETKF, since they do not depend explicitly on the observation order, and they allow for non-diagonal observation error covariance matrices. However, the ENSRF/EAKF might have a good compute performance. |
| | 57 | In our study, Nerger et al. (2015), we discussed that applying localization can lead to stability issues of the ENSRF. The filter performs a loop over all single observations. With localization the assimilation result depends on the order in which these observations are assimilated. This actually leads to the effect that the assimilation result at some grid point does not only depend on the observations within the localization radius **r**, but also on observations further away. This happens if these observation at larger distance are assimilated before the observations within the radius **r** and if they influence the state close to distance **r**. In this case the innovation for the state close to distance **r** is already changed when assimilating the observation within the radius and hence leads to a different result. This effect has implications on the parallelization since keeping the observation order constant over the full model domain leads to a partial serialization of the algorithm. In the implementation of the serial observation processing in PDAF we use the parallelization approach that does not guarantee the same order of the observations. Usually, the differences when changing the observation order are small, but the benefit is a better scaling because the serialization is avoided. Nonetheless, we generally recommend using LESTKF or LETKF, or their global variants ESTKF or ETKF, since they do not depend explicitly on the observation order, and they allow for non-diagonal observation error covariance matrices. Nonetheless, the ENSRF/EAKF methods might have a good compute performance. |
| 112 | | Here |
| 113 | | |
| 114 | | To indicate user-supplied routines we use the prefix `U_`. In the tutorials in `tutorial/` and in the template directory `templates/` these routines exist without the prefix, but with the extension `_pdaf`. The files are named correspondingly. In the section titles below |
| 115 | | |
| 116 | | In the subroutine interfaces |
| | 113 | Here, all user-supplied routines are described that are required in the call to `PDAF_assimilate_ensrf`. For some of the generic routines, we link to the page on [ModifyModelforEnsembleIntegration modifying the model code for the ensemble integration]. |
| | 114 | |
| | 115 | To indicate user-supplied routines we use the prefix `U_`. In the tutorials in `tutorial/` and in the template directory `templates/` these routines exist without the prefix, but with the extension `_pdaf`. The files are named correspondingly. In the section titles below, we provide the name of the template file in parentheses. |
| | 116 | |
| | 117 | In the subroutine interfaces, some variables appear with the suffix `_p`. This suffix indicates that the variable is particular to a model sub-domain, if a domain decomposed model is used. Thus, the value(s) in the variable will be different for different model sub-domains. |
| 145 | | * We recommend to not only determine the size of the observation vector at this point. The routine is a good place to also already gather information about the corresponding indices of the state vector needed later to implement the observation operator. In addition, one can already prepare an array that holds the full observation vector, an array storing the coordinates of the observations and possible an array storing observation error variances (if the observation error covariance matrix is diagonal). The required arrays can be defined in a module like `mod_assimilation`. The information can be used later in `U_localize_covar_serial`. |
| 146 | | || **Note**: PDAF-OMI provides a structured approach for implementing the observation functionality. Because of this we generally recommend to ue the PDAF3 interface that uses PDAF-OMI. || |
| | 146 | * We recommend to not only determine the size of the observation vector at this point. The routine is a good place to also already gather information about the corresponding indices of the state vector needed later to implement the observation operator. In addition, one can already prepare an array that holds the full observation vector, an array storing the coordinates of the observations and possibly an array storing observation error variances (if the observation error covariance matrix is diagonal). The required arrays can be defined in a module like `mod_assimilation`. The information can be used later in `U_localize_covar_serial`. |
| | 147 | || **Note**: PDAF-OMI provides a structured approach for implementing the observation functionality. This approach simpified the implentation. Because of this, we generally recommend to use the [wiki:PDAF3_interface PDAF3 interface] that uses PDAF-OMI. || |
| | 148 | |
| 228 | | The general aspects of this routines have already been described on the [ModifyModelforEnsembleIntegration#U_prepoststepprepoststep_ens_pdaf.F90 page on modifying the model code for the ensemble integration] for the SEIK filter. For completeness, the description is repeated specifically for the EnKF: |
| 229 | | |
| 230 | | The interface of the routine is identical for all filters, but sizes can vary. Also, the particular operations that are performed in the routine can be specific for each filter algorithm. |
| 231 | | |
| 232 | | The interface for this routine is for the LEnKF |
| 233 | | {{{ |
| 234 | | SUBROUTINE |
| | 230 | The general aspects of this routines have already been described on the [ModifyModelforEnsembleIntegration#U_prepoststepprepoststep_ens_pdaf.F90 page on modifying the model code for the ensemble integration] for the SEIK filter. For completeness, the description is repeated specifically for the ENSRF: |
| | 231 | |
| | 232 | The interface of the routine is identical for all filters, but sizes can vary. Additionally, the specific operations performed in the routine can vary for each filter algorithm. |
| | 233 | |
| | 234 | The interface for this routine is for the ENSRF |
| | 235 | {{{ |
| | 236 | SUBROUTINE U_prepoststep(step, dim_p, dim_ens, dim_ens_p, dim_obs_p, & |
| 251 | | The routine `U_prepoststep` is called once at the beginning of the assimilation process. In addition, it is called during the assimilation cycles before the analysis step and after the ensemble transformation. The routine is called by all filter processes (that is `filterpe=1`). |
| 252 | | |
| 253 | | The routine provides for the user the full access to the ensemble of model states. Thus, user-controlled pre- and post-step operations can be performed. For example the forecast and the analysis states and ensemble covariance matrix can be analyzed, e.g. by computing the estimated variances. In addition, the estimates can be written to disk. |
| | 253 | The routine `U_prepoststep` is called once at the beginning of the assimilation process and subsequently during the assimilation cycles before the analysis step and after the ensemble transformation. The routine is called by all filter processes (that is `filterpe=1`). |
| | 254 | |
| | 255 | The routine provides for the user with full access to the ensemble of model states. Thus, user-controlled pre- and post-step operations can be performed. |
| | 256 | For example, the forecast and analysis states, as well as the ensemble covariance matrix, can be analyzed, e.g., by computing the estimated variances. |
| | 257 | In addition, the estimates can be written to disk. |
| 256 | | * If a user considers to perform adjustments to the estimates (e.g. |
| 257 | | * Only for the SEEK filter the state vector (`state_p`) is initialized. For all other filters, the array is allocated, but it can be used freely during the execution of `U_prepoststep`. |
| | 260 | * If a user considers to perform adjustments to the estimates (e.g., for balances), this routine is the right place for it. |
| | 261 | * The state vector (`state_p`) is allocated for all other filters, but it can be used freely during the execution of `U_prepoststep`. One cannot assume that the vector holds the correct state vector. |
| 259 | | * Apart from the size of the array `Uinv`, the interface is identical for all ensemble filters (SEIK/ETKF/EnKF/LSEIK/LETKF/LEnKF). In general it should be possible to use an identical pre/poststep routine for all these filters. |
| 260 | | * The interface through which `U_prepoststep` is called does not include the array of smoothed ensembles. In order to access the smoother ensemble array one has to set a pointer to it using a call to the routine `PDAF_get_smootherens` (see |
| | 263 | * Apart from the size of the array `Uinv`, the interface is identical for all ensemble filters. In general, it should be possible to use an identical pre/poststep routine for all these filters. |
| | 264 | * The interface through which `U_prepoststep` is called does not include the array of smoothed ensembles. In order to access the smoother ensemble array one has to set a pointer to it using a call to the routine `PDAF_get_smootherens` (see the page on [AuxiliaryRoutines auxiliary routines]) |
