About once a year Advent of Code give you a problem that is a gift if you are using a computer algebra system. This year that day was Day 13. The day 13 problem was one of crazy claw machines. Each machine has two buttons that can be pressed to move the claw a given number of X and Y positions and a prize at a given position. A buttons cost 3 tokens to press, and B buttons cost 1 token to press, and we are asked to find the minimum number of tokens needed, IF it is possible to reach the prize. The data presented like so:
Button A: X+94, Y+34
Button B: X+22, Y+67
Prize: X=8400, Y=5400
Button A: X+26, Y+66
Button B: X+67, Y+21
Prize: X=12748, Y=12176
...
Some times the input is harder to parse than the problem is so solve, and this might be one of those cases. I tend to reach for StringTools to take the input apart, but today the right tool is the old school C-style sscanf (after using StringSplit to split at every double linebreak).
machinesL := StringTools:-StringSplit(Trim(input), "\n\n");
machines := map(m->sscanf(m,"Button A: X+%d, Y+%d\nButton B: X+%d, Y+%d\nPrize: X=%d, Y=%d"), machinesL):
Now we have a list of claw machine parameters in the form [A_x, A_y, B_x, B_y, P_x, P_y] and we need to turn those into equations that we can solve. We want the number of A presses a, and B presses b to get the claw to the P_x, P_y position of the claw, it is simple to just write them down:
for m in machines do
eqn := ({m[1]*a+m[3]*b=m[5], m[2]*a+m[4]*b=m[6]});
end do;
Now because of the discrete nature of this problem, we need our variables a and b to be non-negative integers. When solving this, I first reached for isolve like this:
tokens := 0;
for m in machines do
eqn := ({m[1]*a+m[3]*b=m[5], m[2]*a+m[4]*b=m[6]});
sol := isolve(eqn);
if sol <> NULL then
tokens := tokens + eval(3*a+b, sol);
end if;
end do;
Now, sometimes Advent of Code inputs contain a lot of hidden structure. I wrote the code above, it worked on the sample input, so I tried it immediately on my real input (about 300 claw machines like the above) and IT WORKED. But, you might notice that this code does not deal with a couple cases that could have appeared. In particular, it doesn't check that the solutions are positive. It also doesn't handle cases where there is more than one possible solution. The former is easy to check
if sol <> NULL and eval(a,sol) >= 0 and eval(b,sol) >= 0 then
Unfortunately isolve does not handle inequalities, but you could try with solve, but it doesn't save us any checking, because we'd still have to check if the solutions are integers, so we might as well have just solved the equation and then checked if it were a nonnegative integer.
tokens := 0;
for m in machines do
eqn := {m[1]*a+m[3]*b=m[5], m[2]*a+m[4]*b=m[6], a>=0, b>=0};
sol := solve(eqn);
if type(eval(a,sol), integer) and type(eval(b,sol),integer) then
tokens := tokens + eval(3*a+b, sol);
end if;
end do;
ans1 = tokens;
In the multiple solution case we get something like {a = 3 - 2*b, 0 <= b, b <= 3/2} which has some great information in it but might be hard to handle programmatically, so let's see what isolve does with those cases to see if it's easier to deal with
> eqn := { 17*a + 84*b = 7870, 34*a + 168*b = 15740 }:
> constr := { a >= 0, b>= 0 }:
> sol := isolve(eqn);
sol := {a = 458 - 84 _Z1, b = 1 + 17 _Z1}
> constr := eval({ a >= 0, b>= 0 }, sol):
const := {0 <= 1 + 17 _Z1, 0 <= 458 - 84 _Z1}
> obj := eval(3*a+b, sol);
obj := 1375 - 235 _Z1
You can see it's easy to tell if these show up in your input, since your "token" total will have the _Zn variables in it. Now, since everything is simple and linear here, it seems like you could use solve to find the rational value of _Z1 that makes obj=0 and then take the closest integer but it's not so simple, we actually have to deal with the contraints that a and b be positive too. So, it really just makes sense to bring out the big hammer of Optimization:-Minimize which allows us to directly optimize over just the integers. So a full solution looks like this:
tokens := 0:
for m in machines do
eqn := ({m[1]*a+m[3]*b=m[5], m[2]*a+m[4]*b=m[6]});
sol := isolve(eqn);
if sol = NULL then
next;
end if;
constr := eval({ a >= 0, b>= 0 }, sol);
obj := eval(3*a+b, sol);
if not type(obj, constant) then
tokens := tokens + Optimization:-Minimize( obj, constr, assume=integer )[1];
elif andmap(evalb, constr) then
tokens := tokens + obj;
end if;
end do;
But since we're bringing out the big hammer, why not just use Optimization in the first place. The main reason is that Minimize doesn't simply return NULL when it doesn't work, instead it throws an exception, so we need to find all the exceptions that can occur and handle then with a try-catch, thus:
tokens := 0;
for m in machines do
eqn := ({m[1]*a+m[3]*b=m[5], m[2]*a+m[4]*b=m[6]});
try
sol := Optimization:-Minimize(3*a+b, eqn, 'assume'='nonnegint')[1];
tokens := tokens + sol;
catch "no feasible":
end try;
end do;
tokens;
(you can in fact omit the string in the catch: statement, but I can tell you from long experience that that is an excellent way to make your code much much harder to debug)
Alright, so how did people not using Maple solve this problem? The easiest way to solve it, and the one used by all the cheaters scraping the website and using LLM-based code generators that auto-submit solutions to get into the Top 100, was to just check all possible a, b values in 0..100 and take the values than minimize 3*a+b when reaching the prize coordinates. That's only feasible because the problem states the 100 is an upperbound for a and b, but it's also very fast (about 1/10 second in Maple):
tokens := 0:
for m in machines do;
sol := infinity;
for i from 0 to 100 do for j from 0 to 100 do
if i*m[1]+j*m[3]=m[5] and i*m[2]+j*m[4]=m[6] and 3*i+j < sol then
sol := 3*i+j;
end if;
end do; end do;
tokens := tokens + ifelse(sol=infinity,0,sol);
end do;
It does not scale at all to part 2 (which modified everything to be bigger by about 10 trillion), and it seems that foiled all the LLM solvers. So, what solutions scaled in languages without integer equation solvers? Well, the easiest solution is just to solve the general equation using paper and pencil
And you can just hard code that formula in, check that it gives integer values and compute the tokens. As long as you get unique solutions, that looks something like this
solveit := proc(m)
local asol := m[4]*m[5] - m[3]*m[6];
local bsol := m[1]*m[6] - m[2]*m[5];
local deno := m[1]*m[4] - m[2]*m[3];
if deno = 0 then return -2^63; end if; # multiple solution case - not handled
if asol mod deno = 0 and bsol mod deno = 0
and ( ( deno>=0 and asol>=0 and bsol>=0 )
or ( deno<=0 and asol<=0 and bsol<=0 ) )
then
return 3*asol/deno+bsol/deno;
else
return 0;
end if;
end proc:
Which if you have this in Maple, you can impress your friends by auto generating solutions in other languages. Here, for your FORTRAN friends
> CodeGeneration:-Fortran(solveit);
Warning, the following variable name replacements were made: solveit -> cg
integer function cg (m)
doubleprecision m(*)
integer asol
integer bsol
integer deno
asol = int(-m(3) * m(6) + m(4) * m(5))
bsol = int(m(1) * m(6) - m(2) * m(5))
deno = int(m(1) * m(4) - m(2) * m(3))
if (deno .eq. 0) then
cg = -9223372036854775808
return
end if
if (mod(asol, deno) .eq. 0 .and. mod(bsol, deno) .eq. 0 .and. (0
# .le. deno .and. 0 .le. asol .and. 0 .le. bsol .or. deno .le. 0 .a
#nd. asol .le. 0 .and. bsol .le. 0)) then
cg = 3 * asol / deno + bsol / deno
return
else
cg = 0
return
end if
end
Another way that you might solve this without solve is to use a linear algebra library to solve the linear system. It works even if you only have a numeric solver, but you have to be careful about checking for integers:
tokens := 0:
for m in machines do
sol := LinearAlgebra:-LinearSolve(
Matrix(1..2,1..2,[m[[1,3]],m[[2,4]]], datatype=float),
Vector(m[5..6], datatype=float));
if abs(sol[1]-round(sol[1])) < 10^(-8) and abs(sol[2]-round(sol[2])) < 10^(-8)
and sol[1] >= 0 and sol[2] >= 0
then
tokens := tokens + 3*sol[1]+sol[2];
end if;
end do;
Finally, a lot of people solved this sort of thing with the Z3 Theorem prover from Microsoft research which is also way more than you need, but it mostly just uses SMTLIB, which we also have in a library for in Maple, and it can just be used in place of solve
tokens := 0;
for m in machines do
eqn := {m[1]*a+m[3]*b=m[5], m[2]*a+m[4]*b=m[6], a>=0, b>=0};
sol := SMTLIB:-Satisfy(eqn) assuming a::nonnegint, b::nonnegint;
if sol <> NULL and type(eval(a,sol), integer) and type(eval(b,sol),integer) then
tokens := tokens + eval(3*a+b, sol);
end if;
end do;
Notice that Satisfy handled the multiple solution case just by choosing one of the many solutions. It is possible to get SMTLib to optimize but it is slightly more involved, and this post is already too long. This time, I've put all this work in worksheet: Day13-Primes.mw
