CN116894408A - 采用数字化对风洞试验模型堵塞度进行计算的方法 - Google Patents
采用数字化对风洞试验模型堵塞度进行计算的方法 Download PDFInfo
- Publication number
- CN116894408A CN116894408A CN202311161849.2A CN202311161849A CN116894408A CN 116894408 A CN116894408 A CN 116894408A CN 202311161849 A CN202311161849 A CN 202311161849A CN 116894408 A CN116894408 A CN 116894408A
- Authority
- CN
- China
- Prior art keywords
- test model
- grid
- wind tunnel
- projection
- windward
- Prior art date
- Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
- Granted
Links
- 238000012360 testing method Methods 0.000 title claims abstract description 171
- 238000000034 method Methods 0.000 title claims abstract description 43
- 230000000903 blocking effect Effects 0.000 title claims abstract description 40
- 239000007921 spray Substances 0.000 claims abstract description 4
- 238000004364 calculation method Methods 0.000 claims description 54
- 239000013598 vector Substances 0.000 claims description 27
- 230000014509 gene expression Effects 0.000 claims description 7
- 230000005484 gravity Effects 0.000 claims description 5
- 230000008569 process Effects 0.000 abstract description 5
- 239000004973 liquid crystal related substance Substances 0.000 abstract 2
- 230000036544 posture Effects 0.000 description 13
- 230000008859 change Effects 0.000 description 6
- 238000012937 correction Methods 0.000 description 6
- 230000000694 effects Effects 0.000 description 5
- 230000007246 mechanism Effects 0.000 description 5
- 230000005477 standard model Effects 0.000 description 5
- 238000013461 design Methods 0.000 description 4
- 238000012986 modification Methods 0.000 description 3
- 230000004048 modification Effects 0.000 description 3
- 238000012545 processing Methods 0.000 description 3
- 150000001875 compounds Chemical class 0.000 description 2
- 238000011960 computer-aided design Methods 0.000 description 2
- 238000010586 diagram Methods 0.000 description 2
- 238000013507 mapping Methods 0.000 description 2
- 238000005259 measurement Methods 0.000 description 2
- 230000003068 static effect Effects 0.000 description 2
- 238000013519 translation Methods 0.000 description 2
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000009530 blood pressure measurement Methods 0.000 description 1
- 238000006243 chemical reaction Methods 0.000 description 1
- 238000011161 development Methods 0.000 description 1
- 230000006872 improvement Effects 0.000 description 1
- 230000007774 longterm Effects 0.000 description 1
- 238000013139 quantization Methods 0.000 description 1
- 238000011160 research Methods 0.000 description 1
- 238000006467 substitution reaction Methods 0.000 description 1
- 230000009466 transformation Effects 0.000 description 1
Classifications
-
- G—PHYSICS
- G01—MEASURING; TESTING
- G01M—TESTING STATIC OR DYNAMIC BALANCE OF MACHINES OR STRUCTURES; TESTING OF STRUCTURES OR APPARATUS, NOT OTHERWISE PROVIDED FOR
- G01M9/00—Aerodynamic testing; Arrangements in or on wind tunnels
- G01M9/08—Aerodynamic models
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F17/00—Digital computing or data processing equipment or methods, specially adapted for specific functions
- G06F17/10—Complex mathematical operations
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F30/00—Computer-aided design [CAD]
- G06F30/20—Design optimisation, verification or simulation
- G06F30/28—Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2113/00—Details relating to the application field
- G06F2113/08—Fluids
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2119/00—Details relating to the type or aim of the analysis or the optimisation
- G06F2119/14—Force analysis or force optimisation, e.g. static or dynamic forces
-
- Y—GENERAL TAGGING OF NEW TECHNOLOGICAL DEVELOPMENTS; GENERAL TAGGING OF CROSS-SECTIONAL TECHNOLOGIES SPANNING OVER SEVERAL SECTIONS OF THE IPC; TECHNICAL SUBJECTS COVERED BY FORMER USPC CROSS-REFERENCE ART COLLECTIONS [XRACs] AND DIGESTS
- Y02—TECHNOLOGIES OR APPLICATIONS FOR MITIGATION OR ADAPTATION AGAINST CLIMATE CHANGE
- Y02T—CLIMATE CHANGE MITIGATION TECHNOLOGIES RELATED TO TRANSPORTATION
- Y02T90/00—Enabling technologies or technologies with a potential or indirect contribution to GHG emissions mitigation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- General Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Mathematical Physics (AREA)
- General Engineering & Computer Science (AREA)
- Fluid Mechanics (AREA)
- Algebra (AREA)
- Data Mining & Analysis (AREA)
- Mathematical Analysis (AREA)
- Mathematical Optimization (AREA)
- Pure & Applied Mathematics (AREA)
- Computer Hardware Design (AREA)
- Geometry (AREA)
- Evolutionary Computation (AREA)
- Computational Mathematics (AREA)
- Computing Systems (AREA)
- Databases & Information Systems (AREA)
- Software Systems (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Aerodynamic Tests, Hydrodynamic Tests, Wind Tunnels, And Water Tanks (AREA)
Abstract
本发明公开了一种采用数字化对风洞试验模型堵塞度进行计算的方法,涉及风洞试验领域,基于试验模型网格面元和迎风截面平面格子投射的原理,通过在如下公式中引入试验模型的迎风截面积,以对风洞试验模型堵塞度σ进行计算:σ=S W /S FD ;其中,S FD 为风洞试验段横截面积或喷管出口面积,S W 为试验模型迎风截面积。本发明公开了一种采用数字化对风洞试验模型堵塞度进行计算的方法,能在风洞试验过程中,对处于任何姿态的试验模型均可以采用本发明算法进行精确的堵塞度计算。
Description
技术领域
本发明涉及涉及风洞试验领域。更具体地说,本发明涉及一种采用数字化对风洞试验模型堵塞度进行计算的方法。
背景技术
阻塞度是风洞试验模型的最大迎风面积与试验段横截面积(或喷管出口面积)之比。阻塞度可以用来描述风洞阻塞效应的大小。无论低速风洞、高速风洞或超高速风洞,在试验方案设计时,堵塞度都是需要重点考虑的因素之一。风洞试验模型在有限试验段空间里的占位会导致阻塞,当气流流经模型时会产生较为宽阔的侧面绕流和尾流,从而导致阻塞效更为显著。风洞壁面对气流绕流的约束称为“实体阻塞”,对尾流的约束称为“尾流阻塞”,上述两种洞壁干扰即为阻塞效应。
在风洞试验研究中,一般来说,为了得到准确的气动力测量结果,模型的风洞阻塞度不应超过5%。另一方面,为了尽量达到与真实外形的物理相似特别是雷诺数接近,气动噪声测量模型的尺寸需要尽可能大,这就与阻塞度的要求发生了矛盾。为了尽量得到与真实外形雷诺数接近的试验结果,人们采用了各种措施减弱风洞阻塞等洞壁干扰效应,如开口试验段、开槽壁或流线型壁等。但是对洞壁的改进并不能完全消除上述干扰,必须对残存的洞壁效应进行修正。常用的修正方法有映象法、壁压信息法和计算法等。其中,映象法最为简单易用,但阻塞度较大时修正准度难以保证;壁压信息法适应性强,对气动力修正精度高,但存在非定常流动时会受到测量点位置选取和壁压测量结果不确定性的影响。计算法最初在20世纪80年代提出,但一直受到计算量过大、计算结果不确定性大等限制。近十几年来计算技术的长足发展,使计算修正法重新被人们重视。由于可以得到流场细节信息,计算法在机理研究方而有着独到的优势,可以用来揭示不同阻塞度下洞壁干扰产生的机理,并用于建立更准确的修正模型。
除了模型的迎风面积不得超过试验段横截面积的5%,即模型在风洞试验段中的阻塞度不得超过5%,一般而言模型的宽度不得超过试验段宽度的70%,且模型的长度应在风洞静压梯度比较均匀的区域内,以减小风洞试验段内沿风洞轴线方向的静压梯度影响。
无论是风洞试验方案设计还是对试验结果进行堵塞度相关修正,都需要获得堵塞度的量值。对于形状比较规整的试验模型如旋成体之类,通过简单换算基本可以估计堵塞度,但对于复杂外形的试验模型,要对其在任意姿态下的堵塞度精确计算还没有太好的办法。
发明内容
本发明的一个目的是解决至少上述问题和/或缺陷,并提供至少后面将说明的优点。
为了实现本发明的这些目的和其它优点,提供了一种采用数字化对风洞试验模型堵塞度进行计算的方法,基于试验模型网格面元和迎风截面平面格子投射的原理,通过在如下公式中引入试验模型的迎风截面积,以对风洞试验模型堵塞度σ进行计算:
σ=S
W
/S
FD
其中,S FD 为风洞试验段横截面积或喷管出口面积,S W 为试验模型迎风截面积。
优选的是,所述试验模型的迎风截面积获取方式为:
S1、基于试验模型数模并生成表面网格,在风洞试验段内建立基准坐标系,计算试验模型网格面元中心位置、面积及单位法向量;
S2、计算试验模型迎风截面上面元投射格子面积和尺度;
S3、根据试验模型姿态角计算迎风截面投影区域坐标值域,根据试验模型面元投射格子尺度在迎风截面投影值域内划分规范化格子;
S4、遍历投影值域格子及试验模型网格面元记录投射信息,并根据面元-格子投射信息统计以计算迎风截面积。
优选的是,在S1中,所述表面网格采用网格生成软件生成,且表面网格包括四边形结构网格和三角形非结构网格,试验模型基准坐标系基于空气动力学建立;
所述表面网格的面元中心位置(x,y,z)、面积ds及单位法向量基于平面解析几何算法得到。
优选的是,在S2中,试验模型迎风截面上面元投射格子面积和尺度的获取方法为:
S21、试验模型迎风截面上面元投射格子面积为的计算公式如下:
其中,c为试验模型迎风截面上面元投射格子尺度调节系数,N为试验模型面元总数量,ds i 为序号i的面元面积;
S22、试验模型迎风截面上面元投射格子尺度L C 的计算公式如下:
。
优选的是,在S3中,基于试验模型姿态角计算迎风截面投影区域坐标值域的方法如下:
S31、在飞行器处于攻角α和侧滑角β时,任意基准面元纵向位置角变为α 0 +α,侧向位置角变为β 0 +β,则面元中心位置(x,y,z)的计算公式为:
其中,r 0 为面元距坐标原点或重心的距离,其表达式为:
(x 0 ,y 0 ,z 0 )为基准状态下试验模型面元中心位置坐标;
α 0 和β 0 分别为面元基准状态时的纵向位置角和侧向位置角,其表达式为:
其中,面元中心坐标Z分量为负时纵向位置角为正;当面元中心坐标X分量为零,则Y分量为正时侧向位置角为正90度,Y分量为负时侧向位置角为负90度;
S32、试验模型迎风截面投影范围坐标值域分别用(Y max ,Y min )和(Z max ,Z min )表示,则其计算方法如下:
其中,下标j代表所有全部面元对应的坐标值,MAXVAL为取最大值,MINVAL为取最小值;
则试验模型迎风截面所在平面等价于基准状态体轴系坐标系的YZ平面,迎风截面投影范围坐标值域等价于面元中心位置YZ坐标的极大值和极小值。
优选的是,在S3中,基于面元投射格子尺度在迎风截面投影范围内划分规范化格子的方法为:
S33、定义坐标投影范围内在横轴和竖轴上的规范化坐标为(I Y ,I Z ),则一组(I Y ,I Z )坐标对应一个边长为L C 的正方形格子,所述(I Y ,I Z )的计算公式如下:
其中,INT()为最小整数,一组面元坐标(y,z)对应一组规范化整数坐标(I Y ,I Z ),而一组(I Y ,I Z )能对应多组(y,z);
S34、规范化坐标(I Y ,I Z )的值域范围I Y,min、I Y,max、I Z,min、I Z,max由投影范围边界确定,其计算方法如下:
S35、投影范围内规范化格子的数量N YZ 的计算公式如下:
。
优选的是,在S4中,所述迎风截面积计算方式为:
S41、将投影范围内每个格子的面元投影信息通过下述数据结构来表示:
其中,I X ,I Y 为格子规范化坐标,N jXY 为某个格子内投影的面元数量,
I J (j)和K J (j)分别为格子内投影的面元序号和单位法向量状态,i XY 为格子序号;
S42、在进行投射信息记录前,遍历所有格子并初始化,以使所有N jXY 为0,I J (j)和K J (j)为空数组;
S43、遍历所有格子和面元,记录每个格子中投射的面元数量N jXY 、面元序号I J (j)、单位法向量状态K J (j),则迎风截面积为S W 的计算公式为:
其中,N XY 为投影范围内N jXY 不为0的所有格子数量,为每个格子的面积。
优选的是,试验模型任意姿态在风洞试验段基准坐标系中三个维度上的尺度ΔX、ΔY、ΔZ分别通过下式取得:
ΔX=X max -X min
ΔY=Y max -Y min
ΔZ=Z max -Z min
其中,X max 和X min 的计算公式如下:
X max =MAXVAL(x j )
X min =MINVAL(x j )。
本发明至少包括以下有益效果:本发明构建了一种基于面元-格子投射原理的风洞试验模型堵塞度的数字化算法,在模型设计加工和试验方案设计高度数字化条件下,对风洞试验模型堵塞度通过迎风截面积进行快速高效计算,且易于推广至不同应用场景,故风洞试验模型堵塞度可以是任意姿态下的算法结果;此外根据面元随姿态变化坐标变换的结果信息统计,还可以同时得出试验模型在基准坐标系三个维度的值域范围,这些参数都是风洞试验方案需要考虑的一些方面,通过本发明同时提供了数字化的精确量化算法,具体效果如下:
其一,采用完全的数字化方法快速高效地获得风洞试验模型堵塞度及其它关联重要参数。
其二,针对风洞试验过程中,无论试验模型处于任何姿态均可以采用本发明算法精确堵塞度。
其三,根据格子投射面元的单位法向量纵轴分量,可以推断风洞试验中试验模型面元的迎风/背风/顺风状态。
本发明的其它优点、目标和特征将部分通过下面的说明体现,部分还将通过对本发明的研究和实践而为本领域的技术人员所理解。
附图说明
图1为风洞试验模型迎风面积面元-格子投射原理图;
图2为本发明一种风洞试验模型堵塞度的数字化算法流程;
图3为HB2标模数模及基准坐标系示意图;
图4为HB2标模结构化四边形表面网格;
图5为0°攻角、0°侧滑角基准状态下风洞试验模型迎风截面投影情况;
图6为90°侧滑角状态下风洞试验模型迎风截面投影情况;
图7为90°攻角状态下风洞试验模型迎风截面投影情况;
图8为45°攻角、45°侧滑角状态下风洞试验模型迎风截面投影情况。
具体实施方式
下面结合附图对本发明做进一步的详细说明,以令本领域技术人员参照说明书文字能够据以实施。
本发明提供的一种风洞试验模型堵塞度的数字化算法,基于试验模型网格面元-迎风截面平面格子投射原理,如图1所示。在风洞试验方案设计及模型设计加工高度数字化条件下,对风洞试验模型堵塞度通过迎风截面积进行快速高效计算的算法新颖而有效,且易于推广至不同应用场景。本发明对风洞试验模型堵塞度的计算可以是模型处于任意姿态下的算法结果。此外,根据面元随姿态变化坐标变换的结果信息统计,还可以同时得出试验模型在基准坐标系三个维度的值域范围,这些参数都是风洞试验方案需要考虑的一些方面,通过本发明同时提供了数字化的精确量化算法。
本发明提供的一种风洞试验模型堵塞度的数字化算法,如图2流程图所示,具有以下四个步骤:
步骤一、构建或应用试验模型数模并生成表面网格,在风洞试验段内建立基准坐标系,计算试验模型网格面元中心位置、面积及单位法向量
试验模型数模即试验模型设计及交付模型加工中常用的三维计算机辅助设计(CAD:Computer Aided Design)电子文件,风洞试验模型数模的建模软件包括但不限于ProE、UG、Solidworks、CATIA等。
采用网格生成软件生成试验模型表面网格,这些软件包括但不限于Pointwise、GridStar、Gmsh等。
针对试验模型数模生成的表面网格包括四边形结构网格或三角形非结构网格两类;四边形结构网格文件采用PLOT3D格式(行业通用结构化数据格式,源于NASA),三角形非结构网格文件采用UCD格式(Unstructured Cell Data,行业通用非结构网格文件格式)。
根据上述模型装置及其机构在风洞试验段中位置,构建和风洞试验段固连的基准坐标系,具体方式如下:
坐标原点位于试验模型旋转中心位置,记为原点O,坐标为(0,0,0);根据几何关系易知,如旋转中心变化,并不影响迎风面积量值,影响的是迎风截面投影区域在投影平面上的整体平移位置。故本发明算法中设定旋转中心不变,不影响最终算法结果。
坐标X轴为航天器主纵轴,其正方向与气流方向相反,由重心指向航天器头部或前部;
坐标Y轴为航天器横轴,其正向由重心指向航天器右侧;
坐标Z轴根据坐标X轴和坐标Y轴按右手法则确定,具体其正向由重心指向航天器下方。
基准坐标系建立后不再变化,试验模型姿态变化后的面元坐标等信息均在这个基准坐标系中进行描述;
试验模型的迎风截面参数等价于与气流方向垂直的基准坐标系的YZ平面内参数,故YZ平面内参数是迎风面积计算的重要基础。
设四边形面元M 1 M 2 M 3 M 4的坐标为 (x i ,y i ,z i ,),i=1,2,3,4,于是四边形面元中心位置、面积及单位法向量的计算方法如下:
四边形面元的中心点坐标设为(x,y,z,),其坐标值的计算公式为:
其中,面元中心各坐标值即为四边形4个角点相应坐标的平均值;
设面元的单位外法向量为,其计算公式如下:
其中,分别是3个坐标方向上的单位向量;四边形对角线构成的向量/>的表达式为:
设
其中:
于是
于是
设四边形面元的面积为ds,其计算公式如下:
其中,a 1为点M 1到M 2的距离,a 2为点M 1到M 4的距离,b 1为点M 2到M 3的距离,b 2为点M 4到M 3的距离,c为点M 1到M 3的距离。
其中:
对于三角形面元,设三角形的3个角点坐标分别为/>、,于是三角形面元中心位置、面积及单位外法向量的计算方法如下:
三角形面元的中心点坐标设为(x,y,z,),其坐标值的计算公式为:
其中,面元中心各坐标值即为三角形3个角点相应坐标的平均值;
设面元的单位外法向量为,其计算公式如下:
其中,分别是3个坐标方向上的单位向量;三角形两个邻边构成的向量的表达式为:
设三角形面元的面积为ds,其计算公式如下:
步骤二、计算试验模型迎风截面上面元投射格子面积和尺度
根据上述面元参数,计算迎风截面(等价于YZ平面)格子规范化划分基本依据的格子面积和格子尺度的计算方法如下:
定义面元投射格子面积,其计算公式如下:
其中,c定义为面元投射格子尺度系数,用于调节面元投射格子面积和尺度,其值大小决定了一个格子可能投射的面元数量,c值过小可能导致格子无投射面元而丢失,c值过大则显得粗糙可能导致计算精度;根据经验推荐1<<3;
其中,N为试验模型面元总数量,ds i 为序号i的面元面积;
定义试验模型面元投射格子尺度L c ,其计算公式如下:
其中,试验模型面元投射格子尺度的含义即面元投射格子面积相等的正方形边长。
步骤三、根据试验模型姿态角计算迎风截面投影区域坐标值域,根据试验模型面元投射格子尺度在迎风截面投影值域内划分规范化格子
设试验模型处于任意攻角α和侧滑角β时,面元纵向位置角变为α 0 +α,侧向位置角变为β 0 +β,则面元中心位置(x,y,z)的计算公式为:
其中,α和β分别为试验模型随姿态机构绕旋转中心旋转的攻角和侧滑角。
上述式子中,r 0 为面元距坐标原点即旋转中心的距离,即面元矢径,在试验模型姿态变化过程中保持不变,其表达式为:
其中,(x 0 ,y 0 ,z 0 )为基准状态下试验模型面元中心位置坐标。
上述式子中,α 0 和β 0 分别为面元基准状态时的纵向位置角和侧向位置角,其表达式如下:
其中,面元中心坐标Z分量为负时纵向位置角为正;当面元中心坐标X分量为零,则Y分量为正时侧向位置角为正90度,Y分量为负时侧向位置角为负90度。上述式子同样适用于面元处于任意位置处。
设试验模型迎风截面投影区域坐标值域分别用(Y max ,Y min )和(Z max ,Z min )表示,则其计算方法如下:
其中,下标j代表所有全部面元对应的坐标值;上述式子中,表明航天器迎风截面所在平面即等价于基准状态体轴系坐标系的YZ平面,迎风截面投影范围坐标值域即等价于面元中心位置YZ坐标的极大值和极小值。
根据试验模型面元投射格子尺度在迎风截面投影值域内划分规范化格子的方法如下:
定义(I Y ,I Z )为投影值域内在横轴和竖轴上的规范化坐标,一组v坐标对应一个边长为L C 的正方形格子。
面元投影对应的(I Y ,I Z )由面元中心位置在横轴和竖轴上的坐标分量(y,z)转化而来,其计算方法如下:
其中,INT()函数的意义为对括号中取最小整数。一组面元坐标(y,z)对应一组规范化整数坐标v,一组(I Y ,I Z )可能对应多组(y,z)。
投影值域内规范化坐标(I Y ,I Z )的值域范围由投影区域边界确定,计算方法如下:
设计算投影值域内规范化格子的数量为N YZ ,计算公式如下:
由此即完成根据试验模型面元投射格子尺度在迎风截面二维平面坐标系内划分规范化格子的过程。
步骤四、遍历投影值域格子及试验模型网格面元记录投射信息,根据面元-格子投射信息统计计算迎风截面积并结合风洞试验段尺寸计算试验模型堵塞度
首先定义投影值域内每个格子的面元投影信息以数据结构表示,其中I X ,I Y 为格子规范化坐标,N jXY 为某个格子内投影的面元数量,I J (j)和K J (j)分别为格子内投影的面元序号和单位法向量状态,i XY 为格子序号。
进行投射信息记录前,遍历所有格子初始化,所有N jXY 为0,I J (j)和K J (j)为空数组。
遍历所有格子和面元,记录每个格子中投射的面元数量N jXY 、面元序号I J (j)、单位法向量状态K J (j),其中,单位法向量状态取值的规则如下:
若n x <0则K J (j)=1,若n x >0则K J (j)=-1,n x =0则v
其中,上述三种状态分别对应试验模型面元迎风、背风和顺风。
设迎风截面积为S W ,其计算公式为:
其中,N XY 为投影范围内N jXY 不为0的所有格子数量,其值通过数据结构进行逐一判识计数N jXY 不为0状态而得;/>为每个格子的面积即试验模型面元投射格子面积。
设风洞试验模型堵塞度为σ,其计算公式如下:
σ=S
W
/S
FD
其中,S FD 为风洞试验段横截面积或喷管出口面积,S W 为试验模型迎风截面积。
设试验模型任意姿态下,在风洞试验段基准坐标系中三个维度上的占位尺度分别是ΔX、ΔY、ΔZ,其计算公式如下:
ΔX=X max -X min
ΔY=Y max -Y min
ΔZ=Z max -Z min
其中,X max 和X min 的计算公式如下:
X max =MAXVAL(x j )
X min =MINVAL(x j )
其余Y、Z的最大值和最小值见前述公式。
实施例:
本发明的一种风洞试验模型堵塞度的数字化算法,为了更清楚地说明本发明的技术方案,通过实施例进行说明。实施例采用HB-2标模模型作为例子,针对典型的不同风洞试验场景,应用本发明描述的方法对其堵塞度的计算过程进行详细描述。如图2所示,具体算法过程和结果如下。
步骤一、构建或应用试验模型数模并生成表面网格,在风洞试验段内建立基准坐标系,计算试验模型网格面元中心位置、面积及单位法向量
采用建模软件生成的HB2标模模型如图3所示,该模型纵向长度约0.6125米,底部直径约0.2米;模型表面总面积约为0.29146平方米。
在风洞试验段内建立基准状态坐标系,如图3所示。实际风洞试验模型支撑机构运行改变模型姿态时,模型的旋转中心可能会有明显变化,但是根据几何关系易知,即使旋转中心变化,并不影响迎风面积量值,影响的是迎风截面投影区域在投影平面上的整体平移位置。故本发明算法中设定旋转中心不变,不影响最终算法结果。试验模型基准坐标系实际即为试验模型姿态归零情况下的体轴系,因此称之为基准坐标系。该基准状态坐标系建立后不再变化,后续试验模型姿态变化后的面元坐标法向量等信息均在这个基准坐标系中进行描述。
采用网格软件生成的试验模型简化模型表面非结构网格如图4所示,图中给出了一个局部放大图。根据统计,本实例划分的结构化四边形网格面元数量为23814个。
对实施例HB2标模模型的上述23814个面元逐个计算其中心位置、面积和单位法向量。需要注意的是,除了面元面积保持不变,随着模型姿态变化,其每个面元的中心位置和单位法向量都是动态变化的,这也是导致迎风面积随姿态不断变化的几何原因。
步骤二、计算试验模型迎风截面上面元投射格子面积和尺度
根据上述面元的基本信息,计算得出实施例HB2标模模型的面元投射格子面积约为12.239平方毫米,对应面元投射格子尺度约为3.498毫米。在本例计算中,面元投射格子尺度调节系数取c=2;一般情况下要求,建议/>;取值过小可能丢失格子,但取值过大则会影响计算精度。
步骤三、根据试验模型姿态角计算迎风截面投影区域坐标值域,根据试验模型面元投射格子尺度在迎风截面投影值域内划分规范化格子
针对试验模型任意姿态角计算迎风面积的中间步骤之一是获得投影范围的坐标值域,基于空气动力学中攻角和侧滑角的对应关系,如图5-图8,实施例给出了一系列攻角(绕Y轴旋转)和侧滑角(绕Z轴旋转)不同状态下的投影区域坐标值域实际计算结果。具体计算结果如表1所示。
表1
根据上述投影范围值域,结合试验模型面元投射格子尺度在迎风截面投影值域内划分规范化格子,这些格子都是规范化的正方形,它们的规范化坐标可以用整数进行表示。具体计算结果如表2所示,表中给出了规范化格子的值域范围,同时对格子数量进行了统计,需要了解的是这些格子并不是都会被面元进行投射。
表2
步骤四、遍历投影值域格子及试验模型网格面元记录投射信息,根据面元-格子投射信息统计计算迎风截面积并结合风洞试验段尺寸计算试验模型堵塞度
面元-格子投射信息的记录,需要遍历投影范围内的所有格子。首先对属于每个格子的投射信息数据结构中的投影面元数量、序号数组及法向量状态数组初始化,然后再同时遍历格子和面元进行投射匹配同时更新投射信息。
本实施例HB2标模模型的迎风面积和堵塞度在各种姿态条件下计算如表3所示,表中采用了两种典型的风洞口径作为堵塞度的计算基础。实际风洞试验中,对于不同的风洞类型,作为堵塞度计算的风洞尺度参数需根据模型安放方式和安放位置具体考虑,在本发明算法得出试验模型迎风面积量值后,具体的堵塞度只需要经过简单的换算了。本实施例HB2标模模型1万量级的面元数量条件下,用主流PC每个状态的计算时间在毫秒量级,非常快速。
表3
此外,实施例中试验模型在不同姿态下在风洞试验段基准坐标系中三个维度上的占位尺度如表4所示,这些参数都是制定风洞试验方案的主要依据。
表4
以上方案只是一种较佳实例的说明,但并不局限于此。在实施本发明时,可以根据使用者需求进行适当的替换和/或修改。
这里说明的设备数量和处理规模是用来简化本发明的说明的。对本发明的应用、修改和变化对本领域的技术人员来说是显而易见的。
尽管本发明的实施方案已公开如上,但其并不仅仅限于说明书和实施方式中所列运用。它完全可以被适用于各种适合本发明的领域。对于熟悉本领域的人员而言,可容易地实现另外的修改。因此在不背离权利要求及等同范围所限定的一般概念下,本发明并不限于特定的细节和这里示出与描述的图例。
Claims (8)
1.采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,基于试验模型网格面元和迎风截面平面格子投射的原理,通过在如下公式中引入试验模型的迎风截面积,以对风洞试验模型堵塞度σ进行计算:
σ=S
W
/S
FD
其中,S FD 为风洞试验段横截面积或喷管出口面积,S W 为试验模型迎风截面积。
2.如权利要求1所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,所述试验模型的迎风截面积获取方式为:
S1、基于试验模型数模并生成表面网格,在风洞试验段内建立基准坐标系,计算试验模型网格面元中心位置、面积及单位法向量;
S2、计算试验模型迎风截面上面元投射格子面积和尺度;
S3、根据试验模型姿态角计算迎风截面投影区域坐标值域,根据试验模型面元投射格子尺度在迎风截面投影值域内划分规范化格子;
S4、遍历投影值域格子及试验模型网格面元记录投射信息,并根据面元-格子投射信息统计以计算迎风截面积。
3.如权利要求2所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,在S1中,所述表面网格采用网格生成软件生成,且表面网格包括四边形结构网格和三角形非结构网格,试验模型基准坐标系基于空气动力学建立;
所述表面网格的面元中心位置(x,y,z)、面积ds及单位法向量基于平面解析几何算法得到。
4.如权利要求3所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,在S2中,试验模型迎风截面上面元投射格子面积和尺度的获取方法为:
S21、试验模型迎风截面上面元投射格子面积为的计算公式如下:
其中,c为试验模型迎风截面上面元投射格子尺度调节系数,N为试验模型面元总数量,ds i 为序号i的面元面积;
S22、试验模型迎风截面上面元投射格子尺度L C 的计算公式如下:
。
5.如权利要求4所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,在S3中,基于试验模型姿态角计算迎风截面投影区域坐标值域的方法如下:
S31、在飞行器处于攻角α和侧滑角β时,任意基准面元纵向位置角变为α 0 +α,侧向位置角变为β 0 +β,则面元中心位置(x,y,z)的计算公式为:
其中,r 0 为面元距坐标原点或重心的距离,其表达式为:
(x 0 ,y 0 ,z 0 )为基准状态下试验模型面元中心位置坐标;
α 0 和β 0 分别为面元基准状态时的纵向位置角和侧向位置角,其表达式为:
其中,面元中心坐标Z分量为负时纵向位置角为正;当面元中心坐标X分量为零,则Y分量为正时侧向位置角为正90度,Y分量为负时侧向位置角为负90度;
S32、试验模型迎风截面投影范围坐标值域分别用(Y max ,Y min )和(Z max ,Z min )表示,则其计算方法如下:
其中,下标j代表所有全部面元对应的坐标值,MAXVAL为取最大值,MINVAL为取最小值;
则试验模型迎风截面所在平面等价于基准状态体轴系坐标系的YZ平面,迎风截面投影范围坐标值域等价于面元中心位置YZ坐标的极大值和极小值。
6.如权利要求5所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,在S3中,基于面元投射格子尺度在迎风截面投影范围内划分规范化格子的方法为:
S33、定义坐标投影范围内在横轴和竖轴上的规范化坐标为(I Y ,I Z ),则一组(I Y ,I Z )坐标对应一个边长为L C 的正方形格子,所述(I Y ,I Z )的计算公式如下:
其中,INT()为最小整数,一组面元坐标(y,z)对应一组规范化整数坐标(I Y ,I Z ),而一组(I Y ,I Z )能对应多组(y,z);
S34、规范化坐标(I Y ,I Z )的值域范围I Y,min、I Y,max、I Z,min、I Z,max由投影范围边界确定,其计算方法如下:
S35、投影范围内规范化格子的数量N YZ 的计算公式如下:
。
7.如权利要求6所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,在S4中,所述迎风截面积计算方式为:
S41、将投影范围内每个格子的面元投影信息通过下述数据结构来表示:
其中,I X ,I Y 为格子规范化坐标,N jXY 为某个格子内投影的面元数量,
I J (j)和K J (j)分别为格子内投影的面元序号和单位法向量状态,i XY 为格子序号;
S42、在进行投射信息记录前,遍历所有格子并初始化,以使所有N jXY 为0,I J (j)和K J (j)为空数组;
S43、遍历所有格子和面元,记录每个格子中投射的面元数量N jXY 、面元序号I J (j)、单位法向量状态K J (j),则迎风截面积为S W 的计算公式为:
其中,N XY 为投影范围内N jXY 不为0的所有格子数量,为每个格子的面积。
8.如权利要求7所述的采用数字化对风洞试验模型堵塞度进行计算的方法,其特征在于,试验模型任意姿态在风洞试验段基准坐标系中三个维度上的尺度ΔX、ΔY、ΔZ分别通过下式取得:
ΔX=X max -X min
ΔY=Y max -Y min
ΔZ=Z max -Z min
其中,X max 和X min 的计算公式如下:
X max =MAXVAL(x j )
X min =MINVAL(x j )。
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311161849.2A CN116894408B (zh) | 2023-09-11 | 2023-09-11 | 采用数字化对风洞试验模型堵塞度进行计算的方法 |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202311161849.2A CN116894408B (zh) | 2023-09-11 | 2023-09-11 | 采用数字化对风洞试验模型堵塞度进行计算的方法 |
Publications (2)
Publication Number | Publication Date |
---|---|
CN116894408A true CN116894408A (zh) | 2023-10-17 |
CN116894408B CN116894408B (zh) | 2023-12-05 |
Family
ID=88315295
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202311161849.2A Active CN116894408B (zh) | 2023-09-11 | 2023-09-11 | 采用数字化对风洞试验模型堵塞度进行计算的方法 |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN116894408B (zh) |
Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
AU2007229443A1 (en) * | 2006-10-23 | 2008-05-08 | Fujikura Energy Systems Ltd. | Polygonal overhead cable |
CN106524351A (zh) * | 2016-11-15 | 2017-03-22 | 大道网络(上海)股份有限公司 | 一种控制台通风管理系统 |
CN109141805A (zh) * | 2018-09-20 | 2019-01-04 | 北京空间机电研究所 | 一种降落伞气动力系数及力矩系数计算方法及系统 |
CN111272380A (zh) * | 2020-02-24 | 2020-06-12 | 中国空气动力研究与发展中心高速空气动力研究所 | 风洞试验模型位姿视频测量的风轴系自标定方法 |
CN115615655A (zh) * | 2022-11-15 | 2023-01-17 | 中国航空工业集团公司沈阳空气动力研究所 | 一种连续式跨声速风洞试验数据干扰修正方法 |
US20230113577A1 (en) * | 2021-10-08 | 2023-04-13 | Beijing Institute Of Technology | Configuration construction and attitude control method for pyramid deorbit sail |
-
2023
- 2023-09-11 CN CN202311161849.2A patent/CN116894408B/zh active Active
Patent Citations (6)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
AU2007229443A1 (en) * | 2006-10-23 | 2008-05-08 | Fujikura Energy Systems Ltd. | Polygonal overhead cable |
CN106524351A (zh) * | 2016-11-15 | 2017-03-22 | 大道网络(上海)股份有限公司 | 一种控制台通风管理系统 |
CN109141805A (zh) * | 2018-09-20 | 2019-01-04 | 北京空间机电研究所 | 一种降落伞气动力系数及力矩系数计算方法及系统 |
CN111272380A (zh) * | 2020-02-24 | 2020-06-12 | 中国空气动力研究与发展中心高速空气动力研究所 | 风洞试验模型位姿视频测量的风轴系自标定方法 |
US20230113577A1 (en) * | 2021-10-08 | 2023-04-13 | Beijing Institute Of Technology | Configuration construction and attitude control method for pyramid deorbit sail |
CN115615655A (zh) * | 2022-11-15 | 2023-01-17 | 中国航空工业集团公司沈阳空气动力研究所 | 一种连续式跨声速风洞试验数据干扰修正方法 |
Non-Patent Citations (3)
Title |
---|
杨成;李勰;孙军;李志辉;: "类天宫航天器迎风面积建模及变化特性分析", 载人航天, vol. 26, no. 4, pages 429 - 435 * |
王泽江等: "高超声速通气模型直接测力试验", 航空学报, vol. 38, no. 10, pages 30 - 37 * |
马静;贾青;杨志刚;: "基于数值计算的高速列车气动阻力风洞试验缩比模型选取方法", 计算机辅助工程, no. 3, pages 110 - 113 * |
Also Published As
Publication number | Publication date |
---|---|
CN116894408B (zh) | 2023-12-05 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN108986048B (zh) | 基于线激光扫描三维点云快速复合滤波处理方法 | |
CN110389391B (zh) | 一种基于空间域的重磁位场解析延拓方法 | |
Alben et al. | How flexibility induces streamlining in a two-dimensional flow | |
US8833153B2 (en) | Correction of pressure signals measured during supersonic wind tunnel testing | |
Prospathopoulos et al. | Application of a 3D RANS solver on the complex hill of Bolund and assessment of the wind flow predictions | |
Ramírez et al. | New high-resolution-preserving sliding mesh techniques for higher-order finite volume schemes | |
CN115615655A (zh) | 一种连续式跨声速风洞试验数据干扰修正方法 | |
CN110059396B (zh) | 飞行器飞行参数解算方法 | |
CN111652941B (zh) | 基于自适应变异的天牛群优化算法的相机内参校准方法 | |
CN111238489B (zh) | 一种低轨卫星大气阻力摄动建模和计算方法 | |
CN114444216B (zh) | 基于数值模拟的高空条件下飞行器姿态控制方法及系统 | |
Leifsson et al. | Simulation-driven design of low-speed wind tunnel contraction | |
CN115358026A (zh) | 一种基于多元线性回归与曲面拟合的五孔探针数据处理方法 | |
CN115774968A (zh) | 基于递归分解的曲面网格生成方法及系统、计算机设备 | |
CN112163381A (zh) | 一种适用于复杂地形风场流动数值模拟的侧向边界条件设置方法 | |
CN114838906A (zh) | 对风洞来流能谱进行测量的方法及系统、设备、存储介质 | |
Abdul-Rahman et al. | Freeform texture representation and characterisation based on triangular mesh projection techniques | |
CN111754464A (zh) | 一种类pd算法结合icp算法的零件精确找正方法 | |
CN116894408B (zh) | 采用数字化对风洞试验模型堵塞度进行计算的方法 | |
CN111856072B (zh) | 一种气流速度计算方法、系统、设备及一种存储介质 | |
CN114611420A (zh) | 非定常气动力计算精度评估及修正方法 | |
CN110211189A (zh) | ToF相机深度误差建模校正方法及装置 | |
Yu et al. | Mesh optimization using an improved self-organizing mechanism | |
CN111768435B (zh) | 一种应用于零件自动找正的自适应步长点云配正方法 | |
CN113536461B (zh) | 用于高超声速强激波流场气动热预测的湍流模型修正方法 |
Legal Events
Date | Code | Title | Description |
---|---|---|---|
PB01 | Publication | ||
PB01 | Publication | ||
SE01 | Entry into force of request for substantive examination | ||
SE01 | Entry into force of request for substantive examination | ||
GR01 | Patent grant | ||
GR01 | Patent grant |