aboutsummaryrefslogtreecommitdiff
path: root/lmathlib.c
diff options
context:
space:
mode:
authorRoberto Ierusalimschy <roberto@inf.puc-rio.br>2018-03-09 16:23:39 -0300
committerRoberto Ierusalimschy <roberto@inf.puc-rio.br>2018-03-09 16:23:39 -0300
commitdbec41f34ca26f4762bbda3c6e99ba227508d743 (patch)
tree8b3fc30445bc508a604a71197f312f8d22ff7d38 /lmathlib.c
parent0b3db69e4102c148b45a4491ef533ff449b0b927 (diff)
downloadlua-dbec41f34ca26f4762bbda3c6e99ba227508d743.tar.gz
lua-dbec41f34ca26f4762bbda3c6e99ba227508d743.tar.bz2
lua-dbec41f34ca26f4762bbda3c6e99ba227508d743.zip
random floats of different sizes get exactly needed number of random bits
(up to 64)
Diffstat (limited to 'lmathlib.c')
-rw-r--r--lmathlib.c88
1 files changed, 58 insertions, 30 deletions
diff --git a/lmathlib.c b/lmathlib.c
index 062da4a0..2b5fd46d 100644
--- a/lmathlib.c
+++ b/lmathlib.c
@@ -1,5 +1,5 @@
1/* 1/*
2** $Id: lmathlib.c,v 1.121 2018/03/09 14:56:25 roberto Exp roberto $ 2** $Id: lmathlib.c,v 1.122 2018/03/09 15:05:13 roberto Exp roberto $
3** Standard mathematical library 3** Standard mathematical library
4** See Copyright Notice in lua.h 4** See Copyright Notice in lua.h
5*/ 5*/
@@ -10,9 +10,10 @@
10#include "lprefix.h" 10#include "lprefix.h"
11 11
12 12
13#include <float.h>
13#include <limits.h> 14#include <limits.h>
14#include <stdlib.h>
15#include <math.h> 15#include <math.h>
16#include <stdlib.h>
16 17
17#include "lua.h" 18#include "lua.h"
18 19
@@ -250,11 +251,17 @@ static int math_type (lua_State *L) {
250** =================================================================== 251** ===================================================================
251*/ 252*/
252 253
254/* number of binary digits in the mantissa of a float */
255#define FIGS l_mathlim(MANT_DIG)
253 256
254#define twotomin53 (1.0 / 9007199254740992.0) /* 2^-53 */ 257#if FIGS > 64
258/* there are only 64 random bits; use them all */
259#undef FIGS
260#define FIGS 64
261#endif
255 262
256 263
257#if defined(LLONG_MAX) && !defined(LUA_DEBUG) /* { */ 264#if !defined(LUA_USE_C89) && defined(LLONG_MAX) && !defined(LUA_DEBUG) /* { */
258 265
259/* 266/*
260** Assume long long. 267** Assume long long.
@@ -272,15 +279,17 @@ static I xorshift128plus (I *state) {
272 return state[1] + y; 279 return state[1] + y;
273} 280}
274 281
282/* must take care to not shift stuff by more than 63 slots */
275 283
276#define mask53 (~(~0LLU << 53)) 284#define maskFIG (~(~1LLU << (FIGS - 1))) /* use FIGS bits */
285#define shiftFIG (l_mathop(0.5) / (1LLU << (FIGS - 1))) /* 2^(-FIG) */
277 286
278/* 287/*
279** Convert 53 bits from a random integer into a double in the 288** Convert bits from a random integer into a float in the
280** interval [0,1). 289** interval [0,1).
281*/ 290*/
282static double I2d (I x) { 291static lua_Number I2d (I x) {
283 return (x & mask53) * twotomin53; 292 return (lua_Number)(x & maskFIG) * shiftFIG;
284} 293}
285 294
286/* convert an 'I' to a lua_Unsigned */ 295/* convert an 'I' to a lua_Unsigned */
@@ -289,10 +298,10 @@ static double I2d (I x) {
289/* convert a lua_Integer to an 'I' */ 298/* convert a lua_Integer to an 'I' */
290#define Int2I(x) ((I)(x)) 299#define Int2I(x) ((I)(x))
291 300
292#else /* }{ */ 301#else /* no long long }{ */
293 302
294/* 303/*
295** No long long; Use two 32-bit integers to represent a 64-bit quantity. 304** Use two 32-bit integers to represent a 64-bit quantity.
296*/ 305*/
297 306
298#if LUAI_BITSINT >= 32 307#if LUAI_BITSINT >= 32
@@ -303,7 +312,8 @@ typedef unsigned long lu_int32;
303 312
304/* a 64-bit value */ 313/* a 64-bit value */
305typedef struct I { 314typedef struct I {
306 lu_int32 x1, x2; 315 lu_int32 h; /* higher half */
316 lu_int32 l; /* lower half */
307} I; 317} I;
308 318
309 319
@@ -311,31 +321,31 @@ typedef struct I {
311** basic operations on 'I' values 321** basic operations on 'I' values
312*/ 322*/
313 323
314static I pack (int x1, int x2) { 324static I pack (int h, int l) {
315 I result; 325 I result;
316 result.x1 = x1; 326 result.h = h;
317 result.x2 = x2; 327 result.l = l;
318 return result; 328 return result;
319} 329}
320 330
321/* i ^ (i << n) */ 331/* i ^ (i << n) */
322static I Ixorshl (I i, int n) { 332static I Ixorshl (I i, int n) {
323 return pack(i.x1 ^ ((i.x1 << n) | (i.x2 >> (32 - n))), i.x2 ^ (i.x2 << n)); 333 return pack(i.h ^ ((i.h << n) | (i.l >> (32 - n))), i.l ^ (i.l << n));
324} 334}
325 335
326/* i ^ (i >> n) */ 336/* i ^ (i >> n) */
327static I Ixorshr (I i, int n) { 337static I Ixorshr (I i, int n) {
328 return pack(i.x1 ^ (i.x1 >> n), i.x2 ^ ((i.x2 >> n) | (i.x1 << (32 - n)))); 338 return pack(i.h ^ (i.h >> n), i.l ^ ((i.l >> n) | (i.h << (32 - n))));
329} 339}
330 340
331static I Ixor (I i1, I i2) { 341static I Ixor (I i1, I i2) {
332 return pack(i1.x1 ^ i2.x1, i1.x2 ^ i2.x2); 342 return pack(i1.h ^ i2.h, i1.l ^ i2.l);
333} 343}
334 344
335static I Iadd (I i1, I i2) { 345static I Iadd (I i1, I i2) {
336 I result = pack(i1.x1 + i2.x1, i1.x2 + i2.x2); 346 I result = pack(i1.h + i2.h, i1.l + i2.l);
337 if (result.x2 < i1.x2) /* carry? */ 347 if (result.l < i1.l) /* carry? */
338 result.x1++; 348 result.h++;
339 return result; 349 return result;
340} 350}
341 351
@@ -355,28 +365,46 @@ static I xorshift128plus (I *state) {
355 365
356 366
357/* 367/*
358** Converts an 'I' into a double, getting its lower half plus 21 368** Converts an 'I' into a float.
359** (53 - 32) bits from its higher half and joining them into a double.
360*/ 369*/
361 370
362#define mask32 0xffffffff 371#if FIGS <= 32
363#define mask21 (~(~0U << 21)) 372
373/* do not need bits from higher half */
374#define maskHF 0
375#define maskLOW (~(~1U << (FIGS - 1))) /* use FIG bits */
376#define shiftFIG (0.5 / (1U << (FIGS - 1))) /* 2^(-FIG) */
377
378#else /* 32 < FIGS <= 64 */
379
380/* must take care to not shift stuff by more than 31 slots */
381
382/* use FIG - 32 bits from higher half */
383#define maskHF (~(~1U << (FIGS - 33)))
384
385/* use all bits from lower half */
386#define maskLOW (~0)
387
388/* 2^(-FIG) == (1 / 2^33) / 2^(FIG-33) */
389#define shiftFIG ((lua_Number)(1.0 / 8589934592.0) / (1U << (FIGS - 33)))
390
391#endif
364 392
365#define twoto32 4294967296.0 /* 2^32 */ 393#define twoto32 l_mathop(4294967296.0) /* 2^32 */
366 394
367static double I2d (I x) { 395static lua_Number I2d (I x) {
368 return ((x.x1 & mask21) * twoto32 + (x.x2 & mask32)) * twotomin53; 396 return ((x.h & maskHF) * twoto32 + (x.l & maskLOW)) * shiftFIG;
369} 397}
370 398
371static lua_Unsigned I2UInt (I x) { 399static lua_Unsigned I2UInt (I x) {
372 return ((lua_Unsigned)x.x1 << 31 << 1) | x.x2; 400 return ((lua_Unsigned)x.h << 31 << 1) | x.l;
373} 401}
374 402
375static I Int2I (lua_Integer n) { 403static I Int2I (lua_Integer n) {
376 return pack(n, n >> 31 >> 1); 404 return pack(n, n >> 31 >> 1);
377} 405}
378 406
379#endif /* } */ 407#endif /* } */
380 408
381 409
382/* 410/*
@@ -435,7 +463,7 @@ static int math_random (lua_State *L) {
435 I rv = xorshift128plus(state->s); /* next pseudo-random value */ 463 I rv = xorshift128plus(state->s); /* next pseudo-random value */
436 switch (lua_gettop(L)) { /* check number of arguments */ 464 switch (lua_gettop(L)) { /* check number of arguments */
437 case 0: { /* no arguments */ 465 case 0: { /* no arguments */
438 lua_pushnumber(L, (lua_Number)I2d(rv)); /* float between 0 and 1 */ 466 lua_pushnumber(L, I2d(rv)); /* float between 0 and 1 */
439 return 1; 467 return 1;
440 } 468 }
441 case 1: { /* only upper limit */ 469 case 1: { /* only upper limit */