NEL_Z
[physik/nlsop.git] / random.c
1 /*
2  * random functions
3  *
4  * author: hackbard@hackdaworld.dyndns.org
5  *
6  */
7
8 #include <sys/types.h>
9 #include <sys/stat.h>
10 #include <fcntl.h>
11 #include <unistd.h>
12 #include <stdlib.h>
13 #include <math.h>
14
15 #include "random.h"
16
17 static int rand_fd;
18 static u32 *c_ptr,*b_ptr;
19
20 int rand_init(char *rf)
21 {
22  if(rf==NULL)
23  {
24   if((rand_fd=open("/dev/urandom",O_RDONLY))<0)
25   {
26    puts("cannot open /dev/urandom");
27    return -1;
28   }
29  } else
30  {
31   if((rand_fd=open(rf,O_RDONLY))<0)
32   {
33    printf("cannot open %s\n",rf);
34    return -1;
35   }
36  }
37  if((b_ptr=(u32 *)(malloc(BUFSIZE*sizeof(u32))))==NULL)
38  {
39   puts("failed allocating random buffer");
40   return -1;
41  }
42  c_ptr=b_ptr+BUFSIZE;
43  
44  return 1;
45 }
46
47 int rand_close(void)
48 {
49  close(rand_fd);
50  
51  return 1;
52 }
53
54 u32 get_rand(u32 max)
55 {
56  if(c_ptr>=b_ptr+BUFSIZE)
57  {
58 #ifdef MORE_PRINTF
59   printf("getting another %d bytes of random data ...\n",BUFSIZE);
60 #endif
61   if(read(rand_fd,b_ptr,BUFSIZE*sizeof(u32))<BUFSIZE*sizeof(u32))
62   {
63    /* -> assume random file, end reached */
64    puts("warning, will read random file from beginning ...");
65    lseek(rand_fd,0,SEEK_SET);
66    read(rand_fd,b_ptr,BUFSIZE*sizeof(u32));
67   }
68   c_ptr=b_ptr;
69 #ifdef MORE_PRINTF
70   printf("got it!\n");
71 #endif
72  }
73
74  return((u32)(*(c_ptr++)*(max*1.0/((long long unsigned int)URAND_MAX+1))));
75 }
76
77 u32 get_rand_lgp(u32 max,double a,double b)
78 {
79  return((u32)(1.0*max*(-1.0*b+sqrt(b*b+2*a*((b+a/2)*get_rand(URAND_MAX)/((long long unsigned int)URAND_MAX+1))))/a));
80 }
81
82 u32 get_rand_reject(u32 max_x,u32 max_y,u32 *graph)
83 {
84  u32 x,y;
85  unsigned char ok;
86
87  ok=0;
88  while(!ok)
89  {
90   x=get_rand(max_x);
91   y=get_rand(max_y);
92   if(y<=graph[x]) ok=1;
93  }
94
95  return x;
96 }