X-Git-Url: http://www.chiark.greenend.org.uk/ucgi/~ianmdlvl/git?p=matchsticks-search.git;a=blobdiff_plain;f=main.c;h=93dc7ebe5325269db8fda05f5fe954a081489ed3;hp=fd0569be80e712602988a5f79207b77bcb45ff5b;hb=a7a2c255e61d14619fce2e4af2eeee1d2ed6d870;hpb=7c6fcb4307273d71a0f7a5f9cb0e8b7e9d3e26b2 diff --git a/main.c b/main.c index fd0569b..93dc7eb 100644 --- a/main.c +++ b/main.c @@ -34,6 +34,7 @@ #include #include #include +#include #include #include #include @@ -41,6 +42,10 @@ #include +#ifndef VERSION +#define VERSION "(unknown-version)" +#endif + /* * Algorithm. * @@ -53,7 +58,8 @@ * * We search all possible adjacency matrices, and for each one we run * GLPK's simplex solver. We represent the adjacency matrix as an - * array of bitmaps. + * array of bitmaps: one word per input stick, with one bit per output + * stick. * * However, there are a couple of wrinkles: * @@ -71,7 +77,7 @@ * nondecreasing in array order. * * Once we have a solution, we also avoid considering any candidate - * which involves dividing one of the output sticks into so many + * which involves dividing one of the input sticks into so many * fragment that the smallest fragment would necessarily be no bigger * than our best solution. That is, we reject candidates where any of * the hamming weights of the adjacency bitmap words are too large. @@ -86,6 +92,8 @@ typedef uint32_t AdjWord; #define PRADJ "08"PRIx32 +#define FOR_BITS(j,m) for (j=0, j##bit=1; j < (m); j++, j##bit<<=1) + static int n, m, maxhamweight; static AdjWord *adjmatrix; static AdjWord adjall; @@ -94,6 +102,9 @@ static double best; static glp_prob *best_prob; static AdjWord *best_adjmatrix; +static int n_max_frags, m_max_frags; +static int *weight; + static unsigned printcounter; static void iterate(void); @@ -109,6 +120,18 @@ static void progress_eol(void) { static void set_best(double new_best) { best = new_best; + /* + * When computing n_max_frags, we want to set a value that will skip + * anything that won't provide strictly better solutions. So we + * want + * frags < n / best + * _ _ + * <=> frags < | n / best | + * _ _ + * <=> frags <= | n / best | - 1 + */ + n_max_frags = ceil(n / best) - 1; + m_max_frags = ceil(m / best) - 1; } /*----- multicore support -----*/ @@ -150,10 +173,11 @@ typedef struct { pid_t pid; } Worker; static Worker *mc_us; +static bool mc_am_generator; static void multicore_check_for_new_best(void); -#define MAX_NIOVS 3 +#define MAX_NIOVS 4 static AdjWord mc_iter_min; static int mc_niovs; static size_t mc_iovlen; @@ -174,6 +198,7 @@ static void mc_rwvsetup_outer(void) { IOV(maxhamweight, 1); IOV(mc_iter_min, 1); IOV(*adjmatrix, multicore_iteration_boundary); + IOV(*weight, m); } static void mc_rwvsetup_full(void) { @@ -240,8 +265,8 @@ static void mc_iterate_worker(void) { iterate_recurse(mc_org_it_bound, mc_iter_min); multicore_iteration_boundary = mc_org_it_bound; } - LPRINTF("worker %2d reporting",mc_us->w); if (best_adjmatrix) { + LPRINTF("worker %2d reporting",mc_us->w); adjmatrix = best_adjmatrix; mc_rwvsetup_full(); ssize_t r = writev(fileno(mc_us->results), mc_iov, mc_niovs); @@ -284,6 +309,7 @@ static void multicore(void) { genpid = fork(); assert(genpid >= 0); if (!genpid) { + mc_am_generator = 1; LPRINTF("generator running"); iterate(); exit(0); @@ -299,13 +325,15 @@ static void multicore(void) { LPRINTF("reading report from %2d",w); ssize_t sr = preadv(fileno(mc_workers[w].results), mc_iov, mc_niovs, 0); if (!sr) continue; + LPRINTF("got report from %2d",w); maxhamweight = 0; optimise(1); } } static void multicore_check_for_new_best(void) { - if (!ncpus) return; + if (!(mc_us || mc_am_generator)) + return; for (;;) { double msg; @@ -319,7 +347,8 @@ static void multicore_check_for_new_best(void) { } static void multicore_found_new_best(void) { - if (!ncpus) return; + if (!mc_us) + return; if (mc_us /* might be master */) fprintf(stderr," w%-2d ",mc_us->w); ssize_t wrote = write(mc_bus, &best, sizeof(best)); @@ -337,16 +366,22 @@ static void prep(void) { adjmatrix = xalloc_adjmatrix(); glp_term_out(GLP_OFF); setlinebuf(stderr); + weight = calloc(sizeof(*weight), m); assert(weight); + n_max_frags = INT_MAX; + m_max_frags = INT_MAX; } +#if 0 static AdjWord one_adj_bit(int bitnum) { return (AdjWord)1 << bitnum; } +#endif static int count_set_adj_bits(AdjWord w) { - int j, total; - for (j=0, total=0; j best; + return maxhamweight <= m_max_frags; } static bool preconsider_ok(int nwords, bool doprint) { @@ -367,14 +401,13 @@ static bool preconsider_ok(int nwords, bool doprint) { bool had_max = 0; for (i=0, totalfrags=0; i= maxhamweight); - totalfrags += frags; PRINTF("%"PRADJ" ", adjmatrix[i]); - double maxminsize = (double)m / frags; - if (maxminsize <= best) { + if (frags > m_max_frags) { PRINTF(" too fine"); goto out; } + had_max += (frags >= maxhamweight); + totalfrags += frags; } if (!had_max) { /* Skip this candidate as its max hamming weight is lower than @@ -394,6 +427,7 @@ static void optimise(bool doprint) { /* Consider the best answer (if any) for a given adjacency matrix */ glp_prob *prob = 0; int i, j; + AdjWord jbit; /* * Up to a certain point, optimise() can be restarted. We use this @@ -479,8 +513,8 @@ static void optimise(bool doprint) { glp_set_obj_coef(prob, X_minimum, 1); for (i=0; i= n) { printcounter++; optimise(!(printcounter & 0xfff)); @@ -598,8 +635,20 @@ static void iterate_recurse(int i, AdjWord min) { if (i == 0 && (adjmatrix[i] & (1+adjmatrix[i]))) goto again; + FOR_BITS(j,m) + if (adjmatrix[i] & jbit) + weight[j]++; + for (int j = 0; j < m; j++) + if (weight[j] >= n_max_frags) + goto takeout; + iterate_recurse(i+1, adjmatrix[i]); + takeout: + FOR_BITS(j,m) + if (adjmatrix[i] & jbit) + weight[j]--; + again: if (adjmatrix[i] == adjall) return; @@ -631,7 +680,7 @@ static void report(void) { continue; a[x][y] = min + glp_get_col_prim(best_prob, i); } - printf("%d into %d: min fragment %g\n", n, m, min); + printf("%d into %d: min fragment %g [%s]\n", n, m, min, VERSION); for (i = 0; i < n; i++) { for (j = 0; j < m; j++) { if (a[i][j])