Evo mog koda napisanog u C-u, na Linuxu a imam ako treba i za Windowse.
Code:
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <sys/time.h>
#define DIM 512 //dimenzija matrice
#define NISTA 2
#define RANDOM 1
#define NULAMA 0
#define DODAJ 1
#define NOVA 0
void napravi_matricu (int ***M, int dim_matrice, int nacin);
void ispisi_matricu (int **M);
void dealociraj_matricu (int **M, int dim_matrice);
void suma_matrica (int **A, int **B, int **C, int dim_matrice, int nacin); //C = A + B
void suma1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin);
void razlika_matrica (int **A, int **B, int **C, int dim_matrice, int nacin); //C = A - B
void razlika1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin);
void strassen_matrice (int **A, int **B, int **C, int dim_matrice);
int main (void)
{
int **A, **B, **C;
double start, end;
struct timeval vrijeme;
srand( 0 );
napravi_matricu( &A, DIM, RANDOM );
napravi_matricu( &B, DIM, RANDOM );
napravi_matricu( &C, DIM, NULAMA );
gettimeofday( &vrijeme, NULL );
start = vrijeme.tv_sec + (vrijeme.tv_usec/1000000.0);
strassen_matrice( A, B, C, DIM );
gettimeofday( &vrijeme, NULL );
end = vrijeme.tv_sec + (vrijeme.tv_usec/1000000.0);
/* printf( "--- A ---\n" ); ispisi_matricu( A );
printf( "--- B ---\n" ); ispisi_matricu( B );
printf( "--- C ---\n" ); ispisi_matricu( C );
*/ dealociraj_matricu( A, DIM );
dealociraj_matricu( B, DIM );
dealociraj_matricu( C, DIM );
printf( "%d - Vrijeme: %fs\n", DIM, end-start );
return 0;
}
void napravi_matricu (int ***M, int dim_matrice, int nacin)
{
int i, j;
(*M) = (int **)malloc( dim_matrice * sizeof(int *) );
if ((*M) == NULL) {
printf( "Greska: Alokacija memorije nije uspjela!\n" );
exit( -1 );
}
for (i = 0; i < dim_matrice; i++) {
(*M)[i] = (int *)malloc( dim_matrice * sizeof(int) );
if ((*M)[i] == NULL) {
printf( "Greska 2: Alokacija memorije nije uspjela!\n" );
exit( -1 );
}
}
if (nacin == RANDOM)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
(*M)[i][j] = rand() % 101;
else if (nacin == NULAMA)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
(*M)[i][j] = 0;
return;
}
void ispisi_matricu (int **M)
{
int i, j;
for (i = 0; i < DIM; i++) {
for (j = 0; j < DIM; j++)
printf( "%7d", M[i][j] );
printf( "\n" );
}
return;
}
void dealociraj_matricu (int **M, int dim_matrice)
{
int i;
for (i = 0; i < dim_matrice; i++)
free( M[i] );
free( M );
return;
}
void suma_matrica (int **A, int **B, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i][j] + B[i][j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i][j] + B[i][j];
return;
}
void razlika_matrica (int **A, int **B, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i][j] - B[i][j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i][j] - B[i][j];
return;
}
void suma1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i1+i][j1+j] + B[i2+i][j2+j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i1+i][j1+j] + B[i2+i][j2+j];
return;
}
void razlika1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i1+i][j1+j] - B[i2+i][j2+j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i1+i][j1+j] - B[i2+i][j2+j];
return;
}
void strassen_matrice (int **A, int **B, int **C, int dim_matrice)
{
int i, j, k,
dim_submatrice = dim_matrice/2,
**Temp1, **Temp2, **Temp3,
**M1, **M2, **M3, **M4, **M5, **M6, **M7;
if (dim_matrice <= 1) {
C[0][0] = A[0][0] * B[0][0];
return;
}
napravi_matricu( &Temp1, dim_submatrice, NISTA ); napravi_matricu( &Temp2, dim_submatrice, NISTA );
napravi_matricu( &Temp3, dim_submatrice, NISTA );
napravi_matricu( &M1, dim_submatrice, NISTA ); napravi_matricu( &M2, dim_submatrice, NISTA );
napravi_matricu( &M3, dim_submatrice, NISTA ); napravi_matricu( &M4, dim_submatrice, NISTA );
napravi_matricu( &M5, dim_submatrice, NISTA ); napravi_matricu( &M6, dim_submatrice, NISTA );
napravi_matricu( &M7, dim_submatrice, NISTA );
suma1_matrica( A, 0, 0, A, dim_submatrice, dim_submatrice, Temp1, dim_submatrice, NOVA );
suma1_matrica( B, 0, 0, B, dim_submatrice, dim_submatrice, Temp2, dim_submatrice, NOVA );
strassen_matrice( Temp1, Temp2, M1, dim_submatrice );
suma1_matrica( A, dim_submatrice, 0, A, dim_submatrice, dim_submatrice, Temp1, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++)
Temp3[i][j] = B[i][j];
strassen_matrice( Temp1, Temp3, M2, dim_submatrice );
razlika1_matrica( B, 0, dim_submatrice, B, dim_submatrice, dim_submatrice, Temp2, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++)
Temp3[i][j] = A[i][j];
strassen_matrice( Temp3, Temp2, M3, dim_submatrice );
razlika1_matrica( B, dim_submatrice, 0, B, 0, 0, Temp2, dim_submatrice, NOVA );
for (i = dim_submatrice; i < dim_matrice; i++)
for (j = dim_submatrice; j < dim_matrice; j++)
Temp3[(i-dim_submatrice)][(j-dim_submatrice)] = A[i][j];
strassen_matrice( Temp3, Temp2, M4, dim_submatrice );
suma1_matrica( A, 0, 0, A, 0, dim_submatrice, Temp1, dim_submatrice, NOVA );
for (i = dim_submatrice; i < dim_matrice; i++)
for (j = dim_submatrice; j < dim_matrice; j++)
Temp3[(i-dim_submatrice)][(j-dim_submatrice)] = B[i][j];
strassen_matrice( Temp1, Temp3, M5, dim_submatrice );
razlika1_matrica( A, dim_submatrice, 0, A, 0, 0, Temp1, dim_submatrice, NOVA );
suma1_matrica( B, 0, 0, B, 0, dim_submatrice, Temp2, dim_submatrice, NOVA );
strassen_matrice( Temp1, Temp2, M6, dim_submatrice );
razlika1_matrica( A, 0, dim_submatrice, A, dim_submatrice, dim_submatrice, Temp1, dim_submatrice, NOVA );
suma1_matrica( B, dim_submatrice, 0, B, dim_submatrice, dim_submatrice, Temp2, dim_submatrice, NOVA );
strassen_matrice( Temp1, Temp2, M7, dim_submatrice );
razlika_matrica( M4, M5, Temp1, dim_submatrice, NOVA );
suma_matrica( M1, M7, Temp1, dim_submatrice, DODAJ );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[i][j] = Temp1[i][j];
}
suma_matrica( M3, M5, Temp1, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[i][(j+dim_submatrice)] = Temp1[i][j];
}
suma_matrica( M2, M4, Temp1, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[(i+dim_submatrice)][j] = Temp1[i][j];
}
razlika_matrica( M1, M2, Temp1, dim_submatrice, NOVA );
suma_matrica( M3, M6, Temp1, dim_submatrice, DODAJ );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[(i+dim_submatrice)][(j+dim_submatrice)] = Temp1[i][j];
}
dealociraj_matricu( Temp1, dim_submatrice ); dealociraj_matricu( Temp2, dim_submatrice );
dealociraj_matricu( Temp3, dim_submatrice );
dealociraj_matricu( M1, dim_submatrice ); dealociraj_matricu( M2, dim_submatrice );
dealociraj_matricu( M3, dim_submatrice ); dealociraj_matricu( M4, dim_submatrice );
dealociraj_matricu( M5, dim_submatrice ); dealociraj_matricu( M6, dim_submatrice );
dealociraj_matricu( M7, dim_submatrice );
return;
}
#include <stdio.h>
#include <stdlib.h>
#include <time.h>
#include <sys/time.h>
#define DIM 512 //dimenzija matrice
#define NISTA 2
#define RANDOM 1
#define NULAMA 0
#define DODAJ 1
#define NOVA 0
void napravi_matricu (int ***M, int dim_matrice, int nacin);
void ispisi_matricu (int **M);
void dealociraj_matricu (int **M, int dim_matrice);
void suma_matrica (int **A, int **B, int **C, int dim_matrice, int nacin); //C = A + B
void suma1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin);
void razlika_matrica (int **A, int **B, int **C, int dim_matrice, int nacin); //C = A - B
void razlika1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin);
void strassen_matrice (int **A, int **B, int **C, int dim_matrice);
int main (void)
{
int **A, **B, **C;
double start, end;
struct timeval vrijeme;
srand( 0 );
napravi_matricu( &A, DIM, RANDOM );
napravi_matricu( &B, DIM, RANDOM );
napravi_matricu( &C, DIM, NULAMA );
gettimeofday( &vrijeme, NULL );
start = vrijeme.tv_sec + (vrijeme.tv_usec/1000000.0);
strassen_matrice( A, B, C, DIM );
gettimeofday( &vrijeme, NULL );
end = vrijeme.tv_sec + (vrijeme.tv_usec/1000000.0);
/* printf( "--- A ---\n" ); ispisi_matricu( A );
printf( "--- B ---\n" ); ispisi_matricu( B );
printf( "--- C ---\n" ); ispisi_matricu( C );
*/ dealociraj_matricu( A, DIM );
dealociraj_matricu( B, DIM );
dealociraj_matricu( C, DIM );
printf( "%d - Vrijeme: %fs\n", DIM, end-start );
return 0;
}
void napravi_matricu (int ***M, int dim_matrice, int nacin)
{
int i, j;
(*M) = (int **)malloc( dim_matrice * sizeof(int *) );
if ((*M) == NULL) {
printf( "Greska: Alokacija memorije nije uspjela!\n" );
exit( -1 );
}
for (i = 0; i < dim_matrice; i++) {
(*M)[i] = (int *)malloc( dim_matrice * sizeof(int) );
if ((*M)[i] == NULL) {
printf( "Greska 2: Alokacija memorije nije uspjela!\n" );
exit( -1 );
}
}
if (nacin == RANDOM)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
(*M)[i][j] = rand() % 101;
else if (nacin == NULAMA)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
(*M)[i][j] = 0;
return;
}
void ispisi_matricu (int **M)
{
int i, j;
for (i = 0; i < DIM; i++) {
for (j = 0; j < DIM; j++)
printf( "%7d", M[i][j] );
printf( "\n" );
}
return;
}
void dealociraj_matricu (int **M, int dim_matrice)
{
int i;
for (i = 0; i < dim_matrice; i++)
free( M[i] );
free( M );
return;
}
void suma_matrica (int **A, int **B, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i][j] + B[i][j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i][j] + B[i][j];
return;
}
void razlika_matrica (int **A, int **B, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i][j] - B[i][j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i][j] - B[i][j];
return;
}
void suma1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i1+i][j1+j] + B[i2+i][j2+j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i1+i][j1+j] + B[i2+i][j2+j];
return;
}
void razlika1_matrica (int **A, int i1, int j1, int **B, int i2, int j2, int **C, int dim_matrice, int nacin)
{
int i, j;
if (nacin == DODAJ)
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i1+i][j1+j] - B[i2+i][j2+j];
else
for (i = 0; i < dim_matrice; i++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = A[i1+i][j1+j] - B[i2+i][j2+j];
return;
}
void strassen_matrice (int **A, int **B, int **C, int dim_matrice)
{
int i, j, k,
dim_submatrice = dim_matrice/2,
**Temp1, **Temp2, **Temp3,
**M1, **M2, **M3, **M4, **M5, **M6, **M7;
if (dim_matrice <= 1) {
C[0][0] = A[0][0] * B[0][0];
return;
}
napravi_matricu( &Temp1, dim_submatrice, NISTA ); napravi_matricu( &Temp2, dim_submatrice, NISTA );
napravi_matricu( &Temp3, dim_submatrice, NISTA );
napravi_matricu( &M1, dim_submatrice, NISTA ); napravi_matricu( &M2, dim_submatrice, NISTA );
napravi_matricu( &M3, dim_submatrice, NISTA ); napravi_matricu( &M4, dim_submatrice, NISTA );
napravi_matricu( &M5, dim_submatrice, NISTA ); napravi_matricu( &M6, dim_submatrice, NISTA );
napravi_matricu( &M7, dim_submatrice, NISTA );
suma1_matrica( A, 0, 0, A, dim_submatrice, dim_submatrice, Temp1, dim_submatrice, NOVA );
suma1_matrica( B, 0, 0, B, dim_submatrice, dim_submatrice, Temp2, dim_submatrice, NOVA );
strassen_matrice( Temp1, Temp2, M1, dim_submatrice );
suma1_matrica( A, dim_submatrice, 0, A, dim_submatrice, dim_submatrice, Temp1, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++)
Temp3[i][j] = B[i][j];
strassen_matrice( Temp1, Temp3, M2, dim_submatrice );
razlika1_matrica( B, 0, dim_submatrice, B, dim_submatrice, dim_submatrice, Temp2, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++)
Temp3[i][j] = A[i][j];
strassen_matrice( Temp3, Temp2, M3, dim_submatrice );
razlika1_matrica( B, dim_submatrice, 0, B, 0, 0, Temp2, dim_submatrice, NOVA );
for (i = dim_submatrice; i < dim_matrice; i++)
for (j = dim_submatrice; j < dim_matrice; j++)
Temp3[(i-dim_submatrice)][(j-dim_submatrice)] = A[i][j];
strassen_matrice( Temp3, Temp2, M4, dim_submatrice );
suma1_matrica( A, 0, 0, A, 0, dim_submatrice, Temp1, dim_submatrice, NOVA );
for (i = dim_submatrice; i < dim_matrice; i++)
for (j = dim_submatrice; j < dim_matrice; j++)
Temp3[(i-dim_submatrice)][(j-dim_submatrice)] = B[i][j];
strassen_matrice( Temp1, Temp3, M5, dim_submatrice );
razlika1_matrica( A, dim_submatrice, 0, A, 0, 0, Temp1, dim_submatrice, NOVA );
suma1_matrica( B, 0, 0, B, 0, dim_submatrice, Temp2, dim_submatrice, NOVA );
strassen_matrice( Temp1, Temp2, M6, dim_submatrice );
razlika1_matrica( A, 0, dim_submatrice, A, dim_submatrice, dim_submatrice, Temp1, dim_submatrice, NOVA );
suma1_matrica( B, dim_submatrice, 0, B, dim_submatrice, dim_submatrice, Temp2, dim_submatrice, NOVA );
strassen_matrice( Temp1, Temp2, M7, dim_submatrice );
razlika_matrica( M4, M5, Temp1, dim_submatrice, NOVA );
suma_matrica( M1, M7, Temp1, dim_submatrice, DODAJ );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[i][j] = Temp1[i][j];
}
suma_matrica( M3, M5, Temp1, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[i][(j+dim_submatrice)] = Temp1[i][j];
}
suma_matrica( M2, M4, Temp1, dim_submatrice, NOVA );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[(i+dim_submatrice)][j] = Temp1[i][j];
}
razlika_matrica( M1, M2, Temp1, dim_submatrice, NOVA );
suma_matrica( M3, M6, Temp1, dim_submatrice, DODAJ );
for (i = 0; i < dim_submatrice; i++)
for (j = 0; j < dim_submatrice; j++){
C[(i+dim_submatrice)][(j+dim_submatrice)] = Temp1[i][j];
}
dealociraj_matricu( Temp1, dim_submatrice ); dealociraj_matricu( Temp2, dim_submatrice );
dealociraj_matricu( Temp3, dim_submatrice );
dealociraj_matricu( M1, dim_submatrice ); dealociraj_matricu( M2, dim_submatrice );
dealociraj_matricu( M3, dim_submatrice ); dealociraj_matricu( M4, dim_submatrice );
dealociraj_matricu( M5, dim_submatrice ); dealociraj_matricu( M6, dim_submatrice );
dealociraj_matricu( M7, dim_submatrice );
return;
}
E sad on je puno sporiji nego kad napišem običan algoritam za množenje, a trebao bi biti brži.
Znam da bi trebao biti brži za velike matrice, ali čak i kad stavim da mi je red matrice 2048 on je još uvijek sporiji.
Radio sam neka poboljšanja u programu da bude malo brži od ovog al je i taj puno sporiji od običnog. (nisam ga stavio na forum jer je teško čitljiv)
Probao sam maknut pozivanja funkcije razlika1_matrica i suma1_matrica u funkciji strassen_matrice, i stavit taj kod direkno tamo gdje se one pozivaju, al poboljšanje je gotovo nikakvo.
Još sam probao srezat doljnje rekurzije tako da sam ovaj dio koda:
Code:
if (dim_matrice <= 1) {
C[0][0] = A[0][0] * B[0][0];
return;
}
if (dim_matrice <= 1) {
C[0][0] = A[0][0] * B[0][0];
return;
}
zamjenio sa:
Code:
if (dim_matrice <= 64) {
for (i = 0; i < dim_matrice; i++)
for (k = 0; k < dim_matrice; k++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = 0;
for (i = 0; i < dim_matrice; i++)
for (k = 0; k < dim_matrice; k++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i][k] * B[k][j];
return;
}
if (dim_matrice <= 64) {
for (i = 0; i < dim_matrice; i++)
for (k = 0; k < dim_matrice; k++)
for (j = 0; j < dim_matrice; j++)
C[i][j] = 0;
for (i = 0; i < dim_matrice; i++)
for (k = 0; k < dim_matrice; k++)
for (j = 0; j < dim_matrice; j++)
C[i][j] += A[i][k] * B[k][j];
return;
}
i sada je strassen brži (zapravo nije ni pravi strassen sada), al tek kada je red matrice >= 2048. (brži je za 15 sec).
Bio bih zahvalan kad bi mi netko pogledao program i rekao kako da ga ubrzam, ili dao neki svoj program za Strassena napisanog u C/C++, može i nekom drugom sličnom jeziku.