基于分区径向基函数的Level Set界面捕捉研究*

金　秋　李廷秋　林　超

(武汉理工大学交通学院　武汉　430063)



基于分区径向基函数的Level Set界面捕捉研究*

金秋李廷秋林超

(武汉理工大学交通学院武汉430063)

针对传统的Level Set方法的固有缺陷，通过引入径向基函数(radial basis function，RBF)，提出一种分区捕捉动界面的Level Set数值方法，将Level Set偏微分控制方程转化为全微分方程，避免了重新初始化过程，可精确表达界面形状并维持质量守恒.采用径向基函数中的基函数和权重系数来重构符号距离函数，并结合区域划分技术提高计算效率.3个计算算例表明，基于分区径向基函数Level-Set方法，可明显提高界面几何形状的精度和捕捉效率，具有较好的质量守恒性.

水平集；径向基函数；运动界面；分区技术；质量守恒

0　引　　言

动界面捕捉技术研究是模拟海工结构物在复杂海洋环境中运动的必要前提.虽然当前传统Level Set方法已逐步运用于自由液面及动物面捕捉的数值模拟中[1-2]，但因其空间项离散数值误差导致的质量不守恒问题，极大的限制了其应用范围.为了解决传统Level Set方法质量不守恒问题，Cai等[3]构造了一种五阶紧凑CWENO格式，从而提高空间项离散精度，Lsakov等[4]提出基于求解偏微分方程的快速局部Level set算法，通过减少符号距离重新初始化的迭代次数，提高Level Set算法精度.Jemison等[5]提出改进CLSVOF方法，运用Level Set求解相界面法向方向等物理量，结合VOF方法求解流输运问题以保证质量守恒.Enright等[6]提出了粒子Level Set方法(PLSM)，在计算域内设定拉格朗日标记粒子，根据标志粒子修正符号距离函数.

文中提出一种基于分区径向基函数(radial basis function，RBF)的Level Set方法，根据界面实际位置对全局计算域进行分区，突破有限差分求解框架，实现符号距离函数的空间和时间项分离表达，将偏微分方程转化为全微分方程进行求解，有效限制传统方法中空间项离散带来的数值误差，提高计算精度的同时有效维持动界面的质量守恒.

1　理论基础

1.1符号距离函数和Level Set方程

在Level Set方法中，将移动的界面作为零等值面嵌入高一维符号距离函数中，隐式地表达或追踪.构造函数Φ(x,t)，使得在任意时刻，运动界面Γ为Φ(x,t)的零等值面，即

(1)

定义Φ(x,t)初值为t=0时刻，x到界面的符号距离，即

(2)

在已知初始时刻界面位置后，通过对时间求全导数，确定描述界面运动的偏微分形式的控制方程

(3)

式中：v为目标点处的速度，由界面处的法向速度Fn和单位外法向向量Φ表示，则上式可写为

(4)

将界面Γ的演化嵌入至函数Φ中，通过对Φ的求解，可自动获得界面Γ的拓扑结构变化，如界面分裂或界面融合.

传统LevelSet方法中，运用符号距离函数赋初值并用迎风格式的FDM求解运动控制方程，该局部近似的方法产生的数值误差将逐步累积最终影响求解结果，因此需要周期性的重新初始化过程来保证数值求解的稳定性.

1.2径向基函数

给定函数Φ：R+→R，对于数据{Xi,fi}∈Rd⊗R,(i=1,2,…,n)，表示为如下径向基函数的形式

(5)

式中：f(X)为待求物理量(如空间点到界面的距离)；X为计算域中任意空间点；Xi为插值点；λ为权重系数；n为插值点个数；P(X)为一阶多项式.

对于二维问题，其表达式为P(X)=λn+1+λn+2x+λn+3y，φ(‖·‖)为径向基基函数，文中选用MQ基函数

(6)

式中:c为形状参数，取0～1间的常数.

为保证RBF插值有惟一解，权重系数须满足正交性质

(7)

n个插值点组成的径向基函数组写成矩阵形式为

(8)

Aij=φj(xi),i,j=1,2,…,N

Pij=pj(xi),i=1,2,…,N;j=1,2,3.

1.3基于RBF求解Level Set方程

基于径向基函数的散乱点插值功能，对LevelSet函数符号距离函数进行重新描述.在表达形式上实现时间变量和空间变量的分离，当插值点固定不变时，空间项即为已知项，函数中仅有RBF权重系数与时间相关.

(9)

式中：ψ(x)=[φ1(x)…φN(x)1xy]

因此，LevelSet运动方程式(4)转变为

(10)

LevelSet运动偏微分方程由此转变为关于时间的全微分方程，通过更新RBF权重系数，实现更新符号距离函数的目标，继而实现动界面的隐式表达.

结合式(8)、(10)可以改写为

(11)

式中：

采用一阶差分方法对时间项进行离散，即可得到下一时刻权重系数值.

(12)

式中：Δt为时间步长.

1.4区域划分插值技术

由1.2可知，当插值点个数为n时，每一个时间步迭代需要求解一个(n+3)×(n+3)的满秩矩阵.根据实际计算工况对全局计算域划分区域，可减少遍历全局背景网格的区域，减少插值点数量，显著提高矩阵求解速度.同时，在界面区域划分更细的局域网格，有利于更精确的捕捉界面拓扑结构变化.局域区域网格划分示意图见图1.

图1　区域划分示意图

定义一个接近零的足够小数正值δ(一般取为1.5～2倍均匀网格宽度)，根据函数值Φ(x,y)进行判断，在计算域内判断Φ(x,y)与δ的关系，当Φ(x,y)≤δ时记录所有的点集，分别求得xmin{(x,y)|Φ(x,y)≤δ},xmax{(x,y)|Φ(x,y)≤δ},ymin{(x,y)|Φ(x,y)≤δ},ymax{(x,y)|Φ(x,y)≤δ},定位界面位置，保证所划分的局部插值域包含整个界面.引入dx,dy两个安全范围，满足dx>max(ux×Δt),dy>max(vy×Δt)，以此确定局部插值域的四个角点坐标为：(xmin-dx，ymin-dy)，(xmax+dx，ymin-dy)，(xmax+dx，ymax+dy)和(xmin-dx，ymax+dy).由于在计算过程中时间步长根据CFL条件选择，在均匀网格内一般选取dx,dy为2-3倍的网格宽度.

2　计算结果及分析

2.1分区RBF-Level Set与传统Level Set数值实验对比

用上述2种方法分别对典型圆形界面在剪切流中的变形现象进行模拟计算，在相同的计算条件下，将文中提出的分区RBF-Level Set方法与传统Level Set方法进行对比，以验证新方法对质量守恒性的改善效果.

算例1剪切流场.

速度场为

(13)

由图2a)可知，设定计算域为(0,0)×(1,1)，初始时刻圆心位置(0.35,0.5)，半径R=0.25 m，剪切流动旋转中心(0.5,0.5).

图2　分区RBF-Level Set与传统Level Set剪切流场结果对比

本案例中计算了剪切流中圆形界面在不同时刻的运动变形情况.对比可知，分区RBF-Level Set方法能更精确捕捉运动界面的形状与位置.统计不同时间点RBF-Level Set方法和传统Level Set方法界面内质量的相对误差，分别为0.15%/1.05%，0.35%/1.9%和0.75%/2.5%，明显提高了计算结果的精确度.对比验证结果表明，基于径向基函数的Level Set方法可以更好地维持界面形状同时保证质量守恒特性.

2.2分区RBF-Level Set与全局RBF-Level Set方法对比

算例2剪切流场，流场信息同算例2.

图3　分区与全局RBF-Level Set剪切流场结果对比

由图3可知，2种计算方法均能较好的模拟剪切流中界面变形，计算结果基本完全吻合，但在界面变形的细节处，如界面尾端尖点处，分区RBF-LevelSet方法的捕捉更为精确细致.全局方法中每一步的计算遍历整个计算域，计入全局速度场对界面的影响，但消耗了大量计算时间，且单个时间步内远离界面处的速度对界面变形影响较小.文中提出的分区思想随界面移动更新局部区域，划分局部细网格，将插值点集中于界面周围，节省时间的同时更好的模拟界面细节处的变化.

3　分区RBF-Level Set方法界面强迫运动模拟

算例3船球鼻首剖面的横摇运动.

速度场

(14)

设定计算域为计算域为(0,0)×(2,2)，初始时刻球鼻艏处于平衡位置，图3分别表示球鼻艏运动1/4,1/2,3/4,1个周期位置.

图4　RBF-Level球鼻艏横摇运动

利用文中所提出的分区RBF-Level Set方法对船体曲面曲率变化较为复杂的球艏剖面的强迫横摇运动进行捕捉，由图中可以看出，RBF-Level Set界面捕捉结果与模型初始时刻形状十分吻合，尖锐处仍能保持较高分辨率.统计横摇各时刻界面内质量误差，分别为0.08%，0.06%，0.1%，计算结果表明，文中所提出的RBF-Level Set方法也适用于曲率变化较大的界面运动捕捉.

4　结 束 语

提出了一种基于分区RBF-Level Set方程求解方法.数值实验结果表明，运用RBF插值技术分区求解Level Set方程继承了传统方法精确表达运动界面几何形状优点的同时，由于避免了重新初始化过程，明显改善了质量守恒性，具有更好的数值稳定性和精度.

[1]CHENY Y, BOTELLA O. The LS-STAG method: a new immersed boundary/level-set method for the computation of incompressible viscous flows in complex moving geometries with good conservation properties[J]. Computational Physics,2010,229(4):1043-1076.

[2]KRAUSE R, MOHR C. Level set based multi-scale methods for large deformation contact problems[J]. Applied Numerical Mathematics，2011，61:428-442.

[3]CAI L, FENG J H, XIE W X. A CWENO-type central-upwind scheme for ideal MHD equations [J] .Applied Mathmatics and Computation,2005(1):600-612.

[4]ISAKOV V, LEUNG S, QIAN J. A fast local level set method for inverse gravimetry[J]. Communications in Computational Physics,2011,10(4):1044-1070.

[5]JEMISON M, LOCH E, SUSSMAN M,et al. A coupled level set-moment of fluid method for incompressible two-phase flows [J]. Journal of Scientific Computing, 2013,54(3):454-491.

[6]ENRIGH T D, FEDKIW R, FERZIGER J, et al. A hybrid particle level set method for improved interface capturing [J]. Journal of Computational Physics, 2014,183(2):83-116.

Regional RBF based Level Set Method for Interface Capturing

JIN QiuLI TingqiuLIN Chao

(SchoolofTransportation,WuhanUniversityofTechnology,Wuhan430063,China)

The renewal of signed distance function usually relies on reinitialization in conventional Level Set method. An alternative method is proposed in this paper to regionally capture the moving interfaces via the introduction of the Radial Basis Functions (RBF). This approach provides an accuracy and conserved representation of the implicit function, location and shape of interfaces. The signed distance function is rebuilt by based functions and weight coefficients. Its deformation is then considered as an updating of the RBF interpolants, turning the partial differential equation (PDE) problem into an ordinary differential equation (ODE) problem. Reinitialization is found to be no longer necessary and hence it avoids large numerical errors occurring in the conventional method. Three numerical examples show that the proposed regional RBF-Based Level Set method can significantly improve both the accuracy and conservation in capturing moving interfaces.

level set; RBF; moving interfaces; regional technology; mass conservation

2016-07-05

U661.1

10.3963/j.issn.2095-3844.2016.04.013

金秋(1992- )：女，硕士生，主要研究领域为船舶流体力学

*国家自然科学基金项目资助(51579196、51139005)