From 4790fb60acb36159350ee1cda40107e6833ead3f Mon Sep 17 00:00:00 2001 From: Jean-Philippe ARGAUD Date: Wed, 21 Jun 2017 15:45:58 +0200 Subject: [PATCH] Internal tests and warning improvements --- src/daComposant/daAlgorithms/3DVAR.py | 5 +++- src/daComposant/daAlgorithms/4DVAR.py | 6 ++++- src/daComposant/daAlgorithms/AdjointTest.py | 13 +++++---- src/daComposant/daAlgorithms/Blue.py | 5 +++- .../DerivativeFreeOptimization.py | 6 ++++- src/daComposant/daAlgorithms/EnsembleBlue.py | 5 +++- src/daComposant/daAlgorithms/ExtendedBlue.py | 5 +++- .../daAlgorithms/ExtendedKalmanFilter.py | 6 ++++- src/daComposant/daAlgorithms/FunctionTest.py | 5 +++- src/daComposant/daAlgorithms/GradientTest.py | 22 +++++++-------- src/daComposant/daAlgorithms/KalmanFilter.py | 11 ++++---- .../daAlgorithms/LinearLeastSquares.py | 5 +++- src/daComposant/daAlgorithms/LinearityTest.py | 27 +++++++++---------- .../daAlgorithms/NonLinearLeastSquares.py | 5 +++- src/daComposant/daAlgorithms/ObserverTest.py | 2 +- .../daAlgorithms/ParticleSwarmOptimization.py | 5 +++- .../daAlgorithms/QuantileRegression.py | 5 +++- src/daComposant/daAlgorithms/SamplingTest.py | 5 +++- src/daComposant/daAlgorithms/TabuSearch.py | 5 +++- src/daComposant/daAlgorithms/TangentTest.py | 12 +++++---- .../daAlgorithms/UnscentedKalmanFilter.py | 11 ++++---- src/daComposant/daCore/BasicObjects.py | 26 ++++++++++++++++-- 22 files changed, 133 insertions(+), 64 deletions(-) diff --git a/src/daComposant/daAlgorithms/3DVAR.py b/src/daComposant/daAlgorithms/3DVAR.py index 67a060c..1424005 100644 --- a/src/daComposant/daAlgorithms/3DVAR.py +++ b/src/daComposant/daAlgorithms/3DVAR.py @@ -105,9 +105,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # # Correction pour pallier a un bug de TNC sur le retour du Minimum if "Minimizer" in self._parameters and self._parameters["Minimizer"] == "TNC": diff --git a/src/daComposant/daAlgorithms/4DVAR.py b/src/daComposant/daAlgorithms/4DVAR.py index a9f9fb1..ebeea0b 100644 --- a/src/daComposant/daAlgorithms/4DVAR.py +++ b/src/daComposant/daAlgorithms/4DVAR.py @@ -92,9 +92,13 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "EM", "R", "B" ), + optional = ("U", "CM"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # # Correction pour pallier a un bug de TNC sur le retour du Minimum if "Minimizer" in self._parameters and self._parameters["Minimizer"] == "TNC": diff --git a/src/daComposant/daAlgorithms/AdjointTest.py b/src/daComposant/daAlgorithms/AdjointTest.py index d49c5c3..888e1c5 100644 --- a/src/daComposant/daAlgorithms/AdjointTest.py +++ b/src/daComposant/daAlgorithms/AdjointTest.py @@ -76,9 +76,13 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["CurrentState", "Residu", "SimulatedObservationAtCurrentState"] ) + self.requireInputArguments( + mandatory= ("Xb", "HO" ), + optional = ("Y", ), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Direct"].appliedTo Ht = HO["Tangent"].appliedInXTo @@ -118,7 +122,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): Remarque : les nombres inferieurs a %.0e (environ) representent un zero a la precision machine.\n"""%mpr if self._parameters["ResiduFormula"] == "ScalarProduct": - __entete = u" i Alpha ||X|| ||Y|| ||dX|| R(Alpha) " + __entete = u" i Alpha ||X|| ||Y|| ||dX|| R(Alpha)" __msgdoc = u""" On observe le residu qui est la difference de deux produits scalaires : @@ -126,8 +130,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): qui doit rester constamment egal a zero a la precision du calcul. On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - Y doit etre dans l'image de F. S'il n'est pas donne, on prend Y = F(X). - """ + __precision + Y doit etre dans l'image de F. S'il n'est pas donne, on prend Y = F(X).\n""" + __precision # if len(self._parameters["ResultTitle"]) > 0: __rt = unicode(self._parameters["ResultTitle"]) @@ -139,7 +142,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): msgs = u"" msgs += __msgdoc # - __nbtirets = len(__entete) + __nbtirets = len(__entete) + 2 msgs += "\n" + __marge + "-"*__nbtirets msgs += "\n" + __marge + __entete msgs += "\n" + __marge + "-"*__nbtirets diff --git a/src/daComposant/daAlgorithms/Blue.py b/src/daComposant/daAlgorithms/Blue.py index 2d85770..d621343 100644 --- a/src/daComposant/daAlgorithms/Blue.py +++ b/src/daComposant/daAlgorithms/Blue.py @@ -68,9 +68,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Type de simulation pour l'estimation des quantiles", listval = ["Linear", "NonLinear"] ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Tangent"].asMatrix(Xb) Hm = Hm.reshape(Y.size,Xb.size) # ADAO & check shape diff --git a/src/daComposant/daAlgorithms/DerivativeFreeOptimization.py b/src/daComposant/daAlgorithms/DerivativeFreeOptimization.py index eb35174..4781070 100644 --- a/src/daComposant/daAlgorithms/DerivativeFreeOptimization.py +++ b/src/daComposant/daAlgorithms/DerivativeFreeOptimization.py @@ -89,11 +89,15 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # if not PlatformInfo.has_nlopt and not self._parameters["Minimizer"] in ["COBYLA", "POWELL", "SIMPLEX"]: + logging.debug("%s Absence de NLopt, utilisation forcee du minimiseur SIMPLEX"%(self._name,)) self._parameters["Minimizer"] = "SIMPLEX" # # Opérateurs diff --git a/src/daComposant/daAlgorithms/EnsembleBlue.py b/src/daComposant/daAlgorithms/EnsembleBlue.py index 0157578..7c085a3 100644 --- a/src/daComposant/daAlgorithms/EnsembleBlue.py +++ b/src/daComposant/daAlgorithms/EnsembleBlue.py @@ -46,9 +46,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): typecast = numpy.random.seed, message = "Graine fixée pour le générateur aléatoire", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # # Précalcul des inversions de B et R # ---------------------------------- diff --git a/src/daComposant/daAlgorithms/ExtendedBlue.py b/src/daComposant/daAlgorithms/ExtendedBlue.py index 893a099..9271629 100644 --- a/src/daComposant/daAlgorithms/ExtendedBlue.py +++ b/src/daComposant/daAlgorithms/ExtendedBlue.py @@ -68,9 +68,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Type de simulation pour l'estimation des quantiles", listval = ["Linear", "NonLinear"] ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Tangent"].asMatrix(Xb) Hm = Hm.reshape(Y.size,Xb.size) # ADAO & check shape diff --git a/src/daComposant/daAlgorithms/ExtendedKalmanFilter.py b/src/daComposant/daAlgorithms/ExtendedKalmanFilter.py index db3f9ab..66a6fcd 100644 --- a/src/daComposant/daAlgorithms/ExtendedKalmanFilter.py +++ b/src/daComposant/daAlgorithms/ExtendedKalmanFilter.py @@ -59,9 +59,13 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + optional = ("U", "EM", "CM", "Q"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # if self._parameters["EstimationOf"] == "Parameters": self._parameters["StoreInternalVariables"] = True diff --git a/src/daComposant/daAlgorithms/FunctionTest.py b/src/daComposant/daAlgorithms/FunctionTest.py index 6e0f869..4b1adbd 100644 --- a/src/daComposant/daAlgorithms/FunctionTest.py +++ b/src/daComposant/daAlgorithms/FunctionTest.py @@ -65,9 +65,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["CurrentState", "SimulatedObservationAtCurrentState"] ) + self.requireInputArguments( + mandatory= ("Xb", "HO"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Direct"].appliedTo # diff --git a/src/daComposant/daAlgorithms/GradientTest.py b/src/daComposant/daAlgorithms/GradientTest.py index 7be3df4..109088e 100644 --- a/src/daComposant/daAlgorithms/GradientTest.py +++ b/src/daComposant/daAlgorithms/GradientTest.py @@ -102,9 +102,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["CurrentState", "Residu", "SimulatedObservationAtCurrentState"] ) + self.requireInputArguments( + mandatory= ("Xb", "HO"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Direct"].appliedTo if self._parameters["ResiduFormula"] in ["Taylor", "TaylorOnNorm"]: @@ -148,7 +151,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): Remarque : les nombres inferieurs a %.0e (environ) representent un zero a la precision machine.\n"""%mpr if self._parameters["ResiduFormula"] == "Taylor": - __entete = u" i Alpha ||X|| ||F(X)|| ||F(X+dX)|| ||dX|| ||F(X+dX)-F(X)|| ||F(X+dX)-F(X)||/||dX|| R(Alpha) log( R ) " + __entete = u" i Alpha ||X|| ||F(X)|| ||F(X+dX)|| ||dX|| ||F(X+dX)-F(X)|| ||F(X+dX)-F(X)||/||dX|| R(Alpha) log( R )" __msgdoc = u""" On observe le residu issu du developpement de Taylor de la fonction F, normalise par la valeur au point nominal : @@ -166,10 +169,9 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): signifie que F est lineaire et que le residu decroit a partir de l'erreur faite dans le calcul du terme GradientF_X. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision if self._parameters["ResiduFormula"] == "TaylorOnNorm": - __entete = u" i Alpha ||X|| ||F(X)|| ||F(X+dX)|| ||dX|| ||F(X+dX)-F(X)|| ||F(X+dX)-F(X)||/||dX|| R(Alpha) log( R ) " + __entete = u" i Alpha ||X|| ||F(X)|| ||F(X+dX)|| ||dX|| ||F(X+dX)-F(X)|| ||F(X+dX)-F(X)||/||dX|| R(Alpha) log( R )" __msgdoc = u""" On observe le residu issu du developpement de Taylor de la fonction F, rapporte au parametre Alpha au carre : @@ -191,10 +193,9 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): calcul du gradient est correct jusqu'au moment ou le residu est de l'ordre de grandeur de ||F(X)||. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision if self._parameters["ResiduFormula"] == "Norm": - __entete = u" i Alpha ||X|| ||F(X)|| ||F(X+dX)|| ||dX|| ||F(X+dX)-F(X)|| ||F(X+dX)-F(X)||/||dX|| R(Alpha) log( R ) " + __entete = u" i Alpha ||X|| ||F(X)|| ||F(X+dX)|| ||dX|| ||F(X+dX)-F(X)|| ||F(X+dX)-F(X)||/||dX|| R(Alpha) log( R )" __msgdoc = u""" On observe le residu, qui est base sur une approximation du gradient : @@ -204,8 +205,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): qui doit rester constant jusqu'a ce que l'on atteigne la precision du calcul. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision # if len(self._parameters["ResultTitle"]) > 0: __rt = unicode(self._parameters["ResultTitle"]) @@ -217,7 +217,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): msgs = u"" msgs += __msgdoc # - __nbtirets = len(__entete) + __nbtirets = len(__entete) + 2 msgs += "\n" + __marge + "-"*__nbtirets msgs += "\n" + __marge + __entete msgs += "\n" + __marge + "-"*__nbtirets diff --git a/src/daComposant/daAlgorithms/KalmanFilter.py b/src/daComposant/daAlgorithms/KalmanFilter.py index 3ff5e59..50546eb 100644 --- a/src/daComposant/daAlgorithms/KalmanFilter.py +++ b/src/daComposant/daAlgorithms/KalmanFilter.py @@ -48,20 +48,19 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["APosterioriCorrelations", "APosterioriCovariance", "APosterioriStandardDeviations", "APosterioriVariances", "BMA", "CurrentState", "CostFunctionJ", "CostFunctionJb", "CostFunctionJo", "Innovation"] ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + optional = ("U", "EM", "CM", "Q"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # if self._parameters["EstimationOf"] == "Parameters": self._parameters["StoreInternalVariables"] = True # # Opérateurs # ---------- - if B is None: - raise ValueError("Background error covariance matrix has to be properly defined!") - if R is None: - raise ValueError("Observation error covariance matrix has to be properly defined!") - # Ht = HO["Tangent"].asMatrix(Xb) Ha = HO["Adjoint"].asMatrix(Xb) # diff --git a/src/daComposant/daAlgorithms/LinearLeastSquares.py b/src/daComposant/daAlgorithms/LinearLeastSquares.py index 6fc300d..e2525dd 100644 --- a/src/daComposant/daAlgorithms/LinearLeastSquares.py +++ b/src/daComposant/daAlgorithms/LinearLeastSquares.py @@ -41,9 +41,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["OMA", "CurrentState", "CostFunctionJ", "CostFunctionJb", "CostFunctionJo", "SimulatedObservationAtCurrentState", "SimulatedObservationAtOptimum"] ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Tangent"].asMatrix(None) Hm = Hm.reshape(Y.size,-1) # ADAO & check shape diff --git a/src/daComposant/daAlgorithms/LinearityTest.py b/src/daComposant/daAlgorithms/LinearityTest.py index 9733aa9..c2464c2 100644 --- a/src/daComposant/daAlgorithms/LinearityTest.py +++ b/src/daComposant/daAlgorithms/LinearityTest.py @@ -84,9 +84,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["CurrentState", "Residu", "SimulatedObservationAtCurrentState"] ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # def RMS(V1, V2): import math @@ -143,7 +146,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): Remarque : les nombres inferieurs a %.0e (environ) representent un zero a la precision machine.\n"""%mpr if self._parameters["ResiduFormula"] == "CenteredDL": - __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) log10( R ) " + __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) log10( R )" __msgdoc = u""" On observe le residu provenant de la difference centree des valeurs de F au point nominal et aux points perturbes, normalisee par la valeur au @@ -164,10 +167,9 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): cela signifie que le gradient est calculable jusqu'a la precision d'arret de la decroissance quadratique. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision if self._parameters["ResiduFormula"] == "Taylor": - __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) log10( R ) " + __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) log10( R )" __msgdoc = u""" On observe le residu issu du developpement de Taylor de la fonction F, normalisee par la valeur au point nominal : @@ -187,10 +189,9 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): cela signifie que le gradient est bien calcule jusqu'a la precision d'arret de la decroissance quadratique. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision if self._parameters["ResiduFormula"] == "NominalTaylor": - __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) |R-1| en % " + __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) |R-1| en %" __msgdoc = u""" On observe le residu obtenu a partir de deux approximations d'ordre 1 de F(X), normalisees par la valeur au point nominal : @@ -208,10 +209,9 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): l'increment Alpha, c'est sur cette partie que l'hypothese de linearite de F est verifiee. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision if self._parameters["ResiduFormula"] == "NominalTaylorRMS": - __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) |R| en % " + __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) |R| en %" __msgdoc = u""" On observe le residu obtenu a partir de deux approximations d'ordre 1 de F(X), normalisees par la valeur au point nominal : @@ -228,8 +228,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): l'increment Alpha, c'est sur cette partie que l'hypothese de linearite de F est verifiee. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision # if len(self._parameters["ResultTitle"]) > 0: __rt = unicode(self._parameters["ResultTitle"]) @@ -241,7 +240,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): msgs = u"" msgs += __msgdoc # - __nbtirets = len(__entete) + __nbtirets = len(__entete) + 2 msgs += "\n" + __marge + "-"*__nbtirets msgs += "\n" + __marge + __entete msgs += "\n" + __marge + "-"*__nbtirets diff --git a/src/daComposant/daAlgorithms/NonLinearLeastSquares.py b/src/daComposant/daAlgorithms/NonLinearLeastSquares.py index a2e4433..5ea90cf 100644 --- a/src/daComposant/daAlgorithms/NonLinearLeastSquares.py +++ b/src/daComposant/daAlgorithms/NonLinearLeastSquares.py @@ -78,9 +78,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # # Correction pour pallier a un bug de TNC sur le retour du Minimum if "Minimizer" in self._parameters and self._parameters["Minimizer"] == "TNC": diff --git a/src/daComposant/daAlgorithms/ObserverTest.py b/src/daComposant/daAlgorithms/ObserverTest.py index 89e0d5d..4ae19ed 100644 --- a/src/daComposant/daAlgorithms/ObserverTest.py +++ b/src/daComposant/daAlgorithms/ObserverTest.py @@ -30,7 +30,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): BasicObjects.Algorithm.__init__(self, "OBSERVERTEST") def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) print("Results of observer check on all potential variables or commands,") print(" only activated on selected ones by explicit association.") print("") diff --git a/src/daComposant/daAlgorithms/ParticleSwarmOptimization.py b/src/daComposant/daAlgorithms/ParticleSwarmOptimization.py index e288583..6440817 100644 --- a/src/daComposant/daAlgorithms/ParticleSwarmOptimization.py +++ b/src/daComposant/daAlgorithms/ParticleSwarmOptimization.py @@ -97,9 +97,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "BoxBounds", message = "Liste des valeurs de bornes d'incréments de paramètres", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # if ("BoxBounds" in self._parameters) and isinstance(self._parameters["BoxBounds"], (list, tuple)) and (len(self._parameters["BoxBounds"]) > 0): BoxBounds = self._parameters["BoxBounds"] diff --git a/src/daComposant/daAlgorithms/QuantileRegression.py b/src/daComposant/daAlgorithms/QuantileRegression.py index b971f07..0d0cf2a 100644 --- a/src/daComposant/daAlgorithms/QuantileRegression.py +++ b/src/daComposant/daAlgorithms/QuantileRegression.py @@ -73,9 +73,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO" ), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Direct"].appliedTo # diff --git a/src/daComposant/daAlgorithms/SamplingTest.py b/src/daComposant/daAlgorithms/SamplingTest.py index 1d60f73..c94c97a 100644 --- a/src/daComposant/daAlgorithms/SamplingTest.py +++ b/src/daComposant/daAlgorithms/SamplingTest.py @@ -81,9 +81,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): typecast = numpy.random.seed, message = "Graine fixée pour le générateur aléatoire", ) + self.requireInputArguments( + mandatory= ("Xb", "HO"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Direct"].appliedTo # diff --git a/src/daComposant/daAlgorithms/TabuSearch.py b/src/daComposant/daAlgorithms/TabuSearch.py index 76e42b4..49183b0 100644 --- a/src/daComposant/daAlgorithms/TabuSearch.py +++ b/src/daComposant/daAlgorithms/TabuSearch.py @@ -109,9 +109,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # if self._parameters["NoiseDistribution"] == "Uniform": nrange = numpy.ravel(self._parameters["NoiseHalfRange"]) # Vecteur diff --git a/src/daComposant/daAlgorithms/TangentTest.py b/src/daComposant/daAlgorithms/TangentTest.py index cb01499..bff1e73 100644 --- a/src/daComposant/daAlgorithms/TangentTest.py +++ b/src/daComposant/daAlgorithms/TangentTest.py @@ -84,9 +84,12 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): message = "Liste de calculs supplémentaires à stocker et/ou effectuer", listval = ["CurrentState", "Residu", "SimulatedObservationAtCurrentState"] ) + self.requireInputArguments( + mandatory= ("Xb", "HO"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # Hm = HO["Direct"].appliedTo Ht = HO["Tangent"].appliedInXTo @@ -137,7 +140,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): Remarque : les nombres inferieurs a %.0e (environ) representent un zero a la precision machine.\n"""%mpr if self._parameters["ResiduFormula"] == "Taylor": - __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) |R-1|/Alpha " + __entete = u" i Alpha ||X|| ||F(X)|| | R(Alpha) |R-1|/Alpha" __msgdoc = u""" On observe le residu provenant du rapport d'increments utilisant le lineaire tangent : @@ -157,8 +160,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): lineaire ou quasi-lineaire, et le tangent est valide jusqu'a ce que l'on atteigne la precision du calcul. - On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul. - """ + __precision + On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.\n""" + __precision # if len(self._parameters["ResultTitle"]) > 0: __rt = unicode(self._parameters["ResultTitle"]) @@ -170,7 +172,7 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): msgs = u"" msgs += __msgdoc # - __nbtirets = len(__entete) + __nbtirets = len(__entete) + 2 msgs += "\n" + __marge + "-"*__nbtirets msgs += "\n" + __marge + __entete msgs += "\n" + __marge + "-"*__nbtirets diff --git a/src/daComposant/daAlgorithms/UnscentedKalmanFilter.py b/src/daComposant/daAlgorithms/UnscentedKalmanFilter.py index c04fd1d..19bc20d 100644 --- a/src/daComposant/daAlgorithms/UnscentedKalmanFilter.py +++ b/src/daComposant/daAlgorithms/UnscentedKalmanFilter.py @@ -88,9 +88,13 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): name = "Bounds", message = "Liste des valeurs de bornes", ) + self.requireInputArguments( + mandatory= ("Xb", "Y", "HO", "R", "B" ), + optional = ("U", "EM", "CM", "Q"), + ) def run(self, Xb=None, Y=None, U=None, HO=None, EM=None, CM=None, R=None, B=None, Q=None, Parameters=None): - self._pre_run(Parameters) + self._pre_run(Parameters, R, B, Q) # if self._parameters["EstimationOf"] == "Parameters": self._parameters["StoreInternalVariables"] = True @@ -120,11 +124,6 @@ class ElementaryAlgorithm(BasicObjects.Algorithm): # # Opérateurs # ---------- - if B is None: - raise ValueError("Background error covariance matrix has to be properly defined!") - if R is None: - raise ValueError("Observation error covariance matrix has to be properly defined!") - # H = HO["Direct"].appliedControledFormTo # if self._parameters["EstimationOf"] == "State": diff --git a/src/daComposant/daCore/BasicObjects.py b/src/daComposant/daCore/BasicObjects.py index 07f1cc6..4404d2e 100644 --- a/src/daComposant/daCore/BasicObjects.py +++ b/src/daComposant/daCore/BasicObjects.py @@ -492,8 +492,9 @@ class Algorithm(object): self._name = str( name ) self._parameters = {"StoreSupplementaryCalculations":[]} self.__required_parameters = {} - self.StoredVariables = {} + self.__required_inputs = {"RequiredInputValues":{"mandatory":(), "optional":()}} # + self.StoredVariables = {} self.StoredVariables["CostFunctionJ"] = Persistence.OneScalar(name = "CostFunctionJ") self.StoredVariables["CostFunctionJb"] = Persistence.OneScalar(name = "CostFunctionJb") self.StoredVariables["CostFunctionJo"] = Persistence.OneScalar(name = "CostFunctionJo") @@ -526,7 +527,7 @@ class Algorithm(object): self.StoredVariables["SimulationQuantiles"] = Persistence.OneMatrix(name = "SimulationQuantiles") self.StoredVariables["Residu"] = Persistence.OneScalar(name = "Residu") - def _pre_run(self, Parameters ): + def _pre_run(self, Parameters, R=None, B=None, Q=None ): "Pré-calcul" logging.debug("%s Lancement", self._name) logging.debug("%s Taille mémoire utilisée de %.0f Mio", self._name, self._m.getUsedMemory("Mio")) @@ -535,6 +536,20 @@ class Algorithm(object): self.__setParameters(Parameters) # # Corrections et complements + def __test_cvalue( argument, variable, argname): + if argument is None: + if variable in self.__required_inputs["RequiredInputValues"]["mandatory"]: + raise ValueError("%s %s error covariance matrix %s has to be properly defined!"%(self._name,argname,variable)) + elif variable in self.__required_inputs["RequiredInputValues"]["optional"]: + logging.debug("%s %s error covariance matrix %s is not set, but is optional."%(self._name,argname,variable)) + else: + logging.debug("%s %s error covariance matrix %s is not set, but is not required."%(self._name,argname,variable)) + else: + logging.debug("%s %s error covariance matrix %s is set."%(self._name,argname,variable)) + __test_cvalue( R, "R", "Observation" ) + __test_cvalue( B, "B", "Background" ) + __test_cvalue( Q, "Q", "Evolution" ) + # if ("Bounds" in self._parameters) and isinstance(self._parameters["Bounds"], (list, tuple)) and (len(self._parameters["Bounds"]) > 0): logging.debug("%s Prise en compte des bornes effectuee"%(self._name,)) else: @@ -683,6 +698,13 @@ class Algorithm(object): raise ValueError("The value \"%s\" of the parameter named \"%s\" is not allowed, it has to be in the list %s."%( __val, name,listval)) return __val + def requireInputArguments(self, mandatory=(), optional=()): + """ + Permet d'imposer des arguments requises en entrée + """ + self.__required_inputs["RequiredInputValues"]["mandatory"] = tuple( mandatory ) + self.__required_inputs["RequiredInputValues"]["optional"] = tuple( optional ) + def __setParameters(self, fromDico={}): """ Permet de stocker les paramètres reçus dans le dictionnaire interne. -- 2.39.2