~ubuntu-branches/ubuntu/intrepid/ecl/intrepid

« back to all changes in this revision

Viewing changes to src/gmp/mpfr/tests/trandom.c

  • Committer: Bazaar Package Importer
  • Author(s): Peter Van Eynde
  • Date: 2007-04-09 11:51:51 UTC
  • mfrom: (1.1.3 upstream)
  • Revision ID: james.westby@ubuntu.com-20070409115151-ql8cr0kalzx1jmla
Tags: 0.9i-20070324-2
Upload to unstable. 

Show diffs side-by-side

added added

removed removed

Lines of Context:
1
 
/* Test file for the various mpfr_random fonctions.
2
 
 
3
 
Copyright 1999, 2000, 2001, 2002 Free Software Foundation.
4
 
 
5
 
This file is part of the MPFR Library.
6
 
 
7
 
The MPFR Library is free software; you can redistribute it and/or modify
8
 
it under the terms of the GNU Lesser General Public License as published by
9
 
the Free Software Foundation; either version 2.1 of the License, or (at your
10
 
option) any later version.
11
 
 
12
 
The MPFR Library is distributed in the hope that it will be useful, but
13
 
WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY
14
 
or FITNESS FOR A PARTICULAR PURPOSE.  See the GNU Lesser General Public
15
 
License for more details.
16
 
 
17
 
You should have received a copy of the GNU Lesser General Public License
18
 
along with the MPFR Library; see the file COPYING.LIB.  If not, write to
19
 
the Free Software Foundation, Inc., 59 Temple Place - Suite 330, Boston,
20
 
MA 02111-1307, USA. */
21
 
 
22
 
#include <stdio.h>
23
 
#include <stdlib.h>
24
 
#include <math.h>
25
 
#include <time.h>
26
 
#include "gmp.h"
27
 
#include "mpfr.h"
28
 
#include "mpfr-impl.h"
29
 
 
30
 
void test_random _PROTO ((unsigned long, unsigned long, int));
31
 
void test_random2 _PROTO ((unsigned long, unsigned long, int));
32
 
void test_urandomb _PROTO ((unsigned long, unsigned long, int));
33
 
 
34
 
void
35
 
test_random (unsigned long nbtests, unsigned long prec, int verbose)
36
 
{
37
 
  mpfr_t x; 
38
 
  int *tab, size_tab, k; 
39
 
  double d, av = 0, var = 0, chi2 = 0, th; 
40
 
 
41
 
  mpfr_init2(x, prec); 
42
 
 
43
 
  size_tab = (nbtests < 1000 ? nbtests / 50 : 20); 
44
 
  tab = (int *) malloc (size_tab * sizeof(int)); 
45
 
  for (k = 0; k < size_tab; ++k) tab[k] = 0; 
46
 
 
47
 
  for (k = 0; k < nbtests; k++) {
48
 
    mpfr_random(x); 
49
 
    d = mpfr_get_d1 (x); av += d; var += d*d; 
50
 
    tab[(int)(size_tab * d)]++;     
51
 
  }
52
 
 
53
 
  mpfr_clear(x);
54
 
  if (!verbose) { free(tab); return; }
55
 
 
56
 
  av /= nbtests; 
57
 
  var = (var /nbtests) - av*av; 
58
 
 
59
 
  th = (double)nbtests / size_tab; 
60
 
  
61
 
  printf("Average = %.5f\nVariance = %.5f\n", av, var); 
62
 
  printf("Repartition for random. Each integer should be close to %d.\n", 
63
 
         (int)th); 
64
 
 
65
 
  for (k = 0; k < size_tab; k++) {
66
 
    chi2 += (tab[k] - th) * (tab[k] - th) / th; 
67
 
    printf("%d ", tab[k]); 
68
 
    if (((k+1) & 7) == 0) printf("\n"); 
69
 
  }
70
 
 
71
 
  printf("\nChi2 statistics value (with %d degrees of freedom) : %.5f\n\n", 
72
 
         size_tab - 1, chi2); 
73
 
 
74
 
  printf("\n");
75
 
 
76
 
  free(tab);
77
 
  return;
78
 
}
79
 
 
80
 
void
81
 
test_random2 (unsigned long nbtests, unsigned long prec, int verbose)
82
 
{
83
 
  mpfr_t x; 
84
 
  int *tab, size_tab, k; 
85
 
  double d, av = 0, var = 0, chi2 = 0, th; 
86
 
 
87
 
  mpfr_init2(x, prec); 
88
 
 
89
 
  size_tab = (nbtests < 1000 ? nbtests / 50 : 20); 
90
 
  tab = (int *) malloc (size_tab * sizeof(int)); 
91
 
  for (k = 0; k < size_tab; ++k) tab[k] = 0; 
92
 
 
93
 
  for (k = 0; k < nbtests; k++) {
94
 
    mpfr_random2 (x, MPFR_ABSSIZE(x), 0); 
95
 
    d = mpfr_get_d1 (x); av += d; var += d*d; 
96
 
    if (d < 1)
97
 
      tab[(int)(size_tab * d)]++;     
98
 
  }
99
 
 
100
 
  mpfr_clear(x);
101
 
  if (!verbose) { free(tab); return; }
102
 
  
103
 
  av /= nbtests; 
104
 
  var = (var /nbtests) - av*av; 
105
 
 
106
 
  th = (double)nbtests / size_tab; 
107
 
  printf("Average = %.5f\nVariance = %.5f\n", av, var); 
108
 
  printf("Repartition for random2 (taking only values < 1 into account.\n"); 
109
 
 
110
 
  for (k = 0; k < size_tab; k++) {
111
 
    chi2 += (tab[k] - th) * (tab[k] - th) / th; 
112
 
    printf("%d ", tab[k]); 
113
 
    if (((k+1) & 7) == 0) printf("\n"); 
114
 
  }
115
 
 
116
 
  printf("\nChi2 statistics value (with %d degrees of freedom) : %.5f\n\n", 
117
 
         size_tab - 1, chi2); 
118
 
 
119
 
  free(tab); 
120
 
  return;
121
 
}
122
 
 
123
 
void
124
 
test_urandomb (unsigned long nbtests, unsigned long prec, int verbose)
125
 
{
126
 
  mpfr_t x; 
127
 
  int *tab, size_tab, k; 
128
 
  gmp_randstate_t state; 
129
 
  double d, av = 0, var = 0, chi2 = 0, th; 
130
 
 
131
 
  mpfr_init2(x, prec); 
132
 
 
133
 
  size_tab = (nbtests < 1000 ? nbtests / 50 : 20); 
134
 
  tab = (int *) malloc (size_tab * sizeof(int)); 
135
 
  for (k = 0; k < size_tab; ++k) tab[k] = 0; 
136
 
 
137
 
  gmp_randinit (state, GMP_RAND_ALG_LC, 128); 
138
 
  gmp_randseed_ui (state, time(NULL));
139
 
 
140
 
  for (k = 0; k < nbtests; k++) {
141
 
    mpfr_urandomb(x, state); 
142
 
    d = mpfr_get_d1 (x); av += d; var += d*d; 
143
 
    tab[(int)(size_tab * d)]++;     
144
 
  }
145
 
 
146
 
  mpfr_clear(x);
147
 
  gmp_randclear(state);
148
 
  if (!verbose) { free(tab); return; }
149
 
 
150
 
  av /= nbtests; 
151
 
  var = (var /nbtests) - av*av; 
152
 
 
153
 
  th = (double)nbtests / size_tab; 
154
 
  printf("Average = %.5f\nVariance = %.5f\n", av, var); 
155
 
  printf("Repartition for urandomb. Each integer should be close to %d.\n", 
156
 
         (int)th); 
157
 
 
158
 
  for (k = 0; k < size_tab; k++) {
159
 
    chi2 += (tab[k] - th) * (tab[k] - th) / th; 
160
 
    printf("%d ", tab[k]); 
161
 
    if (((k+1) & 7) == 0) printf("\n"); 
162
 
  }
163
 
 
164
 
  printf("\nChi2 statistics value (with %d degrees of freedom) : %.5f\n\n", 
165
 
         size_tab - 1, chi2); 
166
 
 
167
 
  free(tab); 
168
 
  return;
169
 
}
170
 
 
171
 
int
172
 
main (int argc, char *argv[])
173
 
{
174
 
  unsigned long nbtests, prec; int verbose = 0; 
175
 
  
176
 
  if (argc > 1) verbose = 1;
177
 
 
178
 
  if (argc == 1) { nbtests = 10000; } else nbtests = atoi(argv[1]);
179
 
  if (argc <= 2) { prec = 1000; } else prec = atoi(argv[2]); 
180
 
 
181
 
  test_random(nbtests, prec, verbose);
182
 
  test_random2(nbtests, prec, verbose); 
183
 
  test_urandomb(nbtests, prec, verbose);
184
 
  
185
 
  return 0;
186
 
}