Changes between Version 32 and Version 33 of ImplementAnalysislseik


Ignore:
Timestamp:
Sep 6, 2010, 8:52:59 AM (14 years ago)
Author:
lnerger
Comment:

--

Legend:

Unmodified
Added
Removed
Modified

  • ImplementAnalysislseik

    v32 v33  
    6464== User-supplied routines ==
    6565
    66 Here all user-supplied routines are described that are required in the call to `PDAF_put_state_lseik`. For some of the generic routines, we link to the page on [ModifyModelforEnsembleIntegration modifying the model code for the ensemble integration].
     66Here, all user-supplied routines are described that are required in the call to `PDAF_put_state_lseik`. For some of the generic routines, we link to the page on [ModifyModelforEnsembleIntegration modifying the model code for the ensemble integration].
    6767
    6868To indicate user-supplied routines we use the prefix `U_`. In the template directory `templates/` these routines are provided in files with the routines name without this prefix. In the example implementation in `testsuite/src/dummymodel_1D` the routines exist without the prefix, but with the extension `_dummy_D.F90`. In the section titles below we provide the name of the template file in parentheses.
    … …  
    8989
    9090Some hints:
    91  * It can be useful to not only determine the size of the observation vector at this point. One can also already gather information about the location of the observations, which will be used later, e.g. to implement the observation operator. In addition, one can already prepare an array that holds the full observation vector. This can be used later by `U_init_obs_l` to initialize a local vector of observations by selecting the relevant parts of the full observation vector. Arrays for the locations can be defined in a module like `mod_assimilation`.
    92  * The routine is similar to `init_dim_obs` used in the global filters. However, if the global filter is used with a domain-decomposed model, it only initializes the size of the observation vector for the local model sub-domain. This is different for the local filters, as the local analysis also requires observational data from neighboring model sub-domains. Anyway, one can base on an implemented routine `init_dim_obs` to implement `init_dim_obs_full`.
     91 * It can be useful to not only determine the size of the observation vector at this point. One can also already gather information about the location of the observations, which can be used later, e.g. to implement the observation operator. In addition, one can already prepare an array that holds the full observation vector. This can be used later by `U_init_obs_l` to initialize a local vector of observations by selecting the relevant parts of the full observation vector. The required arrays can be defined in a module like `mod_assimilation`.
     92 * The routine is similar to `init_dim_obs` used in the global filters. However, if the global filter is used with a domain-decomposed model, it only initializes the size of the observation vector for the local model sub-domain. This is different for the local filters, as the local analysis also requires observational data from neighboring model sub-domains. Nonetheless, one can base on an implemented routine `init_dim_obs` to implement `init_dim_obs_full`.
    9393
    9494=== `U_obs_op_full` (obs_op_full.F90) ===
    … …  
    9999{{{
    100100SUBROUTINE obs_op_full(step, dim_p, dim_obs_f, state_p, m_state_f)
     101
    101102  INTEGER, INTENT(in) :: step               ! Currrent time step
    102103  INTEGER, INTENT(in) :: dim_p              ! PE-local dimension of state
    … …  
    108109The routine is called during the analysis step, before the loop over the local analysis domain is entered. It has to perform the operation of the observation operator acting on a state vector that is provided as `state_p`. The observed state has to be returned in `m_state_f`. It is the observed state corresponding to the 'full' observation vector.
    109110
    110 Hint:
    111  * If the observation operator involves a global operation, e.g. some global integration, while using domain-decompostion one has to gather the information from the other model domains using MPI communication.
    112  * Analogously to the situation with `init_dim_obs_full`, the routine is similar to `init_dim_obs` used for the global filters. However, with a domain-decompoared model also here `m_state_f` will contain parts of the state vector from neighboring model sub-domains. To make these parts accessible, some parallel communication will be necessary (The state information for a neighboring model sub-domain, will be in the memory of the process that handles that sub-domain). The example implementation in `testsuite/dummymodel_1d` uses the function `MPI_AllGatherV` for this communication.
     111Hint:
     112 * Analogously to the situation with `init_dim_obs_full`, the routine is similar to `init_dim_obs` used for the global filters. However, with a domain-decomposed model `m_state_f` will contain parts of the state vector from neighboring model sub-domains. To make these parts accessible, some parallel communication will be necessary (The state information for a neighboring model sub-domain, will be in the memory of the process that handles that sub-domain). The example implementation in `testsuite/dummymodel_1d` uses the function `MPI_AllGatherV` for this communication.
    113113
    114114=== `U_init_obs_full` (init_obs_full.F90) ===
    … …  
    121121SUBROUTINE init_obs_full(step, dim_obs_f, observation_f)
    122122
    123   INTEGER, INTENT(in) :: step             ! Current time step
    124   INTEGER, INTENT(in) :: dim_obs_f        ! Dimension of full observation vector
     123  INTEGER, INTENT(in) :: step                     ! Current time step
     124  INTEGER, INTENT(in) :: dim_obs_f                ! Dimension of full observation vector
    125125  REAL, INTENT(out)   :: observation_f(dim_obs_f) ! Full observation vector
    126126}}}
    … …  
    129129
    130130Hints:
    131  * As for the other 'full' routines: While the global counterpart of this routine (`init_obs`) has to initialize the observation vector for the local model sub-domain, the 'full' routine has to include observations that spatially belong to neighboring model sub-domains. As an easy choice one can simply initialize a vector of all globally available observations.
     131 * As for the other 'full' routines: While the global counterpart of this routine (`init_obs`) has to initialize the observation vector only for the local model sub-domain, the 'full' routine has to include observations that spatially belong to neighboring model sub-domains. As an easy choice one can simply initialize a vector of all globally available observations.
    132132 * If the adaptive forgetting factor is not used, this routine only has to exist. However, no functionality is required.
    133133
    … …  
    148148
    149149The routine is called during the analysis step during the loop over the local analysis domain.
    150 It has to provide the vector of observations for analysis of the local analysis domain of index `domain_p` in `observation_l` for the current time step.
    151 
    152 Hints:
    153  * For parallel efficiency the LSEIK algorithm is implemented in a way that first the full vectors are initialized. Thus, if `observation_f` has been initialized before `U_init_obs_local` is executed (e.g. by `U_init_dim_obs_full`), the operations performed in this routine will be to select the part of the full observation vector that is relevant for the current local analysis domain.
     150It has to provide the vector of observations for the analysis in the local analysis domain owith index `domain_p` in `observation_l` for the current time step.
     151
     152Hints:
     153 * For parallel efficiency, the LSEIK algorithm is implemented in a way that first the full vectors are initialized. These are then restricted to the local analysis domain when the loop over all local analysis domains is executed. Thus, if `observation_f` has been initialized before `U_init_obs_local` is executed (e.g. by `U_init_dim_obs_full`), the operations performed in this routine will be to select the part of the full observation vector that is relevant for the current local analysis domain.
    154154 * The routine `U_init_dim_obs_local` is executed before this routine. Thus, if that routine already prepares the information which elements of `observation_f` are need for `observation_l` this information can be used efficiently here.
    155155
    … …  
    178178}}}
    179179
    180 The routine is called during the loop over the local analysis domains in the analysis step. In the algorithm the product of the inverse of the observation error covariance matrix with some matrix has to be computed. For the SEIK filter this matrix holds the observed part of the ensemble perturbations for the local analysis domain of index `domain_p`. The matrix is provided as `A_l`. The product has to be given as `C_l`.
     180The routine is called during the loop over the local analysis domains. In the algorithm the product of the inverse of the observation error covariance matrix with some matrix has to be computed. For the SEIK filter this matrix holds the observed part of the ensemble perturbations for the local analysis domain of index `domain_p`. The matrix is provided as `A_l`. The product has to be given as `C_l`.
    181181
    182182This routine is also the place to perform observation localization. To initialize a vector of weights, the routine `PDAF_local_weights` can be called. The procedure is used in the example implementation and also demonstrated in the template routine.
    … …  
    201201
    202202The routine is called during the analysis step before the loop over the local analysis domains is entered.
    203 It has to provide the number of local analysis domains. In case of a domain-decomposed model the number of local analysis domain for the model sub-dmain of the calling process has to be initialized.
     203It has to provide the number of local analysis domains. In case of a domain-decomposed model the number of local analysis domain for the model sub-domain of the calling process has to be initialized.
    204204
    205205Hints:
    … …  
    288288}}}
    289289
    290 The routine is called during the loop over the local analysis domains in the analysis step. It has to provide the local state vector `state_l` that corresponds to the local analysis domain with index `domain_p`. Provided to the routine is the state vector `state_p`. With a domain decomposed model, this is the state for the local model sub-domain.
    291 
    292290The routine is called during the loop over the local analysis domains in the analysis step. It has to initialize the part of the global state vector `state_p` that corresponds to the local analysis domain with index `domain_p`. Provided to the routine is the state vector `state_l` for the local analysis domain.
    293291
    … …  
    312310}}}
    313311
    314 The routine is called during the loop over the local analysis domains in the analysis step. It has to provide a local observation vector `mstate_l` for the observation domain that corresponds to the local analysis domain with index `domain_p`. Provided to the routine is the full observation vector `mstate_f` of which the local part has to be extracted.
     312The routine is called during the loop over the local analysis domains in the analysis step. It has to provide a local observation vector `mstate_l` for the observation domain that corresponds to the local analysis domain with index `domain_p`. Provided to the routine is the full observation vector `mstate_f` from which the local part has to be extracted.
    315313
    316314Hints:
    … …  
    321319=== `U_init_obsvar` (init_obsvar.F90) ===
    322320
    323 This routine is used by the global filter algorithms SEIK and  ETKF as well as the local filters LSEIK and LETKF. The routine is only called if the adaptive forgetting factor is used (`type_forget=1` in the example impementation).
    324 
    325 The interface for this routine is:
    326 {{{
    327 SUBROUTINE init_obsvar(step, dim_obs_p, obs_p, meanvar)
    328 
    329   INTEGER, INTENT(in) :: step          ! Current time step
    330   INTEGER, INTENT(in) :: dim_obs_p     ! PE-local dimension of observation vector
    331   REAL, INTENT(in) :: obs_p(dim_obs_p) ! PE-local observation vector
    332   REAL, INTENT(out)   :: meanvar       ! Mean observation error variance
    333 }}}
    334 
    335 The routine is called in the global filters during the analysis or
    336 by the routine that computes an adaptive forgetting factor (PDAF_set_forget).
    337 The routine has to initialize the mean observation error variance. 
    338 For the global filters this should be the global mean.
     321This routine is used by the global filter algorithms SEIK and  ETKF as well as the local filters LSEIK and LETKF. The routine is only called if the adaptive forgetting factor is used (`type_forget=1` in the example implementation). The difference in this routine between global and local filters is that the global filters use 'global' while the local filters use 'full' quantities.
     322
     323The interface for this routine is:
     324{{{
     325SUBROUTINE init_obsvar(step, dim_obs_f, obs_f, meanvar_f)
     326
     327  INTEGER, INTENT(in) :: step             ! Current time step
     328  INTEGER, INTENT(in) :: dim_obs_f        ! Full dimension of observation vector
     329  REAL, INTENT(in)    :: obs_f(dim_obs_f) ! Full observation vector
     330  REAL, INTENT(out)   :: meanvar_f        ! Mean observation error variance
     331}}}
     332
     333The routine is called in the local filters before the loop over all local analysis domains is entered. The call is by the routine that computes an adaptive forgetting factor (PDAF_set_forget).
     334The routine has to initialize a mean full observation error variance, which should be consistent with the observation vector initialized in `U_init_ob_full`.
     335
    339336
    340337Hints:
    341338 * For a model with domain-decomposition one might use the mean variance for the model sub-domain of the calling process. Alternatively one can compute a mean variance for the full model domain using MPI communication (e.g. the function `MPI_allreduce`).
    342  * The observation vector `obs_p` is provided to the rotine for the case that the observation error variance is relative to the value of the observations.
     339 * The observation vector `obs_p` is provided to the routine for the case that the observation error variance is relative to the value of the observations.
    343340 * If the adaptive forgetting factor is not used, this routine has only to exist for the compilation, but it does not need functionality.
    344341
    … …  
    346343=== `U_init_obsvar_local` (init_obsvar_local.F90) ===
    347344
    348 This routine is used by the local filters LSEIK and LETKF. The routine is only called if the local adaptive forgetting factor is used (`type_forget=2` in the example impementation).
     345This routine is used by the local filters LSEIK and LETKF. The routine is only called if the local adaptive forgetting factor is used (`type_forget=2` in the example implementation).
    349346
    350347The interface for this routine is:
    … …  
    359356}}}
    360357
    361 The routine is called in the local filters during the loop over all local analysis domains by the by the routine that computes a local adaptive forgetting factor (PDAF_set_forget_local). The routine has to initialize a local mean observation error variance for all observations used for the analysis in the current local analysis domain.
     358The routine is called in the local filters during the loop over all local analysis domains by the routine that computes a local adaptive forgetting factor (PDAF_set_forget_local). The routine has to initialize a local mean observation error variance for all observations used for the analysis in the current local analysis domain.
    362359
    363360Hints:
    … …  
    366363== Execution order of user-supplied routines ==
    367364
    368 The user-supplied routines are executed in the order listed below. The order can be important as some routines can perform preparatory work for later routines. For example, `U_init_dim_obs_local` can prepare an index array that provides the information how to localize a 'full' vector of observations. Some hints one this are given with the descriptions of the routine interfaces above.
     365The user-supplied routines are executed in the order listed below. The order can be important as some routines can perform preparatory work for routines executed later on during the analysis. For example, `U_init_dim_obs_local` can prepare an index array that provides the information how to localize a 'full' vector of observations. Some hints one this are given with the descriptions of the routine interfaces above.
    369366
    370367Before the analysis step is called the following is executed:
    371368 1. [#U_collect_statecollect_state.F90 U_collect_state] (called once for each ensemble member)
    372369
    373 When the ensemble integration of the forecast is completed the analysis step is executed. Before the loop over all local analysis domains, the following routines are executed:
     370When the ensemble integration of the forecast is completed, the analysis step is executed. Before the loop over all local analysis domains, the following routines are executed:
    374371 1. [#U_prepoststepprepoststep_seik.F90 U_prepoststep] (call to handle the forecast, called with negative value of the time step)
    375372 1. [#U_init_n_domainsinit_n_domains.F90 U_init_n_domains]
    … …  
    377374 1. [#U_obs_op_fullobs_op_full.F90 U_obs_op_full] (Called `dim_ens` times; once for each ensemble member)
    378375 1. [#U_init_obs_fullinit_obs_full.F90 U_init_obs_full] (Only executed, if the global adaptive forgetting factor is used (`type_forget=1` in the example implemention))
    379  1. [#U_init_obsvarinit_obsvar.F90 U_init_obsvar] (Only execute, if the global adaptive forgetting factor is used (`type_forget=1` in the example implemention))
     376 1. [#U_init_obsvarinit_obsvar.F90 U_init_obsvar] (Only executed, if the global adaptive forgetting factor is used (`type_forget=1` in the example implemention))
    380377
    381378The the loop over all local analysis domains, it is executed:
    … …  
    386383 1. [#U_init_obs_localinit_obs_local.F90 U_init_obs_local]
    387384 1. [#U_global2local_obsglobal2local_obs.F90 U_global2local_obs] (`dim_ens` calls; one call to localize the observed part of each ensemble member)
    388  1. [#U_init_obsvar_localinit_obsvar_local.F90 U_init_obsvar_local] (Only, if the local adaptive forgetting factor is used (`type_forget=2` in the example implementation))
     385 1. [#U_init_obsvar_localinit_obsvar_local.F90 U_init_obsvar_local] (Only called, if the local adaptive forgetting factor is used (`type_forget=2` in the example implementation))
    389386 1. [#U_prodRinvA_localprodrinva_local.F90 U_prodRinvA_local]
    390387 1. [#U_local2global_statelocal2global_state.F90 U_local2global_state] (Called `dim_ens+1` times: Once for each ensemble member and once for the mean state estimate)