htswanalysis/MACS/lib/gsl/gsl-1.11/qrng/test.c

   1 /* Author: G. Jungman (+modifications from O. Teytaud olivier.teytaud@inria.fr)
   2  */
   3 #include <config.h>
   4 #include <stdlib.h>
   5 #include <gsl/gsl_ieee_utils.h>
   6
   7 #include <gsl/gsl_qrng.h>
   8 #include <gsl/gsl_test.h>
   9 #include <math.h>
  10
  11 void test_sobol(void)
  12 {
  13   int status = 0;
  14   double v[3];
  15   /* int i; */
  16
  17   /* test in dimension 2 */
  18   gsl_qrng * g = gsl_qrng_alloc(gsl_qrng_sobol, 2);
  19   gsl_qrng_get(g, v);
  20   gsl_qrng_get(g, v);
  21   gsl_qrng_get(g, v);
  22   status += ( v[0] != 0.25 || v[1] != 0.75 );
  23   gsl_qrng_get(g, v);
  24   status += ( v[0] != 0.375 || v[1] != 0.375 );
  25   gsl_qrng_free(g);
  26
  27   gsl_test (status, "Sobol d=2");
  28
  29   status = 0;
  30   /* test in dimension 3 */
  31   g = gsl_qrng_alloc(gsl_qrng_sobol, 3);
  32   gsl_qrng_get(g, v);
  33   gsl_qrng_get(g, v);
  34   gsl_qrng_get(g, v);
  35   status += ( v[0] != 0.25 || v[1] != 0.75 || v[2] != 0.25 );
  36   gsl_qrng_get(g, v);
  37   status += ( v[0] != 0.375 || v[1] != 0.375 || v[2] != 0.625 );
  38
  39   gsl_test (status, "Sobol d=3");
  40
  41   status = 0;
  42   gsl_qrng_init(g);
  43   gsl_qrng_get(g, v);
  44   gsl_qrng_get(g, v);
  45   gsl_qrng_get(g, v);
  46   status += ( v[0] != 0.25 || v[1] != 0.75 || v[2] != 0.25 );
  47   gsl_qrng_get(g, v);
  48   status += ( v[0] != 0.375 || v[1] != 0.375 || v[2] != 0.625 );
  49   gsl_qrng_free(g);
  50
  51   gsl_test (status, "Sobol d=3 (reinitialized)");
  52 }
  53
  54 void test_halton(void)
  55 {
  56         int status = 0;
  57         double v[1229];
  58         unsigned int i;
  59
  60         /* test in dimension 1229 */
  61
  62         gsl_qrng * g = gsl_qrng_alloc(gsl_qrng_halton, 1229);
  63         for (i=0;i<30;i++)
  64                 gsl_qrng_get(g, v);
  65         gsl_qrng_free(g);
  66
  67         gsl_test (status, "Halton d=1229");
  68         status = 0;
  69
  70         /* test in dimension 2 */
  71         /*should be
  72          * 0.5 0.333333
  73          * 0.25 0.666667
  74          * 0.75 0.111111
  75          * 0.125 0.444444*/
  76         g = gsl_qrng_alloc(gsl_qrng_halton, 2);
  77         gsl_qrng_get(g, v);
  78         gsl_qrng_get(g, v);
  79         gsl_qrng_get(g, v);
  80         gsl_test_rel (v[0], 3.0/4.0, 1e-3, "halton(2) k=2 v[0]");
  81         gsl_test_rel (v[1], 1.0/9.0, 1e-3, "halton(2) k=2 v[1]");
  82
  83         gsl_qrng_get(g, v);
  84         gsl_test_rel (v[0], 1.0/8.0, 1e-3, "halton(2) k=3 v[0]");
  85         gsl_test_rel (v[1], 4.0/9.0, 1e-3, "halton(2) k=3 v[1]");
  86
  87         gsl_qrng_free(g);
  88
  89         /* test in dimension 3 */
  90         g = gsl_qrng_alloc(gsl_qrng_halton, 3);
  91         gsl_qrng_get(g, v);
  92         gsl_qrng_get(g, v);
  93         gsl_qrng_get(g, v);
  94         gsl_test_rel (v[0], 0.75, 1e-3, "halton(3) k=3 v[0]");
  95         gsl_test_rel (v[1], 1.0/9.0, 1e-3, "halton(3) k=3 v[1]");
  96         gsl_test_rel (v[2], 0.6, 1e-3, "halton(3) k=3 v[2]");
  97
  98         gsl_qrng_get(g, v);
  99         gsl_test_rel (v[0], 0.125, 1e-3, "halton(3) k=4 v[0]");
 100         gsl_test_rel (v[1], 4.0/9.0, 1e-3, "halton(3) k=4 v[1]");
 101         gsl_test_rel (v[2], 0.8, 1e-3, "halton(3) k=4 v[2]");
 102
 103         gsl_qrng_init(g);
 104         gsl_qrng_get(g, v);
 105         gsl_qrng_get(g, v);
 106         gsl_qrng_get(g, v);
 107         gsl_test_rel (v[0], 0.75, 1e-3, "halton(3) reinitialized k=3 v[0]");
 108         gsl_test_rel (v[1], 1.0/9.0, 1e-3, "halton(3) reinitialized k=3 v[1]");
 109         gsl_test_rel (v[2], 0.6, 1e-3, "halton(3) reinitialized k=3 v[2]");
 110
 111         gsl_qrng_get(g, v);
 112         gsl_test_rel (v[0], 0.125, 1e-3, "halton(3) reinitialized k=4 v[0]");
 113         gsl_test_rel (v[1], 4.0/9.0, 1e-3, "halton(3) reinitialized k=4 v[1]");
 114         gsl_test_rel (v[2], 0.8, 1e-3, "halton(3) reinitialized k=4 v[2]");
 115
 116         gsl_qrng_free(g);
 117 }
 118
 119 void test_reversehalton(void)
 120 {
 121         int status = 0;
 122         double v[3];
 123
 124         /* test in dimension 2 */
 125         gsl_qrng * g = gsl_qrng_alloc(gsl_qrng_reversehalton, 2);
 126         gsl_qrng_get(g, v);
 127         gsl_qrng_get(g, v);
 128         gsl_qrng_get(g, v);
 129         /* should be
 130          * 0.5 0.666667
 131          * 0.25 0.333333
 132          * 0.75 0.222222
 133          * 0.125 0.888889*/
 134
 135         gsl_test_rel (v[0], 3.0/4.0, 1e-3, "reversehalton(2) k=2 v[0]");
 136         gsl_test_rel (v[1], 2.0/9.0, 1e-3, "reversehalton(2) k=2 v[1]");
 137
 138         gsl_qrng_get(g, v);
 139         gsl_test_rel (v[0], 1.0/8.0, 1e-3, "reversehalton(2) k=2 v[0]");
 140         gsl_test_rel (v[1], 8.0/9.0, 1e-3, "reversehalton(2) k=2 v[1]");
 141
 142         gsl_qrng_free(g);
 143
 144
 145         /* test in dimension 3 */
 146         g = gsl_qrng_alloc(gsl_qrng_reversehalton, 3);
 147         gsl_qrng_get(g, v);
 148         gsl_qrng_get(g, v);
 149         gsl_qrng_get(g, v);
 150         gsl_test_rel (v[0], 0.75, 1e-3, "reversehalton(3) k=3 v[0]");
 151         gsl_test_rel (v[1], 2.0/9.0, 1e-3, "reversehalton(3) k=3 v[1]");
 152         gsl_test_rel (v[2], 0.4, 1e-3, "reversehalton(3) k=3 v[2]");
 153
 154         gsl_qrng_get(g, v);
 155         gsl_test_rel (v[0], 0.125, 1e-3, "reversehalton(3) k=3 v[0]");
 156         gsl_test_rel (v[1], 8.0/9.0, 1e-3, "reversehalton(3) k=3 v[1]");
 157         gsl_test_rel (v[2], 0.2, 1e-3, "reversehalton(3) k=3 v[2]");
 158
 159         status = 0;
 160         gsl_qrng_init(g);
 161         gsl_qrng_get(g, v);
 162         gsl_qrng_get(g, v);
 163         gsl_qrng_get(g, v);
 164         gsl_test_rel (v[0], 0.75, 1e-3, "reversehalton(3) reinitialized k=3 v[0]");
 165         gsl_test_rel (v[1], 2.0/9.0, 1e-3, "reversehalton(3) reinitialized k=3 v[1]");
 166         gsl_test_rel (v[2], 0.4, 1e-3, "reversehalton(3) reinitialized k=3 v[2]");
 167
 168         gsl_qrng_get(g, v);
 169         gsl_test_rel (v[0], 0.125, 1e-3, "reversehalton(3) reinitialized k=3 v[0]");
 170         gsl_test_rel (v[1], 8.0/9.0, 1e-3, "reversehalton(3) reinitialized k=3 v[1]");
 171         gsl_test_rel (v[2], 0.2, 1e-3, "reversehalton(3) reinitialized k=3 v[2]");
 172
 173         gsl_qrng_free(g);
 174 }
 175
 176
 177 void test_nied2(void)
 178 {
 179   int status = 0;
 180   double v[3];
 181   /* int i; */
 182
 183   /* test in dimension 2 */
 184   gsl_qrng * g = gsl_qrng_alloc(gsl_qrng_niederreiter_2, 2);
 185   gsl_qrng_get(g, v);
 186   gsl_qrng_get(g, v);
 187   gsl_qrng_get(g, v);
 188   status += ( v[0] != 0.75 || v[1] != 0.25 );
 189   gsl_qrng_get(g, v);
 190   status += ( v[0] != 0.25 || v[1] != 0.75 );
 191   gsl_qrng_get(g, v);
 192   gsl_qrng_get(g, v);
 193   gsl_qrng_get(g, v);
 194   status += ( v[0] != 0.625 || v[1] != 0.125 );
 195   gsl_qrng_free(g);
 196
 197   gsl_test (status, "Niederreiter d=2");
 198
 199   status = 0;
 200
 201   /* test in dimension 3 */
 202   g = gsl_qrng_alloc(gsl_qrng_niederreiter_2, 3);
 203   gsl_qrng_get(g, v);
 204   gsl_qrng_get(g, v);
 205   gsl_qrng_get(g, v);
 206   status += ( v[0] != 0.75 || v[1] != 0.25 || v[2] != 0.3125 );
 207   gsl_qrng_get(g, v);
 208   status += ( v[0] != 0.25 || v[1] != 0.75 || v[2] != 0.5625 );
 209   gsl_qrng_get(g, v);
 210   gsl_qrng_get(g, v);
 211   gsl_qrng_get(g, v);
 212   status += ( v[0] != 0.625 || v[1] != 0.125 || v[2] != 0.6875 );
 213
 214   gsl_test (status, "Niederreiter d=3");
 215
 216   status = 0;
 217
 218   gsl_qrng_init(g);
 219   gsl_qrng_get(g, v);
 220   gsl_qrng_get(g, v);
 221   gsl_qrng_get(g, v);
 222   status += ( v[0] != 0.75 || v[1] != 0.25 || v[2] != 0.3125 );
 223   gsl_qrng_get(g, v);
 224   status += ( v[0] != 0.25 || v[1] != 0.75 || v[2] != 0.5625 );
 225   gsl_qrng_get(g, v);
 226   gsl_qrng_get(g, v);
 227   gsl_qrng_get(g, v);
 228   status += ( v[0] != 0.625 || v[1] != 0.125 || v[2] != 0.6875 );
 229   gsl_qrng_free(g);
 230
 231
 232   gsl_test (status, "Niederreiter d=3 (reinitialized)");
 233 }
 234
 235
 236 int main()
 237 {
 238
 239   gsl_ieee_env_setup ();
 240
 241   test_sobol();
 242         test_halton();
 243         test_reversehalton();
 244   test_nied2();
 245
 246   exit (gsl_test_summary ());
 247 }