Project Euler 86: Exploring the shortest path from one corner of a cuboid to another.

Project Euler 86: Exploring the shortest path from one corner of a cuboid to another.

Written by on 3 March 2012

Topics: Project Euler

Problem 86 of Project Euler is more of a geometric problem than most of what I have seen here. It requires one insight and then it is fairly straight forward to solve. An insight I am just about  to give you, but let’s start with the problem formulation

A spider, S, sits in one corner of a cuboid room, measuring 6 by 5 by 3, and a fly, F, sits in the opposite corner. By travelling on the surfaces of the room the shortest “straight line” distance from S to F is 10 and the path is shown on the diagram.

However, there are up to three “shortest” path candidates for any given cuboid and the shortest route is not always integer.

By considering all cuboid rooms with integer dimensions, up to a maximum size of M by M by M, there are exactly 2060 cuboids for which the shortest distance is integer when M=100, and this is the least value of M for which the number of solutions first exceeds two thousand; the number of solutions is 1975 when M=99.

Find the least value of M such that the number of solutions first exceeds one million.

The insight we need in order to proceed is how to calculate the shortest path. Once we figure that out the rest should be pretty forward.

Shortest path through a cuboid

For a cuboid with h ≤ w ≤ l we have that the shortest path s is given as . That was not immediately obvious for me when I saw the 3D version of the cuboid, but if you unfold it as the figure to the right. You can see the distance between the two opposite corners marked with a dot is the straight line.

This straight line s, is the hypotenuse in the right angle triangle with the two sides l and w+h which gives us a length  described above.

Finding all integer solutions

What I have done is to find all solutions for a given maximum side length such that l=M and then do this for all side lengths. We have two options, either we check all combinations of 1≤h≤w≤M w and see if it gives us an integer solution. Alternatively we check all combinations of 1<w+h≤2M and then for the integer solutions we figure out how many combinations of w and h that gives us.

You probably figured that I would use the latter as I find it more interesting to investigate. I will start by denoting w+h as wh.Basically we have two cases, either wh ≤ M or wh > M.

In the first case  we get a number of solutions such that h=1 and w=wh-1, h=2 and w = wh-2… and so on until we get that and . Or in other words we get  solutions to the problem.

In the case where wh > M we can only vary w between  and M both included.

Implementing the solution

Now that we have the basics of the solution sorted out we are ready to implement it in C#. I have coded the solution like this

int l = 2;
int count = 0;
int target = 1000000;

while (count < target) {
    l++;
    for (int wh = 3; wh <= 2 * l; wh++) {
        double sqrt = Math.Sqrt(wh * wh + l * l);                    
        if (sqrt == (int)(sqrt)) {
            count += (wh <= l) ? wh / 2 : 1 + (l - (wh+1)/2);
        }
    }
}

This piece of code gives us

M=1818 is the first solution with more than 1.000.000 solutions
Solution took 23,9568 ms

At first I use the square root function twice. After I removed that I almost halved the execution time, so it is pretty obvious that this is what takes the majority of the run time. I don’t think the solution is fast, but I like the simplicity of the coding.

Wrapping up

A fun geometry exercise which was also very much fun to write about since it is easy to illustrate the concepts. I hope you enjoyed the problem as well. If you solved it in another way, or derived the same solution in another way, please tell me about it.

The source code for the problem can be found here.

The fly/spider image used as the featured image today is taken by abznak.

18 Comments For This Post I'd Love to Hear Yours!

  1. Mike says:

    I solved the problem by noticing that (l,wh,s) were Pythagorean triplets. Using the algorithm I described in problem 75, to generate Pythagorean triplets, it was a simple matter to iterate through the triplets.

  2. Kristian says:

    Hi Mike

    That is almost cheating isn’t it? 🙂 So you generated the triplets and then counted how many solutions each triplet would result in, or am I wrong?

    Because it sounds like a more efficient solution than mine.

    /Kristian

  3. Mike says:

    How is it cheating? And yes for a given triplet I did calculate how many ways there were to get wh.

  4. Kristian says:

    Hi Mike

    It is by no means cheating. I always says that with a smile on my face, when someone finds a brilliant shortcut that I missed. So basically all I said is “well done”.

    /Kristian

  5. bill says:

    i cant believe you have do many blogs about project euler. its enough to keep me busy for months. Do u do one every couple of days?

  6. Kristian says:

    I try to do one per week at least.

  7. Rook says:

    Am I stupid …
    I tried both brute force with all combinations and with pythagorean triples and I can`t even get the result from the examples.
    E.g. on C# brute force
    int M = 100;
    int count = 0;
    for (int i = 1; i <= M; i++)
    {
    for (int j = 1; j <= M; j++)
    {
    for (int k = j; k 6^2 + (5+3)^2 = 10^2
    2: 6,7,1 -> 6^2 + (7+1)^2 = 10^2
    3: 6,6,2 -> 6^2 + (6+2)^2 = 10^2
    4: 6,4,4 -> 6^2 + (4+4)^2 = 10^2
    (then the 6 is divided into smaller sides and 5+3=8 is untouched)
    5: 5,1,8 -> 8^2 + (5+1)^2 = 10^2
    6: 4,2,8 -> 8^2 + (2+4)^2 = 10^2
    7: 3,3,8 -> 8^2 + (3+3)^2 = 10^2

    => Which means there are 7 solutions for this case. If calculate this way for all the pythagorean triples I get an answer about double above 2060 …
    Now I am buffed o.0

  8. Kristian says:

    First of all I think your comment got screwed up by the blog. Try enclosing your code in code tags as given below the comment area.

    However, what I think I can see is that you start both i and j from 1, but you should start j from i instead. So you always have that the length is larger than or equal to the width. That should cut the solutions in half.

  9. Gyosko says:

    Hi!First of all great blog!;-)

    I just can’t understand one thing reading this entry anyway:
    Why we have to look for solutions in this range 1

  10. Kristian says:

    That is because both w and h must be smaller than l, but each can be larger than half of l without violating anything.

  11. ozth says:

    Why is floor(wh/2) the number of solutions ?

    Cheers

  12. Kristian says:

    Because each value of h from 1 to floor(wh/2) is a solution, and that means we get floor(wh/2) solutions when you count them.

  13. jeras says:

    I’m surprised that you wouldn’t use pythagorean triples to solve the problem, since the only integral shortest path candidates will come from these triples. The only difficulty then is producing the triples in increasing order of M, where M is the maximum room dimension. For that I used a priority queue (heapq in python).

    It’s important to note that each triple may produce 1 or 2 sets of solutions. For a triple (a, b, c) with a < b < c, you get solutions for M=a if b <= 2a, as well as solutions for M=b.

    For example, the triple (3, 4, 5) produces 2 solutions for M=3 and w+h=4, and 1 solution for M=4 and w+h=3. The triple (5, 12, 13), however, only produces solutions for M=12 and w+h=5, since w+h can never reach 12 if M=5.

    If I translate your solution directly to python, it takes a little over 2.5 seconds. The solution on Dreamshire, which is basically an optimized version of your solution, takes about 1 second. My python solution below takes about 17 ms, so it's orders of magnitude faster.

    #! /usr/bin/env python

    from heapq import heappop, heappush

    def pythagorean_triples():
    def push_triple(a, b, c):
    a, b, c = sorted([a, b, c])
    if b &lt;= 2*a:
    heappush(triples, (a, (a, b, c, 1)))
    heappush(triples, (b, (b, a, c, 1)))

    triples = []
    push_triple(3, 4, 5)
    while True:
    M, (a, b, c, k) = heappop(triples)
    heappush(triples, (M+a, (a, b, c, k+1)))
    if k == 1 and a &gt; b:
    push_triple( a-2*b+2*c, 2*a-b+2*c, 2*a-2*b+3*c)
    push_triple( a+2*b+2*c, 2*a+b+2*c, 2*a+2*b+3*c)
    push_triple(-a+2*b+2*c, -2*a+b+2*c, -2*a+2*b+3*c)
    yield M, (a, b, c, k)

    def solve(target=10**6):
    count = 0
    for M, (a, b, c, k) in pythagorean_triples():
    ka, kb = k*a, k*b
    count += min(ka+1, kb) – (kb+1)//2
    if count &gt; target:
    break
    return M

    if __name__ == '__main__':
    print solve()

  14. jeras says:

    I guess I messed up the code tags. Trying again:

    #! /usr/bin/env python
    
    from heapq import heappop, heappush
    
    def pythagorean_triples():
        def push_triple(a, b, c):
            a, b, c = sorted([a, b, c])
            if b &lt;= 2*a:
                heappush(triples, (a, (a, b, c, 1)))
            heappush(triples, (b, (b, a, c, 1)))
    
        triples = []
        push_triple(3, 4, 5)
        while True:
            M, (a, b, c, k) = heappop(triples)
            heappush(triples, (M+a, (a, b, c, k+1)))
            if k == 1 and a &gt; b:
                push_triple( a-2*b+2*c,  2*a-b+2*c,  2*a-2*b+3*c)
                push_triple( a+2*b+2*c,  2*a+b+2*c,  2*a+2*b+3*c)
                push_triple(-a+2*b+2*c, -2*a+b+2*c, -2*a+2*b+3*c)
            yield M, (a, b, c, k)
    
    def solve(target=10**6):
        count = 0
        for M, (a, b, c, k) in pythagorean_triples():
            ka, kb = k*a, k*b
            count += min(ka+1, kb) - (kb+1)//2
            if count &gt; target:
                break
        return M
    
    if __name__ == '__main__':
        # print solve(2000)
        print solve()
    
  15. Rob says:

    If you use binary search, your solution is very fast (this is using C++ but it shouldn’t be much different for C#):

    int hi = 4000;
    int lo = 1;
    int mid;

    while (lo &lt;= hi) {
    mid = (hi + lo) / 2;

    int count = cuboid_count(mid);

    if (count &lt; 1000000) {
    lo = mid + 1;
    } else {
    hi = mid – 1;
    }
    }
    cout &lt;&lt; cuboid_count(mid) &lt;&lt; endl;

    cout &lt;&lt;&quot;M = &quot; &lt;&lt; mid &lt;&lt; endl;

    Running it:

    $ time ./a.out
    1000457
    M = 1818
    
    real	0m1.026s
    user	0m1.023s
    sys	0m0.003s
    

    This assumes that cuboid_count() is monotonically increasing and I don’t know if it is or not, but this appears to be getting the right answer.

  16. Erebos1337 says:

    Hey Kristian,
    I frequently visit your blog after I’ve solved a problem myself.

    We had a very similar approach and I had the same thoughts about finding the shortest path.

    But I have been missing something in your explanation and thus also in your image.

    There are three possible shortest ways as seen in this picture .

    If you assume that h ≤ w ≤ l, then you are right, that the one, you have drawn is the shortest one. You can prove that it indeed is and this I did.
    I would have liked at least a comment about this, since it seemed like you were ignoring the possibility that one of the other two could be shorter at all, maybe you just knew that from before, though the reader might not.

    About the code I had the same thoughts, that calculating the square root would require much processing time.
    I tried a few things to calculate the square root faster like approximation via newton and such,
    but they didn’t help, it made it even worse in runtime or the solution was wrong.

    But there is one thing that I found interesting:
    The wh in your code (and mine also) are monotonous increasing, thus the length of the hypotenuse is also monotonous increasing, also it is always bigger than l.

    If I want to know if the hypotenuse is an integer, I just need to know if the squared hypotenuse is a square number.

    So what I did is introduce a “pointer” p to the root of l^2 = l.
    Instead of taking the square root of h^2=l*l+wh*wh, I produce increasing square numbers until I reach h^2. If I pass by, then I know h^2 is not a square number and hence h is not an integer. If I find a square number p*p=h^2, then I know, that h is integer and also, I have found the root, which is p (in my code I neeed the root, so huge gain there).

    The best part about this approach is, that since wh is increasing inside the loop, h^2 = l*l+wh*wh is also increasing, so we can reuse the pointer.
    So, for each new wh, we only need to check one or two square numbers.

    To have a comparison I executed your code on my computer (which seems to be way slower) and your original code gave this solution and runtime:

    solution: M=1818: 1000457
    duration: 78 ms

    Then I replaced the part with the square root with my square pointer algorithm and ran the new code with following result:

    solution: M=1818: 1000457
    duration: 21 ms

    So, yes, you were right, calculating the square root takes an immense amount of time.

    Here is your algorithm with the square pointer algorithm:

    int l = 2;
    int count = 0;
    int target = 1000000;

    while (count < target) {
    l++;

    // set pointer for squares to position of l
    // since the length of the hypotenuse is always longer than l
    // and monotonous growing with increasing wh
    int p = l;
    int l2 = l*l;

    for (int wh = 3; wh <= 2 * l; wh++) {

    // save squared length of hypotenuse for convenience
    int h2 = wh*wh + l2;

    // increase pointer until we’ve either found h2 or passed h2
    while(p*p < h2) p++;

    // if we found h2, then we know h2 is a square number, hence h=sqrt(h2) is an integer
    // otherwise p*p holds the first square number bigger than h^2
    // one more advantage: if we needed the square root of h2, we would have found it with p
    if (p*p == h2) {
    count += (wh <= l) ? wh / 2 : 1 + (l – (wh+1)/2);
    // increment p since we’ve found p*p=h2 and the next h2 will be bigger
    p++;
    }
    }
    }

    My approach btw was very similar, except that I first produced the h^2 in such a way that they had to be square numbers. After that I checked if h^2 – l^2 = wh^2 was an integer. Here it was not sufficient to know that wh was an integer, since the value of wh is needed to calculate how many distinct sums w+h=wh exist, luckily my approach gave me the square root (wh) for free 🙂
    My approach was even 3ms faster than the optimized version of yours,
    because there were fewer generated numbers and checks needed,
    but that is only a little gap compared to the other thing.

    Keep the good work up.
    I really like your blog and will continue to visit it always after I’ve finished a new problem.

    Best Regards,
    Thomas

  17. Erebos1337 says:

    Hey Kristian,
    I frequently visit your blog after I’ve solved a problem myself.

    We had a very similar approach and I had the same thoughts about finding the shortest path.

    But I have been missing something in your explanation and thus also in your image.

    There are three possible shortest ways as seen in this picture.

    If you assume that h ≤ w ≤ l, then you are right, that the one, you have drawn is the shortest one. You can prove that it indeed is and this I did.
    I would have liked at least a comment about this, since it seemed like you were ignoring the possibility that one of the other two could be shorter at all, maybe you just knew that from before, though the reader might not.

    About the code I had the same thoughts, that calculating the square root would require much processing time.
    I tried a few things to calculate the square root faster like approximation via newton and such,
    but they didn’t help, it made it even worse in runtime or the solution was wrong.

    But there is one thing that I found interesting:
    The wh in your code (and mine also) are monotonous increasing, thus the length of the hypotenuse is also monotonous increasing, also it is always bigger than l.

    If I want to know if the hypotenuse is an integer, I just need to know if the squared hypotenuse is a square number.

    So what I did is introduce a “pointer” p to the root of l^2 = l.
    Instead of taking the square root of h^2=l*l+wh*wh, I produce increasing square numbers until I reach h^2. If I pass by, then I know h^2 is not a square number and hence h is not an integer. If I find a square number p*p=h^2, then I know, that h is integer and also, I have found the root, which is p (in my code I neeed the root, so huge gain there).

    The best part about this approach is, that since wh is increasing inside the loop, h^2 = l*l+wh*wh is also increasing, so we can reuse the pointer.
    So, for each new wh, we only need to check one or two square numbers.

    To have a comparison I executed your code on my computer (which seems to be way slower) and your original code gave this solution and runtime:

    solution: M=1818: 1000457
    duration: 78 ms

    Then I replaced the part with the square root with my square pointer algorithm and ran the new code with following result:

    solution: M=1818: 1000457
    duration: 21 ms

    So, yes, you were right, calculating the square root takes an immense amount of time.

    Here is your algorithm with the square pointer algorithm:

    int l = 2;
    int count = 0;
    int target = 1000000;
     
    while (count < target) {
        l++;
        
        // set pointer for squares to position of l
        // since the length of the hypotenuse is always longer than l
        // and monotonous growing with increasing wh
        int p = l;
        int l2 = l*l;
        
        for (int wh = 3; wh <= 2 * l; wh++) {
        	
        	// save squared length of hypotenuse for convenience
        	int h2 = wh*wh + l2;
        	
        	// increase pointer until we've either found h2 or passed h2
        	while(p*p < h2) p++;  
        	
        	// if we found h2, then we know h2 is a square number, hence h=sqrt(h2) is an integer
        	// otherwise p*p holds the first square number bigger than h^2
        	// one more advantage: if we needed the square root of h2, we would have found it with p
            if (p*p == h2) {
                count += (wh <= l) ? wh / 2 : 1 + (l - (wh+1)/2);
                // increment p since we've found p*p=h2 and the next h2 will be bigger
                p++;
            }
        }
    }
    

    My approach btw was very similar, except that I first produced the h^2 in such a way that they had to be square numbers. After that I checked if h^2 – l^2 = wh^2 was an integer. Here it was not sufficient to know that wh was an integer, since the value of wh is needed to calculate how many distinct sums w+h=wh exist, luckily my approach gave me the square root (wh) for free 🙂
    My approach was even 3ms faster than the optimized version of yours,
    because there were fewer generated numbers and checks needed,
    but that is only a little gap compared to the other thing.

    Keep the good work up.
    I really like your blog and will continue to visit it always after I’ve finished a new problem.

    Best Regards,
    Thomas

  18. zuggy says:

    If one doesn’t spot the geometric insight, one can get the same result by studying the derivative of the function f, f:[0, c] -> R+, c the length of one of the sides e.g. the longest, x the ‘crossing point’ on that side, f(x) the distance travelled.

Trackbacks For This Post

  1. M : Santa at cuboid world « A Modern Prometheus

Leave a Comment Here's Your Chance to Be Heard!

You can use these html tags: <a href="" title=""> <abbr title=""> <acronym title=""> <b> <blockquote cite=""> <cite> <code> <del datetime=""> <em> <i> <q cite=""> <s> <strike> <strong>

You can use short tags like: [code language="language"][/code]

The code tag supports the following languages: as3, bash, coldfusion, c, cpp, csharp, css, delphi, diff,erlang, groovy, javascript, java, javafx, perl, php, plain, powershell, python, ruby, scala, sql, tex, vb, xml

You can use [latex]your formula[/latex] if you want to format something using latex

This site uses cookies.