kmazouzi / FormationMeso

Compare View

switch
from
...
to
 
Commits (4)

Diff

Showing 4 changed files Inline Diff

lab3/Makefile View file @ 2b4ce20
File was created 1 GCC = gcc
2 CFLAGS = -O3 -fopenmp
3 OMP_FLAG = -fopenmp
4 RM = rm -rf
5
6
7 EXE = omp_heat2D ser_heat2D
8
9 all : $(EXE)
10
11 #.PHONY: all clean purge
12
13
14 pi_ser: ser_heat2D.o
15 $(GCC) $(CFLAGS) -o $@ $^
16
17 pi_task: omp_heat2D.o
lab3/omp_heat2D.c View file @ 2b4ce20
File was created 1 /****************************************************************************
2 * DESCRIPTION:
3 * Serial HEAT2D Example - C Version
4 * This example is based on a simplified
5 * two-dimensional heat equation domain decomposition. The initial
6 * temperature is computed to be high in the middle of the domain and
7 * zero at the boundaries. The boundaries are held at zero throughout
8 * the simulation. During the time-stepping, an array containing two
9 * domains is used; these domains alternate between old data and new data.
10 *
11 * The physical region, and the boundary conditions, are suggested
12 by this diagram;
13
14 u = 0
15 +------------------+
16 | |
17 u = 100 | | u = 100
18 | |
19 +------------------+
20 u = 100
21
22 Interrior point :
23 u[Central] = (1/4) * ( u[North] + u[South] + u[East] + u[West] )
24
25 ****************************************************************************/
26 #include <stdio.h>
27 #include <stdlib.h>
28 #include <math.h>
29 #include <omp.h>
30
31 #define NN 50
32 #define MM 50
33
34 #define ITER_PRINT 100
35 #define PRINT_DATA 1
36
37 #define _EPSILON 0.001
38
39
40 void update(int nx,int ny, float *u, float *unew, float * diff);
41 void inidat(int nx, int ny, float *u, float *unew);
42 void prtdat(int nx, int ny, float *u,const char *fnam);
43
44
45
46
47 int main(int argc, char *argv[])
48 {
49
50 float diff=1.0;
51 float EPSILON=_EPSILON;
52 int N=NN,M=MM;
53
54 if(argc !=3)
55 {
56 fprintf(stderr,"usage %s N EPSILON\n ", argv[0]);
57 fprintf(stderr,"\t\twhere N is GRID size, EPSILON is Tolerance\n");
58 fprintf(stderr,"\t\texample N = 100, EPSILON = 0.1\n");
59 return -1;
60 }
61
62 N = M = atoi(argv[1]);
63 EPSILON = atof(argv[2]);
64
65 float *u = (float *)malloc(N * M * sizeof(float));
66 float *unew = (float *)malloc(N * M * sizeof(float));
67
68 if(u==0 || unew ==0)
69 {
70 perror("Can't allocated data\n");
71 return -1;
72 }
73
74 printf ( "\n" );
75 printf ( "HEATED_PLATE\n" );
76 printf ( " Parallel OpenMP version, using %d Threads\n",omp_get_max_threads() );
77 printf ( " A program to solve for the steady state temperature distribution\n" );
78 printf ( " over a rectangular plate.\n" );
79 printf ( " Spatial grid of %d by %d points.\n\n", M, N );
80
81
82 /* Initialize grid and create input file */
83 printf("Initializing grid\n");
84
85 inidat(N, M,u,unew);
86
87 prtdat(N, M,u, "initial.dat");
88
89 printf("Start computing\n");
90
91 int iter=0;
92
93 /*
94 * iterate until the new solution unew differs from the old solution u
95 * by no more than EPSILON.
96 * */
97
98 while(diff> EPSILON) {
99
100 update(N, M, u, unew,&diff);
101
102 if(iter%ITER_PRINT==0)
103 printf("Iteration %d, diff = %f\n ", iter,diff);
104
105 iter++;
106 }
107
108 prtdat(N, M, u, "final.dat");
109
110 free(u);
111 free(unew);
112 }
113
114
115
116 /****************************************************************************
117 * subroutine update
118 ****************************************************************************/
119 void update(int nx,int ny, float *u, float *unew, float * diff)
120 {
121 int ix, iy;
122 *diff=0.0;
123
124 #pragma omp parallel for shared(nx,ny,u,unew) private (ix,iy)
125 for (ix = 1; ix < nx-1; ix++) {
126 for (iy = 1; iy < ny-1; iy++) {
127 unew[ix*ny+iy] =
128 (u[(ix+1)*ny+iy] + u[(ix-1)*ny+iy] +
129 u[ix*ny+iy+1] + u[ix*ny+iy-1] )/4.0;
130
131 }
132 }
133
134 //compute reduction
135
136
137 float mydiff;
138
139 #pragma omp parallel shared(nx,ny,u,unew, diff) private (ix,iy,mydiff)
140 {
141 mydiff=0.0;
142 #pragma omp for
143 for (ix = 1; ix < nx-1; ix++) {
144 for (iy = 1; iy < ny-1; iy++) {
145 if (mydiff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] ))
146 {
147 mydiff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] );
148 }
149 }
150 }
151
152
153 # pragma omp critical
154 {
155 if (*diff < mydiff )
156 {
157 *diff = mydiff;
158 }
159 }
160
161
162 #pragma omp for
163 for (ix = 1; ix < nx-1; ix++) {
164 for (iy = 1; iy < ny-1; iy++) {
165 u[ix*ny+iy] = unew[ix*ny+iy];
166 }
167 }
168 }
169 }
170
171 /*****************************************************************************
172 * Initialize Data
173 *****************************************************************************/
174 void inidat(int nx, int ny, float *u, float *unew)
175 {
176 int ix, iy;
177
178 /*
179 *Set boundary data and interrior values
180 * */
181 for (ix = 0; ix < nx; ix++)
182 for (iy = 0; iy < ny; iy++) {
183
184 if(ix==0)
185 {
186 u[ix*ny+iy]=0.0;
187 }
188 else
189 if(iy==0 && ix!=0)
190 {
191 u[ix*ny+iy]=100.0;
192 }else
193
194 if(ix==nx-1)
195 {
196 u[ix*ny+iy]=100.0;
197 }else
198
199 if(iy==ny-1 && ix!=0)
200 {
201 u[ix*ny+iy]=100.0;
202 }else
203
204 u[ix*ny+iy]=0.0;
205 }
206 }
207
208 /**************************************************************************
209 * Print Data to files
210 **************************************************************************/
211 void prtdat(int nx, int ny, float *u,const char *fnam)
212 {
213
214 int ix, iy;
215 FILE *fp;
216
217 if(ITER_PRINT==0)return;
218
219 fp = fopen(fnam, "w");
220
221 for (ix = 0 ; ix < nx; ix++) {
222 for (iy =0; iy < ny; iy++) {
223
224 fprintf(fp, "%8.3f", u[ix*ny+iy]);
225
226 if(iy!=ny-1)
227 {
228 fprintf(fp, " ");
229 }else
230 {
231 fprintf(fp, "\n");
232 }
233 }
234 }
235
236 fclose(fp);
237 }
238
239
240
/**************************************************************************** 1 241
* DESCRIPTION: 2
* Serial HEAT2D Example - C Version 3
* This example is based on a simplified 4
* two-dimensional heat equation domain decomposition. The initial 5
* temperature is computed to be high in the middle of the domain and 6
* zero at the boundaries. The boundaries are held at zero throughout 7
* the simulation. During the time-stepping, an array containing two 8
* domains is used; these domains alternate between old data and new data. 9
* 10
* The physical region, and the boundary conditions, are suggested 11
by this diagram; 12
13
u = 0 14
+------------------+ 15
| | 16
u = 100 | | u = 100 17
| | 18
+------------------+ 19
u = 100 20
21
Interrior point : 22
u[Central] = (1/4) * ( u[North] + u[South] + u[East] + u[West] ) 23
24
****************************************************************************/ 25
#include <stdio.h> 26
#include <stdlib.h> 27
#include <math.h> 28
#include <omp.h> 29
30
#define NN 50 31
#define MM 50 32
33
#define ITER_PRINT 100 34
#define PRINT_DATA 1 35
36
#define _EPSILON 0.001 37
38
39
void update(int nx,int ny, float *u, float *unew, float * diff); 40
void inidat(int nx, int ny, float *u, float *unew); 41
void prtdat(int nx, int ny, float *u,const char *fnam); 42
43
44
45
46
int main(int argc, char *argv[]) 47
{ 48
49
float diff=1.0; 50
float EPSILON=_EPSILON; 51
int N=NN,M=MM; 52
53
if(argc !=3) 54
{ 55
fprintf(stderr,"usage %s N EPSILON\n ", argv[0]); 56
fprintf(stderr,"\t\twhere N is GRID size, EPSILON is Tolerance\n"); 57
fprintf(stderr,"\t\texample N = 100, EPSILON = 0.1\n"); 58
return -1; 59
} 60
61
N = M = atoi(argv[1]); 62
EPSILON = atof(argv[2]); 63
64
float *u = (float *)malloc(N * M * sizeof(float)); 65
float *unew = (float *)malloc(N * M * sizeof(float)); 66
67
if(u==0 || unew ==0) 68
{ 69
perror("Can't allocated data\n"); 70
return -1; 71
} 72
73
printf ( "\n" ); 74
printf ( "HEATED_PLATE\n" ); 75
printf ( " Parallel OpenMP version, using %d Threads\n",omp_get_max_threads() ); 76
printf ( " A program to solve for the steady state temperature distribution\n" ); 77
printf ( " over a rectangular plate.\n" ); 78
printf ( " Spatial grid of %d by %d points.\n\n", M, N ); 79
80
81
/* Initialize grid and create input file */ 82
printf("Initializing grid\n"); 83
84
inidat(N, M,u,unew); 85
86
prtdat(N, M,u, "initial.dat"); 87
88
printf("Start computing\n"); 89
90
int iter=0; 91
92
/* 93
* iterate until the new solution unew differs from the old solution u 94
* by no more than EPSILON. 95
* */ 96
97
while(diff> EPSILON) { 98
99
update(N, M, u, unew,&diff); 100
101
if(iter%ITER_PRINT==0) 102
printf("Iteration %d, diff = %f\n ", iter,diff); 103
104
iter++; 105
} 106
107
prtdat(N, M, u, "final.dat"); 108
109
free(u); 110
free(unew); 111
} 112
113
114
115
/**************************************************************************** 116
* subroutine update 117
****************************************************************************/ 118
void update(int nx,int ny, float *u, float *unew, float * diff) 119
{ 120
int ix, iy; 121
*diff=0.0; 122
123
#pragma omp parallel for shared(nx,ny,u,unew) private (ix,iy) 124
for (ix = 1; ix < nx-1; ix++) { 125
for (iy = 1; iy < ny-1; iy++) { 126
unew[ix*ny+iy] = 127
(u[(ix+1)*ny+iy] + u[(ix-1)*ny+iy] + 128
u[ix*ny+iy+1] + u[ix*ny+iy-1] )/4.0; 129
130
} 131
} 132
133
//compute reduction 134
135
136
float mydiff; 137
138
#pragma omp parallel shared(nx,ny,u,unew, diff) private (ix,iy,mydiff) 139
{ 140
mydiff=0.0; 141
#pragma omp for 142
for (ix = 1; ix < nx-1; ix++) { 143
for (iy = 1; iy < ny-1; iy++) { 144
if (mydiff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] )) 145
{ 146
mydiff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] ); 147
} 148
} 149
} 150
151
152
# pragma omp critical 153
{ 154
if (*diff < mydiff ) 155
{ 156
*diff = mydiff; 157
} 158
} 159
160
161
#pragma omp for 162
for (ix = 1; ix < nx-1; ix++) { 163
for (iy = 1; iy < ny-1; iy++) { 164
u[ix*ny+iy] = unew[ix*ny+iy]; 165
} 166
} 167
} 168
} 169
170
/***************************************************************************** 171
* Initialize Data 172
*****************************************************************************/ 173
void inidat(int nx, int ny, float *u, float *unew) 174
{ 175
int ix, iy; 176
177
/* 178
*Set boundary data and interrior values 179
lab3/plot.py View file @ 2b4ce20
File was created 1 #!/usr/bin/python
2 #-*-coding:utf8-*-
3 from sys import argv
4 from pylab import *
5
6 if __name__ == "__main__":
7 if len(argv) != 2:
8 print "Usage: python", argv[0], " matrix.dat"
9 exit (0)
10
11
12 f = open ( argv[1] , 'r')
13 A = [ map(float,line.split()) for line in f ]
14
15 figure(1)
16 imshow(A, interpolation='nearest')
17 grid(True)
18 show()
#!/usr/bin/python 1 19
#-*-coding:utf8-*- 2
from sys import argv 3
from pylab import * 4
5
if __name__ == "__main__": 6
if len(argv) != 2: 7
print "Usage: python", argv[0], " matrix.dat" 8
exit (0) 9
10
11
f = open ( argv[1] , 'r') 12
A = [ map(float,line.split()) for line in f ] 13
14
lab3/ser_heat2D.c View file @ 2b4ce20
File was created 1 /****************************************************************************
2 * DESCRIPTION:
3 * Serial HEAT2D Example - C Version
4 * This example is based on a simplified
5 * two-dimensional heat equation domain decomposition. The initial
6 * temperature is computed to be high in the middle of the domain and
7 * zero at the boundaries. The boundaries are held at zero throughout
8 * the simulation. During the time-stepping, an array containing two
9 * domains is used; these domains alternate between old data and new data.
10 *
11 * The physical region, and the boundary conditions, are suggested
12 by this diagram;
13
14 u = 0
15 +------------------+
16 | |
17 u = 100 | | u = 100
18 | |
19 +------------------+
20 u = 100
21
22 Interrior point :
23 u[Central] = (1/4) * ( u[North] + u[South] + u[East] + u[West] )
24
25 ****************************************************************************/
26 #include <stdio.h>
27 #include <stdlib.h>
28 #include <math.h>
29
30 #define NN 50
31 #define MM 50
32
33 #define ITER_PRINT 100
34 #define PRINT_DATA 1
35
36 #define _EPSILON 0.01
37
38
39 void update(int nx,int ny, float *u, float *unew, float * diff);
40 void inidat(int nx, int ny, float *u, float *unew);
41 void prtdat(int nx, int ny, float *u,const char *fnam);
42
43
44
45
46 int main(int argc, char *argv[])
47 {
48 int N=NN,M=MM;
49
50 float EPSILON=_EPSILON;
51
52 if(argc !=3)
53 {
54 fprintf(stderr,"usage %s N EPSILON\n ", argv[0]);
55 fprintf(stderr,"\t\twhere N is GRID size, EPSILON is Tolerance\n");
56 fprintf(stderr,"\t\texample N = 100, EPSILON = 0.1\n");
57 return -1;
58 }
59
60 N = M = atoi(argv[1]);
61 EPSILON = atof(argv[2]);
62
63 float diff=1.0;
64
65 float *u = (float *)malloc(N * M * sizeof(float));
66 float *unew = (float *)malloc(N * M * sizeof(float));
67
68 if(u==0 || unew ==0)
69 {
70 perror("Can't allocated data\n");
71 return -1;
72 }
73
74 printf ( "\n" );
75 printf ( "HEATED_PLATE\n" );
76 printf ( " Serial version\n" );
77 printf ( " A program to solve for the steady state temperature distribution\n" );
78 printf ( " over a rectangular plate.\n" );
79 printf ( "\n" );
80 printf ( " Spatial grid of %d by %d points.\n", M, N );
81 printf ( " The iteration will end until tolerance <= %f\n\n",EPSILON);
82
83 /* Initialize grid and create input file */
84 printf("Initializing grid\n");
85
86 inidat(N, M,u,unew);
87
88 prtdat(N, M,u, "initial.dat");
89
90 printf("Start computing\n\n");
91
92 int iter=0;
93
94 /*
95 * iterate until the new solution unew differs from the old solution u
96 * by no more than EPSILON.
97 * */
98
99 while(diff> EPSILON) {
100
101 update(N, M, u, unew,&diff);
102
103 if(iter%ITER_PRINT==0)
104
105 printf("Iteration %d, diff = %f\n ", iter,diff);
106
107 iter++;
108 }
109
110 prtdat(N, M, u, "final.dat");
111
112 free(u);
113 free(unew);
114 }
115
116
117
118 /****************************************************************************
119 * subroutine update
120 ****************************************************************************/
121 void update(int nx,int ny, float *u, float *unew, float * diff)
122 {
123 int ix, iy;
124 *diff=0.0;
125
126 for (ix = 1; ix < nx-1; ix++) {
127 for (iy = 1; iy < ny-1; iy++) {
128 unew[ix*ny+iy] = (u[(ix+1)*ny+iy] + u[(ix-1)*ny+iy] + u[ix*ny+iy+1] + u[ix*ny+iy-1] )/4.0
129 ;
130 if (*diff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] ))
131 {
132 *diff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] );
133 }
134 }
135
136 }
137
138
139 for (ix = 1; ix < nx-1; ix++) {
140 for (iy = 1; iy < ny-1; iy++) {
141 u[ix*ny+iy] = unew[ix*ny+iy];
142 }
143 }
144
145 }
146
147 /*****************************************************************************
148 * Initialize Data
149 *****************************************************************************/
150 void inidat(int nx, int ny, float *u, float *unew)
151 {
152 int ix, iy;
153
154 /*
155 *Set boundary data and interrior values
156 * */
157 for (ix = 0; ix < nx; ix++)
158 for (iy = 0; iy < ny; iy++) {
159
160 if(ix==0)
161 {
162 u[ix*ny+iy]=0.0;
163 }
164 else
165 if(iy==0 && ix!=0)
166 {
167 u[ix*ny+iy]=100.0;
168 }else
169
170 if(ix==nx-1)
171 {
172 u[ix*ny+iy]=100.0;
173 }else
174
175 if(iy==ny-1 && ix!=0)
176 {
177 u[ix*ny+iy]=100.0;
178 }else
179
180 u[ix*ny+iy]=0.0;
181 }
182 }
183
184 /**************************************************************************
185 * Print Data to files
186 **************************************************************************/
187 void prtdat(int nx, int ny, float *u,const char *fname)
188 {
189
190 int ix, iy;
191 FILE *fp;
192
193 if(ITER_PRINT==0)return;
194
195 fp = fopen(fname, "w");
196
197 // fprintf ( fp, "%d\n", M );
198 // fprintf ( fp, "%d\n", N );
199
200 for (ix = 0 ; ix < nx; ix++) {
201 for (iy =0; iy < ny; iy++) {
202
203 fprintf(fp, "%6.2f ", u[ix*ny+iy]);
204 }
205 fputc ( '\n', fp);
206 }
207
208 printf (" Data written to the output file %s\n", fname);
209 fclose(fp);
210 }
211
212
213
/**************************************************************************** 1 214
* DESCRIPTION: 2
* Serial HEAT2D Example - C Version 3
* This example is based on a simplified 4
* two-dimensional heat equation domain decomposition. The initial 5
* temperature is computed to be high in the middle of the domain and 6
* zero at the boundaries. The boundaries are held at zero throughout 7
* the simulation. During the time-stepping, an array containing two 8
* domains is used; these domains alternate between old data and new data. 9
* 10
* The physical region, and the boundary conditions, are suggested 11
by this diagram; 12
13
u = 0 14
+------------------+ 15
| | 16
u = 100 | | u = 100 17
| | 18
+------------------+ 19
u = 100 20
21
Interrior point : 22
u[Central] = (1/4) * ( u[North] + u[South] + u[East] + u[West] ) 23
24
****************************************************************************/ 25
#include <stdio.h> 26
#include <stdlib.h> 27
#include <math.h> 28
29
#define NN 50 30
#define MM 50 31
32
#define ITER_PRINT 100 33
#define PRINT_DATA 1 34
35
#define _EPSILON 0.01 36
37
38
void update(int nx,int ny, float *u, float *unew, float * diff); 39
void inidat(int nx, int ny, float *u, float *unew); 40
void prtdat(int nx, int ny, float *u,const char *fnam); 41
42
43
44
45
int main(int argc, char *argv[]) 46
{ 47
int N=NN,M=MM; 48
49
float EPSILON=_EPSILON; 50
51
if(argc !=3) 52
{ 53
fprintf(stderr,"usage %s N EPSILON\n ", argv[0]); 54
fprintf(stderr,"\t\twhere N is GRID size, EPSILON is Tolerance\n"); 55
fprintf(stderr,"\t\texample N = 100, EPSILON = 0.1\n"); 56
return -1; 57
} 58
59
N = M = atoi(argv[1]); 60
EPSILON = atof(argv[2]); 61
62
float diff=1.0; 63
64
float *u = (float *)malloc(N * M * sizeof(float)); 65
float *unew = (float *)malloc(N * M * sizeof(float)); 66
67
if(u==0 || unew ==0) 68
{ 69
perror("Can't allocated data\n"); 70
return -1; 71
} 72
73
printf ( "\n" ); 74
printf ( "HEATED_PLATE\n" ); 75
printf ( " Serial version\n" ); 76
printf ( " A program to solve for the steady state temperature distribution\n" ); 77
printf ( " over a rectangular plate.\n" ); 78
printf ( "\n" ); 79
printf ( " Spatial grid of %d by %d points.\n", M, N ); 80
printf ( " The iteration will end until tolerance <= %f\n\n",EPSILON); 81
82
/* Initialize grid and create input file */ 83
printf("Initializing grid\n"); 84
85
inidat(N, M,u,unew); 86
87
prtdat(N, M,u, "initial.dat"); 88
89
printf("Start computing\n\n"); 90
91
int iter=0; 92
93
/* 94
* iterate until the new solution unew differs from the old solution u 95
* by no more than EPSILON. 96
* */ 97
98
while(diff> EPSILON) { 99
100
update(N, M, u, unew,&diff); 101
102
if(iter%ITER_PRINT==0) 103
104
printf("Iteration %d, diff = %f\n ", iter,diff); 105
106
iter++; 107
} 108
109
prtdat(N, M, u, "final.dat"); 110
111
free(u); 112
free(unew); 113
} 114
115
116
117
/**************************************************************************** 118
* subroutine update 119
****************************************************************************/ 120
void update(int nx,int ny, float *u, float *unew, float * diff) 121
{ 122
int ix, iy; 123
*diff=0.0; 124
125
for (ix = 1; ix < nx-1; ix++) { 126
for (iy = 1; iy < ny-1; iy++) { 127
unew[ix*ny+iy] = (u[(ix+1)*ny+iy] + u[(ix-1)*ny+iy] + u[ix*ny+iy+1] + u[ix*ny+iy-1] )/4.0 128
; 129
if (*diff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] )) 130
{ 131
*diff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] ); 132
} 133
} 134
135
} 136
137
138
for (ix = 1; ix < nx-1; ix++) { 139
for (iy = 1; iy < ny-1; iy++) { 140
u[ix*ny+iy] = unew[ix*ny+iy]; 141
} 142
} 143
144
} 145
146
/***************************************************************************** 147
* Initialize Data 148
*****************************************************************************/ 149
void inidat(int nx, int ny, float *u, float *unew) 150
{ 151
int ix, iy; 152
153
/* 154
*Set boundary data and interrior values 155
* */ 156