Changes between Version 22 and Version 23 of ImplementAnalysisseik


Ignore:
Timestamp:
May 17, 2011, 9:55:08 AM (14 years ago)
Author:
lnerger
Comment:

Some corrections.

Legend:

Unmodified
Added
Removed
Modified

  • ImplementAnalysisseik

    v22 v23  
    2626For completeness we discuss here all user-supplied routines that are specified in the interface to PDAF_put_state_seik. Thus, some of the user-supplied that are explained on the page explaining the modification of the model code for the ensemble integration are repeated here.
    2727
     28The SEIK filter and the ETKF (Ensemble Transform Kalman Filter) are very similar. For this reason, the interface to the user-supplied routines is almost identical. Depending on the implementation it can be possible to use identical routines for the SEIK filter and the ETKF. Differences are marked in the text below.
     29
    2830== `PDAF_put_state_seik` ==
    2931
    30 The general espects of the filter specific routines `PDAF_put_state_*` have been described on the page [ModifyModelforEnsembleIntegration Modification of the model core for the ensemble integration]. Here, we list once more the full interface. Subsequently, the full set of user-supplied routines specified in the call to `PDAF_put_state_seik` is explained.
     32The general espects of the filter specific routines `PDAF_put_state_*` have been described on the page [ModifyModelforEnsembleIntegration Modification of the model code for the ensemble integration]. Here, we list once more the full interface. Subsequently, the full set of user-supplied routines specified in the call to `PDAF_put_state_seik` is explained.
    3133
    3234The interface when using the SEIK filter is the following:
    … …  
    5052Here all user-supplied routines are described that are required in the call to `PDAF_put_state_seik`. For some of the generic routines, we link to the page on [ModifyModelforEnsembleIntegration modifying the model code for the ensemble integration].
    5153
    52 To 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. 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.
     54To indicate user-supplied routines we use the prefix `U_`. In the template directory `templates/` these routines are provided in files with the routine's 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.
     55
     56In 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.
    5357
    5458
    5559=== `U_collect_state` (collect_state.F90) ===
    5660
    57 This routine is independent from the filter algorithm used.
    58 See [ModifyModelforEnsembleIntegration#U_collect_statecollect_state.F90 here] for the description of this routine.
     61This routine is independent of the filter algorithm used.
     62See the page [ModifyModelforEnsembleIntegration#U_collect_statecollect_state.F90 modifying the model code for the ensemble integration] for the description of this routine.
    5963
    6064
    … …  
    7478
    7579Some hints:
    76  * 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 locations of the observations, which will be used later, e.g. to implement the observation operator. An array for the locations can be defined in a module like `mod_assimilation`.
     80 * 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 locations of the observations, which will be used later, e.g. to implement the observation operator. An array for the locations can be defined in a module like `mod_assimilation` of the example implementation.
    7781
    7882
    … …  
    120124
    121125=== `U_prepoststep` (prepoststep_seik.F90) ===
    122  
    123 See [ModifyModelforEnsembleIntegration#U_prepoststepprepoststep_seik.F90 here] for the description of this routine.
     126
     127The routine has already been described on the [ModifyModelforEnsembleIntegration#U_prepoststepprepoststep_seik.F90 page on modifying the model code for the ensemble integration]. For completeness, the description is repeated:
     128
     129The interface of the routine is identical for all filters. However, the particular operations that are performed in the routine can be specific for each filter algorithm. Here, we exemplify the interface on the example of the SEIK filter.
     130
     131The interface for this routine is
     132{{{
     133SUBROUTINE prepoststep(step, dim_p, dim_ens, dim_ens_p, dim_obs_p, &
     134                       state_p, Uinv, ens_p, flag)
     135
     136  INTEGER, INTENT(in) :: step        ! Current time step
     137                         ! (When the routine is called before the analysis -step is provided.)
     138  INTEGER, INTENT(in) :: dim_p       ! PE-local state dimension
     139  INTEGER, INTENT(in) :: dim_ens     ! Size of state ensemble
     140  INTEGER, INTENT(in) :: dim_ens_p   ! PE-local size of ensemble
     141  INTEGER, INTENT(in) :: dim_obs_p   ! PE-local dimension of observation vector
     142  REAL, INTENT(inout) :: state_p(dim_p) ! PE-local forecast/analysis state
     143                                     ! The array 'state_p' is not generally not initialized in the case of SEIK/EnKF/ETKF.
     144                                     ! It can be used freely in this routine.
     145  REAL, INTENT(inout) :: Uinv(dim_ens-1, dim_ens-1) ! Inverse of matrix U
     146  REAL, INTENT(inout) :: ens_p(dim_p, dim_ens)      ! PE-local state ensemble
     147  INTEGER, INTENT(in) :: flag        ! PDAF status flag
     148}}}
     149
     150The 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`).
     151
     152The 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.
     153
     154Hint:
     155 * If a user considers to perform adjustments to the estimates (e.g. for balances), this routine is the right place for it.
     156 * 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`.
     157 * The interface has a difference for ETKF and SEIK: For the ETKF, the array `Uinv` has size `dim_ens` x `dim_ens`. In contrast it has size `dim_ens-1` x `dim_ens-1` for the SEIK filter.
    124158
    125159
    126160=== `U_prodRinvA` (prodrinva.F90) ===
    127161
    128 This routine is used by all filters whose algorithm uses the inverse of the observation error covariance matrix (SEEK, SEIK, and ETKF).
     162This routine is used by all filter algorithms that use the inverse of the observation error covariance matrix (SEEK, SEIK, and ETKF).
    129163
    130164The interface for this routine is:
    … …  
    145179
    146180Hints:
    147  * the routine does not require that the product is implemented as a real matrix-matrix product. Rather, the product can be implemented in its most efficient form. For example, if the observation error covariance matrix is diagonal, only the multiplication of the diagonal with matrix `A_p` has to be implemented.
     181 * The routine does not require that the product is implemented as a real matrix-matrix product. Rather, the product can be implemented in its most efficient form. For example, if the observation error covariance matrix is diagonal, only the multiplication of the diagonal with matrix `A_p` has to be implemented.
    148182 * The observation vector `obs_p` is provided through the interface for cases where the observation error variance is relative to the actual value of the observations.
    149 
     183 * The interface has a difference for SEIK and ETKF: For ETKF the third argument is the ensemble size (`dim_ens`), while for SEIK it is the rank (`rank`) of the covariance matrix (usually ensemble size minus one). In addition, the second dimension of `A_p` and `C_p` has size `dim_ens` for ETKF, while it is `rank` for the SEIK filter. (For most cases, this will be irrelevant, because most usually the ensemble array `ens_p` is used for computations, rather than `Uinv`. However, for the SEIK filter with fixed covariance matrix, `Uinv` is required to compute the estimate analysis error. The fixed covariance matrix mode is not available for the ETKF)
    150184
    151185=== `U_init_obsvar` (init_obsvar.F90) ===
    … …  
    177211For the SEIK filter, the user-supplied routines are essentially executed in the order they are listed in the interface to `PDAF_put_state_seik`. The order can be important as some routines can perform preparatory work for later routines. For example, `U_init_dim_obs` can prepare an index array that provides the information for executing the observation operator in `PDAF_obs_op`.
    178212
    179 Before the analysis step is called the following is executed:
     213Before the analysis step is called the following routine is executed:
    180214 1. [#U_collect_statecollect_state.F90 U_collect_state]
    181215
    182 When the ensemble integration of the forecast is completed the analysis step is executed. During the analysis step the following routines are executed in the given order:
     216The analysis step is executed when the ensemble integration of the forecast is completed. During the analysis step the following routines are executed in the given order:
    183217 1. [#U_prepoststepprepoststep_seik.F90 U_prepoststep] (call to handle the forecast, called with negative value of the time step)
    184218 1. [#U_init_dim_obsinit_dim_obs.F90 U_init_dim_obs]
    185  1. [#U_obs_opobs_op.F90 U_obs_op] (One call to operate on the ensemble mean state)
     219 1. [#U_obs_opobs_op.F90 U_obs_op] (A single call to operate on the ensemble mean state)
    186220 1. [#U_init_obsinit_obs.F90 U_init_obs]
    187221 1. [#U_obs_opobs_op.F90 U_obs_op] (`dim_ens` calls; one call for each ensemble member)