如何用Matplotlib绘制Tupper自指公式的不等式图形？

Tupper自指公式本质是离散二进制像素的判定规则，plt.plot()是用来绘制连续线条的，自然不适用。要画出它的阴影区域，需要用网格渲染类的方法，结合numpy向量化计算提升效率，具体步骤如下：

核心思路 #

公式的逻辑是：对每个(x,y)，判断是否满足不等式，满足的点标记为深色（阴影），不满足的为浅色。其中x的有效范围是0~105，y的范围是N到N+16，这里的N是Tupper公式专用的超大常数（没有这个值画不出自指图形）。

代码实现 #

import numpy as np
import matplotlib.pyplot as plt

# Tupper公式专用的超大常数N
N = 960939379918958884971672962127852754715004339660129306651505519271702802395266424689642842174350718121267153782770623355993237280874144307891325963941337723487857735749823926629715517173716995165232890538221612403238855866184013235585136048828693337902491454229288667081096184496091705183454067827731551705405381627380967602565625016981482083418783163849115590225610003652351370343874461848378737238198224849863465033159410054974700593138339226497249461751545728366702369745461014655997933798537483143786841806593422227898388722980000748404719

# 生成x和y的网格：x从0到105，y从N到N+16
x = np.arange(0, 106)
y = np.arange(N, N+17)
X, Y = np.meshgrid(x, y)

# 向量化计算公式条件，避免循环
floor_y = np.floor(Y)
floor_x = np.floor(X)
mod_y = floor_y % 17
exponent = -17 * floor_x - mod_y
term = (floor_y / 17) * (2 ** exponent)
floor_term = np.floor(term)
result = floor_term % 2
# 得到满足条件的布尔矩阵
mask = result > 0.5

# 绘制图形
plt.figure(figsize=(12, 2))
# 使用pcolormesh渲染网格，cmap用binary实现黑白阴影
plt.pcolormesh(X, Y, mask, cmap='binary', edgecolors='none')
# 反转y轴，因为Tupper公式的y值越大越靠下
plt.gca().invert_yaxis()
# 隐藏坐标轴刻度，让图形更简洁
plt.xticks([])
plt.yticks([])
plt.show()

关键说明 #

• 向量化计算：用numpy的数组操作替代循环，比逐点计算快几个数量级，避免处理106*17=1802个点时的性能问题。
• 坐标轴反转：matplotlib默认y轴向上递增，但Tupper公式的y值越大对应图形的下方，必须反转才能得到正确的自指效果。
• N值的作用：这个超大常数是预先计算好的，它的二进制编码就是Tupper公式的图形本身，没有它只能得到随机的黑白块，无法画出自指图形。

内容的提问来源于stack exchange，提问作者Michael Carpiano