]> SALOME platform Git repositories - modules/adao.git/commitdiff
Salome HOME
Adding tangent test algorithm
authorJean-Philippe ARGAUD <jean-philippe.argaud@edf.fr>
Mon, 5 May 2014 06:42:30 +0000 (08:42 +0200)
committerJean-Philippe ARGAUD <jean-philippe.argaud@edf.fr>
Mon, 5 May 2014 06:42:30 +0000 (08:42 +0200)
doc/en/reference.rst
doc/fr/reference.rst
src/daComposant/daAlgorithms/TangentTest.py [new file with mode: 0644]
src/daSalome/daYacsSchemaCreator/infos_daComposant.py

index 0861af85038c7011eec6a6d5d5091239980afef4..0b89bc0016f52694cd47e056ea575cc0a5f32094 100644 (file)
@@ -759,6 +759,7 @@ Optional and required commands for checking algorithms
 .. index:: single: GradientTest
 .. index:: single: LinearityTest
 .. index:: single: ObserverTest
+.. index:: single: TangentTest
 
 .. index:: single: AlgorithmParameters
 .. index:: single: AmplitudeOfInitialDirection
@@ -924,6 +925,36 @@ for an ADAO checking case`_.
     *"ObserverTest"*, and will not be used. The easiest way is to give "1" as a
     STRING for both, *"ObservationOperator"* having to be of type *Matrix*.
 
+**"TangentTest"**
+
+  *Required commands*
+    *"CheckingPoint",
+    "ObservationOperator"*
+
+  AmplitudeOfInitialDirection
+    This key indicates the scaling of the initial perturbation build as a vector
+    used for the directional derivative around the nominal checking point. The
+    default is 1, that means no scaling.
+
+  EpsilonMinimumExponent
+    This key indicates the minimal exponent value of the power of 10 coefficient
+    to be used to decrease the increment multiplier. The default is -8, and it
+    has to be between 0 and -20. For example, its default value leads to
+    calculate the residue of the scalar product formula with a fixed increment
+    multiplied from 1.e0 to 1.e-8.
+
+  InitialDirection
+    This key indicates the vector direction used for the directional derivative
+    around the nominal checking point. It has to be a vector. If not specified,
+    this direction defaults to a random perturbation around zero of the same
+    vector size than the checking point.
+
+  SetSeed
+    This key allow to give an integer in order to fix the seed of the random
+    generator used to generate the ensemble. A convenient value is for example
+    1000. By default, the seed is left uninitialized, and so use the default
+    initialization from the computer.
+
 Requirements for functions describing an operator
 -------------------------------------------------
 
index 0993e46a55af323e6976bcf319f89e349f9a8859..5fbaf736487ab417a93455f1668c116b9fc63073 100644 (file)
@@ -796,6 +796,7 @@ Commandes optionnelles et requises pour les algorithmes de v
 .. index:: single: GradientTest
 .. index:: single: LinearityTest
 .. index:: single: ObserverTest
+.. index:: single: TangentTest
 
 .. index:: single: AlgorithmParameters
 .. index:: single: AmplitudeOfInitialDirection
@@ -969,6 +970,38 @@ commandes et mots-cl
     plus simple est de donner "1" comme un STRING pour les deux,
     l'*"ObservationOperator"* devant être de type *Matrix*.
 
+**"TangentTest"**
+
+  *Commandes obligatoires*
+    *"CheckingPoint",
+    "ObservationOperator"*
+
+  AmplitudeOfInitialDirection
+    Cette clé indique la mise à l'échelle de la perturbation initiale construite
+    comme un vecteur utilisé pour la dérivée directionnelle autour du point
+    nominal de vérification. La valeur par défaut est de 1, ce qui signifie pas
+    de mise à l'échelle.
+
+  EpsilonMinimumExponent
+    Cette clé indique la valeur de l'exposant minimal du coefficient en
+    puissance de 10 qui doit être utilisé pour faire décroître le multiplicateur
+    de l'incrément. La valeur par défaut est de -8, et elle doit être entre 0 et
+    -20. Par exemple, la valeur par défaut conduit à calculer le résidu de la
+    formule avec un incrément fixe multiplié par 1.e0 jusqu'à 1.e-8.
+
+  InitialDirection
+    Cette clé indique la direction vectorielle utilisée pour la dérivée
+    directionnelle autour du point nominal de vérification. Cela doit être un
+    vecteur. Si elle n'est pas spécifiée, la direction par défaut est une
+    perturbation par défaut autour de zero de la même taille vectorielle que le
+    point de vérification.
+
+  SetSeed
+    Cette clé permet de donner un nombre entier pour fixer la graine du
+    générateur aléatoire utilisé pour générer l'ensemble. Un valeur pratique est
+    par exemple 1000. Par défaut, la graine est laissée non initialisée, et elle
+    utilise ainsi l'initialisation par défaut de l'ordinateur.
+
 Exigences pour les fonctions décrivant un opérateur
 ---------------------------------------------------
 
diff --git a/src/daComposant/daAlgorithms/TangentTest.py b/src/daComposant/daAlgorithms/TangentTest.py
new file mode 100644 (file)
index 0000000..1106e10
--- /dev/null
@@ -0,0 +1,184 @@
+#-*-coding:iso-8859-1-*-
+#
+#  Copyright (C) 2008-2014 EDF R&D
+#
+#  This library is free software; you can redistribute it and/or
+#  modify it under the terms of the GNU Lesser General Public
+#  License as published by the Free Software Foundation; either
+#  version 2.1 of the License.
+#
+#  This library is distributed in the hope that it will be useful,
+#  but WITHOUT ANY WARRANTY; without even the implied warranty of
+#  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU
+#  Lesser General Public License for more details.
+#
+#  You should have received a copy of the GNU Lesser General Public
+#  License along with this library; if not, write to the Free Software
+#  Foundation, Inc., 59 Temple Place, Suite 330, Boston, MA  02111-1307 USA
+#
+#  See http://www.salome-platform.org/ or email : webmaster.salome@opencascade.com
+#
+#  Author: Jean-Philippe Argaud, jean-philippe.argaud@edf.fr, EDF R&D
+
+import logging
+from daCore import BasicObjects
+import numpy, math
+
+# ==============================================================================
+class ElementaryAlgorithm(BasicObjects.Algorithm):
+    def __init__(self):
+        BasicObjects.Algorithm.__init__(self, "TANGENTTEST")
+        self.defineRequiredParameter(
+            name     = "ResiduFormula",
+            default  = "Taylor",
+            typecast = str,
+            message  = "Formule de résidu utilisée",
+            listval  = ["Taylor"],
+            )
+        self.defineRequiredParameter(
+            name     = "EpsilonMinimumExponent",
+            default  = -8,
+            typecast = int,
+            message  = "Exposant minimal en puissance de 10 pour le multiplicateur d'incrément",
+            minval   = -20,
+            maxval   = 0,
+            )
+        self.defineRequiredParameter(
+            name     = "InitialDirection",
+            default  = [],
+            typecast = list,
+            message  = "Direction initiale de la dérivée directionnelle autour du point nominal",
+            )
+        self.defineRequiredParameter(
+            name     = "AmplitudeOfInitialDirection",
+            default  = 1.,
+            typecast = float,
+            message  = "Amplitude de la direction initiale de la dérivée directionnelle autour du point nominal",
+            )
+        self.defineRequiredParameter(
+            name     = "SetSeed",
+            typecast = numpy.random.seed,
+            message  = "Graine fixée pour le générateur aléatoire",
+            )
+        self.defineRequiredParameter(
+            name     = "ResultTitle",
+            default  = "",
+            typecast = str,
+            message  = "Titre du tableau et de la figure",
+            )
+
+    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()
+        #
+        # Paramètres de pilotage
+        # ----------------------
+        self.setParameters(Parameters)
+        #
+        # Opérateurs
+        # ----------
+        Hm = HO["Direct"].appliedTo
+        Ht = HO["Tangent"].appliedInXTo
+        #
+        # Construction des perturbations
+        # ------------------------------
+        Perturbations = [ 10**i for i in xrange(self._parameters["EpsilonMinimumExponent"],1) ]
+        Perturbations.reverse()
+        #
+        # Calcul du point courant
+        # -----------------------
+        Xn      = numpy.asmatrix(numpy.ravel( Xb )).T
+        FX      = numpy.asmatrix(numpy.ravel( Hm( Xn ) )).T
+        NormeX  = numpy.linalg.norm( Xn )
+        NormeFX = numpy.linalg.norm( FX )
+        #
+        # Fabrication de la direction de  l'incrément dX
+        # ----------------------------------------------
+        if len(self._parameters["InitialDirection"]) == 0:
+            dX0 = []
+            for v in Xn.A1:
+                if abs(v) > 1.e-8:
+                    dX0.append( numpy.random.normal(0.,abs(v)) )
+                else:
+                    dX0.append( numpy.random.normal(0.,Xn.mean()) )
+        else:
+            dX0 = numpy.ravel( self._parameters["InitialDirection"] )
+        #
+        dX0 = float(self._parameters["AmplitudeOfInitialDirection"]) * numpy.matrix( dX0 ).T
+        #
+        # Calcul du gradient au point courant X pour l'incrément dX
+        # qui est le tangent en X multiplié par dX
+        # ---------------------------------------------------------
+        GradFxdX = Ht( (Xn, dX0) )
+        GradFxdX = numpy.asmatrix(numpy.ravel( GradFxdX )).T
+        NormeGX  = numpy.linalg.norm( GradFxdX )
+        #
+        # Entete des resultats
+        # --------------------
+        __marge =  12*" "
+        if self._parameters["ResiduFormula"] == "Taylor":
+            __entete = "  i   Alpha     ||X||      ||F(X)||   |     R(Alpha)    |R-1|/Alpha  "
+            __msgdoc = """
+            On observe le résidu provenant du rapport d'incréments utilisant le
+            linéaire tangent :
+
+                          || F(X+Alpha*dX) - F(X) ||
+              R(Alpha) = -----------------------------
+                         || Alpha * TangentF_X * dX ||
+
+            qui doit rester stable en 1+O(Alpha) jusqu'à ce que l'on atteigne la
+            précision du calcul.
+            
+            Lorsque |R-1|/Alpha est inférieur ou égal à une valeur stable
+            lorsque Alpha varie, le tangent est valide, jusqu'à ce que l'on
+            atteigne la précision du calcul.
+            
+            Si |R-1|/Alpha est très faible, le code F est vraisemblablement
+            linéaire ou quasi-linéaire, et le tangent est valide jusqu'à ce que
+            l'on atteigne la précision du calcul.
+
+            On prend dX0 = Normal(0,X) et dX = Alpha*dX0. F est le code de calcul.
+            """
+        #
+        if len(self._parameters["ResultTitle"]) > 0:
+            msgs  = "\n"
+            msgs += __marge + "====" + "="*len(self._parameters["ResultTitle"]) + "====\n"
+            msgs += __marge + "    " + self._parameters["ResultTitle"] + "\n"
+            msgs += __marge + "====" + "="*len(self._parameters["ResultTitle"]) + "====\n"
+        else:
+            msgs  = ""
+        msgs += __msgdoc
+        #
+        __nbtirets = len(__entete)
+        msgs += "\n" + __marge + "-"*__nbtirets
+        msgs += "\n" + __marge + __entete
+        msgs += "\n" + __marge + "-"*__nbtirets
+        #
+        # Boucle sur les perturbations
+        # ----------------------------
+        for i,amplitude in enumerate(Perturbations):
+            dX      = amplitude * dX0
+            #
+            if self._parameters["ResiduFormula"] == "Taylor":
+                FX_plus_dX  = numpy.asmatrix(numpy.ravel( Hm( Xn + dX ) )).T
+                #
+                Residu = numpy.linalg.norm( FX_plus_dX - FX ) / (amplitude * NormeGX)
+                #
+                self.StoredVariables["CostFunctionJ"].store( Residu )
+                msg = "  %2i  %5.0e   %9.3e   %9.3e   |   %11.5e    %5.1e"%(i,amplitude,NormeX,NormeFX,Residu,abs(Residu-1.)/amplitude)
+                msgs += "\n" + __marge + msg
+        #
+        msgs += "\n" + __marge + "-"*__nbtirets
+        msgs += "\n"
+        #
+        # Sorties eventuelles
+        # -------------------
+        print
+        print "Results of tangent check by \"%s\" formula:"%self._parameters["ResiduFormula"]
+        print msgs
+        #
+        self._post_run(HO)
+        return 0
+
+# ==============================================================================
+if __name__ == "__main__":
+    print '\n AUTODIAGNOSTIC \n'
index 326da3827578d4f1db07383c25ceca4758db4349..64431a0209cbd30fd03a0b036aa47a47dd025ec3 100644 (file)
@@ -76,6 +76,7 @@ CheckAlgos = [
     "GradientTest",
     "AdjointTest",
     "ObserverTest",
+    "TangentTest",
     ]
 
 AlgoDataRequirements = {}
@@ -152,6 +153,10 @@ AlgoDataRequirements["AdjointTest"] = [
 AlgoDataRequirements["ObserverTest"] = [
     "Observers",
     ]
+AlgoDataRequirements["TangentTest"] = [
+    "CheckingPoint",
+    "ObservationOperator",
+    ]
 
 AlgoType = {}
 AlgoType["3DVAR"] = "Optim"