目录
文献阅读:求解Richards-Richardson偏微分方程,将新型PINN用于估计绿色基础设施的水力传导率,并作为性能指标
摘要
在本周中,通过阅读文献,了解使用PINN在估计求解Richards-Richardson偏微分方程,将新型PINN用于估计绿色基础设施的水力传导率,并作为绿色基础设施的性能指标,文章使用两个案例进行说明,案例1说明了在均匀土样中的变化,第二个案例为非均匀土柱中的对比。在Fluent中,选用原油输送管道内两相流的CFD分析,进行几何划分和设置求解。理论学习方面,对有luent算法基本知识进行了学习。
Abstract
In the past week, I have read the literature and learned about the use of PINs to solve Richards-Richardson partial differential equations in estimation, and the use of new PINs to estimate the hydraulic conductivity of green infrastructure as a performance indicator of green infrastructure. In Fluent, I use a CFD analysis of two-phase flow in a crude oil pipeline to perform geometric partitions and setups. In terms of theoretical learning, I learned the basic knowledge of Fluent algorithms.
文献阅读:求解Richards-Richardson偏微分方程,将新型PINN用于估计绿色基础设施的水力传导率,并作为性能指标
文献摘要
绿色基础设施(GI)是一种生态知情的雨水管理方法,具有可持续和有效的潜力。基于渗透的GI系统,包括雨水花园、透水人行道、绿色屋顶,渗透地表水和雨水径流,以补给地下水系统。然而由于缺少长期监测,该类系统容易恶化。
但我们知道,这些系统渗透水的能力取决于土壤的非饱和水力传导函数K,这样我们就可以利用体积含水量传感器的数据,通过求解Richards-Richardson偏微分方程(RRE)来估计均匀土柱的K,本文就是使用PINN来进行计算的。
文章介绍并比较了两种不同的深度神经网络架构来解决RRE和估计k:使用HYDRUS-1D模拟了三种类型的土壤水分动力学,并比较了这两种神经网络结构在k值估计方面的结果。数据来源于研究了传感器间放置对k值估计的影响,需要简单的土壤湿度传感器的数据作为输入,这些传感器在地理标志建设时很容易安装甚至改造。
讨论|结论
地理信息系统通过利用植被和渗透性材料将水引入较慢的水文路径,从而管理城市雨水,从而改善水质并提供栖息地。然而,在城市系统中,地理标志的实施受到可用空间的限制,并且需要重新建立地表和地下通道之间的联系。
此外,在安装后立即集中监测的趋势限制了我们理解系统功能演变的能力,并可能掩盖随着时间的推移而发展和恶化的问题。特别是,基于渗透的GI系统,包括雨水花园、透水人行道、绿色屋顶和生物威尔士,很容易堵塞,从而导致系统功能恶化。
基于渗透的地理特征引导地表水和雨水径流补给地下水系统。它们的入渗速率随土壤的水导率(K)而变化,K的降低会限制补给。一般通过测量压力传感器的跌落值或通过水分传感器深度之间土壤水分动态的时间差。然而,这两种方法都不能解决多孔介质的水力特性,因此我们对GI渗透机理的评估能力有限。且在实验室进行模拟的方法并不能完全代表现场状况的条件。
也就是说,水在非饱和多孔介质中的运动由于土壤含水量与导水率之间关系的迟滞而变得复杂。K的推断一般依赖于Richards-Richardson方程(RRE)的数值解。数值方法包括隐式时间离散化方案和空间有限体积或有限元离散化。通常,RRE及其逆的解是使用经验参数模型,将水力导电性定义为基质势W或体积含水量h(土壤含水量)的函数,并定义W与h之间的关系。型依赖于一些土壤性质的知识,如孔隙大小分布和饱和含水量,这些需要估计。但在野外环境中很难估计土壤的初始条件和边界条件,克服这一问题的一种方法是使用高斯过程(GP)回归建模。但是这种方法也有它自己的缺点,主要是它不能很好地处理非线性(在我们的例子中是RRE的非线性参数),并且可能需要局部线性化来规避这一点。(该方法仅限于使用水力参数的预定义模型。此外,GP建模受高斯先验假设的限制)。
文章将所设置的框架与Bandai等人提出的现有框架(架构1)进行了对比,Bandai等人提出的框架本身就是Raissi等人(架构2)描述的PINN模型。这两种架构都使用体积含水量的多深度时间序列数据,并且都使用RRE的残差作为损失函数来表示描述土壤-水动力学的物理约束。
上述架构1和架构2都不需要初始条件或边界条件的知识,这些在GI等现场条件下很难获得;它们不需要任何预定义的水力传导函数形状,也不需要事先的近似;他们只使用从简单的土壤湿度传感器阵列中容易获得的数据,他们不需要矩阵电位测量来找到反解。
文章研究了简单方便的水分传感器阵列设置的使用情况,该阵列使用三个水分传感器通过查找RRE的反解来估计K。三个传感器是本文所述方法应用的绝对最小值(为了估计体积含水量w.r.t.空间的二阶导数以求解RRE),这些传感器是:深度z0的传感器0 (S0),深度z1的传感器1 (S1)和深度z2的传感器2 (S2),如图1所示:
通过实验,第一种架构显示了k的最佳估计,而且,对于传感器放置在较浅的位置,与第二种架构相比,第一种架构具有更低的估计误差。
这两项措施表明,第一种架构可能更适用于具有明显湿度梯度的浅GI应用,例如浅且直接暴露在雨中的绿色屋顶。(这些公式不允许评估这些浅层的蒸散发影响。)除了浅系统之外,第二种架构在所有传感器位置上都有更好的性能,并且与第一种架构相比,对传感器之间的内部距离不敏感,这需要更小的内部距离来给出更好的k估计。因此,第二种架构在传感器放置配置上更加灵活。
模糊的初始条件和边界条件导致了一个不适定的问题。也就是说,传统的方法无法找到RRE的逆解,因为不能保证唯一解。然而,提出的PINN架构和许多其他类似的例子[20,24,25]已经证明,使用时间序列数据(或空间数据)找到微分方程的反解是可能的。
然而,所提出的方法在一些重要的方面是有限的。在制定RRE和网络架构时,没有考虑到一些物理过程。特别是,由于我们假设K与h、h与w之间的关系是唯一的,因此没有考虑迟滞的影响。同时忽略了热水耦合输运和溶质输运的影响。此外,这两种结构在非均匀土壤上的表现都不令人满意,因此,在特定深度K变化很重要的应用中(例如,评估植物根系影响),需要开发一个考虑多层土壤的模型。
理论基础
理查兹-理查森方程(RRE)
首先,支配水在多孔介质中不饱和流动的定律是达西定律:
式中,q为水流速度[长度,时间]-1,K为水力导率[长度,时间]-1,h为水头(总水势)[长度],w为基质水势[长度],Z为标高水头[长度]沿垂直方向Z(正向上)。注意,根据Mualem-Van Genuchten模型[13],我们假设h与w、K与h(或K与w)之间的关系是唯一的,因此没有考虑迟滞。
RRE由Darcy定律和连续性要求(不可压缩流体,不考虑水蒸气)推导而来,定义为:
它包括h和W作为变量。该公式用于第二种PINN方法(第二体系结构)。在本研究中,我们只考虑沿z维的一维流动,因此方程为:
为了得到只依赖于h的公式,我们将水分扩散率D[长度]2[时间]-1定义为:
由式(5)和式(6)可以得到我们提出的方法(第一个架构)中使用的水分公式:
只考虑z维,我们得到::
第一种PINN架构
网络的目标是估计方程(10)右侧的非线性系数,再估计的函数K,然后这些估计系数乘以它们各自的项,最后,将它们相加,以预测体积含水量的导数w.r.t.时间oh的值,即RRE的左侧,如图2a所示:
这些系数的估计是使用一组在三个固定深度h,收集的体积含水量数据来实现的。该网络使用离散时间和h相对于z的导数的有限差分近似,在深度级别z1(即中间传感器S1)进行评估。W表示网络的参数(权重和偏差),X表示其输入。将W、R、T、Z的一阶导数和一阶导数的平方;W、R、T、Z的二阶导数和中间传感器S1测量的体积含水量值作为输入电流和过去的估计值,再经由网络输出Dt,W表示网络的参数(权重和偏差),X表示其输入。
在传感器S1处(深度z1处),使用一阶有限差分近似计算体积含水量w.r.t对空间的一阶和二阶导数。故使用一阶导数w.r.t空间的前后差分作为输入。我们也使用p个过去的估计值:
正向差定义为:
而后向差定义为:
将前向和后向导数的平方值相加作为输入:
传感器S1处二阶导数估计定义为:
使用二阶导数的当前估计值和p过去估计值作为输入的一部分,将电流湿度值传感器S1添加到输入:
网络总结如下,如图2(b)所示:
网络参数W是通过最小化损失函数L得到的,我们将损失函数L定义为RRE残差r的均方误差(MSE)::
其中NB表示训练批的大小。基真一阶导数w.r.t.时间是用以下一阶有限差分近似估计的:
网络f由4层组成(图2b)。前3层是残余层。每个残差层由一个批规范化(批规范)层、一个ReLU激活函数、一个全连接层(带偏置)组成,然后是另一组批规范、ReLU和全连接层。该层的输出与其输入求和,因此称为残差层。然而,对于第一个跳过连接,我们使用大小为16、9、32的完全连接层,以便升级输入特征,使它们具有与层输出相同数量的特征。对于另外两个残差层,输入直接求和到输出。
最后一层(第四层)是由一个完全连接的层组成的。ReLU激活函数仅用于表示扩散率D^的估计的输出。这样做是为了约束它的值为正。对于输出的dk和dD不使用激活函数。用梯形法则对其进行积分,以获得水力传导函数K^的估计。
第二种PINN架构
用的框架是Bandai等人[26,31]提出的框架,它基于Raissi等人[20]的工作中的pin。pinn用于使用上一节中提到的相同的体积含水量数据集来求RRE的逆值,使用的pin由图中所示的三个网络组成:
第一个网络预测基质电位W的值:
WW为网络参数。网络fw由8个完全连接的层组成,每层后使用双曲正切函数tanh作为激活函数。然而,该网络的最后一层使用负指数函数“-exp(x)”作为激活函数,强制矩阵电位的所有预测值为负。完全连接的层由40个神经元组成。将湿势的预测值转换为对数尺度:
将上式输入两个网格,分别预测h和K
当我们使用Mualem-Van Genuchten模型时,该体系结构假设h和w之间以及K和w之间的关系是唯一的(请参阅数据描述部分)。网络fh由3个全连通层组成,网络fK由1个全连通层组成。
两个网络的完全连接层由40个神经元组成,并使用双曲正切(tanh)作为激活函数。通过最小化以下损失函数得到网络参数:
使用自动分化进行评估。在所有测量点对体积含水量的残差和误差进行了评估,由于需要计算二阶导数,使用Tanh为激活函数。
实验设置
体积含水率测量数据采用HYDRUS-1D软件生成。模拟100 cm均质土柱,间隔0.1 cm进行均匀离散。模拟了三种土壤类型:砂壤土、壤土和粉砂壤土。采用MualemVan Genuchten模型对水导率函数和土水保持曲线进行参数化:
三种模拟土壤的函数拟合参数如表1所示:
注意,h和K是w的递增单调函数。所有深度的初始基质势设置为- 1000 cm。采用Neumann边界条件作为底边界条件:
上界条件设为大气条件。为了生成时间序列数据,模拟了3天的上边界水密度通量(表2)。每0.012天采样一次体积含水量。生成场景是为了提供相同数量的训练点,用于估计第一种架构的体积含水量的导数w.r.t.z(中间传感器S1的一阶导数和二阶导数的有限差分估计)。第二个体系结构仅使用场景1中的数据:
使用深度z = - 1、- 3、- 5、- 7、- 9、- 11、- 13、- 15、- 17、- 19 cm (z正向上)的模拟测量数据来评估传感器放置配置的影响,包括深度和传感器间距离。为了比较估算的水力导率函数,我们计算了以下相对误差:
其中Nh为估算K^的点体积含水量值的个数。对于这两种结构,这个误差是在500个均匀间隔的h值下估计的,其中间隔的开始和结束对应于场景1中每个模拟土壤的数据h的最小值和最大值。对于不属于训练数据的h值,对水力传导性值进行线性插值。
结果讨论
均质土柱
显示了描述中的三种不同类型土壤的120种传感器放置配置的估计K函数的中位数,并展示了一个包络线,从估计的水力导率值的第1到第3个四分位数,允许比较两种方法之间的估计方差。表3显示了两种结构对所有土壤的中位数曲线eK (Eq. 38)的误差。我们可以看到,对于所有三种类型的土壤,与建筑1相比,建筑2的中位数曲线的误差较小。而且,架构2的第1到第3个四分位数包络比架构1的四分位数包络要小:
120个传感器位置和第1到第3个四分位包络线的第一和第二结构中位水力导率曲线。
图a、c和e是第一个体系结构的结果。图b、d和f是第二个体系结构的结果:
可以得出结论,架构2在所有传感器放置试验中具有更好的性能。我们还注意到,两种结构的估计误差在较大的水力传导函数下都较小。
当比较图5和表4中每种结构的最佳配置结果时,对于砂壤土类型的土壤,我们注意到第一种结构的误差比第二种结构的误差小,但后者对其他两种土壤的误差更小,还注意到,架构1的传感器配置的平均深度小于架构2:
为了记录传感器间距离对网络性能的影响,我们首先使用平均传感器间距离来比较两种架构的结果,而不考虑它们的深度。在图6中,计算了各配置下相应的水力传导函数的平均误差:
可看出:第二种架构(图6b、d、f)对不同内部距离的影响是不可知的。然而,第一种架构的性能(图6a, c, e)随着内部距离的增大而恶化。这与我们用有限差分近似来描述这个结构是一致的。
结果在所有土壤类型中都是一致的。我们还注意到,第一种结构的误差比第二种结构的误差大一个数量级,特别是当传感器间距离较大时。仅考虑传感器间内距离相等的误差,
可以得出类似的结论。
为了分析传感器深度对性能的影响,在图7中,我们显示了水力传导函数eK的平均相对误差,这取决于传感器放置配置的平均深度,第二种结构的误差随着传感器阵列平均深度的增加而减小(图7b、d、f),说明第二种结构在深度越深时性能越好:
相反,对于第一种架构(图7a, c, e),误差在较浅的深度更低,这意味着第一种架构在具有更锐利湿锋的浅传感器位置上取得了更好的结果。注意,图5中第一种体系结构的最佳传感器配置比第二种体系结构的配置浅,这与最后一个结论一致。然而,对于第一种结构,我们可以看到误差在10-11 cm处达到最大值,然后在更深的深度处再次降低,尽管不像浅深度3-6 cm那么低。仅考虑相对于传感器中间深度z1的误差,即X轴= d1,也可以得出类似的结论。
非均匀土柱
虽然所提出的网络旨在估计均匀土的水力传导性,但我们在模拟的非均匀垂直土柱上测试了它们的性能。这样做是因为我们期望GI的土壤在实践中是非均匀的,长21厘米,由砂壤土、壤土和粉壤土三种土壤组成,从上到下依次排列,每一种长7厘米。使用了相同的120个传感器放置配置,在图8中显示了中位数和第1到第3个四分位数包络图。将估算的水力导电性与三层土的等效垂直水力导电性Kz进行比较,其为:
中位曲线相对于建筑1和建筑2等效水力导率Kz的误差分别为0.984和0.591。预测K值的中位数虽然小于均匀土,但接近于零。此外,估计的方差很大(特别是对于体系结构2)
最佳水力导率估算结果如图9所示。架构1和架构2的误差分别为0.413和0.178:
FLuent实例:原油输送管道内两相流的CFD分析
几何建模部分
使用designmodeler进行模型绘制,首先在XY平面上绘制一个100mm直径的圆,生成草图后,对平面进行旋转,采用挤出命令,生成1000mm长度的管道,后在XY平面上重新生成一个草图平面,绘制一个80mm的圆,使用concept中的草图线命令来将80mm分裂成线条,再通过投影命令来将该圆与1管道分开,其一个为进水口,另为进油口:如下图所示:
网格划分部分
首先,采用自动方法来划分网格,检查几何模型是否可以正常绘制,后使用扫掠来划分网格,算法改为轴对称,再将O形网格分区打开,调试为8,生成网格,采用尺寸划分命令,将进水口处的网格单元缩小,每个单元尺寸为0.001m,后对其进行命名操作进入求解设置:
求解器部分
选择瞬态求解器,打开重力选项,在y上输入-9.8,在材料中添加水和油材料作为我们的两相主体,其中油需要我们主动创建,参数如下:
打开多相流模型,选择VOF,启用隐式体积力选项,因为其中存在重力和表面张力在该水油多相流中。对两相进行名称材料更改,如下图所示:
在相互作用选项中把油和水的表面张力系数填入,打开表面张力建模,打开壁面黏附,这样可以调整靠近壁面的曲率,粘性模型采用经典的K-epsilon,再采用增强壁面函数:
后进行边界条件的调整,油处进口速度改为2m/s,设置水的体积分数为1,在调整面为水,进口速度为2m/s,且进水口处油的体积分数设置为0。其余设置保持默认,壁面选择无滑面。选择二阶迎风格式来进行求解,使用标准初始化,并将原油体积分数改为0,使一开始只有水注入:
结果展示
两相流分析如下:
理论学习部分
总结
本周通过学习Fluent相关的算法,巩固了软件的理解和使用,能更好的理解在什么情况下设置正确的边界条件和如何采取正确的求解方法来进行仿真。阅读文献时会更关注其创新点和实验方式,为了以后的写作中提供更多的灵感。