C++ For Euler Problem 3

Question

Here's my solution for Euler Problem 3:

The prime factors of 13195 are 5, 7, 13 and 29.

What is the largest prime factor of the number 600851475143 ?

I'm using the Euler problems to learn some basics of C++. Are there any C++ "conventions" I'm missing, or any languages features I could use to get a better understanding of C++?

#include <iostream>
#include <vector>
#include <algorithm>

#include "fibonacci.h"
#include <limits>

int get_prime_above(int i);

bool is_prime(int number);

bool has_decimals(float number);


int main() {

    // What number are we getting the prime factors for?
    long target_number = 600851475143;

    // Holds the prime factors.
    std::vector<int> primes;

    // No point dividing by 1, it's a useless factor in this context.
    int current_prime = get_prime_above(1);
    while (target_number > 1) {
        // Does the target divide by the prime evenly?
        long result = target_number % (long) current_prime;
        if (result != 0) {
            // If not, move on to the next prime.
            current_prime = get_prime_above(current_prime);
            continue;
        }

        // If it does, push the prime
        primes.push_back(current_prime);
        // Change the target number to what ever is left to work out
        target_number = (long) target_number / (long) current_prime;
        // And reset the prime.
        current_prime = get_prime_above(1);
    }

    // Sort and reverse the primes,
    // so the highest will be at index 0
    std::sort(primes.begin(), primes.end());
    std::reverse(primes.begin(), primes.end());

    std::cout << "Largest prime: " + std::to_string(primes[0]) << std::endl;
}

/**
 * Get the next prime number above i.
 * @param i
 * @return
 */
int get_prime_above(int i = 1) {
    // Current number we're checking.
    int current_number = 0;
    while (current_number < std::numeric_limits<int>::max()) {

        // If it's lower than what the user wants, just keep incrementing.
        if (current_number <= i) {
            current_number++;
            continue;
        }

        // If it's prime, return it.
        if (is_prime(current_number)) {
            return current_number;
        }

        current_number++;
    }
    throw std::runtime_error("Couldn't find a prime within the int maximum value");
}

/**
 * Is the given number a prime?
 * @param number
 * @return bool
 */
bool is_prime(int number) {
    // Loop all the numbers between (not inclusive of) 1 and number
    // as we know that 1 and number will always have no decimal.
    // If number is 7, it would loop 2..6
    for (int tmp = number - 1; tmp > 1; tmp--) {
        // Divide the number
        float result = ((float) number / (float) tmp);

        // If the result has decimals, it doesn't divide equally
        // so return false.
        if (!has_decimals(result)) {
            return false;
        }
    }

    // All the in between numbers had decimals, so it's a prime.
    return true;
}

/**
 * Is this number a decimal?
 * If converting the float to an int and back again changes the value,
 * then the int conversion removed some decimals.
 * @param number
 * @return bool
 */
bool has_decimals(float number) {
    return number != (float) (int) number;
}

Answer 1

I tried this in VS2019 and I needed to #include <string> in order to compile.

Then I get a warning for the line

long product_of_primes = 600851475143;

which says

truncation from __int64 to long

which means that long is not enough to hold that number. Since it's positive, I changed that to unsigned long long.

For the remaining review, I'll go through your code from top to bottom, trying to understand it. Whenever I find something, I'll point it out.

---

If you change

  // Holds the prime factors.
  std::vector<int> primes;

to

std::vector<int> prime_factors;

you can get rid of the comment and the variable becomes self-explaining.

---

The method

int get_prime_above(int i = 1)

IMHO should not support a default value. As you see in your code, you're calling it with

int current_prime = get_prime_above(1);

and seeing a call with the default value

int current_prime = get_prime_above();

would not make sense to the reader.

Another open question is: which prime above i will the method return? Just any? The next? At this point of understanding and from the comment, I'd propose the name

next_prime_after(int i)

At the same time, I'd change the parameter i to minimum.

---

I stumbled over the line

while (target_number > 1) {

because I was expecting some other number increasing towards that target. However, your code reduces target_number to become 1. Given that target_number is a product of prime numbers, I would call it product_of_primes instead. Again, you can now remove the comment for that variable.

---

Again, you need a comment to explain what a variable means:

// Does the target divide by the prime evenly?
long result = product_of_primes % (long)current_prime;

Note that the comment asks a yes/no question, which is why I would expect a boolean. However, the result is a long. Also, the cast to (long) seems redundant. My proposal:

bool is_divisor = 0 == product_of_primes % current_prime;

This also changes the conditional statement to

if (!is_divisor) {

so the comment here is redundant, so it can be removed:

// If not, move on to the next prime.
current_prime = next_prime_after(current_prime);

---

In order to simplify the branch we have now:

if (!is_divisor) {
    current_prime = next_prime_after(current_prime);
    continue;
}

you can initialize current_prime to 1 and then call next_prime_after() at the beginning of the loop like this:

int current_prime = 1;
while (product_of_primes > 1) {
    current_prime = next_prime_after(current_prime);
    bool is_divisor = 0 == product_of_primes % current_prime;
    if (!is_divisor) continue;
    ...
    current_prime = 1

---

product_of_primes = (long)product_of_primes / (long)current_prime;

This IMHO has redundant casts and can be simplified to

product_of_primes = product_of_primes / current_prime;

This, however, will change the result from 71 to the correct value expected by Project Euler to complete the challenge.

---

current_prime = 1;

IMHO, you don't need to restart at the beginning, since you'll not find prime factors lower than the one you already found.

However, it will not change performance significantly.

---

Instead of reversing

std::reverse(prime_factors.begin(), prime_factors.end());

you can simply get the last element of the vector

std::cout << "Largest prime: " + std::to_string(prime_factors.back()) << std::endl;

or (my preferred option in this case) sort in descending order

std::sort(prime_factors.begin(), prime_factors.end(), std::greater<int>());

---- Dunno, the horizonal rulers stop working in this line

Woah, this is ugly:

if (current_number <= minimum) {
  current_number++;
  continue;

This can simply be replaced by

int current_number = minimum + 1;

And in consequence the loop can be re-written as

int next_prime_after(int minimum) {
    int current_number = minimum + 1;
    while (!is_prime(current_number)) {
        current_number++;
    }
    return current_number;
}

---

In that same method, we probably get a performance benefit of a factor of 2 if you move to odd numbers.

int next_prime_after(int minimum) {
    int current_number = minimum + 1;
    current_number += minimum % 2; // ensure odd number
    while (!is_prime(current_number)) {
        current_number += 2;
    }
    return current_number;
}

---

A had a WTF moment here:

float result = ((float)number / (float)tmp);
// If the result has decimals, it doesn't divide equally
// so return false.
if (!has_decimals(result)) {

Why? Because you already did a proper divisability check somewhere else (cited from the original code):

// Does the target divide by the prime evenly?
long result = target_number % (long) current_prime;

That way, you can get rid of the whole method has_decimals().

---

for (int tmp = number - 1; tmp > 1; tmp--) {

When you learned dividing numbers in school, did you start by dividing an arbiotrary number like 56 by it's predecessor 55 to check if it is a prime factor? Certainly not. So let's not do it in your code as well. Go simple: check if 56 is divisable by 2.

Oh, and BTW: if you can't divide 56 by 2, you can also not divide it by anything > 28, right? Make a huge performance impact of a factor ~10.

bool is_prime(int number) {
    for (int divisor = 2; divisor < number/2; divisor++) {
        bool is_divisor = 0 == number % divisor;
        if (is_divisor) {
            return false;
        }
    }
    return true;
}

This could be reduced to sqrt(2) instead of number/2 if you think about it.

Starting at 3 and testing odd numbers only will increase by another factor of 2:

bool is_prime(int number) {
    if (number % 2 == 0) return false;
    for (int divisor = 3; divisor < number/2; divisor+=2) {
        bool is_divisor = 0 == number % divisor;
        if (is_divisor) {
            return false;
        }
    }
    return true;
}

---

For this Project Euler question, int prime factors see to be good enough. For general purpose factorization, you should consider at least long.

---

My final proposal for you

#include <iostream>
#include <vector>
#include <algorithm>
#include <string>
#include "fibonacci.h"
#include <limits>

int next_prime_after(long minimum);
bool is_prime(long number);

int main() {

    unsigned long long product_of_primes = 600851475143;
    std::vector<int> prime_factors;

    int current_prime = 1;
    while (product_of_primes > 1) {
        current_prime = next_prime_after(current_prime);
        bool is_divisor = 0 == product_of_primes % current_prime;
        if (!is_divisor) continue;

        prime_factors.push_back(current_prime);
        product_of_primes = product_of_primes / current_prime;
        current_prime--;
    }

    std::sort(prime_factors.begin(), prime_factors.end(), std::greater<int>());
    std::cout << "Largest prime: " + std::to_string(prime_factors[0]) << std::endl;
}

/**
 * Get the next prime number above minimum.
 * @param minimum
 * @return
 */
int next_prime_after(long minimum) {
    long current_number = minimum + 1;
    current_number += minimum % 2; // ensure odd number
    while (!is_prime(current_number)) {
        current_number += 2;
    }
    return current_number;
}

/**
 * Is the given number a prime?
 * @param number
 * @return bool
 */
bool is_prime(long number) {
    if (number % 2 == 0) return false;
    for (long divisor = 3; divisor < number/2; divisor+=2) {
        bool is_divisor = 0 == number % divisor;
        if (is_divisor) {
            return false;
        }
    }
    return true;
}

Summary

I hope you agree with my proposals and see that the resulting code is

• easier to read
• bug fixed
• faster by a factor of ~100 (measured using QueryPerformanceCounter())

Answer 2

• while (current_number < std::numeric_limits<int>::max()) { Good that you discovered numeric limits, but there is no point using it here. Just check up to \$\sqrt{n}\$ should be fine.

• Your is_prime function depends on floats, which is not a good idea. This is a simpler implementation, checking whether the modulus is 0 or not (this is done with operator %).

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

• Your algorithm is quite slow, like for the function get_prime_above. Check this answer: Project Euler #3 in Java

Answer 3

    long target_number = 600851475143;

Others commented on the truncation issue here. Just to expand on how you should choose integral types.

1. Use an unsigned integer type if you need
   1. Two's-complement arithmetic or
   2. bit pattern representation
2. Use the smallest signed fixed-size integer type from <cstdint> that will suffice if you are
   1. Storing a lot of integers in a data structure or
   2. Trying to represent a value you wouldn't want to count yourself
3. Inherit types from data structures to avoid mixing signed and unsigned types (std::size_t, std::vector<T, A>::size_type, etc).
4. Use int.

In this case, you've got a really large number you probably would never want to count, so use the smallest signed fixed integer type that can represent it.

    std::int64_t target_number = 600851475143;

---

    int current_prime = get_prime_above(1);
    while (target_number > 1) {

Minimize the scope of your variables.

    for (int current_prime = 2; target_number > 1; ) {

---

    std::sort(primes.begin(), primes.end());
    std::reverse(primes.begin(), primes.end());
    // access primes[0]

You wrote a slower version of std::max_element. Do you even need to keep a container of all prime factors? Would a single element representing the maximum prime factor be sufficient?

---

    long target_number = 600851475143;

    std::vector<int> primes;
    int current_prime = get_prime_above(1);
    while (target_number > 1) {
        long result = target_number % (long) current_prime;
        if (result != 0) {
            current_prime = get_prime_above(current_prime);
            continue;
        }

        primes.push_back(current_prime);
        target_number = (long) target_number / (long) current_prime;
        current_prime = get_prime_above(1);
    }

Looking at your algorithm, you end up duplicating a lot of work because of the reset on current_prime. Consider a target number of 125. You check all values until you find a divisible prime, in this case 2, 3, 4, and 5. You then reduce the target and throw away the work you did. So you check 2, 3, 4, and 5 again, reduce, then trash the work. And again for a third time. This ends up being quite wasteful. The fundamental theorem of arithmetic says that every integer greater than \$1\$ either is a prime or can be represented as the product of prime numbers and it's representation is unique.

$$ 1200 = 2^4 \times 75 = 2^4 \times 3^1 \times 25 = 2^4 \times 3^1 \times 5^2 $$

In the above example, the remaining target will never be divisible by 2 once all factors of 2 are exhausted. You don't need to check it again. The same is true for 3, 4, and 5. Logically, instead of primality testing values, you can shortcut that and just do the divisibility test.

    std::int64_t target_number = 600851475143;
    std::int64_t largest_factor = 0;

    for (std::int64_t candidate = 2; candidate * candidate <= target_number; /* */ ) {
        if (target_number % candidate != 0) {
            ++candidate;
            continue;
        }

        target_number = target_number / candidate;
        largest_factor = candidate;
    }

    if (target_number > largest_factor) {
        largest_factor = target_number;
    }

    std::cout << "Largest prime: " << largest_factor << '\n';

You can, of course, speed this up further using wheels to skip unnecessary division checks. The simplest is skipping all evens.