From 50da396d04714a18fa087ebbd1f2316958dbd6bd Mon Sep 17 00:00:00 2001 From: Tyge Løvset Date: Thu, 21 Jan 2021 08:26:54 +0100 Subject: Reverted namings: crand to crandom, and copt to coption. --- README.md | 2 +- benchmarks/cdeq_benchmark.cpp | 2 +- benchmarks/cmap_benchmark.cpp | 2 +- benchmarks/cmap_benchmark2.cpp | 2 +- benchmarks/cpque_benchmark.cpp | 2 +- benchmarks/crand_benchmark.cpp | 2 +- benchmarks/crand_benchmark2.cpp | 2 +- benchmarks/vector_vs_deque.cpp | 2 +- docs/copt_api.md | 85 ------------------- docs/coption_api.md | 85 +++++++++++++++++++ docs/cpque_api.md | 2 +- docs/crand_api.md | 165 ------------------------------------- docs/crandom_api.md | 165 +++++++++++++++++++++++++++++++++++++ examples/birthday.c | 2 +- examples/cbits_prime.c | 39 +++++++++ examples/crandom_ex.c | 51 ++++++++++++ examples/csmap_ex.c | 2 +- examples/ex_gauss1.c | 2 +- examples/ex_gauss2.c | 2 +- examples/list.c | 2 +- examples/phonebook.c | 72 ++++++++++++++++ examples/priority.c | 2 +- examples/queue.c | 2 +- examples/random.c | 2 +- stc/clist.h | 2 +- stc/copt.h | 176 ---------------------------------------- stc/coption.h | 176 ++++++++++++++++++++++++++++++++++++++++ stc/cpque.h | 2 +- stc/cqueue.h | 2 +- stc/crand.h | 157 ----------------------------------- stc/crandom.h | 157 +++++++++++++++++++++++++++++++++++ 31 files changed, 765 insertions(+), 603 deletions(-) delete mode 100644 docs/copt_api.md create mode 100644 docs/coption_api.md delete mode 100644 docs/crand_api.md create mode 100644 docs/crandom_api.md create mode 100644 examples/cbits_prime.c create mode 100644 examples/crandom_ex.c create mode 100644 examples/phonebook.c delete mode 100644 stc/copt.h create mode 100644 stc/coption.h delete mode 100644 stc/crand.h create mode 100644 stc/crandom.h diff --git a/README.md b/README.md index 11b98ff9..8cbf2b82 100644 --- a/README.md +++ b/README.md @@ -26,7 +26,7 @@ This is a compact headers-only library with the all of the "standard" data conta Others: - [***ccommon*** - General definitions](docs/ccommon_api.md) - [***copt*** - Implements ***copt_get()***, similar to posix **getopt_long()**](docs/copt_api.md) -- [***crand*** - A very efficent modern **pseudo-random number generator**](docs/crand_api.md) +- [***crandom*** - A very efficent modern **pseudo-random number generator**](docs/crandom_api.md) The usage of the containers is similar to the c++ standard containers in STL, so it should be easy if you are familiar with them. All containers mentioned above are generic, except for **cstr** and **cbits**. No casting is used, and containers are therefore diff --git a/benchmarks/cdeq_benchmark.cpp b/benchmarks/cdeq_benchmark.cpp index f44a40a6..ea1c69b3 100644 --- a/benchmarks/cdeq_benchmark.cpp +++ b/benchmarks/cdeq_benchmark.cpp @@ -3,7 +3,7 @@ #include #include #include -#include +#include enum {N = 1000000000, M = 12345, P = 5000, R = 2000}; using_cdeq(i, int); diff --git a/benchmarks/cmap_benchmark.cpp b/benchmarks/cmap_benchmark.cpp index 47f60b00..967c8abe 100644 --- a/benchmarks/cmap_benchmark.cpp +++ b/benchmarks/cmap_benchmark.cpp @@ -1,6 +1,6 @@ #include #include -#include +#include #include #include #include "others/khash.h" diff --git a/benchmarks/cmap_benchmark2.cpp b/benchmarks/cmap_benchmark2.cpp index 33f5c675..54d6e7d9 100644 --- a/benchmarks/cmap_benchmark2.cpp +++ b/benchmarks/cmap_benchmark2.cpp @@ -1,4 +1,4 @@ -#include +#include #include #include #include diff --git a/benchmarks/cpque_benchmark.cpp b/benchmarks/cpque_benchmark.cpp index 396d763d..d3c00803 100644 --- a/benchmarks/cpque_benchmark.cpp +++ b/benchmarks/cpque_benchmark.cpp @@ -1,6 +1,6 @@ #include #include -#include +#include #include #include diff --git a/benchmarks/crand_benchmark.cpp b/benchmarks/crand_benchmark.cpp index 3d42c62b..136dcf30 100644 --- a/benchmarks/crand_benchmark.cpp +++ b/benchmarks/crand_benchmark.cpp @@ -1,7 +1,7 @@ #include #include #include -#include +#include static inline uint64_t rotl64(const uint64_t x, const int k) { return (x << k) | (x >> (64 - k)); } diff --git a/benchmarks/crand_benchmark2.cpp b/benchmarks/crand_benchmark2.cpp index 2692e760..6a517a0b 100644 --- a/benchmarks/crand_benchmark2.cpp +++ b/benchmarks/crand_benchmark2.cpp @@ -1,7 +1,7 @@ #include #include #include -#include "stc/crand.h" +#include "stc/crandom.h" #include "others/pcg_random.hpp" static struct stc32_state { stc64_t rng; uint64_t spare; unsigned n; } stc32_global = diff --git a/benchmarks/vector_vs_deque.cpp b/benchmarks/vector_vs_deque.cpp index c9e612f2..a488b7ab 100644 --- a/benchmarks/vector_vs_deque.cpp +++ b/benchmarks/vector_vs_deque.cpp @@ -8,7 +8,7 @@ #include #include #include -#include +#include typedef struct {const char* first; int second;} Si; using_cvec(si, Si, c_no_compare); diff --git a/docs/copt_api.md b/docs/copt_api.md deleted file mode 100644 index df49f66d..00000000 --- a/docs/copt_api.md +++ /dev/null @@ -1,85 +0,0 @@ -# STC Module [copt](../stc/copt.h): Command line argument parsing - -This describes the API of the *copt_get()* function for command line argument parsing. -See [getopt_long](https://www.freebsd.org/cgi/man.cgi?getopt_long(3)) for a similar posix function. - -## Types - -```c -enum { - copt_no_argument = 0, - copt_required_argument = 1, - copt_optional_argument = 2 -}; -typedef struct { - int ind; /* equivalent to optind */ - int opt; /* equivalent to optopt */ - const char *arg; /* equivalent to optarg */ - const char *faulty; /* points to the faulty option, if any */ - int longindex; /* index of long option; or -1 if short */ - ... -} copt_t; - -typedef struct { - const char *name; - int has_arg; - int val; -} copt_long_t; - -const copt_t copt_inits; -``` - -## Methods - -```c -copt_t copt_init(void); -int copt_get(copt_t *opt, int argc, char *argv[], - const char *shortopts, const copt_long_t *longopts); -``` - -## Example - -```c -#include -#include "stc/copt.h" - -int main(int argc, char *argv[]) { - static copt_long_t long_options[] = { - {"verbose", copt_no_argument, 'V'}, - {"help", copt_no_argument, 'H'}, - {"add", copt_no_argument, 'a'}, - {"append", copt_no_argument, 'b'}, - {"delete", copt_required_argument, 'd'}, - {"create", copt_required_argument, 'c'}, - {"file", copt_required_argument, 'f'}, - {NULL} - }; - copt_t opt = copt_inits; - int c; - while ((c = copt_get(&opt, argc, argv, ":if:lr", long_options)) != -1) { - switch (c) { - case 'V': case 'H': - case 'a': case 'b': - case 'd': case 'c': - case 'i': case 'l': - case 'r': - printf("option: %c\n", c); - break; - case 'f': - printf("filename: %s\n", opt.arg); - break; - case ':': - printf("option %s needs a value\n", opt.faulty); - break; - case '?': - printf("unknown option: %s\n", opt.faulty); - break; - } - } - - for (; opt.ind < argc; ++opt.ind) { - printf("extra arguments: %s\n", argv[opt.ind]); - } - return 0; -} -``` diff --git a/docs/coption_api.md b/docs/coption_api.md new file mode 100644 index 00000000..66038b38 --- /dev/null +++ b/docs/coption_api.md @@ -0,0 +1,85 @@ +# STC Module [copt](../stc/coption.h): Command line argument parsing + +This describes the API of the *coption_get()* function for command line argument parsing. +See [getopt_long](https://www.freebsd.org/cgi/man.cgi?getopt_long(3)) for a similar posix function. + +## Types + +```c +enum { + copt_no_argument = 0, + copt_required_argument = 1, + copt_optional_argument = 2 +}; +typedef struct { + int ind; /* equivalent to posix optind */ + int opt; /* equivalent to posix optopt */ + const char *arg; /* equivalent to posix optarg */ + const char *faulty; /* points to the faulty option, if any */ + int longindex; /* index of long option; or -1 if short */ + ... +} coption; + +typedef struct { + const char *name; + int has_arg; + int val; +} coption_long; + +const coption coption_inits; +``` + +## Methods + +```c +coption coption_init(void); +int coption_get(coption *opt, int argc, char *argv[], + const char *shortopts, const coption_long *longopts); +``` + +## Example + +```c +#include +#include "stc/coption.h" + +int main(int argc, char *argv[]) { + static coption_long long_options[] = { + {"verbose", copt_no_argument, 'V'}, + {"help", copt_no_argument, 'H'}, + {"add", copt_no_argument, 'a'}, + {"append", copt_no_argument, 'b'}, + {"delete", copt_required_argument, 'd'}, + {"create", copt_required_argument, 'c'}, + {"file", copt_required_argument, 'f'}, + {NULL} + }; + coption opt = coption_init(); + int c; + while ((c = coption_get(&opt, argc, argv, ":if:lr", long_options)) != -1) { + switch (c) { + case 'V': case 'H': + case 'a': case 'b': + case 'd': case 'c': + case 'i': case 'l': + case 'r': + printf("option: %c\n", c); + break; + case 'f': + printf("filename: %s\n", opt.arg); + break; + case ':': + printf("option %s needs a value\n", opt.faulty); + break; + case '?': + printf("unknown option: %s\n", opt.faulty); + break; + } + } + + for (; opt.ind < argc; ++opt.ind) { + printf("extra arguments: %s\n", argv[opt.ind]); + } + return 0; +} +``` diff --git a/docs/cpque_api.md b/docs/cpque_api.md index 59203a63..080d580e 100644 --- a/docs/cpque_api.md +++ b/docs/cpque_api.md @@ -60,7 +60,7 @@ cpque_X_value_t cpque_X_value_clone(cpque_X_value_t val); ```c #include #include "stc/cpque.h" -#include "stc/crand.h" +#include "stc/crandom.h" using_cvec(i, int64_t); using_cpque(i, cvec_i, >); // adaptor type, '>' = min-heap diff --git a/docs/crand_api.md b/docs/crand_api.md deleted file mode 100644 index 995c5fa0..00000000 --- a/docs/crand_api.md +++ /dev/null @@ -1,165 +0,0 @@ -# STC Module [crand](../stc/crand.h): Pseudo Random Number Generators -![Random](pics/random.jpg) - -This describes the API of module **crand**. It contains **stc64**, a *64-bit PRNG*, and can generate -bounded uniform and normal distributed random numbers. See [random](https://en.cppreference.com/w/cpp/header/random) -for similar c++ functionality. - -**stc64** is an extremely fast PRNG by Tyge Løvset, suited for parallel usage. It features a -Weyl-sequence as part of the state. It is faster than *sfc64*, *wyhash64*, *pcg64*, and *xoshiro256\*\** -on common platforms. It does not require fast multiplication or 128-bit integer operations. It has a -256 bit state, but updates only 192 bit per generated number. - -There is no *jump function*, but by incrementing the Weyl-increment by 2, it starts a new -unique 2^64 *minimum* length period. Note that for each Weyl-increment (state[3]), the period -length is about 2^126 with a high probability. For a single thread, a minimum period of 2^127 -is generated when the Weyl-increment is incremented by 2 every 2^64 output. - -**stc64** passes *PractRand*, tested up to 8TB output, Vigna's Hamming weight test, and simple -correlation tests, i.e. *n* interleaved streams with only one-bit differences in initial state. - -See the PRNG shootout by Vigna: http://prng.di.unimi.it and the debate between the authors of -xoshiro and pcg (Vigna/O'Neill) PRNGs: https://www.pcg-random.org/posts/on-vignas-pcg-critique.html - -## Types - -| Name | Type definition | Used to represent... | -|:-------------------|:------------------------------------------|:-----------------------------| -| `stc64_t` | `struct {uint64_t state[4];}` | The PRNG engine type | -| `stc64_uniform_t` | `struct {int64_t lower; uint64_t range;}` | Integer uniform distribution | -| `stc64_uniformf_t` | `struct {double lower, range;}` | Real number uniform distr. | -| `stc64_normalf_t` | `struct {double mean, stddev;}` | Normal distribution type | - -## Header file - -All cstr definitions and prototypes may be included in your C source file by including a single header file. -```c -#include "stc/crand.h" -``` - -## Methods - -```c - void stc64_srandom(uint64_t seed); - uint64_t stc64_random(void); - - 1) stc64_t stc64_init(uint64_t seed); - 2) stc64_t stc64_with_seq(uint64_t seed, uint64_t seq); - - 3) uint64_t stc64_rand(stc64_t* rng); - 4) double stc64_randf(stc64_t* rng); - - 5) stc64_uniform_t stc64_uniform_init(int64_t low, int64_t high); - 6) int64_t stc64_uniform(stc64_t* rng, stc64_uniform_t* dist); - 7) stc64_uniformf_t stc64_uniformf_init(double low, double high); - 8) double stc64_uniformf(stc64_t* rng, stc64_uniformf_t* dist); - - 9) stc64_normalf_t stc64_normalf_init(double mean, double stddev); -10) double stc64_normalf(stc64_t* rng, stc64_normalf_t* dist); -``` -`1-2)` PRNG 64-bit engine initializers. `3)` Integer generator, range \[0, 2^64). -`4)` Double RNG with range \[0, 1). `5-6)` Uniform integer RNG with range \[*low*, *high*]. -`7-8)` Uniform double RNG with range \[*low*, *high*). `9-10)` Normal-distributed double -RNG, around 68% of the values fall within the range [*mean* - *stddev*, *mean* + *stddev*]. - -## Example -```c -#include -#include -#include -#include "stc/crand.h" -#include "stc/cstr.h" -#include "stc/cmap.h" -#include "stc/cvec.h" - -// Declare unordered map: int -> int with typetag 'i'. -using_cmap(i, int, size_t); - -// Comparison of map keys. -static int compare(cmap_i_value_t *a, cmap_i_value_t *b) { - return c_default_compare(&a->first, &b->first); -} -// Declare vector of map entries, with comparison function. -using_cvec(e, cmap_i_value_t, compare); - - -int main() -{ - enum {N = 10000000}; - const double Mean = -12.0, StdDev = 6.0, Scale = 74; - - printf("Demo of gaussian / normal distribution of %d random samples\n", N); - - // Setup random engine with normal distribution. - uint64_t seed = time(NULL); - stc64_t rng = stc64_init(seed); - stc64_normalf_t dist = stc64_normalf_init(Mean, StdDev); - - // Create histogram map - cmap_i mhist = cmap_i_init(); - for (size_t i = 0; i < N; ++i) { - int index = (int) round( stc64_normalf(&rng, &dist) ); - cmap_i_emplace(&mhist, index, 0).first->second += 1; - } - - // Transfer map to vec and sort it by map entry keys. - cvec_e vhist = cvec_e_init(); - c_foreach (i, cmap_i, mhist) - cvec_e_push_back(&vhist, *i.ref); - cvec_e_sort(&vhist); - - // Print the gaussian bar chart - cstr_t bar = cstr_init(); - c_foreach (i, cvec_e, vhist) { - size_t n = (size_t) (i.ref->second * StdDev * Scale * 2.5 / N); - if (n > 0) { - cstr_resize(&bar, n, '*'); - printf("%4d %s\n", i.ref->first, bar.str); - } - } - - // Cleanup - cstr_del(&bar); - cmap_i_del(&mhist); - cvec_e_del(&vhist); -} -``` -Output: -``` -Demo of gaussian / normal distribution of 10000000 random samples - -29 * - -28 ** - -27 *** - -26 **** - -25 ******* - -24 ********* - -23 ************* - -22 ****************** - -21 *********************** - -20 ****************************** - -19 ************************************* - -18 ******************************************** - -17 **************************************************** - -16 *********************************************************** - -15 ***************************************************************** - -14 ********************************************************************* - -13 ************************************************************************ - -12 ************************************************************************* - -11 ************************************************************************ - -10 ********************************************************************* - -9 ***************************************************************** - -8 *********************************************************** - -7 **************************************************** - -6 ******************************************** - -5 ************************************* - -4 ****************************** - -3 *********************** - -2 ****************** - -1 ************* - 0 ********* - 1 ******* - 2 **** - 3 *** - 4 ** - 5 * -``` diff --git a/docs/crandom_api.md b/docs/crandom_api.md new file mode 100644 index 00000000..c8ad34ff --- /dev/null +++ b/docs/crandom_api.md @@ -0,0 +1,165 @@ +# STC Module [crandom](../stc/crandom.h): Pseudo Random Number Generators +![Random](pics/random.jpg) + +This describes the API of module **crandom**. It contains **stc64**, a *64-bit PRNG*, and can generate +bounded uniform and normal distributed random numbers. See [random](https://en.cppreference.com/w/cpp/header/random) +for similar c++ functionality. + +**stc64** is an extremely fast PRNG by Tyge Løvset, suited for parallel usage. It features a +Weyl-sequence as part of the state. It is faster than *sfc64*, *wyhash64*, *pcg64*, and *xoshiro256\*\** +on common platforms. It does not require fast multiplication or 128-bit integer operations. It has a +256 bit state, but updates only 192 bit per generated number. + +There is no *jump function*, but by incrementing the Weyl-increment by 2, it starts a new +unique 2^64 *minimum* length period. Note that for each Weyl-increment (state[3]), the period +length is about 2^126 with a high probability. For a single thread, a minimum period of 2^127 +is generated when the Weyl-increment is incremented by 2 every 2^64 output. + +**stc64** passes *PractRand*, tested up to 8TB output, Vigna's Hamming weight test, and simple +correlation tests, i.e. *n* interleaved streams with only one-bit differences in initial state. + +See the PRNG shootout by Vigna: http://prng.di.unimi.it and the debate between the authors of +xoshiro and pcg (Vigna/O'Neill) PRNGs: https://www.pcg-random.org/posts/on-vignas-pcg-critique.html + +## Types + +| Name | Type definition | Used to represent... | +|:-------------------|:------------------------------------------|:-----------------------------| +| `stc64_t` | `struct {uint64_t state[4];}` | The PRNG engine type | +| `stc64_uniform_t` | `struct {int64_t lower; uint64_t range;}` | Integer uniform distribution | +| `stc64_uniformf_t` | `struct {double lower, range;}` | Real number uniform distr. | +| `stc64_normalf_t` | `struct {double mean, stddev;}` | Normal distribution type | + +## Header file + +All cstr definitions and prototypes may be included in your C source file by including a single header file. +```c +#include "stc/crandom.h" +``` + +## Methods + +```c + void stc64_srandom(uint64_t seed); + uint64_t stc64_random(void); + + 1) stc64_t stc64_init(uint64_t seed); + 2) stc64_t stc64_with_seq(uint64_t seed, uint64_t seq); + + 3) uint64_t stc64_rand(stc64_t* rng); + 4) double stc64_randf(stc64_t* rng); + + 5) stc64_uniform_t stc64_uniform_init(int64_t low, int64_t high); + 6) int64_t stc64_uniform(stc64_t* rng, stc64_uniform_t* dist); + 7) stc64_uniformf_t stc64_uniformf_init(double low, double high); + 8) double stc64_uniformf(stc64_t* rng, stc64_uniformf_t* dist); + + 9) stc64_normalf_t stc64_normalf_init(double mean, double stddev); +10) double stc64_normalf(stc64_t* rng, stc64_normalf_t* dist); +``` +`1-2)` PRNG 64-bit engine initializers. `3)` Integer generator, range \[0, 2^64). +`4)` Double RNG with range \[0, 1). `5-6)` Uniform integer RNG with range \[*low*, *high*]. +`7-8)` Uniform double RNG with range \[*low*, *high*). `9-10)` Normal-distributed double +RNG, around 68% of the values fall within the range [*mean* - *stddev*, *mean* + *stddev*]. + +## Example +```c +#include +#include +#include +#include "stc/crandom.h" +#include "stc/cstr.h" +#include "stc/cmap.h" +#include "stc/cvec.h" + +// Declare unordered map: int -> int with typetag 'i'. +using_cmap(i, int, size_t); + +// Comparison of map keys. +static int compare(cmap_i_value_t *a, cmap_i_value_t *b) { + return c_default_compare(&a->first, &b->first); +} +// Declare vector of map entries, with comparison function. +using_cvec(e, cmap_i_value_t, compare); + + +int main() +{ + enum {N = 10000000}; + const double Mean = -12.0, StdDev = 6.0, Scale = 74; + + printf("Demo of gaussian / normal distribution of %d random samples\n", N); + + // Setup random engine with normal distribution. + uint64_t seed = time(NULL); + stc64_t rng = stc64_init(seed); + stc64_normalf_t dist = stc64_normalf_init(Mean, StdDev); + + // Create histogram map + cmap_i mhist = cmap_i_init(); + for (size_t i = 0; i < N; ++i) { + int index = (int) round( stc64_normalf(&rng, &dist) ); + cmap_i_emplace(&mhist, index, 0).first->second += 1; + } + + // Transfer map to vec and sort it by map entry keys. + cvec_e vhist = cvec_e_init(); + c_foreach (i, cmap_i, mhist) + cvec_e_push_back(&vhist, *i.ref); + cvec_e_sort(&vhist); + + // Print the gaussian bar chart + cstr_t bar = cstr_init(); + c_foreach (i, cvec_e, vhist) { + size_t n = (size_t) (i.ref->second * StdDev * Scale * 2.5 / N); + if (n > 0) { + cstr_resize(&bar, n, '*'); + printf("%4d %s\n", i.ref->first, bar.str); + } + } + + // Cleanup + cstr_del(&bar); + cmap_i_del(&mhist); + cvec_e_del(&vhist); +} +``` +Output: +``` +Demo of gaussian / normal distribution of 10000000 random samples + -29 * + -28 ** + -27 *** + -26 **** + -25 ******* + -24 ********* + -23 ************* + -22 ****************** + -21 *********************** + -20 ****************************** + -19 ************************************* + -18 ******************************************** + -17 **************************************************** + -16 *********************************************************** + -15 ***************************************************************** + -14 ********************************************************************* + -13 ************************************************************************ + -12 ************************************************************************* + -11 ************************************************************************ + -10 ********************************************************************* + -9 ***************************************************************** + -8 *********************************************************** + -7 **************************************************** + -6 ******************************************** + -5 ************************************* + -4 ****************************** + -3 *********************** + -2 ****************** + -1 ************* + 0 ********* + 1 ******* + 2 **** + 3 *** + 4 ** + 5 * +``` diff --git a/examples/birthday.c b/examples/birthday.c index dc94fa53..45786fa9 100644 --- a/examples/birthday.c +++ b/examples/birthday.c @@ -2,7 +2,7 @@ #include #include -#include +#include #include using_cmap(ic, uint64_t, uint8_t); diff --git a/examples/cbits_prime.c b/examples/cbits_prime.c new file mode 100644 index 00000000..5c5d3969 --- /dev/null +++ b/examples/cbits_prime.c @@ -0,0 +1,39 @@ +#include +#include + +static inline cbits sieveOfEratosthenes(size_t n) +{ + cbits primes = cbits_with_size(n + 1, true); + cbits_reset(&primes, 0); + cbits_reset(&primes, 1); + + c_forrange (i, size_t, 2, n+1) { + // If primes[i] is not changed, then it is a prime + if (cbits_test(primes, i) && i*i <= n) { + c_forrange (j, size_t, i*i, n+1, i) { + cbits_reset(&primes, j); + } + } + } + return primes; +} + + +int main(void) +{ + int n = 100000000; + printf("computing prime numbers up to %u\n", n); + + cbits primes = sieveOfEratosthenes(n); + puts("done"); + + size_t np = cbits_count(primes); + printf("number of primes: %zu\n", np); + + printf("2 "); + c_forrange (i, int, 3, 1001, 2) { + if (cbits_test(primes, i)) printf("%d ", i); + } + puts(""); + cbits_del(&primes); +} \ No newline at end of file diff --git a/examples/crandom_ex.c b/examples/crandom_ex.c new file mode 100644 index 00000000..c923debd --- /dev/null +++ b/examples/crandom_ex.c @@ -0,0 +1,51 @@ +#include +#include +#include +#include +#include + +int main() +{ + enum {R = 30}; + const size_t N = 1000000000; + uint64_t seed = 1234; // time(NULL); + stc64_t rng = stc64_init(seed); + + uint64_t sum = 0; + + stc64_normalf_t dist2 = stc64_normalf_init(R / 2.0, R / 6.0); + size_t N2 = 10000000; + int hist[R] = {0}; + sum = 0; + c_forrange (N2) { + int n = round((stc64_normalf(&rng, &dist2) + 0.5)); + sum += n; + if (n >= 0 && n < R) ++hist[n]; + } + cstr_t bar = cstr_inits; + c_forrange (i, int, R) { + cstr_resize(&bar, hist[i] * 25ull * R / N2, '*'); + printf("%3d %s\n", i, bar.str); + } + + clock_t diff, before; + + sum = 0; + before = clock(); + c_forrange (N) { + sum += stc64_rand(&rng); + } + diff = clock() - before; + printf("random : %f secs, %zu %f\n", (float) diff / CLOCKS_PER_SEC, N, (double) sum / N); + + stc64_uniform_t dist1 = stc64_uniform_init(0, 1000); + sum = 0; + before = clock(); + c_forrange (N) { + sum += stc64_uniform(&rng, &dist1); + } + diff = clock() - before; + printf("uniform: %f secs, %zu %f\n", (float) diff / CLOCKS_PER_SEC, N, (double) sum / N); + + cstr_del(&bar); +} \ No newline at end of file diff --git a/examples/csmap_ex.c b/examples/csmap_ex.c index aa78028b..ea8fd81c 100644 --- a/examples/csmap_ex.c +++ b/examples/csmap_ex.c @@ -1,6 +1,6 @@ #include #include -#include +#include #include using_csmap(i, int, size_t); diff --git a/examples/ex_gauss1.c b/examples/ex_gauss1.c index 4139210c..85d818d5 100644 --- a/examples/ex_gauss1.c +++ b/examples/ex_gauss1.c @@ -1,7 +1,7 @@ #include #include #include -#include "stc/crand.h" +#include "stc/crandom.h" #include "stc/cstr.h" #include "stc/cmap.h" #include "stc/cvec.h" diff --git a/examples/ex_gauss2.c b/examples/ex_gauss2.c index cff0cef3..434a7e5e 100644 --- a/examples/ex_gauss2.c +++ b/examples/ex_gauss2.c @@ -1,7 +1,7 @@ #include #include #include -#include "stc/crand.h" +#include "stc/crandom.h" #include "stc/cstr.h" #include "stc/csmap.h" diff --git a/examples/list.c b/examples/list.c index 9d6e8d89..394d0b38 100644 --- a/examples/list.c +++ b/examples/list.c @@ -1,7 +1,7 @@ #include #include #include -#include +#include using_clist(fx, double); int main() { diff --git a/examples/phonebook.c b/examples/phonebook.c new file mode 100644 index 00000000..91d01072 --- /dev/null +++ b/examples/phonebook.c @@ -0,0 +1,72 @@ +// The MIT License (MIT) +// Copyright (c) 2018 Maksim Andrianov +// +// Permission is hereby granted, free of charge, to any person obtaining a copy +// of this software and associated documentation files (the "Software"), to +// deal in the Software without restriction, including without limitation the +// rights to use, copy, modify, merge, publish, distribute, sublicense, and/or +// sell copies of the Software, and to permit persons to whom the Software is +// furnished to do so, subject to the following conditions: +// +// The above copyright notice and this permission notice shall be included in +// all copies or substantial portions of the Software. +// +// THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR +// IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, +// FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE +// AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER +// LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING +// FROM, OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS +// IN THE SOFTWARE. + +// Program to emulates the phone book. + +#include +#include +#include + +using_cmap_str(); + +void print_phone_book(cmap_str phone_book) +{ + c_foreach (i, cmap_str, phone_book) + printf("%s\t- %s\n", i.ref->first.str, i.ref->second.str); +} + +int main(int argc, char **argv) +{ + bool erased; + cmap_str phone_book = cmap_inits; + c_push_items(&phone_book, cmap_str, { + {"Lilia Friedman", "(892) 670-4739"}, + {"Tariq Beltran", "(489) 600-7575"}, + {"Laiba Juarez", "(303) 885-5692"}, + {"Elliott Mooney", "(945) 616-4482"}, + }); + + printf("Phone book:\n"); + print_phone_book(phone_book); + + c_try_emplace(&phone_book, cmap_str, "Zak Byers", cstr_from("(551) 396-1880")); + c_try_emplace(&phone_book, cmap_str, "Zak Byers", cstr_from("(551) 396-1990")); + + printf("\nPhone book after adding Zak Byers:\n"); + print_phone_book(phone_book); + + if (cmap_str_find(&phone_book, "Tariq Beltran") != NULL) + printf("\nTariq Beltran is in phone book\n"); + + erased = cmap_str_erase(&phone_book, "Tariq Beltran"); + erased = cmap_str_erase(&phone_book, "Elliott Mooney"); + + printf("\nPhone book after erasing Tariq and Elliott:\n"); + print_phone_book(phone_book); + + cmap_str_insert_or_assign(&phone_book, "Zak Byers", "(555) 396-188"); + + printf("\nPhone book after update phone of Zak Byers:\n"); + print_phone_book(phone_book); + + cmap_str_del(&phone_book); + puts("done"); +} \ No newline at end of file diff --git a/examples/priority.c b/examples/priority.c index 4eceb29f..548cb06d 100644 --- a/examples/priority.c +++ b/examples/priority.c @@ -4,7 +4,7 @@ #include #include #include -#include +#include using_cvec(i, int64_t); using_cpque(i, cvec_i, >); // min-heap (increasing values) diff --git a/examples/queue.c b/examples/queue.c index 8c0a8777..23640a39 100644 --- a/examples/queue.c +++ b/examples/queue.c @@ -1,4 +1,4 @@ -#include +#include #include #include diff --git a/examples/random.c b/examples/random.c index 7325e023..c923debd 100644 --- a/examples/random.c +++ b/examples/random.c @@ -1,7 +1,7 @@ #include #include #include -#include +#include #include int main() diff --git a/stc/clist.h b/stc/clist.h index 89bbc75c..1673d71e 100644 --- a/stc/clist.h +++ b/stc/clist.h @@ -31,7 +31,7 @@ #include #include - #include + #include using_clist(ix, int64_t); int main() { diff --git a/stc/copt.h b/stc/copt.h deleted file mode 100644 index 5fc8fe4d..00000000 --- a/stc/copt.h +++ /dev/null @@ -1,176 +0,0 @@ -/* MIT License - * - * Copyright (c) 2021 Tyge Løvset, NORCE, www.norceresearch.no - * - * Permission is hereby granted, free of charge, to any person obtaining a copy - * of this software and associated documentation files (the "Software"), to deal - * in the Software without restriction, including without limitation the rights - * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell - * copies of the Software, and to permit persons to whom the Software is - * furnished to do so, subject to the following conditions: - * - * The above copyright notice and this permission notice shall be included in all - * copies or substantial portions of the Software. - * - * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR - * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, - * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE - * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER - * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, - * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE - * SOFTWARE. - */ -#ifndef COPT__H__ -#define COPT__H__ - -/* -// Inspired by https://attractivechaos.wordpress.com/2018/08/31/a-survey-of-argument-parsing-libraries-in-c-c -// Fixed major bugs with option arguments (both long and short). -// Added arg->faulty output field, and has a more consistent API. -// -// copt_get() is similar to GNU's getopt_long(). Each call parses one option and -// returns the option name. opt->arg points to the option argument if present. -// The function returns -1 when all command-line arguments are parsed. In this case, -// opt->ind is the index of the first non-option argument. -#include -#include - -int main(int argc, char *argv[]) -{ - copt_long_t longopts[] = { - {"foo", copt_no_argument, 'f'}, - {"bar", copt_required_argument, 'b'}, - {"opt", copt_optional_argument, 'o'}, - {NULL} - }; - const char* optstr = "xy:z::123"; - printf("program -x -y ARG -z [ARG] -1 -2 -3 --foo --bar ARG --opt [ARG] [ARGUMENTS]\n"); - int c; - copt_t opt = copt_init(); - while ((c = copt_get(&opt, argc, argv, optstr, longopts)) != -1) { - switch (c) { - case '?': printf("error: unknown option: %s\n", opt.faulty); break; - case ':': printf("error: missing argument for %s\n", opt.faulty); break; - default: printf("option: %c [%s]\n", c, opt.arg ? opt.arg : ""); break; - } - } - printf("\nNon-option arguments:"); - for (int i = opt.ind; i < argc; ++i) - printf(" %s", argv[i]); - putchar('\n'); - return 0; -} -*/ -#include -#include - -enum { - copt_no_argument = 0, - copt_required_argument = 1, - copt_optional_argument = 2 -}; -typedef struct { - int ind; /* equivalent to optind */ - int opt; /* equivalent to optopt */ - const char *arg; /* equivalent to optarg */ - const char *faulty; /* points to the faulty option */ - int longindex; /* idx of long option; or -1 if short */ - int _i, _pos, _nargs; - char _faulty[4]; -} copt_t; - -typedef struct { - const char *name; - int has_arg; - int val; -} copt_long_t; - -static const copt_t copt_inits = {1, 0, NULL, NULL, -1, 1, 0, 0, {'-', '?', '\0'}}; -static inline copt_t copt_init(void) { return copt_inits; } - -static void _copt_permute(char *argv[], int j, int n) { /* move argv[j] over n elements to the left */ - int k; - char *p = argv[j]; - for (k = 0; k < n; ++k) - argv[j - k] = argv[j - k - 1]; - argv[j - k] = p; -} - -/* @param opt output; must be initialized to copt_init() on first call - * @return ASCII val for a short option; longopt.val for a long option; - * -1 if argv[] is fully processed; '?' for an unknown option or - * an ambiguous long option; ':' if an option argument is missing - */ -static int copt_get(copt_t *opt, int argc, char *argv[], - const char *shortopts, const copt_long_t *longopts) { - int optc = -1, i0, j, posixly_correct = (shortopts[0] == '+'); - if (!posixly_correct) { - while (opt->_i < argc && (argv[opt->_i][0] != '-' || argv[opt->_i][1] == '\0')) - ++opt->_i, ++opt->_nargs; - } - opt->arg = 0, opt->longindex = -1, i0 = opt->_i; - if (opt->_i >= argc || argv[opt->_i][0] != '-' || argv[opt->_i][1] == '\0') { - opt->ind = opt->_i - opt->_nargs; - return -1; - } - if (argv[opt->_i][0] == '-' && argv[opt->_i][1] == '-') { /* "--" or a long option */ - if (argv[opt->_i][2] == '\0') { /* a bare "--" */ - _copt_permute(argv, opt->_i, opt->_nargs); - ++opt->_i, opt->ind = opt->_i - opt->_nargs; - return -1; - } - opt->opt = 0, optc = '?', opt->_pos = -1; - if (longopts) { /* parse long options */ - int k, n_exact = 0, n_partial = 0; - const copt_long_t *o = 0, *o_exact = 0, *o_partial = 0; - for (j = 2; argv[opt->_i][j] != '\0' && argv[opt->_i][j] != '='; ++j) {} /* find the end of the option name */ - for (k = 0; longopts[k].name != 0; ++k) - if (strncmp(&argv[opt->_i][2], longopts[k].name, j - 2) == 0) { - if (longopts[k].name[j - 2] == 0) ++n_exact, o_exact = &longopts[k]; - else ++n_partial, o_partial = &longopts[k]; - } - opt->faulty = argv[opt->_i]; - if (n_exact > 1 || (n_exact == 0 && n_partial > 1)) return '?'; - o = n_exact == 1? o_exact : n_partial == 1? o_partial : 0; - if (o) { - opt->opt = optc = o->val, opt->longindex = o - longopts; - if (o->has_arg != copt_no_argument) { - if (argv[opt->_i][j] == '=') - opt->arg = &argv[opt->_i][j + 1]; - else if (argv[opt->_i][j] == '\0' && opt->_i < argc - 1 && (o->has_arg == copt_required_argument || - argv[opt->_i + 1][0] != '-')) - opt->arg = argv[++opt->_i]; - else if (o->has_arg == copt_required_argument) - optc = ':'; /* missing option argument */ - } - } - } - } else { /* a short option */ - const char *p; - if (opt->_pos == 0) opt->_pos = 1; - optc = opt->opt = argv[opt->_i][opt->_pos++]; - opt->_faulty[1] = optc, opt->faulty = opt->_faulty; - p = strchr((char *) shortopts, optc); - if (p == 0) { - optc = '?'; /* unknown option */ - } else if (p[1] == ':') { - if (argv[opt->_i][opt->_pos] != '\0') - opt->arg = &argv[opt->_i][opt->_pos]; - else if (opt->_i < argc - 1 && (p[2] != ':' || argv[opt->_i + 1][0] != '-')) - opt->arg = argv[++opt->_i]; - else if (p[2] != ':') - optc = ':'; - opt->_pos = -1; - } - } - if (opt->_pos < 0 || argv[opt->_i][opt->_pos] == 0) { - ++opt->_i, opt->_pos = 0; - if (opt->_nargs > 0) /* permute */ - for (j = i0; j < opt->_i; ++j) - _copt_permute(argv, j, opt->_nargs); - } - opt->ind = opt->_i - opt->_nargs; - return optc; -} - -#endif diff --git a/stc/coption.h b/stc/coption.h new file mode 100644 index 00000000..27ec880b --- /dev/null +++ b/stc/coption.h @@ -0,0 +1,176 @@ +/* MIT License + * + * Copyright (c) 2021 Tyge Løvset, NORCE, www.norceresearch.no + * + * Permission is hereby granted, free of charge, to any person obtaining a copy + * of this software and associated documentation files (the "Software"), to deal + * in the Software without restriction, including without limitation the rights + * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell + * copies of the Software, and to permit persons to whom the Software is + * furnished to do so, subject to the following conditions: + * + * The above copyright notice and this permission notice shall be included in all + * copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR + * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE + * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER + * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, + * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ +#ifndef COPTION__H__ +#define COPTION__H__ + +/* +// Inspired by https://attractivechaos.wordpress.com/2018/08/31/a-survey-of-argument-parsing-libraries-in-c-c +// Fixed major bugs with option arguments (both long and short). +// Added arg->faulty output field, and has a more consistent API. +// +// coption_get() is similar to GNU's getopt_long(). Each call parses one option and +// returns the option name. opt->arg points to the option argument if present. +// The function returns -1 when all command-line arguments are parsed. In this case, +// opt->ind is the index of the first non-option argument. +#include +#include + +int main(int argc, char *argv[]) +{ + coption_long longopts[] = { + {"foo", copt_no_argument, 'f'}, + {"bar", copt_required_argument, 'b'}, + {"opt", copt_optional_argument, 'o'}, + {NULL} + }; + const char* optstr = "xy:z::123"; + printf("program -x -y ARG -z [ARG] -1 -2 -3 --foo --bar ARG --opt [ARG] [ARGUMENTS]\n"); + int c; + coption opt = coption_init(); + while ((c = coption_get(&opt, argc, argv, optstr, longopts)) != -1) { + switch (c) { + case '?': printf("error: unknown option: %s\n", opt.faulty); break; + case ':': printf("error: missing argument for %s\n", opt.faulty); break; + default: printf("option: %c [%s]\n", c, opt.arg ? opt.arg : ""); break; + } + } + printf("\nNon-option arguments:"); + for (int i = opt.ind; i < argc; ++i) + printf(" %s", argv[i]); + putchar('\n'); + return 0; +} +*/ +#include +#include + +enum { + copt_no_argument = 0, + copt_required_argument = 1, + copt_optional_argument = 2 +}; +typedef struct { + int ind; /* equivalent to optind */ + int opt; /* equivalent to optopt */ + const char *arg; /* equivalent to optarg */ + const char *faulty; /* points to the faulty option */ + int longindex; /* idx of long option; or -1 if short */ + int _i, _pos, _nargs; + char _faulty[4]; +} coption; + +typedef struct { + const char *name; + int has_arg; + int val; +} coption_long; + +static const coption coption_inits = {1, 0, NULL, NULL, -1, 1, 0, 0, {'-', '?', '\0'}}; +static inline coption coption_init(void) { return coption_inits; } + +static void _copt_permute(char *argv[], int j, int n) { /* move argv[j] over n elements to the left */ + int k; + char *p = argv[j]; + for (k = 0; k < n; ++k) + argv[j - k] = argv[j - k - 1]; + argv[j - k] = p; +} + +/* @param opt output; must be initialized to coption_init() on first call + * @return ASCII val for a short option; longopt.val for a long option; + * -1 if argv[] is fully processed; '?' for an unknown option or + * an ambiguous long option; ':' if an option argument is missing + */ +static int coption_get(coption *opt, int argc, char *argv[], + const char *shortopts, const coption_long *longopts) { + int optc = -1, i0, j, posixly_correct = (shortopts[0] == '+'); + if (!posixly_correct) { + while (opt->_i < argc && (argv[opt->_i][0] != '-' || argv[opt->_i][1] == '\0')) + ++opt->_i, ++opt->_nargs; + } + opt->arg = 0, opt->longindex = -1, i0 = opt->_i; + if (opt->_i >= argc || argv[opt->_i][0] != '-' || argv[opt->_i][1] == '\0') { + opt->ind = opt->_i - opt->_nargs; + return -1; + } + if (argv[opt->_i][0] == '-' && argv[opt->_i][1] == '-') { /* "--" or a long option */ + if (argv[opt->_i][2] == '\0') { /* a bare "--" */ + _copt_permute(argv, opt->_i, opt->_nargs); + ++opt->_i, opt->ind = opt->_i - opt->_nargs; + return -1; + } + opt->opt = 0, optc = '?', opt->_pos = -1; + if (longopts) { /* parse long options */ + int k, n_exact = 0, n_partial = 0; + const coption_long *o = 0, *o_exact = 0, *o_partial = 0; + for (j = 2; argv[opt->_i][j] != '\0' && argv[opt->_i][j] != '='; ++j) {} /* find the end of the option name */ + for (k = 0; longopts[k].name != 0; ++k) + if (strncmp(&argv[opt->_i][2], longopts[k].name, j - 2) == 0) { + if (longopts[k].name[j - 2] == 0) ++n_exact, o_exact = &longopts[k]; + else ++n_partial, o_partial = &longopts[k]; + } + opt->faulty = argv[opt->_i]; + if (n_exact > 1 || (n_exact == 0 && n_partial > 1)) return '?'; + o = n_exact == 1? o_exact : n_partial == 1? o_partial : 0; + if (o) { + opt->opt = optc = o->val, opt->longindex = o - longopts; + if (o->has_arg != copt_no_argument) { + if (argv[opt->_i][j] == '=') + opt->arg = &argv[opt->_i][j + 1]; + else if (argv[opt->_i][j] == '\0' && opt->_i < argc - 1 && (o->has_arg == copt_required_argument || + argv[opt->_i + 1][0] != '-')) + opt->arg = argv[++opt->_i]; + else if (o->has_arg == copt_required_argument) + optc = ':'; /* missing option argument */ + } + } + } + } else { /* a short option */ + const char *p; + if (opt->_pos == 0) opt->_pos = 1; + optc = opt->opt = argv[opt->_i][opt->_pos++]; + opt->_faulty[1] = optc, opt->faulty = opt->_faulty; + p = strchr((char *) shortopts, optc); + if (p == 0) { + optc = '?'; /* unknown option */ + } else if (p[1] == ':') { + if (argv[opt->_i][opt->_pos] != '\0') + opt->arg = &argv[opt->_i][opt->_pos]; + else if (opt->_i < argc - 1 && (p[2] != ':' || argv[opt->_i + 1][0] != '-')) + opt->arg = argv[++opt->_i]; + else if (p[2] != ':') + optc = ':'; + opt->_pos = -1; + } + } + if (opt->_pos < 0 || argv[opt->_i][opt->_pos] == 0) { + ++opt->_i, opt->_pos = 0; + if (opt->_nargs > 0) /* permute */ + for (j = i0; j < opt->_i; ++j) + _copt_permute(argv, j, opt->_nargs); + } + opt->ind = opt->_i - opt->_nargs; + return optc; +} + +#endif diff --git a/stc/cpque.h b/stc/cpque.h index 41113b7a..54ca04a7 100644 --- a/stc/cpque.h +++ b/stc/cpque.h @@ -23,7 +23,7 @@ /* Priority-Queue adapter (implemented as heap), default uses cvec. - #include + #include #include using_cvec(f, float); using_cpque(f, cvec_f, >); // min-heap (increasing values) diff --git a/stc/cqueue.h b/stc/cqueue.h index fe202698..2c49b066 100644 --- a/stc/cqueue.h +++ b/stc/cqueue.h @@ -25,7 +25,7 @@ /* Queue adapter, default uses clist. - #include + #include #include using_cdeq(i, int); using_cqueue(i, cdeq_i); diff --git a/stc/crand.h b/stc/crand.h deleted file mode 100644 index 727c61ac..00000000 --- a/stc/crand.h +++ /dev/null @@ -1,157 +0,0 @@ -/* MIT License - * - * Copyright (c) 2021 Tyge Løvset, NORCE, www.norceresearch.no - * - * Permission is hereby granted, free of charge, to any person obtaining a copy - * of this software and associated documentation files (the "Software"), to deal - * in the Software without restriction, including without limitation the rights - * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell - * copies of the Software, and to permit persons to whom the Software is - * furnished to do so, subject to the following conditions: - * - * The above copyright notice and this permission notice shall be included in all - * copies or substantial portions of the Software. - * - * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR - * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, - * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE - * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER - * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, - * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE - * SOFTWARE. - */ -#ifndef CRANDOM__H__ -#define CRANDOM__H__ - -/* -// crand: Pseudo-random number generator -#include "stc/crand.h" -int main() { - uint64_t seed = 123456789; - stc64_t rng = stc64_init(seed); - stc64_uniform_t dist1 = stc64_uniform_init(1, 6); - stc64_uniformf_t dist2 = stc64_uniformf_init(1.0, 10.0); - stc64_normalf_t dist3 = stc64_normalf_init(1.0, 10.0); - - uint64_t i = stc64_rand(&rng); - int64_t iu = stc64_uniform(&rng, &dist1); - double xu = stc64_uniformf(&rng, &dist2); - double xn = stc64_normalf(&rng, &dist3); -} -*/ -#include "ccommon.h" -#include -#include - -typedef struct {uint64_t state[4];} stc64_t; -typedef struct {int64_t lower; uint64_t range, threshold;} stc64_uniform_t; -typedef struct {double lower, range;} stc64_uniformf_t; -typedef struct {double mean, stddev, next; unsigned has_next;} stc64_normalf_t; - - -/* Stc64: random number generator, range [0, 2^64). PRNG copyright Tyge Løvset, NORCE Research, 2020 */ -STC_API stc64_t stc64_init(uint64_t seed); -STC_API stc64_t stc64_with_seq(uint64_t seed, uint64_t seq); - -STC_INLINE uint64_t stc64_rand(stc64_t* rng) { - uint64_t *s = rng->state; - const uint64_t b = s[1], result = s[0] ^ (s[2] += s[3]|1); - s[0] = (b + (b << 3)) ^ (b >> 11); - s[1] = ((b << 24) | (b >> (64 - 24))) + result; - return result; -} - -/* Global random() */ -static stc64_t stc64_global = {{0x26aa069ea2fb1a4d, 0x70c72c95cd592d04, 0x504f333d3aa0b359, 0x6a09e667a754166b}}; -STC_INLINE void stc64_srandom(uint64_t seed) { stc64_global = stc64_init(seed); } -STC_INLINE uint64_t stc64_random(void) { return stc64_rand(&stc64_global); } - -/* Float64 random number in range [low, high). */ -STC_INLINE double stc64_randf(stc64_t* rng) { - union {uint64_t i; double f;} u = {0x3FF0000000000000ull | (stc64_rand(rng) >> 12)}; - return u.f - 1.0; -} - -/* Int64 uniform distributed RNG, range [low, high]. */ -STC_API stc64_uniform_t stc64_uniform_init(int64_t low, int64_t high); - -/* Float64 uniform distributed RNG, range [low, high). */ -STC_INLINE stc64_uniformf_t stc64_uniformf_init(double low, double high) { - stc64_uniformf_t dist = {low, high - low}; return dist; -} -STC_INLINE double stc64_uniformf(stc64_t* rng, stc64_uniformf_t* dist) { - return stc64_randf(rng)*dist->range + dist->lower; -} - -#if defined(__SIZEOF_INT128__) - #define cmul128(a, b, lo, hi) \ - do { __uint128_t _z = (__uint128_t)(a) * (b); \ - *(lo) = (uint64_t)_z, *(hi) = _z >> 64; } while(0) -#elif defined(_MSC_VER) && defined(_WIN64) - #include - #define cmul128(a, b, lo, hi) (*(lo) = _umul128(a, b, hi), (void)0) -#elif defined(__x86_64__) - #define cmul128(a, b, lo, hi) \ - asm("mulq %[rhs]" : "=a" (*(lo)), "=d" (*(hi)) \ - : [lhs] "0" (a), [rhs] "rm" (b)) -#endif - -/* Unbiased bounded uniform distribution. */ -STC_INLINE int64_t stc64_uniform(stc64_t* rng, stc64_uniform_t* d) { - uint64_t lo, hi; - do { cmul128(stc64_rand(rng), d->range, &lo, &hi); } while (lo < d->threshold); - return d->lower + hi; -} - -/* Normal distributed RNG, Float64. */ -STC_INLINE stc64_normalf_t stc64_normalf_init(double mean, double stddev) { - stc64_normalf_t dist = {mean, stddev, 0.0, 0}; return dist; -} -STC_API double stc64_normalf(stc64_t* rng, stc64_normalf_t* dist); - - -#if !defined(STC_HEADER) || defined(STC_IMPLEMENTATION) - -/* PRNG crand: by Tyge Løvset, NORCE Research, 2020. - * Extremely fast PRNG suited for parallel usage with Weyl-sequence parameter. - * Faster than sfc64, wyhash64, and xoshiro256** on most platforms. - * 256bit state, updates only 192bit per rng. - * 2^63 unique threads with a minimum 2^64 period lengths each. - * 2^127 minimum period length for single thread (double loop). - * Passes PractRand tested up to 8TB output, Vigna's Hamming weight test, - * and simple correlation tests, i.e. interleaved streams with one-bit diff state. - */ - -STC_DEF stc64_t stc64_init(uint64_t seed) { - return stc64_with_seq(seed, seed + 0x3504f333d3aa0b34); -} -STC_DEF stc64_t stc64_with_seq(uint64_t seed, uint64_t seq) { - stc64_t rng = {{seed, seed, seed, (seq << 1u) | 1u}}; - for (int i = 0; i < 8; ++i) stc64_rand(&rng); - return rng; -} - -/* Very fast unbiased uniform int RNG with bounds [low, high] */ -STC_DEF stc64_uniform_t stc64_uniform_init(int64_t low, int64_t high) { - stc64_uniform_t dist = {low, (uint64_t) (high - low + 1)}; - dist.threshold = (uint64_t)(-dist.range) % dist.range; - return dist; -} - -/* Marsaglia polar method for gaussian/normal distribution. */ -STC_DEF double stc64_normalf(stc64_t* rng, stc64_normalf_t* dist) { - double u1, u2, s, m; - if (dist->has_next++ & 1) - return dist->next * dist->stddev + dist->mean; - do { - u1 = 2.0 * stc64_randf(rng) - 1.0; - u2 = 2.0 * stc64_randf(rng) - 1.0; - s = u1*u1 + u2*u2; - } while (s >= 1.0 || s == 0.0); - m = sqrt(-2.0 * log(s) / s); - dist->next = u2 * m; - return (u1 * m) * dist->stddev + dist->mean; -} - -#endif -#endif diff --git a/stc/crandom.h b/stc/crandom.h new file mode 100644 index 00000000..c92a7de6 --- /dev/null +++ b/stc/crandom.h @@ -0,0 +1,157 @@ +/* MIT License + * + * Copyright (c) 2021 Tyge Løvset, NORCE, www.norceresearch.no + * + * Permission is hereby granted, free of charge, to any person obtaining a copy + * of this software and associated documentation files (the "Software"), to deal + * in the Software without restriction, including without limitation the rights + * to use, copy, modify, merge, publish, distribute, sublicense, and/or sell + * copies of the Software, and to permit persons to whom the Software is + * furnished to do so, subject to the following conditions: + * + * The above copyright notice and this permission notice shall be included in all + * copies or substantial portions of the Software. + * + * THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR + * IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY, + * FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE + * AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER + * LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM, + * OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN THE + * SOFTWARE. + */ +#ifndef CRANDOM__H__ +#define CRANDOM__H__ + +/* +// crandom: Pseudo-random number generator +#include "stc/crandom.h" +int main() { + uint64_t seed = 123456789; + stc64_t rng = stc64_init(seed); + stc64_uniform_t dist1 = stc64_uniform_init(1, 6); + stc64_uniformf_t dist2 = stc64_uniformf_init(1.0, 10.0); + stc64_normalf_t dist3 = stc64_normalf_init(1.0, 10.0); + + uint64_t i = stc64_rand(&rng); + int64_t iu = stc64_uniform(&rng, &dist1); + double xu = stc64_uniformf(&rng, &dist2); + double xn = stc64_normalf(&rng, &dist3); +} +*/ +#include "ccommon.h" +#include +#include + +typedef struct {uint64_t state[4];} stc64_t; +typedef struct {int64_t lower; uint64_t range, threshold;} stc64_uniform_t; +typedef struct {double lower, range;} stc64_uniformf_t; +typedef struct {double mean, stddev, next; unsigned has_next;} stc64_normalf_t; + + +/* Stc64: random number generator, range [0, 2^64). PRNG copyright Tyge Løvset, NORCE Research, 2020 */ +STC_API stc64_t stc64_init(uint64_t seed); +STC_API stc64_t stc64_with_seq(uint64_t seed, uint64_t seq); + +STC_INLINE uint64_t stc64_rand(stc64_t* rng) { + uint64_t *s = rng->state; + const uint64_t b = s[1], result = s[0] ^ (s[2] += s[3]|1); + s[0] = (b + (b << 3)) ^ (b >> 11); + s[1] = ((b << 24) | (b >> (64 - 24))) + result; + return result; +} + +/* Global random() */ +static stc64_t stc64_global = {{0x26aa069ea2fb1a4d, 0x70c72c95cd592d04, 0x504f333d3aa0b359, 0x6a09e667a754166b}}; +STC_INLINE void stc64_srandom(uint64_t seed) { stc64_global = stc64_init(seed); } +STC_INLINE uint64_t stc64_random(void) { return stc64_rand(&stc64_global); } + +/* Float64 random number in range [low, high). */ +STC_INLINE double stc64_randf(stc64_t* rng) { + union {uint64_t i; double f;} u = {0x3FF0000000000000ull | (stc64_rand(rng) >> 12)}; + return u.f - 1.0; +} + +/* Int64 uniform distributed RNG, range [low, high]. */ +STC_API stc64_uniform_t stc64_uniform_init(int64_t low, int64_t high); + +/* Float64 uniform distributed RNG, range [low, high). */ +STC_INLINE stc64_uniformf_t stc64_uniformf_init(double low, double high) { + stc64_uniformf_t dist = {low, high - low}; return dist; +} +STC_INLINE double stc64_uniformf(stc64_t* rng, stc64_uniformf_t* dist) { + return stc64_randf(rng)*dist->range + dist->lower; +} + +#if defined(__SIZEOF_INT128__) + #define cmul128(a, b, lo, hi) \ + do { __uint128_t _z = (__uint128_t)(a) * (b); \ + *(lo) = (uint64_t)_z, *(hi) = _z >> 64; } while(0) +#elif defined(_MSC_VER) && defined(_WIN64) + #include + #define cmul128(a, b, lo, hi) (*(lo) = _umul128(a, b, hi), (void)0) +#elif defined(__x86_64__) + #define cmul128(a, b, lo, hi) \ + asm("mulq %[rhs]" : "=a" (*(lo)), "=d" (*(hi)) \ + : [lhs] "0" (a), [rhs] "rm" (b)) +#endif + +/* Unbiased bounded uniform distribution. */ +STC_INLINE int64_t stc64_uniform(stc64_t* rng, stc64_uniform_t* d) { + uint64_t lo, hi; + do { cmul128(stc64_rand(rng), d->range, &lo, &hi); } while (lo < d->threshold); + return d->lower + hi; +} + +/* Normal distributed RNG, Float64. */ +STC_INLINE stc64_normalf_t stc64_normalf_init(double mean, double stddev) { + stc64_normalf_t dist = {mean, stddev, 0.0, 0}; return dist; +} +STC_API double stc64_normalf(stc64_t* rng, stc64_normalf_t* dist); + + +#if !defined(STC_HEADER) || defined(STC_IMPLEMENTATION) + +/* PRNG crandom: by Tyge Løvset, NORCE Research, 2020. + * Extremely fast PRNG suited for parallel usage with Weyl-sequence parameter. + * Faster than sfc64, wyhash64, and xoshiro256** on most platforms. + * 256bit state, updates only 192bit per rng. + * 2^63 unique threads with a minimum 2^64 period lengths each. + * 2^127 minimum period length for single thread (double loop). + * Passes PractRand tested up to 8TB output, Vigna's Hamming weight test, + * and simple correlation tests, i.e. interleaved streams with one-bit diff state. + */ + +STC_DEF stc64_t stc64_init(uint64_t seed) { + return stc64_with_seq(seed, seed + 0x3504f333d3aa0b34); +} +STC_DEF stc64_t stc64_with_seq(uint64_t seed, uint64_t seq) { + stc64_t rng = {{seed, seed, seed, (seq << 1u) | 1u}}; + for (int i = 0; i < 8; ++i) stc64_rand(&rng); + return rng; +} + +/* Very fast unbiased uniform int RNG with bounds [low, high] */ +STC_DEF stc64_uniform_t stc64_uniform_init(int64_t low, int64_t high) { + stc64_uniform_t dist = {low, (uint64_t) (high - low + 1)}; + dist.threshold = (uint64_t)(-dist.range) % dist.range; + return dist; +} + +/* Marsaglia polar method for gaussian/normal distribution. */ +STC_DEF double stc64_normalf(stc64_t* rng, stc64_normalf_t* dist) { + double u1, u2, s, m; + if (dist->has_next++ & 1) + return dist->next * dist->stddev + dist->mean; + do { + u1 = 2.0 * stc64_randf(rng) - 1.0; + u2 = 2.0 * stc64_randf(rng) - 1.0; + s = u1*u1 + u2*u2; + } while (s >= 1.0 || s == 0.0); + m = sqrt(-2.0 * log(s) / s); + dist->next = u2 * m; + return (u1 * m) * dist->stddev + dist->mean; +} + +#endif +#endif -- cgit v1.2.3