离开计算器可以走多远: x + ln x =0 2015 年秋季研究生课程《高等应用数学》期末测验的一道附加题: 求方程 x + ln x = 0 的实根,要求给出分析和计算过程,结果的有效数字越多越好,且不能使用计算器。 在上一年度的期末测验中也有一道类似的估算题,难度稍微要大一些。168 名研究生参加测验,几尽无视附加题,也或许是瞄一眼整个人都不好了,直接放弃。仅有十余名学生动笔写了点东西,而仅有一名学生稍稍走得远一些,但也还是差临门一脚。此次稍稍降低了难度,在确保思想性和技巧性不差的情况下适当地降低了所需手工计算的难度。142 名研究生和2 名本科生参加测验,虽然学生数少了,但约半数学生动了笔,不过都走得不是太远,相比较离目标最近的还是去年的那名学生。考虑到考试时间的限制,而且没有任何提示,能走出一两步也已经很不容易了,但即使时间不限,不提示,不讨论,不会的也还是很可能不会。 本年度首次有本科生选修该研究生课程,自然对这两名本科生的表现较为期待,就谈谈他和她如何处理这道附加题。他先用作图法找到根所在的区间 (0,1) ,然后采用牛顿迭代格式 x k +1 = x k - ( x k + ln x k )/(1 + 1/ x k ), k = 0, 1, 2, ... 给出 x 0 = 1, x 1 = 0.5, x 2 = ? ,然而由于不能使用计算器,他没能走得更远。她给的第一种思路也是牛顿迭代法,并给出了 x 0 = 0.6 和 x 1 = 0.55 ,但未详细说明 0.55 是如何手算出来的。她显然意识到在这个思路下没有计算器是走不了更远的,而且牛顿迭代法是本科掌握的方法,在大学数学复习课上有提到,但不是该课程考察的知识。进而她在第二种思路中试图构造一个小参数 0.1 以采用级数法,然而经过一些尝试发现未能找到合适的构造就放弃了。事实上,在后一种思路中,她已经跨出了大大的一步,若稍加提示,必将天堑变通途。 作为附加题,题目可以较为开放,任何方法都可以尝试,目的也是希望学生能够灵活运用课程中教授的方法。如果不限制使用计算器,牛顿迭代法无疑是可选的最佳方案,迭代效率高。通常的数学软件都采用了牛顿迭代法,如Mathematica 软件的FindRoot 函数或Matlab 软件的fzero 函数。这个时代的人啊,敲几下键盘立马就可以获得精度极高的结果,例如用Mathematica 软件输入FindRoot , { x ,1}] ,可得{ x -0.567143} 。然而,重要的并不是这个附加题的结果。 该课程以林家翘先生和西格尔的《自然科学中确定性问题的应用数学》第一、二卷为主要教材,旨在通过案例讲授“应用数学过程”的思想,贯穿始终介绍了各种各样的摄动方法:级数法、参数微商法、迭代法、庞加莱方法、分部积分法、拉普拉斯方法、傅里叶分析、小波分析、匹配法、WKB 方法、多重尺度法、均匀化方法。解答该附加题也自然可以只限定在这些方法中尝试,容易猜测级数法和迭代法等正则摄动法最为可能可行。 摄动法的思想是将精确解作微小的扰动变成近似解,使得它满足一些可解的方程用以代替原来很难或不能精确求解的方程。因此,对解有初步的了解是很重要的。若把 x + ln x 看作函数,可以采用零值定理找出解所在的大概区间,如 0 x 1 ,当然也可以采用作图法。可是,居然还有更直接的方法,考虑到 x 和 ln x 必须一正一负,立即可得 0 x 1 。这么直接而清晰的一个方法是一名研究生写下的,虽然他写下这点后就没有再走多远。由解所在的区间说明可以将 x 视为小量,那么借助泰勒展开令人讨厌的 ln x 似乎就可以不那么讨厌了。倘若将 ln x 在 x = 0 附近作泰勒展开,即麦克劳林展开,那是行不通的。换个思路,1 - x 也是小量,由泰勒展开(见附A ),可得: ln x = ln ~ - (1 - x ) - (1 - x ) 2 /2 - (1 - x ) 3 /3 - (1 - x ) 4 /4 - (1 - x ) 5 /5 - ..., 至此,可以得到一个能够手工计算的迭代式子,如 x k +1 = 1/2+1/2* , k = 0, 1, 2, ..., 但算起来很辛苦。例如,若保留方括号中两项进行迭代,有: x 0 = 0.5, x 1 = 0.5833, x 2 = 0.5555, x 3 = 0.5640, x 4 = 0.5613, x 5 = 0.5622, ... 辛辛苦苦算来算去也仅能精确到小数点后两位。为提高精度必须保留更高阶的部分,计算难度越发大了去。如果把方程改写为 x *exp( x ) = 1 或 exp( - x ) = x ,将左端作麦克劳林展开再构造迭代式,这样的思路也可行,但也还是会算得很辛苦。至此说明可以通过构造特殊的迭代格式进行计算,与牛顿迭代法相比,效率低一些,但至少是能够手工计算的。实际上,“有效数字越多越好”这个要求的信息量很大,暗示不能蛮干。嗯,还是可以再走得更远一些的。 迭代法通过不断迭代来提高结果的精度,但每次迭代所需的计算量越来越大,而且还有一个严重的缺点,每次迭代都要重复计算出已经精确到的部分。这个缺点是有办法克服的,比方改用级数法,这在课程教学中已经说明过。然而,整个教科书几乎都是围绕着包含小参数的问题在讨论,该附加题所给的方程可不包含小量,那似乎就没辙了。嗯,得先构造一个小量。迭代法可以不用明确小参数,但是可以帮助构造小参数。如果先将方程改写为 x = exp( - x ) ,再运用麦克劳林展开,可得: x = 1 - x + x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... 然后基于迭代法思想将方程进一步改写为: 2x -1 = x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... 等式左边的每一项至少是 O (x) ,而右边的每一项都是 o ( x ) 。若略去右边的部分可得 x = 1/2 ,可将其取为小量,记为 a = 1/2 。方程可以进一步改写为: 2 x -2 a = x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... 进而,可以采用级数法求解得到: x = a + 1/4* a 2 + 1/24* a 3 - 1/192* a 4 - 13/1920* a 5 +... 考虑到 a = 1/2 ,上式的每一项都可以较轻松地通过手工计算完成,各项分别为0.5, 0.0625, 0.0052083, -0.0003255, -0.0002115 等等。因此,若保留前五项可得结果为 0.567171 ,精确到了小数点后四位;项数取得越多精度越高,但收敛速度越来越慢。实际上,也可取小量 c = 1/10 = 0.1 ,正如那名本科生所期待的。进而将方程改写为: 2 x -10 c = x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... 其级数解中每一项的数量级将非常清晰,见附B 。然而,对于这个题嘛,其实际的计算量没有太多差别,因为收敛速度是一样的。 还能走多远呢?路漫漫其修远兮 ~~~ 郑志军 2016 年1 月7 日 附 A :麦克劳林展开 ln(1 + u ) = u - u 2 /2 + u 3 /3 - u 4 /4 + u 5 /5 - u 6 /6 + ... = Sum exp( u ) = 1 + u + u 2 /2! + u 3 /3! + u 4 /4! + u 5 /5! + u 6 /6! + ... = Sum 附 B :解答小结 1. 方程 x + ln x = 0 的根满足 0 x 1 ,说明 x 可以作为小量 2. 由麦克劳林展开有: x = exp( - x ) = 1 - x + x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... 3. 由迭代法思想有:2 x -1 = x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... ,右边的每一项都是 o ( x ) 4. 取小量 c = 1/10 ,将方程进一步改写为:2 x -10 c = x 2 /2! - x 3 /3! + x 4 /4! - x 5 /5! + ... 5. 由级数法可得: x = 5 c + 25/4* c 2 + 125/24* c 3 - 625/192* c 4 - 8125/384* c 5 - 146875/4608* c 6 + 1140625/64512* c 7 + 191171875/1032192* c 8 + ... 6. 手工计算可得: x = 5 c + 6.25 c 2 + 5.2083 c 3 - 3.255 c 4 - 21.15 c 5 - 31.8 c 6 + 17. c 7 + 185. c 8 +... = 0.567143...
第一颗原子弹爆炸当量的估算 费米对世界上第一颗原子弹的当量的估算,被传为美谈。 有人质疑费米是否作出过这种估算,认为是以讹传讹,下面给出比较可信的证据: http://www.nuclearfiles.org/menu/key-issues/nuclear-weapons/history/pre-cold-war/manhattan-project/trinity/index.htm http://www.nuclearfiles.org/menu/key-issues/nuclear-weapons/history/pre-cold-war/manhattan-project/trinity/eyewitness-enrico-fermi_1945-07-16.htm ” Trinity Test, July 16, 1945 Eyewitness Report by Enrico Fermi Observations During the Explosion at Trinity on July 16, 1945 On the morning of the 16th of July, I was stationed at the Base Camp at Trinity in a position about ten miles from the site of the explosion. The explosion took place at about 5:30 A.M. I had my face protected by a large board in which a piece of dark welding glass had been inserted. My first impression of the explosion was the very intense flash of light, and a sensation of heat on the parts of my body that were exposed. Although I did not look directly towards the object, I had the impression that suddenly the countryside became brighter than in full daylight. I subsequently looked in the direction of the explosion through the dark glass and could see something that looked like a conglomeration of flames that promptly started rising. After a few seconds the rising flames lost their brightness and appeared as a huge pillar of smoke with an expanded head like a gigantic mushroom that rose rapidly beyond the clouds probably to a height of 30,000 feet. After reaching its full height, the smoke stayed stationary for a while before the wind started dissipating it. About 40 seconds after the explosion the air blast reached me. I tried to estimate its strength by dropping from about six feet small pieces of paper before, during, and after the passage of the blast wave. Since, at the time, there was no wind I could observe very distinctly and actually measure the displacement of the pieces of paper that were in the process of falling while the blast was passing. The shift was about 2 1/2 meters, which, at the time, I estimated to correspond to the blast that would be produced by ten thousand tons of T.N.T. Source: U.S. National Archives, Record Group 227, OSRD-S1 Committee, Box 82 folder 6, Trinity. Transcription: Thank you Gene Dannen for transcribing this document. “ 下面是译言网的中文翻译,我借用一下: http://article.yeeyan.org/view/mjysci/121826 ” 译言网的中文翻译如下: “ 三位一体测试目击者报告 1945年7月16日 恩里科·费米 对位于三位一体试验场爆炸的观测 在7月16日的早晨,我被派驻在距爆炸现场10 英里的三位一体基地。 爆炸发生在上午大约5时 30分,我的面部由一块中间焊了黑玻璃的大板所保护。我的第一印象是爆炸形成了非常强烈的闪光,身体的暴露部分感到一股热浪袭来。虽然我并未直视爆炸点,但当时的印象是此处突然变得比白天还亮。接着我通过黑玻璃看向爆炸中心,可以看到一些火焰开始聚集,并立刻上升。几秒钟后上升的火焰亮度降低,烟的顶部扩大像一个巨大的蘑菇,很快就超越大概是3万英尺高空的云层。在达到其高度的极限时,烟柱维持了一段时间才开始被风吹散。 在爆炸发生后约40秒冲击波才到达我这里。我尝试通过这个方法估计其爆炸强度:在爆炸开始前,爆炸中,冲击波通过后分别在大约六英尺高撒落小纸片。因为在爆炸前,我清晰地观察到确实没有风造成飘落的纸片偏移。而冲击波到来时,纸片偏移了大约2.5米,当时,我就估计爆炸产生的TNT当量是10000吨。 资料来源:美国国家档案馆,记录组227,OSRD-S1委员会,82箱第6文件夹,“三位一体”。 抄本:感谢 Gene Dannen 抄写本文件。 “ 我不知道费米是怎么估算这一当量 - 10000吨 TNT - 的,所以我自己来估算了一下,方法如下,数据全部仅仅来自费米的原始报告: 1、费米距离爆心10英里,大约为16公里,则在此半径范围内,一千米高度的空气体积接近800立方千米,空气密度为1.29,则空气总质量约为10^12千克,即10亿吨。 2、纸片初始高度为2米左右,落地需要1秒钟,水平偏移了2.5米,即纸片水平速度约为2.5米每秒,将此速度粗略作为冲击波到达此处时候的速度。 3、假设原子弹爆炸时将这一范围内的空气全部加热为这一速度,则所需能量从动能公式计算出约为,3×10^12焦耳。 4、一千克TNT 的热值约为4*10^6焦耳,则上述的空气总能量约为 0.75 千吨的TNT 的能量,取整为一千吨TNT 当量。 5、考虑到前面计算空气体积的时候,是计算了一个厚度为一千米的圆饼而不是一个半径为16千米的半球, 如果按照半球计算,则体积会增加一个数量级,费米描述中就有(3万英尺)也是接近一万米的高度,但是空气密度随高度而减小,所以折中一下,取一个因子 5, 再考虑到冲击波只是原子弹爆炸释放能量的形式的一种,再取一个因子2。 6、这样,就正好得到了 一万吨 TNT 的 总当量。 很有意思的是, 我完成了上述估算以后,检索有没有其它的估算方法时,发现科学松鼠会上的一篇文章: http://songshuhui.net/archives/40461 费米问题,理科生的脑筋急转弯 Comments 候戏 发表于 2010-07-25 09:39 作者的方案是: ” 1)假设纸片做自由落体运动,初速度相当于气浪的速度,这个计算在初中物理习题中常见。 2)假设原子弹爆炸能量全部转化为空气的动能,爆炸之后,气浪形成球面向各个方向扩展,扩展到费米所在地时,球体内总的空气质量可以通过空气密度乘以球体体积算出来。 3)总能量等于气浪速度平方乘上空气质量。然后转换为 TNT当量单位 ,完成。 “ 思路和方法跟我完全一样。 难道就这么一条路? 我反而郁闷了。 ^_^。。。。。。。。。。。。。。 顺便问问,编辑 MM 不给本文一朵 小红花吗 ?
上次提到: 1945 年 7 月 16 日第一颗原子弹在新墨西哥州的沙漠引爆成功 40 秒后,试爆现场附近的科学家感受到冲击波的气浪。费米伸手向空中撒了一把碎纸片,落在他身后 2.5 米的地方。根据这个“实验数据”费米估算出的爆炸威力相当于一万吨 TNT ,而后来仪器测到的数值大约是两万吨。 他是如何估算的?我们以后再说。 我们先来看另外一个更精确估算的例子: 1947 年英国的物理学家 G. I. Taylor 根据美国公布的一组爆炸火球的照片(见下图),更准确地估算出爆炸当量是 1.7 万吨 TNT 。 他又是如何估算的呢? 我们不妨看一下我们知道什么:从照片里我们可以得到不同时间( t )火球的半径( R );而我们想知道爆炸的能量( E )。怎么才能用已知的物理量 R 和 t 表示出能量 E 呢? 能量的量纲是: 2 = 2 -2 。用 R 和 t 容易得到速度的量纲 ,但是显然我们缺少含质量量纲的物理量。 我们知道爆炸火球是在大气中膨胀的,其膨胀速度一定与大气的性质有关。因为是热膨胀,那么有关的物理量就是大气的温度、比热、和质量密度。所以我们选择质量密度 r , = -3 。 这样我们可以用这四个物理量组成一个无量纲量 P =E/ r R 5 t -2 即: E= P r R 5 t -2 ,或者写成 R 5/2 =(E/ P r ) 1/2 t 。 Taylor 就 是从上面这组美国公布的爆炸火球的照片,得到火球半径随时间变化的曲线(见下图): lnR=(2/5)lnt+(1/5)ln(E/ P r ) 。 这里截距 (1/5)ln(E/ P r ) 可以从图上直接得到!取 P ~O(1) ,则可以直接计算能量当量 。 在国际单位制下,代入 r =1.25kg/m 3 ,假设 无量纲量 P ~O(1) ,就估算出爆炸能量 E~8x10 13 J~2 万吨 TNT ( 1 吨 TNT~4x10 9 J )。 【文章中两张图都来自:】 G. I. Taylor, “The Formation of a Blast Wave by a Very Intense Explosion. II”, Proc. Roy Soc. London A201 , 175 (1950).