文章快速检索     高级检索
  空气动力学学报  2022, Vol. 40 Issue (1): 26-32  DOI: 10.7638/kqdlxxb-2021.0117

引用本文  

陈子玉, 周楷文, 刘应征, 等. 基于压缩感知的高频响流场重构方法及其应用[J]. 空气动力学学报, 2022, 40(1): 26-32.
CHEN Z, ZHOU K, LIU Y, et al. High-frequency flow field reconstruction method based on compressed sensing and its applications[J]. Acta Aerodynamica Sinica, 2022, 40(1): 26-32.

基金项目

国家自然科学基金(12072196);上海市科学技术委员会“科技创新行动计划”(19JC1412900)

作者简介

陈子玉,江苏扬州人,女,硕士生,研究方向:流动控制,实验测量方法. E-mail:chen_z_y@sjtu.edu.cn

文章历史

收稿日期:2021-07-08
修订日期:2021-08-19
优先出版时间:2021-11-22
基于压缩感知的高频响流场重构方法及其应用
陈子玉 , 周楷文 , 刘应征 , 温新     
上海交通大学 机械与动力工程学院,上海 200240
摘要:粒子图像测速(PIV)方法具有高空间分辨率的优势,但是往往受到采样频率的限制(一般在15 Hz以下),难以完成高频响测量。压缩感知(CS)能够基于稀疏采样数据获得高频信息,但如果直接应用于所有的数据点则计算量过大。基于亚采样(sub-Nyquist)PIV数据,本文提出了基于压缩感知和本征正交分解(POD)的高频响流场重构方法。首先采用POD对数据降维,同时获取空间模态和相应的亚采样时间系数,将亚采样时间系数作为观测值,选取适当的稀疏基,通过求解基追踪问题来计算高频响的模态系数。结合空间模态和所得到的时间分辨模态系数,可以重构高频响的非定常流场。利用该方法分别对周期性的振荡器流场和非周期性的不同直径圆柱流场进行重构,检测该方法的适应性。结果表明,压缩感知方法无需侵入式的辅助测量,可以为周期性流场提供准确的重构,重构误差低于3%,而对于非周期性复杂流场,则出现较大的高频噪声。因此,本文所提出的方法可以应用到周期性流场中以提高测量数据的时间分辨率。
关键词粒子图像测速    本征正交分解    压缩感知    流场重构    亚采样    
High-frequency flow field reconstruction method based on compressed sensing and its applications
CHEN Ziyu , ZHOU Kaiwen , LIU Yingzheng , WEN Xin     
School of Mechanical Engineering, Shanghai Jiao Tong University, Shanghai 200240, China
Abstract: Though particle image velocimetry (PIV) has a high spatial resolution, it is often limited by the temporal resolution, which is typically below 15 Hz and difficult to acquire high-frequency flow information. Compressive sensing (CS) is able to capture the high-frequency information based on sparse sampled data, but if directly employed to all spatial points, the massive data poses a great computational cost. In this study, a high-frequency response flow field reconstruction method via coupling POD and CS is proposed to solve the above difficulty. POD is used first to reduce the dimensionality of the sampled PIV data, and the spatial modes and corresponding time coefficients are calculated at the same time. By taking the time coefficients at the sub-Nyquist sampling points in CS as the observable and selecting an appropriate sparse basis, the high-frequency mode coefficients can be computed via the basic pursuit method. The unsteady flow field with high-frequency temporal information can then be reconstructed via combining the spatial modes and the corresponding resolved time coefficients. Two cases, i.e. the periodic oscillator flow and the non-periodic double-cylinder wake flow, are used to test the adaptivity of the proposed method. The results show that without the aid from contact measurement, the proposed method can accurately reconstruct the original flow field for period flow and the reconstruction error is less than 3%, while for nonperiodic complex flow, the error is significantly increased due to high-frequency noise. Overall, the proposed method can be used for periodic flow to improve the temporal resolution of the measured experimental data.
Keywords: particle image velocimetry    proper orthogonal decomposition    compressive sensing    flow field reconstruction    sub-Nyquist sampling    
0 引 言

流场的时空变化包含了重要的流动信息,但实际测量却难以同时得到较高时间分辨率和空间分辨率的流场。粒子图像测速(PIV)[1-4]诞生于20世纪80年代,是一种可获得高空间分辨率的光学流场测试技术,但其时间分辨率通常限于15 Hz以下。尽管利用高速相机和激光器可以获得相对较高时间分辨率的PIV流场信息,但价格昂贵,且有时受实验条件限制,而实际应用中高频流场非常常见,因此有必要发展高频PIV流场重构方法。

为了提高PIV测量流场的时间分辨率,国内外诸多学者进行了相关研究。Druault等[5]对所测得的PIV数据进行本征正交分解,然后直接对模态系数进行样条插值,从而得到高时间分辨率的流场。这种方法鲁棒性好且易于实施,但这种方法本身依赖于高时间分辨率采样,频率提高有限。之后,许多研究者借用具有高时间分辨率的局部点测量技术,如热线、麦克风等,来提高PIV的时间分辨率。He等[6]在使用PIV 测量技术的同时,附加了34个局部点测量的高频采样传感器,结合线性随机估计和数据同化,重构出时间分辨率的湍流流场,但这种方法需要解控制方程,计算量较大,且在某些工况时,无法安装大量的高频采样传感器。Jin[7]等使用双向循环神经网络,并借助于高频局部点测量和PIV测量,恢复出时间分辨率的圆柱绕流流场,但这种方法仍需要消耗较大计算量,且加入的局部测量点同样会扰乱流场。黄路[8]针对非定常流场,结合互相关算法和卡尔曼滤波,利用仿真验证了时间分辨率实时自适应调整算法的有效性,但这种方法依赖于纳秒级分幅成像的硬件装置。因此,目前针对该问题的研究存在许多缺陷,如仅使用PIV测量的时间分辨率重构依赖于自身的高时间分辨率测量或依赖于昂贵的硬件,而结合侵入式点测量的重构在某些工况下甚至无法安装点测量仪器或会扰乱流场。

考虑上述问题,本文提出了一种基于压缩感知(CS)从亚采样的PIV数据中重构高时间分辨率流场的方法。该方法无需借助侵入式的高频点测量,使用压缩感知对亚采样的POD系数进行高分辨率重构,可实现12倍的频率放大,且不依赖于高时间分辨率的PIV测量,利用常用的PIV测量装置即可实现。该方法将随机采样得到的PIV数据进行本征正交分解,将所得到的时间系数利用压缩感知进行重构,继而得到高时间分辨率的系数,结合分解得到的空间模态,重构出时间分辨率的流场。本方法首先在较为简单的周期性实验数据中进行验证,并且为探索该方法的适用边界,将该方法进一步应用至非周期的复杂流场(不同直径的双圆柱尾涡流场)中。研究结果表明,该方法可以对周期性流场实现较好的重构,而对于非周期流场则表现较差。

1 流场重构原理 1.1 压缩感知概述

压缩感知理论[9-10]由Candès、陶哲轩及Donoho等提出。该理论指出,只要信号在某个变换域是稀疏的,那么就可以用一个与变换基不相关的观测矩阵将变换所得的高维信号投影到一个低维空间上,然后通过求解一个优化问题就可以从这些少量的投影中以高概率重构出原信号[11]。由于压缩感知方法基于信号的随机采样,观测矩阵 ${\boldsymbol{ \phi }} $ 定义了在高时间分辨率全流场F上的随机采样规则,基于这样的随机采样规则,亚采样的PIV数据G得以产生,可以定义为:

$ {\boldsymbol{G}}={{{\boldsymbol{\phi}} }}^{\rm{{T}}}{\boldsymbol{F}} $ (1)

如式(3)所示,对随机采样的PIV数据G进行本征正交分解,可以获得空间模态 $ \pi $ 和相应于数据G的亚采样时间系数 $ {A}^{sa} $ 。在本文中,由于时间系数在离散余弦变换的变换空间中是稀疏的[12],因此可将离散余弦变换(DCT)作为所选取的变换域。通过压缩感知,可以从亚采样时间系数 $ {A}^{sa} $ 中重构出系数A,重构过程求解最优化问题:

$ \begin{array}{l} {\rm{Min}}\;\;\;{\left\| S \right\|_0}\\ {\rm{s.t.}}\;\;\;\;\; {A^{sa}} = {{\boldsymbol{\phi}} ^{\rm{T}}}A = {{\boldsymbol{\phi}} ^{\rm{T}}}\psi S \end{array} $ (2)

其中, $ \psi $ 为离散余弦变换,S为时间分辨率系数A在离散余弦变换中的系数。但是,求解上述问题是一个NP难问题[13]。Chen、Donoho和Saunders[14]指出,求解一个更加简单的 $ {L}_{1} $ 优化问题会产生同等的解,使得原问题简化成一个线性规划问题。本文中,选择典型的基追踪算法,使用Matlab中的linprog函数进行求解。结合所得到的空间模态,时间分辨率的流场重建完成。更多细节可见Baraniuk等[13]、Zhao等[15]、Stankovic等[16]对压缩感知理论的论述。

1.2 本征正交分解概述

本征正交分解[17-18]是一种对流场降维,获取流场中主要流动成分的方法,已经得到了广泛应用。记g为瞬态流场数据集G在某时刻的数据,g可以表示为:

$ g=\bar{g}+\sum _{i=1}^{w}{\pi }_{i}{a}_{i}^{sa} $ (3)

其中: $ \bar{g} $ 表示平均流场, $ {\pi }_{i} $ 表示空间模态, $ w $ 表示模态总数, $ {a}_{i}^{sa} $ 表示相应于某一模态的时间系数。由于在该方法中,数据G为PIV采样数据的集合,因此 $ {a}_{i}^{sa} $ 所记录的时刻点相应于数据G的采样时间。

时间系数 $ {a}^{sa} $ 可通过奇异值分解计算得到:

$ ({a}_{{\rm{norm}}}^{sa},\lambda )=\mathrm{s}\mathrm{v}\mathrm{d}\left({\boldsymbol{C}}\right) $ (4)

其中,C表示对流动数据内积 $ \left({\boldsymbol{G}},{\boldsymbol{G}}'\right) $ 进行空间互相关计算得到的矩阵。将数据G映射到式(4)计算所得的时间系数 $ {a}^{sa} $ 上并进行归一化即可得到空间模态 $ \pi $ 。特征值 $ \lambda $ 是逐渐递减的序列,表示模态 $ \pi $ 所占有的能量。模态所占有的能量比例通常被记为 $\dfrac{{\lambda }_{i}}{\displaystyle\sum {\lambda }_{i}} {\text{×}} 100{\%}$

1.3 流场重构流程

图1所示,首先对低时间分辨率随机采样的PIV数据G进行本征正交分解,得到空间模态及相应的时间系数,由于该时间系数与PIV数据G遵循相同的采样规则,因此式(1)中的采样矩阵即为时间系数的观测矩阵,选取离散余弦变换作为稀疏基,从而重构得到高时间分辨率的模态系数,结合前面所得到的空间模态,即可完成流场重构。


图 1 流场重构流程 Fig.1 Flowchat of the flow field reconstruction method
2 算例验证与分析 2.1 振荡器流场的验证与分析

振荡器无需任何移动部件即可产生周期性的射流。在流动控制、换热等应用领域中,这种振荡器射流具有良好的鲁棒性[19],但其振荡频率远高于PIV采样频率。针对这个问题,该时间分辨率流场重构方法首先应用于振荡器实验数据,来验证该方法的可靠性。PIV拍摄区域如图2所示。为验证重构准确性,使用TR-PIV获取真实流场信息,有关该实验的更多详细信息可见Wen等[20]的论文。在本文中,PIV采样数据为间隔采样时间超过0.444 s的随机采样,振荡器频率为2.687 Hz,而重构所得的高时间分辨率全流场数据采样频率为27 Hz。


图 2 振荡器PIV拍摄区域示意图 Fig.2 Schematic of the PIV measurement region of the oscillator flow field

为了选择适当数量的模态,对PIV数据的POD模态累积能量谱进行分析,由图3所示,可以发现,在一共7208个模态中,随着模态数的增加,每一模态所占据的能量比例逐渐降低,前7个模态占据了振荡器射流约72%的能量,而前100个模态占据了约88%的能量,因此,在下面的分析中,首先选取前7个模态用于重构,继而选择前100个模态用于对比。


图 3 振荡器流场POD模态能量谱 Fig.3 Energy spectrum of the POD modes of the oscillator flow field

图4对比了该模态系数在时间域上的真实值及重构值。低频采样数据G所得的模态系数用黑色菱形标出,可以看出该重构得到的模态系数与真实值有着较好的拟合。尽管相应于其他模态的系数这里未给出,但重构效果基本类似,不再赘述。


图 4 振荡器流场时间域上重构模态系数与真实值对比 Fig.4 Comparison of reconstructed POD coefficients with ground truth data in the temporal domain for the oscillator flow field

基于图4重构得到的模态系数,结合POD分解所得空间模态,完成高时间分辨率流场重构。随机选取某周期中未被采样的三张快照观察其重构效果(图5)。可见,7阶模态重构流场准确获取了真实流场的主模态相关信息,但由于所选模态较少,部分细节信息有所丢失。针对该问题,进一步增加重构阶数至100阶。当使用100阶模态进行重构时,重构流场的细节信息更为丰富。


图 5 振荡器流场非采样点重构与真实值对比 Fig.5 Comparison of the reconstructed flow field with the ground truth data at a non-sampling point for the flow oscillator

为量化重构准确率,误差由Ruscher等[21]的分析定义为:

$ {e}_{i}=\frac{{\Big|u}_{i}^{{\rm{recon}}}-{u}_{i}^{{\rm{truth}}}\Big|}{{u}_{{\rm{max}}}^{{\rm{truth}}}-{u}_{{\rm{min}}}^{{\rm{truth}}}} {\text{×}} 100{{\text{%}}} $ (5)

其中, ${u}_{i}^{{\rm{recon}}}$ 表示在一张PIV快照中空间 $ i $ 位置的重构值, $ {u}_{i}^{{\rm{truth}}} $ 表示空间 $ i $ 位置的真实参考值, $ {u}_{{\rm{max}}}^{{\rm{truth}}} $ $ {u}_{{\rm{min}}}^{{\rm{truth}}} $ 分别表示真实PIV数据F在时间域和空间域的最大值和最小值。最终每张快照的误差计算值对上式所有所得 $ {e}_{i} $ 进行平均,继而排除采样点位置处的快照数据,对所得每张快照的误差进行平均。所得误差曲线如图6所示。当选取模态阶数由3阶升至20阶时,重构误差从2.8%降至2.7%,而随着重构阶数进一步增多至100阶时,误差稍微增加。经过分析,认为该误差的增加是由高阶模态包含更多噪声所导致的。


图 6 振荡器流场重构误差随选取模态数量的变化 Fig.6 Error incurred in data reconstruction vs the number of POD modes employed of oscillator flow field
2.2 双圆柱流场的验证与分析

为了探索提出的方法在非周期流场中的重构效果,对不同大小的双圆柱尾涡流场原始PIV数据进行了处理分析(图7)。实验在循环水槽[22]中进行,实验的具体细节参见文献[23]。圆柱采样频率为250 Hz,因此本文用于重构的数据间隔采样大于0.048 s。


图 7 双圆柱PIV拍摄区域示意图 Fig.7 Schematic of the PIV measurement region of double-cylinder wake flow

图8为 POD模态累积能量谱,可以发现,第一个模态占据了超过80%的能量,前7个模态占据了约92%的能量,因此,在本案例中,仍然选用前7个模态进行重构。


图 8 双圆柱流场POD模态能量谱 Fig.8 Energy spectrum of the POD modes of the double-cylinder wake flow

相应于第一模态,图9对比了该模态系数的真实值及重构值,如图9(a)所示,由低频采样数据 $ {\boldsymbol{G}} $ 所得的模态系数用蓝色五角星标出。在非采样点,重构得到的模态系数偏离于真实参考值,存在较多的高频噪声,而结合图9(b)所示,尽管在非采样点重构值偏离于真实值,但也是在真实值附近波动。这是由于 $ {L}_{1} $ 范数无法区分稀疏系数尺度的位置,所以存在低尺度能量搬移到了高尺度的现象,从而容易出现一些人工效应。即,如图9(a)所示的一维信号高频振荡现象。整体上来说,重构信号在欧氏距离上仍然逼近原信号。相应于其他模态的系数未在图中示出,但重构效果基本类似,所以不再赘述。


图 9 双圆柱流场时间域上重构模态系数与真实值对比 Fig.9 Comparison of the reconstructed POD coefficients with the ground truth data in the temporal domain for the double-cylinder wake flow

随机选取1张在非采样点的重构快照,见图10。可以发现虽然压缩感知方法可以重构流场主要特征,但是难以捕捉高频流场细节。相比于周期性的振荡器射流流场,误差显著增大。这说明,该方法对非周期性流场重构适应性较低。


图 10



图 10 双圆柱流场非采样点7阶模态重构流场与真实流场对比 Fig.10 Comparison of the reconstrcucted flow field using leading 7 POD modes with the ground truth data at a non-sampling point for the double-cylinder wake flow
3 结 论

研究提出了一种基于压缩感知和POD的PIV流场高时间分辨率重构方法,对低频随机采样数据进行本征正交分解,然后对时间系数进行压缩感知重构,结合分解所得的空间模态,重构出时间分辨率的流场。利用该方法处理了周期性的亚采样振荡器射流实验流场和非周期性的双圆柱实验流场。结果表明,对于周期性流场,该重构方法的重构误差在3%以下,能够准确恢复时间分辨率的流场。而对于非周期性复杂流场,则出现较大的高频噪声。

该方法为一种纯数据驱动的方法,可以实现任意选择的高阶或低阶重构,而且无需任何接触式的传感器,不会扰乱流场,且计算量较低。该方法不依赖于高时间分辨率的PIV测量,利用常用的PIV测量装置即可实现,且可实现12倍的频率放大。因此,该方法对于周期性流场的时间分辨率重构是一种行之有效的方案。但目前该方法对非周期流场的重构适应性较低,需要做进一步研究。

参考文献
[1]
许联锋, 陈刚, 李建中, 等. 粒子图像测速技术研究进展[J]. 力学进展, 2003, 33(4): 533-540.
XU L F, CHEN G, LI J Z, et al. Reserch progress of particle image velocimetry[J]. Advances in Mechanics, 2003, 33(4): 533-540. DOI:10.3321/j.issn:1000-0992.2003.04.010 (in Chinese)
[2]
WESTERWEEL J. Fundamentals of digital particle image velocimetry[J]. Measurement Science and Technology, 1997, 8(12): 1379-1392. DOI:10.1088/0957-0233/8/12/002
[3]
徐玉明, 迟卫, 莫立新. PIV测试技术及其应用[J]. 舰船科学技术, 2007, 29(3): 101-105.
XU Y M, CHI W, MO L X. PIV measurement technique and its application[J]. Ship Science and Technology, 2007, 29(3): 101-105. (in Chinese)
[4]
盛森芝, 徐月亭, 袁辉靖. 近十年来流动测量技术的新发展[J]. 力学与实践, 2002, 24(5): 1-14.
SHENG S Z, XU Y T, YUAN H J. New development in the technology of flow measurement over the last decade[J]. Mechanics in Engineering, 2002, 24(5): 1-14. DOI:10.3969/j.issn.1000-0879.2002.05.001 (in Chinese)
[5]
DRUAULT P, GUIBERT P, ALIZON F. Use of proper orthogonal decomposition for time interpolation from PIV data[J]. Experiments in Fluids, 2005, 39(6): 1009-1023. DOI:10.1007/s00348-005-0035-3
[6]
HE C X, LIU Y Z. Time-resolved reconstruction of turbulent flows using linear stochastic estimation and sequential data assimilation[J]. Physics of Fluids, 2020, 32(7): 075106. DOI:10.1063/5.0014249
[7]
JIN X W, LAIMA S J, CHEN W L, et al. Time-resolved reconstruction of flow field around a circular cylinder by recurrent neural networks based on non-time-resolved particle image velocimetry measurements[J]. Experiments in Fluids, 2020, 61(4): 1-23. DOI:10.1007/s00348-020-2928-6
[8]
黄路. 时间分辨率实时自适应PIV测量技术研究[D]. 武汉: 华中科技大学, 2018.
[9]
CANDÈS E J, ROMBERG J, TAO T. Robust uncertainty principles: exact signal reconstruction from highly incomplete frequency information[J]. IEEE Transactions on Information Theory, 2006, 52(2): 489-509. DOI:10.1109/TIT.2005.862083
[10]
DONOHO D L. Compressed sensing[J]. IEEE Transactions on Information Theory, 2006, 52(4): 1289-1306. DOI:10.1109/TIT.2006.871582
[11]
石光明, 刘丹华, 高大化, 等. 压缩感知理论及其研究进展[J]. 电子学报, 2009, 37(5): 1070-1081.
SHI G M, LIU D H, GAO D H, et al. Advances in theory and application of compressed sensing[J]. Acta Electronica Sinica, 2009, 37(5): 1070-1081. DOI:10.3321/j.issn:0372-2112.2009.05.028 (in Chinese)
[12]
DO T T, GAN L, NGUYEN N, et al. Sparsity adaptive matching pursuit algorithm for practical compressed sensing[C]//42nd Asilomar Conference on Signals, Systems and Computers, Pacific Grove, CA, USA. IEEE, 2008: 581-587. doi: 10.1109/ACSSC.2008.5074472.
[13]
BARANIUK R G. Compressive sensing[J]. IEEE Signal Processing Magazine, 2007, 24(4): 118-121. DOI:10.1109/MSP.2007.4286571
[14]
CHEN S S, DONOHO D L, SAUNDERS M A. Atomic decomposition by basis pursuit[J]. SIAM Journal on Scientific Computing, 1998, 20(1): 33-61. DOI:10.1137/s1064827596304010
[15]
ZHAO C, MA S W, GAO W. Image compressive-sensing recovery using structured Laplacian sparsity in DCT domain and multi-hypothesis prediction[C]//2014 IEEE International Conference on Multimedia and Expo (ICME), Chengdu, China. IEEE, 2014: 1-6. doi: 10.1109/ICME.2014.6890254.
[16]
STANKOVIĆ L, BRAJOVIĆ M. Analysis of the reconstruction of sparse signals in the DCT domain applied to audio signals[J]. IEEE/ACM Transactions on Audio, Speech, and Language Processing, 2018, 26(7): 1220-1235. DOI:10.1109/TASLP.2018.2819819
[17]
SANGHI S, HASAN N. Proper orthogonal decomposition and its applications[J]. Asia-Pacific Journal of Chemical Engineering, 2011, 6(1): 120-128. DOI:10.1002/apj.481
[18]
SIROVICH L, KIRBY M. Low-dimensional procedure for the characterization of human faces[J]. Josa A, 1987, 4(3): 519-524. DOI:10.1364/JOSAA.4.000519
[19]
GREGORY J, TOMAC M N. A review of fluidic oscillator development and application for flow control[C]//43rd Fluid Dynamics Conference, San Diego, CA. Reston, Virginia: AIAA, 2013. doi: 10.2514/6.2013-2474.
[20]
WEN X, LI Z Y, ZHOU L L, et al. Flow dynamics of a fluidic oscillator with internal geometry variations[J]. Physics of Fluids, 2020, 32(7): 075111. DOI:10.1063/5.0012471
[21]
RUSCHER C J, DANNENHOFFER J F, GLAUSER M N. Repairing occluded data for a Mach 0.6 jet via data fusion[J]. AIAA Journal, 2016, 55(1): 255-264. DOI:10.2514/1.J054785
[22]
LIU Y Z, SHI L L, YU J. TR-PIV measurement of the wake behind a grooved cylinder at low Reynolds number[J]. Journal of Fluids and Structures, 2011, 27(3): 394-407. DOI:10.1016/j.jfluidstructs.2010.11.013
[23]
ZHANG Q S, LIU Y Z, WANG S F. The identification of coherent structures using proper orthogonal decomposition and dynamic mode decomposition[J]. Journal of Fluids and Structures, 2014, 49: 53-72. DOI:10.1016/j.jfluidstructs.2014.04.002