1 //                                               -*- C++ -*-
2 /**
3  *  @brief Implementation of the Sobol sequence
4  *
5  *  Copyright 2005-2021 Airbus-EDF-IMACS-ONERA-Phimeca
6  *
7  *  This library is free software: you can redistribute it and/or modify
8  *  it under the terms of the GNU Lesser General Public License as published by
9  *  the Free Software Foundation, either version 3 of the License, or
10  *  (at your option) any later version.
11  *
12  *  This library is distributed in the hope that it will be useful,
13  *  but WITHOUT ANY WARRANTY; without even the implied warranty of
14  *  MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
15  *  GNU Lesser General Public License for more details.
16  *
17  *  You should have received a copy of the GNU Lesser General Public License
18  *  along with this library.  If not, see <http://www.gnu.org/licenses/>.
19  *
20  */
21 
22 #ifndef OPENTURNS_SOBOLSEQUENCE_HXX
23 #define OPENTURNS_SOBOLSEQUENCE_HXX
24 
25 #include "openturns/LowDiscrepancySequenceImplementation.hxx"
26 
27 BEGIN_NAMESPACE_OPENTURNS
28 /**
29  * @class SobolSequence
30  */
31 
32 class OT_API SobolSequence :
33   public LowDiscrepancySequenceImplementation
34 {
35   CLASSNAME
36 
37 public:
38   // this implementation supports dimensions up to 1111
39   static const UnsignedInteger MaximumDimension;
40 
41   // this implementation has a cycle of 2^62 = ~5e18, thanks to the use of 64 bits integers (Unsigned64BitsInteger)
42   static const UnsignedInteger MaximumBase2Logarithm;
43 
44   // this value is 2^-MaximumBase2Logarithm, precomputed to speed up generation
45   static const Scalar Epsilon;
46 
47   // maximum number of columns in InitialBase array = 8
48   static const UnsignedInteger MaximumInitialDegree;
49 
50   // numbers used to generate the coefficients of directionNumber_[][] each row corresponds to a component (dimension)
51   static const UnsignedInteger InitialBase[];
52 
53   // a primitive polynomial used to generate the sequence
54   static const Unsigned64BitsInteger PrimitivePolynomial[];
55 
56 public:
57 
58   /** Constructor with parameters */
59   explicit SobolSequence(const UnsignedInteger dimension = 1);
60 
61   /** Virtual constructor */
62   SobolSequence * clone() const override;
63 
64   /** Initialize the sequence */
65   void initialize(const UnsignedInteger dimension) override;
66 
67   /** Generate a quasi-random vector of numbers uniformly distributed over [0, 1[ */
68   using LowDiscrepancySequenceImplementation::generate;
69   Point generate() const override;
70 
71   /** String converter */
72   String __repr__() const override;
73 
74   /** Method save() stores the object through the StorageManager */
75   void save(Advocate & adv) const override;
76 
77   /** Method load() reloads the object from the StorageManager */
78   void load(Advocate & adv) override;
79 
80 protected:
81   /** The numbers used to generate the sequence */
82   Unsigned64BitsIntegerPersistentCollection base_;
83   mutable Unsigned64BitsIntegerPersistentCollection coefficients_;
84 
85   /** Current seed */
86   mutable Unsigned64BitsInteger seed_;
87 
88 private:
89   /** return 2^n */
90   static Unsigned64BitsInteger inline power2(const UnsignedInteger n);
91 
92   /** Returns the position of the lowest '0' in the binary representation of an unsigned integer */
93   static UnsignedInteger computePositionOfLowest0Bit(const Unsigned64BitsInteger number);
94 
95 }; /* class SobolSequence */
96 
97 END_NAMESPACE_OPENTURNS
98 
99 #endif /* OPENTURNS_SOBOLSEQUENCE_HXX */
100