Exercise 6, Sheet 1

mainEx6.cpp

@@ -0,0 +1,111 @@
+#include <iostream>
+#include <vector>
+#include <chrono>
+#include <cmath>
+#include <algorithm>
+#include <utility>  
+#include <iomanip> 
+#include "mayer_primes.h"
+using namespace std;
+
+#include "mayer_primes.h"
+
+int single_goldbach(int k, const vector<int>& primes, const vector<char>& is_prime) //how many pairs (p,q) satisfy p+q=k 
+{
+    int count = 0;
+    for (size_t i = 0; i < primes.size(); ++i) 
+    {
+        int p = primes[i];
+        if (p > k / 2) break; 
+        int q = k - p; //compute the complementary number 
+        if (q >= 0 && q < (int)is_prime.size() && is_prime[q]) //check if the q is prime too 
+            ++count; 
+    }
+    return count;
+}
+
+vector<pair<int,int>> count_goldbach(int n, const vector<int>& primes, const vector<char>& is_prime)
+{
+    vector<pair<int,int>> results;
+    for (int k = 4; k <= n; k += 2) //loop over evn numbers 
+    {
+        int c = single_goldbach(k, primes, is_prime); // how many pairs do we have for each k 
+        results.push_back({k, c}); //append the result for each k 
+    }
+    return results;
+}
+
+vector<vector<pair<int,int>>> all_decompositions(int n, const vector<int>& primes, const vector<char>& is_prime)
+{
+    vector<vector<pair<int,int>>> out; //container for all the decopomistions of one k 
+    for (int k = 4; k <= n; k += 2) //loop over even numbers 
+    {
+        vector<pair<int,int>> decomp; // temporary vector to accumulate pairs of that particular k 
+        for (size_t i = 0; i < primes.size(); ++i)
+        {
+            int p = primes[i];
+            if (p > k / 2) break;
+            int q = k - p;
+            if (q >= 0 && q < (int)is_prime.size() && is_prime[q])
+                decomp.push_back({p, q}); //same as earlier, check if q and p are both prime, if yes, saves the tuple 
+        }
+        out.push_back(std::move(decomp)); //add all of the tuples to the vector out for this k
+    }
+    return out;
+}
+
+vector<char> make_is_prime(int n, const vector<int>& primes)
+{
+    vector<char> is_prime(n + 1, 0); //makes a vector 0...n with =1 if x is prime, 0 if not. 
+    for (int p : primes)
+    {
+        if (p <= n) is_prime[p] = 1;
+    }
+    return is_prime;
+}
+
+int main()
+{
+    //vector<int> test_ns = {10000, 100000, 400000, 1000000, 2000000};
+    vector<int> test_ns = {10000, 100000, 400000};
+    for (int n : test_ns)
+    {
+        cout << "Computing up to n = " << n << " even numbers from 4 to n \n";
+        auto t0 = chrono::system_clock::now();
+
+        vector<int> primes = get_primes<int>(n); //generates all prime
+        auto t1 = chrono::system_clock::now();
+        chrono::duration<double> gen_time = t1 - t0;
+        
+        cout << "Generated " << primes.size() << " primes up to " << n 
+             << " in " << gen_time.count() << " s\n"; 
+
+        vector<char> is_prime = make_is_prime(n, primes); //creates the boolean vector for that particular n
+
+        auto start = chrono::system_clock::now();
+        vector<pair<int,int>> results = count_goldbach(n, primes, is_prime); //computes for all even k up to n 
+        auto end = chrono::system_clock::now();
+        chrono::duration<double> elapsed = end - start;
+
+        int max_k = 0;
+        int max_count = 0;
+        for (size_t i = 0; i < results.size(); ++i) //find the k with most decomposition
+        {
+            int k = results[i].first;
+            int count = results[i].second;
+            if (count > max_count)
+            {
+                max_count = count;
+                max_k = k;
+            }
+        }
+        cout << "Max decompositions: k = " << max_k
+            << " with " << max_count << " decompositions\n";
+
+        cout << "Time to count decompositions (count_goldbach): "
+             << elapsed.count() << " s\n";
+        cout << "--- \n";
+    }
+
+    return 0;
+}