chiark / gitweb /
wip lp, results, compiles
[matchsticks-search.git] / main.c
diff --git a/main.c b/main.c
index e246752ef153b4272fb333f13e511fd05aa38252..bdcb18bfbb2c2fa54582aa14b383539e92a28e49 100644 (file)
--- a/main.c
+++ b/main.c
@@ -2,9 +2,12 @@
 #include <stdio.h>
 #include <stdint.h>
 #include <stdlib.h>
+#include <string.h>
+#include <assert.h>
 #include <inttypes.h>
 
 #include <publib.h>
+#include <glpk.h>
 
 typedef uint32_t AdjWord;
 #define PRADJ "08"PRIx32
@@ -14,14 +17,20 @@ static AdjWord *adjmatrix;
 static AdjWord adjall;
 
 static double best;
+static glp_prob *best_prob;
+static AdjWord *best_adjmatrix;
+
+static AdjWord *xalloc_adjmatrix(void) {
+  return xmalloc(sizeof(*adjmatrix)*n);
+}
 
 static void prep(void) {
   adjall = ~((~(AdjWord)0) << m);
-  adjmatrix = xmalloc(sizeof(*adjmatrix)*n);
+  adjmatrix = xalloc_adjmatrix();
 }
 
 static AdjWord one_adj_bit(int bitnum) {
-  return (AdjWord)1 << j;
+  return (AdjWord)1 << bitnum;
 }
 
 static int count_set_adj_bits(AdjWord w) {
@@ -32,7 +41,8 @@ static int count_set_adj_bits(AdjWord w) {
 }
 
 static void optimise(void) {
-  int i, totalfrags;
+  int i, j, totalfrags;
+
   for (i=0, totalfrags=0; i<n; i++) {
     int frags = count_set_adj_bits(adjmatrix[i]);
     totalfrags += frags;
@@ -70,19 +80,17 @@ static void optimise(void) {
 
   glp_prob *prob = glp_create_prob();
 
-  int Y_totals_i = glp_add_rows(prob, i);
-  int Y_totals_j = glp_add_rows(prob, j);
+  int Y_totals_i = glp_add_rows(prob, n);
+  int Y_totals_j = glp_add_rows(prob, m);
   int X_minimum = glp_add_cols(prob, 1);
-  int rows = glp_get_num_rows(prob);
-  int cols = glp_get_num_rows(cols);
 
   int next_matrix_entry = 1; /* wtf GLPK! */
-  int matrix_entries_size = next_matrix_entry + i + j + totalfrags*2;
+  int matrix_entries_size = next_matrix_entry + n + m + totalfrags*2;
   double matrix_entries[matrix_entries_size];
   int matrix_entries_XY[2][matrix_entries_size];
 
 #define ADD_MATRIX_ENTRY(Y,X) ({                       \
-      assert(matrix_entries_size < next_matrix_entry); \
+      assert(next_matrix_entry < matrix_entries_size); \
       matrix_entries_XY[0][next_matrix_entry] = X;     \
       matrix_entries_XY[1][next_matrix_entry] = Y;     \
       matrix_entries[next_matrix_entry] = 0;           \
@@ -97,11 +105,11 @@ static void optimise(void) {
 
   /* \forall_i x_totals_i = m */
   /* \forall_i x_totals_j = n */
-  for (i=0; i<n; i++) glp_set_row_bounds(prob, Y_totals_i+i, GLP_FX, m,m);
-  for (j=0; j<m; j++) glp_set_row_bounds(prob, Y_totals_j+j, GLP_FX, n,n);
+  for (i=0; i<n; i++) glp_set_row_bnds(prob, Y_totals_i+i, GLP_FX, m,m);
+  for (j=0; j<m; j++) glp_set_row_bnds(prob, Y_totals_j+j, GLP_FX, n,n);
 
   /* x_minimum >= 0 */
-  glp_set_col_bounds(prob, X_minimum, GLP_LB, 0, 0);
+  glp_set_col_bnds(prob, X_minimum, GLP_LO, 0, 0);
 
   /* objective is maximising x_minimum */
   glp_set_obj_dir(prob, GLP_MAX);
@@ -131,12 +139,54 @@ static void optimise(void) {
 
   assert(next_matrix_entry == matrix_entries_size);
 
-  for (row=1; row<=rows; row++) {
-    glp_load_matrix(prob, next_matrix_entry,
-                   matrix_entries_XY[1], matrix_entries_XY[0],
-                   matrix_entries);
+  glp_load_matrix(prob, next_matrix_entry,
+                 matrix_entries_XY[1], matrix_entries_XY[0],
+                 matrix_entries);
+
+  int r = glp_simplex(prob, NULL);
+  switch (r) {
+  case 0:;
+
+    double got = glp_get_obj_val(prob);
+    printf("%g", got);
+    if (got <= best) {
+      printf("\n");
+      break;
+    }
+
+    best = got;
+
+    if (best_prob) glp_delete_prob(best_prob);
+    best_prob = prob;
+
+    free(best_adjmatrix);
+    best_adjmatrix = xalloc_adjmatrix();
+    memcpy(best_adjmatrix, adjmatrix, sizeof(*adjmatrix)*n);
+
+    printf(" <--\n");
+    return;
+
+#define OKERR(e) \
+  case e: printf(" " #e "\n"); break;
+  OKERR(GLP_ESING);
+  OKERR(GLP_ECOND);
+  OKERR(GLP_EBOUND);
+  OKERR(GLP_EFAIL);
+  OKERR(GLP_ENOPFS);
+  OKERR(GLP_ENODFS);
+
+#define BADERR(e) \
+  case e: printf(" " #e " CRASHING\n"); exit(-1);
+  BADERR(GLP_EBADB);
+  BADERR(GLP_EOBJLL);
+  BADERR(GLP_EOBJUL);
+  BADERR(GLP_EITLIM);
+  BADERR(GLP_ETMLIM);
+
+  default: printf(" r=%d! CRASHING\n", r); exit(-1);
+  }
 
-  printf("nyi\n");
+  glp_delete_prob(prob);
 }
 
 static void iterate_recurse(int i, AdjWord min) {