#86 - Cuboid route

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.

cubeImg

However, there are up to three “shortest” path candidates for any given cuboid and the shortest route doesn’t always have integer length.

It can be shown that there are exactly 2060 distinct cuboids, ignoring rotations, with integer dimensions, up to a maximum size of $M$ by $M$ by $M$, for which the shortest route has integer length when $M=100$. This is the least value of $M$ for which the number of solutions first exceeds two thousand; the number of solutions when $M=99$ is 1975.

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


Since this is an optimization problem, the quickest way to solve this is single-variabled calculus. Next, I also want to show the three possible candidates for the shortest paths that the problem talked about. Like above, the spider would have to go over an edge, and so, **each path corresponds to each edge spider goes over **. She can either go over the 6 edge, the 5 edge, or the 3 edge. I’ve displayed them below.

routeImg

Which edge corresponds to the actual shortest path from $S$ to $F$? Intuitively, it is the longest edge, because the spider is cutting across it through both legs (before going over the edge, and after). If the spider cross any of the shorter edges, she would have cover that longest edge in just of the two legs, which would lead a slightly longer distance.

I’ll assume a cuboid measuring $a\times b\times c$, with $a\geq b\geq c$. Our assumption is that the spider goes over edge $a$, and so the variable we would be optimizing is where she does it. Let’s call this point $x$. The picture below helps with visualizing.

xpointImg

Each of the two parts are hypotenuses of two right triangles. The triangle on the $ab$ face, has sides $x$ and $b$, and thus has a hyponetuse of $\sqrt{x^2+b^2}$. The other triangle has sides $c$ and $a-x$. Therefore, our distance function is \[f(x)=\sqrt{x^2+b^2}+\sqrt{(a-x)^2+c^2}\]

To minimize this, we take the derivative and set it equal to 0. This will tell us the crossing point, which we plug into $f(x)$ to find the distance. Lots of algebra and simplying is needed: \[\begin{aligned} \frac{df}{dx}&=\frac{1}{2}(x^2+b^2)^{-\frac{1}{2}}(2x)+ \frac{1}{2}(x^2-2ax+a^2+c^2)^{-\frac{1}{2}} (2x-2a) \\ 0 &= \frac{x}{\sqrt{x^2+b^2}}+\frac{x-a} {\sqrt{x^2-2ax+a^2+c^2}} \\ \frac{a-x}{\sqrt{x^2-2ax+a^2+c^2}} &= \frac{x}{\sqrt{x^2+b^2}} \\ (a-x)\sqrt{x^2+b^2} &= x\sqrt{x^2-2ax+a^2+c^2} \\ (a-x)^2(x^2+b^2) &= x^2(x^2-2ax+a^2+c^2) \\ (x^2-2ax+a^2)(x^2+b^2) &= x^4-2ax^3+a^2x^2+c^2x^2 \\ x^4+b^2x^2-2ax^3-2ab^2x+a^2x^2+a^2b^2 &= x^4-2ax^3+a^2x^2+c^2x^2 \\ (b^2-c^2)x^2-2ab^2x+a^2b^2 &= 0 \end{aligned}\]

The last line is a quadratic equation in $x$: \[\begin{aligned} x &= \frac{2ab^2\pm\sqrt{4a^2b^4-4(b^2-c^2)a^2b^2}} {2(b^2-c^2)} \\ &= \frac{2ab^2\pm\sqrt{4a^2b^4-4a^2b^4+4a^2b^2c^2}} {2b^2-2c^2} \\ &= \frac{2ab^2\pm\sqrt{4a^2b^2c^2}}{2b^2-2c^2} \\ &= \frac{2ab^2\pm 2abc}{2b^2-2c^2} \\ &= \frac{ab^2\pm abc}{b^2-c^2} \end{aligned}\]

We choose the negative sign, since $x$ has to be less than $a$. \[\frac{ab^2-abc}{b^2-c^2} = \frac{ab(b-c)}{(b-c)(b+c)}=\boxed{\frac{ab}{b+c}}\]

Now we have the minimizing point, which we plug into our function to find the distance. \[\begin{aligned} f\left(\frac{ab}{b+c}\right) &= \sqrt{\left(\frac{ab}{b+c}\right)^2 + b^2} + \sqrt{\left(a-\frac{ab}{b+c}\right)^2+c^2} \\ &= \sqrt{\frac{(ab)^2}{(b+c)^2}+\frac{b^2(b+c)^2}{(b+c)^2}} + \sqrt{\left(\frac{a(b+c)} {b+c}-\frac{ab}{b+c}\right)^2 + c^2} \\ &= \frac{b}{b+c}\sqrt{a^2+(b+c)^2} + \frac{c}{b+c}\sqrt{a^2+(b+c)^2} \\ &= \left(\frac{b}{b+c}+\frac{c}{b+c}\right) \sqrt{a^2+(b+c)^2} \\ &= \boxed{\sqrt{a^2+(b+c)^2}} \end{aligned}\]

An extra step you can do is assume $x$ is on the smaller sides $b$ or $c$, and convince yourself those distances are larger than this one.

Let the minimum distance be $D$. We want $D$ to be an integer, which menas the expression inside the square root needs to be a perfect square: \[a^2+(b+c)^2=D^2\]

This looks exactly like the Pythagorean formula, which means $\mathbf{a}$, $\mathbf{b+c}$, and $\mathbf{D}$ all have to form a Pythagorean triple. For example, having $a=6$, $b=5$, and $c=3$ corresponds to the triple ${6,8,10}$ and indeed 10 is the shortest path, as the example in the problem shows. We have looped through Pythagorean triples before in #39 - Integer right triangles. But there are some cases we must consider…

When we loop through triples $\alpha, \beta, \gamma$ (with $\alpha\leq\beta\leq\gamma$), the hypotenuse corresponds to $D$ directly. But we have choices for the other two sides. Either our largest side $a=\alpha$, or $a=\beta$, with differing consequences.

$a=\alpha$

If we set our largest side to be $\alpha$, then we must split up $\beta$ into two values $b$ and $c$ such that $a$ are greater than both. This will only work if $\beta \leq 2$. Assuming this is true, how many solutions are there for a fixed triple in this case?

Without any constraints, there are $\lfloor \beta/2\rfloor$ distinct pairs of sums. Of these, exactly $\beta-\alpha$ of them have a number which is greater than $\alpha$, so we need to subtract these out. However, we also have to add back in the sum with an operand that is equal to $\alpha$, since that is also allowed. Therefore, the total number of solutions in this case are $\boxed{\lfloor \beta/2\rfloor - (\beta-\alpha)+1}$.

$a=\beta$

Since $\beta\geq\alpha$, all partitions of $\alpha$ are valid solutions. Thus, the number of solutions in this case is $\lfloor \alpha/2\rfloor$.

Implementation

We can use our recursive function from Problem 39, but when do we stop the generation? Well, $a$ is our biggest side, and this can’t exceed $M$. Based on our two cases above, we stop when either $\alpha > M$, or when both $\beta >M$ and $\beta >2\alpha$, as this prevents us from partitioning $\beta$.

This method is quick for a specific value of $M$. However, the problem requires the first $M$ whose solutions exceed one million, which implies we need to compute solutions for all values until that point. This can get really slow. To combat this, I increment $M$ by 100 in my loop. Eventually when we overshoot the value, we can slowly decrement the value back down to find the exact $M$ we need.

# file: "problem086.py"
def genPythagoreanTriples(triple, maxSide, A, B, C):
    # Find min and second largest
    # value. Never the last value...
    if triple[0] == min(triple):
        minVal, secVal = triple[0], triple[1]
    else:
        minVal, secVal = triple[1], triple[0]
    # Base case...
    if min(triple) > maxSide or (secVal > maxSide and secVal > 2 * minVal):
        return
    yield np.sort(triple)
    # Multiply each matrix
    for matrix in [A, B, C]:
        for multTriple in genPythagoreanTriples(np.dot(matrix, triple), maxSide, A, B, C):
            yield multTriple
def calculateSols(triple, maxSide, A, B, C):
    s = 0
    for pythTriple in genPythagoreanTriples(triple, maxSide, A, B, C):
        print(pythTriple)
        # Check to see if the second value is
        # less than twice the lowest value.
        # Otherwise, we can't break up the second value
        # into two sides less than the lowest...
        if pythTriple[1] <= 2 * pythTriple[0]:
            # Loop through, keeping the minimum side
            # constant...
            # Number of solutions is secVal // 2 - (secVal - firstVal) + 1
            for i in range(1, maxSide // pythTriple[0] + 1):
                multTriple = i * pythTriple
                v1 = multTriple[1] // 2 - (multTriple[1] - multTriple[0]) + 1
                s += v1
        # Keep the second value constant,
        # breaking the lowest value into all
        # possible unique partitions...
        for i in range(1, maxSide // pythTriple[1] + 1):
            v2 = int((i * pythTriple)[0] // 2)
            s += v2
    return s
A = np.array([
    [1, -2, 2],
    [2, -1, 2],
    [2, -2, 3]
], dtype=object)
B = np.array([
    [1, 2, 2],
    [2, 1, 2],
    [2, 2, 3]
], dtype=object)
C = np.array([
    [-1, 2, 2],
    [-2, 1, 2],
    [-2, 2, 3]
], dtype=object)

triple = np.array([3,4,5], dtype=object)
maxSide = 100
maxSolutions = 1000000
numSolutions = 1
while numSolutions < maxSolutions:
    numSolutions = calculateSols(triple, maxSide, A, B, C)
    # Save time by increasing the side by 100
    maxSide += 100
# Now slowly decrease the maxSide until we're less
while numSolutions > maxSolutions:
    maxSide -= 1
    numSolutions = calculateSols(triple, maxSide, A, B, C)
# Add one because we went one
# over with the while loop.
print('M = {} with {} solutions.'.format(maxSide + 1, calculateSols(triple, maxSide + 1, A, B, C)))

Running it all together, we get an output of,

M = 1818 with 1000457 solutions.
3.3048678999766707 seconds.

Thus, the first side length which has the number of solutions exceed one million is 1818.