fortran常用算法的程序设计举例

1第八章常用算法的程序设计举例第一章算法第二章计算机和计算机程序第四章逻辑运算和选择结构第五章循环结构的实现第六章Fortran的数据结构第七章数据的输入、输出第三章Fortran语言程序设计初步2一、数值积分几何意义：badxxf)(3近似求小曲边梯形面积的方法：（1）用小矩形代替小曲边梯形；（2）用小梯形代替小曲边梯形；（3）在小区间范围内，用一条抛物线代替该区间内的f(x)。4read(*,*)a,b,nx=ah=(b-a)/nf0=exp(x)s=0.0do10i=1,nsi=f0*hs=s+six=x+hf0=exp(x)10continuewrite(*,100)a,b,nwrite(*,200)s100format(1x,'a=',f10.3,3x,'b=',$f10.3,3x,'n=',i4)200format(1x,'s=',f15.8)end1.矩形法))1((:hiafhsii个小矩形面积第10:dxex求5read(*,*)a,b,nx=ah=(b-a)/ns=0.0do10i=1,nsi=(sin(x+(i-1)*h)+$sin(x+i*h))*h/2.0s=s+si10continuewrite(*,100)a,b,nwrite(*,200)s100format(1x,'a=',f10.3,3x,$'b=',f10.3,3x,'n=',i4)200format(1x,'s=',f15.8)end2.梯形法hhiafihafsii2))1(()(:个小梯形面积第10sin:xdx求6其他几种程序变形...f1=sin(a)...do10i=1,nf2=sin(a+i*h)si=(f1+f2)*h/2.0s=s+sif1=f210continue......x2=a...do10i=1,nx1=x2x2=x2+hsi=(sin(x1)+$sin(x2))*h/2.0s=s+si10continue......f0=sin(a)h=(b-a)/ns=f0do10i=1,nf=sin(a+i*h)s=s+2.0*f10continues=(s-sin(b))*h/2.0...012()2())2bnnahfxdxfffff73.Sinpson法取a，b中点c—((a+b)/2，0)，通过f(a)，f(b)，f(c)三点可作唯一一条抛物线f1(x)。根据抛物线定积分求值公式，有：2)]()(4)([3)(1abhbfcfafhdxxfba其中如果将(a，b)分成两个小区间(a，c)和(c，b):12(){()()4[()(3)]2(2)}3bafxdxsshfafbfahfahfah2ach其中8如果将(a，b)分成四个小区间:(){()()4[()(3)(5)3(7)]2[(2)(4)(6)}bahfxdxfafbfahfahfahfahfahfahfah42abh其中如果将(a，b)分成n个小区间:(){()()4[()(3)((21))]32[(2)(4)((22))}bahfxdxfafbfahfahfanhfahfahfanhnabh2其中9read(*,*)a,b,nh=(b-a)/(2.0*n)s=0.0fa=1.0/(1.0+a)fb=1.0/(1.0+b)x=a+hf2=0.0f4=1.0/(1.0+x)do10i=1,n-1x=x+hf2=f2+1.0/(1.0+x)x=x+hf4=f4+1.0/(1.0+x)10continues=h/3.0*(fa+fb+4.0*f4+2.0*f2)write(*,100)a,b,nwrite(*,150)s100format(1x,'a=',f8.2,2x,'b=',f8.2,$2x,'n=',i4)150format(1x,'s=',f16.7)end))22(()4()2(2))12(()3()(4)()(3hnafhafhafhnafhafhafbfafhs小区间面积101:xdx求三种求定积分的方法中，矩形法的误差较大，梯形法次之，辛普生法最好。10二、解一元方程（解非线性函数）1.直接迭代法read(*,*)x,mdo10i=1,mx1=(-x**3-2.0*x*x-2.0)/2.0write(*,100)i,x1if(abs(x-x1).gt.1e-6)thenx=x1elsestopendif10continuewrite(*,200)m100format(1x,'i=',i3,5x,'x1=',f15.7)200format(1x,'computationhasnot',$'convergedafter',i4,'iteration')end022223xxx2)22(23xxx因有收敛问题，要设最大循环次数。11有的g(x)是收敛的，而有的g(x)是不收敛的。同一个g(x)，对某些x0是收敛的，对有的x0则是不收敛的。如果g(x)具有一阶导数连续，且对于所有的x，若|g’(x)|≤q1（q为一个定数），那么x=g(x)对于任意的x0均收敛，且q愈小，收敛速度愈快。如果不满足对所有的x存在|g’(x)|≤q1，则可能对有的x0收敛，对有的x0不收敛。因此要恰当的选择g(x)形式和初值x0。122.牛顿迭代法read(*,*)xn=110x1=xf=x1**3-2.0*x1**2+4.0*x1+1.0f1=3.0*x1**2-4.0*x1+4.0x=x1-f/f1write(*,100)n,x1,xn=n+1if(abs(x-x1).gt.1e-6)goto10100format(1x,'n=',i3,3x,'x1=',f15.7,$3x,'x=',f15.7)end322()241=00()344fxxxxxfxxx求在附近的一个实根1()()nnnnfxxxfx1112()()fxfxxx133.二分法5read(*,*)x1,x2f1=x1**3-6.0*x1-1.0f2=x2**3-6.0*x2-1.0if(sign(f1,f2).eq.f1)goto510x=(x1+x2)/2.0f=x**3-6.0*x-1.0if(sign(f,f1).eq.f)thenx1=xf1=felsex2=xf2=fendifif((abs(x1-x2).gt.1e-5).and.$abs(f).gt.1e-6)goto10if(abs(f).gt.1e-6)x=(x1+x2)/2.0write(*,100)x100format(1x,'x=',f15.7)end3()6102fxxxx用二分法求在附近的一个实根144.弦截法（割线法）5read(*,*)x1,x2f1=x1**3-2.0*x1**2+7.0*x1+4.0f2=x2**3-2.0*x2**2+7.0*x2+4.0if(sign(f1,f2).eq.f1)goto5f=1.020if((abs(x1-x2).gt.1e-5).and.$abs(f).gt.1e-6)thenx=x2-(x2-x1)/(f2-f1)*f2f=x**3-2.0*x**2+7.0*x+4.0if(sign(f,f1).eq.f)thenx1=xf1=felsex2=xf2=fendifgoto20endifif(abs(f).gt.1e-6)x=(x1+x2)/2.0write(*,100)x100format(1x,'x=',f15.7)end)()()(212122xfxfxfxxxx047223xxx15以上方法都是近似求根，得到不是准确值而是近似值。但只要给定的误差足够小，就可以认为它们之间足够近似。事实上，只有少数的方程才能用解析的方法求出准确的根值。计算机可以解任何有实根的一元方程，它采用的基本方法就是迭代，经过多次迭代，使近似根逐渐趋近于真实根。迭代可以用循环实现，正好充分发挥计算机快速运算的特点。16三、求函数极值的极小值求54)(:2xxxfFibonacci搜索算法，或0.618法，或黄金值搜索法。17f1=x1*x1-4.0*x1+5.0f2=x2*x2-4.0*x2+5.0if(f1.gt.f2)thenf=f2x=x2elsef=f1x=x1endifwrite(*,204)x,f200format(12x,'x1',14x,'f1',$13x,'x2',13x,'f2'/)202format(1x,4f15.7)204format('0','x=',f10.6,5x,$'f(x)=',f10.7)endreallow,high,x1,x2read(*,*)low,highwrite(*,200)x1=low+0.618*(high-low)x2=high-0.618*(high-low)10if(high-low.gt.1e-4)thenf1=x1*x1-4.0*x1+5.0f2=x2*x2-4.0*x2+5.0write(*,202)x1,f1,x2,f2if(f1.gt.f2)thenhigh=x1x1=x2x2=high-0.618*(high-low)elselow=x2x2=x1x1=low+0.618*(high-low)endifgoto10endif18五、计算机模拟计算机模拟（ComputerSimulation）,又称“仿真”：用计算机模仿实物系统进行测试，从测试的结果获得期望的资料。根据模拟对象的不同特点，可分为：确定性模拟（DeterministicMode）；随机性模拟（StochasticMode）。19小球以10m/s沿45°斜抛，落地反弹方向同前，速度减小10%，求前三次周期轨迹。20if(abs(y).gt.0.001)thent=t-dtdt=0.1*dtt=t+dtx=v*t+x0y=v*t-0.5*g*t**2elsewrite(*,100)t,x,yflag=.false.endifendifgoto5endifv=0.9*vt=0.0x0=xwrite(*,*)10continue100format(1x,'t=',f8.4,3x,'x=',$f12.6,3x,'y=',f12.6)endlogicalflagparameter(g=9.8)read(*,*)v0,dt=0.0x0=0.0v=v0*cos(3.1415926/4.0)do10i=1,3dt=dflag=.true.x=v*t+x0y=v*t-0.5*g*t**25if(flag)thenif(y.ge.0.0)thenif(dt.ge.d)$write(*,100)t,x,yt=t+dtx=v*t+x0y=v*t-0.5*g*t**2else21上机目的：1.掌握数值积分的矩形法、梯形法和Simpson法；2.掌握求解一元方程的迭代法、牛顿迭代法、二分法和弦截法；3.掌握计算机模拟的一般方法。上机内容：1.调试课本中的所有程序；2.习题第5、10题。