两步迭代法线性方程组求解(共10页)

1、精选优质文档-倾情为你奉上MATLAB程序设计实践课程考核一 编程实现以下科学计算算法，并举一例应用之（参考书籍精通MATLAB科学计算，王正林等著，电子工业出版社，2009年）“两步迭代法线性方程组求解”（1）算法说明两步迭代法的一般公式为： yn=g（xn） xn+1=h(yn)它比较常用两种形式：（1） yn=xn-f(xn)f'(xn) ; xn+1=yn-f(yn)f'(xn)（2） yn=xn-f(xn)f'(xn) ; xn+1=yn-f(yn)(yn-xn)2fyn-f(xn)在MATLAB中编程实现的两步迭代法TwoStep。功能：用两步迭代法求函数在

2、某个区间上的一个零点。调用格式：root=TwoStep(f,a,b,type,eps)。其中，f为函数名； a为区间左端点； b为区间右端点； type为两步迭代法的形式； eps为根的精度； root为求出的函数零点。（2）程序代码function root= TwoStep(f,a,b,type,eps)%两步迭代法求函数f在区间a,b上的一个零点%函数名：f%区间左端点：a%区间右端点：b%两步迭代法的形式：type（可取1,2）%根的精度：eps%求出的函数零点：rootif(nargin=4) eps=1.0e-4;endf1=subs(sym(f),findsym(sym(f),

3、a);f2=subs(sym(f),findsym(sym(f),b);if(f1=0) root=a;endif(f2=0) root=b;endif(f1*f2>0) disp('两端点函数值乘积大于0！'); return;else tol=1; fun=diff(sym(f); fa=subs(sym(f),findsym(sym(f),a); fb=subs(sym(f),findsym(sym(f),b);dfa=subs(sym(fun),findsym(sym(fun),a);dfb=subs(sym(fun),findsym(sym(fun),b); i

4、f(dfa>dfb) root=a; else root=b; end while(tol>eps) %第一种形式的两步迭代法 if(type=1) r1=root; fx1=subs(sym(f),findsym(sym(f),r1); dfx=subs(sym(fun),findsym(sym(fun),r1); r2=r1-fx1/dfx; %第一步迭代 fx2=subs(sym(f),findsym(sym(f),r2); root=r2-fx2/dfx; %第二步迭代 tol=abs(root-r1); else %第二种形式的两步迭代法 r1=root; fx1subs

5、(sym(f),findsym(sym(f),r1); dfx=subs(sym(fun),findsym(sym(fun),r1); r2=r1-fx1/dfx; %第一步迭代 fx2=subs(sym(f),findsym(sym(f),r2); root=r2-fx2*(r2-r1)/(2*fx2-fx1); %第二步迭代 tol=abs(root-r1); end endend（3）流程图（可放大）（4）举例说明两步迭代法求解非线性方程应用实例。利用两步迭代法求解方程logx-sin(x)+1=0在区间0.1,3的一个根。解：在MATLAB命令框中输入：r=TwoStep('(

6、log(x)-sin(x)+1)',0.1,3,1)；r=TwoStep('(log(x)-sin(x)+1)',0.1,3,2)；从上例可以看出两种迭代结果是一样的。二、编程解决以下科学计算问题。1、如下图所示为一个双电感并联单调谐网络，求回路的通频带B及满足回路阻抗大于50k的频率范围。28.2KRS R1 2 R2 3 1000PF0.75mH L1 0.25mH L2（1）算法说明先把回路变换为一个等效单电感谐振回路，把信号源的内阻Ra变为并接在该单电感回路上的等效内阻Rse，如下图所示。设      

7、  则       其它两支路的等效阻抗分别为      总阻抗是这3个支路阻抗的并联     由此求出频率范围。（2）流程图输入等效电路各电阻R以及感抗L,电容C根据R，L，C算出通频带B，谐振频率F与品质因数Q根据F，Q，B，求出等效阻抗Z1，Z2以及总阻抗Z。输出通频带B和频率范围F1利用FIND函数求出符合条件的频率范围F1结束开始（3）程序代码r1=2;r2=3;L1=0.75e-3;L2=0.25e-3;

8、C=1000e-12;rs=28200;L=L1+L2;r=r1+r2;rse=rs*(L/L1)2;%折算内阻f0=1/(2*pi*sqrt(C*L) %谐振频率Q0=sqrt(L/C)/r,r0=L/C/r; %空载（即不接信号源时）的回路Q0值re=r0*rse/(r0+rse), %折算内阻与回路电阻的并联Q=Q0*re/r0,B=f0/Q, %实际Q值和通带s=log10(f0);f=logspace(s-.1,s+.1,501); w=2*pi*f   %设定计算的频率范围及数组%等效单回路中两个电抗支路的阻抗z1e=r1+j*w*L;z2e=r2+1

9、./(j*w*C);ze=1./(1./z1e+1./z2e+1./rse); %等效单回路中三个支路的并联阻抗%画对数幅频特性曲线subplot(2,1,1),loglog(w,abs(ze),xlabel('频率'),ylabel('幅度'),title('幅频特性曲线'),grid%画相频特性曲线axis(min(w),max(w),0.9*min(abs(ze),1.1*max(abs(ze)subplot(2,1,2),semilogx(w,angle(ze)*180/pi)axis(min(w),max(w),-100,100),xl

10、abel('频率'),ylabel('相位'),title('相频特性曲线'),gridfh=w(find(abs(1./(1./z1e+1./z2e)>5e4)/2/pi; %幅频特性大于50k的频带fhmin=min(fh),fhmax=max(fh),（4）程序运行结果f0 =  1.5915e+005Q0 =   200rse =  4.0085e+004Q =   40.0853B =  3.9704e+003w

11、 = 1.0e+006fhmin =  1.5770e+005fhmax =  1.6063e+005即电路中的各个参数为：谐振频率f0159.15kHz空载品质因数Q0200等效信号源内阻rse50.133k考虑内阻后的品质因数Q40.0853通频带B3.97kHz回路阻抗大于50k的频率范围fhmin157.7kHzfhmax160.63kHz谐振频率附近的幅频和相频特征曲线如下图所示：从以上运行结果可知本题所求的通频带B为3.97kHz，满足回路阻抗大于50k的频率范围为（157.7，160.63）kHz。2、求压紧三次样条曲线，经过点（-3,2），（-2,0），（1,3），（4,1