Algoritma Strassen
Algoritma Strassen dalam matematika, khususnya aljabar linear adalah sebuah algoritma yang dinamakan oleh Volker Strassen yang merupakan sebuah algoritma yang digunakan untuk perkalian matriks yang secara asimtot lebih cepat dari pada algoritma perkalian matriks standar dan sangat berguna dalam penggunaanya untuk matriks yang berukuran besar.
Sejarah
Volker Strassen mempublikasikan algoritma Strassen tahun 1969. Meskipun algoritma ini hanya sedikit lebih cepat daripada algoritma standar untuk perkalian matriks, dialah yang pertama menjelaskan bahwa eliminasi Gauss adalah tidak optimal. Dalam tulisannya, dia memulai penelitian untuk melengkapi algoritma-algoritma yang lebih cepat seperti algoritma Winograd dari Shmuel Winograd pada 1980, dan yang lebih kompleks algoritma Coppersmith-Winograd dipublikasikan pada 1987.
Algoritma
Misalkan A, B dua matriks persegi pada ring R. Kita ingin menghitung produk matriks C sebagai
Jika matriks A, B bukan bertipe 2n x 2n kita isi baris-baris dan kolom-kolom yang kosong dengan nol.
Kita partisi A, B dan C kedalam matriks blok yang berukuran sama.
dengan
lalu
Dengan konstruksi ini kita tidak mengurangi jumlah dari perkalian-perkalian. Kita masih memerlukan 8 perkalian-perkalian untuk menghitung matriks-matriks Ci,j , dengan jumlah perkalian yang sama kita perlukan ketika menggunakan matriks perkalian standar.
Sekarang sampai pada bagian terpenting. Kita tetapkan matriks baru
Yang kemudian digunakan untuk mengekspresikan Ci,j dalam bentuk Mk. Karena kita telah mendefenisikan Mk kita bisa mengeliminasi satu perkalian matriks dan mengurangi jumlah perkalian-perkalian menjadi 7 (satu perkalian matriks untuk tiap Mk) dan ekspresi Ci,j sebagai
Kita iterasikan bagian diatas ke-n kali proses sampai submatriks-submatriks menjadi angka-angka.
Algoritma Strassen pada penerapannya mengubah metode standar dari perkalian matriks agar submatriks-submatriks yang cukup kecil menjadi lebih efisien. Fakta-fakta agar algoritma Strassen lebih efisien bergantung pada implementasi khusus dan hardware.
Analisi Numerik
Perkalian matriks standar melakukan
perkalian-perkalian dari elemen-elemen dalam ring R. Kita anggap penjumlahan-penjumlahan diperlukan karena bergantung pada R, yang bisa jauh lebih cepat daripada perkalian-perkalian dalam implementasi pada komputer terutama jika ukuran dari entri matriks melebihi ukuran kata dari mesin.
Dengan algoritma Strassen kita bisa mengurangi jumlah perkalian-perkalian
- .
Pengurangan dalam jumlah perkalian bagaimanapun akan sampai saat pilihan dari sedikit pengurangan kestabilan numerik.
Contoh program sederhana pada Matlab
function c = strass(a,b) nmin = 2; %misalkan matriks a dan b berukuran 2 x 2 [n,n] = size(a); if n <= nmin; c = a*b; else %entri matriks a dan b berukuran n x n; n=2^k; k=2,3,... %misalkan entri matriks a dan b berukuran n=2^2 atau 4 x 4 a11=a(1:2,1:2); a12=a(1:2,3:4); a21=a(3:4,1:2); a22=a(3:4,3:4); b11=b(1:2,1:2); b12=b(1:2,3:4); b21=b(3:4,1:2); b22=b(3:4,3:4); p1 = (a11+a22)*(b11+b22); p2 = (a21+a22)*b11; p3 = a11*(b12-b22); p4 = a22*(b21-b11); p5 = (a11+a12)*b22; p6 = (a21-a11)*(b11+b12); p7 = (a12-a22)*(b21+b22); c = [p1+p4-p5+p7 p3+p5; p2+p4 p1-p2+p3+p6]; end
Catatan: program diatas hanya untuk matriks berukuran 1x1, 2x2, 4x4. Untuk matriks yang berukuran lebih besar, masih diperlukan penyempurnaan. Agar programnya bisa berjalan.
Contoh program dalam bahasa C++
Program dalam bahasa C++ di bawah ini telah teruji dan memberikan hasil yang diinginkan. Untuk mempermudah pemahaman, program dipecah menjadi fungsi dan prosedur yang lebih kecil. Program ini dikompilasi di lingkungan linux Ubuntu 9.04 (Intrepid Ibex) dengan compiler g++. Kalau mau dickompilasi di lingkungan windows, silahkan hilangkan baris gettimeofday, fungsi ini adalah fungsi yang khusus di linux. Berikut ini adalah source codenya.
Untuk membandingkan algoritma strassen dan algoritma perkalian matrix konvensional, beberapa algoritma perkalian matrix secara konvensional (dalam istilah lain, algoritma Brute Force) juga disertakan.
/* Zulfikar Hakim (NIM 13507085) * Program Studi Teknik Informatika * Sekolah Teknik Elektro dan Informatika * Institut Teknologi Bandung * zulfikar.hakim@gmail.com * * * Tugas Kuliah Strategi Algoritma * IF3051 */
#include <iostream> #include <math.h> #include <cstdlib> #include <sys/time.h> #include <vector>
using namespace std;
static int n;
int randomizeData() { //membangkitkan sebuah elemen matrix secara random //mengembalikan nilai random untuk digunakan di elemen matrix. return rand()%3; }
void initMatrix(vector <vector <int> > *MatrixInput) /* Menginisiasi setiap elemen dari matrix dengan * hasil random data dari fungsi randomizeData * */ { for(int i=0;i<n;i++) { for(int j=0;j<n;j++) { (*MatrixInput)[i][j]=randomizeData(); } } }
void printMatrix(vector < vector<int> > MatrixInput,int panjangMatrix) /* menampilkan matrix ke layar jika panjang matrix * kurang dari 25 */ { if (panjangMatrix<=25) { for(int i=0;i<panjangMatrix;i++) { for(int j=0;j<panjangMatrix;j++) { cout<<MatrixInput[i][j]<<" "; if(MatrixInput[i][j]<10) { cout<<" "; } } cout<<endl; } } }
void samakanMatrix(vector <vector <int> > matrixBrute, vector <vector <int> > *matrixStrassen) /* Menyamakan isi dari sebuah matrixBruteForce dengan matrixStrassen * Hal ini diperlukan karen matrixStrassen memiliki banyak elemen N x N * di mana N adalah 1 atau bilangan kelipatan 2*/ { for(int i=0;i<n;i++) { for(int j=0;j<n;j++) { (*matrixStrassen)[i][j]=matrixBrute[i][j]; } } }
void bruteMultiplication(vector< vector<int> > matrix1, vector< vector<int> > matrix2, vector< vector<int> > *matrixHasil, int *kali, int *jumlah, struct timeval *output) /* Melakuka perkalian antar matrix dengan cara konvensional, * dan menghasilkan jumlah perkalian, jumlah operasi penjumlahan, * danjuga waktu yang dibutuhkan dalam timeval.*/ { timeval mulai; timeval selesai; gettimeofday(&mulai,NULL); for(int i=0;i<n;i++) { for(int j=0;j<n;j++) { (*matrixHasil)[i][j]=0; for(int k=0;k<n;k++) { (*matrixHasil)[i][j]+=matrix1[i][k]*matrix2[k][j];; (*kali)++; (*jumlah)++; } } } gettimeofday(&selesai,NULL); (*output).tv_sec=selesai.tv_sec-mulai.tv_sec; (*output).tv_usec=selesai.tv_usec-mulai.tv_usec; }
int isKelipatanDua(int inputBilangan, int *pnjgSeharusnya) /* Memeriksa apakah inputBilangan adalah bilangan kelipatan dua * Jika ya mengembalikan 1 dan memberikan pnjgSeharusnya yang merupakan * bilangan kelipatan dua terdekat dari inpuBilangan * Jika tidak mengembalikan 0*/ { int iterator; iterator=1; while(iterator<inputBilangan) { iterator=iterator*2; } if (iterator==inputBilangan) { *pnjgSeharusnya=iterator; return 1; } else /*if iterator > inputBilangan*/ { *pnjgSeharusnya=iterator; return 0; } }
vector <vector <int> > strassenAddition(vector <vector <int> > matrix1, vector <vector <int> > matrix2, int *operasiJumlah) /* Melakukan perkalian dua buah matrix dengan algoritma strassen dan * dilakukan dengan cara rekursif*/ { unsigned int ukuranMatrix=matrix1.size(); vector <vector <int> > retval(ukuranMatrix,vector<int>(ukuranMatrix)); for(unsigned int i=0;i<ukuranMatrix;i++) { for(unsigned int j=0;j<ukuranMatrix;j++) { retval[i][j]=matrix1[i][j]+matrix2[i][j]; (*operasiJumlah)++; } } return retval; }
vector <vector <int> > strassenSubtraction(vector <vector <int> > matrix1, vector <vector <int> > matrix2, int *operasiKurang) /* Melakukan pengurangan dua buah matrix. Fungsi ini digunakan dalam fungsi * strassenMultiplication */ { unsigned int ukuranMatrix=matrix1.size(); vector <vector <int> > retval(ukuranMatrix,vector<int>(ukuranMatrix)); for(unsigned int i=0;i<ukuranMatrix;i++) { for(unsigned int j=0;j<ukuranMatrix;j++) { retval[i][j]=matrix1[i][j]-matrix2[i][j]; (*operasiKurang)++; } } return retval; }
vector <vector <int> > getMatrix11(vector <vector <int> > matrixInput,unsigned int ukuranInput) /* Mendapatkan matrix kiri atas dari sebuah matrixBesar yang harus dipartisi*/ { vector <vector <int> > retval(ukuranInput,(vector<int>(ukuranInput))); for(unsigned int i=0;i<ukuranInput;i++) { for(unsigned int j=0;j<ukuranInput;j++) { retval[i][j]=matrixInput[i][j]; } } return retval; }
vector<vector <int> > getMatrix12(vector <vector <int> > matrixInput,unsigned int ukuranInput) /* Mendapatkan matrix kanan atas dari matrix besar yang harus dipartisi*/ { vector <vector <int> > retval(ukuranInput,(vector<int>(ukuranInput))); for(unsigned int i=0;i<ukuranInput;i++) { for(unsigned int j=ukuranInput;j<ukuranInput*2;j++) { retval[i][j-ukuranInput]=matrixInput[i][j]; } } return retval; }
vector <vector <int> > getMatrix21(vector <vector <int> > matrixInput,unsigned int ukuranInput) /* Mendapatkan matrix kiri bawah dari matrix besar yang harus dipartisi*/ { vector <vector <int> > retval(ukuranInput,(vector<int>(ukuranInput))); for(unsigned int i=ukuranInput;i<matrixInput.size();i++) { for(unsigned int j=0;j<ukuranInput;j++) { retval[i-ukuranInput][j]=matrixInput[i][j]; } } return retval; }
vector <vector <int> > getMatrix22(vector <vector <int> > matrixInput,unsigned int ukuranInput) /* Mendapatkan matrix kiri bawah dari matrix besar yang harus dipartisi*/ { vector <vector <int> > retval(ukuranInput,(vector<int>(ukuranInput))); for(unsigned int i=ukuranInput;i<ukuranInput*2;i++) { for(unsigned int j=ukuranInput;j<ukuranInput*2;j++) { retval[i-ukuranInput][j-ukuranInput]=matrixInput[i][j]; } } return retval; }
int get11(vector <vector <int> > matrixInput) /* Mendapatkan elemen matrix dengan index 0,0 jika matrix berukuran 2x2*/ { return matrixInput[0][0]; }
int get12(vector <vector <int> > matrixInput) /* Mendapatkan elemen matrix dengan index 0,1 jika matrix berukuran 2x2*/ { return matrixInput[0][1]; }
int get21(vector <vector <int> > matrixInput) /* Mendapatkan elemen matrix dengan index 1,0 jika matrix berukuran 2x2*/ { return matrixInput[1][0]; }
int get22(vector <vector <int> > matrixInput) /* Mendapatkan elemen matrix dengan index 1,1 jika matrix berukuran 2x2*/ { return matrixInput[1][1]; }
vector <vector <int> > JoinMatrix(vector <vector <int> > C11, vector <vector <int> > C12, vector <vector <int> > C21, vector <vector <int> > C22) /* Menggabungkan 4 buah matrix menjadi matrix yang lebih besar*/ { int retvalSize=C11.size()*2; int subMatrixSize=C11.size(); vector <vector <int> > retval(retvalSize, vector<int>(retvalSize)); //mengopikan C11 for(int i=0;i<subMatrixSize;i++) { for(int j=0;j<subMatrixSize;j++) { retval[i][j] =C11[i][j]; retval[i][j+subMatrixSize] =C12[i][j]; retval[i+subMatrixSize][j] =C21[i][j]; retval[i+subMatrixSize][j+subMatrixSize]=C22[i][j]; } } return retval; }
void strassenMultiplication(vector <vector <int> > matrixA, vector <vector <int> > matrixB, int *jumlahKali, int *jumlahTambah, vector <vector <int> > *matrixHasil) /* Menjalankan perkalian matrix strassen dan mengembalikan * Jumlah operasi kali dalam jumlahKali * Jumlah operasi tambah dalam jumlahTambah * serta mengembalikan matrixHasil yang telah dikalikan*/ { int ukuranSubMatrix=matrixA.size()/2; if (matrixA.size()==1) { (*matrixHasil)[0][0]=matrixA[0][0]*matrixB[0][0]; (*jumlahTambah)++; } else if (ukuranSubMatrix==1) { //kalikan rumusnya dengan cara biasa int A11,A21,A12,A22; int B11,B21,B12,B22; int M1,M2,M3,M4,M5,M6,M7; int C11,C21,C12,C22; A11=get11(matrixA); A21=get21(matrixA); A12=get12(matrixA); A22=get22(matrixA); B11=get11(matrixB); B21=get21(matrixB); B12=get12(matrixB); B22=get22(matrixB); M1=(A12-A22)*(B21+B22); M2=(A11+A22)*(B11+B22); M3=(A11-A21)*(B11+B12); M4=(A11+A12)*B22; M5=A11*(B12-B22); M6=A22*(B21-B11); M7=(A21+A22)*B11; (*jumlahKali)+=7; (*jumlahTambah)+=10; C11=M1+M2-M4+M6; C12=M4+M5; C21=M6+M7; C22=M2-M3+M5-M7; (*jumlahTambah)+=8; (*matrixHasil)[0][0]=C11; (*matrixHasil)[0][1]=C12; (*matrixHasil)[1][0]=C21; (*matrixHasil)[1][1]=C22; } else { //rekursif vector< vector<int> > A11(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > A12(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > A21(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > A22(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > B11(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > B12(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > B21(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > B22(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > C11(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > C12(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > C21(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > C22(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M1(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M2(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M3(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M4(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M5(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M6(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); vector< vector<int> > M7(ukuranSubMatrix, vector<int>(ukuranSubMatrix)); A11=getMatrix11(matrixA,ukuranSubMatrix); A21=getMatrix21(matrixA,ukuranSubMatrix); A12=getMatrix12(matrixA,ukuranSubMatrix); A22=getMatrix22(matrixA,ukuranSubMatrix); B11=getMatrix11(matrixB,ukuranSubMatrix); B21=getMatrix21(matrixB,ukuranSubMatrix); B12=getMatrix12(matrixB,ukuranSubMatrix); B22=getMatrix22(matrixB,ukuranSubMatrix); strassenMultiplication(strassenSubtraction(A12,A22,jumlahTambah),strassenAddition(B21,B22,jumlahTambah),jumlahKali,jumlahTambah,&M1); strassenMultiplication(strassenAddition(A11,A22,jumlahTambah),strassenAddition(B11,B22,jumlahTambah),jumlahKali,jumlahTambah,&M2); strassenMultiplication(strassenSubtraction(A11,A21,jumlahTambah),strassenAddition(B11,B12,jumlahTambah),jumlahKali,jumlahTambah,&M3); strassenMultiplication(strassenAddition(A11,A12,jumlahTambah),B22,jumlahKali,jumlahTambah,&M4); strassenMultiplication(A11,strassenSubtraction(B12,B22,jumlahTambah),jumlahKali,jumlahTambah,&M5); strassenMultiplication(A22,strassenSubtraction(B21,B11,jumlahTambah),jumlahKali,jumlahTambah,&M6); strassenMultiplication(strassenAddition(A21,A22,jumlahTambah),B11,jumlahKali,jumlahTambah,&M7); C11=strassenAddition(strassenSubtraction(strassenAddition(M1,M2,jumlahTambah),M4,jumlahTambah),M6,jumlahTambah); C12=strassenAddition(M4,M5,jumlahTambah); C21=strassenAddition(M6,M7,jumlahTambah); C22=strassenSubtraction(strassenAddition(strassenSubtraction(M2,M3,jumlahTambah),M5,jumlahTambah),M7,jumlahTambah); (*matrixHasil)=JoinMatrix(C11,C12,C21,C22); } }
int main(){ //Deklarasi data cout<<"____________________________________________________________"<<endl<<endl; cout<<"----------------------DIVIDE AND CONQUER--------------------"<<endl; cout<<"____________________________________________________________"<<endl<<endl; cout<<"Masukkan ukuran matrix!\n"; cin>>n; cout<<endl<<"Memulai perkalian dengan brute force"<<endl; srand(time(0)); /*Memulai perkalian matrix konvensional*/ int jumlah_perkalian, jumlah_penjumlahan; jumlah_perkalian=0; jumlah_penjumlahan=0; vector< vector<int> > MatrixA(n, vector<int>(n,0)); vector< vector<int> > MatrixB(n, vector<int>(n,0)); vector< vector<int> > hasil(n, vector<int>(n,0)); initMatrix(&MatrixA); initMatrix(&MatrixB); timeval zul; bruteMultiplication(MatrixA,MatrixB,&hasil,&jumlah_perkalian,&jumlah_penjumlahan,&zul); if(n<=25) { cout<<"matrix pertama yang dibangkitkan :"<<endl; printMatrix(MatrixA,n); cout<<endl; cout<<"matrix kedua yang dibangkitkan :"<<endl; printMatrix(MatrixB,n); cout<<endl; cout<<"Matrix hasil perkalian :"<<endl; printMatrix(hasil,n); } cout<<"Jumlah Operasi Penjumlahan="<<jumlah_penjumlahan<<endl; cout<<"Jumlah Operasi Perkalian ="<<jumlah_perkalian<<endl; cout<<"Waktu yang dibutuhkan untuk bruteforce= "<<zul.tv_sec<<" sec "<<endl; cout<<endl<<"Memulai perkalian dengan algoritma strassen"<<endl; /*Memulai perkalian strassen*/ int panjangSeharusnya; int apakahKelipatanDua; apakahKelipatanDua=isKelipatanDua(n,&panjangSeharusnya); vector< vector<int> > MatrixC(panjangSeharusnya, vector<int>(panjangSeharusnya,0)); vector< vector<int> > MatrixD(panjangSeharusnya, vector<int>(panjangSeharusnya,0)); vector< vector<int> > MatrixHasil(panjangSeharusnya,vector<int>(panjangSeharusnya)); jumlah_perkalian=0; jumlah_penjumlahan=0; samakanMatrix(MatrixA,&MatrixC); samakanMatrix(MatrixB,&MatrixD); cout<<endl; timeval mulaiStrassen,selesaiStrassen; gettimeofday(&mulaiStrassen,NULL); strassenMultiplication(MatrixC,MatrixD,&jumlah_perkalian,&jumlah_penjumlahan, &MatrixHasil); gettimeofday(&selesaiStrassen,NULL); timeval elapsedStrassen; elapsedStrassen.tv_sec=selesaiStrassen.tv_sec-mulaiStrassen.tv_sec; elapsedStrassen.tv_usec=selesaiStrassen.tv_usec-mulaiStrassen.tv_usec; cout<<"Hasil perkalian dengan algoritma strassen"<<endl; printMatrix(MatrixHasil,n); cout<<"Jumlah Operasi Penjumlahan="<<jumlah_penjumlahan<<endl; cout<<"Jumlah Operasi Perkalian ="<<jumlah_perkalian<<endl; cout<<"Waktu yang dibutuhkan ="<<elapsedStrassen.tv_sec<<" sec "<<endl; cout<<"SELESAI"<<endl<<endl; return 0; }
Referensi
- Strassen, Volker, Gaussian Elimination is not Optimal, Numer. Math. 13, p. 354-356, 1969
- Thomas H. Cormen, Charles E. Leiserson, Ronald L. Rivest, and Clifford Stein. Introduction to Algorithms, Second Edition. MIT Press and McGraw-Hill, 2001. ISBN 0-262-03293-7. Chapter 28: Section 28.2: Strassen's algorithm for matrix multiplication, pp.735–741.
Pranala luar
- (Inggris) Weisstein, Eric W. "Strassen's Formulas". MathWorld. (also includes formulas for fast matrix inversion)