弦截法非线性方程求解.doc

《MATLAB程序设计实践》课程考核 一、,编程实现以下科学计算算法，并举一例应用之。（参考书籍《精通MATLAB科学计算》，王正林等著，电子工业出版社，2009年） “弦截法非线性方程求解” 1、算法说明: 弦截法的算法过程如下： 过两点(a,f(a)),(b,f(b))作一直线，它与x轴有一个交点，记为x1.. 如果f(a)f(x1)0,过两点(a,f(a)),(x1,f(x1))作一直线，它与x轴的交点记为 x2，否则过两点 (b,f(b)),(x1,f(x1))作一直线，它与x轴的交点记为x2 ； 如此下去，直到|xn-xn-1|，就可以认为xn为f(x)=0在区间[a,b]上的一个根。 Xk的递推公式为：且 在MATLAB中编程实现的弦截法的函数为：Secant. 功能：用弦截法求函数在某个区间的一个零点。 调用格式：root=Secant(f,a,b,eps). 其中，f为函数名； a为区间左端点； b为区间右端点； eps为根的精度； root为求出的函数零点。 2、流程图: 3、M文件: function root=Secant(f,a,b,eps) %弦截法求函数f在区间[a,b]上的一个零点 %函数名：f %区间左端点：a %区间右端点：b %根的精度：eps %求出的函数零点：root if(nargin==3) eps=1.0e-4; end f1=subs(sym(f),findsym(sym(f)),a); f2=subs(sym(f),findsym(sym(f)),b); if(f1==0) root=a; end if(f2==0) root=b; end if(f1*f20) disp(两端点函数值大于0!); return; else tol=1; fa=subs(sym(f),findsym(sym(f)),a); fb=subs(sym(f),findsym(sym(f)),b); root=a-(b-a)*fa./(fb-fa); %迭代初始值 while (toleps) r1=root; fx=subs(sym(f),findsym(sym(f)),r1); s=fx*fa; if(s==0) root=r1; else if(s0) root=b-(r1-b)*fb/(fx-fb); %用递推公式2 else root=a-(r1-a)*fa/(fx-fa); %用递推公式1 end end tol=abs(root-r1) end end 4、程序应用举例: 采用弦截法求方程lgx+＝2在区间[1，4]上的一个根。 解：在MATLAB命令窗口中输入： r=Secant(sqrt(x)+log(x)-2,1,4) 输出计算结果为： r = 1.8773 由计算结果知方程在lgx+＝2在区间[1，4]上的一个根为1.8773 5、运行流程图: 6、运行结果: 二.分析单自由度阻尼系统的阻尼系数对其固有振动模态的影响。 1、算法说明: 根据题目意思可理解为在一定范围内取不同的阻尼系数值，根据单自由度阻尼系统的动力学方程，则可以分别算出对应的振动曲线，即振动规律曲线。用matlab进行编程计算，通过调用subplot函数分别绘制x-t二维和三维图形，实现其振动规律曲线，并通过曲线对阻尼振动规律进行分析。 2、流程图 3、M文件（U.m）： wn=10; %公共参数 x0=1; %初始位置 v0=0; %初始速度 tf=2; %终点时间 for j=1:10 zeta(j)=0.1*j; %阻尼系数 wd(j)=wn*sqrt(1-zeta(j)^2); A=sqrt((wn*x0*zeta(j)+v0)^2+(x0*wd(j))^2)/wd(j); %振幅A a=atan2(wd(j)*x0,v0+zeta(j)*wn*x0); %相位角 t=0:tf/1000:tf x(j,:)=A*exp(-zeta(j)*wn*t).*sin(wd(j)*t+a); end subplot(1,2,1) %绘制x,t二维图象 for j=1:2:10 plot(t,x(j,:)) text(0.3,x(j,151),[zeta=,num2str(