.. warning::
- this change in execution mode is extremely powerful and flexible. It is
+ This change in execution mode is extremely powerful and flexible. It is
therefore recommended that the user both use it, and at the same time be
attentive to the interaction of the different choices he makes, to avoid, for
example, an unintended performance deterioration, or computer conflicts that
.. warning::
- in case of doubt, it is recommended NOT TO ACTIVATE this parallelism.
+ In case of doubt, it is recommended NOT TO ACTIVATE this parallelism.
It is also recalled that one have to choose the type "*multi*" for the default
container in order to launch the scheme, to allow a really parallel execution.
execution of this "*observer*" is simply never activated.
.. warning::
+
If not using the default available templates, it is up to the user to make
carefully established function scripts or external programs that do not
crash before being registered as an "*observer*" function. The debugging
.. index:: single: TangentOperator
.. index:: single: AdjointOperator
-**In general, it is recommended to use the first functional form rather than
-the second one. A small performance improvement is not a good reason to use a
-detailed implementation as this second functional form.**
+.. warning::
+
+ In general, it is recommended to use the first functional form rather than
+ the second one. A small performance improvement is not a good reason to use a
+ detailed implementation as this second functional form.
The second one consist in providing directly the three associated operators
:math:`O`, :math:`\mathbf{O}` and :math:`\mathbf{O}^*`. This is done by using
where :math:`\mathbf{u}` is the control over one state increment. In fact, the
direct operator has to be applied to a pair of variables :math:`(X,U)`.
-Schematically, the operator has to be set as::
+Schematically, the operator :math:`O` has to be set up as a function applicable
+on a pair :math:`\mathbf{(X, U)}` as follows::
def DirectOperator( pair = (X, U) ):
""" Direct non-linear simulation operator """
.. warning::
- it is strongly recommended not to use this explicit "multiple" functions
+ It is strongly recommended not to use this explicit "multiple" functions
definition without a very strong computing justification. This treatment is
already done by default in ADAO to increase performances. Only the very
experienced user, seeking to manage particularly difficult cases, can be
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.. index:: single: Template
-.. index:: single: AnalysisPrinter
-.. index:: single: AnalysisSaver
-.. index:: single: AnalysisPrinterAndSaver
+.. index:: single: AnalysisPrinter (UserPostAnalysis)
+.. index:: single: AnalysisSaver (UserPostAnalysis)
+.. index:: single: AnalysisPrinterAndSaver (UserPostAnalysis)
These examples present Python commands or scripts which allow to obtain or to
treat the output of an algorithm run. To help the user, they are directly
be build from these simple examples, the main objective being to help the user
to elaborate the exact procedure he needs in output.
-The first example (named "*AnalysisPrinter*" in the inputs of type
-"*Template*") consists in printing, in the standard log output, the value of the
-analysis or the optimal state, noted as :math:`\mathbf{x}^a` in the section
-:ref:`section_theory`. It is realized by the commands::
+The first example (named "*AnalysisPrinter*" in the inputs of type "*Template*"
+for "*UserPostAnalysis*" section) consists in printing, in the standard log
+output, the value of the analysis or the optimal state, noted as
+:math:`\mathbf{x}^a` in the section :ref:`section_theory`. It is realized by
+the commands::
import numpy
xa=numpy.ravel(ADD.get('Analysis')[-1])
contain a real unidimensional vector, whatever the previous computing choices
are.
-A second example (named "*AnalysisSaver*" in the inputs of type "*Template*")
-consists in saving on file the value of the analysis or the optimal state
-:math:`\mathbf{x}^a`. It is realized by the commands::
+A second example (named "*AnalysisSaver*" in the inputs of type "*Template*"
+for "*UserPostAnalysis*" section) consists in saving on file the value of the
+analysis or the optimal state :math:`\mathbf{x}^a`. It is realized by the
+commands::
import numpy
xa=numpy.ravel(ADD.get('Analysis')[-1])
The chosen recording file is a text one named ``/tmp/analysis.txt``.
It is easy to combine these two examples by building a third one (named
-"*AnalysisPrinterAndSaver*" in the inputs of type "*Template*"). It consists in
-simultaneously printing in the standard log output and in saving on file the
-value of :math:`\mathbf{x}^a`. It is realized by the commands::
+"*AnalysisPrinterAndSaver*" in the inputs of type "*Template*" for
+"*UserPostAnalysis*" section). It consists in simultaneously printing in the
+standard log output and in saving on file the value of :math:`\mathbf{x}^a`. It
+is realized by the commands::
import numpy
xa=numpy.ravel(ADD.get('Analysis')[-1])
ADD.setUserPostAnalysis(Template = "AnalysisPrinter")
.. warning::
+
If not using the default available templates, it is up to the user to make
carefully established function scripts or external programs that do not
crash before being registered as an "*observer*" function. The debugging
.. warning::
- in its present version, this command is experimental, and therefore remains
+ In its present version, this command is experimental, and therefore remains
subject to changes in future versions.
.. warning::
- in this particular version, this algorithm or some of its variants are
+ In this particular version, this algorithm or some of its variants are
experimental, and therefore remain subject to change in future versions.
.. math:: \forall n \in \mathbb{N}, \mathbf{x}_{n+1} = M(\mathbf{x}_{n},\mathbf{u}_{n},t_n\rightarrow t_{n+1})
for an indexing :math:`t_n` of discrete times with :math:`n\in\mathbf{N}`.
-:math:`M` is the discrete evolution operator obtained from :math:`\mathcal{D}`.
-Usually, we omit the time notation in the evolution operator :math:`M`.
-Approximating the :math:`\mathcal{D}` operator by :math:`M` introduces (or
-adds, if it already exists) a :math:`\epsilon` model error.
+:math:`M` is the discrete evolution operator, symbolically obtained from
+:math:`\mathcal{D}` by the discretization scheme. Usually, we omit the time
+notation in the evolution operator :math:`M`. Approximating the
+:math:`\mathcal{D}` operator by :math:`M` introduces (or adds, if it already
+exists) a :math:`\epsilon` model error.
We can then characterize two types of estimates in dynamics, which we describe
hereafter on the discrete time dynamical system: `State estimation in
where :math:`\mathbf{x}` is the system state to be estimated,
:math:`\mathbf{x}_{n}` and :math:`\mathbf{y}_{n}` are respectively the
-(computed) unobserved and (measured) observed state of the system, :math:`M`
+computed (unobserved) and measured (observed) state of the system, :math:`M`
and :math:`H` are the incremental evolution and observation operators,
respectively, :math:`\mathbf{\epsilon}_{n}` and :math:`\mathbf{\nu}_{n}` are
the evolution and observation noise or error, respectively, and
.. math:: \mathbf{y}_{n} = G(\mathbf{x}_{n},\mathbf{a}_{n}) + \mathbf{\nu}_{n}
-where, this time, the choices of the evolution and observation error models
+where, this time, the choice of the evolution and observation error models
:math:`\mathbf{\epsilon}_{n}` and :math:`\mathbf{\nu}_{n}` condition the
-performance of convergence and observation tracking. The estimation of the
-parameters :math:`\mathbf{a}` is done by using pairs
-:math:`(\mathbf{x}_{n},\mathbf{y}_{n})` of corresponding inputs and outputs.
+performance of convergence and observation tracking (while the error
+representations come from the behavior of the physics in the case of state
+estimation). The estimation of the parameters :math:`\mathbf{a}` is done by
+using pairs :math:`(\mathbf{x}_{n},\mathbf{y}_{n})` of corresponding inputs and
+outputs.
In this case of parameter estimation, in order to apply data assimilation
methods, we therefore impose the hypothesis that the evolution operator is the
A special case concerns the joint estimation of state and parameters used in a
dynamic system. One seeks to jointly estimate the state :math:`\mathbf{x}`
-(which depends on time) and the parameters :math:`\mathbf{a}` (which does not
-depend on time). There are several ways to deal with this problem, but the most
-general one is to use a state vector augmented by the parameters, and to extend
-the operators accordingly.
+(which depends on time) and the parameters :math:`\mathbf{a}` (which here does
+not depend on time). There are several ways to deal with this problem, but the
+most general one is to use a state vector augmented by the parameters, and to
+extend the operators accordingly.
To do this, using the notations of the previous two subsections, we define the
auxiliary variable :math:`\mathbf{w}` such that:
.. math:: \mathbf{y}_{n} = \tilde{H}(\mathbf{w}_{n}) + \mathbf{\nu}_{n}
-The incremental evolution and observation operators are therefore respectively
-the augmented operators :math:`\tilde{M}` and :math:`\tilde{H}`, and are
-directly usable in data assimilation with ADAO.
+where :math:`\mathbf{w}_{n}=[\mathbf{x}_n~~\mathbf{a}_n]^T`. The incremental
+evolution and observation operators are therefore respectively the augmented
+operators :math:`\tilde{M}` and :math:`\tilde{H}`, and are directly usable in
+data assimilation with ADAO.
Conceptual scheme for data assimilation in dynamics
+++++++++++++++++++++++++++++++++++++++++++++++++++
simulation is often called "*high fidelity simulation*" (or "*full scale
simulation*").
-In a generic way, **reduction methods thus aim at reducing the computational
-cost of the optimization while controlling as much as possible the numerical
-error implied by this reduction**.
+In a generic way, **different strategies to reduce the cost of the optimization
+calculation exist, and some of them also allow to control the numerical error
+implied by this reduction**.
To establish this, one seeks to reduce at least one of the ingredients that
make up the data assimilation or optimization problem. One can thus classify
time. ADAO integrates some of the most robust methods, but this aspect is
still largely the subject of research and development.
-One can end this quick overview of reduction methods by pointing out that their
+One can end this quick overview of reduction methods highlighting that their
use is ubiquitous in real applications and in numerical tools, and that ADAO
allows to use proven methods without even knowing it.
.. warning::
- ce changement de mode d'exécution est extrêmement puissant et souple. Il est
+ Ce changement de mode d'exécution est extrêmement puissant et souple. Il est
donc recommandé à l'utilisateur à la fois de l'utiliser, et en même temps
d'être attentif à l'interaction des différents choix qu'il effectue, pour
éviter par exemple une dégradation involontaire des performances, ou des
.. warning::
- en cas de doute, il est recommandé de NE PAS ACTIVER ce parallélisme.
+ En cas de doute, il est recommandé de NE PAS ACTIVER ce parallélisme.
On rappelle aussi qu'il faut choisir dans YACS un container par défaut de type
"*multi*" pour le lancement du schéma, pour permettre une exécution
simplement jamais activée.
.. warning::
+
Si les modèles disponibles par défaut ne sont pas utilisés, il revient à
l'utilisateur de faire des scripts de fonctions soigneusement établis ou
des programmes externes qui ne se plantent pas avant d'être enregistrés
.. warning::
- en général, il est recommandé d'utiliser la première forme fonctionnelle
+ En général, il est recommandé d'utiliser la première forme fonctionnelle
plutôt que la seconde. Un petit accroissement de performances n'est pas une
bonne raison pour utiliser l'implémentation détaillée de cette seconde forme
fonctionnelle.
.. warning::
- il est recommandé de ne pas utiliser cette troisième forme fonctionnelle sans
+ Il est recommandé de ne pas utiliser cette troisième forme fonctionnelle sans
une solide raison numérique ou physique. Un accroissement de performances
n'est pas une bonne raison pour utiliser la complexité de cette troisième
forme fonctionnelle. Seule une impossibilité à utiliser les première ou
fichier script ou un code externe nommé ici
"*Physical_simulation_functions.py*", contenant trois fonctions nommées
"*DirectOperator*", "*TangentOperator*" et "*AdjointOperator*" comme
-précédemment. Voici le squelette d'aiguillage::
+précédemment. Voici le squelette d'aiguillage:
+::
import Physical_simulation_functions
import numpy, logging, codecs, pickle
Dans certains cas, l'opérateur d'évolution ou d'observation doit être contrôlé
par un contrôle d'entrée externe, qui est donné *a priori*. Dans ce cas, la
forme générique du modèle incrémental :math:`O` est légèrement modifiée comme
-suit:
+suit :
.. math:: \mathbf{y} = O( \mathbf{x}, \mathbf{u})
où :math:`\mathbf{u}` est le contrôle sur l'incrément d'état. En effet,
l'opérateur direct doit être appliqué à une paire de variables :math:`(X,U)`.
Schématiquement, l'opérateur :math:`O` doit être construit comme une fonction
-applicable sur une paire:math:`\mathbf{(X, U)}` suit::
+applicable sur une paire :math:`\mathbf{(X, U)}` comme suit :
+::
def DirectOperator( paire = (X, U) ):
""" Opérateur direct de simulation non-linéaire """
utilise le plus souvent comme référence les valeurs par défaut
:math:`\mathbf{x}^b` (ébauche, ou valeur nominale). Pourvu que chaque
composante de :math:`\mathbf{x}^b` soit non nulle, on peut ensuite procéder par
-correction multiplicative. Pour cela, on peut par exemple poser:
+correction multiplicative. Pour cela, on peut par exemple poser :
.. math:: \mathbf{x} = \mathbf{\alpha}\mathbf{x}^b
et optimiser ensuite le paramètre multiplicatif :math:`\mathbf{\alpha}`. Ce
paramètre a pour valeur par défaut (ou pour ébauche) un vecteur de 1. De
manière similaire, on peut procéder par correction additive si c'est plus
-judicieux pour la physique sous-jacente. Ainsi, dans ce cas, on peut poser:
+judicieux pour la physique sous-jacente. Ainsi, dans ce cas, on peut poser :
.. math:: \mathbf{x} =\mathbf{x}^b + \mathbf{\alpha}
.. warning::
- il est fortement recommandé de ne pas utiliser cette gestion explicite de
+ Il est fortement recommandé de ne pas utiliser cette gestion explicite de
fonctions "multiples" sans une très solide raison informatique pour le faire.
Cette gestion est déjà effectuée par défaut dans ADAO pour l'amélioration des
performances. Seul l'utilisateur très averti, cherchant à gérer des cas
série d'arguments, pour restituer en sortie la série des valeurs
correspondantes. En pseudo-code, la fonction "multiple", ici nommée
``MultiFunctionO``, représentant l'opérateur classique :math:`O` nommé
-"*DirectOperator*", effectue::
+"*DirectOperator*", effectue :
+::
def MultiFunctionO( Inputs ):
""" Multiple ! """
++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
.. index:: single: Template
-.. index:: single: AnalysisPrinter
-.. index:: single: AnalysisSaver
-.. index:: single: AnalysisPrinterAndSaver
+.. index:: single: AnalysisPrinter (UserPostAnalysis)
+.. index:: single: AnalysisSaver (UserPostAnalysis)
+.. index:: single: AnalysisPrinterAndSaver (UserPostAnalysis)
Ces exemples présentent des commandes ou scripts Python qui permettent d'obtenir
ou de traiter les sorties d'une exécution d'algorithme. Pour aider
dont il a besoin en sortie.
Le premier exemple (appelé "*AnalysisPrinter*" dans les entrées de type
-"*Template*") consiste à afficher, dans la sortie standard d'exécution, la
-valeur de l'analyse ou de l'état optimal, noté :math:`\mathbf{x}^a` dans la
-partie :ref:`section_theory`. Cela se réalise par les commandes::
+"*Template*" pour la section "*UserPostAnalysis*") consiste à afficher, dans la
+sortie standard d'exécution, la valeur de l'analyse ou de l'état optimal, noté
+:math:`\mathbf{x}^a` dans la partie :ref:`section_theory`. Cela se réalise par
+les commandes::
import numpy
xa=numpy.ravel(ADD.get('Analysis')[-1])
précédents.
Un second exemple (appelé "*AnalysisSaver*" dans les entrées de type
-"*Template*") consiste à enregistrer sur fichier la valeur de l'analyse ou de
-l'état optimal :math:`\mathbf{x}^a`. Cela se réalise par les commandes::
+"*Template*" pour la section "*UserPostAnalysis*") consiste à enregistrer sur
+fichier la valeur de l'analyse ou de l'état optimal :math:`\mathbf{x}^a`. Cela
+se réalise par les commandes::
import numpy
xa=numpy.ravel(ADD.get('Analysis')[-1])
Le fichier d'enregistrement choisi est un fichier texte ``/tmp/analysis.txt``.
Il est aisé de combiner ces deux exemples pour en construire un troisième
-(appelé "*AnalysisPrinterAndSaver*" dans les entrées de type "*Template*"). Il
-consiste à simultanément afficher dans la sortie standard d'exécution et à
-enregistrer sur fichier la valeur de :math:`\mathbf{x}^a`. Cela se réalise par
-les commandes::
+(appelé "*AnalysisPrinterAndSaver*" dans les entrées de type "*Template*" pour
+la section "*UserPostAnalysis*"). Il consiste à simultanément afficher dans la
+sortie standard d'exécution et à enregistrer sur fichier la valeur de
+:math:`\mathbf{x}^a`. Cela se réalise par les commandes::
import numpy
xa=numpy.ravel(ADD.get('Analysis')[-1])
ADD.setUserPostAnalysis(Template = "AnalysisPrinter")
.. warning::
+
Si les modèles disponibles par défaut ne sont pas utilisés, il revient à
l'utilisateur de faire des scripts soigneusement établis et vérifiés, ou
des programmes externes qui ne se plantent pas, avant d'être enregistrés
.. warning::
- dans sa présente version, cet commande est expérimentale, et reste donc
+ Dans sa présente version, cet commande est expérimentale, et reste donc
susceptible de changements dans les prochaines versions.
.. warning::
- dans la présente version, cet algorithme ou certaines de ses variantes sont
+ Dans la présente version, cet algorithme ou certaines de ses variantes sont
expérimentaux, et restent donc susceptibles de changements dans les
prochaines versions.
.. math:: \forall n \in \mathbb{N}, \mathbf{x}_{n+1} = M(\mathbf{x}_{n},\mathbf{u}_{n},t_n\rightarrow t_{n+1})
pour une indexation :math:`t_n` des temps discrets avec :math:`n\in\mathbb{N}`.
-:math:`M` est l'opérateur d'évolution discret issu de :math:`\mathcal{D}`.
-Usuellement, on omet la notation du temps dans l'opérateur d'évolution
-:math:`M`. L'approximation de l'opérateur :math:`\mathcal{D}` par :math:`M`
-introduit (ou ajoute, si elle existe déjà) une erreur de modèle
-:math:`\epsilon`.
+:math:`M` est l'opérateur d'évolution discret, issu symboliquement de
+:math:`\mathcal{D}` par le schéma de discrétisation. Usuellement, on omet la
+notation du temps dans l'opérateur d'évolution :math:`M`. L'approximation de
+l'opérateur :math:`\mathcal{D}` par :math:`M` introduit (ou ajoute, si elle
+existe déjà) une erreur de modèle :math:`\epsilon`.
-On peut alors caractériser deux types d'estimations en dynamique, que l'on
+On peut alors caractériser deux types d'estimation en dynamique, que l'on
décrit ci-après sur le système dynamique en temps discret : `Estimation d'état
en dynamique`_ et `Estimation de paramètres en dynamique`_. Combinés, les deux
types peuvent permettre de faire une `Estimation conjointe d'état et de
.. math:: \mathbf{y}_{n} = H(\mathbf{x}_{n}) + \mathbf{\nu}_{n}
où :math:`\mathbf{x}` est l'état à estimer du système, :math:`\mathbf{x}_{n}`
-et :math:`\mathbf{y}_{n}` sont respectivement l'état (calculé) non observé et
-(mesuré) observé du système, :math:`M` et :math:`H` sont respectivement les
+et :math:`\mathbf{y}_{n}` sont respectivement l'état calculé (non observé) et
+mesuré (observé) du système, :math:`M` et :math:`H` sont respectivement les
opérateurs d'évolution incrémentale et d'observation,
:math:`\mathbf{\epsilon}_{n}` et :math:`\mathbf{\nu}_{n}` sont respectivement
les bruits ou erreurs d'évolution et d'observation, et :math:`\mathbf{u}_{n}`
.. math:: \mathbf{y}_{n} = G(\mathbf{x}_{n},\mathbf{a}_{n}) + \mathbf{\nu}_{n}
-où, cette fois, les choix des modèles d'erreurs d'évolution et d'observation
-:math:`\mathbf{\epsilon}_{n}` et :math:`\mathbf{\nu}_{n}` conditionnent la
-performance de la convergence et du suivi des observations. L'estimation des
-paramètres :math:`\mathbf{a}` se fait par utilisation de paires
-:math:`(\mathbf{x}_{n},\mathbf{y}_{n})` d'entrées et de sorties
-correspondantes.
+où, cette fois, le choix des modèles d'erreurs d'évolution et d'observation
+:math:`\mathbf{\epsilon}_{n}` et :math:`\mathbf{\nu}_{n}` conditionne la
+performance de la convergence et du suivi des observations (alors que les
+représentations d'erreurs proviennent du comportement de la physique dans le
+cas de l'estimation d'état). L'estimation des paramètres :math:`\mathbf{a}` se
+fait par utilisation de paires :math:`(\mathbf{x}_{n},\mathbf{y}_{n})`
+d'entrées et de sorties correspondantes.
Dans ce cas de l'estimation de paramètres, pour appliquer les méthodes
d'assimilation de données, on impose donc l'hypothèse que l'opérateur
Un cas spécial concerne l'estimation conjointe d'état et de paramètres utilisés
dans un système dynamique. On cherche à estimer conjointement l'état
:math:`\mathbf{x}` (qui dépend du temps) et les paramètres :math:`\mathbf{a}`
-(qui ne dépendent pas du temps). Il existe plusieurs manières de traiter ce
+(qui ici ne dépendent pas du temps). Il existe plusieurs manières de traiter ce
problème, mais la plus générale consiste à utiliser un vecteur d'état augmenté
par les paramètres, et à étendre les opérateurs en conséquence.
.. math:: \mathbf{y}_{n} = \tilde{H}(\mathbf{w}_{n}) + \mathbf{\nu}_{n}
-Les opérateurs d'évolution incrémentale et d'observation sont donc
-respectivement les opérateurs augmentés :math:`\tilde{M}` et :math:`\tilde{H}`,
-et sont directement utilisables en assimilation de données avec ADAO.
+avec :math:`\mathbf{w}_{n}=[\mathbf{x}_n~~\mathbf{a}_n]^T`. Les opérateurs
+d'évolution incrémentale et d'observation sont donc respectivement les
+opérateurs augmentés :math:`\tilde{M}` et :math:`\tilde{H}`, et sont
+directement utilisables en assimilation de données avec ADAO.
Schéma conceptuel pour l'assimilation de données en dynamique
+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++
simulation physique complète est souvent appelée "*simulation haute fidélité*"
(ou "*full scale simulation*").
-De manière générale, **les méthodes de réduction visent donc à réduire le coût de
-calcul de l'optimisation tout en contrôlant au mieux l'erreur numérique
-impliquée par cette réduction**.
+De manière générale, **différentes stratégies de réduction du coût du calcul
+d'optimisation existent, et certaines permettent également de contrôler au
+mieux l'erreur numérique impliquée par cette réduction**.
Pour établir cela, on cherche à réduire au moins l'un des ingrédients qui
composent le problème d'assimilation de données ou d'optimisation. On peut
certaines méthodes parmi les plus robustes, mais cet aspect fait toujours
largement l'objet de recherches et d'évolutions.
-On peut terminer ce rapide tour d'horizon des méthodes de réduction par le fait
-que leur usage est omni-présent dans les applications réelles et dans les
-outils numériques, et qu'ADAO permet d'utiliser les méthodes éprouvées sans
-même le savoir.
+On peut terminer ce rapide tour d'horizon des méthodes de réduction en
+soulignant que leur usage est omni-présent dans les applications réelles et
+dans les outils numériques, et qu'ADAO permet d'utiliser des méthodes éprouvées
+sans même le savoir.
Xb = Xt + normal(0, 20%*Xt)
-Pour décrire la matrice des covariances d'erreur d'ébauche math:`\mathbf{B}`,
+Pour décrire la matrice des covariances d'erreur d'ébauche :math:`\mathbf{B}`,
on fait comme précédemment l'hypothèse d'erreurs décorrélées (c'est-à-dire, une
matrice diagonale, de taille 3x3 parce-que :math:`\mathbf{x}^b` est de taille
3) et d'avoir la même variance de 0,1 pour toutes les variables. On obtient :
