简介面向科学计算与机器学习交叉领域研究者的 PINN 物理信息网络源码包围绕离散时间识别、离散时间推理、连续时间识别、连续时间推理四类典型任务提供可运行的 Python 实现与配套数据。压缩包共 559 个文件其中 txt 文本文件多达 503 个主要保存训练日志、配置说明与结果记录13 个 py 脚本为核心代码完整封装模型训练、推理和误差评估流程12 个 csv 误差表记录不同参数下的 error_table 与 lambda 误差便于横向对比精度另有 eps 与 pdf 格式图表、mat 数据文件等辅助分析材料整体约 475.95MB直接解压即可查阅。已有 1546 人学习下载。通过源码可掌握 PINN 在时间维度上识别与推理的建模思路利用 csv 表中的误差指标和可视化图表能快速定位调参方向适合用于课程设计、论文复现或在此基础上扩展新的物理约束模型。1. PINN 四种任务模式的定位与选型先明确一点这套资源不是单个脚本而是把物理信息神经网络PINN按时间域划分成离散时间识别、离散时间推理、连续时间识别、连续时间推理四条技术路线每条路线配套独立可跑的 Python 源码和 CSV 误差记录表。很多人拿到代码后第一反应是“四个文件差不多”实际上这四者解决的是两类完全不同的问题识别identification是已知系统部分结构、反推未知参数推理inference是给定观测数据、重构整个动态系统的状态轨迹。时间域离散与连续的区别则决定了网络如何构造时间导数项——离散时间用 Runge-Kutta 等积分步进结构连续时间则依赖自动微分直接计算微分算子。适用于三类读者正在复现 PINN 基线实验的研究生、需要把物理约束嵌入神经网络的算法工程师、以及做数字孪生或动力学系统建模的从业者。如果只是为了跑通 PINN 做数据拟合这四个代码可能显得“重”但一旦涉及方程形式未知、参数反演或噪声观测下的状态估计离散与连续、识别与推理的区分就是绕不开的起点。下面从时间域的数学设定差异切入逐步拆解每一类代码的模型结构、损失函数和调试要点最后说明如何用附带的 CSV 误差表做收敛性验证。2. 连续时间推理与识别自动微分驱动的物理约束2.1 连续时间推理的方程约束构造连续时间场景假设系统的状态由常微分方程组ODE控制已知方程形式但未知部分参数或初值。PINN 推理任务的核心是构造一个神经网络 u_θ(t)使其输出同时满足两个条件一是逼近观测数据二是让方程残差尽可能小。这里的方程残差定义为R(t) du_θ/dt - f(u_θ, t, λ)其中 λ 是待学习的物理参数。在 PyTorch 中实现时时间导数 du_θ/dt 直接通过torch.autograd.grad对输入时间求一阶导得到这正是 PINN 相比传统数值方法的根本差异——不需要离散化微分算子而是利用自动微分链式法则把物理约束嵌入损失函数。import torch import torch.nn as nn class ContinuousInferencePINN(nn.Module): def __init__(self, n_input1, n_hidden50, n_output2, n_layers4): super().__init__() layers [nn.Linear(n_input, n_hidden), nn.Tanh()] for _ in range(n_layers - 2): layers [nn.Linear(n_hidden, n_hidden), nn.Tanh()] layers [nn.Linear(n_hidden, n_output)] self.net nn.Sequential(*layers) def forward(self, t): return self.net(t) def residual(self, t, lambda_param): t t.clone().requires_grad_(True) u self.forward(t) du_dt torch.autograd.grad( u, t, grad_outputstorch.ones_like(u), create_graphTrue )[0] # 以 Lotka-Volterra 两变量系统为例 x, y u[:, 0:1], u[:, 1:2] fx lambda_param[0] * x - lambda_param[1] * x * y fy -lambda_param[2] * y lambda_param[3] * x * y return du_dt - torch.cat([fx, fy], dim1)这段代码的关键有三处。第一requires_grad_(True)必须在残差计算之前重新设置因为每次前向传播的时间张量可能来自一个不需要梯度的 DataLoader。第二grad_outputstorch.ones_like(u)表示对输出逐元素求梯度这一参数漏写是常见的报错点——autograd 要求显式提供与输出同形状的 grad_outputs否则无法计算雅可比向量积。第三残差返回的是 du_dt 与方程右侧的差值维度必须与网络输出严格一致如果系统是三维状态空间而输出维度误设为 2后续损失拼接时会抛出形状不匹配错误。2.1.1 损失函数的多目标加权推理任务的损失函数由三部分构成观测数据拟合损失、方程残差损失、初始条件损失。权重 λ_data、λ_res 和 λ_ic 的选择直接影响训练结果。一般经验值是 λ_res 取 1.0、λ_data 为 0.1 到 1.0、λ_ic 为 1.0。原因在于方程残差是物理约束需要优先满足观测数据若噪声较大过高的 λ_data 会把噪声完全拟合成动力学模式反而使物理推理失真。def loss_function(model, t_data, u_data, t_col, t_ic, u_ic, lambda_param): u_pred model(t_data) loss_data torch.mean((u_pred - u_data) ** 2) residual model.residual(t_col, lambda_param) loss_res torch.mean(residual ** 2) u_ic_pred model(t_ic) loss_ic torch.mean((u_ic_pred - u_ic) ** 2) return 0.5 * loss_data 1.0 * loss_res 1.0 * loss_ic这里的配置点是配点 t_col 的采样方式。常见做法是在时间区间内均匀采样 2000 个点如果系统存在高频振荡可以在预测值变化剧烈的区域加密采样。另一个实操细节是 λ_param 必须注册为nn.Parameter而非普通 Tensor否则优化器不会更新它。2.2 连续时间识别的参数反演策略识别与推理的区别在于推理只重建状态轨迹识别则要在重建轨迹的同时反推出方程中的未知参数 λ。实现上两者共用同一个网络差异体现在优化器的参数列表中——识别任务必须把 λ 传入优化器。这里有个隐蔽的边界如果初始 λ 与实际值偏离过大梯度会被残差项主导导致 λ 的更新方向受方程右侧非线性项支配而不稳定。缓解手段是分阶段训练前 500 轮固定 λ 只优化网络权重待残差降到一定阈值后再联合优化。optimizer torch.optim.Adam([ {params: model.parameters(), lr: 1e-3}, {params: [lambda_param], lr: 5e-3} ]) for epoch in range(total_epochs): if epoch 500: # 解锁 lambda 的优化前面固定可以避免初期发散 lambda_param.requires_grad_(True) optimizer.zero_grad() loss loss_function(model, t_data, u_data, t_col, t_ic, u_ic, lambda_param) loss.backward() optimizer.step()注意lambda_param.requires_grad_(True)的时机。严格说如果定义时没有设置 requires_gradAdam 会直接跳过这个参数所以不建议用 if 分支控制 requires_grad 开关更稳妥的做法是用两个优化器——前 500 轮只用optimizer_model之后换用包含两个参数组的联合优化器。另外识别的评估指标不只是最终 λ 的绝对误差还要看 λ 在训练末期的波动幅度——如果连续 200 轮的变化量小于 1e-4说明收敛如果一直振荡说明观测数据的信息量不足以唯一确定参数需要增加配点密度或引入正则化项。3. 离散时间推理Runge-Kutta 积分器与观测融合3.1 从连续残差到离散步进结构当观测数据只在若干时间步上获得且相邻步之间没有高频率采样连续时间模型会面临一个可观察性问题方程残差在稀疏区间内可能对应多条轨迹网络没有足够约束选出正确的一条。离散时间推理的做法是利用 Runge-Kutta 积分器把状态按固定步长在时间上推进用一步积分后的预测值与下一时刻观测值做比较这样每一次时间步的更新都是一个“预测-校正”闭环观测信息被逐步注入而非一次性约束。3.1.1 经典四阶 Runge-Kutta 的 PINN 化改造以四阶 Runge-KuttaRK4为例在 PINN 中构造步进过程如下给定当前时刻的状态 u_n通过 RK4 计算 u_{n1} 的预测值损失函数比较该预测值与观测点 u_{n1}^{obs} 的偏差。RK4 的积分步长 Δt 属于超参数需要结合观测数据的时间间隔来设定。以附着动力学模型中常见的步长选择为例当观测间隔为 0.1 时Δt 通常设置在 0.01 到 0.05 之间步长过大会引入积分误差过小则计算量成倍增加。def rk4_step(model, t_n, u_n, dt, lambda_param): k1 model.rhs(t_n, u_n, lambda_param) k2 model.rhs(t_n dt / 2, u_n dt * k1 / 2, lambda_param) k3 model.rhs(t_n dt / 2, u_n dt * k2 / 2, lambda_param) k4 model.rhs(t_n dt, u_n dt * k3, lambda_param) return u_n dt / 6 * (k1 2 * k2 2 * k3 k4)这里的model.rhs是用神经网络表示的动力学函数输入的是当前时刻的状态和时间输出状态的导数。与连续时间模型中直接对时间自动微分不同离散时间法的核心在网络内部不再计算时间导数而是把网络整体当作一个 ODE 右侧函数时间的推进由 RK4 公式完成。这种设计的优势是训练时可以批量处理多个时间步——每个时刻的输出都经过整个 RK4 链路梯度会通过 K 系数回传到网络的每一层相当于在时间轴上做了反向传播。3.1.2 稀疏观测下的损失权重分配离散时间推理的损失函数按时间步加权越靠后的时间步预测误差会被 RK4 积分过程累积放大因此权重应随步数增加而递减。具体操作中步进索引 i 对应的损失权重取 1/(i1) 或 0.9^i。前几个时间步的观测对模型参数更新的贡献更大因为误差尚未经过长时间积分污染。def discrete_inference_loss(model, t_seq, u_seq, dt, lambda_param, gamma0.9): total_loss 0.0 u_current u_seq[0] for i in range(1, len(t_seq)): u_pred_next rk4_step(model, t_seq[i-1], u_current, dt, lambda_param) obs_next u_seq[i] weight gamma ** i total_loss weight * torch.mean((u_pred_next - obs_next) ** 2) # 这里要注意训练时下一时刻的输入应使用预测值还是观测值 u_current u_pred_next.detach() # 阻断梯度流动避免误差累积爆炸 return total_loss代码中的u_current u_pred_next.detach()是决定训练稳定性的关键。如果直接把预测值传入下一步运算梯度会沿着时间维度回传很多步导致梯度消失或爆炸断开后每一步的 RK4 结构相对独立虽然牺牲了长期依赖但稳定性显著提升。这是与序列模型用 Teacher Forcing 和 Free Running 进行取舍时的常见做法在可观察性足够时先用观测值做引导后续再切换为全预测模式微调。3.2 初始化与积分步长的联动关系离散时间模型对初始化比较敏感。如果网络的随机初始化使 RK4 的 K 系数计算值过大u 会在最初几步就偏离合理范围丢失有效信息。常见初始化方案是将网络最后一层权重置为较小的均值为 0、标准差 0.01且偏置初始化为满足先验的平均变化率。另一个做法是把 RK4 内部的 dt 用一个可学习的缩放因子修正——网络初始预测的导数尺度与实际系统偏差较大时可通过这个因子让模型先学到正确的变化率再逐步调整。成功经验是步长先取观测间隔的 1/10 保证积分稳定性训练 300 轮后逐步把 dt 提到 1/5 以覆盖更长的动力学特征。4. 离散时间识别离散观测下的参数反演实现离散时间识别是四条路线中最“工程化”的一条逻辑上需要在 RK4 步进结构下加入参数 λ 的梯度更新实现上要注意 λ 的梯度经过多条 RK4 路径后出现的量级失衡问题——如果 λ 初始值偏大方程右侧的乘积项会让导数预测值成倍放大K 系数随之膨胀。误差表error_lambda_1_table_2.csv中记录的正是不同 λ 初始值下识别误差的对比数据下面先用代码说明实现再基于误差表解释参数-初始值的敏感性。def discrete_identification_loss(model, t_seq, u_seq, dt, lambda_param): u u_seq[0] loss 0.0 all_pred [u] for i in range(len(t_seq) - 1): u_pred rk4_step(model, t_seq[i], u, dt, lambda_param) loss torch.mean((u_pred - u_seq[i 1]) ** 2) u u_pred # 这里不 detach让 lambda 能获得跨步梯度 all_pred.append(u_pred) return loss, all_pred与推理代码的显著区别是在整个时间序列上保留计算图不执行 detach实现“检测参数在长期积分效果上的影响”。但这会造成训练负担成倍增加——时间步数为 N 时反向传播的复杂度为 O(N × 网络层数)。对于一个 200 步的时间序列GPU 显存占用大约是推理模式的 20 倍以上这是限制离散时间识别可扩展性的主要瓶颈。缓解方案是分段优化——把序列切成数段每段保留计算图段间用 detach 切割相当于把长期依赖折中在数个子时段内。交替使用两种策略前 100 轮用分段模式快速逼近大致参数范围后续切换为全程模式精调。def segment_identification_loss(model, t_seq, u_seq, dt, lambda_param, seg_len20): total_loss 0.0 n_segments (len(t_seq) - 1) // seg_len # 每个段的起点直接使用观测值保证段内轨迹误差是真实偏差 for seg in range(n_segments): start seg * seg_len u u_seq[start] for i in range(start, start seg_len): u_pred rk4_step(model, t_seq[i], u, dt, lambda_param) total_loss torch.mean((u_pred - u_seq[i 1]) ** 2) u u_pred # 段末不将梯度传到下一段防止段间误差累积 u u.detach() return total_loss分段模式中每个段长度为 20 步显存占用约为全程模式的 1/10。观察配套误差表error_lambda_2_table_1.csv当 λ2 的初始值偏离真值超过 50% 时分段模式的识别误差比全程模式高出约 30%说明分段虽然省显存但精调的精度下降当初始值偏离在 20% 以内时两者误差几乎没有差异这就是工程上“分阶段先粗后精”的依据。另一个值得注意的点是识别误差表里的误差不是每个时间步的均方误差而是 λ 估计值与真实值的相对误差。因此训练结束后需单独打印 λ 的最终收敛值警惕“总损失下降但参数没收敛”的假收敛情况——优化器很容易把网络权重调得足够扭曲以拟合错误的 λ 值这种现象在识别任务中被称为“参数补偿效应”。4.1 多参数联合识别的正则化与可辨识性当系统存在多个未知参数时离散时间识别面临可辨识性边界——不同参数组合可能在同一组观测数据下产生几乎相同的轨迹。此时单纯的最小化损失无法统一确定参数需要添加先验正则化项。工程实现中针对性处理方式为在损失函数中加入 λ 与实际值估计先验的偏差项 λ_prior使用 L2 正则回归目标为初始猜测 λ0def identification_loss_with_prior(model, t_seq, u_seq, dt, lambda_param, lambda_prior, prior_weight0.05): loss, _ discrete_identification_loss(model, t_seq, u_seq, dt, lambda_param) prior_loss prior_weight * torch.sum((lambda_param - lambda_prior) ** 2) return loss prior_loss需要注意的是prior_weight 不应设置过大。若超过 0.1λ 会被强制拉向先验值而网络为了降低数据拟合误差会倾向于“把参数承担的作用转嫁给其他可调参数”反而扭曲了物理含义。比较合理的做法是先不做先验训练一轮观察 λ 的收敛趋势若明显向某一方向偏移再根据趋势设置一个偏向性的先验这是连续与离散识别任务中共通的传统技巧。5. 从 CSV 误差表反推训练诊断收敛、稳定性与数值病态这组资源中最容易被忽略但价值最高的是error_table_1.csv、error_lambda_1_table_1.csv、error_lambda_2_table_1.csv及对应 table_2 版本共六个误差表记录不同实验配置下的误差演化或分布。这部分并不是训练的必备输入而是用来定位模型失败模式的“病历表”——通过读表的思路能快速判断一个 PINN 模型是在数据拟合阶段失败、在物理约束阶段发散还是在参数更新阶段陷入病态。5.1 读懂误差表的三个关键字段字段名含义诊断指向iteration当前训练轮次或时间步曲线太陡或振荡说明学习率过高data_loss观测数据拟合损失持续高位说明网络容量不足或观测有噪声residual_loss方程残差损失指数下降后停滞说明需要增加配点或调整权重lambda_error参数估计相对误差数值跳动说明病态单调不降说明不可辨识table_1 与 table_2 的结构差异体现了识别与推理两种任务下误差指标的关注点不同。连续时间推理的 error_table 通常记录每轮训练后的总损失和验证集残差离散时间识别的 lambda 误差表则专门记录每个候选 λ 在当前 epoch 的估计值变化。若 lambda 误差先下降后上升典型原因是优化器在残差损失收敛后继续优化数据损失时破坏了参数稳定性——常见做法是冻结 λ单独对网络权重进行若干轮微调再解冻参数重新联合优化。5.2 用误差表定位 PINN 典型失败场景场景一residual_loss 持续在 1e-1 量级不下降。原因指向配点 t_col 在相空间覆盖面不够——系统若有多个吸引子或振荡模态均匀采样可能漏掉关键区域。调整方向是改用自适应配点从当前残差较大的位置周围加密采样迭代更新配点集合。场景二残余损失下降但 lambda_error 不降反升。这是参数补偿效应的直接体现另一条排查思路是路径检查训练过程中方程右侧的每一项对损失的梯度贡献对 lotka-volterra 系统如果 λ1 对应的线性项梯度比非线性交互项的小很多损失函数会将参数误差推给高梯度项此时需要对 λ 进行归一化或给低梯度项更大的权重。场景三离散时间推理的验证损失在某个 epoch 后突然陡增。这通常是 RK4 步长 dt 过大导致积分路径发散——误差表上会表现为该 epoch 之后的数据记录中断或残余损失出现数量级跳变。# 一个实际可用的自适应配点方法 def adaptive_resample(model, t_old, t_range, n_sample2000): # 先评估当前残差 with torch.no_grad(): residual model.residual(t_old, current_lambda) # 残差极大值附近加密取分位数作为阈值 thresh torch.quantile(residual.abs(), 0.9) hard_points t_old[residual.abs() thresh] random_points torch.rand(n_sample - len(hard_points), 1) * (t_range[1] - t_range[0]) t_range[0] return torch.cat([hard_points, random_points]).requires_grad_(True)5.3 误差表对比从单次运行到系统性调参资源和报告中往往有能力较强的训练策略版本差异如下同一 λ 初始值table_1 的记录间隔是每 10 轮记录一次table_2 是每 50 轮记录一次。若要做训练曲线对比需要对误差表执行重采样统一记录点。其次error_lambda_2_table_1.csv与error_lambda_2_table_2.csv的区别在于后者加入了 L2 先验约束观察两组误差表的 lambda_error 列能直接判断先验约束是否有效。如果加了先验后 lambda 的收敛速度显著加快且没有出现后期的发散说明该参数的观测信息不足先验是必要的如果与原模型的收敛速度几乎一致说明该参数已经被观测数据充分约束删掉先验不影响精度却能减少超参数的调试成本。# 快速分析误差表的收敛速度 python - EOF import pandas as pd import numpy as np df1 pd.read_csv(error_lambda_2_table_1.csv) df2 pd.read_csv(error_lambda_2_table_2.csv) # 看最后 100 轮的参数估计波动 var1 np.var(df1[lambda_estimate][-100:]) var2 np.var(df2[lambda_estimate][-100:]) print(fTable1 last100 variance: {var1:.6f}) print(fTable2 last100 variance: {var2:.6f}) # 若 var2 var1说明先验正则有效 # 若两者都很大说明参数不可辨识需要调整模型结构或增加观测 EOF注意此处的 lambda_estimate 字段名需要和实际 CSV 列名匹配资源中给出的 CSV 可能是英文列名如estimated_lambda或lambda_1等。读取之前先用df.columns确认存在列名不一致不要强行匹配。6. 进阶解耦训练与参数估计的双优化器收敛技巧把四类代码综合到一个工程框架中后训练不稳定性的本源大多数时候不是网络结构而是网络参数与物理参数的更新速率不匹配。这里给出一个在离散和连续、识别与推理场景通用的双优化器技巧将网络权重和物理参数 λ 分配不同的学习率并按照特定 epoch 动态调整两阶段的学习率同时对 optimizer 的配置做有节奏的调整。import torch.optim as optim model ContinuousInferencePINN(n_input1, n_hidden64, n_output2, n_layers5) lambda_param torch.nn.Parameter(torch.tensor([1.2, 0.8, 0.5, 0.3], dtypetorch.float32)) optimizer_net optim.Adam(model.parameters(), lr1e-3) optimizer_lam optim.Adam([lambda_param], lr5e-3) scheduler_net optim.lr_scheduler.CosineAnnealingLR(optimizer_net, T_max500, eta_min1e-5) scheduler_lam optim.lr_scheduler.StepLR(optimizer_lam, step_size200, gamma0.8) for epoch in range(2000): optimizer_net.zero_grad() optimizer_lam.zero_grad() loss loss_function(model, t_data, u_data, t_col, t_ic, u_ic, lambda_param) loss.backward() optimizer_net.step() optimizer_lam.step() scheduler_net.step() scheduler_lam.step()双优化器机制下如果 λ 的学习率是网络权重的 5 倍λ 会快速逼近真实范围但网络还没来得及对新的方程形式建立足够的表征能力残差可能暂时上升。这是正常现象。在误差表上表现为前 200 轮 lambda_error 陡降但 residual_loss 上升200 轮后网络权重逐渐适应残差下降λ 同步微调。观测到这个模式说明优化器节奏整体正常若 lambda_error 下降过慢可把optimizer_lam的 lr 提升到 1e-2但高于 2e-2 时参数发散风险显著增加。另一种实用技巧是对 λ 额外施加软约束。在识别任务中每 50 轮打印一次 λ 的标准差若超过初始值的 30%使用指数移动平均EMA平滑 λ 更新。另外连续时间与离散时间模型对学习率的容忍度差异很大——RK4 结构梯度逆变会放大某些维度因此离散时间识别中 λ 学习率通常要调低至连续时间的 1/2 到 1/3。如果想从 CSV 表直接判断当前学习率是否合理观察 error_table 中数据损失下降到前期最小值的十分之一时所消耗的 epoch 数与后续波动次数即可波动不超过十次的配置往往效果更好。最终将这四套代码串联成一套完整工作流时我的建议路径是先跑连续时间推理确认模型具备基本拟合能力再切换为连续识别反推参数然后改用离散推理处理稀疏观测最后落到离散识别进行联合反演验证——这样每一层新增的不确定性都有了基线对照排查问题时的备选方案也更充分。本文还有配套的精品资源点击获取
PINN四种任务模式：离散/连续时间识别与推理的代码实现
简介面向科学计算与机器学习交叉领域研究者的 PINN 物理信息网络源码包围绕离散时间识别、离散时间推理、连续时间识别、连续时间推理四类典型任务提供可运行的 Python 实现与配套数据。压缩包共 559 个文件其中 txt 文本文件多达 503 个主要保存训练日志、配置说明与结果记录13 个 py 脚本为核心代码完整封装模型训练、推理和误差评估流程12 个 csv 误差表记录不同参数下的 error_table 与 lambda 误差便于横向对比精度另有 eps 与 pdf 格式图表、mat 数据文件等辅助分析材料整体约 475.95MB直接解压即可查阅。已有 1546 人学习下载。通过源码可掌握 PINN 在时间维度上识别与推理的建模思路利用 csv 表中的误差指标和可视化图表能快速定位调参方向适合用于课程设计、论文复现或在此基础上扩展新的物理约束模型。1. PINN 四种任务模式的定位与选型先明确一点这套资源不是单个脚本而是把物理信息神经网络PINN按时间域划分成离散时间识别、离散时间推理、连续时间识别、连续时间推理四条技术路线每条路线配套独立可跑的 Python 源码和 CSV 误差记录表。很多人拿到代码后第一反应是“四个文件差不多”实际上这四者解决的是两类完全不同的问题识别identification是已知系统部分结构、反推未知参数推理inference是给定观测数据、重构整个动态系统的状态轨迹。时间域离散与连续的区别则决定了网络如何构造时间导数项——离散时间用 Runge-Kutta 等积分步进结构连续时间则依赖自动微分直接计算微分算子。适用于三类读者正在复现 PINN 基线实验的研究生、需要把物理约束嵌入神经网络的算法工程师、以及做数字孪生或动力学系统建模的从业者。如果只是为了跑通 PINN 做数据拟合这四个代码可能显得“重”但一旦涉及方程形式未知、参数反演或噪声观测下的状态估计离散与连续、识别与推理的区分就是绕不开的起点。下面从时间域的数学设定差异切入逐步拆解每一类代码的模型结构、损失函数和调试要点最后说明如何用附带的 CSV 误差表做收敛性验证。2. 连续时间推理与识别自动微分驱动的物理约束2.1 连续时间推理的方程约束构造连续时间场景假设系统的状态由常微分方程组ODE控制已知方程形式但未知部分参数或初值。PINN 推理任务的核心是构造一个神经网络 u_θ(t)使其输出同时满足两个条件一是逼近观测数据二是让方程残差尽可能小。这里的方程残差定义为R(t) du_θ/dt - f(u_θ, t, λ)其中 λ 是待学习的物理参数。在 PyTorch 中实现时时间导数 du_θ/dt 直接通过torch.autograd.grad对输入时间求一阶导得到这正是 PINN 相比传统数值方法的根本差异——不需要离散化微分算子而是利用自动微分链式法则把物理约束嵌入损失函数。import torch import torch.nn as nn class ContinuousInferencePINN(nn.Module): def __init__(self, n_input1, n_hidden50, n_output2, n_layers4): super().__init__() layers [nn.Linear(n_input, n_hidden), nn.Tanh()] for _ in range(n_layers - 2): layers [nn.Linear(n_hidden, n_hidden), nn.Tanh()] layers [nn.Linear(n_hidden, n_output)] self.net nn.Sequential(*layers) def forward(self, t): return self.net(t) def residual(self, t, lambda_param): t t.clone().requires_grad_(True) u self.forward(t) du_dt torch.autograd.grad( u, t, grad_outputstorch.ones_like(u), create_graphTrue )[0] # 以 Lotka-Volterra 两变量系统为例 x, y u[:, 0:1], u[:, 1:2] fx lambda_param[0] * x - lambda_param[1] * x * y fy -lambda_param[2] * y lambda_param[3] * x * y return du_dt - torch.cat([fx, fy], dim1)这段代码的关键有三处。第一requires_grad_(True)必须在残差计算之前重新设置因为每次前向传播的时间张量可能来自一个不需要梯度的 DataLoader。第二grad_outputstorch.ones_like(u)表示对输出逐元素求梯度这一参数漏写是常见的报错点——autograd 要求显式提供与输出同形状的 grad_outputs否则无法计算雅可比向量积。第三残差返回的是 du_dt 与方程右侧的差值维度必须与网络输出严格一致如果系统是三维状态空间而输出维度误设为 2后续损失拼接时会抛出形状不匹配错误。2.1.1 损失函数的多目标加权推理任务的损失函数由三部分构成观测数据拟合损失、方程残差损失、初始条件损失。权重 λ_data、λ_res 和 λ_ic 的选择直接影响训练结果。一般经验值是 λ_res 取 1.0、λ_data 为 0.1 到 1.0、λ_ic 为 1.0。原因在于方程残差是物理约束需要优先满足观测数据若噪声较大过高的 λ_data 会把噪声完全拟合成动力学模式反而使物理推理失真。def loss_function(model, t_data, u_data, t_col, t_ic, u_ic, lambda_param): u_pred model(t_data) loss_data torch.mean((u_pred - u_data) ** 2) residual model.residual(t_col, lambda_param) loss_res torch.mean(residual ** 2) u_ic_pred model(t_ic) loss_ic torch.mean((u_ic_pred - u_ic) ** 2) return 0.5 * loss_data 1.0 * loss_res 1.0 * loss_ic这里的配置点是配点 t_col 的采样方式。常见做法是在时间区间内均匀采样 2000 个点如果系统存在高频振荡可以在预测值变化剧烈的区域加密采样。另一个实操细节是 λ_param 必须注册为nn.Parameter而非普通 Tensor否则优化器不会更新它。2.2 连续时间识别的参数反演策略识别与推理的区别在于推理只重建状态轨迹识别则要在重建轨迹的同时反推出方程中的未知参数 λ。实现上两者共用同一个网络差异体现在优化器的参数列表中——识别任务必须把 λ 传入优化器。这里有个隐蔽的边界如果初始 λ 与实际值偏离过大梯度会被残差项主导导致 λ 的更新方向受方程右侧非线性项支配而不稳定。缓解手段是分阶段训练前 500 轮固定 λ 只优化网络权重待残差降到一定阈值后再联合优化。optimizer torch.optim.Adam([ {params: model.parameters(), lr: 1e-3}, {params: [lambda_param], lr: 5e-3} ]) for epoch in range(total_epochs): if epoch 500: # 解锁 lambda 的优化前面固定可以避免初期发散 lambda_param.requires_grad_(True) optimizer.zero_grad() loss loss_function(model, t_data, u_data, t_col, t_ic, u_ic, lambda_param) loss.backward() optimizer.step()注意lambda_param.requires_grad_(True)的时机。严格说如果定义时没有设置 requires_gradAdam 会直接跳过这个参数所以不建议用 if 分支控制 requires_grad 开关更稳妥的做法是用两个优化器——前 500 轮只用optimizer_model之后换用包含两个参数组的联合优化器。另外识别的评估指标不只是最终 λ 的绝对误差还要看 λ 在训练末期的波动幅度——如果连续 200 轮的变化量小于 1e-4说明收敛如果一直振荡说明观测数据的信息量不足以唯一确定参数需要增加配点密度或引入正则化项。3. 离散时间推理Runge-Kutta 积分器与观测融合3.1 从连续残差到离散步进结构当观测数据只在若干时间步上获得且相邻步之间没有高频率采样连续时间模型会面临一个可观察性问题方程残差在稀疏区间内可能对应多条轨迹网络没有足够约束选出正确的一条。离散时间推理的做法是利用 Runge-Kutta 积分器把状态按固定步长在时间上推进用一步积分后的预测值与下一时刻观测值做比较这样每一次时间步的更新都是一个“预测-校正”闭环观测信息被逐步注入而非一次性约束。3.1.1 经典四阶 Runge-Kutta 的 PINN 化改造以四阶 Runge-KuttaRK4为例在 PINN 中构造步进过程如下给定当前时刻的状态 u_n通过 RK4 计算 u_{n1} 的预测值损失函数比较该预测值与观测点 u_{n1}^{obs} 的偏差。RK4 的积分步长 Δt 属于超参数需要结合观测数据的时间间隔来设定。以附着动力学模型中常见的步长选择为例当观测间隔为 0.1 时Δt 通常设置在 0.01 到 0.05 之间步长过大会引入积分误差过小则计算量成倍增加。def rk4_step(model, t_n, u_n, dt, lambda_param): k1 model.rhs(t_n, u_n, lambda_param) k2 model.rhs(t_n dt / 2, u_n dt * k1 / 2, lambda_param) k3 model.rhs(t_n dt / 2, u_n dt * k2 / 2, lambda_param) k4 model.rhs(t_n dt, u_n dt * k3, lambda_param) return u_n dt / 6 * (k1 2 * k2 2 * k3 k4)这里的model.rhs是用神经网络表示的动力学函数输入的是当前时刻的状态和时间输出状态的导数。与连续时间模型中直接对时间自动微分不同离散时间法的核心在网络内部不再计算时间导数而是把网络整体当作一个 ODE 右侧函数时间的推进由 RK4 公式完成。这种设计的优势是训练时可以批量处理多个时间步——每个时刻的输出都经过整个 RK4 链路梯度会通过 K 系数回传到网络的每一层相当于在时间轴上做了反向传播。3.1.2 稀疏观测下的损失权重分配离散时间推理的损失函数按时间步加权越靠后的时间步预测误差会被 RK4 积分过程累积放大因此权重应随步数增加而递减。具体操作中步进索引 i 对应的损失权重取 1/(i1) 或 0.9^i。前几个时间步的观测对模型参数更新的贡献更大因为误差尚未经过长时间积分污染。def discrete_inference_loss(model, t_seq, u_seq, dt, lambda_param, gamma0.9): total_loss 0.0 u_current u_seq[0] for i in range(1, len(t_seq)): u_pred_next rk4_step(model, t_seq[i-1], u_current, dt, lambda_param) obs_next u_seq[i] weight gamma ** i total_loss weight * torch.mean((u_pred_next - obs_next) ** 2) # 这里要注意训练时下一时刻的输入应使用预测值还是观测值 u_current u_pred_next.detach() # 阻断梯度流动避免误差累积爆炸 return total_loss代码中的u_current u_pred_next.detach()是决定训练稳定性的关键。如果直接把预测值传入下一步运算梯度会沿着时间维度回传很多步导致梯度消失或爆炸断开后每一步的 RK4 结构相对独立虽然牺牲了长期依赖但稳定性显著提升。这是与序列模型用 Teacher Forcing 和 Free Running 进行取舍时的常见做法在可观察性足够时先用观测值做引导后续再切换为全预测模式微调。3.2 初始化与积分步长的联动关系离散时间模型对初始化比较敏感。如果网络的随机初始化使 RK4 的 K 系数计算值过大u 会在最初几步就偏离合理范围丢失有效信息。常见初始化方案是将网络最后一层权重置为较小的均值为 0、标准差 0.01且偏置初始化为满足先验的平均变化率。另一个做法是把 RK4 内部的 dt 用一个可学习的缩放因子修正——网络初始预测的导数尺度与实际系统偏差较大时可通过这个因子让模型先学到正确的变化率再逐步调整。成功经验是步长先取观测间隔的 1/10 保证积分稳定性训练 300 轮后逐步把 dt 提到 1/5 以覆盖更长的动力学特征。4. 离散时间识别离散观测下的参数反演实现离散时间识别是四条路线中最“工程化”的一条逻辑上需要在 RK4 步进结构下加入参数 λ 的梯度更新实现上要注意 λ 的梯度经过多条 RK4 路径后出现的量级失衡问题——如果 λ 初始值偏大方程右侧的乘积项会让导数预测值成倍放大K 系数随之膨胀。误差表error_lambda_1_table_2.csv中记录的正是不同 λ 初始值下识别误差的对比数据下面先用代码说明实现再基于误差表解释参数-初始值的敏感性。def discrete_identification_loss(model, t_seq, u_seq, dt, lambda_param): u u_seq[0] loss 0.0 all_pred [u] for i in range(len(t_seq) - 1): u_pred rk4_step(model, t_seq[i], u, dt, lambda_param) loss torch.mean((u_pred - u_seq[i 1]) ** 2) u u_pred # 这里不 detach让 lambda 能获得跨步梯度 all_pred.append(u_pred) return loss, all_pred与推理代码的显著区别是在整个时间序列上保留计算图不执行 detach实现“检测参数在长期积分效果上的影响”。但这会造成训练负担成倍增加——时间步数为 N 时反向传播的复杂度为 O(N × 网络层数)。对于一个 200 步的时间序列GPU 显存占用大约是推理模式的 20 倍以上这是限制离散时间识别可扩展性的主要瓶颈。缓解方案是分段优化——把序列切成数段每段保留计算图段间用 detach 切割相当于把长期依赖折中在数个子时段内。交替使用两种策略前 100 轮用分段模式快速逼近大致参数范围后续切换为全程模式精调。def segment_identification_loss(model, t_seq, u_seq, dt, lambda_param, seg_len20): total_loss 0.0 n_segments (len(t_seq) - 1) // seg_len # 每个段的起点直接使用观测值保证段内轨迹误差是真实偏差 for seg in range(n_segments): start seg * seg_len u u_seq[start] for i in range(start, start seg_len): u_pred rk4_step(model, t_seq[i], u, dt, lambda_param) total_loss torch.mean((u_pred - u_seq[i 1]) ** 2) u u_pred # 段末不将梯度传到下一段防止段间误差累积 u u.detach() return total_loss分段模式中每个段长度为 20 步显存占用约为全程模式的 1/10。观察配套误差表error_lambda_2_table_1.csv当 λ2 的初始值偏离真值超过 50% 时分段模式的识别误差比全程模式高出约 30%说明分段虽然省显存但精调的精度下降当初始值偏离在 20% 以内时两者误差几乎没有差异这就是工程上“分阶段先粗后精”的依据。另一个值得注意的点是识别误差表里的误差不是每个时间步的均方误差而是 λ 估计值与真实值的相对误差。因此训练结束后需单独打印 λ 的最终收敛值警惕“总损失下降但参数没收敛”的假收敛情况——优化器很容易把网络权重调得足够扭曲以拟合错误的 λ 值这种现象在识别任务中被称为“参数补偿效应”。4.1 多参数联合识别的正则化与可辨识性当系统存在多个未知参数时离散时间识别面临可辨识性边界——不同参数组合可能在同一组观测数据下产生几乎相同的轨迹。此时单纯的最小化损失无法统一确定参数需要添加先验正则化项。工程实现中针对性处理方式为在损失函数中加入 λ 与实际值估计先验的偏差项 λ_prior使用 L2 正则回归目标为初始猜测 λ0def identification_loss_with_prior(model, t_seq, u_seq, dt, lambda_param, lambda_prior, prior_weight0.05): loss, _ discrete_identification_loss(model, t_seq, u_seq, dt, lambda_param) prior_loss prior_weight * torch.sum((lambda_param - lambda_prior) ** 2) return loss prior_loss需要注意的是prior_weight 不应设置过大。若超过 0.1λ 会被强制拉向先验值而网络为了降低数据拟合误差会倾向于“把参数承担的作用转嫁给其他可调参数”反而扭曲了物理含义。比较合理的做法是先不做先验训练一轮观察 λ 的收敛趋势若明显向某一方向偏移再根据趋势设置一个偏向性的先验这是连续与离散识别任务中共通的传统技巧。5. 从 CSV 误差表反推训练诊断收敛、稳定性与数值病态这组资源中最容易被忽略但价值最高的是error_table_1.csv、error_lambda_1_table_1.csv、error_lambda_2_table_1.csv及对应 table_2 版本共六个误差表记录不同实验配置下的误差演化或分布。这部分并不是训练的必备输入而是用来定位模型失败模式的“病历表”——通过读表的思路能快速判断一个 PINN 模型是在数据拟合阶段失败、在物理约束阶段发散还是在参数更新阶段陷入病态。5.1 读懂误差表的三个关键字段字段名含义诊断指向iteration当前训练轮次或时间步曲线太陡或振荡说明学习率过高data_loss观测数据拟合损失持续高位说明网络容量不足或观测有噪声residual_loss方程残差损失指数下降后停滞说明需要增加配点或调整权重lambda_error参数估计相对误差数值跳动说明病态单调不降说明不可辨识table_1 与 table_2 的结构差异体现了识别与推理两种任务下误差指标的关注点不同。连续时间推理的 error_table 通常记录每轮训练后的总损失和验证集残差离散时间识别的 lambda 误差表则专门记录每个候选 λ 在当前 epoch 的估计值变化。若 lambda 误差先下降后上升典型原因是优化器在残差损失收敛后继续优化数据损失时破坏了参数稳定性——常见做法是冻结 λ单独对网络权重进行若干轮微调再解冻参数重新联合优化。5.2 用误差表定位 PINN 典型失败场景场景一residual_loss 持续在 1e-1 量级不下降。原因指向配点 t_col 在相空间覆盖面不够——系统若有多个吸引子或振荡模态均匀采样可能漏掉关键区域。调整方向是改用自适应配点从当前残差较大的位置周围加密采样迭代更新配点集合。场景二残余损失下降但 lambda_error 不降反升。这是参数补偿效应的直接体现另一条排查思路是路径检查训练过程中方程右侧的每一项对损失的梯度贡献对 lotka-volterra 系统如果 λ1 对应的线性项梯度比非线性交互项的小很多损失函数会将参数误差推给高梯度项此时需要对 λ 进行归一化或给低梯度项更大的权重。场景三离散时间推理的验证损失在某个 epoch 后突然陡增。这通常是 RK4 步长 dt 过大导致积分路径发散——误差表上会表现为该 epoch 之后的数据记录中断或残余损失出现数量级跳变。# 一个实际可用的自适应配点方法 def adaptive_resample(model, t_old, t_range, n_sample2000): # 先评估当前残差 with torch.no_grad(): residual model.residual(t_old, current_lambda) # 残差极大值附近加密取分位数作为阈值 thresh torch.quantile(residual.abs(), 0.9) hard_points t_old[residual.abs() thresh] random_points torch.rand(n_sample - len(hard_points), 1) * (t_range[1] - t_range[0]) t_range[0] return torch.cat([hard_points, random_points]).requires_grad_(True)5.3 误差表对比从单次运行到系统性调参资源和报告中往往有能力较强的训练策略版本差异如下同一 λ 初始值table_1 的记录间隔是每 10 轮记录一次table_2 是每 50 轮记录一次。若要做训练曲线对比需要对误差表执行重采样统一记录点。其次error_lambda_2_table_1.csv与error_lambda_2_table_2.csv的区别在于后者加入了 L2 先验约束观察两组误差表的 lambda_error 列能直接判断先验约束是否有效。如果加了先验后 lambda 的收敛速度显著加快且没有出现后期的发散说明该参数的观测信息不足先验是必要的如果与原模型的收敛速度几乎一致说明该参数已经被观测数据充分约束删掉先验不影响精度却能减少超参数的调试成本。# 快速分析误差表的收敛速度 python - EOF import pandas as pd import numpy as np df1 pd.read_csv(error_lambda_2_table_1.csv) df2 pd.read_csv(error_lambda_2_table_2.csv) # 看最后 100 轮的参数估计波动 var1 np.var(df1[lambda_estimate][-100:]) var2 np.var(df2[lambda_estimate][-100:]) print(fTable1 last100 variance: {var1:.6f}) print(fTable2 last100 variance: {var2:.6f}) # 若 var2 var1说明先验正则有效 # 若两者都很大说明参数不可辨识需要调整模型结构或增加观测 EOF注意此处的 lambda_estimate 字段名需要和实际 CSV 列名匹配资源中给出的 CSV 可能是英文列名如estimated_lambda或lambda_1等。读取之前先用df.columns确认存在列名不一致不要强行匹配。6. 进阶解耦训练与参数估计的双优化器收敛技巧把四类代码综合到一个工程框架中后训练不稳定性的本源大多数时候不是网络结构而是网络参数与物理参数的更新速率不匹配。这里给出一个在离散和连续、识别与推理场景通用的双优化器技巧将网络权重和物理参数 λ 分配不同的学习率并按照特定 epoch 动态调整两阶段的学习率同时对 optimizer 的配置做有节奏的调整。import torch.optim as optim model ContinuousInferencePINN(n_input1, n_hidden64, n_output2, n_layers5) lambda_param torch.nn.Parameter(torch.tensor([1.2, 0.8, 0.5, 0.3], dtypetorch.float32)) optimizer_net optim.Adam(model.parameters(), lr1e-3) optimizer_lam optim.Adam([lambda_param], lr5e-3) scheduler_net optim.lr_scheduler.CosineAnnealingLR(optimizer_net, T_max500, eta_min1e-5) scheduler_lam optim.lr_scheduler.StepLR(optimizer_lam, step_size200, gamma0.8) for epoch in range(2000): optimizer_net.zero_grad() optimizer_lam.zero_grad() loss loss_function(model, t_data, u_data, t_col, t_ic, u_ic, lambda_param) loss.backward() optimizer_net.step() optimizer_lam.step() scheduler_net.step() scheduler_lam.step()双优化器机制下如果 λ 的学习率是网络权重的 5 倍λ 会快速逼近真实范围但网络还没来得及对新的方程形式建立足够的表征能力残差可能暂时上升。这是正常现象。在误差表上表现为前 200 轮 lambda_error 陡降但 residual_loss 上升200 轮后网络权重逐渐适应残差下降λ 同步微调。观测到这个模式说明优化器节奏整体正常若 lambda_error 下降过慢可把optimizer_lam的 lr 提升到 1e-2但高于 2e-2 时参数发散风险显著增加。另一种实用技巧是对 λ 额外施加软约束。在识别任务中每 50 轮打印一次 λ 的标准差若超过初始值的 30%使用指数移动平均EMA平滑 λ 更新。另外连续时间与离散时间模型对学习率的容忍度差异很大——RK4 结构梯度逆变会放大某些维度因此离散时间识别中 λ 学习率通常要调低至连续时间的 1/2 到 1/3。如果想从 CSV 表直接判断当前学习率是否合理观察 error_table 中数据损失下降到前期最小值的十分之一时所消耗的 epoch 数与后续波动次数即可波动不超过十次的配置往往效果更好。最终将这四套代码串联成一套完整工作流时我的建议路径是先跑连续时间推理确认模型具备基本拟合能力再切换为连续识别反推参数然后改用离散推理处理稀疏观测最后落到离散识别进行联合反演验证——这样每一层新增的不确定性都有了基线对照排查问题时的备选方案也更充分。本文还有配套的精品资源点击获取