Dazbo's Advent of Code solutions, written in Python

Work in progress…

SymPy is so cool! I came across when I was trying to solve some equations as part of 2023 AoC puzzles. It is computer alegebra systems (CAS): a library that allows us to perfrom algebra and mathematical computations.

You can use it for things like:

- Simplifying, rearranging, factoring and expanding algebraic expressions and equations.
- Solving equations of various types, including linear, polynomial, and differential.
- Performing calculus, including differentiation and integration.
- Plotting graphs.
- Generating LaTeX output.

We need to define variables before we can use them in SymPy. We do this by creating `symbols`

.

One interesting observation here is that when you evaluate a SymPy expression a Jupyter notebook, it is rendered in *pretty-printed* (mathematical) format. But when we print using `print()`

, it is printed in conventional Python style.

If you want to explicitly pretty-print from a Python notebook, you can use `IPython.display`

.

Alternatively, you can optionally run `init_printing()`

to enable the best printer for your environment. Or run `init_session()`

to automatically import everything in SymPy, create some common Symbols, setup plotting, and run `init_printing()`

.

```
import sympy
from sympy import latex
from IPython.display import display, Markdown
sympy.init_session() # set up default behavour and printing - not essential
a, b, x, y = sympy.symbols("a b x y") # define symbols (variables)
a = (x + 1)**2 # set a to be an expresssion
b = x**2 + 2*x + 1 # set b to be an equivalent expression
display(Markdown(f"$a = {latex(a)}$")) # pretty-print a
print(f"{latex('a = ')}{latex(a)}") # If we just want the raw latex
display(Markdown(f"$b = {latex(b)}$")) # pretty-print b
# Expand a...
display(Markdown("Expanding $a$..."))
expanded = sympy.expand(a)
display(Markdown(f"$a = {latex(expanded)}$"))
# Show a-b without simplifying...
a_minus_b = a-b
display(Markdown(f"$a-b = {latex(a_minus_b)}$"))
# Simplify...
display(Markdown("Simplifying..."))
display(Markdown(f"$a-b$ $= {latex(sympy.simplify(a_minus_b))}$"))
```

We can also test if two expressions are the same, e.g.

```
a.equals(b)
```

You can’t simply assign values to the Python variable. If you did so, they would cease to be a symbol. Instead, you need to use `subs()`

to associate a value with a SymPy symbol.

```
x, y = sympy.symbols("x y")
expr = x+y
result = expr.subs({x: 2, y: 5})
print(f"{expr=}, {result=}")
```

```
import sympy
x = sympy.symbols("x")
expr = x**2
deriv = sympy.diff(expr)
print(deriv)
integral = sympy.integrate(deriv)
print(integral)
```

Here, I use `solve()`

to determine the values of $x$. I’ve specified `dict=True`

such that $x$ can be retrieved by passing its symbol as a dictionary key to the returned `solutions`

dictionary.

```
expr = x**2 + 3*x - 10
solutions = sympy.solve(expr, x, dict=True)
solutions
for solution in solutions:
print(f"x={solution[x]}")
```

Note that the expressions still need to be written in valid Python. So we can write `3*x`

, but we can’t write `3x`

.

What if you want to ignore *complex* solutions and only include *real* solutions? In this case, you can use `solveset()`

, and specify the domain as `domain=sympy.S.Reals`

. (Whereas `domain=S.Complexes`

is the default.)

```
expressions = []
expressions.append(x**2 - 9) # there are only real solutions
expressions.append(x**2 + 9) # there are only complex solutions
for expr in expressions:
display(expr)
solutions = sympy.solveset(expr, x)
print(f"There are {len(solutions)} solutions for x...")
display(solutions)
display(expressions[-1])
solutions = sympy.solveset(expressions[-1], x, domain=sympy.S.Reals)
print(f"There are {len(solutions)} real solutions for x.")
```

We can evaluate to obtain the float value of a symbol, and display to an arbitrary level of precision:

```
expr = sympy.sqrt(8)
display(expr)
display(expr.evalf(4)) # to 4 digits of precisions
```

In this example taken from 2023 Day 6, I’m solving a quadratic:

\[h^2 - th + d = 0\]```
def solve_part2_sympy_quadratic(data):
""" h^2 - th + d = 0 """
race_duration = int("".join(x for x in data[0].split(":")[1].split()))
distance = int("".join(x for x in data[1].split(":")[1].split()))
logger.debug(f"{race_duration=}, {distance=}")
# solve using quadratic with SymPy
h = sympy.symbols("h", real=True)
equation = sympy.Eq(h**2 - race_duration*h + distance, 0)
solutions = sympy.solve(equation, h, dict=True)
answers = [solution[h].evalf() for solution in solutions] # there should be two
```

In this example taken from 2023 Day 24, I have a series of equations that are necessary to find several unknown variables. Specifically:

\[\begin{align} t &= \frac{x_{r} - x_{h}}{v_{x_{h}} - v_{x_{r}}} = \frac{y_{r} - y_{h}}{v_{y_{h}} - v_{y_{r}}} = \frac{z_{r} - z_{h}}{v_{z_{h}} - v_{z_{r}}} \\ \notag \\ (x_{r} - x_{h})(v_{y_{h}} - v_{y_{r}}) &= (y_{r} - y_{h})(v_{x_{h}} - v_{x_{r}}) \\ (y_{r} - y_{h})(v_{z_{h}} - v_{z_{r}}) &= (z_{r} - z_{h})(v_{y_{h}} - v_{y_{r}}) \\ (z_{r} - z_{h})(v_{x_{h}} - v_{x_{r}}) &= (x_{r} - x_{h})(v_{z_{h}} - v_{z_{r}}) \\ \end{align}\]And here’s the code:

```
def solve_part2(data: list[str]):
"""
Determine the sum of the rock's (x,y,z) coordinate at t=0, for a rock that will hit every hailstone
in our input data. The rock has constant velocity and is not affected by collisions.
"""
stones = parse_stones(data)
logger.debug(f"We have {len(stones)} stones.")
# define SymPy rock symbols - these are our unknowns representing:
# initial rock location (xr, yr, zr)
# rock velocity (vxr, vyr, vzr)
xr, yr, zr, vxr, vyr, vzr = sympy.symbols("xr yr zr vxr vyr vzr")
equations = [] # we assemble a set of equations that must be true
for stone in stones[:10]: # we don't need ALL the stones to find a solution. We need just enough.
x, y, z = stone.posn
vx, vy, vz = stone.velocity
equations.append(sympy.Eq((xr-x)*(vy-vyr), (yr-y)*(vx-vxr)))
equations.append(sympy.Eq((yr-y)*(vz-vzr), (zr-z)*(vy-vyr)))
solutions = sympy.solve(equations, dict=True) # SymPy does the hard work
if solutions:
solution = solutions[0]
logger.info(solution)
return sum([solution[xr], solution[yr], solution[zr]])
logger.info("No solutions found.")
return None
```