src/daComposant/daCore/BasicObjects.py

   1 #-*-coding:iso-8859-1-*-
   2 #
   3 #  Copyright (C) 2008-2014 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
  23 __doc__ = """
  24     Définit les outils généraux élémentaires.
  25
  26     Ce module est destiné à etre appelée par AssimilationStudy pour constituer
  27     les objets élémentaires de l'algorithme.
  28 """
  29 __author__ = "Jean-Philippe ARGAUD"
  30
  31 import logging
  32 import numpy
  33 import Persistence
  34
  35 # ==============================================================================
  36 class Operator:
  37     """
  38     Classe générale d'interface de type opérateur
  39     """
  40     NbCallsAsMatrix = 0
  41     NbCallsAsMethod = 0
  42     #
  43     def __init__(self, fromMethod=None, fromMatrix=None):
  44         """
  45         On construit un objet de ce type en fournissant à l'aide de l'un des
  46         deux mots-clé, soit une fonction python, soit matrice.
  47         Arguments :
  48         - fromMethod : argument de type fonction Python
  49         - fromMatrix : argument adapté au constructeur numpy.matrix
  50         """
  51         self.__NbCallsAsMatrix, self.__NbCallsAsMethod = 0, 0
  52         if   fromMethod is not None:
  53             self.__Method = fromMethod
  54             self.__Matrix = None
  55             self.__Type   = "Method"
  56         elif fromMatrix is not None:
  57             self.__Method = None
  58             self.__Matrix = numpy.matrix( fromMatrix, numpy.float )
  59             self.__Type   = "Matrix"
  60         else:
  61             self.__Method = None
  62             self.__Matrix = None
  63             self.__Type   = None
  64
  65     def isType(self):
  66         return self.__Type
  67
  68     def appliedTo(self, xValue):
  69         """
  70         Permet de restituer le résultat de l'application de l'opérateur à un
  71         argument xValue. Cette méthode se contente d'appliquer, son argument
  72         devant a priori être du bon type.
  73         Arguments :
  74         - xValue : argument adapté pour appliquer l'opérateur
  75         """
  76         if self.__Matrix is not None:
  77             self.__addOneMatrixCall()
  78             return self.__Matrix * xValue
  79         else:
  80             self.__addOneMethodCall()
  81             return self.__Method( xValue )
  82
  83     def appliedControledFormTo(self, (xValue, uValue) ):
  84         """
  85         Permet de restituer le résultat de l'application de l'opérateur à une
  86         paire (xValue, uValue). Cette méthode se contente d'appliquer, son
  87         argument devant a priori être du bon type. Si la uValue est None,
  88         on suppose que l'opérateur ne s'applique qu'à xValue.
  89         Arguments :
  90         - xValue : argument X adapté pour appliquer l'opérateur
  91         - uValue : argument U adapté pour appliquer l'opérateur
  92         """
  93         if self.__Matrix is not None:
  94             self.__addOneMatrixCall()
  95             return self.__Matrix * xValue
  96         elif uValue is not None:
  97             self.__addOneMethodCall()
  98             return self.__Method( (xValue, uValue) )
  99         else:
 100             self.__addOneMethodCall()
 101             return self.__Method( xValue )
 102
 103     def appliedInXTo(self, (xNominal, xValue) ):
 104         """
 105         Permet de restituer le résultat de l'application de l'opérateur à un
 106         argument xValue, sachant que l'opérateur est valable en xNominal.
 107         Cette méthode se contente d'appliquer, son argument devant a priori
 108         être du bon type. Si l'opérateur est linéaire car c'est une matrice,
 109         alors il est valable en tout point nominal et il n'est pas nécessaire
 110         d'utiliser xNominal.
 111         Arguments : une liste contenant
 112         - xNominal : argument permettant de donner le point où l'opérateur
 113           est construit pour etre ensuite appliqué
 114         - xValue : argument adapté pour appliquer l'opérateur
 115         """
 116         if self.__Matrix is not None:
 117             self.__addOneMatrixCall()
 118             return self.__Matrix * xValue
 119         else:
 120             self.__addOneMethodCall()
 121             return self.__Method( (xNominal, xValue) )
 122
 123     def asMatrix(self, ValueForMethodForm = "UnknownVoidValue"):
 124         """
 125         Permet de renvoyer l'opérateur sous la forme d'une matrice
 126         """
 127         if self.__Matrix is not None:
 128             self.__addOneMatrixCall()
 129             return self.__Matrix
 130         elif ValueForMethodForm is not "UnknownVoidValue": # Ne pas utiliser "None"
 131             self.__addOneMethodCall()
 132             return numpy.matrix( self.__Method( (ValueForMethodForm, None) ) )
 133         else:
 134             raise ValueError("Matrix form of the operator defined as a function/method requires to give an operating point.")
 135
 136     def shape(self):
 137         """
 138         Renvoie la taille sous forme numpy si l'opérateur est disponible sous
 139         la forme d'une matrice
 140         """
 141         if self.__Matrix is not None:
 142             return self.__Matrix.shape
 143         else:
 144             raise ValueError("Matrix form of the operator is not available, nor the shape")
 145
 146     def nbcalls(self, which=None):
 147         """
 148         Renvoie les nombres d'évaluations de l'opérateur
 149         """
 150         __nbcalls = (
 151             self.__NbCallsAsMatrix+self.__NbCallsAsMethod,
 152             self.__NbCallsAsMatrix,
 153             self.__NbCallsAsMethod,
 154             Operator.NbCallsAsMatrix+Operator.NbCallsAsMethod,
 155             Operator.NbCallsAsMatrix,
 156             Operator.NbCallsAsMethod,
 157             )
 158         if which is None: return __nbcalls
 159         else:             return __nbcalls[which]
 160
 161     def __addOneMatrixCall(self):
 162         self.__NbCallsAsMatrix   += 1 # Decompte local
 163         Operator.NbCallsAsMatrix += 1 # Decompte global
 164
 165     def __addOneMethodCall(self):
 166         self.__NbCallsAsMethod   += 1 # Decompte local
 167         Operator.NbCallsAsMethod += 1 # Decompte global
 168
 169 # ==============================================================================
 170 class Algorithm:
 171     """
 172     Classe générale d'interface de type algorithme
 173
 174     Elle donne un cadre pour l'écriture d'une classe élémentaire d'algorithme
 175     d'assimilation, en fournissant un container (dictionnaire) de variables
 176     persistantes initialisées, et des méthodes d'accès à ces variables stockées.
 177
 178     Une classe élémentaire d'algorithme doit implémenter la méthode "run".
 179     """
 180     def __init__(self, name):
 181         """
 182         L'initialisation présente permet de fabriquer des variables de stockage
 183         disponibles de manière générique dans les algorithmes élémentaires. Ces
 184         variables de stockage sont ensuite conservées dans un dictionnaire
 185         interne à l'objet, mais auquel on accède par la méthode "get".
 186
 187         Les variables prévues sont :
 188             - CostFunctionJ  : fonction-cout globale, somme des deux parties suivantes
 189             - CostFunctionJb : partie ébauche ou background de la fonction-cout
 190             - CostFunctionJo : partie observations de la fonction-cout
 191             - GradientOfCostFunctionJ  : gradient de la fonction-cout globale
 192             - GradientOfCostFunctionJb : gradient de la partie ébauche de la fonction-cout
 193             - GradientOfCostFunctionJo : gradient de la partie observations de la fonction-cout
 194             - CurrentState : état courant lors d'itérations
 195             - Analysis : l'analyse
 196             - Innovation : l'innovation : d = Y - H Xb
 197             - SigmaObs2 : indicateur de correction optimale des erreurs d'observation
 198             - SigmaBck2 : indicateur de correction optimale des erreurs d'ébauche
 199             - MahalanobisConsistency : indicateur de consistance des covariances
 200             - OMA : Observation moins Analysis : Y - Xa
 201             - OMB : Observation moins Background : Y - Xb
 202             - AMB : Analysis moins Background : Xa - Xb
 203             - APosterioriCovariance : matrice A
 204         On peut rajouter des variables à stocker dans l'initialisation de
 205         l'algorithme élémentaire qui va hériter de cette classe
 206         """
 207         logging.debug("%s Initialisation"%str(name))
 208         self._name = str( name )
 209         self._parameters = {}
 210         self.__required_parameters = {}
 211         self.StoredVariables = {}
 212         #
 213         self.StoredVariables["CostFunctionJ"]            = Persistence.OneScalar(name = "CostFunctionJ")
 214         self.StoredVariables["CostFunctionJb"]           = Persistence.OneScalar(name = "CostFunctionJb")
 215         self.StoredVariables["CostFunctionJo"]           = Persistence.OneScalar(name = "CostFunctionJo")
 216         self.StoredVariables["GradientOfCostFunctionJ"]  = Persistence.OneVector(name = "GradientOfCostFunctionJ")
 217         self.StoredVariables["GradientOfCostFunctionJb"] = Persistence.OneVector(name = "GradientOfCostFunctionJb")
 218         self.StoredVariables["GradientOfCostFunctionJo"] = Persistence.OneVector(name = "GradientOfCostFunctionJo")
 219         self.StoredVariables["CurrentState"]             = Persistence.OneVector(name = "CurrentState")
 220         self.StoredVariables["Analysis"]                 = Persistence.OneVector(name = "Analysis")
 221         self.StoredVariables["Innovation"]               = Persistence.OneVector(name = "Innovation")
 222         self.StoredVariables["SigmaObs2"]                = Persistence.OneScalar(name = "SigmaObs2")
 223         self.StoredVariables["SigmaBck2"]                = Persistence.OneScalar(name = "SigmaBck2")
 224         self.StoredVariables["MahalanobisConsistency"]   = Persistence.OneScalar(name = "MahalanobisConsistency")
 225         self.StoredVariables["OMA"]                      = Persistence.OneVector(name = "OMA")
 226         self.StoredVariables["OMB"]                      = Persistence.OneVector(name = "OMB")
 227         self.StoredVariables["BMA"]                      = Persistence.OneVector(name = "BMA")
 228         self.StoredVariables["APosterioriCovariance"]    = Persistence.OneMatrix(name = "APosterioriCovariance")
 229         self.StoredVariables["SimulationQuantiles"]      = Persistence.OneMatrix(name = "SimulationQuantiles")
 230
 231     def get(self, key=None):
 232         """
 233         Renvoie l'une des variables stockées identifiée par la clé, ou le
 234         dictionnaire de l'ensemble des variables disponibles en l'absence de
 235         clé. Ce sont directement les variables sous forme objet qui sont
 236         renvoyées, donc les méthodes d'accès à l'objet individuel sont celles
 237         des classes de persistance.
 238         """
 239         if key is not None:
 240             return self.StoredVariables[key]
 241         else:
 242             return self.StoredVariables
 243
 244     def has_key(self, key=None):
 245         """
 246         Vérifie si l'une des variables stockées est identifiée par la clé.
 247         """
 248         return self.StoredVariables.has_key(key)
 249
 250     def keys(self):
 251         """
 252         Renvoie la liste des clés de variables stockées.
 253         """
 254         return self.StoredVariables.keys()
 255
 256     def run(self, Xb=None, Y=None, H=None, M=None, R=None, B=None, Q=None, Parameters=None):
 257         """
 258         Doit implémenter l'opération élémentaire de calcul d'assimilation sous
 259         sa forme mathématique la plus naturelle possible.
 260         """
 261         raise NotImplementedError("Mathematical assimilation calculation has not been implemented!")
 262
 263     def defineRequiredParameter(self, name = None, default = None, typecast = None, message = None, minval = None, maxval = None, listval = None):
 264         """
 265         Permet de définir dans l'algorithme des paramètres requis et leurs
 266         caractéristiques par défaut.
 267         """
 268         if name is None:
 269             raise ValueError("A name is mandatory to define a required parameter.")
 270         #
 271         self.__required_parameters[name] = {
 272             "default"  : default,
 273             "typecast" : typecast,
 274             "minval"   : minval,
 275             "maxval"   : maxval,
 276             "listval"  : listval,
 277             "message"  : message,
 278             }
 279         logging.debug("%s %s (valeur par défaut = %s)"%(self._name, message, self.setParameterValue(name)))
 280
 281     def getRequiredParameters(self, noDetails=True):
 282         """
 283         Renvoie la liste des noms de paramètres requis ou directement le
 284         dictionnaire des paramètres requis.
 285         """
 286         if noDetails:
 287             ks = self.__required_parameters.keys()
 288             ks.sort()
 289             return ks
 290         else:
 291             return self.__required_parameters
 292
 293     def setParameterValue(self, name=None, value=None):
 294         """
 295         Renvoie la valeur d'un paramètre requis de manière contrôlée
 296         """
 297         default  = self.__required_parameters[name]["default"]
 298         typecast = self.__required_parameters[name]["typecast"]
 299         minval   = self.__required_parameters[name]["minval"]
 300         maxval   = self.__required_parameters[name]["maxval"]
 301         listval  = self.__required_parameters[name]["listval"]
 302         #
 303         if value is None and default is None:
 304             __val = None
 305         elif value is None and default is not None:
 306             if typecast is None: __val = default
 307             else:                __val = typecast( default )
 308         else:
 309             if typecast is None: __val = value
 310             else:                __val = typecast( value )
 311         #
 312         if minval is not None and __val < minval:
 313             raise ValueError("The parameter named \"%s\" of value \"%s\" can not be less than %s."%(name, __val, minval))
 314         if maxval is not None and __val > maxval:
 315             raise ValueError("The parameter named \"%s\" of value \"%s\" can not be greater than %s."%(name, __val, maxval))
 316         if listval is not None:
 317             if typecast is list or typecast is tuple or type(__val) is list or type(__val) is tuple:
 318                 for v in __val:
 319                     if v not in listval:
 320                         raise ValueError("The value \"%s\" of the parameter named \"%s\" is not allowed, it has to be in the list %s."%(v, name, listval))
 321             elif __val not in listval:
 322                 raise ValueError("The value \"%s\" of the parameter named \"%s\" is not allowed, it has to be in the list %s."%( __val, name,listval))
 323         return __val
 324
 325     def setParameters(self, fromDico={}):
 326         """
 327         Permet de stocker les paramètres reçus dans le dictionnaire interne.
 328         """
 329         self._parameters.update( fromDico )
 330         for k in self.__required_parameters.keys():
 331             if k in fromDico.keys():
 332                 self._parameters[k] = self.setParameterValue(k,fromDico[k])
 333             else:
 334                 self._parameters[k] = self.setParameterValue(k)
 335             logging.debug("%s %s : %s"%(self._name, self.__required_parameters[k]["message"], self._parameters[k]))
 336
 337 # ==============================================================================
 338 class Diagnostic:
 339     """
 340     Classe générale d'interface de type diagnostic
 341
 342     Ce template s'utilise de la manière suivante : il sert de classe "patron" en
 343     même temps que l'une des classes de persistance, comme "OneScalar" par
 344     exemple.
 345
 346     Une classe élémentaire de diagnostic doit implémenter ses deux méthodes, la
 347     méthode "_formula" pour écrire explicitement et proprement la formule pour
 348     l'écriture mathématique du calcul du diagnostic (méthode interne non
 349     publique), et "calculate" pour activer la précédente tout en ayant vérifié
 350     et préparé les données, et pour stocker les résultats à chaque pas (méthode
 351     externe d'activation).
 352     """
 353     def __init__(self, name = "", parameters = {}):
 354         self.name       = str(name)
 355         self.parameters = dict( parameters )
 356
 357     def _formula(self, *args):
 358         """
 359         Doit implémenter l'opération élémentaire de diagnostic sous sa forme
 360         mathématique la plus naturelle possible.
 361         """
 362         raise NotImplementedError("Diagnostic mathematical formula has not been implemented!")
 363
 364     def calculate(self, *args):
 365         """
 366         Active la formule de calcul avec les arguments correctement rangés
 367         """
 368         raise NotImplementedError("Diagnostic activation method has not been implemented!")
 369
 370 # ==============================================================================
 371 class Covariance:
 372     """
 373     Classe générale d'interface de type covariance
 374     """
 375     def __init__(self,
 376             name          = "GenericCovariance",
 377             asCovariance  = None,
 378             asEyeByScalar = None,
 379             asEyeByVector = None,
 380             ):
 381         """
 382         Permet de définir une covariance :
 383         - asCovariance : entrée des données, comme une matrice compatible avec
 384           le constructeur de numpy.matrix
 385         - asEyeByScalar : entrée des données comme un seul scalaire de variance,
 386           multiplicatif d'une matrice de corrélation identité, aucune matrice
 387           n'étant donc explicitement à donner
 388         - asEyeByVector : entrée des données comme un seul vecteur de variance,
 389           à mettre sur la diagonale d'une matrice de corrélation, aucune matrice
 390           n'étant donc explicitement à donner
 391         """
 392         self.__name       = str(name)
 393         #
 394         self.__B          = None
 395         self.__is_scalar  = False
 396         self.__is_vector  = False
 397         self.__is_matrix  = False
 398         if asEyeByScalar is not None:
 399             self.__is_scalar = True
 400             self.__B         = numpy.abs( float(asEyeByScalar) )
 401             self.shape       = (0,0)
 402             self.size        = 0
 403         elif asEyeByVector is not None:
 404             self.__is_vector = True
 405             self.__B         = numpy.abs( numpy.array( numpy.ravel( asEyeByVector ), float ) )
 406             self.shape       = (self.__B.size,self.__B.size)
 407             self.size        = self.__B.size**2
 408         elif asCovariance is not None:
 409             self.__is_matrix = True
 410             self.__B         = numpy.matrix( asCovariance, float )
 411             self.shape       = self.__B.shape
 412             self.size        = self.__B.size
 413         else:
 414             pass
 415             # raise ValueError("The %s covariance matrix has to be specified either as a matrix, a vector for its diagonal or a scalar multiplying an identity matrix."%self.__name)
 416         #
 417         self.__validate()
 418
 419     def __validate(self):
 420         if self.ismatrix() and min(self.shape) != max(self.shape):
 421             raise ValueError("The given matrix for %s is not a square one, its shape is %s. Please check your matrix input."%(self.__name,self.shape))
 422         if self.isscalar() and self.__B <= 0:
 423             raise ValueError("The %s covariance matrix is not positive-definite. Please check your scalar input %s."%(self.__name,self.__B_scalar))
 424         if self.isvector() and (self.__B <= 0).any():
 425             raise ValueError("The %s covariance matrix is not positive-definite. Please check your vector input."%(self.__name,))
 426         if self.ismatrix() and logging.getLogger().level < logging.WARNING: # La verification n'a lieu qu'en debug
 427             try:
 428                 L = numpy.linalg.cholesky( self.__B )
 429             except:
 430                 raise ValueError("The %s covariance matrix is not symmetric positive-definite. Please check your matrix input."%(self.__name,))
 431
 432     def isscalar(self):
 433         return self.__is_scalar
 434
 435     def isvector(self):
 436         return self.__is_vector
 437
 438     def ismatrix(self):
 439         return self.__is_matrix
 440
 441     def getI(self):
 442         if   self.ismatrix():
 443             return Covariance(self.__name+"I", asCovariance  = self.__B.I )
 444         elif self.isvector():
 445             return Covariance(self.__name+"I", asEyeByVector = 1. / self.__B )
 446         elif self.isscalar():
 447             return Covariance(self.__name+"I", asEyeByScalar = 1. / self.__B )
 448         else:
 449             return None
 450
 451     def getT(self):
 452         if   self.ismatrix():
 453             return Covariance(self.__name+"T", asCovariance  = self.__B.T )
 454         elif self.isvector():
 455             return Covariance(self.__name+"T", asEyeByVector = self.__B )
 456         elif self.isscalar():
 457             return Covariance(self.__name+"T", asEyeByScalar = self.__B )
 458
 459     def cholesky(self):
 460         if   self.ismatrix():
 461             return Covariance(self.__name+"C", asCovariance  = numpy.linalg.cholesky(self.__B) )
 462         elif self.isvector():
 463             return Covariance(self.__name+"C", asEyeByVector = numpy.sqrt( self.__B ) )
 464         elif self.isscalar():
 465             return Covariance(self.__name+"C", asEyeByScalar = numpy.sqrt( self.__B ) )
 466
 467     def choleskyI(self):
 468         if   self.ismatrix():
 469             return Covariance(self.__name+"H", asCovariance  = numpy.linalg.cholesky(self.__B).I )
 470         elif self.isvector():
 471             return Covariance(self.__name+"H", asEyeByVector = 1.0 / numpy.sqrt( self.__B ) )
 472         elif self.isscalar():
 473             return Covariance(self.__name+"H", asEyeByScalar = 1.0 / numpy.sqrt( self.__B ) )
 474
 475     def diag(self, msize=None):
 476         if   self.ismatrix():
 477             return numpy.diag(self.__B)
 478         elif self.isvector():
 479             return self.__B
 480         elif self.isscalar():
 481             if msize is None:
 482                 raise ValueError("the size of the %s covariance matrix has to be given in case of definition as a scalar over the diagonal."%(self.__name,))
 483             else:
 484                 return self.__B * numpy.ones(int(msize))
 485
 486     def asfullmatrix(self, msize=None):
 487         if   self.ismatrix():
 488             return self.__B
 489         elif self.isvector():
 490             return numpy.matrix( numpy.diag(self.__B), float )
 491         elif self.isscalar():
 492             if msize is None:
 493                 raise ValueError("the size of the %s covariance matrix has to be given in case of definition as a scalar over the diagonal."%(self.__name,))
 494             else:
 495                 return numpy.matrix( self.__B * numpy.eye(int(msize)), float )
 496
 497     def trace(self, msize=None):
 498         if   self.ismatrix():
 499             return numpy.trace(self.__B)
 500         elif self.isvector():
 501             return float(numpy.sum(self.__B))
 502         elif self.isscalar():
 503             if msize is None:
 504                 raise ValueError("the size of the %s covariance matrix has to be given in case of definition as a scalar over the diagonal."%(self.__name,))
 505             else:
 506                 return self.__B * int(msize)
 507
 508     def __repr__(self):
 509         return repr(self.__B)
 510
 511     def __str__(self):
 512         return str(self.__B)
 513
 514     def __add__(self, other):
 515         if   self.ismatrix():
 516             return self.__B + numpy.asmatrix(other)
 517         elif self.isvector() or self.isscalar():
 518             _A = numpy.asarray(other)
 519             _A.reshape(_A.size)[::_A.shape[1]+1] += self.__B
 520             return numpy.asmatrix(_A)
 521
 522     def __radd__(self, other):
 523         raise NotImplementedError("%s covariance matrix __radd__ method not available for %s type!"%(self.__name,type(other)))
 524
 525     def __sub__(self, other):
 526         if   self.ismatrix():
 527             return self.__B - numpy.asmatrix(other)
 528         elif self.isvector() or self.isscalar():
 529             _A = numpy.asarray(other)
 530             _A.reshape(_A.size)[::_A.shape[1]+1] = self.__B - _A.reshape(_A.size)[::_A.shape[1]+1]
 531             return numpy.asmatrix(_A)
 532
 533     def __rsub__(self, other):
 534         raise NotImplementedError("%s covariance matrix __rsub__ method not available for %s type!"%(self.__name,type(other)))
 535
 536     def __neg__(self):
 537         return - self.__B
 538
 539     def __mul__(self, other):
 540         if   self.ismatrix() and isinstance(other,numpy.matrix):
 541             return self.__B * other
 542         elif self.ismatrix() and (isinstance(other,numpy.ndarray) \
 543                                or isinstance(other,list) \
 544                                or isinstance(other,tuple)):
 545             if numpy.ravel(other).size == self.shape[1]: # Vecteur
 546                 return self.__B * numpy.asmatrix(numpy.ravel(other)).T
 547             elif numpy.asmatrix(other).shape[0] == self.shape[1]: # Matrice
 548                 return self.__B * numpy.asmatrix(other)
 549             else:
 550                 raise ValueError("operands could not be broadcast together with shapes %s %s in %s matrix"%(self.shape,numpy.asmatrix(other).shape,self.__name))
 551         elif self.isvector() and (isinstance(other,numpy.matrix) \
 552                                or isinstance(other,numpy.ndarray) \
 553                                or isinstance(other,list) \
 554                                or isinstance(other,tuple)):
 555             if numpy.ravel(other).size == self.shape[1]: # Vecteur
 556                 return numpy.asmatrix(self.__B * numpy.ravel(other)).T
 557             elif numpy.asmatrix(other).shape[0] == self.shape[1]: # Matrice
 558                 return numpy.asmatrix((self.__B * (numpy.asarray(other).transpose())).transpose())
 559             else:
 560                 raise ValueError("operands could not be broadcast together with shapes %s %s in %s matrix"%(self.shape,numpy.ravel(other).shape,self.__name))
 561         elif self.isscalar() and isinstance(other,numpy.matrix):
 562             return self.__B * other
 563         elif self.isscalar() and (isinstance(other,numpy.ndarray) \
 564                                or isinstance(other,list) \
 565                                or isinstance(other,tuple)):
 566             if len(numpy.asarray(other).shape) == 1 or numpy.asarray(other).shape[1] == 1 or numpy.asarray(other).shape[0] == 1:
 567                 return self.__B * numpy.asmatrix(numpy.ravel(other)).T
 568             else:
 569                 return self.__B * numpy.asmatrix(other)
 570         else:
 571             raise NotImplementedError("%s covariance matrix __mul__ method not available for %s type!"%(self.__name,type(other)))
 572
 573     def __rmul__(self, other):
 574         if self.ismatrix() and isinstance(other,numpy.matrix):
 575             return other * self.__B
 576         elif self.isvector() and isinstance(other,numpy.matrix):
 577             if numpy.ravel(other).size == self.shape[0]: # Vecteur
 578                 return numpy.asmatrix(numpy.ravel(other) * self.__B)
 579             elif numpy.asmatrix(other).shape[1] == self.shape[0]: # Matrice
 580                 return numpy.asmatrix(numpy.array(other) * self.__B)
 581             else:
 582                 raise ValueError("operands could not be broadcast together with shapes %s %s in %s matrix"%(self.shape,numpy.ravel(other).shape,self.__name))
 583         elif self.isscalar() and isinstance(other,numpy.matrix):
 584             return other * self.__B
 585         else:
 586             raise NotImplementedError("%s covariance matrix __rmul__ method not available for %s type!"%(self.__name,type(other)))
 587
 588     def __len__(self):
 589         return self.shape[0]
 590
 591 # ==============================================================================
 592 if __name__ == "__main__":
 593     print '\n AUTODIAGNOSTIC \n'