流体模拟《Fluid Simulation for Computer Graphics 》——从入门到放弃(二)模拟视角/Navier-Stokes流体方程/MAC法

从入门到放弃系列

深刻认识到了好好学数值分析的重要性,,所以人生就是【被劝】-【不听】-【踩坑】-【劝人】……

哪里理解错误啥的请指出3Q


1 Simulation Viewpoints(模拟视角)

跟踪运动的流体或可变形固体(连续体)的运动有两种途径:欧拉视角,拉格朗日视角。

  1. 拉格朗日视角

在这种视角下,流体由空间中的无限多个点组成,每个点都标记为一个单独的颗粒,有各自的位置向量,速度向量。

  1. 欧拉视角

在空间上保持固定点,并且查看这些点所处的流体(速度/温度等)如何随时间变化,这相当于使用固定的网格。在这种视角下,更容易在固定的欧拉网格上对一些流体量的空间导数进行分析及对这些导数进行数值逼近。

2 Fluid Equations(流体方程)

大多数流体由不可压流体的Navier-Stokes方程(纳维一斯托克斯方程)决定:

2.1 Symbols

向量u表示流体速度,ρ表示流体密度,p表示压力(流体在任何物体上施加的单位面积的力),向量g表示重力加速度,(三维空间里表示为(0,-9.81,0)m/s^2),ν表示运动粘度(用于测量流体的粘性)。

2.2 动量方程

把上面的(1)式分为三个部分:

  1. Material derivative(物质导数)

其中q表示通过流体传输的模拟量之一(比如温度),t表示时间。∂q/∂t就表示了q在固定点变化的速度,这是一种欧拉度量。第二项u⋅∇q纠正多少变化是由于流体流过的差异。

比如在一个模拟空间中,温度在原点处为0,温度变化为T(x)=20x(即向右变热),流体以恒定速度c流动,水分子的温度不变,那么根据拉格朗日视角,温度的物质导数为:

也就是说在空间固定处,温度以-20c的恒定速率变化,如果水流速度为0,那么固定点温度不变,如果水以恒定速度c=2向右流动,那么固定点的温度将以-40的速度下降。

我们可以知道在拉格朗日观点下每个水分子的温度都没有变化,因此导数为0,而欧拉视角下导数值取决于流体流动的速度和方向。

当上面说的q是速度时,称为自对流,因为此时速度有两种不同的作用,即流体运动的速度和随着流体平移而变的速度,这使我们得出与动量方程相同的方程:

(论文中的对流指由于流体在特定方向上的运动而导致的物质通过流体的传输机制,在流体模拟领域中用于描述模拟量q的物质导数为0)

  1. 作用在流体上的外力

重力加速度g

  1. 作用在流体上的内力

内力是影响流体粒子与附近其他粒子相互作用的两种流体力:

  • 第一种:1\ρ ∇p体现了压力对流体的影响。压力可以使流体体积保持稳定且无速度耗散。高压流体区域推动低压流体区域,我们可以通过负压力梯度-∇p来测量粒子处的压力不平衡(因为梯度运算符代表了最陡的上升方向,而高压区域向低压区域推进,因此要在每个点取负压力梯度)。
  • 第二种:第二种内力是流体粘度引起的。粘性流体试图抵抗变形,这种内力试图使这个粒子以周围粒子的平均速度运动,也就是使邻近粒子之间的速度差异最小化。所以∇⋅∇u就是用拉普拉斯算子衡量一个量与周围平均数的差。

2.3 不可压缩条件

真实的流体,甚至是水这样的液体体积也会改变,但通常变化不大。这样的研究对象称为可压缩流,研究起来很复杂,而宏观上的影响又不大,所以就应用层面来说没必要考虑。

所以在流体模拟中一个前提条件是流体体积变化率为0(即不可压缩),也就是前文的(2)式。

3 粘度

在大多数情况下粘度只有很小的影响,所以也可以不考虑,而归入误差范围。不考虑粘度的Navier-Stokes方程称为欧拉方程,理想的没有粘度的流体称为inviscid,方程是:

还有两个边界情况,实心墙和自由表面:

  • 实心墙:流体不能流进或流出实心墙,所以速度的法向分量为0,即)。如果是移动的固体的话那么就是固体速度的法向分量
  • 自由表面:自由表面是停止对流体建模的地方,在模拟区域外可能存在另一种流体(比如空气)。因为只有压力的差异(理想条件不可压缩体积不变),所以自由表面p=0

4 数值模拟

4.1 流体算法

把上面的理想情况(不可压缩无粘度)的Navier-Stokes方程(7)进行进一步拆分,对每一部分分别求解:(拆分依据就是2.2中的物质导数/外力/内力,其中(9)式就是物质导数为0,(10)式那个就是重力方向(y方向)求偏导,(11)式那个就是压力反方向求偏导,不可压缩)

基础流体算法:

前三步说的是,从初始无散度速度0开始;然后时间n=0,1,2,……,也就是t_0,t_1,t_2,...;然后确定一个Δt,也就是t_n到t_(n+1),后三步就,,0.0应该说的是一个速度,然后算上重力影响后的速度,然后根据时间变化的速度

4.2 MAC法(标记网格法)

为了使空间离散化,就引入一个叫MAC法的好东西。

引入MAC具体可以看这篇论文:【F. Harlow and J. Welch. “Numerical Calculation of Time-Dependent Viscous Incompressible Flow of Fluid with Free Surface.” Phys. Fluids 8 (1965), 2182–2189.】

,,我是看不动了QWQ,,就看看百度出来的各种介绍吧,,

MAC法:

  • marked-and-cell,标记网格法,一种偏微分方程数值解法
  • 用差分法和标记点相结合求解不可压缩自由表面流动问题
  • 用Navier-Stokes方程
  • 计算区域按照欧拉矩形网格单元划分
  • 压力和速度作为基本未知量
  • 压力定义在格子中心,也就是(i,j)
  • r方向速度分量定义在格子左右边界中点,也就是(i,j)(i+1,j)的速度是u_(i+1/2,j)
  • y方向速度分量定义在格子上下边界中点,同理,(i,j)(i,j+1)的速度是v_(i,j+1/2)
  • 多维也是同理的
  • 用差分方法由动量方程和泊松方程分别计算出速度和压力
  • 在网格中布置适量的无质量的标记点,这些点来确定自由表面
  • 位置由拉格朗日坐标表示(欧拉坐标就是固定在空间中的参考坐标,我们一般理解的那种;拉格朗日坐标是嵌在物体质点上,随物体一起运动和变形的坐标,这时候的坐标值就相当于一个质点的姓名,只和那个质点有关,坐标系性质可能会变但是那个质点的坐标值是不随物体的运动或者变形改变的。拉格朗日坐标的优点就是推导公式方便,求导的时候坐标值不变只要对时间求导就好了,用欧拉坐标的话坐标值也会变就得考虑质点的运动)
  • 使用双变量线性插值计算标记点的速度
  • 计算过程中跟踪每个标记点,从而确定自由表面的形状/位置/演变清况

既然看到了MAC(标记网格法),就顺便提一下另一个PIC(质点网格法),这个是计算二维非定常可压缩理想流动问题的欧拉-拉格朗日混合方法,特别适用于计算具有多种介质和大变形流动的问题,具体可以自己百度……


先结个尾,下一篇继续

打赏
  • 版权声明: 本博客所有文章除特别声明外,著作权归作者所有。转载请注明出处!
  • Copyrights © 2018-2021 LeFlacon

奶茶一杯 快乐起飞

支付宝
微信