jump to navigation

Combine matrices June 28, 2011

Posted by viboon in : Matlab , add a comment
horzcat(1, [2 3], [4 5 6]) %horizontally concatenates matrices
vertcat(1, [2; 3], [4; 5; 6]) %vertically concatenates matrices

This is actually equivalent to:

v1=1;
v2=[2; 3];
v3=[4; 5; 6];
v=[v1; v2; v3]

Tags: , ,

Related posts:

การหา inverse matrix ด้วย LU Factorization และ MATLAB September 25, 2010

Posted by viboon in : Mathematics, Matlab , add a comment

กำหนดให้ [A][x]=[B] โดยที่ [A] เป็น matrix ขนาด nxn และ [B] ขนาด nx1 โดยตัวอย่างนี้ เราต้องการหา matrix x ที่มีขนาด nx1

ปกติแล้วเรามักคุ้นเคยกับการหา A-1 ซึ่งก็หาได้จาก

A-1 = (1/det(A))*adj(A)T

หรือใน MATLAB จะใช้คำสั่ง

1
inv(A)

หรือ

1
A^(-1)

ก็ได้

แต่การหา x โดยการใช้วิธีนี้ มักมี error เกิดขึ้นอันเนื่องมาจากการปัดเศษ โดยเฉพาะในการคำนวณในลักษณะ numerical computation ซึ่ง error ที่เกิดขึ้นจะสะสม ทำให้คำตอบสุดท้ายผิดเพี้ยนไปได้ ซึ่งหลักการนี้สามารถเห็นได้จาก numerical methods ประเภท implicit ทั่วไป เช่น backward euler หรือ Crank-Nicolson เป็นต้น

วิธีที่จะนำเสนอคือ LU Factorization ซึ่งถือว่าเป็น numerical method ที่เก่าแก่วิธีการหนึ่ง เกิดจากการคิดค้นวิจัยในช่วงปี 1955 ถึง 1965 โดยชื่อวิธีการตั้งตาม C.F. Gauss

การทำ LU Factorization ซึ่งใช้หลักการ Gauss elimination จะแยกการคำนวณออกเป็น 2 ส่วนคือ Forward substitution และ Backward substitution ตามลำดับ

เริ่มจาก กำหนดให้ [A]=[L][U] โดยที่ [L] คือ lower matrix และ [U] คือ upper matrix ซึ่งเมื่อนำเข้าไปแทนในสมการตั้งต้นจะได้ [L][U][x]=[B]

สิ่งที่เราต้องการหาก็คือ [L] และ [U] ซึ่งหาได้จากการใช้ row operation จัดรูป matrix A ให้อยู่ในรูปของ upper matrix ก่อน
[A|I] -> [U|R]
โดยที่ [U]=[R][A]
จากสมการนี้ สังเกตได้ว่า [R] ก็คือ [L]-1 นั่นเอง หรือ [U]=[L]-1[A]

ในขั้นตอนถัดไปคือการทำ Forward substitution โดยกำหนดให้ [z]=[U][x] ดังนั้น เราจะได้ [L][z]=[B] ซึ่งเราสามารถทำการแก้สมการเพื่อหา [z] ได้โดย [z]=[L]-1[B] หรือก็คือ [z]=[R][B]

ต่อมาเป็นการทำ Backward substitution โดยแทน [z] เข้าไปในสมการก่อนหน้า เพื่อใช้ในการหา [x] โดย [U][x]=[z] ขั้นตอนนี้จะต้องหา inverse U ซึ่งเราสามารถหาได้ง่ายเพราะเป็น triangular matrix โดยที่ [x]=[U]-1[z]

นอกจากจะใช้ row operation ดังที่แสดงไว้ก่อนหน้านี้ การหา [U] และ [L] เรายังสามารถใช้การหาค่าใน matrix ได้โดยตรงหากเป็น matrix ขนาดไม่ใหญ่มาก เช่น

สำหรับการใช้ MATLAB เราสามารถหา [L] และ [U] ได้จาก

1
[L,U]=lu(A)

และหากต้องการหา [x] ก็ใช้คำสั่ง

1
x=inv(U)*inv(L)*B

นอกจากนี้ MATLAB ยังมี operator ที่คิดขึ้นมาเฉพาะ MATLAB เพื่อจัดกการกับ error ที่เกิดขึ้นในการหา inverse แบบ inv(A) คือ การใช้ backward slash “\”

หากพิจารณาปัญหาตามตัวอย่างที่กล่าวมาแล้ว [x] สามารถหาได้ด้วยคำสั่ง

1
x=A\B

error ที่เกิดจากการปัดเศษ หรือ roundoff error ที่เกิดขึ้นในขั้นตอนการหา inverse matrix แบบดั้งเดิมนั้นอาจดูไม่เป็นเรื่องสำคัญ หากเป็นการแก้ปัญหาที่ไม่ซับซ้อนมาก แต่สำหรับการคำนวณที่ต้องการความถูกต้องมากๆ อย่างเช่นการแก้ปัญหาด้วย numerical methods แบบต่างๆ เช่น finite difference method, finite element method หรือ finite volume method การสะสมของ roundoff error นำไปสู่การประมาณค่าสุดท้ายที่ผิดพลาดอย่างมากได้

รูปตัวอย่างจาก: http://www.alkires.com/103/6_lu.pdf


Tags: , , , , , ,

Related posts:

Checking the difference of matrices July 15, 2010

Posted by viboon in : Matlab , add a comment
1
2
3
4
A=rand(3,3)
B=A; A(1)=A(1)+5;
A == B
isequal(A,B)

Tags: ,

Related posts:

Shift the entire row/column of a matrix October 10, 2009

Posted by viboon in : Matlab , add a comment

This method can help in matrix operations without any loop or many steps of matrix expansion and delete. Row and column can be shifted together within a one-line command.

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
m=[1 2 3;4 5 6;7 8 9]
m =
     1     2     3
     4     5     6
     7     8     9
circshift(m, [1 0])
ans =
     7     8     9
     1     2     3
     4     5     6
circshift(m, [2 0])
ans =
     4     5     6
     7     8     9
     1     2     3
circshift(m, [0 1])
ans =
     3     1     2
     6     4     5
     9     7     8

Tags: ,

Related posts:

Rotate and transpose the matrix October 10, 2009

Posted by viboon in : Matlab , add a comment

See how they are different

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
m=[1 2 3;4 5 6;7 8 9]
m =
     1     2     3
     4     5     6
     7     8     9
m'
ans =
     1     4     7
     2     5     8
     3     6     9
rot90(m) %rotate the matrix by 90 degrees
ans =
     3     6     9
     2     5     8
     1     4     7

Tags: ,

Related posts:

Delete the entire row/column in the matrix October 10, 2009

Posted by viboon in : Matlab , add a comment
1
2
3
4
5
6
7
8
9
10
11
m=rand(4)
m =
    0.7094    0.6551    0.9597    0.7513
    0.7547    0.1626    0.3404    0.2551
    0.2760    0.1190    0.5853    0.5060
    0.6797    0.4984    0.2238    0.6991
m(3,:)=[]
m =
    0.7094    0.6551    0.9597    0.7513
    0.7547    0.1626    0.3404    0.2551
    0.6797    0.4984    0.2238    0.6991

Tags: ,

Related posts:

Find matrix size in Matlab September 30, 2009

Posted by viboon in : Matlab , add a comment
1
2
3
4
5
m=rand(3,4)
m =
    0.2769    0.8235    0.9502    0.3816
    0.0462    0.6948    0.0344    0.7655
    0.0971    0.3171    0.4387    0.7952
1
2
3
size(m)
ans =
     3     4
1
2
3
size(m,1)
ans =
     3
1
2
3
size(m,2)
ans =
     4
1
2
3
numel(m)
ans =
    12

Tags: ,

Related posts:

Useful matrix functions in Matlab August 21, 2009

Posted by viboon in : Matlab , add a comment

Diagonal matrix (with optional shifting from the diagonal) — 1D matrix as an input.

1
2
3
4
5
6
7
X = diag([1 2 3 4],2); disp(M);
     0     0     1     0     0     0
     0     0     0     2     0     0
     0     0     0     0     3     0
     0     0     0     0     0     4
     0     0     0     0     0     0
     0     0     0     0     0     0

Build a diagonal matrix — number, 1D and 2D matrices are applicable.

1
2
3
4
5
6
7
blkdiag(1, [2 2; 2 2])
 
ans =
 
     1     0     0
     0     2     2
     0     2     2

Repeat matrix

1
2
3
4
5
6
repmat([0 1], 2,3)
 
ans =
 
     0     1     0     1     0     1
     0     1     0     1     0     1

Tags: ,

Related posts: