Working version

I finally implemented the algorithm. I have no idea if it's implemented
correctly, but that's a problem for another time. Oh, yeah. I should
free the benchmark data after use.

Signed-off-by: wmb <wmb@teknik.io>

benchmark-data.h

@@ -9,4 +9,6 @@ extern vector *benchmark_data;
 extern const size_t N_VECTORS;
 extern const size_t VECTOR_SIZE;
 
+void init_benchmark_data(void);
+
 #endif /* BENCHMARK_DATA_H */

main.c

@@ -10,8 +10,6 @@
 
 #include "benchmark-data.h"
 
-/* TODO: Use typedefs for solution and vector instead of using int (*s)[230]
- * and double (*v)[42]. */
 
 /* Types */
 typedef int *solution;
@@ -24,23 +22,25 @@ bool solution_valide(const solution s);
 void print_solution(const solution s, const char *name);
 vector *centres_gravite(const solution s);
 int solution_voisine(solution *dest, const solution src);
-#if 0
 double fonction_objective(const solution s);
 
-static double calculer_inter(const solution s, vector *cgs);
-#endif /* 0 */
+static double calculer_inter(vector *cgs);
 static double calculer_intra(const solution s, vector *cgs);
 
+
+/* Variables */
 solution solution_opt;
 solution new_solution;
 
 int n_clusters = 5;
 int n_iterations = 1000;
 
+
 int main(int argc, char *argv[])
 {
 	int c;
 	int i;
+	double f0;
 
 	while ((c = getopt(argc, argv, ":c:i:h")) != -1) {
 		int k;
@@ -66,6 +66,9 @@ int main(int argc, char *argv[])
 
 	printf("Using %d clusters\n", n_clusters);
 
+	/* Init benchmark_data */
+	init_benchmark_data();
+
 	srand(time(NULL));
 
 	if (!(solution_opt = calloc(N_VECTORS, sizeof(solution_opt[0]))))
@@ -75,16 +78,28 @@ int main(int argc, char *argv[])
 	while (!trouver_solution_initiale())
 		continue; /* do it again */
 
+	f0 = fonction_objective(solution_opt);
+
 	print_solution(solution_opt, "initiale");
+	printf("F(solution) = %f\n", f0);
 
 	for (i = 0; i < n_iterations; i++) {
+		double f1;
 
 		while (!solution_voisine(&new_solution, solution_opt))
 			continue; /* do it again */
 
 		/* TODO: Calculer les centres de gravités, les distances moyennes inter-/intra-cluster,
 		 * et faire le reste de l'algorithme */
-		print_solution(new_solution, "new solution");
+
+		if ((f1 = fonction_objective(new_solution)) > f0) {
+			print_solution(new_solution, "optimale trouvée");
+			printf("F(solution) = %f\n", f1);
+			free(solution_opt);
+			solution_opt = new_solution;
+			new_solution = NULL;
+			f0 = f1;
+		}
 	}
 
 	exit(EXIT_SUCCESS);
@@ -95,6 +110,11 @@ double distance(const vector *u, const vector *v)
 	size_t i = 0;
 	double sum = 0.0;
 
+	if (!u || !v) {
+		fputs("distance(): u or v is null.\n", stderr);
+		abort();
+	}
+
 	/* d(u, v) = sqrt( (u1 - v1)^2 + (u2 - v2)^2 + ... + (un - vn)^2 ) */
 
 	for (; i < VECTOR_SIZE; i++)
@@ -162,19 +182,23 @@ vector *centres_gravite(const solution s)
 
 	int cluster;
 	size_t i = 0;
-	vector *cgs = calloc(n_clusters, sizeof(*cgs));
-	int *vectors_in_cluster = calloc(n_clusters, sizeof(*vectors_in_cluster));
+
+	/* cgs[n_clusters] est le centre de gravité global */
+	vector *cgs = calloc(n_clusters + 1, sizeof(*cgs));
+
+	int *vectors_in_cluster = calloc(n_clusters + 1, sizeof(*vectors_in_cluster));
 
 	if (!cgs || !vectors_in_cluster)
 		abort();
 
-	for (cluster = 0; cluster < n_clusters; cluster++) {
+	for (cluster = 0; cluster < n_clusters + 1; cluster++) {
 		if (!(cgs[cluster] = calloc(VECTOR_SIZE, sizeof(*cgs[cluster]))))
 			abort();
 
 		/* Initialiser à 0 */
 		while (i < VECTOR_SIZE)
 			cgs[cluster][i++] = 0;
+		if (cluster < n_clusters)
 		vectors_in_cluster[cluster] = 0;
 	}
 
@@ -184,11 +208,15 @@ vector *centres_gravite(const solution s)
 		cluster = s[i];
 		++vectors_in_cluster[cluster];
 
-		for (j = 0; j < VECTOR_SIZE; j++)
+		for (j = 0; j < VECTOR_SIZE; j++) {
 			cgs[cluster][j] += benchmark_data[i][j];
+			cgs[n_clusters][j] += benchmark_data[i][j];
+		}
 	}
 
-	for (cluster = 0; cluster < n_clusters; cluster++)
+	vectors_in_cluster[n_clusters] = N_VECTORS;
+
+	for (cluster = 0; cluster < n_clusters + 1; cluster++)
 		for (i = 0; i < VECTOR_SIZE; i++)
 			cgs[cluster][i] /= (double) vectors_in_cluster[cluster];
 
@@ -230,14 +258,10 @@ int solution_voisine(solution *dest, const solution src)
 	return valid;
 }
 
-#if 0
 double fonction_objective(const solution s)
 {
-	/* NOTE: this is very ineffecient on processor time, because
-	 * each of these functions re-calculates the same centre_gravite
-	 * twice. */
-
-	double **cgs = calloc(n_clusters, sizeof *cgs);
+	int i;
+	double **cgs = centres_gravite(s);
 	double inter;
 	double intra;
 
@@ -245,26 +269,27 @@ double fonction_objective(const solution s)
 	if (!cgs)
 		abort();
 
-	inter = calculer_inter(s, cgs);
-	inter = calculer_inter(s, cgs);
+	inter = calculer_inter(cgs);
+	intra = calculer_intra(s, cgs);
+
+	for (i = 0; i < n_clusters + 1; i++)
+		free(cgs[i]);
+
+	free(cgs);
 
 	return inter - intra;
 }
 
-static double calculer_inter(const solution s, vector *cgs)
+static double calculer_inter(vector *cgs)
 {
-	const size_t n = sizeof(*s) / sizeof((*s)[0]);
 	double inter = 0.0;
 	int cluster;
 
-	for (cluster = 0; cluster < n_clusters; cluster++) {
-		double intra_for_this_cluster = 0.0;
-		int vectors_in_this_cluster = 0;
-		size_t i;
-
-		centre_gravite(cluster, s, &cg);
-	}
+	for (cluster = 0; cluster < n_clusters; cluster++)
+		inter += distance((const vector *) cgs[n_clusters],
+						  (const vector *) cgs[cluster]);
 
+	return inter / (double) n_clusters;
 }
 
 static double calculer_intra(const solution s, vector *cgs)
@@ -283,7 +308,7 @@ static double calculer_intra(const solution s, vector *cgs)
 
 			++vectors_in_this_cluster;
 			intra_for_this_cluster += distance((const vector *) cgs[cluster],
-					(const vector *) benchmark_data[i]);
+											   (const vector *) benchmark_data[i]);
 		}
 
 		intra += intra_for_this_cluster / (double) vectors_in_this_cluster;
@@ -291,4 +316,3 @@ static double calculer_intra(const solution s, vector *cgs)
 
 	return intra / (double) n_clusters;
 }
-#endif /* 0 */