CN115081121A - 一种考虑溪流现象的结冰模拟方法和存储介质 - Google Patents

一种考虑溪流现象的结冰模拟方法和存储介质 Download PDF

Info

Publication number
CN115081121A
CN115081121A CN202211003148.1A CN202211003148A CN115081121A CN 115081121 A CN115081121 A CN 115081121A CN 202211003148 A CN202211003148 A CN 202211003148A CN 115081121 A CN115081121 A CN 115081121A
Authority
CN
China
Prior art keywords
water film
flow rate
mass flow
icing
moment
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
Application number
CN202211003148.1A
Other languages
English (en)
Other versions
CN115081121B (zh
Inventor
周志宏
吴主龙
熊华杰
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Sichuan University
Original Assignee
Sichuan University
Priority date (The priority date 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 date listed.)
Filing date
Publication date
Application filed by Sichuan University filed Critical Sichuan University
Priority to CN202211003148.1A priority Critical patent/CN115081121B/zh
Publication of CN115081121A publication Critical patent/CN115081121A/zh
Application granted granted Critical
Publication of CN115081121B publication Critical patent/CN115081121B/zh
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Classifications

    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/10Geometric CAD
    • G06F30/15Vehicle, aircraft or watercraft design
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F30/00Computer-aided design [CAD]
    • G06F30/20Design optimisation, verification or simulation
    • G06F30/28Design optimisation, verification or simulation using fluid dynamics, e.g. using Navier-Stokes equations or computational fluid dynamics [CFD]
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2111/00Details relating to CAD techniques
    • G06F2111/10Numerical modelling
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2113/00Details relating to the application field
    • G06F2113/08Fluids
    • GPHYSICS
    • G06COMPUTING; CALCULATING OR COUNTING
    • G06FELECTRIC DIGITAL DATA PROCESSING
    • G06F2119/00Details relating to the type or aim of the analysis or the optimisation
    • G06F2119/14Force analysis or force optimisation, e.g. static or dynamic forces

Landscapes

  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • General Physics & Mathematics (AREA)
  • Theoretical Computer Science (AREA)
  • Geometry (AREA)
  • Mathematical Optimization (AREA)
  • General Engineering & Computer Science (AREA)
  • Mathematical Analysis (AREA)
  • Evolutionary Computation (AREA)
  • Pure & Applied Mathematics (AREA)
  • Computer Hardware Design (AREA)
  • Aviation & Aerospace Engineering (AREA)
  • Computational Mathematics (AREA)
  • Automation & Control Theory (AREA)
  • Algebra (AREA)
  • Computing Systems (AREA)
  • Fluid Mechanics (AREA)
  • Mathematical Physics (AREA)
  • Management, Administration, Business Operations System, And Electronic Commerce (AREA)

Abstract

本申请涉及结冰模拟领域,提供了一种考虑溪流现象的结冰模拟方法和存储介质。该方法包括:以水膜平均速度、水膜厚度、溪流现象作用力、结冰质量流率和蒸发质量流率为未知数建立欧拉水膜模型的动量方程和质量方程;以结冰质量流率和蒸发质量流率为未知数建立Messinger模型热力学平衡方程;根据上一时刻的水膜厚度,采用假设修正法获取当前时刻的结冰质量流率和蒸发质量流率;获取当前时刻的溪流现象作用力,并根据当前时刻的结冰质量流率、蒸发质量流率、溪流现象作用力、动量方程和质量方程,获取当前时刻的水膜平均速度、水膜厚度,通过迭代计算获取每个时刻的结冰质量流率。通过上述方法,可以有效解决结冰模拟数值误差较大、精度较低的问题。

Description

一种考虑溪流现象的结冰模拟方法和存储介质
技术领域
本申请涉及结冰模拟领域,更具体地,涉及一种考虑溪流现象的结冰模拟方法和存储介质。
背景技术
飞机穿越云层时,过冷水滴撞击到机体后,可能发生相变并导致结冰现象。结冰会改变飞机的外形与绕流流场,破坏气动性能,降低操纵性与稳定性,威胁飞行安全,严重时导致空难事故。因此对飞机结冰模拟的研究非常重要。
目前,现有技术在对飞机结冰时的水膜流动进行模拟时,存在结冰模拟数值误差较大、精度较低的问题。
发明内容
本申请发明人在通过长期实践发现,在物面水膜厚度差异较大时,结冰模拟数值在这个区域的误差较为明显,精度较低。在经过发明人的长期思考后,发现这是由于溪流现象造成的。溪流现象是指,在水膜较厚的位置,水膜的表面张力小,因此水膜更容易流动,在水膜较薄的位置,水膜的表面张力大,因此水膜更难流动。具体地,在实际中,水膜流动多在水膜较厚的位置,在水膜厚度较小的位置,水膜流动较少,因此水膜流动非常不均匀,但是在现有技术的结冰模型中,得到的水膜流动是较为均匀的,因此在物面水膜厚度差异较大的区域,结冰模拟数值的误差较大、精度较低。
基于此,本申请提出了一种考虑溪流现象的结冰模拟方法,以水膜平均速度
Figure 879605DEST_PATH_IMAGE001
、水膜厚度h和溪流现象作用力
Figure 468849DEST_PATH_IMAGE002
作为未知数建立欧拉水膜模型的动量方程,并以水膜平均速度
Figure 983007DEST_PATH_IMAGE001
,水膜厚度h、结冰质量流率
Figure 847058DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 333534DEST_PATH_IMAGE004
作为未知数建立欧拉水膜模型的质量方程;以结冰质量流率
Figure 890417DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 829554DEST_PATH_IMAGE004
作为未知数建立Messinger模型的热力学平衡方程;根据上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的结冰质量流率
Figure 497296DEST_PATH_IMAGE003
和蒸发质量流率
Figure 900595DEST_PATH_IMAGE004
;获取当前时刻的溪流现象作用力
Figure 362801DEST_PATH_IMAGE002
,并根据当前时刻的结冰质量流率
Figure 792164DEST_PATH_IMAGE003
、蒸发质量流率
Figure 263596DEST_PATH_IMAGE004
和溪流现象作用力
Figure 255823DEST_PATH_IMAGE002
,结合欧拉水膜模型的动量方程和质量方程,获取当前时刻的水膜平均速度
Figure 357771DEST_PATH_IMAGE001
、水膜厚度h,通过迭代计算获取每个时刻的结冰质量流率
Figure 68238DEST_PATH_IMAGE003
。如此,可以有效解决现有技术在对飞机结冰时的水膜流动进行模拟时,存在的结冰模拟数值误差较大、精度较低的问题。
第一方面,提供了一种考虑溪流现象的结冰模拟方法,该方法包括:以水膜平均速度
Figure 281045DEST_PATH_IMAGE001
、水膜厚度h和溪流现象作用力
Figure 393357DEST_PATH_IMAGE002
作为未知数建立欧拉水膜模型的动量方程,并以所述水膜平均速度
Figure 462945DEST_PATH_IMAGE001
,所述水膜厚度h、结冰质量流率
Figure 863970DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 411626DEST_PATH_IMAGE004
作为未知数建立欧拉水膜模型的质量方程;以所述结冰质量流率
Figure 644024DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 884513DEST_PATH_IMAGE004
作为未知数建立Messinger模型的热力学平衡方程;根据上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的结冰质量流率
Figure 504325DEST_PATH_IMAGE003
和蒸发质量流率
Figure 121252DEST_PATH_IMAGE004
;获取当前时刻的溪流现象作用力
Figure 942577DEST_PATH_IMAGE002
,并根据所述当前时刻的结冰质量流率
Figure 557229DEST_PATH_IMAGE003
、蒸发质量流率
Figure 464005DEST_PATH_IMAGE004
和溪流现象作用力
Figure 619043DEST_PATH_IMAGE002
,结合所述欧拉水膜模型的动量方程和质量方程,获取当前时刻的水膜平均速度
Figure 763717DEST_PATH_IMAGE001
、水膜厚度h,通过迭代计算获取每个时刻的结冰质量流率
Figure 346008DEST_PATH_IMAGE003
,其中,初始时刻的水膜厚度h为0。
第二方面,还提供了一种计算机可读存储介质,所述计算机可读存储介质中存储有程序代码,所述程序代码可被处理器调用执行上述方法。
综上所述,本申请至少具有如下技术效果:
1.本申请提供的其中一种考虑溪流现象的结冰模拟方法,在欧拉水膜模型中将溪流现象作用力作为未知数,在进行结冰模拟时考虑溪流现象对结冰模型的影响,相比现有技术,结冰模拟数值误差较小、精度较高。
2.本申请提供的其中一种考虑溪流现象的结冰模拟方法,考虑了水膜流动对结冰的影响,建立的结冰模型更加精确,以及,还考虑了空气对水膜的压力梯度、水膜表面张力梯度,计算得到的结冰量更精确。
因此,本申请提供的方案可以有效解决现有技术在对飞机结冰时的水膜流动进行模拟时,存在的结冰模拟数值误差较大、精度较低的问题。
附图说明
为了更清楚地说明本申请实施例中的技术方案,下面将对实施例描述中所需要使用的附图作简单地介绍,显而易见地,下面描述中的附图仅仅是本申请的一些实施例,对于本领域技术人员来讲,在不付出创造性劳动的前提下,还可以根据这些附图获得其他的附图。
图1示出了本申请实施例1提供的一种考虑溪流现象的结冰模拟方法的流程示意图;
图2示出了本申请实施例1提供的溪流现象的示意图;
图3示出了本申请实施例1提供的物面位置与水膜厚度的示意图;
图4示出了本申请实施例1提供的物面位置与结冰率的示意图;
图5示出了本申请实施例2提供的一种计算机可读存储介质的结构框图。
具体实施方式
为了使本技术领域的人员更好地理解本申请方案,下面将结合本申请实施例中的附图,对本申请实施例中的技术方案进行清楚、完整地描述,显然,所描述的实施例仅仅是本申请一部分实施例,而不是全部的实施例。基于本申请中的实施例,本领域普通技术人员在没有作出创造性劳动前提下所获得的所有其他实施例,都属于本申请保护的范围。
目前,现有技术在对飞机结冰时的水膜流动进行模拟时,存在结冰模拟数值误差较大、精度较低的问题。
具体地,在物面水膜厚度差异较大时,结冰模拟数值在这个区域的误差较为明显,精度较低。具体地,在实际中,水膜流动多在水膜较厚的位置,在水膜厚度较小的位置,水膜流动较少,因此水膜流动非常不均匀,但是在现有技术的结冰模型中,得到的水膜流动是较为均匀的,因此在物面水膜厚度差异较大的区域,结冰模拟数值的误差较大、精度较低。
因此,为了解决上述缺陷,本申请实施例提供了一种考虑溪流现象的结冰模拟方法,以水膜平均速度
Figure 5659DEST_PATH_IMAGE001
、水膜厚度h和溪流现象作用力
Figure 902071DEST_PATH_IMAGE002
作为未知数建立欧拉水膜模型的动量方程,并以水膜平均速度
Figure 432409DEST_PATH_IMAGE001
,水膜厚度h、结冰质量流率
Figure 451181DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 332549DEST_PATH_IMAGE004
作为未知数建立欧拉水膜模型的质量方程;以结冰质量流率
Figure 23863DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 939866DEST_PATH_IMAGE004
作为未知数建立Messinger模型的热力学平衡方程;根据上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的结冰质量流率
Figure 598381DEST_PATH_IMAGE003
和蒸发质量流率
Figure 170308DEST_PATH_IMAGE004
;获取当前时刻的溪流现象作用力
Figure 470839DEST_PATH_IMAGE002
,并根据当前时刻的结冰质量流率
Figure 241349DEST_PATH_IMAGE003
、蒸发质量流率
Figure 274027DEST_PATH_IMAGE004
和溪流现象作用力
Figure 129988DEST_PATH_IMAGE002
,结合欧拉水膜模型的动量方程和质量方程,获取当前时刻的水膜平均速度
Figure 234210DEST_PATH_IMAGE001
、水膜厚度h,通过迭代计算获取每个时刻的结冰质量流率
Figure 796909DEST_PATH_IMAGE003
。如此,可以有效解决现有技术在对飞机结冰时的水膜流动进行模拟时,存在的结冰模拟数值误差较大、精度较低的问题。
下面对本申请所涉及到的考虑溪流现象的结冰模拟方法进行介绍。
实施例1
请参照图1,图1为本申请实施例1提供的一种考虑溪流现象的结冰模拟方法的流程示意图。应说明的是:本申请方法步骤的标号并非为了限制其顺序,而是为了区分不同的步骤。
该考虑溪流现象的结冰模拟方法可以包括以下步骤:
步骤S110:以水膜平均速度
Figure 62806DEST_PATH_IMAGE001
、水膜厚度h和溪流现象作用力
Figure 406062DEST_PATH_IMAGE002
作为未知数建立欧拉水膜模型的动量方程,并以所述水膜平均速度
Figure 248729DEST_PATH_IMAGE001
,所述水膜厚度h、结冰质量流率
Figure 728252DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 430628DEST_PATH_IMAGE004
作为未知数建立欧拉水膜模型的质量方程。
在绝对光滑的物面上,沿着物面流动的水膜应该是均匀、平齐的,但是由于水膜接触角的作用,实际上水膜流动是参差不齐的,如图2所示,图2为溪流现象的示意图,210为连续的水膜的流动方向,220为水膜流动的分离点,230为水膜分离后形成的溪流,240为没有水膜的干区。而现有技术按照理想的均匀、平齐的水膜流动来构建结冰模型,得到的结冰模拟效果较差,精度较低。在分析现有技术的缺陷的基础上,发明人发现了是溪流现象造成了结冰模型精度较低,因此,发明人在构建结冰模型时,考虑了溪流现象作用力,将溪流现象作用力
Figure 995602DEST_PATH_IMAGE002
作为未知数建立欧拉水膜方程。
在示例性实施例中,溪流现象作用力
Figure 379310DEST_PATH_IMAGE002
的表达式可以是:
Figure 978918DEST_PATH_IMAGE005
,其中,q为差异参数,
Figure 586617DEST_PATH_IMAGE006
为水膜接触角,w为湿面积分数。
作为一种可选实施方式,湿面积分数w可以反应水膜厚度对水膜流动的影响。
在示例性实施例中,所述欧拉水膜模型的动量方程为:
Figure 842149DEST_PATH_IMAGE007
Figure 826286DEST_PATH_IMAGE008
,且
Figure 280401DEST_PATH_IMAGE009
所述欧拉水膜模型的质量方程为:
Figure 262263DEST_PATH_IMAGE010
其中,
Figure 67408DEST_PATH_IMAGE011
为水的密度,
Figure 795849DEST_PATH_IMAGE012
为水膜平均速度采用二次剖面形式带来的对流修正项,且
Figure 104470DEST_PATH_IMAGE012
Figure 319551DEST_PATH_IMAGE001
的函数,P为水膜相关作用力,
Figure 284096DEST_PATH_IMAGE013
为空气对水膜的压力,
Figure 875614DEST_PATH_IMAGE014
为水膜表面张力,
Figure 38742DEST_PATH_IMAGE015
为空气剪切力,
Figure 362407DEST_PATH_IMAGE016
为撞击质量流率,
Figure 142145DEST_PATH_IMAGE017
为水的动力粘性系数。
本申请实施例在结冰模型中加入欧拉水膜模型,考虑了水膜流动对结冰的影响,建立的结冰模型更加精确,以及,在欧拉水膜模型中还考虑了空气对水膜的压力梯度、水膜表面张力梯度,计算得到的结冰量更精确。
步骤S120:以所述结冰质量流率
Figure 6195DEST_PATH_IMAGE003
以及蒸发质量流率
Figure 227092DEST_PATH_IMAGE004
作为未知数建立Messinger模型的热力学平衡方程。
在示例性实施例中,所述Messinger模型的热力学平衡方程包括:
Figure 783975DEST_PATH_IMAGE018
Figure 785430DEST_PATH_IMAGE019
其中,
Figure 653504DEST_PATH_IMAGE020
为结冰放热,
Figure 791224DEST_PATH_IMAGE021
为气动热,
Figure 519009DEST_PATH_IMAGE022
为水滴动能转化热,
Figure 945442DEST_PATH_IMAGE023
为水膜流入流出所带走的能量,
Figure 416875DEST_PATH_IMAGE024
为气流与物面的对流换热,
Figure 612364DEST_PATH_IMAGE025
蒸发吸热,
Figure 511050DEST_PATH_IMAGE026
为水滴自身能量,且
Figure 221517DEST_PATH_IMAGE027
Figure 434323DEST_PATH_IMAGE028
Figure 546636DEST_PATH_IMAGE029
Figure 616223DEST_PATH_IMAGE030
Figure 31897DEST_PATH_IMAGE031
Figure 845132DEST_PATH_IMAGE032
Figure 811951DEST_PATH_IMAGE033
Figure 990123DEST_PATH_IMAGE034
Figure 675182DEST_PATH_IMAGE035
Figure 292108DEST_PATH_IMAGE036
其中,Lf为融解潜热,htc为气流与物面的对流换热系数,
Figure 316696DEST_PATH_IMAGE037
为恢复因子,
Figure 728086DEST_PATH_IMAGE038
为气流速度,
Figure 900441DEST_PATH_IMAGE039
为空气比热容,
Figure 993162DEST_PATH_IMAGE016
为撞击质量流率,V为水滴相对物面的速度,
Figure 934573DEST_PATH_IMAGE040
为当前微元流出的质量流率,
Figure 516864DEST_PATH_IMAGE041
为当前微元的下游温度,
Figure 111269DEST_PATH_IMAGE042
为当前微元流入的质量流率,
Figure 804419DEST_PATH_IMAGE043
为当前微元的上游温度,
Figure 865915DEST_PATH_IMAGE044
为空气来流温度,Ld为蒸发潜热,
Figure 822370DEST_PATH_IMAGE045
为水的比热容,
Figure 703739DEST_PATH_IMAGE046
为空气密度,
Figure 200579DEST_PATH_IMAGE047
为水蒸气气体常数,
Figure 851003DEST_PATH_IMAGE048
为传热传质类比准则数,
Figure 978359DEST_PATH_IMAGE049
为物面温度。
气流与物面的对流换热系数htc可以通过等效粗糙度模型修正后的SA(Spalart-Allmaras)湍流模型计算得到。恢复因子
Figure 612603DEST_PATH_IMAGE050
可以取0.9。撞击质量流率
Figure 850817DEST_PATH_IMAGE016
可以由水滴流场计算得到。水蒸气气体常数
Figure 355748DEST_PATH_IMAGE047
可以取461.4。传热传质类比准则数
Figure 716322DEST_PATH_IMAGE048
可以取1。
步骤S130:根据上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的结冰质量流率
Figure 572283DEST_PATH_IMAGE003
和蒸发质量流率
Figure 351538DEST_PATH_IMAGE004
在示例性实施例中,步骤S130还可以包括:根据所述上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的物面温度
Figure 242134DEST_PATH_IMAGE051
、结冰质量流率
Figure 508030DEST_PATH_IMAGE003
和蒸发质量流率
Figure 788970DEST_PATH_IMAGE004
具体地,可以是:
根据上一时刻的水膜厚度h,获取当前时刻的微元水膜提供的质量流率
Figure 696883DEST_PATH_IMAGE052
,且
Figure 379668DEST_PATH_IMAGE053
Figure 285308DEST_PATH_IMAGE054
为时间步长。
在本申请实施例中,微元水膜可以为物面上划分的一个网格单元的水膜,时间步长可以为每2个时刻之间的时间间隔。
在本申请实施例中,可以通过Messinger模型的热力学平衡方程看出,
Figure 115860DEST_PATH_IMAGE055
Figure 233989DEST_PATH_IMAGE051
的函数,因此,该热力学平衡方程可以看做关于
Figure 99177DEST_PATH_IMAGE056
Figure 907208DEST_PATH_IMAGE051
的二元方程。并且,
Figure 225057DEST_PATH_IMAGE056
Figure 146877DEST_PATH_IMAGE051
需要同时满足下述关系:
当物面温度
Figure 866571DEST_PATH_IMAGE057
,撞击到物面的水滴不结冰,结冰质量流率为最小值,即
Figure 848434DEST_PATH_IMAGE058
当物面温度
Figure 653579DEST_PATH_IMAGE059
,撞击到物面的水滴可能部分结冰,另一部分形成水膜,也可能不结冰,还可能全部结冰,即
Figure 175827DEST_PATH_IMAGE060
当物面温度
Figure 687711DEST_PATH_IMAGE061
,撞击到物面的水滴会全部结冰,结冰质量流率到达最大值,即
Figure 902791DEST_PATH_IMAGE062
,且此时不再有液态水蒸发吸热,因此蒸发质量流率
Figure 867336DEST_PATH_IMAGE063
根据该热力学平衡方程和上述
Figure 458855DEST_PATH_IMAGE056
Figure 816456DEST_PATH_IMAGE051
的关系,可以使用假设修正法得到
Figure 202438DEST_PATH_IMAGE056
Figure 919858DEST_PATH_IMAGE051
所述Messinger模型的假设修正法为:
1.将当前时刻的物面温度的假设值
Figure 783909DEST_PATH_IMAGE064
设置为273.15K,并根据所述Messinger模型的热力学平衡方程,获取当前时刻的结冰质量流率的假设值
Figure 67123DEST_PATH_IMAGE065
在本申请实施例中,将当前时刻的物面温度假设为
Figure 827268DEST_PATH_IMAGE064
Figure 828723DEST_PATH_IMAGE064
只是计算过程中的物面温度的假设值,实际要得到的是修正后的物面温度
Figure 434147DEST_PATH_IMAGE051
的值。
将根据假设值
Figure 837447DEST_PATH_IMAGE064
得到的当前时刻的结冰质量流率的假设值记为
Figure 768494DEST_PATH_IMAGE065
Figure 726085DEST_PATH_IMAGE065
只是计算过程中的结冰质量流率的假设值,实际要得到的是修正后的结冰质量流率
Figure 463097DEST_PATH_IMAGE056
的值。
具体地,当知道当前时刻的物面温度
Figure 655657DEST_PATH_IMAGE056
的假设值
Figure 554343DEST_PATH_IMAGE066
时,可以求解热力学平衡方程,即求解已知一个未知数的二元方程,从而得到当前时刻的结冰质量流率的假设值
Figure 468072DEST_PATH_IMAGE065
2.若
Figure 743196DEST_PATH_IMAGE067
,则在当前时刻,
Figure 793191DEST_PATH_IMAGE068
,并得到当前时刻的
Figure 862778DEST_PATH_IMAGE069
,根据所述Messinger模型的热力学平衡方程、当前时刻的
Figure 60541DEST_PATH_IMAGE056
的值,和当前时刻的
Figure 77039DEST_PATH_IMAGE055
的值,得到当前时刻的
Figure 309437DEST_PATH_IMAGE070
的值。
具体地,若当前时刻的结冰质量流率的假设值
Figure 487609DEST_PATH_IMAGE065
大于当前时刻的
Figure 438247DEST_PATH_IMAGE016
与当前时刻的
Figure 995786DEST_PATH_IMAGE052
之和,则说明:将物面温度假设为273.15K时,根据前述
Figure 82691DEST_PATH_IMAGE056
Figure 431764DEST_PATH_IMAGE051
的关系,实际上结冰质量流率最大值为
Figure 604119DEST_PATH_IMAGE071
,即实际上的结冰放热的最大值为
Figure 962419DEST_PATH_IMAGE072
,但实际上的结冰质量流率的最大值小于通过假设计算得出的
Figure 903830DEST_PATH_IMAGE065
,即物面温度达到273.15K所需要的结冰放热
Figure 689384DEST_PATH_IMAGE073
会由于实际上的结冰质量流率的限制而供给不足,这种情况下,实际上的结冰质量流率只能等于
Figure 349035DEST_PATH_IMAGE071
,从而推导出实际上的物面温度小于273.15K。另一方面,也可以从Messinger模型的热力学平衡方程看出,结冰质量流率和物面温度是正相关关系,当
Figure 42185DEST_PATH_IMAGE065
过于大,以至于大于实际上的结冰质量流率的最大值,在修正结冰质量流率时,需要将结冰质量流率相对于其假设值
Figure 41365DEST_PATH_IMAGE065
减小,此时修正后的物面温度相对于其假设值
Figure 60136DEST_PATH_IMAGE064
也会减小,即
Figure 876258DEST_PATH_IMAGE074
Figure 373099DEST_PATH_IMAGE075
。因此,将当前时刻的结冰质量流率的值修正为
Figure 289102DEST_PATH_IMAGE071
,即修正后的当前时刻的结冰质量流率
Figure 416458DEST_PATH_IMAGE075
,此时,撞击到物面的水滴全部结冰,不再有液态水蒸发吸热,可以得到当前时刻的蒸发质量流率
Figure 50702DEST_PATH_IMAGE069
,此时,又可以根据修正后的当前时刻的结冰质量流率
Figure 85654DEST_PATH_IMAGE076
的值,求解Messinger模型的热力学平衡方程,即求解已知一个未知数的二元方程,从而得到修正后的当前时刻的物面温度
Figure 793847DEST_PATH_IMAGE070
的值。
3.若
Figure 154421DEST_PATH_IMAGE077
,则在当前时刻,
Figure 948065DEST_PATH_IMAGE078
Figure DEST_PATH_IMAGE079
,并根据所述Messinger模型的热力学平衡方程、当前时刻的
Figure 521128DEST_PATH_IMAGE080
的值,和当前时刻的
Figure 364056DEST_PATH_IMAGE081
的值,得到当前时刻的
Figure 895531DEST_PATH_IMAGE082
的值。
具体地,若当前时刻的结冰质量流率
Figure 176471DEST_PATH_IMAGE080
的假设值大于等于0,且当前时刻的结冰质量流率
Figure 818805DEST_PATH_IMAGE080
的假设值小于等于当前时刻的
Figure 501590DEST_PATH_IMAGE083
与当前时刻的
Figure 203967DEST_PATH_IMAGE084
之和,则说明:此时达到热力学平衡,无需对结冰质量流率和物面温度进行修正。因此,将当前时刻的结冰质量流率的假设值
Figure 768940DEST_PATH_IMAGE085
作为当前时刻的结冰质量流率的值,即当前时刻的结冰质量流率
Figure 152648DEST_PATH_IMAGE086
,将当前时刻的物面温度的假设值
Figure 17836DEST_PATH_IMAGE087
作为当前时刻的物面温度的值,即当前时刻的物面温度
Figure 828797DEST_PATH_IMAGE088
。然后,再根据Messinger模型的热力学平衡方程得到当前时刻的蒸发质量流率
Figure 615488DEST_PATH_IMAGE082
的值。
4.若
Figure 865203DEST_PATH_IMAGE089
,则在当前时刻,
Figure 519651DEST_PATH_IMAGE090
,根据所述Messinger模型的热力学平衡方程、当前时刻的
Figure 298251DEST_PATH_IMAGE080
的值,得到当前时刻的
Figure 837817DEST_PATH_IMAGE081
的值和当前时刻的
Figure 563328DEST_PATH_IMAGE082
的值。
具体地,若当前时刻的结冰质量流率
Figure 137528DEST_PATH_IMAGE080
的假设值小于0,则说明:将物面温度假设为273.15K时,根据前述
Figure 87030DEST_PATH_IMAGE080
Figure 317154DEST_PATH_IMAGE081
的关系,实际上结冰质量流率的最小值为0,即实际上的结冰放热的最小值为
Figure 643093DEST_PATH_IMAGE091
,但实际上的结冰质量流率的最小值大于通过假设计算得出的
Figure 71800DEST_PATH_IMAGE085
,即物面温度达到273.15K所需要的结冰放热
Figure 395465DEST_PATH_IMAGE092
供给过剩,这种情况下,实际上的结冰质量流率只能等于0,从而推导出实际上的物面温度大于273.15K。另一方面,也可以从Messinger模型的热力学平衡方程看出,结冰质量流率和物面温度是正相关关系,当
Figure 909623DEST_PATH_IMAGE085
过于小,以至于小于实际上的结冰质量流率的最小值,在修正结冰质量流率时,需要将结冰质量流率相对于其假设值
Figure 773674DEST_PATH_IMAGE085
增大,此时修正后的物面温度相对于其假设值
Figure 263080DEST_PATH_IMAGE087
也会增大,即
Figure 819963DEST_PATH_IMAGE093
Figure 759100DEST_PATH_IMAGE094
。因此,将当前时刻的结冰质量流率的值修正为
Figure 426842DEST_PATH_IMAGE095
,即修正后的当前时刻的结冰质量流率
Figure 830142DEST_PATH_IMAGE094
,此时,又可以根据修正后的当前时刻的结冰质量流率
Figure 495609DEST_PATH_IMAGE080
的值,求解Messinger模型的热力学平衡方程,即求解已知一个未知数的二元方程,从而得到修正后的当前时刻的物面温度
Figure 984360DEST_PATH_IMAGE081
的值,并得到当前时刻的蒸发质量流率
Figure 659055DEST_PATH_IMAGE082
的值。
步骤S140:获取当前时刻的溪流现象作用力
Figure 916861DEST_PATH_IMAGE096
,并根据所述当前时刻的结冰质量流率
Figure 753229DEST_PATH_IMAGE097
、蒸发质量流率
Figure 463697DEST_PATH_IMAGE082
和溪流现象作用力
Figure 735890DEST_PATH_IMAGE096
,结合所述欧拉水膜模型的动量方程和质量方程,获取当前时刻的水膜平均速度
Figure 785886DEST_PATH_IMAGE098
、水膜厚度h,通过迭代计算获取每个时刻的结冰质量流率
Figure 855473DEST_PATH_IMAGE097
,其中,初始时刻的水膜厚度h为0。
在示例性实施例中,若湿面积分数w反应水膜厚度对水膜流动的影响,则步骤S140可以包括子步骤S141:
子步骤S141:判断所述上一时刻的水膜厚度h是否大于预设厚度标准
Figure 53236DEST_PATH_IMAGE099
;若大于,则所述当前时刻溪流现象作用力
Figure 69734DEST_PATH_IMAGE100
,且w=1;若不大于,则所述当前时刻溪流现象作用力
Figure 770974DEST_PATH_IMAGE100
,且w=0。
示例性地,若上一时刻的水膜厚度h大于预设厚度标准
Figure 277041DEST_PATH_IMAGE099
,则说明水膜处于润湿状态,水膜更容易流动,根据结冰实验经验,w取1;若上一时刻的水膜厚度h小于预设厚度标准
Figure 696521DEST_PATH_IMAGE101
,则说明水膜处于干燥状态,水膜更难流动,根据结冰实验经验,w取0。
示例性地,若物面为机翼表面,如图3所示,图3为物面位置与水膜厚度的示意图,其中,横坐标为:物面位置,具体指物面位置距离机翼前缘厚度方向的长度,单位为米;纵坐标为:水膜厚度,单位为米;曲线310指:水膜接触角为120°时,水膜厚度与物面位置的关系;曲线320指:水膜接触角为0°时,水膜厚度与物面位置的关系。如图4所示,图4为物面位置与结冰率的示意图,其中,横坐标为:物面位置,具体指物面位置距离机翼前缘厚度方向的长度,单位为米;纵坐标为:结冰率,单位为千克每平方米每秒;曲线410指:水膜接触角为120°时,水膜厚度与结冰率的关系;曲线420指:水膜接触角为0°时,水膜厚度与结冰率的关系。
本申请实施例提供的考虑溪流现象的结冰模拟方法,在欧拉水膜模型中将溪流现象作用力作为未知数,在进行结冰模拟时考虑溪流现象对结冰模型的影响,相比现有技术,结冰模拟数值误差较小、精度较高。
实施例2
请参考图5,图5示出了本申请实施例2提供的一种计算机可读存储介质的结构框图。该计算机可读存储介质500中存储有程序代码510,所述程序代码510可被处理器调用执行上述方法实施例中所描述的方法。
计算机可读存储介质500可以是诸如闪存、EEPROM(电可擦除可编程只读存储器)、EPROM(可擦除可编程只读存储器)、硬盘或者ROM之类的电子存储器。可选地,计算机可读存储介质500包括非易失性计算机可读存储介质(non-transitory computer-readablestorage medium)。计算机可读存储介质500具有执行上述方法中的任何方法步骤的程序代码510的存储空间。这些程序代码可以从一个或者多个计算机程序产品中读取或者写入到这一个或者多个计算机程序产品中。程序代码510可以例如以适当形式进行压缩。
最后应说明的是:以上实施例仅用以说明本申请的技术方案,而非对其限制;尽管参照前述实施例对本申请进行了详细的说明,本领域的普通技术人员当理解:其依然可以对前述各实施例所记载的技术方案进行修改,或者对其中部分技术特征进行等同替换;而这些修改或者替换,并不驱使相应技术方案的本质脱离本申请各实施例技术方案的精神和范围。

Claims (6)

1.一种考虑溪流现象的结冰模拟方法,其特征在于,所述方法包括:
S110.以水膜平均速度
Figure DEST_PATH_IMAGE001
、水膜厚度h和溪流现象作用力
Figure 330330DEST_PATH_IMAGE002
作为未知数建立欧拉水膜模型的动量方程,并以所述水膜平均速度
Figure DEST_PATH_IMAGE003
,所述水膜厚度h、结冰质量流率
Figure 855990DEST_PATH_IMAGE004
以及蒸发质量流率
Figure DEST_PATH_IMAGE005
作为未知数建立欧拉水膜模型的质量方程;
S120.以所述结冰质量流率
Figure 731542DEST_PATH_IMAGE004
以及蒸发质量流率
Figure 557415DEST_PATH_IMAGE005
作为未知数建立Messinger模型的热力学平衡方程;
S130.根据上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的结冰质量流率
Figure 922538DEST_PATH_IMAGE004
和蒸发质量流率
Figure 342018DEST_PATH_IMAGE005
S140.获取当前时刻的溪流现象作用力
Figure 552419DEST_PATH_IMAGE002
,并根据所述当前时刻的结冰质量流率
Figure 373745DEST_PATH_IMAGE004
、蒸发质量流率
Figure 909768DEST_PATH_IMAGE005
和溪流现象作用力
Figure 675599DEST_PATH_IMAGE002
,结合所述欧拉水膜模型的动量方程和质量方程,获取当前时刻的水膜平均速度
Figure 565058DEST_PATH_IMAGE001
、水膜厚度h,通过迭代计算获取每个时刻的结冰质量流率
Figure 631103DEST_PATH_IMAGE004
,其中,初始时刻的水膜厚度h为0。
2.根据权利要求1所述的考虑溪流现象的结冰模拟方法,其特征在于,所述溪流现象作用力
Figure 947814DEST_PATH_IMAGE006
,其中,q为差异参数,
Figure DEST_PATH_IMAGE007
为水膜接触角,w为湿面积分数。
3.根据权利要求2所述的考虑溪流现象的结冰模拟方法,其特征在于,所述获取当前时刻的溪流现象作用力
Figure 947081DEST_PATH_IMAGE002
包括:
判断所述上一时刻的水膜厚度h是否大于预设厚度标准
Figure 499285DEST_PATH_IMAGE008
若大于,则所述当前时刻溪流现象作用力
Figure DEST_PATH_IMAGE009
,且w=1;
若不大于,则所述当前时刻溪流现象作用力
Figure 623099DEST_PATH_IMAGE009
,且w=0。
4.根据权利要求3所述的考虑溪流现象的结冰模拟方法,其特征在于,所述欧拉水膜模型的动量方程为:
Figure 235345DEST_PATH_IMAGE010
Figure DEST_PATH_IMAGE011
,且
Figure 444610DEST_PATH_IMAGE012
所述欧拉水膜模型的质量方程为:
Figure DEST_PATH_IMAGE013
其中,
Figure 534926DEST_PATH_IMAGE014
为水的密度,
Figure DEST_PATH_IMAGE015
为水膜平均速度采用二次剖面形式带来的对流修正项,且
Figure 513246DEST_PATH_IMAGE015
Figure 296394DEST_PATH_IMAGE016
的函数,P为水膜相关作用力,
Figure DEST_PATH_IMAGE017
为空气对水膜的压力,
Figure 992955DEST_PATH_IMAGE018
为水膜表面张力,
Figure DEST_PATH_IMAGE019
为空气剪切力,
Figure 355803DEST_PATH_IMAGE020
为撞击质量流率,
Figure DEST_PATH_IMAGE021
为水的动力粘性系数。
5.根据权利要求4所述的考虑溪流现象的结冰模拟方法,其特征在于,所述Messinger模型的热力学平衡方程包括:
Figure 578843DEST_PATH_IMAGE022
Figure DEST_PATH_IMAGE023
其中,
Figure 1734DEST_PATH_IMAGE024
为结冰放热,
Figure DEST_PATH_IMAGE025
为气动热,
Figure 185591DEST_PATH_IMAGE026
为水滴动能转化热,
Figure DEST_PATH_IMAGE027
为水膜流入流出所带走的能量,
Figure 352130DEST_PATH_IMAGE028
为气流与物面的对流换热,
Figure DEST_PATH_IMAGE029
蒸发吸热,
Figure 305042DEST_PATH_IMAGE030
为水滴自身能量,且
Figure DEST_PATH_IMAGE031
Figure 898835DEST_PATH_IMAGE032
Figure DEST_PATH_IMAGE033
Figure 429042DEST_PATH_IMAGE034
Figure DEST_PATH_IMAGE035
Figure 399272DEST_PATH_IMAGE036
Figure DEST_PATH_IMAGE037
Figure 472270DEST_PATH_IMAGE038
Figure DEST_PATH_IMAGE039
Figure 236964DEST_PATH_IMAGE040
其中,Lf为融解潜热,htc为气流与物面的对流换热系数,
Figure DEST_PATH_IMAGE041
为恢复因子,
Figure 395413DEST_PATH_IMAGE042
为气流速度,
Figure 434913DEST_PATH_IMAGE043
为空气比热容,
Figure 768942DEST_PATH_IMAGE020
为撞击质量流率,V为水滴相对物面的速度,
Figure 235696DEST_PATH_IMAGE044
为当前微元流出的质量流率,
Figure DEST_PATH_IMAGE045
为当前微元的下游温度,
Figure 881441DEST_PATH_IMAGE046
为当前微元流入的质量流率,
Figure DEST_PATH_IMAGE047
为当前微元的上游温度,
Figure 927894DEST_PATH_IMAGE048
为空气来流温度,Ld为蒸发潜热,
Figure DEST_PATH_IMAGE049
为水的比热容,
Figure 834539DEST_PATH_IMAGE050
为空气密度,
Figure DEST_PATH_IMAGE051
为水蒸气气体常数,
Figure 206615DEST_PATH_IMAGE052
为传热传质类比准则数,
Figure DEST_PATH_IMAGE053
为物面温度;
所述步骤S130包括:
根据所述上一时刻的水膜厚度h,采用所述Messinger模型的假设修正法,获取当前时刻的物面温度
Figure 74076DEST_PATH_IMAGE053
、结冰质量流率
Figure 455379DEST_PATH_IMAGE004
和蒸发质量流率
Figure 764001DEST_PATH_IMAGE005
6.一种计算机可读存储介质,其特征在于,所述计算机可读存储介质中存储有程序代码,所述程序代码可被处理器调用执行所述权利要求1-5任一项所述方法。
CN202211003148.1A 2022-08-22 2022-08-22 一种考虑溪流现象的结冰模拟方法和存储介质 Active CN115081121B (zh)

Priority Applications (1)

Application Number Priority Date Filing Date Title
CN202211003148.1A CN115081121B (zh) 2022-08-22 2022-08-22 一种考虑溪流现象的结冰模拟方法和存储介质

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
CN202211003148.1A CN115081121B (zh) 2022-08-22 2022-08-22 一种考虑溪流现象的结冰模拟方法和存储介质

Publications (2)

Publication Number Publication Date
CN115081121A true CN115081121A (zh) 2022-09-20
CN115081121B CN115081121B (zh) 2022-11-01

Family

ID=83244692

Family Applications (1)

Application Number Title Priority Date Filing Date
CN202211003148.1A Active CN115081121B (zh) 2022-08-22 2022-08-22 一种考虑溪流现象的结冰模拟方法和存储介质

Country Status (1)

Country Link
CN (1) CN115081121B (zh)

Cited By (1)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116562192A (zh) * 2023-07-06 2023-08-08 中国空气动力研究与发展中心计算空气动力研究所 一种飞机结冰冰形预测方法、装置、设备及存储介质

Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20080190853A1 (en) * 2004-01-27 2008-08-14 Alberta Research Council Inc. Method and Apparatus For Separating Liquid Droplets From a Gas Stream
CN102663238A (zh) * 2012-03-21 2012-09-12 南京航空航天大学 基于液态水分布的结冰表面粗糙度衡量方法
CN102682144A (zh) * 2011-11-30 2012-09-19 天津空中代码工程应用软件开发有限公司 直升机旋翼飞行结冰的数值模拟方法
US20140257770A1 (en) * 2011-11-30 2014-09-11 Ming Lu Numerical simulation method for the flight-icing of helicopter rotary-wings
US20140257771A1 (en) * 2011-11-30 2014-09-11 Ming Lu Numerical simulation method for aircrasft flight-icing
CN109060295A (zh) * 2018-08-29 2018-12-21 南京航空航天大学 多功能水膜发生实验装置及方法

Patent Citations (6)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20080190853A1 (en) * 2004-01-27 2008-08-14 Alberta Research Council Inc. Method and Apparatus For Separating Liquid Droplets From a Gas Stream
CN102682144A (zh) * 2011-11-30 2012-09-19 天津空中代码工程应用软件开发有限公司 直升机旋翼飞行结冰的数值模拟方法
US20140257770A1 (en) * 2011-11-30 2014-09-11 Ming Lu Numerical simulation method for the flight-icing of helicopter rotary-wings
US20140257771A1 (en) * 2011-11-30 2014-09-11 Ming Lu Numerical simulation method for aircrasft flight-icing
CN102663238A (zh) * 2012-03-21 2012-09-12 南京航空航天大学 基于液态水分布的结冰表面粗糙度衡量方法
CN109060295A (zh) * 2018-08-29 2018-12-21 南京航空航天大学 多功能水膜发生实验装置及方法

Non-Patent Citations (6)

* Cited by examiner, † Cited by third party
Title
CHENGXIANG ZHU: "3D ICE ACCRETION SIMULATION FOR COMPLEX CONFIGURATION BASING ON IMPROVED MESSINGER MODEL", 《WORLD SCIENTIFIC》 *
ZHIHONG ZHOU ET AL: "Icing Numerical Simulation for Single and Multi-Element Airfoils", 《28TH AIAA APPLIED AERODYNAMICS CONFERENCE》 *
廖佩琳等: "流域源区溪流CO_2的来源与扩散过程研究综述", 《山地学报》 *
易贤等: "考虑相变时间效应的结冰试验相似参数", 《实验流体力学》 *
王方年等: "GASFLOW程序液膜模型开发及初步验证", 《原子能科学技术》 *
郑梅等: "来流速度对防冰表面溢流水流动换热的影响", 《空气动力学学报》 *

Cited By (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
CN116562192A (zh) * 2023-07-06 2023-08-08 中国空气动力研究与发展中心计算空气动力研究所 一种飞机结冰冰形预测方法、装置、设备及存储介质
CN116562192B (zh) * 2023-07-06 2023-09-12 中国空气动力研究与发展中心计算空气动力研究所 一种飞机结冰冰形预测方法、装置、设备及存储介质

Also Published As

Publication number Publication date
CN115081121B (zh) 2022-11-01

Similar Documents

Publication Publication Date Title
Fortin et al. Heat and mass transfer during ice accretion on aircraft wings with an improved roughness model
Croce et al. FENSAP-ICE: Analytical model for spatial and temporal evolution of in-flight icing roughness
Fortin et al. New roughness computation method and geometric accretion model for airfoil icing
US11161629B2 (en) System for numerical simulation and test verification of icing characteristics of an aerostat
CN115081121B (zh) 一种考虑溪流现象的结冰模拟方法和存储介质
CN112989727B (zh) 一种防冰系统的壁面温度模拟方法
Ozcer et al. FENSAP-ICE: Numerical prediction of ice roughness evolution, and its effects on ice shapes
CN113255067B (zh) 考虑挥舞作用的旋翼结冰冰形计算方法、存储介质及终端
KR20140115835A (ko) 공기 흡입구의 표면 결빙 시뮬레이션 방법
CN112985347A (zh) 一种飞机结冰表面粗糙度计算方法
Fujiwara et al. Comparison of computational and experimental ice accretions of large swept wings
CN115292806A (zh) 考虑周期性边界的三维热气防冰系统表面温度计算方法
Homola et al. Turbine size and temperature dependence of icing on wind turbine blades
CN115099048A (zh) 一种基于Messinger模型的结冰数值模拟方法、系统、设备及介质
Dong et al. Experimental investigation on anti-icing performance of an engine inlet strut
Löwe et al. Inception of ice accretion by ice crystal impact
Knobbe-Eschen et al. Numerical and experimental investigations of wind-turbine blade aerodynamics in the presence of ice accretion
WO2019186151A1 (en) Methods and apparatus for simulating liquid collection on aerodynamic components
Addy et al. A study of the effects of altitude on thermal ice protection system performance
Ghenai et al. Verification and validation of NASA LEWICE 2.2 icing software code
Anderson et al. Overview of icing physics relevant to scaling
Prikhod’ko et al. Numerical Simulation of the Processes of Icing on Airfoils with Formation of a “Barrier” Ice
Alegre Full-scale runback ice: accretion and aerodynamic study
Addy Jr et al. Altitude Effects on Thermal Ice Protection System Performance; a Study of an Alternative Approach
Abbadi et al. Preliminary analysis of ice accretion prediction on wind turbine blades

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