diff --git a/examples/more/00_example_template/example b/examples/more/00_example_template/example index 5944427..820272d 100755 Binary files a/examples/more/00_example_template/example and b/examples/more/00_example_template/example differ diff --git a/examples/more/01_sample_from_cdf/example b/examples/more/01_sample_from_cdf/example index 7e0ceb8..093bbd9 100755 Binary files a/examples/more/01_sample_from_cdf/example and b/examples/more/01_sample_from_cdf/example differ diff --git a/examples/more/02_sample_from_cdf_beta/example b/examples/more/02_sample_from_cdf_beta/example index 7e44064..7ed3878 100755 Binary files a/examples/more/02_sample_from_cdf_beta/example and b/examples/more/02_sample_from_cdf_beta/example differ diff --git a/examples/more/03_ci_beta/example b/examples/more/03_ci_beta/example index e1192aa..c014e3d 100755 Binary files a/examples/more/03_ci_beta/example and b/examples/more/03_ci_beta/example differ diff --git a/examples/more/04_nuclear_war/example b/examples/more/04_nuclear_war/example index 5e8be10..2a62491 100755 Binary files a/examples/more/04_nuclear_war/example and b/examples/more/04_nuclear_war/example differ diff --git a/examples/more/05_burn_10kg_fat/example b/examples/more/05_burn_10kg_fat/example index 9123b3c..89e5e34 100755 Binary files a/examples/more/05_burn_10kg_fat/example and b/examples/more/05_burn_10kg_fat/example differ diff --git a/examples/more/06_nuclear_recovery/example b/examples/more/06_nuclear_recovery/example index 40ef062..7844d96 100755 Binary files a/examples/more/06_nuclear_recovery/example and b/examples/more/06_nuclear_recovery/example differ diff --git a/examples/more/07_algebra/example b/examples/more/07_algebra/example index 0154034..c40ca97 100755 Binary files a/examples/more/07_algebra/example and b/examples/more/07_algebra/example differ diff --git a/examples/more/08_algebra_and_conversion/example b/examples/more/08_algebra_and_conversion/example index f304f6d..ca7ab7b 100755 Binary files a/examples/more/08_algebra_and_conversion/example and b/examples/more/08_algebra_and_conversion/example differ diff --git a/examples/more/09_ergonomic_algebra/example b/examples/more/09_ergonomic_algebra/example index 4b9a5c2..f72590a 100755 Binary files a/examples/more/09_ergonomic_algebra/example and b/examples/more/09_ergonomic_algebra/example differ diff --git a/examples/more/10_twitter_thread_example/example b/examples/more/10_twitter_thread_example/example index d8d5747..422c236 100755 Binary files a/examples/more/10_twitter_thread_example/example and b/examples/more/10_twitter_thread_example/example differ diff --git a/examples/more/11_billion_lognormals_paralell/example b/examples/more/11_billion_lognormals_paralell/example index 64d7522..7ca1dcb 100755 Binary files a/examples/more/11_billion_lognormals_paralell/example and b/examples/more/11_billion_lognormals_paralell/example differ diff --git a/examples/more/11_billion_lognormals_paralell/example.c b/examples/more/11_billion_lognormals_paralell/example.c index c9d49d7..1988157 100644 --- a/examples/more/11_billion_lognormals_paralell/example.c +++ b/examples/more/11_billion_lognormals_paralell/example.c @@ -21,7 +21,7 @@ int main() sampler_parallel(sampler, results, n_threads, n_samples); double avg = array_sum(results, n_samples) / n_samples; - printf("Average of 1B lognormal(0,10): %f", avg); + printf("Average of 1B lognormal(0,10): %f\n", avg); free(results); diff --git a/examples/more/12_time_to_botec_parallel/example b/examples/more/12_time_to_botec_parallel/example index eb823ae..cf2afa4 100755 Binary files a/examples/more/12_time_to_botec_parallel/example and b/examples/more/12_time_to_botec_parallel/example differ diff --git a/examples/more/13_parallelize_min/example b/examples/more/13_parallelize_min/example index 4daa170..e8ee480 100755 Binary files a/examples/more/13_parallelize_min/example and b/examples/more/13_parallelize_min/example differ diff --git a/examples/more/14_check_confidence_interval/example b/examples/more/14_check_confidence_interval/example index 7d0b50d..f6fc748 100755 Binary files a/examples/more/14_check_confidence_interval/example and b/examples/more/14_check_confidence_interval/example differ diff --git a/squiggle_more.c b/squiggle_more.c index 1b48ca3..872f55b 100644 --- a/squiggle_more.c +++ b/squiggle_more.c @@ -10,14 +10,40 @@ /* Parallel sampler */ void sampler_parallel(double (*sampler)(uint64_t* seed), double* results, int n_threads, int n_samples) { - if ((n_samples % n_threads) != 0) { - fprintf(stderr, "Number of samples isn't divisible by number of threads, aborting\n"); - exit(1); - } + + // Division terminology: + // a = b * quotient + reminder + // a = (a/b)*b + (a%b) + // dividend: a + // divisor: b + // quotient = a / b + // reminder = a % b + // "divisor's multiple" := (a/b)*b + + // now, we have n_samples and n_threads + // to make our life easy, each thread will have a number of samples of: a/b (quotient) + // and we'll compute the remainder of samples separately + // to possibly do by Jorge: improve so that the remainder is included in the threads + + int quotient = n_samples / n_threads; + int remainder = n_samples % n_threads; + int divisor_multiple = quotient * n_threads; + uint64_t** seeds = malloc(n_threads * sizeof(uint64_t*)); + // printf("UINT64_MAX: %lu\n", UINT64_MAX); + srand(1); for (uint64_t i = 0; i < n_threads; i++) { seeds[i] = malloc(sizeof(uint64_t)); - *seeds[i] = i + 1; // xorshift can't start with 0 + // Constraints: + // - xorshift can't start with 0 + // - the seeds should be reasonably separated and not correlated + *seeds[i] = (uint64_t)rand() * (UINT64_MAX / RAND_MAX); + // printf("#%ld: %lu\n",i, *seeds[i]); + + // Other initializations tried: + // *seeds[i] = 1 + i; + // *seeds[i] = (i + 0.5)*(UINT64_MAX/n_threads); + // *seeds[i] = (i + 0.5)*(UINT64_MAX/n_threads) + constant * i; } int i; @@ -25,14 +51,18 @@ void sampler_parallel(double (*sampler)(uint64_t* seed), double* results, int n_ { #pragma omp for for (i = 0; i < n_threads; i++) { - int lower_bound = i * (n_samples / n_threads); - int upper_bound = ((i + 1) * (n_samples / n_threads)) - 1; + int lower_bound_inclusive = i * quotient; + int upper_bound_not_inclusive = ((i + 1) * quotient); // note the < in the for loop below, // printf("Lower bound: %d, upper bound: %d\n", lower_bound, upper_bound); - for (int j = lower_bound; j < upper_bound; j++) { + for (int j = lower_bound_inclusive; j < upper_bound_not_inclusive; j++) { results[j] = sampler(seeds[i]); } } } + for (int j = divisor_multiple; j < n_samples; j++) { + results[j] = sampler(seeds[0]); + // we can just reuse a seed, this isn't problematic because we are not doing multithreading + } for (uint64_t i = 0; i < n_threads; i++) { free(seeds[i]);