I'm trying to learn some multithreading and this seems like a good example, want to apply this to frustum culling later. The goal is to find prime numbers in randomly generated vector of numbers. Single threaded solution is simple, multithreaded way is to divide vector of random numbers in two parts, solve each one in the same time and merge results.

Here is the time it took to solve:

mt:

solved in: 17185.4 ms

solved in: 15832.2 ms

solved in: 16012.5 ms

solved in: 16354.9 ms

st:

solved in: 26942 ms

solved in: 30577.9 ms

solved in: 26403.4 ms

solved in: 29625.7 ms

It obviously works, but is the code correct?

This is from main.cpp:

 #include "primeNumbers.h"
int main(){
numbers.Solve(solver);
delete solver;
}


This is the rest of the code:

//primeNumbers.h
#include <mutex>
#include <vector>
#include <iostream>
public:
bool isPrime(int n) {
for (int i = 2; i < n/2; ++i) {
if (n % i == 0)
return false;
}

return true;
}

virtual void solve(const std::vector<int>& source, std::vector<int>& result) = 0;

};

private:
virtual void solve(const std::vector<int>& source, std::vector<int>& result) override{
for (unsigned i = 0; i < source.size(); ++i) {
if (isPrime(source[i]))
result.push_back(source[i]);
}
}
};

private:
virtual void solve(const std::vector<int>& source, std::vector<int>& result) override {
std::vector<int> r1;
std::vector<int> r2;

//split source
int half = source.size() / 2;
int end = source.size();

std::mutex m;

for (int i = 0; i < half; ++i) {
//std::lock_guard<std::mutex> lock(m);  //do i need to lock mutex before accessing source
if (isPrime(source[i]))
r1.push_back(source[i]);
}
});
for (int i = half; i < end; ++i) {
//std::lock_guard<std::mutex> lock(m);  //do i need to lock mutex before accessing source
if (isPrime(source[i]))
r2.push_back(source[i]);
}
});

if(t1.joinable()) t1.join();
if(t2.joinable()) t2.join();

//merge results
result.insert(result.end(), r1.begin(), r1.end());
result.insert(result.end(), r2.begin(), r2.end());
}
};

public:
PrimeNumbers(int size, unsigned min, unsigned max) {
std::mt19937 rng;
for (int i = 0; i < size; ++i) {
rng.seed(std::random_device()());
std::uniform_int_distribution<std::mt19937::result_type> dist6(min, max);
randoms.push_back(dist6(rng));
}
std::cout << "Numbers are generated" << std::endl;
}

double time = glfwGetTime();
solver->solve(randoms, primes);
std::cout << "solved in: " << (glfwGetTime() - time) * 1000.0 << " ms" << std::endl;
}

void printPrimes() {
for (unsigned i = 0; i < primes.size(); ++i)
std::cout << primes[i] << std::endl;
}

std::vector<int> randoms;
std::vector<int> primes;
};

• "It obviously works, but is the code correct?" Have you tested the output? Does it work the way you expected it to?
– Mast
Aug 30, 2018 at 9:25
• Output is tested and time in displayed in ms. I am not sure whether or not to use mutex lock before accessing source vector in thread's lambda function and maybe i missed something else. Aug 30, 2018 at 9:41
• Concurrently accessing or modifying different elements is safe -cplusplus.com/reference/vector/vector/operator[]
– user52292
Aug 30, 2018 at 12:56

This is a poor way to test primality:

bool isPrime(int n) {
for (int i = 2; i < n/2; ++i) {
if (n % i == 0)
return false;
}
return true;
}


Even sticking with trial division, we can stop at √n, which is generally much less than n/2. It's also well worth treating 2 as a special case, to reduce the iterations by half (i.e. i += 2 instead of ++i).

Why are we using signed integers, given that we don't need any negative values?

You say the "output is tested", but don't show any of the test code. It's useful to see how the test works, to identify tests that are missing or erroneously passing. I mention this because it's notoriously difficult to create reliable tests for multithreaded code.

In PrimeNumberMultiThread::solve, we have opportunity to reduce duplication, and we don't need to hard-code the number of threads to two. Instead, we can use as many threads as we have cores available, using a for loop to start them and another for loop to join them all.

The mutex m and the commented-out lock_guard should just be removed: these aren't necessary, because only one thread accesses each vector between the worker threads starting and finishing (which points are memory barriers).

Instead of result.insert(), the general case might benefit from std::move algorithm - in any case, it's well worth a reserve before adding the results.

With OpenMP, much of the work can be automated for you - you'll just need a custom reduction to combine vectors - I think that's available from OpenMP 4.0 onwards.

An alternative to using vectors is to use an array of flags to indicate which values are to be included in the result (the flags could even be in the input vector: overwrite the non-primes with zero to discard them). Then there's no reallocation within either thread, and you can use OpenMP of older versions because we no longer need a custom reduction.

Here's one version that uses the OpenMP approach I've just described:

#include <algorithm>
#include <iostream>
#include <iterator>
#include <random>
#include <vector>

bool isPrime(unsigned n) {
if (n < 2) return false;
if (n < 4) return true;
if (n%2 == 0) return false;
for (auto i = 3ul;  i*i < n;  i += 2) {
if (n % i == 0)
return false;
}

return true;
}

std::vector<unsigned> solve(std::vector<int>& source)
{
auto const len = source.size();
std::size_t result_size = 0;

#pragma omp parallel for reduction(+:result_size)
for (std::size_t i = 0;  i < len;  ++i) {
if (isPrime(source[i])) {
++result_size;
} else {
source[i] = 0;
}
}

std::vector<unsigned> result;
result.reserve(result_size);

std::copy_if(source.begin(), source.end(),
std::back_inserter(result),
[](unsigned i){ return i; });

return result;
}

int main(){
const std::size_t size = 1000000u;
const unsigned min = 0u;
const unsigned max = 5000000u;

std::vector<int> randoms;
randoms.reserve(size);

{
std::uniform_int_distribution<std::mt19937::result_type> dist(min, max);
std::mt19937 rng;
rng.seed(0);            // for a reproducible test
for (std::size_t i = 0;  i < size;  ++i) {
randoms.push_back(dist(rng));
}
}

auto primes = solve(randoms);
std::clog << "Found " << primes.size() << " primes "
<< " from " << randoms.size() << " inputs\n";
}


Compile with g++ -fopenmp, or the equivalent for your compiler.

To test a single-threaded version, we just have to remove (manually, or with #ifdef) the #pragma omp line.

Note that since C++17, we're able to specify an execution policy for std::copy_if, so we could reduce that code even further, especially if we're not fussy about the order of the results:

#include <execution>

std::vector<unsigned> solve(std::vector<int>& source)
{
std::vector<unsigned> result;

std::copy_if(std::execution::par_unseq
source.begin(), source.end(),
std::back_inserter(result),
[](unsigned i){ return isPrime(i); });

return result;
}


I haven't been able to test this, as I'm still working with an old standard library that doesn't have <execution>.

• Thanks for the answer. Its not about prime numbers, its about multithreading. Primes are just an example. This is what i was looking for "Concurrently accessing or modifying different elements is safe -". Aug 30, 2018 at 13:03
• I've reached that bit of the code now. Aug 30, 2018 at 13:17
• Concurrently reading/writing different elements is safe; concurrently updating different vectors is safe, concurrently modifying a single vector would be unsafe. Aug 30, 2018 at 13:23