Tác giả: Nguyễn RR Thành Trung, Nguyễn Mạnh Quân
Thông thường, để đạt được độ phức tạp thuật toán như mong muốn, cách làm thường là tìm ra một thuật toán ban đầu làm cơ sở, rồi từ đó dùng các kỹ năng để giảm độ phức tạp của thuật toán. Trong bài viết này, tôi xin giới thiệu với bạn đọc một kỹ năng khá thông dụng: nhân ma trận.
Trước khi đọc bài viết này, nếu bạn chưa có khái niệm gì về ma trận, bạn có thể tham khảo định nghĩa về ma trận trong một tài liệu khác.
Trước hết, tôi xin nhắc lại tóm tắt khái niệm về phép nhân ma trận:
Để thực hiện phép nhân ma trận trên máy tính, ta có thể thực hiện thuật toán với độ phức tạp $\mathcal{O}(M * N * P)$ như sau:
for i:=1 to M do
for j:=1 to P do
begin
C[i,j]:=0;
for k:=1 to N do
C[i,j]:=C[i,j]+A[i,k] * B[k,j];
end;
Đối với phép nhân các ma trận vuông kích thước $N * N$, có thuật toán nhân ma trận Strassen với độ phức tạp $\mathcal{O}(N^{\log{7}})$ theo tư tưởng chia nhỏ ma trận (tương tự cách nhân nhanh 2 số lớn)). Tuy nhiên cài đặt rất phức tạp và trên thực tế với giá trị $N$ thường gặp, cách này không chạy nhanh hơn nhân ma trận thông thường $\mathcal{O}(N^3)$.
Cần chú ý thêm là phép nhân ma trận không có tính giao hoán (do có thể thực hiện nhân 2 ma trận $A$ kích thước $M * N$ và ma trận $B$ kích thước $N * P$ nhưng không thể thực hiện phép nhân $B * A$ nếu $P \ne M$). Tuy nhiên, nhân ma trận lại có tính kết hợp:
$(A * B) * C = A * (B * C)$
Chúng ta hãy cùng xem xét một ví dụ kinh điển nhất trong ứng dụng của phép nhân ma trận.
Bài toán
Dãy Fibonacci được định nghĩa như sau:
F(0) = 1
F(1) = 1
...
F(i) = F(i-1) + F(i-2), i >= 2
Yêu cầu: Cho $N$ $(N \le 10^9)$, tính $F(N)$.
Bạn có thể nộp bài thử ở VNOI - LATGACH4.
Phân tích
Hiển nhiên cách làm thông thường là tính lần lượt các $F(j)$. Tuy nhiên, cách làm này hoàn toàn không hiệu quả với $N$ lên đến $10^9$, và ta cần một cách tiếp cận khác.
Ta xét các lớp số:
Ta hình dung mỗi lớp là một ma trận 1x2. Tiếp đó, ta sẽ biến đổi từ lớp $i-1$ đến lớp $i$. Sau mỗi lần biến đổi như vậy, ta tính thêm được một giá trị $F(i+1)$. Để thực hiện phép biến đổi này, chú ý là các số ở lớp sau chỉ phụ thuộc vào lớp ngay trước nó theo các phép cộng, ta tìm được cách biến đổi bằng nhân ma trận:
Chắc hẳn đọc đến đây bạn đọc sẽ thắc mắc, làm thế nào để tìm được ma trận trên? Để tìm được ma trận này, ta làm như sau:
Ta có:
Bây giờ ta sẽ cần tìm cách tăng tốc việc tính $[0, 1; 1, 1]^T$ ( * ). Việc tính nhanh ( * ) cũng hoàn toàn tương tự việc ta tính $a^T$ với $a$ là số nguyên. Sau đây là đoạn code minh hoạ. Trong đoạn code này, để bạn đọc dễ hiểu, tôi bỏ qua yếu tố về tính toán số lớn, và thực hiện các phép tính với kiểu số 32-bit.
type
matrix=array[0..1,0..1] of longint;
const
a: matrix=((0,1),(1,1));
//Định nghĩa phép nhân 2 ma trận
operator * (a,b:matrix) c:matrix;
var
i,j,k:longint;
begin
fillchar(c,sizeof(c),0);
for i:=0 to 1 do
for j:=0 to 1 do
for k:=0 to 1 do
c[i,j]:=c[i,j]+a[i,k]*b[k,j];
end;
//Tính a^n
function power(n:longint):matrix;
var
temp:matrix;
begin
if n=1 then exit(a);
temp:=power(n div 2);
temp:=temp*temp;
if n mod 2=1 then temp:=temp*a;
exit(temp);
end;
Bây giờ chúng ta sẽ cùng xem xét một ví dụ tổng quát hơn của ví dụ 1.
Bài toán: SPOJ - SEQ
Cho số nguyên $N$ $(N \le 100)$ và 2 dãy số độ dài $N$: $a_1, a_2, …, a_N$; $b_1, b_2, …, b_N$. Dãy số $c$ được định nghĩa như sau:
Yêu cầu: Tính $c_k$ với $k \le 10^9$.
Phân tích
Cũng như trong ví dụ 1, ta xét các lớp số:
Ta cũng sẽ áp dụng phép nhân ma trận để biến đổi từ lớp $i$ sang lớp $i+1$ như sau:
Để xây dựng ma trận vuông như trên, ta thực hiện tương tự như trong ví dụ trước: Phân tích $a_{i+1}$ đến $a_{i+N}$ dưới dạng $a_i, …, a_{i+N-1}$:
Từ đó, ta thu được cách làm như trong ví dụ 1. Cài đặt cụ thể xin nhường lại cho bạn đọc.
Chú ý rằng ta hoàn toàn có thể thay thế phép nhân và phép cộng trong định nghĩa phép nhân ma trận, chỉ cần đảm bảo giữ nguyên tính chất kết hợp. Cụ thể hơn, thay vì $C(i,j) = \sum{A(i,k) * B(k,j)}$, ta có thể định nghĩa phép nhân ma trận: $C(i,j) = min(A(i,k) + B(k,j))$. Từ đó, ta có thể thu được một lớp các bài toán khác.
Sau đây là một ví dụ minh hoạ cho nhóm các bài toán này
Bài toán
Cho đồ thị có hướng $N$ đỉnh $(N \le 100)$. Tính ma trận $C(k)$ kích thước $N * N$, với $C(k) [i,j]$ là độ dài đường đi ngắn nhất từ $i$ đến $j$ đi qua đúng $k$ cạnh
Phân tích
Xét ma trận $A$ là ma trận kề của đồ thị đã cho. Ta có:
Như vậy, nếu ta thay phép nhân và phép cộng trong việc nhân ma trận thông thường lần lượt bởi phép cộng và phép lấy min, ta thu được một phép ”nhân ma trận” mới, tạm dùng ký hiệu x
, thì:
C1 = A
C2 = C1 x C1 = A x C1
C3 = C1 x C2 = A x C2
C4 = C1 x C3 = A x C3
...
Ck = C1 x C(k-1) = A x C(k-1)
Do đó, $C(k) = A^k$
Như vậy, bài toán được đưa về bài toán tính lũy thừa của một ma trận, ta hoàn toàn có thể giải tương tự các ví dụ trước. Cài đặt phép nhân ma trận mới này hoàn toàn không phức tạp hơn cài đặt phép nhân ma trận thông thường. Việc cài đặt xin nhường lại cho bạn đọc.
Bài toán
Người ta mới tìm ra một loại vi khuẩn mới. Chúng sống thành $N$ bầy $(N \le 100)$, đánh số từ 0 đến $N-1$. Ban đầu, mỗi bầy này chỉ có một con vi khuẩn. Tuy nhiên, mỗi giây, số lượng vi khuẩn trong các bầy lại có sự thay đổi. Ví dụ:
Các thay đổi này tuân theo một số quy luật cho trước. Tại mỗi giây chỉ xảy ra đúng một quy luật. Các quy luật này được thực hiện nối tiếp nhau và theo chu kỳ. Có nghĩa là, nếu đánh số các quy luật từ 0 đến $M-1$, tại giây thứ $S$ thì quy luật được áp dụng sẽ là $(S-1) \space mod \space M$ $(M \le 1000)$
Nhiệm vụ của bạn là tìm xem, với một bộ các quy luật cho trước, sau $T$ đơn vị thời gian $(T \le 10^{18})$, mỗi bầy có bao nhiêu vi khuẩn.
Các loại quy luật có thể có:
A i 0
: Tất cả các vi khuẩn thuộc bầy $i$ chết.B i k
: Số vi khuẩn trong bầy $i$ tăng lên $k$ lần.C i j
: số vi khuẩn bầy thứ $i$ tăng lên một số lượng bằng với số vi khuẩn bầy $j$.D i j
: Các vi khuẩn thuộc bầy $j$ di chuyển toàn bộ sang bầy $i$.E i j
: Các vi khuẩn thuộc bầy $i$ và bầy $j$ đổi vị trí cho nhau.F 0 0
: Vị trí các vi khuẩn di chuyển trên vòng tròn.Phân tích
Cách làm đơn giản nhất là chúng ta mô phỏng lại số lượng vi khuẩn trong mỗi bầy qua từng đơn vị thời gian. Cách làm này có độ phức tạp $\mathcal{O}(T * N * k)$ với $\mathcal{O}(k)$ là độ phức tạp cho xử lý số lớn. Cách này không thể chạy được với $T$ lớn.
Ta hình dung số lượng vi khuẩn trong mỗi bầy trong một đơn vị thời gian là một dãy số. Như vậy, mỗi quy luật cho trước thực chất là một phép biến đổi từ một dãy số thành một dãy số mới, và ta hoàn toàn có thể thực hiện biến đổi này bằng một phép nhân ma trận.
Cụ thể hơn, ta coi số lượng vi khuẩn trong $N$ bầy tại một thời điểm xác định là một ma trận $1 * N$, và mỗi phép biến đổi là một ma trận $N * N$. Khi áp dụng mỗi phép biến đổi, ta nhân hai ma trận nói trên với nhau.
Bây giờ, xét trường hợp $N = 4$, tôi xin lần lượt mô tả các ma trận tương ứng với các phép biến đổi:
Biến đổi: A 2 0
1 0 0 0
0 1 0 0
0 0 0 0
0 0 0 1
Biến đổi: B 2 6
1 0 0 0
0 1 0 0
0 0 6 0
0 0 0 1
Biến đổi: C 1 3
1 0 0 0
0 1 0 0
0 0 1 0
0 1 0 1
Biến đổi: D 1 3
1 0 0 0
0 1 0 0
0 0 1 0
0 1 0 0
Biến đổi: E 1 3
1 0 0 0
0 0 0 1
0 0 1 0
0 1 0 0
Biến đổi: F 0 0
0 1 0 0
0 0 1 0
0 0 0 1
1 0 0 0
Cũng như các bài toán trước, ta sẽ cố gắng áp dụng việc tính toán lũy thừa, kết hợp với phép nhân ma trận để giảm độ phức tạp từ $T$ xuống $\log{T}$. Tuy nhiên, có thể thấy việc sử dụng phép lũy thừa trong bài toán này phần nào phức tạp hơn bởi các ma trận được cho không giống nhau. Để giải quyết vấn đề này, ta làm như sau:
Gọi $X_1, X_2, …, X_m$ là các ma trận tương ứng với các phép biến đổi được cho.
Đặt $X = X_1 * X_2 * … * X_m$.
Đặt $S = [1, 1, …, 1]$ (dãy số lượng vi khuẩn tại thời điểm đầu tiên).
Như vậy, $Y = S * X^t * X_1 * X_2 * … * X_r$ là ma trận thể hiện số lượng vi khuẩn tại thời điểm $M * t + r$.
Như vậy, thuật toán đến đây đã rõ. Ta phân tích $T = M * t + r$, nhờ đó, ta có thể giải quyết bài toán trong $\mathcal{O}(N^3 * M)$ cho bước tính ma trận $X$ và $\mathcal{O}(N^3 * (\log{T/M} + M)$ cho bước tính $Y$. Bài toán được giải quyết.
Trong phần Mở Đầu ta đã có thuật toán nhân hai ma trận $A$ kích cỡ $M * N$ và $B$ kích cỡ $N * P$ cần độ phức tạp $O(M * N * P)$. Giả sử ta có thêm ma trận $C$ có kích cỡ $P * Q$ và ta cần tính tích $A * B * C$. Xét hai cách thực hiện phép nhân này:
Như vậy là hai cách thực hiện khác nhau cần hai độ phức tạp khác nhau. Chọn $M = N = 500, P = 1000, Q = 2$, cách 1 sẽ cần tới $500 * 1000 * (500 + 2) = 251 * 10^6$ phép tính, trong khi cách 2 chỉ cần $500 * 2 * (500 + 1000) = 1.5 * 10^6$ phép tính, nghĩa là cách 1 chậm hơn cách 2 tới gần 200 lần.
Khi độ dài của dãy ma trận tăng lên, sự khác biệt có thể còn lớn hơn nữa. Ví dụ trên đã cho thấy rằng trong một số trường hợp thứ tự thực hiện phép nhân ma trận có ý nghĩa rất lớn đối với việc tìm lời giải của các bài toán.
Trong thực tế, bài toán xác định thứ tự nhân ma trận hiệu quả nhất là một bài toán rất phổ biến, bạn có thể tìm đọc chi tiết thêm ở Mục 3.5 Phép Nhân Tổ Hợp dãy Ma Trận trong ebook của thầy Lê Minh Hoàng.
Giải thuật Freivalds là một ví dụ điển hình về việc áp dụng thứ tự thực hiện phép nhân ma trận để giảm độ phức tạp tính toán của phép nhân một dãy ma trận. Bài toán đặt ra là cho ba ma trận vuông $A, B, C$ có kích cỡ $N * N$ với $N \le 1000$. Ta cần kiểm tra xem $C$ có phải là tích của $A$ và $B$, nói cách khác ta cần kiểm tra $A*B = C$ có phải là mệnh đề đúng hay không (đây chính là bài VMATRIX - VNOI Marathon 2014).
Phân tích
Cách làm thông thường là nhân trực tiếp hai ma trận $A, B$ rồi so sánh kết quả với $C$. Như phân tích trong phần Mở Đầu độ phức tạp của cách làm này là $O(N^3)$, với $N = 1000$ thì cách làm này không đủ nhanh. Giải thuật Freivalds thực hiện việc kiểm tra thông qua thuật toán xác suất kiểu Monte Carlo với $k$ lần thử cho xác suất kết luận sai là xấp xỉ $1 / 2^k$, mỗi lần thử có độ phức tạp $O(N^2)$. Các bước cơ bản của một phép thử Freivalds như sau:
True
nếu $P$ chỉ gồm phần tử $0$ (bằng với vector $0$) và False
nếu ngược lại.Ta thực hiện $k$ lần thử, nếu gặp phép thử trả về False
thì ta kết luận là $A * B \neq C$. Ngược lại nếu sau $k$ phép thử mà luôn thấy True
thì ta kết luận $A * B = C$. Vì xác suất lỗi giảm theo hàm mũ của $k$ nên thông thường chỉ cần chọn $k$ vừa đủ là sẽ thu được xác suất đúng rất cao ($k = 5$ với bài VMATRIX ở trên). Một nhận xét quan trọng khác là cận trên của đánh giá xác suất kiểm tra lỗi không phụ thuộc vào kích cỡ $N$ của ma trận được cho mà chỉ phụ thuộc vào số lần thực hiện phép thử.
Xét bước thứ 2, ta thấy rằng phép thử Freivalds chỉ có ý nghĩa nếu như ta có thể thực hiện phép nhân $A * B * v$ trong thời gian $O(N^2)$ (vì phép nhân $C * v$ đã đạt sẵn $O(N^2)$ rồi). Thay vì thực hiện tuần tự từ trái qua phải sẽ cần $O(N^3)$, ta thực hiện theo thứ tự $A * (B * v)$. Vì kết quả của phép nhân $B$ và $v$ là một ma trận $N * 1$ nên độ phức tạp tổng cộng sẽ là $O(N^2)$. Trên tất cả các phép thử, độ phức tạp là $O(k * N^2)$.