工程科学学报，第40卷，第1期：99-107,2018年1月 Chinese Journal of Engineering,Vol.40,No.1:99-107,January 2018 DOI:10.13374/j.issn2095-9389.2018.01.013;http://journals.ustb.edu.cn 基于单相LBM模拟大平板反重力充型过程 包羽冲，刘林，骈松，张照，李日⑧ 河北工业大学材料科学与工程学院.天津300130 ☒通信作者，E-mail:sdzllr(@163.com 摘要采用单相格子Boltzmann方法研究大平板的反重力充型过程，该模型不需考虑气相格子的变化，从而提高了计算效 率.针对该方法，本文新提出了一种权重系数重新分配的方法来处理格点中的液相排出及分配问题.首先用该模型计算了单 浇口条件下的大平板型腔反重力充填过程，以相同参数下的高速相机成像的水力充填实验为参照，数值模拟的流场特征及流 体形态与实验结果吻合良好.另外，还采用线速度分布云图，并提出了自由表面的高度差判据来分析充型过程中的流场区域 特征和流体平稳程度.在此基础上，继续用该模型研究了双浇口和圆柱扰流条件下的大平板反重力充型过程.双浇口条件下 由于浇道间的互相影响，流场中形成的漩涡多于单浇口：圆柱扰流条件下的充填方式会降低流体的晃动程度，提高充型的稳 定性. 关键词反重力充型：格子Boltzmann方法：自由表面：流场：数值模拟 分类号TG249.9 Simulation of large plate castings in counter-gravity mould filling process based on single-phase LBM BAO Yu-chong,LIU Lin,PIAN Song,ZHANG Zhao,LI Ri School of Material Science and Engineering,Hebei University of Technology,Tianjin 300130,China Corresponding author,E-mail:sdzllr@163.com ABSTRACT The single-phase free surface lattice Boltzmann model(LBM)was used to investigate the counter-gravity mould filling process of large plate castings.Computational efficiency can be improved by ignoring the variation of the gas lattice.In this paper,a weight coefficient redistribution method was proposed for dealing with the liquid phase discharge problem and distribution in a grid. First,the model was used to calculate the anti-gravity filling process of a large plate cavity with a single gate,with the same parame- ters,the water filling experiment of high speed camera is taken as a reference,the characteristics of the flow field and the fluid mor- phology of numerical simulation are in good agreement with the experimental results.Furthermore,the velocity profile of the flow field is used to analyze the characteristics of the flow field,and a criterion called "the altitude difference of free surface"is introduced to e- valuate the stability of the fluid.Subsequently,the counter-gravity filling process of a plate with a double gate and the same plate with one gate and a cylindrical obstruction to turbulent flow are investigated.Due to the mutual influence between the two gates,the vortex formed by the double gates is greater than that formed by the single gate.However,the obstacles under the conditions of cylindrical tur- bulent flow obstruction could reduce the degree of fluid sloshing and improve the stability of the filling. KEY WORDS counter-gravity mould filling;lattice Boltzmann model;free surface;flow field;numerical simulation 充型过程在铸件生产过程中起着重要作用，许 多铸造缺陷如卷气、夹杂、缩孔、冷隔等与其有关 收稿日期：2017-04-28 基金项目：国家自然科学基金资助项目(51475138)
工程科学学报,第 40 卷,第 1 期:99鄄鄄107,2018 年 1 月 Chinese Journal of Engineering, Vol. 40, No. 1: 99鄄鄄107, January 2018 DOI: 10. 13374 / j. issn2095鄄鄄9389. 2018. 01. 013; http: / / journals. ustb. edu. cn 基于单相 LBM 模拟大平板反重力充型过程 包羽冲, 刘 林, 骈 松, 张 照, 李 日苣 河北工业大学材料科学与工程学院, 天津 300130 苣通信作者, E鄄mail: sdzllr@ 163. com 摘 要 采用单相格子 Boltzmann 方法研究大平板的反重力充型过程,该模型不需考虑气相格子的变化,从而提高了计算效 率. 针对该方法,本文新提出了一种权重系数重新分配的方法来处理格点中的液相排出及分配问题. 首先用该模型计算了单 浇口条件下的大平板型腔反重力充填过程,以相同参数下的高速相机成像的水力充填实验为参照,数值模拟的流场特征及流 体形态与实验结果吻合良好. 另外,还采用线速度分布云图,并提出了自由表面的高度差判据来分析充型过程中的流场区域 特征和流体平稳程度. 在此基础上,继续用该模型研究了双浇口和圆柱扰流条件下的大平板反重力充型过程. 双浇口条件下 由于浇道间的互相影响,流场中形成的漩涡多于单浇口;圆柱扰流条件下的充填方式会降低流体的晃动程度,提高充型的稳 定性. 关键词 反重力充型; 格子 Boltzmann 方法; 自由表面; 流场; 数值模拟 分类号 TG249郾 9 收稿日期: 2017鄄鄄04鄄鄄28 基金项目: 国家自然科学基金资助项目(51475138) Simulation of large plate castings in counter鄄gravity mould filling process based on single鄄phase LBM BAO Yu鄄chong, LIU Lin, PIAN Song, ZHANG Zhao, LI Ri 苣 School of Material Science and Engineering, Hebei University of Technology, Tianjin 300130, China 苣Corresponding author, E鄄mail: sdzllr@ 163. com ABSTRACT The single鄄phase free surface lattice Boltzmann model (LBM) was used to investigate the counter鄄gravity mould filling process of large plate castings. Computational efficiency can be improved by ignoring the variation of the gas lattice. In this paper, a weight coefficient redistribution method was proposed for dealing with the liquid phase discharge problem and distribution in a grid. First, the model was used to calculate the anti鄄gravity filling process of a large plate cavity with a single gate, with the same parame鄄 ters, the water filling experiment of high speed camera is taken as a reference, the characteristics of the flow field and the fluid mor鄄 phology of numerical simulation are in good agreement with the experimental results. Furthermore, the velocity profile of the flow field is used to analyze the characteristics of the flow field, and a criterion called “the altitude difference of free surface冶 is introduced to e鄄 valuate the stability of the fluid. Subsequently, the counter鄄gravity filling process of a plate with a double gate and the same plate with one gate and a cylindrical obstruction to turbulent flow are investigated. Due to the mutual influence between the two gates, the vortex formed by the double gates is greater than that formed by the single gate. However, the obstacles under the conditions of cylindrical tur鄄 bulent flow obstruction could reduce the degree of fluid sloshing and improve the stability of the filling. KEY WORDS counter鄄gravity mould filling; lattice Boltzmann model; free surface; flow field; numerical simulation 充型过程在铸件生产过程中起着重要作用,许 多铸造缺陷如卷气、夹杂、缩孔、冷隔等与其有关
·100· 工程科学学报，第40卷，第1期 铸件充型过程中的数值模拟是材料加工领域的研究 过求解粒子分布函数f的演化得到流体宏观密度 难点，国内外学者做了大量研究].反重力充型是 p= ∑f和速度u-∑ef/p.演化方程为o: 一些特种铸造工艺，如低压铸造、真空吸铸法等的充 f(x,t+△)=(1-ω)f(x,t+△t)+(1) 型方式，研究其充型过程特征对获得相应的合格优 其中：w=1/行，r为量纲一的松弛时间：f(x,t+△1)= 质铸件意义重大，采用数值模拟技术是现代最为先 f(x+△e1,t)称为迁移步；e:为粒子在i方向的速 进的研究反重力充型过程的方法 度：I表示i的反向：△t是时间步长：平衡分布函数 目前对于反重力充填的模拟研究，常用方法是 为： 通过求解宏观Navier--Stokes方程来模拟流体运动， 用体积函数法(volume of fluid,VOF)或水平集(level w+-+] (2) st)等方法来捕捉界面[4，如1997年吕衣礼等[s]用 当i=0时，10：=4/9，当i=1,2,3,4时0：=1/9，当 SOLA-VOF solution algorithm volume of fluid,SO- i=5,6,7,8时0，=1/36.固体壁面边界条件采用完 LA-VOF)方法，2000年李梅娥等[6用CG-VOF方 全反弹格式，边界处的分布函数： conjugate gradient volume of fluid,CG-VOF), f(xo,t+△t)=f(x1,t) (3) 2003年李日等[]用压力修正法(SIMPLE方法)等. 其中I表示i的反向且指向壁面，x。表示固体边界， 这些方法计算繁琐复杂，难于收敛，计算效率低，并 x1表示靠近固体边界的第一层流体格点 且上述反重力数值模拟研究仅仅做了较为简单的大 1.2单相自由表面模型 平板单浇口空腔反重力充填方式的研究，没有能够 大平板的俦造充型过程是气液两相流问题，本 涉及到更加复杂的反重力充填方式如多浇口和复杂 文采用Thurey提出的单相LBM方法[)模拟了带自 内腔型式的反重力充填研究 由表面流体的运动.如图1所示，它在标准的单相 近年来，不同于求解连续N-S方程，格子Boltz- 模型的基础上将每个格子标记为气相G(空格子)、 mann方法(lattice Boltzmann method,LBM)迅速发 液相F(满格子)或自由界面f(界面格子).满格子 展，LBM方法从分子运动论、统计力学的观点和理 区域演化采用上述标准的LBM模型，气相格子不求 论出发，以微观粒子尺度为基础，建立离散的速度模 解，大大提高了求解和计算效率.自由界面格子区 型.在满足质量、动量和能量守恒的条件下，得出粒 域做特殊处理，动态充型过程通过自由界面格子的 子分布函数，然后对粒子分布函数进行统计计算，得 移动来表现. 到压力、流速等宏观变量.该方法易于并行处理，非 常适合大规模的计算，具有高效率、高精确度等优 点.LBM方法应用于流体研究方面，成功的模拟出 气相 流体的交融、飞溅、圆柱绕流[8-)等复杂流体的运 液相 动，适合处理多相流问题及自由表面的运动[2-15] 但用LBM方法对材料加工领域的铸造反重力充型 过程的研究极少有文献报道. 图1格子类型的划分 Fig.1 Grid type classification 本文采用一种单相LBM方法进行反重力充填 数值模拟研究，对自由表面运动时的界面格点多余 1.2.1界面格子的质量交换 液相的处理提出了一种权重系数重新分配的新方 界面格子之间通过迁移步完成与非空格子之间 法.除对传统的反重力充填大平板单浇口型式进行 的质量交换，质量变化为： 研究外，还对大平板反重力充填的双浇口及圆柱扰 △m:(x,t+△t)= 流两种型式进行了研究，阐述了不同条件下的大平 (+Ate)+(+e)-f(x)] 板反重力充型过程流体内部涡流特征、自由表面的 波动、速度场的分布以及其形成原因，这是以往的大 (4) 平板反重力充填所未涉及的研究内容 其中，定义每个格子的质量分数s=m/p,e的值在 0~1之间，其中气相ε=0，液相ε=1.m为格子质 1模型描述 量，流体为不可压缩流体，P为液体密度. 1.1LBM模型 1.2.2界面分布函数的重构 本文采用基于单步松弛时间的D2Q9模型，通 由于气相格子区域不进行计算，从气相迁移到
工程科学学报,第 40 卷,第 1 期 铸件充型过程中的数值模拟是材料加工领域的研究 难点,国内外学者做了大量研究[1鄄鄄3] . 反重力充型是 一些特种铸造工艺,如低压铸造、真空吸铸法等的充 型方式,研究其充型过程特征对获得相应的合格优 质铸件意义重大,采用数值模拟技术是现代最为先 进的研究反重力充型过程的方法. 目前对于反重力充填的模拟研究,常用方法是 通过求解宏观 Navier鄄鄄 Stokes 方程来模拟流体运动, 用体积函数法(volume of fluid,VOF)或水平集(level set)等方法来捕捉界面[4] ,如 1997 年吕衣礼等[5] 用 SOLA鄄鄄VOF ( solution algorithm volume of fluid, SO鄄 LA鄄鄄VOF)方法,2000 年李梅娥等[6] 用 CG鄄鄄 VOF 方 法( conjugate gradient volume of fluid, CG鄄鄄 VOF), 2003 年李日等[7] 用压力修正法( SIMPLE 方法)等. 这些方法计算繁琐复杂,难于收敛,计算效率低,并 且上述反重力数值模拟研究仅仅做了较为简单的大 平板单浇口空腔反重力充填方式的研究,没有能够 涉及到更加复杂的反重力充填方式如多浇口和复杂 内腔型式的反重力充填研究. 近年来,不同于求解连续 N鄄鄄 S 方程,格子 Boltz鄄 mann 方法( lattice Boltzmann method, LBM )迅速发 展,LBM 方法从分子运动论、统计力学的观点和理 论出发,以微观粒子尺度为基础,建立离散的速度模 型. 在满足质量、动量和能量守恒的条件下,得出粒 子分布函数,然后对粒子分布函数进行统计计算,得 到压力、流速等宏观变量. 该方法易于并行处理,非 常适合大规模的计算,具有高效率、高精确度等优 点. LBM 方法应用于流体研究方面,成功的模拟出 流体的交融、飞溅、圆柱绕流[8鄄鄄11] 等复杂流体的运 动,适合处理多相流问题及自由表面的运动[12鄄鄄15] . 但用 LBM 方法对材料加工领域的铸造反重力充型 过程的研究极少有文献报道. 本文采用一种单相 LBM 方法进行反重力充填 数值模拟研究,对自由表面运动时的界面格点多余 液相的处理提出了一种权重系数重新分配的新方 法. 除对传统的反重力充填大平板单浇口型式进行 研究外,还对大平板反重力充填的双浇口及圆柱扰 流两种型式进行了研究,阐述了不同条件下的大平 板反重力充型过程流体内部涡流特征、自由表面的 波动、速度场的分布以及其形成原因,这是以往的大 平板反重力充填所未涉及的研究内容. 1 模型描述 1郾 1 LBM 模型 本文采用基于单步松弛时间的 D2Q9 模型,通 过求解粒子分布函数 f i 的演化得到流体宏观密度 籽 = 移 f i 和速度 u = 移 ei f i / 籽. 演化方程为[16] : f i(x,t + 驻t) = (1 - 棕)f * i (x,t + 驻t) + 棕f eq i (1) 其中:棕 =1/ 子,子 为量纲一的松弛时间;f * i (x,t + 驻t) = f i(x + 驻teI,t)称为迁移步;ei 为粒子在 i 方向的速 度;I 表示 i 的反向;驻t 是时间步长;平衡分布函数 f eq i 为: f eq i = wi 籽 [1 + 3eiu c 2 - 3 2 u 2 c 2 + 9 2 (eiu) 2 c 4 ] (2) 当 i = 0 时,wi = 4 / 9,当 i = 1,2,3,4 时 wi = 1 / 9,当 i = 5,6,7,8 时 wi = 1 / 36. 固体壁面边界条件采用完 全反弹格式,边界处的分布函数: f i(x0 ,t + 驻t) = f I(x1 ,t) (3) 其中 I 表示 i 的反向且指向壁面,x0 表示固体边界, x1 表示靠近固体边界的第一层流体格点. 1郾 2 单相自由表面模型 大平板的铸造充型过程是气液两相流问题,本 文采用 Thurey 提出的单相 LBM 方法[17]模拟了带自 由表面流体的运动. 如图 1 所示,它在标准的单相 模型的基础上将每个格子标记为气相 G(空格子)、 液相 F(满格子)或自由界面 If(界面格子). 满格子 区域演化采用上述标准的 LBM 模型,气相格子不求 解,大大提高了求解和计算效率. 自由界面格子区 域做特殊处理,动态充型过程通过自由界面格子的 移动来表现. 图 1 格子类型的划分 Fig. 1 Grid type classification 1郾 2郾 1 界面格子的质量交换 界面格子之间通过迁移步完成与非空格子之间 的质量交换,质量变化为: 驻mi(x,t + 驻t) = 着(x + 驻tei,t) + 着(x,t) 2 [f I(x + 驻tei,t) - f i(x,t)] (4) 其中,定义每个格子的质量分数 着 = m / 籽,着 的值在 0 ~ 1 之间,其中气相 着 = 0,液相 着 = 1. m 为格子质 量,流体为不可压缩流体,籽 为液体密度. 1郾 2郾 2 界面分布函数的重构 由于气相格子区域不进行计算,从气相迁移到 ·100·
包羽冲等：基于单相LBM模拟大平板反重力充型过程 ·101· 界面格子的分布函数缺失，为了使界面处分布函数 相的分布去向.相对于前两种方法的计算更为简便 保持完整，那些来源于气相格点最终应该流入界面 且易于理解，同时也更容易编程.模拟结果显示该 层格点的分布函数就必须根据宏观边界条件构造. 方法能很好的显示反重力充型的流体形态，适用于 在没有表面张力的情况下，宏观边界条件即为界面 低速流动. 处的压力边界条件.从气相迁移过来的分布函 如图2所示，0为过满或过空格点，周围格点定 数为： 义权重系数，当i=1,2,3,4时0：=1/9，当i=5,6, fi(x,t+At)=f(p:,u)+f(p:,u)-f(x,t) 7,8时心，=1/36.对于所有周围界面格点进行权重 (5) 系数加和W=∑ω：，每个界面格点分配到的质量 其中P为当地大气压力密度，“为界面格子速度， △t为一个计算步长.(5)式的物理意义为同一界面 △m= ·m“,其中m“是所有界面格点分配到的总 格子中气相和液相的总动量在充型过程中守恒. 质量. 1.2.3界面重构 自由表面的迁移是由于整体格子类型的重新初 始化实现，上一步中的所有格子类型在计算步完成 后重新定义.气体与液体之间不能直接变换，所以 界面的重构只能以界面格点变为气体或者界面格点 变为液体来进行.重构由质量与密度的关系来决 定，在1.2.1中质量交换完成后，每个界面格子的质 量为： 图2中心格子的液相向周围排出 m(x,t+△)=m(x,)+∑△m, (6) Fig.2 Liquid distribution in central grid 每个格子在充满状态时E为1，在空状态时ε为0. 1.2.5外力项 当界面格子的质量分数ε大于1时为过满状态，小 本文中大平板的反重力充型受到重力作用，式 于0时为过空状态，在这两种状态下的界面格子类 (1)右端加入外力项为： 型发生了更新.为了不让界面点在相邻时间步内变 f(x,t+△t)=(1-ω)f后(x,t+△t）+f+F: 换过于剧烈，引入常量k=103.分为如表1两种 (7) 情况. 根据动量及质量守恒，外力项表示为F:=t(e,F) 表1界面格点的演变 其中6=4/3，.234=1/3,5.67.8=1/12.重力F= Table 1 Evolution of interface grid Pg,相应的带有重力项的平衡分布函数为： 质量分数 格点状态 转变状态 排出液相质量 3e43w9(e,)23(e8)1 6>1+k 过满 液相 m(x,t+At)-p =wp1+2-2。+2c 2c 8<-k 过空 气相 m(x,t) (8) 对于过满格点周围的所有气相格点重新定义 宏观速度修正为“=（名c)/p+宁 为界面格点，这些新界面格点的宏观量为周围非 1.3浇道入口条件和边界条件 空格点宏观量的均值，分布函数为此宏观量下的 初始条件为：浇道入口速度恒定为(x,y,t= 平衡分布函数.过空界面格点周围的所有液相格 0)=0.2ms1,方向竖直向上.初始分布函数为平 点定义为界面格点.为保证质量守恒，过满和过空 衡分布函数. 界面格点排出的液相需分配到周围的所有界面格 边界条件：计算区域四周为固体壁面，采用完全 点中. 反弹格式[19-20】 1.2.4权重系数再分配 对于1.2.3中的多余液相排出，一般采用沿界 2大平板反重力充填模拟结果与水力实验 面法向量分配和沿速度方向分配的方法[].本文 结论对比 中提出了一种权重系数再分配的新方法，利用 2.1大平板反重力充填验证算例及水力实验结论 D2Q9模型中的分布函数的权重系数来处理多余液 本文选取西北工业大学国家凝固技术重点实验
包羽冲等: 基于单相 LBM 模拟大平板反重力充型过程 界面格子的分布函数缺失,为了使界面处分布函数 保持完整,那些来源于气相格点最终应该流入界面 层格点的分布函数就必须根据宏观边界条件构造. 在没有表面张力的情况下,宏观边界条件即为界面 处的压力边界条件. 从气相迁移过来的分布函 数为[18] : f I(x,t + 驻t) = f eq i (籽g,u) + f eq I (籽g,u) - f i(x,t) (5) 其中,籽g 为当地大气压力密度,u 为界面格子速度, 驻t 为一个计算步长. (5)式的物理意义为同一界面 格子中气相和液相的总动量在充型过程中守恒. 1郾 2郾 3 界面重构 自由表面的迁移是由于整体格子类型的重新初 始化实现,上一步中的所有格子类型在计算步完成 后重新定义. 气体与液体之间不能直接变换,所以 界面的重构只能以界面格点变为气体或者界面格点 变为液体来进行. 重构由质量与密度的关系来决 定,在 1郾 2郾 1 中质量交换完成后,每个界面格子的质 量为: m(x,t + 驻t) = m(x,t) + 移i 驻mi (6) 每个格子在充满状态时 着 为 1,在空状态时 着 为 0. 当界面格子的质量分数 着 大于 1 时为过满状态,小 于 0 时为过空状态,在这两种状态下的界面格子类 型发生了更新. 为了不让界面点在相邻时间步内变 换过于剧烈,引入常量 k = 10 - 3 . 分为如表 1 两种 情况. 表 1 界面格点的演变 Table 1 Evolution of interface grid 质量分数 格点状态 转变状态 排出液相质量 着 > 1 + k 过满 液相 m(x,t + 驻t) - 籽 着 < - k 过空 气相 m(x,t) 对于过满格点周围的所有气相格点重新定义 为界面格点,这些新界面格点的宏观量为周围非 空格点宏观量的均值,分布函数为此宏观量下的 平衡分布函数. 过空界面格点周围的所有液相格 点定义为界面格点. 为保证质量守恒,过满和过空 界面格点排出的液相需分配到周围的所有界面格 点中. 1郾 2郾 4 权重系数再分配 对于 1郾 2郾 3 中的多余液相排出,一般采用沿界 面法向量分配和沿速度方向分配的方法[17] . 本文 中提出了一种权重系数再分配的新方法, 利用 D2Q9 模型中的分布函数的权重系数来处理多余液 相的分布去向. 相对于前两种方法的计算更为简便 且易于理解,同时也更容易编程. 模拟结果显示该 方法能很好的显示反重力充型的流体形态,适用于 低速流动. 如图 2 所示,0 为过满或过空格点,周围格点定 义权重系数,当 i = 1,2,3,4 时 wi = 1 / 9,当 i = 5,6, 7,8 时 wi = 1 / 36. 对于所有周围界面格点进行权重 系数加和 W = 移i 棕i,每个界面格点分配到的质量 驻m * = 棕i W ·m ex ,其中 m ex是所有界面格点分配到的总 质量. 图 2 中心格子的液相向周围排出 Fig. 2 Liquid distribution in central grid 1郾 2郾 5 外力项 本文中大平板的反重力充型受到重力作用,式 (1)右端加入外力项为[1] : f i(x,t + 驻t) = (1 - 棕)f * i (x,t + 驻t) + 棕f eq i + Fi (7) 根据动量及质量守恒,外力项表示为 Fi = t * i (eiF). 其中 t * 0 = 4 / 3,t * 1,2,3,4 = 1 / 3,t * 5,6,7,8 = 1 / 12. 重力 F = 籽0 g,相应的带有重力项的平衡分布函数为: f eq i = wi 籽 [1 + 3eiu c 2 - 3 2 u 2 c 2 + 9 2 (eiu) 2 c 4 - 3(eig) 2 ] c (8) 宏观速度修正为 u = ( 移 8 i = 0 ei f i ) 籽 + 1 2 g. 1郾 3 浇道入口条件和边界条件 初始条件为:浇道入口速度恒定为 u( x,y,t = 0) = 0郾 2 m·s - 1 ,方向竖直向上. 初始分布函数为平 衡分布函数. 边界条件:计算区域四周为固体壁面,采用完全 反弹格式[19鄄鄄20] . 2 大平板反重力充填模拟结果与水力实验 结论对比 2郾 1 大平板反重力充填验证算例及水力实验结论 本文选取西北工业大学国家凝固技术重点实验 ·101·
·102· 工程科学学报，第40卷，第1期 室的单浇口薄壁大平板反重力水充填实验作为所建 (b) 立模型的验证算例2)， 如图3(a)所示，大平板高600mm、宽500mm、 厚5mm,浇口处于大平板底部中心，宽度为100mm. 充填介质为水，实验时水中放置白色示踪粒子，采 用高速摄影技术记录充填过程.图3(b)给出了充 填过程中的某一时刻充填状态，中间处的自由表 面较两侧高，内部存在明显的涡流区域，左右两侧 基本对称.图3(c)是图3(b)的流场特征图，可 图4实验和模拟的涡流特征.()单浇口大平板实验涡流特征： 见，充型过程中会形成射流区(J区)、主涡流区(J2 (b)单浇口大平板数值模拟涡流特征 Fig.4 Experimental and simulation vortex characteristics:(a)ex- 区)、前涡流区(J3区)和后涡流区(J区)四个特征 perimental characteristics of single gate large plate vortex;(b)nu- 区域 merical simulation characteristics of single gate large plate vortex 500 流从而形成前涡流区J3·由图5(b1),当t=6s时， a (c) 将右下侧拐角处放大，可看出在两个下侧角部出现 了涡旋区J4,其旋转方向与主涡流区相反，但强度 较弱.随着充型的继续，如图5(c1),当t=9s时，流 体呈现平稳向上移动的态势，自由表面中间高，两侧 低.射流区不断增大，四个侧涡流区消失，仅存在两 100 个主涡流区 图3单浇口大平板反重力充填实验.(a)计算区域（单位： a mm):(b)充填实验：(c)涡流特征示意图 Fig.3 Counter gravity filling insingle gate of large plate:(a)calcu- lation area (unit:mm);(b)filling experiment;(c)characteristic distribution of vortex 2.2对大平板单浇口反重力充型过程的计算验证 2.2.1对比实验结论和数值模拟的流场特征 数值模拟中采用与实验相同的充填物质水，相 (al) b1) c1) 同的充填速度0.2m·s-1,都采用图3(a)中的大平 板单浇口反重力充填计算区域.图4(a)为现有单 浇口大平板反重力充填实验结论给出的流场特征. 图4(b)是单相LBM的数值模拟得到的流场特征. 射流区及各部位涡流区位置和方向基本一致，证明 了本文模型的正确性. 2.2.2对比实验结论和数值模拟的充填形态 图5(al)、(bl)和(cl)为大平板单浇口反重力 充型过程的LBM数值计算结果.模拟结果的流体 形态与大平板单浇口反重力充型过程的水力实验图 图5实验和数值模拟的流体形态.(a)t=3s:(b)t=6s:(c)t=9 5(a)、(b)和(c)基本一致，结合图4的实验验证， s;(al)t=3s;(bl)t=6s;(cl)t=9s 2.2.1进一步证明了本文模型的正确性. Fig.5 Fluid morphology of experiment and numerical calculation: 在计算区域内，流动呈现左右对称的状态.图5 (a)t=3s;(b)t=6s;(c)t=9s;(al)t=3s;(bl)t=6s; (c1)t=9s (al),t=3s,为充型开始，随着射流区J1的扩大，主 涡流区J,对称出现在射流区两侧，此时自由表面中 2.2.3流场区域线速度的分布 间部位明显高于两侧，当流体撞击左右壁面后，两侧 对整体流动区域进行线速度分析有利于对 流体在两侧壁的挤压作用下迅速向上移动并产生回 卷气和夹杂成因的理解，以图5(b1)作为线速度
工程科学学报,第 40 卷,第 1 期 室的单浇口薄壁大平板反重力水充填实验作为所建 立模型的验证算例[21] . 如图 3(a)所示,大平板高 600 mm、宽 500 mm、 厚 5 mm,浇口处于大平板底部中心,宽度为 100 mm. 充填介质为水,实验时水中放置白色示踪粒子,采 用高速摄影技术记录充填过程. 图 3( b)给出了充 填过程中的某一时刻充填状态,中间处的自由表 面较两侧高,内部存在明显的涡流区域,左右两侧 基本对称. 图 3 ( c) 是图 3 ( b) 的流场特征图,可 见,充型过程中会形成射流区( J1区) 、主涡流区( J2 区) 、前涡流区( J3区)和后涡流区( J4区)四个特征 区域. 图 3 单浇口大平板反重力充填实验. ( a) 计算区域( 单位: mm);(b)充填实验;(c)涡流特征示意图 Fig. 3 Counter gravity filling insingle gate of large plate:( a) calcu鄄 lation area (unit: mm); ( b) filling experiment; ( c) characteristic distribution of vortex 2郾 2 对大平板单浇口反重力充型过程的计算验证 2郾 2郾 1 对比实验结论和数值模拟的流场特征 数值模拟中采用与实验相同的充填物质水,相 同的充填速度 0郾 2 m·s - 1 ,都采用图 3( a)中的大平 板单浇口反重力充填计算区域. 图 4( a)为现有单 浇口大平板反重力充填实验结论给出的流场特征. 图 4(b)是单相 LBM 的数值模拟得到的流场特征. 射流区及各部位涡流区位置和方向基本一致,证明 了本文模型的正确性. 2郾 2郾 2 对比实验结论和数值模拟的充填形态 图 5(a1)、(b1)和(c1)为大平板单浇口反重力 充型过程的 LBM 数值计算结果. 模拟结果的流体 形态与大平板单浇口反重力充型过程的水力实验图 5(a)、(b) 和( c) 基本一致,结合图 4 的实验验证, 2郾 2郾 1 进一步证明了本文模型的正确性. 在计算区域内,流动呈现左右对称的状态. 图 5 (a1),t = 3 s,为充型开始,随着射流区 J1 的扩大,主 涡流区 J2 对称出现在射流区两侧,此时自由表面中 间部位明显高于两侧,当流体撞击左右壁面后,两侧 流体在两侧壁的挤压作用下迅速向上移动并产生回 图 4 实验和模拟的涡流特征. (a)单浇口大平板实验涡流特征; (b)单浇口大平板数值模拟涡流特征 Fig. 4 Experimental and simulation vortex characteristics: ( a) ex鄄 perimental characteristics of single gate large plate vortex; ( b) nu鄄 merical simulation characteristics of single gate large plate vortex 流从而形成前涡流区 J3 . 由图 5(b1),当 t = 6 s 时, 将右下侧拐角处放大,可看出在两个下侧角部出现 了涡旋区 J4,其旋转方向与主涡流区相反,但强度 较弱. 随着充型的继续,如图 5(c1),当 t = 9 s 时,流 体呈现平稳向上移动的态势,自由表面中间高,两侧 低. 射流区不断增大,四个侧涡流区消失,仅存在两 个主涡流区. 图 5 实验和数值模拟的流体形态. (a)t = 3 s;(b) t = 6 s;(c) t = 9 s;(a1) t = 3 s;(b1) t = 6 s;(c1) t = 9 s Fig. 5 Fluid morphology of experiment and numerical calculation: (a) t = 3 s; (b) t = 6 s; ( c) t = 9 s; ( a1) t = 3 s; ( b1) t = 6 s; (c1) t = 9 s 2郾 2郾 3 流场区域线速度的分布 对整体流动区域进行线速度分析有利于对 卷气和夹杂成因的理解,以图 5 ( b1 ) 作为线速度 ·102·
包羽冲等：基于单相LBM模拟大平板反重力充型过程 ·103· 分析对象.图6为依据图5(b1)做出的线速度分 差持续减小，这个阶段的自由表面最不平整，流体 布云图. 运动最为活跃.当高度差下降到一定值后，在t= 图6的1处为靠近平板壁面处的流体，由于壁 3.6s时再次出现增大趋势，这是因为两侧流体的 面处的阻滞作用，流体处于静止状态，如果此处有夹 撞击壁面后的回流超过了中间的高度，这会增加 杂和气体，则不易排出：2处的射流区线速度最大， 卷气和夹杂的危险性.接着，在短暂上升之后，高 但随着高度的增加，在重力的作用下，线速度呈现减 度差再次减小.在7.2s后，自由表面高度差一直 小的趋势：3处区域处于前涡流区和主涡流区之间 处于10mm左右的动态平衡，这是因为充型过程 的交汇处，流体运动较为活跃，线速度较大：4和5 中运动的流体与两侧壁面的撞击而产生的晃动效 处分别为前涡流区和主涡流区，这两处涡流区线速 果，但是总的来说高度差趋于平稳，卷气和夹杂的 度较小且中间部位线速度几乎为零，若气体和夹杂 危险性较低. 进入该区域，则会滞留于此不易排出：6处为倾斜度 80 较大的自由表面区域，流体速度几乎与浇道口速度 70 接近，流动活跃. 60 速度ms 目50 0.22 0.20 0.18 30 0.151184 20 0.121951 0.1 10 0.08 0.0587435 4 10 0.0340444 0.02 时间/s -6 图7单浇口反重力充填自由表面高度差值随时间的变化 Fig.7 Variation of free surface altitude difference with time in sin- gle-gate counter gravity filling 2.3对大平板双浇口反重力充型过程的计算分析 1一壁面：2一射流区：3一涡流交汇区域：4一前涡流区：5一主涡流 2.3.1流体运动的特征 区：6一自由表面 图8(a)为双浇口反重力充填的大平板计算区 图6t=6s时单浇口流体线速度分布 域，大平板仍采用图3(a)的尺寸，两浇口宽度均为 Fig.6 Fluid line velocity distribution of single gate in six seconds 50mm,两浇口中心分别处于底边上的两个三等分 2.2.4自由表面高度差 点.图8(b)、(c)、(d)和(e)为大平板双浇口反重 在反重力充填金属液面的上升过程中，如果自 力充型过程的数值计算动态过程.在双浇口条件 由表面能够平稳上升，波动较小，不被后续上升的金 下，在两个浇口的两侧分别形成各自的主涡流区. 属液冲击破裂，则可以防止表面氧化皮及夹杂的卷 不同于单浇口，在两个浇道之间的区域，线速度的水 入，从而保证铸件内部金属液体的纯净度. 平分量方向相反，在相互冲击后，形成图8(b)中的 为了表征自由表面波动的剧烈程度，本文引 中间部位的两个涡流，本文称之为“内涡流区”.当 入了自由表面的“高度差”作为判断剧烈程度的特 t=6s时，从图8(c)可看出，内涡流区被拉长且分 征值，即自由表面上的最高点与最低点的高度差 离，从原来的两个涡流区演化为4个涡流区.在t= 值.通过计量模拟图像显示的网格数，来统计自由 9s时，从图8(d)中可以观察到，靠近两侧壁的主涡 表面上各点高度，然后计算出高度差，以此度量流 流区也出现类似的分离现象.这是由于流体与两侧 体在型腔内的波动起伏程度以及因晃动造成的形 壁面撞击后，上部流体与下部流体速度对称反向，形 状变化. 成了后续出现的涡流分离的现象.在t=2s时，如 如图7所示，在t=1.7s时高度差为最大值70 图8()所示，两侧壁的四个涡流区分别合并为两个 mm,此时流体由中间浇口射入型腔后，撞击到左右 主涡流区，而中心区域的四个小涡流区则完全消失， 壁面，在壁面的挤压作用下迅速上升，此时高度差最 只剩下中间的射流区 大，也是卷入气体和夹杂的可能性最大的阶段.在 2.3.2流场区域线速度的分布 其后的充填过程中，由于重力作用，上升速度持续 以图8（©）作为线速度分析对象.图9为依据 减小，使得从t=1.7s到t=3.2s这一阶段的高度 图8（©）做出的线速度分布云图.从图9中可看出
包羽冲等: 基于单相 LBM 模拟大平板反重力充型过程 分析对象. 图 6 为依据图 5( b1)做出的线速度分 布云图. 图 6 的 1 处为靠近平板壁面处的流体,由于壁 面处的阻滞作用,流体处于静止状态,如果此处有夹 杂和气体,则不易排出;2 处的射流区线速度最大, 但随着高度的增加,在重力的作用下,线速度呈现减 小的趋势;3 处区域处于前涡流区和主涡流区之间 的交汇处,流体运动较为活跃,线速度较大;4 和 5 处分别为前涡流区和主涡流区,这两处涡流区线速 度较小且中间部位线速度几乎为零,若气体和夹杂 进入该区域,则会滞留于此不易排出;6 处为倾斜度 较大的自由表面区域,流体速度几乎与浇道口速度 接近,流动活跃. 1—壁面;2—射流区;3—涡流交汇区域;4—前涡流区;5—主涡流 区;6—自由表面 图 6 t = 6 s 时单浇口流体线速度分布 Fig. 6 Fluid line velocity distribution of single gate in six seconds 2郾 2郾 4 自由表面高度差 在反重力充填金属液面的上升过程中,如果自 由表面能够平稳上升,波动较小,不被后续上升的金 属液冲击破裂,则可以防止表面氧化皮及夹杂的卷 入,从而保证铸件内部金属液体的纯净度. 为了表征自由表面波动的剧烈程度,本文引 入了自由表面的“高度差冶作为判断剧烈程度的特 征值,即自由表面上的最高点与最低点的高度差 值. 通过计量模拟图像显示的网格数,来统计自由 表面上各点高度,然后计算出高度差,以此度量流 体在型腔内的波动起伏程度以及因晃动造成的形 状变化. 如图 7 所示,在 t = 1郾 7 s 时高度差为最大值 70 mm,此时流体由中间浇口射入型腔后,撞击到左右 壁面,在壁面的挤压作用下迅速上升,此时高度差最 大,也是卷入气体和夹杂的可能性最大的阶段. 在 其后的充填过程中,由于重力作用,上升速度持续 减小,使得从 t = 1郾 7 s 到 t = 3郾 2 s 这一阶段的高度 差持续减小,这个阶段的自由表面最不平整,流体 运动最为活跃. 当高度差下降到一定值后,在 t = 3郾 6 s 时再次出现增大趋势,这是因为两侧流体的 撞击壁面后的回流超过了中间的高度,这会增加 卷气和夹杂的危险性. 接着,在短暂上升之后,高 度差再次减小. 在 7郾 2 s 后,自由表面高度差一直 处于 10 mm 左右的动态平衡,这是因为充型过程 中运动的流体与两侧壁面的撞击而产生的晃动效 果,但是总的来说高度差趋于平稳,卷气和夹杂的 危险性较低. 图 7 单浇口反重力充填自由表面高度差值随时间的变化 Fig. 7 Variation of free surface altitude difference with time in sin鄄 gle鄄gate counter gravity filling 2郾 3 对大平板双浇口反重力充型过程的计算分析 2郾 3郾 1 流体运动的特征 图 8(a)为双浇口反重力充填的大平板计算区 域,大平板仍采用图 3(a)的尺寸,两浇口宽度均为 50 mm,两浇口中心分别处于底边上的两个三等分 点. 图 8(b)、( c)、( d)和( e)为大平板双浇口反重 力充型过程的数值计算动态过程. 在双浇口条件 下,在两个浇口的两侧分别形成各自的主涡流区. 不同于单浇口,在两个浇道之间的区域,线速度的水 平分量方向相反,在相互冲击后,形成图 8( b)中的 中间部位的两个涡流,本文称之为“内涡流区冶. 当 t = 6 s 时,从图 8( c) 可看出,内涡流区被拉长且分 离,从原来的两个涡流区演化为 4 个涡流区. 在 t = 9 s 时,从图 8(d)中可以观察到,靠近两侧壁的主涡 流区也出现类似的分离现象. 这是由于流体与两侧 壁面撞击后,上部流体与下部流体速度对称反向,形 成了后续出现的涡流分离的现象. 在 t = 12 s 时,如 图 8(e)所示,两侧壁的四个涡流区分别合并为两个 主涡流区,而中心区域的四个小涡流区则完全消失, 只剩下中间的射流区. 2郾 3郾 2 流场区域线速度的分布 以图 8(c)作为线速度分析对象. 图 9 为依据 图 8(c)做出的线速度分布云图. 从图 9 中可看出, ·103·
·104 工程科学学报，第40卷，第1期 500 速度ms 0.19 (a) 0.18 0.17 0.16 0.15 0.0979416 0.08 0.05 0.03 0.01 5050L -6 1一壁面：2一主涡流区：3一自由表面：4一射流区：5一下部内涡流 区：6一上部内涡流区 图9t=6s时双浇口流体线速度分布 Fig.9 Fluid line velocity distribution of double gate in six seconds 2.4对大平板单浇口带圆柱扰流的反重力充型问 题的计算分析 2.4.1流体运动的特征 图10(a)为单浇口带圆柱扰流的反重力充型 大平板计算区域，大平板仍采用图3(a)的尺寸. 浇口宽度为100mm,处于大平板底部中间.圆柱 图8双浇口充填的数值计算结果.(a)计算区域（单位：mm): (b)t=3s(c)t=6s:(d)t=9s:(e)t=12s 障碍体穿过平板，半径为30mm,圆心距离平板底 Fig.8 Results of numerical calculation in double gate:(a)calcula- 部120mm.图10(b)、(c)和(d)为大平板带圆柱 tion area (unit:mm);(b)t=3s;(c)t=6s;(d)t=9s;(e)t= 扰流的反重力充型过程的数值计算结果.充型开 12s 始时（如图10(b)),t=3s,与单浇口充型情况类 1区域靠近平板壁面处的流体由于壁面的阻力基本 似，也会形成主涡流区、前涡流区和后涡流区；二 处于静止：2区域主涡流区和单浇口的情况一样，周 者的区别在于，射流区将会从圆柱周围绕过，之后 围线速度小，中间基本接近于零：3区域自由表面斜 在圆柱体上方有交汇趋势，向中间聚集，形成回旋 度越大，线速度越大，形状波动就越大：4区域射流 区域.在图10(c)中，当t=6s时，前涡流区已变得 区，速度最大，且随着高度的增加速度呈现减小趋 很小，圆柱体左上方和右上方的回旋区域演变成 势：5区域的下部内涡流由于处于流体内部，线速度 两个涡流区，且沿着型腔中轴线左右对称，称之为 基本接近于零，两内涡流速度相反，速度互相抵消：6 “障碍涡流”，该涡流区域随着充型过程的延续不 区域的上部内涡流在自由表面附近，是由于两侧流 断扩大.由图10(d)中可见，当t=9s时，流体中 体向中间靠拢形成，在挤压的作用下内侧速度很小， 的主涡流区明显小于无圆柱体的主涡流区，这是 越往外侧，线速度越大.可知双浇口除了壁面处和 因为圆柱体的阻碍作用使得流体运动方向发生改 主涡流区外，在两浇道之间的涡旋区域线速度也接 变，在圆柱体上方产生的障碍涡流区抑制了主涡 近于零，对比单浇口，双浇口的充填方式出现卷气和 流区的扩散. 夹杂的可能性会大大增加. 2.4.2流场区域线速度的分布 双浇口自由表面“高度差”与单浇口变化情况 以图10(c)作为线速度分析对象，图11为依据 一致，充型开始过程中波动较大，最后由于重力作用 图10(c)做出的线速度分布云图.如图11所示，1 高度差波动较小，由于流体在平板内的晃动效果，同 区域靠近平板壁面处的流体由于壁面的阻力基本处 样呈现高度差小于10mm的动态平衡.但是双浇口 于静止：2区域主涡流区线速度处于静止状态：3处 由于内部涡旋区多于单浇口的涡旋区，所以内部卷 的前涡流区中心速度接近于零，外侧处于自由表面 气夹杂的可能要大于单浇口的情况 附近，离中心区域越远速度越大：4区域的自由表面
工程科学学报,第 40 卷,第 1 期 图 8 双浇口充填的数值计算结果. ( a) 计算区域(单位:mm); (b) t = 3 s (c) t = 6 s;(d) t = 9 s;(e) t = 12 s Fig. 8 Results of numerical calculation in double gate: (a) calcula鄄 tion area (unit: mm); (b) t = 3 s; (c) t = 6 s; (d) t = 9 s; (e) t = 12 s 1 区域靠近平板壁面处的流体由于壁面的阻力基本 处于静止;2 区域主涡流区和单浇口的情况一样,周 围线速度小,中间基本接近于零;3 区域自由表面斜 度越大,线速度越大,形状波动就越大;4 区域射流 区,速度最大,且随着高度的增加速度呈现减小趋 势;5 区域的下部内涡流由于处于流体内部,线速度 基本接近于零,两内涡流速度相反,速度互相抵消;6 区域的上部内涡流在自由表面附近,是由于两侧流 体向中间靠拢形成,在挤压的作用下内侧速度很小, 越往外侧,线速度越大. 可知双浇口除了壁面处和 主涡流区外,在两浇道之间的涡旋区域线速度也接 近于零,对比单浇口,双浇口的充填方式出现卷气和 夹杂的可能性会大大增加. 双浇口自由表面“高度差冶 与单浇口变化情况 一致,充型开始过程中波动较大,最后由于重力作用 高度差波动较小,由于流体在平板内的晃动效果,同 样呈现高度差小于 10 mm 的动态平衡. 但是双浇口 由于内部涡旋区多于单浇口的涡旋区,所以内部卷 气夹杂的可能要大于单浇口的情况. 1—壁面;2—主涡流区;3—自由表面;4—射流区;5—下部内涡流 区;6—上部内涡流区 图 9 t = 6 s 时双浇口流体线速度分布 Fig. 9 Fluid line velocity distribution of double gate in six seconds 2郾 4 对大平板单浇口带圆柱扰流的反重力充型问 题的计算分析 2郾 4郾 1 流体运动的特征 图 10( a)为单浇口带圆柱扰流的反重力充型 大平板计算区域,大平板仍采用图 3 ( a) 的尺寸. 浇口宽度为 100 mm,处于大平板底部中间. 圆柱 障碍体穿过平板,半径为 30 mm,圆心距离平板底 部 120 mm. 图 10( b) 、( c)和( d) 为大平板带圆柱 扰流的反重力充型过程的数值计算结果. 充型开 始时(如图 10 ( b) ) ,t = 3 s,与单浇口充型情况类 似,也会形成主涡流区、前涡流区和后涡流区;二 者的区别在于,射流区将会从圆柱周围绕过,之后 在圆柱体上方有交汇趋势,向中间聚集,形成回旋 区域. 在图 10( c)中,当 t = 6 s 时,前涡流区已变得 很小,圆柱体左上方和右上方的回旋区域演变成 两个涡流区,且沿着型腔中轴线左右对称,称之为 “障碍涡流冶 ,该涡流区域随着充型过程的延续不 断扩大. 由图 10 ( d) 中可见,当 t = 9 s 时,流体中 的主涡流区明显小于无圆柱体的主涡流区,这是 因为圆柱体的阻碍作用使得流体运动方向发生改 变,在圆柱体上方产生的障碍涡流区抑制了主涡 流区的扩散. 2郾 4郾 2 流场区域线速度的分布 以图 10(c)作为线速度分析对象,图 11 为依据 图 10(c)做出的线速度分布云图. 如图 11 所示,1 区域靠近平板壁面处的流体由于壁面的阻力基本处 于静止;2 区域主涡流区线速度处于静止状态;3 处 的前涡流区中心速度接近于零,外侧处于自由表面 附近,离中心区域越远速度越大;4 区域的自由表面 ·104·
包羽冲等：基于单相LBM模拟大平板反重力充型过程 ·105· 500 域比主涡流区线速度大，运动更为活跃，结合图10 (©)进一步说明了新形成的涡流区抑制了主涡流区 a 的扩散 0.24 0.22 0.20 0.18 230 0.16 0.14 0.12 0.06 .100 0.0524027 0.0322975 d 0.02 0.0179708 -6 1一壁面：2一主涡流区：3一前涡流区：4一自由表面：5一射流区； 6一障碍涡流区 图11t=3s时圆柱扰流流体线速度分布 图10带圆柱扰流的充填数值计算结果.(a)计算区域（单位： Fig.11 Fluid line velocity distribution of cylindrical obstacle in three mm);(b)t=3s;(c)t=6s;(d)t=9s seconds Fig.10 Results of numerical calculation in cylindrical obstacle:(a) calculation area unit:mm);(b)t=3s;(c)t=6s;(d)t=9s 2.4.3有无圆柱扰流条件下的自由表面高度差对比 图12(a)显示圆柱扰流条件下自由表面的高度 斜度最大，线速度最大，与无圆柱体充型对比，整体 差最大值为55mm.在t=2.6s后，由于两侧流体撞 自由表面区域线速度均值明显较后者小，充型较为 击壁面后的挤压作用，流体迅速向上涨，而中间部分 平稳，说明圆柱体障碍物对充型过程自由表面的运 受到圆柱障碍物的阻流作用导致高度差的短暂增 动影响较大，有稳定流动的作用：5处射流区与无障 大.由于圆柱周围流体运动速度慢，两侧流体上涨 碍对比，越接近圆柱障碍射流区速度越小，在撞击圆 速度快，导致圆柱正上方区域形成凹陷，使得在t= 柱后射流区被一分为二，向左上方和右上方延展，在 5.2s后出现高度差的第二次短暂增大.在之后的 这两个延展方向线速度呈现减小梯度，有明显的绕 充型过程中由于晃动作用，自由表面高度差一直处 流现象，与现有圆柱绕流结论定性一致；6区域的流 于小范围波动小于6mm的动态平衡状态.从图12 体回旋向中间聚集，线速度较大，运动较为活跃，在 (b)中将无圆柱障碍物与圆柱扰流条件下自由表面 之后的运动过程中逐渐演变成障碍涡流区域，该区 高度差进行对比，可看出由于圆柱的阻流作用，自由 80 80 ·一圆柱扰流充型 b ·一无障碍充型 70 70 ·圆柱扰流充型 60 60 50 50 40 40 30 20 20 ◆ 10 10 6 0 时间s 时间s 图12自由表面的高度差.(a)圆柱扰流自由表面的高度差值随时间的变化：(b)单浇口与圆柱扰流的自由表面高度差对比 Fig.12 Altitude difference of free surface:(a)free surface variation of altitude difference with time in cylindrical obstacle;(b)comparison of free surface altitude difference between single gate and cylindrical obstacle
包羽冲等: 基于单相 LBM 模拟大平板反重力充型过程 图 10 带圆柱扰流的充填数值计算结果. ( a) 计算区域(单位: mm);(b) t = 3 s;(c) t = 6 s;(d) t = 9 s Fig. 10 Results of numerical calculation in cylindrical obstacle: (a) calculation area (unit: mm); (b) t = 3 s; (c) t = 6 s; (d) t = 9 s 斜度最大,线速度最大,与无圆柱体充型对比,整体 自由表面区域线速度均值明显较后者小,充型较为 平稳,说明圆柱体障碍物对充型过程自由表面的运 动影响较大,有稳定流动的作用;5 处射流区与无障 碍对比,越接近圆柱障碍射流区速度越小,在撞击圆 柱后射流区被一分为二,向左上方和右上方延展,在 这两个延展方向线速度呈现减小梯度,有明显的绕 图 12 自由表面的高度差. (a)圆柱扰流自由表面的高度差值随时间的变化;(b)单浇口与圆柱扰流的自由表面高度差对比 Fig. 12 Altitude difference of free surface: (a) free surface variation of altitude difference with time in cylindrical obstacle; (b) comparison of free surface altitude difference between single gate and cylindrical obstacle 流现象,与现有圆柱绕流结论定性一致;6 区域的流 体回旋向中间聚集,线速度较大,运动较为活跃,在 之后的运动过程中逐渐演变成障碍涡流区域,该区 域比主涡流区线速度大,运动更为活跃,结合图 10 (c)进一步说明了新形成的涡流区抑制了主涡流区 的扩散. 1—壁面;2—主涡流区;3—前涡流区;4—自由表面;5—射流区; 6—障碍涡流区 图 11 t = 3 s 时圆柱扰流流体线速度分布 Fig. 11 Fluid line velocity distribution of cylindrical obstacle in three seconds 2郾 4郾 3 有无圆柱扰流条件下的自由表面高度差对比 图 12(a)显示圆柱扰流条件下自由表面的高度 差最大值为 55 mm. 在 t = 2郾 6 s 后,由于两侧流体撞 击壁面后的挤压作用,流体迅速向上涨,而中间部分 受到圆柱障碍物的阻流作用导致高度差的短暂增 大. 由于圆柱周围流体运动速度慢,两侧流体上涨 速度快,导致圆柱正上方区域形成凹陷,使得在 t = 5郾 2 s 后出现高度差的第二次短暂增大. 在之后的 充型过程中由于晃动作用,自由表面高度差一直处 于小范围波动小于 6 mm 的动态平衡状态. 从图 12 (b)中将无圆柱障碍物与圆柱扰流条件下自由表面 高度差进行对比,可看出由于圆柱的阻流作用,自由 ·105·
·106· 工程科学学报，第40卷，第1期 表面的高度差最大值明显小于无障碍充型，两次短 filling and heat transfer of anti-gravity thin-walled casting.J 暂增大幅度后者也小于前者.在t=6s后的动态平 Northuestern Polytech Univ,1997,15(3):327 衡状态中，前者自由表面一直处于10mm左右的高 (吕衣礼，杨根仓，周尧和.薄壁铸件反重力铸造充型和传热 度差值，而后者该值一直处于5mm左右.结合 的数值模拟研究.西北工业大学学报.1997,15(3)：327) [6]Li M E,Wang Y X,Lii Y L,et al.Numerical simulation of 3.3.2中无障碍和3.5.2中有障碍的流体线速度对 counter gravity filling pattems for thin plate casting.Spec Cast 比分析，前者自由表面线速度比后者大，高度差值进 Nonferrous Alloy,2000,12(2):45 一步说明了当有圆柱障碍物扰流时自由面更加平 (李梅娥，王友序，吕衣礼，等.薄板件反重力充型形态的数 整，说明圆柱障碍物在充型中有使流体更加平稳的 值模拟研究.特种铸造及有色合金，2000,12(2)：45) 作用. [7]Li R.Wang Y X,Yang G C.Numerical simulation for large thin- walled plate castings in counter-gravity mould filling.Spec Cast 3结论 Nonferrous Alloy,2003(3):12 (1)本文采用格子Boltzmann单相自由表面模 (李日，王友序，杨根仓.薄壁大平板铸件的反重力充型数值 模拟.特种铸造及有色合金，2003(3)：12) 型进行流场计算.针对此方法，提出了一种权重系 [8]Janssen C,Krafczyk M.A lattice Boltzmann approach for free-sur- 数重新分配的方法来处理液相分配.并用所建立的 face-flow simulations on non-uniform block structured grids.Com- 模型模拟了大平板反重力充型过程的流体流动过 put Math Appl,2010,59(7):2215 程.数值模拟结果与水力模拟实验结论的流场特征 [9]Tang B,Li J F,Wang T S.Single-phase lattice Boltzmann model 和流体形态吻合很好，证明了该模型可用于计算反 for free-surface flow.J Tsinghua Unir Sci Tech,2008,48(11 ) 重力充型的过程 2017 (2)对单浇口、双浇口和圆柱扰流条件下大平 (汤波，李俊峰，王天舒.带自由面流体运动的单相格子Bot- 板反重力充型流动过程的模拟结果表明，双浇口浇 mam方法模拟.清华大学学报（自然科学版），2008,48 (11):2017) 注时两浇道之间的互相影响会导致流体内部运动更 [10]Zhou Y,Guo J H,Chen H X.Lattice Boltzmann simulation of 加复杂，在两浇口之间形成了两个反向对称内涡流 tow liquid droplets impacting on a liquid film.Chin J Comput 区，且该涡流区和主涡流区都会被拉长并且一分为 Phs,2010,27(1):31 二.圆柱扰流时，在圆柱体上侧会形成障碍涡流区， (周轶，郭加宏，陈红勋.格子Boltzman方法模拟双液滴撞 会抑制主涡流区的扩散.有圆柱扰流与无圆柱扰流 击液膜的流动过程.计算物理，2010,27(1)：31) 相比，前者整体线速度更加均匀且自由表面高度差 [11]Han WJ,Yan K,Chu X S,et al.Influence factors to numerical 小，说明圆柱障碍物能明显改善充型流体的稳定性， simulation of splash flow based on LBM.J Ship Mech,2015,19 降低流体与壁面冲击产生的晃动幅度，使充型更加 (4):356 平稳. (韩文骥，颜开，褚学森，等.LBM方法数值模拟飞溅流动 的影响因素研究.船舶力学，2015,19(4)：356) 参考文献 [12]Attar E,Korner C.Lattice Boltzmann model for thermal free sur- face flows with liquid-solid phase transition.Int I Heat Fluid [1]Ginzburg I,Steiner K.Lattice Boltzmann model for free-surface Floc,2011,32(1):156 flow and its application to filling process in casting.J Comput [13]Schreiber M.Neumann P,Zimmer S,et al.Free-Surface Lat- Phys,2003,185(6):61 tice-Boltzmann simulation on many-core architectures.Proc Com- [2]Peng S Y,Chen LL,Zhang M Y,et al.Numerical simulation two put Sci,2011,4:984 phase flows of casting filling process using SOLA particle level set [14]Janssen C,Krafczyk M.A lattice Boltzmann approach for free- method.Appl Math Model,2010,34(12):4106 surface-flow simulations on non-unifor block-structured grids. [3]Zhang SJ.Hou H,Yang J,et al.Numerical simulation technolo- Comput Math Appl,2010,59(7):2215 gy of the filling process.Mater Rev,2007,21(3):104 [15]Janssen C,Krafezyk M.Free surface flow simulations on GPG- (张舒娟，侯华，杨品，等.充型过程的数值模拟技术。材料 PUs using the LBM.Comput Math Appl,2011,61(12):3549 导报，2007,21(3)：104) [16]He Y L,Wang Y,Li Q.Lattice Boltzmann Method:Theory and [4]Zhou J X.Current status and development trend of casting numeri- Applications.Beijing:Science Press,2009 cal simulation technology.Foundry,2012,61(10):1105 (何雅玲，王勇，李庆.格子Boltzmann方法的理论及应用. (周建新.铸造计算机模拟仿真技术现状及发展趋势.铸造， 北京：科学出版社，2009) 2012,61(10):1105) [17]Thurey N.Physically based Animation of Free Surface Flows with [5]Lu YL,YangG C.Zhou Y H.On numerical simulation of mould the Lattice Bolizmann Method Dissertation ]Freistaat Bayem:
工程科学学报,第 40 卷,第 1 期 表面的高度差最大值明显小于无障碍充型,两次短 暂增大幅度后者也小于前者. 在 t = 6 s 后的动态平 衡状态中,前者自由表面一直处于 10 mm 左右的高 度差值,而后者该值一直处于 5 mm 左右. 结合 3郾 3郾 2 中无障碍和 3郾 5郾 2 中有障碍的流体线速度对 比分析,前者自由表面线速度比后者大,高度差值进 一步说明了当有圆柱障碍物扰流时自由面更加平 整,说明圆柱障碍物在充型中有使流体更加平稳的 作用. 3 结论 (1)本文采用格子 Boltzmann 单相自由表面模 型进行流场计算. 针对此方法,提出了一种权重系 数重新分配的方法来处理液相分配. 并用所建立的 模型模拟了大平板反重力充型过程的流体流动过 程. 数值模拟结果与水力模拟实验结论的流场特征 和流体形态吻合很好,证明了该模型可用于计算反 重力充型的过程. (2)对单浇口、双浇口和圆柱扰流条件下大平 板反重力充型流动过程的模拟结果表明,双浇口浇 注时两浇道之间的互相影响会导致流体内部运动更 加复杂,在两浇口之间形成了两个反向对称内涡流 区,且该涡流区和主涡流区都会被拉长并且一分为 二. 圆柱扰流时,在圆柱体上侧会形成障碍涡流区, 会抑制主涡流区的扩散. 有圆柱扰流与无圆柱扰流 相比,前者整体线速度更加均匀且自由表面高度差 小,说明圆柱障碍物能明显改善充型流体的稳定性, 降低流体与壁面冲击产生的晃动幅度,使充型更加 平稳. 参 考 文 献 [1] Ginzburg I, Steiner K. Lattice Boltzmann model for free鄄surface flow and its application to filling process in casting. J Comput Phys, 2003, 185(6): 61 [2] Peng S Y, Chen L L, Zhang M Y, et al. Numerical simulation two phase flows of casting filling process using SOLA particle level set method. Appl Math Model, 2010, 34(12): 4106 [3] Zhang S J, Hou H, Yang J, et al. Numerical simulation technolo鄄 gy of the filling process. Mater Rev, 2007, 21(3): 104 (张舒娟, 侯华, 杨晶, 等. 充型过程的数值模拟技术. 材料 导报, 2007, 21(3): 104) [4] Zhou J X. Current status and development trend of casting numeri鄄 cal simulation technology. Foundry, 2012, 61(10): 1105 (周建新. 铸造计算机模拟仿真技术现状及发展趋势. 铸造, 2012, 61(10): 1105) [5] L俟 Y L, Yang G C, Zhou Y H. On numerical simulation of mould filling and heat transfer of anti鄄gravity thin鄄walled casting. J Northwestern Polytech Univ, 1997, 15(3): 327 (吕衣礼, 杨根仓, 周尧和. 薄壁铸件反重力铸造充型和传热 的数值模拟研究. 西北工业大学学报, 1997, 15(3): 327) [6] Li M E, Wang Y X, L俟 Y L, et al. Numerical simulation of counter gravity filling patterns for thin plate casting. Spec Cast Nonferrous Alloy, 2000, 12(2): 45 (李梅娥, 王友序, 吕衣礼, 等. 薄板件反重力充型形态的数 值模拟研究. 特种铸造及有色合金, 2000, 12(2): 45) [7] Li R, Wang Y X, Yang G C. Numerical simulation for large thin鄄 walled plate castings in counter鄄gravity mould filling. Spec Cast Nonferrous Alloy, 2003(3): 12 (李日, 王友序, 杨根仓. 薄壁大平板铸件的反重力充型数值 模拟. 特种铸造及有色合金, 2003(3): 12) [8] Janssen C, Krafczyk M. A lattice Boltzmann approach for free鄄sur鄄 face鄄flow simulations on non鄄uniform block structured grids. Com鄄 put Math Appl, 2010, 59(7): 2215 [9] Tang B, Li J F, Wang T S. Single鄄phase lattice Boltzmann model for free鄄surface flow. J Tsinghua Univ Sci Tech, 2008, 48(11): 2017 (汤波, 李俊峰, 王天舒. 带自由面流体运动的单相格子 Bolt鄄 zmann 方法模拟. 清华大学学报( 自然科学版), 2008, 48 (11): 2017) [10] Zhou Y, Guo J H, Chen H X. Lattice Boltzmann simulation of tow liquid droplets impacting on a liquid film. Chin J Comput Phys, 2010, 27(1): 31 (周轶, 郭加宏, 陈红勋. 格子 Boltzman 方法模拟双液滴撞 击液膜的流动过程. 计算物理, 2010, 27(1): 31) [11] Han W J, Yan K, Chu X S, et al. Influence factors to numerical simulation of splash flow based on LBM. J Ship Mech, 2015, 19 (4): 356 (韩文骥, 颜开, 褚学森, 等. LBM 方法数值模拟飞溅流动 的影响因素研究. 船舶力学, 2015, 19(4): 356) [12] Attar E, K觟rner C. Lattice Boltzmann model for thermal free sur鄄 face flows with liquid鄄solid phase transition. Int J Heat Fluid Flow, 2011, 32(1): 156 [13] Schreiber M, Neumann P, Zimmer S, et al. Free鄄Surface Lat鄄 tice鄄Boltzmann simulation on many鄄core architectures. Proc Com鄄 put Sci, 2011, 4: 984 [14] Janssen C, Krafczyk M. A lattice Boltzmann approach for free鄄 surface鄄flow simulations on non鄄 uniform block鄄structured grids. Comput Math Appl, 2010, 59(7): 2215 [15] Janssen C, Krafczyk M. Free surface flow simulations on GPG鄄 PUs using the LBM. Comput Math Appl, 2011, 61(12): 3549 [16] He Y L, Wang Y, Li Q. Lattice Boltzmann Method: Theory and Applications. Beijing: Science Press, 2009 (何雅玲, 王勇, 李庆. 格子 Boltzmann 方法的理论及应用. 北京: 科学出版社, 2009) [17] Thurey N. Physically based Animation of Free Surface Flows with the Lattice Boltzmann Method [Dissertation]. Freistaat Bayern: ·106·
包羽冲等：基于单相LBM模拟大平板反重力充型过程 ·107· University of Erlangen-Nuremberg,2007 [20]Nie D M,Lin J Z.Boundary conditions in lattice Boltzmann [18]Komer C.Thies M,Hofmann T,et al.Lattice Boltzmann model method.Chin J Comput Phys,2004,21(1):21 for free surface flow for modeling foaming.J Stat Phys,2005, (聂德明，林建忠.格子Boltzmann方法中的边界条件.计算 121(1-2):179 物理.2004,21(1)：21) [19]Guo Z L,Zheng C G.Li Q,et al.Theory and Applications of [21]Li Y L,Zhang Q X,Zhou Y H.Flow morphology of thin-walled Lattice Boltzmann Method.Wuhan:Hubei Science and Technolo- castings during anti-gravity mould filling.Spec Cast Nonferrous gy Press,2002 Aloy,1994(4):1 (郭照立，郑楚光，李青，等.流体动力学的格子Boltzmann (吕衣礼，张启勋，周尧和.薄壁件反重力铸造中的充填流 方法.武汉：湖北科学技术出版社，2002) 体形态.特种铸造及有色合金，1994(4)：1)
包羽冲等: 基于单相 LBM 模拟大平板反重力充型过程 University of Erlangen鄄Nuremberg, 2007 [18] Korner C, Thies M, Hofmann T, et al. Lattice Boltzmann model for free surface flow for modeling foaming. J Stat Phys, 2005, 121(1鄄2): 179 [19] Guo Z L, Zheng C G, Li Q, et al. Theory and Applications of Lattice Boltzmann Method. Wuhan: Hubei Science and Technolo鄄 gy Press, 2002 (郭照立, 郑楚光, 李青, 等. 流体动力学的格子 Boltzmann 方法. 武汉: 湖北科学技术出版社, 2002) [20] Nie D M, Lin J Z. Boundary conditions in lattice Boltzmann method. Chin J Comput Phys, 2004, 21(1): 21 (聂德明, 林建忠. 格子 Boltzmann 方法中的边界条件. 计算 物理, 2004, 21(1): 21) [21] L俟 Y L, Zhang Q X, Zhou Y H. Flow morphology of thin鄄walled castings during anti鄄gravity mould filling. Spec Cast Nonferrous Alloy, 1994(4): 1 (吕衣礼, 张启勋, 周尧和. 薄壁件反重力铸造中的充填流 体形态. 特种铸造及有色合金, 1994(4): 1) ·107·
