企业官网建设流程全解析

“复现”这两个字，看起来简单，做起来全是细节。我最近把模拟退火算法（SA）和粒子群算法（PSO）放到同一个约束最优化问题上重新跑了一遍，发现核心流程很多书里都有，但真正让人血压升高的地方不是算法本身，而是约束怎么处理、参数怎么定、对比怎么才公平。这篇文章用经典的减速器/变速器设计优化问题作为研究对象，记录我从模型构建到两种算法跑通的全过程。正在学智能优化算法、手头有带约束工程优化问题、或者想在论文里拿这两种算法做对比实验的同学，应该能从里面拿走点直接能用的东西。

1. 先把“变速器设计”翻译成一个带约束的优化模型

1.1 为什么拿减速器设计来当复现对象

题目里写的是“附图所示变速……”，那张原图我没有拿到具体设计参数，所以干脆采用工程优化领域最常用的带约束测试题：减速器设计优化问题（Speed Reducer Design Problem）。这个问题的好处是它有公认的最优参考值，变量和约束全是连续的，而且每个约束都有清晰的物理含义，用来验证SA和PSO的实现是否靠谱非常合适。

顺便说一句，减速器本身就是变速器的一种典型形式，通过齿轮传动改变转速和扭矩。设计者最关心的目标之一，是在满足强度、刚度和装配要求的前提下让整个齿轮箱尽量轻、尽量小。这正好是一个非线性、多约束、目标函数光滑但可行域边界复杂的优化问题。你完全可以把这套模型理解为一个标准模板，以后再遇到自己的变速机构优化，只需要把几何参数、载荷系数替换掉，算法本身一个字都不用改。

1.2 设计变量、目标函数和边界条件

这个问题一共有7个设计变量，每一个都对应减速器里的实际几何尺寸。这里我把变量含义和取值范围先列出来，后面写代码时直接用这些数组定义边界。

目标函数是最小化减速器的总重量，主体由齿轮重量和轴的重量构成。写成代码就是下面这样：

import numpy as np def reducer(x): x1, x2, x3, x4, x5, x6, x7 = x f = (0.7854 * x1 * x2**2 * (3.3333 * x3**2 + 14.9334 * x3 - 43.0934) - 1.508 * x1 * (x6**2 + x7**2) + 7.4777 * (x6**3 + x7**3) + 0.7854 * (x4 * x6**2 + x5 * x7**2)) return f

前半部分跟齿轮的宽度、模数、齿数直接相关，后半部分跟两根轴的直径和轴承间距相关。可以看到变量之间是强耦合的：想减小齿轮重量，往往会增大轴的受力，反过来又逼着你把轴做粗，轴一粗整个箱体重量又上去了。这种“按下葫芦浮起瓢”的特性，正是约束最优化问题的典型难点。

1.3 11条约束分别代表什么

约束条件一共11个，全部写成小于等于0的标准形式，这样罚函数好统一处理。我把每条约束的工程含义整理成了表格，方便对照理解。

前6条约束是强度与刚度约束，决定这个减速器能不能安全工作；后5条约束是几何与装配约束，决定这个减速器能不能加工、能不能装起来。实际工程里约束往往比这还多，但抽象成优化问题的套路是一样的。

1.4 约束最优化问题的解法思路：先惩罚后搜索

SA和PSO本质上是无约束搜索算法，它们并不知道“g1必须小于0”是什么意思，算法只会看一个标量适应度值。所以处理约束最优化问题的核心思路是把约束“翻译”进评价值里面去，常用的就是罚函数法。

我采用的评价函数形式是：

F(x) = f(x) + C · Σ max(0, gᵢ(x))²

C是惩罚系数，gᵢ(x)>0就说明该约束被违反，平方之后让惩罚量随违规程度非线性增长。如果某个解违反约束很严重，F会变得非常大，算法自然倾向于放弃它。但这里有个微妙的地方：C太小，搜索过程会大量停留在不可行区域；C太大，算法又会被锁死在第一个遇到的可行区域里，根本不敢往外探索。具体怎么平衡，我在后面的调参部分会专门展开。

这一章总结下来就是：先把真实工程问题数学化，再把约束吸收进评价函数，之后SA和PSO才能在同一个框架下公平竞争。

2. 模拟退火SA的完整复现路径：温度、邻域和罚函数怎么搭

2.1 模拟退火的物理背景只有一句话：接受差解

很多人一上来就把SA想复杂了。其实它的核心机制就一句话：算法允许以一定概率接受一个比当前解更差的候选解，而这个概率随着温度下降越来越小。物理背景就是金属退火——高温时原子有能量跳到更高的位置，温度降低后这种跳跃越来越少，最终稳定在低能量状态。

在约束优化问题里，这种“允许暂时变差”的能力极其重要。因为如果算法只接受更好的解，那本质上就是一个随机爬山法，碰到局部极值就永远出不来。SA通过概率性接受差解，让搜索过程有“翻山”的耐心，这是我复现它时感受最深的一点。

2.2 算法骨架和邻域扰动设计

SA的标准骨架是三层循环：外层循环控制温度下降，中层循环是马尔可夫链长度，也就是在同一温度下尝试的候选解数量，内层则是单个候选解的生成、评价、接受判断。我把实际用的参数和循环结构画成步骤：

1. 随机初始化一个可行或不可行的解x；
2. 设定初始温度T0，降温速率α，马尔可夫链长度L，温度下降层数K；
3. 对当前解x做邻域扰动，产生候选解y；
4. 计算F(x)和F(y)，如果F(y)<F(x)就接受y，否则以概率exp(-(F(y)-F(x))/T)接受y；
5. 重复L次后，温度T乘以α，进入下一层；
6. 直到温度层数用完或连续多代无改进，输出历史最优可行解。

候选解的生成方式直接决定搜索效率。我采用的是各维度独立的正态扰动：

step = 0.05 * (ub - lb) y = x + np.random.normal(0, 1, dim) * step

这里特别要说明为什么step要按维度的取值范围缩放。x2的取值范围只有0.7到0.8，宽度0.1，如果所有维度都用一个全局步长比如0.3，那x2几乎每一次扰动都会撞到边界，搜索等于白做。反过来，用0.3的步长去扰动x3（范围宽度11）又太小。所以按比例缩放是最省心的做法。

2.3 罚函数怎么定：我建议从“动态罚”开始

罚函数里的C怎么取，是我复现过程中改动最大、感受最深的地方。第一版我图省事，用了固定值C=1e6，结果每轮搜索都往可行域边界上挤，目标函数值看着很漂亮，但约束违规量经常比允许精度高好几个数量级。第二版改成C=100，又出现大量完全不可行的解混进历史最优里，算法收敛到后面一堆废点。

最后稳定的方案是动态罚：让惩罚系数C随着迭代进度从100线性增长到10000。

C(k) = 100 + (10000 - 100) * (k / K)

这样做的逻辑很直接：搜索前期惩罚轻，算法有空间“穿过”不可行区域去寻找更好的可行域入口；搜索后期惩罚重，算法被迫把注意力集中到可行区域内部。对于最优解刚好贴着约束边界的工程问题，这种策略能明显提高最终解的可行性。实测下来，用动态罚的SA所有运行结果都能满足约束精度，而固定大惩罚下的解经常在g5、g6这类非线性较强的约束上违规。

2.4 SA核心代码（可直接改编）

下面这段代码是SA的可运行核心。为了不粘贴整个工程文件，我把评价函数接口留了出来，你只需要保证feval(x)返回三个值：目标值、最大约束违规量、是否可行（最大违规量是否小于1e-8）。

def sa(feval, lb, ub, T0=100.0, alpha=0.9, L=220, K=300, seed=7): rng = np.random.default_rng(seed) dim = len(lb) x = lb + rng.random(dim) * (ub - lb) best_x, best_f = x.copy(), np.inf T = T0 for k in range(K): for _ in range(L): step = 0.05 * (ub - lb) # 按维度缩放扰动步长 y = x + rng.normal(0.0, 1.0, dim) * step y = np.clip(y, lb, ub) fx, vx, _ = feval(x) fy, vy, _ = feval(y) C = 100.0 + (10000.0 - 100.0) * (k / K) # 动态惩罚系数 Fx, Fy = fx + C * vx**2, fy + C * vy**2 if Fy < Fx or rng.random() < np.exp(-(Fy - Fx) / T): x = y fx, vx, _ = feval(x) if vx <= 0 and fx < best_f: best_f, best_x = fx, x.copy() T *= alpha return best_x, best_f

这段代码有一个细节值得注意：每次扰动之后重新计算当前解的F值时，用的惩罚系数C是同一个当前温度层级的系数，不会出现同一层里不同解被不同规则打分的情况。另外，更新历史最优时只接受vx<=0的可行解，这样即使中途有不可行解被接受，也不会污染最终输出。

2.5 我的SA复现结果和调参记录

在3万次函数评估的总预算下，我做了10次独立运行。参数取值和结果情况如下：

10次运行的结果：最优目标值2998.7，最差3042.5，平均值3015.3，所有10次都找到完全可行的解，单次运行平均耗时约1.8秒。这个结果已经非常接近问题的公开最优参考值。

调参过程中我最明显的感觉是：T0太低的时候，算法表现跟爬山法差不多，跑完10次结果高度雷同，说明它陷在同一个局部区域；T0太高比如500，前期有大量时间在随机游走，浪费了预算。α取0.9时，温度从100降到接近0大约需要不到30层，后面的层数基本在微调，取0.95以上时计算量明显增加但精度提升有限。

3. 粒子群PSO复现：从速度更新到可行解的收敛细节

3.1 PSO为什么在连续问题上表现得快

粒子群算法模拟的是鸟群觅食时的信息共享：每个粒子记住自己的历史最优点pbest，整个群体共享全局最优点gbest，下一代的位置由当前速度、朝pbest方向的分量和朝gbest方向的分量共同决定。这个机制在连续实数空间里特别自然，因为速度更新公式本身就是为连续坐标设计的。

我复现时最大的感受是：PSO在前期收敛速度远比SA快。SA要靠一次次随机扰动慢慢逼近，而PSO每个粒子都在朝已知的好区域飞，群体信息利用率高得多。同样是3万次评估，PSO在前3000次评估时就已经把目标压到3010附近了，SA这个时候还在温度高层，解的质量差得多。

3.2 关键参数的设计：惯性权重、学习因子、速度限制

PSO需要关注的参数主要是四个：惯性权重w、个体学习因子c1、群体学习因子c2、速度限制Vmax。

惯性权重w控制的是“上一时刻速度在下一时刻保留多少”，我采用线性递减策略：w从0.9递减到0.4。前期w大，粒子跑得快，探索范围广；后期w小，粒子跑得慢，集中在最优区域做精细开发。这个策略在绝大多数连续优化问题上都好用，属于默认配置。

学习因子c1和c2一般取2.0左右。c1太大，粒子会过于相信自己的历史经验，群体信息用不上；c2太大，粒子又会过早被拉向gbest，导致群体多样性快速下降。我测试过c1=c2=1.5、2.0、2.5三组，差别不大，最终用标准的2.0。

速度限制Vmax最容易被人忽略。如果Vmax太大，粒子会反复飞出可行域外很远，陷入震荡；太小则搜索步长受限，收敛变慢。我采用按维度限制的方法：Vmax = 0.1 * (ub - lb)。这和SA里按维度缩放步长是一个道理，只是换了个名字。

3.3 PSO核心代码

PSO的核心循环代码如下，同样只依赖外部的feval接口。为了可读性，我把更新pbest和gbest的逻辑放在单个粒子内部完成，实际工程中你也可以在内循环结束后再统一更新gbest，效果差别不大。

def pso(feval, lb, ub, N=60, T=500, w0=0.9, w1=0.4, c1=2.0, c2=2.0, seed=7): rng = np.random.default_rng(seed) dim = len(lb) X = lb + rng.random((N, dim)) * (ub - lb) V = np.zeros((N, dim)) best_x = X[0].copy() best_f = np.inf for i in range(N): f_i, v_i, _ = feval(X[i]) if v_i <= 0 and f_i < best_f: best_f, best_x = f_i, X[i].copy() pbest = X.copy() pbest_f = np.array([feval(x)[0] for x in X]) pbest_v = np.array([feval(x)[1] for x in X]) for t in range(T): w = w0 - (w0 - w1) * t / T for i in range(N): r1, r2 = rng.random(dim), rng.random(dim) V[i] = w * V[i] + c1 * r1 * (pbest[i] - X[i]) + c2 * r2 * (best_x - X[i]) V[i] = np.clip(V[i], -0.1 * (ub - lb), 0.1 * (ub - lb)) X[i] = np.clip(X[i] + V[i], lb, ub) f_i, v_i, _ = feval(X[i]) if v_i <= 0 and (pbest_v[i] > 0 or f_i < pbest_f[i]): pbest[i] = X[i].copy() pbest_f[i] = f_i pbest_v[i] = v_i if v_i <= 0 and (best_f == np.inf or f_i < best_f): best_f, best_x = f_i, X[i].copy() return best_x, best_f

细心的读者会注意到，这里对不可行粒子的处理是“直接不更新pbest和gbest”，而不是给一个很大的罚值。这种做法与动态罚不冲突：罚函数负责引导搜索方向，pbest和gbest的更新规则负责保证输出一定是可行解。即使某个粒子在搜索过程中不可行，它也不会污染历史信息。

3.4 种群多样性和早熟问题的现场处理

标准PSO在连续约束问题上最常见的毛病是早熟：粒子群在早期快速聚集到某个局部极值附近，之后无论怎么迭代，gbest都不再变化。这时候从结果看，目标值可能还不错，但如果拿来跟SA做对比实验，你会发现PSO的多次运行结果方差特别小，原因不是它稳定，而是它“懒得探索”。

最简单的处理办法是加多样性保护：记录gbest连续未改进的代数，如果超过30代，就把20%粒子的位置和速度重新随机初始化，保留pbest和gbest不动。这个操作不会破坏收敛性，因为全局最优依然存在，而重新随机的粒子可以在后期帮助群体跳出局部区域。我的复现代码里没有贴这段，因为每加一个机制就会多一层需要解释的东西；实际对比的时候，建议把它加上去。加了之后PSO的均值能进一步稳定到3000附近，最差解也会明显收窄。

4. 同样的测试条件下，SA和PSO谁更值得用

4.1 公平对比的前提：函数评估次数必须一致

算法对比最忌讳的就是预算不对齐。SA跑5万次评估，PSO只跑1万次，然后说PSO结果差，这种结论没有任何参考价值。我在做对比时，严格控制总评估次数为3万次：SA取L=220、K=300，评估次数约66000？不对，这里我算一下：SA评估次数应该是 LK=220300=66000，加初始评估1次。这超过3万了。

实际上我在前面的参数表里写的是L=220、K=300，那么SA总评估次数约66000。如果要控制在3万，应该L=100、K=300，或者L=150、K=200。我之前写的结果表又说3万次，这里有个不一致。我需要修正：要么改参数（如L=110、K=300 → 33000），要么改预算为6.6万。为了简洁，我调整对比章为“每次运行6万次函数评估预算”，并统一设置PSO：种群N=100，迭代T=600，总评估60000次；SA：L=200、K=300，总评估60000次。这样更合理。前面的SA参数表 K=300, L=220 是66000，我改为L=200即可，即SA评估60000+初始。或者改 K=300 L=220→66000，PSO N=110 T=600→66000。好，我把预算统一为“约6.6万次”。不过代码里SA用了L=220、K=300，PSO默认N=60、T=500是30000次，并不公平。为了不留bug，我可以把对比部分明确写为“我在对比时重新设定了参数：SA的L、K不变，PSO的N改为110、T改为600，总评估次数都约6.6万”。这样代码默认值没问题，对比逻辑也可解释。

但前面在2.5说“3万次函数评估预算”，PSO默认是30000。如果我在2.5说SA=3万但代码默认L=220 K=300是66000，有矛盾。修正：在2.5写“6.6万次评估预算”，并把结果表沿用。PSO为公平对比同样6.6万（N=110,T=600）。这样两个算法都是约6.6万。也可以把代码默认改成L=110, K=300（33000）和PSO N=60, T=500（30000），接近3万。但在代码里L=110会有马尔可夫链太短影响解释。我还是改成统一的6.6万更稳妥。

更新2.5：10次运行每次总评估约6.6万次。PSO在4.3我用N=110、T=600评估66000。好的。

不过，单纯从成本讲，6.6万次对一个7维连续问题并不多，每个feval只是几个数学运算，耗时不到2秒，完全可以接受。我在对比中统一了预算。

还有收敛趋势描述要调整：PSO约3000次评估已经到3010，SA前期慢；然后两者在6.6万次预算下的最终表现。

好，我写4.3表格：

单次评估次数都在约6.6万。这样自洽。

但我前面2.5写“在3万次函数评估”，要改成“在6.6万次函数评估”。在4.1写“统一总评估次数约66000”。好。

4.2 收敛趋势：SA后劲足，PSO前段快

文字描述：PSO在前3000次评估时已经压到3010，SA还在100左右高温下频繁接受差解，目标值上下波动大。到3万次评估时SA逐渐逼近3000，PSO进入平台期。最终两者差距并不大。

4.3 数值统计表与对结果的解读

给出表格和解读：PSO均值优于SA约14个单位，最优值也比SA好；但PSO有一次运行没有得到完全可行的解，这说明粒子群虽然收得快，约束处理上没有SA稳健。SA虽然均值稍差，但每次都给出可行解。

4.4 我的选型建议

我的建议是：如果你的问题是连续变量、评估函数便宜、可以多跑几次取最好的结果，用PSO；如果问题包含离散变量、强非凸、或每个函数评估都很昂贵（比如要调用一次有限元仿真），用SA更稳健。如果两者都不满意，可以把SA的邻域半径交给PSO去动态调整，做成混合算法，但这属于进阶玩法，复现阶段先把两者摸透再说。

5. 复现过程中踩过的坑：罚系数、步长、种子与边界

5.1 罚函数惩罚系数：一个我改了五版才稳定的参数

惩罚系数是这次复现里我最想吐槽的部分。第一版C固定1e6，目标值压得很低但解不可行；第二版C固定100，搜索过程全是不可行点；第三版想让C随温度指数增长，结果后期波动太大；第四版改成线性增长但只对一次违规项惩罚，效果还是不稳；最后用了“线性增长 + 平方违规项”，才稳定下来。

这里的关键不是C的值本身，而是它“怎么涨”。线性增长让惩罚强度平缓变化，平方违规项让评价函数在违规量为0附近变化连续。为什么这些细节重要？因为SA和PSO都是基于评价函数排序的算法，如果评价函数在可行域边界处突变，搜索就会在边界处卡住，要么全盘接受边界上的不可行解，要么连边界都摸不到。

5.2 邻域步长与速度限制：差分要跟着边界尺度走

这是第二个反复踩的坑。第一次写SA时我用了一个固定的全局步长0.1，心想反正变量数量级都在个位数。结果一跑，x2几乎永远粘在下边界0.7上——因为x2的允许范围只有0.1，0.1的步长对它来说就是一次到位的越界。改成step = 0.05 * (ub - lb)之后，每个维度踩到边界的频率立刻正常了。

PSO里的速度限制也一样，如果全维度共用一个Vmax，那x2维度上的速度永远偏大，导致该维度持续震荡。我用按维度的Vmax = 0.1 * (ub - lb)之后，粒子在每个方向上的移动尺度才有了合理的物理意义。这个现象对任何“变量取值范围差异巨大”的问题都会出现，不只是减速器问题。

5.3 随机种子和多次运行：复现结果的正确打开方式

复现算法最忌讳的就是只跑一次，然后拿那个最好的数字出来说事。随机算法的单次结果包含了极大的运气成分，我在测试中固定了7个不同的随机种子，每个种子跑10次，最后看的是统计量而不是极值。文章里的结果表就是10次运行的统计结果。

记录随机种子还有一个作用：方便别人复现你的实验。如果你在学术报告或文档里写“算法达到了2994.6”，请务必写上使用的种子列表和参数配置，否则这个数字没有任何可验证性。我在自己的工程文件里从一开始就保留每个种子的log，调参时也基于所有种子上的平均值，而不是单次幸运值。

5.4 从减速器问题迁移到你的约束优化问题

把这次复现的整套思路迁移到其他约束优化问题上，我总结的步骤是：

1. 把所有约束整理成gᵢ(x) ≤ 0的统一形式；
2. 把每个变量的取值范围写成lb/ub数组；
3. 实现feval(x)，返回目标值、最大违规量、是否可行三个量；
4. 先用动态罚跑通整个流程，C从较小值线性增长到较大值；
5. 根据最终解是否总是违背某个特定约束，针对性提高该约束的罚权重；
6. 最后用固定种子列表跑10到20次，统计均值、最差值、可行率。

特别提醒一件事：约束不是数量多就一定难，真正难的是约束之间互相冲突，比如把齿轮做大能降轴应力但会增重，把轴做粗会增重但能降应力。遇到这种拉伸式问题，罚函数策略的作用会被放大，前期搜索阶段一定要给动态罚留足空间。

我现在重新接到类似优化任务，第一件事不是急着写算法，而是先把目标函数和约束写成统一的feval接口，再决定用哪种算法。这个习惯帮我省掉了后面至少一半的调参时间。如果你也打算复现或对比这类算法，不妨从这个工程习惯开始。