Skip to content

S
supercomputer_lab
Commit 099699b9 authored by Gaurav Kukreja's avatar Gaurav Kukreja
Browse files
Options

Sequential Optimization of Gauss 

 * residual interleaving
Signed-off-by: 's avatarGaurav Kukreja <mailme.gaurav@gmail.com>
parent c083da69
Hide whitespace changes
Inline Side-by-side
Showing with 152 additions and 141 deletions
heat-seq/Makefile
View file @ 099699b9
# Intel compiler
CC = icc
CFLAGS = -g
CFLAGS = -g -O3
MPICC = mpicc
......
heat-seq/heat.c
View file @ 099699b9
......@@ -16,162 +16,164 @@
void usage( char *s )
{
fprintf(stderr,
"Usage: %s <input file> [result file]\n\n", s);
fprintf(stderr,
"Usage: %s <input file> [result file]\n\n", s);
}
int main( int argc, char *argv[] )
{
unsigned iter;
FILE *infile, *resfile;
char *resfilename;
unsigned iter;
FILE *infile, *resfile;
char *resfilename;
int Interleaving_Count;
double *tmp;
// algorithmic parameters
algoparam_t param;
int np;
double runtime, flop;
double residual;
// check arguments
if( argc < 2 )
{
usage( argv[0] );
return 1;
}
// check input file
if( !(infile=fopen(argv[1], "r")) )
{
fprintf(stderr,
"\nError: Cannot open \"%s\" for reading.\n\n", argv[1]);
usage(argv[0]);
return 1;
}
// check result file
resfilename= (argc>=3) ? argv[2]:"heat.ppm";
if( !(resfile=fopen(resfilename, "w")) )
{
fprintf(stderr,
"\nError: Cannot open \"%s\" for writing.\n\n",
resfilename);
usage(argv[0]);
return 1;
}
// check input
if( !read_input(infile, &param) )
{
fprintf(stderr, "\nError: Error parsing input file.\n\n");
usage(argv[0]);
return 1;
}
// algorithmic parameters
algoparam_t param;
int np;
print_params(&param);
double runtime, flop;
double residual;
// set the visualization resolution
param.visres = param.max_res;
// check arguments
if( argc < 2 )
{
usage( argv[0] );
return 1;
}
param.u = 0;
param.uhelp = 0;
param.uvis = 0;
// check input file
if( !(infile=fopen(argv[1], "r")) )
{
fprintf(stderr,
"\nError: Cannot open \"%s\" for reading.\n\n", argv[1]);
param.act_res = param.initial_res;
usage(argv[0]);
return 1;
}
// loop over different resolutions
while(1) {
// free allocated memory of previous experiment
if (param.u != 0)
finalize(&param);
// check result file
resfilename= (argc>=3) ? argv[2]:"heat.ppm";
if( !initialize(&param) )
if( !(resfile=fopen(resfilename, "w")) )
{
fprintf(stderr, "Error in Jacobi initialization.\n\n");
fprintf(stderr,
"\nError: Cannot open \"%s\" for writing.\n\n",
resfilename);
usage(argv[0]);
usage(argv[0]);
return 1;
}
fprintf(stderr, "Resolution: %5u\r", param.act_res);
// check input
if( !read_input(infile, &param) )
{
fprintf(stderr, "\nError: Error parsing input file.\n\n");
usage(argv[0]);
return 1;
}
print_params(&param);
// set the visualization resolution
param.visres = param.max_res;
param.u = 0;
param.uhelp = 0;
param.uvis = 0;
// full size (param.act_res are only the inner points)
np = param.act_res + 2;
// starting time
runtime = wtime();
residual = 999999999;
param.act_res = param.initial_res;
iter = 0;
// loop over different resolutions
while(1) {
switch( param.algorithm ) {
// free allocated memory of previous experiment
if (param.u != 0)
finalize(&param);
case 0: // JACOBI
Interleaving_Count = INTERLEAVING_COUNT;
if( !initialize(&param) )
{
fprintf(stderr, "Error in Jacobi initialization.\n\n");
residual = relax_jacobi_return_residual(param.u, param.uhelp, np, np, Interleaving_Count);
tmp=param.u;
param.u=param.uhelp;
param.uhelp=tmp;
break;
usage(argv[0]);
}
case 1: // GAUSS
fprintf(stderr, "Resolution: %5u\r", param.act_res);
relax_gauss(param.u, np, np);
residual = residual_gauss( param.u, param.uhelp, np, np);
break;
}
iter = iter + Interleaving_Count;
// full size (param.act_res are only the inner points)
np = param.act_res + 2;
// solution good enough ?
if (residual < 0.000005) break;
// starting time
runtime = wtime();
residual = 999999999;
// max. iteration reached ? (no limit with maxiter=0)
if (param.maxiter>0 && iter>=param.maxiter) break;
if (iter % 100 == 0)
fprintf(stderr, "residual %f, %d iterations\n", residual, iter);
}
iter = 0;
while(1) {
switch( param.algorithm ) {
case 0: // JACOBI
Interleaving_Count = INTERLEAVING_COUNT;
residual = relax_jacobi_return_residual(param.u, param.uhelp, np, np, Interleaving_Count);
tmp=param.u;
param.u=param.uhelp;
param.uhelp=tmp;
break;
case 1: // GAUSS
// Flop count after <i> iterations
flop = iter * 11.0 * param.act_res * param.act_res;
// stopping time
runtime = wtime() - runtime;
residual = relax_gauss_return_residual(param.u, np, np);
// residual = residual_gauss( param.u, param.uhelp, np, np);
break;
}
if(param.algorithm == 0)
iter = iter + Interleaving_Count;
else
iter++;
fprintf(stderr, "Resolution: %5u, ", param.act_res);
fprintf(stderr, "Time: %04.3f ", runtime);
fprintf(stderr, "(%3.3f GFlop => %6.2f MFlop/s, ",
flop/1000000000.0,
flop/runtime/1000000);
fprintf(stderr, "residual %f, %d iterations)\n", residual, iter);
// solution good enough ?
if (residual < 0.000005) break;
// for plot...
printf("%5d %f\n", param.act_res, flop/runtime/1000000);
// max. iteration reached ? (no limit with maxiter=0)
if (param.maxiter>0 && iter>=param.maxiter) break;
if (iter % 100 == 0)
fprintf(stderr, "residual %f, %d iterations\n", residual, iter);
}
// Flop count after <i> iterations
flop = iter * 7.0 * param.act_res * param.act_res;
// stopping time
runtime = wtime() - runtime;
fprintf(stderr, "Resolution: %5u, ", param.act_res);
fprintf(stderr, "Time: %04.3f ", runtime);
fprintf(stderr, "(%3.3f GFlop => %6.2f MFlop/s, ",
flop/1000000000.0,
flop/runtime/1000000);
fprintf(stderr, "residual %f, %d iterations)\n", residual, iter);
// for plot...
printf("%5d %f\n", param.act_res, flop/runtime/1000000);
if (param.act_res + param.res_step_size > param.max_res) break;
param.act_res += param.res_step_size;
}
if (param.act_res + param.res_step_size > param.max_res) break;
param.act_res += param.res_step_size;
}
coarsen( param.u, np, np,
param.uvis, param.visres+2, param.visres+2 );
coarsen( param.u, np, np,
param.uvis, param.visres+2, param.visres+2 );
write_image( resfile, param.uvis,
param.visres+2,
param.visres+2 );
write_image( resfile, param.uvis,
param.visres+2,
param.visres+2 );
finalize( &param );
finalize( &param );
return 0;
return 0;
}
heat-seq/heat.h
View file @ 099699b9
......@@ -55,12 +55,14 @@ int coarsen(double *uold, unsigned oldx, unsigned oldy ,
double residual_gauss( double *u, double *utmp,
unsigned sizex, unsigned sizey );
#endif
void relax_gauss( double *u,
double relax_gauss_return_residual( double *u,
unsigned sizex, unsigned sizey );
// Jacobi: relax_jacobi.c
#if 0
double residual_jacobi( double *u,
unsigned sizex, unsigned sizey );
#endif
double relax_jacobi_return_residual( double *u, double *utmp,
unsigned sizex, unsigned sizey, int Interleaving_Count );
......
heat-seq/relax_gauss.c
View file @ 099699b9
......@@ -7,6 +7,8 @@
#include "heat.h"
#if 0
/*
* Residual (length of error vector)
* between current solution and next after a Gauss-Seidel step
......@@ -47,27 +49,32 @@ double residual_gauss( double *u, double *utmp,
return sum;
}
#endif
/*
* One Gauss-Seidel iteration step
*
* Flop count in inner body is 4
*/
void relax_gauss( double *u,
unsigned sizex, unsigned sizey )
double relax_gauss_return_residual( double *u,
unsigned sizex, unsigned sizey )
{
unsigned i, j;
unsigned i, j;
double temp, diff, residual;
residual =0;
for( i=1; i<sizey-1; i++ )
{
for( j=1; j<sizex-1; j++ )
for( i=1; i<sizey-1; i++ )
{
u[i*sizex+j]= 0.25 * (u[ i*sizex + (j-1) ]+
u[ i*sizex + (j+1) ]+
u[ (i-1)*sizex + j ]+
u[ (i+1)*sizex + j ]);
for( j=1; j<sizex-1; j++ )
{
temp = u[i*sizex+j];
u[i*sizex+j]= 0.25 * (u[ i*sizex + (j-1) ]+
u[ i*sizex + (j+1) ]+
u[ (i-1)*sizex + j ]+
u[ (i+1)*sizex + j ]);
diff = u[i*sizex+j] - temp;
residual += diff*diff;
}
}
}
return residual;
}
heat-seq/test.dat
View file @ 099699b9
50 # iterations
100 # initial resolution
2900 # max resolution (spatial resolution)
200 # resolution step size
0 # Algorithm 0=Jacobi 1=Gauss
1024 # iterations
256 # initial resolution
4096 # max resolution (spatial resolution)
3840 # resolution step size
1 # Algorithm 0=Jacobi 1=Gauss
2 # number of heat sources
0.0 0.0 1.0 1.0 # (x,y), size temperature
1.0 1.0 1.0 0.5
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment