moebiusband/CG-Bench
1
0
Fork 0
Code Issues Pull Requests Actions Packages Projects Releases Wiki Activity

Finish mm format file reader

Browse Source
This commit is contained in:
Jan Eitzinger 2025-01-03 18:50:57 +01:00
parent d355c4bbdb
commit d267a6d550
Signed by: moebiusband
GPG Key ID: 2574BA29B90D6DD5
7 changed files with 403 additions and 53 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

1
.gitignore vendored
View File

@ -1,2 +1,3 @@
/build /build
.clangd .clangd
exe-*

300
data/matrix_band_klein.mtx Normal file
View File

@ -0,0 +1,300 @@
%%MatrixMarket matrix coordinate real general
100 100 298
1 1 2
2 2 3
3 3 3
4 4 3
5 5 3
6 6 3
7 7 3
8 8 3
9 9 3
10 10 3
11 11 3
12 12 3
13 13 3
14 14 3
15 15 3
16 16 3
17 17 3
18 18 3
19 19 3
20 20 3
21 21 3
22 22 3
23 23 3
24 24 3
25 25 3
26 26 3
27 27 3
28 28 3
29 29 3
30 30 3
31 31 3
32 32 3
33 33 3
34 34 3
35 35 3
36 36 3
37 37 3
38 38 3
39 39 3
40 40 3
41 41 3
42 42 3
43 43 3
44 44 3
45 45 3
46 46 3
47 47 3
48 48 3
49 49 3
50 50 3
51 51 3
52 52 3
53 53 3
54 54 3
55 55 3
56 56 3
57 57 3
58 58 3
59 59 3
60 60 3
61 61 3
62 62 3
63 63 3
64 64 3
65 65 3
66 66 3
67 67 3
68 68 3
69 69 3
70 70 3
71 71 3
72 72 3
73 73 3
74 74 3
75 75 3
76 76 3
77 77 3
78 78 3
79 79 3
80 80 3
81 81 3
82 82 3
83 83 3
84 84 3
85 85 3
86 86 3
87 87 3
88 88 3
89 89 3
90 90 3
91 91 3
92 92 3
93 93 3
94 94 3
95 95 3
96 96 3
97 97 3
98 98 3
99 99 3
100 100 2
1 2 -1
2 3 -1
3 4 -1
4 5 -1
5 6 -1
6 7 -1
7 8 -1
8 9 -1
9 10 -1
10 11 -1
11 12 -1
12 13 -1
13 14 -1
14 15 -1
15 16 -1
16 17 -1
17 18 -1
18 19 -1
19 20 -1
20 21 -1
21 22 -1
22 23 -1
23 24 -1
24 25 -1
25 26 -1
26 27 -1
27 28 -1
28 29 -1
29 30 -1
30 31 -1
31 32 -1
32 33 -1
33 34 -1
34 35 -1
35 36 -1
36 37 -1
37 38 -1
38 39 -1
39 40 -1
40 41 -1
41 42 -1
42 43 -1
43 44 -1
44 45 -1
45 46 -1
46 47 -1
47 48 -1
48 49 -1
49 50 -1
50 51 -1
51 52 -1
52 53 -1
53 54 -1
54 55 -1
55 56 -1
56 57 -1
57 58 -1
58 59 -1
59 60 -1
60 61 -1
61 62 -1
62 63 -1
63 64 -1
64 65 -1
65 66 -1
66 67 -1
67 68 -1
68 69 -1
69 70 -1
70 71 -1
71 72 -1
72 73 -1
73 74 -1
74 75 -1
75 76 -1
76 77 -1
77 78 -1
78 79 -1
79 80 -1
80 81 -1
81 82 -1
82 83 -1
83 84 -1
84 85 -1
85 86 -1
86 87 -1
87 88 -1
88 89 -1
89 90 -1
90 91 -1
91 92 -1
92 93 -1
93 94 -1
94 95 -1
95 96 -1
96 97 -1
97 98 -1
98 99 -1
99 100 -1
2 1 -1
3 2 -1
4 3 -1
5 4 -1
6 5 -1
7 6 -1
8 7 -1
9 8 -1
10 9 -1
11 10 -1
12 11 -1
13 12 -1
14 13 -1
15 14 -1
16 15 -1
17 16 -1
18 17 -1
19 18 -1
20 19 -1
21 20 -1
22 21 -1
23 22 -1
24 23 -1
25 24 -1
26 25 -1
27 26 -1
28 27 -1
29 28 -1
30 29 -1
31 30 -1
32 31 -1
33 32 -1
34 33 -1
35 34 -1
36 35 -1
37 36 -1
38 37 -1
39 38 -1
40 39 -1
41 40 -1
42 41 -1
43 42 -1
44 43 -1
45 44 -1
46 45 -1
47 46 -1
48 47 -1
49 48 -1
50 49 -1
51 50 -1
52 51 -1
53 52 -1
54 53 -1
55 54 -1
56 55 -1
57 56 -1
58 57 -1
59 58 -1
60 59 -1
61 60 -1
62 61 -1
63 62 -1
64 63 -1
65 64 -1
66 65 -1
67 66 -1
68 67 -1
69 68 -1
70 69 -1
71 70 -1
72 71 -1
73 72 -1
74 73 -1
75 74 -1
76 75 -1
77 76 -1
78 77 -1
79 78 -1
80 79 -1
81 80 -1
82 81 -1
83 82 -1
84 83 -1
85 84 -1
86 85 -1
87 86 -1
88 87 -1
89 88 -1
90 89 -1
91 90 -1
92 91 -1
93 92 -1
94 93 -1
95 94 -1
96 95 -1
97 96 -1
98 97 -1
99 98 -1
100 99 -1

52
src/main.c
View File

@ -7,67 +7,42 @@
#include <unistd.h> #include <unistd.h>
#include "comm.h" #include "comm.h"
#include "matrix.h"
#include "progress.h" #include "progress.h"
#include "solver.h" #include "solver.h"
#include "timing.h" #include "timing.h"
static FILE* initResidualWriter()
{
FILE* fp;
fp = fopen("residual.dat", "w");
if (fp == NULL) {
printf("Error!\n");
exit(EXIT_FAILURE);
}
return fp;
}
static void writeResidual(FILE* fp, double ts, double res)
{
fprintf(fp, "%f, %f\n", ts, res);
}
int main(int argc, char** argv) int main(int argc, char** argv)
{ {
int rank; int rank;
double timeStart, timeStop; double timeStart, timeStop;
Parameter p; Matrix m;
Solver s; Solver s;
commInit(s.comm, argc, argv); // commInit(s.comm, argc, argv);
initParameter(&p);
FILE* fp; // FILE* fp;
if (commIsMaster(s.comm)) fp = initResidualWriter(); // if (commIsMaster(s.comm)) fp = initResidualWriter();
if (argc != 2) { if (argc != 2) {
printf("Usage: %s <configFile>\n", argv[0]); printf("Usage: %s <configFile>\n", argv[0]);
exit(EXIT_SUCCESS); exit(EXIT_SUCCESS);
} }
readParameter(&p, argv[1]); matrixRead(&m, argv[1]);
commPartition(s.comm, p.jmax, p.imax); matrixDump(&m);
if (commIsMaster(s.comm)) { exit(EXIT_SUCCESS);
printParameter(&p);
}
initSolver(&s, &p); double te = 1000;
// #ifndef VERBOSE
// initProgress(s.te);
// #endif
double te = p.te;
double t = 0.0; double t = 0.0;
double res = 0.0; double res = 0.0;
timeStart = getTimeStamp(); timeStart = getTimeStamp();
while (t <= te) { while (t <= te) {
if (commIsMaster(s.comm)) writeResidual(fp, t, res); // if (commIsMaster(s.comm)) writeResidual(fp, t, res);
t += p.dt; t += 0.5;
#ifdef VERBOSE #ifdef VERBOSE
if (commIsMaster(s.comm)) { if (commIsMaster(s.comm)) {
@ -78,13 +53,10 @@ int main(int argc, char** argv)
#endif #endif
} }
timeStop = getTimeStamp(); timeStop = getTimeStamp();
#ifndef VERBOSE
stopProgress();
#endif
if (commIsMaster(s.comm)) { if (commIsMaster(s.comm)) {
printf("Solution took %.2fs\n", timeStop - timeStart); printf("Solution took %.2fs\n", timeStop - timeStart);
} }
if (commIsMaster(s.comm)) fclose(fp); // if (commIsMaster(s.comm)) fclose(fp);
commFinalize(s.comm); commFinalize(s.comm);
return EXIT_SUCCESS; return EXIT_SUCCESS;
} }

58
src/matrix.c
View File

@ -31,7 +31,6 @@ static inline int compareRow(const void* a, const void* b)
void matrixRead(Matrix* m, char* filename) void matrixRead(Matrix* m, char* filename)
{ {
int ret_code;
MM_typecode matcode; MM_typecode matcode;
FILE* f; FILE* f;
int M, N, nz; int M, N, nz;
@ -81,6 +80,8 @@ void matrixRead(Matrix* m, char* filename)
exit(EXIT_FAILURE); exit(EXIT_FAILURE);
} }
printf("Read matrix %s with %d non zeroes and %d rows", filename, nz, M);
m->nr = M; m->nr = M;
m->nnz = nz; m->nnz = nz;
Entry* mm; Entry* mm;
@ -121,10 +122,59 @@ void matrixRead(Matrix* m, char* filename)
} }
fclose(f); fclose(f);
size_t mms = cursor;
// sort by column // sort by column
qsort(mm, cursor, sizeof(Entry), compareColumn); qsort(mm, mms, sizeof(Entry), compareColumn);
// sort by row // sort by row
qsort(mm, cursor, sizeof(Entry), compareRow); qsort(mm, mms, sizeof(Entry), compareRow);
m->rowPtr = (CG_UINT*)allocate(64, (m->nr + 1) * sizeof(CG_UINT));
m->colInd = (CG_UINT*)allocate(64, m->nnz * sizeof(CG_UINT));
m->val = (CG_FLOAT*)allocate(64, m->nnz * sizeof(CG_FLOAT));
int* valsPerRow = (int*)allocate(64, m->nr * sizeof(int));
for (int i = 0; i < m->nr; i++) {
valsPerRow[i] = 0;
}
for (int i = 0; i < mms; i++) {
valsPerRow[mm[i].row]++;
}
int* offsets = (int*)allocate(64, (m->nr + 1) * sizeof(int));
m->rowPtr[0] = 0;
// convert to CRS format
for (int rowID = 0; rowID < m->nr; rowID++) {
m->rowPtr[rowID + 1] = m->rowPtr[rowID] + valsPerRow[rowID];
// loop over all elements in Row
for (int id = m->rowPtr[rowID]; id < m->rowPtr[rowID + 1]; id++) {
m->val[id] = (CG_UINT)mm[id].val;
m->colInd[id] = (CG_FLOAT)mm[id].col;
}
}
}
void matrixDump(Matrix* m)
{
CG_UINT numRows = m->nr;
CG_UINT* rowPtr = m->rowPtr;
CG_UINT* colInd = m->colInd;
CG_FLOAT* val = m->val;
printf("Matrix: %lld non zeroes, number of rows %lld\n", numRows, m->nnz);
for (int rowID = 0; rowID < numRows; rowID++) {
printf("Row [%d]: ", rowID);
for (size_t rowEntry = rowPtr[rowID]; rowEntry < rowPtr[rowID + 1]; rowEntry++) {
printf("[%lld]:%.2f ", colInd[rowEntry], val[rowEntry]);
}
printf("\n");
}
} }

1
src/matrix.h
View File

@ -20,5 +20,6 @@ typedef struct {
} Matrix; } Matrix;
extern void matrixRead(Matrix* m, char* filename); extern void matrixRead(Matrix* m, char* filename);
extern void matrixDump(Matrix* m);
#endif // __MATRIX_H_ #endif // __MATRIX_H_

28
src/solver.c Normal file
View File

@ -0,0 +1,28 @@
/* Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
* All rights reserved.
* Use of this source code is governed by a MIT-style
* license that can be found in the LICENSE file.*/
#include <stdlib.h>
#include "matrix.h"
#include "solver.h"
#include "util.h"
void spMVM(Matrix* m, double* x, double* y)
{
CG_UINT numRows = m->nr;
CG_UINT* rowPtr = m->rowPtr;
CG_UINT* colInd = m->colInd;
CG_FLOAT* val = m->val;
for (size_t rowID = 0; rowID < numRows; rowID++) {
CG_FLOAT tmp = y[rowID];
// loop over all elements in row
for (size_t rowEntry = rowPtr[rowID]; rowEntry < rowPtr[rowID + 1]; rowEntry++) {
tmp += val[rowEntry] * x[colInd[rowEntry]];
}
y[rowID] = tmp;
}
}

16
src/solver.h
View File

@ -1,26 +1,24 @@
/* /* Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
* Copyright (C) NHR@FAU, University Erlangen-Nuremberg.
* All rights reserved. This file is part of nusif-solver. * All rights reserved. This file is part of nusif-solver.
* Use of this source code is governed by a MIT style * Use of this source code is governed by a MIT style
* license that can be found in the LICENSE file. * license that can be found in the LICENSE file. */
*/
#ifndef __SOLVER_H_ #ifndef __SOLVER_H_
#define __SOLVER_H_ #define __SOLVER_H_
#include "comm.h" #include "comm.h"
#include "mpi.h" #include "matrix.h"
#include "parameter.h" #include "parameter.h"
#include "util.h"
typedef struct { typedef struct {
/* geometry and grid information */ /* geometry and grid information */
/* parameters */ /* parameters */
double eps, omega; double eps, omega;
int itermax; int itermax;
int levels, presmooth, postsmooth;
double **r, **e;
/* communication */ /* communication */
Comm* comm; Comm* comm;
} Solver; } Solver;
void initSolver(Solver*, Parameter*); void initSolver(Solver*, Parameter*);
double solve(Solver*, double*, double*); void spMVM(Matrix* m, CG_FLOAT* x, CG_FLOAT* y);
#endif #endif // __SOLVER_H_