9.22 更新【Ti example 投稿】taichi实现等离子体模拟——以双流不稳定性为例

9.22 更新

之前的投稿内容过于繁琐,这里重新提交一份简单易懂的~

使用Taichi实现了一套简短的等离子体particle-in-cell code。
代码链接:
JiaoLuhuai/pic88: a short pic code (github.com)
这里是以“双流不稳定性”现象为例进行演示的。(简单来说,双流不稳定性是这样一种现象,等离子体中有两束有明显速度差异的粒子,在满足一定条件下会产生不稳定的现象然后渐渐趋于稳定。这里无意讨论物理机制,仅作演示。)
需要说明的是为了简化计算,模拟中的空间是1维的。

效果如下
(图中横坐标为粒子位置,纵坐标为粒子速度。)

初始状态:

出现不稳定性的状态:

逐渐趋于稳定的状态:

总的动态图:
pic88

–分割线–

以下是初版:

使用taichi实现了一套电磁的等离子体模拟particle-in-cell code,将taichi的应用拓展到等离子体模拟领域。这里以“双流不稳定性”这种物理现象为例进行展示。(由于等离子体模拟很难有直观酷炫的视觉效果,双流不稳定性已经是矮子里面拔高个儿了。)(双流不稳定性简单来说就是等离子体中有两束速度有明显差异的粒子,在满足一定条件下会产生不稳定现象激发出波动然后又趋于稳定。这里无意讨论物理机制,仅用模拟中的粒子状态来演示。)
(实现思路主要参考京都大学开发的KEMPO1,见 space.rish.kyoto-u.ac.jp/software/
(需要说明的是为了简化计算,模拟中空间是1维的,只有x坐标,速度是3维的。这也是等离子体PIC模拟中一种常见的处理手段。)

演示双流不稳定性的32位精度代码:
taichi-PIC-code/two_stream_instability.py at main · JiaoLuhuai/taichi-PIC-code (github.com)
这里给出32位精度代码一方面是考虑支持不同的GPU驱动,另一方面是控制行数在300行以内。(在下面有给出64位精度的代码,比300行多一点点,但泛用性更高,可以模拟多种等离子体物理过程,如果能用于评example就太好了。)

效果如下:
初始状态:

出现不稳定性的状态:

逐渐趋于稳定的状态:


(需要注意的是图中横坐标代表粒子的x坐标,纵坐标代表粒子的速度vx,不同颜色代表速度不同的两束粒子。)

支持多种等离子体物理过程模拟的64位精度代码(默认设置为模拟双流不稳定性):
taichi-PIC-code/kempo1_ti.py at main · JiaoLuhuai/taichi-PIC-code (github.com)
这里给出的64位精度代码用于与KEMPO1(MATLAB)对比。

taichi版本64位精度模拟code与KEMPO1(MATLAB)对比效果如下:
taichi版本与matlab版本效果对比
taichi版本与matlab版本效果对比_2
taichi版本与matlab版本效果对比_3
可以看出二者看不出差异,证明了taichi代码的可靠性。

更详细的64位精度模拟code与KEMPO1(MATLAB)的模拟效果GIF对比见 taichi-PIC-code/taichi版本与matlab版本效果对比.gif at main · JiaoLuhuai/taichi-PIC-code (github.com)

5 个赞

欢迎投稿~~ :star_struck:

1 个赞

哇!这个很有新意!!

很棒!可否请你调整一下代码格式 (比如使用 yapf 等自动格式化工具),让代码风格更规范?

顺便,你的代码里面有很多 (ux, uy, uz),ren? 之类名字和含义相似的变量,大量使用这种全局变量比较影响代码的可读性,是否可以考虑把相似的变量 (假设它们是一个向量在不同方向的分量) 放在一个 ti.math.vec3 中?

感谢您的建议,这十分有帮助!

1.代码格式已经借助yapf进行调整。
2. (ux,uy,uz)含义是u=gamma*v,其中gamma是洛伦兹因子(相对论因子)。它是作为速度求解中的一个中间变量,并非全局变量。目前已添加相关注释。
3. 代码中renorm部分大量出现"ren?"格式的变量名,这是代表对应物理量的归一化系数。它们并非全局变量。只在进行初始化的python scope中使用了它们,taichi kernel中并未用到。

代码文件已在GitHub中更新,再次感谢您的建议。

请问楼主,改变 kempo1_ti.py的哪些输入,还能有其他的模拟么?

比如,改变iex=1(电磁),与iex=2 ,似乎效果上没有区别
(iex = 1 #1电磁,2静电,0静磁)

因为背景磁场是0,背景磁场由第一种粒子的回旋频率wc决定,wc=0,背景磁场也为0.

有,具体参考http://space.rish.kyoto-u.ac.jp/software/
里面有个PDF,有介绍原理和一些算例

1 个赞

有熟悉taichi的GUI,能帮忙把楼主的代码,加两句(比如gui.set_image(pixels)+ ti.VideoManager),做成gif或mp4么?在下无能,看楼主的代码后,也没摸着头脑,如何把楼主的图像变成视频,有请高手添把柴,谢谢了~~

还是不行,我把wc=20000,效果和wc=0一样一样,不知原因

楼主,网友zemora的意思是ux、uy、uz(假如是一个向量的三个分量),他建议用ti.math.vec3来存放ux,uy,uz,这样代码可读性会比较好,重点不在是否是全局变量。

您好,改的是哪个文件里的wc,如果是two_stream_instability里的那确实不会有区别,因为那个文件没写电场磁场的旋度相关的方程,因为静电模式下用不到。

好的,好的,我有时间会去看看,谢谢您和zemora的建议。

我改的是 kempo1_ti.py,改了的是:
wc=20000,且iex=1,

运行效果和双流不稳定性(i32)几乎一样~~

另外,说个感觉,这个双流不稳定性,非常像中国传统文化的,阴阳鱼,也就是taichi图形的logo,阴阳鱼的形成,非常像~~~

1 个赞

1.我给出的kempo1_ti设置里vpe也就是垂直于磁场的热速度是0,因此磁场还是不会对粒子产生作用。
2.另外,在京都大学开发的KEMPO(MATLAB)的3.4版本及之前,双流不稳定性算例是以电磁模式演示的,后来为了简化计算改为了静电模式。所以电磁模式下也是可以得到类似结果的。

可以参考这个文档 Export Your Results | Taichi Docs

或者有同学整理的 taichi 导出 gif 视频的方法~ taichi导出:图片、视频、gif - 知乎

谢谢指教,另外。我意思是能不能从gui.show这个封装函数把视频缓冲区里读出(类似这两句,
pixels_img = pixels.to_numpy()
video_manager.write_frame(pixels_img)),
用在内存里直接操作(制作视频文件)而不用先写到磁盘为png等图片文件,然后再制作视频ti video。。。。

额,看了几遍,没找到楼主你说的pdf哎