微槽湍流减阻模拟的三维稳态数值模拟方法

本发明涉及流体动力学领域，具体涉及微槽湍流减阻模拟的三维稳态数值模拟方法。

背景技术：

1、目前，公知的微槽湍流减阻数值模拟方法主要分为两类。

2、一类是基于全尺寸计算域的湍流减阻模拟方法。该类方法为使微槽表面进入完全湍流状态，需在计算域前端设置足够长度的湍流过渡段，然而微槽几何尺寸（通常为微米级）远小于过渡段长度（通常为毫米至厘米级，甚至米量级），二者尺度相差数个数量级，导致几何建模与网格划分存在极大困难。例如在comsol仿真软件中进行几何建模时，软件易因尺度差异自动忽略微槽细节，致使无法开展微槽减阻模拟；即便通过精细设置保留微槽几何特征，微槽区域所需的高密度网格与过渡区域的大尺度网格之间存在悬殊差异，又会引发整体网格数量剧增、计算量过于庞大的问题，严重制约仿真效率。

3、另一类是基于微尺度计算域的减阻模拟方法，为当前主流技术路线。该类方法通过截取包含若干微槽的局部区域作为计算域，并结合周期性边界条件来模拟无限大流场中的微槽减阻特性，可显著降低计算规模。然而，现有微尺度计算域方法均针对规则截面微槽（如v形、梯形、矩形、半圆形等）展开，通过展向阵列排布与流向周期边界条件实现模拟，其中流向周期性边界的适用前提是计算域入口面与出口面的几何形貌完全一致。上述条件仅对规则几何截面的理想微槽有效，对于激光刻蚀制备的实际微槽表面则难以适用。受激光加工过程中热力耦合作用的影响，实际刻蚀得到的微槽底部在展向和流向上并非严格平直一致，微槽截面形状沿程亦存在一定偏差，导致微尺度计算域的入口面与出口面几何形貌不完全相同，进而使得传统基于周期边界条件的微尺度计算域模拟方法失效。此外，目前该类微尺度减阻模拟大多借助fluent流体分析软件完成，在comsol平台下针对非规则微槽的高效高精度建模方法尚属空白。

4、综上，现有微槽减阻模拟技术存在以下技术缺陷：全尺寸计算域方法存在尺度差异悬殊、建模困难、计算量巨大的问题；常规微尺度计算域方法仅适用于规则截面微槽，无法处理激光刻蚀等工艺制备的非规则形貌微槽表面，适用范围受限。

技术实现思路

1、为了克服现有微槽湍流减阻模拟技术中，全尺寸计算域存在尺度悬殊导致建模困难与计算效率低下，以及常规微尺度计算域仅适用于规则截面微槽、无法处理激光刻蚀非规则形貌微槽表面的不足，本发明提供微槽湍流减阻模拟的三维稳态数值模拟方法，该方法通过构建适用于非规则微槽形貌的微尺度计算域，采用充分发展流动入口与压力出口替代传统周期边界条件，结合comsol物理场控制网格自动划分技术与密度修正的精度验证方法，该修正模式基于普朗特平板湍流边界层经验公式建立微尺度计算域与真机特征尺度之间的等效关联，以空气密度为修正参量，通过初次仿真求光滑壁面剪应力、真机尺寸特征参数下理论计算求光滑壁面剪应力、密度修正后仿真重算的三步流程完成模型校准，实现对激光刻蚀等非规则微槽表面减阻特性的高效、高精度数值模拟，为复杂形貌微槽的减阻效果评估与机理研究提供可靠的数字化分析工具。

2、为了实现上述目的，本发明提供如下技术方案：微槽湍流减阻模拟的三维稳态数值模拟方法，包括以下步骤：

3、步骤一：基于激光刻蚀非规则微槽几何数据以及约束设计，采用充分发展入口与压力出口的微尺度模拟流向周期性边界条件，以构建三维计算域几何模型；

4、步骤二：物理场控制网格通过识别三维计算域几何模型的流体域体积、各壁面边界，结合sst湍流模型和低雷诺数壁面处理的物理场约束进行网格划分，自动生成体网格与满足无量纲壁面距离y+ < 1的近壁边界层网格的离散网格；

5、步骤三：完成网格划分后，进入先后两次调用的cfd稳态求解，两次调用均采用同一套cfd稳态求解设置；

6、第一次初次仿真调用，将体网格与满足无量纲壁面距离y+ < 1的近壁边界层网格的离散网格、初始未修正的流体物性，以及sst湍流模型、离散格式、两步amg-gmres求解器配置和收敛判据作为输入进行求解，求解完成后输出初次仿真的流场数据，并从中提取光滑壁面剪应力，将光滑壁面剪应力作为密度校准计算的输入参数；得到修正后的空气密度，触发第二次调用，即校准后重算；

7、第二次调用直接复用步骤二输出的离散网格，仅把流体物性替换为修正后空气密度，重新开展cfd稳态求解，输出校准完成后的流场以及光滑壁面剪应力；

8、步骤四：基于第一次初次仿真求解得到的光滑壁面剪应力、真机机翼特征弦长、来流速度以及普朗特平板湍流经验公式，计算真机尺度雷诺数、理论摩擦阻力系数、理论剪应力、修正系数，得到修正后空气密度；基于修正后的空气密度，触发cfd重求解，即步骤三的第二次调用；cfd模型重求解完成后，得到校准后的摩擦阻力系数和剪应力。

9、优选的，步骤一的计算域内同时布置光滑表面与微槽表面，微槽表面采用非规则实际或模拟形貌，无需基于规则截面阵列排布，不强制微槽底部共线，不要求计算域入口面与出口面几何形貌一致；入口边界采用充分发展流动条件，出口采用压力出口条件，展向两侧为对称边界，壁面为无滑移壁面条件；计算域设置10道微槽，计算域高度不小于最大微槽深度的10倍，计算域长度为100μm。

10、优选的，步骤二具体为：定义三维计算域的流体域体积为ω，三维计算域流体域体积ω的所有壁面边界集合={光滑壁面，微槽壁面，入口面，出口面， 展向侧面}；采用sst湍流模型与低雷诺数壁面处理约束（comsol中“湍流，sst”设置中的壁处理类型选择“低雷诺数”），对流体域ω进行物理场控制网格划分（comsol中网格设置下的序列类型选为“物理场控制网格”，其单元大小选择“常规”），自动生成体网格单元集合t与近壁边界层网格单元集合b；

11、近壁第一层网格的无量纲壁面距离满足约束：

12、；

13、式中，yp为近壁第一层网格形心到光滑壁面或微槽壁面的2倍法向距离，ρ0为流体密度，即初始空气密度，μ为流体动力粘度；uτ为壁面摩擦速度，壁面摩擦速度uτ表示为：

14、；

15、式中，τ2为壁面剪应力，表示为：

16、；

17、式中，u0为入口平均速度，理论摩擦阻力系数cf，theo，表示为：

18、；

19、式中，rel为真机尺度雷诺数，表示为：

20、；

21、式中，l为真机机翼特征弦长；

22、优选的，步骤三中的cfd稳态求解设置如下：

23、设置流动为稳态不可压缩流动，湍流模型基于sst模型通过旋转-曲率校正，适配近壁曲率流动；纳维-斯托克斯方程与湍流输运方程的对流项采用流线扩散与侧风扩散格式（comsol内置的自洽稳定性选项“流线扩散”、“侧风扩散”），配合动态亚格子时间尺度（comsol内置的数值求解选项，针对微槽流场近壁区与主流区特征时间尺度差异，依据局部流场湍流特征自适应分配迭代步长），在保障近壁壁面剪应力求解精度的同时提升稳态迭代收敛效率；

24、双次调用cfd稳态求解，两次调用采用同一套cfd稳态求解配置，求解域、网格、湍流模型、离散格式、求解器、收敛判据均保持一致，仅流体密度参数替换为修正值；

25、分两步迭代求解：

26、输入参数集合包括体网格单元集合t，近壁边界层网格单元集合b，初始空气密度ρ，流体动力粘度μ，sst模型，离散格式，amg-gmres求解器，以及收敛残差阈值；

27、第一步为壁距离场初始化迭代：以体网格单元集合t、近壁边界层网格单元集合b、流体物性参数、sst湍流模型参数、离散格式为输入，采用代数多重网格法amg预处理器结合广义最小残差法gmres迭代求解器，求解壁距离输运方程，生成全流场壁距离场，为低雷诺数壁面处理提供近壁距离基准；迭代过程监测壁距离变量的相对残差，达到收敛残差阈值1×10-6时判定该步迭代收敛；

28、第二步为流动-湍流耦合稳态迭代：复用第一步得到的壁距离场，以对应调用阶段的空气密度、动力粘度μ为流体物性输入，采用amg与gmres耦合迭代求解稳态不可压缩纳维-斯托克斯方程与sst湍流输运方程组，迭代过程同步更新速度场、压力场、湍动能场、耗散率场及壁面剪应力；迭代过程监测壁面阻力以及各流动变量和湍流变量的相对残差，达到收敛残差阈值0.001或壁面阻力趋于稳定时判定该步迭代收敛；

29、第一次初次仿真调用，以初始空气密度ρ0作为流体物性输入完成上述两步迭代后，输出初次仿真流场数据，提取光滑壁面的平均壁面剪应力τ1；

30、第二次校准重算调用，复用网格集合{t，b}，仅将流体物性中的空气密度替换为修正后空气密度ρ′，其余输入参数与第一次调用完全一致，重复上述两步迭代后，输出校准后的速度场、压力场以及光滑壁面平均剪应力τnum。

31、优选的，基于普朗特平板理论的密度修正与精度验证，修正与验证流程如下依次执行：

32、(1) 计算真机尺度雷诺数：

33、；

34、(2) 由普朗特平板湍流边界层摩擦阻力系数经验公式，计算理论摩擦阻力系数：

35、；

36、(3) 由壁面剪应力公式计算理论壁面剪应力：

37、；

38、(4) 定义密度修正系数k为仿真剪应力与理论剪应力的比值：

39、；

40、(5) 计算修正后的空气密度：

41、；

42、(6) 将ρ′代入步骤三的第二次调用，得到修正后的光滑壁面平均剪应力τnum，进而计算模拟摩擦阻力系数：

43、；

44、(7) 计算模拟值与理论值的相对误差：

45、；

46、当相对误差δ ≤ δset时，判定高精度模型建立完成，δset为预设误差阈值。

47、本发明具有如下优点：

48、本发明通过构建适用于非规则形貌的微尺度计算域并采用充分发展入口与压力出口替代传统周期边界条件，突破了常规微尺度计算域方法仅适用于规则截面微槽的限制，实现了对激光刻蚀等工艺制备的实际非规则微槽表面的减阻特性模拟。与全尺寸计算域方法相比，本发明采用百微米级微尺度计算域，从根本上解决了微槽尺寸与湍流过渡段尺度悬殊导致的建模困难与计算量巨大的问题，显著提升了仿真效率。在网格划分方面，利用comsol物理场控制网格自动生成满足低雷诺数要求的高质量边界层网格，无需手动配置边界层参数，简化了前处理流程，同时保证了近壁区求解精度。此外，本发明提出的基于密度修正的两步法精度验证策略，通过普朗特经验公式校准模型，有效保证了微尺度计算域下壁面摩擦阻力计算的准确性，使得微尺度计算域的模拟结果能够等效反映真实机翼尺度的流动特性。

49、综上，本发明为复杂形貌微槽的减阻效果评估与机理研究提供了一种高效、高精度且适用性广的数值模拟方法，对推动微槽减阻技术在航空等领域的工程化应用具有重要意义。