-
Notifications
You must be signed in to change notification settings - Fork 0
/
Matrix_multi.c
132 lines (126 loc) · 2.94 KB
/
Matrix_multi.c
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
#include <stdlib.h>
#include <stdio.h>
#include <math.h>
#include <omp.h>
#define NUM_THREADS 3
double r8_uniform_01(int *seed)
{
int k;
double r;
k = *seed / 127773;
*seed = 16807 * (*seed - k * 127773) - k * 2836;
if (*seed < 0)
{
*seed = *seed + 2147483647;
}
r = (double)(*seed) * 4.656612875E-10;
return r;
}
void r8_mxm(int l, int m, int n)
{
double *a;
double *b;
double *c;
int i;
int j;
int k;
int ops;
double rate;
int seed;
double time_begin;
double funcTime;
double time_elapsed;
double time_stop;
/*
Allocate the matrices.
*/
a = (double *)malloc(l * n * sizeof(double));
b = (double *)malloc(l * m * sizeof(double));
c = (double *)malloc(m * n * sizeof(double));
/*
Assign values to the B and C matrices.
*/
seed = 123456789;
for (k = 0; k < l * m; k++)
{
b[k] = r8_uniform_01(&seed);
}
for (k = 0; k < m * n; k++)
{
c[k] = r8_uniform_01(&seed);
}
/*
Compute A = B * C.
*/
funcTime = omp_get_wtime();
time_begin = omp_get_wtime();
omp_set_num_threads(NUM_THREADS);
int nthreads, id;
#pragma omp parallel shared(a, b, c, l, m, n) private(i, j, k)
{
#pragma omp for
for (j = 0; j < n; j++)
{
for (i = 0; i < l; i++)
{
a[i + j * l] = 0.0;
for (k = 0; k < m; k++)
{
a[i + j * l] = a[i + j * l] + b[i + k * l] * c[k + j * m];
}
}
}
/*id = omp_get_thread_num();
if (id == 0)
nthreads = omp_get_num_threads();*/
}
time_stop = omp_get_wtime();
/*printf(" Working number of threads = %d\n", nthreads);*/
/*
Report.
*/
ops = l * n * (2 * m);
time_elapsed = time_stop - time_begin -(time_begin - funcTime);
rate = (double)(ops) / time_elapsed / 1000000.0;
printf("\n");
printf("R8_MXM matrix multiplication timing.\n");
printf(" A(LxN) = B(LxM) * C(MxN).\n");
printf(" L = %d\n", l);
printf(" M = %d\n", m);
printf(" N = %d\n", n);
printf(" Floating point OPS roughly %d\n", ops);
printf(" Elapsed time dT = %f\n", time_elapsed);
printf(" Rate = MegaOPS/dT = %f\n", rate);
free(a);
free(b);
free(c);
return;
}
int main(int argc, char *argv[])
{
int id;
int l;
int m;
int n;
printf("\n");
printf("MXM\n");
printf(" C/OpenMP version.\n");
printf("\n");
printf(" Matrix multiplication tests.\n");
printf("\n");
printf(" Number of processors available = %d\n",
omp_get_num_procs());
printf(" Max number of threads = %d\n",
omp_get_max_threads());
l = 200;
m = 200;
n = 200;
r8_mxm(l, m, n);
/*
Terminate.
*/
printf("\n");
printf("MXM:\n");
printf(" Normal end of execution.\n");
return 0;
}