# Largest prime factor for a given number

I created a program to find the largest prime factor for a given number. The program works but unfortunately takes very long time to compute a huge number like 600851475143.

How can I optimise this program to run faster?

I don't know what does this multicore thing is, but I think it is a possible way to improve this program.

// A program to find the largest prime factor for a given number

#include <iostream>
#include <math.h>

using namespace std;

bool IsPrime (int n) {
if (n < 2)      return false;
if (n < 4)      return true;
if (n % 2 == 0) return false;

for (int i = 3; i <= int(sqrt(n)) + 1; i += 2)
if (n % i == 0)
return false;
return true;
}

int main () {

long long x;
int lpf = 0; //Largest Prime Factor

cin >> x;

for (int n = 0; n <= x; n++)
if (IsPrime(n) && x % n == 0)
lpf = n;

cout << lpf;
}

• You have optimised your loop by jumping by 2. There is another small improvement that you increment alternatively by 2 and 4. This is the equivalent of removing multiples of 2 and 3 from all tests. – Martin York Jan 6 '13 at 15:36
• Consider improving question title. It doesn't apply to any techniques that's related with "for" loop optimization (nothing related with rewinding or skipping proved conditions with inlining or anything else). My suggestion is "Largest prime factor algorithm review" – ony Jan 6 '13 at 19:16

## 6 Answers

You have done the optimization of removing multiples of 2 (by incrementing by 2).

There is another slight improvement that removes multiples of 2 and 3 from the test.

bool IsPrime (int n) {
if (n == 2)     return true;
if (n == 3)     return true;
if (n < 5)      return false;
if (n % 2 == 0) return false;
if (n % 3 == 0) return false;

int inc = 4;
int dec = -2;
for (int i = 5; i <= int(sqrt(n)) + 1; i += inc)
{
if (n % i == 0)
return false;
inc = inc + dec;  // makes inc go 2/4/2/4/2/4/.....
dec = dec * -1;
}
return true;
}


Another optimization for IsPrime() is that you only need to test using primes. ie no pointing in testing with 6 (its not prime) as all values that are divisible by 6 are also divisible by 3 (is prime). So if your loop only tested by using other primes:

// Assumes you are calling this function with monotimically inreassing
// values of n so that the foundPrimes vector will be constructed in order.
bool IsPrime (int n) {

static std::vector<int>   foundPrimes = { 2, 3, 5, 7, 11 };

for(std::size_t loop = 0; loop <  foundPrimes.size(); ++loop)
{
if (n == foundPrimes[loop])     { return true;}
if (n % foundPrimes[loop] == 0) { return false;}
}
start = foundPrimes.last() + 2;

for (int i = start; i <= int(sqrt(n)) + 1; i += 2)
{
if (n % i == 0)
return false;
}
foundPrimes.push_back(n);
return true;
}


The next optimization is called the Sieve of Eratosthenes. This needs a minor re-write but is very efficient. You would use the sieve to calculate all the primes you need first. Then just look up the values be checking if they are in the sieve.

• Last one have error. Consider that IsPrime(6857); IsPrime(12);. After storing prime number 6857 you'll start from 6857+5 and thus everything beetween 11 and 6857 will be considered as a prime numbers. As for the first one - I don't know such approach at all and have no idea how it works. Isn't it much easier to just jump in +6 (skipping all divisors of 2 and 3). It would be great to hear how that works. – ony Jan 6 '13 at 19:08
• If first one is variant of wheel factorization please confirm. – ony Jan 6 '13 at 19:14
• @ony: There is no error in the second one. You forgot to read the comments (Assumes you are calling this function with monotonically increasing values). In the first one skipping +6 is not the same. Its a standard optimization for prime: This might explaining it better codereview.stackexchange.com/a/7342/507 – Martin York Jan 7 '13 at 3:28
• @ony: I don't think it is wheel factorization. As that requires you to use another method to build a set of primes first. – Martin York Jan 7 '13 at 3:34
• not only monotonically but also without gaps. Consider then IsPrime(6857) == true; IsPrime(28561 /* 13^4 */) == true; (6857, 28561 is monotonically increasing and 13^3 < 6857). Look at the algorithm of wheel factorization it will for analyzing only sequence of 5, 7, 11, 13, 17, 19, 23, 25, 29, 31 - exactly the same as that loop does. You simply change the way you represent wheel addition your wheel is 6 and the only 2 sectors in that wheel may contain primes starting from 5. Distances between those sectors are +2 and +4 which is what that loop does (adds 5 +2, +4, +2, +4). – ony Jan 7 '13 at 9:24

As it was already pointed out, in C++ you should be including cmath instead of math.h.

The following two checks are bad for performance:

if (n < 2) return false;
if (n < 4) return true;


They will incur in two comparisons that would only return if there was very very little work to do. So its optimizing the already fast cases, and adding work to the already lengthy cases. On the other hand, this check is ok:

if (n % 2 == 0) return false;


Moving to your loop, there is no need to calculate sqrt(n) each iteration. You could instead do:

for (int i = 3, top = int(sqrt(n)) + 1; i <= top; i += 2)
...

• For the two if statements, I don't have much to do with. I edited my for loop as you pointed out. Thank you very much. – Eissa Jan 5 '13 at 18:34
• Pulling the sqrt(n) out of the test will be unlikely to help. compilers can do that optimization quite well already. – Martin York Jan 6 '13 at 15:35
• @Loki Astari: The fact that the compiler can optimize it is no reason to be lazy and express the incorrect intention. Plus, what makes you think this typical homework assignment is being compiled on anything modern? I wouldn't dare go around telling people compilers can figure out externally declared functions produce no side-effects so they can reorder your code and do once what you asked it to do N-times... it doesn't hold for the general case. Eventually this could be changed to work with a 3rd party library for very big numbers. – K-ballo Jan 6 '13 at 15:48
• @K-ballo: I disagree that it is the lazy option. I put it to you that this is the most intuitive version to read and thus should be the default way to write it. Peephole optimizations (done by a human) like this should only be done when you can show that it is sub-optimal and that mangling the code does actual get you an improvement. Also when you do it you should also write a huge comment explaining why you have made the code look messey otherwise somedy will come along and tidy it up after you thus destroying the optimization. – Martin York Jan 6 '13 at 16:11
• @LokiAstari: I share your opinion in general, but not for this particular case. Those are two different codes, asking the compiler to do two different things, with different complexity guarantees. They are not conceptually the same thing, even if the optimizer generates the same code for both. – K-ballo Jan 7 '13 at 16:41

First, a couple of suggestions:

• If you're using <iostream>, you might as well use #include <cmath> instead of #include <math.h>, for consistency if nothing else.
• avoid using namespace std;, it's not worth the loss of readability for saving a couple of instances of std::.
• You already test for 3 (< 4), so why not start your loop from 5?
• Use {} for your loops, even if the content is a single statement. I'm happy enough with the IsPrime ifs seeing as they're essentially early-exit preconditions, but anything else should have them as well.

In terms of your algorithm, the biggest issue is that you "forget" which are prime numbers for each iteration of the (outer) loop, meaning that you have to work them out every time.

What you could do is to build up a (descending) list of primes previously seen, then check whether each list number is a factor of the current value - if the list is descending, then a match is automatically the LPF. That will give you a much smaller searchspace.

• Thank you for your suggestions. I found that changing the main function to the following helped speeding up the program execution: int main(){ long long x; cin >> x; long long lpf = 0;//Largest Prime Factor for (long long n = 3; n*n <= x; n++) if (IsPrime(n) && x%n == 0) lpf = n; cout << lpf; } – Eissa Jan 5 '13 at 19:35

You could use different approach of searching LPF in a single pass (with O(sqrt(n))).

• assume your current x as a potential last prime number
• loop n through number 2 up while your current x is greater
• if you found some n that is factor of x divide your x at it as much as you can.
• if by some reason everything is left from x is 1 than that n is the largest prime factor
• otherwise you've got your new assumption about largest prime number and continue reducing it through cycle
• if you've reached situation when n is larger than x (which is greter than 1) then your x is actually the largest prime factor of your original x

Sample:

#include <iostream>
#include <cstdlib>
long long lpf(long long x) {
long long n;
// since x changes anyway I think it's better to use n*n instead of sqrt
for(n = 2; (n*n) < x ; ++n)
{
if ((x%n) == 0)
{
x /= n;
while((x%n) == 0) x /= n;
if (x == 1) return n;
}
}
return x;
}
int main(int argc, char *argv[]) {
if (argc < 2) return 1;
long long x = atoll(argv);
std::cout << "x = " << x << std::endl;
std::cout << "LPF(x) = " << lpf(x) << std::endl;
return 0;
}


This is variation of Trivial division algorithm without prime numbers at all. See Factoring algorithms for some more.

• It's strange that this answer - which is the only valid answer in this topic (apart from the one by PassingStranger which came much later) - received so little attention (upvotes) from punters. All the other respondents 'optimised' irrelevant fluff and completely missed the boat... – DarthGizka Feb 13 '16 at 21:27

Hmm... Am I the only one who sees complexity of original algorithm nearly O(n^2) even with caching prime numbers in IsPrime? If you still want to stick with prime numbers for factoring consider building lazy primes list and walk over it rather than over integers (it also can be used for algorithm I've posted earlier).

Long time ago I used this sinppet of code (with some fresh updates):

class Primes {
public:
typedef unsigned int Num;
typedef std::list<Num> list;
typedef list::const_iterator iterator;

static void init(size_t pagesz = 10)
{
page_size = pagesz;
primes.clear();
primes.push_back(2);
primes.push_back(3);
}

static iterator begin() { return primes.begin(); }

static Num next(iterator &i)
{
iterator j = i++;
if (i == primes.end()) {
Num x = primes.back();
size_t k;
for(x+=2,k=0;k<page_size;x+=2) {
if (isPrime(x)) {
primes.push_back(x);
++k;
}
}
// move i through end()
i = j;
++i;
}
return (*j);
}

private:
static list primes;
static size_t page_size;

static bool isPrime(Num v)
{
for(iterator i = begin(); i != primes.end(); ++i) {
const Num p = *i;
if (p > v) return true;
if ((v%p) == 0) return false;
}
return true;
}
};


P.S. wheel factorization with more factors than just 2 can also be used here.

Quick answer: Don't optimize IsPrime, optimize main.

First short optimization is to swap

if (IsPrime(n) && x % n == 0)


for

if (x % n == 0 && IsPrime(n))


That way you won't even run IsPrime most of the time.

As for the for-loop, every optimization you are doing in IsPrime you should do to that loop.
That means that you start with n = 2 as a special case and then only go through odd numbers.
If you also used the method described in ony's answer, you could just copy the relevant parts from IsPrime and not call the function at all.

Edit: Here's what it would look like if you put it all in main:

int main () {

long long x;
int lpf = 1; //Largest Prime Factor

cin >> x;

// Doing division by two separately.
if (x % 2 == 0) {
lpf = 2;
do {
x /= 2;
} while (x % 2 == 0)
}

// The squareroot here allows us to stop when we are down to a single prime.
for (int n = 3, end = (int)sqrt(x); n <= end; n += 2) {
if (x % n == 0) {
lpf = n;
do {
// This part will also guarantee that you won't get any composite divisors.
x /= n;
} while (x % n == 0)
end = (int)sqrt(x);
}
}

// This snippet will give you the largest prime unless it was caught earlier.
if (x > 1)
lpf = x;

cout << lpf;
}


In short, don't just use IsPrime.
Use the methods in it on the actual problem at hand.