物理 备战数学建模赛打卡帖
这个帖子先用来规划写论文的进度
那个建模能力赛就是北师大组织的那个
顺便问一下,坛上还有没有谁参加了的啊
说一个扎心的事实,留给我的时间比别人少了将近一半
只是因为我们的县城中学没有通知我们,导致最近才知道……
import numpy as np
from scipy.integrate import solve_ivp
import matplotlib.pyplot as plt
# -------------------------------
# 定义 SEIR 模型
def seir_model(t, y, beta, sigma, gamma):
S, E, I, R = y
dSdt = -beta * S * I / N
dEdt = beta * S * I / N - sigma * E
dIdt = sigma * E - gamma * I
dRdt = gamma * I
return [dSdt, dEdt, dIdt, dRdt]
# -------------------------------
# 适应度函数:最小化 MSE
def fitness(params, data):
beta, sigma, gamma = params
sol = solve_ivp(seir_model, [0, T], [S0, E0, I0, R0], args=(beta, sigma, gamma), t_eval=t_points)
I_pred = sol.y[2]
mse = np.mean((I_pred - data)**2)
return 1 / (mse + 1e-8) # 避免除零
-----------------
# 遗传算法主循环
def genetic_algorithm():
population_size = 50
generations = 100
mutation_rate = 0.1
crossover_rate = 0.7
# 初始化种群
pop = np.random.uniform(0.01, 0.5, (population_size, 3)) # beta, sigma, gamma
for gen in range(generations):
fitnesses = [fitness(ind, real_data) for ind in pop]
# 选择
parents = []
for _ in range(population_size // 2):
idx1, idx2 = np.random.choice(len(pop), 2, replace=False)
if fitnesses[idx1] > fitnesses[idx2]:
parents.append(pop[idx1])
else:
parents.append(pop[idx2])
# 交叉
new_pop = []
for i in range(0, len(parents), 2):
p1, p2 = parents[i], parents[i+1]
if np.random.random() < crossover_rate:
point = np.random.randint(1, 3)
child1 = np.concatenate([p1[:point], p2[point:]])
child2 = np.concatenate([p2[:point], p1[point:]])
new_pop.extend([child1, child2])
else:
new_pop.extend([p1, p2])
# 变异
for i in range(len(new_pop)):
if np.random.random() < mutation_rate:
j = np.random.randint(3)
new_pop[i][j] += np.random.normal(0, 0.05)
new_pop[i][j] = np.clip(new_pop[i][j], 0.01, 0.5)
pop = np.array(new_pop)
best_idx = np.argmax(fitnesses)
best_params = pop[best_idx]
print("最佳参数:", best_params)
return best_params