7
\$\begingroup\$

I am always super interested in both science and art, science and art are two different ways to describe the objective reality, they are two sides of the same coin, in many areas they overlap, and mathematics is the most prominent one of such areas.

So I used the last three days to write a bunch of functions that generate artistic geometrical shapes using mathematics and visualize the shapes using matplotlib.

Previously I had absolutely zero experience with matplotlib, in the past I have used GeoGebra 5 to draw shapes, but doing it manually is really slow and prone to error, and I don't use brushes because I want my images to be precise, so I have decided to draw pictures with code and here I am.

In this gigantic script I have combined:

Hailstone Sequence

Fibonacci Sequence

Prime Factorization

Dice Average Formula

The following equation is used to calculate the average value of a set of dice.

AV = ( (M + 1 )/ 2 ) * N

Where AV is the average dice value
M is the max value of all dice
N is the total number of die

source

Full Saturation Spectrum

Explanations by me: the RGB color model uses one byte for each channel to represent colors, one byte is eight bits, 2 to the power of 8 is 256, so each channel can have a value from 0 to 255, and the total number of colors possible in this model is 16777216 (2 to the power of 24).

I don't know about much about color theory, but in my observation, fully saturated colors have only one or two bytes that are set, and at least one of the bytes is fully set. So the total number of fully saturated colors is 6 * 255 = 1530.

And of course, programming, I have written over 400 lines of code, spent more than 12 hours on this project, expressed math in code, just to create, something with absolutely no practical value whatsoever—ART.


Code

from functools import lru_cache
from matplotlib.patches import Arc, Circle, Rectangle
from itertools import combinations
import matplotlib.pyplot as plt
import random
import math
import re
FIBONACCI = [2, 3, 5, 8, 13, 21, 34, 55, 89, 144, 233]
DIRECTIONS = [('right', 'up'), ('right', 'down'), ('left', 'down'), ('left', 'up')]
ANGLES = [(90, 180), (0, 90), (-90, 0), (180, 270)]
CCW_ANGLES = [(0, 90), (90, 180), (180, 270), (270, 360)]
cycles = [
 (4, 2, 1),
 (-2, -1),
 (-14, -7, -20, -10, -5),
 (
 -50, -25, -74,
 -37, -110, -55,
 -164, -82, -41,
 -122, -61, -182,
 -91, -272, -136,
 -68, -34, -17
 )
]
loops = {term: len(cycle) for cycle in cycles for term in cycle}
stops = {term: cycle[-1] for cycle in cycles for term in cycle}
@lru_cache(maxsize=None)
def fibonacci(n):
 s = []
 a, b = 0, 1
 for i in range(n):
 s.append(a)
 a, b = b, a + b
 return s
@lru_cache(maxsize=None)
def spectrum_position(n):
 if not isinstance(n, int):
 raise TypeError('`n` should be an integer')
 if n < 0:
 raise ValueError('`n` must be non-negative')
 n %= 1530
 if 0 <= n < 255:
 return f'#ff{n:02x}00'
 elif 255 <= n < 510:
 return f'#{510-n:02x}ff00'
 elif 510 <= n < 765:
 return f'#00ff{n-510:02x}'
 elif 765 <= n < 1020:
 return f'#00{1020-n:02x}ff'
 elif 1020 <= n < 1275:
 return f'#{n-1020:02x}00ff'
 elif 1275 <= n < 1530:
 return f'#ff00{1530-n:02x}'
@lru_cache(maxsize=None)
def factorize(n):
 assert isinstance(n, int)
 assert n >= 2
 factors = []
 while not n % 2:
 factors.append(2)
 n /= 2
 
 for i in range(3, int(n**.5+1), 2):
 while not n % i:
 factors.append(i)
 n /= i
 
 if n > 1:
 factors.append(int(n))
 return factors
@lru_cache(maxsize=None)
def dice_average(n):
 factors = factorize(n)
 if len(factors) == 1:
 return (1, n * 2 - 1)
 else:
 maxf = max(factors)
 factors.remove(maxf)
 return (math.prod(factors), maxf * 2 - 1)
@lru_cache(maxsize=None)
def hailstone(n):
 if not isinstance(n, int):
 raise TypeError('n should be an integer')
 
 if n == 0:
 raise ValueError('n should not be 0')
 
 sequence = [n]
 loop_iter = 0
 looping = False
 while True:
 if not looping and n in loops:
 loop_terms = loops[n]
 looping = True
 
 if looping:
 loop_iter += 1
 if loop_iter == loop_terms:
 break
 
 if n % 2:
 n = 3 * n + 1
 else:
 n /= 2
 
 sequence.append(int(n))
 
 return sequence
def sin(x): return math.sin(math.radians(x))
def cos(x): return math.cos(math.radians(x))
def tan(d): return math.tan(math.radians(d))
@lru_cache(maxsize=None)
def negate(ls): return [-e for e in ls]
@lru_cache(maxsize=None)
def rotate(point, deg):
 x, y = point
 x1 = x * cos(deg) - y * sin(deg)
 y1 = y * cos(deg) + x * sin(deg)
 return (x1, y1)
@lru_cache(maxsize=None)
def square_vertices(init_pos, side, direction):
 x0, y0 = init_pos
 hd, vd = direction
 x1 = x0 + side if hd == 'right' else x0 - side
 y1 = y0 + side if vd == 'up' else y0 - side
 return [(x0, y0), (x1, y0), (x1, y1), (x0, y1)]
def make_square(init_pos, side, direction, color=None):
 vertices = square_vertices(init_pos, side, direction)
 x, y = zip(*vertices)
 minx = min(x)
 miny = min(y)
 if not color:
 color = '#'+random.randbytes(3).hex()
 return Rectangle((minx, miny), side, side, color=color, alpha=0.5), vertices
@lru_cache(maxsize=None)
def first_level(number):
 assert number > 2
 half = 180 / number
 x1 = (16 / (1 + tan(half) ** 2)) ** .5
 y1 = (16 - x1 ** 2) ** .5
 k = y1 - tan(half - 90) * x1
 x2 = -k / tan(half - 90)
 r = ((x1 - x2) ** 2 + y1 ** 2) ** .5
 return (x2, r)
@lru_cache(maxsize=None)
def next_level(distance, radius):
 base = (distance - radius) / (distance + radius)
 new_radius = radius * base
 new_distance = distance - radius - new_radius
 return (new_distance, new_radius)
@lru_cache(maxsize=None)
def fade_color(color, strength):
 assert re.match(r'^#[0-9a-f]{6}$', color)
 assert 0 < strength <= 1
 color = int(color[1:], 16)
 r, g, b = [round(((color>>i)%256)*strength) for i in (16, 8, 0)]
 return f'#{r:02x}{g:02x}{b:02x}'
def fibonacci_spiral(n, c=12, bsquare=True, barc=True, reverse=False):
 assert barc or bsquare
 colors = [spectrum_position(int(1530*i/c)) for i in range(c)]
 sequence = fibonacci(n+1)
 init_pos = (0, 0)
 vertexes = []
 shapes = []
 for i, n in enumerate(sequence[1:]):
 direction = DIRECTIONS[i % 4] if not reverse else DIRECTIONS[::-1][i % 4]
 angle = ANGLES[i % 4] if not reverse else CCW_ANGLES[i % 4]
 square, vertices = make_square(init_pos, n, direction, colors[i%c])
 center = vertices[[1, 3][i % 2]]
 a, b = angle
 arc = Arc(center, 2*n, 2*n, theta1=a, theta2=b, edgecolor='k', lw=3.0)
 if bsquare: shapes.append(square)
 if barc: shapes.append(arc)
 init_pos = vertices[-2]
 vertexes.extend(vertices)
 
 return shapes, vertexes
def fibonacci_arcs(n, degree, color='k' ,reverse=False, alpha=0.75, lw=2):
 sequence = fibonacci(n)
 init_pos = (0, 0)
 vertexes = []
 shapes = []
 for i, n in enumerate(sequence[1:]):
 direction = DIRECTIONS[i % 4] if not reverse else DIRECTIONS[::-1][i % 4]
 angle = ANGLES[i % 4] if not reverse else CCW_ANGLES[i % 4]
 vertices = square_vertices(init_pos, n, direction)
 rotated_vertices = [rotate(pos, degree) for pos in vertices]
 index = [1, 3][i % 2]
 center = rotated_vertices[index]
 a, b = [(i + degree) % 360 for i in angle]
 arc = Arc(center, 2*n, 2*n, theta1=a, theta2=b, edgecolor=color, lw=lw, alpha=alpha)
 init_pos = vertices[-2]
 shapes.append(arc)
 vertexes.extend(rotated_vertices)
 
 return shapes, vertexes
################################################################
# plotting functions
def plot_hailstone(n, real=False):
 assert n > 1
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 numbers = list(range(-n, 0)) + list(range(1, n))
 for i in numbers:
 ax.plot(hailstone(i))
 
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 if real:
 plt.axis('scaled')
 plt.box(False)
 plt.set_cmap('rainbow')
 plt.show()
def plot_spiral(n, c, barc=True, bsquare=True, reverse=False):
 shapes, vertexes = fibonacci_spiral(n, c, barc=barc, bsquare=bsquare, reverse=reverse)
 x, y = zip(*vertexes)
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 for shape in shapes:
 ax.add_patch(shape)
 
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 plt.axis('scaled')
 plt.box(False)
 ax = plt.gca()
 ax.set_xlim([min(x), max(x)])
 ax.set_ylim([min(y), max(y)])
 plt.show()
def plot_lolipop(n, arms, lw=1):
 arcs = []
 points = []
 u = 360/arms
 for i in range(arms):
 shapes, vertices = fibonacci_arcs(n, i*u, color='k', alpha=1, lw=lw)
 arcs.extend(shapes)
 points.extend(vertices)
 
 x, y = zip(*points)
 a, b = points[-2]
 radius = (a*a+b*b)**.5
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 for arc in arcs:
 ax.add_patch(arc)
 
 ax.add_patch(Arc((0, 0), 2*radius, 2*radius, theta1=0, theta2=360, edgecolor='k', lw=lw))
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 plt.axis('scaled')
 plt.box(False)
 ax = plt.gca()
 ax.set_xlim([min(x), max(x)])
 ax.set_ylim([min(y), max(y)])
 plt.show()
def plot_sunflower(n, cw, color=True, alpha=0.75):
 assert cw in FIBONACCI and cw <= 144
 i = FIBONACCI.index(cw)
 ccw = FIBONACCI[i+1]
 arcs = []
 points = []
 cwu = 360/cw
 ccwu = 360/ccw
 color_cwu = 1530/cw
 color_ccwu = 1530/ccw
 for i in range(cw):
 color1 = spectrum_position(round(i*color_cwu)) if color else 'k'
 shapes, vertices = fibonacci_arcs(n, i*cwu, color=color1, alpha=alpha)
 arcs.extend(shapes)
 points.extend(vertices)
 
 for i in range(ccw):
 color2 = spectrum_position(round(1530-i*color_ccwu)) if color else 'k'
 shapes, vertices = fibonacci_arcs(n, (90+i*ccwu)%360, color=color2, reverse=True, alpha=alpha)
 arcs.extend(shapes)
 points.extend(vertices)
 
 x, y = zip(*points)
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 for arc in arcs:
 ax.add_patch(arc)
 
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 plt.axis('scaled')
 plt.box(False)
 ax = plt.gca()
 ax.set_xlim([min(x), max(x)])
 ax.set_ylim([min(y), max(y)])
 plt.show()
def recur_plot_circle(division, levels, fade=0.9):
 assert division > 2
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 rotation_diff = 90 / division
 unit_rotation = 360 / division
 color_increment = 1530 / division
 distance, radius = first_level(division)
 max_distance = distance + radius
 for level in range(levels):
 for index in range(division):
 center = rotate((0, distance), level * rotation_diff + index * unit_rotation)
 color = spectrum_position(round(index * color_increment))
 color = fade_color(color, fade**level)
 ax.add_patch(Circle(center, radius, color=color))
 distance, radius = next_level(distance, radius)
 
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 plt.axis('scaled')
 plt.box(False)
 ax = plt.gca()
 ax.set_xlim([-max_distance, max_distance])
 ax.set_ylim([-max_distance, max_distance])
 plt.show()
def plot_star(number):
 assert number > 2
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 unit_rotation = 360 / number
 points = []
 for i in range(number):
 points.append(rotate((0, 4), unit_rotation*i))
 
 for a, b in combinations(points, r=2):
 x, y = zip(*[a, b])
 ax.plot(x, y)
 
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 plt.axis('scaled')
 plt.box(False)
 ax = plt.gca()
 ax.set_xlim([-4, 4])
 ax.set_ylim([-4, 4])
 plt.set_cmap('rainbow')
 plt.show()
def plot_tetragram(n):
 assert n > 2 and not (n - 2) % 3
 averages = [dice_average(i) for i in range(2, n)]
 fig = plt.figure(frameon=False)
 ax = fig.add_subplot(111)
 ax.set_axis_off()
 for a, b, c in zip(averages[::3], averages[1::3], averages[2::3]):
 x, y = zip(*[a, b, c])
 ax.triplot(x, y)
 ax.triplot(negate(x), y)
 ax.triplot(negate(x), negate(y))
 ax.triplot(x, negate(y))
 ax.triplot(y, x)
 ax.triplot(negate(y), x)
 ax.triplot(negate(y), negate(x))
 ax.triplot(y, negate(x))
 
 fig.subplots_adjust(left=0, bottom=0, right=1, top=1, wspace=0, hspace=0)
 plt.axis('scaled')
 plt.box(False)
 plt.set_cmap('rainbow')
 plt.show()
if __name__ == '__main__':
 plot_hailstone(128)
 plot_hailstone(12, 1)
 plot_spiral(24, 6)
 plot_spiral(24, 6, reverse=True)
 plot_lolipop(25, 12)
 for i in (2, 3, 5, 8, 13):
 plot_sunflower(24, i)
 
 for i in (3, 6, 12, 24):
 recur_plot_circle(i, round(i/2))
 
 for i in range(4, 13):
 plot_star(i)
 
 for i in range(14, 32, 3):
 plot_tetragram(i)

All examples:

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here

enter image description here


How can the script be improved? I am really new to all this stuff.

toolic
15k5 gold badges29 silver badges208 bronze badges
asked May 14, 2022 at 8:31
\$\endgroup\$

1 Answer 1

1
\$\begingroup\$

Thanks for posting this code; all the cool images really caught my eye.

This is just a mini review since I don't have much to offer in the math arena.

DRY

In the spectrum_position function, some of the same numbers are repeated in the elif lines.

if 0 <= n < 255:
 ...
elif 255 <= n < 510:
 ...
elif 510 <= n < 765:

When the if condition is false, we know that n must be greater than or equal to 255. Therefore, the lower bound check is not needed. We can simplify:

elif 255 <= n < 510:

as:

elif n < 510:

The same is true for the remaining lower bound checks. This makes the code a little easier to read, and it may be slightly more efficient performing fewer checks.

In the plot_tetragram function, these two expressions are repeated several times:

negate(x)
negate(y)

Consider assigning each to a variable:

x, y = zip(*[a, b, c])
neg_x = negate(x)
neg_y = negate(y)
ax.triplot(x, y)
ax.triplot(neg_x, y)

Readability

In lines like this:

radius = (a*a+b*b)**.5

it might be easier to understand the code using the sqrt function:

radius = math.sqrt(a*a+b*b)

This is part of the math library, which you are already importing.

Documentation

It would be nice to have a docstring at the top of the code to summarize its purpose, as well as some sprinkled throughout the code for some of the more complex functions.

Magic number

The number 1530 is used in a few places throughout the code. Consider making it a named constant, for example:

NUM_COLORS = 6*255
answered Jun 6 at 17:02
\$\endgroup\$
2
  • 1
    \$\begingroup\$ I would use math.hypot(a, b) in place of math.sqrt(a*a+b*b) for simplicity. \$\endgroup\$ Commented Jun 9 at 11:12
  • \$\begingroup\$ @TobySpeight: Absolutely! I had a feeling there was some other function, but I wasn't aware of hypot. Thanks! \$\endgroup\$ Commented Jun 9 at 11:14

Your Answer

Draft saved
Draft discarded

Sign up or log in

Sign up using Google
Sign up using Email and Password

Post as a guest

Required, but never shown

Post as a guest

Required, but never shown

By clicking "Post Your Answer", you agree to our terms of service and acknowledge you have read our privacy policy.

Start asking to get answers

Find the answer to your question by asking.

Ask question

Explore related questions

See similar questions with these tags.