← Back to branch summary

~ubuntu-branches/debian/stretch/openbabel/stretch

« back to all changes in this revision

Viewing changes to src/rand.cpp

  • Committer: Bazaar Package Importer
  • Author(s): Michael Banck
  • Date: 2008-07-22 23:54:58 UTC
  • mfrom: (3.1.10 intrepid)
  • Revision ID: james.westby@ubuntu.com-20080722235458-3o606czluviz4akx
Tags: 2.2.0-2
* Upload to unstable.
* debian/control: Updated descriptions.
* debian/patches/gauss_cube_format.patch: New patch, makes the 
  gaussian cube format available again.
* debian/rules (DEB_DH_MAKESHLIBS_ARGS_libopenbabel3): Removed.
* debian/rules (DEB_CONFIGURE_EXTRA_FLAGS): Likewise.
* debian/libopenbabel3.install: Adjust formats directory to single 
  version hierarchy.

Show diffs side-by-side

added added

removed removed

Lines of Context:
src/rand.cpp
2
2
rand.cpp - Pseudo random number generator.
3
3
 
4
4
Copyright (C) 1998-2001 by OpenEye Scientific Software, Inc.
5
 
Some portions Copyright (C) 2001-2005 by Geoffrey R. Hutchison
 
5
Some portions Copyright (C) 2001-2006 by Geoffrey R. Hutchison
6
6
 
7
7
This file is part of the Open Babel project.
8
8
For more information, see <http://openbabel.sourceforge.net/>
22
22
 * Roger Sayle, Metaphorics LLC
23
23
 * Version 1.2, September 1998
24
24
 */
 
25
#include <openbabel/babelconfig.h>
25
26
 
26
27
#include <stdlib.h>
27
28
#include <stdio.h>
28
29
#include <math.h>
29
 
#include "mol.h"
30
 
#include "math/vector3.h"
31
 
#include "obutil.h"
 
30
 
 
31
#include <openbabel/rand.h>
 
32
 
 
33
#if TIME_WITH_SYS_TIME
 
34
#include <sys/time.h>
 
35
#include <time.h>
 
36
#else
 
37
#if HAVE_SYS_TIME_H
 
38
#include <sys/time.h>
 
39
#else
 
40
#include <time.h>
 
41
#endif
 
42
#endif
32
43
 
33
44
#ifndef True
34
45
#define True   1
51
62
namespace OpenBabel
52
63
{
53
64
#define MAXPRIMES  256
54
 
static int primes[MAXPRIMES] = {
55
 
                                   1,    2,    3,    5,    7,   11,   13,   17,   19,   23,
56
 
                                   29,   31,   37,   41,   43,   47,   53,   59,   61,   67,
57
 
                                   71,   73,   79,   83,   89,   97,  101,  103,  107,  109,
58
 
                                   113,  127,  131,  137,  139,  149,  151,  157,  163,  167,
59
 
                                   173,  179,  181,  191,  193,  197,  199,  211,  223,  227,
60
 
                                   229,  233,  239,  241,  251,  257,  263,  269,  271,  277,
61
 
                                   281,  283,  293,  307,  311,  313,  317,  331,  337,  347,
62
 
                                   349,  353,  359,  367,  373,  379,  383,  389,  397,  401,
63
 
                                   409,  419,  421,  431,  433,  439,  443,  449,  457,  461,
64
 
                                   463,  467,  479,  487,  491,  499,  503,  509,  521,  523,
65
 
                                   541,  547,  557,  563,  569,  571,  577,  587,  593,  599,
66
 
                                   601,  607,  613,  617,  619,  631,  641,  643,  647,  653,
67
 
                                   659,  661,  673,  677,  683,  691,  701,  709,  719,  727,
68
 
                                   733,  739,  743,  751,  757,  761,  769,  773,  787,  797,
69
 
                                   809,  811,  821,  823,  827,  829,  839,  853,  857,  859,
70
 
                                   863,  877,  881,  883,  887,  907,  911,  919,  929,  937,
71
 
                                   941,  947,  953,  967,  971,  977,  983,  991,  997, 1009,
72
 
                                   1013, 1019, 1021, 1031, 1033, 1039, 1049, 1051, 1061, 1063,
73
 
                                   1069, 1087, 1091, 1093, 1097, 1103, 1109, 1117, 1123, 1129,
74
 
                                   1151, 1153, 1163, 1171, 1181, 1187, 1193, 1201, 1213, 1217,
75
 
                                   1223, 1229, 1231, 1237, 1249, 1259, 1277, 1279, 1283, 1289,
76
 
                                   1291, 1297, 1301, 1303, 1307, 1319, 1321, 1327, 1361, 1367,
77
 
                                   1373, 1381, 1399, 1409, 1423, 1427, 1429, 1433, 1439, 1447,
78
 
                                   1451, 1453, 1459, 1471, 1481, 1483, 1487, 1489, 1493, 1499,
79
 
                                   1511, 1523, 1531, 1543, 1549, 1553, 1559, 1567, 1571, 1579,
80
 
                                   1583, 1597, 1601, 1607, 1609, 1613
81
 
                               };
 
65
  static int primes[MAXPRIMES] = {
 
66
    1,    2,    3,    5,    7,   11,   13,   17,   19,   23,
 
67
    29,   31,   37,   41,   43,   47,   53,   59,   61,   67,
 
68
    71,   73,   79,   83,   89,   97,  101,  103,  107,  109,
 
69
    113,  127,  131,  137,  139,  149,  151,  157,  163,  167,
 
70
    173,  179,  181,  191,  193,  197,  199,  211,  223,  227,
 
71
    229,  233,  239,  241,  251,  257,  263,  269,  271,  277,
 
72
    281,  283,  293,  307,  311,  313,  317,  331,  337,  347,
 
73
    349,  353,  359,  367,  373,  379,  383,  389,  397,  401,
 
74
    409,  419,  421,  431,  433,  439,  443,  449,  457,  461,
 
75
    463,  467,  479,  487,  491,  499,  503,  509,  521,  523,
 
76
    541,  547,  557,  563,  569,  571,  577,  587,  593,  599,
 
77
    601,  607,  613,  617,  619,  631,  641,  643,  647,  653,
 
78
    659,  661,  673,  677,  683,  691,  701,  709,  719,  727,
 
79
    733,  739,  743,  751,  757,  761,  769,  773,  787,  797,
 
80
    809,  811,  821,  823,  827,  829,  839,  853,  857,  859,
 
81
    863,  877,  881,  883,  887,  907,  911,  919,  929,  937,
 
82
    941,  947,  953,  967,  971,  977,  983,  991,  997, 1009,
 
83
    1013, 1019, 1021, 1031, 1033, 1039, 1049, 1051, 1061, 1063,
 
84
    1069, 1087, 1091, 1093, 1097, 1103, 1109, 1117, 1123, 1129,
 
85
    1151, 1153, 1163, 1171, 1181, 1187, 1193, 1201, 1213, 1217,
 
86
    1223, 1229, 1231, 1237, 1249, 1259, 1277, 1279, 1283, 1289,
 
87
    1291, 1297, 1301, 1303, 1307, 1319, 1321, 1327, 1361, 1367,
 
88
    1373, 1381, 1399, 1409, 1423, 1427, 1429, 1433, 1439, 1447,
 
89
    1451, 1453, 1459, 1471, 1481, 1483, 1487, 1489, 1493, 1499,
 
90
    1511, 1523, 1531, 1543, 1549, 1553, 1559, 1567, 1571, 1579,
 
91
    1583, 1597, 1601, 1607, 1609, 1613
 
92
  };
82
93
 
83
 
static unsigned int isqrt( unsigned int val )
84
 
{
 
94
  static unsigned int isqrt( unsigned int val )
 
95
  {
85
96
    register unsigned int temp;
86
97
    register unsigned int rem;
87
98
    register int i,result;
88
99
 
89
100
    i = 16;
90
101
    while( !(val&((unsigned int)3<<30)) && i )
91
 
    {
 
102
      {
92
103
        val <<= 2;
93
104
        i--;
94
 
    }
 
105
      }
95
106
 
96
107
    if( i )
97
 
    {
 
108
      {
98
109
        rem = (val>>30)-1;
99
110
        val <<= 2;
100
111
        result = 1;
101
112
        i--;
102
113
 
103
114
        while( i )
104
 
        {
 
115
          {
105
116
            rem = (rem<<2) | (val>>30);
106
117
            result <<= 1;
107
118
            val <<= 2;
108
119
 
109
120
            temp = result<<1;
110
121
            if( rem > temp )
111
 
            {
 
122
              {
112
123
                rem -= temp|1;
113
124
                result |= 1;
114
 
            }
 
125
              }
115
126
            i--;
116
 
        }
 
127
          }
117
128
        return result;
118
 
    }
 
129
      }
119
130
    else
120
 
        return 0;
121
 
}
 
131
      return 0;
 
132
  }
122
133
 
123
 
static int IsOddPrime( unsigned int x )
124
 
{
 
134
  static int IsOddPrime( unsigned int x )
 
135
  {
125
136
    register unsigned int root;
126
137
    register unsigned int i;
127
138
 
128
139
    root = isqrt(x);
129
140
    for( i=2; i<MAXPRIMES-1; i++ )
130
 
    {
 
141
      {
131
142
        if( (x%primes[i]) == 0 )
132
 
            return False;
 
143
          return False;
133
144
        if( (unsigned int) primes[i] >= root )
134
 
            return True;
135
 
    }
 
145
          return True;
 
146
      }
136
147
 
137
148
    for( i=primes[MAXPRIMES-1]; i<=root; i+=2 )
138
 
        if( (x%i) == 0 )
139
 
            return False;
 
149
      if( (x%i) == 0 )
 
150
        return False;
140
151
    return True;
141
 
}
 
152
  }
142
153
 
143
 
static int RelativelyPrime( unsigned int x, unsigned int y )
144
 
{
 
154
  static int RelativelyPrime( unsigned int x, unsigned int y )
 
155
  {
145
156
    if( BothEven(x,y) )
146
 
        return False;
 
157
      return False;
147
158
 
148
159
    if( IsEven(x) )
149
 
    {
 
160
      {
150
161
        do
151
 
        {
 
162
          {
152
163
            x >>= 1;
153
 
        }
 
164
          }
154
165
        while( IsEven(x) );
155
 
    }
 
166
      }
156
167
    else
157
 
        while( IsEven(y) )
158
 
            y >>= 1;
 
168
      while( IsEven(y) )
 
169
        y >>= 1;
159
170
 
160
171
    while( x != y )
161
 
    {
 
172
      {
162
173
        if( x > y )
163
 
        {
 
174
          {
164
175
            x -= y;
165
176
            do
166
 
            {
 
177
              {
167
178
                x >>= 1;
168
 
            }
 
179
              }
169
180
            while( IsEven(x) );
170
 
        }
 
181
          }
171
182
        else if( x < y )
172
 
        {
 
183
          {
173
184
            y -= x;
174
185
            do
175
 
            {
 
186
              {
176
187
                y >>= 1;
177
 
            }
 
188
              }
178
189
            while( IsEven(y) );
179
 
        }
180
 
    }
 
190
          }
 
191
      }
181
192
    return( x == 1 );
182
 
}
 
193
  }
183
194
 
184
 
void DoubleAdd( DoubleType *x, unsigned int y )
185
 
{
 
195
  void DoubleAdd( DoubleType *x, unsigned int y )
 
196
  {
186
197
    x->lo += y;
187
198
    if( x->lo < y )
188
 
        x->hi++;
189
 
}
 
199
      x->hi++;
 
200
  }
190
201
 
191
 
void DoubleMultiply( unsigned int x, unsigned int y, DoubleType *z )
192
 
{
 
202
  void DoubleMultiply( unsigned int x, unsigned int y, DoubleType *z )
 
203
  {
193
204
    register unsigned int x0, x1, x2, x3;
194
205
    register unsigned int hx, lx;
195
206
    register unsigned int hy, ly;
207
218
    x1 += HiPart(x0);
208
219
    x1 += x2;
209
220
    if( x1 < x2 )
210
 
        x3 += (1<<16);
 
221
      x3 += (1<<16);
211
222
 
212
223
    z->hi = HiPart(x1) + x3;
213
224
    z->lo = (LoPart(x1)<<16) + LoPart(x0);
214
 
}
 
225
  }
215
226
 
216
 
static int LeadingZeros( unsigned int x )
217
 
{
 
227
  static int LeadingZeros( unsigned int x )
 
228
  {
218
229
    static int table[256] = {
219
 
                                0,1,2,2,3,3,3,3,4,4,4,4,4,4,4,4,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
220
 
                                6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
221
 
                                7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
222
 
                                7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
223
 
                                8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,
224
 
                                8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,
225
 
                                8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,
226
 
                                8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8
227
 
                            };
 
230
      0,1,2,2,3,3,3,3,4,4,4,4,4,4,4,4,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,5,
 
231
      6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,6,
 
232
      7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
 
233
      7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,7,
 
234
      8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,
 
235
      8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,
 
236
      8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,
 
237
      8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8,8
 
238
    };
228
239
 
229
240
    if( x >= (1<<16) )
230
 
    {
 
241
      {
231
242
        if( x >= (1<<24) )
232
 
        {
 
243
          {
233
244
            return  8-table[x>>24];
234
 
        }
 
245
          }
235
246
        else
236
 
            return 16-table[x>>16];
237
 
    }
 
247
          return 16-table[x>>16];
 
248
      }
238
249
    else if( x >= (1<<8) )
239
 
    {
 
250
      {
240
251
        return 24-table[x>>8];
241
 
    }
 
252
      }
242
253
    else
243
 
        return 32-table[x];
244
 
}
 
254
      return 32-table[x];
 
255
  }
245
256
 
246
 
unsigned int DoubleModulus( DoubleType *n,  unsigned int d )
247
 
{
 
257
  unsigned int DoubleModulus( DoubleType *n,  unsigned int d )
 
258
  {
248
259
    register unsigned int d1, d0;
249
260
    register unsigned int r1, r0;
250
261
    register unsigned int m,s;
251
262
 
252
263
    s = LeadingZeros(d);
253
264
    if( s > 0 )
254
 
    {
 
265
      {
255
266
        d = d<<s;
256
267
        n->hi = (n->hi<<s) | (n->lo>>(32-s));
257
268
        n->lo = n->lo << s;
258
 
    }
 
269
      }
259
270
 
260
271
    d1 = HiPart(d);
261
272
    d0 = LoPart(d);
263
274
    m = ((unsigned int)(n->hi/d1)) * d0;
264
275
    r1 = ((n->hi%d1)<<16) + HiPart(n->lo);
265
276
    if( r1 < m )
266
 
    {
 
277
      {
267
278
        r1 += d;
268
279
        if( (r1>=d) && (r1<m) )
269
 
            r1 += d;
270
 
    }
 
280
          r1 += d;
 
281
      }
271
282
    r1 -= m;
272
283
 
273
284
    m = ((unsigned int)(r1/d1)) * d0;
274
285
    r0 = ((r1%d1)<<16) + LoPart(n->lo);
275
286
    if( r0 < m )
276
 
    {
 
287
      {
277
288
        r0 += d;
278
289
        if( (r0>=d) && (r0<m) )
279
 
            r0 += d;
280
 
    }
 
290
          r0 += d;
 
291
      }
281
292
    r0 -= m;
282
293
 
283
294
    return r0 >> s;
284
 
}
 
295
  }
285
296
 
286
 
static int DeterminePotency( unsigned int m, unsigned int a )
287
 
{
 
297
  static int DeterminePotency( unsigned int m, unsigned int a )
 
298
  {
288
299
    auto DoubleType d;
289
300
    register unsigned int k;
290
301
    register unsigned int b;
294
305
    k = b;
295
306
    s = 1;
296
307
    while( (k!=0) && (s<100) )
297
 
    {
 
308
      {
298
309
        DoubleMultiply(k,b,&d);
299
310
        k = DoubleModulus(&d,m);
300
311
        s++;
301
 
    }
 
312
      }
302
313
    return s;
303
 
}
 
314
  }
304
315
 
305
 
static int DetermineFactors( unsigned int x, unsigned int *factors )
306
 
{
 
316
  static int DetermineFactors( unsigned int x, unsigned int *factors )
 
317
  {
307
318
    register unsigned int *ptr;
308
319
    register unsigned int half;
309
320
    register unsigned int i;
311
322
    half = x/2;
312
323
    ptr = factors;
313
324
    for( i=1; i<MAXPRIMES; i++ )
314
 
    {
 
325
      {
315
326
        if( (x%primes[i]) == 0 )
316
 
            *ptr++ = primes[i];
 
327
          *ptr++ = primes[i];
317
328
        if( (unsigned int)(primes[i]) >= half )
318
 
            return ptr-factors;
319
 
    }
 
329
          return ptr-factors;
 
330
      }
320
331
 
321
332
    for( i=primes[MAXPRIMES-1]+2; i<=half; i+=2 )
322
 
        if( IsOddPrime(i) && ((x%i)==0) )
323
 
            *ptr++ = i;
 
333
      if( IsOddPrime(i) && ((x%i)==0) )
 
334
        *ptr++ = i;
324
335
    return ptr-factors;
325
 
}
 
336
  }
326
337
 
327
 
static unsigned int DetermineIncrement( unsigned int m )
328
 
{
 
338
  static unsigned int DetermineIncrement( unsigned int m )
 
339
  {
329
340
    register unsigned int hi,lo;
330
341
    register unsigned int half;
331
342
    register unsigned int i;
333
344
    /* 1/2 + sqrt(3)/6 */
334
345
    hi = (int)floor(0.7886751345948*m+0.5);
335
346
    if( RelativelyPrime(m,hi) )
336
 
        return hi;
 
347
      return hi;
337
348
 
338
349
    /* 1/2 - sqrt(3)/6 */
339
350
    lo = (int)floor(0.2113248654052*m+0.5);
340
351
    if( RelativelyPrime(m,lo) )
341
 
        return lo;
 
352
      return lo;
342
353
 
343
354
    half = m/2;
344
355
    for( i=1; i<half; i++ )
345
 
    {
 
356
      {
346
357
        if( RelativelyPrime(m,hi+i) )
347
 
            return hi+i;
 
358
          return hi+i;
348
359
        if( RelativelyPrime(m,hi-i) )
349
 
            return hi-i;
 
360
          return hi-i;
350
361
        if( RelativelyPrime(m,lo+i) )
351
 
            return lo+i;
 
362
          return lo+i;
352
363
        if( RelativelyPrime(m,lo-i) )
353
 
            return lo-i;
354
 
    }
 
364
          return lo-i;
 
365
      }
355
366
    return 1;
356
 
}
 
367
  }
357
368
 
358
 
int DetermineSequence( unsigned int m, unsigned int *pm,
359
 
                       unsigned int *pa,
360
 
                       unsigned int *pc )
361
 
{
 
369
  int DetermineSequence( unsigned int m, unsigned int *pm,
 
370
                         unsigned int *pa,
 
371
                         unsigned int *pc )
 
372
  {
362
373
    auto unsigned int fact[MAXFACT];
363
374
    register unsigned int a=0, c;
364
375
    register unsigned int b;
368
379
    register int i;
369
380
 
370
381
    do
371
 
    {
 
382
      {
372
383
        best = 0;
373
384
 
374
385
        count = DetermineFactors(m,fact);
375
386
        if( (m&3) == 0 )
376
 
            fact[0] = 4;
377
 
 
378
 
#ifdef DEBUG
379
 
 
380
 
        printf("factors(%d):",m);
381
 
#endif
 
387
          fact[0] = 4;
382
388
 
383
389
        if( count )
384
 
        {
385
 
#ifdef DEBUG
386
 
            printf("%u",fact[0]);
387
 
            for( i=1; i<count; i++ )
388
 
                printf(",%u",fact[i]);
389
 
#endif
390
 
 
 
390
          {
391
391
            for( b=m-2; b>0; b-- )
392
 
            {
 
392
              {
393
393
                flag = True;
394
394
                for( i=0; i<count; i++ )
395
 
                    if( b%fact[i] )
 
395
                  if( b%fact[i] )
396
396
                    {
397
 
                        flag = False;
398
 
                        break;
 
397
                      flag = False;
 
398
                      break;
399
399
                    }
400
400
 
401
401
                if( flag )
402
 
                {
 
402
                  {
403
403
                    pot = DeterminePotency(m,b+1);
404
404
                    if( pot > best )
405
 
                    {
 
405
                      {
406
406
                        best = pot;
407
407
                        a = b+1;
408
 
                    }
409
 
                }
410
 
            }
411
 
 
412
 
#ifdef DEBUG
413
 
            if( best )
414
 
                printf(" [%d]",best);
415
 
#endif
416
 
 
417
 
        }
418
 
#ifdef DEBUG
419
 
        printf("\n");
420
 
#endif
 
408
                      }
 
409
                  }
 
410
              }
 
411
 
 
412
          }
421
413
 
422
414
        m++;
423
 
    }
 
415
      }
424
416
    while( best < 3 );
425
417
    m--;
426
418
 
430
422
    *pa = a;
431
423
    *pc = c;
432
424
    return best;
433
 
}
434
 
 
435
 
 
436
 
void GenerateSequence( unsigned int p, unsigned int m,
437
 
                       unsigned int a, unsigned int c )
438
 
{
 
425
  }
 
426
 
 
427
 
 
428
  void GenerateSequence( unsigned int p, unsigned int m,
 
429
                         unsigned int a, unsigned int c )
 
430
  {
439
431
    register unsigned int i;
440
432
    register unsigned int x;
441
433
    auto DoubleType d;
442
434
 
443
435
    x = 0;  /* seed */
444
436
    for( i=0; i<p; i++ )
445
 
    {
446
 
        printf("%u\n",x);
 
437
      {
 
438
        //        printf("%u\n",x);
447
439
 
448
440
        do
449
 
        {
 
441
          {
450
442
            DoubleMultiply(a,x,&d);
451
443
            DoubleAdd(&d,c);
452
444
            x = DoubleModulus(&d,m);
453
 
        }
 
445
          }
454
446
        while( x >= p );
455
 
    }
456
 
}
457
 
 
458
 
//**********************************************
459
 
//***** Member functions from Random class *****
460
 
//**********************************************
461
 
 
462
 
OBRandom::OBRandom(bool useSysRand)
463
 
{
 
447
      }
 
448
  }
 
449
 
 
450
  //**********************************************
 
451
  //***** Member functions from Random class *****
 
452
  //**********************************************
 
453
 
 
454
  OBRandom::OBRandom(bool useSysRand)
 
455
  {
464
456
 
465
457
    this->OBRandomUseSysRand= useSysRand;
466
458
    p = 70092;
467
459
    DetermineSequence(p,&m,&a,&c);
468
460
    x = 0;  /* seed */
469
 
}
 
461
  }
470
462
 
471
 
int OBRandom::NextInt()
472
 
{
 
463
  int OBRandom::NextInt()
 
464
  {
473
465
    if (OBRandomUseSysRand)
474
 
    {
 
466
      {
475
467
        return(rand());
476
 
    }
 
468
      }
477
469
    do
478
 
    {
 
470
      {
479
471
        DoubleMultiply(a,x,&d);
480
472
        DoubleAdd(&d,c);
481
473
        x = DoubleModulus(&d,m);
482
 
    }
 
474
      }
483
475
    while( x >= p );
484
476
 
485
477
    return(x);
486
 
}
 
478
  }
487
479
 
488
 
double OBRandom::NextFloat()
489
 
{
 
480
  double OBRandom::NextFloat()
 
481
  {
490
482
 
491
483
    if (OBRandomUseSysRand)
492
 
    {
 
484
      {
493
485
        return(double(rand())/double(RAND_MAX));
494
 
    }
 
486
      }
495
487
    do
496
 
    {
 
488
      {
497
489
        DoubleMultiply(a,x,&d);
498
490
        DoubleAdd(&d,c);
499
491
        x = DoubleModulus(&d,m);
500
 
    }
 
492
      }
501
493
    while( x >= p );
502
494
 
503
495
    return((double)x/p);
504
 
}
 
496
  }
505
497
 
506
 
void OBRandom::TimeSeed()
507
 
{
 
498
  void OBRandom::TimeSeed()
 
499
  {
508
500
#ifdef WIN32
509
501
    // for VC++ do it this way
510
502
    time_t ltime;
523
515
#endif
524
516
 
525
517
#endif
526
 
}
 
518
  }
527
519
 
528
520
} //end namespace OpenBabel
529
521