for i=1:n
p1=A(i,i);%差商矩阵对角线上的元素就是Newdon插值多项式的系数 for j=1:i-1
p1=conv(p1,[1 -x(j)]);%计算Newdon插值多项式的基项 end
p1=[zeros(1,n-i),p1];%向量相加,维数必须相同。把向量的元素补齐 p=p+p1; end
if nargin==3
f=polyval(p,x0);%计算插值多项式在x0处的值 else
f=poly2str(p,'x');%把插值多项式的向量形式转化为插值多项式的符号形式 end
5、基本Guass消去法求解线性方程组
function x=EqtsBasicGuass(A,b)
%基本Guass消去法求解线性方程组Ax=b %x=EqtsBasicGuass(A,b) %x:解向量,列向量 %A:线性方程组的矩阵 %b:列向量
%
%应用举例:
%A=[2 2 3;4 7 7;-2 4 5]; b=[3;1;-7]; %x=EqtsBasicGuass(A,b)
%检查输入参数 if size(A,1) ~= size(b,1) disp('输入参数有误!'); x=' '; return; end %(A|b) A=[A b];
%消去过程 n=size(A,1); l=zeros(n); for k=1:n-1 for i=k+1:n
l(i,k)=A(i,k)/A(k,k); end for i=k+1:n for j=k+1:n+1
A(i,j)=A(i,j)-l(i,k)*A(k,j);
end for j=1:k A(i,j)=0; end end end
%回代过程 x=zeros(n,1); x(n)=A(n,n+1)/A(n,n); for i=n-1:-1:1 y=0; for j=i+1:n y=y+A(i,j)*x(j); end
x(i)=(A(i,n+1)-y)/A(i,i); end return;
6、三角分解法求解线性方程组
function x=EqtsDoolittleLU(A,b) %Doolittle分解法求解线性方程组Ax=b %x=EqtsDoolittleLU(A,b)
%x:解向量,列向量 %A:线性方程组的矩阵 %b:列向量 %
%应用举例:
%A=[6 2 1 -1;2 4 1 0;1 1 4 -1;-1 0 -1 3];b=[6;-1;5;-5]; %x=EqtsDoolittleLU(A,b)
%检查输入参数 if size(A,1) ~= size(b,1) disp('输入参数有误!'); x=' '; return; end %分解
%把L和U的元素存储在A的相应位置上 n=length(b); for k=1:n for j=k:n z=0; for r=1:k-1
z=z+A(k,r)*A(r,j); end
A(k,j)=A(k,j)-z;
end for i=k+1:n z=0; for r=1:k-1
z=z+A(i,r)*A(r,k); end
A(i,k)=(A(i,k)-z)/A(k,k); end end %求解
x=zeros(size(b)); for i=1:n z=0; for k=1:i-1 z=z+A(i,k)*x(k); end x(i)=b(i)-z; end for i=n:-1:1 z=0; for k=i+1:n z=z+A(i,k)*x(k); end
x(i)=(x(i)-z)/A(i,i);