一个自由液面的瞬态问题：Level set还是动网格？

soliton

看来大家有点跑题了，我再来探讨下此问题解决的思路。
首先，使用的是pde系数型和相场法两相流 ，因为我觉得相场法比较稳定，水平集法也可这样参考。
设常数b1=1，b2=-1，标量表达式b=b1+(b2-b1)*Vf2_chns；很明显 ...
facklaoer 发表于 2010-9-9 22:03

最近，液滴在纯粹表面张力下回复-变形的例子，尽管还不完全真实，但有解决的苗头。

开始试一下facklaoer大侠这个想法。

初步的结果显示有一点点问题：

在纯粹的两相中，phi=0和phi=1都没有问题。问题就是界面上的邻域（过渡层）上，b有值，但既不等于1也不等于-1，此时PDE方程中系数c/a/f都有值。就这一层没有满足上面界面上T=0的条件，解出来数值差了好远。

试了直接赋值，效果也不好。还在思考中。

soliton

看来大家有点跑题了，我再来探讨下此问题解决的思路。
首先，使用的是pde系数型和相场法两相流 ，因为我觉得相场法比较稳定，水平集法也可这样参考。
设常数b1=1，b2=-1，标量表达式b=b1+(b2-b1)*Vf2_chns；很明显 ...
facklaoer 发表于 2010-9-9 22:03

试了两天，不得不承认，fachlaoer大侠所提出的这个想法，似乎是不行了。

我前两天在逻辑上，对这个想法深信不疑。以为是通过PDE方程系数的形式，巧妙地把那个“移动的边界”所要满足的边界条件施加进去的一个有效的方法。但计算的结果完全不对。

没错，在a等于零，c=f=1（或者-1）时，计算域都是满足△T=-1的。问题是在界面上，即使a=1,c=f=0，理所应当得出T=0。但这里有个陷阱：就是事实上，作为控制方程，这里整个计算域上都会算出T=0，而不单单是界面上！

所以，问题依然存在：怎么把这个边界条件施加到一个演化的界面上？

求高手指点！

soliton

给facklaoer大侠写的几句话： 一眨眼，这个问题都做了两个半月了。而且还不一定能够解决。有一点点沮丧。但这个问题确实有它好玩的地方，能够让人放弃一段时间的主业，而不断地玩一下。 我只能说大概解决了这个问题的50%。一个相关的帖子是：http://forum.simwe.com/thread-948772-1-6.html。那里是初始形状为椭圆的情况：没有所谓“电磁力”的时候，一个椭圆振荡着最后变成圆，加上这个帖子一开始所列的那个电磁力之后，还是震荡着变成圆，只不过振荡的幅值要大，时间要长。那个帖子中给了动画。 关键，现在我还不能弄出由“某个初始形状”由电磁力变成各种花样的计算。 退了一步，弄了个初始就为“五边形”的计算。我有动画，但是那是通过POSTMOVIE在matlab里对N个fem后处理之后拼出来的一个.avi的格式，比较大不好传上来。（我用了动网格中stop condition_remesh_restart这一技巧） 截几个图：

soliton

理论上，这种模式五的周期是0.51秒不到。上面的几张图给的是几个过渡状态：0.25s，0.26s到0.27秒速度的反向(箭头)，然后到0.50秒之后恢复到初始的形状。（目前在模式三时，是完全类似的结果。模式四计算起来比较困难点，主要是收敛的问题）

但是，上面的结果是没有那个poisson方程（也就是没有“电磁力”）的结果！

加上电磁力之后的表现是：由于N-S方程中多出来一项体积力，计算的收敛性变得很差。总是不能顺利地算下去。我做了一个极端的例子，在那个F的公式中，将磁场强度设为零（也就是没有电磁力了！),而在模型中保留那个poisson方程，但计算仍然不收敛！我真是百思不得其解了。

也就是说：
在我现在的这个模型里面有三个物理场：（1）动网格；（2）N-S方程；（3）POISSON方程。（3）所算出的结果以某个体积力在流场中起作用。也就是（3）对（1）（2）单向起作用。

我在没有（3）的时候，能算出结果来。在将（3）所求出的解设置为“零”的时候（也就是（3）完全不起作用的时候），却算不出结果了！

为什么？

soliton

把这个帖子翻出来，给自己提提神！

从上个星期开始，这个问题有彻底解决的可能了。在某些情况下，由十分简单的初始形状，能“自发”地在电磁力和表面张力作用下，“涌现”出实验中观察到的一些“花样”。至少最近十天所得到的结果都和理论的、实验的预测相当，其误差是实际三维简化为二维之后的合理结果。

尽管还没有到最终解决的程度，但可以初步总结几点：1）动网格；2）表面张力和电磁力通过弱解形式写入边界；3）网格的划分要不停试，要适应速度场的大小；4）中间以网格质量为标准停下来，再重新生成计算域再restart（这个还没有自动）；5）电磁力的耦合，没有facklaoer担心的那么复杂，因为是二维的而已。

三个月了，确实不容易！可能还要弄三个月。太艰苦了。

*午夜流星*

看了半天，竟然不知道你想表达什么意思。

“就是在Level set中，怎么把每时每刻的那个边界条件T=0写进去？因为这个边界是不断变化的，并且在Comsol中，这个界面是内部的边界，且和解出来的函数Phi有关。”

Level set本来就是追踪界面变形的，不知道你所谓的“边界条件T=0”什么意思。

另外你这个物理过程有点像分形里面的“指进”现象，就是本来应该均匀变化的，由于局部振荡，或者压力或者粘度等的不均匀，形成“指进”，生成比较漂亮的图形，一般随机性比较强，模拟起来不容易。因为不太了解具体你的物理内涵，也不好多说，可能我的看法跟你的研究完全不是一回事情。

如果真是“指进”，那模拟起来是很有难度的，除非你人为的增加振荡或者由于数值误差产生局部的不均匀，并且这种差异逐渐放大。

soliton

看了半天，竟然不知道你想表达什么意思。

“就是在Level set中，怎么把每时每刻的那个边界条件T=0写进去？因为这个边界是不断变化的，并且在Comsol中，这个界面是内部的边界，且和解出来的函数Phi有关。”

Lev ...
*午夜流星* 发表于 2010-12-15 17:56

谢谢午版关注。其实那个所谓“边界”t=0的条件，是针对液滴和空气之间那个不断变化的“界面”而言的。水平集的好处是实际有界面，而“心中”无界面（由标量场phi来定），而偏偏那个表征电磁力的关于T的方程，就需要一个真实的物理的界面！矛盾在这里，不知道我说清楚了没有。Facklaoer大侠第一时间就给了一个解决方案，我仔细学习实践后，发现是有漏洞的（上面写到了）。

两个月前，我就放弃level set 做这个问题了。难点有二：（1）计算的收敛非常差。事实上我后来找到一个电磁场的方程来取代那个Delta^2(T)=-1的标量场方程，进而可以算这个椭圆抖动的过程，但速度奇慢，不能忍受。（2）水平集存在的相泄漏，导致液滴越算越小。幸亏当初放弃了，因为现在看来，有时候液滴要抖动、变形几十秒。用水平集算这个问题，一个参数就要一个月吧。

至于你说的“指进”，我觉得有相通之处，但不完全相同。我看过有些分形混沌的书，知道你说的那个意思。

我说相通，就是在我这个系统内，存在“正反馈”。垂直于这个金属液滴施加低频的交流磁场后，由于电磁感应，感应电流会在液滴平面内“环形流动”，电流的散度为零，这决定了“凹”下去的地方，电力线的密度会大一些，所以此处电磁力比别处要大（这个力在界面上某些时刻指向外，某些时刻指向内）。如果凹的地方（它有回复为圆形的冲动），这时候电磁力指向内部，还比较大，那就会更凹，结果电磁力更大，更凹。。。这个就是正反馈了。

本质上讲，这是个稳定性分析的问题，就是什么条件下是稳定的（圆形），什么时候是不稳定（花样）。回答这个问题，要用“确定性”的数值模拟，方法论上我也觉得有问题。这也就是我这几天，机器不停地算，想要干地事情：遍历参数空间，勾勒稳定性那个所谓“map”，确实是有难度。

这个里面，有没有混沌，我直觉可能有。但我现在这个算法，可能算不到那个参数里面去：动网格不允许有结构”拓扑“上的变化。但是那些花样中的一些，我确实是得到了。

*午夜流星*

个人感觉，你似乎对LEVEL SET的算法没有理解。 LEVEL SET 本来就是用一个虚拟的薄体积来代替面，实际也就是具有了厚度，从物理意义来说，还是具有“真实性”的，因为所谓的多相流的界面，实际上也不绝对就是理想的面，也是有过渡厚度的，这个暂且不谈。
LEVEL SET 中物理量实际是在厚度内过渡的。所谓“边界”t=0的条件，你说是真实的物理的条件，我不知道现在的计算方法是否可以计算，如果完全理想的话，好像没有办法计算，因为表面张力实际也是个面应力，也是按照体积力来计算的，好像基本上流体软件都是这么处理的。因为数值不连续，突然间断的话，很难处理。
我个人认为，“边界”t=0的条件，实际还是应该可以处理成过渡条件的，本身来说应该不影响你的物理意义。
我这几年用COMSOL计算介电泳力，也是一个面应力，都是这么处理的。

可能我还是没有理解你的物理模型，不过只要你结果能计算的不错，就行了。不过COMSOL计算，有时候不正确的设置会产生很“奇妙”的结果，这个要反复验证。

soliton

个人感觉，你似乎对LEVEL SET的算法没有理解。 LEVEL SET 本来就是用一个虚拟的薄体积来代替面，实际也就是具有了厚度，从物理意义来说，还是具有“真实性”的，因为所谓的多相流的界面，实际上也不绝对就是理想的面 ...
*午夜流星* 发表于 2010-12-16 08:17

午版真是高手中的高手！

我完全理解你的意思。“LEVEL SET 本来就是用一个虚拟的薄体积来代替面，实际也就是具有了厚度”。其实我说：“水平集的好处是实际有界面，而“心中”无界面（由标量场phi来定）”，这里的实际是物理实际，而“心中”就是指它的算法上，界面是个“虚”的东西，是按照phi来确定的（通过这个phi的值来定厚度的）。大概是你说的这个意思。

但问题的难点不在这里。在于我在主贴中写的关于T的那个方程，“只针对液态金属部分”成立，所谓T=0的边界，也只是在金属液“边界”上成立，和周围的部分一点关系也没有。这下好，如果还要用水平集的话，这个时刻变化的边界的边界条件如何加？难点在这里。

我后来没办法抛弃了这个PDE，通过“参数化几何”，由phi来定金属和空气的电导率（界面上还是要光滑），直接解电磁场的方程，才把这个电磁场＋水平集的问题弄得差不多。但问题是计算速度和相泄漏的问题，才逼我去用动网格的。

再说动网格上表面张力的处理。是借助于弱解形式，把表面张力写到方程系统的边界上（不是直接施加到边界条件上的）来处理的。（我有专门的帖子请教大家这个问题）我数学底子不好，其实这里还是有点问题没有完全想明白。但这个思路，已经和水平集没有任何关系了。

soliton

后面的一段时间，重心在动网格和弱解上了。

那个表面张力弱解处理的帖子在这里：http://forum.simwe.com/thread-955092-1-1.html。一直没有人来答疑解惑。我的问题，目前完全就只剩下这个帖子里提出来的了。欢迎高手来指点一二。

soliton

秀一下最终的计算结果，同时另开贴请教一些问题。 这个问题确实神奇。我的机器没日没夜算了两个月，确实是，磁场的频率和磁场强度决定了会出来什么花样！下面的图初始条件都是一样的，从一个略微变形的椭圆开始，在这个简单的电磁场作用下，我几乎得到了所有想要的2-3-4-5-6-7种模式。

但为什么会出来这些模式，是一个太专门的问题，我刚刚才明白了一个大概。用数值模拟的方法来做一个稳定性分析问题，确实不是个好的做法。

shoney

我最近也在琢磨comsol中的这个level set，我的是跟一个heat transfer耦合的，想在移动界面上施加一个温度的条件，即T=T1(T1是一个常数)，另外想得到穿过这个界面的热流差。一直没有好方法解决这两个问题。

soliton

32# shoney

如果用level set，想在移动的边界上加进去这个温度和热流，我个人的感觉，真的有困难。不如尝试动网格？

如果Level Set 可以，千万上来介绍一下经验！谢谢！