All Classes Namespaces Files Functions Variables Typedefs Enumerations Enumerator Friends Macros Modules Pages
kaldi-math.cc
Go to the documentation of this file.
1 // base/kaldi-math.cc
2 
3 // Copyright 2009-2011 Microsoft Corporation; Yanmin Qian;
4 // Saarland University; Jan Silovsky
5 
6 // See ../../COPYING for clarification regarding multiple authors
7 //
8 // Licensed under the Apache License, Version 2.0 (the "License");
9 // you may not use this file except in compliance with the License.
10 // You may obtain a copy of the License at
11 //
12 // http://www.apache.org/licenses/LICENSE-2.0
13 //
14 // THIS CODE IS PROVIDED *AS IS* BASIS, WITHOUT WARRANTIES OR CONDITIONS OF ANY
15 // KIND, EITHER EXPRESS OR IMPLIED, INCLUDING WITHOUT LIMITATION ANY IMPLIED
16 // WARRANTIES OR CONDITIONS OF TITLE, FITNESS FOR A PARTICULAR PURPOSE,
17 // MERCHANTABLITY OR NON-INFRINGEMENT.
18 // See the Apache 2 License for the specific language governing permissions and
19 // limitations under the License.
20 
21 #include "base/kaldi-math.h"
22 #ifndef _MSC_VER
23 #include <pthread.h>
24 #include <stdlib.h>
25 #endif
26 #include <string>
27 
28 namespace kaldi {
29 // These routines are tested in matrix/matrix-test.cc
30 
32  KALDI_ASSERT(n > 0);
33  n--;
34  n |= n >> 1;
35  n |= n >> 2;
36  n |= n >> 4;
37  n |= n >> 8;
38  n |= n >> 16;
39  return n+1;
40 }
41 
42 #ifndef _MSC_VER
43 static pthread_mutex_t _RandMutex = PTHREAD_MUTEX_INITIALIZER;
44 #endif
45 
46 int Rand(struct RandomState* state) {
47 #ifdef _MSC_VER
48  // On Windows, just call Rand()
49  return rand();
50 #else
51  if (state) {
52  return rand_r(&(state->seed));
53  } else {
54  int rs = pthread_mutex_lock(&_RandMutex);
55  KALDI_ASSERT(rs == 0);
56  int val = rand();
57  rs = pthread_mutex_unlock(&_RandMutex);
58  KALDI_ASSERT(rs == 0);
59  return val;
60  }
61 #endif
62 }
63 
65  // we initialize it as Rand() + 27437 instead of just Rand(), because on some
66  // systems, e.g. at the very least Mac OSX Yosemite and later, it seems to be
67  // the case that rand_r when initialized with rand() will give you the exact
68  // same sequence of numbers that rand() will give if you keep calling rand()
69  // after that initial call. This can cause problems with repeated sequences.
70  // For example if you initialize two RandomState structs one after the other
71  // without calling rand() in between, they would give you the same sequence
72  // offset by one (if we didn't have the "+ 27437" in the code). 27437 is just
73  // a randomly chosen prime number.
74  seed = Rand() + 27437;
75 }
76 
77 bool WithProb(BaseFloat prob, struct RandomState* state) {
78  KALDI_ASSERT(prob >= 0 && prob <= 1.1); // prob should be <= 1.0,
79  // but we allow slightly larger values that could arise from roundoff in
80  // previous calculations.
81  KALDI_COMPILE_TIME_ASSERT(RAND_MAX > 128 * 128);
82  if (prob == 0) return false;
83  else if (prob == 1.0) return true;
84  else if (prob * RAND_MAX < 128.0) {
85  // prob is very small but nonzero, and the "main algorithm"
86  // wouldn't work that well. So: with probability 1/128, we
87  // return WithProb (prob * 128), else return false.
88  if (Rand(state) < RAND_MAX / 128) { // with probability 128...
89  // Note: we know that prob * 128.0 < 1.0, because
90  // we asserted RAND_MAX > 128 * 128.
91  return WithProb(prob * 128.0);
92  } else {
93  return false;
94  }
95  } else {
96  return (Rand(state) < ((RAND_MAX + static_cast<BaseFloat>(1.0)) * prob));
97  }
98 }
99 
100 int32 RandInt(int32 min_val, int32 max_val, struct RandomState* state) {
101  // This is not exact.
102  KALDI_ASSERT(max_val >= min_val);
103  if (max_val == min_val) return min_val;
104 
105 #ifdef _MSC_VER
106  // RAND_MAX is quite small on Windows -> may need to handle larger numbers.
107  if (RAND_MAX > (max_val-min_val)*8) {
108  // *8 to avoid large inaccuracies in probability, from the modulus...
109  return min_val +
110  ((unsigned int)Rand(state) % (unsigned int)(max_val+1-min_val));
111  } else {
112  if ((unsigned int)(RAND_MAX*RAND_MAX) >
113  (unsigned int)((max_val+1-min_val)*8)) {
114  // *8 to avoid inaccuracies in probability, from the modulus...
115  return min_val + ( (unsigned int)( (Rand(state)+RAND_MAX*Rand(state)))
116  % (unsigned int)(max_val+1-min_val));
117  } else {
118  throw std::runtime_error(std::string()
119  +"rand_int failed because we do not support "
120  +"such large random numbers. "
121  +"(Extend this function).");
122  }
123  }
124 #else
125  return min_val +
126  (static_cast<int32>(Rand(state)) % static_cast<int32>(max_val+1-min_val));
127 #endif
128 }
129 
130 // Returns poisson-distributed random number.
131 // Take care: this takes time proportinal
132 // to lambda. Faster algorithms exist but are more complex.
133 int32 RandPoisson(float lambda, struct RandomState* state) {
134  // Knuth's algorithm.
135  KALDI_ASSERT(lambda >= 0);
136  float L = expf(-lambda), p = 1.0;
137  int32 k = 0;
138  do {
139  k++;
140  float u = RandUniform(state);
141  p *= u;
142  } while (p > L);
143  return k-1;
144 }
145 
146 void RandGauss2(float *a, float *b, RandomState *state) {
147  KALDI_ASSERT(a);
148  KALDI_ASSERT(b);
149  float u1 = RandUniform(state);
150  float u2 = RandUniform(state);
151  u1 = sqrtf(-2.0f * logf(u1));
152  u2 = 2.0f * M_PI * u2;
153  *a = u1 * cosf(u2);
154  *b = u1 * sinf(u2);
155 }
156 
157 void RandGauss2(double *a, double *b, RandomState *state) {
158  KALDI_ASSERT(a);
159  KALDI_ASSERT(b);
160  float a_float, b_float;
161  // Just because we're using doubles doesn't mean we need super-high-quality
162  // random numbers, so we just use the floating-point version internally.
163  RandGauss2(&a_float, &b_float, state);
164  *a = a_float;
165  *b = b_float;
166 }
167 
168 
169 } // end namespace kaldi
170 
171 
Relabels neural network egs with the read pdf-id alignments.
Definition: chain.dox:20
float RandUniform(struct RandomState *state=NULL)
Returns a random number strictly between 0 and 1.
Definition: kaldi-math.h:151
#define M_PI
Definition: kaldi-math.h:44
bool WithProb(BaseFloat prob, struct RandomState *state)
Definition: kaldi-math.cc:77
static pthread_mutex_t _RandMutex
Definition: kaldi-math.cc:43
int32 RoundUpToNearestPowerOfTwo(int32 n)
Definition: kaldi-math.cc:31
int32 RandPoisson(float lambda, struct RandomState *state)
Definition: kaldi-math.cc:133
float BaseFloat
Definition: kaldi-types.h:29
void RandGauss2(float *a, float *b, RandomState *state)
Definition: kaldi-math.cc:146
struct rnnlm::@11::@12 n
int Rand(struct RandomState *state)
Definition: kaldi-math.cc:46
#define KALDI_ASSERT(cond)
Definition: kaldi-error.h:169
int32 RandInt(int32 min_val, int32 max_val, struct RandomState *state)
Definition: kaldi-math.cc:100
#define KALDI_COMPILE_TIME_ASSERT(b)
Definition: kaldi-utils.h:128