一种基于sph算法的流体热运动仿真系统的制作方法
【专利摘要】本发明涉及一种基于SPH算法的流体热运动仿真系统，属于热力学模拟与可视化流动领域。该方法包括如下步骤：首先根据计算实例和计算要求对流体进行建模，然后基于SPH算法求解流体控制方程和修正方程，实现流体拉格朗日热运动的仿真模拟，最后应用可视化处理技术以图形及动画的方式直观的展示流体的热运动过程。与现有网格化CFD仿真系统相比，在模拟某些高度非线性问题(如射流问题、海洋运动问题、溃坝问题、动态边界问题)时，会得到更为合理的仿真结果；与现有无网格CFD仿真系统相比，提升了计算稳定性和计算精度，并且填补了流体热运动仿真和湍流仿真的空白。
【专利说明】一种基于SPH算法的流体热运动仿真系统
【技术领域】
[0001] 本发明为一种基于光滑粒子流体动力学算法（SPH)的流体热运动仿真系统，属于 热力学模拟与可视化流动领域。
【背景技术】
[0002] 流体的热运动仿真在实际工程中有着广泛的应用背景，例如建筑物隔热、核反应 堆安全设计、微电子设备的冷却等。
[0003] 光滑粒子流体动力学方法是一种基于粒子的无网格数值模拟方法.该方法被应 用于解决各种流体力学和固体力学问题，包括流固耦合问题、凝固问题、射流问题、多相流 问题、热传导问题等.光滑粒子动力学方法作为最经典的无网格方法，其算法被不断完善， 多种修正技术被提出并应用，并且发展了相应的前处理与后处理技术，其已经成为较为完 善的粒子型无网格方法.
[0004] 近年来，使用传统的网格化计算流体动力学（CFD)仿真系统在解决实际流体工程 问题时遇到了许多难点，其模拟高度非线性问题（如射流问题、海洋运动问题、溃坝问题、 动态边界问题）时往往会受到网格变形的制约，无法实现某些高度非线性问题的模拟，这 些高度非线性问题往往需要通过采用无网格CFD方法来解决；而现有的无网格CFD仿真系 统又无法解决热力学问题。
【发明内容】
[0005] 本发明目的是为了解决现有网格化CFD仿真系统在解决某些高度非线性问题时 因为网格的变形导致结果不够准确以及现有无网格CFD仿真系统无法对流体的热运动进 行仿真的问题，提供一种基于SPH算法的流体热运动仿真系统。
[0006] 本发明技术方案的思想是：根据计算实例和计算要求，首先对流体进行建模，然后 基于SPH算法求解流体控制方程和修正方程，实现流体拉格朗日热运动的仿真模拟，最后 应用可视化处理技术以图形及动画的方式直观的展示流体的热运动过程。
[0007] 本发明是通过以下技术方案实现的：
[0008] -种基于SPH算法的流体热运动仿真系统，包括前处理模块、计算模块和后处理 模块，整个仿真系统的模拟流程如图1所示。
[0009] 1.前处理模块
[0010] 所述前处理模块，用于根据实际模拟目标建立CFD仿真模型，一般CFD仿真模型需 具备以下内容：离散方程、离散方法和离散模型，本仿真系统中离散方程为流体的流动控制 方程（质量守恒方程、动量守恒方程、能量守恒方程），尚散方法为SPH算法，尚散模型为根 据用户要求构建的粒子模型；因此，系统将前处理模块分为以下两个子模块：SPH算法设置 和粒子模型建模；
[0011] 1)SPH算法设置
[0012] 前处理模块中的SPH算法设置将决定如何使用SPH算法，其包括以下四部分内容： 计算方式设置、SPH算法格式设置、修正方法设置和SPH算法参数设置；
[0013] a)计算方式设置决定了该仿真系统的计算方式，以充分利用计算资源，其可选类 型包括：中央处理器（CPU)并行、图形处理器（GPU)并行；其中CPU并行通过标准消息传递 界面（MPI)技术实现，GPU并行通过通用并行计算架构（CUDA)实现。
[0014] b)若要使用SPH算法解算流体控制方程，则需要根据模拟对象来确定如何使用 SPH算法，即需要设置SPH算法的格式；SPH算法的格式包括：时间积分格式、核函数格式、 粘性解算格式；其中时间积分格式决定了模拟时间推进的方式，可选格式包括：蛙跳法、韦 尔莱积分法、辛积分法；核函数为SPH算法的最基本函数，其决定了 SPH算法如何离散流动 控制方程，可选格式包括：三次样条核函数、五次样条核函数；粘性解算格式为动量守恒方 程中流体粘性的计算方式，可选格式包括：层流粘性、人工粘性、亚粒子模型湍流粘性（湍 流模型）；
[0015] C)由于SPH算法易于产生震荡误差，则需要在算法中引入一些修正方法，修正 方法包括：密度重构技术，人工位移技术（Xu R，Stansby P, Laurence D.Accuracy and stability in incompressible SPH(ISPH)based on the projection method and a new approach [J]· Journal of Computational Physics, 2009, 228 (18) :6703-6725);其中密度 重构技术用于修正密度计算的震荡误差，可选类型包括：薛帕德滤波、S -SPH ;人工位移技 术用于修正粒子位置的计算误差。
[0016] d)SPH算法参数包括：粒子间距离、时间步长、模拟时长、核函数的影响半径（光滑 长度）、修正方程参数、人工声速（人工声速为人工状态方程的参数，人工状态方程用于压 力求解）。
[0017] 2)粒子模型建模
[0018] 粒子模型建模模块用于对粒子模型进行建模，包括粒子坐标、粒子属性（流体或 者固体）、边界条件和初始条件；
[0019] 所述粒子坐标的生成通过输入相应外形参数实现，根据外形参数的几何特征及给 定的粒子间距离，便可求得该几何形状内所有粒子的坐标，例如：按顺序输入矩形的四个顶 点坐标，即可生成矩形粒子域；
[0020] 所述粒子属性设置即是将边界粒子属性设置为固体，其他粒子属性设置为流体的 过程；
[0021] 所述边界条件确定了边界固体粒子对内部流体粒子的作用方式，可设置以下两方 面：固体壁面边界条件、温度边界条件；
[0022] 计算域边界部分粒子的属性为固体粒子，边界粒子位置不变，并对靠近的流体粒 子产生一定作用力，作用力施加方式可选类型包括：边界力法、镜像粒子法；其中边界力法 通过施加排斥力实现固体边界；镜像粒子法通过施加排斥动量实现固体边界；
[0023] 当存在具有一定温度特性的边界时，该部分边界粒子需要设置温度边界条件，以 使其温度特性影响其附近粒子，温度边界条件可选类型包括：恒温边界条件、绝热边界条 件；恒温边界条件情况下，该部分固体粒子温度恒定，持续与周围粒子进行热交换；绝热边 界条件情况下，该部分固体粒子与其他粒子不发生任何热交换。
[0024] 所述初始条件为计算域内粒子的初始物理参数，包括：速度、压力、温度、参考温 度、密度、运动粘度、重力、瑞利数、普朗特数。
[0025] 当完成了前处理模块中的SPH算法设置和粒子模型建模，便可将数据导入计算模 块，进行流体热运动的解算。
[0026] 2.计算模块
[0027] 所述计算模块用于使用SPH算法对流体的流动控制方程和修正方程进行求解，得 到各时间步的各物理参数的数值模拟结果。该模块首先读入前处理模块所建立的仿真模 型，从而获得SPH算法设置信息和粒子模型信息，进而开始SPH计算循环，在每个计算循环 中，首先SPH算法会将流体的流动控制方程（质量守恒方程、动量守恒方程、能量守恒）施 加到计算域内所有的粒子上；然后所有粒子的物理参数会按照SPH算法离散后的控制方程 产生变化；最后当所有粒子完成一次计算后，求解前处理模块中选取的修正方程，使需要修 正的物理参数产生修正变化，由此完成所有粒子物理参数的更新，该时间步完成；当完成粒 子状态更新后，按照时间积分格式更新时间步，并同时输出包含所有粒子物理参数的数据 文件，然后进入下一循环，直至时间步迭代到模拟时长结束。
[0028] 3.后处理模块
[0029] 所述后处理模块，用于对计算模块输出的仿真数据进行动态可视化处理，以直观 展现流体热运动过程中流体的运动过程及各个参数的变化过程。该模块的处理过程如下： 首先读入计算模块输出的数据文件，进而根据每个时间步仿真数据中的位置、速度、温度、 密度等物理参数，将粒子位置分布以坐标点的形式表示其分布，得到整个计算域内粒子分 布图，将其他物理参数用连续变化的颜色表示其分布，得到整个计算域内物理参数分布云 图，最后利用所有时间步的粒子分布图及物理参数分布云图生成整个模拟过程的粒子运动 动画和各物理参数变化动画。
[0030] 有益效果
[0031] 本无网格CFD仿真系统由于采用了 SPH算法，较之已有的网格化CFD仿真系统，在 模拟某些高度非线性问题（如射流问题、海洋运动问题、溃坝问题、动态边界问题）时，会得 到更为合理的仿真结果；
[0032] 本无网格CFD仿真系统较之已有的无网格CFD仿真系统存在如下优点：
[0033] 1)添加了能量守恒方程的求解，从而使本仿真系统得以解决热力学问题；
[0034] 2)添加了多种修正方程的求解，从而提升了系统的计算稳定性和计算精度；
[0035] 3)添加了湍流模型，从而使本仿真系统可以解决流体热运动湍流问题。
【专利附图】
【附图说明】
[0036] 图1是本发明系统的组成结构以及数据流向示意图；
[0037] 图2为本发明实施例的封闭方腔自然对流物理模型；
[0038] 图3为本发明实施例的粒子模型；
[0039] 图4为本发明实施例某一时间步所生成的流场温度分布云图。
【具体实施方式】
[0040] 为了使本发明的目的、技术方案和优点更加清晰明白，下面结合实施例和附图，对 本发明实施例做进一步的详细说明。
[0041] 下面以封闭方腔自然对流模拟为例介绍本发明基于SPH算法的流体热运动仿真 系统的实施过程。
[0042] 按图1所示流程实现该实施例的仿真模拟，具体实现步骤如下：
[0043] 1.前处理模块包括：
[0044] 1) SPH算法设置
[0045] 计算方式：GPU并行。
[0046] SPH算法格式设置如下，时间积分格式：娃跳法；核函数格式：五次样条核函数；粘 性结算格式：层流粘性。
[0047] 修正方法：人工位移技术。
[0048] SPH算法参数：粒子间距离dx = 0. 01、时间步长10'模拟时长0. 5、核函数的影响 半径（光滑长度）h = 1. 5dx = 0. 015、人工位移技术强度参数0. 04、人工声速400。
[0049] 2)粒子模型建模
[0050] 图2为封闭方腔自然对流物理模型，其中L(L= 1)表示物理模型边长，其中 $： = 〇表示该部分粒子采用了恒温边界条件，τΗ和t分别表示恒高温和恒低温，g表示重 力的方向。
[0051] 粒子坐标和粒子属性：以图2中正方形的四个顶点坐标为控制点，建立如图3所示 的粒子模型，内部正方形域设置为流体粒子，正方形域外围设置三层边界（固体）粒子。
[0052] 边界条件：四个壁面均采用固体壁面边界条件，采用边界力法；左右壁面采用恒 温温度边界条件，上下壁面采用绝热温度边界条件。
[0053] 初始条件：速度0、压力0、左侧壁面温度为1，右侧壁面温度为0、参考温度为0. 5、 密度1. 225、运动粘度v = 10'重力g = 9. 81、瑞利数Ra = 105、普朗特数Pr = 0· 71 (以 上均为无量纲形式）。
[0054] 2.计算模块进行以下处理：
[0055] 1)读入前处理模块输入的信息；
[0056] 2)根据前处理模块所输入的信息，设置SPH算法与粒子模型；
[0057] 3)开始循环，根据所选计算方法与参数利用SPH算法解算质量守恒方程、动量守 恒方程、能量守恒方程和人工位移技术方程；
[0058] 4)根据解算结果，更新粒子信息，包括：速度、位置、温度、密度；
[0059] 5)输出该时间步的数据文件；
[0060] 6)通过时间更新进入下一循环，直至时间步迭代到模拟时长循环结束。
[0061] 3.后处理模块进行以下处理：
[0062] 1)读入计算模块生成的数据；
[0063] 2)生成单个时间步的粒子分布图及物理参数分布云图，如图4所示为根据某一时 间步数据生成的流场温度分布云图，可以看出所得温度分布较为合理；
[0064] 3)用所有时间步的粒子分布图及物理参数分布云图生成整个模拟过程的粒子运 动动画和各物理参数变化动画，以动态展示流场的变化情况。
[0065] 以上所述仅为本发明的具体实施例，并不用于限定本发明的保护范围，凡在本发 明的精神和原则之内，所做的任何修改、等同替换、改进等，均应包含在本发明的保护范围 之内。
【权利要求】
1. 一种基于SPH算法的流体热运动仿真系统，其特征在于，包括前处理模块、计算模 块和后处理模块；所述前处理模块用于根据实际模拟目标建立CFD仿真模型，所述计算模 块用于使用SPH算法对流体的流动控制方程和修正方程进行求解，得到各时间步的仿真数 据，所述后处理模块用于对计算模块输出的仿真数据进行动态可视化处理，以直观展现流 体热运动过程中流体的运动过程及各个参数的变化过程。
2. 根据权利要求1所述的一种基于SPH算法的流体热运动仿真系统，其特征在于，前处 理模块包含SPH算法设置和粒子模型建模两个子模块； 所述SPH算法设置模块决定计算模块如何使用SPH算法，包括以下内容：计算方式设 置、SPH算法格式设置、修正方法设置和SPH算法参数设置； a) 计算方式设置决定了该仿真系统的计算方式，包括中央处理器（CPU)并行和图形处 理器（GPU)并行两种计算方式；其中CPU并行通过标准消息传递界面（MPI)技术实现，GPU 并行通过通用并行计算架构（CUDA)实现； b) SPH算法格式设置包括时间积分格式、核函数格式和粘性解算格式三项内容；其中 时间积分格式决定了模拟时间推进的方式，可选格式包括：蛙跳法、韦尔莱积分法、辛积分 法；核函数为SPH算法的最基本函数，其决定了 SPH算法如何离散流动控制方程，可选格式 包括：三次样条核函数、五次样条核函数；粘性解算格式为动量守恒方程中流体粘性的计 算方式，可选格式包括：层流粘性、人工粘性、亚粒子模型湍流粘性； c) 由于SPH算法易于产生震荡误差，则需要在算法中引入一些修正方法，修正方法包 括：密度重构技术，人工位移技术；其中密度重构技术用于修正密度计算的震荡误差，可选 类型包括：薛帕德滤波、S -SPH ;人工位移技术用于修正粒子位置的计算误差； d) SPH算法参数包括：粒子间距离、时间步长、模拟时长、核函数的影响半径、修正方程 参数、和人工声速； 所述粒子模型建模模块用于对粒子模型进行建模，其包括以下四部分内容：粒子坐标、 粒子属性、边界条件和初始条件； a) 粒子坐标的生成通过输入相应外形参数实现，根据外形参数的几何特征及给定的粒 子间距离，便可求得该几何形状内所有粒子的坐标； b) 粒子属性设置即是将边界粒子属性设置为固体，其他粒子属性设置为流体的过程； c) 边界条件确定了边界固体粒子对内部流体粒子的作用方式，包括两方面内容：固体 壁面边界条件、温度边界条件； 固体壁面边界条件即是计算域边界部分粒子的属性为固体粒子，边界粒子位置不变， 并对靠近的流体粒子产生一定作用力，作用力施加方式包括边界力法和镜像粒子法两种方 式；其中边界力法通过施加排斥力实现固体边界；镜像粒子法通过施加排斥动量实现固体 边界； 温度边界条件即是当存在具有一定温度特性的边界时，该部分边界粒子需要设置温度 边界条件，以使其温度特性影响其附近流体粒子，温度边界条件包括恒温边界条件和绝热 边界条件两种；恒温边界条件情况下，该部分固体粒子温度恒定，持续与周围粒子进行热交 换；绝热边界条件情况下，该部分固体粒子与其他粒子不发生任何热交换； d) 初始条件为计算域内粒子的初始物理参数，包括：速度、压力、温度、参考温度、密 度、运动粘度、重力、瑞利数、普朗特数。
【文档编号】G06F17/50GK104143027SQ201410377151
【公开日】2014年11月12日 申请日期:2014年8月1日 优先权日:2014年8月1日
【发明者】雷娟棉, 杨浩, 黄灿, 王锁柱 申请人:北京理工大学