]> SALOME platform Git repositories - modules/adao.git/commitdiff
Salome HOME
Minor documentation and code review corrections (21)
authorJean-Philippe ARGAUD <jean-philippe.argaud@edf.fr>
Thu, 10 Feb 2022 20:29:07 +0000 (21:29 +0100)
committerJean-Philippe ARGAUD <jean-philippe.argaud@edf.fr>
Thu, 10 Feb 2022 20:29:07 +0000 (21:29 +0100)
doc/en/snippets/CurrentIterationNumber.rst
doc/en/snippets/Observation.rst
doc/fr/snippets/Observation.rst
src/daComposant/daAlgorithms/3DVAR.py
src/daComposant/daCore/BasicObjects.py
src/daComposant/daCore/NumericObjects.py

index d1eea0079e46b55842c9c7d27dcf32d04a330bc8..840942180d7f053910746c701daab2253461dc8d 100644 (file)
@@ -1,8 +1,8 @@
 .. index:: single: CurrentIterationNumber
 
 CurrentIterationNumber
-  *List of integers*. Each element is the iteration index at the current step during the
-  iterative algorithm procedure.
+  *List of integers*. Each element is the iteration index at the current step
+  during the iterative algorithm procedure.
 
   Example:
   ``i = ADD.get("CurrentIterationNumber")[-1]``
index 53c095c9efa01cd62f2bc49528f478e185faf29d..49518fd4ef5fed5e2491dbc2efc85e0fd3d435bd 100644 (file)
@@ -1,8 +1,9 @@
 .. index:: single: Observation
 
 Observation
-  *Vector*. The variable indicates the observation vector used for data
-  assimilation or optimization, usually noted as :math:`\mathbf{y}^o`. Its
-  value is defined as a "*Vector*" or "*VectorSerie*" type object. Its
-  availability in output is conditioned by the boolean "*Stored*" associated
-  with input.
+  *List of vectors*. The variable indicates the observation vector used for
+  data assimilation or optimization, and usually noted :math:`\mathbf{y}^o`.
+  Its value is defined as an object of type "*Vector*" if it is a single
+  observation (temporal or not) or "*VectorSeries*" if it is a succession of
+  observations. Its availability in output is conditioned by the boolean
+  "*Stored*" associated in input.
index 81bcb3b741df53859dd8aa0bd1b9284242d5640a..da59d99e72978d16025ab298b33cb8a0e373bf31 100644 (file)
@@ -1,8 +1,9 @@
 .. index:: single: Observation
 
 Observation
-  *Vecteur*. La variable désigne le vecteur d'observation utilisé en
+  *Liste de vecteurs*. La variable désigne le vecteur d'observation utilisé en
   assimilation de données ou en optimisation, et usuellement noté
   :math:`\mathbf{y}^o`. Sa valeur est définie comme un objet de type "*Vector*"
-  ou "*VectorSerie*". Sa disponibilité en sortie est conditionnée par le
-  booléen "*Stored*" associé en entrée.
+  si c'est une unique observation (temporelle ou pas) ou "*VectorSerie*" si
+  c'est une succession d'observations. Sa disponibilité en sortie est
+  conditionnée par le booléen "*Stored*" associé en entrée.
index 28ae3bc4e4dd97c0c22b4d99818659034d902238..1500028e2577ae4566635bd72526b2f9de06578a 100644 (file)
@@ -43,7 +43,6 @@ class ElementaryAlgorithm(BasicObjects.Algorithm):
             listadv  = [
                 "3DVAR-Std",
                 "Incr3DVAR",
-                "OneCorrection3DVAR-Std",
                 ],
             )
         self.defineRequiredParameter(
index 47b8c371f78d70b8a6284fb577352a0fa8e2cdd0..09f4bdb35f05980eafe0af9c88e1be38a2bcf805 100644 (file)
@@ -625,6 +625,7 @@ class Algorithm(object):
             - CurrentIterationNumber : numéro courant d'itération dans les algorithmes itératifs, à partir de 0
             - CurrentOptimum : état optimal courant lors d'itérations
             - CurrentState : état courant lors d'itérations
+            - CurrentStepNumber : numéro courant de pas de mesure dans les algorithmes temporels
             - GradientOfCostFunctionJ  : gradient de la fonction-coût globale
             - GradientOfCostFunctionJb : gradient de la partie ébauche de la fonction-coût
             - GradientOfCostFunctionJo : gradient de la partie observations de la fonction-coût
@@ -684,6 +685,7 @@ class Algorithm(object):
         self.StoredVariables["CurrentIterationNumber"]               = Persistence.OneIndex(name  = "CurrentIterationNumber")
         self.StoredVariables["CurrentOptimum"]                       = Persistence.OneVector(name = "CurrentOptimum")
         self.StoredVariables["CurrentState"]                         = Persistence.OneVector(name = "CurrentState")
+        self.StoredVariables["CurrentStepNumber"]                    = Persistence.OneIndex(name  = "CurrentStepNumber")
         self.StoredVariables["ForecastCovariance"]                   = Persistence.OneMatrix(name = "ForecastCovariance")
         self.StoredVariables["ForecastState"]                        = Persistence.OneVector(name = "ForecastState")
         self.StoredVariables["GradientOfCostFunctionJ"]              = Persistence.OneVector(name = "GradientOfCostFunctionJ")
@@ -1094,6 +1096,7 @@ class PartialAlgorithm(object):
         self.StoredVariables["CostFunctionJb"]                       = Persistence.OneScalar(name = "CostFunctionJb")
         self.StoredVariables["CostFunctionJo"]                       = Persistence.OneScalar(name = "CostFunctionJo")
         self.StoredVariables["CurrentIterationNumber"]               = Persistence.OneIndex(name  = "CurrentIterationNumber")
+        self.StoredVariables["CurrentStepNumber"]                    = Persistence.OneIndex(name  = "CurrentStepNumber")
         #
         self.__canonical_stored_name = {}
         for k in self.StoredVariables:
@@ -2099,7 +2102,10 @@ class Covariance(object):
             return self.__C + numpy.asmatrix(other)
         elif self.isvector() or self.isscalar():
             _A = numpy.asarray(other)
-            _A.reshape(_A.size)[::_A.shape[1]+1] += self.__C
+            if len(_A.shape) == 1:
+                _A.reshape((-1,1))[::2] += self.__C
+            else:
+                _A.reshape(_A.size)[::_A.shape[1]+1] += self.__C
             return numpy.asmatrix(_A)
 
     def __radd__(self, other):
index c1feb221981cc671808b5293eb96286c7c7e3de8..d1fbf96029312f9ce614ca09ef1693d515cbb9aa 100644 (file)
@@ -477,24 +477,24 @@ class FDApproximation(object):
             else:                return _HaY
 
 # ==============================================================================
-def EnsembleOfCenteredPerturbations( __bgcenter, __bgcovariance, __nbmembers ):
-    "Génération d'un ensemble de taille __nbmembers-1 d'états aléatoires centrés"
+def EnsembleOfCenteredPerturbations( __bgCenter, __bgCovariance, __nbMembers ):
+    "Génération d'un ensemble de taille __nbMembers-1 d'états aléatoires centrés"
     #
-    __bgcenter = numpy.ravel(__bgcenter)[:,None]
-    if __nbmembers < 1:
-        raise ValueError("Number of members has to be strictly more than 1 (given number: %s)."%(str(__nbmembers),))
+    __bgCenter = numpy.ravel(__bgCenter)[:,None]
+    if __nbMembers < 1:
+        raise ValueError("Number of members has to be strictly more than 1 (given number: %s)."%(str(__nbMembers),))
     #
-    if __bgcovariance is None:
-        _Perturbations = numpy.tile( __bgcenter, __nbmembers)
+    if __bgCovariance is None:
+        _Perturbations = numpy.tile( __bgCenter, __nbMembers)
     else:
-        _Z = numpy.random.multivariate_normal(numpy.zeros(__bgcenter.size), __bgcovariance, size=__nbmembers).T
-        _Perturbations = numpy.tile( __bgcenter, __nbmembers) + _Z
+        _Z = numpy.random.multivariate_normal(numpy.zeros(__bgCenter.size), __bgCovariance, size=__nbMembers).T
+        _Perturbations = numpy.tile( __bgCenter, __nbMembers) + _Z
     #
     return _Perturbations
 
 # ==============================================================================
-def EnsembleOfBackgroundPerturbations( __bgcenter, __bgcovariance, __nbmembers, __withSVD = True):
-    "Génération d'un ensemble de taille __nbmembers-1 d'états aléatoires centrés"
+def EnsembleOfBackgroundPerturbations( __bgCenter, __bgCovariance, __nbMembers, __withSVD = True):
+    "Génération d'un ensemble de taille __nbMembers-1 d'états aléatoires centrés"
     def __CenteredRandomAnomalies(Zr, N):
         """
         Génère une matrice de N anomalies aléatoires centrées sur Zr selon les
@@ -508,31 +508,31 @@ def EnsembleOfBackgroundPerturbations( __bgcenter, __bgcovariance, __nbmembers,
         Zr = numpy.dot(Q,Zr)
         return Zr.T
     #
-    __bgcenter = numpy.ravel(__bgcenter).reshape((-1,1))
-    if __nbmembers < 1:
-        raise ValueError("Number of members has to be strictly more than 1 (given number: %s)."%(str(__nbmembers),))
-    if __bgcovariance is None:
-        _Perturbations = numpy.tile( __bgcenter, __nbmembers)
+    __bgCenter = numpy.ravel(__bgCenter).reshape((-1,1))
+    if __nbMembers < 1:
+        raise ValueError("Number of members has to be strictly more than 1 (given number: %s)."%(str(__nbMembers),))
+    if __bgCovariance is None:
+        _Perturbations = numpy.tile( __bgCenter, __nbMembers)
     else:
         if __withSVD:
-            _U, _s, _V = numpy.linalg.svd(__bgcovariance, full_matrices=False)
-            _nbctl = __bgcenter.size
-            if __nbmembers > _nbctl:
+            _U, _s, _V = numpy.linalg.svd(__bgCovariance, full_matrices=False)
+            _nbctl = __bgCenter.size
+            if __nbMembers > _nbctl:
                 _Z = numpy.concatenate((numpy.dot(
                     numpy.diag(numpy.sqrt(_s[:_nbctl])), _V[:_nbctl]),
-                    numpy.random.multivariate_normal(numpy.zeros(_nbctl),__bgcovariance,__nbmembers-1-_nbctl)), axis = 0)
+                    numpy.random.multivariate_normal(numpy.zeros(_nbctl),__bgCovariance,__nbMembers-1-_nbctl)), axis = 0)
             else:
-                _Z = numpy.dot(numpy.diag(numpy.sqrt(_s[:__nbmembers-1])), _V[:__nbmembers-1])
-            _Zca = __CenteredRandomAnomalies(_Z, __nbmembers)
-            _Perturbations = __bgcenter + _Zca
+                _Z = numpy.dot(numpy.diag(numpy.sqrt(_s[:__nbMembers-1])), _V[:__nbMembers-1])
+            _Zca = __CenteredRandomAnomalies(_Z, __nbMembers)
+            _Perturbations = __bgCenter + _Zca
         else:
-            if max(abs(__bgcovariance.flatten())) > 0:
-                _nbctl = __bgcenter.size
-                _Z = numpy.random.multivariate_normal(numpy.zeros(_nbctl),__bgcovariance,__nbmembers-1)
-                _Zca = __CenteredRandomAnomalies(_Z, __nbmembers)
-                _Perturbations = __bgcenter + _Zca
+            if max(abs(__bgCovariance.flatten())) > 0:
+                _nbctl = __bgCenter.size
+                _Z = numpy.random.multivariate_normal(numpy.zeros(_nbctl),__bgCovariance,__nbMembers-1)
+                _Zca = __CenteredRandomAnomalies(_Z, __nbMembers)
+                _Perturbations = __bgCenter + _Zca
             else:
-                _Perturbations = numpy.tile( __bgcenter, __nbmembers)
+                _Perturbations = numpy.tile( __bgCenter, __nbMembers)
     #
     return _Perturbations
 
@@ -552,9 +552,9 @@ def EnsembleOfAnomalies( __Ensemble, __OptMean = None, __Normalisation = 1.):
     return __Normalisation * (numpy.asarray( __Ensemble ) - __Em)
 
 # ==============================================================================
-def EnsembleErrorCovariance( __Ensemble, __quick = False ):
+def EnsembleErrorCovariance( __Ensemble, __Quick = False ):
     "Renvoie l'estimation empirique de la covariance d'ensemble"
-    if __quick:
+    if __Quick:
         # Covariance rapide mais rarement définie positive
         __Covariance = numpy.cov( __Ensemble )
     else:
@@ -638,6 +638,9 @@ def CovarianceInflation(
     else:
         __InflationFactor = float(__InflationFactor)
     #
+    __InputCovOrEns = numpy.asarray(__InputCovOrEns)
+    if __InputCovOrEns.size == 0: return __InputCovOrEns
+    #
     if __InflationType in ["MultiplicativeOnAnalysisCovariance", "MultiplicativeOnBackgroundCovariance"]:
         if __InflationFactor < 1.:
             raise ValueError("Inflation factor for multiplicative inflation has to be greater or equal than 1.")
@@ -659,17 +662,20 @@ def CovarianceInflation(
             raise ValueError("Inflation factor for additive inflation has to be greater or equal than 0.")
         if __InflationFactor < mpr:
             return __InputCovOrEns
-        __n, __m = numpy.asarray(InputCovOrEns).shape
+        __n, __m = __InputCovOrEns.shape
         if __n != __m:
             raise ValueError("Additive inflation can only be applied to squared (covariance) matrix.")
-        __OutputCovOrEns = (1. - __InflationFactor) * __InputCovOrEns + __InflationFactor * numpy.identity(__n)
+        __tr = __InputCovOrEns.trace()/__n
+        if __InflationFactor > __tr:
+            raise ValueError("Inflation factor for additive inflation has to be small over %.0e."%__tr)
+        __OutputCovOrEns = (1. - __InflationFactor)*__InputCovOrEns + __InflationFactor * numpy.identity(__n)
     #
     elif __InflationType == "HybridOnBackgroundCovariance":
         if __InflationFactor < 0.:
             raise ValueError("Inflation factor for hybrid inflation has to be greater or equal than 0.")
         if __InflationFactor < mpr:
             return __InputCovOrEns
-        __n, __m = numpy.asarray(InputCovOrEns).shape
+        __n, __m = __InputCovOrEns.shape
         if __n != __m:
             raise ValueError("Additive inflation can only be applied to squared (covariance) matrix.")
         if __BackgroundCov is None:
@@ -679,10 +685,10 @@ def CovarianceInflation(
         __OutputCovOrEns = (1. - __InflationFactor) * __InputCovOrEns + __InflationFactor * __BackgroundCov
     #
     elif __InflationType == "Relaxation":
-        raise NotImplementedError("InflationType Relaxation")
+        raise NotImplementedError("Relaxation inflation type not implemented")
     #
     else:
-        raise ValueError("Error in inflation type, '%s' is not a valid keyword."%InflationType)
+        raise ValueError("Error in inflation type, '%s' is not a valid keyword."%__InflationType)
     #
     return __OutputCovOrEns
 
@@ -853,7 +859,7 @@ def Apply3DVarRecentringOnEnsemble(__EnXn, __EnXf, __Ynpu, __HO, __R, __B, __Bet
 # ==============================================================================
 def multiXOsteps(selfA, Xb, Y, U, HO, EM, CM, R, B, Q, oneCycle):
     """
-    Prévision multi-pas avec une correction par pas en X (multi-méthodes)
+    Prévision multi-pas avec une correction par pas (multi-méthodes)
     """
     #
     # Initialisation