企业官网建设流程全解析

简介：一份面向数学建模竞赛与课程学习的幻灯片课件，围绕动物群体的常微分方程模型，系统讲解如何用微分方程刻画种群动态、分析平衡点稳定性，并以ACM-85试题A为实例，推导有限资源环境下最优捕捞策略与最大净利润条件。资源包仅含1个演示文稿文件，大小约1.05MB，内容高度凝练，适合备赛学生、建模爱好者及相关课程师生快速查阅。目前已有113人学习，属于小而精的建模方法资料。课件从单种群开发模型讲起，覆盖收获率与最大可承受产量的关系、稳定与不稳定平衡点的判别，并引入弱肉强食的Volterra模型，展示狐兔数量交替波动的生态平衡过程；同时结合价格、捕捞成本和增长率，给出净利润最大化时捕捞方案的设计思路。通过这份幻灯片，读者可以快速掌握将实际问题抽象为常微分方程、做定性分析并指导决策的完整链条，为处理同类种群管理或资源优化题目提供可直接借鉴的建模步骤。

1. 为什么数学建模题里的动物群体常微分方程模型值得单独练

常微分方程建模在数学建模竞赛里的出现频率极高，而动物群体是最容易入门的载体：状态变量直接是种群数量，方程结构直观，结果又能落到“保护生态还是控制虫害”这类实际决策上。建模课件里讲动物群体的常微分方程模型时，公式一般列得很清楚，真正卡住大部分人的是后续那几步——参数改了不分析稳定性，初值换了不观察数值解是否发散，最后画出的图对不上题目要求。这篇文章按“先推方程，再用 Python 求解器跑通，最后用一个完整案例把参数影响、稳定性分析和调试方法串起来”的顺序展开。正在备赛的学生，以及刚接触动态系统建模的开发者，都可以顺着这条路径把一个动物群体建模题跑出有依据的相图与平衡点结论。

2. 从单种群增长到捕食关联：动物群体ODE模型的建立过程

2.1 Malthus 模型为什么只能当初始版本

对动物群体建模，第一步是确定状态变量与时间尺度。对单一种群，设第 t 时刻的种群数量为 N(t)，把这个量对时间求导就得到方程。最原始的 Malthus 模型写作

dN/dt = r * N

其中 r 称为内禀增长率，量纲是“单位时间内单个个体对种群增长的贡献”。这个一阶线性方程的解是N(t) = N(0) * exp(r*t)，即指数增长。它适合描述资源充足的阶段性增长，但真实生态系统中，种群会受食物、生存空间等条件限制。引入环境承载力 K，得到 Logistic 方程：

dN/dt = r * N * (1 - N/K)

非线性项(1 - N/K)体现了密度制约。当 N 远小于 K 时，方程接近指数增长；当 N 接近 K 时增长率趋零；一旦 N 超过 K，增长率变为负值，种群回落。参数 r 决定趋近 K 的快慢，K 决定长期平衡值，这两个参数是后续参数估计的对象。这里不要匆匆带过。不少竞赛题的区分点就藏在追问里：K 为什么是常数？如果环境波动明显，K 需要看成时间函数，模型便从自治系统变为非自治系统，后续相平面分析思路要跟着换。

2.2 Lotka-Volterra 方程：捕食、竞争、互惠共用一套符号规则

两物种模型最经典的是 Lotka-Volterra 捕食-被捕食模型。设 x 为被捕食者数量，y 为捕食者数量，方程为：

dx/dt = α*x - β*x*y dy/dt = δ*x*y - γ*y

参数含义：α 是被捕食者在没有捕食者时的净增长率；β 是捕食行为造成的猎物损失系数；δ 描述猎物转化为捕食者生物量的效率；γ 是捕食者在没有猎物时的死亡率。注意-β*x*y与+δ*x*y来自同一次捕食事件，但两个系数通常不等，因为能量在营养级之间传递时有损耗。

“状态变量 + 相互作用项”的写法可以覆盖其他物种关系，整理成表：

写竞争或互惠模型时，最容易被忽略的是环境承载力如何分配。两个物种共享同一种食物和生存空间，那么 x 的承载力项应为1 - (x + c1*y)/K1，其中 c1 表示单位 y 个体对 x 资源的占用比例。c1=1 表示完全重叠，c1=0 表示完全分离。先把这个写清楚再做参数标定，可以避免后期模拟出现“两物种之间没有关系、曲线却一抬一落”的伪相关。

2.3 建模时先过的两关：量纲一致与非负约束

写方程时我先做两个检查，避免模型偏离实际。第一，量纲一致性：方程右端每一项都必须是“数量/时间”。β*x*y看起来是数量的平方，但 β 本身带着1/(数量·时间)的量纲。写代码时不声明单位，但参数标定必须保证同一套时间单位，混用“天”和“月”是数值爆炸的最常见来源。第二，状态变量非负：真实种群数量不会为负，ODE 解算器并不天然知道这个约束。参数或初值组合不合适时，求解结果可能冒出负值，此时要改参数区间，而不是简单地把负值截断为零。截断动作会破坏解算器内部的连续性判断，后面的轨迹会出现解释不通的拐角。

3. 用SciPy把动物群体ODE跑通：最小代码与求解器参数

3.1 模型函数与求解器分离：一套代码可复用

SciPy 生态里，积分接口推荐使用integrate.solve_ivp。早期常用的odeint在新的 SciPy 版本里已经处于遗留状态，不放进新代码。我习惯把模型参数全部通过args传入求解器，而不是把参数常量写在函数体内部，这样后面的参数扫描可以直接复用同一个模型函数：

import numpy as np from scipy.integrate import solve_ivp import matplotlib.pyplot as plt def lotka_volterra(t, z, alpha, beta, delta, gamma): x, y = z # z[0] 为猎物，z[1] 为捕食者 dx = alpha * x - beta * x * y dy = delta * x * y - gamma * y return [dx, dy] params = (1.1, 0.4, 0.1, 0.4) z0 = [10, 3] t_span = (0, 60) t_eval = np.linspace(0, 60, 1200) sol = solve_ivp( lotka_volterra, t_span, z0, args=params, method="RK45", rtol=1e-6, atol=1e-9, t_eval=t_eval, ) print(sol.success, sol.message)

模型函数的参数顺序是时间 t、状态数组 z 以及四个模型参数。函数返回的dx、dy顺序必须与z[0]、z[1]对应，这里先写猎物再写捕食者。solve_ivp的输入里，t_span是积分区间，t_eval只决定输出点的采样密度，不影响内部自适应步长，因此可以适当加密采样点让绘图更平滑。args的传入顺序与模型函数形参顺序一致，后续用functools.partial固定一部分参数也很方便。

3.2 返回值里先看这三个字段

sol.success为 False 时说明积分提前终止，sol.message会给出终止原因。sol.t是时间序列，sol.y的形状为(状态数, 采样点数)。对两物种模型，sol.y[0]是猎物曲线，sol.y[1]是捕食者曲线。一个比较常见的错误是直接拿sol.y整体画图，那样会得到两条形状相近的曲线，容易看混。相平面图则要单独取sol.y[0]作为横坐标、sol.y[1]作为纵坐标，这样轨道闭合性才看得清楚。

3.3 求解器怎么选：刚性问题先于“默认即可”

求解器的选择不复杂，按问题特性对号入座即可：

判断是否刚性，最省事的做法是保持同一组参数，把 method 换成 LSODA 或 Radau，比较计算量。如果 RK45 内部推进需要几千步，而 LSODA 只用几百步，说明原始问题是刚性的，应当坚持用隐式方法。动物群体模型中的刚性通常来自增长率的时间尺度差异，比如猎物以“天”为单位增长，捕食者存活期以“月”为单位，两边特征时间相差很大。

容差参数的影响也值得单独说。rtol是相对误差容限，atol是绝对误差容限。动物群体的状态变量在 0 附近徘徊时，atol过大会容忍负值，结果出现类似物种灭绝又莫名回到正值的情况。经验是把atol设在初始状态量级小 3 到 5 个数量级的位置，初值在几十的量级时，atol取 1e-9 到 1e-12 都比较合适。

提示：参数都写死在模型函数里不算错，但会拖慢后续参数扫描。把参数扫描逻辑放在solve_ivp调用层，用循环遍历参数网格，灵敏度分析的代码结构会清晰很多。

4. 兔与狐狸：一个动物群体ODE建模案例的参数影响与输出检验

4.1 参数来源与量纲换算

竞赛题里，参数通常需要从题目给的生态背景数据估算。以兔-狐系统为例：假设兔子在没有狐狸时每周净增长率为 0.8，即周率 α=0.8；每只狐狸每周对兔子群体造成的捕食压力系数 β=0.4；狐狸从每单位兔子生物量中获得增长的效率 δ=0.1；没有兔子时狐狸每周死亡率 γ=0.4。代码里统一以“天”作为时间单位，这几个参数都除以 7，换算成日率后再传给solve_ivp。

换算错误很隐蔽，因为模型本身不会报警，只会让输出曲线的时间尺度整体漂移。确认量纲是否统一的方法是把周率和日率两组参数各跑一遍，比较平衡点位置：平衡点应当落在同一位置，只是到达平衡的快慢不同。若平衡点变了，说明参数单位混用，需要回到换算环节排查。

4.2 基准模拟与相平面图

继续使用第 3 章的模型函数，完成模拟后绘制时间序列与相平面：

fig, axes = plt.subplots(1, 2, figsize=(10, 3.5)) axes[0].plot(sol.t, sol.y[0], label="rabbit") axes[0].plot(sol.t, sol.y[1], label="fox") axes[0].set_xlabel("time (days)") axes[0].set_ylabel("population") axes[0].legend() axes[1].plot(sol.y[0], sol.y[1]) axes[1].axhline(alpha / beta, color="gray", linestyle="--") axes[1].axvline(gamma / delta, color="gray", linestyle="--") axes[1].set_xlabel("rabbit") axes[1].set_ylabel("fox") plt.tight_layout()

相平面里的竖直虚线是猎物零增长线x* = gamma/delta = 4，水平虚线是捕食者零增长线y* = alpha/beta = 2.75，交点就是平衡点。注意两条线的方向：猎物没有净变化时，得到的是捕食者数量的固定值，对应水平方向的 y 值；捕食者没有净变化时，得到的是猎物数量的固定值，对应竖直方向的 x 值。把这两条线画反之后，后续所有稳定性讨论都会偏离。

4.3 单参数扰动：周期比稳态更值得记录

保持其他参数不变，只把 β 从 0.4 提高到 0.6，捕食效率上升后，时间序列的振荡周期会变短，振幅也会增大。再把 γ 从 0.4 降到 0.2，捕食者死亡率降低，平衡点向右移动。单参数扫描需要记录两个量：平衡点位置与振荡主周期。主周期可用scipy.signal.find_peaks对sol.y[0]做峰值检测后取相邻峰间隔的平均值：

from scipy.signal import find_peaks peaks, _ = find_peaks(sol.y[0]) # 只对猎物时间序列检测峰值 periods = np.diff(sol.t[peaks]) print("mean period:", periods.mean(), "days")

如果数值噪声干扰了峰值检测，可以先对sol.y[0]做一次滑动平均再交给find_peaks。这段代码本身简单，但说明一个建模习惯：ODE 模型的可观察量不只包含稳态均值，振荡周期和相位关系同样有分析价值。竞赛报告里写出“系统周期随捕食效率上升而变短”这类结论，比只贴一张图更有信息含量。

4.4 初值改变与参数改变要分开讨论

初值对 ODE 演化的影响需要分模型区别对待。对单一物种 Logistic 模型，初值只影响趋近 K 的路径，不影响最终平衡值；对 Lotka-Volterra 捕食模型，不同初值落在相平面不同的闭合轨道上，参数决定平衡点的位置，初值决定环绕中心走哪条轨道。报告中每写一个动力学结论，都要先确认它是由参数驱动还是初值驱动。切换初值后如果系统形态发生质变，比如从闭合轨道变成发散的螺旋，说明模型里可能还有另一个不稳定平衡点，需要进入下一章的稳定性分析去确认。

5. 稳定性分析与数值排错：动物群体ODE模型发散排查

5.1 平衡点与雅可比矩阵：从零增长线到特征值判断

相平面里零增长线的交点就是平衡点，但平衡点是否稳定，要看平衡点处雅可比矩阵的特征值。对 Lotka-Volterra 模型做扰动展开，雅可比矩阵为

J = [[α - β*y, -β*x], [δ*y, δ*x - γ]]

在平衡点(x*, y*) = (γ/δ, α/β)处，矩阵变成

J = [[0, -β*γ/δ], [δ*α/β, 0]]

特征值是纯虚数，系统在平衡点附近做周期振荡，不收敛也不发散。用代码验证：

import numpy as np alpha, beta, delta, gamma = 1.1, 0.4, 0.1, 0.4 x_star = gamma / delta y_star = alpha / beta J = np.array([ [0, -beta * gamma / delta], [delta * alpha / beta, 0], ]) eigvals = np.linalg.eigvals(J) print("特征值:", eigvals)

输出是一对共轭纯虚根。纯虚根对应的运动形态是环绕平衡点的闭合轨道，而不是向内收敛的螺旋。这个差别在相平面图上肉眼可见：闭合轨道意味着系统对初值有记忆，同一组参数下不同初值对应不同振幅的振荡。现实中纯虚根很难长期成立，环境随机波动会把轨道推离原闭合曲线，因此实际建模常把猎物方程加上密度制约项，改成dx/dt = α*x*(1 - x/Kx) - β*x*y，运动特征随之从中性稳定变为阻尼振荡。

5.2 数值发散时的三步排查

“模拟曲线飞到天上去”的原因，多数不是模型写错，而是求解器设置与模型时间尺度不匹配。按顺序排查效率最高。

第一步，显式设置max_step。solve_ivp默认最大步长可能偏大，当参数数量级很小时，这一步会跨过高频振荡周期，轨迹直接跳到系统范围之外，画出来就是高速跳变线。设max_step=0.1后观察结果是否稳定。

第二步，切换求解器。同一组参数同时用 RK45 和 LSODA 跑一遍，比较sol.y的最大值。如果 RK45 的结果到了 1e6，LSODA 的结果在几百这个量级，基本可判定是刚性问题引起的数值不稳定，不是方程本身发散。

第三步，查时间单位。周率与日率的混用是低级却高发的错误。把整组参数都除以 7 后重跑，如果平衡点位置不变、只有趋近速度变化，说明换算正确；如果平衡点位置漂移，就从每个参数进入方程的系数开始查，不要只盯模型函数里的那几行。

5.3 用断言提前暴露问题，不靠肉眼扫图

参数扫描过程中，建议把单次求解封装起来，在内层加断言：

def solve_lv(params, z0, t_span): sol = solve_ivp(lotka_volterra, t_span, z0, args=params, method="RK45") assert sol.success, f"solver failed at params {params}: {sol.message}" return sol

只要某一组参数让积分器失败，异常信息会直接指出是哪一组参数引起的。sol.status字段同样有用：0 表示正常，1 表示事件触发提前停止，-1 表示积分失败。没有自定义事件时 status 几乎总是 0，一旦出现非 0 状态，优先读sol.message的文本，而不是反复调整画图范围去猜测原因。

6. 用零增长线图快速完成动物群体ODE的参数区间设计

6.1 零增长线是建模学习阶段的调试器

零增长线图在第 4 章已经画过，它的价值在于比时间序列更容易看出参数变化的方向。把猎物方程置零得到一条关于 y 的直线，把捕食者方程置零得到一条关于 x 的直线。参数改变时这两条线在相平面里平移，交点的移动就是平衡位置的变化。调参之前先用几何工具把可行域框出来，能省下大量无效的积分运算。

6.2 用代码批量过滤无意义参数组合

参数区间设计不需要只靠手工点图。写一组循环遍历 α、β、δ、γ 的候选区间，把平衡点不落在第一象限的组合过滤掉：

gamma, delta = 0.4, 0.1 alpha_range = np.linspace(0.5, 2.0, 50) beta_range = np.linspace(0.2, 0.8, 50) usable = [] for alpha in alpha_range: for beta in beta_range: x_star = gamma / delta y_star = alpha / beta if x_star > 0 and y_star > 0: usable.append((alpha, beta))

这里y_star = alpha / beta的正性要求其实由参数本身为正自动满足，但这套写法可以扩展出更严格的条件，比如要求x_star与y_star都落在某个生态观测区间内。把零增长线的代数条件写进筛选逻辑后，整个参数网格扫描就变成了预报步骤，而不是事后整理结果。

6.3 三步闭合成一套固定流程

常见的做法是把流程固定成三步：先用零增长线筛选平衡点位置，给出参数可行区间；再用solve_ivp对区间内有代表性的点计算时间序列；最后回到雅可比矩阵判断稳定类型。每一步都能验证前一步的结果，写报告时只需保留最后一步的图。绘图调试有一个小技巧：检查零增长线位置时用灰白底的单色图，不要加颜色填充，颜色会干扰平衡点附近微小偏移的判断。

本文还有配套的精品资源，点击获取