视线扫过结论部分时,他整个人猛地在椅子上坐直了。
petschek在原稿里,白纸黑字地写了这样一段话。
【上述模型的成立,完全依赖于在中心扩散区存在某种局域的反常电阻机制。本文无法提供这种机制的微观物理推导。作者坦率承认,这一假设是本模型最为脆弱的部分。如果未来的工作不能在微观动力学尺度上证明这种反常电阻机制的存在,那么本模型只能被视为一种几何拓扑上的可能性,而非一种物理上的必然。】
江临把这段话反反复复读了三遍。
每读一遍,教材的公信力就在心里崩塌一层。
教材永远不会引用这一段,它只会不遗余力地展示petschek那个漂亮绝顶的x型几何,画出那四条威风凛凛的慢模激波。
但petschek自己,在模型诞生之初,就已经把刀架在了自己的脖子上。
他清楚地知道,那四条华丽的激波,那个名叫反常电阻的微观机制,是建立在一个没有被证明的沙盒之上的。
后来的几十年物理史,简直就像是针对petschek这句预言的精准打击。
超级计算机发展起来后,无数的数值模拟反复表明。
在简单的均匀电阻d条件下,petschek的x型结构根本站不住脚。
只要稍微跑长一点时间,那个漂亮的x就会坍缩,退化回Sweet-parker那条又长又无聊的面条。
如果不人为地在代码中心区域塞进去一个局域的反常电阻项,petschek的快速重联几何就像是用纸牌搭的城堡,一阵微风就能让它垮塌。
而教材把这段惨烈的历史轻描淡写地写成了,petschek模型在后来的高精度数值模拟中遇到了一些挑战。
挑战?
那更像是物理条件在拒绝这种几何的无条件成立。
petschek本人,在一九六四年就已经预见到了这一切。
他在原稿里明明白白地把自己的软肋暴露给了世界。
江临第三次起身,在东墙的空白处画了一个巨大的x。
在x的左边,写【几何成立】。
在x的右边,写【微观机制未证】。
在x的下方,写【petschek知道,petschek写了,教材没写。】
写完,他热血难凉,走出石屋。
风机二号的尾舵在低风中轻轻摆动。
看上去就像是风机在追着风转。
如果有人现在问他,风机为什么能追着风转?
如果用最简单的教材语言来回答,那就是一句话。
【风机的尾舵面积产生了空气动力学力矩,使得机舱迎风。】
这句话对吗?
对。
它完美地描述了宏观现象。
但这句话有什么用吗?
没用。
因为真正决定风机能不能追风,追得快不快,会不会在湍流中被撕裂的,根本不是那个宏观的尾舵面积。
而是底下那些极其繁琐复杂,毫不起眼的参数的耦合:阻尼弹簧的刚度,偏航轴承的摩擦系数,叶轮旋转时的陀螺推力偏置,控制系统的阵风频率响应……
只要其中一个微观参数崩溃,那个宏观上看起来天经地义的追风动作,立刻就会变成一场灾难。
江临坐在沙地上,脑子里突然闪过petschek论文里的那个x型几何。
下一刻,他冲回石屋。
站在那个大大的x的下方,添上了一句话。
【漂亮的宏观几何背后,永远站着一个没被证明的微观机制。】
读完Sweet、parker、petschek的原稿之后,江临意识到一个严重的问题。
自己一直在纸面上做评判,从来没有真正弄脏过自己的手。
理论物理的尽头,如果不能用实验去验证,就只能靠数值模拟去逼近。
教材上说petschek模型在简单的均匀电阻d模拟中退化成了Sweet-parker,这只是一个结论。
江临决定摸一遍这背后的泥。
他要亲手把petschek的模型放进代码里,看着它崩溃。
他选了数值模拟历史上非常经典的一篇论文。
biska,d.(1986).ic reconnection via current sheets.
这篇论文的核心工作,就是用直接数值模拟研究二维电阻d中电流片的演化,并给出了petschek结构无法长期维持的实锤。
江临在工作站里新建了一个项目。
【Replicate_biska_1986】
他要用c++,从零开始,手敲每一行求解器代码。
二维笛卡尔网格,均匀电阻假设,反平行磁场初始条件,空间离散采用中心差分,时间推进采用半隐式的预测-校正法。
他试图严格按照论文里的描述去设置Lundquist数、初始扰动幅度、边界条件和网格尺度。
但真正开始写代码时,他立刻发现,老论文里的文字描述,根本不足以完整复现一个数值实验。
论文里写,我们在边界上采用了自由流出条件。
但自由流出是一阶外推还是二阶外推?
有没有做特征波分解处理无反射边界?
如果是简单的零梯度外推,在非线性演化后期绝对会产生数值反射回波,污染中心扩散区。
论文里写,我们引入了一个小幅度的磁通量扰动。
但扰动函数是什么?
是高斯型还是正弦型?
截断范围怎么选?
扰动的高频分量有没有在初始化时做平滑处理?
……
这些细节,在论文里统统没有。
或者说,当年受限于篇幅,作者不可能把几万行代码的细节逐一写进文章里。
但对于复现者来说,每一个缺失的细节,都是一个致命的变量。
江临没有把这些不确定性糊弄过去了事。
他在项目的REAd文件的第一行写——
# 非严格复现。尽可能逼近。所有未明实现项,将单独列出参数空间进行扫描。
第一次正式运行,工作站硬盘灯连续亮了十七天。
第十七天下午,终端上的时间步循环终于停下。
江临颤抖着手,用python的tplotlib库调出了最后一个时间步的磁场流函数等值线图。
屏幕上出现了一条缓慢演化的电流片。
他调出时间序列动画。
在模拟的早期,也就是刚加入初始扰动的时候,系统确实在一瞬间出现过petschek那种华丽的x型几何。
那是一个令人心跳加速的瞬间。
但随着时间的推移,中心扩散区无法将堆积的物质及时排出。
那个x的夹角越来越小,重联点被慢慢拉长,变薄。
最后,在动画的结尾,那个x完全坍缩,最后退回那条熟悉得令人不舒服的长电流片。
和biska论文里的截图几乎一模一样。
很漂亮。
江临却看得脊背发凉。
太顺了。
他决定做压力测试。
回到配置文件,把Lundquist数提高一个量级。
但这一次,他做了一个违背常规的操作,故意没有同步加密网格。
他想知道,当前这套网格的承受边界到底在哪里。
系统是真的因为物理规律演化成这样,还是因为数值格式把它限制成了这样?
代码重新编译,运行。
这一次,程序跑到第八天深夜,工作站发出一声尖锐的报错蜂鸣。
ERRoR: Floating point exception (core dued).
跑崩了。
江临立刻调出崩溃前几个时间步的快照。
不是因为更高的Lundquist数更容易算。
而是因为电阻减小后,物理上的电流片扩散区厚度变得极其薄,被压到了接近甚至小于单个网格尺度的地步。
当物理尺度小于网格尺度时,代码里的数值结构就强行接管了物理演化。
江临在屏幕上看到了可怕的一幕。
x型几何刚刚尝试形成,电流片的中央就开始出现非物理的高频振荡。
红蓝相间的色块像病毒一样在网格间蔓延,几个时间步之后,整幅精密的磁场图散成了一片毫无意义的马赛克。
也就是所谓的数值爆炸。
江临感觉这片马赛克在嘲笑他。
于是关掉可视化窗口,一头扎进几万行的c++代码里。
开始排查
第一天,他怀疑是数值黏性不够,压不住高频振荡,于是引入了人工耗散项。
结果振荡没了,但重联速率被人工耗散接管,变成了假数据。
第二天,他怀疑是网格尺度太大,局部做了A。
结果交界面上出现了虚假的反射波。
第三天,他怀疑是边界反射,把流出边界改成了完全匹配层。
第四天到第二十天……
他像个在黑暗中乱撞的困兽,在这三件事之间反复来回。
每修掉一个问题,另一个更诡异的数值伪影就跳出来。
有时候是质
petschek在原稿里,白纸黑字地写了这样一段话。
【上述模型的成立,完全依赖于在中心扩散区存在某种局域的反常电阻机制。本文无法提供这种机制的微观物理推导。作者坦率承认,这一假设是本模型最为脆弱的部分。如果未来的工作不能在微观动力学尺度上证明这种反常电阻机制的存在,那么本模型只能被视为一种几何拓扑上的可能性,而非一种物理上的必然。】
江临把这段话反反复复读了三遍。
每读一遍,教材的公信力就在心里崩塌一层。
教材永远不会引用这一段,它只会不遗余力地展示petschek那个漂亮绝顶的x型几何,画出那四条威风凛凛的慢模激波。
但petschek自己,在模型诞生之初,就已经把刀架在了自己的脖子上。
他清楚地知道,那四条华丽的激波,那个名叫反常电阻的微观机制,是建立在一个没有被证明的沙盒之上的。
后来的几十年物理史,简直就像是针对petschek这句预言的精准打击。
超级计算机发展起来后,无数的数值模拟反复表明。
在简单的均匀电阻d条件下,petschek的x型结构根本站不住脚。
只要稍微跑长一点时间,那个漂亮的x就会坍缩,退化回Sweet-parker那条又长又无聊的面条。
如果不人为地在代码中心区域塞进去一个局域的反常电阻项,petschek的快速重联几何就像是用纸牌搭的城堡,一阵微风就能让它垮塌。
而教材把这段惨烈的历史轻描淡写地写成了,petschek模型在后来的高精度数值模拟中遇到了一些挑战。
挑战?
那更像是物理条件在拒绝这种几何的无条件成立。
petschek本人,在一九六四年就已经预见到了这一切。
他在原稿里明明白白地把自己的软肋暴露给了世界。
江临第三次起身,在东墙的空白处画了一个巨大的x。
在x的左边,写【几何成立】。
在x的右边,写【微观机制未证】。
在x的下方,写【petschek知道,petschek写了,教材没写。】
写完,他热血难凉,走出石屋。
风机二号的尾舵在低风中轻轻摆动。
看上去就像是风机在追着风转。
如果有人现在问他,风机为什么能追着风转?
如果用最简单的教材语言来回答,那就是一句话。
【风机的尾舵面积产生了空气动力学力矩,使得机舱迎风。】
这句话对吗?
对。
它完美地描述了宏观现象。
但这句话有什么用吗?
没用。
因为真正决定风机能不能追风,追得快不快,会不会在湍流中被撕裂的,根本不是那个宏观的尾舵面积。
而是底下那些极其繁琐复杂,毫不起眼的参数的耦合:阻尼弹簧的刚度,偏航轴承的摩擦系数,叶轮旋转时的陀螺推力偏置,控制系统的阵风频率响应……
只要其中一个微观参数崩溃,那个宏观上看起来天经地义的追风动作,立刻就会变成一场灾难。
江临坐在沙地上,脑子里突然闪过petschek论文里的那个x型几何。
下一刻,他冲回石屋。
站在那个大大的x的下方,添上了一句话。
【漂亮的宏观几何背后,永远站着一个没被证明的微观机制。】
读完Sweet、parker、petschek的原稿之后,江临意识到一个严重的问题。
自己一直在纸面上做评判,从来没有真正弄脏过自己的手。
理论物理的尽头,如果不能用实验去验证,就只能靠数值模拟去逼近。
教材上说petschek模型在简单的均匀电阻d模拟中退化成了Sweet-parker,这只是一个结论。
江临决定摸一遍这背后的泥。
他要亲手把petschek的模型放进代码里,看着它崩溃。
他选了数值模拟历史上非常经典的一篇论文。
biska,d.(1986).ic reconnection via current sheets.
这篇论文的核心工作,就是用直接数值模拟研究二维电阻d中电流片的演化,并给出了petschek结构无法长期维持的实锤。
江临在工作站里新建了一个项目。
【Replicate_biska_1986】
他要用c++,从零开始,手敲每一行求解器代码。
二维笛卡尔网格,均匀电阻假设,反平行磁场初始条件,空间离散采用中心差分,时间推进采用半隐式的预测-校正法。
他试图严格按照论文里的描述去设置Lundquist数、初始扰动幅度、边界条件和网格尺度。
但真正开始写代码时,他立刻发现,老论文里的文字描述,根本不足以完整复现一个数值实验。
论文里写,我们在边界上采用了自由流出条件。
但自由流出是一阶外推还是二阶外推?
有没有做特征波分解处理无反射边界?
如果是简单的零梯度外推,在非线性演化后期绝对会产生数值反射回波,污染中心扩散区。
论文里写,我们引入了一个小幅度的磁通量扰动。
但扰动函数是什么?
是高斯型还是正弦型?
截断范围怎么选?
扰动的高频分量有没有在初始化时做平滑处理?
……
这些细节,在论文里统统没有。
或者说,当年受限于篇幅,作者不可能把几万行代码的细节逐一写进文章里。
但对于复现者来说,每一个缺失的细节,都是一个致命的变量。
江临没有把这些不确定性糊弄过去了事。
他在项目的REAd文件的第一行写——
# 非严格复现。尽可能逼近。所有未明实现项,将单独列出参数空间进行扫描。
第一次正式运行,工作站硬盘灯连续亮了十七天。
第十七天下午,终端上的时间步循环终于停下。
江临颤抖着手,用python的tplotlib库调出了最后一个时间步的磁场流函数等值线图。
屏幕上出现了一条缓慢演化的电流片。
他调出时间序列动画。
在模拟的早期,也就是刚加入初始扰动的时候,系统确实在一瞬间出现过petschek那种华丽的x型几何。
那是一个令人心跳加速的瞬间。
但随着时间的推移,中心扩散区无法将堆积的物质及时排出。
那个x的夹角越来越小,重联点被慢慢拉长,变薄。
最后,在动画的结尾,那个x完全坍缩,最后退回那条熟悉得令人不舒服的长电流片。
和biska论文里的截图几乎一模一样。
很漂亮。
江临却看得脊背发凉。
太顺了。
他决定做压力测试。
回到配置文件,把Lundquist数提高一个量级。
但这一次,他做了一个违背常规的操作,故意没有同步加密网格。
他想知道,当前这套网格的承受边界到底在哪里。
系统是真的因为物理规律演化成这样,还是因为数值格式把它限制成了这样?
代码重新编译,运行。
这一次,程序跑到第八天深夜,工作站发出一声尖锐的报错蜂鸣。
ERRoR: Floating point exception (core dued).
跑崩了。
江临立刻调出崩溃前几个时间步的快照。
不是因为更高的Lundquist数更容易算。
而是因为电阻减小后,物理上的电流片扩散区厚度变得极其薄,被压到了接近甚至小于单个网格尺度的地步。
当物理尺度小于网格尺度时,代码里的数值结构就强行接管了物理演化。
江临在屏幕上看到了可怕的一幕。
x型几何刚刚尝试形成,电流片的中央就开始出现非物理的高频振荡。
红蓝相间的色块像病毒一样在网格间蔓延,几个时间步之后,整幅精密的磁场图散成了一片毫无意义的马赛克。
也就是所谓的数值爆炸。
江临感觉这片马赛克在嘲笑他。
于是关掉可视化窗口,一头扎进几万行的c++代码里。
开始排查
第一天,他怀疑是数值黏性不够,压不住高频振荡,于是引入了人工耗散项。
结果振荡没了,但重联速率被人工耗散接管,变成了假数据。
第二天,他怀疑是网格尺度太大,局部做了A。
结果交界面上出现了虚假的反射波。
第三天,他怀疑是边界反射,把流出边界改成了完全匹配层。
第四天到第二十天……
他像个在黑暗中乱撞的困兽,在这三件事之间反复来回。
每修掉一个问题,另一个更诡异的数值伪影就跳出来。
有时候是质