採用GAUSS列主消元法求解線性方程組(MATLAB)
阿新 • • 發佈:2019-02-12
%%求解任意線性方程組的解
clc;
clear all;
format long e
disp('線性方程組求解,請輸入引數');
n=input('維數n=');
A=input('矩陣A=');
b=input('右端項b=');
eps=input('控制精度eps=');
b=b'; %%變為列向量
A=[A b]; %%矩陣增廣
for k=1:n-1
B=A(k:n,k); %%先將第k列可能作為主元的元素取出方至矩陣B
P=max(abs(B)); %%選主元P
if(P<eps) %%控制小主元
disp('無解');
break;
else
u=find((abs(B))==P); %%計算主元所在行相對與k行的位置
if(u~=1)
A([k,u],:)=A([u,k],:); %%換行
end
m=A(k+1:n,k)/A(k,k); %%求出各行行乘數並放至矩陣m
for i=1:length(m)
A(k+i,k:n+1)=A(k+i,k:n+1)-m(i)*A(k,k:n+1); %%消元按行進行
end
end
end
if A(n,n)==0
disp('無解'); %%若矩陣A不滿秩,則無解
else
x(n)=A(n,n+1)/A(n,n); %%由最後一行首先求出方程組的第一個解x(n)
for i=n-1:-1:1 %%計算第i個解x(i)
for j=1:1:n-i %%利用回代思想
A(i,n+1)=A(i,n+1)-A(i,i+j)*x(i+j); %%減去已知部分
end
x(i)=A(i,n+1)/A(i,i);
end
end
disp('方程組的解');
x=x' %%輸出方程組的解
clc;
clear all;
format long e
disp('線性方程組求解,請輸入引數');
n=input('維數n=');
A=input('矩陣A=');
b=input('右端項b=');
eps=input('控制精度eps=');
b=b'; %%變為列向量
A=[A b]; %%矩陣增廣
for k=1:n-1
B=A(k:n,k); %%先將第k列可能作為主元的元素取出方至矩陣B
P=max(abs(B)); %%選主元P
if(P<eps) %%控制小主元
disp('無解');
break;
else
u=find((abs(B))==P); %%計算主元所在行相對與k行的位置
if(u~=1)
A([k,u],:)=A([u,k],:); %%換行
end
m=A(k+1:n,k)/A(k,k); %%求出各行行乘數並放至矩陣m
for i=1:length(m)
A(k+i,k:n+1)=A(k+i,k:n+1)-m(i)*A(k,k:n+1); %%消元按行進行
end
end
end
if A(n,n)==0
disp('無解'); %%若矩陣A不滿秩,則無解
else
x(n)=A(n,n+1)/A(n,n); %%由最後一行首先求出方程組的第一個解x(n)
for i=n-1:-1:1 %%計算第i個解x(i)
for j=1:1:n-i %%利用回代思想
A(i,n+1)=A(i,n+1)-A(i,i+j)*x(i+j); %%減去已知部分
end
x(i)=A(i,n+1)/A(i,i);
end
end
disp('方程組的解');
x=x' %%輸出方程組的解