methods, for those that require such a technique. Inflation can be applied in
various ways, according to the following options: multiplicative or additive
by the specified inflation factor, applied on the background or on the
- analysis, applied on covariances or on anomalies. Only one type of inflation
- is applied at the same time, and the default value is
+ analysis, applied on covariances or on anomalies. The *multiplicative
+ inflation on anomalies*, that are obtained by subtracting the ensemble mean,
+ is elaborated by multiplying these anomalies by the inflation factor, then by
+ rebuilding the ensemble members by adding the previously evaluated mean. Only
+ one type of inflation is applied at the same time, and the default value is
"MultiplicativeOnAnalysisAnomalies". The possible names are in the following
list: [
"MultiplicativeOnAnalysisAnomalies",
its method "*New()*" as illustrated in the following lines (the ``case`` object
name being let free to the user choice)::
- from numpy import array, matrix
+ from numpy import array
from adao import adaoBuilder
case = adaoBuilder.New()
Diagonal matrix, with either 1 or a given vector on the diagonal
"""
if diagonal is not None:
- S = numpy.diag( diagonal )
+ S = numpy.diagflat( diagonal )
else:
- S = numpy.matrix(numpy.identity(int(size)))
+ S = numpy.identity(int(size))
return S
We can then define the background state :math:`\mathbf{x}^b` as a random
""" Direct non-linear simulation operator """
#
# --------------------------------------> EXAMPLE TO BE REMOVED
- if type(XX) is type(numpy.matrix([])): # EXAMPLE TO BE REMOVED
- HX = XX.A1.tolist() # EXAMPLE TO BE REMOVED
- elif type(XX) is type(numpy.array([])): # EXAMPLE TO BE REMOVED
- HX = numpy.matrix(XX).A1.tolist() # EXAMPLE TO BE REMOVED
- else: # EXAMPLE TO BE REMOVED
- HX = XX # EXAMPLE TO BE REMOVED
+ HX = 1. * numpy.ravel( XX ) # EXAMPLE TO BE REMOVED
# --------------------------------------> EXAMPLE TO BE REMOVED
#
- return numpy.array( HX )
+ return HX
We does not need the linear companion operators ``"TangentOperator"`` and
``"AdjointOperator"`` because they will be approximated using ADAO
L'inflation peut être appliquée de diverses manières, selon les options
suivantes : multiplicative ou additive du facteur d'inflation spécifié,
appliquée sur l'ébauche ou sur l'analyse, appliquée sur les covariances ou
- sur les anomalies. Un seul type d'inflation est appliqué à la fois, et la
- valeur par défaut est "MultiplicativeOnAnalysisAnomalies". Les noms possibles
- sont dans la liste suivante : [
+ sur les anomalies. L'*inflation multiplicative sur les anomalies*, qui sont
+ obtenues en retranchant la moyenne d'ensemble, est effectuée en multipliant
+ ces anomalies par le facteur d'inflation, puis en reconstituant les membres
+ de l'ensemble par ajout de la moyenne préalablement calculée. Un seul type
+ d'inflation est appliqué à la fois, et la valeur par défaut est
+ "MultiplicativeOnAnalysisAnomalies". Les noms possibles sont dans la liste
+ suivante : [
"MultiplicativeOnAnalysisAnomalies",
"MultiplicativeOnBackgroundAnomalies",
].
méthode "*New()*" comme illustré dans les quelques lignes suivantes (le nom
``case`` de l'objet étant quelconque, au choix de l'utilisateur)::
- from numpy import array, matrix
+ from numpy import array
from adao import adaoBuilder
case = adaoBuilder.New()
Diagonal matrix, with either 1 or a given vector on the diagonal
"""
if diagonal is not None:
- S = numpy.diag( diagonal )
+ S = numpy.diagflat( diagonal )
else:
- S = numpy.matrix(numpy.identity(int(size)))
+ S = numpy.identity(int(size))
return S
On définit ensuite l'état d'ébauche :math:`\mathbf{x}^b` comme une perturbation
""" Direct non-linear simulation operator """
#
# --------------------------------------> EXAMPLE TO BE REMOVED
- if type(XX) is type(numpy.matrix([])): # EXAMPLE TO BE REMOVED
- HX = XX.A1.tolist() # EXAMPLE TO BE REMOVED
- elif type(XX) is type(numpy.array([])): # EXAMPLE TO BE REMOVED
- HX = numpy.matrix(XX).A1.tolist() # EXAMPLE TO BE REMOVED
- else: # EXAMPLE TO BE REMOVED
- HX = XX # EXAMPLE TO BE REMOVED
+ HX = 1. * numpy.ravel( XX ) # EXAMPLE TO BE REMOVED
# --------------------------------------> EXAMPLE TO BE REMOVED
#
- return numpy.array( HX )
+ return HX
On n'a pas besoin des opérateurs linéaires associés ``"TangentOperator"`` et
``"AdjointOperator"`` car ils vont être approximés en utilisant les capacités
# Verbosité et logging
if logging.getLogger().level < logging.WARNING:
self._parameters["optiprint"], self._parameters["optdisp"] = 1, 1
- if PlatformInfo.has_scipy:
- import scipy.optimize
- self._parameters["optmessages"] = scipy.optimize.tnc.MSG_ALL
- else:
- self._parameters["optmessages"] = 15
+ self._parameters["optmessages"] = 15
else:
self._parameters["optiprint"], self._parameters["optdisp"] = -1, 0
- if PlatformInfo.has_scipy:
- import scipy.optimize
- self._parameters["optmessages"] = scipy.optimize.tnc.MSG_NONE
- else:
- self._parameters["optmessages"] = 15
+ self._parameters["optmessages"] = 0
#
return 0
elif __Series is not None:
self.__is_series = True
if isinstance(__Series, (tuple, list, numpy.ndarray, numpy.matrix, str)):
- #~ self.__V = Persistence.OneVector(self.__name, basetype=numpy.matrix)
self.__V = Persistence.OneVector(self.__name)
if isinstance(__Series, str):
__Series = PlatformInfo.strmatrix2liststr(__Series)
# logging.debug("MULTF Internal multifonction calculations end")
return __multiHX
-# ==============================================================================
-def CostFunction3D(_x,
- _Hm = None, # Pour simuler Hm(x) : HO["Direct"].appliedTo
- _HmX = None, # Simulation déjà faite de Hm(x)
- _arg = None, # Arguments supplementaires pour Hm, sous la forme d'un tuple
- _BI = None,
- _RI = None,
- _Xb = None,
- _Y = None,
- _SIV = False, # A résorber pour la 8.0
- _SSC = [], # self._parameters["StoreSupplementaryCalculations"]
- _nPS = 0, # nbPreviousSteps
- _QM = "DA", # QualityMeasure
- _SSV = {}, # Entrée et/ou sortie : self.StoredVariables
- _fRt = False, # Restitue ou pas la sortie étendue
- _sSc = True, # Stocke ou pas les SSC
- ):
- """
- Fonction-coût générale utile pour les algorithmes statiques/3D : 3DVAR, BLUE
- et dérivés, Kalman et dérivés, LeastSquares, SamplingTest, PSO, SA, Tabu,
- DFO, QuantileRegression
- """
- if not _sSc:
- _SIV = False
- _SSC = {}
- else:
- for k in ["CostFunctionJ",
- "CostFunctionJb",
- "CostFunctionJo",
- "CurrentOptimum",
- "CurrentState",
- "IndexOfOptimum",
- "SimulatedObservationAtCurrentOptimum",
- "SimulatedObservationAtCurrentState",
- ]:
- if k not in _SSV:
- _SSV[k] = []
- if hasattr(_SSV[k],"store"):
- _SSV[k].append = _SSV[k].store # Pour utiliser "append" au lieu de "store"
- #
- _X = numpy.asmatrix(numpy.ravel( _x )).T
- if _SIV or "CurrentState" in _SSC or "CurrentOptimum" in _SSC:
- _SSV["CurrentState"].append( _X )
- #
- if _HmX is not None:
- _HX = _HmX
- else:
- if _Hm is None:
- raise ValueError("COSTFUNCTION3D Operator has to be defined.")
- if _arg is None:
- _HX = _Hm( _X )
- else:
- _HX = _Hm( _X, *_arg )
- _HX = numpy.asmatrix(numpy.ravel( _HX )).T
- #
- if "SimulatedObservationAtCurrentState" in _SSC or \
- "SimulatedObservationAtCurrentOptimum" in _SSC:
- _SSV["SimulatedObservationAtCurrentState"].append( _HX )
- #
- if numpy.any(numpy.isnan(_HX)):
- Jb, Jo, J = numpy.nan, numpy.nan, numpy.nan
- else:
- _Y = numpy.asmatrix(numpy.ravel( _Y )).T
- if _QM in ["AugmentedWeightedLeastSquares", "AWLS", "AugmentedPonderatedLeastSquares", "APLS", "DA"]:
- if _BI is None or _RI is None:
- raise ValueError("Background and Observation error covariance matrix has to be properly defined!")
- _Xb = numpy.asmatrix(numpy.ravel( _Xb )).T
- Jb = 0.5 * (_X - _Xb).T * _BI * (_X - _Xb)
- Jo = 0.5 * (_Y - _HX).T * _RI * (_Y - _HX)
- elif _QM in ["WeightedLeastSquares", "WLS", "PonderatedLeastSquares", "PLS"]:
- if _RI is None:
- raise ValueError("Observation error covariance matrix has to be properly defined!")
- Jb = 0.
- Jo = 0.5 * (_Y - _HX).T * _RI * (_Y - _HX)
- elif _QM in ["LeastSquares", "LS", "L2"]:
- Jb = 0.
- Jo = 0.5 * (_Y - _HX).T * (_Y - _HX)
- elif _QM in ["AbsoluteValue", "L1"]:
- Jb = 0.
- Jo = numpy.sum( numpy.abs(_Y - _HX) )
- elif _QM in ["MaximumError", "ME"]:
- Jb = 0.
- Jo = numpy.max( numpy.abs(_Y - _HX) )
- elif _QM in ["QR", "Null"]:
- Jb = 0.
- Jo = 0.
- else:
- raise ValueError("Unknown asked quality measure!")
- #
- J = float( Jb ) + float( Jo )
- #
- if _sSc:
- _SSV["CostFunctionJb"].append( Jb )
- _SSV["CostFunctionJo"].append( Jo )
- _SSV["CostFunctionJ" ].append( J )
- #
- if "IndexOfOptimum" in _SSC or \
- "CurrentOptimum" in _SSC or \
- "SimulatedObservationAtCurrentOptimum" in _SSC:
- IndexMin = numpy.argmin( _SSV["CostFunctionJ"][_nPS:] ) + _nPS
- if "IndexOfOptimum" in _SSC:
- _SSV["IndexOfOptimum"].append( IndexMin )
- if "CurrentOptimum" in _SSC:
- _SSV["CurrentOptimum"].append( _SSV["CurrentState"][IndexMin] )
- if "SimulatedObservationAtCurrentOptimum" in _SSC:
- _SSV["SimulatedObservationAtCurrentOptimum"].append( _SSV["SimulatedObservationAtCurrentState"][IndexMin] )
- #
- if _fRt:
- return _SSV
- else:
- if _QM in ["QR"]: # Pour le QuantileRegression
- return _HX
- else:
- return J
-
# ==============================================================================
if __name__ == "__main__":
print('\n AUTODIAGNOSTIC\n')
if __InflationType in ["MultiplicativeOnAnalysisCovariance", "MultiplicativeOnBackgroundCovariance"]:
if __InflationFactor < 1.:
raise ValueError("Inflation factor for multiplicative inflation has to be greater or equal than 1.")
- if __InflationFactor < 1.+mpr:
+ if __InflationFactor < 1.+mpr: # No inflation = 1
return __InputCovOrEns
__OutputCovOrEns = __InflationFactor**2 * __InputCovOrEns
#
elif __InflationType in ["MultiplicativeOnAnalysisAnomalies", "MultiplicativeOnBackgroundAnomalies"]:
if __InflationFactor < 1.:
raise ValueError("Inflation factor for multiplicative inflation has to be greater or equal than 1.")
- if __InflationFactor < 1.+mpr:
+ if __InflationFactor < 1.+mpr: # No inflation = 1
return __InputCovOrEns
__InputCovOrEnsMean = __InputCovOrEns.mean(axis=1, dtype=mfp).astype('float')
__OutputCovOrEns = __InputCovOrEnsMean[:,numpy.newaxis] \
elif __InflationType in ["AdditiveOnAnalysisCovariance", "AdditiveOnBackgroundCovariance"]:
if __InflationFactor < 0.:
raise ValueError("Inflation factor for additive inflation has to be greater or equal than 0.")
- if __InflationFactor < mpr:
+ if __InflationFactor < mpr: # No inflation = 0
return __InputCovOrEns
__n, __m = __InputCovOrEns.shape
if __n != __m:
elif __InflationType == "HybridOnBackgroundCovariance":
if __InflationFactor < 0.:
raise ValueError("Inflation factor for hybrid inflation has to be greater or equal than 0.")
- if __InflationFactor < mpr:
+ if __InflationFactor < mpr: # No inflation = 0
return __InputCovOrEns
__n, __m = __InputCovOrEns.shape
if __n != __m: