1 2 /* 3 * Copyright 2013 Google Inc. 4 * 5 * Use of this source code is governed by a BSD-style license that can be 6 * found in the LICENSE file. 7 */ 8 9 #include "Test.h" 10 #include "SkRandom.h" 11 #include "SkTSort.h" 12 13 static bool anderson_darling_test(double p[32]) { 14 // Min and max Anderson-Darling values allowable for k=32 15 const double kADMin32 = 0.202; // p-value of ~0.1 16 const double kADMax32 = 3.89; // p-value of ~0.99 17 18 // sort p values 19 SkTQSort<double>(p, p + 31); 20 21 // and compute Anderson-Darling statistic to ensure these are uniform 22 double s = 0.0; 23 for(int k = 0; k < 32; k++) { 24 double v = p[k]*(1.0 - p[31-k]); 25 if (v < 1.0e-30) { 26 v = 1.0e-30; 27 } 28 s += (2.0*(k+1)-1.0)*log(v); 29 } 30 double a2 = -32.0 - 0.03125*s; 31 32 return (kADMin32 < a2 && a2 < kADMax32); 33 } 34 35 static bool chi_square_test(int bins[256], int e) { 36 // Min and max chisquare values allowable 37 const double kChiSqMin256 = 206.3179; // probability of chance = 0.99 with k=256 38 const double kChiSqMax256 = 311.5603; // probability of chance = 0.01 with k=256 39 40 // compute chi-square 41 double chi2 = 0.0; 42 for (int j = 0; j < 256; ++j) { 43 double delta = bins[j] - e; 44 chi2 += delta*delta/e; 45 } 46 47 return (kChiSqMin256 < chi2 && chi2 < kChiSqMax256); 48 } 49 50 // Approximation to the normal distribution CDF 51 // From Waissi and Rossin, 1996 52 static double normal_cdf(double z) { 53 double t = ((-0.0004406*z*z* + 0.0418198)*z*z + 0.9)*z; 54 t *= -1.77245385091; // -sqrt(PI) 55 double p = 1.0/(1.0 + exp(t)); 56 57 return p; 58 } 59 60 static void test_random_byte(skiatest::Reporter* reporter, int shift) { 61 int bins[256]; 62 memset(bins, 0, sizeof(int)*256); 63 64 SkMWCRandom rand; 65 for (int i = 0; i < 256*10000; ++i) { 66 bins[(rand.nextU() >> shift) & 0xff]++; 67 } 68 69 REPORTER_ASSERT(reporter, chi_square_test(bins, 10000)); 70 } 71 72 static void test_random_float(skiatest::Reporter* reporter) { 73 int bins[256]; 74 memset(bins, 0, sizeof(int)*256); 75 76 SkMWCRandom rand; 77 for (int i = 0; i < 256*10000; ++i) { 78 float f = rand.nextF(); 79 REPORTER_ASSERT(reporter, 0.0f <= f && f < 1.0f); 80 bins[(int)(f*256.f)]++; 81 } 82 REPORTER_ASSERT(reporter, chi_square_test(bins, 10000)); 83 84 double p[32]; 85 for (int j = 0; j < 32; ++j) { 86 float f = rand.nextF(); 87 REPORTER_ASSERT(reporter, 0.0f <= f && f < 1.0f); 88 p[j] = f; 89 } 90 REPORTER_ASSERT(reporter, anderson_darling_test(p)); 91 } 92 93 // This is a test taken from tuftests by Marsaglia and Tsang. The idea here is that 94 // we are using the random bit generated from a single shift position to generate 95 // "strings" of 16 bits in length, shifting the string and adding a new bit with each 96 // iteration. We track the numbers generated. The ones that we don't generate will 97 // have a normal distribution with mean ~24108 and standard deviation ~127. By 98 // creating a z-score (# of deviations from the mean) for one iteration of this step 99 // we can determine its probability. 100 // 101 // The original test used 26 bit strings, but is somewhat slow. This version uses 16 102 // bits which is less rigorous but much faster to generate. 103 static double test_single_gorilla(skiatest::Reporter* reporter, int shift) { 104 const int kWordWidth = 16; 105 const double kMean = 24108.0; 106 const double kStandardDeviation = 127.0; 107 const int kN = (1 << kWordWidth); 108 const int kNumEntries = kN >> 5; // dividing by 32 109 unsigned int entries[kNumEntries]; 110 111 SkMWCRandom rand; 112 memset(entries, 0, sizeof(unsigned int)*kNumEntries); 113 // pre-seed our string value 114 int value = 0; 115 for (int i = 0; i < kWordWidth-1; ++i) { 116 value <<= 1; 117 unsigned int rnd = rand.nextU(); 118 value |= ((rnd >> shift) & 0x1); 119 } 120 121 // now make some strings and track them 122 for (int i = 0; i < kN; ++i) { 123 value <<= 1; 124 unsigned int rnd = rand.nextU(); 125 value |= ((rnd >> shift) & 0x1); 126 127 int index = value & (kNumEntries-1); 128 SkASSERT(index < kNumEntries); 129 int entry_shift = (value >> (kWordWidth-5)) & 0x1f; 130 entries[index] |= (0x1 << entry_shift); 131 } 132 133 // count entries 134 int total = 0; 135 for (int i = 0; i < kNumEntries; ++i) { 136 unsigned int entry = entries[i]; 137 while (entry) { 138 total += (entry & 0x1); 139 entry >>= 1; 140 } 141 } 142 143 // convert counts to normal distribution z-score 144 double z = ((kN-total)-kMean)/kStandardDeviation; 145 146 // compute probability from normal distibution CDF 147 double p = normal_cdf(z); 148 149 REPORTER_ASSERT(reporter, 0.01 < p && p < 0.99); 150 return p; 151 } 152 153 static void test_gorilla(skiatest::Reporter* reporter) { 154 155 double p[32]; 156 for (int bit_position = 0; bit_position < 32; ++bit_position) { 157 p[bit_position] = test_single_gorilla(reporter, bit_position); 158 } 159 160 REPORTER_ASSERT(reporter, anderson_darling_test(p)); 161 } 162 163 static void test_range(skiatest::Reporter* reporter) { 164 SkMWCRandom rand; 165 166 // just to make sure we don't crash in this case 167 (void) rand.nextRangeU(0, 0xffffffff); 168 169 // check a case to see if it's uniform 170 int bins[256]; 171 memset(bins, 0, sizeof(int)*256); 172 for (int i = 0; i < 256*10000; ++i) { 173 unsigned int u = rand.nextRangeU(17, 17+255); 174 REPORTER_ASSERT(reporter, 17 <= u && u <= 17+255); 175 bins[u - 17]++; 176 } 177 178 REPORTER_ASSERT(reporter, chi_square_test(bins, 10000)); 179 } 180 181 static void TestRandom(skiatest::Reporter* reporter) { 182 // check uniform distributions of each byte in 32-bit word 183 test_random_byte(reporter, 0); 184 test_random_byte(reporter, 8); 185 test_random_byte(reporter, 16); 186 test_random_byte(reporter, 24); 187 188 test_random_float(reporter); 189 190 test_gorilla(reporter); 191 192 test_range(reporter); 193 } 194 195 #include "TestClassDef.h" 196 DEFINE_TESTCLASS("Random", RandomTestClass, TestRandom) 197
