changes

README.md

+# Heatmap
+
+Das Programm verwendet zum Parallelisieren SPMD und Threads.
+
+## Datendekomposition
+
+Jeder Bildpunkt braucht seine Nachbarn zur iterativen Berechnung, also müssen die Strukturen bei dem "verteilten Speicher" Modell möglichst gut erhalten bleiben. Bei der threadparallelen Verarbeitung
+
+## SPMD
+
+Das Programm läuft mehrfach. Die Prozesse kommunizieren dabei über Nachrichten. Jede Instanz kann sich selbst eindeutig identifizieren. Dies kann genutzt werden und Datenmengen (gleichmäßig) zuzuordnen. 
\ No newline at end of file

args.c

+
+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+
+void process_args(int argc, char **argv, int *m, int *n, double *eps, double *delta_t) {
+	if ((argc < 3) || (argc > 5)) {
+		fprintf(stderr, "Usage: heat <m> <n> [<epsilon> [delta_t]]!\n\n");
+		fprintf(stderr, "   <m>       -- X-Size of matrix\n");
+		fprintf(stderr, "   <n>       -- Y-Size of matrix\n");
+		fprintf(stderr, "   <epsilon> -- accuracy parameter\n");
+		fprintf(stderr, "   <delta_t> -- time step\n");
+		exit(1);
+	}
+
+	/*
+	** Erstes und zweites Argument: Größe der Matrix
+	*/
+	*m = atoi(argv[1]);
+	if ((*m < 3) || (*m > 6000)) {
+		fprintf(stderr, "Error: size m out of range [3 .. 6000] !\n");
+		exit(1);
+	}
+	*n = atoi(argv[2]);
+	if ((*n < 3) || (*n > 6000)) {
+		fprintf(stderr, "Error: size n out of range [3 .. 6000] !\n");
+		exit(1);
+	}
+	/*
+	** Drittes (optionales) Argument: "Genauigkeitsfaktor" eps
+	*/
+	if (argc >= 4) {
+		*eps = atof(argv[3]);
+	}
+	if ((*eps <= 0) || (*eps > 0.01)){
+		fprintf(stderr, "Error: epsilon must be between 0.0000000001 and 0.01 )\n");
+		exit(1);
+	}
+	/*
+	** Viertes (optionales) Argument: delta_t
+	*/
+	if (argc >= 5) {
+		*delta_t = atof(argv[4]);
+	}
+	if ((*delta_t <= 0) || (*delta_t > 1.0)){
+		fprintf(stderr, "Error: delta_t (%f) must be between 0.0000000001 and 1.0\n", *delta_t);
+		exit(1);
+	}
+}
\ No newline at end of file

matrix.c

+#include <stdio.h>
+#include <stdlib.h>
+#include <math.h>
+#include <mpi.h> 
+#include <sys/time.h>
+
+/* HILFSFUNKTIONEN ******************************************************* */
+
+/*
+** Allokiere ein 2D-Array mit m * n Elementen (m Zeilen mit n Spalten).
+** Also so etwas wie ``new double[m][n]'' ...
+*/
+double ** New_Matrix(int m, int n)
+{
+	double **res;
+	double *array;
+	int i;
+
+	/*
+	** Ergebnis ist ein 1D-Array mit Zeigern auf die Zeilen des 2D-Arrays
+	*/
+	res = (double **)malloc(sizeof(double *) * m);
+	if (res == NULL)
+		return NULL;
+
+	/*
+	** Das eigentliche Array wird aus Effizienzgr�nden in einem St�ck
+	** allokiert. Das entspricht dem Speicher-Layout, den es auch bei
+	** statischer Allokation h�tte.
+	*/
+	array = (double *)malloc(sizeof(double) * m * n);
+	if (array == NULL)
+		return NULL;
+
+	/*
+	** Jeder Zeilen-Zeiger wird nun mit dem richtigen Wert besetzt.
+	** (Zeile i beginnt bei Eintrag i*n)
+	*/
+	for (i = 0; i < m; i++) {
+		res[i] = &array[i * n];
+	}
+	return res;
+}
+
+/*
+** Freigabe eines 2D-Arrays.
+*/
+void Delete_Matrix(double **matrix)
+{
+	free(matrix[0]);
+	free(matrix);
+}
+
+/*
+** Gibt die mxn-Matrix 'a' in die Datei 'Matrix.txt' aus.
+*/
+void Write_Matrix(double **a, int m, int n)
+{
+	int i, j;
+	FILE *f;
+
+	/* Datei zum Schreiben �ffnen */
+	f = fopen("Matrix.txt", "w");
+	if (f == NULL) {
+		fprintf(stderr, "Cannot open file 'Matrix.txt' for writing!");
+		exit(1);
+	}
+
+	/* Gr��e der Matrix schreiben */
+	fprintf(f, "%d\n", m);	
+	fprintf(f, "%d\n\n", n);
+
+	/* Matrix-Elemente zeilenweise in die Datei schreiben */
+	for (i = 0; i < m; i++) {
+		for (j = 0; j < n; j++) {
+			fprintf(f, "%.10lf\n", a[i][j]);
+		}
+		fprintf(f, "\n");
+	}
+
+	/* Datei schlie�en */
+	fclose(f);
+}
+
+
+/*
+** Matrix initialisieren
+*/
+void Init_Matrix(double **a, int m, int n, int init_value) {
+	int i,j;
+	for (i=0; i<m; i++) {
+		for (j=0; j<n; j++) {
+			a[i][j] = init_value;
+		}
+	}
+}
+
+void Insert_Matrix(double **a, int a_dim0, int a_dim1, int pos_dim0, int pos_dim1, double **b, int b_dim0, int b_dim1, int offset_dim0_left, int offset_dim0_right, int offset_dim1_left, int offset_dim1_right) {
+	int i,j;
+	for(i=0; i < b_dim0 - offset_dim0_left - offset_dim0_right; i++) {
+		for(j=0; j < b_dim1 - offset_dim1_left - offset_dim1_right; j++) {
+			a[pos_dim0 + i][pos_dim0 + j] = b[offset_dim0_left + i][offset_dim1_right + j];
+		}
+	}
+}
\ No newline at end of file

proj.c

@@ -4,72 +4,35 @@
 #include <math.h>
 #include <mpi.h> 
 #include <sys/time.h>
+#include "project.h"
 
 double eps = 0.001;
 double delta_t = 0.000001;
 double alpha = 1;
 
-typedef struct {
-	int start_m;	// first tile to process in X-axis
-	int start_n;	// first tile to process in Y-axis
-	int end_m;	// last tile to process in X-axis
-	int end_n;	// last tile to process in Y-axis
-} t_process_info;
-
 int main(int argc, char **argv)
 {
+
 	MPI_Init(&argc, &argv);
-	int i, j, m, n, err;
-	double **a;
+	int m, n;
+	double **partial_field;
 	double start, end;
 	t_process_info pi;
 
-	if ((argc < 3) || (argc > 5)) {
-		fprintf(stderr, "Usage: heat <m> <n> [<epsilon> [delta_t]]!\n\n");
-		fprintf(stderr, "   <m>       -- X-Size of matrix\n");
-		fprintf(stderr, "   <n>       -- Y-Size of matrix\n");
-		fprintf(stderr, "   <epsilon> -- accuracy parameter\n");
-		fprintf(stderr, "   <delta_t> -- time step\n");
-		exit(1);
-	}
+	process_args(argc, argv, &m, &n, &eps, &delta_t);
 
-	/*
-	** Erstes und zweites Argument: Größe der Matrix
-	*/
-	m = atoi(argv[1]);
-	if ((m < 3) || (m > 6000)) {
-		fprintf(stderr, "Error: size m out of range [3 .. 6000] !\n");
-		exit(1);
-	}
-	n = atoi(argv[2]);
-	if ((n < 3) || (n > 6000)) {
-		fprintf(stderr, "Error: size n out of range [3 .. 6000] !\n");
-		exit(1);
-	}
-	/*
-	** Drittes (optionales) Argument: "Genauigkeitsfaktor" eps
-	*/
-	if (argc >= 4) {
-		eps = atof(argv[3]);
-	}
-	if ((eps <= 0) || (eps > 0.01)){
-		fprintf(stderr, "Error: epsilon must be between 0.0000000001 and 0.01 )\n");
-		exit(1);
-	}
-	/*
-	** Viertes (optionales) Argument: delta_t
-	*/
-	if (argc >= 5) {
-		delta_t = atof(argv[4]);
-	}
-	if ((delta_t <= 0) || (delta_t > 1.0)){
-		fprintf(stderr, "Error: delta_t (%f) must be between 0.0000000001 and 1.0\n", delta_t);
+/*
+	a = New_Matrix(m, n);
+	if (a == NULL) {
+		fprintf(stderr, "Can't allocate matrix !\n");
 		exit(1);
 	}
 
+	Init_Matrix(a, m, n, 0);
+*/
 	int pid, num_p;
-	if((err = MPI_Comm_rank(MPI_COMM_WORLD, &pid))) {
-		fprintf(stderr, "Cannot fetch PID due to %d\n", err);
+	if(MPI_Comm_rank(MPI_COMM_WORLD, &pid)) {
+		fprintf(stderr, "Cannot fetch PID\n");
 		exit(1);
 	}
 
@@ -122,6 +85,12 @@ int main(int argc, char **argv)
 	pi.end_m = (coord[0]+1)*m_per_pro -1;
 	pi.start_n = coord[1] * n_per_pro;
 	pi.end_n = (coord[0]+1) * n_per_pro -1;
+	if(pi.end_m > m - 1) {
+		pi.end_m = m - 1;
+	}
+	if(pi.end_n > n - 1) {
+		pi.end_n = n - 1;
+	}
 
 	// print own info
 	printf("pid(%d,%d) from (%d, %d) to (%d,%d)\n",
@@ -133,10 +102,47 @@ int main(int argc, char **argv)
 		pi.end_m
 	);
 
+	double delta_a;
+	partial_field = New_Matrix(pi.end_m - pi.start_m + 2, pi.end_n - pi.start_n + 2);
+	double **partial_field_tmp = New_Matrix(pi.end_m - pi.start_m + 2, pi.end_n - pi.start_n + 2);
+	int i,j;
+	double **swap;
+	double hx = 1.0/(double)m;
+	double hy = 1.0/(double)n;
+	double hx_square = hx * hx;
+	double hy_square = hy * hy;
+
+	double max_delta_t = 0.25*((min(hx,hy))*(min(hx,hy)))/alpha;  /* minimaler Wert f�r Konvergenz */
+	if (delta_t > max_delta_t) { 
+		delta_t = max_delta_t;
+		printf ("Info: delta_t set to %.10lf.\n", delta_t);
+	}
+	for(i = 1; i < pi.end_m - pi.start_m + 1; i++) {		// catch edges
+		for(j = 1; j < pi.end_n - pi.start_n + 1; j++) {	// catch edges
+			delta_a = alpha * 
+				    ( (partial_field[i+1][j] + partial_field[i-1][j] - 2.0 * partial_field[i][j]) / (hy_square)
+					 +(partial_field[i][j-1] + partial_field[i][j+1] - 2.0 * partial_field[i][j]) / (hx_square) );
+			delta_a = delta_a * delta_t;
+			partial_field_tmp[i][j] = partial_field[i][j] + delta_a;
+
+			if(delta_a > maxdiff)
+				maxdiff = delta_a;
+		}
+	}
+	swap = partial_field_tmp;
+	partial_field_tmp = partial_field;
+	partial_field = swap;
+	
+
+	if(pid == 0) {
+		Write_Matrix(partial_field, pi.end_m - pi.start_m + 2, pi.end_n - pi.start_n + 2);
+	}
+
 	// write neighbar_com function
 
+
+	//MPI_Gather()
 	MPI_Finalize();
 
 	return 0;
-}
-
+}
\ No newline at end of file

project.h

+double ** New_Matrix(int m, int n);
+void Delete_Matrix(double **matrix);
+void Write_Matrix(double **a, int m, int n);
+void Init_Matrix(double **a, int m, int n, int init_value);
+void Insert_Matrix(double **a, int a_dim0, int a_dim1, int pos_dim0, int pos_dim1, double **b, int b_dim0, int b_dim1, int offset_top, int offset_right, int offset_bottom, int offset_left);
+void process_args(int argc, char **argv, int *m, int *n, double *eps, double *delta_t);
+
+typedef struct {
+	int start_m;	// first tile to process in X-axis
+	int start_n;	// first tile to process in Y-axis
+	int end_m;	// last tile to process in X-axis
+	int end_n;	// last tile to process in Y-axis
+} t_process_info;
\ No newline at end of file