Add Pi sample

Signed-off-by: Dmitri Mokhov <[email protected]>

examples/CMakeLists.txt

@@ -66,6 +66,7 @@ tbb_add_example(parallel_for_each parallel_preorder)
 tbb_add_example(parallel_pipeline square)
 
 tbb_add_example(parallel_reduce convex_hull)
+tbb_add_example(parallel_reduce pi)
 tbb_add_example(parallel_reduce primes)
 
 tbb_add_example(task_arena fractal)

examples/README.md

@@ -19,6 +19,7 @@ This directory contains example usages of oneAPI Threading Building Blocks.
 | parallel_for_each/parallel_preorder | Parallel preorder traversal of a graph.
 | parallel_pipeline/square | Another string transformation example that squares numbers read from a file.
 | parallel_reduce/convex_hull | Parallel version of convex hull algorithm (quick hull).
+| parallel_reduce/pi | Parallel version of calculating π by numerical integration.
 | parallel_reduce/primes | Parallel version of the Sieve of Eratosthenes.
 | task_arena/fractal |The example calculates two classical Mandelbrot fractals with different concurrency limits.
 | task_group/sudoku | Compute all solutions for a Sudoku board.

examples/parallel_reduce/README.md

@@ -4,4 +4,5 @@ Examples using `parallel_reduce` algorithm.
 | Code sample name | Description
 |:--- |:---
 | convex_hull | Parallel version of convex hull algorithm (quick hull).
+| pi | Parallel version of calculating π by numerical integration.
 | primes | Parallel version of the Sieve of Eratosthenes.

examples/parallel_reduce/pi/CMakeLists.txt

@@ -0,0 +1,33 @@
+# Copyright (c) 2023 Intel Corporation
+#
+# Licensed under the Apache License, Version 2.0 (the "License");
+# you may not use this file except in compliance with the License.
+# You may obtain a copy of the License at
+#
+#     http://www.apache.org/licenses/LICENSE-2.0
+#
+# Unless required by applicable law or agreed to in writing, software
+# distributed under the License is distributed on an "AS IS" BASIS,
+# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+# See the License for the specific language governing permissions and
+# limitations under the License.
+
+cmake_minimum_required(VERSION 3.1)
+
+project(pi CXX)
+
+include(../../common/cmake/common.cmake)
+
+set_common_project_settings(tbb)
+
+add_executable(pi main.cpp pi.cpp)
+
+target_link_libraries(pi TBB::tbb Threads::Threads)
+target_compile_options(pi PRIVATE ${TBB_CXX_STD_FLAG})
+
+set(EXECUTABLE "$<TARGET_FILE:pi>")
+set(ARGS "")
+set(PERF_ARGS auto 100000000000)
+
+add_execution_target(run_pi pi ${EXECUTABLE} "${ARGS}")
+add_execution_target(perf_run_pi pi ${EXECUTABLE} "${PERF_ARGS}")

examples/parallel_reduce/pi/README.md

@@ -0,0 +1,23 @@
+# Pi sample
+Parallel version of calculating π by numerical integration.
+
+## Building the example
+```
+cmake <path_to_example>
+cmake --build .
+```
+
+## Running the sample
+### Predefined make targets
+* `make run_pi` - executes the example with predefined parameters
+* `make perf_run_pi` - executes the example with suggested parameters to measure the oneTBB performance
+
+### Application parameters
+Usage:
+```
+pi [n-of-threads=value] [n-of-intervals=value] [silent] [-h] [n-of-threads [n-of-intervals]]
+```
+* `-h` - prints the help for command line options.
+* `n-of-threads` - the number of threads to use; a range of the form low\[:high\], where low and optional high are non-negative integers or `auto` for a platform-specific default number.
+* `n-of-intervals` - the number of intervals to subdivide into, must be a positive integer.
+* `silent` - no output except elapsed time.

examples/parallel_reduce/pi/common.h

@@ -0,0 +1,48 @@
+/*
+    Copyright (c) 2023 Intel Corporation
+
+    Licensed under the Apache License, Version 2.0 (the "License");
+    you may not use this file except in compliance with the License.
+    You may obtain a copy of the License at
+
+        http://www.apache.org/licenses/LICENSE-2.0
+
+    Unless required by applicable law or agreed to in writing, software
+    distributed under the License is distributed on an "AS IS" BASIS,
+    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+    See the License for the specific language governing permissions and
+    limitations under the License.
+*/
+
+#ifndef TBB_examples_pi_H
+#define TBB_examples_pi_H
+
+#include <cstdlib>
+
+typedef std::size_t number_t;
+typedef double pi_t;
+
+extern const number_t chunk_size;
+extern number_t num_intervals;
+extern pi_t inv_intervals;
+
+extern bool silent;
+
+inline pi_t pi_kernel(number_t i) {
+    pi_t dx = (pi_t(i) + pi_t(0.5)) * inv_intervals;
+    return pi_t(4.0) / (pi_t(1.0) + dx * dx);
+}
+
+inline double pi_slice_kernel(number_t slice, number_t slice_size = chunk_size) {
+    pi_t pi = pi_t(0.0);
+    for (number_t i = slice; i < slice + slice_size; ++i) {
+        pi += pi_kernel(i);
+    }
+    return pi;
+}
+
+void init_threading(int p);
+void destroy_threading();
+double compute_pi_parallel();
+
+#endif //  TBB_examples_pi_H

examples/parallel_reduce/pi/main.cpp

@@ -0,0 +1,101 @@
+/*
+    Copyright (c) 2023 Intel Corporation
+
+    Licensed under the Apache License, Version 2.0 (the "License");
+    you may not use this file except in compliance with the License.
+    You may obtain a copy of the License at
+
+        http://www.apache.org/licenses/LICENSE-2.0
+
+    Unless required by applicable law or agreed to in writing, software
+    distributed under the License is distributed on an "AS IS" BASIS,
+    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+    See the License for the specific language governing permissions and
+    limitations under the License.
+*/
+
+#include "oneapi/tbb/tick_count.h"
+
+#include "common/utility/get_default_num_threads.hpp"
+#include "common/utility/utility.hpp"
+
+#include "common.h"
+
+const number_t chunk_size = 4096; // Multiple of 16, to fit float datatype to a vector register.
+
+// number of intervals
+number_t num_intervals = 1000000000;
+pi_t inv_intervals = pi_t(0.0);
+
+bool silent = false;
+
+double compute_pi_serial() {
+    double ret = 0;
+
+    inv_intervals = pi_t(1.0) / num_intervals;
+
+    number_t tail = num_intervals % chunk_size;
+    number_t last = num_intervals - tail;
+
+    for (number_t slice = 0; slice < last; slice += chunk_size) {
+        ret += pi_slice_kernel(slice);
+    }
+    ret += pi_slice_kernel(last, tail);
+    ret *= inv_intervals;
+
+    return ret;
+}
+
+int main(int argc, char* argv[]) {
+    try {
+        tbb::tick_count main_start_time = tbb::tick_count::now();
+        // zero number of threads means to run serial version
+        utility::thread_number_range threads(utility::get_default_num_threads, 0);
+
+        utility::parse_cli_arguments(
+            argc,
+            argv,
+            utility::cli_argument_pack()
+                //"-h" option for for displaying help is present implicitly
+                .positional_arg(threads, "n-of-threads", utility::thread_number_range_desc)
+                .positional_arg(num_intervals, "n-of-intervals", "number of intervals")
+                .arg(silent, "silent", "no output except time elapsed"));
+
+        for (int p = threads.first; p <= threads.last; p = threads.step(p)) {
+            pi_t pi;
+            double compute_time;
+            if (p == 0) {
+                //run a serial version
+                tbb::tick_count compute_start_time = tbb::tick_count::now();
+                pi = compute_pi_serial();
+                compute_time = (tbb::tick_count::now() - compute_start_time).seconds();
+            }
+            else {
+                //run a parallel version
+                init_threading(p);
+                tbb::tick_count compute_start_time = tbb::tick_count::now();
+                pi = compute_pi_parallel();
+                compute_time = (tbb::tick_count::now() - compute_start_time).seconds();
+                destroy_threading();
+            }
+
+            if (!silent) {
+                if (p == 0) {
+                    std::cout << "Serial run:\tpi = " << pi << "\tcompute time = " << compute_time
+                              << " sec\n";
+                }
+                else {
+                    std::cout << "Parallel run:\tpi = " << pi << "\tcompute time = " << compute_time
+                              << " sec\t on " << p << " threads\n";
+                }
+            }
+        }
+
+        utility::report_elapsed_time((tbb::tick_count::now() - main_start_time).seconds());
+        return 0;
+    }
+    catch (std::exception& e) {
+        std::cerr << "error occurred. error text is :\"" << e.what() << "\"\n";
+        return 1;
+    }
+}

examples/parallel_reduce/pi/pi.cpp

@@ -0,0 +1,55 @@
+/*
+    Copyright (c) 2023 Intel Corporation
+
+    Licensed under the Apache License, Version 2.0 (the "License");
+    you may not use this file except in compliance with the License.
+    You may obtain a copy of the License at
+
+        http://www.apache.org/licenses/LICENSE-2.0
+
+    Unless required by applicable law or agreed to in writing, software
+    distributed under the License is distributed on an "AS IS" BASIS,
+    WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
+    See the License for the specific language governing permissions and
+    limitations under the License.
+*/
+
+#include "common.h"
+#include "oneapi/tbb/blocked_range.h"
+#include "oneapi/tbb/global_control.h"
+#include "oneapi/tbb/parallel_reduce.h"
+
+struct reduce_body {
+    double my_pi;
+    reduce_body() : my_pi(0) {}
+    reduce_body(reduce_body& x, tbb::split) : my_pi(0) {}
+    void operator()(const tbb::blocked_range<number_t>& r) {
+        my_pi += pi_slice_kernel(r.begin(), r.size());
+    }
+    void join(const reduce_body& y) {
+        my_pi += y.my_pi;
+    }
+};
+
+double compute_pi_parallel() {
+    inv_intervals = pi_t(1.0) / num_intervals;
+
+    double ret = 0.0;
+
+    reduce_body body;
+    tbb::parallel_reduce(tbb::blocked_range<number_t>(0, num_intervals), body);
+
+    ret = body.my_pi * inv_intervals;
+
+    return ret;
+}
+
+static std::unique_ptr<tbb::global_control> gc;
+
+void init_threading(int p) {
+    gc.reset(new tbb::global_control(tbb::global_control::max_allowed_parallelism, p));
+}
+
+void destroy_threading() {
+    gc.reset();
+}