从硬球碰撞数据中学习H函数:用DeepSets与孪生网络发现时间箭头
1. 项目概述与核心动机在统计物理和热力学的核心地带有一个困扰了物理学家们一个多世纪的基本问题宏观世界那清晰、单向的时间箭头究竟是如何从微观粒子那些完全可逆的碰撞规则中“涌现”出来的这就是著名的“可逆性悖论”。玻尔兹曼在1872年提出的H定理为这个谜题提供了一个优雅的数学解答。它通过定义一个称为H函数的量并证明其在孤立系统演化过程中永不增加dH/dt ≤ 0为热力学第二定律的熵增原理奠定了微观动力学基础。然而H定理的证明依赖于“分子混沌”假设这本身就是一个关于系统初始条件和演化的统计性断言其普适性至今仍是理论物理中一个活跃的讨论话题。作为一名长期关注计算物理与机器学习交叉领域的研究者我一直在思考我们能否不依赖于先验的理论假设而是直接从系统演化的“数据”中让机器自己“发现”这个标志时间方向的H函数这不仅仅是一个理论游戏更是一种全新的物理规律发现范式。如果成功意味着我们可以将这种“从数据中学习物理不变性”的方法推广到那些理论尚未完善或过于复杂的系统比如活性物质、复杂流体或多原子分子气体。最近一项结合了DeepSets架构和孪生神经网络思想的工作在这个方向上迈出了令人兴奋的一步。它利用硬球碰撞这一经典的、微观可逆的模型生成了大量的动力学轨迹数据然后训练一个神经网络去学习一个与H函数行为一致的标量函数。这个想法极具启发性我们不直接告诉模型物理定律是什么而是通过设计一个巧妙的损失函数告诉模型“时间箭头”的方向——即系统演化的“前”与“后”。模型需要自己从混沌的粒子速度数据中提炼出那个能够区分时序、单调变化的“隐秩序”而这个秩序理论上就应该对应着H函数或其仿射变换。在接下来的内容里我将详细拆解这项工作的技术路径、实现细节并分享我在复现和思考过程中获得的一些实操心得与扩展思路。2. 理论基础从玻尔兹曼方程到H定理要理解机器学习方法的目标我们必须先回到经典的物理图景中。考虑一个由大量全同粒子组成的稀薄气体系统。玻尔兹曼的伟大洞见在于他不再追踪每一个粒子的精确轨迹这在数学和计算上都是不可能的而是转向描述粒子在相空间中的概率分布函数 f(r,v, t)。这个函数告诉我们在时刻t位于空间位置r附近速度在v附近的粒子数密度是多少。2.1 玻尔兹曼方程的精髓系统的演化由玻尔兹曼方程描述。在无外场且忽略空间不均匀性的简化情况下这也是后续模拟的设定方程简化为 ∂f(v, t)/∂t (∂f/∂t)_coll等式右边被称为“碰撞积分”它复杂地描述了粒子间碰撞如何改变速度分布。其具体形式对于硬球或硬盘模型是已知的它包含了所有可能的碰撞对和碰撞参数核心是“碰撞前”分布乘积与“碰撞后”分布乘积的差值。这个方程的建立引入了一个关键的、也是备受争议的假设分子混沌假设。它假定在碰撞发生前两个粒子的速度是互不相关的。正是这个假设向系统中注入了时间不对称性打破了微观动力学的可逆对称。注意分子混沌假设并非逻辑必然它是为了从复杂的多体关联中“封闭”方程而做的统计性近似。朗福德定理在一定极限如玻尔兹曼-格拉德极限和短时标下证明了其合理性但对于一般情况和长时间演化它仍然是一个有效的“信仰之跃”。我们机器学习方法的有趣之处在于它试图绕过这个先验假设直接从数据中“嗅出”时间的方向。2.2 H定理与时间箭头玻尔兹曼进一步定义了一个泛函即H函数 H(t) ∫ f(v, t) ln f(v, t) dv通过对玻尔兹曼方程进行一番巧妙的运算可以证明 dH/dt ≤ 0这个不等式就是H定理的核心。它意味着只要分子混沌假设成立H函数就会随时间单调递减或不变直到系统达到平衡态此时分布函数变为著名的麦克斯韦-玻尔兹曼分布且dH/dt 0。这里的物理图像极其深刻H函数在平衡态取最小值这与热力学熵在平衡态取最大值对于孤立系统是等价的实际上H -S/k_B 常数。因此H函数的单调递减为宏观世界观测到的熵增时间箭头提供了一个微观的、动力学的解释。系统从一种有序非平衡H值较高的状态不可逆地演化到最无序平衡H值最低的状态。我们机器学习模型要学习的正是这样一个函数。模型输入是某一时刻所有粒子的速度集合 {v_i(t)}输出是一个标量h(t)。我们希望训练得到的h(t)能够满足对于系统演化的真实轨迹当 t1 t2 时有 h(t1) ≥ h(t2)。换句话说模型要学会从一堆看似混乱的速度数据中抽取出那个隐含的、单向变化的“时间之箭”。3. 数据生成硬球碰撞模拟的实操要点理论是优美的但我们需要数据来喂养模型。研究选择了二维硬盘系统作为 playground这是统计物理中最经典、研究最透彻的模型之一。它的规则简单纯粹粒子是刚性的圆盘除了发生弹性碰撞时瞬间交换动量其余时间均做匀速直线运动。这种简单性使得其微观动力学是完全确定性和时间可逆的但大量粒子的集体行为却展现出不可逆的统计趋势是验证我们想法的完美沙盒。3.1 模拟参数与初始化策略我参考了原始研究并进行了复现以下是关键的模拟参数和初始化细节系统规模N 1000个全同硬盘。这个数量足够产生显著的统计行为又不会让模拟计算成本过高。空间设置粒子被放置在一个边长为1的二维正方形盒子中采用周期性边界条件。这相当于一个无限大的平铺空间避免了边界效应带来的复杂影响是处理大量粒子系统的标准做法。粒子属性每个盘质量m1直径d 6×10^{-3}。较小的直径相对于盒子保证了系统的低密度更接近玻尔兹曼方程适用的“稀薄气体” regime。初始化这是整个实验设计的精妙之处。为了清晰地观察到H函数的变化我们需要系统从一个明显的非平衡态出发。位置随机均匀分布在盒子内确保无重叠。速度方向随机均匀分布。速度大小所有粒子被赋予相同的速率v 0.12。这是一个高度有序的状态速度分布是一个δ函数远离平衡态的麦克斯韦-玻尔兹曼分布因此初始H值会很高。实操心得初始化速度大小均匀分布而非平衡的麦克斯韦分布至关重要。如果直接从平衡分布初始化H函数从一开始就接近最小值其随时间下降的趋势会非常微弱甚至被数值噪声淹没导致模型难以学习。这种“从有序到无序”的鲜明对比为模型提供了强烈的学习信号。3.2 碰撞检测与动力学推进硬球模型的模拟核心是事件驱动Event-Driven或时间驱动Time-Driven的碰撞处理。本研究采用了时间驱动法即设定一个固定的时间步长dt逐步推进系统。运动在每个时间步根据粒子当前速度更新其位置r_i(tdt) r_i(t) v_i(t) * dt。碰撞检测检查所有粒子对(i, j)。计算相对位移 Δrr_j -r_i考虑周期性边界条件的最小镜像约定。如果 |Δr| d则判定为重叠即发生了碰撞。碰撞处理对于发生碰撞的粒子对进行二维弹性碰撞计算。动量仅在两球心连线方向交换。计算公式为v_i v_i - (Δv·n)nv_j v_j (Δv·n)n其中 Δvv_i -v_jn Δr/ |Δr| 是单位法向量。同时必须施加碰撞条件仅在相对速度在法向的分量指向内时即 Δv·n 0才应用上述规则否则粒子正在分离不应发生碰撞。忽略此条件会导致粒子“粘”在一起的数值错误。数据记录记录下每个时间步所有粒子的速度向量 {v_i(t)}。总共模拟 N_run 500 条独立轨迹每条轨迹包含 N_t 500 个时间步。关于时间步长dt的选择这是一个需要权衡的参数。dt太大会错过许多碰撞事件见图2的模拟结果dt太小计算成本激增。通过分析检测到的碰撞次数随帧率1/dt的变化发现当帧率高于约60 Hz即 dt ~0.0167后碰撞次数趋于稳定。因此选择 dt 0.016 作为模拟步长在精度和效率间取得了良好平衡。模拟完成后我们可以通过统计速度分布来验证系统是否趋向平衡。如图3所示随着时间演化初始的均匀速率分布一个圆环逐渐扩散最终与理论预测的二维麦克斯韦-玻尔兹曼分布吻合得很好。这直观地展示了系统从非平衡态向平衡态的弛豫过程也是H函数下降的宏观体现。4. 模型架构DeepSets与置换不变性现在我们有了数据对于500条轨迹每条有500个时间快照每个快照是1000个粒子的二维速度向量。如何设计一个模型来处理这样的数据并学习我们想要的标量函数h(t)这里面临一个核心挑战粒子的标签索引i是任意的、无物理意义的。交换两个粒子的数据系统的物理状态完全没有改变。因此我们学习的函数h必须对输入粒子集合的顺序具有置换不变性。即对于任何排列π都有 h({v_1,v_2, ...,v_N}) h({v_π(1),v_π(2), ...,v_π(N)})。4.1 DeepSets 架构解析DeepSets 是处理此类集合数据的经典神经网络架构。其核心思想非常直观且数学上优美任何关于集合的置换不变函数都可以分解为“对每个元素进行变换”和“对变换结果进行聚合”两个步骤。具体到我们的模型函数 h(V_t) 被构造为 h(V_t) ρ ( Σ_{i1 to N} φ(v_i) )其中φ: 一个共享的多层感知机MLP称为“编码器”。它独立地处理每一个粒子的速度向量v_i将其映射到一个高维特征空间。所有粒子共享同一套φ网络的权重。Σ (求和): 聚合函数。这里采用了最简单的元素求和sum pooling。求和操作天然就是置换不变的——无论怎么打乱顺序总和不变。其他如均值、最大值等对称函数也可用但求和能保留更多的信息例如粒子总数。ρ: 另一个MLP称为“解码器”。它接收聚合后的全局特征向量即所有φ(v_i)的和并输出最终的标量值 h。这种结构强制模型满足置换不变性。模型无法通过记忆粒子的顺序来“作弊”它必须学习基于每个粒子个体属性速度及其集体统计特性的表示。4.2 网络实现细节在我的PyTorch复现中φ和ρ都实现为简单的全连接网络。具体超参数参考了原研究的设置并进行了微调输入层接收一个粒子的二维速度向量 (v_x, v_y)。φ网络包含两个隐藏层每层神经元数量为nhidden一个关键的超参数研究尝试了16, 32, 64。使用ReLU激活函数引入非线性。聚合对批次中所有粒子的φ输出进行求和。ρ网络同样包含两个隐藏层神经元数也为nhidden使用ReLU激活。输出层为线性层输出单个标量值。通过调整nhidden我们可以控制模型的容量。容量太小可能欠拟合无法捕捉复杂模式容量太大可能过拟合或导致训练不稳定。注意事项在实现DeepSets时需要特别注意批处理batching的维度。通常我们会将数据组织为(batch_size, num_particles, feature_dim)的张量。φ网络需要独立地应用于每个粒子的特征上这可以通过view和reshape操作或者使用torch.nn.Linear配合适当的维度变换来实现。聚合操作求和则沿着num_particles维度进行。确保这些操作不影响批次维度是关键。5. 损失函数设计从单调性约束到稳定学习模型架构决定了如何表示函数而损失函数则定义了什么是好的函数。我们的目标很明确学到的函数h(t)应该随时间单调不增。如何将这一物理约束转化为可优化的损失函数5.1 基础损失函数及其问题最直接的思路是采用对比学习或排序学习的范式。对于同一条轨迹中相邻的两个时间步 t 和 t1我们希望模型输出满足 h(V_t) ≥ h(V_{t1})。一个简单的损失函数可以设计为 L^{(0)}t h(V_t) - h(V{t1}) λ [h(V_{t1}) - h(V_t)]^2让我们拆解这个损失第一项 h(V_t) - h(V_{t1})这是单调性约束的核心。如果 h(V_t) h(V_{t1})即函数值随时间增加这一项为负产生一个正的损失梯度推动模型调整参数以减少这种情况。如果 h(V_t) ≥ h(V_{t1})这一项非负对损失没有正贡献理想情况。第二项 λ [h(V_{t1}) - h(V_t)]^2这是一个正则化项。λ 是一个超参数研究中设为0.2。它的作用是防止模型走向一个平凡的、无信息的解。试想如果模型简单地将所有输入都输出为同一个常数那么第一项恒为0损失看起来很小但这显然不是我们想要的H函数。第二项惩罚输出值的变化鼓励模型在满足单调性的前提下让相邻时间步的输出尽可能接近从而学习到一个平滑变化的函数而不是剧烈振荡的。然而在实际训练中发现使用这个损失函数学到的 h(t) 曲线振荡非常剧烈。模型似乎找到了一个“取巧”的策略它让 h(t) 在某些步长大幅下降在另一些步长又大幅上升但总体趋势可能是下降的从而在平均意义上降低损失。这不符合我们对H函数平滑、单调递减的物理预期。5.2 改进的Leaky ReLU损失函数问题的根源在于基础损失函数对“正确行为”h下降和“错误行为”h上升的惩罚是对称的。但根据H定理我们只要求函数不增加≤ 0并没有要求它必须严格递减。也就是说函数在一段时间内保持恒定斜率为0是完全允许的甚至是平衡态的特征。因此一个更符合物理的改进是引入非对称惩罚。研究采用了带泄露的ReLULeaky ReLU函数来重构损失项 L^{α}t LeakyReLU( h(V_t) - h(V{t1}), α ) λ [h(V_{t1}) - h(V_t)]^2其中LeakyReLU(x, α) x if x 0 else α * x。这里α是一个介于0和1之间的参数。这个修改带来了关键的变化当 h(V_t) - h(V_{t1}) 0 时正确函数下降损失为正值x模型到惩罚。等等为什么正确还要受罚这里需要仔细理解我们的目标是让这个差值尽可能为0非正即可而不是尽可能大。过大的下降虽然满足不等式但并不是H函数的典型行为H函数通常是指数弛豫到平衡。所以轻微的惩罚可以防止函数下降过快。当 h(V_t) - h(V_{t1}) 0 时错误函数上升损失为 α * x一个负值。由于x是负的α*x就是一个正数但它的绝对值是 |αx|。而正确情况下的惩罚是 |x|。因为 0 ≤ α 1所以对于相同的差值绝对值函数上升错误所受到的惩罚 |αx|要小于函数下降正确所受到的惩罚 |x|。这种非对称性使得模型更倾向于让差值停留在0附近即函数不变或者缓慢下降而不是剧烈振荡。这更贴合H函数在接近平衡时变化平缓的物理图像。研究中还采用了动态调整α的策略随着训练周期epoche从0增加到总周期Ne让 α(e) 从0线性增加到接近1如0.99。这样在训练初期非对称性很强鼓励模型快速找到一个大致单调的解后期α接近1损失函数趋近于对称进行微调。6. 训练流程与模型评估6.1 训练设置与优化我们将所有500条轨迹的时序数据组织成样本对 (V_t, V_{t1})。模型训练采用标准的迷你批mini-batch随机梯度下降使用Adam优化器学习率设置为4e-4。训练过程可以理解为让一个孪生神经网络学习排序。如图4所示两个共享权重的DeepSets模型即“孪生”分别处理时刻t和t1的粒子速度集合输出h(V_t)和h(V_{t1})。损失函数计算这对输出的差异并按照上述规则进行优化。通过在所有时间步和所有轨迹上平均损失模型逐渐学会赋予更晚时间状态以更小或相等的h值。6.2 临界周期与训练动态一个有趣的现象是随着α在训练中逐渐增大损失函数会经历一个拐点。在训练早期由于α较小非对称惩罚明显模型很容易找到一个使损失为负的解因为错误上升的惩罚很轻。随着α增大至接近1损失函数越来越对称那个使损失为负的“舒适区”逐渐消失。最终损失会由负转正。研究将损失首次变为正的训练周期称为“临界周期” e_c。表II展示了不同隐藏层维度nhidden模型的最终损失和临界周期。可以看到模型容量越大nhidden64临界周期来得越晚e_c111说明模型有更强的能力去拟合一个在更严格的对称损失下仍然表现良好的函数。最终所有模型在α接近1时其损失 ⟨L^{α_final}⟩ 都非常接近于0这表明模型成功学习到了一个几乎完全满足 h(V_t) ≈ h(V_{t1}) 的函数即一个近乎常数的函数。这听起来像是失败了别急这恰恰是成功的第一步——模型学会了识别平衡态H函数不变。关键在于在达到平衡之前它必须经历一个下降的过程。6.3 与真实H函数的对标模型输出的h(t)与真实的H(t)之间存在两个固有的自由度平移不变性损失函数只关心h(t)的相对差值因此给所有输出加上一个常数b损失不变。缩放不变性如公式(14)所示将h(t)乘以一个常数c同时将正则化系数λ除以c损失函数整体乘以c优化动态不变。因此我们只能期望模型学到与H函数相差一个仿射变换的函数˜h(t) a * h(t) b。为了比较我们需要为每条轨迹n单独拟合出最优的缩放因子a^[n]和平移量b^[n]通过最小化 ˜h^ n 与真实 H^ n 之间的均方误差MSE来实现。真实H函数的计算由于我们只有离散的粒子速度样本无法得到连续的分布函数f(v, t)。需要采用核密度估计或直方图法来近似。研究中采用了直方图法将速度空间划分为细小的网格bin统计每个bin中的粒子数归一化后得到离散的概率分布然后代入H(t) Σ p_i ln p_i 的离散版本即香农熵的负值进行计算。这里bin大小的选择需要谨慎太粗会丢失信息太细会因样本不足而产生噪声。6.4 结果分析与讨论图5展示了经过仿射变换后的模型输出 ˜h(t) 与真实H(t)的对比曲线。表III则给出了不同容量模型的拟合误差和仿射参数统计。模型容量影响从 nhidden16 增加到 32拟合的MSE显著下降从0.0429降至0.0155说明模型容量提升带来了更好的学习能力。但从32增加到64MSE略有上升到0.0171可能出现了轻微的过拟合或者需要更精细的超参数调优和更长的训练。参数一致性值得注意的是对于 nhidden64 的模型其仿射参数a和b在不同轨迹间的变异系数CV明显小于较小模型。这表明更大容量的模型学习到的内部表示在不同初始条件的轨迹间更具一致性它捕捉到的是更普适的、与具体轨迹细节无关的H函数本质这是一个非常积极的信号。符号与趋势所有模型的平均缩放因子 ⟨a⟩ 都是负数。这很好理解模型被训练为输出随时间增加的值因为损失函数鼓励 h(t) ≥ h(t1)而真实的H函数是随时间减少的。所以需要一个负的a来翻转趋势使其与H函数一致。实操心得与扩展思考数据标准化在将速度数据输入网络前进行适当的标准化如减去均值、除以标准差可以加速训练并提高稳定性。但要注意如果标准化是基于整个数据集计算的可能会泄露“全局时间信息”更好的做法是基于每个独立轨迹或每个批次进行标准化。超越硬球模型这项方法的魅力在于其通用性。我们可以尝试将其应用于更复杂的相互作用系统例如多原子分子气体粒子具有内部自由度如转动、振动碰撞动力学更复杂。活性物质如自驱动粒子鸟群、细菌其相互作用可能导致集体运动flocking其“熵”或有序度的定义本身就是一个有趣的问题。 在这些系统中我们可能没有像H定理这样清晰的理论目标函数但模型仍然可以学习到一个单调变化的“序参量”这可能帮助我们发现新的非平衡稳态或动态相变。解释性与可视化一个值得深入的方向是理解模型究竟学到了什么。我们可以通过可视化φ网络编码器将单个粒子速度映射到的特征空间或者分析ρ网络解码器的权重来尝试解释模型是如何从微观速度中聚合出宏观不可逆信息的。这有助于桥接机器学习模型与物理直觉。7. 常见问题与排查技巧实录在复现和实验过程中我遇到了一些典型问题以下是排查思路和解决方案的总结。问题现象可能原因排查与解决思路模型输出h(t)剧烈振荡不单调1. 使用了基础的对称损失函数 L^{(0)}。2. 学习率设置过高。3. 批处理大小batch size太小梯度估计噪声大。1.切换到Leaky ReLU损失函数并设置α 1如0.1引入非对称惩罚。这是最有效的解决方案。2. 逐步降低学习率尝试 1e-3, 5e-4, 1e-4。3. 增大batch size例如从32增加到128或256使梯度更新更稳定。损失函数很快降为零但模型输出为常数1. 正则化系数λ过大过度惩罚了输出变化。2. 模型容量nhidden过小表达能力不足。3. 数据本身H函数变化不明显如已接近平衡。1.减小λ值尝试0.1, 0.05, 0.01。λ的作用是防止平解但过大会扼杀任何变化。2. 增加隐藏层维度或层数提升模型容量。3. 检查模拟数据确保初始化在远离平衡的状态如均匀速度大小并观察真实H(t)曲线是否有明显下降段。训练后期损失突然上升或变得不稳定1. 动态调整的α(e)增长过快导致损失函数性质突变。2. 没有使用梯度裁剪gradient clipping梯度爆炸。1. 调整α(e)的调度策略如采用更平滑的S型sigmoid增长或更慢的线性增长减小斜率。2. 在优化器更新步骤前添加梯度裁剪torch.nn.utils.clip_grad_norm_将梯度范数限制在一个阈值内如1.0。不同训练轮次得到的模型差异很大1. 随机种子未固定。2. 数据划分或采样存在随机性且未控制。3. 模型初始化方差过大。1.固定所有随机种子Python, NumPy, PyTorch的随机种子。这是可复现性的基础。2. 确保数据加载的顺序是确定的或保存固定的训练/验证集索引。3. 使用PyTorch中更稳定的初始化方法如对线性层使用init.kaiming_normal_。计算真实H(t)时得到负值或异常值1. 直方图bin划分不当太多或太少。2. 速度数据未进行适当的范围处理导致有些bin为空。3. 在计算 p * log(p) 时对p0的情况未处理。1.自适应选择bin数量例如使用Freedman-Diaconis规则。可以尝试多个bin大小观察H(t)曲线是否平滑收敛。2. 确保直方图覆盖所有速度数据范围可以基于所有时间步的全局速度范围来定义bin边缘。3. 在计算熵时添加一个极小值ε如1e-10到概率p上避免log(0)的出现H -np.sum(p * np.log(p epsilon))。模型在验证集上过拟合1. 模型过于复杂层数多、神经元多。2. 训练数据量不足。3. 没有使用任何正则化技术。1. 首先尝试减小模型容量减少nhidden。2. 增加模拟的轨迹数量N_run。3. 在损失函数中除了λ项可以引入权重衰减L2正则化或在网络中添加Dropout层注意在ρ网络中添加比在φ网络更有效因为φ网络需要为每个粒子生成稳定的特征。一个关键的调试技巧始终同步可视化三样东西1) 模型输出的原始h(t)曲线2) 经过仿射拟合后的˜h(t)与真实H(t)的对比曲线3) 损失函数在训练过程中的变化曲线。这三者结合能迅速定位问题是出在模型架构、损失函数、还是数据本身。8. 总结与未来展望通过这个项目我们实践了一种“让机器从动力学数据中自行发现物理规律”的范式。我们并没有将玻尔兹曼方程或H定理的形式硬编码到模型中而是仅仅通过一个蕴含了“时间箭头”方向的损失函数以及一个尊重系统基本对称性置换不变性的神经网络架构就成功地让模型从硬球碰撞的混沌轨迹中抽取出与H函数行为高度一致的宏观量。这项工作的价值不仅在于它验证了H定理在简单系统中的有效性更在于它提供了一套可扩展的方法论框架。DeepSets处理集合数据的能力使其可以轻松应用于粒子数可变或关注整体性质而非个体身份的系统。孪生网络与对比学习的思路可以将“时间方向”这种抽象约束转化为具体的、可优化的机器学习任务。在我个人的复现和探索中最深的体会是损失函数设计是连接物理先验与机器学习模型的桥梁。Leaky ReLU损失中对“不增加”和“减少”的非对称处理完美地捕捉了H定理不等式dH/dt ≤ 0的精髓——它允许平衡态的平稳但禁止热力学第二定律的违反。这种对物理原理的细致翻译是成功的关键。未来一个令我兴奋的方向是将此方法应用于非平衡稳态的研究。例如在持续能量注入的活性物质系统中系统会达到一个非平衡稳态。我们能否训练模型学习一个类似H的函数该函数在弛豫到稳态的过程中单调变化并在稳态达到极值这样的函数可能成为描述这些复杂系统的新颖序参量。此外探索更复杂的网络架构如图神经网络它还能利用粒子的空间位置信息或将这种方法与符号回归结合以发现H函数可能的解析形式都是极具潜力的研究方向。这个项目就像是一把钥匙它打开了一扇门让我们看到机器学习如何以一种数据驱动、约束引导的方式帮助我们发现甚至重新定义那些支配复杂系统演化的基本规律。它需要的不仅是编程技巧更是对物理问题的深刻理解以及将这种理解转化为算法约束的创造力。