迭代公式:
)x(0初始向量x(k?1)?B0xk(?f)
function [x,n]=Jacobi_Solve(A,b,x0,dalt) % Jacobi跌代法解线性方程组P213 %[x,n]=Jacobi_Solve(A,b,x0,dalt) % A 方程组系数 % b 常数项(列向量) % x0 初始值,默认为0 % dalt 精度,默认为 10 % x 返回跌代结果 % n 返回跌代次数 e=1; i=0; r=size(b); a=b; if nargin<4 dalt=1e-8; end if nargin<3 x=zeros(r); else x=x0; end r=r(1); for t=1:r a(t)=A (t,t); A (t,t)=0; A (t,:)=A (t,:)/a(t); end b=b./a; while e>=dalt Y=b-A*x; e=max(abs(Y-x)); x=Y; i=i+1; end if nargout>1 n=i; end
11
例1 对不同初值用 Jacobi 迭代法解例 1 中方程组AX=b,比较结果.
>>x0=[1.5 2.0 3 1 1.5 1.8]';
>> [x,n]=Jacobi_Solve(A,f',x0) >> [x,n]=Jacobi_Solve(A,f')
(2) 超松弛跌代法解线性方程组
超松弛跌代法是Gauss-Seidel跌代法的一种加速,是解大型稀疏矩阵方程组的有效方法之一. 但需要选择好的加速因子(最佳松弛因子).
function [x,n]=SOR_Solve(A,b,w,x0,dalt) % 超松弛跌代法解线性方程组P213 %[x,n]=SOR_Solve(A,b,w,x0,dalt) % A 方程组系数 % b 常数项(列向量) % w 松弛因子 % x0 初始值,默认为0 % dalt 精度,默认为 10 % x 返回跌代结果 % n 返回跌代次数 r=size(b); a=b; if nargin<5 dalt=1e-8; end if nargin<4 x=zeros(r); else x=x0; end
r=r(1); m=0; e=1; for t=1:r a(t)=A (t,t); A (t,t)=0; A (t,:)=A (t,:)/a(t); end b=b./a; while e>dalt
12
e=0; for i=1:r t=x(i);
x(i)=(1-w)*x(i)+w*(b(i)-A (i,:)*x); t=abs(x(i)-t); if t>e e=t; end end m=m+1; end if nargout>1 n=m; end
例2 对不同松弛因子用SOR解例 1 中方程组A X =b,并与例 2 比较结果. 解:
>> [x,n]=SOR_Solve(A,f',1.42) >> [x,n]=SOR_Solve(A,f',1)
数值分析实验三
第四章 函数的数值逼近
一. 实验目的
1) Matlab中多项式的表示及多项式运算 2) 用Matlab实现拉格朗日及牛顿插值法 3) 用多项式插值法拟合数据
二. 实验内容
1. 多项式运算 (1) 多项式表示 对多项式
p?a0xn?a1xn?1?an?1x?an，用以下的行向量表示：
p??a0,a1,,an?1,an?。这样把多项式转化为向量运算。
① >> p=[1,-5,6,-33];
13
>> poly2sym(p) %符号表示多项式 ans =
x^3-5*x^2+6*x-33 ② >> a=[1,2,3;2,3,4;3,4,5];
>> p1=poly(a) %矩阵a的特征多项式 p1 =
1.0000 -9.0000 -6.0000 -0.0000 ③ >> root=[-5 -3+4i -3-4i];
>> p=poly(root) p =
1 11 55 125 >> poly2sym(p) ans =
x^3+11*x^2+55*x+125
(2) 多项式运算
① >> p=[1,11,55,125]; >> a=1:2:10 ; >> b=[1,2;2,1]; >> polyval (p,a) ans =
192 416 800 >> polyval(p,b) ans = 192 287 287 192
>> polyvalm(p,b) ans = 248 168 168 248
>> roots(p) ans =
-5.0000 -3.0000 + 4.0000i -3.0000 - 4.0000i
%生成以 root中元素为根的多项式 %建立多项式p
%向量a %方阵b
%计算a对应元素的多项式值 1392 2240
%矩阵b对应元素多项式值 %按矩阵多项式计算 %多项式p的所有根 14
② >> p=[2 -5 6 -1 9]; >> d=[3 -90 -18];
>> pd=conv(p,d) %多项式p,d相乘 pd =
6 -195 432 -453 9 -792 -162
>> p1=deconv(pd,d) %p1 为多项式pd 除 d(等于p) p1 =
2 -5 6 -1 9 >> p=[1 2 3 4 5 6]; >> q=[3 2 1];
>> [s,r]=deconv(p,q) %s 为商,r余项 s =
0.3333 0.4444 0.5926 0.7901 r =
0 0 0 0 2.8272 5.2099
>> Dp=polyder(p) %对p求导所得多项式 Dp =
5 8 9 8 5
2．拉格朗日插值法 (1) 算法实现
function [p]=Lag_polyfit(X,Y) %Matlab函数文件 Lag_polyfit.m %[p]=Lag_polyfit(X,Y) % X 拟合自变量 % Y 拟合函数值
% p 所得的拟合多项式系数 if size(X)~=size(Y) error('变量不匹配');
end %如果要拟合的函数值与自变量维数不一样,则退出报错 tic %开始记时
format long g %设置最合适的数字格式 r=size(Y); n=r(2); %n为要拟合的数据长度 p=zeros(1,n); %保存所得多项式系数
15