首页/文章/ 详情

PNAS:基于机器学习从湍流结构预测统计数据:以循环流动模式作为二维湍流的基础

4天前浏览23
 
0. 研究意义    
湍流领域的一个长期挑战是将个别的相干结构与流动的更为众所周知的统计特性联系起来。在本文中,我们通过将二维湍流表示为在流动中瞬时实现的精确不稳定周期轨道之间的马尔可夫链,建立了这一联系。为了找到动态相关的解,我们开发了一种基于标量损失函数优化的方法,这种方法克服了以往算法的限制,并在高雷诺数下具有良好的效果。随后,我们使用神经网络根据最近的不稳定解对湍流快照进行标记,从而实现马尔可夫表示,其不变测度可以重现混沌流动的完整概率密度函数(PDFs)。  
1. 摘要  

动力系统方法将湍流视为在高维状态空间中的轨迹运动[Hopf, Commun. Appl. Maths 1, 303 (1948)]。混沌动力学由填充惯性流形的不稳定简单不变解所塑造。研究的希望是将这种图景转化为一种预测框架,在其中流动的统计特性可以通过各个简单不变解的统计量的加权和来得到。然而,两个突出的障碍阻碍了这一目标的实现:1)已知解的匮乏,以及2)缺乏预测所需权重的合理理论。在本文中,我们描述了一种方法,可以实质性地解决这些问题,从而提供有力的证据,证明可以通过一组不稳定的周期轨道来重建完全发展的湍流流动的概率密度函数(PDFs)。我们寻找解的方法使用了自动微分技术,通过最小化与轨迹相关的损失函数构建高质量的初始猜测。我们利用这种方法在湍流的二维Kolmogorov流中找到了数百个解。随后,通过将湍流轨迹转化为马尔可夫链来学习权重,并利用深度卷积自编码器确定给定快照的最近解,从而计算出稳健的统计预测。在本研究中,我们成功地用一组简单的不变状态重现了时空混沌系统的PDFs,并且我们提供了一个自维持动力学过程与湍流更为人熟知的统计特性之间的有趣联系。

2. 绪论  

对于湍流的广为人接受的观点最初由Hopf倡导,是将湍流视为在非常高维状态空间中的轨道。湍流被认为是在这个轨道“弹跳”于不稳定简单不变解之间时的瞬时实现。这一观点对机械理解(可以在个别解的动力学中找到)以及量化个别动力学事件与长期统计之间的关系都很有吸引力。近几十年来,实现这一方法的尝试极大地改善了我们对过渡剪切流的理解:例如,管道中湍流的出现已经被与在某些临界雷诺数Re之后鞍结分岔中出现的有限振幅行波解的出现联系起来,而在所谓“最小”湍流配置中发现的不稳定周期轨道(UPOs)则揭示了在壁面湍流中起作用的一些自维持机制。然而,由于计算方法的局限性,无论是在寻找UPOs方面还是在使用它们标记完全发展的湍流实现方面,这些想法在高Re流动中的扩展都受到了限制。

对于完全发展的湍流,动力系统的观点设想状态空间中布满了简单的不变解(平衡态、行波、不稳定周期轨道),其复杂交织的稳定和不稳定流形为湍流的“弹球”运动构建了支架。这一图景提出了一种基于适当加权的不变解特性的预测理论,其中权重反映了在该解附近停留的相对时间。然而,找到足够多的重要解,更不用说识别出适当的权重,已经被证明是极其昂贵的。特别是在高雷诺数(Re)下,不变解的集 合急剧增加,并且每个解变得越来越不稳定,阻碍了它们的识别。

一个核心且广为人知的问题是,牛顿-拉夫森根寻找算法在高维问题中对初始猜测质量的敏感性。过去的工作主要依赖于“循环流动分析”来生成足够好的不稳定周期轨道(UPO)猜测,其中要求湍流轨道至少在一个完整周期内与UPO重合。在实际操作中,这对UPO的不稳定性设定了一个限制,随着雷诺数(Re)的增加,极大地阻碍了搜索。第二个问题是识别流动何时几乎重现的策略。通常这只是通过初始和最终状态之间差异的欧几里得范数来完成,因此,近似重现的阈值必须设定得相当高。最后,第三个缺陷是许多动态相关的周期轨道尚未被发现,特别是那些耗散率高于平均水平的轨道。这些困难促使了一些替代方法的出现,用于生成猜测,甚至是UPO的识别,但即使在最简单的稳态湍流模型问题中,这些方法也未能充分证明动力学与统计之间的联系。

即使拥有完整的UPO集 合,一个深刻的理论问题是如何预测每个UPO在湍流流动的“UPO展开”中应被计入的权重。周期轨道理论提供了一种在低维动力系统中有效的理论方案。然而,将这一理论应用于Navier–Stokes方程仍然具有挑战性。早期将该理论应用于二维Navier–Stokes方程的尝试显示,其效果并不优于使用等权重的对照实验,尽管可用的解集显然太小。在这一点上,甚至确定是否有可能用UPOs来展开任意湍流流动,仍然是流体力学中的一个未解问题。

在本文中,我们提出了解决这一问题的计算方法,这些方法克服了许多早期的限制,介绍了用于寻找和收敛不稳定周期轨道(UPOs)以及通过标记湍流数据来定义权重的方法,该标记基于在状态空间中哪个解最接近。与早期工作相比,我们的UPO检测方法不需要精心构建初始猜测,并且能够产生大量动态相关的UPOs。为此,我们改进了最近开发的完全可微分流动求解器,通过对涉及整个解轨迹的损失函数执行梯度下降,我们能够找到UPO的高质量猜测。这使得我们能够显式搜索具有某些特性(例如高耗散率)的周期轨道,并且能够成功地从任意的湍流快照开始收敛大量的UPOs。随后,我们训练了深度神经网络,以学习湍流的精确低阶表示,我们可以利用这些表示来测量在任何时刻湍流轨道最接近哪个UPO。结果是一个马尔可夫湍流动力学,它不仅允许我们通过链的不变测度定义UPOs的权重,还提供了对极端事件路径的洞察。从不变测度中发现的权重可以稳健地再现湍流吸引子的统计特性,包括完整的耗散概率密度函数(PDF),实现了Hopf的原始设想。虽然我们在模型问题中开发了这些方法,但其底层方法论有望改变我们对更高雷诺数(Re)下典型流动的理解。

3. 周期轨道搜索策略  

3.1 二维湍流

我们在广泛研究的湍流“Kolmogorov”流中展示了我们的不稳定周期轨道(UPO)搜索方法:在这里,我们研究的是在周期性域中的二维单色强迫湍流。其控制方程为:

该问题通过垂直方向上的基波数进行无量纲化,我们使用常用的强迫波长n=4。我们还经常使用面外涡量,无论是在收敛精确解的公式中还是在用于标记发现解的神经网络训练中。我们考虑了两个雷诺数Re = 40和Re = 100,在这些情况下观察到自维持湍流,其中Re = 100的情况明显处于渐近区(16)。在Re = 40时,之前已发现大约50个UPOs,且它们的平均耗散率都较低,而在Re = 100时仅收敛了9个UPOs。

控制方程[1]在水平方向的连续平移下具有等变性。

因此,许多简单的不变解是相对平衡态(行波)或周期轨道。同时,还存在离散平移反射对称性和旋转对称性。在此,我们寻找具有某一周期和位移的相对不稳定周期轨道(RPO),其满足以下条件:


3.2 自动微分用于周期轨道

我们使用JAX-CFD求解方程1(以及材料与方法中描述的等效速度-涡量形式),JAX-CFD是一个完全可微的流体求解器,可以通过自动微分计算时间前向映射,相对于初始条件的梯度,并达到机器精度。这种能力构成了周期轨道搜索策略的基础。我们使用了JAX-CFD的“标准”有限差分原始变量公式和光谱涡量版本。前者用于构建稳健的周期轨道猜测,原因将在下文讨论,而后者用于在牛顿求解器中的最终收敛以及与先前报道的结果进行比较(所有这些结果都是通过光谱代码获得的)。

与早期通过在时间序列中识别“近似重现”并将其直接输入到根求解器中来寻找UPO的尝试相比,我们的方法是通过标量损失函数的基于梯度的优化进行搜索——无需任何显式的初始条件选择。这个损失函数只是方程2的一个缩放范数:

并且它依赖于初始条件u、未知的位移,α和周期T。使用JAX库及其扩展可以高效地计算所有这些变量的梯度。通常情况下,我们认为那些可以将损失降低到L≤0.015的猜测适合传递给牛顿求解器进行进一步处理——直接使用优化器收敛的速度太慢,因此它被用作牛顿方法猜测的有效预处理器。

在Re = 40时,我们还将明确针对某些周期T*,并尝试找到平均耗散率高于某些阈值D*的周期轨道(此类UPO在先前的结果中大多缺失)。我们通过在损失函数中添加适当的项来实现这一目标:

JAX-CFD 的原始变量公式允许存在恒定的背景垂直速度。上文描述的基本“Kolmogorov”流动中,但原则上该问题允许任意一均匀的跨流速度,一旦设定后,其在时间上保持恒定。添加这一效应将改变周期轨道——其中许多轨道可能通过同伦方法回到v0=0——并且还可能在鞍节点分岔中引入与原始配置无关的新解。在此,我们发现添加这一效应总体上是积极的,因为它可以防止优化器陷入浅局部极小值——这在直接使用JAX-CFD的光谱涡量公式进行搜索时是一个常见问题。最近在一种变分公式中尝试寻找周期轨道时也有类似的观察,其中使用了非散度速度场作为初始猜测。在基于梯度的优化后,使用JAX-CFD的光谱版本进行UPO的牛顿收敛时,不可能存在恒定的背景流,因为我们通过Δψ=-ω解诱导速度,假设流函数ψ是周期性的,由此得到u和v。

对于短周期(T≈10,Re=40),弱背景垂直速度对UPO的影响较小,光谱牛顿求解能够在几步内收敛到附近的v0 = 0解。对于较长周期,弱垂直流的影响更大,但通常可以通过附加的优化过程有效地消除,这一过程会对垂直速度进行惩罚:

在这种情况下,我们设置μ = 103并使用较小的学习率(通常为η = 10-2)来小心地将近闭合回路变形为垂直速度接近于零的状态。

4. 不稳定周期轨道    

4.1状态密度在Re = 40时的表现

我们首先在研究较多的Kolmogorov流动问题中展示我们方法的威力,雷诺数为Re = 40。尽管许多前人已经研究过这一配置,但我们仍然对其状态密度随周期ρ(T)的变化知之甚少,也未能识别出任何局部高耗散的不稳定周期轨道(UPOs)。受此启发,我们使用损失函数 [4a] 对特定周期的解进行了全面搜索。鉴于周期轨道理论中“素周期”(prime cycles)的重要性,我们首先聚焦于短周期。我们在该范围内以 0.5 为步长递增目标,并在每个步长处进行 50 次优化计算。每次计算均由湍流吸引子的随机快照初始化。此前已知在此范围内存在四个解,周期为T∈ {2.83, 2.92, 5.38, 6.72}。

我们还使用损失函数 [4b] 初始化了对高耗散解的独立搜索。我们进行了三次计算,搜索平均耗散率高于阈值的解。在变量v0公式中,大量解被收敛,尽管在目标T*上成功率不一。找到这些常见解是预期之内的,但我们也发现了许多之前未知的大量解。事实上,对于T < 8的情况,我们收敛了38个独特的解,其中包括许多以前搜索方法无法访问的高耗散状态。

值得注意的是,这些短周期解几乎涵盖了整体流动中产生和耗散事件的全部范围(见图1右侧面板)。实际缺失的是低耗散事件——这些事件与较慢的动力学相关,UPO通常具有较长的周期。为找到这些状态,我们还搜索了目标周期为T*∈{12.5, 15, 17.5, 20}的解。通常情况下,优化器在面对这些较长周期的轨道时表现困难,因为从随机初始条件开始,fT(u)很可能在解流形上离得很远,梯度信息并不特别有用。然而,我们确实获得了几个较长周期的解,它们全都是低耗散的。这部分工作将从更仔细的初始条件选择方法中获益,以确保fT(u)在解流形上“接近”,而不必依赖近似重现。

图1. (左图)所有在 Re=40时周期 T<8 的收敛轨道的耗散与周期的关系(平均耗散率以白色方框表示)。参考文献16中列出的四个已知解以蓝色标识,平均耗散率以圆圈表示,周期为 T∈{2.83,2.92,5.38,6.72}(周期为 T=7.16 的解也可能在参考文献27中被发现)。用十字标识的解稍后将在图2中展示。(右图)在 Re=40 时所有收敛解的能量产生与耗散的关系,包括左图未显示的一些较长周期的轨道(所有收敛的UPOs在附录中列出)。湍流状态的概率密度(从长时间计算运行 tlong=2.5×105 得出)以灰色显示。等高线级数以对数间隔,最小值为 106。所有数值均以层流值 Dl=Re/(2n2) 归一化。

我们在图2中检查了一些𝑅𝑒=40的不稳定周期轨道(UPOs),并在轨道上的四个点报告了面外涡量的快照。低耗散解(图中显示了𝑇=2.83的UPO)都具有一个共同的结构,即涡旋结构位于一对倾斜的涡量条纹之上(这让人联想到该配置中的第一个非平凡平衡解)。在𝑇=2.83的情况下,周期对应于右侧面板中一个椭圆形涡旋块的完整旋转;周期为𝑇=5.38的轨道类似,但涉及的是同号涡旋对的旋转。

图2. 面外涡量在𝑅𝑒=40下的三个不稳定周期轨道(UPOs)上以时间等间隔提取的四个点处的快照。(顶部)在𝑅𝑒=40时最短的轨道——已知的𝑇=2.83解,平均耗散率⟨𝐷/𝐷𝑙⟩=0.095。(中间)一个新的高耗散UPO,周期为𝑇=2.90,平均耗散率⟨𝐷/𝐷𝑙⟩=0.246。(底部)一个新的高耗散UPO,周期为𝑇=3.27,平均耗散率⟨𝐷/𝐷𝑙⟩=0.285。在所有情况下,等高线的范围为−10到10。这些解都在图1中用十字标识。

相比之下,高耗散结构表现出更加多样化的动力学。例如,图2中报告的两个“爆发”UPO之一具有晶体状结构,其周期内保持着四个大振幅的涡核。图2中的另一个高耗散解则表现为一个强大的偶极子结构,从左到右穿过域移动。不出所料,许多新解具有比先前记录的解更高维度的不稳定流形(详见附录中的完整细节)。


4.2 雷诺数Re=100时的周期轨道

我们现在将注意力转向雷诺数Re=100下的强湍流情况,此前通过任何方法仅获得了少数几个解。在该Re值下先前获得的九个解均为短周期,T<5。我们在此报告的大多数解都是使用基本损失函数[3]计算的,未添加其他物理条件(例如,方程4a和4b中包含的内容)。

我们初始化了大量计算,起始周期从T0∈{1,2,3,4,5}中选择,另外也对T0=2.5进行了计算,这是受T0=2和T0=3计算成功率的启发。当一批大约100次的猜测未能产生新的、未包含在我们已收敛的UPO集 合中的解时,我们停止了这些计算。我们还使用损失函数[4b]进行了高耗散搜索,寻找平均耗散率⟨D/Dl⟩≥0.03的周期轨道,但发现这远不如在Re=40时的等效计算有效。相比之下,标准损失函数[3]无需扩展便返回了广泛的周期轨道,而不是相同的一组解。总体而言,我们观察到解收敛的成功率在5%到15%之间,取决于起始周期的选择。需要强调的是,过程中的起点是来自湍流计算的随机快照,这与其他方法形成对比,后者的成功率接近于零,即使猜测是经过更仔细构建的。

到目前为止,我们的计算在Re=100时共产生了151个独特的UPO,这些状态总结在图3中(完整列表见附录)。在这里,我们看到大多数UPO在相空间中表现出高度的局部化,几乎都以二维投影中的小闭合环出现。当从动能E的角度可视化时,这些状态表现得更加局部化(见附录,图S4),并且为了全面覆盖相关的湍流概率密度函数(PDFs),需要比在Re=40时的等效计算更多的状态。

图3. 在Re=100时能量产生率与耗散率的关系。灰色背景是从长时间湍流计算tlong=2.5×105得出的概率密度函数(PDF)。等高线级数以对数间隔,最小值为 10−6。闭合环是151个收敛UPO的二维投影。所有数值均以层流值 Dl=Re/(2n2)归一化。所有周期轨道及相关属性(周期、位移、Floquet指数)均列在附录中。

我们在图4中展示了一些收敛的UPO。解表现出多种不同的动力学行为。通常,这些状态主要由两个大涡块主导[由于逆向级联所致]——见图4中间两排。在这种类型的解中存在很大差异,除了图中所示的行为(例如,一个强大的静止涡旋和一个共转涡旋对),还存在具有三个或更多同号共转涡旋的状态——详见附录。更强耗散的状态(图4的顶部和底部)显示出更加多样的涡量动力学,需要进一步的工作来评估我们方法所发现的各种UPO“类别”。

图4. 面外涡量在 Re=100下的四个UPOs上以时间等间隔提取的四个点处的快照。从上到下,UPOs的周期和平均耗散率分别为:(T,⟨D/Dl⟩)=(1.424,0.057)、(1.794,0.053)、(4.212,0.027)和(1.164,0.078)。涡量等高线级数在 ±10之间。

5. 从周期轨道影子中构建马尔可夫链    

5.1 识别UPO影子

鉴于我们找到的众多收敛UPO及其在 Re=40 时对I–D平面的特别良好覆盖,我们现在尝试通过标记湍流时间序列快照,找出哪个UPO“最接近”,以验证并可视化Hopf的原始猜想。我们的目标是使用我们收集的UPO来划分状态空间,并将湍流轨道转化为离散时间马尔可夫过程的实现。

为了准确测量到最近UPO的距离,我们在“DenseNet”配置中训练了高度精确的深度卷积自动编码器(有关完整的架构和训练细节,见材料与方法),并将使用基于这些网络中的潜在表示的可观测量,而不是物理空间中快照之间的距离。这些模型的准确性已在我们最近的工作中展示过,覆盖了广泛的Re范围——即使在罕见的最高耗散事件中,准确性也得到了保持。

这些自动编码器由编码器 E:RNx×Ny→Rm(其中 m=128在 Re=40 时,m=512在 Re=100时)和解码器 D:Rm→RNx×Ny组成,使得[DE](ω)≈ω。给定一个编码后的快照E(ω),我们首先通过将E投影到所谓的“潜在傅里叶模态”上来构建一个流向平移不变的可观测量,这些模态是离散平移算子 TαE(ω):=E(Tαω)(对于某些固定的流向平移α)的特征向量。然后,我们使用这些投影构建一个向量可观测量ψ(ω),其具有性质ψ(Tsω)=ψ(ω)对于所有s∈R(详细信息见材料与方法)。最后,我们计算每个周期轨道以及它的15个离散对称 品的周期平均值 {ψ(SmRqftj))T:0≤m≤7,q∈{0,1}}。

最近的周期轨道与快照 ω\omegaω 的对应关系根据以下方式确定:

其中1≤j≤Np(Np是在给定Re下我们库中周期轨道的总数),并且我们在每个时间点上搜索离散对称性。由于我们的大多数UPOs都是短周期且在相空间中局部化的,因此与UPO嵌入的时间平均值进行比较是稳健的,尽管更复杂的方法也可以在时间方向上进行搜索。

我们在Re=40时构建了长度为T=2.5×104的长轨迹,在Re=100时构建了长度为T=104的长轨迹,其中快照在前者中以δt=1的间隔分隔,在后者中以δt=0.25的间隔分隔。每个案例中的快照间隔都是基于我们库中UPO的典型周期(例如,在Re=40时常见的周期为T5,而在Re=100时,许多解的周期为1≤T≤2)以及观察周期解“影子”动力学发生的动机。湍流时间序列被转换为形式为 POi→POi→POi→POk→POj→POj→的标签序列。


5.2 离散时间马尔可夫链和统计预测

图5展示了上述UPO标记协议的一个示例,我们在此为了说明目的使用了更精细的δt。图中的示例显示了一个扩展的高耗散爆发事件,该事件在大约t30时恢复到较为平静的低耗散动力学。曲线根据方程6确定的最接近的UPO进行着色,这表明流动在长时间内保持接近某个特定的UPO,并且这一特定序列可以通过少量的精确解很好地描述。事实上,这个示例轨迹超过一半的时间都在接近仅四个UPO的区域内,高耗散事件多次采样相同的解(图5中的深蓝色曲线,另见图例)。我们将该示例轨迹沿途的快照与图右侧面板中最接近的UPO的快照进行比较,可以看到湍流中的许多相似的定性流动特征在周期解中得到了再现。同时也清楚地表明,匹配还可以进一步改善——更稳健的标记协议也可以在所有UPO的时间方向上进行搜索,尽管这样做会显著增加计算开销。

图5. (A) 在 Re=40 时耗散率的示例时间演化,颜色根据方程6确定的最近UPO进行标记。在此区间内最频繁访问的五个UPO在图例中突出显示,并以其周期 TTT 进行标注。(B) 涡量的快照(顶部)来自用于生成(A)中耗散图的直接数值模拟 (DNS),快照对应的时间在(A)中用方框/数字标识,同时也展示了最接近的UPO的快照(底部),其中我们选择了水平位移s和轨道上的时间τ,以最小||ω−Tsfτ(ωPO||2。

我们现在使用上一节末尾描述的长时间序列来构建转移概率矩阵P,其元素为 Pij:=P(POi→POj)。通过简单地计数状态之间的转移,然后对行进行归一化∑jPij=1i,可以轻松计算出这一量。图6A展示了Re=40时的转移矩阵,其中状态按耗散排序(低耗散在顶部,高耗散在底部)。我们还展示了通过πTP=πT获得的不变测度。

图6. (A)在Re=40时的不变测度π和转移矩阵P(显示为转移概率的对数,快照间隔为 δt=1)。状态按平均耗散率从低到高排序(最低在顶部/最左侧)。红色水平线标识了文中讨论的网 关或“混合”状态。红色垂直线强调了具有有限转移概率进入该混合区域的广泛状态范围,涵盖了高耗散和低耗散事件,而虚线垂直线表示文中讨论的 D/Dl=0.15的“阈值” 。(B) 展开式 [7] 中的权重(也是马尔可夫链的不变测度wj=πj)与每个UPO的增长Floquet指数之和(取其实部)∑jσj,σj>0的关系图。

转移矩阵中存在许多有趣的特征,值得进一步讨论。转移矩阵通常在其对角线上显示出最大的概率,这意味着对于大多数状态而言,最有可能的结果是在该特定UPO附近停留另一个时间片刻(此处为δt=1在Re=40时)。这与湍流轨道“影子”单个循环解的行为一致。非对角线上的非零概率还表明,转移倾向于发生在具有相似耗散率的状态之间。在Re=40时,混沌动力学可以分为低耗散的“静止”状态和较为罕见的高耗散爆发事件(大致归一化耗散 D/Dl0.15),这一划分在图6的转移矩阵中很明显,其中图中的状态按平均耗散排序,值D/Dl=0.15用黑色虚线标出;在高耗散状态之间有多条路径,然而这些爆发事件的转移似乎通过少数几个低耗散的网 关UPO(大约4个)发生,这些UPO也能将系统抛回到非常低的耗散状态(见图6中的红色线条)。

湍流的马尔可夫视角可以扩展为以周期轨道理论为基础进行统计预测。为此,我们寻求一组固定的权重{wj}Npj=1,其中Np是找到的UPO的总数,使得任何可观测量γ:M→Rn在状态空间 M上的空间平均可以构建为UPO统计量的线性叠加。

其中Γ是正在考虑的平均值,而Γj是针对周期轨道j计算的相同统计量。根据我们的转移矩阵,我们可以简单地将权重定义为wj≡πj,其中∑jwj=1(根据定义)。这些权重在图6的下方面板B中作为基础UPOs不稳定性的函数进行了分析(通过Floquet指数的增长率之和来衡量)。权重与不稳定性水平呈负相关,其中高度不稳定的状态在重构中显得不太重要。这反映了这样一个事实:高度耗散的UPOs(它们对分布中非常弱的右尾贡献)往往更加不稳定。

为了展示上述UPO展开式[7]的表现,图7A报告了在Re=40时的长时间范围内(tlong=2.5×105)计算的耗散率、能量产生率和动能的PDF,并叠加了通过表达式(如方程7)计算的基于UPO的PDF。无论是耗散还是产生的再现都非常完整——尽管我们缺少一些与较长轨道相关的非常低的值(如前文所述)。在这两种情况下,最值得注意的是高 DDD/高 III 尾部的精确重构,这在所有先前的尝试中都未能实现。动能E也被重现得相当高标准[特别是与以往尝试应用周期轨道理论的结果相比]。缺失的状态在E/El0.35和E/El0.45处,可能与缺少最低耗散的轨道有关。在附录中对此进行了进一步探讨,我们还展示了解在能量-耗散平面上的分布。

图7. 基于UPO的统计预测,这些预测是根据马尔可夫链的不变测度计算的,定义了UPO展开中的权重(方程7),适用于(A) Re=40 和(B) Re=100(转移矩阵在附录中报告)。从左到右的前三个面板显示了耗散率、产生率和能量的概率密度函数(PDF):虚线蓝色曲线是通过长时间湍流计算(tlong=2.5×105)获得的“真实”PDF,而填充的灰色曲线是基于UPO的重构。最后两个面板比较了平均速度剖面和均方根速度波动(左侧为u,右侧为v))——在流向、离散对称性和时间上取平均值。蓝色和橙色曲线分别表示u和v的DNS真实值,黑色曲线表示基于UPO的重构。

图7A还测试了UPO预测的平均速度剖面和均方根速度波动(同样使用相同的权重),这些结果在离散对称性、流向坐标和时间上取平均值。所有三个预测值都接近真实的时间平均值,这再次表明UPO基于预测的表现相较于以往的最先进方法有了显著提升。

在图7B中,我们使用在 Re=100 时收集的大量151个UPO进行了类似的分析。对应的转移矩阵包含在附录中,表明存在影子效应,但在平均耗散上可能的转移大范围分布——尽管显然这一图景尚不完整,因为我们缺少许多重要状态,这一点从之前的结果中可以看出(例如图3)。尽管如此,图7B中报告的统计数据是令人鼓舞的,代表了即使在更低 Re 时也无法实现的显著进步。特别是,平均剖面的预测几乎完美,而均方根速度的误差仅为 O(10%)。PDF显示,缺失的状态与更大的耗散和能量值相关,这里概述的UPO搜索过程为通过进一步计算填补这些空白提供了明确的策略。

6. 讨论    

在本文中,我们汇集了令人信服的证据,支持湍流是一种在不稳定简单不变集 合之间弹跳的高维“弹球”的观点(这一观点通常归功于Hopf)。为此,我们设计了一个方法来找到大量动态相关的不稳定周期轨道(UPO)——这是该领域长期存在的局限性——并提出了一种方法,通过最近的UPO对湍流快照进行准确标记,其中距离在深度卷积自动编码器的潜在空间中测量。结果是湍流的马尔可夫图景,这在适中的雷诺数下不仅提供了对动态路径和极端事件路径的洞察,还提供了对混沌动力学的稳健统计预测。

UPO搜索策略被构建为一个基于梯度的优化问题,并在一个完全可微的流动求解器中实现。基于损失的这种方法允许针对具有特定特征的解(例如,高耗散、高能量)进行定向搜索,在适中的雷诺数 Re=40 下应用此方法揭示了大量先前无法检测到的短周期解——包括高耗散和低耗散。该方法在更高的雷诺数 Re=100 下仍然有效,我们再次发现了大量短周期解。在高雷诺数下的状态在状态空间中表现出高度局部化,并显示出丰富有趣的涡旋动力学。

为了根据最近的UPO对涡量快照进行标记,我们训练了高度精确的深度卷积自动编码器,并使用这些网络的低维潜在空间中的可观测量来衡量相似性。然后,我们能够将长时间的湍流时间序列视为马尔可夫链,其中每个UPO都是一个独立的状态,构建转移矩阵,并使用其不变测度进行统计预测。在 Re=40 下,该方法特别有效,因为新的UPO库几乎涵盖了完全湍流状态下看到的产生和耗散事件的全范围,结果是统计预测非常稳健。即使在 Re=100 下状态集不完整,速度矩统计的预测也相当稳健,并且相较于以前的方法,即使在更低的雷诺数下也是一个显著的改进。

尽管收敛了大量新解,仍然有一些重要的状态在 Re=40 时缺失,而在 Re=100 时则缺失了大量解。这些空白中的一些(例如 Re=40 时的低耗散事件)可能会从改进的快照选择过程受益。在当前配置中,改进可能涉及诸如仅选择返回到状态空间相似区域的快照——不是近似重现,而是更弱的要求,例如通过这里训练的自动编码器确定——或者通过搜索离散对称性。这在研究更复杂的三维流动时可能是一个重要的考虑因素。

尽管如此,在 Re=100 时状态空间覆盖中的大范围空白,结合组装解时不断增加的计算成本,自然会引起一些关于“UPO-弹球”概念适用性的担忧,以分解二维和三维多尺度湍流的高雷诺数流动的统计数据。例如,在 Re=100 时使用 T≈2.5 的新解收敛,网格点数为 Nx×Ny=512×512 时,平均需要约24小时的GPU计算时间(使用Nvidia V100显卡,此计算假设从随机湍流快照开始的成功率为5%)。此外,还有用于标记快照的自动编码器训练成本(在约105个快照的训练数据集上需要约48小时的GPU时间)以及构建稳健的马尔可夫链所需的长参考轨迹成本——这本身就产生了准确的统计数据,UPOs用于重构这些数据。

即使完全的统计覆盖具有挑战性,寻找具有特定物理特性的精确循环流动的能力也非常有用,如果目标是理解个别动态事件在诸如能量级联等众所周知的统计现象中的重要性。同样值得注意的是,三维情况可能与本文研究的二维流动有很大不同:在二维情况下,我们预计在高雷诺数下获得的一些解可能会在 Re→∞ 时连接到欧拉方程的解,这些解预计会形成充满域的涡旋,而有壁湍流预计需要更复杂的多尺度描述,具有附着于边界的不同大小的涡旋。早期对三维最小流动单元中弱瞬态湍流的研究表明,个别UPOs的统计数据有时可能与湍流本身非常相似,这为三维情况提供了谨慎的乐观态度。

本文概述的UPO搜索方法使我们能够明确针对我们迄今为止构建的PDF中的空白进行搜索,我们相信这一能力意味着尽管目前所需的计算时间似乎很极端,但仍有可能通过UPOs实现稳健的统计覆盖。方法的许多方面都有改进的空间,这些改进可能会加速搜索过程,包括i) 更加考虑的初始快照选择和ii) 修改优化(例如通过改变范数)以提高收敛率。这种表示的潜力远远超出简单的统计重构:UPO表示允许计算混沌系统的敏感性,而通常统计量的梯度是不可计算的,而最大的希望是在一个雷诺数下的构建可以用来评估个别动态事件在更高雷诺数下的统计中的作用,通过延续解及其权重。

7. 材料与方法    

7.1 模拟

使用JAX-CFD的有限差分版本进行的模拟在大小为 Nx×Ny=256×256 的网格上进行,雷诺数为 Re=40,以及在大小为 Nx×Ny=512×512 的网格上进行,雷诺数为 Re=100。对于算法中的牛顿求解组件,我们使用了JAX-CFD的光谱版本,并在光谱模拟中将分辨率与前期有限差分优化中的分辨率匹配。时间步长由基于速度估计的CFL条件决定,该速度估计是层流基本剖面的两倍,即 2×Re/(2n2),这一速度通常比湍流状态下观察到的速度大得多。

代码的光谱版本在涡量-速度形式下求解Navier–Stokes方程,

其中 ω:=(×u)z(与方程1相比)。与有限差分的原始变量公式不同,没有可能的背景恒定流动。每个时间步的速度场都是由涡量引起的,并通过求解泊松方程 Δψ=ω 找到,其中流函数 ψ 通过以下关系导出诱导的速度分量:u=yψ,v=−∂xψ。

8. 神经网络和距离度量    

8.1 架构

卷积自动编码器是编码器 Em 和解码器 Dm 的组合,使得

通过编码器Em:RNx×Ny→Rm进行降维,其中我们在Re=40的数据中将m固定为128,在Re=100时固定为512。

编码器是一个全卷积架构,降维通过在卷积“密集块”后的最大池化来完成。每个密集块由三个连续的卷积层组成,每个层接收前一卷积层的输出并与同一密集块中的所有上游层的输出连接。每个密集块内的卷积都会创建额外的32个特征图。每个密集块后,我们应用最大池化,接着是一个卷积层,将特征图的数量减少到32。整个编码器由六个密集块/最大池化组合组成。在最内部的表示中,编码器生成一个形状为4×8×M的图像,其中 M=m/32(m是最内部潜在表示的指定维度)。Kolmogorov流中的离散平移反射对称性确定了垂直方向的值“8”。

在整个网络中,我们使用“GELU”激活函数,除了解码器输出使用tanh(输入数据被归一化ω→ω/ωnorm使得max|ω(x,y)|≤1)。解码器模块的结构与编码器相似,但顺序相反,使用上采样代替最大池化。


8.2 训练

我们使用修改后的“均方误差”作为损失函数:

其中额外的项(对涡量平方的均方误差)旨在鼓励网络学习高效地表示更为罕见的高耗散事件(而不是改变基础数据的分布),因为该项在整体损失中占据越来越重要的地位。

我们在每个 Re 下使用 N=105 样本进行训练,这些样本由1000个独立轨迹组成,快照间隔为一个对流时间单位。我们应用了数据增强:在x和y方向上随机平移,并应用旋转对称性。所有模型都使用Adam优化器,学习率为η=5×104,并训练500个周期,批量大小为64。


8.3 潜在傅里叶分析

潜在傅里叶分析是一种自动编码器可解释性技术,它利用了控制方程/边界条件中的连续对称性。在潜在傅里叶分析中,我们寻求一个操作符来对涡量场的嵌入进行连续平移,即

在实际操作中,我们使用 “动态模态分解”算法来确定一个近似Tα,它为测试数据集的嵌入及其α平移对应项之间的映射提供了最佳(最小二乘意义上)的操作符。经验上,我们观察到随着α的逐渐减小,只有少量非零的潜在波数会达到饱和(在Re=40时lmax=3,在Re=100时lmax=7),这表明网络嵌入的模式具有由l设置的基本周期性。每个潜在波数可能是高度简并的。

然后,我们可以写出一个快照嵌入的表示,受到任意平移 s∈R 的影响:

这与标准的傅里叶变换有着明确的联系。这里,Pl是与波数l相关的简并特征空间上的投影算符:

每个d(l)简并方向上的投影算符定义为,其中ξk(l)是子空间l中的特征向量,匕首符号表示伴随特征向量。

为了定义我们的平移不变可观测量ψ(ω),我们对与l≤3相关的单个特征空间上的投影进行SVD(大多数能量都包含在这些模式中),通过特征值|λ|>0.9的数量确定简并度d(l)。然后使用每个子空间中SVD的左奇异向量确定平移不变的可观测量:

通过将对l>0的投影取绝对值,确保了ψ(ω)≈ψ(Tsω)s∈R。

与在物理空间中测量涡量场之间的距离相比,使用ψ作为可观测量显著提高了UPO检测中的循环流动分析的效果。



翻译转载自:《PNAS》: "Recurrent flow patterns as a basis for two-dimensional turbulence: Predicting statistics from structures"  


 


来源:多相流在线
湍流理论材料多尺度控制管道
著作权归作者所有,欢迎分享,未经许可,不得转载
首次发布时间:2024-10-26
最近编辑:4天前
积鼎科技
联系我们13162025768
获赞 108粉丝 104文章 295课程 0
点赞
收藏
作者推荐

Comput Geotech: 用ChatGPT写的MATLAB代码能做CFD?

摘要:ChatGPT作为大语言模型(LLMs)的代表,最近在我们的社会中引发了革命性的变化,其在各种应用中的有效性也得到了越来越多的报道。本研究旨在探讨在岩土工程领域中,利用ChatGPT对对话提示的响应所驱动的编程性能的潜力。测试的示例包括渗流分析、边坡稳定性分析以及部分饱和砂的X射线计算机断层图像处理。在每个案例中,首先通过叙述性解释提供问题的属性,如几何形状、初始条件和边界条件,以生成MATLAB代码,并执行以评估其正确性和功能性。对于出现的任何错误和意外结果,通过附加提示进行进一步优化,直到获得正确的结果。当提供的提示基于对给定问题的全面理解时,ChatGPT能够生成具有相当水平的数值代码,并表现出在优化过程中值得信赖的意识。尽管ChatGPT可能无法替代整个编程过程,但它可以帮助减少粗心的语法错误,并协助设计逻辑编程的基本框架。SIMPOP1. 引言 ChatGPT(由OpenAI开发)是基于Transformer架构的最先进的大语言模型(LLM),用于生成类似人类的文本响应。ChatGPT能够接收对话提示,进行处理,并生成相关的文本响应。虽然还有其他的大语言模型,如BERT和T5(Devlin等人,2018年;Raffel等人,2020年),但ChatGPT通过其对上下文的细腻理解、对各种查询的适应性以及生成连贯响应的卓越能力脱颖而出。正因为如此,ChatGPT开创了人工智能领域的新篇章,并成为现代大语言模型(LLMs)所能实现的目标和对社会带来影响的代表。自2019年11月发布完整的GPT-2模型以来,关于ChatGPT在各个领域的应用已有广泛的报道和讨论,这些内容在表1中进行了简要总结。随着对ChatGPT的优缺点、潜在风险以及伦理问题的讨论,ChatGPT的应用领域似乎在迅速扩展,尽管其生成响应的可信度和真实性尚未被正式报道和公布。表1. ChatGPT的已报道应用人们认为,现在是时候设想、讨论并分享如何在科学和工程领域采用ChatGPT的强大功能,尤其是在这些领域中已有的文献参考数量有限的情况下。Bahrini等人(2023年)通过让ChatGPT生成系统与信息工程专业的初级考试问卷,评估了其在教育中的有效性。Tsai等人(2023年)提出了利用ChatGPT辅助的解决问题程序,用于开发化学工程教育课程,尤其是在建立蒸汽动力装置、相反应和传质模型的情况下。Kashefi和Mukerji(2023年)探讨了ChatGPT在生成和改进不同编程语言的数值算法代码方面的能力,应用于诸如泊松方程、扩散方程、不可压缩Navier-Stokes方程和可压缩无粘流动等各种数学问题。在建筑领域,ChatGPT不仅能够生成满足项目要求的施工进度表(Prieto等人,2023年),还能够促进危险识别和安全培训(Uddin等人,2023年)。对先前研究的回顾得出的结论是,进一步探索ChatGPT在数值建模中的潜在应用是值得的。因此,本研究旨在调查ChatGPT在编写编程代码方面的能力,这是ChatGPT提供的最可靠功能之一,用于解决诸如渗流和边坡稳定性等著名且基础的岩土工程问题。同时,还尝试了对部分饱和砂的X射线CT图像进行图像处理。对于每个示例,首先定义了给定问题的概述,并通过逻辑提示向ChatGPT请求生成MATLAB编程代码。通过改进提示反复检查结果,直到获得正确的结果。通过这些用于自动化编写编程代码的指定示例,展示了将ChatGPT集成到解决传统岩土工程问题中的潜在应用。2. 问题描述与解决过程 选择了岩土工程中的三个基础问题,以探讨ChatGPT-4.0在通过提示生成MATLAB代码时的性能。以下是问题定义及后续优化过程的描述。2.1. 二维渗流分析 图1a中渗流的配置显示,土层高度为10米,宽度为25米,一道宽度为1米的板桩插入土层中,深度为5米。初始水位设置在板桩的左侧,以在两侧施加水力梯度。假设土壤是均质且各向同性的,土层底部被定义为不透水边界。板桩作为不透水屏障,允许渗流通过土层。(a) 板桩周围的二维渗流分析(b) 使用Fellenius切片法的二维边坡稳定性分析(c) 通过图像处理对部分饱和砂的XCT图像进行相位分割图1. 岩土工程中问题场景的示意图拉普拉斯方程用于描述土壤域内水力头的分布。可以应用有限差分法(FDM)来获得近似解。为了实施FDM,将土层等效离散为在x方向上间隔为dx、y方向上间隔为dy的正方形网格,如图1a右侧网格化框中所示。通过中心差分法来近似拉普拉斯算子,第i行j列节点的水力头可以通过以下公式获得:基于水位、初始条件和边界条件,可以迭代计算每个节点的水力头,直到收敛为止。最终水力头值的分布可以用于描绘等势线和流线。使用ChatGPT生成的代码所得结果与使用商业软件GeoStudio SEEP/W得到的结果进行了比较。2.2. 边坡稳定性分析边坡稳定性分析的关键参数包括边坡几何形状,如高度(h)和边坡倾角(β),以及土壤性质,包括黏聚力(c)、单位重度(γ)和摩擦角(ϕ)。边坡的稳定性通过安全系数(FS)来量化,FS定义为可用抗剪强度与作用在潜在破坏面上的剪应力之比。需要注意的是,最小化FS的潜在破坏弧及其位置指示了最脆弱的区域,称为“临界破坏面”。Fellenius切片法,又称普通切片法,被用于给定边坡几何形状的边坡稳定性分析(Fellenius, 1936)。该方法假设潜在破坏面为圆形滑动面,并将可能滑动的土体分割为多个垂直切片以评估FS。它还假设切片间的作用力平行于每个切片的基底,且作用在侧面的力在破坏面上的法向为零。因此,FS可以仅基于力矩平衡推导出来,从而进行确定性的边坡稳定性分析。在此分析中,假设孔隙水压力为零。潜在破坏面假设为如图1b中红色弧线所示的圆弧,其中心坐标为C (x, y),半径为r。作用在每个切片基底(如蓝色阴影梯形所示)上的法向力(Nn)和剪切力(Tn)为:其中,Wn 是第 n 个切片的重量,αn 是第 n 个切片基底的倾角。基于力矩平衡的安全系数(FS)定义为相对于中心 C 的每个切片的驱动力矩(Md)之和与抗力矩(Mr)之和的比值。第 n 个切片的驱动力矩、抗力矩和 FS 可以通过以下公式确定:通过改变潜在破坏弧的中心和半径,可以计算出不同情况下的安全系数(FS),最终确定最小的FS及其对应的几何形状。ChatGPT实现了上述过程,并使用商业软件GeoStudio SLOPE/W对结果进行了验证。2.3. 相位分割的图像分析 部分饱和土壤包括三种相态,即固体颗粒、水和空隙(空气),其切片通过二维X射线计算机断层成像在8位灰度下可视化,如图1c-i所示。A-B段的CT值轮廓线如图1c-ii所示,可以大致勾勒出相邻相态之间逐渐变化的边界。图1c-iii中的强度水平(即CT值)直方图清楚地表明,空气-水界面和水-固体界面的边界值分别位于大约91和113附近。使用传统的分割方法时,由于水的CT值与界面处逐渐变化的CT值相似,水不可避免地会在固-气界面附近被误分割。因此,从灰度图像中分割多个相位不仅需要采用传统的阈值分割方法,还需要使用先进的图像处理技术,如迭代膨胀-侵蚀法(Adams和Bischof, 1994; Serra, 1982)。在本研究的范围内,目标是仅基于可感知的提示逐步开发由ChatGPT驱动的数值代码,而无需明确提供图像处理算法。2.4. 通过使用ChatGPT生成代码解决问题的过程 图2概述了问题解决的工作流程。对于每个案例,通过描述几何形状、初始条件、边界条件和计算方法来清晰地定义问题。需要注意的是,提示由叙述性描述组成,不包含任何数学公式,从而允许ChatGPT根据其理解自由生成数值代码。在请求ChatGPT生成MATLAB代码后,生成的代码被执行。一旦代码运行无误,用户会仔细审查代码和后续的结果以进行验证。根据此审查,用户会提供额外的提示来纠正代码,直到它完全没有错误。这些提示会不断被细化和优化,以确保ChatGPT生成的结果与预期的正确答案一致。图2. 通过使用ChatGPT进行自动代码生成的解决问题流程图3. 用户提示、ChatGPT回答及结果 3.1. 二维渗流分析 图3展示了首次尝试为二维渗流分析生成MATLAB代码的过程。图3的上框中包含了全面的问题定义,如几何细节、初始和边界条件、空间域网格配置的规格以及迭代计算的收敛标准。需要注意的是,虽然提示明确要求使用有限差分法(FDM),但并未详细说明该方法的具体方程或操作原理。图3下部分显示了由ChatGPT生成的结构合理的MATLAB代码。为了清晰起见,生成的代码中每个部分基于用户提供的提示被用不同颜色标注。图3. 第一次向ChatGPT提供的渗流分析提示(上框)及ChatGPT生成的相应MATLAB代码(下框)在第一次尝试中,用户定义的特定尺寸变量正确地反映在第21行。第24到29行(以绿色高亮)建立了土层的分析域、初始条件以及有限差分法(FDM)应用的收敛标准。土层根据用户定义的间隔dx和dy被划分为网格,每个网格点被分配给矩阵变量H,表示总水头值。该矩阵H根据初始条件初始化,反映了板桩两侧的平均水位。随后的代码部分(第32-34行,以黄色高亮)定义了板桩在土壤域中的位置,并通过为两侧的水位分配常量水头值H1和H2来建立边界条件。利用用户指定的阈值作为FDM的收敛标准,ChatGPT自主构建了第37-52行(以蓝色高亮)的FDM迭代计算代码,该代码集成了while循环、for循环和if语句。值得注意的是,在这个迭代部分中,第43-45行包含了一个排除板桩内计算的条件语句。为了简洁起见,绘制流网的部分被省略。第一次尝试生成的MATLAB代码运行无误,ChatGPT生成的流网如图4a所示。GeoStudio SEEP/W生成的流网结果也显示在图4b中。对比这两个图可以发现,关键差异源于第32行中与水位相关的不正确边界条件语句。在矩阵H的第一行中,板桩左侧应指定为高水位H1 = 10米,而右侧应为低水位H2 = 0米,作为常量水头值。图4. 通过ChatGPT生成的MATLAB代码在第一次尝试中分析的流网(a)与在相同几何和边界条件下由GeoStudio SEEP/W分析的基准结果(b)的比较第二次尝试的提示如图5所示,旨在纠正矩阵H的边界条件,代码也在第14-15行通过正确分配常量水头值进行了修正。然而,这一修正并未产生正确的流网,如图6a所示。原因是“越界”情况下缺少了用于特定节点(如AF、FE、ED、CH、HG和GB线上)的有限差分法(FDM)计算所需的四个相邻值。因此,第三次尝试的提示详细解释了常用的“镜像”技术,以用于那些缺少相邻值的节点的FDM计算,这部分内容在图5中以蓝色高亮显示。首先在板桩区域分配了“NaN”值(第18-19行),然后检查“NaN”值并正确实施镜像技术(第30-33行)。图5. 为解决错误向ChatGPT提供的第二次和第三次尝试提示(上方两个框),以及ChatGPT生成的相应修订后的MATLAB代码(下方框)(a) 第二次尝试结果 (b) 第三次尝试结果图6. 通过ChatGPT生成的MATLAB代码分析的流网通过迭代改进提示,最终生成的流网(图6b)使用ChatGPT自动生成的MATLAB代码,与基准流网准确匹配。尽管SEEP/W基于有限元法(FEM),但两种方法生成的矩阵H中的水头值比较后显示出差异极小。3.2. 边坡稳定性分析第一次尝试的提示指定使用Fellenius方法来计算安全系数(FS),并包括变量的声明,如边坡几何形状、土壤性质、坐标和表示潜在破坏面的弧的半径(图7)。需要注意的是,虽然建议使用Fellenius方法,但其基本原理并未解释。ChatGPT生成的代码以用户自定义函数“Fellenius”的形式显示在图7的下部分。从第6行到第16行,生成了滑动体的坐标、滑动体切片的代码,以及每个切片重量的计算。第19行到第21行定义了驱动力和抗力的基本方程,以计算FS。图7. 第一次向ChatGPT提供的用于创建计算FS的MATLAB函数的提示(上框)及ChatGPT生成的相应MATLAB代码(下框)初始代码虽然没有语法错误,但包含逻辑错误,如图8所示。计算边坡的安全系数(FS)需要识别滑动体、将其分为N个切片,然后基于每个切片的重量计算驱动力和抗力。错误情况1表明滑动体识别不正确,具体是由于边坡长度范围设置错误(红色箭头显示的错误定义长度和黑色箭头显示的正确长度)。在错误情况2中,分成N个部分的切片被视为矩形,且每个切片的高度错误地设置为第n个y坐标值(红色箭头所示)。这应通过考虑边坡的y坐标与切片基底的y坐标之间的差异进行修正,如黑色箭头所示。错误情况3中,用于计算FS评估的切片基底倾角的变量错误地设置为与边坡相关的倾角。图8. (a) 三种不同破坏类型的分析域示意图,以及使用ChatGPT在第一次尝试中生成的MATLAB代码计算的FS与GeoStudio SLOPE/W结果的准确性评估比较。(b) ChatGPT在第一次尝试中创建的MATLAB函数中的三个错误导致了不正确的FS计算结果。为了纠正这些错误,向ChatGPT提供了额外的提示以指导代码改进,修订后的MATLAB函数代码如图9所示。错误1被修正为确保所有部分都包含在边坡长度内(第9行)。为纠正错误2,ChatGPT添加了一个for循环,根据每个切片在边坡长度内的位置不同来计算其高度(第16-24行),并考虑到梯形的形状。对于错误3,替代了在确定边坡几何形状时提供的倾角变量(beta),在第29行引入了一个新定义的角度(alpha),以考虑未定义的第n个切片的倾斜度。这个变量描述了从破坏圆心到第n个切片基底的线与垂直线之间形成的角度。图9. 为解决三种不同错误而向ChatGPT提供的第二次尝试提示(上框)及ChatGPT生成的相应修订后的MATLAB代码(下框)该代码在第30-32行使用角度变量alpha计算驱动力和抗力,从而得出FS值。为了简洁起见,主函数(使用for循环和MATLAB函数代码来计算潜在破坏圆每个坐标的FS值)被省略。对于具有不同半径的潜在破坏圆心范围,ChatGPT得到了最小的FS值为1.630(图10a)。热图表示了给定坐标下计算得到的FS值的分布。这个结果与GeoStudio SLOPE/W的结果相同(图10b)。在这两种情况下,趾部破坏类型的临界破坏圆的坐标x、y和半径r分别确定为-0.30、14.78和14.77,对应的最小FS值计算为1.630。图10. 使用热图可视化潜在破坏面的FS计算结果,并展示临界破坏面。结果分别通过(a) ChatGPT创建的MATLAB函数和(b) GeoStudio SLOPE/W获得。(边坡尺寸:h = 10米,β = 40°,土壤性质:Ф = 35°,c = 12 kN/m²,γ = 17 kN/m³)。3.3. 相位分割的图像分析 图11展示了对部分饱和砂的X射线CT图像中三种相态(即颗粒、水和空气)进行分割的第一次尝试。每个相态的CT值范围被大致提供,并在提示中指定了分割相态的颜色。通过使用名为Link Reader的ChatGPT插件,上传了8位灰度图像。ChatGPT生成的代码与提示中的输入几乎相同,首次输出的图像显示在图11底部。需要注意的是,为了简洁起见,加载图像的代码被省略。由于部分体积效应(PVE),使用第一次尝试中编码的传统阈值方法无法避免在颗粒和空气边界周围捕获到的不必要的环状水相。由于空间分辨率的限制,PVE往往会导致相位边界处不同物体的模糊和融合(Glover和Pelc, 1980; Herman, 1980)。图11. 向ChatGPT提供的用于生成执行相位分割的MATLAB代码的第一次尝试提示(上框)及ChatGPT生成的相应MATLAB代码(下框)因此,提供了第二次提示,告知ChatGPT部分体积效应(PVE)的起源以及所需的输出图像,而未指定相关算法。图12展示了第二次提示、ChatGPT生成的MATLAB代码以及输出图像。ChatGPT自愿提出了腐蚀和膨胀的方法,这确实已被广泛应用于克服PVE(Haralick等,1987;Serra,1982)以去除伪影,正如ChatGPT在代码之前的回答中准确解释的那样。ChatGPT最初设置的 `strel` 函数和跨度值 `3` 可以由用户调整以产生所需的结果。然而,这种方法往往过度去除了环状水相,导致饱和度从59.03%降低到51.66%,如图12底部的输出图像所示。第一次尝试中捕获的被去除水相的轮廓在图12底部接触区的高亮图像中以红色呈现。图12. 向ChatGPT提供的用于生成执行相位分割的MATLAB代码的第二次尝试提示(上框)及ChatGPT生成的相应MATLAB代码(下框)在第三次尝试中,提示建议使用替代策略来对抗PVE,例如“边缘保留平滑”或“平滑约束”,这些方法已被用于在减少噪声的同时保留PVE引入伪影的重要结构细节(Rudin等, 1992;Smith和Brady, 1997;Tomasi和Manduchi, 1998)。图13展示了第三次提示、ChatGPT生成的代码以及结果。这导致ChatGPT使用了MATLAB的 `imbilatfilt` 函数(第15行),该函数基于像素的接近度和强度相似性执行非线性和边缘保留平滑。由于这种方法特别补偿了粒子接触区域附近的水相,饱和度略微回升至52.83%。图13. 向ChatGPT提供的用于生成执行相位分割的MATLAB代码的第三次尝试提示(上框)及ChatGPT生成的相应MATLAB代码(下框)4. 结论与讨论 本文探讨了在岩土工程中利用ChatGPT辅助编程来解决渗流、边坡稳定性和X射线图像处理等基本但重要的问题的可行性。对于每个案例,通过对话式提示向ChatGPT提供了问题定义、初始条件和边界条件,而无需直接输入数值或数学定义与方程。对初始提示生成的MATLAB代码进行了校对并执行,以验证结果的准确性。通过为需要进一步修改的代码提供附加提示,代码被逐步更新,直到能够生成准确且可接受的结果。提出的用于自动化编程任务的ChatGPT-MATLAB框架引出了以下几点观察:ChatGPT能够逻辑性地建立常见编程序列中的一系列过程,包括变量和域的定义、控制方程的制定、迭代操作与收敛(如有必要)以及结果的可视化。由于ChatGPT在处理复杂数学问题或高级逻辑任务时并不总是准确(Plevris等, 2023),因此有必要通过清晰描述的附加提示来识别和改进初始代码中生成不佳的部分,直到生成预期的结果。从ChatGPT得出的最终输出与商业软件获得的结果一致,这验证了ChatGPT的编程功能表现。ChatGPT能够显著减少手动编码中常见的语法错误,这一点非常有利。然而,本研究强调了在使用ChatGPT解决给定的岩土工程问题时,人类专业知识的重要性。尽管ChatGPT为半自动化计算任务提供了有前途的途径,但其有效性取决于用户提供详细且准确提示的能力。研究结果进一步表明,ChatGPT并不能替代人类的判断力和专业知识,而是可以加速和简化编程过程的工具。用户必须具备基础的学科知识和编程原理,才能有效地利用ChatGPT。本质上,ChatGPT能够加速代码生成的复杂过程,但并不能完全替代它,仍然需要人类的监督和伦理考量,以确保结果的准确性和相关性。在本研究中展示的借助ChatGPT进行数值编码的案例,仅是生成式人工智能在工程、研究和教育领域广泛潜在应用中的一个例子。预计这一技术的能力将扩展到处理更复杂的工程问题,例如多过程耦合建模、设计任务以及大数据分析。生成式人工智能的潜在应用及其在科学和工程领域的作用,是应该由我们社会共同讨论的问题。更广泛的社区应参与这些可能性的讨论,因为集体的经验和知识将提升生成式人工智能在研究和教育中的应用效果。翻译转载自:《Computers and Geotechnics》: &quot;A ChatGPT-MATLAB framework for numerical modeling in geotechnical engineering applications&quot;来源:多相流在线

未登录
还没有评论
课程
培训
服务
行家
VIP会员 学习 福利任务 兑换礼品
下载APP
联系我们
帮助与反馈