จงหาคำตอบของระบบสมการ - มหาวิทยาลัยเทคโนโลยีสุรนารี

Report
ปฏิบตั ิการคอมพิวเตอร์
เรื่ อง ระบบสมการหลายตัวแปร
อ.ดร.ชโลธร ธรรมแท้
สาขาวิชาวิศวกรรมเครื่ องกล มหาวิทยาลัยเทคโนโลยีสุรนารี
Linear Systems: GAUSS Elimination
ปั ญหา n สมการ - n ตัวแปร
หรื อเขียนในรูป [A]{x}={b}
GAUSS Elimination
•กำจัดไปข้ ำงหน้ ำ
•แปลง Matrix ให้ เป็ น
Matrix สำมเหลี่ยมบน
•จะเหลือสมกำรสุดท้ ำยเป็ น 1 สมกำร
1 ตัวแปร
•จะได้ คำตอบของตัวแปรสุดท้ ำย
•จำกนันแทนค่
้
ำกลับขึ ้นมำทีละสมกำร
Forward Elimination
แถวล่ำง = แถวล่ำง - (สปส ล่ำง / สปส บน )แถวบน
2
row 2  row 2  
 row 1
 4 
1
row 3  row 3    row 1
4
  1 .5 
row 3  row 3  
 row 2
 3 
Back Substitution
จำก Forward Elimination
ได้ Matrix สำมเหลี่ยมบน
แทนค่ำ
กลับ
หรื อ
MATLAB: GaussElim.m
% ขันตอน
้
Forward Elimination %
factor = A(i,j)/A(j,j);
% factor = สปส ล่ำง / สปส บน
A(i,j:n) = A(i,j:n) - factor*A(j,j:n); %แถวล่ำง-factor*แถวบน
b(i)= b(i) - factor*b(j);
%แถวล่ำง-factor*แถวบน
% ขันตอน
้
Back Substitution%
for k = n-1:-1:1
% loop Back sub
x(k) = (b(k) - A(k,k+1:n)*x(k+1:n))/A(k,k)
end
Download ได้ จาก http://eng.sut.ac.th/me/subject.html
ปัญหาของวิธี GAUSS Elimination
สมกำรในรูปต่อไปนี ้ จะไม่สำมำรถทำ Gauss Elimination ได้
•เนื่องจำก a11 เป็ น 0 ทำให้ ขนตอนกำรก
ั้
ำจัด row2 เกิดกำรหำรด้ วย 0
•แก้ ปัญหำโดยสลับ Row3 กับ Row1 แล้ วจึงเริ่มคำนวณ
•ระหว่ำงกำรคำนวณอำจเกิด สปส ของ column ซ้ ำยเป็ นศูนย์ ให้ ข้ำมกำร
ทำงำนในแถวนันไป
้ ใน MATLAB คือบรรทัด: if A(i,j)~=0
•เรี ยกว่ำกำรทำ Pivoting
ปัญหาของวิธี GAUSS Elimination
The exact solution is
ทำ Forward Elimination
0 . 0003 x1
 3 . 000 x 2  2 . 0001
 9999 x 2  6666
หำกปั ดเศษ x2 จะได้ คำตอบ x1 ต่ำงกันอย่ำงมำก ดังนี ้
สมการชุดเดียวกัน หากสลับแถว
ทำ Forward Elimination
1 . 0000 x1
 1 . 0000 x 2  1 . 0000
2 . 9997 x 2  1 . 9998
หำกปั ดเศษ x2 จะได้ คำตอบ x1 ต่ำงกันเล็กน้ อย : Partial Pivoting
แบบฝึ กหัด
จงหาคาตอบของระบบสมการ
1.
2.
Linear Systems: LU factorization
[ A ]  [ L ][ U ]
u 11  a11
u 12  a12
u 13  a13
l 21 u 11  a 21
l 21 u 12  a 22
l 21 u 13  a 23
l 31 u 11  a 31
l 31 u 12  l 32 u 22  a 32
l 31 u 13  l 32 u 23  u 33  a 32
เมื่อ [A] เป็ นค่ำที่ทรำบ แก้ สมกำรทัง้ 9 ทีละสมกำร จะได้ [L] และ [U]
LU factorization
[ A ]{ x }  {b }
[ L ][ U ]{ x }  {b }
ให้ [U]{x}={y}
แก้ [L]{y}={b}
ให้ [U]{x}={y}
1

l
 21
 l 31
 u 11

0

 0
ทำให้ [L]{y}={b}
0
1
l 32
u 12
u 22
0
0   y1   b1 
   
0  y 2    b2 

1   y 3   b3 
หำคำตอบของ {y} ได้ ง่ำยจำก
matrix สำมเหลี่ยมล่ำง
u 13   x1   y1 
   
u 23  x 2    y 2 

u 33   x 3   y 3 
หำคำตอบของ {x} ได้ ง่ำยจำก
matrix สำมเหลี่ยมบน
MATLAB function: lu()
[L,U]=lu(A)
สมกำร
MATLAB function: \
• คำสัง่ \ ใน MATLAB ใช้ ในกำรหำคำตอบของสมกำร
[A]{x}={b}
• ต้ องกำรหำคำตอบ [L]{y}={b} พิมพ์
>> y=L\b
• ต้ องกำรหำคำตอบ [U]{x}={y} พิมพ์
>> x=U\y
คาถาม: A\b ไห้ ผลลัพธ์ เป็ น อะไร
MATLAB function: inv()
• inv() ใช้ ในกำรหำ inverse matrix เช่น inv(A)
1
จำกคุณสมบัติ
[ A] [ A]  [ I ]
คูณเข้ ำทำงซ้ ำย
[ A ] [ A ]{ x }  [ A ] {b }
ทำให้ ได้ {x}
1
1
1
[ I ]{ x }  [ A ] {b }
• ดังนันหำค
้
ำตอบของสมกำร [A]{x}={b} ได้ จำกคำสัง่
>> x=inv(A)*b
Linear Systems: : Jacobi Iteration
สมการ
จัดรูป
เดาค่ าเริ่มต้ น x1,x2,x3 เช่ น = 0 แล้ วแทนค่ ากลับลงไปในสมการที่จัดรู ปแล้ ว
x1  [ 7 . 85  0 . 1( 0 )  0 . 2 ( 0 )] / 3  2 . 6167
Iteration 1:
x 2  [  19 . 3  0 . 1( 0 )  0 . 3 ( 0 )] / 7   2 . 7571
x 3  [ 71 . 4  0 . 3 ( 0 )  0 . 2 ( 0 )] / 10  7 . 140
Linear Systems: : Jacobi Iteration
แทนค่ า x1,x2,x3 จาก iteration 1 ในสมการที่จัดรู ปแล้ ว
x1  [ 7 . 85  0 . 1(  2 . 757 )  0 . 2 ( 7 . 14 )] / 3  2 . 8783
Iteration 2: x  [  19 . 3  0 . 1( 2 . 6167 )  0 . 3 ( 7 . 14 )] / 7   2 . 6824
2
x 3  [ 71 . 4  0 . 3 ( 2 . 6717 )  0 . 2 (  2 . 7571 )] / 10  7 . 1138
ทาซา้ แทนค่ า x1,x2,x3 จาก iteration 2 และทาไปเรื่ อยๆ
Iteration 3 : X1=2.904500, X2=-2.674905, X3=7.111217
Iteration 4 : X1=2.907117, X2=-2.674157, X3=7.110955
x1 , 2 , 3  x1 , 2 , 3
new
ES 
old
x
new
1, 2 , 3
 100 %
Jacobi.m
X1old=0; X2old=0; X3old=0; % กำหนดค่ำเริ่ มต้ น
for i=1:10
X1new = (7.85 + 0.1*X2old + 0.2*X3old) /3 ; %สมกำร 1
X2new = (-19.3 - 0.1*X1old + 0.3*X3old) /7; %สมกำร 2
X3new = (71.4 - 0.3*X1old + 0.2*X2old) /10; %สมกำร 3
X1old=X1new; X2old=X1new; X3old=X1new; %เก็บค่ำเก่ำ
fprintf ('iter=%d X1=%f X2=%f X3=%f \n', i,X1new,X2new,X3new);
end
Download ได้ จาก http://eng.sut.ac.th/me/subject.html
แบบฝึ กหัด
• จงหำผลเฉลยของสมกำร โดยใช้ วิธี Jacobi Iteration
• เปรียบเทียบผลที่ได้ กับวิธีอ่ นื ๆ
Nonlinear Systems
สมการไม่ เชิงเส้ นดังในบท root of equations สามารถมีหลายสมการ
หลายตัวแปร
คาตอบของสมการ
คือจุดตัดกันของ
เส้ นกราฟทัง้ 2
Nonlinear Systems: Jacobi Iteration
สามารถประยุกต์ ใช้ Jacobi Iteration ในการหาคาตอบของ
Nonlinear systems
Case1
Case2
แบบฝึ กหัด
• จงหำผลเฉลยของสมกำร โดยใช้ วิธี Jacobi Iteration
• เปรียบเทียบระหว่ างการจัดรูปสมการ Case1 กับ Case2 ให้
ผลลัพธ์ ต่างกันอย่ างไร
• กำหนดค่ำเริ่มต้ น X1=2, X2=3 ทำซ ้ำ 5 ครัง้

similar documents