/src/quantlib/ql/math/randomnumbers/haltonrsg.cpp
Line | Count | Source |
1 | | /* -*- mode: c++; tab-width: 4; indent-tabs-mode: nil; c-basic-offset: 4 -*- */ |
2 | | |
3 | | /* |
4 | | Copyright (C) 2003 Ferdinando Ametrano |
5 | | |
6 | | This file is part of QuantLib, a free-software/open-source library |
7 | | for financial quantitative analysts and developers - http://quantlib.org/ |
8 | | |
9 | | QuantLib is free software: you can redistribute it and/or modify it |
10 | | under the terms of the QuantLib license. You should have received a |
11 | | copy of the license along with this program; if not, please email |
12 | | <quantlib-dev@lists.sf.net>. The license is also available online at |
13 | | <https://www.quantlib.org/license.shtml>. |
14 | | |
15 | | This program is distributed in the hope that it will be useful, but WITHOUT |
16 | | ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS |
17 | | FOR A PARTICULAR PURPOSE. See the license for more details. |
18 | | */ |
19 | | |
20 | | // =========================================================================== |
21 | | // NOTE: The following copyright notice applies to the original code, |
22 | | // |
23 | | // Copyright (C) 2002 Peter Jäckel "Monte Carlo Methods in Finance". |
24 | | // All rights reserved. |
25 | | // |
26 | | // Permission to use, copy, modify, and distribute this software is freely |
27 | | // granted, provided that this notice is preserved. |
28 | | // =========================================================================== |
29 | | |
30 | | #include <ql/math/randomnumbers/haltonrsg.hpp> |
31 | | #include <ql/math/randomnumbers/rngtraits.hpp> |
32 | | #include <ql/math/primenumbers.hpp> |
33 | | |
34 | | namespace QuantLib { |
35 | | |
36 | | HaltonRsg::HaltonRsg(Size dimensionality, |
37 | | unsigned long seed, |
38 | | bool randomStart, |
39 | | bool randomShift) |
40 | 0 | : dimensionality_(dimensionality), sequence_(std::vector<Real>(dimensionality), 1.0), |
41 | 0 | randomStart_(dimensionality, 0UL), randomShift_(dimensionality, 0.0) { |
42 | |
|
43 | 0 | QL_REQUIRE(dimensionality>0, |
44 | 0 | "dimensionality must be greater than 0"); |
45 | | |
46 | 0 | if (randomStart || randomShift) { |
47 | 0 | RandomSequenceGenerator<MersenneTwisterUniformRng> |
48 | 0 | uniformRsg(dimensionality_, seed); |
49 | 0 | if (randomStart) |
50 | 0 | randomStart_ = uniformRsg.nextInt32Sequence(); |
51 | 0 | if (randomShift) |
52 | 0 | randomShift_ = uniformRsg.nextSequence().value; |
53 | 0 | } |
54 | 0 | } |
55 | | |
56 | 0 | const HaltonRsg::sample_type& HaltonRsg::nextSequence() const { |
57 | 0 | ++sequenceCounter_; |
58 | 0 | for (Size i=0; i<dimensionality_; ++i) { |
59 | 0 | double h = 0.0; |
60 | 0 | unsigned long b = PrimeNumbers::get(i); |
61 | 0 | double f = 1.0; |
62 | 0 | unsigned long k = sequenceCounter_+randomStart_[i]; |
63 | 0 | while (k != 0U) { |
64 | 0 | f /= b; |
65 | 0 | h += (k%b)*f; |
66 | 0 | k /= b; |
67 | 0 | } |
68 | 0 | sequence_.value[i] = h+randomShift_[i]; |
69 | 0 | sequence_.value[i] -= long(sequence_.value[i]); |
70 | 0 | } |
71 | 0 | return sequence_; |
72 | 0 | } |
73 | | |
74 | | } |
75 | | |