test/test6903/Verification_des_mono_et_multi_fonctions_E.py

   1 # -*- coding: utf-8 -*-
   2 #
   3 # Copyright (C) 2008-2022 EDF R&D
   4 #
   5 # This library is free software; you can redistribute it and/or
   6 # modify it under the terms of the GNU Lesser General Public
   7 # License as published by the Free Software Foundation; either
   8 # version 2.1 of the License.
   9 #
  10 # This library is distributed in the hope that it will be useful,
  11 # but WITHOUT ANY WARRANTY; without even the implied warranty of
  12 # MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
  13 # Lesser General Public License for more details.
  14 #
  15 # You should have received a copy of the GNU Lesser General Public
  16 # License along with this library; if not, write to the Free Software
  17 # Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
  18 #
  19 # See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
  20 #
  21 # Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
  22 "Verification du fonctionnement correct d'entrees en mono ou multi-fonctions"
  23
  24 import sys
  25 import unittest
  26 import numpy
  27 from adao import adaoBuilder
  28
  29 # ==============================================================================
  30
  31 def ElementaryFunction01( InputArgument ):
  32     """
  33     Exemple de fonction non-lineaire et non-carree
  34
  35     L'argument en entree est un vecteur au sens mathematique, c'est-a-dire
  36     une suite ordonnee de valeurs reelles. Au sens informatique, c'est tout
  37     objet qui peut se transformer en une serie continue unidimensionnelle de
  38     valeurs reelles. L'argument en sortie est un vecteur Numpy 1D.
  39     """
  40     #
  41     if isinstance( InputArgument, (numpy.ndarray, numpy.matrix, list, tuple) ) or type(InputArgument).__name__ in ('generator','range'):
  42         _subX = numpy.ravel(InputArgument)
  43     else:
  44         raise ValueError("ElementaryFunction01 unkown input type: %s"%(type(InputArgument).__name__,))
  45     #
  46     _OutputArgument = []
  47     _OutputArgument.extend( (-1 + _subX).tolist() )
  48     _OutputArgument.extend( numpy.cos(_subX/2).tolist() )
  49     _OutputArgument.extend( numpy.exp((3.14 * _subX)).tolist() )
  50     #
  51     return numpy.ravel( _OutputArgument )
  52
  53 def MultiFonction01( xSerie ):
  54     """
  55     Exemple de multi-fonction
  56
  57     Pour une liste ordonnee de vecteurs en entree, renvoie en sortie la liste
  58     correspondante de valeurs de la fonction en argument
  59     """
  60     if not (isinstance( xSerie, (list, tuple) ) or type(xSerie).__name__ in ('generator','range')):
  61         raise ValueError("MultiFonction01 unkown input type: %s"%(type(xSerie),))
  62     #
  63     _ySerie = []
  64     for _subX in xSerie:
  65         _ySerie.append( ElementaryFunction01( _subX ) )
  66     #
  67     return _ySerie
  68
  69 # ==============================================================================
  70 class Test_Adao(unittest.TestCase):
  71     def test1(self):
  72         """
  73         Verification du fonctionnement identique pour les algorithmes temporels
  74         en utilisant une fonction non-lineaire et non-carree
  75         """
  76         print("\n        "+self.test1.__doc__.strip()+"\n")
  77         Xa = {}
  78         #
  79         for algo in ("ExtendedKalmanFilter", "KalmanFilter", "EnsembleKalmanFilter", "UnscentedKalmanFilter", "4DVAR"):
  80             print("")
  81             msg = "Algorithme en test en MonoFonction : %s"%algo
  82             print(msg+"\n"+"-"*len(msg))
  83             #
  84             adaopy = adaoBuilder.New()
  85             adaopy.setAlgorithmParameters(Algorithm=algo, Parameters={"EpsilonMinimumExponent":-10, "SetSeed":1000})
  86             adaopy.setBackground         (Vector = [0,1,2])
  87             adaopy.setBackgroundError    (ScalarSparseMatrix = 1.)
  88             adaopy.setObservation        (Vector = [0.5,1.5,2.5,0.5,1.5,2.5,0.5,1.5,2.5])
  89             adaopy.setObservationError   (DiagonalSparseMatrix = "1 1 1 1 1 1 1 1 1")
  90             adaopy.setObservationOperator(OneFunction = ElementaryFunction01)
  91             adaopy.setEvolutionError     (ScalarSparseMatrix = 1.)
  92             adaopy.setEvolutionModel     (Matrix = "1 0 0;0 1 0;0 0 1")
  93             adaopy.setObserver("Analysis",Template="ValuePrinter")
  94             adaopy.execute()
  95             Xa["Mono/"+algo] = adaopy.get("Analysis")[-1]
  96             del adaopy
  97         #
  98         for algo in ("ExtendedKalmanFilter", "KalmanFilter", "EnsembleKalmanFilter", "UnscentedKalmanFilter", "4DVAR"):
  99             print("")
 100             msg = "Algorithme en test en MultiFonction : %s"%algo
 101             print(msg+"\n"+"-"*len(msg))
 102             #
 103             adaopy = adaoBuilder.New()
 104             adaopy.setAlgorithmParameters(Algorithm=algo, Parameters={"EpsilonMinimumExponent":-10, "SetSeed":1000})
 105             adaopy.setBackground         (Vector = [0,1,2])
 106             adaopy.setBackgroundError    (ScalarSparseMatrix = 1.)
 107             adaopy.setObservation        (Vector = [0.5,1.5,2.5,0.5,1.5,2.5,0.5,1.5,2.5])
 108             adaopy.setObservationError   (DiagonalSparseMatrix = "1 1 1 1 1 1 1 1 1")
 109             adaopy.setObservationOperator(OneFunction = MultiFonction01, InputFunctionAsMulti = True)
 110             adaopy.setEvolutionError     (ScalarSparseMatrix = 1.)
 111             adaopy.setEvolutionModel     (Matrix = "1 0 0;0 1 0;0 0 1")
 112             adaopy.setObserver("Analysis",Template="ValuePrinter")
 113             adaopy.execute()
 114             Xa["Multi/"+algo] = adaopy.get("Analysis")[-1]
 115             del adaopy
 116         #
 117         print("")
 118         msg = "Tests des ecarts attendus :"
 119         print(msg+"\n"+"="*len(msg))
 120         for algo in ("ExtendedKalmanFilter", "KalmanFilter", "EnsembleKalmanFilter", "UnscentedKalmanFilter", "4DVAR"):
 121             verify_similarity_of_algo_results(("Multi/"+algo, "Mono/"+algo), Xa, 1.e-20)
 122         print("  Les resultats obtenus sont corrects.")
 123         print("")
 124         #
 125         return 0
 126
 127 # ==============================================================================
 128 def almost_equal_vectors(v1, v2, precision = 1.e-15, msg = ""):
 129     """Comparaison de deux vecteurs"""
 130     print("    Difference maximale %s: %.2e"%(msg, max(abs(v2 - v1))))
 131     return max(abs(v2 - v1)) < precision
 132
 133 def verify_similarity_of_algo_results(serie = [], Xa = {}, precision = 1.e-15):
 134     print("  Comparaisons :")
 135     for algo1 in serie:
 136         for algo2 in serie:
 137             if algo1 is algo2: break
 138             assert almost_equal_vectors( Xa[algo1], Xa[algo2], precision, "entre %s et %s "%(algo1, algo2) )
 139     print("  Algorithmes dont les resultats sont similaires a %.0e : %s\n"%(precision, serie,))
 140     sys.stdout.flush()
 141
 142 #===============================================================================
 143 if __name__ == "__main__":
 144     print("\nAUTODIAGNOSTIC\n==============")
 145     sys.stderr = sys.stdout
 146     unittest.main(verbosity=2)