commit 843d9afdb48f1344bd1ab1473ee39be0cccb23ca
parent b776f967e46827d2936be82466ce7cdfd8cfc0d1
Author: Sebastiano Tronto <sebastiano.tronto@gmail.com>
Date: Wed, 26 May 2021 11:26:53 +0200
First commit
Diffstat:
| A | README.md | | | 4 | ++++ |
| A | billiard3d.py | | | 197 | +++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++ |
2 files changed, 201 insertions(+), 0 deletions(-)
diff --git a/README.md b/README.md
@@ -0,0 +1,4 @@
+See [Wikipedia: Arithmetic Billiards](https://en.wikipedia.org/wiki/Arithmetic_billiards).
+
+This Python script uses matplotlib to draw pictures for a 3d arithmetic
+billiard, either as a 2d projection or as an axonometry.
diff --git a/billiard3d.py b/billiard3d.py
@@ -0,0 +1,197 @@
+import matplotlib.pyplot as plt
+from mpl_toolkits.mplot3d import Axes3D
+from math import pi, sin, cos, degrees, gcd #, lcm (python3.9+ only)
+
+# Python 3.8 or lower
+def lcm(a, b):
+ return a*b // gcd(a,b)
+
+###############################################################################
+# This part is generic for any size
+###############################################################################
+def is_corner(n, sizes, point):
+ for i in range(n):
+ if point[i] % sizes[i] != 0:
+ return False
+ return True
+
+def is_bouncing(n, sizes, point):
+ for i in range(n):
+ if point[i] % sizes[i] == 0:
+ return True
+ return False
+
+def lcm_list(L):
+ ret = 1
+ for e in L:
+ ret = lcm(ret, e)
+ return ret
+
+# Generic n-dimensional path.
+def get_path(n, sizes, r):
+ if n <= 1 or len(sizes) != n or len(r) != n:
+ print("Invalid sizes or starting point:", n, sizes, r)
+ return
+
+ path = [r]
+ p = r
+ d = [1] * n
+ steps = 0
+ while steps < 2*lcm_list(sizes):
+ p = [p[i] + d[i] for i in range(n)]
+ path.append(p)
+ steps = steps+1
+ d = [d[i] if p[i] % sizes[i] != 0 else -d[i] for i in range(n)]
+ return path
+
+###############################################################################
+# This part is specific for 2d drawings
+###############################################################################
+class Transformation:
+ def __init__(self, angle=0.0, mirror=False, shift=(0,0)):
+ self.angle = angle
+ self.mirror = mirror
+ self.shift = shift
+
+ def rotate(self, point):
+ c, s = cos(self.angle), sin(self.angle)
+ return [c*point[0] - s*point[1], s*point[0] + c*point[1]]
+
+ def reflect(self, point):
+ return [point[0], -point[1] if self.mirror else point[1]]
+
+ def translate(self, point):
+ return [point[i] + self.shift[i] for i in range(2)]
+
+ def apply_to(self, point):
+ return self.translate(self.rotate(self.reflect(point)))
+
+def draw_line_2d(p1, p2, c="black", w=1):
+ plt.plot([p1[0], p2[0]], [p1[1], p2[1]], color=c, linewidth=w)
+
+def draw_grid_and_rectangle(sizes, t):
+ for i in range(0, sizes[0]+1):
+ p1 = t.apply_to([i,0])
+ p2 = t.apply_to([i,sizes[1]])
+ color, width = ("red", 1) if i % sizes[0] == 0 else ("grey", 0.5)
+ draw_line_2d(p1, p2, c=color, w=width)
+ for i in range(0, sizes[1]+1):
+ p1 = t.apply_to([0,i])
+ p2 = t.apply_to([sizes[0],i])
+ color, width = ("red", 1) if i % sizes[1] == 0 else ("grey", 0.5)
+ draw_line_2d(p1, p2, c=color, w=width)
+
+def draw_path_2d(sizes, r, transformation):
+ if len(sizes) != 2:
+ print("Cannot draw non-2d path")
+ return
+
+ plt.axis("off")
+ plt.axis("equal")
+ draw_grid_and_rectangle(sizes, transformation)
+
+ path = [transformation.apply_to(p) for p in get_path(2, sizes, r)]
+ plt.plot([p[0] for p in path], [p[1] for p in path], color="blue")
+
+###############################################################################
+# This part is specific for drawing the 2d projections of 3d billiards
+###############################################################################
+class Face:
+ def __init__(self, fixed_coordinate, value, transformation):
+ self.f = fixed_coordinate
+ self.v = value
+ self.t = transformation
+
+ def contains(self, point):
+ return point[self.f] == self.v
+
+ def proj(self, point):
+ return [point[i] for i in range(3) if i != self.f]
+
+def draw_bouncing_points_3d2d(sizes, r, face):
+ path_3d = get_path(3, sizes, r)
+ bp = [face.t.apply_to(face.proj(p))
+ for p in path_3d if is_bouncing(3, sizes, p) and face.contains(p)]
+
+ plt.scatter([p[0] for p in bp], [p[1] for p in bp], color="green")
+
+def draw_3d_projections(sizes, r):
+ a, b, c = sizes
+
+ bottom = Face(2, 0, Transformation())
+ top = Face(2, c, Transformation(shift=(a+c, -(b+c))))
+ front = Face(1, 0, Transformation(mirror=True))
+ back = Face(1, b, Transformation(shift=(0,b)))
+ left = Face(0, 0, Transformation(angle=pi/2))
+ right = Face(0, a, Transformation(angle=pi/2, mirror=True, shift=(a,0)))
+
+ for face in [bottom, top, front, back, left, right]:
+ draw_path_2d(face.proj(sizes), face.proj(r), face.t)
+ draw_bouncing_points_3d2d(sizes, r, face)
+ plt.show()
+
+###############################################################################
+# This part is for drawing 3d pictures
+###############################################################################
+def draw_line_3d(ax, p1, p2, c="black", w=1):
+ ax.plot([p1[0], p2[0]], [p1[1], p2[1]], [p1[2], p2[2]], color=c, lw=w)
+
+def draw_box_3d(ax, s):
+ color, width = ("red", 1)
+
+ draw_line_3d(ax, ( 0, 0, 0), (s[0], 0, 0), color, width)
+ draw_line_3d(ax, ( 0, s[1], 0), (s[0], s[1], 0), color, width)
+ draw_line_3d(ax, ( 0, 0, s[2]), (s[0], 0, s[2]), color, width)
+ draw_line_3d(ax, ( 0, s[1], s[2]), (s[0], s[1], s[2]), color, width)
+
+ draw_line_3d(ax, ( 0, 0, 0), ( 0, s[1], 0), color, width)
+ draw_line_3d(ax, (s[0], 0, 0), (s[0], s[1], 0), color, width)
+ draw_line_3d(ax, ( 0, 0, s[2]), ( 0, s[1], s[2]), color, width)
+ draw_line_3d(ax, (s[0], 0, s[2]), (s[0], s[1], s[2]), color, width)
+
+ draw_line_3d(ax, ( 0, 0, 0), ( 0, 0, s[2]), color, width)
+ draw_line_3d(ax, (s[0], 0, 0), (s[0], 0, s[2]), color, width)
+ draw_line_3d(ax, ( 0, s[1], 0), ( 0, s[1], s[2]), color, width)
+ draw_line_3d(ax, (s[0], s[1], 0), (s[0], s[1], s[2]), color, width)
+
+def draw_3d_picture(sizes, r):
+ ax = plt.figure().add_subplot(111, projection="3d")
+
+ plt.axis("off")
+ ax.set_xlim(0, max(sizes))
+ ax.set_ylim(0, max(sizes))
+ ax.set_zlim(0, max(sizes))
+
+ draw_box_3d(ax, sizes)
+
+ path = get_path(3, sizes, r)
+ plt.plot([p[0] for p in path], [p[1] for p in path], [p[2] for p in path],
+ color="blue")
+
+ plt.show()
+
+###############################################################################
+# Billiard data / user input
+###############################################################################
+def user_input():
+ a = int(input("Value for a: "))
+ b = int(input("Value for b: "))
+ c = int(input("Value for c: "))
+ ra = int(input("a-coordinate of r: "))
+ rb = int(input("b-coordinate of r: "))
+ rc = int(input("c-coordinate of r: "))
+ pic = input("Choose p for projections or 3 for 3d (empty = both): ")
+ return [a,b,c], [ra,rb,rc], pic
+
+# You can choose to write your input here or to get it interactively
+#sizes, r, pic = [15, 9, 7], [2, 0, 0], "p"
+#sizes, r, pic = [2, 3, 4], [0, 0, 0], "3"
+sizes, r, pic = user_input()
+
+if pic == "p":
+ draw_3d_projections(sizes, r)
+elif pic == "3":
+ draw_3d_picture(sizes, r)
+else:
+ draw_3d_projections(sizes, r)
+ draw_3d_picture(sizes, r)
