1/* Test file for mpfr_urandomb 2 3Copyright 1999, 2000, 2001, 2002, 2003, 2004, 2006, 2007, 2008, 2009, 2010, 2011 Free Software Foundation, Inc. 4Contributed by the Arenaire and Cacao projects, INRIA. 5 6This file is part of the GNU MPFR Library. 7 8The GNU MPFR Library is free software; you can redistribute it and/or modify 9it under the terms of the GNU Lesser General Public License as published by 10the Free Software Foundation; either version 3 of the License, or (at your 11option) any later version. 12 13The GNU MPFR Library is distributed in the hope that it will be useful, but 14WITHOUT ANY WARRANTY; without even the implied warranty of MERCHANTABILITY 15or FITNESS FOR A PARTICULAR PURPOSE. See the GNU Lesser General Public 16License for more details. 17 18You should have received a copy of the GNU Lesser General Public License 19along with the GNU MPFR Library; see the file COPYING.LESSER. If not, see 20http://www.gnu.org/licenses/ or write to the Free Software Foundation, Inc., 2151 Franklin St, Fifth Floor, Boston, MA 02110-1301, USA. */ 22 23#include <stdio.h> 24#include <stdlib.h> 25 26#include "mpfr-test.h" 27 28 29static void 30test_urandomb (long nbtests, mpfr_prec_t prec, int verbose) 31{ 32 mpfr_t x; 33 int *tab, size_tab, k, sh, xn; 34 double d, av = 0, var = 0, chi2 = 0, th; 35 mpfr_exp_t emin; 36 37 size_tab = (nbtests >= 1000 ? nbtests / 50 : 20); 38 tab = (int *) calloc (size_tab, sizeof(int)); 39 if (tab == NULL) 40 { 41 fprintf (stderr, "trandom: can't allocate memory in test_urandomb\n"); 42 exit (1); 43 } 44 45 mpfr_init2 (x, prec); 46 xn = 1 + (prec - 1) / mp_bits_per_limb; 47 sh = xn * mp_bits_per_limb - prec; 48 49 for (k = 0; k < nbtests; k++) 50 { 51 mpfr_urandomb (x, RANDS); 52 /* check that lower bits are zero */ 53 if (MPFR_MANT(x)[0] & MPFR_LIMB_MASK(sh)) 54 { 55 printf ("Error: mpfr_urandomb() returns invalid numbers:\n"); 56 mpfr_print_binary (x); puts (""); 57 exit (1); 58 } 59 d = mpfr_get_d1 (x); av += d; var += d*d; 60 tab[(int)(size_tab * d)]++; 61 } 62 63 /* coverage test */ 64 emin = mpfr_get_emin (); 65 set_emin (1); /* the generated number in [0,1[ is not in the exponent 66 range, except if it is zero */ 67 k = mpfr_urandomb (x, RANDS); 68 if (MPFR_IS_ZERO(x) == 0 && (k == 0 || mpfr_nan_p (x) == 0)) 69 { 70 printf ("Error in mpfr_urandomb, expected NaN, got "); 71 mpfr_dump (x); 72 exit (1); 73 } 74 set_emin (emin); 75 76 mpfr_clear (x); 77 if (!verbose) 78 { 79 free(tab); 80 return; 81 } 82 83 av /= nbtests; 84 var = (var / nbtests) - av * av; 85 86 th = (double)nbtests / size_tab; 87 printf("Average = %.5f\nVariance = %.5f\n", av, var); 88 printf("Repartition for urandomb. Each integer should be close to %d.\n", 89 (int)th); 90 91 for (k = 0; k < size_tab; k++) 92 { 93 chi2 += (tab[k] - th) * (tab[k] - th) / th; 94 printf("%d ", tab[k]); 95 if (((k+1) & 7) == 0) 96 printf("\n"); 97 } 98 99 printf("\nChi2 statistics value (with %d degrees of freedom) : %.5f\n\n", 100 size_tab - 1, chi2); 101 102 free(tab); 103 return; 104} 105 106int 107main (int argc, char *argv[]) 108{ 109 long nbtests; 110 mpfr_prec_t prec; 111 int verbose = 0; 112 113 tests_start_mpfr (); 114 115 if (argc > 1) 116 verbose = 1; 117 118 nbtests = 10000; 119 if (argc > 1) 120 { 121 long a = atol(argv[1]); 122 if (a != 0) 123 nbtests = a; 124 } 125 126 if (argc <= 2) 127 prec = 1000; 128 else 129 prec = atol(argv[2]); 130 131 test_urandomb (nbtests, prec, verbose); 132 133 if (argc == 1) /* check also small precision */ 134 { 135 test_urandomb (nbtests, 2, 0); 136 } 137 138 tests_end_mpfr (); 139 return 0; 140} 141