diff --git a/lab3/omp_heat2D.c b/lab3/omp_heat2D.c index b4d9ed9..2744257 100644 --- a/lab3/omp_heat2D.c +++ b/lab3/omp_heat2D.c @@ -28,13 +28,13 @@ Interrior point : #include #include -#define N 500 -#define M 500 +#define NN 50 +#define MM 50 #define ITER_PRINT 100 #define PRINT_DATA 1 -#define EPSILON 1e-1 +#define _EPSILON 0.001 void update(int nx,int ny, float *u, float *unew, float * diff); @@ -48,32 +48,44 @@ int main(int argc, char *argv[]) { float diff=1.0; + float EPSILON=_EPSILON; + int N=NN,M=MM; + + if(argc !=3) + { + fprintf(stderr,"usage %s N EPSILON\n ", argv[0]); + fprintf(stderr,"\t\twhere N is GRID size, EPSILON is Tolerance\n"); + fprintf(stderr,"\t\texample N = 100, EPSILON = 0.1\n"); + return -1; + } + + N = M = atoi(argv[1]); + EPSILON = atof(argv[2]); float *u = (float *)malloc(N * M * sizeof(float)); float *unew = (float *)malloc(N * M * sizeof(float)); - + if(u==0 || unew ==0) { perror("Can't allocated data\n"); return -1; } - printf ( "\n" ); - printf ( "HEATED_PLATE\n" ); - printf ( " Parallel OpenMP version, using %d Threads\n",omp_get_max_threads() ); - printf ( " A program to solve for the steady state temperature distribution\n" ); - printf ( " over a rectangular plate.\n" ); - printf ( " Spatial grid of %d by %d points.\n\n", M, N ); + printf ( "\n" ); + printf ( "HEATED_PLATE\n" ); + printf ( " Parallel OpenMP version, using %d Threads\n",omp_get_max_threads() ); + printf ( " A program to solve for the steady state temperature distribution\n" ); + printf ( " over a rectangular plate.\n" ); + printf ( " Spatial grid of %d by %d points.\n\n", M, N ); /* Initialize grid and create input file */ printf("Initializing grid\n"); - + inidat(N, M,u,unew); prtdat(N, M,u, "initial.dat"); - - + printf("Start computing\n"); int iter=0; @@ -82,19 +94,19 @@ int main(int argc, char *argv[]) * iterate until the new solution unew differs from the old solution u * by no more than EPSILON. * */ - + while(diff> EPSILON) { update(N, M, u, unew,&diff); - + if(iter%ITER_PRINT==0) - printf("Iteration %d, diff = %f\n ", iter,diff); + printf("Iteration %d, diff = %f\n ", iter,diff); iter++; } prtdat(N, M, u, "final.dat"); - + free(u); free(unew); } @@ -115,49 +127,45 @@ void update(int nx,int ny, float *u, float *unew, float * diff) 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; - + } } -//compute reduction + //compute reduction - float mydiff; - + float mydiff; + #pragma omp parallel shared(nx,ny,u,unew, diff) private (ix,iy,mydiff) - { -mydiff=0.0; + { + mydiff=0.0; #pragma omp for - for (ix = 1; ix < nx-1; ix++) { - for (iy = 1; iy < ny-1; iy++) { - if (mydiff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] )) - { - mydiff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] ); + for (ix = 1; ix < nx-1; ix++) { + for (iy = 1; iy < ny-1; iy++) { + if (mydiff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] )) + { + mydiff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] ); + } } } - } # pragma omp critical - { - if (*diff < mydiff ) - { - *diff = mydiff; - } - } - - -#pragma omp for - for (ix = 1; ix < nx-1; ix++) { - for (iy = 1; iy < ny-1; iy++) { - u[ix*ny+iy] = unew[ix*ny+iy]; + { + if (*diff < mydiff ) + { + *diff = mydiff; + } } - } - - - } +#pragma omp for + for (ix = 1; ix < nx-1; ix++) { + for (iy = 1; iy < ny-1; iy++) { + u[ix*ny+iy] = unew[ix*ny+iy]; + } + } + } } /***************************************************************************** @@ -178,22 +186,22 @@ void inidat(int nx, int ny, float *u, float *unew) u[ix*ny+iy]=0.0; } else - if(iy==0 && ix!=0) - { - u[ix*ny+iy]=100.0; - }else - - if(ix==nx-1) - { - u[ix*ny+iy]=100.0; - }else - - if(iy==ny-1 && ix!=0) - { - u[ix*ny+iy]=100.0; - }else - - u[ix*ny+iy]=( float ) ( 2 * nx + 2 * ny - 4 ); + if(iy==0 && ix!=0) + { + u[ix*ny+iy]=100.0; + }else + + if(ix==nx-1) + { + u[ix*ny+iy]=100.0; + }else + + if(iy==ny-1 && ix!=0) + { + u[ix*ny+iy]=100.0; + }else + + u[ix*ny+iy]=0.0; } } @@ -202,12 +210,12 @@ void inidat(int nx, int ny, float *u, float *unew) **************************************************************************/ void prtdat(int nx, int ny, float *u,const char *fnam) { - + int ix, iy; FILE *fp; if(ITER_PRINT==0)return; - + fp = fopen(fnam, "w"); for (ix = 0 ; ix < nx; ix++) { diff --git a/lab3/ser_heat2D.c b/lab3/ser_heat2D.c index bb14733..176bd1c 100644 --- a/lab3/ser_heat2D.c +++ b/lab3/ser_heat2D.c @@ -27,13 +27,13 @@ Interrior point : #include #include -#define N 200 -#define M 200 +#define NN 50 +#define MM 50 #define ITER_PRINT 100 #define PRINT_DATA 1 -#define EPSILON 1e-1 +#define _EPSILON 0.01 void update(int nx,int ny, float *u, float *unew, float * diff); @@ -45,6 +45,20 @@ void prtdat(int nx, int ny, float *u,const char *fnam); int main(int argc, char *argv[]) { + int N=NN,M=MM; + + float EPSILON=_EPSILON; + + if(argc !=3) + { + fprintf(stderr,"usage %s N EPSILON\n ", argv[0]); + fprintf(stderr,"\t\twhere N is GRID size, EPSILON is Tolerance\n"); + fprintf(stderr,"\t\texample N = 100, EPSILON = 0.1\n"); + return -1; + } + + N = M = atoi(argv[1]); + EPSILON = atof(argv[2]); float diff=1.0; @@ -63,8 +77,8 @@ int main(int argc, char *argv[]) printf ( " A program to solve for the steady state temperature distribution\n" ); printf ( " over a rectangular plate.\n" ); printf ( "\n" ); - printf ( " Spatial grid of %d by %d points.\n\n", M, N ); - + printf ( " Spatial grid of %d by %d points.\n", M, N ); + printf ( " The iteration will end until tolerance <= %f\n\n",EPSILON); /* Initialize grid and create input file */ printf("Initializing grid\n"); @@ -72,9 +86,8 @@ int main(int argc, char *argv[]) inidat(N, M,u,unew); prtdat(N, M,u, "initial.dat"); - - printf("Start computing\n"); + printf("Start computing\n\n"); int iter=0; @@ -88,7 +101,8 @@ int main(int argc, char *argv[]) update(N, M, u, unew,&diff); if(iter%ITER_PRINT==0) - printf("Iteration %d, diff = %f\n ", iter,diff); + + printf("Iteration %d, diff = %f\n ", iter,diff); iter++; } @@ -111,10 +125,8 @@ void update(int nx,int ny, float *u, float *unew, float * diff) for (ix = 1; ix < nx-1; ix++) { for (iy = 1; iy < ny-1; iy++) { - 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; - + 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 + ; if (*diff < fabs (unew[ix*ny+iy] - u[ix*ny+iy] )) { *diff = fabs ( unew[ix*ny+iy] - u[ix*ny+iy] ); @@ -152,27 +164,27 @@ void inidat(int nx, int ny, float *u, float *unew) else if(iy==0 && ix!=0) { - u[ix*ny+iy]=100.0; + u[ix*ny+iy]=100.0; }else if(ix==nx-1) { - u[ix*ny+iy]=100.0; + u[ix*ny+iy]=100.0; }else if(iy==ny-1 && ix!=0) { - u[ix*ny+iy]=100.0; + u[ix*ny+iy]=100.0; }else - u[ix*ny+iy]=( float ) ( 2 * nx + 2 * ny - 4 ); + u[ix*ny+iy]=0.0; } } /************************************************************************** * Print Data to files **************************************************************************/ -void prtdat(int nx, int ny, float *u,const char *fnam) +void prtdat(int nx, int ny, float *u,const char *fname) { int ix, iy; @@ -180,23 +192,20 @@ void prtdat(int nx, int ny, float *u,const char *fnam) if(ITER_PRINT==0)return; - fp = fopen(fnam, "w"); + fp = fopen(fname, "w"); + + // fprintf ( fp, "%d\n", M ); + // fprintf ( fp, "%d\n", N ); for (ix = 0 ; ix < nx; ix++) { for (iy =0; iy < ny; iy++) { - fprintf(fp, "%8.3f", u[ix*ny+iy]); - - if(iy!=ny-1) - { - fprintf(fp, " "); - }else - { - fprintf(fp, "\n"); - } + fprintf(fp, "%6.2f ", u[ix*ny+iy]); } + fputc ( '\n', fp); } + printf (" Data written to the output file %s\n", fname); fclose(fp); }