Decreasing the runtime of this non-abundant sums solution

Question

I have written a solution to problem 23 of Project Euler, and after doing some testing and editing, I've managed to get the run time down from 1.7443 seconds to 1.56616 seconds. The problem is

Find the sum of all the positive integers which cannot be written as the sum of two abundant numbers.

What are some other ways that I could improve this program? I'm pasting it entirely, and it should be noted that I'm compiling with g++ -o3.

#include <chrono>
#include <cmath>
#include <iostream>
#include <ctime>
#include <ratio>

int main()
{
//Value given in the question
const int LIMIT = 28124;
//Used to measure the runtime
std::chrono::high_resolution_clock::time_point t1 = std::chrono::high_resolution_clock::now();

//Sieve of Erasthotenes to remove a significant number of iterations in the next step
int rootMax = (int)sqrt((double)LIMIT);

bool set[LIMIT+1] = {0};

for(int i = 2; i <= rootMax; ++i)
{
    if(!set[i])
    {
        for(int j = (i*i); j <= LIMIT; j += i)
        {
            set[j] = true;
        }
    }
}

bool isSumOfAb[LIMIT] = {0};

int abundant[6965];
int count = 0;
for(int i = 12; i < LIMIT; ++i) //This can start at 12, as we know that 12 is the 1st abundant number
{
    if(!set[i]) // Prime numbers won't be abundant, so next number
    {
        continue;
    }

    int sum = 0;
    for(int j = 1; j < i; ++j)
    {
        if(i % j == 0)
        {
            sum += j;   
        }
    } 

    if(sum > i)
    {
        abundant[count] = i;
        ++count;
    }
}   

//Lots of work to do here.
//Establish all of the sums and compare
int a, c;
for(int i = 0; i < 6965; ++i)
{
    a = abundant[i];

    for(int j = i; j < 6965; ++j)
    {
        c = a + abundant[j];
        //Big speed increase when made like this as opposed to `if(c <= LIMIT) isSumOfAb[c] = true;`
        if(c > LIMIT)
        {
            break;  
        }else
        {
            isSumOfAb[c] = true;
        }
    }
}

int sum = 0;
for(int i = 24; i < LIMIT; ++i) //can start at 24, as that is the 1st sum of 2 abundant numbers
{
    if(!isSumOfAb[i])
    {
        sum += i;
    }
}

std::cout << "Sum: " <<  sum;
std::chrono::high_resolution_clock::time_point t2 = std::chrono::high_resolution_clock::now();
std::chrono::duration<double> time_span = std::chrono::duration_cast<std::chrono::duration<double>>(t2 - t1);
std::cout << "\nIt took me " << time_span.count() << " seconds.";
}

Answer 1

After some tests (adding more time-points), we can find that the double for-loop used to find all abundant numbers is the slowest. Now we can think about speeding it up, by shortening the inner loop. We start with sum=1 and try every number we find can divide i upto sqrt(i). We immediatelly know that the product will of course divide the number as well, but need to take special care for numbers of the form i=j*j.

//  stop is the biggest integer such that stop*stop <= i
//  the reasoning behind using it is similar to the reasoning for using it in first phase
//  when marking all primes, but we need to take care for div=i/j when i%j == 0
    int sum = 1, stop = sqrt(i);
    for(int j = 2;;)
    {
        int div = i/j;         // divison and modulo will probably be optimized
        int mod = i%j;         // as one instruction (we could use some divmod built-in)

        if(mod == 0)           // when this condition is met we know that:
        {                      // i = j*div
            if(j == div)
            {                  // but can be i = j*j when j=div
                sum += j;      // we need to sum it only once in such case (e.g. 16=4*4)
                break;         // and can exit the loop, because j=stop here
            }
            sum += j+div;      // this is the tricky part when we sum both j and div
        }                      // e.g. for 15: 1 + 3+5 (j=3, div=5, stop=3 here)
        if(++j > stop)
            break;             // moving the condition here is another optimization
    }                          // ...and could be rewritten as do-while loop

The original code needed more than 8.8s on my old PC, but the current version needs only 0.11s. Full code:

#include <chrono>
#include <iostream>
#include <algorithm>

using the_clock = std::chrono::high_resolution_clock;
using time_point = the_clock::time_point;
using time_span = std::chrono::duration<double>;

// all integers greater than 28123 can be written as the sum of two abundant numbers
static const int LIMIT = 28124;
// smallest abundant number is 12
static const int minAbundant = 12;

// all prime numbers will have false in this array
static bool not_prime[LIMIT+1]; // zero initialized
static const int maxPrimeRoot = (int)sqrt(LIMIT);
static inline bool is_prime(int x) {
    return !not_prime[x];
}

// all abundant numbers
static const int abundantSz = 6965; // we know that, see final output
static int abundant[abundantSz];

// all numbers that are sum of abundant numbers will be marked true here
static bool isSumOfAb[LIMIT];

int main()
{
//  1. Mark all primes:
    time_point t1 = the_clock::now();

    for(int i = 2; i <= maxPrimeRoot; i++)
    {
        if(is_prime(i))
        {
            for(int j = (i*i); j <= LIMIT; j += i)
                not_prime[j] = true;
        }
    }

//  2. Find all abundant numbers
    time_point t2 = the_clock::now();

    abundant[0] = minAbundant;
    int count = 1;
    for(int i = minAbundant+1; i < LIMIT; ++i)
    {
        if(is_prime(i)) // Prime numbers won't be abundant, so next number
            continue;

    //  BIG CHANGE HERE
    //  we sum both the number we find can divide i
    //  and the product + stop after sqrt(i)
    //  we need to take special care for numbers of the form i*i (e.g. 16=4*4)
        int sum = 1, stop = sqrt(i);
        for(int j = 2;;)
        {
            int div = i/j;
            int mod = i%j;
            if(mod == 0) {
                if(j == div) {
                    sum += j;
                    break; }
                sum += j+div;
            }
            if(++j > stop)
                break;
        }

        if(sum > i)
        {
            abundant[count] = i;
            ++count;
        }
    }

//  3. Mark are numbers that are sum of abundant numbers
    time_point t3 = the_clock::now();

    for(int i = 0; i < count; ++i)
    {
        int a = abundant[i];

        for(int j = i; j < count; j++)
        {
            int c = a + abundant[j];
            if(c > LIMIT)
                break;
            isSumOfAb[c] = true;
        }
    }

//  4. Sum them all
    time_point t4 = the_clock::now();

    int sum = 0;
    for(int i = 24; i < LIMIT; ++i) //can start at 24, as that is the 1st sum of 2 abundant numbers
    {
        if(!isSumOfAb[i])
            sum += i;
    }

//  Output
    time_point t5 = the_clock::now();
    using std::cout;
    cout << "Sum: " <<  sum;
    cout << "\nIt took me " << time_span(t5-t1).count() << " seconds.";
    cout << "\nTimes: " << time_span(t2 - t1).count()
        << ", " << time_span(t3 - t2).count()
        << ", " << time_span(t4 - t3).count()
        << ", " << time_span(t5 - t4).count();
    cout << "\nAbundant count: " << count; // just a check

    if(sum != 4179595)
        cout << "\nERROR - bad sum\n";
}

You can see I have moved few variables outside main() and added few using to reduce the code a bit. But I am really not the one to talk about code-styling :D