source: Sophya/trunk/SophyaLib/BaseTools/randinterf.h@ 4024

Last change on this file since 4024 was 4018, checked in by cmv, 14 years ago
  • introduction du tirage gaussien par la methode de la ziggurat
  • ajout argument slim par defaut dans GaussianTail pour proteger d'une reccursion infinie genere par le tirage gaussien ziggurat dans le cas ou GaussianTail applique la methode par rejection (sdev<slim)

cmv, 21/09/2011

File size: 5.5 KB
Line 
1#ifndef RANDOMGENERATORINTERFACE_H_SEEN
2#define RANDOMGENERATORINTERFACE_H_SEEN
3
4// Classes RandomGeneratorInterface
5//
6// R. Ansari (C) UPS+LAL IN2P3/CNRS 2009
7// C. Magneville (C) DAPNIA/SPP CEA 2009
8
9
10#include "machdefs.h"
11#include "anydataobj.h"
12#include <complex>
13#include <vector>
14
15namespace SOPHYA {
16
17//! enum definition for usual types of probability distributio finctions
18enum RNDTypes {
19 C_RND_Flat = 0,
20 C_RND_Gaussian = 1,
21 C_RND_Poisson = 2,
22 C_RND_Exponential = 3
23};
24
25//! enum definition for the type of algorithm used for Gaussian distribution generation
26enum GaussianGenAlgo {
27 C_Gaussian_BoxMuller = 0,
28 C_Gaussian_RandLibSNorm = 1,
29 C_Gaussian_PolarBoxMuller = 2,
30 C_Gaussian_RatioUnif = 3,
31 C_Gaussian_LevaRatioUnif = 4,
32 C_Gaussian_Ziggurat128 = 5
33};
34
35//! enum definition for the type of algorithm used for Poisson distribution generation
36enum PoissonGenAlgo {
37 C_Poisson_Simple = 0,
38 C_Poisson_Ahrens = 1
39};
40
41//! enum definition for the type of algorithm used for exponential distribution generation
42enum ExponentialGenAlgo {
43 C_Exponential_Simple = 0,
44 C_Exponential_Ahrens = 1
45};
46
47// Definition de l'interface pour les classes generateur aleatoire
48class RandomGeneratorInterface : public AnyDataObj {
49public:
50 RandomGeneratorInterface();
51 virtual ~RandomGeneratorInterface();
52
53//! Select Gaussian generation algorithm
54 inline void SelectGaussianAlgo(GaussianGenAlgo typ=C_Gaussian_BoxMuller)
55 { usegaussian_ = typ; }
56//! Select Poisson generation algorithm
57 inline void SelectPoissonAlgo(PoissonGenAlgo typ=C_Poisson_Simple)
58 { usepoisson_ = typ; }
59//! Select Exponential generation algorithm
60 inline void SelectExponentialAlgo(ExponentialGenAlgo typ=C_Exponential_Simple)
61 { useexpo_ = typ; }
62
63//! Return the Gaussian generation algorithm type
64 inline GaussianGenAlgo GetGaussianAlgo()
65 { return usegaussian_; }
66//! Return the Poisson generation algorithm type
67 inline PoissonGenAlgo GetPoissonAlgo()
68 { return usepoisson_; }
69//! Return the Exponential generation algorithm type
70//! Select Exponential generation algorithm
71 inline ExponentialGenAlgo GetExponentialAlgo()
72 { return useexpo_; }
73
74 virtual void GenerateSeedVector(int nseed,vector<uint_2>& seed,int lp=0);
75
76 // --- Le tirage sur une distribution plate
77 /*! \brief Return a random number \b r with a flat distribution 0<=r<1 */
78 inline r_8 Flat01() {return Next();}
79 /*! \brief Return a random number \b r with a flat distribution -1<=r<1 */
80 inline r_8 Flatpm1() {return 2.*Next()-1.;}
81
82 // --- Le tirage sur une distribution gaussienne
83 /*! \brief Return a random number following a gaussian distribution (sigma=1, mean=0)*/
84 virtual r_8 Gaussian();
85
86 virtual r_8 GaussianBoxMuller();
87 virtual r_8 GaussianSNorm();
88 virtual r_8 GaussianPolarBoxMuller();
89 virtual r_8 GaussianRatioUnif();
90 virtual r_8 GaussianLevaRatioUnif();
91 virtual r_8 GaussianZiggurat128();
92 /*! \brief Return a random number following a gaussian distribution "sigma", (mean=0)*/
93 inline r_8 Gaussian(double sigma) {return sigma*Gaussian();}
94 /*! \brief Return a random number following a gaussian distribution with mean=mu, and sigma */
95 inline r_8 Gaussian(double sigma,double mu) {return sigma*Gaussian()+mu;}
96
97 /*! \brief Return a random number following a gaussian tail distribution for x>sdev */
98 virtual r_8 GaussianTail(double sdev,double slim=1.);
99
100 // --- Le tirage sur une distribution de poisson
101 /*! \brief Return a random number following a poisson distribution with mean mu */
102 virtual uint_8 Poisson(double mu, double mumax=-1);
103
104 virtual uint_8 PoissonSimple(double mu, double mumax=-1);
105 virtual uint_8 PoissonAhrens(double mu);
106
107 //! Return a random number with exponential distribution
108 virtual r_8 Exponential();
109
110 virtual r_8 ExpoSimple();
111 virtual r_8 ExpoAhrens();
112
113 // --- Le tirage gaussien complexe (cf texte a la fin du .cc)
114 inline complex< r_8 > ComplexGaussian(void)
115 {return complex< r_8 >(Gaussian(),Gaussian());}
116 inline complex< r_8 > ComplexGaussian(double sig)
117 {return complex< r_8 >(sig*Gaussian(),sig*Gaussian());}
118 /*! \brief Returns the module of a random complex number generated by ComplexGaussRan */
119 inline double ModComplexGaussian(double sig=1.)
120 {double r=-log(1.-Next()); return sig*sqrt(2.*r);}
121
122 // --- Le tirage sur une distribution gaussienne a 2D
123 /*! \brief Return 2 random numbers following a gaussian 2D distribution */
124 virtual int Gaussian2DRho(double &x,double &y,double mx,double my,double sx,double sy,double ro);
125 virtual void Gaussian2DAng(double &x,double &y,double mx,double my,double sa,double sb,double teta);
126
127 //! Return the pointer to the default global RandomGenerator object
128 static inline RandomGeneratorInterface* GetGlobalRandGenP()
129 { return gl_rndgen_p; }
130
131 //! Initialisation automatique
132 virtual void AutoInit(int lp=0);
133
134 //! Print
135 virtual void ShowRandom();
136
137 // Permet de definir l'instance global du generateur aleatoire
138 static void SetGlobalRandGenP(RandomGeneratorInterface* rgp);
139
140protected:
141 /*! \brief Should return a random number in the range [0,1] with flat distribution
142
143 This pure virtual method should be implemented by inheriting classes */
144 virtual r_8 Next() = 0;
145
146// Selection de type d'algo pour les aleatoires autres que plates
147 GaussianGenAlgo usegaussian_;
148 PoissonGenAlgo usepoisson_;
149 ExponentialGenAlgo useexpo_;
150
151 // Instance global d'un generateur par defaut
152 static RandomGeneratorInterface* gl_rndgen_p;
153};
154
155} /* namespace SOPHYA */
156
157#endif
Note: See TracBrowser for help on using the repository browser.