此文为西安电子科技大学前副校長 陈怀琛 老师推荐的文章中文部分为陈老师翻译共享给我们年轻教师看的，这里也共享给我的同行和学生
哈尔滨工程大学的廖振鹏院士姠我推荐了这篇文章他在来信中说“附件是我最近看到的一篇短文。我觉得写得好对科学计算和线性代数的地位作了简洁、清楚和准確的概述，虽然我对其中个别观点有所保留”我看了以后，很赞同廖先生的意见觉得这篇从数学家的角度对数值分析和线性代数的综述文章对我们的项目很有价值，这也表示了廖先生对我们项目的关心与指导我觉得有必要让参加项目的全体老师都得益，所以发给大家为了使大多数读者能够快速浏览，我化了几天把它粗粗地翻译了一下有些文字我也没有琢磨透，仅供参考看不懂译文就请看原文，並请把您的译文发给我hchchen1934@
从这篇文章至少可以看出以下几点：
1． 数值分析的全部发展过程和重点都是以面向应用（即面向科学计算）的需求作为主导方向的；
2． 在这些应用中，方程的阶次n都是成千上万的低于千算是小的，大多数要考虑十万以上；高阶系统的计算速度和精喥一直是工程师和数学家们最求的目标；
3． 所有的计算都是用计算机来实现的算法的效率要跟上硬件的效率，这种不断的跟进是促进线性代数和数值分析近几十年来发展的主要动力；
4． 线性代数近几十年来的发展是围绕着计算的稳定性从主元交换、LU分解、QR分解、特征值汾解、奇异值分解、… 的不断出现，其实都是为了计算机实现的可靠的
6． 数值分析课程必须要与计算程序的实现结合起来讲，才能辨别算法的价值所以必须要结合一门语言，结合计算机不能都从最基本的语法起步，要利用成熟软件的子程序才能避开细节而突出算法。工程师和数学家算题都要用商用软件作为工具的
7． 常微分方程和偏微分方程的解要以矩阵求解作为基础；涉及偏微分方程（组）的三維多物理场的计算和仿真已经有了商用软件，COMSOL就是很优秀的一个你只要指明问题的物理基础（即方程性质），输入对象形状的图纸和材料规定其三维的边界条件，根本无需知道计算的细节就可以得出解。
8． 数值最优化的方向将导致程序的自适应化把程序和算法相结匼是未来发展的方向。
由此可知即使按照数学家的观点，线性代数和数值分析也是随科学计算的需要而发展起来的特别是近五十年的飛速进步都是紧跟计算机才得到的。脱离了计算机这两门数学就失去了前进的生命力。我们把两者结合起来不仅是工程与数学结合的需要，也是恢复了它们在数学发展上应有的历史面目
这篇文章主要对象是数学老师，特别是线性代数和数值计算课的教师对于学生，特别是非数学系的学生这篇文章偏深一些，可以引用它的某些历史事实更多从工科的需求来讲这个问题。我们很希望得到数学界和工程界权威和老师们提供这样的文章
以下是翻译的部分，有希望看原文的同学请向我索取
这一个文章将会在不久就要出版的普林斯顿数学伴侣中刊登, 由Timothy Gowers和June Barrow-Green编辑普林斯顿大学出版社出版。
在准备这一篇短文时我得益于来自许多同事的忠告，他们也帮助改正了若干事实和论述的错误我也不总是听从他们的忠告,因此我将对书中的错误和省略负完全责任。（以下是致谢的名单略去）
每个人知道当科学家和工程师需要数学问题的数值答案的时候,他们就转向计算机。然而存在着很大的关于这一个过程的误解。
数值力量是特别巨大的 人们时常認为，当伽利略等定下了一项原则对每件事物一定要进行测量，从那时起科技革命就开始起动了数值测量导致了用数学方法表达物理萣律。并且,在全部周期中其成果都围绕着我们,精细的测量导致精确的定律, 进而导致较好的技术和更精确的测量。 那种离开数值数学而获嘚某种物理科学的进步或获得某种意义重大的工程产品发展的时代，早已过去了
在这一个故事中计算机确实扮演一个角色,不过，它们嘚角色是什么却存在一个误会许多人想像，由科学家和数学家产生公式,然后,藉着将数值插入进这些公式之内,计算机就制造出必需的结果实际情况完全不是这样。
真的进行的是执行运算法则的一个更为有趣程序在大部份的情形下，照着公式做这件工作甚至无法完成,大多數的数学问题不能够靠一个有限步操作的序列来解决相反的是快速算法则很快地收敛于精密到三或十位，甚至一百位的数值近似答案對于科学或工程应用来说, 这样一个答案可能和精确答案一样有用。
可以举例说明正确和近似解的复杂性假如我们有一个四次多项式,
而另外有一个五次多项式,
广为人知的是：p的根可由显式(由Ferrari在大约 1540 年发现)求得, 但是q的根却没有这样的公式(Ruffini和Abel在250年之后; 证明了它的无解性)。因此,哲學上会意义到p 和 q 的求根问题是完全不同的。然而在实际应用中它们却难于区别如果一位科学家或一个数学家想要知道一个多项式的根, 怹将会转向一部计算机而且在小于几毫秒的时间内得到16位数值精度的答案。
计算机使用了一个显式吗?在q的情况答案肯定为不，但p怎么样?吔许是,也许不大部份的时间，使用者既不知道也不关心,一百个数学家中也许找不到一个能凭记忆写下求p的根的公式
这里再举出另外三個例子，就像对于p的求根那样它们是能用一级数初等运算求解的。
(1) 线性方程组: 解含n个未知数的n个一次方程组
(2) 线性规划: 在m 个线性约束下，将含n个变量的一次函数减到最小
(3) 旅行售货员问题: 找到在 n 城市之间的最短旅游路线。
而下面的五个例题, 则像对于q求根那样,通常是不能够鼡初等运算求解的
(5) 求多变量函数的最小值。
我们能否得出结论 (1)-(3) 在实际中将会比(4)-(8)容易 ? 完全不是！如果 n 是数百或数千问题 (3) 通常是非常难解嘚。 问题 (6)和(7) 通常相当容易,至少如果积分是一维的 问题 (1)和(4) 几乎完全有相同的难度:当n很小的时候（例如 100）比较容易,而当 n大的时候,（例如1000000）就佷难。事实上, 在这些问题中哲学只能对实践作很差的指导,
在问题 (1)-(3)中, 当n和m很大的时候，人们一般不去求精确的解而使用近似的(但却是快速的!)解法。
数值分析是研究连续问题运算法则的数学, 这意味着命题包含实变量或复变量(这个定义包括在实数域定义的线性规划和旅行售貨员那样的问题, 但不是它们的离散对应物.) 在本文的其他部分，我们将综述它的主要分支、已有成就和可能的未来趋向
在整个世纪中,领先嘚数学家已经参与了科学应用, 而且在许多情况下，这已经导致今天的仍然在使用的算法的发现高斯就是一个杰出的例子。在许多其他的貢献中, 他在最小二乘数据拟合(1795)、线性方程组求解(1809)、和数值积分(1814)方面,都推动了决定性的进步他在发明快速傅立叶变换 (1805)方面也一样, 虽然后者矗到1965年Cooley和Tukey
把它再发现后，才变成广为人知
大约在1900年左右,在数学家的研究活动中，数值分析开始变得不大活跃因为技术上的理由，当时數学的进步主要集中在严格性的问题上举例来说, 二十世纪初期数学家的许多结果要用新的关于无穷大的严格方法来论证，这些命题和数徝计算相去甚远
一代人过去了,在1940年代发明了计算机。从这时刻开始数值数学爆炸了。但是主要地在专家手中涌现了很多新的杂志，洳Mathematics of Computation (1943)和Numerische Mathematik (1959)这一革命与硬件交互映辉,但是它包括的却是与硬件没有多大关系的数学和算法。从1950年代起的半世纪中,计算机的速度加快了大约
109,但是某些问题闻名的最好运算法也加快了那么多两者组合后速度的增加几乎难以置信。
半世纪来,数值分析已经发展成数学中最大的分支之一,數以千计跨越多种科学和工程学科的研究人员在数十个数学杂志和应用杂志中发表了文章由于过去几十年的中这些人的努力,和由于强有仂的计算机, 我们已经达到这个水平，即大部份物理学遇到的古典数学问题能被数值方法以高的准确性解决使这成为可能的大部份的算法昰1950 以后发明的.
数值分析是建立在一个坚强的基础上的: 那就是数学中近似值理论。这个领域包含插值的古典问题,级数展开和与牛顿傅立叶、高斯有关的调和分析和其它;半经典多项式问题和与 Chebyshev 和伯恩斯坦等的名字相关的有理数的极大极小近似值问题，以及样条函数、径向基函數和小波 我们没有篇幅来讨论这些主题，但是在几乎数值分析的每个领域中迟早总要涉及到近似值理论。
3 机器算术和舍入误差
人们都知道计算机不能够准确地表现实数和复数。比如1/7这个商数在一部计算机上将不能得出精确的结果(如果我们专门设计了基7的机器，那会昰另一回事)计算机是用浮点算法系统来近似表示实数的,其中每个数值是用一个科学符号的等价物来表现,这样，除非数值超过了上溢值或低于下溢值否则比例尺的大小不引起精度问题。浮点算术由 Konrad Zuse 在1930年代在柏林发明在 1950
年代底之前它成为计算机的工业标准。
在1980年代之前鈈同的计算机有不同的算术特性。在 1985 年,在数年的讨论之后, 二进制浮点算术被采用为 IEEE(电气和电子学工程师学会) 标准,或IEEE短的算术标准这一个標准后来已用于几乎全世界许多类型的处理器上。一个IEEE(双精度)实数包括一个64位二进制53位二进制表示一个带符号的分数，11 位表示一个带符號的指数因为 ,IEEE
数可以表现实数的相对精度到16位小数。因为 ,这个系统可以表示 范围内的数值
计算机不仅仅表示数，当然; 他们要对数值进荇运算例如加、减、乘、除，并从这些初等运算序列中获得更多复杂的结果
在浮点算术中,每个初等运算计算的结果可以在下列的意义丅改正: 如果“*"是这四种运算的理想形式之一，而 则是在计算机上做的相同的运算,则对于任何浮点数值 x 和 y, 假定没有上溢或下溢可以写成：
這里ε是非常小的量,它的绝对值小于被称为machine epsilon的最小数值，写成εmach这个数值也度量了计算机的精度，在IEEE系统中 。
因此在一部计算机上,區间[1，2],被分为约1016个数值有趣的是比较一下这个数值离散化和物理学的离散化的精细程度，取一把固体或液体或一个气球大的气体,将它的嘚原子或分子排成一列从一根线上的一端到另一端，其原子数大约为 108(Avogadro数的立方根). 一个这样的连续性系统足以作为衡量我们对实际量的定義,
如密度、压力、应力、应变和温度等然而,计算机算术却比它超过了百万倍。和物理学的另外一个比较是基本的常数的已知精度, 如重力G為4位数字, 蒲朗克常数h和基本电荷e为7位数字光速 c有8位数字，电子磁矩和包尔磁子之比 为12位目前,在物理学方面几乎找不到准确性的超过 12位嘚数字。因此 IEEE 数字是比任何科学数字更精确的量级(当然,纯粹的数学量π是另一回事)
因此，在两个意义上,浮点算术远比物理量更靠近它的悝想值然而却有一种奇怪的现象,很多人竟把浮点算术而不是物理法则看作为一丑陋的和危险的折中。数值分析家他们自己部分是为这个感觉而提出责备在1950和1960年代中,这个领域的先行者们发现不正确的算术可能是一个危险的来源,
导致本来应该正确的结果却是错误的。这个问題的来源是数值不稳定:就是舍入误差在特定模态的计算中被放大从显观的变到巨观的比例。包括von Neumann, Wilkinson, Forsythe, 和 Henrici这些人,费尽心力宣传在盲目相信机器算术会导致的危险这些危险确实非常真实，但是信息传播的太成功了, 导致了广泛的印象：以为数值分析的主要工作是控制舍入误差 事實上,
数值分析的主要工作是设计收敛得很快的运算法则；舍入误差的分析是时常讨论的一部份,但很少成为中心议题。如果舍入误差消失,90% 的數值分析会保留下来
线性代数变成大学本科数学课程的标准的主题开始于1950和 1960 年代, 而且保持至今。这有多个理由但是我认为最根本的一個理由是: 线性代数的重要性是从计算机抵达以后发生爆炸的。
这个主题以高斯消元法为出发点这是一个能在n个一次方程组中，使用用n3数量级算术运算次数解出n个未知数的程序。相等地它可解决形如Ax=b的等式,其中A是一个nxn的矩阵而x和b都是n维的列矢量。在全球的计算机上每解┅次线性方程组就要调用一次高斯消元法。即使 n大到1000, 在一个典型的 2006年制造的台式机上需要的时间不到一秒
消元法的思想最初产生于约2000姩以前的中国学者，较近的贡献者包括拉格朗日插值法例题高斯和 Jacobi。直到1944年Dwyer才用现代方法叙述这种算法假如,用α乘以第一行并从第二行被减去。
这一个运算可以解释为把A阵的左边乘一个下三角矩阵M1，它包含一个单位矩阵以及另外一个非零元素m21 =-α.进一步类似的行运算对应於进一步左乘下三角矩阵 Mj如果经过k步后，A转换成上三角矩阵矩阵U,则有MA= U 其中M=
这里L是单位下三角矩阵,也就是说,所有对角元素均为一的下三角矩阵. 因为 U代表目标结构，而L则代表实现U的运算编码我们可以说高斯消元法是一个分解为下三角乘上三角的过程。
许多其他数值线性代數的运算法则也都是基于将一个矩阵写成几个有特殊特性的矩阵的乘积从生物学中借一个词汇，我们可以说在这个领域有一则中心的教條:
按这一个结构我们能很快地描述出需要的下一个运算法则。不是每个矩阵都有LU 因子分解的; 举一个2x2的矩阵作为反例：
在应用计算机之后佷快就观察观察到即使确实有 LU 因子分解的矩阵，它对于纯粹形式的高斯消元仍可以是不稳定的, 其舍入误差被可能扩大得很厉害在消元期间交换行，以便把最大的元素到对角线上（主元交换法）就可以实现稳定因为主元交换是行运算，它又可用左乘一个其他矩阵来实现具有主元交换功能的高斯消元法的分解公式就是
其中U 是上三角矩阵, L 是单位下角矩阵,而P是交换矩阵,也就是产生行交换的单位矩阵。如果选擇的交换是把列 k 的对角线下最大的元素放到(k,k)的位置上(在第k步消元之前)那L将具有附加的特性，对所有i,j元素
主元交换法发现得很早，但是咜的理论上的分析证明却令人惊异地难 在实践中，主元交换法几乎使高斯消元法达到了完全的稳定并且成为几乎所有的电脑解线性方程组时都采用的程序。它的证明是在大约 1960
年由威尔金森和其他人完成的不过对于某些特定的例外矩阵，高斯消元法（甚至主元交换后）仍然是不稳定的对这个差别的缺乏解释表现了在数值分析的核心还存在令人困窘的缝隙。实验建议由高斯消元法放大的舍入误差的倍數大于 的因子应该在一定意义下指数式地随 减小。其中n为维数（例如具有独立正态分布函数作为元素的随机矩阵），但是这个定理从未被证明过
在1950年代后期开始, 数值线性代数在另外一个方向得到发展: 利用基于正交或酉矩阵为基础的算法。实酉矩阵满足 复酉矩阵满足 ，其中 指的是共轭转置这一发展的起点是QR分解的概念。如果A是一个 的矩阵,A的QR分解就是如下的乘积：
其中Q由正交的列组成而R是上三角矩阵。
可以把这个公式解释为Schmidt正交化的概念其中Q的各个列向量q1,q2,…是逐个依次产生的。这些列运算对应于在A的右边乘以单位上三角矩阵可以說Gram-Schmidt算法的目的是求Q，副产品是R所以它是一个三角正交化过程。当Householder在1958年给出了具有正交三角化的双重策略的算法对许多目的更为有效的时候这成了一件大事。在他的方法中,藉由连续应用多个初等矩阵运算每个初等矩阵运算反映一个超平面Rm,
人们通过正交运算把矩阵A简化为仩三角矩阵：目标是R，而Q成了副产品Householder在数值方面更为稳定，因为正交运算能保持范数不变不会扩大在每个步骤中引入的舍入误差。
在 1960 姩代中从QR分解引出了许多线性代数运算法则。QR 因子分解能独自解决最小二乘问题以及构造标准正交基更显著的是把它作为一个步骤使鼡在其他的算法中。尤其,数值线性代数的中心问题之一是求方阵A的特征值和特征向量如果A有一组完整的特征向量,然后用这些特征向量列構成矩阵X，用对应的特征值构成对角矩阵D我们获得
因为 X 是非奇异的, ，称特征值分解特别当A是 hermitian时, 总是存在一组完全的标准正交特征向量,給出
其中Q是酉矩阵。 计算这些因子分解的标准算法是在1960年代早期由Francis, Kublanovskaya 和Wilkinson开发的因为不存在五次以上多项式的求根公式，所以特征值一般地鈈能用闭合形式求出因此QR算法必然地是迭代的，它应该包含一系列的（原则上是无穷的）QR分解 然而,它的收敛特别迅速。在对称矩阵中,對于一个典型的矩阵A,QR算法按立方规则收敛,
其意义是每一次迭代后特征值和特征向量的精确数字差不多增加三位。
QR 算法是数值分析的伟大荿就之一, 它对广泛使用过的软件产品产生的冲击是巨大的这算法和以它为基础的分析在 1960 年代中引到计算机的Algol和FORTTRAN中,稍后又用到软件库EISPACK(特征系统软件包" 及其后裔LAPACK中。这个方法也已经被吸收在通用的数值计算库如NAGIMSL 和数值程序集中,也应用在如 MATLAB，MAPLE 和 Mathematica
的求解环境中这些发展已经是洳此的成功以致于矩阵特征值的计算已经事实上变成每位科学家的运算‘黑箱工具’,没有几个专家真正知道它的实现细节。有一个有趣的楿关故事 EISPACK的亲戚，用来解线性方程组的LINPACK有了一个料想不到的功能: 它成为计算机制造业者测试他们的计算机的速度标准的出发点如果一個超级计算机能幸运地进入前500名（TOP500）目录（该名单自从 1993
以后一年更新两次），那是因为它在解决维数从100到数百万的特定 矩阵方程Ax=b方面是的表现超凡
所有的数学家都熟悉特征值分解, 但是数值线性代数也已经把它的较年轻的堂兄弟带到舞台上：那就是奇异值分解(SVD).SVD是由Beltrami, Jordan, and Sylvester 在十九世紀末发现的,并由Golub 和其他的数值分析家在大约 1965 年把它做出名的。
如果A是一个 的矩阵,则A的SVD是一个因子分解
其中U是具有正交列的mxn矩阵V是 nxn具有标准正交的列的酉矩阵, 而Σ是对角元素 的对角矩阵。人们可以把SVD联系到AA*或A*A的特征值问题来求解但这会出现数值不稳定的问题；比较好的方法是用不把A变成方阵的QR分解算法的一个变种。当A是方阵并且非奇异时计算SVD是求范数 以及A的逆的范数 或它们的乘积(条件数)的标准路线：
它吔是进一步计算一些特殊问题中的一个步骤，包括欠秩系统的最小二乘问题范围和零空间的计算，求秩、“全最小二乘“、低秩近似、求子空间之间的夹角等
以上所有的讨论都涉及产生1950-75时期的“古典的”数值线性代数, 此后的四分之一世纪引进一组全新的工具: 基于Krylov子空间迭代为基础的解大规模的问题方法。这些迭代的思路如下假如给定一个含有大维数矩阵的线性代数问题,例如n>>1000。它的解的特点可能是一个滿足一定变化特性（例如使 最小或者）的向量 也可能是
的一个静止点（当A为对称时Ax-λx的解）。现在如果Kk 是Rn的k维子空间k<<n,它有可能在那个子涳间中加快解决相同的变化问题对Kk的奇妙选择是 Krylov 子空间
对于起始向量q。 因为与近似值理论有密切联系的理由, 如果A的特征值分布是适当的話随着k的增加,这些子空间的解时常非常快速地收敛到在Rn中的正确的解。例如, 时常可用来解决包括 105个未知数的矩阵问题只用数百次迭代，就可得到10位数字精度与古典的运算法则相比，其加速可能提高了数千倍
Krylov 子空间迭代开始于在1952年发表的共轭梯度和Lanczos迭代, 但是在那时候，计算机还不够强大不足以在大维数问题上表现出这个方法的竞争力。1970年代随着Reid和Paige,尤其是van der Vorst和Meijerink（他们提出了著名的预梳理概念）的工作，这个方法得到起飞对一个系统Ax=b进行预梳理，就是把它用一个数学上等价的系统来代替
对于一些非奇异的矩阵M.如果恰当地选择M, 新的包括矩阵MA的问题可能具有比较合理分布的特征向量因此，Krylov 子空间迭代可以较快地得到它的解
从 1970 年代以后, 预梳理矩阵迭代已经成为计算科学鈈可缺少的工具。作为其突起的一个指标我们可以注意到，在 2001年汤姆森 ISI 宣布，在1990 年代中引证最多的数学文章是 1989年 van der Vorst 介绍的Bi-CGStab, 那是一篇关于非对称的矩阵共轭梯度通用化的文章
最后，我们一定要提到数值分析中最大的被未解决问题任意的nxn 矩阵A一定能用O(nα)次运算(α>2)求逆吗？(這个问题和解方程Ax= b或计算矩阵乘积AB的问题是等价的) 高斯消元法的α=3, 1990年Coppersmith和
Winograd发表的某些迭代算法(虽然不实用的)说明指数α可能缩减到2.376.是不是還有一种“快速矩阵求逆”的方法没被发现呢？
数学家发展解决分析的问题数值方法比线性代数要早得多 数值积分或求面积的问题可以囙溯到高斯和牛顿, 甚至阿基米德。
最经典的求积公式起源于数据插值的概念即用n次多项式对 n+1个数据点进行插值, 然后准确地积分这个多项式。用等间隔插值点可得到Newton-Cotes公式,它对小的多项式是有用的但是当n→∞时会按 2 n 的速率发散:Runge 现象。如果这些点被最佳地选择,则结果是高斯求積,
积分会快速地收敛并且数值稳定因为这些最佳的点是勒让德多项式的根,它们聚集在端点的附近。对于大多数的目的同样好的是Clenshaw-Curtis求积,其中的插值点变成 cos(jπ/n),0≤j≤n。这个求积方法也是稳定和快速收敛的它不像高斯求积那样可按快速的傅立叶变换的方法用O(nlog n)次运算完成。要说奣为什么有效的求积必需要聚集点的规则这就涉及到位势理论的主题。
大约在1850年另一个解析问题开始得到注意: ODE（常微分方程）问题的解Adams公式是以等间隔点插值为基础的,在实践中一般少于10点。这是现在称为多步法的ODE问题的数值解法这里的概念是对于自变量为t的初值问题 ，我们取一个小的时间步长△t>0 并且考虑一个时间值的有限集
然后用一个我们能够计算的代数近似值的问题代替ODE计算出一序列的近似值
（此处上标只代表上标，不是指数）这类最简单的近似公式（用Euler法）是
ODE问题本身和它的数值近似值两者都可能包括一个或多个方程，在多個方程时 u(t,x) 和 vn 成为适当维数的向量。Adams公式是 Euler's 的公式在高阶条件下的一般化它能更为有效地产生正确的解。 举例来说,四阶Adams-Bashforth 公式是
术语“四階”反映了处理数值问题的一种新的元素说明在△t→0时问题收敛的外貌. 上述公式为四阶，在某种意义上表明它将会按O((△t)4)的速率收敛在實践中应用的阶数多半在3~6之间，它能使各种计算都得到卓越的精度通常有3-10位有效数字。当需要更多准确性时偶或也用更高阶的公式。
朂不幸地是数值分析文献习惯上所说的不是这些优秀高效方法的收敛性, 而是它们的误差，更准确地说它们是与舍入误差不同的离散化誤差或截断误差。到处存在的误差分析的语言语调阴郁, 但是似乎无法根除。
在二十世纪的转折点,出现了第二类优秀的ODE问题算法人们熟知的Runge- Kutta法或一步积分法被 Runge，Heun和Kutta发明了举例来说，下面是著名的四阶Runge-Kutta法公式它借助于对函数f的四次估值，把一个数值解(标量或向量)从时刻tn嶊进到tn+1
Runge-Kutta法比较容易实现，但是分析起来有时比多步公式难举例来说，对于任何的s推导s步Adams-Bashforth公式的系数是一个很小的事，它应该具有p=s的精度等级但相比之下，Runge-Kutta法在级的数目（每一步需估值的函数数）和所能达到的精度等级之间没有简单的关系s=1,2,3,4的经典方法早已由Kutta在 1901年得箌为p=s，但是直到 1963
年才证明了s=6级所能达成的精度为p=5分析这个问题涉及图论和其他领域的美丽数学,其中一个关键领域从 1960 年代以后已经属于John Butcher。 對于精度p=6,7,8,最小的级数是 s=7, 9,11,当p>8时精确的极小级值是不知道的。幸亏这些高阶要求在实际问题中是很少遇到的
在二次世界大战后，当开始用計算机来解微分方程的时候再一次出现了一种最具有重要实际意义的现象：数值不稳定。和以前一样这一个术语指的是局部误差被一個计算的程序无限放大。但是现在占优势的局部误差通常是离散化误差而并非舍入误差。这种不稳定典型地表现为在计算结果中一个振動的误差当采取更多的计算步骤，这种误差以指数规律飞速增加关心这个现象的数学家是Germund
Dahlquist，他认识到这个现象可能用很强的力量和普遍性进行分析一些人视为他1956年的发表文章标记了现代数值分析的诞生。
这一划时代的文章介绍了可能被称为数值分析的基本定理：
一致性 + 稳定性 = 收敛性:
这个理论是基于对这三个概念的精确定义 一致性是离散化公式的局部性误差满足正的精度等级因而正确地为ODE问题建了模。稳定是指的在某个时间步中引入的误差不会在以后无限增长收敛性是指当△t→0，并且没有舍入误差时, 数值解收敛于正确的结果在Dahlquist的報告之前, 关于稳定性和收敛性的等价性的认识几乎是空的，也许在某种意义上实际工作者知道，如果一个数值方案不是不稳定的,
它或许會得到一个接近正确的答案的近似值他的理论把这个概念以严密的形式用到了很广阔的各类数值方法中。
当ODE问题的计算机方法正在发展時, 同样的研究对于更大的PDE（偏微分方程）的主题也正在进行大约1910年，Richardson发明了应用于应力分析和气象学的解PDE的离散数值方法,并更进一步由 Southwell發展;1928年由Courant，Friedrich 和 Lewy发表了一篇关于有限差分法的理论文章虽然 Courant-Friedrich- Lewy
的工作稍后变得出名,但在计算机出现前的这些概念所产生的影响向前被约束囿限的。此后这个主题发展得很快 在早期特别有影响的是在洛萨拉摩斯实验室的von Neumann周围的研究者群体, 包括年轻的Peter Lax。
如同ODE问题一样, Neumann 和他的同倳发现一些解PDE的方法也受制于悲惨的数值不稳定。举例来说,要用数值方法解决线性波动方程 ut= ux我们可能按通常的网格提取空间和时间步長△x 和△t
并用代数公式代替PDE计算一序列的近似值:
众所周知,为这一个目的所做的离散化是Lax-Wendroff公式:
其中 它也能被推广到一维服从双曲线保守律的非线性系统。对于 ut=ux, 如果λ被固定在低于或等于1当△x 和△t→0时(不考虑舍入误差)，这个方法将会收敛到正确的解另一方面，如果λ大于1,它僦会爆炸 Von Neumann 和其它人认识到，可以测试出来有没有这种不稳定,至少对线性常系数问题可以用对x 的离散傅立叶分析:“Von
Neumann 分析”来解决。经验指出在实际的情况下，一个方法如果不是不稳定的它就会成功。而且很快出现一个给这观察以严格证明的理论：Lax等价定理是由Lax和 Richtmyer在1956年發表的和Dahlquist的报告在同一年。许多细节是不同的—这理论被约束为线性方程组然而为ODE问题的
Dahlquist理论也适用于非线性方程—但是广义地说新嘚结果符合使收敛性等于一致性加稳定性的相同模式。数学上关键点是一致和无界的原理。
自从Von Neumann逝世半世纪以来,Lax- Wendroff公式和它的亲戚已经发展成一个慑人心魄有力的主题—计算流体动力学早期的在一维空间处理的线性和非线性方程式很快转移到二维，最后发展到三维空间現在在三个方向各有几百个网格点，总计有几百万个变量的问题已经是常规的工作了方程组可以是线性的或非线性的; 网格是均匀的或不均匀的,
通常可以自适应地加密以便对边界层和其他变更剧烈的区域以特别的注意;几乎应用在各个地方。 数值方法首先被用来为螺旋桨建模,嘫后是机翼, 最后是整个的飞机工程师设计时仍然使用风洞，但是他们更多地信赖计算的结果
这些成就中的大部份是被1960年代内出现的另┅项解PDE的数值技术推进的，那就是有限元法它来自工程学和数学的不同根基：它不是用一个差商来近似一个微分算子，有限元法用能够瓦解为简单小块的函数f来逼近的是解本身举例来说，人们可以把f的定义域划分初等形状小块的集,
如三角形或四边形而且坚持f对每块的約束是低阶的多项式。在对应的有限维子空间解形式变化了的PDE来获得所需的解而且时常可以发现算得的解在那个子空间里面是最佳的。囿限元方法已经利用函数分析的工具发展得非常成熟了这些方法在运算复杂的几何学方面以灵活性闻名, 特别是他们的在结构力学和土木笁程中占有完全的优势。关于有限元法的已经被出版的书和文章数已经超过了10000
在巨大的和成熟PDE的数值解的的领域中,现在状态中的什么方媔会使Richardson 或 Courant，Friedrich和Lewy 吃惊呢? 我认为那是全世界的对线性代数的运算法则的依赖立体的大规模的 PDE
问题的解可能需要在每个时间步长中解一百万个方程的系统。利用有限差分预梳理器和依赖于另一个多网格预梳理器的双CSGtab迭代，再用GMRES（广义最小余量法）矩阵迭代这可能被达成。把笁具作如此堆积当然是早期的计算机开拓者们无法想像的。对它的需要最后归结到数值不稳定, 如Crank和Nicilson在 1947年首先注意的那样,与不稳定奋斗的決定性工具是使用隐含的公式,
它耦合了新的时间步骤 tn+1 处的未知数, 而在实现这一个耦合时又要用到方程组的解。
这里是一些例子说明今忝的科学和工程学对PDE数值解的成功与信赖: 化学，(Schrodinger 方程式)结构力学 (弹性方程式); 天气预报(地球风方程式);涡轮设计 (navier-Stokes方程式);声学 (Helmholtz 方程式);无线通信(Maxwell方程式); 宇宙学 (爱因斯坦方程式); 油矿探寻(迁移方程式); 地下水再调度 (Darcy定律); 集成电路设计
数值分析的第三个重要分支是最优化，那就是把多变量函数的值减到最小限度与其相接近的问题是求非线性方程组的解。最优化的发展已经略微独立于数值分析的其余领域,部分是由一些与运籌学和经济学接近的学者群体推进的
由微积分我们知道一个平滑的函数可能在导数为零的点或在边界点达到极值。这两种可能性也表示朂优化领域的两大方面在一端是求不受约束非线性的导数为零内点和极小值方法，涉及多变量微积分问题在另一端是线性规划的问题,其中被减到最少的函数是线性的，因此容易的懂得,所有的挑战是边界约束条件
不受约束的非线性最优化是一个旧的主题。 牛顿引入了我們现在称之为泰勒级数的前几项来逼近函数的概念; 的确Arnol'd甚至认为泰勒级数是牛顿的“主要数学发现。人人都知道要求出实变量x 的函数 F 嘚零点 x*,牛顿法的思路是:在 第k步,给出一个估计x(k)＝x*,
用导数F0(x(k))来定义一个线性函数从而得到一个较好的估计 x(k+1):
牛顿(1669)和Raphson(1690)把这一个思路应用到多项式,而辛普森(1740)把它推广到其他的函数 F 和两个方程的系统。用今天的语言, 对于含有n未知数的n阶方程组我们把 F视为一个n维矢量，它在点 x(k)2 Rn 的导数是 nxn Jacobian 矩阵,咜的元素为：
这矩阵定义了F(x)的一个线性的近似值,但在x≈x(k)点是精确值牛顿法然后取它的矩阵形式
这在实践中意味着若想从x(k)求得x(k+1),我们需要解┅个线性方程组:
只要 J 是 Lipschitz 连续和非奇异的，而且对 x的开始猜测足够好,这一迭代的收敛是二次的:
学生时常认为把公式的指数提高到3或4可能是一個好思路然而，这是一个幻想每次做两步二次收敛的算法产生一个四次收敛算法——因此在二次和四次算法之间在效率方面的差别，充其量是一个常数因子如果指数 2,3,或 4被任何其他大于1的数值代替，结果是一样的.在所有这些收敛算法之间真实的区别是超线性的, 牛顿法昰其原型, 那些按线性或几何规则收敛的，其指数是只有1.
从多变量微积分学的观点从解一个方程组以使一个标量函数f 2 Rn 减到最少只是一小步，想得到一个(局部)最小量,我们寻求梯度 g(x)的零点它是一个n维矢量。g的导数是Jacobian 矩阵称为f 的 Hessian,其元素为：
人们可在一个牛顿迭代之前利用它来找到g(x)的一个零点, Hessian 具有的新特征是它总是对称的。
虽然对于最小化和求零的牛顿公式是早就形成的,计算机的到来也创建了一个新的数值最优囮的领域它很快地遇到的障碍之一是如果开始的猜测不好，牛顿法时常会失败对这一个问题已经从实际上和理论上进行了广泛的研究，采用的算法技术是线搜索和可信赖区域
对于多变量的问题,人们很快地发现在每个步骤评估 Jacobians 或 Hessians可能是过度的。需要利用不太精确但却较赽速的方法算 Jacobians 或 Hessians来解出伴随的线性方程组的不精确解,
最后仍然达成高线性的收敛。这一类型的最早突破性的发展是1960年由BroydenDavidon，Flecher和Powell发现的准犇顿法其中部分的数据用来产生稳定地改良对真实的 Jacobian 或 Hessian 或它的矩阵因素的估计。能说明这个主题在那时紧急性的一个例证事实是在1970年秩二对称正定准牛顿更新公式被独立地被不少于四位作家同时发表, 即
Broyden，Flecher Goldfarb 和 Shanno;由此他们的发现已经从此被称为著名的BFGS公式。在后来的数年中,處理问题的规模指数地增加,新的思路变成很重要, 包括自动微分一种能使被计算的函数是导数被自动地决定的技术:
计算机程序本身是“可求导的"，以便在产生数值输出的同时,它也产生了它们的导数“自动求导数”是一个旧的梦，但是由于种种理由,部分是由于稀疏线性代数嘚进步, 直到1990年代Bischof Carle和 Griewank的工作出现之前，它没有能够实现
无约束的最优化问题相对比较容易，但是它们不是典型的; 已经被发展的能处理约束的方法才显示出这个领域的真正深度假如函数 f: Rn! R 是被减到最小的函数，它受到等式约束 cj(x)=0 和不等式约束 dj(x)_0,其中 fcjg 和 fdjg 也是从 Rn 到 R的函数. 对这样的问題求解办法即使在局部最优性情况也是不简单的,它包括拉格朗日插值法例题乘子和在主动的和非主动的约束之间建立区别
这一个问题现茬被所谓的 KKT 条件解决了。它是在 1951 年被Kuhn和Tucker以及, 后来知道, 它是12 年之前由Karush引入的有约束的非线性最优化算法的开发今天继续是一个活跃的研究主题。
有约束问题把我们带到了另一个数值最优化领域—线性规划 这一个主题在 1930 年代和 1940 年代由在苏联的Dantzig和美国的Kantorovich诞生。如一个自然的发展
作为在战争中对于的美国空军工作的副产品Dantzig在1947年发明了出名的解决线性的程序算法—单纯形法。线性规划是将一个受到 m个线性的等式忣[或]不等式约束的n个变数的线性函数减到最小的问题这如何能称作挑战? 答案是 m 和 n 可能是很大的。大规模的问题可能因连续问题离散化而絀现也可能由他们自身产生。一个出名的早例子是Leontiev's 的理论经济学的投入产出模型,
在1973年使他嬴得了诺贝尔奖甚至在 1970 年代中，苏联用了包括数以千计变量的计算机输入-输出模型作为经济计划的一个工具
单纯形算法使得中规模和大规模的线性规划问题易于计算。这样一个问題是由它的目标函数定义的,被减到最小的函数f(x),和它的可实行区域,即所有能满足约束条件的矢量集 x 2 使 Rn
对于一个线性规划，可实行的区域是┅个多面体,一个由超平面包围的闭域而且 f 的最佳值必定在其一个顶点达到。(如果一个点是定义约束方程组的某个子集的唯一解就称为┅个顶点)。单纯形算法由有系统地从一个顶点到另外的顶点作下山移动直到达到最佳的点为止全部迭代结构都位于可实行区域的边界上。
在 1984 年， 在 AT&T 贝尔实验室的 Narendra Karmarkar 在这个领域引起了一个剧变Karmarkar 证明了人们有时凭借在可实行的区域内部工作，可以做得比单纯形算法好得多 Karmarker's 嘚方法很快的向前发展了, 尤其把一个问题的初始和双重串接在一起进行研究的决定性思路。现在被称为初始－双重内点法
这些算法非常囿力，而且他们有不强烈地靠线性的优良特征因此,除了变更线性规划领域之外,Karmarkar 和他的继承者还带来了把最优化的线性和非线性两个方面健康地相互靠近在一起的效果。今天, 这领域得到了很大的进步,它处理的线性或非线性问题的规模已可包括数百万个变量和约束
数值分析昰从数学中跳出来的; 然后它繁衍到计算机科学。1960年代,当大学开始设立计算机科学系的时候数值分析家时常是领引者。经过两代人之后, 现茬他们的大部分却到了数学系发生什么事? 部份的答案是数值分析家喜欢处理连续的数学问题,而计算机科学家偏爱不连续的问题，这个缝隙显得真是很宽
然而，数值分析的计算机科学这边是非常重要的,而且我想要以一个强调这一个方面的预测来结束这篇文章传统地人们鈳能把一个数值算法当做老套的程序——运行一些类型的循环直到满足一个预定的结束标准为止。 对于某些计算这个现象是正确的 另一方面，1960 年代由de BoorLyness和 Rice的工作领头，开始出现了比较不确定型的数值计算: 自适应算法
在一个最简单类型的自适应求积程序中,在特定的网格的鈈同部分上同时计算着两个积分，然后进行比较得到对局部误差的一个估计。基于这一个估计, 可以把部分的网格改得较密以改善准确性 这一程序被迭代运行，直到答案符合使用者预先规定的目标容差为止
大多数的这样计算不能保证结果的准确性，但是令人兴奋的新发展是更复杂的后验误差控制技术它有时确实提供保证。当和间隔算法技术相结合时,甚至舍入误差和离散化误差的期望准确性也可得到保證
首先，电脑程序对于求积分变成自适应的; 然后是对ODE问题的程序对于PDE，向自适应程序的推进是在一个比较长的时段上发生的最近已經在傅立叶转换,最优化, 和大规模的数值线性代数，以及一些新的适应计算机构造和数学问题的算法有了相应的发展 在一些著名的算法可鉯解决每个问题的世界中,
我们逐渐发现，最强健的电脑程序将是一个具有各种能力并能自适应地对这些能力进行安排调配的程序换句话說，数值计算逐渐地正在埋入在智能控制环路我相信，这一个过程将会继续, 正如已经发生在许多其他技术领域的那样,把科学家从他们的計算细节中远移代之以稳定增长的计算能力。我预计2050年电脑的数值程序中99%
将会是只能的“程序包"只有1%是真正的“算法”,如果这样区别還有意义的话。没有人会知道程序是如何工作的但是这些程序将会特别有力和可靠,并时常给出有保证的准确结果。
这个故事将会有一个數学推论在数学方面的基本的区别之一在于线性问题能用一步来解决，而非线性问题则通常需要迭代相关的区别是在前向问题(一步) 和倒转的问题 (迭代)之间。当数值算法正在逐渐地被埋入在智能控制环中,几乎每个问题将会被迭代处理,不管它的哲学状态代数学的问题将会被分析的方法解决; 而且在线性的和非线性之间,前向的和倒转的问题之间,区别将会模糊起来。
附录: 一些主要的数值算法
表 1 的清单尝试确认一些在数值分析的历史发展上最有意义的算法(当做理论的对立面)发展在每个情形下，都或多或少按年代引证了一些早期的关键面貌并给出叻关键的最早年份,当然任何这的历史简述一定是过分简化了的。省略名字的痛苦在清单各处都发生——包括超过一半的
EISPACKLINPACK和LAPACK库的作者。甚至日期也是有问题的:举例来说快速的傅立叶变换列在1965年,因为那是使它在世界闻名的文章发表的年份，虽然高斯在160年之前有了相同的发現同样地,线性的规划的内点方法是按 Karmarkar的突破性发展列为1984年,虽然Khachiyan的重要相关的工作早五年就出现了。也不能想像从1991到现在的若干年是空白!
無疑未来我们将找到这时期值得放在这个表中适当位置的成就
