CN108717493A - 一种面向结构化四边网格生成的二维区域自动分解方法 - Google Patents

Abstract

本发明公开了一种面向结构化四边网格生成的二维区域自动分解方法。当前，针对复杂几何外形的结构化网格生成，主要依赖人工将复杂区域手动分解为规则的子区域，增加了出错的可能性。本发明步骤：建立分析实体对象的二维几何模型，并离散成三角形网格模型；根据几何边界约束建立二维几何模型内部区域光滑矢量场的控制偏微分方程，并采用边界元数值方法进行求解；建立与光滑矢量场对应的光滑标架场，并通过分析光滑标架场的奇异结构实现对二维几何模型的自动分解。本发明可以自动快速地获得任意复杂二维几何模型的高质量分区结果，为后续多块结构化四边形网格生成提供输入。

Description

一种面向结构化四边网格生成的二维区域自动分解方法

技术领域

本发明涉及数值模拟领域前处理中的分块结构化四边形网格生成的区域分解过程，具体涉及一种面向结构化四边网格生成的二维区域自动分解方法。

背景技术

网格生成是有限元法、有限体积法和有限差分法等数值模拟技术中的前处理过程，该过程将连续的几何区域剖分成有限个基本几何形体的组合，这些基本几何形体被称为网格单元。常用的单元类型有三角形单元、四边形单元、四面体单元和六面体单元等；根据网格单元的拓扑结构是否规则可将网格分为结构化网格和非结构化网格。在数值模拟中，结构化网格由于具有比非结构化网格更高的求解精度和更少的求解自由度而备受青睐。

相对其它类型的网格，结构网格具有很强的拓扑约束。因而，针对复杂几何模型，生成一套高质量的结构网格是一项非常费时且困难的工作。通过将模型区域分解为多个子域，并在各子域中生成结构化网格，而各子域之间的拓扑不要求保持结构化特性，将得到分块结构化网格，这将在一定程度上降低网格生成难度，也能避免产生完全结构化网格中可能出现的扭曲单元。正因如此，在实际工程中，分块结构化网格比完全结构化网格应用范围更广。

对于复杂几何模型，对其进行区域分解具有和生成其结构化网格一样的难度，自动、可控、高质量的区域分解至今仍是一个开放性难题。当前，针对复杂几何外形的结构化网格生成，主要依赖人工将复杂区域手动分解为规则的子区域，再在各子区域里生成规则的四边形网格。人工环节的引入一方面要求工程人员具有较强的工程经验，另一方面非常耗时且增加了出错的可能性。这些缺点制约了自动化分块结构化网格生成工具的诞生。

发明内容

为避免手工区域分解存在的缺点，本发明提供一种面向结构化四边网格生成的二维区域自动分解方法，能够针对任意复杂二维几何模型快速自动获得高质量的区域分解结果，从而提升分块结构化网格生成效率和分块结构化网格质量。

本发明采用的技术方案是：建立分析实体对象的二维几何模型，并离散成三角形网格模型；根据几何边界约束建立二维几何模型内部区域光滑矢量场的控制偏微分方程，并采用边界元数值方法进行求解；建立与光滑矢量场对应的光滑标架场，并通过分析光滑标架场的奇异结构实现对二维几何模型的自动分解，具体包括以下步骤：

步骤1、建立分析实体对象的二维几何模型，在二维几何模型区域上定义矢量场与标架场之间的映射关系，具体如下：

建立分析实体对象的二维几何模型，将二维几何模型区域Ω离散成三角形网格T_h，并在T_h上定义矢量场和标架场，分别由u(p₁)和v_k(p₁)表示，0≤k≤3，v₀(p₁)垂直于v₁(p₁)且与v₂(p₁)反向，v₃(p₁)与v₁(p₁)反向，p₁为网格顶点；记θ_d(p₁)为x轴沿逆时针旋转到与u(p₁)方向相同时所转过的最小角度；θ_c,k(p₁)为x轴沿逆时针旋转到达与v_k(p₁)方向相同时所转过的最小角度，其中θ_c,k(p₁)，0≤k≤3为四个互不相等的角度值，矢量场和标架场的映射关系用式(1)表示：

θ_d(p₁)＝4θ,θ＝min{θ_c,0(p₁),θ_c,1(p₁),θ_c,2(p₁),θ_c,3(p₁)} (1)

步骤2、在二维几何模型区域上建立光滑矢量场最优化模型，具体如下：

以二维几何模型区域边界上矢量向区域内推进时满足最小变化为目标，以二维几何模型区域边界的几何特征约束为约束条件建立光滑矢量场最优化模型，即目标函数为约束条件为：

其中，为u_i(p)的拉普拉斯算子，p为二维几何模型区域内的点，i为点p的坐标分量x或y，u_i(p)为u(p)在x或y方向的分量，和分别为狄利克雷边界条件和诺依曼边界条件，n为边界上的点的外法向量；u_i,0(p)为二维几何模型区域边界上的矢量沿i方向的分量，为u_i,0(p)的法向导数。光滑矢量场最优化模型的最优解等价为如下拉普拉斯边值问题的解：

式(3)中，Δu_i(p)为u_i(p)的拉普拉斯算子的另一种表达形式。

为使二维几何模型区域分解后的区域分割线垂直于区域分解后的边界，三角形网格边界点的标架满足如下条件：

式中，p₀为三角形网格边界上的网格顶点，θ_h(p₀)和n(p₀)分别为三角形网格边界上以p₀点为交点的两条边朝向区域内的夹角值和该夹角平分线的单位向量，t(p₀)为垂直于n(p₀)的单位向量，v_k(p₀)为p₀处的标架，0≤k≤3。

通过公式(4)求得三角形网格边界上的标架信息后，通过步骤1定义的矢量场与标架场之间的映射关系求解三角形网格边界上的矢量边界条件，其中，三角形网格边界矢量的模为1。

步骤3、求解矢量拉普拉斯方程

根据公式(3)，二维几何模型区域内部点p₂处的矢量沿i方向的分量通过如下边界积分公式求解：

其中，二维几何模型区域边界上的点ds_p表示二维几何模型区域边界上p点邻域的微分，G(p₃,p₂)为Laplace方程的基本解，即 n3为p₃的外法向量。由公式(3)，在二维几何模型区域边界上每个点的u_i(p₃)和中有一个值是已知的，并且另一个值通过式(6)得到：

其中，二维几何模型区域边界上的点且p₄与p₃为不同点，由于在点p₃处光滑，则n4为p₄的外法向量。采用边界元数值方法求解公式(6)，将划分成N段，记为k1在1～N中取值，N＞2，得到如下离散形式：

其中，q_k1为位于上的点，q_j为位于上的点，nq_j为q_j的外法向量。公式(7)表示一个由N个未知数和N个方程组成的方程组，通过求解该方程组，得到位于二维几何模型区域边界点上的边值u_i(p)或随后，通过公式(5)求解位于二维几何模型区域内的内部点处的矢量，最终获得整个二维几何模型区域内的矢量场。

步骤4、建立光滑标架场并利用其实现区域分解

通过步骤1定义的矢量场与标架场之间的映射关系，求得步骤3中获得的矢量场对应的标架场。通过分析该标架场获得奇异点的位置，并以该标架场为输入数据，利用龙格-库塔法寻找以奇异点为起始点和终止点的流线。最后，以这些流线将二维几何模型区域分解为多个四边结构的子区域。

本发明的有益效果是，本发明通过建立覆盖几何模型区域的矢量场物理模型，并采用边界元法进行数值求解，随后建立与之一一映射的标架场，有效地增强了获得标架场的精度和速度，最终通过分析标架场的奇异结构，将二维几何模型区域精确且自动分解成多个四边形结构的子区域。该方法实现了任意复杂二维几何模型区域的自动快速分解，提高了分块结构化网格生成效率。

本发明面向分块结构化网格生成需求，提出一种针对任意复杂二维几何模型的自动区域分解方法，该发明可以分析二维飞机翼形外流场区域结构并将其自动分解为多个具有四边结构的子区域，并且该发明也适合其他复杂二维几何模型的区域分解。

附图说明

图1a和图1b分别为位于p₁点处的方向向量和标架；

图2a为本实施例中二维几何模型区域边界示意图；

图2b图为图2a中二维几何模型区域离散成三角形网格的示意图；

图2c为本发明方法针对图2a中二维几何模型区域形成的标架场边界图；

图2d为本发明方法针对图2a中二维几何模型区域形成的矢量场边界图；

图3a为本发明方法针对图2a中二维几何模型区域形成的矢量场图；

图3b为本发明方法针对图2a中二维几何模型区域形成的标架场图；

图3c为本发明方法针对图2a中二维几何模型区域形成的多个四边结构子区域分解图；

图4为本实施例中的二维飞机翼形外流场区域示意图；

图5为本发明方法针对图4中二维飞机翼形外流场区域形成的多个四边结构子区域分解图；

图6a为本发明方法针对图4中二维飞机翼形外流场区域形成的全局四边形网格；

图6b为图6a中A部分的局部放大图；

图6c为图6a中B部分的局部放大图。

具体实施方式

下面结合附图对本发明作进一步说明。

一种面向结构化四边网格生成的二维区域自动分解方法，具体步骤如下：

1、在二维几何模型区域上定义矢量场与标架场之间的映射关系

将如图2a所示的二维几何模型区域Ω，离散成三角形网格T_h(如图2b)，并在T_h上定义矢量场和标架场，分别由u(p₁)和v_k(p₁)，0≤k≤3表示，v₀(p₁)垂直于v₁(p₁)且与v₂(p₁)反向，v₃(p₁)与v₁(p₁)反向，p₁为网格顶点；如图1a和图1b所示。记θ_d(p₁)为x轴沿逆时针旋转到与u(p₁)方向相同时所转过的最小角度；θ_c,k(p₁)为x轴沿逆时针旋转到达与v_k(p₁)方向相同时所转过的最小角度，其中θ_c,k(p₁)，0≤k≤3为四个互不相等的角度值，矢量场和标架场的映射关系用式(1)表示：

θ_d(p₁)＝4θ,θ＝min{θ_c,0(p₁),θ_c,1(p₁),θ_c,2(p₁),θ_c,3(p₁)} (1)

2、在二维几何模型区域上建立光滑矢量场最优化模型

以二维几何模型区域边界上矢量向区域内推进时满足最小变化为目标，以二维几何模型区域边界的几何特征约束为约束条件建立光滑矢量场最优化模型，即目标函数为约束条件为：

其中，p为二维几何模型区域内(内部或边界)的点,i为点p的坐标分量x或y，u_i(p)为u(p)在x或y方向的分量，和分别为狄利克雷(Dirichlet)边界条件和诺依曼(Neumann)边界条件；u_i,0(p)为二维几何模型区域边界上的矢量沿i方向的分量，为u_i,0(p)的法向导数。光滑矢量场最优化模型的最优解等价为如下拉普拉斯(Laplace)边值问题的解：

为使二维几何模型区域分解后的区域分割线垂直于区域分解后的边界，三角形网格边界点的标架满足如下条件：

式中，p₀为三角形网格边界上的网格顶点，θ_b(p₀)和n(p₀)分别为三角形网格边界上以该点为交点的两条边朝向区域内的夹角值和该夹角平分线的单位向量，t(p₀)为垂直于n(p₀)的单位向量，v_k(p₀)为p₀处的标架，0≤k≤3。

通过公式(4)求得三角形网格边界上的标架信息(如图2c)后，通过步骤1定义的矢量场与标架场之间的映射关系求解三角形网格边界上的矢量边界条件，其中，三角形网格边界矢量的模为1，如图2d所示。

3、求解矢量拉普拉斯方程

根据公式(3)，二维几何模型区域内部点p₂(除边界以外的点)处的矢量分量通过如下边界积分公式求解：

其中，二维几何模型区域边界上的点ds_p表示二维几何模型区域边界上p点邻域的微分，G(p₃,p₂)为Laplace方程的基本解，即 由边界积分公式(5)可知，求解该方程的关键在于获得u_i(p₃)和而由公式(3)可知，在二维几何模型区域边界上每个点的u_i(p₃)和中有一个值是已知的，并且另一个值通过式(6)得到：

其中，二维几何模型区域边界上的点且p₄与p₃为不同点，由于在点p₃处光滑，则采用边界元数值方法求解公式(6)，将划分成N段(记为)，得到如下离散形式：

其中，q_k1(k1＝1～N)为位于上的点。可知，公式(7)表示一个由N个未知数和N个方程组成的方程组。通过求解该方程组，求解位于二维几何模型区域边界点上未知的边值u_i(p)或随后，通过公式(5)求解位于二维几何模型区域内的内部点处的矢量，最终获得整个二维几何模型区域内的矢量场(如图3a)。

4、建立光滑标架场并利用其实现区域分解

通过步骤1定义的矢量场与标架场之间的映射关系，求得步骤3中获得的矢量场对应的标架场(如图3b)。通过分析该标架场获得奇异点(矢量的模为0的点)的位置，并以该标架场为输入数据，利用龙格-库塔法(Runge-Kutta)寻找以奇异点为起始点和终止点的流线。最后，以这些流线将二维几何模型区域分解为多个四边结构的子区域(如图3c)。

实施例

图4为二维飞机翼形形状(包含中间的主翼及两端的襟翼共三段翼形结构)及其外流场计算域，图5展示了使用本发明对某二维翼形计算域的自动区域分解结果。本发明可以自动实现对该计算域的自动分解，得到多块高质量的四边结构子域，从而提升分块结构化网格(如图6a、6b和6c所示)生成效率和分块结构化网格质量。

Claims (1)

1.一种面向结构化四边网格生成的二维区域自动分解方法，其特征在于：该方法具体包括以下步骤：

步骤1、建立分析实体对象的二维几何模型，在二维几何模型区域上定义矢量场与标架场之间的映射关系，具体如下：

建立分析实体对象的二维几何模型，将二维几何模型区域Ω离散成三角形网格T_h，并在T_h上定义矢量场和标架场，分别由u(p₁)和v_k(p₁)表示，0≤k≤3，v₀(p₁)垂直于v₁(p₁)且与v₂(p₁)反向，v₃(p₁)与v₁(p₁)反向，p₁为网格顶点；记θ_d(p₁)为x轴沿逆时针旋转到与u(p₁)方向相同时所转过的最小角度；θ_c，k(p₁)为x轴沿逆时针旋转到达与v_k(p₁)方向相同时所转过的最小角度，其中θ_c，k(p₁)，0≤k≤3为四个互不相等的角度值，矢量场和标架场的映射关系用式(1)表示：

θ_d(p₁)＝4θ，θ＝min{θ_c，0(p₁)，θ_c，1(p₁)，θ_c，2(p₁)，θ_c，3(p₁)} (1)

步骤2、在二维几何模型区域上建立光滑矢量场最优化模型，具体如下：

以二维几何模型区域边界上矢量向区域内推进时满足最小变化为目标，以二维几何模型区域边界的几何特征约束为约束条件建立光滑矢量场最优化模型，即目标函数为约束条件为：

其中，为u_i(p)的拉普拉斯算子，p为二维几何模型区域内的点，i为点p的坐标分量x或y，u_i(p)为u(p)在x或y方向的分量，和分别为狄利克雷边界条件和诺依曼边界条件，n为边界上的点的外法向量；u_i，0(p)为二维几何模型区域边界上的矢量沿i方向的分量，为u_i，0(p)的法向导数；光滑矢量场最优化模型的最优解等价为如下拉普拉斯边值问题的解：

式(3)中，Δu_i(p)为u_i(p)的拉普拉斯算子的另一种表达形式；

为使二维几何模型区域分解后的区域分割线垂直于区域分解后的边界，三角形网格边界点的标架满足如下条件：

式中，p₀为三角形网格边界上的网格顶点，θ_b(p₀)和n(p₀)分别为三角形网格边界上以p₀点为交点的两条边朝向区域内的夹角值和该夹角平分线的单位向量，t(p₀)为垂直于n(p₀)的单位向量，v_k(p₀)为p₀处的标架，0≤k≤3；

通过公式(4)求得三角形网格边界上的标架信息后，通过步骤1定义的矢量场与标架场之间的映射关系求解三角形网格边界上的矢量边界条件，其中，三角形网格边界矢量的模为1；

步骤3、求解矢量拉普拉斯方程

根据公式(3)，二维几何模型区域内部点p₂处的矢量沿i方向的分量通过如下边界积分公式求解：

其中，二维几何模型区域边界上的点ds_p表示二维几何模型区域边界上p点邻域的微分，G(p₃，p₂)为Laplace方程的基本解，即 n3为p₃的外法向量；由公式(3)，在二维几何模型区域边界上每个点的u_i(p₃)和中有一个值是已知的，并且另一个值通过式(6)得到：

其中，二维几何模型区域边界上的点且p₄与p₃为不同点，由于在点p₃处光滑，则n4为p₄的外法向量；采用边界元数值方法求解公式(6)，将划分成N段，记为k1在1～N中取值，N＞2，得到如下离散形式：

其中，q_k1为位于上的点，q_j为位于上的点，nq_j为q_j的外法向量；公式(7)表示一个由N个未知数和N个方程组成的方程组，通过求解该方程组，得到位于二维几何模型区域边界点上的边值u_i(p)或随后，通过公式(5)求解位于二维几何模型区域内的内部点处的矢量，最终获得整个二维几何模型区域内的矢量场；

步骤4、建立光滑标架场并利用其实现区域分解

通过步骤1定义的矢量场与标架场之间的映射关系，求得步骤3中获得的矢量场对应的标架场；通过分析该标架场获得奇异点的位置，并以该标架场为输入数据，利用龙格-库塔法寻找以奇异点为起始点和终止点的流线；最后，以这些流线将二维几何模型区域分解为多个四边结构的子区域。