一种涉及热传导与动态黏度的真实感流体仿真方法技术

本发明专利技术公开了一种涉及热传导与动态黏度的真实感流体仿真方法，其步骤为：1）基于光滑粒子流体动力学（SPH）模型，对流体内部及流体与外界的热传导过程进行离散建模，并根据相变焓对相变温度的影响以模拟相变过程；2）将动态黏度的计算引入以展示流体运动过程中的细节；3）调用PCISPH算法完成剩余的流体运动模拟过程；4）利用GPU加速的算法，在通用并行计算架构CUDA上将热传导、相变、流体黏度变化等过程并行处理，实现真实感相变流体的快速仿真。应用本发明专利技术能够真实、高效地模拟不同流体与外界的热传导及黏度变化过程，增强了现有方法中的模拟细节，提升了流体仿真的真实感。

【技术实现步骤摘要】

本专利技术属于计算机图形学领域，具体地说是一种涉及热传导与动态黏度的真实感流体仿真方法，其部分技术包括SPH流体模拟方法，热传导与相变原理以及GPU并行加速等。
技术介绍
流体运动作为自然界普遍的物理现象，其仿真技术一直以来都是计算机图形学领域一个重要的研究方向，并且在影视制作、游戏开发、虚拟现实等领域都有着广泛的应用前景。基于物理的流体模拟主要通过求解流体力学中经典的纳维-斯托克斯(N-S)方程来进行，目前较流行的两种方法分别为拉格朗日法及欧拉法。其中拉格朗日法使用粒子系统模拟流体，概念易于理解，并且在处理细节及较大形变等方面表现突出，是目前应用较为广泛的一种模拟方法。光滑粒子流体动力学(SPH)方法作为一种拉格朗日方法常用来模拟流体的运动规律。考虑到多数流体的不可压缩性，SPH的许多改进算法中逐渐完善了这一问题，如WCSPH、PCISPH方法等。其中PCISPH方法是目前应用较多的一种方法，其利用了预测-校正技术能够取得时间步长较大且每步迭代开销较低的优点，大幅提升了模拟效率，是一种非常适合模拟不可压缩流体的方法。然而随着对真实感的进一步需求，仅仅考察流体的速度、位移已满足不了大众的视觉体验。因而流体的温度、状态等其他物理属性也必须考虑，以展示流体运动的细节。同时，当某些流体的温度发生改变时，其自身黏度往往也随之改变，尤其在诸如火山爆发、蜂蜜溶解、巧克力熔化等场景中体现得尤为明显。所以如何真实地模拟同种流体在不同温度下的黏性差异也是影响上述场景真实感仿真的一大因素，如果能很好地改善上述流体模拟过程中的细节，则对大众的视觉体验来说将会是一个巨大的提升。
技术实现思路
本专利技术的目的在于将现实生活中存在于流体中的热传导过程及黏度变化引入到基于物理的流体仿真过程中，提出了一种涉及热传导与动态黏度的真实感流体仿真方法，该方法引入了温度对流体黏度的影响，对流体热传导过程中的细节进行完善，获得了更加逼真的视觉效果。实现本专利技术目的的具体技术方案是：一种涉及热传导与动态黏度的真实感流体仿真方法，包括以下步骤：a)基于光滑粒子流体动力学(SPH)模型模拟流体及热传导过程，具体包括：ⅰ)邻居查找首先将所有流体及固体用离散粒子表示，并查找得到与每个粒子i距离在光滑核半径h内的邻居粒子集合Ni；ⅱ)各粒子导热系数计算利用SPH算法的插值方式，考虑距离不同的邻居粒子对某个粒子i的影响，计算不同状态下粒子i的导热系数ki，其具体公式为：其中Phasei表示粒子i的当前状态，mj、ρj、kj分别表示粒子i的邻居粒子j的质量、密度、导热系数，Wij为形如的光滑核函数，它是关于粒子i与j的位置xi与xj之间距离||xi-xj||及光滑核半径h的函数；kfluid与ksolid分别为该种流体在液态及固态下实际的导热系数；iii)各粒子温度变化率计算首先计算流体内部的热传导，对每个粒子i，其温度随时间变化率计算公式为： dT i d t = Δ ( k i T i ) ]]>其中Ti及ki分别为粒子i的温度及导热系数，Δ为拉普拉斯算子；同时Δ(kiTi)的计算可通过公式Δ(kiTi)＝及得到，其中表示梯度；其次计算流体与外界的热传导，对每个粒子i，其温度随时间变化率计算公式为： ( dT i d t ) e x t = Σ j k i m J R i 2 ρ i ρ j ( T b - T i ) W i j ]]>其中Tb为固体粒子温度，Ri可通过计算得到，ρ0为流体的静态密度；然后更新粒子i的温度Ti，具体公式为： T i ( t + Δ t ) = T i ( t ) + [ ( dT i d t ) e x t + dT i d t ] Δ t ]]>其中Δt为时间步长，Ti(t)与Ti(t+Δt)分别为粒子i在时刻t与t+Δt的温度值；iv)相变处理由于相变焓的存在，分别设置两个温度阈值Tmelt和Tsolid；当粒子温度Ti大于Tmelt时，粒子转化为液体粒子，当Ti小于Tsolid时，粒子转化为固体粒子；当Ti处于Tsolid与Tmelt之间时，则将其看作为临界状态的流体粒子；b)动态黏度计算的引入，具体包括：在上述过程得到每个粒子i的最新时刻温度Ti后，由于不同温度对流体黏性系数存在一定的影响，在SPH模型中引入动态黏度的计算，实现流体温度变化过程中细节的逼真模拟，其中温度Ti与流体黏性系数μi的变化关系具体为：log(log(μi+γ))＝q-y log(Ti)其中γ为本文档来自技高网...

【技术保护点】
一种涉及热传导与动态黏度的真实感流体仿真方法，其特征在于该方法包括以下步骤：a)基于光滑粒子流体动力学(SPH)模型模拟流体及热传导过程，具体包括：ⅰ)邻居查找首先将所有流体及固体用离散粒子表示，并查找得到与每个粒子i距离在光滑核半径h内的邻居粒子集合Ni；ⅱ)各粒子导热系数计算利用SPH算法的插值方式，考虑距离不同的邻居粒子对某个粒子i的影响，计算不同状态下粒子i的导热系数ki，其具体公式为：其中Phasei表示粒子i的当前状态，mj、ρj、kj分别表示粒子i的邻居粒子j的质量、密度、导热系数，Wij为形如的光滑核函数，它是关于粒子i与j的位置xi与xj之间距离||xi‑xj||及光滑核半径h的函数；kfluid与ksolid分别为该种流体在液态及固态下实际的导热系数；ⅲ)各粒子温度变化率计算首先计算流体内部的热传导，对每个粒子i，其温度随时间变化率计算公式为：dTidt=Δ(kiTi)]]>其中Ti及ki分别为粒子i的温度及导热系数，Δ为拉普拉斯算子；同时Δ(kiTi)的计算可通过公式及得到，其中表示梯度；其次计算流体与外界的热传导，对每个粒子i，其温度随时间变化率计算公式为：(dTidt)ext=ΣjkimjRi2ρiρj(Tb-Ti)Wij]]>其中Tb为固体粒子温度，Ri可通过计算得到，ρ0为流体的静态密度；然后更新粒子i的温度Ti，具体公式为：Ti(t+Δt)=Ti(t)+[(dTidt)ext+dTidt]Δt]]>其中Δt为时间步长，Ti(t)与Ti(t+Δt)分别为粒子i在时刻t与t+Δt的温度值；ⅳ)相变处理由于相变焓的存在，分别设置两个温度阈值Tmelt和Tsolid；当粒子温度Ti大于Tmelt时，粒子转化为液体粒子，当Ti小于Tsolid时，粒子转化为固体粒子；当Ti处于Tsolid与Tmelt之间时，则将其看作为临界状态的流体粒子；b)动态黏度计算的引入，具体包括：在上述过程得到每个粒子i的最新时刻温度Ti后，由于不同温度对流体黏性系数存在一定的影响，在SPH模型中引入动态黏度的计算，实现流体温度变化过程中细节的逼真模拟，其中温度Ti与流体黏性系数μi的变化关系具体为：log(log(μi+γ))＝q‑y log(Ti)其中γ为常数，通常取值在0.6到0.9之间，q与y则为与流体自身有关的特定参数，为定值；c)PCISPH算法调用利用PCISPH算法，对每个粒子i分别计算其所受到压力、黏性力以及外力而产生的加速度aipre、aivis及aiother，其计算公式分别为：aipre=-1miΣj(mjpi+pj2ρj▿Wij)]]>aivis=12miΣj(mj(μi+μj)vj-viρj▿2Wij)]]>aiother＝g其中，其中m、p、ρ、μ、v分别表示粒子的质量、压强、密度、黏性系数及速度，g为重力加速度；然后通过粒子i所受的加速度计算得到其在下一时刻的位置以模拟流体运动变化规律；d)GPU并行加速，具体包括：利用通用并行计算架构CUDA，将所有粒子的物理属性传输至GPU，并对粒子进行排序，加速邻居查找过程；对于每个粒子，为其开辟一个独立的线程以计算其密度、温度、导热系数、加速度，并更新流体粒子的速度与位置，从而并行地提升流体仿真的效率。...

【技术特征摘要】
1.一种涉及热传导与动态黏度的真实感流体仿真方法，其特征在于该方法包括以下步骤：a)基于光滑粒子流体动力学(SPH)模型模拟流体及热传导过程，具体包括：ⅰ)邻居查找首先将所有流体及固体用离散粒子表示，并查找得到与每个粒子i距离在光滑核半径h内的邻居粒子集合Ni；ⅱ)各粒子导热系数计算利用SPH算法的插值方式，考虑距离不同的邻居粒子对某个粒子i的影响，计算不同状态下粒子i的导热系数ki，其具体公式为：其中Phasei表示粒子i的当前状态，mj、ρj、kj分别表示粒子i的邻居粒子j的质量、密度、导热系数，Wij为形如的光滑核函数，它是关于粒子i与j的位置xi与xj之间距离||xi-xj||及光滑核半径h的函数；kfluid与ksolid分别为该种流体在液态及固态下实际的导热系数；ⅲ)各粒子温度变化率计算首先计算流体内部的热传导，对每个粒子i，其温度随时间变化率计算公式为： dT i d t = Δ ( k i T i ) ]]>其中Ti及ki分别为粒子i的温度及导热系数，Δ为拉普拉斯算子；同时Δ(kiTi)的计算可通过公式及得到，其中表示梯度；其次计算流体与外界的热传导，对每个粒子i，其温度随时间变化率计算公式为： ( dT i d t ) e x t = Σ j k i m j R i 2 ρ i ρ j ( T b - T i ) W i j ]]>其中Tb为固体粒子温度，Ri可通过计算得到，ρ0为流体的静态密度；然后更新粒子i的温度Ti，具体公式为： T i ( t + Δ t ) = T i ( t ) + [ ( dT i d t ) e x t + dT i d t ...

【专利技术属性】
技术研发人员：王长波，张申帆，孔凡龙，李晨，
申请(专利权)人：华东师范大学，
类型：发明
国别省市：上海;31