Matlab中Gauss全主元消元法求解线性方程组的程序

Matlab中Gauss全主元消元法求解线性方程组的程序

function [x,qa]=gaussq(a,b) %输出想x为解，qa为全主元变换后的a矩阵d=[a b];

RA=rank(a);RD=rank(a);L=length(b);n=size(a);pos=1:n(1);

if RA~=RD

fprintf('无解')

else

if RA~=L

fprintf('有无数多个解')

else

fprintf('有唯一解')

for q=1:n

big=max(max(abs(a(q:n,q:n))));

for r=q:n

for t=q:n

if big==abs(a(r,t))

zhuh=r;

zhul=t;

end

end

end

p=a(q,:); %换主行

a(q,:)=a(zhuh,:);

a(zhuh,:)=p;

bb=b(q);

b(q)=b(zhuh);

b(zhuh)=bb;

p=a(:,q); %换主列

a(:,q)=a(:,zhul);

a(:,zhul)=p;

p=pos(q); %记录由于换主列而造成的解的位置的变化 pos(q)=pos(zhul);

pos(zhul)=p;

end

c=[a b];

for j=1:L-1 %化为上三角阵

for i=(j+1):L

m=c(i,j)/c(j,j);

c(i,:)=c(i,:)-c(j,:)*m;

end

end

x(L,1)=c(L,L+1)/c(L,L);

for k=L-1:-1:1 %求解x

x(k,1)=(c(k,L+1)-c(k,k+1:L)*x(k+1:L))/c(k,k);

end

y=[1:n(1)]'; %交换被列调换时打乱的解的位置

for w=1:n

for v=1:n

if (pos(v)==w)

y(w)=x(v);

end

end

end

x=y;

end

end qa=a;