Báo cáo bài tập lớn môn Phương pháp tính

Đại học Quốc gia Thành phố Hồ Chí Minh
Trường đại học Bách Khoa
Báo cáo bài tập lớn môn Phương pháp tính
Nhóm 14
TP.HCM, tháng 12, 2010
Thành viên:
Nguyễn Quốc Thái
Lê Việt Anh
Hồ Thiện Ân
Nguyễn Duy Anh
Trần Xuân An
Vũ Đức Anh
GVHD: thầy Nguyễn Hồng Lộc
Đề bài
Chứng minh: Ma trận vuông A đối xứng, xác định dương phân tích được dưới dạng [A]=[B][B]T với [B] là ma trận tam giác dưới.
Áp dụng: cho ma trận [A] như sau:
 với α là số thỏa điều kiện xác định dương.
Bài làm
Xét ma trận [A] là ma trận vuông cấp 3 thỏa hai điều kiện là [A] đối xứng và [A] xác định dương.
Giả sử phân tích được [A]=[B][B]T với [B] là ma trận tam giác dưới như sau:
Sau khi thực hiện phép nhân ta sẽ có:
Vế phải là một ma trận đối xứng. Cân bằng các phần tử của 2 ma trận ta có:
Vậy với trường hợp n=3 ta đã tìm được ma trận [B] thỏa mãn yêu cầu.
Tổng quát, với ma trận [A] vuông cấp n ta sẽ có:
 với 
Cân bằng với phần tử của ma trận [A] ta có:
 với 
Do ma trận [B] là ma trận tam giác dưới nên đối với cột thứ nhất ta có:
Đối với cột khác, rút bij ra khỏi tổng ta có:
Nếu i=j (phần tử trên đường chéo) thì:
 với j=2,3,4,,n
Và phần tử nằm ngoài đường chéo là:
 với 
Vậy ta đã tìm được ma trận [B] thỏa mãn yêu cầu.
Áp dụng với ma trận [A] đề cho:
Ta sẽ trình bày lời giải bằng chương trình Matlab.
Giải thuật gồm có 2 bước chính:
Nhập α từ bàn phím và kiểm tra xem α vừa nhập có thỏa điều kiện xác định dương hay không (điều kiện đối xứng đã thỏa với mọi α). Nếu thỏa chuyển qua bước tiếp theo. Nếu không thỏa, yêu cầu nhập lại α.
Khi α đã thỏa mãn điều kiện, ta bắt đầu tính [B] dựa vào công thức đã chứng minh ở câu 1, sau đó xuất kết quả [B] ra màn hình.
Để tạo 1 file chương trình Matlab, ta thực hiện như sau: chọn menu File à New à M-File. Sau đó tạo 2 file nhapa.m và ptcholeski.m với nội dung như bên dưới.
File nhapa.m
a=input('Nhap a: ');
A=[1,2,3,4,5;2,5,6,8,10;3,6,10,9,11;4,8,9,26,a;5,10,11,a,450];
disp('Ma tran A la:'),disp(A);
xdd=true;
for i=1:length(A)
 if (det(A(1:i,1:i))<=0)
 xdd=false;
 break;
 end;
end;
f=xdd;
if f==1
 disp('Ma tran A da thoa dieu kien xac dinh duong.');
end;
while f==0
 disp('Ma tran A khong xac dinh duong, vui long nhap lai a!!!');
 nhapa;
end
File ptcholeski.m
function [B] = c3choleski(N,A)
if nargin<2, error('Ham co toi thieu 2 doi so');end;
%Doan code kiem tra tinh doi xung cua ma tran A
for i=1:N
 for j=1:N
 if A(i,j)~=A(j,i)
 error('Ma tran khong doi xung.');
 end;
 end;
end;
%Doan code nay co the bo trong bai toan nay (do ma tran A da doi xung)
B=zeros(N);%tao 1 ma tran B la ma tran vuong cap N co cac phan tu la 0
B(1,1)=sqrt(A(1,1));
for i=2:N, B(i,1)=A(i,1)/B(1,1);end;
for k=2:N
 Ak=0;
 for j=1:k-1
 Ak=Ak+B(k,j)*B(k,j);
 end;
 B(k,k)=sqrt(A(k,k)-Ak);
 for i=k+1:N
 Ak=0;
 for j=1:k-1
 Ak=Ak+B(i,j)*B(k,j);
 end;
 B(i,k)=(A(i,k)-Ak)/B(k,k);
 end;
end;
disp('Ma tran B can tim la:'); disp(B);%xuat [B] ra man hinh
Ta lưu 2 file này vào thư mục chính của Matlab.
Cách chạy chương trình:
Ta tạo 1 file btlppt.m với nội dung như sau:
clc;% xoa man hinh hien tai cho trong
nhapa;
ptcholeski(5,A);%N=5
Sau đó, tại dấu nhắc lệnh của Matlab, ta gõ dòng lệnh btlppt; rồi gõ Enter
Khi đó, chương trình sẽ yêu cầu nhập vào α:
Nếu α không thỏa (trong hình α=9) thì màn hình sẽ như sau:
Ta phải nhập lại α cho tới khi nào α thỏa mãn.
Ví dụ, khi ta nhập lại α=25 thì sẽ thỏa mãn điều kiện và màn hình sẽ như sau: 
Nếu α thỏa điều kiện ngay từ lần nhập đầu (như trong hình α=50), ta sẽ có: 
Vậy từ chương trình, với α thỏa điều kiện, ta sẽ xác định được ma trận [B] cần tìm.