读书笔记:A Model for Soap Film Dynamics with Evolving Thickness

3 SOAP FILM DYNAMICS WITH VARYING THICKNESS

3.1 Definitions and Setting

流体在薄膜上的流动遵循带有表面张力项的无粘NS方程:

这个x上面一点就是速度的意思。(1)右端第一项是表面张力。里面的u代表气-液界面,l代表液-气界面。H是平均曲率,n是法向。这里u和l得名于“upper”和“lower”.

然后\(\delta_{S^u},\delta_{S^l}\)是狄拉克函数,在\(S^u,S^l\)之外消逝。

压强p在除边界之外的地方二阶可导,在每个泡泡里面的空气区域中取恒定值。 S是surface:

这里面A是\(R^2\)的子集。然后在一个点可以确定局部坐标系:

这里\(S_a\)和\(S_b\)代表偏导。

厚度场为h,因此薄膜上下表面为:

然后它定义了一个\(\Gamma\):

这个\(\Gamma\)实际上是一个“柱”,也就是以S(a,b)表示的这个表面点为中心,向上下到薄膜表面的一个线段。

然后X就是这个“柱”的重心:

文中将用X代表表面S。

3.2 Normal Acceleration

因为厚度h足够小,所以可以忽略其上的变化,用表面上点的速度代表这个“柱”的速度:

单看法向加速度的话,根据NS方程,就是:

最后一个等号是根据狄拉克函数的性质,以及如下等式:

这里面,\(\nabla^2_s\)是表面S上的Laplace-Beltrami算子。我们认为薄膜两边的body force一样,因为它足够薄。然后也假设\(\nabla_{S^u}^2S^u+\nabla_{S^l}^2S^l=2\nabla_S^2S\).

然后考虑压强项。在薄膜表面,也就是“柱”\(\Gamma\)的两端,根据Ishida et al’s [2017]的近似:

也就是说,压强差乘以法向。然后把中间的压强z方向梯度积分起来,其实也等于两端压强差:

我们认为法向在整个“柱”上都是一样的,也就是\(n^u=n^l=n\),结果就是:

其实就是两端压强差,这么点事……

这样的话,合起来的加速度公式就是:

影响这个法向运动的主要是:表面张力(曲率)、两端压强差、外力。

3.3 Tangential Acceleration

切向加速度就是:

这里面我们认为压强p是二阶可导的,同样\(n^u=n^l=n\).

然后由于法向的x,y分量为零,也就是:

这里面的\(\nabla_S\)是表面梯度(surface gradient),P是一个“柱”的压强积分:

注意:切向没有表面张力。也就是说,文章只算毛细力,忽略了Marangoni力。

3.4 Thickness Evolution

首先回顾厚度的定义:

然后我们看,在\(x^u\)处的法向力为:

类似地

把这两个合到一起算出来h的变化率:

这里把上下表面的表面梯度之差写成了\(\nabla_S^2h\)。然后

可以发现,如果忽略掉Q,h的方程就是一个advected wave equation.然后Q的作用是“挤压”这个film,保持不可压性,就像流体求解器里面的压强场一样。

3.5 Pressure

压差\(\nabla p\)实际上是每个泡泡体积守恒产生的一个修正项。我们使用Ishida et al. [2017]的方法做这个。

然后我们看方程(20)里的P。根据质量守恒,厚度h满足:

使用一个类似投影法的方法。首先对(20)做对流,然后加上外力,然后剩下:

这里带波浪线的是对流&加外力之后的速度场。把它和(26)联立可以得到:

这个就相当于N-S求解器里的投影了。这个和拿triangle mesh做的二维不可压流体求解器几乎一致,但是它允许质量在切向速度和法向的h变化率之间转意。所以,速度场本身是可压的,散度会被h的变化率吸收掉。

本文的算法没有真正地计算Q。它是用Q去修正数值误差导致的体积变化。每个时间步结束之后,用Q去scale up厚度项h,保证薄膜的质量守恒。

3.6 Evaporation and Bursting Bubbles

蒸发就是每步把h削薄一点:

破裂就是当某个地方的h太小的时候,删掉一块表面。把这个bubble的质量分给周围其他bubble.

4 IMPLEMENTATION

大致可以总结为:

1、切向advection.
2、根据切向速度和h算P。
3、根据算出来的P更新切向速度。
4、Advect厚度。
5、给h加表面张力。
6、用Q对h做体积修正。
7、计算法向速度,其中有表面张力,需要用到h。
8、做空气体积守恒。

发表评论

电子邮件地址不会被公开。 必填项已用*标注