简介这是一份面向智能优化课程设计与多目标优化学习者的完整代码方案使用Python语言实现非支配排序遗传算法第二代NSGA-II并专门对接CEC-2021竞赛中的典型多目标优化问题既可作为课程设计提交的参考工程也可作为算法入门的实战模板。整个压缩包共包含185个文件体积约969KB其中12个Python脚本负责主程序、算法核心、问题定义与工具函数100个mat文件存放迭代过程中的种群和结果数据50个txt文件记录运行输出另有11张png图片用于可视化展示以及少量xml配置、markdown说明和辅助的m评估脚本目录组织清晰方便按需检索。目前该代码包已有125人学习下载。借助这份工程读者可以对照代码体会选择、交叉、变异、非支配排序、拥挤距离计算与精英保留等关键步骤也能直接运行CEC-2021赛题进行测试再利用超体积HV等评估脚本分析解的收敛性与分布性适合在校本科生、研究生以及备战优化竞赛的开发者。1. 一个课设 zip 里的 NSGA-Ⅱ凭什么能解决 CEC-2021如果你正在准备智能优化课设大概率会遇到这类包名Python 实现 NSGA-Ⅱ题目是 CEC-2021 竞赛里的多目标问题。它不是一个玩具 demo而是把 NSGA-Ⅱ 的完整流程——非支配排序、拥挤度距离、精英保留、SBX 交叉、多项式变异——用 Python 和 NumPy 落成可以跑出实验图的代码。多数课设交的内容是一段能出结果的主程序、几个测试函数定义、一组收敛曲线和帕累托前沿图再加一份讲清楚原理和参数的报告。这个包适合正在修进化计算课程、需要交代码和实验数据的学生也适合想拿现成多目标基准快速验证自己想法的新手。它解决的问题很具体不依赖 MATLAB 工具箱用最直白的方式看到 NSGA-Ⅱ 在 CEC-2021 问题上的收敛行为和多样性表现。2. NSGA-Ⅱ 的三块基石把“哪些解该留下”从算法上拆清楚NSGA-Ⅱ 本质上在反复回答一个问题当前种群里有那么多解哪些应该保留到下一代。答案由三个机制拼出来快速非支配排序决定解的等级拥挤度距离决定同一等级里保留谁精英保留机制保证父代里真正好的解不会被子代覆盖。三者叠起来就是“先分好坏再分稀疏程度好的与更有特点的都留下”。如果你只是把别人的代码跑通就交差一旦被问到“为什么这段代码要先算 rank 再算 crowding_distance”很容易卡壳下面把每一块都拆开讲。2.1 快速非支配排序把种群按“支配等级”切层多目标优化里最基本的比较是支配。对全部目标都做最小化的问题解 A 支配解 B 的条件是A 在所有目标上不劣于 B且至少在一个目标上严格优于 B。用 Python 写这个判断只需要几行import numpy as np def dominates(a, b): 目标向量 a 是否支配 b全部维度 且至少一个维度严格 。 return np.all(a b) and np.any(a b)快速非支配排序的核心思想是“逐层剥洋葱”先把当前集合里不被任何解支配的个体全部拿出来标为第 1 层从集合里剔除它们再找剩余个体中不被支配的个体标为第 2 层循环直到所有个体都有了层级编号。注意“不被任何解支配”和“不被其他解支配”的区别前者指的是全局非支配个体。为什么不能只保留第 1 层因为第 1 层里可能挤着几十上百个非支配个体它们彼此之间无法用支配关系区分必须交给拥挤度距离排序。层级编号本身还会影响后续选择rank 小的个体在二元锦标赛里永远优先被选中这是 NSGA-Ⅱ 收敛性的来源之一。这个步骤复杂度最坏是 O(MN²)M 是目标数N 是种群大小。课设里常见 M2 或 3N100~300直接双重循环完全能接受。不要为了常数级优化把代码改成复杂的数据结构能看懂、能复现比微优化更重要。2.2 拥挤度距离同一层内保留“最有特点”的解当某一层级的人数比下一代剩余名额多时必须从这一层淘汰一部分解。淘汰规则可以直观理解为在目标空间里距离相邻个体越远的解越“稀疏”越值得保留挤成一团的解则优先淘汰。具体算法是按每个目标分别排序并累加距离。def crowding_distance(fitness, front_indices): 计算给定前沿层个体的拥挤度距离返回与 front_indices 对应的数组。 n_obj fitness.shape[1] dist np.zeros(len(front_indices)) for m in range(n_obj): idx sorted(front_indices, keylambda i: fitness[i, m]) dist[idx[0]] np.inf dist[idx[-1]] np.inf if len(idx) 2: continue f_min fitness[idx[0], m] f_max fitness[idx[-1], m] span f_max - f_min if span 0: continue for k in range(1, len(idx) - 1): dist[idx[k]] (fitness[idx[k 1], m] - fitness[idx[k - 1], m]) / span return dist边界个体人为赋 np.inf是所有 NSGA-Ⅱ 实现里必须一致认同的约定。如果边界个体不是无穷大经过多代筛选前沿两端会被不断淘汰最终帕累托前沿的两个端点永远不会被找到。这也是一个高频踩坑点后面避坑清单会专门讲。另一个容易忽视的细节是每个目标上的距离要做归一化。不同目标的数值范围可能差几个数量级不除以 span 的话量纲大的目标会主导拥挤度的计算。上面代码里(fmax - fmin)做了归一化如果没有这段多目标问题里某个目标范围特别大时拥挤度会退化成只看那一个目标。2.3 精英保留与二元锦标赛这一代如何衔接下一代选择父代时用的是二元锦标赛从当前种群中随机抽两个个体优先选 rank 较小者rank 相同选拥挤度更大者。拥挤度大的个体周围更空旷选它有利于维持探索性。精英保留发生在子代生成之后。比较有代表性的做法是把父代和子代合并成规模 2N 的集合按 rank 从小到大的顺序填充下一代先填满第 1 层再填第 2 层直到某一层人数超过剩余名额就在这一层内按拥挤度距离从大到小截断。这么做的好处是父代里最优秀的前沿个体不会因为随机育种而整体丢失——这就是 NSGA-Ⅱ 被称为“精英保留”策略的原因。一个需要留意的实现细节填充子代时用来排序的 fitness 必须在合并后的集合上重新计算一次如果你沿用旧个体对象里的 rank 属性合并后会乱套。正确做法是每代对 2N 个个体重新做非支配排序和拥挤度计算再截断。2.4 CEC-2021 的多模态设定把多样性从加分项变成硬指标CEC-2021 的多目标竞赛题目里有不少是多模态问题典型特征是一个目标向量对应多组决策变量帕累托前沿由相互分离的多个子前沿组成。算法如果只守住一条前沿段可能 IGD 值看着还行但实际漏掉了一半的全局最优段。对这类问题拥挤度距离只能保证“同一段前沿内分布均匀”保证不了“不同子前沿都被找到”。所以实验上必须有多次独立运行每次只找到其中一段平均指标也会暴露这一点。很多人跑完单次实验看到一组漂亮曲线就下结论这是多模态问题里最典型的误判后面第 5 章会把复现实验的具体要求展开。3. 把 CEC-2021 问题接进 Python选函数、写函数、配参数3.1 CEC-2021 这个问题集课设应该怎么选CEC-2021 不是一张函数表而是一个按赛道划分的竞赛集合。课设中最常被引用的多目标赛道是带复杂形态的问题组其中包含多模态多目标问题 MMF 系列MMF1~MMF12 这类编号具体以当年公布为准以及带约束的 CMOP 系列。它们的共同特点是目标函数可解析、单次评估开销低、决策变量维度不高适合在课设周期内反复跑实验。我一般会建议选 4 个左右问题构成实验集而不是只跑一个一个结构简单的基准问题做 sanity check两个 MMF 类多模态问题做多样性验证有余力再选一个带约束的问题。这样报告里可以写“算法在不同问题形态下的表现差异”比单跑一个问题的说服力强太多。3.2 用纯 NumPy 写目标函数结构与边界CEC-2021 官方题目有精确方程组但课设里你可以先写一个与 MMF 结构相似的教学版函数跑通流程后再换成官方定义。下面这个双目标函数用于验证算法在“多个帕累托子前沿”上的表现决策变量 x1 和 x2 都在 [0, 1]def mmf_like_problem(x): 教学版多模态双目标问题结构与 MMF 系列相似。 x: 长度为 2 的决策向量两个目标都做最小化。 x1, x2 x f1 x1 # 多模态项x2 的周期性会让同一个 f1 对应不同质量的 f2 g 2.0 - np.abs(np.sin(4.0 * np.pi * x2)) ** 0.5 f2 g * (1.0 - np.sqrt(f1)) return np.array([f1, f2])注意这段代码只是教学示意不是 CEC-2021 官方公式。真正的竞赛题有精确的变量边界、偏移参数和约束条件写进报告时必须引用官方定义。用这个函数做开发调试的好处是它明显存在多个子前沿NSGA-Ⅱ 如果实现正确最终前沿会呈现出多个分离的段如果实现有问题很容易看出多样性丢失。写目标函数时另一个常见误用是把所有目标当最大化处理。NSGA-Ⅱ 的标准形式通常假定全部目标最小化如果你的题目里某个目标是利润、收益这类越大越好的量要么取负号要么把支配判断里的比较符号统一改掉。混用是算法表现莫名其妙的头号原因。3.3 参数配置哪些值能直接抄哪些必须自己试课程设计不追求极致调参但合理的默认参数能让结果稳定很多。下表是一组经过大量课设验证的基准值参数推荐取值说明种群大小100~300100 起步MMF 类问题 200 更稳最大代数500~1000以目标函数评估次数为预算更科学交叉概率 pc0.9SBX 交叉概率高不等于所有维度都交叉变异概率 pm1 / n_var每个维度独立的变异概率SBX 分布指数20控制子代与父代的接近程度多项式变异分布指数20控制变异步长独立运行次数20~30统计平均和方差单次结果不算数随机种子固定一个数保证结果可复现交叉概率 0.9 是 NSGA-Ⅱ 文献里的经典值它对很多问题都足够稳定变异概率取 1/n_var 保证平均每个个体约有一个维度发生变异这是保持探索性的底线。分布指数 eta 越大产生的子代越贴近父代收敛快但探索弱eta 越小搜索越跳跃。如果你希望把“评估次数”写进报告可以设定总预算为种群大小乘以代数例如 200 个个体跑 500 代就是 100000 次评估。CEC-2021 这类竞赛问题允许的评估预算通常是公开的按官方预算执行会让你后面写可比性论述时更站得住脚。4. 从零拼出一个可运行的 NSGA-Ⅱ 最小系统4.1 个体、初始种群与评价函数直接一步到位的做法是用 numpy 二维数组表示种群每一行是一个个体每一列是一个决策维。这种方法代码量最少也方便向量化def initialize_population(pop_size, lb, ub): 生成初始种群。 lb, ub: 各决策变量的下界和上界数组。 n_var len(lb) rng np.random.default_rng(seed20240401) pop lb rng.random((pop_size, n_var)) * (ub - lb) return pop def evaluate_population(pop, problem): 对种群中每个个体逐一评估目标函数。 fitness np.array([problem(ind) for ind in pop]) return fitness这里把随机数生成器单独拿出来是有意为之。自己实现 NSGA-Ⅱ 时不要到处调用np.random.random()而是统一使用一个rng实例。这样只要改一个 seed整条实验链都能复现否则每次运行结果不同报告里的图根本没法解释。4.2 合并、排序、截断完整实现非支配排序与拥挤度第 2 章已经给了 dominates 和 crowding_distance 的片段这里把它们串成一个完整的“合并选下一代”函数。它接受合并后的 fitness 数组返回下一代个体索引def select_next_generation(fitness): 对合并种群执行快速非支配排序再按层级填充。 返回下一代在合并种群中的索引列表。 n_pop fitness.shape[0] dominated_count np.zeros(n_pop, dtypeint) dominated_list [[] for _ in range(n_pop)] # 第一轮统计每个个体被谁支配 for i in range(n_pop): for j in range(n_pop): if i j: continue if dominates(fitness[i], fitness[j]): dominated_list[i].append(j) elif dominates(fitness[j], fitness[i]): dominated_count[i] 1 # 找出第 1 层 front [] rank np.zeros(n_pop, dtypeint) for i in range(n_pop): if dominated_count[i] 0: rank[i] 1 front.append(i) # 逐层剥离 all_fronts [] while front: all_fronts.append(front) next_front [] for i in front: for j in dominated_list[i]: dominated_count[j] - 1 if dominated_count[j] 0: rank[j] len(all_fronts) 1 next_front.append(j) front next_front # 按层级填下一代层内按拥挤度截断 next_indices [] for front in all_fronts: if len(next_indices) len(front) n_pop // 2: next_indices.extend(front) else: remain n_pop // 2 - len(next_indices) if remain 0: break dist crowding_distance(fitness, front) # 在 front 内按拥挤度从大到小取 remain 个 order np.argsort([dist[i] for i in front])[::-1] for k in range(remain): next_indices.append(front[order[k]]) break return np.array(next_indices)这个实现有几个点需要解释。第一n_pop // 2假设父代和子代各占一半合并后共 2N 个个体选回 N 个。第二从第 1 层开始逐层填充跨层时一旦剩余名额不足先算整层的拥挤度再按拥挤度降序取前 remain 个。第三传入的 fitness 必须是合并种群的目标值而不是某个个体对象上残留的旧值。如果你觉得这段代码的嵌套循环看着复杂说明你的方向感是对的。这个函数的正确性直接决定整个算法质量建议写完以后用一个已知结果的小例子做单元测试比如构造 4 个双目标个体手算一遍支配关系再和函数输出对照。4.3 SBX 交叉与多项式变异两个必须亲手写的算子NSGA-Ⅱ 文献里默认使用模拟二进制交叉 SBX 和多项式变异。下面是可供课设直接使用的实现逐维处理保证每个决策变量都受到边界约束def sbx_crossover(p1, p2, lb, ub, eta_c20): 模拟二进制交叉p1、p2 是两个父代个体。 n len(p1) c1, c2 np.empty(n), np.empty(n) for i in range(n): if np.random.random() 0.5: c1[i], c2[i] p1[i], p2[i] continue if abs(p1[i] - p2[i]) 1e-12: c1[i], c2[i] p1[i], p2[i] continue if p1[i] p2[i]: y1, y2 p2[i], p1[i] else: y1, y2 p1[i], p2[i] beta 1.0 2.0 * min(y1 - lb[i], ub[i] - y2) / (y2 - y1 1e-10) alpha 2.0 - beta ** (-(eta_c 1.0)) rand np.random.random() if rand 1.0 / alpha: beta_q (rand * alpha) ** (1.0 / (eta_c 1.0)) else: beta_q (1.0 / (2.0 - rand * alpha)) ** (1.0 / (eta_c 1.0)) c1[i] 0.5 * (y1 y2 - beta_q * (y2 - y1)) c2[i] 0.5 * (y1 y2 beta_q * (y2 - y1)) c1[i] np.clip(c1[i], lb[i], ub[i]) c2[i] np.clip(c2[i], lb[i], ub[i]) return c1, c2这段交叉代码里最关键的是beta和beta_q的换算。SBX 不是简单地把两个父代取平均而是按照分布指数 eta_c 生成一个与父代距离相关的概率分布子代可能离父代很远也可能很近。eta_c 越大子代越接近父代。最后那行np.clip是保险丝——所有算子变异后的个体都要重新裁剪到变量边界内这是保证后续评估不出 NaN 的基本功。def polynomial_mutation(ind, lb, ub, eta_m20): 对单个个体执行多项式变异。 n len(ind) mut ind.copy() for i in range(n): if np.random.random() 1.0 / n: continue delta1 (mut[i] - lb[i]) / (ub[i] - lb[i]) delta2 (ub[i] - mut[i]) / (ub[i] - lb[i]) rand np.random.random() if rand 0.5: mul (2.0 * rand) ** (1.0 / (eta_m 1.0)) - 1.0 else: mul 1.0 - (2.0 * (1.0 - rand)) ** (1.0 / (eta_m 1.0)) mut[i] mut[i] mul * (ub[i] - lb[i]) return np.clip(mut, lb, ub)变异概率 1/n 决定平均每个个体只有一个维度变异。很多翻车现场是把这个概率当成整个个体的变异概率结果每个体被改得面目全非算法永远不收敛。如果你想调整探索强度优先调 eta_m 而不是把变异概率拉到 0.5。4.4 主循环把上面所有碎片拼起来有了初始化、选择、交叉、变异、环境选择之后主循环其实很短def run_nsga2(problem, lb, ub, pop_size200, max_gen500, pc0.9): pop initialize_population(pop_size, lb, ub) fitness evaluate_population(pop, problem) for gen in range(max_gen): # 1. 二元锦标赛选父代 offspring [] while len(offspring) pop_size: idx np.random.choice(pop_size, 2, replaceFalse) # 简化这里用随机选两个父代直接交叉完整版应按 rank 和拥挤度做锦标赛。 p1, p2 pop[idx[0]], pop[idx[1]] if np.random.random() pc: c1, c2 sbx_crossover(p1, p2, lb, ub) else: c1, c2 p1.copy(), p2.copy() offspring.append(polynomial_mutation(c1, lb, ub)) if len(offspring) pop_size: offspring.append(polynomial_mutation(c2, lb, ub)) offspring np.array(offspring) # 2. 合并父代和子代 combined_pop np.vstack([pop, offspring]) combined_fitness evaluate_population(combined_pop, problem) # 3. 环境选择 keep_idx select_next_generation(combined_fitness) pop combined_pop[keep_idx] fitness combined_fitness[keep_idx] return pop, fitness这段主循环里刻意做了简化父代选择直接随机抽没有严格按 rank 拥挤度做锦标赛。原因是主循环已经较长真正完整版本应当在 while 循环里先计算当前父代每个个体的 rank 和拥挤度再执行二元锦标赛。你可以在自己项目里把np.random.choice那段换成基于 rank 和拥挤度的选择函数。这个流程跑完以后fitness 里 rank 为 1 的个体就是最终找到的帕累托前沿。把它们按目标值保存下来是报告里所有配图的数据来源。4.5 输出结果保存前沿数据和前沿图课程设计的报告需要图和表这一步别省import matplotlib.pyplot as plt # 提取第 1 层个体 rank1_mask np.zeros(len(fitness), dtypebool) # 简便做法自己写了个简化标记正式项目中用 select_next_generation 的返回值 # 这里直接按“非支配”筛选 final_front [] for i in range(len(fitness)): if not any(dominates(fitness[j], fitness[i]) for j in range(len(fitness)) if j ! i): final_front.append(fitness[i]) final_front np.array(final_front) np.savetxt(pareto_front.txt, final_front, headerf1 f2, fmt%.6f) plt.scatter(final_front[:, 0], final_front[:, 1], s12) plt.xlabel(f1) plt.ylabel(f2) plt.title(NSGA-II on MMF-like problem) plt.savefig(pareto_front.png, dpi200)这里附加一段筛选前沿的循环是因为 select_next_generation 返回的是索引最后一代按非支配关系重新筛一遍更直观。保存 txt 的意义是绘图、IGD 计算、HV 计算都要基于同一份数据避免以后要从图上反向读取数据点再算指标。5. 课设避坑清单现象、原因、解决5.1 现象跑了 500 代帕累托前沿只有一条细线所有点挤在一起原因拥挤度距离的边界个体没有赋 np.inf。没有这个 inf前沿两端在每次截断时被最先淘汰种群逐渐收缩到中间区域多样性丢失。解决在 crowding_distance 函数里对每一维目标排序后的第一个和最后一个个体直接赋予 inf并且累加过程不要把这些边界个体计数进去。改完之后重新观察前沿图两端应该出现明显散开的点。5.2 现象支配判断结果完全混乱rank 分层和手算对不上原因把最小化目标和最大化目标混在同一个 fitness 数组里或者把np.all(a b) and np.any(a b)写成了np.all(a b)。后者会漏掉大量本应被支配的个体导致第 1 层异常庞大。解决统一在 problem 函数里把所有目标转成最小化测试阶段构造两个已知目标向量手算一遍 dominates 结果打印出来核对。等号的处理尤其要注意a b里的等号是支配判断成立的必要条件不能省。5.3 现象子代个体出现 NaN目标函数报错原因交叉或变异后新个体超出变量边界带到目标函数里计算出非法值。SBX 生成的子代理论上可能越界多项式变异同理。解决每个算子的最后一步都调用np.clip(x, lb, ub)。这是投入最小、收益最大的保险丝。另外在 evaluate_population 里加一个np.isnan(fitness).any()的检查一旦发现直接打印出问题的个体方便回溯是哪个算子造成的。5.4 现象结果每次都变两张图完全不一样报告没法写原因全程使用全局随机状态没有固定随机种子。np.random.seed 设置了但被打断或者用多线程并行实验时随机状态互相污染。解决在程序入口固定np.random.seed(20240401)并且只用一个rng np.random.default_rng(seed)实例贯穿所有随机操作每次独立运行前重新设置种子。课程设计要跑 20 次取平均我一般是为每次运行配置不同的种子号比如 1 到 20并把种子号写进文件名比如pareto_front_seed05.txt这样数据可追溯。5.5 现象算法看似收敛但换一个问题集就彻底失效原因只在一个问题上调通了参数参数过度拟合。例如把变异概率调大让多模态问题刚好能跳出局部段换到另一个更平滑的问题上反而震荡不收敛。解决实验设计里至少选 4 个问题一个平滑、一个多模态、一个带约束统一用同一组默认参数跑一遍。报告里写明“本实验统一采用种群 200、代数 500、pc0.9、pm1/n”而不是针对每个问题单独调一套参数这样结论才有泛化性。6. 进阶验证用 IGD/HV 把“表现好”写进报告课程设计只放几张前沿图评审老师很容易追问“好在哪里”。IGD 和 HV 是两个最常用的量化指标。IGD 衡量近似前沿与真实前沿之间的平均距离越小越好HV 衡量近似前沿在目标空间里覆盖参考点的体积越大越好。IGD 的计算需要真实帕累托前沿参考点。CEC-2021 官方问题大多公布了参考前端也可以自己用高精度网格采样来生成。代码里用 scipy 的距离矩阵即可from scipy.spatial import distance def igd(approx_front, true_front): IGD 真实前沿每个点到近似前沿最近距离的平均值。 approx_front: 算法找出的前沿点集 true_front: 真实前沿采样点集 d distance.cdist(true_front, approx_front) return float(np.mean(np.min(d, axis1)))注意矩阵的行列方向cdist(true_front, approx_front) 的每一行对应一个真实前沿点对每一行取最小值就是该真实前沿点到近似集合的距离。算出来数值越小说明近似前沿越接近完整前沿。这里的常见错误是把矩阵写反得到的是另一种指标数值含义完全不同。HV 的计算可以用蒙特卡洛近似在目标空间里以参考点为顶角生成大量随机点统计落在近似前沿之下的比例再乘以包围盒体积。参考点一般选取各目标最大值稍大的点。HV 的优点是不需要知道真实前沿适合竞赛题没有官方参考解的场合。缺点是它对目标值尺度敏感同一组数据换参考点结果差异很大报告里必须写清楚参考点取了多少。实验设计方面我每次跑 20 到 30 次独立实验记录每一代的 IGD 均值和中位数画收敛曲线时用均值加阴影表示方差区间。只把最后一次迭代的前沿导出来画帕累托图而把每一代的指标数据单独存 CSV。这样做的好处是最后写报告时随时能反查“第 100 代到底收敛到了什么程度”。如果你想在课设里超出平均水平可以多做一步把 NSGA-Ⅱ 跑出来的结果和随机搜索、单目标遗传算法各跑一遍做对比。哪怕对比结果不如 NSGA-Ⅱ也说明了你理解“多目标问题为什么需要多目标算法”。这个对比不需要额外写代码改一下主循环里的环境选择部分就能实现。做完这一步课程设计的深度基本上就够看了。我的个人习惯是代码全部跑通之后新建一个README.md写明运行环境、依赖版本、一键复现的命令连 Python 版本和 NumPy 版本都标注清楚。课设时间越紧越需要这份“后悔药”——两周后重新打开代码你绝不会记得当初用什么环境跑的。希望这份拆解能帮你少踩几个坑。本文还有配套的精品资源点击获取