奇阶幻方公式法全解析：从罗伯法到模运算与Python实现

如果你第一次接触“幻方”大概是在某个数学兴趣小组或者奥数题里一个3×3的方格把1到9填进去让每行、每列、两条对角线的和都相等。很多人会下意识地开始试数试几分钟能试出3阶等到了5阶、7阶靠碰运气基本就走不通了。我第一次见到“奇阶幻方公式法”的时候脑子里冒出来的问题是这种“填数游戏”居然能有一个通用公式让任意奇数阶幻方都机械式地填出来后来认真推了一遍才发现它不仅是个公式背后还有一套非常漂亮的模运算原理。这篇文章就把奇阶幻方公式法彻底拆开讲从最经典的手算口诀到数学上的递推式、闭式公式再到可以直接复制运行的Python代码。无论你是刚接触幻方的数学爱好者还是想在项目里用幻方做个数据排列、生成测试用例的开发者都能从里面拿到可以直接用的东西。1. 奇阶幻方先搞懂我们到底在填什么1.1 幻方的定义与核心参数幻方也叫魔方、纵横图是一个n×n的方阵里面填入1到n²这些不同的正整数使得每一行、每一列以及两条主对角线上的数字之和都相等。这个相等的和叫“幻和”它其实可以提前算出来因为1到n²所有数字的总和是固定的S 1 2 ... n² n²(n²1)/2既然方阵有n行每行和相等那么幻和就是M S / n n(n²1)/2所以3阶幻方的幻和是155阶是657阶是175。这个数字本身可以作为验证工具你构造完一个幻方之后随便抽一行求和如果不是这个数那一定填错了。这里有一个很容易被忽略的小知识不同阶数的幻方构造难度差异非常大。奇数阶、双偶数阶能被4整除、单偶数阶形如4k2比如6、10、14各有各的构造路线。其中奇数阶恰好是公式法最成熟、最统一的一类本文要讲的“奇阶幻方公式法”就是专门解决这个问题的。1.2 为什么奇数阶值得单独拿出来讲从历史上看幻方最早的记载就是三阶幻方“洛书”也就是九宫图。古人把1到9排成三行三列得到的就是一个标准奇阶幻方。这个3阶幻方的神奇之处在于它数字很小但排列的规律已经隐含了“奇数阶幻方”最重要的特征——中心数字是整个序列的中位数5四个角上都是偶数四边中间都是奇数。这些性质在更高阶的奇阶幻方里也有体现比如中心位置永远是(n²1)/2。奇数阶幻方之所以能形成一套“公式法”核心原因是奇数阶下“斜向移动”可以形成完美的循环路径而偶数阶会在这个循环上提前“撞车”。更直观地说你从第一行中间出发按照“右上一步”的规则走奇数阶情况下绕n步就能回到起点偶数阶情况下这个周期性会被破坏。这一点在第3节我会用模运算严格展开。1.3 公式法解决的是“从无到有”的问题幻方问题可以分成两类验证一个给定方阵是不是幻方这个简单暴力求和就行但“从无到有”构造一个n阶幻方在没有公式的时候难度陡增。如果你写一个回溯搜索程序去尝试所有排列n3时9!种排列n5时25!种大约1.55×10²⁵种即使剪枝也不可能在合理时间内完成。奇阶幻方公式法解决的就是这个“从无到有”的构造问题。它把“搜索”降维成了“按规则执行”每填一个数字下一步落在哪里完全由当前坐标和一个固定规则决定不需要判断、不需要回溯按部就班就能填完。这也是组合数学里一个非常经典的“构造性算法”案例。2. 罗伯法三步口诀从3阶一路填到任意奇阶2.1 操作口诀在西方数学界最经典的奇阶幻方构造法叫Siamese method中文资料里通常叫“罗伯法”或“连续摆数法”。它最早可以追溯到17世纪但思路和几百年前法国外交官de la Loubère在暹罗今泰国记录的玩法一脉相承。操作规则只有三条第一步把数字1放在第一行的正中间。第二步每个数字都放在前一个数字的“右上方”也就是行号减1、列号加1。第三步两个补充规则。如果“右上方”越出了方阵边界就把这个位置映射到方阵的另一侧类似于世界地图从左边出国从右边回来的效果。如果“右上方”已经被数字占据了就把当前数字放在前一个数字的正下方一格。这三条规则看上去简单但里面每一步都有讲究。起点为什么在第一行中间因为这样整个路径最对称。为什么越界要回绕因为在模n的世界里方阵的左右边界是接壤的上下边界也是接壤的。为什么占位就往下这一步正是保证所有n²个格子都能被走遍的关键也是整个公式法的“气眼”。2.2 五阶幻方完整手填过程规则光看不够我建议你拿笔亲自填一遍5阶。下面是前几步1放在第0行第2列5阶方阵中间列是第2列。2放在1的右上方行-1变成第4行回绕列1变成第3列所以2在第4行第3列。3放在2的右上方第3行第4列。4放在3的右上方第2行第5列越界回绕到第0列所以4在第2行第0列。5放在4的右上方第1行第1列。6的“右上方”是第0行第2列那个位置已经被1占据了所以按照规则6要放在5的正下方也就是第2行第1列。继续这个流程填完25个数字后得到如下方阵17241815235714164613202210121921311182529可以验证一下每一行的和都是65每一列的和也都是65两条对角线的和同样是65。从1到25全部出现且没有重复。这个方阵还能看出几个特征中心位置是13恰好是(251)/2任意关于中心对称的两个数相加都等于26比如17和9、24和2、1和25。这个“中心对称互补”的性质在罗伯法生成的方阵里是天然的。2.3 口诀里最容易误解的一句话我在给朋友讲这个口诀的时候发现十个人里有三四个会在“被占住就往下”这里卡住。有人理解成“右上方格子的正下方”还有人理解成“随便找下面一个空位”。正确理解是在准备填新数字时先算出“前一个数字的右上位置”如果这个位置被占了新数字应该放在“前一个数字自己”的正下方而不是右上方格子的正下方。举个例子5阶幻方里填6的时候前一个数字5在第1行第1列它的右上方是第0行第2列已经被1占据。这时候6放在谁的正下方是5的正下方也就是第2行第1列。如果你误放成右上角那个格子的正下方也就是第1行第2列整个路径就会崩掉后面连续好几个格子都会错位。这个细节同样会影响代码实现。很多人写循环的时候先算出next_r, next_c (r-1)%n, (c1)%n然后发现square[next_r][next_c]不为0就写r (next_r1)%n这其实是错的。正确写法应该是r (r1)%n保持列不变因为新数字要放在“当前数字”的正下方。3. 把口诀翻译成模运算公式法的数学内核3.1 坐标递推式如果只满足于“会填”那口诀就够了。但要想搞清楚这个公式法“为什么永远成立”就得把口诀翻译成数学语言。这里用(r, c)表示行、列坐标所有坐标都在0到n-1之间并约定下标从0开始。罗伯法的每一步可以写成下面的递推式初始条件r0 0c0 n//2第一行正中间。对第k步k 0, 1, ..., n²-1对应数字k1先计算右上方向的候选位置nr (r - 1) % nnc (c 1) % n。如果(nr, nc)还没被填过就移动到(nr, nc)。否则移动到(r 1) % n列保持为c。这里% n就是越界回绕的数学化表达行从0再往上就回到n-1列从n-1再往右就回到0。这个递推式可以很方便地写成程序也是后面推导闭式公式的起点。3.2 碰撞为什么恰好每n步发生一次罗伯法最神奇的地方在于你会发现“被占位”这个例外不是随机出现的它非常有规律。以5阶为例碰撞发生在填6、填11、填16、填21的时候正好是每5步一次。3阶也是碰撞发生在填4、填7的时候正好每3步一次。这不是巧合。假如从头到尾没有碰撞所有数字都严格走“右上”方向那么从起点(0, c0)出发第t步的位置就是(r, c) ((-t) % n, (c0 t) % n)当t是n的倍数时显然(-t) % n 0、(c0 t) % n c0也就是回到了起点(0, c0)。起点已经被数字1占据了所以第n1个数字也就是tn的那一步必然发生碰撞。于是碰撞严格发生在t n, 2n, 3n, ...这些位置。这个规律是整个闭式公式能够推导出来的基石。3.3 为什么这样填出来的数字不重复要证明罗伯法有效首先要排除一个风险会不会填着填着某个候选位置被占向下一格之后下一步又撞到别的已填格子如果出现这种情况规则就失灵了。证明不重复其实不复杂思路可以这样看把所有t步分成若干段每段长度n。每段内部除了段首那一步是“向下”移动其余n-1步都是“右上”移动。这样每走完一段相当于整体行进的净效果是向右上走了n-1步再向下走1步。换算成坐标系净效果其实等价于行号增加2、列号增加n-1在模n意义下列号等于减1。也就是说每一段结束后的起始位置会整体平移。由于n是奇数这个“每段一次向下”的修正量不会与已有的路径产生周期重叠最终恰好把整个n×n网格遍历一次。严格证明需要用到群论里的循环群生成元但这篇文章不展开。直观上你也可以相信既然每一步都有唯一且可执行的去处而且总步数和格子数完全一致那么只要过程中没有任何一步走到死胡同就必然形成一个覆盖全图的回路。罗伯法的构造保证了死胡同不会出现。3.4 行列对角线为何必然等和这是整个公式法最核心的问题就算所有格子都填满了凭什么每行每列的和都相等这里有一个非常漂亮的证明方法理解它比死记口诀有用得多。把数字k1拆成高位和低位两部分k a * n b其中a k // n取值范围0到n-1相当于“第几组”b k % n取值范围也是0到n-1相当于“组内第几个”。数字本身是a * n b 1。注意幻方要求1到n²每个数出现一次而a和b的各n种取值组合恰好产生了n²个不同数字。目标就变成了证明在任意一行里n个数字对应的a值互不相同b值也互不相同。因为如果每行里a遍历0到n-1、b也遍历0到n-1那么这一行的数字和就是Σ(a*n b 1) n * Σa Σb n n * n(n-1)/2 n(n-1)/2 n n(n²1)/2恰好等于幻和。列和对角线同理。接下来要用到闭式公式里a和b与坐标的关系。先说结论对于罗伯法生成的方阵位置(r, c)上的数字满足一组模关系其中一个关键性质是对固定的行r当b从0变到n-1时a会以2a ≡ rb (mod n)的方式随之变化。因为2和n互质所以a和b是一一对应的b遍历0到n-1时a也就遍历0到n-1。列的证明同理。这条性质保证了行列和的恒等性不需要逐个数去求和一行公式就能算完。4. 闭式公式与代码实现从“一步一步走”到“直接算出坐标”4.1 序号到坐标的闭式公式有了递推式我们其实已经可以写程序了。但如果能更进一步把“第t步填到哪个格子”直接用一个公式算出来那就更方便了不需要维护“当前坐标”也不需要判断碰撞直接用序号算坐标。推导思路是这样从起点到第t步在没有碰撞的情况下应该是(-t, c0t)。而碰撞恰好发生在t n, 2n, ...所以第t步之前包括第t步这一步的碰撞次数是q t // n每一次碰撞相当于把原本“右上移动”的那一步改成了“向下移动”。比较一下两种移动右上移动是(-1, 1)向下移动是(1, 0)差值为(2, -1)。所以每发生一次碰撞最终位置就比无碰撞的理论位置多(2, -1)的修正量。设c0 n//2可以得到r (-t 2q) % nc (c0 t - q) % n其中q t // n。这个公式可以用5阶方阵验证几个点。比如数字1对应t0q0r0c2正确数字6对应t5q1r(-52)%52c(25-1)%51正确数字17对应t16q3r(-166)%50c(216-3)%50正确。多验几个点和手填结果完全一致。这个闭式公式的价值在于它可以随机访问任意一个数字的坐标而不需要从1开始一路填到它。比如你想知道99阶幻方里数字5000在哪个格子直接用这个公式一算就出来了。4.2 坐标到数字的反查公式反过来给定一个坐标(r, c)能不能直接算出这个格子里是什么数字这个在验证幻方时特别有用。推导思路是把t写成t q*n s其中q t // ns t % n。把上一节的两个公式展开r ≡ (-s 2q) mod nc ≡ (c0 s - q) mod n两式相加消掉s可以得到q ≡ (r c - c0) mod n有了q之后再代回去求ss ≡ (2q - r) mod n最后数字就是t1 q*n s 1。这一组公式可以直接用来做“坐标查值”。举个例子5阶方阵里坐标(1, 4)的数字c02q(14-2)%53s(2*3-1)%50t15数字是16和手填结果一致。4.3 三份可以直接用的Python代码理论讲完上代码。我写项目的时候比较喜欢把构造和验证分开这样后续改任何一边另一边也能独立测试。第一份代码是标准的递推版罗伯法def magic_square_siamese(n): if n % 2 0: raise ValueError(n must be odd) square [[0] * n for _ in range(n)] r, c 0, n // 2 for num in range(1, n * n 1): square[r][c] num nr, nc (r - 1) % n, (c 1) % n if square[nr][nc] ! 0: r (r 1) % n else: r, c nr, nc return square第二份是闭式公式版直接按坐标填数字def magic_square_formula(n): if n % 2 0: raise ValueError(n must be odd) c0 n // 2 square [[0] * n for _ in range(n)] for t in range(n * n): q t // n r (-t 2 * q) % n c (c0 t - q) % n square[r][c] t 1 return square第三份是验证函数检查行、列、两条对角线的和是否都等于幻和同时检查1到n²是否恰好出现一次def verify(square): n len(square) target n * (n * n 1) // 2 flat [x for row in square for x in row] if sorted(flat) ! list(range(1, n * n 1)): return False if any(sum(row) ! target for row in square): return False if any(sum(square[i][j] for i in range(n)) ! target for j in range(n)): return False if sum(square[i][i] for i in range(n)) ! target: return False if sum(square[i][n - 1 - i] for i in range(n)) ! target: return False return True如果想把坐标反查公式也用上可以这样写def cell_value_formula(n, row, col): c0 n // 2 q (row col - c0) % n s (2 * q - row) % n return q * n s 1用这个函数可以从方阵里任取一个坐标直接得到格子里的数字不需要遍历整个方阵。在生成大阶幻方后做抽样检查时这个函数能省不少时间。4.4 代码中的模运算陷阱闭式公式里有好几处模运算在不同编程语言里踩坑的可能性完全不一样。Python的%永远返回非负值所以(-5) % 5结果是0这符合数学上的模运算定义。但在C、C、Java里负数做%的结果是负数或0(-5) % 5在C语言里是0还算好但(-6) % 5在Python里是4在C语言里却是-1。如果你把公式直接翻译成C语言r (-t 2*q) % n在t较大时会因为负数取模问题直接算错。解决方案有两个一是先加若干个n再取模比如r ((-t 2*q) % n n) % n二是老老实实用无符号整数运算并确保中间结果非负。我在实际开发里吃过一次这个亏当时生成一个31阶幻方前30阶都正常第31阶对角线校验一直失败查了半天才发现是C语言的负数取模问题。这个问题非常隐蔽因为不是所有坐标都会出问题只有-t 2*q为负数的那些格子才会错。5. 变形、边界与实战避坑5.1 等价变换换起点、换方向、旋转镜像罗伯法生成的幻方不是唯一的。我们可以在不破坏性质的前提下对标准幻方做各种变换。最安全的一类是旋转和镜像把整个方阵旋转90度、180度、270度或者沿水平轴、垂直轴、对角线做镜像得到的新方阵依然是幻方。原因是旋转和镜像保持了“行、列、对角线”的对应关系只是把它们交换成了彼此的镜像。另一类等价变换是改变初始移动方向。标准罗伯法从第一行中间开始向右上走。你也可以改成向左上走只要每一步对称碰撞规则不变得到的仍然是幻方。甚至可以从最中间一列的下方开始从下往上走本质上只是把方阵旋转了180度。这些变体生成的结果通常能覆盖一部分“看起来完全不同”的幻方。但要注意不要随意把起点放在第一行的非中间列。起点改变后碰撞周期和路径覆盖性都会发生变化很可能走着走着提前撞墙最后行列和也不再相等。我曾经试过把起点放在第一行第一列3阶和5阶都无法生成完整幻方。5.2 另一种奇阶公式法斜排法阶梯法在聊奇阶幻方公式法的时候不能不提另一个历史悠久的构造方式——“斜排法”又叫阶梯法。南宋数学家杨辉在《续古摘奇算法》里记载的“九子斜排上下对易左右相更四维挺出”就是3阶幻方的斜排构造。这个方法推广到任意奇数阶的思路是先把1到n²按自然顺序从左上到右下斜向排成一个菱形这个菱形的宽和高都是2n-1中间重叠部分是n×n方阵然后把方阵边界外凸出的“小三角”整体平移到对侧的空位里得到的方阵就是幻方。这个方法和罗伯法在数学上有非常深的联系。有些资料称它们是同一个方法的不同呈现方式因为从最终结果看斜排法生成的幻方可以通过旋转、镜像与罗伯法生成的幻方建立对应关系。对于想要从数学上理解幻方结构的读者斜排法的几何直观更强它把“为什么数字要斜着走”直接画在了图里。5.3 常见的三个坑第一个坑是起点放错。很多人记住了“1放中间”但忘记是“第一行中间”随手放在第一列的中间。结果后面的碰撞周期完全乱掉。正确做法是用0-indexed坐标起点固定在(0, n//2)。如果是手算可以先用笔在纸上把第一行中间一格标出来。第二个坑是碰撞后处理方向搞混。新数字永远放在“当前数字”的正下方不是“右上角那个格子”的正下方也不是随机找一个空位。这个规则在代码里就是r (r 1) % n列保持不变。写程序时如果图省事写成r (nr 1) % n多半会得到错误方阵。第三个坑是在验证时只检查行和列漏掉对角线。幻方的定义里两条对角线也是必须等和的。罗伯法生成的方阵天然满足对角线等和但如果你的代码实现里某个细节错了一点很可能行列和依然碰巧相等唯独对角线差那么一两个数。我建议验证函数里把对角线和普通行的检查分开写这样出错时能更快定位。5.4 为什么偶数阶不能直接套这组公式奇阶幻方公式法看着很通用但不要试图把它直接套在偶数阶上。原因可以回到碰撞规律那一节碰撞在标准罗伯法里严格每n步发生一次这个结论依赖于“无碰撞情况下第t步回到起点”这个事实。偶数阶时从第一行中间出发右上移动n步之后不一定回到起点因为第一行中间的位置在偶数阶里没有“正中”的格子而是中间偏左或偏右。起步位置和回绕周期的配合关系被打破路径会在中途提前冲突导致整个构造失败。偶数阶幻方有自己的一套方法双偶数阶4、8、12等常用“对称交换法”单偶数阶6、10、14等常用“斯特雷奇法”这两套方法又是完全不同的逻辑不适合硬套罗伯法。所以拿到一个n时第一件事永远是判断奇偶性再决定用哪套公式。最后分享一个我实际用下来的经验。如果只是手算一个奇数阶幻方闭式公式反而不如递推口诀直观因为手算时连续取模很容易算错但如果你是在写程序、做大阶幻方闭式公式就非常爽了它支持随机访问任意坐标的数字生成后可以直接抽样验证不用等整个方阵填完。我现在的习惯是写代码时用闭式公式生成然后写一个独立的验证函数再随机抽几个坐标用反查公式交叉核对。这套组合拳我用了很多年几乎没出过错。