﻿ 一种基于聚类分析的二维激波模式识别算法
 文章快速检索 高级检索

A two-dimensional shock wave pattern recognition algorithm based on cluster analysis
CHANG Siyuan, BAI Xiaozheng, LIU Jun
School of Aeronautics and Astronautics, Dalian University of Technology, Dalian 116024, China
Abstract: Based on compressible and inviscid flow solutions computed by shock-capturing solvers, a new technique for detecting and recognizing two-dimensional shock wave interaction patterns is proposed. The implementation process of this algorithm is illustrated from three aspects. Firstly, using a traditional shock wave detection approach based on local flow parameters, a series of grid-cells near the shock waves are identified. Next, the shock cells are divided into various clusters by means of a classical K-means clustering algorithm, and the category of each cluster is defined according to its adjacent information. Finally, a criterion is introduced to merge related adjacent clusters and to further determine the locations of shock interaction points. The clusters contained in each shock wave are recorded, and then all the fitting shock lines can be obtained by the Bézier curve algorithm. Numerical experiments show that this newly developed technique can be used in different types of mesh, The new technique produces fitted shock lines with high quality and positional accuracy. Meanwhile, the multiple shock wave interaction patterns are clearly recognized and provide a good visualization method for analyzing the motion and evolution of shock waves in complex, unsteady flow.
Keywords: shock-capturing    shock wave detection    cluster analysis    computational fluid dynamics    shock interaction    unsteady flow

1 二维激波模式识别算法

 图 1 定常激波反射：几何模型与压强云图 Fig. 1 Steady shock-wall reflections: geometric model and pressure contours
1.1 辨识激波单元

 图 2 定常激波反射：初始辨识的激波单元 Fig. 2 Steady shock-wall reflections: initially detected shock-cells

1) 计算所有单元格心压强梯度▽p在速度矢量v方向上的分量，即流向压强梯度：

 $\delta p = \frac{\mathit{\boldsymbol{v}}}{{\left\| \mathit{\boldsymbol{v}} \right\|}} \cdot \nabla p$ （1）

2) 由于数值计算难免会存在耗散和误差，会造成一些非激波区域的δp值为一个较小的正量，故过滤掉满足式(2)的单元：

 $\delta p ＜ {\xi _1}\delta {p_{\max }}$ （2）

3) 然而几乎不可能只通过设置满意的ξ1，既能完整地辨识多激波又能去除全部噪声。此外，因全局滤波强度不大，经第2步辨识的激波带可能会出现一些伪激波单元，为了在一定程度上提高辨识方法的普适性，考虑局部特性对满足式(3)的单元进行二次去噪：

 $\delta p ＜ {\xi _2}\delta p_{_{\max }}^i$ （3）

 图 3 定常激波反射：激波带的优化 Fig. 3 Steady shock-wall reflections: optimization of detected shock bands
1.2 激波单元的聚类分析

 图 4 K-means聚类算法示意图 Fig. 4 Schematic diagram of K-means clustering algorithm

1) 选择K个初始“聚类中心”。在图 4(b)中，初始化选择了两个聚类中心，即图中红、蓝“×”形点。

2) 对任意一个数据点，求其到K个聚类中心的“距离”(取欧氏距离)，将该数据点归到距离其最近的聚类中心所在的“簇”。如图 4(c)所示，所有数据点在经过第一轮遍历后被分为红、蓝标识的两个簇。

3) 利用均值等方法更新K个簇的中心值。此时对当前标记为红色和蓝色的数据点分别求出其新的质心，作为新的聚类中心，如图 4(d)所示，红色和蓝色聚类中心的位置发生了变动。

4) 重复步骤2)~3)，若所有聚类中心在迭代更新后位置都保持不变，则迭代结束；否则继续迭代。最终得到的两个簇及其聚类中心如图 4(f)所示。

K-means聚类算法操作简单，易于实现，且十分高效，一般仅需5~10次迭代即可收敛。但是研究表明，该算法对初始聚类中心的位置选取比较敏感，不同的随机初始聚类中心得到的聚类结果可能差异显著。因此，在对激波单元聚类时，合理的初始聚类中心位置显得至关重要，其既不能偏离激波带过远，也不能分布的过密或过稀。下面介绍下针对激波带初始聚类中心的选取方法。

 图 5 定常激波反射：确定初始聚类中心 Fig. 5 Steady shock-wall reflections: determination of initial clustering centers

 图 6 定常激波反射：K-means聚类结果 Fig. 6 Steady shock-wall reflections: results of K-means clustering

1) 终止簇。若流场内部某簇仅有1个相邻簇，则该簇为终止簇；其一般位于流场内部激波起始/终止处，注意该簇不与边界相邻。

2) 普通簇。若流场内部某簇仅有2个相邻簇，则该簇为普通簇；其一般呈“串状”普遍分布于激波带中，是数量最多的簇。

3) 分叉簇。若流场内部某簇有3个及以上相邻簇，则该簇为分叉簇；其一般分布于厚度较大的激波带处，如三波点、四波点和弱激波处等。

4) 边界簇。若某簇直接与流场边界相邻，则该簇为边界簇；其通常存在于激波-壁面反射处、壁面拐角处以及流场进出口边界处等。

1.3 激波的模式识别与拟合

1) 判断是否存在需要进行合并的簇。具体分两方面，若某分叉簇与其他分叉簇相邻，则将这些分叉簇合并；若某边界簇与其他边界簇(注意它们必须对应同一个边界)相邻，则将这些边界簇合并。

2) 若进行了簇合并的操作，需要重新对簇进行分类，并计算新簇的聚类中心。如图 6(a)无需进行簇合并，而图 6(b)经簇合并后变成图 7，即原4个相邻分叉簇合并后形成1个新的普通簇，一同紧靠下壁面的2个相邻边界簇合并后形成1个新的边界簇。

 图 7 定常激波反射：簇合并后 Fig. 7 Steady shock-wall reflections: after cluster-merge

3) 识别关键激波点。规定若某分叉簇与周边3个簇(见图 6(a))或4个簇(见图 23)相邻，则该分叉簇对应的聚类中心即为三波点或四波点；若某分叉簇仅与1个簇相邻，则其聚类中心即为激波终止/起始点(见图 25)；若某边界簇与2个簇相邻(图 7下边界处)，则取该边界簇对应的边界中心点作为激波-壁面反射点；同理，若某边界簇仅与1个簇相邻(图 7右边界处)，则边界中心点即为激波-边界干扰点。通过以上近邻分析，便可识别出激波-激波干扰或激波-边界干扰等模式。

4) 记录各条激波所包含的簇。从第3)步得到的某一关键簇开始搜索，依次记录下该分支上两两相邻的所有普通簇，直至搜索到另一关键簇停止，此时该一系列簇即对应一条激波；继续搜索确定下一激波分支上的簇，直到遍历完流场中所有的簇。

5) 拟合激波线。依次连接激波包含的一系列簇所对应的聚类中心，便可得到各条激波线；然而此时得到的激波线可能比较曲折(图 7中的白色实线)，为此，本文进一步采用著名的Bézier曲线拟合算法[22]获取更加光滑的激波线，下面对其作简要介绍。

 ${\mathit{\boldsymbol{B}}_n}\left( \lambda \right) = \sum\limits_{i = 0}^n {\frac{{n!}}{{i!\left( {n - i} \right)!}}{\mathit{\boldsymbol{P}}_i}} {\left( {1 - \lambda } \right)^{n - i}}{\lambda ^i}{\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} {\kern 1pt} \lambda \in \left[ {0, 1} \right]$ （4）

 $\begin{array}{l} {\mathit{\boldsymbol{B}}_3}\left( \lambda \right) = {\mathit{\boldsymbol{P}}_0}{\left( {1 - \lambda } \right)^3} + 3{\mathit{\boldsymbol{P}}_1}\lambda {\left( {1 - \lambda } \right)^3} + \\ \begin{array}{*{20}{c}} {}&{3{\mathit{\boldsymbol{P}}_2}{\lambda ^2}\left( {1 - \lambda } \right) + } \end{array}{\mathit{\boldsymbol{P}}_3}{\lambda ^3}\begin{array}{*{20}{c}} {}&{\lambda \in \left[ {0, 1} \right]} \end{array} \end{array}$ （5）

 图 8 定常激波反射：拟合的激波线 Fig. 8 Steady shock-wall reflections: fitting shock lines

 图 9 定常激波反射：激波线拟合结果与云图对比 Fig. 9 Steady shock-wall reflections: comparison of shock line fitting results and contours
1.4 算法小结

 图 10 算法流程图 Fig. 10 Flowchart of algorithm

2 算例验证

2.1 定常斜激波

 图 11 定常斜激波：聚类结果 Fig. 11 Steady oblique shock wave: clustering results
 图 12 定常斜激波：激波线拟合值与理论值对比 Fig. 12 Steady oblique shock wave: comparison of shock line fitting results and theoretical values

 来流马赫数 楔面角/(°) 斜激波的激波角/(°) 相对偏差|(βF－βT)/βT| 理论值βT 拟合值βF 1.469 3 46.502 3 46.438 3 0.14% 1.959 20 54.998 1 55.323 7 0.59%
2.2 定常激波反射

 图 13 定常激波反射：不同种类网格下的压强云图 Fig. 13 Steady shock-wall reflections: pressure contours of different types of meshes

 图 14 定常激波反射：基于结构网格的聚类结果 Fig. 14 Steady shock-wall reflections: clustering results based on structured meshes
 图 15 定常激波反射：基于混合网格的聚类结果 Fig. 15 Steady shock-wall reflections: clustering results based on hybrid meshes

 图 16 定常激波反射：非均匀网格下的压强云图 Fig. 16 Steady shock-wall reflections: pressure contours of non-uniform meshes
 图 17 定常激波反射：基于非均匀网格的聚类结果 Fig. 17 Steady shock-wall reflections: clustering results based on non-uniform meshes

 图 18 定常激波反射：激波线拟合结果对比 Fig. 18 Steady shock-wall reflections: comparison of fitting shock lines
2.3 同侧和异侧激波干扰

 图 19 激波干扰：几何模型与压强云图 Fig. 19 Shock-shock interactions: geometric model and pressure contours

 图 20 激波干扰：辨识出的激波带 Fig. 20 Shock-shock interactions: detected shock bands
 图 21 激波干扰：聚类与合并 Fig. 21 Shock-shock interactions: clustering and merging

 图 22 激波干扰：拟合的激波线 Fig. 22 Shock-shock interactions: fitting shock lines

 图 23 激波干扰：不同网格尺寸下的聚类结果 Fig. 23 Shock-shock interactions: clustering results at different mesh sizes
 图 24 激波干扰：3种网格尺寸下拟合的激波线对比 Fig. 24 Shock-shock interactions: comparison of fitting shock lines at three different mesh sizes

2.4 前向台阶非定常流动

 图 25 前向台阶绕流：不同时刻激波带的聚类结果 Fig. 25 Forward step flow: clustering results of shock bands at different times

 图 26 前向台阶绕流：拟合的激波线与压强云图对比 Fig. 26 Forward step flow: comparisons of fitting shock lines with pressure contours

 图 27 前向台阶绕流：4个时刻拟合的激波线 Fig. 27 Forward step flow: fitting shock lines at four different times
3 结论

1) 该算法在传统基于流场当地参数进行激波探测的方法基础上，对提取的激波单元进一步处理；先后结合K-means聚类算法和近邻分析成功实现了激波干扰点的准确定位；最后采用Bézier曲线算法拟合得到各条质量较高的激波线，数值试验表明激波线的大部分位置与形状受流场网格疏密影响不大。

2) 对于同侧激波干扰，干扰点的位置受流场本身分辨率的影响较大，密网格下一般会获得更为精准的结果；而对于异侧激波干扰，其干扰点的识别结果通常比较准确，且网格依赖性不强。

3) 该算法自动化程度较高，基本不需要人工干预，且适用于任意网格类型，对复杂非定常流场中的运动激波辨识也具有良好的可靠性和准确性。

 [1] BUNING P, STEGER J. Graphics and flow visualization in computational fluid dynamics[C]//7th Computational Physics Conference, 1985: 1507. [2] DARMOFAL D L. Hierarchal visualization of three-dimensional vortical flow calculations[R]. Cambridge: Massachusetts Institute of Technology Cambridge Computational Fluid Dynamics Lab, 1991: 42-43. [3] LIOU S P, SINGH A, MEHLIG S, et al. An image analysis based approach to shock identification in CFD[C]//33rd Aerospace Sciences Meeting and Exhibit, 1995: 117. [4] LOVELY D, HAIMES R. Shock detection from computational fluid dynamics results[C]//14th Computational Fluid Dynamics Conference, 1999: 3285. [5] PAGENDARM H G, SEITZ B. An algorithm for detection and visualization of discontinuities in scientific data fields applied to flow data with shock waves[J]. Scientific Visualization:Advanced Software Techniques, 1993, 161-177. Click to display the text [6] MA K L, VAN ROSENDALE J, VERMEER W. 3D shock wave visualization on unstructured grids[C]//Proceedings of 1996 Symposium on Volume Visualization, 1996: 87-94. [7] VAN ROSENDALE J. Floating shock fitting via Lagrangian adaptive meshes[R]. 1994: 89-94. [8] KANAMORI M, SUZUKI K. Shock wave detection in two-dimensional flow based on the theory of characteristics from CFD data[J]. Journal of Computational Physics, 2011, 230(8): 3085-3092. Click to display the text [9] KANAMORI M, SUZUKI K. Three-dimensional shock wave detection based on the theory of characteristics[J]. AIAA Journal, 2013, 51(9): 2126-2132. Click to display the text [10] AKHLAGHI H, DALIRI A, SOLTANI M R. Shock-wave-detection technique for high-speed rarefied-gas flows[J]. AIAA Journal, 2017, 55(11): 3747-3756. Click to display the text [11] MONFORT M, LUCIANI T, KOMPERDA J, et al. A deep learning approach to identifying shock locations in turbulent combustion tensor fields[M]//Modeling, Analysis, and Visualization of Anisotropy, 2017: 375-392. [12] LIU Y, LU Y T, WANG Y Q, et al. A CNN-based shock detection method in flow visualization[J]. Computers & Fluids, 2019, 184: 1-9. Click to display the text [13] 马千里, 李思昆, 曾亮. 基于两级采样的非结构化网格流场多激波特征可视化方法[J]. 计算机研究与发展, 2012, 49(7): 1450-1459. MA Q L, LI S K, ZENG L. Visualization of multi-shock features for unstructured-grid flows based on two-level sampling[J]. Journal of Computer Research and Development, 2012, 49(7): 1450-1459. (in Chinese) Cited By in Cnki | Click to display the text [14] WU Z N, XU Y Z, WANG W B, et al. Review of shock wave detection method in CFD post-processing[J]. Chinese Journal of Aeronautics, 2013, 26(3): 501-513. Click to display the text [15] PACIORRI R, BONFIGLIOLI A. Recognition of shock-wave patterns from shock-capturing solutions[C]//Computational Modelling of Objects Represented in Images, 2012: 91-96. [16] PACIORRI R, BONFIGLIOLI A. A shock-fitting technique for 2D unstructured grids[J]. Computers & Fluids, 2009, 38(3): 715-726. Click to display the text [17] 刘君, 邹东阳, 董海波. 动态间断装配法模拟斜激波壁面反射[J]. 航空学报, 2016, 37(3): 836-846. LIU J, ZOU D Y, DONG H B. A moving discontinuity fitting technique to simulate shock waves impinged on a straight wall[J]. Acta Aeronautica et Astronautica Sinica, 2016, 37(3): 836-846. (in Chinese) Cited By in Cnki (8) | Click to display the text [18] ZOU D Y, XU C G, DONG H B, et al. A shock-fitting technique for cell-centered finite volume methods on unstructured dynamic meshes[J]. Journal of Computational Physics, 2017, 345: 866-882. Click to display the text [19] CHANG S Y, BAI X Z, ZOU D Y, et al. An adaptive discontinuity fitting technique on unstructured dynamic grids[J]. Shock Waves, 2019, 29(8): 1103-1115. Click to display the text [20] 吴夙慧, 成颖, 郑彦宁, 等. K-means算法研究综述[J]. 数据分析与知识发现, 2011, 27(5): 28-35. WU S H, CHENG Y, ZHENG Y N, et al. Survey on K-means algorithm[J]. New Technology of Library and Information Service, 2011, 27(5): 28-35. (in Chinese) Cited By in Cnki (324) | Click to display the text [21] MACQUEEN J. Some methods for classification and analysis of multivariate observations[C]//Proceedings of the Fifth Berkeley Symposium on Mathematical Statistics and Probability, 1967, 1(14): 281-297. [22] BÉZIER P E. Numerical control-mathematics and applications[M]. 1972: 256. [23] 刘君, 徐春光, 白晓征. 有限体积法和非结构动网格[M]. 北京: 科学出版社, 2016: 42-72. LIU J, XU C G, BAI X Z. Finite volume methods and unstructured dynamic grids technique[M]. Beijing: Science Press, 2016: 42-72. (in Chinese) [24] WOODWARD P, COLELLA P. The numerical simulation of two-dimensional fluid flow with strong shocks[J]. Journal of Computational Physics, 1984, 54(1): 115-173. Click to display the text [25] ALAUZET F, LOSEILLE A. A decade of progress on anisotropic mesh adaptation for computational fluid dynamics[J]. Computer-Aided Design, 2016, 72: 13-39. Click to display the text [26] KALLINDERIS Y, LYMPEROPOULOU E M, ANTONELLIS P. Flow feature detection for grid adaptation and flow visualization[J]. Journal of Computational Physics, 2017, 341: 182-207. Click to display the text
http://dx.doi.org/10.7527/S1000-6893.2019.23626

0

#### 文章信息

CHANG Siyuan, BAI Xiaozheng, LIU Jun

A two-dimensional shock wave pattern recognition algorithm based on cluster analysis

Acta Aeronautica et Astronautica Sinica, 2020, 41(8): 123626.
http://dx.doi.org/10.7527/S1000-6893.2019.23626