ADD file via upload

NSGA2.py

@@ -0,0 +1,145 @@
+import numpy as np
+import matplotlib.pyplot as plt
+from mpl_toolkits.mplot3d import Axes3D
+# 定义目标函数
+def objective(x):
+    return np.array([x[0]**2 + x[1]**2, (x[0]-2)**2 + x[1]**2])
+
+# 计算拥挤度
+def crowding_distance(F):
+    n, m = F.shape
+    d = np.zeros(n)
+    for i in range(m):
+        idx = np.argsort(F[:, i])
+        d[idx[0]] = np.inf
+        d[idx[-1]] = np.inf
+        for j in range(1, n-1):
+            d[idx[j]] += (F[idx[j+1], i] - F[idx[j-1], i]) / (F[idx[-1], i] - F[idx[0], i])
+    return d
+
+# 非支配排序
+def non_dominated_sorting(F):
+    n = F.shape[0]
+    S = [[] for i in range(n)]
+    fronts = []
+    N = np.zeros(n)
+    F1 = []
+    for p in range(n):
+        S[p] = []
+        N[p] = 0
+        for q in range(n):
+            if np.all(F[p] <= F[q]) and np.any(F[p] < F[q]):
+                if q not in S[p]:
+                    S[p].append(q)
+            elif np.all(F[q] <= F[p]) and np.any(F[q] < F[p]):
+                N[p] += 1
+        if N[p] == 0:
+            F1.append(p)
+    fronts.append(F1)
+    i = 0
+    while len(fronts[i]) > 0:
+        Q = []
+        for p in fronts[i]:
+            for q in S[p]:
+                N[q] -= 1
+                if N[q] == 0:
+                    Q.append(q)
+        i += 1
+        fronts.append(Q)
+    return fronts[:-1]
+
+# 锦标赛选择
+def tournament_selection(F, I, K):
+    idx = np.random.choice(I, K, replace=False)
+    rank = np.argsort(np.argsort(F[idx, 0])) + np.argsort(np.argsort(F[idx, 1]))
+    return idx[np.argmin(rank)]
+
+# 选择操作
+def selection(F, I, K):
+    Xp = np.zeros((K, len(xl)))
+    for i in range(K):
+        p1 = tournament_selection(F, I, 2)
+        p2 = tournament_selection(F, I, 2)
+        if np.random.random() < 0.5:
+            Xp[i] = X[p1]
+        else:
+            Xp[i] = X[p2]
+    return Xp
+
+# 交叉操作
+def crossover(x1, x2, pc):
+    if np.random.random() < pc:
+        alpha = np.random.uniform(-0.5, 1.5, len(x1))
+        c1 = alpha*x1 + (1-alpha)*x2
+        c2 = alpha*x2 + (1-alpha)*x1
+        return c1, c2
+    else:
+        return x1, x2
+
+# 变异操作
+def mutation(x, pm, xl, xu):
+    for i in range(len(x)):
+        if np.random.random() < pm:
+            x[i] = np.random.uniform(xl[i], xu[i])
+    return x
+
+# NSGA-II算法
+def nsga2(objective, xl, xu, n, T, pc, pm):
+    # 初始化种
+    X = np.random.uniform(xl, xu, (n, len(xl)))
+    # 迭代 T 次
+    for t in range(T):
+        # 计算目标函数值和帕累托前沿
+        F = np.zeros((n, 2))
+        for i in range(n):
+            F[i] = objective(X[i])
+        fronts = non_dominated_sorting(F)
+
+        # 计算拥挤度
+        CD = np.zeros(n)
+        for i, front in enumerate(fronts):
+            FD = F[front]
+            CD[front] = crowding_distance(FD)
+
+        # 选择操作
+        I = np.concatenate(fronts)
+        K = len(I)
+        Xp = selection(F, I, K)
+
+        # 交叉操作
+        for i in range(0, K, 2):
+            Xp[i], Xp[i+1] = crossover(Xp[i], Xp[i+1], pc)
+
+        # 变异操作
+        for i in range(K):
+            Xp[i] = mutation(Xp[i], pm, xl, xu)
+
+        # 生成新种群
+        X = Xp
+
+    # 计算最终帕累托前沿和解集
+    F = np.zeros((n, 2))
+    for i in range(n):
+        F[i] = objective(X[i])
+    fronts = non_dominated_sorting(F)
+    # 生成网格点
+    x = np.linspace(-5, 5, 100)
+    y = np.linspace(-5, 5, 100)
+    X, Y = np.meshgrid(x, y)
+    Z1 = (X - 2) ** 2 + (Y + 1) ** 2
+    Z2 = 4 * (X + 3) ** 2 + (Y - 1) ** 2
+    # 绘制等高线图
+    fig = plt.figure()
+    ax = fig.add_subplot(111, projection='3d')
+    ax.plot_surface(X, Y, Z1, cmap="Blues", alpha=0.5)
+    ax.plot_surface(X, Y, Z2, cmap="Oranges", alpha=0.5)
+    # 绘制 Pareto 前沿
+    for front in fronts:
+        FD = F[front]
+        ax.plot(FD[:, 0], FD[:, 1], "ro-", label="Pareto Front")
+
+    ax.set_xlabel("x1")
+    ax.set_ylabel("x2")
+    ax.set_zlabel("Objective Values")
+    plt.legend()
+    plt.show()