src/vppinfra/random_isaac.c

   1 /*
   2  * Copyright (c) 2015 Cisco and/or its affiliates.
   3  * Licensed under the Apache License, Version 2.0 (the "License");
   4  * you may not use this file except in compliance with the License.
   5  * You may obtain a copy of the License at:
   6  *
   7  *     http://www.apache.org/licenses/LICENSE-2.0
   8  *
   9  * Unless required by applicable law or agreed to in writing, software
  10  * distributed under the License is distributed on an "AS IS" BASIS,
  11  * WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
  12  * See the License for the specific language governing permissions and
  13  * limitations under the License.
  14  */
  15 /*
  16   ------------------------------------------------------------------------------
  17   By Bob Jenkins, 1996, Public Domain
  18   MODIFIED:
  19   960327: Creation (addition of randinit, really)
  20   970719: use context, not global variables, for internal state
  21   980324: renamed seed to flag
  22   980605: recommend ISAAC_LOG2_SIZE=4 for noncryptography.
  23   010626: note this is public domain
  24   ------------------------------------------------------------------------------
  25
  26   Modified for CLIB by Eliot Dresselhaus.
  27   Dear Bob, Thanks for all the great work. - Eliot
  28
  29   modifications copyright (c) 2003 Eliot Dresselhaus
  30
  31   Permission is hereby granted, free of charge, to any person obtaining
  32   a copy of this software and associated documentation files (the
  33   "Software"), to deal in the Software without restriction, including
  34   without limitation the rights to use, copy, modify, merge, publish,
  35   distribute, sublicense, and/or sell copies of the Software, and to
  36   permit persons to whom the Software is furnished to do so, subject to
  37   the following conditions:
  38
  39   The above copyright notice and this permission notice shall be
  40   included in all copies or substantial portions of the Software.
  41
  42   THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND,
  43   EXPRESS OR IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF
  44   MERCHANTABILITY, FITNESS FOR A PARTICULAR PURPOSE AND
  45   NONINFRINGEMENT. IN NO EVENT SHALL THE AUTHORS OR COPYRIGHT HOLDERS BE
  46   LIABLE FOR ANY CLAIM, DAMAGES OR OTHER LIABILITY, WHETHER IN AN ACTION
  47   OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, OUT OF OR IN CONNECTION
  48   WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE SOFTWARE.
  49 */
  50
  51 /* ISAAC is Bob Jenkins' random number generator.
  52    http://burtleburtle.net/bob/rand/isaacafa.html */
  53
  54 #include <vppinfra/random_isaac.h>
  55
  56 #if uword_bits != 32 && uword_bits != 64
  57 #error "isaac only works for 32 or 64 bit words"
  58 #endif
  59
  60 #if uword_bits == 32
  61
  62 #define ind32(mm,x)  (*(u32 *)((u8 *)(mm) + ((x) & ((ISAAC_SIZE-1)<<2))))
  63 #define rngstep32(mix,a,b,mm,m,m2,r,x,y)                \
  64 {                                                       \
  65   x = *m;                                               \
  66   a = (a^(mix)) + *(m2++);                              \
  67   *(m++) = y = ind32(mm,x) + a + b;                     \
  68   *(r++) = b = ind32(mm,y>>ISAAC_LOG2_SIZE) + x;        \
  69 }
  70
  71 void
  72 isaac (isaac_t * ctx, uword * results)
  73 {
  74   u32 a, b, c, x, y, *m, *mm, *m2, *r, *mend;
  75
  76   mm = ctx->memory;
  77   r = results;
  78   a = ctx->a;
  79   b = ctx->b;
  80   c = ctx->c;
  81
  82   b += ++c;
  83   mend = m2 = mm + ARRAY_LEN (ctx->memory) / 2;
  84   m = mm;
  85   while (m < mend)
  86     {
  87       rngstep32 (a << 13, a, b, mm, m, m2, r, x, y);
  88       rngstep32 (a >> 6, a, b, mm, m, m2, r, x, y);
  89       rngstep32 (a << 2, a, b, mm, m, m2, r, x, y);
  90       rngstep32 (a >> 16, a, b, mm, m, m2, r, x, y);
  91     }
  92
  93   m2 = mm;
  94   while (m2 < mend)
  95     {
  96       rngstep32 (a << 13, a, b, mm, m, m2, r, x, y);
  97       rngstep32 (a >> 6, a, b, mm, m, m2, r, x, y);
  98       rngstep32 (a << 2, a, b, mm, m, m2, r, x, y);
  99       rngstep32 (a >> 16, a, b, mm, m, m2, r, x, y);
 100     }
 101
 102   ctx->a = a;
 103   ctx->b = b;
 104   ctx->c = c;
 105 }
 106
 107 /* Perform 2 isaac runs with different contexts simultaneously. */
 108 void
 109 isaac2 (isaac_t * ctx, uword * results)
 110 {
 111 #define _(n) \
 112   u32 a##n, b##n, c##n, x##n, y##n, * m##n, * mm##n, * m2##n, * r##n, * mend##n
 113
 114   _(0);
 115   _(1);
 116   (void) mend1;                 /* "set but unused variable" error on mend1 with gcc 4.9  */
 117 #undef _
 118
 119 #define _(n)                                                    \
 120 do {                                                            \
 121   mm##n = ctx[(n)].memory;                                      \
 122   r##n = results + (n) * ISAAC_SIZE;                            \
 123   a##n = ctx[(n)].a;                                            \
 124   b##n = ctx[(n)].b;                                            \
 125   c##n = ctx[(n)].c;                                            \
 126   b##n += ++c##n;                                               \
 127   mend##n = m2##n = mm##n + ARRAY_LEN (ctx[(n)].memory) / 2;    \
 128   m##n = mm##n;                                                 \
 129 } while (0)
 130
 131   _(0);
 132   _(1);
 133
 134 #undef _
 135
 136   while (m0 < mend0)
 137     {
 138       rngstep32 (a0 << 13, a0, b0, mm0, m0, m20, r0, x0, y0);
 139       rngstep32 (a1 << 13, a1, b1, mm1, m1, m21, r1, x1, y1);
 140       rngstep32 (a0 >> 6, a0, b0, mm0, m0, m20, r0, x0, y0);
 141       rngstep32 (a1 >> 6, a1, b1, mm1, m1, m21, r1, x1, y1);
 142       rngstep32 (a0 << 2, a0, b0, mm0, m0, m20, r0, x0, y0);
 143       rngstep32 (a1 << 2, a1, b1, mm1, m1, m21, r1, x1, y1);
 144       rngstep32 (a0 >> 16, a0, b0, mm0, m0, m20, r0, x0, y0);
 145       rngstep32 (a1 >> 16, a1, b1, mm1, m1, m21, r1, x1, y1);
 146     }
 147
 148   m20 = mm0;
 149   m21 = mm1;
 150   while (m20 < mend0)
 151     {
 152       rngstep32 (a0 << 13, a0, b0, mm0, m0, m20, r0, x0, y0);
 153       rngstep32 (a1 << 13, a1, b1, mm1, m1, m21, r1, x1, y1);
 154       rngstep32 (a0 >> 6, a0, b0, mm0, m0, m20, r0, x0, y0);
 155       rngstep32 (a1 >> 6, a1, b1, mm1, m1, m21, r1, x1, y1);
 156       rngstep32 (a0 << 2, a0, b0, mm0, m0, m20, r0, x0, y0);
 157       rngstep32 (a1 << 2, a1, b1, mm1, m1, m21, r1, x1, y1);
 158       rngstep32 (a0 >> 16, a0, b0, mm0, m0, m20, r0, x0, y0);
 159       rngstep32 (a1 >> 16, a1, b1, mm1, m1, m21, r1, x1, y1);
 160     }
 161
 162   ctx[0].a = a0;
 163   ctx[0].b = b0;
 164   ctx[0].c = c0;
 165   ctx[1].a = a1;
 166   ctx[1].b = b1;
 167   ctx[1].c = c1;
 168 }
 169
 170 #define mix32(a,b,c,d,e,f,g,h)                  \
 171 {                                               \
 172    a^=b<<11; d+=a; b+=c;                        \
 173    b^=c>>2;  e+=b; c+=d;                        \
 174    c^=d<<8;  f+=c; d+=e;                        \
 175    d^=e>>16; g+=d; e+=f;                        \
 176    e^=f<<10; h+=e; f+=g;                        \
 177    f^=g>>4;  a+=f; g+=h;                        \
 178    g^=h<<8;  b+=g; h+=a;                        \
 179    h^=a>>9;  c+=h; a+=b;                        \
 180 }
 181
 182 void
 183 isaac_init (isaac_t * ctx, uword * seeds)
 184 {
 185   word i;
 186   u32 a, b, c, d, e, f, g, h, *m, *r;
 187
 188   ctx->a = ctx->b = ctx->c = 0;
 189   m = ctx->memory;
 190   r = seeds;
 191
 192   a = b = c = d = e = f = g = h = 0x9e3779b9;   /* the golden ratio */
 193
 194   for (i = 0; i < 4; ++i)       /* scramble it */
 195     mix32 (a, b, c, d, e, f, g, h);
 196
 197   /* initialize using the contents of r[] as the seed */
 198   for (i = 0; i < ISAAC_SIZE; i += 8)
 199     {
 200       a += r[i];
 201       b += r[i + 1];
 202       c += r[i + 2];
 203       d += r[i + 3];
 204       e += r[i + 4];
 205       f += r[i + 5];
 206       g += r[i + 6];
 207       h += r[i + 7];
 208       mix32 (a, b, c, d, e, f, g, h);
 209       m[i] = a;
 210       m[i + 1] = b;
 211       m[i + 2] = c;
 212       m[i + 3] = d;
 213       m[i + 4] = e;
 214       m[i + 5] = f;
 215       m[i + 6] = g;
 216       m[i + 7] = h;
 217     }
 218
 219   /* do a second pass to make all of the seed affect all of m */
 220   for (i = 0; i < ISAAC_SIZE; i += 8)
 221     {
 222       a += m[i];
 223       b += m[i + 1];
 224       c += m[i + 2];
 225       d += m[i + 3];
 226       e += m[i + 4];
 227       f += m[i + 5];
 228       g += m[i + 6];
 229       h += m[i + 7];
 230       mix32 (a, b, c, d, e, f, g, h);
 231       m[i] = a;
 232       m[i + 1] = b;
 233       m[i + 2] = c;
 234       m[i + 3] = d;
 235       m[i + 4] = e;
 236       m[i + 5] = f;
 237       m[i + 6] = g;
 238       m[i + 7] = h;
 239     }
 240 }
 241 #endif /* uword_bits == 32 */
 242
 243 #if uword_bits == 64
 244
 245 #define ind64(mm,x)  (*(u64 *)((u8 *)(mm) + ((x) & ((ISAAC_SIZE-1)<<3))))
 246 #define rngstep64(mix,a,b,mm,m,m2,r,x,y)                \
 247 {                                                       \
 248   x = *m;                                               \
 249   a = (mix) + *(m2++);                                  \
 250   *(m++) = y = ind64(mm,x) + a + b;                     \
 251   *(r++) = b = ind64(mm,y>>ISAAC_LOG2_SIZE) + x;        \
 252 }
 253
 254 void
 255 isaac (isaac_t * ctx, uword * results)
 256 {
 257   u64 a, b, c, x, y, *m, *mm, *m2, *r, *mend;
 258
 259   mm = ctx->memory;
 260   r = results;
 261   a = ctx->a;
 262   b = ctx->b;
 263   c = ctx->c;
 264
 265   b += ++c;
 266   mend = m2 = mm + ARRAY_LEN (ctx->memory) / 2;
 267   m = mm;
 268   while (m < mend)
 269     {
 270       rngstep64 (~(a ^ (a << 21)), a, b, mm, m, m2, r, x, y);
 271       rngstep64 (a ^ (a >> 5), a, b, mm, m, m2, r, x, y);
 272       rngstep64 (a ^ (a << 12), a, b, mm, m, m2, r, x, y);
 273       rngstep64 (a ^ (a >> 33), a, b, mm, m, m2, r, x, y);
 274     }
 275
 276   m2 = mm;
 277   while (m2 < mend)
 278     {
 279       rngstep64 (~(a ^ (a << 21)), a, b, mm, m, m2, r, x, y);
 280       rngstep64 (a ^ (a >> 5), a, b, mm, m, m2, r, x, y);
 281       rngstep64 (a ^ (a << 12), a, b, mm, m, m2, r, x, y);
 282       rngstep64 (a ^ (a >> 33), a, b, mm, m, m2, r, x, y);
 283     }
 284
 285   ctx->a = a;
 286   ctx->b = b;
 287   ctx->c = c;
 288 }
 289
 290 /* Perform 2 isaac runs with different contexts simultaneously. */
 291 void
 292 isaac2 (isaac_t * ctx, uword * results)
 293 {
 294 #define _(n) \
 295   u64 a##n, b##n, c##n, x##n, y##n, * m##n, * mm##n, * m2##n, * r##n, * mend##n
 296
 297   _(0);
 298   _(1);
 299
 300 #undef _
 301
 302 #define _(n)                                                    \
 303 do {                                                            \
 304   mm##n = ctx[(n)].memory;                                      \
 305   r##n = results + (n) * ISAAC_SIZE;                            \
 306   a##n = ctx[(n)].a;                                            \
 307   b##n = ctx[(n)].b;                                            \
 308   c##n = ctx[(n)].c;                                            \
 309   b##n += ++c##n;                                               \
 310   mend##n = m2##n = mm##n + ARRAY_LEN (ctx[(n)].memory) / 2;    \
 311   m##n = mm##n;                                                 \
 312 } while (0)
 313
 314   _(0);
 315   _(1);
 316
 317 #undef _
 318
 319   (void) mend1;                 /* compiler warning */
 320
 321   while (m0 < mend0)
 322     {
 323       rngstep64 (~(a0 ^ (a0 << 21)), a0, b0, mm0, m0, m20, r0, x0, y0);
 324       rngstep64 (~(a1 ^ (a1 << 21)), a1, b1, mm1, m1, m21, r1, x1, y1);
 325       rngstep64 (a0 ^ (a0 >> 5), a0, b0, mm0, m0, m20, r0, x0, y0);
 326       rngstep64 (a1 ^ (a1 >> 5), a1, b1, mm1, m1, m21, r1, x1, y1);
 327       rngstep64 (a0 ^ (a0 << 12), a0, b0, mm0, m0, m20, r0, x0, y0);
 328       rngstep64 (a1 ^ (a1 << 12), a1, b1, mm1, m1, m21, r1, x1, y1);
 329       rngstep64 (a0 ^ (a0 >> 33), a0, b0, mm0, m0, m20, r0, x0, y0);
 330       rngstep64 (a1 ^ (a1 >> 33), a1, b1, mm1, m1, m21, r1, x1, y1);
 331     }
 332
 333   m20 = mm0;
 334   m21 = mm1;
 335   while (m20 < mend0)
 336     {
 337       rngstep64 (~(a0 ^ (a0 << 21)), a0, b0, mm0, m0, m20, r0, x0, y0);
 338       rngstep64 (~(a1 ^ (a1 << 21)), a1, b1, mm1, m1, m21, r1, x1, y1);
 339       rngstep64 (a0 ^ (a0 >> 5), a0, b0, mm0, m0, m20, r0, x0, y0);
 340       rngstep64 (a1 ^ (a1 >> 5), a1, b1, mm1, m1, m21, r1, x1, y1);
 341       rngstep64 (a0 ^ (a0 << 12), a0, b0, mm0, m0, m20, r0, x0, y0);
 342       rngstep64 (a1 ^ (a1 << 12), a1, b1, mm1, m1, m21, r1, x1, y1);
 343       rngstep64 (a0 ^ (a0 >> 33), a0, b0, mm0, m0, m20, r0, x0, y0);
 344       rngstep64 (a1 ^ (a1 >> 33), a1, b1, mm1, m1, m21, r1, x1, y1);
 345     }
 346
 347   ctx[0].a = a0;
 348   ctx[0].b = b0;
 349   ctx[0].c = c0;
 350   ctx[1].a = a1;
 351   ctx[1].b = b1;
 352   ctx[1].c = c1;
 353 }
 354
 355 #define mix64(a,b,c,d,e,f,g,h)                  \
 356 {                                               \
 357    a-=e; f^=h>>9;  h+=a;                        \
 358    b-=f; g^=a<<9;  a+=b;                        \
 359    c-=g; h^=b>>23; b+=c;                        \
 360    d-=h; a^=c<<15; c+=d;                        \
 361    e-=a; b^=d>>14; d+=e;                        \
 362    f-=b; c^=e<<20; e+=f;                        \
 363    g-=c; d^=f>>17; f+=g;                        \
 364    h-=d; e^=g<<14; g+=h;                        \
 365 }
 366
 367 void
 368 isaac_init (isaac_t * ctx, uword * seeds)
 369 {
 370   word i;
 371   u64 a, b, c, d, e, f, g, h, *m, *r;
 372
 373   ctx->a = ctx->b = ctx->c = 0;
 374   m = ctx->memory;
 375   r = seeds;
 376
 377   a = b = c = d = e = f = g = h = 0x9e3779b97f4a7c13LL; /* the golden ratio */
 378
 379   for (i = 0; i < 4; ++i)       /* scramble it */
 380     mix64 (a, b, c, d, e, f, g, h);
 381
 382   for (i = 0; i < ISAAC_SIZE; i += 8)   /* fill in mm[] with messy stuff */
 383     {
 384       a += r[i];
 385       b += r[i + 1];
 386       c += r[i + 2];
 387       d += r[i + 3];
 388       e += r[i + 4];
 389       f += r[i + 5];
 390       g += r[i + 6];
 391       h += r[i + 7];
 392       mix64 (a, b, c, d, e, f, g, h);
 393       m[i] = a;
 394       m[i + 1] = b;
 395       m[i + 2] = c;
 396       m[i + 3] = d;
 397       m[i + 4] = e;
 398       m[i + 5] = f;
 399       m[i + 6] = g;
 400       m[i + 7] = h;
 401     }
 402
 403   /* do a second pass to make all of the seed affect all of mm */
 404   for (i = 0; i < ISAAC_SIZE; i += 8)
 405     {
 406       a += m[i];
 407       b += m[i + 1];
 408       c += m[i + 2];
 409       d += m[i + 3];
 410       e += m[i + 4];
 411       f += m[i + 5];
 412       g += m[i + 6];
 413       h += m[i + 7];
 414       mix64 (a, b, c, d, e, f, g, h);
 415       m[i] = a;
 416       m[i + 1] = b;
 417       m[i + 2] = c;
 418       m[i + 3] = d;
 419       m[i + 4] = e;
 420       m[i + 5] = f;
 421       m[i + 6] = g;
 422       m[i + 7] = h;
 423     }
 424 }
 425 #endif /* uword_bits == 64 */
 426
 427
 428 /*
 429  * fd.io coding-style-patch-verification: ON
 430  *
 431  * Local Variables:
 432  * eval: (c-set-style "gnu")
 433  * End:
 434  */