#include <stdio.h>
#include <cstdlib>
#include "mmio.h"
#include <cstring>
#include <time.h>
//-------------------------------------------------------------------------------------------
#ifndef MM_IO_H
#define MM_IO_H
#define MM_MAX_LINE_LENGTH 1025
#define MatrixMarketBanner "%%MatrixMarket"
#define MM_MAX_TOKEN_LENGTH 64
typedef char MM_typecode[4];
char *mm_typecode_to_str(MM_typecode matcode);
int mm_read_banner(FILE *f, MM_typecode *matcode);
int mm_read_mtx_crd_size(FILE *f, int *M, int *N, int *nz);
int mm_read_mtx_array_size(FILE *f, int *M, int *N);
int mm_write_banner(FILE *f, MM_typecode matcode);
int mm_write_mtx_crd_size(FILE *f, int M, int N, int nz);
int mm_write_mtx_array_size(FILE *f, int M, int N);
/********************* MM_typecode query fucntions ***************************/
#define mm_is_matrix(typecode) ((typecode)[0]=='M')
#define mm_is_sparse(typecode) ((typecode)[1]=='C')
#define mm_is_coordinate(typecode)((typecode)[1]=='C')
#define mm_is_dense(typecode) ((typecode)[1]=='A')
#define mm_is_array(typecode) ((typecode)[1]=='A')
#define mm_is_complex(typecode) ((typecode)[2]=='C')
#define mm_is_real(typecode) ((typecode)[2]=='R')
#define mm_is_pattern(typecode) ((typecode)[2]=='P')
#define mm_is_integer(typecode) ((typecode)[2]=='I')
#define mm_is_symmetric(typecode)((typecode)[3]=='S')
#define mm_is_general(typecode) ((typecode)[3]=='G')
#define mm_is_skew(typecode) ((typecode)[3]=='K')
#define mm_is_hermitian(typecode)((typecode)[3]=='H')
int mm_is_valid(MM_typecode matcode); /* too complex for a macro */
/********************* MM_typecode modify fucntions ***************************/
#define mm_set_matrix(typecode) ((*typecode)[0]='M')
#define mm_set_coordinate(typecode) ((*typecode)[1]='C')
#define mm_set_array(typecode) ((*typecode)[1]='A')
#define mm_set_dense(typecode) mm_set_array(typecode)
#define mm_set_sparse(typecode) mm_set_coordinate(typecode)
#define mm_set_complex(typecode)((*typecode)[2]='C')
#define mm_set_real(typecode) ((*typecode)[2]='R')
#define mm_set_pattern(typecode)((*typecode)[2]='P')
#define mm_set_integer(typecode)((*typecode)[2]='I')
#define mm_set_symmetric(typecode)((*typecode)[3]='S')
#define mm_set_general(typecode)((*typecode)[3]='G')
#define mm_set_skew(typecode) ((*typecode)[3]='K')
#define mm_set_hermitian(typecode)((*typecode)[3]='H')
#define mm_clear_typecode(typecode) ((*typecode)[0]=(*typecode)[1]= \
(*typecode)[2]=' ',(*typecode)[3]='G')
#define mm_initialize_typecode(typecode) mm_clear_typecode(typecode)
/********************* Matrix Market error codes ***************************/
#define MM_COULD_NOT_READ_FILE 11
#define MM_PREMATURE_EOF 12
#define MM_NOT_MTX 13
#define MM_NO_HEADER 14
#define MM_UNSUPPORTED_TYPE 15
#define MM_LINE_TOO_LONG 16
#define MM_COULD_NOT_WRITE_FILE 17
/******************** Matrix Market internal definitions ********************
MM_matrix_typecode: 4-character sequence
ojbect sparse/ data storage
dense type scheme
string position: [0] [1] [2] [3]
Matrix typecode: M(atrix) C(oord) R(eal) G(eneral)
A(array) C(omplex) H(ermitian)
P(attern) S(ymmetric)
I(nteger) K(kew)
***********************************************************************/
#define MM_MTX_STR "matrix"
#define MM_ARRAY_STR "array"
#define MM_DENSE_STR "array"
#define MM_COORDINATE_STR "coordinate"
#define MM_SPARSE_STR "coordinate"
#define MM_COMPLEX_STR "complex"
#define MM_REAL_STR "real"
#define MM_INT_STR "integer"
#define MM_GENERAL_STR "general"
#define MM_SYMM_STR "symmetric"
#define MM_HERM_STR "hermitian"
#define MM_SKEW_STR "skew-symmetric"
#define MM_PATTERN_STR "pattern"
/* high level routines */
int mm_write_mtx_crd(char fname[], int M, int N, int nz, int I[], int J[],
double val[], MM_typecode matcode);
int mm_read_mtx_crd_data(FILE *f, int M, int N, int nz, int I[], int J[],
double val[], MM_typecode matcode);
int mm_read_mtx_crd_entry(FILE *f, int *I, int *J, double *real, double *img,
MM_typecode matcode);
int mm_read_unsymmetric_sparse(const char *fname, int *M_, int *N_, int *nz_,
double **val_, int **I_, int **J_);
#endif
//-------------------------------------------------------------------------------------------
//#define RAND_MAX = 0x7FFF;
#define maksymalny_random 100
#define minimalny_random 0
void quicksort(int * Wiersze, int * Kolumny, double * Wartosci, int Lo, int Hi)
{
int i,j;
int x;
x = Wiersze[(Lo+Hi)/2];
i = Lo;
j = Hi;
do
{
while (Wiersze[i] < x) ++i;
while (Wiersze[j] > x) --j;
if (i<=j)
{
int tmp = Wiersze[i]; //zamiana wierszy, kolumn i wartosci
Wiersze[i] = Wiersze[j];
Wiersze[j] = tmp;
tmp = Kolumny[i];
Kolumny[i] = Kolumny[j];
Kolumny[j] = tmp;
double tmpd = Wartosci[i];
Wartosci[i] = Wartosci[j];
Wartosci[j] = tmpd;
++i;
--j;
}
} while(i < j);
if (Lo < j) quicksort(Wiersze, Kolumny, Wartosci, Lo, j);
if (Hi > i) quicksort(Wiersze, Kolumny, Wartosci, i, Hi);
}
int main(int argc, char *argv[])
{
int ret_code;
MM_typecode matcode;
FILE *f; //plik wejsciowy
int M, N, nz; //wiersze / kolumny / ilosc elementow na liscie
int i, *I, *J; //iterator / wiersze / kolumny
double *val; //wartosci dla kolejnych komorek macierzy
double *wektor; //wektor liczb dla (wiersz*kolumna)
double *wynik; //wektor wynikowy
//dane do losowania wektora wejsciowego:
// const double minimalny_random=-20; //zakres dolny randoma
// const double maksymalny_random=20; //zakres gorny
// const int RAND_MAX = 0x7FFF; // z <stdlib.h>
//CSR - zajmiemy dużo miejsca na wektor kolumnowy. Póki co nie będę pisał samorozszerzajšcej się tablicy
int * indeksy_kolumn;
int * indeksy_wierszy;
double * wartosci;
//SEKCJA MAGICZNEGO WCZYTYWANIA
if (argc < 2)
{
fprintf(stderr, "Usage: %s [martix-market-filename]\n", argv[0]);
exit(1);
}
else
{
if ((f = fopen(argv[1], "r")) == NULL)
exit(1);
}
printf("Wczytywanie macierzy: %s...",argv[1]);
if (mm_read_banner(f, &matcode) != 0)
{
printf("Could not process Matrix Market banner.\n");
exit(1);
}
if (mm_is_complex(matcode) && mm_is_matrix(matcode) && mm_is_sparse(matcode) )
{
printf("Sorry, this application does not support ");
printf("Market Market type: [%s]\n", mm_typecode_to_str(matcode));
exit(1);
}
if ((ret_code = mm_read_mtx_crd_size(f, &M, &N, &nz)) !=0)
exit(1);
I = (int *) malloc(nz * sizeof(int));
J = (int *) malloc(nz * sizeof(int));
val = (double *) malloc(nz * sizeof(double));
for (i=0; i<nz; i++)
{
fscanf(f, "%d %d %lg\n", &I[i], &J[i], &val[i]);
I[i]--; /* adjust from 1-based to 0-based */
J[i]--;
}
if (f !=stdin) fclose(f);
//KONIEC MAGICZNEGO WCZYTWANIA
printf(" OK!\n\nPodsumowanie: \n ");
mm_write_banner(stdout, matcode);
printf(" wierszy | kolumn | wpisow:\n ");
mm_write_mtx_crd_size(stdout, M, N, nz);
printf("\nQuicksort wzgledem wierszy...");
quicksort(I,J,val,0,nz-1);
printf(" OK!\n");
//niepotrzebne wyswietlanie wszystkich wczytanych elementow macierzy
//!
/*printf("\nwiersz | kolumna | wartosc\n");
for (i=0; i<nz; i++)
fprintf(stdout, "%d %d %20.19g\n", I[i]+1, J[i]+1, val[i]); //*/
printf("Konwersja do CSR...");
indeksy_wierszy = (int*)malloc((M+1)*sizeof(int)); //maksymalny (potencjalnie) potrzebny rozmiar
indeksy_kolumn=J;
int indeks=-1;
for (i=0; i<nz; i++)
if (indeks!=I[i])
{
indeks=I[i];
indeksy_wierszy[indeks]=i;
}
indeksy_wierszy[indeks+1]=nz; //na ostatniej pozycji dodajemy numer+1 ostatniego elementu z tablicy
free(I);
printf(" OK!\n");
//!
/*for (i=0; i<=M; i++)
printf("%d/",indeksy_wierszy[i]); //*/
wektor = (double*)malloc(N*sizeof(double)); //ilosc elementow w wektorze wejsciowym
//przygotowanie nazwy wektora do wczytania
int dl=strlen(argv[1]);
char * wektor_name = (char*)malloc((dl+8)*sizeof(char));
memcpy(wektor_name,argv[1],dl-4);
memmove(wektor_name+dl-4,"_wektor.txt",12);
if ((f = fopen(wektor_name, "r")) == NULL) //probujemy wczytac wektor z pliku. Jesli nie istnieje - tworzymy go
{
printf("Brak pliku wektora! -> Generuje wektor: %s...",wektor_name);
srand(time(NULL));
for (i=0; i<N; i++)
wektor[i] = (double)rand() / (((double)RAND_MAX + 1) / (maksymalny_random-minimalny_random)) + minimalny_random;
f = fopen(wektor_name, "w");
for (i=0; i<N; i++)
fprintf(f, "%lf ", wektor[i]);
} else {
printf("Wczytywanie wektora: %s...",wektor_name);
for (i=0; i<N; i++)
fscanf(f, "%lf", &wektor[i]);
}
fclose(f);
printf(" OK!\nMnozenie...");
wynik = (double*)malloc(M*sizeof(double)); /*zapiszemy caly wektor wynikowy (mało optymalne, ale na teraz musi starczyc)
M, poniewaz: M*N x N*1 = M*1 */
for (i=0; i<M; i++) //wyczyscmy pamiec..
wynik[i]=0;
/*for (i=0; i<nz; i++) //teraz obsluga indeksow I=wiersze, J=kolumny | M=wiesze, N=kolumny
wynik[I[i]]+=val[i]*wektor[J[i]]; //wzorek na */
int wiersz=-1; //mnożenie z wykorzystaniem konstrukcji CSR
while (indeksy_wierszy[++wiersz]!=nz)
for (i=indeksy_wierszy[wiersz]; i<indeksy_wierszy[wiersz+1]; i++)
wynik[wiersz]+=val[i]*wektor[indeksy_kolumn[i]];
printf(" OK!\n");
//przygotowanie nazwy wyjsciowej
dl=strlen(argv[1]);
char * wynik_name = (char*)malloc((dl+7)*sizeof(char));
memcpy(wynik_name,argv[1],dl-4);
memmove(wynik_name+dl-4,"_wynik.txt",11);
printf("Zapisywanie: %s...",wynik_name);
f = fopen(wynik_name, "w"); //wpisujemy do pliku wyjsciowego dane
for (i=0; i<M; i++) //standardowe zapisanie cišgu elementow wektora
fprintf(f, "%lf ", wynik[i]);
fclose(f);
printf(" OK!\n");
//czyszczenie pamieci:
free(wynik_name);
free(wektor_name);
free(wektor);
free(indeksy_wierszy);
free(J);
free(val);
printf("\nKoniec, wcisnij [enter] by wyjsc.");
getchar();
return 0;
}