From c9ca66bd402470d380fdbca166f62f21251ccfca Mon Sep 17 00:00:00 2001 From: Mason Kramer Date: Tue, 5 Jul 2022 14:03:40 -0400 Subject: [PATCH 1/4] Deterministic Rayon monte carlo --- Cargo.toml | 1 + examples/rayon-monte-carlo.rs | 67 +++++++++++++++++++++++++++++++++++ 2 files changed, 68 insertions(+) create mode 100755 examples/rayon-monte-carlo.rs diff --git a/Cargo.toml b/Cargo.toml index 98ba373c68f..b013057427b 100644 --- a/Cargo.toml +++ b/Cargo.toml @@ -83,3 +83,4 @@ libc = { version = "0.2.22", optional = true, default-features = false } rand_pcg = { path = "rand_pcg", version = "0.3.0" } # Only to test serde1 bincode = "1.2.1" +rayon = "1.5.3" diff --git a/examples/rayon-monte-carlo.rs b/examples/rayon-monte-carlo.rs new file mode 100755 index 00000000000..a79c44cd90f --- /dev/null +++ b/examples/rayon-monte-carlo.rs @@ -0,0 +1,67 @@ +// Copyright 2018 Developers of the Rand project. +// Copyright 2013-2018 The Rust Project Developers. +// +// Licensed under the Apache License, Version 2.0 or the MIT license +// , at your +// option. This file may not be copied, modified, or distributed +// except according to those terms. + +//! # Monte Carlo estimation of π with a chosen seed and rayon for parallelism +//! +//! Imagine that we have a square with sides of length 2 and a unit circle +//! (radius = 1), both centered at the origin. The areas are: +//! +//! ```text +//! area of circle = πr² = π * r * r = π +//! area of square = 2² = 4 +//! ``` +//! +//! The circle is entirely within the square, so if we sample many points +//! randomly from the square, roughly π / 4 of them should be inside the circle. +//! +//! We can use the above fact to estimate the value of π: pick many points in +//! the square at random, calculate the fraction that fall within the circle, +//! and multiply this fraction by 4. + +#![cfg(all(feature = "std", feature = "std_rng"))] + +use rand::distributions::{Distribution, Uniform}; +use rand_chacha::rand_core::SeedableRng; +use rand_chacha::ChaCha8Rng; +use rayon::prelude::*; + +static SEED: u64 = 0; + +fn main() { + let range = Uniform::new(-1.0f64, 1.0); + + let total = 1_000_000; + let cha_cha = ChaCha8Rng::seed_from_u64(SEED); + + let in_circle: usize = (0..total) + .into_par_iter() + .map_init( + || cha_cha.clone(), + |rng, i| { + // ChaCha supports indexing into its stream. We need this because due to work-stealing, + // Rayon does not guarantee that the same work items will be run in the same order + // between runs of the program. We can only guarantee determinism by setting the stream + // according to the work number. + rng.set_word_pos(i); + let a = range.sample(rng); + let b = range.sample(rng); + if a * a + b * b <= 1.0 { + 1 + } else { + 0 + } + }, + ) + .reduce(|| 0, |a, b| a + b); + // prints something close to 3.14159... + println!( + "π is approximately {}", + 4. * (in_circle as f64) / (total as f64) + ); +} From c4d3821a3c6682e2ed990fd5e0abd9d8ce3485d6 Mon Sep 17 00:00:00 2001 From: Mason Kramer Date: Wed, 6 Jul 2022 10:15:20 -0400 Subject: [PATCH 2/4] Update deterministic mt with a batching example --- examples/rayon-monte-carlo.rs | 46 +++++++++++++++++------------------ 1 file changed, 23 insertions(+), 23 deletions(-) diff --git a/examples/rayon-monte-carlo.rs b/examples/rayon-monte-carlo.rs index a79c44cd90f..859e9ad4f5d 100755 --- a/examples/rayon-monte-carlo.rs +++ b/examples/rayon-monte-carlo.rs @@ -27,41 +27,41 @@ #![cfg(all(feature = "std", feature = "std_rng"))] use rand::distributions::{Distribution, Uniform}; -use rand_chacha::rand_core::SeedableRng; -use rand_chacha::ChaCha8Rng; +use rand_chacha::{rand_core::SeedableRng, ChaCha8Rng}; use rayon::prelude::*; static SEED: u64 = 0; +static BATCH_SIZE: u64 = 10_000; +static BATCHES: u64 = 1000; fn main() { let range = Uniform::new(-1.0f64, 1.0); - let total = 1_000_000; - let cha_cha = ChaCha8Rng::seed_from_u64(SEED); - - let in_circle: usize = (0..total) + // We have to manually batch the work so that we can get Rng construction + // out of the hot loop. Due to work-stealing, ParallelIterator::map_init + // doesn't guarantee the same # and order of work items per run. + let in_circle = (0..BATCHES) .into_par_iter() - .map_init( - || cha_cha.clone(), - |rng, i| { - // ChaCha supports indexing into its stream. We need this because due to work-stealing, - // Rayon does not guarantee that the same work items will be run in the same order - // between runs of the program. We can only guarantee determinism by setting the stream - // according to the work number. - rng.set_word_pos(i); - let a = range.sample(rng); - let b = range.sample(rng); + .map(|i| { + let mut rng = ChaCha8Rng::seed_from_u64(SEED); + // using set_stream() here force every batch to use a different part of the rng stream + rng.set_stream(i); + let mut count = 0; + for _ in 0..BATCH_SIZE { + let a = range.sample(&mut rng); + let b = range.sample(&mut rng); if a * a + b * b <= 1.0 { - 1 - } else { - 0 + count += 1; } - }, - ) - .reduce(|| 0, |a, b| a + b); + } + count + }) + .reduce(|| 0usize, |a, b| a + b); + + // prints something close to 3.14159... println!( "π is approximately {}", - 4. * (in_circle as f64) / (total as f64) + 4. * (in_circle as f64) / ((BATCH_SIZE * BATCHES) as f64) ); } From 8742a3cefa08565db99ccc96ab18c27894f39b78 Mon Sep 17 00:00:00 2001 From: Mason Kramer Date: Wed, 6 Jul 2022 10:59:30 -0400 Subject: [PATCH 3/4] discuss determinism in the context of rayon + rand --- examples/rayon-monte-carlo.rs | 22 +++++++++++++++++----- 1 file changed, 17 insertions(+), 5 deletions(-) diff --git a/examples/rayon-monte-carlo.rs b/examples/rayon-monte-carlo.rs index 859e9ad4f5d..72fa0b266c4 100755 --- a/examples/rayon-monte-carlo.rs +++ b/examples/rayon-monte-carlo.rs @@ -23,6 +23,20 @@ //! We can use the above fact to estimate the value of π: pick many points in //! the square at random, calculate the fraction that fall within the circle, //! and multiply this fraction by 4. +//! +//! Note on determinism: +//! It's slightly tricky to build a parallel simulation using Rayon +//! which is both efficient *and* reproducible. +//! +//! Rayon's ParallelIterator api does not guarantee that the work will be +//! batched into identical batches on every run, so we can't simply use +//! map_init to construct one RNG per Rayon batch. +//! +//! Instead, we do our own batching, so that a Rayon work item becomes a +//! batch. Then we can fix our rng stream to the batched work item. This +//! Then we amortizes the cost of constructing the Rng over BATCH_SIZE trials. +//! Manually batching also turns out to be faster for the nondeterministic +//! version of this program as well. #![cfg(all(feature = "std", feature = "std_rng"))] @@ -37,14 +51,13 @@ static BATCHES: u64 = 1000; fn main() { let range = Uniform::new(-1.0f64, 1.0); - // We have to manually batch the work so that we can get Rng construction - // out of the hot loop. Due to work-stealing, ParallelIterator::map_init - // doesn't guarantee the same # and order of work items per run. let in_circle = (0..BATCHES) .into_par_iter() .map(|i| { let mut rng = ChaCha8Rng::seed_from_u64(SEED); - // using set_stream() here force every batch to use a different part of the rng stream + // We chose ChaCha because it's fast, has suitable statical properties for simulation, + // and because it supports this set_stream() api, which lets us chose a different stream + // per work item. ChaCha supports 2^64 independent streams. rng.set_stream(i); let mut count = 0; for _ in 0..BATCH_SIZE { @@ -58,7 +71,6 @@ fn main() { }) .reduce(|| 0usize, |a, b| a + b); - // prints something close to 3.14159... println!( "π is approximately {}", From 36b7be28ac0a9db6fb39d18da37767a1cb18c9dd Mon Sep 17 00:00:00 2001 From: Mason Kramer Date: Wed, 6 Jul 2022 11:02:13 -0400 Subject: [PATCH 4/4] reword the discussion --- examples/rayon-monte-carlo.rs | 8 ++++---- 1 file changed, 4 insertions(+), 4 deletions(-) diff --git a/examples/rayon-monte-carlo.rs b/examples/rayon-monte-carlo.rs index 72fa0b266c4..041f7379c48 100755 --- a/examples/rayon-monte-carlo.rs +++ b/examples/rayon-monte-carlo.rs @@ -33,10 +33,10 @@ //! map_init to construct one RNG per Rayon batch. //! //! Instead, we do our own batching, so that a Rayon work item becomes a -//! batch. Then we can fix our rng stream to the batched work item. This -//! Then we amortizes the cost of constructing the Rng over BATCH_SIZE trials. -//! Manually batching also turns out to be faster for the nondeterministic -//! version of this program as well. +//! batch. Then we can fix our rng stream to the batched work item. +//! Batching amortizes the cost of constructing the Rng from a fixed seed +//! over BATCH_SIZE trials. Manually batching also turns out to be faster +//! for the nondeterministic version of this program as well. #![cfg(all(feature = "std", feature = "std_rng"))]