数学建模社区-数学中国
标题:
高斯消元法解线性方程组
[打印本页]
作者:
2744557306
时间:
2023-12-31 17:28
标题:
高斯消元法解线性方程组
这段代码是用于解线性方程组的高斯消元法（Gaussian elimination）。以下是代码的主要步骤：
8 [2 e6 `/ n8 i2 U
) H$ k* Z3 @! d o3 S) V
1.初始化： 定义系数矩阵 a 和常数向量 b，以及一个排列矩阵 L，用于记录行的交换顺序。
1 ]& `8 ]3 d1 D. k$ `/ o# l3 ~4 r
p% L+ [( U3 W g
a=[1,-1,2,-1;2,-2,3,-3;1,1,1,0;1,-1,4,3]; % 系数矩阵 a
0 C' h$ t0 F5 Z0 @3 L% `
b=[-8,-20,-2,4]'; % 常数向量 b
+ }6 g: c0 s8 z0 V
L=[1,2,3,4]; % 排列矩阵 L
/ z" @4 [+ E( G# e) G6 G
n=length(b);
# t8 P* r- ^6 \+ r4 Z. G8 O- f
9 P& D1 F. G* o1 M! z0 N# Y
: p2 O: I4 c- i7 k& r
2.高斯消元： 通过一系列行变换将系数矩阵转化为上三角矩阵，并相应地更新常数向量。这里使用了列主元素法，即每次选取绝对值最大的元素所在的行与当前行进行交换。
0 v2 n* N/ V) L' ^2 |/ y3 F: X* x7 C
( I# X* |2 C. B7 V# o3 k
for k=1:n-1
) K: e4 E+ H' ~& T$ n
[p,q]=find(abs(a)==max(max(abs(a(k:n,k:n)))));
$ u- T% f/ v0 `
5 Z% b9 W9 T8 S) b* h/ H1 U
if(p~=k | q~=k)
; v m) N4 E, k( e+ ^' e ~: p
t=a(k,
;
9 ~2 z( B5 P7 C9 v
a(k,
=a(p,
;
' n$ U: {0 P1 v# z6 a! C: }! @2 e
a(p,
=t;
* U5 q+ `' P/ j: n
r=a(:,k);
3 R- a' x- Q6 E, t! _9 u7 L& h
a(:,k)=a(:,q);
0 A1 ~$ v9 {' ]
a(:,q)=r;
' P; @3 Q/ I _3 z0 R& T; g
t=L(k);
: ]. p# l* d0 J" v
L(k)=L(q);
' X$ g* M5 z4 G- ^8 o5 S
L(q)=t;
; \% b. J, e. p+ m! ?2 O. X2 Z7 i
u=b(k);
3 b3 N% w O% f5 |, o
b(k)=b(p);
% n' Z! J7 w/ B# Y# n S, S
b(p)=u;
: \- w, M' v; {3 @* C) v
end
8 T) ^8 Y" y# J. F5 J- }; ]
m(k+1:n,k)=a(k+1:n,k)./a(k,k);
$ a* L/ ~! ?4 ]
a(k+1:n,k:n)=a(k+1:n,k:n)-m(k+1:n,k)*a(k,k:n);
8 @! y2 u- N9 P/ D& ?3 K
b(k+1:n)=b(k+1:n)-m(k+1:n,k)*b(k);
3 P% t( w( l9 i5 o' L5 u) E, E9 z1 H
end
3 ]2 P" E% q- B# j' e
" _$ y! X' {2 k! U3 r8 E
; r( n' I! G6 V* T! R
3.回代： 通过回代过程求解方程组。从最后一行开始，逐步计算未知数的值。
$ u5 g# p6 A# S. E# [$ p
& l0 P2 C5 T7 b: o- l$ G
y(n)=b(n)/a(n,n);
1 a# P$ b* n* O9 k5 T! A( ]4 r
for i=n-1:-1:1
7 S8 b5 `1 V; F" `5 s/ K; n
sum=0;
0 _8 G8 R$ [1 S. b- w1 ~
for j=i+1:n
1 W) f- F0 a. ~
sum=sum+a(i,j)*y(j);
* |. {' O. O, R' A7 r. D
end
) w7 x0 k0 v' Y6 k5 B
y(i)=(b(i)-sum)/a(i,i);
0 L: Y$ Z* N3 r7 _
end
+ c8 s/ w3 X8 [6 ~8 N, y6 j% |/ u
$ o- s! p9 u; @1 ]+ l
7 X! i8 J! O3 K5 R [' R
4.输出结果： 将解存储在 x 中，并输出结果。
% E# U8 \0 b) R [% S3 n
( {/ B: q4 j" Q' j4 B
x(L(n))=y(n);
! B7 r; Q; |3 `$ b
x(L(1:n-1))=y(1:n-1);
1 W' l" |3 q# k. m& z4 B
jie=x'
6 Y3 d: r2 h8 A2 J. H/ _' }7 X
+ i: ]1 Y0 V" _2 B0 `
最后，解向量 jie 包含了线性方程组的解。请注意，这段代码在求解之前进行了列主元素的行交换，以提高数值稳定性。
* g' J. ^1 D4 R9 Z. a
# ?5 m% z: b+ N7 b3 D9 N
5 }4 w- q( v. S# [, P" V d* S3 g9 R+ H
欢迎光临 数学建模社区-数学中国 (http://www.madio.net/)
Powered by Discuz! X2.5
;=a(p,;=t;* U5 q+ `' P/ j: n