Gappy POD 方法重构湍流数据的研究1)

李天一 Buzzicotti Michele Biferale Luca 万敏平 , 陈十一

* (南方科技大学力学与航空航天工程系,广东深圳 518055)

† (罗马第二大学物理系和意大利核物理研究所,意大利罗马 00133)

引言

本征正交分解(proper orthogonal decomposition,POD)是一种统计分析、简化数据集的方法.对于一个函数系综,例如在实验中获得的数据,POD 可以给出一组模态分解的基,具有很多应用中所需要的良好的性质.其中,最引人注目的就是最优性(optimality),即POD 给出了采用有限模态来刻画一个无限维过程主要成分的最有效的方式.

POD 方法由Lumley 首次引入到湍流问题的背景中[1].在大部分应用中,POD 主要用来分析实验或计算数据以获取其中的主要模态[2-5],即拟序结构(coherent structures).此外,很多研究采用POD 提供的基函数构成的低维子空间来对原始问题建立降维模型[6-10].

Everson 和Sirovich[11]于1995 年首次提出了采用POD 来重构存在缺失的(gappy)数据的方法,即gappy POD.此后很多学者针对不同的问题使用并发展了gappy POD 重构方法[12-15].Venturi 和Karniadakis[16]对于圆柱绕流的直接数值模拟数据进行了gappy POD 重构研究.通过引入适当修正,他们提升了gappy POD 的鲁棒性和精度.同时,研究了流场缺失率大小对gappy POD 重构的影响,并给出了相关解释.Gunes 等[17]对比了gappy POD 和克里格(Kriging)插值对于非稳态圆柱绕流数据的重构效果,发现当流场的时间分辨率不高或者空间缺失率较大时,克里格插值比gappy POD 重构更加有效.

在前人的工作中,以下两个因素对gappy POD的影响均未被系统研究.第一,数据本身的复杂程度.当雷诺数Re较低时,流动处于层流状态,此时描述流动的POD 模态较少[16-17];当Re很高时,流动为复杂的湍流状态,需要用来描述流动的POD 模态往往很多.给定流场的分辨率和破损区域,流动的POD模态数量越多,重构时所需要确定的未知系数就越多.第二,破损区域的几何形状.很多关于gappy POD 的研究考虑了不同缺失率(缺失面积)对重构效果的影响[15-16],但是其中的缺失区域均是随机生成的,并没有对其几何形状的约束.在破损区域面积相同时,破损区域的几何形状不同会导致流场损失的相干信息不同,从而可能会对gappy POD 的重构效果产生影响.

在本文中,选取旋转湍流流场来进行gappy POD 重构的研究.旋转湍流是一个具有丰富物理现象的湍流的典型范例.其中,通过外力输入的能量不仅导致了大尺度的气旋和反气旋结构,也产生了小尺度的、间歇的、均匀各向同性的高度非高斯的扰动,这使得旋转湍流表现为分布在不同空间尺度和时间尺度的混乱流动[18-24].此外,旋转湍流对于很多地球物理现象也有着重要的意义[25-29].

本文通过旋转湍流场数据,研究了不同流场复杂程度以及破损区域的面积大小和几何形状对gappy POD 重构结果的影响.通过更加严格地表述gappy POD 重构过程,其推导了gappy POD 重构误差的公式,并解释了它与流场复杂程度以及破损区域的关系.本工作期待为gappy POD 的应用和发展提供参考.

1 POD 理论与gappy POD 重构

1.1 POD 简介

为了在1.2 节中系统、清楚地描述gappy POD重构过程,本节简要介绍POD 分解的过程.在进行gappy POD 重构之前,需要由Nc个已知的完整训练数据来求得POD 模态.

对于第c个训练数据uc(x) ,可计算其在位置x处相对所有训练数据系综平均的偏移

1.2 Gappy POD 重构

本节介绍gappy POD 的重构过程,相比Everson和Sirovich[11]中的描述,这里进一步在数学上严格分析了不同因素对重构效果的影响,包括数据复杂程度、采用的POD 模态个数以及缺失区域.

假定考虑的流动数据最多可以由Nflow个POD模态描述,那么Nflow的大小反映了流动的自由度或复杂程度.当仅用前Ncomp个POD 分量描述系统时,对任意一个数据样本 φ (x),有

式中的截断残差r取决于Ncomp,为

一个存在部分缺失的样本可以被表示为

其中在缺失的位置x处有m=0,而在其他位置处有m=1.对于该样本,有

Everson 和Sirovich[11]中通过最小二乘拟合,即最小化如下误差

2 旋转湍流流场重构

2.1 湍流数据简介和重构问题描述

本工作选用旋转湍流流场来进行gappy POD 重构的研究.TURB-Rot 是一个公开的旋转湍流数据集[30],包含105 600 张垂直于旋转方向的二维流场切片,切片的分辨率为 6 4×64,所有切片均来自同一个三维旋转湍流的直接数值模拟.由于数据集中的旋转湍流场高度二维化,所有切片均可看作来自同一系综.数值模拟中采用了周期边界条件,因此可对二维流场切片进行如下傅里叶谱展开

由于旋转流场的各向异性,gappy POD 对u的不同分量的重构效果并不相同.不失一般性,本文考虑速度模u=‖u‖2的gappy POD 重构.在由速度u计算速度模u时,需要注意消除混淆误差[31],使得u可以写为

式中kf为32.可以看到,kf的大小确定了u所包含的傅里叶模态个数,从而确定了流场的复杂程度.因此可以通过仅保留 ‖k‖＜kf=8 的模态来获得一个低复杂度的流场数据.值得指出的是,对于现在所考虑的均匀湍流系统,可以证明其傅里叶模态就是POD模态[1].

根据式(1)～ 式(5),选取Nc=84 480 个速度模的二维切片(训练集),可以得到kf=8和kf=32 时POD 模态的特征值曲线如图1 所示.图1 表明对于kf=8和kf=32 的数据分别有193 和3205 个特征值非零的POD 模态,这与相应的傅里叶模态数目相等.

图1 POD 模态的特征值曲线Fig.1 Eigenvalues of the POD modes

获得系统的POD 模态后,选取20 480 个非训练集的数据作为测试集,来研究不同破损区域对gappy POD 重构效果的影响.图2(a)～ 图2(e) 给出了拥有相同破损面积Agap=32×32 和不同几何的破损区域示意图.其中,图2(a)～ 图2(d) 中的破损区域为均匀排布的方形区域,但是不同图片中的方形边长和数量不同;图2(e) 中的破损区域为随机分布的像素点.真实的完整流场如图2(f)所示.

图2 (a)-(e) 破损面积相同但破损区域几何不同的部分缺失流场,(f) 真实的完整流场Fig.2 (a)-(e) Damaged flow fields with gaps of the same area but different geometries.(f) The complete flow field

图2 (a)-(e) 破损面积相同但破损区域几何不同的部分缺失流场,(f) 真实的完整流场 (续)Fig.2 (a)-(e) Damaged flow fields with gaps of the same area but different geometries.(f) The complete flow field (continued)

2.2 湍流场复杂程度的影响

本节研究当分辨率相同时,湍流场数据的复杂程度对gappy POD 重构的影响.下面以Nflow代表湍流场数据的复杂程度,即特征值非零的POD 模态数目;同时以Ncomp代表gappy POD 重构采用的POD模态数目.实际计算中往往会选取Ncomp≤Nflow,以在获得一定重构精度的同时降低计算量.下面对于复杂程度不同的数据,考虑采用不同的Ncomp对gappy POD 重构误差的影响.

对于2.1 节中kf=8和kf=32 的数据,分别有Nflow=193和Nflow=3205 .为了定量描述重构效果,选取Nc=20 480 个测试数据,它们与训练数据来自数值模拟中时间相关性很弱的不同时间段.定义破损区域中的均方重构误差为

其中积分域为无破损的完整区域,A代表完整区域面积,u0为平均速度

对于一个位于流场中心的 8 ×8 的正方形破损(图2(a)),图3对于kf=8 和kf=32的数据给出了破损区域中的均方重构误差MSEgap关于NcompNflow的变化曲线.注意图中Ncomp的取值范围为 [ 1,Nflow-1] .

图3 破损区域中的均方重构误差 M SEgap 关于Ncomp/Nflow的变化曲线Fig.3 Normalized mean square error in thegap,M SEgap,as a function of Ncomp/Nflow

图4 的最小奇异值 σ min 关于 N comp / N flow 的变化曲线Fig.4 The minimum singular value of ,σ min,as a function of N comp /Nflow

本节讨论了湍流场的复杂程度对gappy POD 重构的影响.湍流场的复杂程度由对应的POD 模态决定,它们在已知点上的值构成了矩阵(式(13)).当破损区域使得为列满秩时,重构误差大小约等于流场POD 展开的截断误差除以的最小奇异值σmin(式(27)).当流场复杂度较低时,σmin的值并不十分小,重构误差主要由截断误差所影响,可以根据对重构精度和计算量的要求采用较少的POD 模态来进行重构.当流场复杂度较高时,由于截断误差会被极小的 σmin放大,此时必须采用全部POD 模态来进行有效的重构.

2.3 破损区域面积大小和几何形状的影响

本节研究破损区域大小和几何形状对gappy POD 重构的影响.考虑kf=32 的高复杂度的数据,根据2.2 节中的讨论,始终选取Ncomp=Nflow来进行重构.

给定破损区域面积Agap,定义破损尺寸(gap size)为图5 给出了破损区域中的均方重构误差MS Egap关于破损尺寸的变化曲线,其中的误差棒代表重构误差在Nc=20 480 个测试数据上的标准差.可以看到,对任何一个破损形状,总是存在一个临界的破损尺寸,使得破损小于等于该尺寸时重构误差很小,而大于该尺寸时重构误差很大..对于图2(a)～ 图2(e) 中所示的破损形状,

图5 对于不同几何形状的破损,破损区域中的均方重构误差MS Egap关于破损尺寸的变化曲线Fig.5 Normalized mean square error in the gap,M S Egap,as a function of the gap size for different gap geometries

由图5 和图6 还可以看出,对于不同的几何形状,相应的临界破损尺寸也不同.例如,对于图2(a)所示的一个方形的破损,临界破损尺寸为12;而对于图2(e) 所示的随机点破损,临界破损尺寸为30.这表明对于同样的破损尺寸,当破损形状不同时,gappy POD 的重构效果不同.图7 给出了破损尺寸为16 时,不同破损几何的重构结果图.其中,前3 列分别为破损区域、重构结果和原始区域的示意图;第4 列为红色参考线上重构结果与真实结果的分布;第5 列显示了在固定的x1处的平均重构误差

图6 对于不同几何形状的破损,N flow-r 关于破损尺寸的变化曲线Fig.6 N flow-r as a function of the gap size for different gap geometries

图7 相同破损面积,不同破损几何的gappy POD 重构结果(对应不同行).第1 列:破损的流场;第2 列:重构的流场;第3 列:原始流场;第4 列:第一列中红色参考线上重构结果(虚线)与真实结果(实线)的分布;第5 列:平均重构误差(红色),Δ u(x1),及其在测试集上的平均(黑色),〈 Δu(x1)〉 . x1 和 x2 分别代表水平和竖直方向Fig.7 Gappy POD reconstruction results for the same gap area and different gap geometries (one for each row).1st column:damaged image in input.2nd column:image generated in output.3rd column:ground truth.4th column:generated (dashed) and ground truth (solid) profiles along the vertical line shown in the 1st column.5th column:mean reconstruction error,Δ u(x1) ,for each image (red line) and the average error,〈 Δu(x1)〉,over the test images (black curve).Note that x1 and x2 denotes the horizontal and vertical directions,respectively

式中的积分区域为存在破损的x2位置,Lgap为相应的破损长度. 〈 Δu(x1)〉 代表 Δu(x1) 在测试数据上的平均.可以看到,除了一个大的正方形破损,gappy POD 对于其他几何形状的破损均可精确重构流场.虽然破损区域的总面积相同,相比于其他几何形状,一个大的正方形破损会导致损失的信息具有更大的相关性,因此导致gappy POD 无法重构流场.

本节通过考虑kf=32 的高复杂度的湍流数据,并选取Ncomp=Nflow,研究了破损区域大小和几何形状对gappy POD 重构的影响.对于任何一种几何形状的破损,都对应着一个保证重构成功的最大尺寸,称为临界破损尺寸.破损的几何形状影响着损失的信息具有的相关性,从而影响着上述临界破损尺寸的大小.

3 结论

本文研究了gappy POD 在湍流数据重构上的应用,主要考虑了湍流场复杂程度以及破损区域的大小和几何形状的影响.

通过更严格地表述了gappy POD 重构的过程,本文指出重构误差由两部分构成.第一部分来自流场POD 展开的截断误差,经由POD 基函数在已知点上组成的矩阵的最小特征值放大.当湍流场的复杂度较低时,这一项随着采用的POD 模态数目增大而减小,因此实际应用中可根据精度和计算量要求合适选取POD 模态数目.然而,当湍流场复杂度较高时,POD 基函数在已知点上组成的矩阵的最小特征值非常小,此时即使有很小的POD 截断误差也会被极大地放大.因此,为了对高复杂度湍流场进行gappy POD 重构,必须采用其所有POD 模态使得截断误差为零.

Gappy POD 重构误差的第二部分来自POD 基函数在已知点上组成的矩阵的非列满秩性,它主要取决于破损区域的面积大小和几何形状.对于固定几何形状的破损,当破损区域增大时,流场损失的信息逐渐增多,使得上述矩阵逐渐变为非列满秩的,从而使得重构失效.而对于相同的破损面积,当破损导致的损失信息包含的相关性越大(如一个大的方形破损),上述矩阵越容易变为非列满秩的,从而使得重构失效;而当破损导致的损失信息包含的相关性不大时(如随机点组成的破损),上述矩阵易于维持列满秩性质,此时流场可以被精确重构.