Function which finds the sum of factors from one to N

Question

After optimizing my function which I plan to use in Project Euler 95, I've run into a wall concerning further improvements. This code is destined for my library of project Euler functions, so while readability is nice, it definitely takes a backseat to performance.

I have verified that the output is correct by non-rigorous means (Pasting 10+ non consecutive random numbers for many ranges of output into Wolfram Alpha, no discrepancies), but it remains possible, (but unlikely) that there may be errors. Current run time is between 470-490ms for input size 10⁷. Algorithmic or stylistic review is greatly appreciated.

public class sumDivisors2 {

    /**
     * @param args the command line arguments
     */
    public static void main(String[] args) {
        int[] divisors = divisors((int)Math.pow(10, 7)); 
    }
    /**
     * 
     * @param limit
     * @returns an array of integers, such that n[i] = the sum of divisors of i + 1
     * That is the 0th index is 1, 1st index is 2, etc. 
     * All divisors of i, including i and 1, are included in the sum. 
     * (The second index of the array will contain the value three)
     */
    public static int[] divisors(int limit){
        /**
         * The general idea behind the algorithm is to compute the sum of divisors for 
         * each successive prime, raised to a given power, and multiply them out over the
         * entire sieve.
         */
        //Setup
        long start = System.currentTimeMillis();
        int[] primes = generatePrimes(limit);
        long pStop = System.currentTimeMillis();
        int[] data = new int[limit];
        long fStart = System.currentTimeMillis();
        fill(data, 1);
        long fStop = System.currentTimeMillis();
        int index = 0;

        //Use bitshifting to quickly multiply out all values of two
        int max_pow_two = (int) Math.floor(Math.log(limit) / Math.log(2));
        for (int i = 1; i <= max_pow_two; i++) {
            for (int j = 1 << i; j <= limit; j = j + (1 << i + 1)) {
                data[j - 1] = (data[j - 1] << i + 1) - 1;
            }
        }

        for (index = 1; index < primes.length; index++) {
            int maxPow = (int) Math.floor(Math.log(limit) / Math.log(primes[index]));
            if (maxPow == 1) {
                //Break when prime is greater than sqrt(n)
                break;
            }

            int skipConstant = primes[index] - 1;
            int[] powers = new int[maxPow];
            for (int j = powers.length - 1; j > -1; j--) {
                int skipCount = 0;
                powers[j] = ((int) Math.pow(primes[index], j + 1));
                if (j < powers.length - 1) {
                    for (int k = powers[j]; k <= limit; k += powers[j]) {
                        data[k - 1] *= (powers[j] * primes[index] - 1) / (primes[index] - 1);
                        skipCount++;
                        if (skipCount == skipConstant) { //Skips every nth multiple to  multiplications                                                      //that were performed by an earlier iteration
                            k += powers[j];
                            skipCount = 0;
                        }
                    }
                } else {
                    for (int k = powers[j]; k <= limit; k += powers[j]) {
                        data[k - 1] *= (powers[j] * primes[index] - 1) / (primes[index] - 1);
                    }
                }
            }
        }
        //Iterate over the rest of the primes, which will never have more than a   first power
        for (int i = index; i < primes.length; i++) {
            for (int j = primes[i]; j <= limit; j += primes[i]) {
                data[j - 1] *= primes[i] + 1;
            }
        }
        long stop = System.currentTimeMillis();
    //        for (int i = 0; i < data.length; i++) {
    //            System.out.printf("The sum of the factors of %d is %d\n", i +          1, data[i]);
    //        }
        System.out.println("Computation time: " + (stop - start) + " ms");
        System.out.println("\tPrime generation time: " + (pStop - start) + " ms");
        System.out.println("\tArray fill time: " + (fStop - fStart) + " ms");
        System.out.println("Verification for 10 ^ 7: 24902280 = " + data[limit - 1]);
        return data;
    }
    /**
     * Simple implementation of the Sieve of Eratosthenes 
     * @param max the integer to generate primes up to
     * @returns an array consisting of all the primes 2 <= n <= max
     */
    public static int[] generatePrimes(int max) {
        boolean[] isComposite = new boolean[max + 1];
        int numPrimes = 0;
        for (int i = 2; i * i <= max; i++) {
            if (!isComposite[i]) {
                for (int j = i; i * j <= max; j++) {
                    isComposite[i * j] = true;
                }
            }
        }
        for (int i = 2; i <= max; i++) {
            if (!isComposite[i]) {
                numPrimes++;
            }
        }
        //System.out.println(numPrimes);
        int[] primes = new int[numPrimes];
        int index = 0;
        for (int i = 2; i <= max; i++) {
            if (!isComposite[i]) {
                primes[index++] = i;
            }
        }
        return primes;
    }
    /**
     * Quickly fill an int array with a given value 
     * @param array the int array   
     * @param value the value
     */
    private static void fill(int[] array, int value) {

        int len = array.length;

        if (len > 0) {
            array[0] = value;
        }

        for (int i = 1; i < len; i += i) {
            System.arraycopy(array, 0, array, i, ((len - i) < i) ? (len - i) : i);
        }

    }

}

Answer 1

Your methods are quite long, which is problematic. They're doing too much, and they're too complicated to follow.

Let's boil divisors down to just the control flow statements...

divisors {
    for {
        for {}
    }

    for {
        if {}
        for {
            if {
                for {
                    if {}
                }
            }
            else {
                for {}
            }
        }
    }

    for {
        for {}
    }
}

This is way too complicated for a single method. We need some refactoring!

The code is confusing. The comments are sparse and it certainly doesn't help that, as an example, this particular method has three different i, three different j, and two k as a result of all the nesting and multiple loops.

We've also cluttered timing code all over the place.

Comments mentioning "quickly" doing things aren't very helpful either. Should we assume the other methods are not "quick" as they don't have comments explicitly mentioning they are quick?

So... let's get some focus on what can be done more specifically...

---

/**
 * Quickly fill an int array with a given value 
 * @param array the int array   
 * @param value the value
 */
private static void fill(int[] array, int value) {

    int len = array.length;

    if (len > 0) {
        array[0] = value;
    }

    for (int i = 1; i < len; i += i) {
        System.arraycopy(array, 0, array, i, ((len - i) < i) ? (len - i) : i);
    }

}

Given the comment above this method, the actual code in this method isn't even remotely what I'd expect to see. I'd expect this to simply be a wrapper around this:

for (int i = 0; i < array.length; ++i) {
    array[i] = value;
}

If it is anything different than this, then here's what I expect to see...

1. Unit tests verifying it produces the same result as the above code.
2. Unit tests verifying it is faster than the above code.
3. A pile of comments explaining the exact logic of this non-obvious approach (or as a minimum, a link to somewhere explaining why this is faster/better than the straight forward approach).

---

You public static int[] generatePrimes(int max) seems less than ideal, and I think it's doing more work than necessary...

How about this:

public static int[] generatePrimes(int max) {
    int[] primes = new int[max];
    int primeCount = 0;

    // special cases
    if (max > 2) { primes[primeCount++] = 2; }
    if (max > 3) { primes[primeCount++] = 3; }

    // special cases take care of everything less than 5
    // always increment i by 2 or 4
    // this skips all multiples of 2, 3
    for (int i = 5, w = 2; i < max; i += w, w = 6 - w) {
        boolean isPrime = true;

        // we want to skip the same multiples our outer loop skips
        // if our outerloop skipped the multiples, 
        // there is no use checking if it's divisible by those multiples...
        for (int j = 5, x = 2; j * j < i; j += x, x = 6 -x) {
            if (i % j == 0) {
                isPrime = false;
                break;
            }
        }

        if (isPrime) {
            primes[primeCount++] = i;
        }
    }

    int[] returnArray = new int[primeCount];
    System.arraycopy(primes, 0, returnArray, 0, primeCount); 

    return returnArray;
}

This method eliminates several of your loops. It eliminates several unnecessary checks. And only does what it really needs to do.

And actually, I'm going to argue that the inner loop should be replaced with a simple call to an isPrime method.

public static boolean isPrime(int value) {
    // special cases
    if (value == 2) { return true; }
    if (value == 3) { return true; }
    if (value % 2 == 0) { return false; }
    if (value % 3 == 0) { return false; }

    // special cases takes care of everything less than 5
    // special cases takes care of all multiples of 2 or 3
    // always increment i by 2 or 4 to skip all multiples of 2, 3
    for (int i = 5, w = 2; i * i < value; i += w, w = 6 - w) {
        if (value % i == 0) {
            return false;
        }
    }

    return true;
}

Now our generatePrimes is that much simpler:

public static int[] generatePrimes(int max) {
    int[] primes = new int[max];
    int primeCount = 0;

    // special cases
    if (max > 2) { primes[primeCount++] = 2; }
    if (max > 3) { primes[primeCount++] = 3; }

    // special cases take care of everything less than 5
    // always increment i by 2 or 4
    // this skips all multiples of 2, 3
    for (int i = 5, w = 2; i < max; i += w, w = 6 - w) {
        if (isPrime(i)) {
            primes[primeCount++] = i;
        }
    }

    int[] returnArray = new int[primeCount];
    System.arraycopy(primes, 0, returnArray, 0, primeCount); 

    return returnArray;
}

Answer 2

Just a little remark as an addition to the other excellent answers: Don't repeat yourself, always try to to find an abstraction for common tasks. E.g. If you need to time methods, don't do it "manually". Of course the "right" solution depends on your needs, but I would suggest something very simple as:

public class StopWatch {
    private Map<String, Long> times = new HashMap<>();

    public <T> T time(String name, Supplier<T> supplier) {
        long start = System.nanoTime();
        T result = supplier.get();
        long stop = System.nanoTime();
        times.put(name, stop - start);
        return result;
    }

    public Long getTime(String name) {
        return times.get(name);
    }
}

Usage:

StopWatch watch = new StopWatch();
String s = watch.time("version", () -> System.getProperty("java.version"));
System.out.println("needed " + watch.getTime("version") + " ns for result " + s);

If you don't use such tools, you clutter your code, and after some code changes you might even end up measuring something you didn't intend. Further it is much easier to change behavior afterwards (e.g. printing the timings immediately, or log the timings instead) if you have only one piece of code dedicated to that task.

Answer 3

Floating point rounding errors

You use this expression more than once in your program:

        int maxPow = (int) Math.floor(Math.log(limit) / Math.log(primes[index]));

It's supposed to find the maximum power that the prime number can be raised to without exceeding the limit. So for example, if limit=16 and prime=2, then maxPow should be 4. While your expression is mathematically correct, the problem is that floating point rounding errors could cause your expression to be off by one.

As an example, if you use limit=243 and prime=3, this should result in an answer of 5 because \$3^5 = 243\$. However, using your expression, I got an answer of 4. This is because on my computer, Math.log(243)/Math.log(3) resulted in 4.999999999999999, and the floor operation rounded that down to 4.

I suggest finding maxPow using all integer arithmetic instead.

Improve sieve

Your current sieve:

    for (int i = 2; i * i <= max; i++) {
        if (!isComposite[i]) {
            for (int j = i; i * j <= max; j++) {
                isComposite[i * j] = true;
            }
        }
    }
    for (int i = 2; i <= max; i++) {
        if (!isComposite[i]) {
            numPrimes++;
        }
    }

could be improved to:

    for (int i = 2; i * i <= max; i++) {
        if (!isComposite[i]) {
            int increment = i + i;
            for (int j = i * i; j <= max; j += increment) {
                isComposite[j] = true;
            }
            numPrimes++;
        }
    }

Notice I replaced the multiplication in the inner loop with addition. Also, I merged the two loops into one.

Unnecessary array

You have an array called powers, but you only ever use the array element you just computed, powers[j]. Therefore you could eliminate the array and just use a single int to hold the current power.