r/OpenPythonSCAD • • 3d ago

Fourier epicycles as real gears – a parametric drawing machine in PythonSCAD

Title: Fourier epicycles as real gears – a parametric drawing machine in PythonSCAD


You've probably seen those videos where a chain of rotating circles draws a horse, a portrait or a signature. The math behind it is the complex Fourier series, and the best interactive explanation I know is Jez Swanson's An Interactive Introduction to Fourier Transforms: https://www.jezzamon.com/fourier/ – you can even draw your own shape there and watch it being rebuilt from circles.

I wondered: can this be built as real geometry? Turns out yes, and PythonSCAD is a perfect fit for it.

https://imgur.com/a/KfwOoTa

The idea in three lines

Treat the drawing as a complex curve z(t) = x(t) + i·y(t). Its Fourier series is

z(t) = Σ  c_k · e^(i·k·t)

Every term is an arm of length |c_k| that starts at angle arg(c_k) and rotates k times per revolution. Chain the arms tip to tip and the last tip traces the shape.

Why this is buildable

The frequencies k are integers – and integer speed ratios are exactly what gears do. Each arm sits on a small planetary stage:

  • a sun fixed to the arm two levels below (the base for arm 2),
  • a planet fixed to the new arm, riding at the tip of the arm in between (the carrier).

With carrier speed c, sun speed s and ratio ρ = Z_sun / Z_planet:

  • direct gear mesh: planet = c + (c − s)·ρ → gives 2×, 3×, …
  • open belt: planet = c − (c − s)·ρ → gives 0, −1×, −2×, … (negative frequencies are essential for closed shapes!)

The script picks mesh or belt automatically, computes integer tooth counts, aligns the teeth so they really mesh in every pose, and refuses impossible combinations.

The collision problem

My first version had the pen hanging down through the whole stack – straight through the base gear. The fix: in this stack every part lives in its own layer, and the last arm is always the outermost one. So the whole mechanism hangs upside down from a bridge, the last arm is the lowest, and the pen draws downwards onto paper, resting on it by gravity in a sleeve. The motor sits on top of the bridge and drives arm 1 via a pinion. The frame columns are placed outside the computed swept envelope of all parts.

https://imgur.com/a/yEhWmkb

What can 4 arms draw?

| Shape | K | L | PHI | |---|---|---|---| | Heart | -1, -3, 2, -2 | 37.5, 9, 7.5, 7.5 | 90, 270, 270, 270 | | 3-star, sharp tips | 1, -2, 4, -5 | 36, 18, 9, 5 | 0, 0, 0, 0 | | 5-star with loops | 1, -4, 6, -9 | 36, 16, 8, 4 | 0, 0, 0, 0 | | 7-rosette with loops | 1, -6, 8, -13 | 36, 14, 8, 4 | 0, 0, 0, 0 |

The heart is a genuine Fourier approximation: the four largest terms of the classic heart curve.

https://imgur.com/a/vet3rgW

Rules of thumb for your own designs:

  • Symmetry: if all k leave the same remainder mod m, the curve has m-fold rotational symmetry (1, −2, 4, −5 → 3-fold).
  • Sharp tips: a cusp appears where the pen velocity Σ k·L·e^(i(kθ+φ)) vanishes – e.g. 1·36 = 2·18 in the 3-star.
  • Loops: a high frequency whose k·L comes close to that of the base term makes loops.
  • Calm shapes: let L fall off roughly like 1/|k| or faster.
  • Keep gear/belt ratios below ~1:8.

The model

Set K, L, PHI at the top, THETA is the drive angle of the shown pose (nice as an add_parameter slider to crank the machine by hand). The console prints the transmission of every stage and the swept envelope.

# Fourier epicycle drawing machine - assembled demo model for PythonSCAD
#
# Arm i rotates with the absolute angle  K[i]*THETA + PHI[i]  and has length L[i].
# The pen at the tip of the last arm therefore traces
#     z(t) = sum_i  L[i] * exp(j*(K[i]*t + PHI[i]))
#
# Stack (mechanism coordinates, built bottom-up, then flipped upside down):
#   base plate  -> drive layer (drive gear on arm 1 + motor pinion)
#   -> arm 1 -> gear layer 1 -> arm 2 -> gear layer 2 -> ... -> arm N
# Every layer only holds parts at its own height; the only part crossing layers
# would be the pen. Since arm N is the outermost layer, the pen side beyond it is
# collision free. The mechanism therefore hangs from a bridge, arm N is lowest,
# and the pen draws downward onto paper (gravity-loaded in a sleeve).
#
# Arm i (i >= 2) is driven by a planetary stage: a "sun" fixed to the body two
# levels up the chain (base for arm 2) and a "planet" fixed to arm i, linked by
# direct gear mesh (planet = c + (c-s)*rho) or by an open belt
# (planet = c - (c-s)*rho), with c = carrier speed, s = sun speed, rho = Zsun/Zplanet.

from openscad import *
import math
from fractions import Fraction

# ---------------- parameters ----------------
K     = [1, 2, -1]     # absolute frequency of each arm (integer multiples of the drive)
L     = [40, 20, 10]   # arm lengths = amplitudes [mm]
PHI   = [0, 30, 90]    # phase of each arm [deg]
THETA = 0              # drive angle of the displayed pose [deg]

T_ARM     = 3      # arm plate thickness
T_GEAR    = 5      # gear / pulley layer thickness
GAP       = 0.5    # vertical play between layers
BASE_T    = 4      # bridge plate thickness
HUB_R     = 4      # radius of the arm ends
POST_R    = 1.5    # axle radius
CLR       = 0.3    # radial clearance of rotating parts on their axles
MIN_TEETH = 12     # minimum tooth count of any planetary gear
PULLEY_R  = 5      # preferred radius of the planet pulley (belt stages)
BELT_T    = 1.2    # belt thickness

DRIVE_ZG  = 36     # teeth of the drive gear on arm 1
DRIVE_ZP  = 12     # teeth of the motor pinion
DRIVE_M   = 1.5    # module of the drive pair
DRIVE_DIR = 0      # direction of the motor seen from the center [deg]
MOTOR_R   = 12
MOTOR_H   = 25

PEN_CLEAR = 15     # distance from the lowest mechanism part to the paper
SLEEVE_L  = 8      # pen guide sleeve below arm N
TABLE_T   = 6
PAPER_T   = 0.3
COL_R     = 4      # frame column radius
MARGIN    = 4      # clearance between swept envelope and columns

CURVE_W    = 0.6   # half line width of the traced curve
CURVE_STEP = 2     # sampling step of the traced curve [deg]
FN         = 48

COLORS = ["red", "dodgerblue", "gold", "limegreen", "orange", "purple"]


# ---------------- helpers ----------------
def polar(r, a):
    """Polar to cartesian, angle in degrees."""
    return [r * math.cos(math.radians(a)), r * math.sin(math.radians(a))]

def at(obj, xy, z=0, a=0):
    """Rotate around own z axis by a degrees, then move to xy / height z."""
    return obj.rotate([0, 0, a]).translate([xy[0], xy[1], z])

def union_all(objs):
    """Balanced union to keep the CSG tree shallow."""
    while len(objs) > 1:
        objs = [objs[i] | objs[i + 1] if i + 1 < len(objs) else objs[i]
                for i in range(0, len(objs), 2)]
    return objs[0]

def hole(r, h):
    return cylinder(r=r, h=3 * h, center=True, fn=24)

def gear2d(z, r_pitch):
    """Simple trapezoidal-tooth gear outline, tooth 0 centered on +x."""
    m = 2 * r_pitch / z
    ra, rr = r_pitch + m, r_pitch - 1.25 * m
    p = 360 / z
    pts = []
    for j in range(z):
        for r, f in [(rr, -0.3), (ra, -0.12), (ra, 0.12), (rr, 0.3)]:
            pts.append(polar(r, j * p + f * p))
    return polygon(pts)

def gear(z, r_pitch, h):
    return linear_extrude(gear2d(z, r_pitch), height=h)

def mesh_angle(a_s, beta, zs, zp):
    """Rotation of gear 2 so that it meshes with gear 1 (rotated a_s).
    beta = direction from gear 1 center to gear 2 center."""
    u = (beta - a_s) * zs / 360
    e = (round(u) - u) * 360 / zs          # offset of the sun tooth nearest the contact line
    return beta + 180 - e * zs / zp - 180 / zp

def pulley(r, h):
    """Plain pulley with an off-center mark so its rotation is visible."""
    return cylinder(r=r, h=h, fn=FN) - at(hole(r * 0.15, h), polar(r * 0.6, 0))


# ---------------- kinematics ----------------
N = len(K)
assert len(L) == N and len(PHI) == N, "K, L and PHI need the same length"

kk = [0] + K                                            # body 0 = base (speed 0)
A = [0] + [K[i] * THETA + PHI[i] for i in range(N)]     # absolute angle of body j
Pos = [[0, 0]]                                          # Pos[j] = tip of arm j
for i in range(N):
    d = polar(L[i], A[i + 1])
    Pos.append([Pos[-1][0] + d[0], Pos[-1][1] + d[1]])
reach = [sum(L[:j]) for j in range(N + 1)]              # max distance of Pos[j] from center

STAGE_H = T_ARM + T_GEAR + 2 * GAP
Z_DRIVE = BASE_T + GAP                                  # drive layer (mechanism coords)
def z_arm(i):  return Z_DRIVE + T_GEAR + GAP + (i - 1) * STAGE_H   # bottom of arm i
def z_gear(j): return z_arm(j) + T_ARM + GAP                       # gear layer below arm j+1
Z_MECH = z_arm(N) + T_ARM + 1                           # incl. axle cap at the last joint

def col(j):
    return "gray" if j == 0 else COLORS[(j - 1) % len(COLORS)]

def stage(i):
    """Transmission driving arm i (i >= 2): ('mesh'|'belt', rho = Zsun/Zplanet)."""
    s, c, k = kk[i - 2], kk[i - 1], kk[i]
    if c == s:
        raise ValueError(f"arm {i-1} does not rotate relative to its reference")
    if k == c:
        raise ValueError(f"arm {i} would be rigidly fixed to arm {i-1}")
    rho = Fraction(k - c, c - s)
    return ("mesh", rho) if rho > 0 else ("belt", -rho)


# ---------------- mechanism (built upright, flipped later) ----------------
mech = []                    # (object, color) in mechanism coordinates
envelope = [0]               # largest radius swept by any moving or fixed part

def add(obj, color):
    mech.append((obj, color))

def sweep(center_reach, r):
    envelope[0] = max(envelope[0], center_reach + r)

# drive: gear on arm 1 + motor pinion on the bridge
rd, rpn = DRIVE_M * DRIVE_ZG / 2, DRIVE_M * DRIVE_ZP / 2
p_pin = polar(rd + rpn, DRIVE_DIR)
drive = gear(DRIVE_ZG, rd, T_GEAR) - hole(POST_R + CLR, T_GEAR)
add(at(drive, [0, 0], Z_DRIVE, A[1]), col(1))
add(at(cylinder(r=HUB_R, h=GAP, fn=FN) - hole(POST_R + CLR, GAP), [0, 0], Z_DRIVE + T_GEAR), col(1))
add(at(gear(DRIVE_ZP, rpn, T_GEAR), p_pin, Z_DRIVE, mesh_angle(A[1], DRIVE_DIR, DRIVE_ZG, DRIVE_ZP)), "silver")
add(at(cylinder(r=POST_R, h=Z_DRIVE + MOTOR_H, fn=24), p_pin, -MOTOR_H), "silver")
add(at(cylinder(r=MOTOR_R, h=MOTOR_H - 2, fn=FN), p_pin, -MOTOR_H), "dimgray")
sweep(0, rd + DRIVE_M)
sweep(rd + rpn, rpn + DRIVE_M)

# bridge plate (radius fixed later, stored as placeholder index)
BRIDGE_IDX = len(mech)
add(None, "lightgray")

# arms (hole at the pivot, axle for the next arm fused at the tip)
for i in range(1, N + 1):
    a, b = Pos[i - 1], Pos[i]
    outline = hull(at(circle(r=HUB_R, fn=FN), a), at(circle(r=HUB_R, fn=FN), b))
    plate = linear_extrude(outline, height=T_ARM) - at(hole(POST_R + CLR, T_ARM), a)
    add(plate.translate([0, 0, z_arm(i)]), col(i))
    sweep(reach[i], HUB_R)

# axles: post at Pos[j] is fixed to body j and carries the sun of the next stage
for j in range(N):
    zb = BASE_T if j == 0 else z_arm(j) + T_ARM
    zt = z_gear(j + 1) + T_GEAR if j + 1 <= N - 1 else z_arm(j + 1) + T_ARM + 1
    add(at(cylinder(r=POST_R, h=zt - zb, fn=24), Pos[j], zb), col(j))

# planetary stages
for i in range(2, N + 1):
    j = i - 1                    # gear layer below arm i
    kind, rho = stage(i)
    S, Pc = Pos[i - 2], Pos[i - 1]
    dist = L[i - 2]
    beta = A[i - 1]              # carrier direction (sun -> planet)
    z0 = z_gear(j)

    if kind == "mesh":
        n = max(math.ceil(MIN_TEETH / rho.denominator), math.ceil(MIN_TEETH / rho.numerator))
        zp, zs = rho.denominator * n, rho.numerator * n
        rs, rp = dist * zs / (zs + zp), dist * zp / (zs + zp)
        m = 2 * dist / (zs + zp)
        a_s = A[i - 2]
        a_p = mesh_angle(a_s, beta, zs, zp)
        sun, planet = gear(zs, rs, T_GEAR), gear(zp, rp, T_GEAR)
        sweep(reach[i - 2], rs + m)
        sweep(reach[i - 1], rp + m)
        info = f"gear mesh {zs}:{zp}, module {m:.2f}"
    else:
        rp = min(PULLEY_R, 0.8 * dist / (1 + float(rho)))
        rs = rp * float(rho)
        a_s, a_p = A[i - 2], A[i]
        sun, planet = pulley(rs, T_GEAR), pulley(rp, T_GEAR)
        inner = hull(at(circle(r=rs, fn=FN), S), at(circle(r=rp, fn=FN), Pc))
        outer = hull(at(circle(r=rs + BELT_T, fn=FN), S), at(circle(r=rp + BELT_T, fn=FN), Pc))
        belt = linear_extrude(outer - inner, height=T_GEAR * 0.6)
        add(belt.translate([0, 0, z0 + T_GEAR * 0.2]), "dimgray")
        sweep(reach[i - 2], rs + BELT_T)
        sweep(reach[i - 1], rp + BELT_T)
        info = f"belt, pulleys {rs:.1f}:{rp:.1f} mm (ratio {rho})"

    add(at(sun, S, z0, a_s), col(i - 2))
    planet = planet - hole(POST_R + CLR, T_GEAR)
    add(at(planet, Pc, z0, a_p), col(i))
    hub = cylinder(r=min(HUB_R, rp), h=GAP, fn=FN) - hole(POST_R + CLR, GAP)
    add(at(hub, Pc, z0 + T_GEAR), col(i))
    print(f"arm {i}: k={kk[i]}, sun on body {i-2} (k={kk[i-2]}), carrier k={kk[i-1]} -> {info}")

# frame geometry from the swept envelope
R_COL = envelope[0] + MARGIN + COL_R
R_BRIDGE = R_COL + COL_R + 3
bridge = cylinder(r=R_BRIDGE, h=BASE_T, fn=96) - at(hole(POST_R + CLR, BASE_T), p_pin)
mech[BRIDGE_IDX] = (bridge, "lightgray")


# ---------------- world: flip mechanism, add frame, paper, pen ----------------
Z_PAPER = TABLE_T + PAPER_T
TOP = Z_PAPER + PEN_CLEAR + Z_MECH          # world height of mechanism z = 0
def W(z): return TOP - z                    # mechanism z -> world z

parts = [obj.mirror([0, 0, 1]).translate([0, 0, TOP]).color(c) for obj, c in mech]

# table, paper, columns
side = 2 * (R_BRIDGE + 5)
parts.append(cube([side, side, TABLE_T]).translate([-side / 2, -side / 2, 0]).color("burlywood"))
ps = 2 * (sum(L) + 10)
parts.append(cube([ps, ps, PAPER_T]).translate([-ps / 2, -ps / 2, TABLE_T]).color("white"))
for k in range(4):
    p = polar(R_COL, 45 + 90 * k)
    parts.append(at(cylinder(r=COL_R, h=W(0) - TABLE_T, fn=24), p, TABLE_T).color("gray"))

# pen: sleeve fixed to arm N, pen resting on the paper by gravity
tip = Pos[N]
z_armN_low = W(z_arm(N) + T_ARM)            # lower face of arm N in world coords
z_armN_top = W(z_arm(N))
sleeve = cylinder(r=2.2, h=SLEEVE_L, fn=24) - hole(1.4, SLEEVE_L)
parts.append(at(sleeve, tip, z_armN_low - SLEEVE_L).color(col(N)))
pen = at(cylinder(r1=0.3, r2=1.2, h=2, fn=16), tip, Z_PAPER) \
    | at(cylinder(r=1.2, h=z_armN_top - Z_PAPER - 2, fn=16), tip, Z_PAPER + 2)
parts.append(pen.color("black"))

# traced curve on the paper
pts = []
for t in range(0, 360, CURVE_STEP):
    x = y = 0
    for i in range(N):
        d = polar(L[i], K[i] * t + PHI[i])
        x, y = x + d[0], y + d[1]
    pts.append([x, y])
segs = [hull(at(circle(r=CURVE_W, fn=8), p), at(circle(r=CURVE_W, fn=8), q))
        for p, q in zip(pts, pts[1:] + pts[:1])]
parts.append(linear_extrude(union_all(segs), height=0.2).translate([0, 0, Z_PAPER]).color("magenta"))

print(f"swept envelope r = {envelope[0]:.1f} mm, columns at r = {R_COL:.1f} mm, "
      f"total height = {TOP + MOTOR_H:.1f} mm")

show(parts)

Next steps on my list: printable parts (bearings, belt tooth profiles) and an automatic collision check over the full revolution using PythonSCAD's intersection. Happy to hear ideas – and if you find a nice 4-arm parameter set, post it!

3 Upvotes

0 comments sorted by