personal/squiggle.c
1
0
Fork 1
Code Issues Pull Requests Packages Projects Releases Wiki Activity

formatting pass.

Browse Source
This commit is contained in:
NunoSempere 2023-07-23 16:30:42 +02:00
parent c8fd237bbf
commit d531d5571f
1 changed files with 10 additions and 10 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

20
squiggle.c
View File

@ -35,11 +35,11 @@ uint64_t xorshift64(uint64_t* seed)
// https://en.wikipedia.org/wiki/Xorshift // https://en.wikipedia.org/wiki/Xorshift
// Also some drama: <https://www.pcg-random.org/posts/on-vignas-pcg-critique.html>, <https://prng.di.unimi.it/> // Also some drama: <https://www.pcg-random.org/posts/on-vignas-pcg-critique.html>, <https://prng.di.unimi.it/>
uint64_t x = *seed; uint64_t x = *seed;
x ^= x << 13; x ^= x << 13;
x ^= x >> 7; x ^= x >> 7;
x ^= x << 17; x ^= x << 17;
return *seed = x; return *seed = x;
} }
// Distribution & sampling functions // Distribution & sampling functions
@ -70,9 +70,9 @@ double sample_normal(double mean, double sigma, uint64_t* seed)
return (mean + sigma * sample_unit_normal(seed)); return (mean + sigma * sample_unit_normal(seed));
} }
double sample_lognormal(double logmean, double logsigma, uint64_t* seed) double sample_lognormal(double logmean, double logstd, uint64_t* seed)
{ {
return expf(sample_normal(logmean, logsigma, seed)); return exp(sample_normal(logmean, logstd, seed));
} }
double sample_to(double low, double high, uint64_t* seed) double sample_to(double low, double high, uint64_t* seed)
@ -84,8 +84,8 @@ double sample_to(double low, double high, uint64_t* seed)
double loglow = logf(low); double loglow = logf(low);
double loghigh = logf(high); double loghigh = logf(high);
double logmean = (loglow + loghigh) / 2; double logmean = (loglow + loghigh) / 2;
double logsigma = (loghigh - loglow) / (2.0 * NORMAL95CONFIDENCE); double logstd = (loghigh - loglow) / (2.0 * NORMAL95CONFIDENCE);
return sample_lognormal(logmean, logsigma, seed); return sample_lognormal(logmean, logstd, seed);
} }
double sample_gamma(double alpha, uint64_t* seed) double sample_gamma(double alpha, uint64_t* seed)
@ -105,7 +105,7 @@ double sample_gamma(double alpha, uint64_t* seed)
v = 1.0 + c * x; v = 1.0 + c * x;
} while (v <= 0.0); } while (v <= 0.0);
v = v * v * v; v = v * v * v;
u = sample_unit_uniform(seed); u = sample_unit_uniform(seed);
if (u < 1.0 - 0.0331 * (x * x * x * x)) { // Condition 1 if (u < 1.0 - 0.0331 * (x * x * x * x)) { // Condition 1
// the 0.0331 doesn't inspire much confidence // the 0.0331 doesn't inspire much confidence