CN114528781A - CFD numerical simulation method combining actual measurement wind-induced snow drift - Google Patents
CFD numerical simulation method combining actual measurement wind-induced snow drift Download PDFInfo
- Publication number
- CN114528781A CN114528781A CN202210155363.7A CN202210155363A CN114528781A CN 114528781 A CN114528781 A CN 114528781A CN 202210155363 A CN202210155363 A CN 202210155363A CN 114528781 A CN114528781 A CN 114528781A
- Authority
- CN
- China
- Prior art keywords
- snow
- phase
- wind
- model
- numerical simulation
- 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
- 238000004088 simulation Methods 0.000 title claims abstract description 70
- 238000000034 method Methods 0.000 title claims abstract description 48
- 241001609370 Puschkinia scilloides Species 0.000 title claims abstract description 32
- 238000005259 measurement Methods 0.000 title abstract description 13
- 239000002245 particle Substances 0.000 claims abstract description 41
- 238000012805 post-processing Methods 0.000 claims abstract description 12
- 230000000007 visual effect Effects 0.000 claims abstract description 10
- 238000004364 calculation method Methods 0.000 claims description 43
- 239000000463 material Substances 0.000 claims description 16
- 230000003628 erosive effect Effects 0.000 claims description 14
- 230000008021 deposition Effects 0.000 claims description 13
- 239000012530 fluid Substances 0.000 claims description 12
- 238000011161 development Methods 0.000 claims description 11
- 230000008859 change Effects 0.000 claims description 8
- 230000001052 transient effect Effects 0.000 claims description 7
- 230000001133 acceleration Effects 0.000 claims description 6
- 230000006870 function Effects 0.000 claims description 6
- 230000005514 two-phase flow Effects 0.000 claims description 6
- 238000009825 accumulation Methods 0.000 claims description 5
- 238000004422 calculation algorithm Methods 0.000 claims description 4
- 230000005484 gravity Effects 0.000 claims description 4
- 238000012821 model calculation Methods 0.000 claims description 4
- 230000021715 photosynthesis, light harvesting Effects 0.000 claims description 4
- 238000004062 sedimentation Methods 0.000 claims description 4
- 239000000725 suspension Substances 0.000 claims description 4
- 230000009191 jumping Effects 0.000 claims description 3
- 230000000903 blocking effect Effects 0.000 claims description 2
- 238000010187 selection method Methods 0.000 claims description 2
- 238000010008 shearing Methods 0.000 claims 5
- 230000001276 controlling effect Effects 0.000 claims 3
- 230000001105 regulatory effect Effects 0.000 claims 3
- 239000000243 solution Substances 0.000 claims 3
- 239000006104 solid solution Substances 0.000 claims 1
- 238000013461 design Methods 0.000 abstract description 2
- 230000004907 flux Effects 0.000 description 7
- 239000000203 mixture Substances 0.000 description 6
- 239000007787 solid Substances 0.000 description 6
- 238000002474 experimental method Methods 0.000 description 5
- 238000011160 research Methods 0.000 description 4
- 230000007704 transition Effects 0.000 description 4
- 238000010586 diagram Methods 0.000 description 3
- 230000008676 import Effects 0.000 description 3
- 230000008569 process Effects 0.000 description 3
- 238000005516 engineering process Methods 0.000 description 2
- 230000002730 additional effect Effects 0.000 description 1
- 238000004458 analytical method Methods 0.000 description 1
- 230000009286 beneficial effect Effects 0.000 description 1
- 238000004891 communication Methods 0.000 description 1
- 230000008878 coupling Effects 0.000 description 1
- 238000010168 coupling process Methods 0.000 description 1
- 238000005859 coupling reaction Methods 0.000 description 1
- 238000009792 diffusion process Methods 0.000 description 1
- 230000000694 effects Effects 0.000 description 1
- 238000004215 lattice model Methods 0.000 description 1
- 238000004519 manufacturing process Methods 0.000 description 1
- 230000000704 physical effect Effects 0.000 description 1
- 230000002265 prevention Effects 0.000 description 1
- 230000003746 surface roughness Effects 0.000 description 1
Images
Classifications
-
- 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
- G06F30/00—Computer-aided design [CAD]
- G06F30/10—Geometric CAD
- G06F30/13—Architectural design, e.g. computer-aided architectural design [CAAD] related to design of buildings, bridges, landscapes, production plants or roads
-
- 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/23—Design optimisation, verification or simulation using finite element methods [FEM] or finite difference methods [FDM]
-
- G—PHYSICS
- G06—COMPUTING; CALCULATING OR COUNTING
- G06F—ELECTRIC DIGITAL DATA PROCESSING
- G06F2111/00—Details relating to CAD techniques
- G06F2111/10—Numerical modelling
-
- 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
- Y02A—TECHNOLOGIES FOR ADAPTATION TO CLIMATE CHANGE
- Y02A90/00—Technologies having an indirect contribution to adaptation to climate change
- Y02A90/10—Information and communication technologies [ICT] supporting adaptation to climate change, e.g. for weather forecasting or climate simulation
Landscapes
- Engineering & Computer Science (AREA)
- Physics & Mathematics (AREA)
- Theoretical Computer Science (AREA)
- Geometry (AREA)
- General Physics & Mathematics (AREA)
- Computer Hardware Design (AREA)
- Evolutionary Computation (AREA)
- General Engineering & Computer Science (AREA)
- Mathematical Analysis (AREA)
- Pure & Applied Mathematics (AREA)
- Mathematical Optimization (AREA)
- Fluid Mechanics (AREA)
- Computing Systems (AREA)
- Mathematical Physics (AREA)
- Algebra (AREA)
- Architecture (AREA)
- Civil Engineering (AREA)
- Structural Engineering (AREA)
- Computational Mathematics (AREA)
- Management, Administration, Business Operations System, And Electronic Commerce (AREA)
- Aerodynamic Tests, Hydrodynamic Tests, Wind Tunnels, And Water Tanks (AREA)
Abstract
Description
技术领域technical field
本发明属于土木工程领域,涉及风致雪漂移的CFD数值模拟方法。The invention belongs to the field of civil engineering and relates to a CFD numerical simulation method of wind-induced snow drift.
背景技术Background technique
近些年来,在寒冷地区暴雪给人们的出行,通讯及生产生活带来了极大的不便,同时也极大地威胁了房屋以及公共建筑的安全,其中风致雪漂移容易导致建筑局部位置的荷载过大,而这种建筑上的不均匀雪荷载往往会对结构产生额外的作用,使得结构更容易失稳。因而对风雪的研究以及预防变得越来越重要。In recent years, blizzards in cold areas have brought great inconvenience to people's travel, communication, production and life, and also greatly threatened the safety of houses and public buildings. Among them, wind-induced snow drifts can easily lead to excessive loads at local locations of buildings. However, uneven snow loads on such buildings tend to have additional effects on the structure, making it more susceptible to instability. Therefore, the research and prevention of wind and snow has become more and more important.
在风雪灾害的早期研究中,囿于计算机的条件,国内外的学者们以实测以及实验研究为主,随着社会的不断进步和科学技术的迅猛发展,以计算流体力学为基础的数值模拟研究逐渐被大众所采纳。目前的风致雪漂移数值模拟研究采用单向耦合的假设,参数设置采用经验值,单向耦合即为假设不考虑雪相对空气相的影响,将雪相视为被动运输的固体介质,而雪相的沉积与侵蚀在实际中会与对积雪边界条件产生影响,从而进一步影响流场分布,而该方法本质上为稳定的风场下对雪相的扩散进行计算,这与实际情况大相径庭,降低了模拟的准确性。In the early research on wind and snow disasters, limited by the computer conditions, scholars at home and abroad focused on actual measurement and experimental research. With the continuous progress of society and the rapid development of science and technology, numerical simulation based on computational fluid dynamics Research is gradually being adopted by the general public. The current numerical simulation study of wind-induced snow drift adopts the assumption of one-way coupling, and the parameter setting adopts the empirical value. In practice, the deposition and erosion will affect the snow boundary conditions, thereby further affecting the flow field distribution, and this method is essentially the calculation of the snow phase diffusion in a stable wind field, which is very different from the actual situation and reduces the the accuracy of the simulation.
发明内容SUMMARY OF THE INVENTION
本发明是要解决现有的数值模拟方法不考虑风雪之间的相对运动、准确性差的技术问题,而提供一种结合实测的风致雪漂移的CFD数值模拟方法。该方法考虑风雪相对运动,从风雪两相流的建模,到参数设置,模拟计算等一系列过程进行了优化,得到风致雪漂移的CFD模拟数值。The present invention aims to solve the technical problems that the relative motion between wind and snow is not considered and the accuracy is poor in the existing numerical simulation method, and provides a CFD numerical simulation method combined with the measured wind-induced snow drift. The method considers the relative motion of wind and snow, and optimizes a series of processes from the modeling of wind and snow two-phase flow, to parameter setting, and simulation calculation, and obtains the CFD simulation value of wind-induced snow drift.
本发明的结合实测的风致雪漂移的CFD数值模拟方法,按以下步骤进行:The CFD numerical simulation method of the present invention combined with the measured wind-induced snow drift is carried out according to the following steps:
步骤一:收集数值模拟需要的基本信息:收集具体研究的建筑物的几何特征、该建筑所处的地理位置、地貌参数、参考高度H0处平均风速V0、风剖面规范、湍动能规范、湍动强度规范、积分湍流尺度规范、湍能耗散率规范;同时测得建筑物所在地积雪的雪颗粒密度、雪颗粒直径ds、休止角、风速、风向及建筑物的初始积雪深度,计算阈值摩擦速度、跃移雪粒速度、地面粗糙度、建筑物墙面粗糙度、雪面粗糙度和雪粒沉降速度wf;Step 1: Collect the basic information required for numerical simulation: collect the geometric characteristics of the building under study, the geographic location of the building, the topographic parameters, the average wind speed V 0 at the reference height H 0 , the wind profile specification, the turbulent kinetic energy specification, Turbulence intensity specification, integrated turbulence scale specification, turbulent energy dissipation rate specification; meanwhile, the snow particle density, snow particle diameter d s , angle of repose, wind speed, wind direction and the initial snow depth of the building are measured. , calculate the threshold friction velocity, the snow particle velocity of jumping, the ground roughness, the building wall roughness, the snow surface roughness and the snow particle settling velocity w f ;
步骤二:建立数值模拟需要的网格模型:利用步骤一所收集的建筑物几何信息和初始积雪分布形态进行模型的建立以及网格的划分,得到网格文件;Step 2: Establish a grid model required for numerical simulation: use the building geometric information and the initial snow distribution pattern collected in step 1 to establish a model and divide the grid to obtain a grid file;
步骤三:编写数值模拟用的空气项边界条件:利用步骤一所收集的建筑物所在地的地理位置信息、参考高度平均风速、风剖面、湍动能、湍动强度、积分湍流尺度、湍流耗散率规范信息编写空气相的边界条件,再在Visual Studio中以C语言、Fluent二次开发接口为基础,编写指数型风剖面、湍动能、湍流耗散率程序,保存C语言文件,并将其与步骤二网格模型文件放在同一文件夹中;Step 3: Write the air term boundary conditions for numerical simulation: use the geographic location information of the building location collected in Step 1, the average wind speed at the reference height, the wind profile, the turbulent kinetic energy, the turbulent intensity, the integrated turbulent scale, and the turbulent dissipation rate. The boundary conditions of the air phase are written in the normative information, and then based on the C language and the Fluent secondary development interface in Visual Studio, the exponential wind profile, turbulent kinetic energy, and turbulent dissipation rate programs are written, and the C language file is saved. Step 2: Grid model files are placed in the same folder;
步骤四:编写数值模拟用的雪相边界条件:利用步骤一所收集的雪颗粒密度、建筑物所在地的重力加速度、阈值摩擦速度,采用雪浓度经验公式及其公式中所采用的雪相运动速度、雪颗粒摩擦速度、跃移和悬移层临界高度信息,在Visual Studio中以C语言、Fluent二次开发接口为基础编写雪相体积分数公式,作为雪项的入口边界条件,保存C语言文件,并将其与步骤二模型网格文件、步骤三的C语言文件放在同一文件夹中;Step 4: Write the snow boundary conditions for numerical simulation: use the snow particle density collected in step 1, the gravitational acceleration at the location of the building, and the threshold friction speed, and use the empirical formula of snow concentration and the snow movement speed used in the formula. , Snow particle friction velocity, jump and critical height information of the suspension layer, write the snow phase volume fraction formula based on C language and Fluent secondary development interface in Visual Studio, as the inlet boundary condition of the snow term, save the C language file , and put it in the same folder as the model grid file in
步骤五:数值模拟需要的程序接入,具体操作如下:打开Ansys-Fluent界面选择2D/3D、双精度、并行,在General Options-Working Directory下选择包含步骤二的网格文件、步骤三的C语言文件、步骤四的C语言文件的文件夹目录,选择环境下Set upCompilation Environment for UDF接口打开软件;选择Ansys-Fluent主菜单下的UserDefined–Functions–Compiled,选择步骤三、四的C语言文件进行加载编译,得到空气相及雪相的入口速度及体积分数;Step 5: Program access required for numerical simulation. The specific operations are as follows: Open the Ansys-Fluent interface and select 2D/3D, double precision, and parallelism. Under General Options-Working Directory, select the grid
步骤六:采用Ansys-Fluent软件,先设置模型、材料、求解器参数,再导入步骤二所得的网格模型、步骤五所得的空气相及雪相的入口速度及体积分数,进行风雪两相流模拟,得到壁面剪切力、空气相及雪相的速度及体积分数的收敛文件;Step 6: Using Ansys-Fluent software, first set the model, material, and solver parameters, and then import the grid model obtained in
步骤七:利用步骤六计算出的收敛文件,进行壁面剪切力和速度的后处理:将Ansys-Fluent计算出的收敛文件以Tecplot支持格式保存,用Tecplot进行后处理计算,通过提取计算的计算域三维空间中x、y、z方向的速度、压力、壁面剪切力的分量,用Tecplot自定义函数Specify Equations,做出计算域范围内的压力、速度、壁面剪切力的云图、等值线图;Step 7: Use the convergence file calculated in Step 6 to perform post-processing of the wall shear force and velocity: Save the convergence file calculated by Ansys-Fluent in the format supported by Tecplot, use Tecplot for post-processing calculation, and calculate by extracting the calculation The components of the velocity, pressure, and wall shear force in the x, y, and z directions in the domain three-dimensional space. Use the Tecplot custom function Specify Equations to make the cloud map and equivalent value of the pressure, velocity, and wall shear force within the calculation domain. line graph;
步骤八:数值模拟雪沉积侵蚀的后处理计算:基于步骤七计算获得的壁面剪切力,通过Tecplot提取模型底面或底面横、纵方向中线上的壁面剪切力,通过公式把壁面剪切力转换为壁面的摩擦速度,其中τ为壁面剪切力,ρ为空气密度;对壁面摩擦速度与阈值磨擦速度进行比较,当壁面摩擦速度大于阈值磨擦速度时,积雪发生侵蚀;当壁面摩擦速度小于阈值磨擦速度时,积雪发生沉积;分别计算侵蚀量、沉积量,并且计算积雪改变厚度,最后将积雪改变厚度与建筑物的初始积雪深度进行叠加,得到积雪的最终分布形态,完成结合实测的风致雪漂移的CFD数值模拟。Step 8: Post-processing calculation of numerical simulation of snow deposition erosion: Based on the wall shear force calculated in Step 7, extract the wall shear force on the bottom surface of the model or the horizontal and vertical midline of the bottom surface through Tecplot, and use the formula Convert the wall shear force to the wall friction velocity, where τ is the wall shear force, ρ is the air density; compare the wall friction velocity with the threshold friction velocity, when the wall friction velocity is greater than the threshold friction velocity, the snow will erode ; When the friction velocity of the wall is less than the threshold friction velocity, the snow is deposited; the erosion amount and deposition amount are calculated respectively, and the snow change thickness is calculated. Finally, the snow change thickness and the building's initial snow depth are superimposed to obtain the The final distribution of snow, complete the CFD numerical simulation combined with the measured wind-induced snow drift.
更进一步地,步骤二中所述的网格文件的求解过程是:在Ansys-Icem中按照建筑物几何特征确定模型计算域且计算域设置满足阻塞率小于3%,接着在Icem中依据模型建立网格模型,为保证计算顺利收敛及模拟的准确性,依据模型大小设置合理的网格尺寸,调控网格偏斜率(Skewness)三角形与四面体不大于0.95,调控网格纵横比(Aspect Ratio),调控网格压扁程度(Squish)在0-1之间,最后以Ansys-Fluent格式输出求解的网格文件。Further, the solution process of the grid file described in the
更进一步地,其中合理的网格尺寸是指:网格模型控制结构化网格质量在0.6~1,非结构化网格质量在0.3~1;网格纵横比为:流动核心区内网格纵横比小于5:1,边界层网格纵横比小于10:1。Further, the reasonable grid size refers to: the grid model controls the structured grid quality to be 0.6-1, and the unstructured grid quality is 0.3-1; the grid aspect ratio is: the grid in the flow core area. The aspect ratio is less than 5:1, and the aspect ratio of the boundary layer mesh is less than 10:1.
更进一步地,步骤四中所述的雪浓度经验公式为1990版的Pomeroy&Gray、1992版的Pomeroy&Male所提出的雪浓度经验公式。Further, the snow concentration empirical formula described in
更进一步地,步骤六中所述的风雪两相流模拟,具体的过程如下:Further, the simulation of the two-phase flow of wind and snow described in step 6, the specific process is as follows:
1)采用Ansys-Fluent软件,在General下Time一栏选择瞬态Transient;1) Using Ansys-Fluent software, select Transient Transient in the Time column under General;
2)在模型Models下多项流Multiphase一栏选择Mixture模型,在模型Models下Viscous一栏选择k-epsilon两方程模型下的Realizable k-epsilon模型,在模型Models下多项流Multiphase一栏选择Phase Interaction-slip,设置滑动速度代数方程;2) Select the Mixture model in the Multiphase column under Models, select the Realizable k-epsilon model under the k-epsilon two-equation model in the Viscous column under Models, and select Phase in the Multiphase column under Models Interaction-slip, set the sliding speed algebraic equation;
3)在材料Materials下选择流体Fluid,添加空气相,设置Viscosity粘度、Density密度并且新建材料雪相,设置雪相Viscosity粘度与空气相一致并依据步骤一设置Density密度、粒径;3) Select Fluid under Material Materials, add air phase, set Viscosity viscosity, Density density and create a new material snow phase, set snow phase Viscosity viscosity to be consistent with air and set Density density and particle size according to step 1;
4)在体网格Cell Zone Conditions下,把网格计算域内设置为所需要的流体fluid、固体Solid,在Boundary Conditions设置边界条件:入口边界为速度入口、计算域两侧及顶部为对称边界、地面及模型表面为无滑移壁面、出口边界为自由出口;4) Under the volume grid Cell Zone Conditions, set the required fluid fluid and solid solid in the grid computational domain, and set the boundary conditions in the Boundary Conditions: the inlet boundary is the velocity inlet, the two sides and top of the computational domain are symmetrical boundaries, The ground and model surfaces are non-slip walls, and the exit boundary is a free exit;
5)在Boundary Conditions-inlet(入口)接入步骤三的空气相边界条件:设置风剖面、湍动能及耗散率,在Boundary Conditions-inlet(入口)-snow雪相-Multiphase下Volume Fraction接入步骤四的雪相体积分数;Fluent-Mixture模型采用有限体积法计算,由于风雪流只有两相,且两相的积分分数和为1,因此只需给出雪相的体积分数,空气相由软件计算无需设置;5) In the Boundary Conditions-inlet (inlet), access the air phase boundary conditions of step 3: set the wind profile, turbulent kinetic energy and dissipation rate, and access the Volume Fraction in the Boundary Conditions-inlet (inlet)-snow snow phase-Multiphase The volume fraction of the snow phase in
6)在求解器设置Solution下的选择方法Methods中,选用SIMPLE算法;6) In the selection method Methods under the solver setting Solution, select the SIMPLE algorithm;
7)在求解器设置Solution下的选择Monitors-Residual中,设置连续性方程、动量方程以及k-epsilon的收敛精度,其中k-epsilon的收敛精度设置为1e-6精度;7) In the Select Monitors-Residual under the Solver setting Solution, set the Continuity Equation, Momentum Equation and the Convergence Accuracy of k-epsilon, where the Convergence Accuracy of k-epsilon is set to 1e- 6 Accuracy;
8)在求解器设置Solution下的初始化Initialization中,选用标准初始化Standard Initialization,在入口初始化;8) In the initialization Initialization under the solver setting Solution, select the standard initialization Standard Initialization, and initialize at the entrance;
9)最后在Run Calculation设置计算的时间步长、计算步数,开始计算,得到了壁面剪切力、空气相及雪相的速度及体积分数的收敛文件。9) Finally, set the calculation time step and the number of calculation steps in Run Calculation, start the calculation, and obtain the convergence files of the wall shear force, the velocity and volume fraction of the air phase and the snow phase.
更进一步地,步骤八中,侵蚀量的计算公式为:其中Aero为常系数,u*为入口摩擦速度,u*t为阈值速度。Further, in
更进一步地,步骤八中,沉积量的计算公式为:其中wf中为雪的沉降速度,φs为雪的质量浓度,φs=ρsf,f为雪的体积分数,ρs为雪颗粒密度。Further, in
更进一步地,步骤八中,积雪改变厚度的计算公式为:其中qs为qero或qdep,γ为积雪的最大体积分数,ρs为雪颗粒密度。Furthermore, in
本发明的方法运用Icem CFD软件进行模型的建立,将实测所测得的立方体及周边积雪分布运用于模型的建立中,考虑风致雪漂移的初始积雪分布,并采用Ansys Fluent进行数值模拟分析,本方法基于Euler-Euler模型的两相流理论,结合运用Visual Studio改进的湍流模型的数值计算开发分析程序,将实测所测得的风速、风向、雪颗粒的物理特性(形状、密度、尺寸)、风雪漂移的阈值速度、雪颗粒的降落速度、风雪漂移的运动模式、休止角等数据运用于模拟中,对建筑的积雪分布进行模拟,完成结合实测的风致雪漂移的CFD数值模拟。将模拟结果与风致雪漂移经典实验进行对照,并同未与实测相结合的方法进行对照,验证了该模拟方法的准确性。本发明考虑了风雪之间的相对运动并与实测进行的结合,提高了风致雪漂移数值模拟的准确性。The method of the invention uses Icem CFD software to establish a model, applies the measured cube and surrounding snow distribution in the establishment of the model, considers the initial snow distribution of wind-induced snow drift, and uses Ansys Fluent to carry out numerical simulation analysis , This method is based on the two-phase flow theory of the Euler-Euler model, combined with the numerical calculation and analysis program of the improved turbulence model using Visual Studio, and the measured wind speed, wind direction, and physical properties of snow particles (shape, density, size, etc. ), the threshold speed of wind and snow drift, the falling speed of snow particles, the movement mode of wind and snow drift, the angle of repose and other data are used in the simulation to simulate the snow distribution of the building, and complete the CFD value combined with the measured wind-induced snow drift. simulation. The simulation results were compared with the classical experiment of wind-induced snow drift, and with the method that was not combined with the actual measurement, to verify the accuracy of the simulation method. The invention considers the relative motion between wind and snow and combines it with actual measurement, thereby improving the accuracy of numerical simulation of wind-induced snow drift.
本发明的结合实测的风致雪漂移的CFD数值模拟方法的整体流程示意如图1所示,本发明的结合实测的风致雪漂移的数值模拟方法具有下列优点:The overall flow diagram of the CFD numerical simulation method combined with the measured wind-induced snow drift of the present invention is shown in Figure 1. The numerical simulation method of the present invention combined with the measured wind-induced snow drift has the following advantages:
1、考虑了风相与雪相之间的相对滑移速度,通过求解两相之间的相对速度来考虑两相之间的作用。1. The relative slip velocity between the wind phase and the snow phase is considered, and the effect between the two phases is considered by solving the relative velocity between the two phases.
2、将实测获得的风吹雪的积雪分布形态运用在模型的建立中,考虑了雪相的沉积与侵蚀对计算域内流场分布的影响。2. The snow distribution pattern of wind-blown snow obtained from the actual measurement is used in the establishment of the model, and the influence of snow deposition and erosion on the flow field distribution in the calculation domain is considered.
3、运用实验工具测得了各项基于当地条件的数据并将其运用在数值模拟中,使得模拟结果与实际更为贴合,更为准确。3. Use experimental tools to measure various data based on local conditions and apply them in numerical simulation, so that the simulation results are more in line with reality and more accurate.
本发明的方法有助于认识风雪运动的规律,有助于预防建筑倒塌、减少人员伤亡及经济损失。该方法可用于寒冷地区的建筑结构设计领域。The method of the present invention helps to recognize the law of wind and snow movement, helps prevent the collapse of buildings, and reduces casualties and economic losses. This method can be used in the field of building structure design in cold regions.
附图说明Description of drawings
图1是本发明的结合实测的风致雪漂移的CFD数值模拟方法的整体流程示意图;1 is a schematic diagram of the overall flow of the CFD numerical simulation method of the present invention in conjunction with the measured wind-induced snow drift;
图2是实施例1的网格划分以及计算域的图;Fig. 2 is the grid division of embodiment 1 and the figure of computation domain;
图3是实施例1的经步骤七得到的积雪分布形式图;Fig. 3 is the snow distribution form figure that obtains through step 7 of embodiment 1;
图4是实施例1得到的单位时间内的积雪通量、经典实验得到的单位时间内的积雪通量及实施例2即未与实测相结合的方法得到的单位时间内的积雪通量的对照图。Fig. 4 shows the snow flux per unit time obtained in Example 1, the snow flux per unit time obtained by the classical experiment, and the snow flux per unit time obtained by the method not combined with actual measurement in Example 2 Quantitative comparison chart.
具体实施方式Detailed ways
用下面的实施例验证本发明的有益效果:Verify the beneficial effects of the present invention with the following examples:
实施例1:本实施例为结合实测的风致雪漂移的数值模拟方法,按以下步骤进行:Embodiment 1: This embodiment is a numerical simulation method combined with the measured wind-induced snow drift, and is carried out according to the following steps:
步骤一:收集数值模拟需要的基本信息:收集具体研究的建筑物的几何特征:边长为1m的立方体模型,该建筑所处的地理位置地貌参数为0.16、参考高度10m处平均风速5m/s,风剖面采用公式z0与v0分别表示参考点的高度(1m)及对应的参考平均风速(5m/s),z表示任一高度;α表示地貌参数,对于B类地貌,α取值为0.16;湍动能采用公式湍能耗散率采用公式其中,Cμ为经验常数,取0.09;I和l分别为湍流强度和湍流长度尺度,其中:Step 1: Collect the basic information required for numerical simulation: Collect the geometric characteristics of the building to be studied: a cube model with a side length of 1m, the geographic location of the building, and the geomorphological parameters of 0.16. The average wind speed at a reference height of 10m is 5m/s , the wind profile adopts the formula z 0 and v 0 represent the height of the reference point (1m) and the corresponding reference average wind speed (5m/s), respectively, z represents any height; α represents the landform parameter, for B-type landform, α is 0.16; turbulent kinetic energy using the formula Turbulent energy dissipation rate using the formula Among them, C μ is an empirical constant, which is taken as 0.09; I and l are the turbulent intensity and turbulent length scales, respectively, where:
湍动能强度采用公式积分湍流尺度采用公式其中,zG、zb、α'、均为地貌参数,在日本AIJ规范中,B类地貌对应Ⅱ类平地,以上三个地貌参数取值分别为5m、350m和0.15。The turbulent kinetic energy intensity adopts the formula The integral turbulence scale adopts the formula Among them, z G , z b , α', are all geomorphological parameters. In the Japanese AIJ code, the B-type landform corresponds to the II-level flat land, and the values of the above three geomorphic parameters are 5m, 350m and 0.15 respectively.
通过采用分样筛、锥形漏斗、Snow Fork以及HOBO小型气象站仪器测得当地的雪颗粒直径0.00015m、休止角45°、雪密度145kg/m3、计算阈值摩擦速度为0.2m/s、跃移雪粒速度0.2m/s、地面粗糙度系数为0.25、雪粒沉降速度0.2m/s。The local snow particle diameter is 0.00015m, the angle of repose is 45°, the snow density is 145kg/m 3 , the calculated threshold friction speed is 0.2m/s, The jumping snow particle velocity is 0.2m/s, the ground roughness coefficient is 0.25, and the snow particle settling velocity is 0.2m/s.
步骤二:建立数值模拟需要的网格模型:在Ansys-Icem CFD中按照建筑物在步骤一所收集的建筑物几何信息和初始积雪分布形态进行模型的建立以及网格的划分,确定模型计算域,长、宽、高分别为20m、10m、8m,且计算域设置阻塞率为1.25%,接着在Icem CFD中依据模型建立网格模型,依据模型大小网格数为360万。网格模型控制结构化网格质量在0.8,非结构化网格质量在0.8,最后以Ansys-Fluent格式输出求解的网格文件。网格划分以及计算域的图如图2所示;图2中1为对称边界,2为自由发展边界,3为速度入口,4为无滑移壁面。Step 2: Establish the grid model required for numerical simulation: In Ansys-Icem CFD, the model is established and the grid is divided according to the building geometric information and the initial snow distribution shape collected in step 1, and the model calculation is determined. The length, width and height are 20m, 10m and 8m respectively, and the blockage rate of the computational domain is set to 1.25%. Then, a grid model is established based on the model in Icem CFD, and the number of grids is 3.6 million according to the size of the model. The mesh model controls the structured mesh quality at 0.8, the unstructured mesh quality at 0.8, and finally outputs the solved mesh file in Ansys-Fluent format. The graph of mesh division and computational domain is shown in Fig. 2; in Fig. 2, 1 is the symmetric boundary, 2 is the free development boundary, 3 is the velocity inlet, and 4 is the non-slip wall.
步骤三:编写数值模拟用的空气项边界条件:利用步骤一所收集的建筑物所在地理位置信息、参考高度平均风速、风剖面、湍动能、湍动强度、积分湍流尺度、湍流耗散率规范信息编写空气相的边界条件,在Visual Studio中以C语言、Fluent二次开发接口为基础,编写指数型风剖面、湍动能、湍流耗散率,保存C语言文件,并将其与步骤二网格模型文件放在同一文件夹中;Step 3: Write the boundary conditions of the air term for numerical simulation: use the geographic location information of the building collected in Step 1, the average wind speed at the reference height, the wind profile, the turbulent kinetic energy, the turbulent intensity, the integral turbulence scale, and the turbulent dissipation rate specification Write the boundary conditions of the air phase based on the C language and the Fluent secondary development interface in Visual Studio, write the exponential wind profile, turbulent kinetic energy, and turbulent dissipation rate, save the C language file, and compare it with the second step. Lattice model files are placed in the same folder;
步骤四:编写数值模拟用的雪相边界条件:利用步骤一所收集的雪颗粒密度145kg/m3、当地重力加速度9.81m/s2、阈值磨擦速度0.2m/s,采用Pomeroy&Gray(1990)、Pomeroy&Male(1992)所提出的雪浓度经验公式,即公式其中ρs表示雪颗粒密度,g为重力加速度,取9.81m/s2,up为跃移层雪粒平均速度,up=2.8u*t,u*为入口摩擦速度;根据Beyers等实测拟合得到的摩擦速度与10m高度处风速u10的关系确定,雪相运动速度取5m/s、雪颗粒摩擦速度取0.15m/s、u*t为阈值速度;取0.15m/s;hs为跃移层和悬移层的临界高度,雪相运动速度5m/s、雪颗粒摩擦速度0.15m/s、在Visual Studio中以C语言、Fluent二次开发接口为基础编写雪相体积分数公式作为雪项的入口边界条件,保存C语言文件,并将其与步骤二模型网格文件、步骤三的C语言文件放在同一文件夹中;Step 4: Write the snow boundary conditions for numerical simulation: using the density of snow particles collected in step 1 of 145kg/m 3 , the local acceleration of gravity 9.81m/s 2 , and the threshold friction velocity of 0.2m/s, Pomeroy & Gray (1990), The empirical formula for snow concentration proposed by Pomeroy & Male (1992) is the formula where ρ s is the density of snow particles, g is the acceleration of gravity, which is 9.81m/s 2 , u p is the average speed of snow particles in the transition layer, u p = 2.8u *t , u * is the inlet friction speed; according to the actual measurement of Beyers et al. The relationship between the friction speed obtained by fitting and the wind speed u 10 at a height of 10m is determined, The speed of snow phase movement is 5m/s, the friction speed of snow particles is 0.15m/s, u *t is the threshold speed; 0.15m/s is taken; h s is the critical height of the transition layer and the suspension layer, Snow phase motion speed is 5m/s, snow particle friction speed is 0.15m/s, and the snow phase volume fraction formula is written in Visual Studio based on C language and Fluent secondary development interface as the inlet boundary condition of the snow term, and the C language file is saved , and put it in the same folder as the model grid file in
步骤五:数值模拟需要的程序接入,具体操作如下:打开Ansys-Fluent界面选择3D、双精度、并行。在General Options-Working Directory下选择包含步骤二的网格文件、步骤三、四的自定义程序的文件夹目录,选择环境下Set up Compilation Environmentfor UDF接口打开软件。选择Ansys-Fluent主菜单下User-Defined-Functions-Compiled下选择步骤三、四的C语言文件进行加载编译,得到空气相及雪相的入口速度及体积分数;Step 5: Program access required for numerical simulation. The specific operations are as follows: Open the Ansys-Fluent interface and select 3D, double precision, and parallel. Under General Options-Working Directory, select the folder directory containing the grid file in
步骤六:采用Ansys-Fluent软件,先设置模型、材料、求解器参数,再导入步骤二所得的网格模型、步骤五所得的空气相及雪相的入口速度及体积分数,进行风雪两相流模拟,得到壁面剪切力、空气相及雪相的速度及体积分数的收敛文件,具体如下:Step 6: Using Ansys-Fluent software, first set the model, material, and solver parameters, and then import the grid model obtained in
1)在General下Time一栏选择瞬态Transient;1) Select transient Transient in the Time column under General;
2)在模型Models下多项流Multiphase一栏选择Mixture模型、在模型Models下Viscous一栏选择k-epsilon两方程模型下的可实现Realizable k-epsilon模型,在模型Models下多项流Multiphase一栏选择Phase Interaction-slip,设置相对滑动速度0.25m/s;2) Select the Mixture model in the Multiphase column under Models, select the Realizable k-epsilon model under the k-epsilon two-equation model in the Viscous column under Models, and select the Multiphase multiphase model under Models. Select Phase Interaction-slip and set the relative sliding speed to 0.25m/s;
3)在材料Materials下选择流体Fluid,添加空气相,设置Viscosity粘度、Density密度,并且新建材料雪相,设置雪相Viscosity粘度与空气一致,依据步骤一设置Density密度145kg/m3(由Snow Fork雪特性分析仪测得)、粒径0.00015m(由分样筛测得);3) Under the material Materials, select Fluid, add the air phase, set the Viscosity viscosity and Density density, and create a new material snow phase, set the snow phase Viscosity viscosity to be the same as air, and set the Density density to 145kg/m 3 according to step 1 (by Snow Fork Snow characteristic analyzer), particle size 0.00015m (measured by sample sieve);
4)在体网格Cell Zone Conditions下的网格计算域内设置需要的流体fluid、固体Solid,在Boundary Conditions设置边界条件:入口边界为速度入口、计算域两侧及顶部为对称边界、地面及模型表面为无滑移壁面、出口边界为自由出口;4) Set the required fluid fluid and solid solid in the grid computing domain under the volume grid Cell Zone Conditions, and set the boundary conditions in the Boundary Conditions: the inlet boundary is the velocity inlet, the sides and top of the computational domain are symmetrical boundaries, the ground and the model The surface is a non-slip wall, and the outlet boundary is a free outlet;
5)在Boundary Conditions-inlet(入口)接入步骤三的空气相边界条件:设置风剖面、湍动能及耗散率,在Boundary Conditions-inlet(入口)-snow雪相-Multiphase下Volume Fraction接入步骤四的雪相体积分数;Fluent-Mixture模型采用有限体积法计算,由于风雪流只有两相,且两相的积分分数和为1,因此只需给出雪相的体积分数,空气相由软件计算无需设置。5) In the Boundary Conditions-inlet (inlet), access the air phase boundary conditions of step 3: set the wind profile, turbulent kinetic energy and dissipation rate, and access the Volume Fraction in the Boundary Conditions-inlet (inlet)-snow snow phase-Multiphase The volume fraction of the snow phase in
6)在求解器设置Solution下选择方法Methods中选用SIMPLE算法。6) Select the SIMPLE algorithm in Methods under the solver setting Solution.
7)在求解器设置Solution下选择Monitors-Residual中设置连续性方程、动量方程以及k-epsilon的收敛精度,本实施例设置收敛精度为1e-6精度。7) Under the solver setting Solution, select Monitors-Residual to set the convergence precision of the continuity equation, the momentum equation and the k-epsilon. In this embodiment, the convergence precision is set to 1e- 6 precision.
8)在求解器设置Solution下初始化Initialization中选用标准初始化StandardInitialization在入口初始化。8) Select the standard initialization StandardInitialization in the initialization initialization under the solver setting Solution to initialize at the entrance.
9)最后计算在Run Calculation设置计算的步数为1000步,精度为10-6,开始计算,得到壁面剪切力、空气相及雪相的速度及体积分数的收敛文件。9) Final calculation Set the number of calculation steps in Run Calculation to 1000 steps and the accuracy to 10 -6 , start the calculation, and obtain the convergence files of wall shear force, air phase and snow phase velocity and volume fraction.
步骤七:利用步骤六计算出的收敛文件,进行壁面剪切力和速度的后处理:将Ansys-Fluent计算出的收敛文件以Tecplot支持格式保存,用Tecplot进行后处理计算,通过提取计算的计算域三维空间中x、y、z方向的速度、压力、壁面剪切力的分量,用Tecplot自定义函数Specify Equations,做出计算域范围内的压力、速度、壁面剪切力的云图、等值线图;其中积雪分布形式图如图3所示。Step 7: Use the convergence file calculated in Step 6 to perform post-processing of the wall shear force and velocity: Save the convergence file calculated by Ansys-Fluent in the format supported by Tecplot, use Tecplot for post-processing calculation, and calculate by extracting the calculation The components of the velocity, pressure, and wall shear force in the x, y, and z directions in the domain three-dimensional space. Use the Tecplot custom function Specify Equations to make the cloud map and equivalent value of the pressure, velocity, and wall shear force within the calculation domain. Line diagram; the distribution of snow cover is shown in Figure 3.
步骤八:数值模拟雪沉积侵蚀的后处理计算:基于步骤七计算获得的壁面剪切力:基于步骤七计算获得的壁面剪切力,通过Tecplot提取模型底面的壁面剪切力,通过自定义函数把壁面剪切应力转换为壁面的摩擦速度,将壁面摩擦速度与阈值磨擦速度进行比较,当壁面摩擦速度大于阈值磨擦速度时为侵蚀,采用公式计算侵蚀量,其中Aero为常系数,取7.0×10-4;当壁面摩擦速度小于阈值磨擦速度时沉积,采用公式计算沉积量,其中wf为雪的沉降速度,φs为雪的质量浓度,φs=ρsf,f为雪的体积分数,ρs为雪颗粒密度;再采用公式计算积雪改变厚度,其中qs为qero或qdep,γ为积雪的最大体积分数,取0.62;最后将积雪改变厚度与步骤一的建筑物的初始积雪深度进行叠加,得到积雪的最终分布形态,完成结合实测的风致雪漂移的CFD数值模拟。Step 8: Post-processing calculation for numerical simulation of snow deposition erosion: Based on the wall shear force calculated in step 7: Based on the wall shear force calculated in step 7, extract the wall shear force on the bottom of the model through Tecplot, and use a custom function Convert the wall shear stress to the friction velocity of the wall, and compare the wall friction velocity with the threshold friction velocity. When the wall friction velocity is greater than the threshold friction velocity, it is erosion, using the formula Calculate the amount of erosion, where A ero is a constant coefficient, taking 7.0×10 -4 ; when the wall friction velocity is less than the threshold friction velocity, the formula is used for deposition Calculate the amount of deposition, where w f is the sedimentation velocity of snow, φ s is the mass concentration of snow, φ s = ρ s f, f is the volume fraction of snow, and ρ s is the density of snow particles; then use the formula Calculate the changed thickness of snow cover, where q s is q ero or q dep , γ is the maximum volume fraction of snow cover, which is taken as 0.62; finally, superimpose the changed thickness of snow cover with the initial snow cover depth of the building in step 1 to obtain the The final distribution of snow, complete the CFD numerical simulation combined with the measured wind-induced snow drift.
实施例2:本实施例在未结合实测、不考虑风相和雪相的相对滑移速度条件下对风致雪漂移的数值进行模拟,具体的方法按以下步骤进行:Embodiment 2: This embodiment simulates the numerical value of wind-induced snow drift without combining actual measurement and without considering the relative slip velocity of wind phase and snow phase. The specific method is carried out according to the following steps:
步骤一:收集数值模拟需要的基本信息:具体研究的建筑物的几何特征:边长为1m的立方体模型,该建筑所处的地理位置地貌参数为0.16、参考高度10m处平均风速5m/s,风剖面采用公式z与v0分别表示参考点的高度(1m)及对应的参考平均风速(5m/s);α表示地貌参数,对于B类地貌,α取值为0.16;湍动能采用公式湍能耗散率采用公式其中,Cμ为经验常数,取0.09;I和l分别为湍流强度和湍流长度尺度,湍动能强度采用公式积分湍流尺度采用公式其中,zG、zb、α'、均为地貌参数,在日本AIJ规范中,B类地貌对应Ⅱ类平地,以上三个地貌参数取值分别为5m、350m和0.15。Step 1: Collect the basic information required for numerical simulation: the geometric characteristics of the building under study: a cube model with a side length of 1m, the geographic location of the building where the geomorphological parameters are 0.16, the average wind speed at a reference height of 10m is 5m/s, The wind profile adopts the formula z and v 0 represent the height of the reference point (1m) and the corresponding reference average wind speed (5m/s), respectively; α represents the geomorphic parameter, for Type B landform, the value of α is 0.16; the turbulent kinetic energy adopts the formula Turbulent energy dissipation rate using the formula Among them, C μ is an empirical constant, which is taken as 0.09; I and l are the turbulent flow intensity and turbulent flow length scale, respectively, and the turbulent kinetic energy intensity adopts the formula The integral turbulence scale adopts the formula Among them, z G , z b , α', are all geomorphological parameters. In the Japanese AIJ code, the B-type landform corresponds to the II-level flat land, and the values of the above three geomorphic parameters are 5m, 350m and 0.15 respectively.
通过经验采用相关数值:雪颗粒直径0.00015m、不考虑休止角、雪密度150kg/m3、计算阈值摩擦速度为0.2m/s、跃移雪粒速度0.2m/s、地面粗糙度系数为0.25、雪粒沉降速度0.2m/s。Relevant numerical values are adopted through experience: the diameter of snow particles is 0.00015m, the angle of repose is not considered, the snow density is 150kg/m 3 , the calculated threshold friction speed is 0.2m/s, the snow particle speed is 0.2m/s, and the ground roughness coefficient is 0.25 , Snow particle settling velocity 0.2m/s.
步骤二:建立数值模拟需要的网格模型:在Ansys-Icem CFD中按照建筑物在步骤一所收集的建筑物几何信息和初始积雪分布形态进行模型的建立以及网格的划分,确定模型计算域,长、宽、高分别为20m、10m、8m,且计算域设置阻塞率为1.25%,接着在Icem CFD中依据模型建立网格模型,依据模型大小网格数为360万。网格模型控制结构化网格质量在0.8,非结构化网格质量在0.8,最后以Ansys-Fluent格式输出求解的网格文件。Step 2: Establish the grid model required for numerical simulation: In Ansys-Icem CFD, the model is established and the grid is divided according to the building geometric information and the initial snow distribution shape collected in step 1, and the model calculation is determined. The length, width and height are 20m, 10m and 8m respectively, and the blockage rate of the computational domain is set to 1.25%. Then, a grid model is established based on the model in Icem CFD, and the number of grids is 3.6 million according to the size of the model. The mesh model controls the structured mesh quality at 0.8, the unstructured mesh quality at 0.8, and finally outputs the solved mesh file in Ansys-Fluent format.
步骤三:编写数值模拟用的空气项边界条件:利用步骤一所收集的建筑物所在地理位置信息、参考高度平均风速、风剖面、湍动能、湍动强度、积分湍流尺度、湍流耗散率规范信息编写空气相的边界条件,以现有技术在Visual Studio中以C语言、Fluent二次开发接口为基础,编写指数型风剖面、湍动能、湍流耗散率,保存C语言文件,并将其与步骤二网格模型文件放在同一文件夹中;Step 3: Write the boundary conditions of the air term for numerical simulation: use the geographic location information of the building collected in Step 1, the average wind speed at the reference height, the wind profile, the turbulent kinetic energy, the turbulent intensity, the integral turbulence scale, and the turbulent dissipation rate specification Write the boundary conditions of the air phase based on the existing technology in Visual Studio, based on the C language and the Fluent secondary development interface, write the exponential wind profile, turbulent kinetic energy, and turbulent dissipation rate, save the C language file, and put it Put it in the same folder as the
步骤四:编写数值模拟用的雪相边界条件:利用步骤一所收集的雪颗粒密度150kg/m3、当地重力加速度9.81m/s2、阈值磨擦速度0.2m/s,采用Pomeroy&Gray(1990)、Pomeroy&Male(1992)所提出的雪浓度经验公式,即公式其中ρs表示雪颗粒密度,g为重力加速度,取9.81m/s2,up为跃移层雪粒平均速度,up=2.8u*t,u*为入口摩擦速度,根据Beyers等实测拟合得到的摩擦速度与10m高度处风速u10的关系确定,雪相运动速度取5m/s、雪颗粒摩擦速度取0.15m/s、u*t为阈值速度;取0.15m/s;hs为跃移层和悬移层的临界高度,雪相运动速度5m/s、雪颗粒摩擦速度0.15m/s、在Visual Studio中以C语言、Fluent二次开发接口为基础编写雪相体积分数公式作为雪项的入口边界条件,保存C语言文件,并将其与步骤二模型网格文件、步骤三的C语言文件放在同一文件夹中;Step 4: Write the snow phase boundary conditions for numerical simulation: using the snow particle density of 150kg/m 3 collected in step 1, the local gravitational acceleration 9.81m/s 2 , and the threshold friction velocity of 0.2m/s, using Pomeroy & Gray (1990), The empirical formula for snow concentration proposed by Pomeroy & Male (1992) is the formula where ρ s is the density of snow particles, g is the acceleration of gravity, 9.81m/s 2 is taken, u p is the average speed of snow particles in the transition layer, u p =2.8u *t , u * is the friction velocity at the entrance, according to the actual measurement of Beyers et al. The relationship between the friction speed obtained by fitting and the wind speed u 10 at a height of 10m is determined, The speed of snow phase movement is 5m/s, the friction speed of snow particles is 0.15m/s, u *t is the threshold speed; 0.15m/s is taken; h s is the critical height of the transition layer and the suspension layer, Snow phase motion speed is 5m/s, snow particle friction speed is 0.15m/s, and the snow phase volume fraction formula is written in Visual Studio based on C language and Fluent secondary development interface as the inlet boundary condition of the snow term, and the C language file is saved , and put it in the same folder as the model grid file in
步骤五:数值模拟需要的程序接入,具体操作如下:将步骤二的网格文件、步骤三、四的自定义程序存放相同文件夹,打开Ansys-Fluent界面选择3D、双精度、并行。在GeneralOptions-Working Directory下选择上述文件夹目录,选择环境下Set up CompilationEnvironment for UDF接口打开软件。选择Ansys-Fluent主菜单下User-Defined-FunctionsStep 5: Program access required for numerical simulation. The specific operations are as follows: save the grid file in
-Compiled下选择步骤三、四的C语言文件进行加载编译,得到空气相及雪相的入口速度及体积分数;- Under Compiled, select the C language files in
步骤六:将数值模拟软件的模型、材料、求解器参数进行具体设置,采用Ansys-Fluent导入网格模型、基于步骤五导入的二次开发编程进行风雪两相流模拟,具体如下:Step 6: Set the model, material and solver parameters of the numerical simulation software in detail, use Ansys-Fluent to import the grid model, and simulate the two-phase flow of wind and snow based on the secondary development programming imported in
1)在General下Time一栏选择瞬态Transient;1) Select transient Transient in the Time column under General;
2)在模型Models下多项流Multiphase一栏选择Mixture模型、在模型Models下Viscous一栏选择k-epsilon两方程模型下的可实现Realizable k-epsilon模型,在模型Models下多项流Multiphase一栏选择Phase Interaction-slip,此时不考虑风相与雪相的相对滑移速度;2) Select the Mixture model in the Multiphase column under Models, select the Realizable k-epsilon model under the k-epsilon two-equation model in the Viscous column under Models, and select the Multiphase multiphase model under Models. When Phase Interaction-slip is selected, the relative slip velocity of the wind phase and the snow phase is not considered at this time;
3)在材料Materials下选择流体Fluid,添加空气相,设置Viscosity粘度、Density密度,并且新建材料雪相,设置雪相Viscosity粘度与空气一致,依据步骤一设置Density密度150kg/m3、粒径0.00015m;3) Under the material Materials, select Fluid, add the air phase, set the Viscosity viscosity and Density density, and create a new material snow phase, set the snow phase Viscosity viscosity to be the same as air, and set the Density density to 150kg/m 3 and the particle size to 0.00015 according to step 1. m;
4)在体网格Cell Zone Conditions下的网格计算域内设置需要的流体fluid、固体Solid,在Boundary Conditions设置边界条件:入口边界为速度入口、计算域两侧及顶部为对称边界、地面及模型表面为无滑移壁面、出口边界为自由出口;4) Set the required fluid fluid and solid solid in the grid computing domain under the volume grid Cell Zone Conditions, and set the boundary conditions in the Boundary Conditions: the inlet boundary is the velocity inlet, the sides and top of the computational domain are symmetrical boundaries, the ground and the model The surface is a non-slip wall, and the outlet boundary is a free outlet;
5)在Boundary Conditions-inlet(入口)接入步骤三的空气相边界条件:设置风剖面、湍动能及耗散率,在Boundary Conditions-inlet(入口)-snow雪相-Multiphase下Volume Fraction接入步骤四的雪相体积分数;Fluent-Mixture模型采用有限体积法计算,由于风雪流只有两相,且两相的积分分数和为1,因此只需给出雪相的体积分数,空气相由软件计算无需设置;5) In the Boundary Conditions-inlet (inlet), access the air phase boundary conditions of step 3: set the wind profile, turbulent kinetic energy and dissipation rate, and access the Volume Fraction in the Boundary Conditions-inlet (inlet)-snow snow phase-Multiphase The volume fraction of the snow phase in
6)在求解器设置Solution下选择方法Methods中选用SIMPLE算法。6) Select the SIMPLE algorithm in Methods under the solver setting Solution.
7)在求解器设置Solution下选择Monitors-Residual中设置连续性方程、动量方程以及k-epsilon的收敛精度,本实施例设置收敛精度为1e-6精度。7) Under the solver setting Solution, select Monitors-Residual to set the convergence precision of the continuity equation, the momentum equation and the k-epsilon. In this embodiment, the convergence precision is set to 1e- 6 precision.
8)在求解器设置Solution下初始化Initialization中选用标准初始化StandardInitialization在入口初始化。8) Select the standard initialization StandardInitialization in the initialization initialization under the solver setting Solution to initialize at the entrance.
9)最后计算在Run Calculation设置计算的计算步数步数为1000步,精度为10-6,开始计算,得到壁面剪切力、空气相及雪相的速度及体积分数的收敛文件。9) Final calculation Set the number of calculation steps in Run Calculation to 1000 steps, and the accuracy is 10 -6 , start the calculation, and obtain the convergence files of the wall shear force, the velocity and volume fraction of the air phase and the snow phase.
步骤七:利用步骤六计算出的收敛文件,进行壁面剪切力和速度的后处理:在Ansys-Fluent计算收敛后文件以Tecplot支持格式保存。用Tecplot后处理计算,通过提取计算的x、y、z的速度、压力、壁面剪切力的分量,用Tecplot自定义函数Specify Equations,做出计算域范围为内的压力、速度、壁面剪切力的云图、等值线图。Step 7: Use the convergence file calculated in Step 6 to post-process the wall shear force and velocity: After the Ansys-Fluent calculation converges, the file is saved in the Tecplot supported format. Post-processing calculation with Tecplot, by extracting the calculated components of x, y, z velocity, pressure, and wall shear force, and use Tecplot's custom function Specify Equations to make pressure, velocity, and wall shear within the computational domain. Force cloud map, contour map.
步骤八:数值模拟雪沉积侵蚀的后处理计算:基于步骤七计算获得的壁面剪切力,通过Tecplot提取模型底面上壁面剪切力,通过自定义函数把壁面剪切应力转换为壁面的摩擦速度,对壁面摩擦速度与阈值磨擦速度进行比较,当壁面摩擦速度大于阈值磨擦速度时侵蚀,采用公式其中Aero为常系数,取7.0×10-4;当壁面摩擦速度小于阈值磨擦速度时沉积,采用公式其中wf中为雪的沉降速度,φs为雪的质量浓度(φs=ρsf,f为雪的体积分数,ρs为雪颗粒密度)计算侵蚀量、沉积量;基于步骤一的初始积雪深度,采用公式计算积雪改变厚度、叠加积雪改变量(其中,统一用qs代替qero和qdep,γ为积雪的最大体积分数,取0.62),得到积雪的最终分布形态,完成风致雪漂移的数值模拟。Step 8: Post-processing calculation for numerical simulation of snow deposition erosion: Based on the wall shear force calculated in Step 7, extract the wall shear force on the bottom of the model through Tecplot, and convert the wall shear stress into the friction velocity of the wall through a custom function , the wall friction velocity is compared with the threshold friction velocity. When the wall friction velocity is greater than the threshold friction velocity, the erosion is carried out using the formula where A ero is a constant coefficient, which is taken as 7.0×10 -4 ; when the wall friction velocity is less than the threshold friction velocity, the deposition is made using the formula Wherein w f is the sedimentation velocity of snow, φ s is the mass concentration of snow (φ s = ρ s f, f is the volume fraction of snow, ρ s is the density of snow particles) to calculate the amount of erosion and deposition; Initial snow depth, using the formula Calculate the change thickness of snow cover and the change amount of superimposed snow cover (where q s is used to replace q ero and q dep , and γ is the maximum volume fraction of snow cover, which is taken as 0.62) to obtain the final distribution form of snow cover, and complete the wind-induced snow drift numerical simulation.
将实施例1的模拟结果显示的单位时间内的积雪通量、实施例2在未与实测相结合、不考虑风相和雪相的相对滑移速度的条件下得到的单位时间内的积雪通量及经典实验的单位时间内的积雪通量曲线图进行比较,如图4所示,从图4可以看出,实施例2在未与实测相结合、不考虑风相和雪相的相对滑移速度的条件下得到的单位时间内的积雪通量与经典实验方法得到的曲线相差较大,而实施例1的模拟结果与经典试验的结果能很好地吻合,说明实施例1的模拟方法的准确性高。The snow flux per unit time shown by the simulation results of Example 1, and the snow accumulation per unit time obtained in Example 2 without combining with the actual measurement and ignoring the relative slip velocities of the wind and snow phases The snow flux and the snow flux curve per unit time of the classic experiment are compared, as shown in Figure 4. It can be seen from Figure 4 that Example 2 is not combined with the actual measurement, and does not consider the wind phase and snow phase. The snow flux per unit time obtained under the condition of relative slip velocity is quite different from the curve obtained by the classical experimental method, and the simulation results of Example 1 are in good agreement with the results of the classical experiment. The accuracy of the simulation method of 1 is high.
Claims (8)
Priority Applications (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210155363.7A CN114528781B (en) | 2022-02-21 | 2022-02-21 | A CFD numerical simulation method for wind-induced snow drift combined with actual measurements |
Applications Claiming Priority (1)
Application Number | Priority Date | Filing Date | Title |
---|---|---|---|
CN202210155363.7A CN114528781B (en) | 2022-02-21 | 2022-02-21 | A CFD numerical simulation method for wind-induced snow drift combined with actual measurements |
Publications (2)
Publication Number | Publication Date |
---|---|
CN114528781A true CN114528781A (en) | 2022-05-24 |
CN114528781B CN114528781B (en) | 2025-01-14 |
Family
ID=81625510
Family Applications (1)
Application Number | Title | Priority Date | Filing Date |
---|---|---|---|
CN202210155363.7A Active CN114528781B (en) | 2022-02-21 | 2022-02-21 | A CFD numerical simulation method for wind-induced snow drift combined with actual measurements |
Country Status (1)
Country | Link |
---|---|
CN (1) | CN114528781B (en) |
Cited By (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115186569A (en) * | 2022-07-07 | 2022-10-14 | 中科三清科技有限公司 | Floating object drift simulation method and device, storage medium and electronic equipment |
CN118261020A (en) * | 2024-05-31 | 2024-06-28 | 北京航空航天大学 | A method, device, medium and product for predicting the shedding of high-altitude ice crystal accretion |
CN118364554A (en) * | 2024-05-17 | 2024-07-19 | 哈尔滨工业大学 | Integrated wind-induced snow drift numerical simulation method for early stage of building design |
Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20070219766A1 (en) * | 2006-03-17 | 2007-09-20 | Andrew Duggleby | Computational fluid dynamics (CFD) coprocessor-enhanced system and method |
CN109406741A (en) * | 2018-10-24 | 2019-03-01 | 黑龙江大学 | A kind of measurement method of snow load |
CN111008425A (en) * | 2019-12-09 | 2020-04-14 | 桂林理工大学 | Wind and snow resistance design method for large-span hyperbolic roof |
-
2022
- 2022-02-21 CN CN202210155363.7A patent/CN114528781B/en active Active
Patent Citations (3)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
US20070219766A1 (en) * | 2006-03-17 | 2007-09-20 | Andrew Duggleby | Computational fluid dynamics (CFD) coprocessor-enhanced system and method |
CN109406741A (en) * | 2018-10-24 | 2019-03-01 | 黑龙江大学 | A kind of measurement method of snow load |
CN111008425A (en) * | 2019-12-09 | 2020-04-14 | 桂林理工大学 | Wind and snow resistance design method for large-span hyperbolic roof |
Cited By (5)
Publication number | Priority date | Publication date | Assignee | Title |
---|---|---|---|---|
CN115186569A (en) * | 2022-07-07 | 2022-10-14 | 中科三清科技有限公司 | Floating object drift simulation method and device, storage medium and electronic equipment |
CN115186569B (en) * | 2022-07-07 | 2023-01-31 | 中科三清科技有限公司 | Floating object drift simulation method and device, storage medium and electronic equipment |
CN118364554A (en) * | 2024-05-17 | 2024-07-19 | 哈尔滨工业大学 | Integrated wind-induced snow drift numerical simulation method for early stage of building design |
CN118261020A (en) * | 2024-05-31 | 2024-06-28 | 北京航空航天大学 | A method, device, medium and product for predicting the shedding of high-altitude ice crystal accretion |
CN118261020B (en) * | 2024-05-31 | 2024-07-23 | 北京航空航天大学 | Method, device, medium and product for predicting falling off of ice crystal ice deposit in high altitude |
Also Published As
Publication number | Publication date |
---|---|
CN114528781B (en) | 2025-01-14 |
Similar Documents
Publication | Publication Date | Title |
---|---|---|
CN114528781A (en) | CFD numerical simulation method combining actual measurement wind-induced snow drift | |
Shao et al. | PIGNN-CFD: A physics-informed graph neural network for rapid predicting urban wind field defined on unstructured mesh | |
Beyers et al. | Modeling transient snowdrift development around complex three-dimensional structures | |
Liu et al. | LES study on the turbulent flow fields over complex terrain covered by vegetation canopy | |
Haider et al. | Urban flood modelling using computational fluid dynamics | |
Koblitz et al. | Computational Fluid Dynamics model of stratified atmospheric boundary‐layer flow | |
Kocaer et al. | Experimental and numerical investigation of flow over ogee spillway | |
Gorlé et al. | Dispersion in the wake of a rectangular building: validation of two Reynolds-averaged Navier–Stokes modelling approaches | |
CN117172152A (en) | A CFD dynamic grid numerical simulation method based on measured wind-induced snow drift | |
Moreira et al. | Numerical study of the neutral atmospheric boundary layer over complex terrain | |
Chen et al. | DriftScalarDyFoam: an OpenFOAM-based multistage solver for drifting snow and its distribution around buildings | |
CN117952032A (en) | A method for numerical simulation of wind-snow flow on roof using dynamic grids considering wind-snow coupling | |
CN110990926B (en) | Urban surface building hydrodynamic simulation method based on area correction rate | |
CN114925624A (en) | Natural river channel three-dimensional water flow numerical simulation method | |
KR102392067B1 (en) | 3-dimensihonal wind flow analyzing system by stages of development in step-up street canyons using commutational fluid dynamics, and analyzing method using the same | |
Freire et al. | A one-dimensional stochastic model of turbulence within and above plant canopies | |
Wright et al. | Environmental applications of computational fluid dynamics | |
Rethore et al. | A CFD model of the wake of an offshore wind farm: using a prescribed wake inflow | |
Bechmann et al. | Atmospheric flow over terrain using hybrid RANS/LES | |
Basirat et al. | Eulerian–Eulerian model application to simulate scouring downstream of sluice gate | |
Ansari | Boundary shear stress distribution and flow structures in trapezoidal channels | |
Hu et al. | Research on KMC-based evolution simulation of dust particles in virtual campus environment | |
Sridhar | Numerical prediction of wind flow over complex terrain with shallow and steep hills | |
Silva et al. | Analysis and validation of a CFD simulation of the wind through a horizontal axis wind turbine as an actuator disc with a porous jump condition | |
Weng et al. | Study on CFD Simulation of Wind Filed in the Near-Earth Boundary Layer |
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 |