Project Euler Problem 530: sum of sum of greatest common divisors

Question

Here is a description of Project Euler problem 530:

Every divisor \$d\$ of a number \$n\$ has a complementary divisor \$n/d\$.

Let \$f(n)\$ be the sum of the greatest common divisor of \$d\$ and \$n/d\$ over all positive divisors \$d\$ of \$n\$, that is \$\displaystyle f(n)=\sum_{d/n}\gcd\left(d, \frac{n}{d}\right)\$.

Let \$F\$ be the summatory function of \$f\$, that is \$\displaystyle F(k) = \sum_{n=1}^{k}f(n)\$.

You are given that \$F(10)=32\$ and \$F(1000)=12776\$.

Find \$F(10^{15})\$.

My solution for Euler Problem 530 is as follows:

import java.math.BigInteger;

public class problem530 {

public static void main(String[] args) {

    double start = System.nanoTime();

    BigInteger result = BigInteger.ZERO;
    BigInteger one = BigInteger.ONE;

    for (long i = 1; i <= 1000; i++){
            if ( i == 1) {
                result = result.add(one);
            }
            if (PrimeTest.testLong(i) == true) {
                result = result.add(one).add(one);
            }
            else if (i != 1 && PrimeTest.testLong(i) == false) {
                result = result.add(one);
                for(long j = 1; j*2 <= i; j++) {
                    if (i%j == 0) {
                    long sub = i / j;
                    result = result.add(BigInteger.valueOf(GCD(sub, j)));
                }
            }
        }
    }

    System.out.println(result);

    double duration = (System.nanoTime() - start)/1000000000;

    System.out.println("Your code took " + duration + " seconds to execute.");

}
}

I have PrimeTest and GCD written in libraries. The code works (it returns the correct answer for the test cases). It takes unimaginably long to get the answer for \$F(10^{15})\$ though (essentially, I won't get the answer). I'm looking for ways to optimise it.

At this link, @alan wrote something there about a sieve (presumably in C++?) which I don't quite understand. Does anyone have any ideas for optimisation?

• I perform GCD operations using the Eulerian algorithm.

• Prime Testing is performed for:

  n = 1; n*n <= (number being tested); n++
  

To give an idea of how long this code takes to execute:

• It ran F(10)=32 in 0.008 seconds.
• It ran F(1000)=12766 in 0.019 seconds.
• It ran F(100,000)=2907546 in 29.724 seconds.

This was the optimisation suggested in another post on codereview.

It seems to me that the alternative summation suggested there is wrong! F(10) would give you 22, not 32. I'm not sure how the alternative summation was derived.

Answer 1

Special Cases and Branching

If you look at your prime checking, your code is:

if (prime(i) == true) {
    // ...
}
else if (prime(i) == false) {
    // ...
}

This ends up with you checking composite numbers for primality twice! Once in each branch! The second part could have just been else. Additionally, explicitly checking for == true is an anti-pattern. Prefer:

if (prime(i)) {
    // ...
} else {
    // ...
}

That said, why is being prime a special case? Your non-prime logic would've worked just fine on primes too! And for 1! So just drop all the branches completely and simply iterate.

Need a better algorithm

10^15 is simply too big to brute force the answer. Even dropping all the unnecessary extra work you're doing is likely insufficient. You'll have to try to come up with some optimal substructure this problem has so that you can use dynamic programming.

Answer 2

I don't usually use java, but as to the algorithm, to make the suggestion of the real optimization work, you need to consider 2 things.

First, if you want to calculate F(10), this algorithm will go from d=1 to 3 in the first loop. Then from 1 to n*d <= 10. However, with this solution you are doing the gcd of (2,5) but never doing the one for (5,2).

To solve this, you could just add the gcd result*2. However, in special cases like (1,1),(2,2),(3,3) you only want to add it once because they are square numbers. So you need to add gcd*1 if d == n, and gcd*2 if d != n

This gets us to the second consideration. Problem with this is that we will do gcd*2 for (1,2) but then we will do gcd*2 with (2,1) so we will add some cases 4 times instead of 2 as we wanted.

The fix for this is pretty easy, instead of going from n=1 to nd <=k, go from n=d to nd <= k. Why? Because if n < d, we have already passed over those pairs before ((2,1) pair was already seen when d=1 and n=2)

To make all this more clear, here is a short python implementation: (I used d as x and n as y)

def solve(n):
    solution = 0
    for x in xrange(1,int(math.sqrt(n))+1): # from 1 to sqrt(n)
        y = x # to avoid repetitions, from x to x*y <= n
        while x*y<=n :
            solution+= gcd(x,y)*(2 if x != y else 1) # this adds gcd*1 if x ==y and gcd*2 if x != y
            y+=1 
    print solution   

This runs pretty fast as we are not even doing prime checks (which would be faster with a sieve, but not necessary with this solution). Python is pretty slow so it would take a lot to solve 10^15, but you could try to code it in java and check the speed improvement

In python for n = 500000, it takes a bit less than 2 seconds. I did try the same in java, but it seems to take 2.8 seconds for n = 500000, even slower than python. Might be because I had never used BigInteger before, so there might be some improvements to make, here is the code:

    double start = System.nanoTime();

    BigInteger n = new BigInteger("500000");

    BigInteger x = BigInteger.ONE;

    BigInteger solution = BigInteger.ZERO;

    while(x.multiply(x).compareTo(n) != 1){


        BigInteger y = x;


        while(y.multiply(x).compareTo(n) != 1){

            BigInteger gc = gcd(x,y);

            if(x.compareTo(y) == 0){
                solution = solution.add(gc);
            }else{
                solution = solution.add(gc).add(gc);

            }

            y = y.add(BigInteger.ONE);

        }

        x = x.add(BigInteger.ONE);

    }

    double duration = (System.nanoTime() - start)/1000000000;

    System.out.println(solution);
    System.out.println("Your code took " + duration + " seconds to execute.");

For further optimizations, you could memorize some gcd calls as we are clearly repeating a lot of them