How the algorithm works
We saw on the 2D exploration page that the smaller sums lie close to the contours where z-x is a constant. The contour page implements an algorithm to follow these contours. We start where the line V(x,y)=x^3+y^3-(x+n)^3=0 intersects the line x=y, which gives us 2x^3=(x+n)^3, or x=n/(2^[1/3]-1) (approx n/0.26), which we round up to the next integer. Once we have a point on the line we can look for the next point by checking the adjacent points for which has the smallest V. Because we have chosen points in the domain x>y>0 we only have to check two adjacent points, (x+1,y) and (x+1,y+1).
To do this we split the function V(x,y) into two components, V(x,y)=Vy(y)-Vx(x), where Vx(x)=3nx^2+3n^2x+n^3 and Vy(y)=y^3. From these we can derive a set deltas (similar to derivatives):
- ΔF(x)=F(x+1)-F(x)
- Δ^2F(x)=ΔF(x+1)-ΔF(x)
- Δ^3F(x)=Δ^2F(x+1)-Δ^2F(x)
- Δ^4F(x)=Δ^3F(x+1)-Δ^3F(x)
- Etc...
When we want to calculate the next member in the sequence we can work these equations backwards. As stated above there will be an infinite sequence of these, but for polynomials we eventually get to Δ^nF(x)=constant. For cubic equations n will be 3. If we try an example, F(x)=x^3:
| x | 1 | 2 | 3 | 4 | 5 |
|---|---|---|---|---|---|
| F(x) | 1 | 8 | 27 | 64 | 125 |
| ΔF(x) | 7 | 19 | 37 | 61 | |
| Δ^2F(x) | 12 | 18 | 24 | ||
| Δ^3F(x) | 6 | 6 | |||
So to calculate 6^3 we perform the following steps:
- Δ^2F(x+1)=Δ^2F(x)+6 (24+6=30)
- ΔF(x+1)=ΔF(x)+Δ^2F(x+1) (61+30=91)
- F(x+1)=F(x)+ΔF(x+1) (125+91=216)
We can preform a similar option to move to smaller x but we need to update the deltas in the reverse order.
Vy(y)
Vy(y) is essentially the same as F(x)=x^3, which we have already looked at. We can intialise the deltas more generically, starting at any y values as:
| Vy | y^3 | y^3+3y^2+3y+1 | y^3+6y^2+12y+8 | y^3+9y^2+27y+27 | y^3+12y^2+48y+64 |
|---|---|---|---|---|---|
| ΔVy | 3y^2+3y+1 | 3y^2+9y+7 | 3y^2+15y+19 | 3y^2+21y+37 | |
| Δ^2Vy | 6y+6 | 6y+12 | 6y+18 | ||
| Δ^3Vy | 6 | 6 |
Vx(x)
Vx(x) is a little more interesting, V(x)=x^3-(x+n)^3=3nx^2+3n^2x+n^3. To make this a little clearer, both here and in the implementation, we can define: a=3n, b=3n^2 and c=n^3, and construct a delta table for Vx. This time we reach the constant term on row 3 (2a=6n).
| Vx | ax^2+bx+c | a(x+1)^2+b(x+1)+c | a(x+2)^2+b(x+2)+c | a(x+3)^2+b(x+3)+c |
|---|---|---|---|---|
| ΔVx | 2ax+a+b | 2ax+3a+b | 2ax+5a+b | |
| Δ^2Vx | 2a | 2a |
Contour following
We have chosen x>y>0 so that we are in an octant where, as we follow the contour, x and y are both increasing, with x increasing faster than y. If we start on a point close to the line there are only two potential next points to consider, (x+1,y) and (x+1,y+1). As V=Vy(y)-Vx(x) the values of these next points can be calculated using:
V(x+1,y)=V(x,y)-ΔVx and
V(x+1,y+1)=V(x+1,y)+ΔVy.
The next point on the line is whichever of these has the smaller value of |V|. This graph shows some points near the start of the n=2 contour.
You can see that because the gradient is less than one the line is formed of horizontal sequences of dots. Because these dots all have the same y value V is large and positive on the left and large and negative on the right. The smallest V values will be somewhere near the middle. As the value of X increases the horizontal legs become longer and it becomes feasible to skip a lot of intermediate values and just look at the two where V transits from positive to negative. This is the basis of the second and third scanning options. The nice thing about this is that because the number of points that can be skipped increases as x increases the algorithm can scan larger and larger x ranges in the same amount of time.