Line data Source code
1 : //* This file is part of the MOOSE framework
2 : //* https://mooseframework.inl.gov
3 : //*
4 : //* All rights reserved, see COPYRIGHT for full restrictions
5 : //* https://github.com/idaholab/moose/blob/master/COPYRIGHT
6 : //*
7 : //* Licensed under LGPL 2.1, please see LICENSE for details
8 : //* https://www.gnu.org/licenses/lgpl-2.1.html
9 :
10 : #pragma once
11 :
12 : // MOOSE includes
13 : #include "MooseError.h"
14 : #include "DataIO.h"
15 :
16 : #include <unordered_map>
17 :
18 : // External library includes
19 : #include "randistrs.h"
20 :
21 : /**
22 : * This class encapsulates a useful, consistent, cross-platform random number generator
23 : * with multiple utilities.
24 : *
25 : * 1. SIMPLE INTERFACE:
26 : * There are three static functions that are suitable as a drop in replacement for the
27 : * random number capabilities available in the standard C++ library.
28 : * This interface DOES NOT SUPPORT restart/recover
29 : *
30 : * 2. ADVANCED INTERFACE:
31 : * When creating an instance of this class, one can maintain an arbitrary number of
32 : * multiple independent streams of random numbers. Furthermore, the state of these
33 : * generators can be saved and restored for all streams by using the "saveState" and
34 : * "restoreState" methods. Finally, this class uses a fast hash map so that indexes
35 : * for the generators are not required to be contiguous.
36 : */
37 : class MooseRandom
38 : {
39 : public:
40 : /**
41 : * The method seeds the random number generator
42 : * @param seed the seed number
43 : */
44 51840 : static inline void seed(unsigned int seed) { mt_seed32new(seed); }
45 :
46 : /**
47 : * This method returns the next random number (Real format) from the generator
48 : * @return the next random number in the range [0,1) with 64-bit precision
49 : */
50 13404 : static inline Real rand() { return mt_ldrand(); }
51 :
52 : /**
53 : * This method returns the next random number (Real format) from the generator,
54 : * drawn from a normal distribution centered around mean, with a width of sigma
55 : * @param mean center of the random number distribution
56 : * @param sigma width of the random number distribution
57 : * @return the next random number following a normal distribution of width sigma around mean
58 : * with 64-bit precision
59 : */
60 2 : static inline Real randNormal(Real mean, Real sigma) { return rd_normal(mean, sigma); }
61 :
62 : /**
63 : * Return next random number drawn from a standard distribution.
64 : */
65 1 : static inline Real randNormal() { return randNormal(0.0, 1.0); }
66 :
67 : /**
68 : * This method returns the next random number (long format) from the generator
69 : * @return the next random number in the range [0,max(uinit32_t)) with 32-bit number
70 : */
71 993 : static inline uint32_t randl() { return mt_lrand(); }
72 :
73 : /**
74 : * The method seeds one of the independent random number generators
75 : * @param i the index of the generator
76 : * @param seed the seed number
77 : */
78 402027 : inline void seed(std::size_t i, unsigned int seed) { mts_seed32new(&(_states[i].first), seed); }
79 :
80 : /**
81 : * This method returns the next random number (Real format) from the specified generator
82 : * @param i the index of the generator
83 : * @return the next random number in the range [0,1) with 64-bit precision
84 : */
85 6347263 : inline Real rand(std::size_t i)
86 : {
87 : // mooseAssert(_states.find(i) != _states.end(), "No random state initialized for id: " << i);
88 6347263 : return mts_ldrand(&(_states[i].first));
89 : }
90 :
91 : /**
92 : * This method returns the next random number (Real format) from the specified generator,
93 : * drawn from a normal distribution centered around mean, with a width of sigma
94 : * @param i the index of the generator
95 : * @param mean center of the random number distribution
96 : * @param sigma width of the random number distribution
97 : * @return the next random number following a normal distribution of width sigma around mean
98 : * with 64-bit precision
99 : */
100 : inline Real randNormal(std::size_t i, Real mean, Real sigma)
101 : {
102 : mooseAssert(_states.find(i) != _states.end(), "No random state initialized for id: " << i);
103 : return rds_normal(&(_states[i].first), mean, sigma);
104 : }
105 :
106 : /**
107 : * Return next random number drawn from a standard distribution.
108 : */
109 : inline Real randNormal(std::size_t i) { return randNormal(i, 0.0, 1.0); }
110 :
111 : /**
112 : * This method returns the next random number (long format) from the specified generator
113 : * @param i the index of the generator
114 : * @return the next random number in the range [0,max(uinit32_t)) with 32-bit number
115 : */
116 3093768 : inline uint32_t randl(std::size_t i)
117 : {
118 : mooseAssert(_states.find(i) != _states.end(), "No random state initialized for id: " << i);
119 3093768 : return mts_lrand(&(_states[i].first));
120 : }
121 :
122 : /**
123 : * This method returns the next random number (long format) from the specified generator
124 : * within the supplied range.
125 : *
126 : * @param lower lower bounds of value
127 : * @param upper upper bounds of value
128 : * @param i the index of the generator
129 : * @return the next random number in the range [0,max(uinit32_t)) with 32-bit number
130 : */
131 248 : inline uint32_t randl(std::size_t i, uint32_t lower, uint32_t upper)
132 : {
133 : mooseAssert(_states.find(i) != _states.end(), "No random state initialized for id: " << i);
134 248 : return rds_iuniform(&(_states[i].first), lower, upper);
135 : }
136 :
137 : /**
138 : * This method saves the current state of all generators which can be restored at a later time
139 : * (i.e. re-generate the same sequence of random numbers of this generator
140 : */
141 3233 : void saveState()
142 : {
143 3233 : _saved = true;
144 3233 : std::for_each(_states.begin(),
145 : _states.end(),
146 401171 : [](std::pair<const std::size_t, std::pair<mt_state, mt_state>> & pair)
147 401171 : { pair.second.second = pair.second.first; });
148 3233 : }
149 :
150 : /**
151 : * This method restores the last saved generator state
152 : */
153 13258 : void restoreState()
154 : {
155 : mooseAssert(_saved, "saveState() must be called prior to restoreState().");
156 13258 : std::for_each(_states.begin(),
157 : _states.end(),
158 1010437 : [](std::pair<const std::size_t, std::pair<mt_state, mt_state>> & pair)
159 1010437 : { pair.second.first = pair.second.second; });
160 13258 : }
161 :
162 : /**
163 : * Return the number of states.
164 : */
165 862 : inline std::size_t size() { return _states.size(); }
166 :
167 : private:
168 : /**
169 : * We store a pair of states in this map. The first one is the active state, the
170 : * second is the backup state. It is used to restore state at a later time
171 : * to the active state.
172 : */
173 : std::unordered_map<std::size_t, std::pair<mt_state, mt_state>> _states;
174 :
175 : // for restart capability
176 : friend void dataStore<MooseRandom>(std::ostream & stream, MooseRandom & v, void * context);
177 : friend void dataLoad<MooseRandom>(std::istream & stream, MooseRandom & v, void * context);
178 :
179 : /// Flag to make certain that saveState is called prior to restoreState
180 : bool _saved = false;
181 : };
182 :
183 : template <>
184 : inline void
185 55 : dataStore(std::ostream & stream, MooseRandom & v, void * context)
186 : {
187 55 : storeHelper(stream, v._states, context);
188 55 : }
189 : template <>
190 : inline void
191 11 : dataLoad(std::istream & stream, MooseRandom & v, void * context)
192 : {
193 11 : loadHelper(stream, v._states, context);
194 11 : }
|