中学时代，没有什么比一张数学图更可怕的了。许多问题没有现成的解决方案，而其他问题是可以解决的，但需要大量的计算工作。在这种情况下，您需要算法来提供帮助。希望大家不要被算法吐槽。万一你真的吐了，往好的方面看：另一顿午餐:]在本教程中，你将学习算法的概念以及如何使用它们来解决难以分析的问题。使用playgrounds，可以很容易地可视化解决方案以便于审查。即使您不是数学爱好者，对物理或计算机科学不感兴趣，您仍然会在本教程中找到一些有趣的东西。您所需要的只是一些微积分和基础物理知识。您将学习如何使用数值算法解决两个难题。但为了清楚起见，这两个问题也可以解析求解。当分析方法不起作用时，算法是理想的，即使这样，通过比较这两种方法也更容易理解它们的本质。什么是数值算法？简单来说，数值算法是一类不依赖于封闭形式解析解的解决数学问题的方法。问题来了，什么是闭式解析解？如果有一个公式，我们可以利用已知数，通过有限的运算步骤得到一个准确的结果，这就叫做闭式解析解。更简单地理解一下：如果你能用代数找到一个表达式来解决一个未知的问题，并且你可以通过代入一些已知的数来得到结果，那么就意味着你得到了一个封闭形式的解析方法。何时使用数值算法？许多问题没有现成的解决方案，而其他问题是可以解决的，但需要大量的计算工作。在这种情况下，您需要算法来提供帮助。例如，假设您需要编写一个物理引擎来计算许多对象在有限时间内的行为。此时你可以使用数值算法来更快地得到结果。这也有副作用，更快的计算意味着结果不会很准确。然而，在大多数情况下，这样的结果就足够了。天气预报得益于数值计算。天气变化的速度和影响它的因素数量惊人。这是一个非常复杂的系统，只有数值模拟才能完成预测未来的任务。可能是因为缺少这些算法，你的iPhone一直告诉你要下雨，即使外面仍然阳光明媚。作为热身，让我们玩一个游戏，您可以找到给定数字的平方根。两种方法都将使用二分法。令人惊讶的是，您可能已经知道这种方法，但不知道它的名字。回想一下，我们小时候可能玩过一个游戏，其中选出一个不超过100的数字，然后其他人来猜。你只会提示他猜数字是大是小。游戏开始。小明让小强开始猜，小强说我猜1，小明说太小了，小强又猜了2，小明说太小了。小强接着猜了3、4……最后选了5，终于猜对了。5步就猜到了，还不错。但是如果小明选择78，要猜对还需要一些时间。这个游戏如果用二分法来玩，猜的速度会快很多。二分法我们知道数字在1到100之间，所以除了一一猜，还是乱猜。我们将数字分为两个区间：a=[1,50]，b=[51,100]。接下来我们判断这个数在哪个区间，这个很简单，问这个数是不是50就可以了，如果小于50，那么就不考虑区间b了。接下来，我们继续将区间a分成两半，重复这一步。例如：数字为60。区间为：a=[1,50],b=[51,100]。第一步，小强说50（也就是区间a的上限），小明说太小了。这时，小强知道这个数一定在b区间。第二步，将区间b分解为：c=[51,75],d=[76,100]。还是选择c区间上限75，小强告诉他太大了。这意味着该数字必须在c区间内。继续分解。..使用这种方法，可以在7次试验中得到正确答案，每次试验需要进行60次试验。1.50->小2.75->大3.62->大4.56->小5.59->小6.61->大7.60->对！要计算x的平方根，过程类似。数字x的平方根介于0和x之间。也可以记为：(0,x)。若数x大于等于1，则记为：[1,x]。对区间进行分解得到a=(0,x/2],b=(x/2,x)。如果x为9，则区间为[1,9]，分解区间为a=[1,5],b=(5,9]，则中间值m为(1+9)/2=5，接下来检查m*m-x是否大于期望精度，这里检查m*m是否大于或小于等于x，如果更大，我们继续使用(0,x/2]区间检查，否则，使用(x/2,x]区间。看实际步骤，初始化m=5，期望准确值为0.1:1。计算m*m-x:5*5-9=112。检查结果是否小于或等于预期值：25-11<=0.1?3。显然如果不满足，继续下一个区间：m*m是否大于9？是的。接下来使用区间[1,5]，测试值为(1+5)/2=3。4.计算m*m-x:3*3-9=05。检查是否满足预期值：9-9<=0.1?6。完毕。注意：您可能想知道括号是什么意思？区间有固定的格式（下限和上限）。不同的括号代表不同的区间边界。括号表示边界不在区间内（即开区间），方括号表示边界在区间内（闭区间）。如(0,a]包含a但不包含0。上例中：区间(0,a/2]包含a/2但不包含0；区间(a/2,a]表示都大于a/2，一个数小于a。在Playground上使用二分法现在，是时候应用我们学到的理论了。让我们自己实现二分算法。创建一个新的playground文件并添加以下代码：funcbisection(x:Double)->Double{//1varlower=1.0//2varupper=x//3varm=(lower+upper)/2varepsilon=1e-10//4while(fabs(m*m-x)>epsilon){//5m=(lower+upper)/2ifm*m>x{upper=m}else{lower=m}}returnm}查看代码各部分含义：1.设置区间下限为12。设置区间上限为x3找到中间值，定义想要的精度度4.检查运算是否能满足精度5.如果不能，找到新的中间值，定义新的区间，继续搜索.添加以下代码来测试函数：letbis=bisection(2.5）在这一行右边的m=(lower+upper)/2中，可以看到代码执行了35次，也就是说需要35步才能找到结果。#p#立即我们可以看到playground的一个可爱功能的好处：可以查看值的完整历史记录。二分法可以成功逼近实际结果。通过数值历史图，您可以看到数值算法是如何逼近正确解的。按option+cmd+enter打开辅助编辑器，点击m=(lower+upper)/2行代码右边的圆圈按钮，在辅助编辑器中添加一条历史记录。你会看到方法一点一点跳转到正确的结果。古典数学也决定了一个算法需要追溯到远古时代。它起源于公元前1750年的巴比伦，在100年前亚历山大城的希耶罗所著的《度量论》（Metrica）一书中有所描述。这就是为什么它被称为“英雄方法”。您可以通过此链接了解更多关于Hiro的信息。这个方法用到了函数，这里要计算a的平方根。如果你能找到函数曲线的“零点”，函数的值为0，那么找到x的值将得到a的平方根。为此，我们首先选择一个任意的x值作为起始值，计算对应于该值的点的切线（函数的切线），然后查看切线上的0点（函数的交点）线和x轴）。然后我们再次使用这个0点作为起始值，重复多次直到满足精度。因为每一次计算正切值都会更接近真实的0值，这个过程会逐渐接近真实的结果。下图为求解a=9时a的平方根，起始值为1。起始点x0=1，生成红色切线，然后生成x1点，生成紫??色线；x2分，连接蓝线，终于找到正确答案了。当古典数学遇到Playground，我们就有了很多古代巴比伦人没有的好东西，比如：playground。在操场上添加以下代码并检查：funcheron(x:Double)->Double{//1varxOld=0.0varxNew=(x+1.0)/2.0varepsilon=1e-10//2while(fabs(xNew-xOld)>epsilon){//3xOld=xNewxNew=(xOld+x/xOld)/2}returnxNew}这段代码做了什么？1.创建一个变量来存储当前结果。xOld是上次计算的结果，xNew是实际结果。找到一种方法来初始化xNew是一个很好的起点epsilon是所需的精度在这个例子中，我们计算平方根到小数点后10位。2.使用while循环查看是否达到了所需的精度。3.如果达不到精度，则将xNew设置为xOld，继续下一次迭代。使用下面的代码来验证这个函数的作用：lether=heron(2.5)Hero方法只需要迭代5次就可以找到正确的结果。点击xNew=(xOld+x/xOld)/2行右侧的圆角按钮，添加值历史记录，可以看到第一次迭代可以找到一个很好的近似值。模拟谐波振荡器运动在本节中，我们将研究如何使用数值积分算法来模拟简谐系统的运动——一个基本的动力系统。该系统可以描述许多现象，例如钟摆的摆动和弹簧的振动。特别是，它可以描述物体在一段时间内移动一定偏移量的场景。对于此示例，假设质量附加到弹簧。为了简化问题，我们忽略阻力和重力。所以唯一起作用的力是弹簧的弹力，它把物体拉回原来的位置。在这样的假设下，你只需要用到两条力学定律就可以解决问题：牛顿第二运动定律，它描述了物体所受的力与加速度之间的关系。胡克定律，描述弹力与物体偏移量之间的比例关系。这里k是弹簧系数，x是物体的偏移量。简单的振动方程因为弹力是作用在物体上的唯一力，所以我们整理以上两个方程：简写：也记为，即共振频率的平方。更精确的方程式如下：其中A是振幅，这里是物体的偏移量，称为相位差。这两个值最初都设置为常量值。如果设时间t=0，那么，,可以简单的算出幅值和相位差：看个例子，我们有一个质量为2kg的物体，连接着一个弹簧，弹簧的力系数为k=196牛/米。在初始时间t=0，弹簧移动了0.1米。我们可以用公式计算振幅、相位差和谐振频率。在Playground上做实验：在使用这个公式计算任意时间点对应的值之前，我们需要添加一些代码。回到playground文件，在末尾添加如下代码：//1typealiasSolver=(Double,Double,Double,Double,Double)->Void//2structHarmonicOscillator{varkSpring=0.0varmass=0.0varphase=0.0varamplitude=0.0vardeltaT=0.0init(kSpring:Double,mass:Double,phase:Double,amplitude:Double,deltaT:Double){self.kSpring=kSpringself.mass=massself.phase=phaseself.amplitude=amplitudeself.deltaT=deltaT}//3funcsolveUsingSolver(solver:Solver){solver(kSpring,mass,phase,amplitude,deltaT)}}这些代码块做了什么？1.定义一个函数类型的类型别名。该函数有5个Double类型的参数并返回空。2.创建一个结构来定义谐波振荡器。3.这个函数只是简单的创建了Clousure，用来解决振动问题。（没有真正的计算代码）我们来看看精确求解代码如下：funcsolveExact(amplitude:Double,phase:Double,kSpring:Double,mass:Double,t:Double){varx=0.0//1letomega=sqrt(kSpring/mass)vari=0.0whilei<100.0{//2x=amplitude*sin(omega*i+phase)i+=t}}该方法包含求解运动方程所需的所有参数。1.计算共振频率2.在while循环中计算物体当前位置，i作为下一次计算的自增变量。添加以下测试代码：letosci=HarmonicOscillator(kSpring:0.5,mass:10,phase:10,amplitude:50,deltaT:0.1)osci.solveUsingSolver(solveExact)这个解决方案中的方法有点奇怪：它传递了参数in，但没有返回数据，也没有显示任何内容。#p#这样做有什么好处？使用该函数的目的是在while循环中计算振动过程中的具体动态值。在Playground中，可以观察这些动态值的历史。在x=amplitude*sin(omega*i+phase)处添加值历史，我们可以看到运动的轨迹。既然已经验证了第一个精确算法，那么我们来看看数值计算方案。欧拉法欧拉法是一种数值积分方法。该方法于1768年由莱昂哈德·欧拉(LeonhardEuler)在InstitutionesCalculiIntegralis《积分学原理》一书中首次提出。要查看该方法的详细信息，请参考：http://en.wikipedia.org/wiki/Euler_method欧拉方法的核心思想是利用折线来逼近一条曲线。具体方法是计算给定点的斜率，然后绘制一条具有相同斜率的折线。在这条多段线的末端，继续计算斜率，绘制另一条线。如您所见，算法的准确性取决于折线的长度。你想知道deltaT是做什么的吗？该数值算法中使用步长，这对算法的准确性很重要。大步长导致精度低，但执行块更快。反之，精度会提高，速度会降低。deltaT变量表示步长。我们将这个值初始化为0.1，这意味着我们每0.1秒计算一次对象的位置。在欧拉算法中，表示折线在x轴上的投影长度为0.1。在写代码之前，需要再看看这个公式：二阶微分方程可以转化为两个一阶微分方程。可以写成and我们可以用差商来做换算：又得到：有了这些方程，我们就可以直接实现欧拉法了。在solveExact方法后添加代码：funcsolveEuler(amplitude:Double,phase:Double,kSpring:Double,mass:Double,t:Double){//1varx=amplitude*sin(phase)letomega=sqrt(kSpring/mass)vari=0.0//2varv=amplitude*omega*cos(phase);varvold=vvarxoldEuler=xwhile<100{//3v-=omega*omega*x*t//4x+=vold*txoldEuler=xvold=vi+=t}}内容:1.设置当前位置和omega值。2.计算当前速度。3.在while循环中，新的速度是通过一阶函数计算出来的。4.使用速度计算新位置，最后，将i的值增加步长t。现在，添加以下代码进行测试：osci.solveUsingSolver(solveEuler)在xoldEuler=x位置添加历史值，并查看。您会看到此方法显示一个振幅递增的正弦函数。因为欧拉法并不是一个精确的算法，这里0.1的步长过大也会导致计算结果不准确。这是另一个步长为0.01的函数图像，效果明显更好。所以，在使用欧拉法的时候，你要想到步长越小，结果越好。但是通过折衷的步长，它更容易实现。VelocityVerlet最后讨论的算法称为VelocityVerlet。其思路与欧拉法相同，只是计算新位置的方式略有不同。欧拉在计算新位置时忽略了实际加速度，公式为：，而速度Verlet算法在计算时考虑了加速度：。当步长相等时，结果更好。在solveEuler方法后添加新代码：funcsolveVerlet(amplitude:Double,phase:Double,kSpring:Double,mass:Double,t:Double){//1varx=amplitude*sin(phase)varxoldVerlet=xletomega=sqrt(kSpring/mass)varv=amplitude*omega*cos(phase)varvold=vvari=0.0while<100{//2x=xoldVerlet+v*t+0.5*omega*omega*t*tv-=omega*omega*x*txoldVerlet=xi+=t}}代码内容：1、循环前所有代码与欧拉法中相同。2、根据速度计算新的位置，增加i的值。在Playground中测试代码：osci.solveUsingSolver(solveVerlet)不要忘记在xoldVerlet=x代码行上添加值历史记录：下一步做什么？您可以通过此链接下载完整的项目。希望您在这次数字世界旅行中玩得开心。了解一些算法，即使只是一些有趣的经典数学历史，也会对您有所帮助。
