= Online Mode: Initializing PDAF = {{{ #!html

Online Mode: Implementation Guide

  1. Main page
  2. Adapting the parallelization
  3. Initializing PDAF
  4. Modifications for ensemble integration
  5. Implementing the analysis step
  6. Memory and timing information
}}} [[PageOutline(2-3,Contents of this page)]] == Overview == The PDAF release provides example code for the online mode in `tutorial/online_2D_parallelmodel` and `tutorial/online_2D_serialmodel`. We refer to these codes to use it as a basis. The initialization of PDAF is done in the subroutine `init_pdaf`. This routine is called from the model code and sets all variables required for the initialization. Then the subroutine `PDAF_init` is called, which performs the actual internal initialization of PDAF. Thus, it sets the internal parameters for the data assimilation, chooses the data assimilation method, allocates the internal ensemble array, and initializes the ensemble. By using `init_pdaf`, only a single additional subroutine call has to be inserted to the model source code for the initialization. The file `templates/online/init_pdaf.F90` provides a commented template for this routine, which can be used as the basis of the implementation. `PDAF_init` itself calls a user-supplied call-back routine to initialize the ensemble of model states. In the example and templates, this routine can be found in `init_ens_pdaf.F90`. To initialize the ensemble forecasting, also the subroutine `PDAF_init_forecast` is called in `init_pdaf`. == Inserting `init_pdaf` into the model code == The right place to insert the routine `init_pdaf` into the model code is in between the initialization part of the model and the time stepping loop. At this point, the regular model initialization is completed, which allows PDAF to initialize the ensemble of model states. == Using `init_pdaf` == In `init_pdaf` a number of variables are defined that are used in the call to `PDAF_init` as described below. There are also a few variables that are initialized in `init_pdaf` for later use in user-supplied call-back routies. For the tutorial example, these variables are described below in the section '[#Othervariablesfortheassimilation Other variables for the assimilation]'. The template file `templates/online/init_pdaf.F90` provide some more functionality which can be adapted. The example implementation and the template code allow to specify all options at run time using a command line parser. These options are specified as the combination `-VARIABLE VALUE`. This method provides a convenient way to define an experiment and could also be used for other models. The parser module is provided by the file `tutorial/offline_2D_parallel/parser_mpi.F90`. An alterantive could be to use a configuration file. == Calling `PDAF_init` == In the tutorial codes and the template, the call to `PDAF_init` is fully implemented. Here, we provide an overview of the arguments that are set in the call to `PDAF_init`. The call to `PDAF_init` has the following structure: {{{ CALL PDAF_init(filtertype, subtype, step_null, & filter_param_i, length_filter_param_i, & filter_param_r, length_filter_param_r, & COMM_model, COMM_filter, COMM_couple, & task_id, n_modeltasks, filterpe, & init_ens_offline, screen, status_pdaf) }}} The required arguments are described below. In the list, we mark those variables bold, which one might like to change, like the type of the DA method. The other variables are required, but usually not changed by the user. * **filtertype**:[[BR]] An integer defining the type of the DA method. (See the [#Noteonavailableoptions Note on Available Options]) * **subtype**:[[BR]] An integer defining the sub-type of the filter algorithm. (See the [#Noteonavailableoptions Note on Available Options]) * `step_null`: An integer defining the initial time step. For some cases it can use useful to set `step_null` larger to 0. * **filter_param_i**:[[BR]] Integer array collecting options for PDAF. The first two variables are mandatory and equal for all filters. Further variables are optional (See the [#Noteonavailableoptions Note on Available Options]). The mandatory variables are in the following order: 1. The size of the local state vector for the current process. (see [#Definingthestatevector information on defining the state vector]) 1. The ensemble size for all ensemble-based filters * **length_filter_param_i**:[[BR]] An integer defining the length of the array `filter_param_i`. The entries in the array are parsed up to this index. * **filter_param_r**:[[BR]] Array collecting real-valued options for PDAF. The first value is mandatory and equal for all filters. Further variables are optional (See the [#Noteonavailableoptions Note on Available Options]). The mandatory variable is: 1. The value of the forgetting factor controlling inflation (required to be larger than zero) * **length_filter_param_r**:[[BR]] An integer defining the length of the array `filter_param_r`. The entries in the array are parsed up to this index. * `COMM_model`:[[BR]] The communicator variable `COMM_model` as initialized by `init_parallel_pdaf`. If the model-communicator is named differently in the actual program, the name has to be adapted * `COMM_filter`:[[BR]] The communicator variable `COMM_filter` as initialized by `init_parallel_pdaf`. * `COMM_couple`:[[BR]] The communicator variable `COMM_couple` as initialized by `init_parallel_pdaf`. * `task_id`:[[BR]] The index of the model tasks as initialized by `init_parallel_pdaf`. Always 1 for the offline mode * `n_modeltasks`:[[BR]] The number of model tasks as defined before the call to `init_parallel_pdaf`. * `filterpe`:[[BR]] The flag showing if a process belongs to `COMM_filter` as initialized by `init_parallel_pdaf`. * **init_ens_pdaf**:[[BR]] The name of the user-supplied routine that is called by `PDAF_init` to initialize the ensemble of model states. (See below: '[#User-suppliedroutineinit_ens_pdaf User-supplied routine init_ens_pdaf]' * `screen`:[[BR]] An integer defining whether information output is written to the screen (i.e. standard output). The following choices are available: * 0: quite mode - no information is displayed. * 1: Display standard information (recommended) * 2: as 1 plus display of timing information during the assimilation process * `status_pdaf`:[[BR]] An integer used as status flag of PDAF. If `status_pdaf` is zero upon exit from `PDAF_init`, the initialization was successful. An error occurred for non-zero values. (The error codes are documented in the routine PDAF_init.) PDAF uses two arrays **filter_param_i** and **filter_param_r** to respectively specify integer and real-valued options for PDAF. As described above, 2 integer values (state vector size, ensemble size) and 1 real value (forgetting factor) are mandatory. Additional options can be set by specifying a larger array and setting the corresponding size value (`length_filter_param_i`, `length_filter_param_r`). However, with PDAF V3.0 it can be more convenient to use the subroutines `PDAF_set_iparam` and `PDAF_set_rparam`, which are explained further below. We recommended to check the value of `status_pdaf` in the program after PDAF_init (and potentially `PDAF_set_iparam and `PDAF_set_rparam`) are executed. Only if its value is 0, the initialization was successful. === Note on available options === A **list of available values of `filtertype`** as well as an **overview of available integer and real-valued options** for each DA method can be found on the page [wiki:AvailableOptionsforInitPDAF Available options for the different DA methods]. The **available options for a specific DA method** can also be displayed by running the assimilation program for the selected DA method setting `subtype = -1`. (In the tutorial and template codes one can set `-subtype -1` on the command line). Generally, available options and valid settings are also listed in `mod_assimilation.F90` of the tutorials and template codes, but this might not be up-to-date in all cases. == Other variables for the assimilation == The routine `init_pdaf` in the example also initializes several variables that are not used to call `PDAF_init`. These variables control some functionality of the user-supplied routines for the data assimilation system and are shared with these routines through the Fortran module `mod_assimilation`. These variables are for example: * `delt_obs`: An integer specifying the number of time steps between two analysis steps (used later in `next_observation_pdaf`) * `cradius`: Localization cut-off radius in grid points for the observation domain * `sradius`: support radius, if observation errors are weighted (i.e. `locweight>0`) * `locweight`: Type of localizing weight (see further below) The latter three variables are localization parameters that are used later in the subroutines called by `init_dim_obs_l_pdafomi`. We recommend to define such configuration variables at this central position, so that all configuration is done at one place. Of course, their definition should be adapted to the particular model used. The setting of `locweight` influences the weight function for the localization. With the PDAF3 interface (and generally with PDAF-OMI), the choices are standardized as follows ||= '''locweight''' =||= '''0''' =||= '''1''' =||= '''2''' =||= '''3''' =||= '''4''' =|| ||= function =|| unit weight || exponential || 5-th order[[BR]]polynomial || 5-th order[[BR]]polynomial || 5-th order[[BR]]polynomial || ||= regulation =|| - || - || - || regulation using[[BR]]mean variance || regulation using variance[[BR]]of single observation point || ||= '''cradius''' =|||||||||||| weight=0 if distance > cradius || ||= '''sradius''' =|| no impact || weight = exp(-d / sradius) |||||||| weight = 0 if d >= sradius[[BR]] else[[BR]] weight = f(sradius, distance) || Here, 'regulation' refers to the regulated localization introduced in Nerger, L., Janjić, T., Schröter, J., Hiller, W. (2012). A regulated localization scheme for ensemble-based Kalman filters. Quarterly Journal of the Royal Meteorological Society, 138, 802-812. ​[https://doi.org/10.1002/qj.945 doi:10.1002/qj.945]. Apart from the generic variables for localization, we also specify variables that are specific for each observation type, for example in the tutorial code, we specify * `assim_A`: Flag whether to assimialtion observations of type A * `rms_obs_A`: Assumed observation error standard deviation of observation type A == User-supplied routine `init_ens_pdaf` == The user-supplied routine that we named `init_ens_pdaf` here, is the call-back routine that is called by PDAF through the defined interface described below. The routine is called by all MPI processes that compute the filter analysis step (i.e. those for which 'filterpe' is set to true. In the standard configuration of `init_parallel_pdaf` these are all processes of the first model task, i.e. task_id=1.) `init_ens_pdaf` is only called by `PDAF_init` if no error occurred before; thus the status flag is zero. The interface details can be looked up in the template and tutorial codes. It is the following: {{{ SUBROUTINE init_ens_pdaf(filtertype, dim_p, dim_ens, & state_p, Ainv, ens_p, flag) }}} with * `filtertype`, `integer, intent(in)`:[[BR]]The type of filter algorithm as given in the call to `PDAF_init` * `dim_p`, `integer, intent(in)`:[[BR]] The size of the state dimension for the calling process as specified in the call to `PDAF_init` * `dim_ens`, `integer, intent(in)`:[[BR]]The size of the ensemble * `state_p`, `real, dimension(dim_p), intent(inout)`:[[BR]]Array for the local model state of the calling process (can be used freely for ensemble-based methods) * `Ainv`, `real, dimension(dim_ens-1, dim_ens-1), intent(inout)`:[[BR]]A possible weight matrix (Not relevant for ensemble-based methods) * `ens_p`, `real, dimension(dim_p, dim_ens), intent(inout)`:[[BR]] The ensemble array, which has to be filled with the ensemble of model states. * `flag`, `integer, intent(inout)`:[[BR]]Status flag for PDAF. It is 0 upon entry and can be set by in the user-supplied routine, depending on the success of the ensemble initialization. Preferably, values above 102 should be used for failures to avoid conflicts with the error codes defined within PDAF_init. === Defining the state vector === The ensemble initialization routine `init_ens_pdaf` is the first location at which the user has to fill a state vector (or array of state vectors). A state vector is the collection of all model fields that are handled in the analysis step of the assimilation procedure into a single vector. Usually one concatenates the different model fields as complete fields. Thus, the vector could contain a full 3-dimensional temperature field, followed by the salinity field (in case of an ocean model), and then followed by the 3 fields of the velocity components. The logical definition of the state vector will also be utilized in several other user-supplied routines. E.g. in routines that fill model fields from a state vector or in the routine providing the observation operator. The actual setup of the state vector should be done in `init_pdaf`. The tutorial example `tutorial/online_2D_serialmodel_2fields` demonstrates a possible setup of the state vector for 2 fields. Here, one defines the number of fields, the dimension of each included field as well as the offset of each field in the state vector. === Initialization for ensemble-based filters === For the ensemble-based filters and the ensemble/hybrid 3D-Var methods, only the array `ens_p` needs to be initialized by the ensemble of model states. If a parallel model with domain decomposition is used, the full ensemble for the local sub-domain has to be initialized. The arrays `state_p` and `Ainv` are allocated to their correct sizes because they are used during the assimilation cycles. They are not yet initialized and it is allowed to use these arrays in the initialization. An exception from this is EnKF for which `Ainv` is allocated only with size (1,1), because `Ainv` is not used for EnKF. == Setting additional options == In PDAF V3.0 we added the possibility to set options for PDAF after the call to `PDAF_init`. For these there are the subroutines {{{ PDAF_set_iparam(id, value, status_pdaf) }}} to set interger parameters and {{{ PDAF_set_rparam(id, value, status_pdaf) }}} to set REAL (floating point) parameters. The arguments are * `id`:[[BR]] The index value of a parameter * `value`:[[BR]] The value of the parameter with index `id` * `status_pdaf`:[[BR]] Status flag for PDAF. Both routines increment in the input value. The increment is 0 for no error (this allows to check `flag` once after all calls to `PDAF_init`, `PDAF_set_iparam`, and `PDAF_set_rparam`.) The tutorial code uses these routines for a few settings while the template code include an extended set of calls specific for different DA methods. An overview of available integer and real-valued options for each DA method can be found on the page [wiki:AvailableOptionsforInitPDAF Available options for the different DA methods]. The available options for a specific DA method can also be displayed by running the assimilation program for the selected DA method setting `subtype = -1`. (In the tutorial and template codes one can set `-subtype -1` on the command line). == Calling `PDAF_init_forecast` == In `init_pdaf` also the routine `PDAF_init_forecast` is called. This initializes the information, whether further model integrations have to be performed and how many time steps have to be computed in the next forecast aphse. In addition, the ensemble state vectors are written into the model fields so that they can be propagates. In the tutorial codes and the template, the call to `PDAF_init_forecast` is fully implemented. If one bases on the template code, the required files for the call-back routines in the call to `PDAF_init_forecast` are present, but without functionality. This allows to compile and test the code for the initialization step. We explain the arguments of `PDAF_init_forecast` on the following page on [wiki:ModifyModelforEnsembleIntegration Modifying the the model code for the ensemble integration]. == Testing the PDAF initialization == The PDAF initialization can be tested by compiling the assimilation program (one can out-comment the call to `PDAF3_assim_offline` if one likes to focus on the initialization) and executing it. Standard output from PDAF_init looks like the following: {{{ PDAF ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ PDAF +++ PDAF +++ PDAF +++ Parallel Data Assimilation Framework +++ PDAF +++ +++ PDAF +++ Version 3.0beta +++ PDAF +++ +++ PDAF +++ Please cite +++ PDAF +++ L. Nerger and W. Hiller, Computers and Geosciences +++ PDAF +++ 2013, 55, 110-118, doi:10.1016/j.cageo.2012.03.026 +++ PDAF +++ when publishing work resulting from using PDAF +++ PDAF +++ +++ PDAF +++ PDAF itself can also be cited as +++ PDAF +++ L. Nerger. Parallel Data Assimilation Framework +++ PDAF +++ (PDAF). Zenodo. 2024. doi:10.5281/zenodo.7861812 +++ PDAF ++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ PDAF: Initialize filter PDAF ++++++++++++++++++++++++++++++++++++++++++++++++++++++ PDAF +++ Error Subspace Transform Kalman Filter (ESTKF) +++ PDAF +++ +++ PDAF +++ Nerger et al., Mon. Wea. Rev. 140 (2012) 2335 +++ PDAF +++ doi:10.1175/MWR-D-11-00102.1 +++ PDAF ++++++++++++++++++++++++++++++++++++++++++++++++++++++ PDAF: Initialize Parallelization PDAF Parallelization - Filter on model PEs: PDAF Total number of PEs: 6 PDAF Number of parallel model tasks: 6 PDAF PEs for Filter: 1 PDAF # PEs per ensemble task and local ensemble sizes: PDAF Task 1 2 3 4 5 6 PDAF #PEs 1 1 1 1 1 1 PDAF N 1 1 1 1 1 1 PDAF: Call ensemble initialization Initialize state ensemble --- read ensemble from files --- Ensemble size: 6 PDAF: Initialization completed }}} The correctness of the ensemble initialization in `init_ens_pdaf` should be checked by the user.