Chris@10: /*
Chris@10:  * Copyright (c) 2000 Matteo Frigo
Chris@10:  * Copyright (c) 2000 Massachusetts Institute of Technology
Chris@10:  *
Chris@10:  * This program is free software; you can redistribute it and/or modify
Chris@10:  * it under the terms of the GNU General Public License as published by
Chris@10:  * the Free Software Foundation; either version 2 of the License, or
Chris@10:  * (at your option) any later version.
Chris@10:  *
Chris@10:  * This program is distributed in the hope that it will be useful,
Chris@10:  * but WITHOUT ANY WARRANTY; without even the implied warranty of
Chris@10:  * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
Chris@10:  * GNU General Public License for more details.
Chris@10:  *
Chris@10:  * You should have received a copy of the GNU General Public License
Chris@10:  * along with this program; if not, write to the Free Software
Chris@10:  * Foundation, Inc., 51 Franklin Street, Fifth Floor, Boston, MA  02110-1301  USA
Chris@10:  *
Chris@10:  */
Chris@10: 
Chris@10: #include "ifftw.h"
Chris@10: #include <string.h>
Chris@10: 
Chris@10: /* GNU Coding Standards, Sec. 5.2: "Please write the comments in a GNU
Chris@10:    program in English, because English is the one language that nearly
Chris@10:    all programmers in all countries can read."
Chris@10: 
Chris@10:                     ingemisco tanquam reus
Chris@10: 		    culpa rubet vultus meus
Chris@10: 		    supplicanti parce [rms]
Chris@10: */
Chris@10: 
Chris@10: #define VALIDP(solution) ((solution)->flags.hash_info & H_VALID)
Chris@10: #define LIVEP(solution) ((solution)->flags.hash_info & H_LIVE)
Chris@10: #define SLVNDX(solution) ((solution)->flags.slvndx)
Chris@10: #define BLISS(flags) (((flags).hash_info) & BLESSING)
Chris@10: #define INFEASIBLE_SLVNDX ((1U<<BITS_FOR_SLVNDX)-1)
Chris@10: 
Chris@10: 
Chris@10: #define MAXNAM 64  /* maximum length of registrar's name.
Chris@10: 		      Used for reading wisdom.  There is no point
Chris@10: 		      in doing this right */
Chris@10: 
Chris@10: 
Chris@10: #ifdef FFTW_DEBUG
Chris@10: static void check(hashtab *ht);
Chris@10: #endif
Chris@10: 
Chris@10: /* x <= y */
Chris@10: #define LEQ(x, y) (((x) & (y)) == (x))
Chris@10: 
Chris@10: /* A subsumes B */
Chris@10: static int subsumes(const flags_t *a, unsigned slvndx_a, const flags_t *b)
Chris@10: {
Chris@10:      if (slvndx_a != INFEASIBLE_SLVNDX) {
Chris@10: 	  A(a->timelimit_impatience == 0);
Chris@10: 	  return (LEQ(a->u, b->u) && LEQ(b->l, a->l));
Chris@10:      } else {
Chris@10: 	  return (LEQ(a->l, b->l) 
Chris@10: 		  && a->timelimit_impatience <= b->timelimit_impatience);
Chris@10:      }
Chris@10: }
Chris@10: 
Chris@10: static unsigned addmod(unsigned a, unsigned b, unsigned p)
Chris@10: {
Chris@10:      /* gcc-2.95/sparc produces incorrect code for the fast version below. */
Chris@10: #if defined(__sparc__) && defined(__GNUC__)
Chris@10:      /* slow version  */
Chris@10:      return (a + b) % p;
Chris@10: #else
Chris@10:      /* faster version */
Chris@10:      unsigned c = a + b;
Chris@10:      return c >= p ? c - p : c;
Chris@10: #endif
Chris@10: }
Chris@10: 
Chris@10: /*
Chris@10:   slvdesc management:
Chris@10: */
Chris@10: static void sgrow(planner *ego)
Chris@10: {
Chris@10:      unsigned osiz = ego->slvdescsiz, nsiz = 1 + osiz + osiz / 4;
Chris@10:      slvdesc *ntab = (slvdesc *)MALLOC(nsiz * sizeof(slvdesc), SLVDESCS);
Chris@10:      slvdesc *otab = ego->slvdescs;
Chris@10:      unsigned i;
Chris@10: 
Chris@10:      ego->slvdescs = ntab;
Chris@10:      ego->slvdescsiz = nsiz;
Chris@10:      for (i = 0; i < osiz; ++i)
Chris@10: 	  ntab[i] = otab[i];
Chris@10:      X(ifree0)(otab);
Chris@10: }
Chris@10: 
Chris@10: static void register_solver(planner *ego, solver *s)
Chris@10: {
Chris@10:      slvdesc *n;
Chris@10:      int kind;
Chris@10: 
Chris@10:      if (s) { /* add s to solver list */
Chris@10: 	  X(solver_use)(s);
Chris@10: 
Chris@10: 	  A(ego->nslvdesc < INFEASIBLE_SLVNDX);
Chris@10: 	  if (ego->nslvdesc >= ego->slvdescsiz)
Chris@10: 	       sgrow(ego);
Chris@10: 
Chris@10: 	  n = ego->slvdescs + ego->nslvdesc;
Chris@10: 
Chris@10: 	  n->slv = s;
Chris@10: 	  n->reg_nam = ego->cur_reg_nam;
Chris@10: 	  n->reg_id = ego->cur_reg_id++;
Chris@10: 	  
Chris@10: 	  A(strlen(n->reg_nam) < MAXNAM);
Chris@10: 	  n->nam_hash = X(hash)(n->reg_nam);
Chris@10: 
Chris@10: 	  kind = s->adt->problem_kind;
Chris@10: 	  n->next_for_same_problem_kind = ego->slvdescs_for_problem_kind[kind];
Chris@10: 	  ego->slvdescs_for_problem_kind[kind] = ego->nslvdesc;
Chris@10: 
Chris@10: 	  ego->nslvdesc++;
Chris@10:      }
Chris@10: }
Chris@10: 
Chris@10: static unsigned slookup(planner *ego, char *nam, int id)
Chris@10: {
Chris@10:      unsigned h = X(hash)(nam); /* used to avoid strcmp in the common case */
Chris@10:      FORALL_SOLVERS(ego, s, sp, {
Chris@10: 	  UNUSED(s);
Chris@10: 	  if (sp->reg_id == id && sp->nam_hash == h
Chris@10: 	      && !strcmp(sp->reg_nam, nam))
Chris@10: 	       return sp - ego->slvdescs;
Chris@10:      });
Chris@10:      return INFEASIBLE_SLVNDX;
Chris@10: }
Chris@10: 
Chris@10: /* Compute a MD5 hash of the configuration of the planner.
Chris@10:    We store it into the wisdom file to make absolutely sure that
Chris@10:    we are reading wisdom that is applicable */
Chris@10: static void signature_of_configuration(md5 *m, planner *ego)
Chris@10: {
Chris@10:      X(md5begin)(m);
Chris@10:      X(md5unsigned)(m, sizeof(R)); /* so we don't mix different precisions */
Chris@10:      FORALL_SOLVERS(ego, s, sp, {
Chris@10: 	  UNUSED(s);
Chris@10: 	  X(md5int)(m, sp->reg_id);
Chris@10: 	  X(md5puts)(m, sp->reg_nam);
Chris@10:      });
Chris@10:      X(md5end)(m);
Chris@10: }
Chris@10: 
Chris@10: /*
Chris@10:   md5-related stuff:
Chris@10: */
Chris@10: 
Chris@10: /* first hash function */
Chris@10: static unsigned h1(const hashtab *ht, const md5sig s)
Chris@10: {
Chris@10:      unsigned h = s[0] % ht->hashsiz;
Chris@10:      A(h == (s[0] % ht->hashsiz));
Chris@10:      return h;
Chris@10: }
Chris@10: 
Chris@10: /* second hash function (for double hashing) */
Chris@10: static unsigned h2(const hashtab *ht, const md5sig s)
Chris@10: {
Chris@10:      unsigned h = 1U + s[1] % (ht->hashsiz - 1);
Chris@10:      A(h == (1U + s[1] % (ht->hashsiz - 1)));
Chris@10:      return h;
Chris@10: }
Chris@10: 
Chris@10: static void md5hash(md5 *m, const problem *p, const planner *plnr)
Chris@10: {
Chris@10:      X(md5begin)(m);
Chris@10:      X(md5unsigned)(m, sizeof(R)); /* so we don't mix different precisions */
Chris@10:      X(md5int)(m, plnr->nthr);
Chris@10:      p->adt->hash(p, m);
Chris@10:      X(md5end)(m);
Chris@10: }
Chris@10: 
Chris@10: static int md5eq(const md5sig a, const md5sig b)
Chris@10: {
Chris@10:      return a[0] == b[0] && a[1] == b[1] && a[2] == b[2] && a[3] == b[3];
Chris@10: }
Chris@10: 
Chris@10: static void sigcpy(const md5sig a, md5sig b)
Chris@10: {
Chris@10:      b[0] = a[0]; b[1] = a[1]; b[2] = a[2]; b[3] = a[3];
Chris@10: }
Chris@10: 
Chris@10: /*
Chris@10:   memoization routines :
Chris@10: */
Chris@10: 
Chris@10: /*
Chris@10:    liber scriptus proferetur
Chris@10:    in quo totum continetur
Chris@10:    unde mundus iudicetur
Chris@10: */
Chris@10: struct solution_s {
Chris@10:      md5sig s;
Chris@10:      flags_t flags;
Chris@10: };
Chris@10: 
Chris@10: static solution *htab_lookup(hashtab *ht, const md5sig s, 
Chris@10: 			     const flags_t *flagsp)
Chris@10: {
Chris@10:      unsigned g, h = h1(ht, s), d = h2(ht, s);
Chris@10:      solution *best = 0;
Chris@10: 
Chris@10:      ++ht->lookup;
Chris@10: 
Chris@10:      /* search all entries that match; select the one with
Chris@10: 	the lowest flags.u */
Chris@10:      /* This loop may potentially traverse the whole table, since at
Chris@10: 	least one element is guaranteed to be !LIVEP, but all elements
Chris@10: 	may be VALIDP.  Hence, we stop after at the first invalid
Chris@10: 	element or after traversing the whole table. */
Chris@10:      g = h;
Chris@10:      do {
Chris@10: 	  solution *l = ht->solutions + g;
Chris@10: 	  ++ht->lookup_iter;
Chris@10: 	  if (VALIDP(l)) {
Chris@10: 	       if (LIVEP(l)
Chris@10: 		   && md5eq(s, l->s)
Chris@10: 		   && subsumes(&l->flags, SLVNDX(l), flagsp) ) { 
Chris@10: 		    if (!best || LEQ(l->flags.u, best->flags.u))
Chris@10: 			 best = l;
Chris@10: 	       }
Chris@10: 	  } else 
Chris@10: 	       break;
Chris@10: 
Chris@10: 	  g = addmod(g, d, ht->hashsiz);
Chris@10:      } while (g != h);
Chris@10: 
Chris@10:      if (best) 
Chris@10: 	  ++ht->succ_lookup;
Chris@10:      return best;
Chris@10: }
Chris@10: 
Chris@10: static solution *hlookup(planner *ego, const md5sig s, 
Chris@10: 			 const flags_t *flagsp)
Chris@10: {
Chris@10:      solution *sol = htab_lookup(&ego->htab_blessed, s, flagsp);
Chris@10:      if (!sol) sol = htab_lookup(&ego->htab_unblessed, s, flagsp);
Chris@10:      return sol;
Chris@10: }
Chris@10: 
Chris@10: static void fill_slot(hashtab *ht, const md5sig s, const flags_t *flagsp,
Chris@10: 		      unsigned slvndx, solution *slot)
Chris@10: {
Chris@10:      ++ht->insert;
Chris@10:      ++ht->nelem;
Chris@10:      A(!LIVEP(slot));
Chris@10:      slot->flags.u = flagsp->u;
Chris@10:      slot->flags.l = flagsp->l;
Chris@10:      slot->flags.timelimit_impatience = flagsp->timelimit_impatience;
Chris@10:      slot->flags.hash_info |= H_VALID | H_LIVE;
Chris@10:      SLVNDX(slot) = slvndx;
Chris@10: 
Chris@10:      /* keep this check enabled in case we add so many solvers
Chris@10: 	that the bitfield overflows */
Chris@10:      CK(SLVNDX(slot) == slvndx);     
Chris@10:      sigcpy(s, slot->s);
Chris@10: }
Chris@10: 
Chris@10: static void kill_slot(hashtab *ht, solution *slot)
Chris@10: {
Chris@10:      A(LIVEP(slot)); /* ==> */ A(VALIDP(slot));
Chris@10: 
Chris@10:      --ht->nelem;
Chris@10:      slot->flags.hash_info = H_VALID;
Chris@10: }
Chris@10: 
Chris@10: static void hinsert0(hashtab *ht, const md5sig s, const flags_t *flagsp, 
Chris@10: 		     unsigned slvndx)
Chris@10: {
Chris@10:      solution *l;
Chris@10:      unsigned g, h = h1(ht, s), d = h2(ht, s); 
Chris@10: 
Chris@10:      ++ht->insert_unknown;
Chris@10: 
Chris@10:      /* search for nonfull slot */
Chris@10:      for (g = h; ; g = addmod(g, d, ht->hashsiz)) {
Chris@10: 	  ++ht->insert_iter;
Chris@10: 	  l = ht->solutions + g;
Chris@10: 	  if (!LIVEP(l)) break;
Chris@10: 	  A((g + d) % ht->hashsiz != h);
Chris@10:      }
Chris@10: 
Chris@10:      fill_slot(ht, s, flagsp, slvndx, l);
Chris@10: }
Chris@10: 
Chris@10: static void rehash(hashtab *ht, unsigned nsiz)
Chris@10: {
Chris@10:      unsigned osiz = ht->hashsiz, h;
Chris@10:      solution *osol = ht->solutions, *nsol;
Chris@10: 
Chris@10:      nsiz = (unsigned)X(next_prime)((INT)nsiz);
Chris@10:      nsol = (solution *)MALLOC(nsiz * sizeof(solution), HASHT);
Chris@10:      ++ht->nrehash;
Chris@10: 
Chris@10:      /* init new table */
Chris@10:      for (h = 0; h < nsiz; ++h) 
Chris@10: 	  nsol[h].flags.hash_info = 0;
Chris@10: 
Chris@10:      /* install new table */
Chris@10:      ht->hashsiz = nsiz;
Chris@10:      ht->solutions = nsol;
Chris@10:      ht->nelem = 0;
Chris@10: 
Chris@10:      /* copy table */
Chris@10:      for (h = 0; h < osiz; ++h) {
Chris@10: 	  solution *l = osol + h;
Chris@10: 	  if (LIVEP(l))
Chris@10: 	       hinsert0(ht, l->s, &l->flags, SLVNDX(l));
Chris@10:      }
Chris@10: 
Chris@10:      X(ifree0)(osol);
Chris@10: }
Chris@10: 
Chris@10: static unsigned minsz(unsigned nelem)
Chris@10: {
Chris@10:      return 1U + nelem + nelem / 8U;
Chris@10: }
Chris@10: 
Chris@10: static unsigned nextsz(unsigned nelem)
Chris@10: {
Chris@10:      return minsz(minsz(nelem));
Chris@10: }
Chris@10: 
Chris@10: static void hgrow(hashtab *ht)
Chris@10: {
Chris@10:      unsigned nelem = ht->nelem;
Chris@10:      if (minsz(nelem) >= ht->hashsiz)
Chris@10: 	  rehash(ht, nextsz(nelem));
Chris@10: }
Chris@10: 
Chris@10: #if 0
Chris@10: /* shrink the hash table, never used */
Chris@10: static void hshrink(hashtab *ht)
Chris@10: {
Chris@10:      unsigned nelem = ht->nelem;
Chris@10:      /* always rehash after deletions */
Chris@10:      rehash(ht, nextsz(nelem));
Chris@10: }
Chris@10: #endif
Chris@10: 
Chris@10: static void htab_insert(hashtab *ht, const md5sig s, const flags_t *flagsp,
Chris@10: 			unsigned slvndx)
Chris@10: {
Chris@10:      unsigned g, h = h1(ht, s), d = h2(ht, s);
Chris@10:      solution *first = 0;
Chris@10: 
Chris@10:      /* Remove all entries that are subsumed by the new one.  */
Chris@10:      /* This loop may potentially traverse the whole table, since at
Chris@10: 	least one element is guaranteed to be !LIVEP, but all elements
Chris@10: 	may be VALIDP.  Hence, we stop after at the first invalid
Chris@10: 	element or after traversing the whole table. */
Chris@10:      g = h;
Chris@10:      do {
Chris@10: 	  solution *l = ht->solutions + g;
Chris@10: 	  ++ht->insert_iter;
Chris@10: 	  if (VALIDP(l)) {
Chris@10: 	       if (LIVEP(l) && md5eq(s, l->s)) {
Chris@10: 		    if (subsumes(flagsp, slvndx, &l->flags)) {
Chris@10: 			 if (!first) first = l;
Chris@10: 			 kill_slot(ht, l);
Chris@10: 		    } else {
Chris@10: 			 /* It is an error to insert an element that
Chris@10: 			    is subsumed by an existing entry. */
Chris@10: 			 A(!subsumes(&l->flags, SLVNDX(l), flagsp));
Chris@10: 		    }
Chris@10: 	       }
Chris@10: 	  } else 
Chris@10: 	       break;
Chris@10: 
Chris@10: 	  g = addmod(g, d, ht->hashsiz);
Chris@10:      } while (g != h);
Chris@10: 
Chris@10:      if (first) {
Chris@10: 	  /* overwrite FIRST */
Chris@10: 	  fill_slot(ht, s, flagsp, slvndx, first);
Chris@10:      } else {
Chris@10: 	  /* create a new entry */
Chris@10:  	  hgrow(ht);
Chris@10: 	  hinsert0(ht, s, flagsp, slvndx);
Chris@10:      }
Chris@10: }
Chris@10: 
Chris@10: static void hinsert(planner *ego, const md5sig s, const flags_t *flagsp, 
Chris@10: 		    unsigned slvndx)
Chris@10: {
Chris@10:      htab_insert(BLISS(*flagsp) ? &ego->htab_blessed : &ego->htab_unblessed,
Chris@10: 		 s, flagsp, slvndx );
Chris@10: }
Chris@10: 
Chris@10: 
Chris@10: static void invoke_hook(planner *ego, plan *pln, const problem *p, 
Chris@10: 			int optimalp)
Chris@10: {
Chris@10:      if (ego->hook)
Chris@10: 	  ego->hook(ego, pln, p, optimalp);
Chris@10: }
Chris@10: 
Chris@10: #ifdef FFTW_RANDOM_ESTIMATOR
Chris@10: /* a "random" estimate, used for debugging to generate "random"
Chris@10:    plans, albeit from a deterministic seed. */
Chris@10: 
Chris@10: unsigned X(random_estimate_seed) = 0;
Chris@10: 
Chris@10: static double random_estimate(const planner *ego, const plan *pln, 
Chris@10: 			      const problem *p)
Chris@10: {
Chris@10:      md5 m;
Chris@10:      X(md5begin)(&m);
Chris@10:      X(md5unsigned)(&m, X(random_estimate_seed));
Chris@10:      X(md5int)(&m, ego->nthr);
Chris@10:      p->adt->hash(p, &m);
Chris@10:      X(md5putb)(&m, &pln->ops, sizeof(pln->ops));
Chris@10:      X(md5putb)(&m, &pln->adt, sizeof(pln->adt));
Chris@10:      X(md5end)(&m);
Chris@10:      return ego->cost_hook ? ego->cost_hook(p, m.s[0], COST_MAX) : m.s[0];
Chris@10: }
Chris@10: 
Chris@10: #endif
Chris@10: 
Chris@10: double X(iestimate_cost)(const planner *ego, const plan *pln, const problem *p)
Chris@10: {
Chris@10:      double cost =
Chris@10: 	  + pln->ops.add
Chris@10: 	  + pln->ops.mul
Chris@10: 	  
Chris@10: #if HAVE_FMA
Chris@10: 	  + pln->ops.fma
Chris@10: #else
Chris@10: 	  + 2 * pln->ops.fma
Chris@10: #endif
Chris@10: 	  
Chris@10: 	  + pln->ops.other;
Chris@10:      if (ego->cost_hook)
Chris@10: 	  cost = ego->cost_hook(p, cost, COST_MAX);
Chris@10:      return cost;
Chris@10: }
Chris@10: 
Chris@10: static void evaluate_plan(planner *ego, plan *pln, const problem *p)
Chris@10: {
Chris@10:      if (ESTIMATEP(ego) || !BELIEVE_PCOSTP(ego) || pln->pcost == 0.0) {
Chris@10: 	  ego->nplan++;
Chris@10: 
Chris@10: 	  if (ESTIMATEP(ego)) {
Chris@10: 	  estimate:
Chris@10: 	       /* heuristic */
Chris@10: #ifdef FFTW_RANDOM_ESTIMATOR
Chris@10: 	       pln->pcost = random_estimate(ego, pln, p);
Chris@10: 	       ego->epcost += X(iestimate_cost)(ego, pln, p);
Chris@10: #else
Chris@10: 	       pln->pcost = X(iestimate_cost)(ego, pln, p);
Chris@10: 	       ego->epcost += pln->pcost;
Chris@10: #endif
Chris@10: 	  } else {
Chris@10: 	       double t = X(measure_execution_time)(ego, pln, p);
Chris@10: 	       
Chris@10: 	       if (t < 0) {  /* unavailable cycle counter */
Chris@10: 		    /* Real programmers can write FORTRAN in any language */
Chris@10: 		    goto estimate;
Chris@10: 	       }
Chris@10: 
Chris@10: 	       pln->pcost = t;
Chris@10: 	       ego->pcost += t;
Chris@10: 	       ego->need_timeout_check = 1;
Chris@10: 	  }
Chris@10:      }
Chris@10:      
Chris@10:      invoke_hook(ego, pln, p, 0);
Chris@10: }
Chris@10: 
Chris@10: /* maintain dynamic scoping of flags, nthr: */
Chris@10: static plan *invoke_solver(planner *ego, const problem *p, solver *s, 
Chris@10: 			   const flags_t *nflags)
Chris@10: {
Chris@10:      flags_t flags = ego->flags;
Chris@10:      int nthr = ego->nthr;
Chris@10:      plan *pln;
Chris@10:      ego->flags = *nflags;
Chris@10:      PLNR_TIMELIMIT_IMPATIENCE(ego) = 0;
Chris@10:      A(p->adt->problem_kind == s->adt->problem_kind);
Chris@10:      pln = s->adt->mkplan(s, p, ego);
Chris@10:      ego->nthr = nthr;
Chris@10:      ego->flags = flags;
Chris@10:      return pln;
Chris@10: }
Chris@10: 
Chris@10: /* maintain the invariant TIMED_OUT ==> NEED_TIMEOUT_CHECK */
Chris@10: static int timeout_p(planner *ego, const problem *p)
Chris@10: {
Chris@10:      /* do not timeout when estimating.  First, the estimator is the
Chris@10: 	planner of last resort.  Second, calling X(elapsed_since)() is
Chris@10: 	slower than estimating */
Chris@10:      if (!ESTIMATEP(ego)) {
Chris@10: 	  /* do not assume that X(elapsed_since)() is monotonic */
Chris@10: 	  if (ego->timed_out) {
Chris@10: 	       A(ego->need_timeout_check);
Chris@10: 	       return 1;
Chris@10: 	  }
Chris@10: 
Chris@10: 	  if (ego->timelimit >= 0 &&
Chris@10: 	      X(elapsed_since)(ego, p, ego->start_time) >= ego->timelimit) {
Chris@10: 	       ego->timed_out = 1;
Chris@10: 	       ego->need_timeout_check = 1;
Chris@10: 	       return 1;
Chris@10: 	  }
Chris@10:      }
Chris@10: 
Chris@10:      A(!ego->timed_out);
Chris@10:      ego->need_timeout_check = 0;
Chris@10:      return 0;
Chris@10: }
Chris@10: 
Chris@10: static plan *search0(planner *ego, const problem *p, unsigned *slvndx, 
Chris@10: 		     const flags_t *flagsp)
Chris@10: {
Chris@10:      plan *best = 0;
Chris@10:      int best_not_yet_timed = 1;
Chris@10: 
Chris@10:      /* Do not start a search if the planner timed out. This check is
Chris@10: 	necessary, lest the relaxation mechanism kick in */
Chris@10:      if (timeout_p(ego, p))
Chris@10: 	  return 0;
Chris@10: 
Chris@10:      FORALL_SOLVERS_OF_KIND(p->adt->problem_kind, ego, s, sp, {
Chris@10: 	  plan *pln;
Chris@10: 
Chris@10: 	  pln = invoke_solver(ego, p, s, flagsp);
Chris@10: 
Chris@10: 	  if (ego->need_timeout_check) 
Chris@10: 	       if (timeout_p(ego, p)) {
Chris@10: 		    X(plan_destroy_internal)(pln);
Chris@10: 		    X(plan_destroy_internal)(best);
Chris@10: 		    return 0;
Chris@10: 	       }
Chris@10: 
Chris@10: 	  if (pln) {
Chris@10: 	       /* read COULD_PRUNE_NOW_P because PLN may be destroyed
Chris@10: 		  before we use COULD_PRUNE_NOW_P */
Chris@10: 	       int could_prune_now_p = pln->could_prune_now_p;
Chris@10: 
Chris@10: 	       if (best) {
Chris@10: 		    if (best_not_yet_timed) {
Chris@10: 			 evaluate_plan(ego, best, p);
Chris@10: 			 best_not_yet_timed = 0;
Chris@10: 		    }
Chris@10: 		    evaluate_plan(ego, pln, p);
Chris@10: 		    if (pln->pcost < best->pcost) {
Chris@10: 			 X(plan_destroy_internal)(best);
Chris@10: 			 best = pln;
Chris@10: 			 *slvndx = sp - ego->slvdescs;
Chris@10: 		    } else {
Chris@10: 			 X(plan_destroy_internal)(pln);
Chris@10: 		    }
Chris@10: 	       } else {
Chris@10: 		    best = pln;
Chris@10: 		    *slvndx = sp - ego->slvdescs;
Chris@10: 	       }
Chris@10: 
Chris@10: 	       if (ALLOW_PRUNINGP(ego) && could_prune_now_p) 
Chris@10: 		    break;
Chris@10: 	  }
Chris@10:      });
Chris@10: 
Chris@10:      return best;
Chris@10: }
Chris@10: 
Chris@10: static plan *search(planner *ego, const problem *p, unsigned *slvndx, 
Chris@10: 		    flags_t *flagsp)
Chris@10: {
Chris@10:      plan *pln = 0;
Chris@10:      unsigned i;
Chris@10: 
Chris@10:      /* relax impatience in this order: */
Chris@10:      static const unsigned relax_tab[] = {
Chris@10: 	  0, /* relax nothing */
Chris@10: 	  NO_VRECURSE,
Chris@10: 	  NO_FIXED_RADIX_LARGE_N,
Chris@10: 	  NO_SLOW,
Chris@10: 	  NO_UGLY
Chris@10:      };
Chris@10: 
Chris@10:      unsigned l_orig = flagsp->l;
Chris@10:      unsigned x = flagsp->u;
Chris@10: 
Chris@10:      /* guaranteed to be different from X */
Chris@10:      unsigned last_x = ~x; 
Chris@10: 
Chris@10:      for (i = 0; i < sizeof(relax_tab) / sizeof(relax_tab[0]); ++i) {
Chris@10: 	  if (LEQ(l_orig, x & ~relax_tab[i]))
Chris@10: 	       x = x & ~relax_tab[i];
Chris@10: 
Chris@10: 	  if (x != last_x) {
Chris@10: 	       last_x = x;
Chris@10: 	       flagsp->l = x;
Chris@10: 	       pln = search0(ego, p, slvndx, flagsp);
Chris@10: 	       if (pln) break;
Chris@10: 	  }
Chris@10:      }
Chris@10: 
Chris@10:      if (!pln) {
Chris@10: 	  /* search [L_ORIG, U] */
Chris@10: 	  if (l_orig != last_x) {
Chris@10: 	       last_x = l_orig;
Chris@10: 	       flagsp->l = l_orig;
Chris@10: 	       pln = search0(ego, p, slvndx, flagsp);
Chris@10: 	  }
Chris@10:      }
Chris@10: 
Chris@10:      return pln;
Chris@10: }
Chris@10: 
Chris@10: #define CHECK_FOR_BOGOSITY						\
Chris@10:      if ((ego->bogosity_hook ?						\
Chris@10: 	  (ego->wisdom_state = ego->bogosity_hook(ego->wisdom_state, p)) \
Chris@10: 	  : ego->wisdom_state) == WISDOM_IS_BOGUS)			\
Chris@10: 	  goto wisdom_is_bogus;
Chris@10: 
Chris@10: static plan *mkplan(planner *ego, const problem *p)
Chris@10: {
Chris@10:      plan *pln;
Chris@10:      md5 m;
Chris@10:      unsigned slvndx;
Chris@10:      flags_t flags_of_solution;
Chris@10:      solution *sol;
Chris@10:      solver *s;
Chris@10: 
Chris@10:      ASSERT_ALIGNED_DOUBLE;
Chris@10:      A(LEQ(PLNR_L(ego), PLNR_U(ego)));
Chris@10: 
Chris@10:      if (ESTIMATEP(ego))
Chris@10: 	  PLNR_TIMELIMIT_IMPATIENCE(ego) = 0; /* canonical form */
Chris@10: 
Chris@10: 
Chris@10: #ifdef FFTW_DEBUG
Chris@10:      check(&ego->htab_blessed);
Chris@10:      check(&ego->htab_unblessed);
Chris@10: #endif
Chris@10: 
Chris@10:      pln = 0;
Chris@10: 
Chris@10:      CHECK_FOR_BOGOSITY;
Chris@10: 
Chris@10:      ego->timed_out = 0;
Chris@10: 
Chris@10:      ++ego->nprob;
Chris@10:      md5hash(&m, p, ego);
Chris@10: 
Chris@10:      flags_of_solution = ego->flags;
Chris@10: 
Chris@10:      if (ego->wisdom_state != WISDOM_IGNORE_ALL) {
Chris@10: 	  if ((sol = hlookup(ego, m.s, &flags_of_solution))) { 
Chris@10: 	       /* wisdom is acceptable */
Chris@10: 	       wisdom_state_t owisdom_state = ego->wisdom_state;
Chris@10: 	       
Chris@10: 	       /* this hook is mainly for MPI, to make sure that
Chris@10: 		  wisdom is in sync across all processes for MPI problems */
Chris@10: 	       if (ego->wisdom_ok_hook && !ego->wisdom_ok_hook(p, sol->flags))
Chris@10: 		    goto do_search; /* ignore not-ok wisdom */
Chris@10: 	       
Chris@10: 	       slvndx = SLVNDX(sol);
Chris@10: 	       
Chris@10: 	       if (slvndx == INFEASIBLE_SLVNDX) {
Chris@10: 		    if (ego->wisdom_state == WISDOM_IGNORE_INFEASIBLE)
Chris@10: 			 goto do_search;
Chris@10: 		    else
Chris@10: 			 return 0;   /* known to be infeasible */
Chris@10: 	       }
Chris@10: 	       
Chris@10: 	       flags_of_solution = sol->flags;
Chris@10: 	       
Chris@10: 	       /* inherit blessing either from wisdom
Chris@10: 		  or from the planner */
Chris@10: 	       flags_of_solution.hash_info |= BLISS(ego->flags);
Chris@10: 	       
Chris@10: 	       ego->wisdom_state = WISDOM_ONLY;
Chris@10: 	       
Chris@10: 	       s = ego->slvdescs[slvndx].slv;
Chris@10: 	       if (p->adt->problem_kind != s->adt->problem_kind)
Chris@10: 		    goto wisdom_is_bogus;
Chris@10: 	       
Chris@10: 	       pln = invoke_solver(ego, p, s, &flags_of_solution);
Chris@10: 	       
Chris@10: 	       CHECK_FOR_BOGOSITY; 	  /* catch error in child solvers */
Chris@10: 	       
Chris@10: 	       sol = 0; /* Paranoia: SOL may be dangling after
Chris@10: 			   invoke_solver(); make sure we don't accidentally
Chris@10: 			   reuse it. */
Chris@10: 	       
Chris@10: 	       if (!pln)
Chris@10: 		    goto wisdom_is_bogus;
Chris@10: 	       
Chris@10: 	       ego->wisdom_state = owisdom_state;
Chris@10: 	       
Chris@10: 	       goto skip_search;
Chris@10: 	  }
Chris@10: 	  else if (ego->nowisdom_hook) /* for MPI, make sure lack of wisdom */
Chris@10: 	       ego->nowisdom_hook(p);  /*   is in sync across all processes */
Chris@10:      }
Chris@10: 
Chris@10:  do_search:
Chris@10:      /* cannot search in WISDOM_ONLY mode */
Chris@10:      if (ego->wisdom_state == WISDOM_ONLY)
Chris@10: 	  goto wisdom_is_bogus;
Chris@10: 
Chris@10:      flags_of_solution = ego->flags;
Chris@10:      pln = search(ego, p, &slvndx, &flags_of_solution);
Chris@10:      CHECK_FOR_BOGOSITY; 	  /* catch error in child solvers */
Chris@10: 
Chris@10:      if (ego->timed_out) {
Chris@10: 	  A(!pln);
Chris@10: 	  if (PLNR_TIMELIMIT_IMPATIENCE(ego) != 0) {
Chris@10: 	       /* record (below) that this plan has failed because of
Chris@10: 		  timeout */
Chris@10: 	       flags_of_solution.hash_info |= BLESSING;
Chris@10: 	  } else {
Chris@10: 	       /* this is not the top-level problem or timeout is not
Chris@10: 		  active: record no wisdom. */
Chris@10: 	       return 0;
Chris@10: 	  }
Chris@10:      } else {
Chris@10: 	  /* canonicalize to infinite timeout */
Chris@10: 	  flags_of_solution.timelimit_impatience = 0;
Chris@10:      }
Chris@10: 
Chris@10:  skip_search:
Chris@10:      if (ego->wisdom_state == WISDOM_NORMAL ||
Chris@10: 	 ego->wisdom_state == WISDOM_ONLY) {
Chris@10: 	  if (pln) {
Chris@10: 	       hinsert(ego, m.s, &flags_of_solution, slvndx);
Chris@10: 	       invoke_hook(ego, pln, p, 1);
Chris@10: 	  } else {
Chris@10: 	       hinsert(ego, m.s, &flags_of_solution, INFEASIBLE_SLVNDX);
Chris@10: 	  }
Chris@10:      }
Chris@10: 
Chris@10:      return pln;
Chris@10: 
Chris@10:  wisdom_is_bogus:
Chris@10:      X(plan_destroy_internal)(pln);
Chris@10:      ego->wisdom_state = WISDOM_IS_BOGUS;
Chris@10:      return 0;
Chris@10: }
Chris@10: 
Chris@10: static void htab_destroy(hashtab *ht)
Chris@10: {
Chris@10:      X(ifree)(ht->solutions);
Chris@10:      ht->solutions = 0;
Chris@10:      ht->nelem = 0U;
Chris@10: }
Chris@10: 
Chris@10: static void mkhashtab(hashtab *ht)
Chris@10: {
Chris@10:      ht->nrehash = 0;
Chris@10:      ht->succ_lookup = ht->lookup = ht->lookup_iter = 0;
Chris@10:      ht->insert = ht->insert_iter = ht->insert_unknown = 0;
Chris@10: 
Chris@10:      ht->solutions = 0;
Chris@10:      ht->hashsiz = ht->nelem = 0U;
Chris@10:      hgrow(ht);			/* so that hashsiz > 0 */
Chris@10: }
Chris@10: 
Chris@10: /* destroy hash table entries.  If FORGET_EVERYTHING, destroy the whole
Chris@10:    table.  If FORGET_ACCURSED, then destroy entries that are not blessed. */
Chris@10: static void forget(planner *ego, amnesia a)
Chris@10: {
Chris@10:      switch (a) {
Chris@10: 	 case FORGET_EVERYTHING:
Chris@10: 	      htab_destroy(&ego->htab_blessed);
Chris@10: 	      mkhashtab(&ego->htab_blessed);
Chris@10: 	      /* fall through */
Chris@10: 	 case FORGET_ACCURSED:
Chris@10: 	      htab_destroy(&ego->htab_unblessed);
Chris@10: 	      mkhashtab(&ego->htab_unblessed);
Chris@10: 	      break;
Chris@10: 	 default:
Chris@10: 	      break;
Chris@10:      }
Chris@10: }
Chris@10: 
Chris@10: /* FIXME: what sort of version information should we write? */
Chris@10: #define WISDOM_PREAMBLE PACKAGE "-" VERSION " " STRINGIZE(X(wisdom))
Chris@10: static const char stimeout[] = "TIMEOUT";
Chris@10: 
Chris@10: /* tantus labor non sit cassus */
Chris@10: static void exprt(planner *ego, printer *p)
Chris@10: {
Chris@10:      unsigned h;
Chris@10:      hashtab *ht = &ego->htab_blessed;
Chris@10:      md5 m;
Chris@10: 
Chris@10:      signature_of_configuration(&m, ego);
Chris@10: 
Chris@10:      p->print(p, 
Chris@10: 	      "(" WISDOM_PREAMBLE " #x%M #x%M #x%M #x%M\n",
Chris@10: 	      m.s[0], m.s[1], m.s[2], m.s[3]);
Chris@10: 
Chris@10:      for (h = 0; h < ht->hashsiz; ++h) {
Chris@10: 	  solution *l = ht->solutions + h;
Chris@10: 	  if (LIVEP(l)) {
Chris@10: 	       const char *reg_nam;
Chris@10: 	       int reg_id;
Chris@10: 
Chris@10: 	       if (SLVNDX(l) == INFEASIBLE_SLVNDX) {
Chris@10: 		    reg_nam = stimeout;
Chris@10: 		    reg_id = 0;
Chris@10: 	       } else {
Chris@10: 		    slvdesc *sp = ego->slvdescs + SLVNDX(l);
Chris@10: 		    reg_nam = sp->reg_nam;
Chris@10: 		    reg_id = sp->reg_id;
Chris@10: 	       }
Chris@10: 
Chris@10: 	       /* qui salvandos salvas gratis
Chris@10: 		  salva me fons pietatis */
Chris@10: 	       p->print(p, "  (%s %d #x%x #x%x #x%x #x%M #x%M #x%M #x%M)\n",
Chris@10: 			reg_nam, reg_id, 
Chris@10: 			l->flags.l, l->flags.u, l->flags.timelimit_impatience, 
Chris@10: 			l->s[0], l->s[1], l->s[2], l->s[3]);
Chris@10: 	  }
Chris@10:      }
Chris@10:      p->print(p, ")\n");
Chris@10: }
Chris@10: 
Chris@10: /* mors stupebit et natura
Chris@10:    cum resurget creatura */
Chris@10: static int imprt(planner *ego, scanner *sc)
Chris@10: {
Chris@10:      char buf[MAXNAM + 1];
Chris@10:      md5uint sig[4];
Chris@10:      unsigned l, u, timelimit_impatience;
Chris@10:      flags_t flags;
Chris@10:      int reg_id;
Chris@10:      unsigned slvndx;
Chris@10:      hashtab *ht = &ego->htab_blessed;
Chris@10:      hashtab old;
Chris@10:      md5 m;
Chris@10: 
Chris@10:      if (!sc->scan(sc, 
Chris@10: 		   "(" WISDOM_PREAMBLE " #x%M #x%M #x%M #x%M\n",
Chris@10: 		   sig + 0, sig + 1, sig + 2, sig + 3))
Chris@10: 	  return 0; /* don't need to restore hashtable */
Chris@10: 
Chris@10:      signature_of_configuration(&m, ego);
Chris@10:      if (m.s[0] != sig[0] || m.s[1] != sig[1] ||
Chris@10: 	 m.s[2] != sig[2] || m.s[3] != sig[3]) {
Chris@10: 	  /* invalid configuration */
Chris@10: 	  return 0;
Chris@10:      }
Chris@10:      
Chris@10:      /* make a backup copy of the hash table (cache the hash) */
Chris@10:      {
Chris@10: 	  unsigned h, hsiz = ht->hashsiz;
Chris@10: 	  old = *ht;
Chris@10: 	  old.solutions = (solution *)MALLOC(hsiz * sizeof(solution), HASHT);
Chris@10: 	  for (h = 0; h < hsiz; ++h)
Chris@10: 	       old.solutions[h] = ht->solutions[h];
Chris@10:      }
Chris@10: 
Chris@10:      while (1) {
Chris@10: 	  if (sc->scan(sc, ")"))
Chris@10: 	       break;
Chris@10: 
Chris@10: 	  /* qua resurget ex favilla */
Chris@10: 	  if (!sc->scan(sc, "(%*s %d #x%x #x%x #x%x #x%M #x%M #x%M #x%M)",
Chris@10: 			MAXNAM, buf, &reg_id, &l, &u, &timelimit_impatience,
Chris@10: 			sig + 0, sig + 1, sig + 2, sig + 3))
Chris@10: 	       goto bad;
Chris@10: 
Chris@10: 	  if (!strcmp(buf, stimeout) && reg_id == 0) {
Chris@10: 	       slvndx = INFEASIBLE_SLVNDX;
Chris@10: 	  } else {
Chris@10: 	       if (timelimit_impatience != 0)
Chris@10: 		    goto bad;
Chris@10: 
Chris@10: 	       slvndx = slookup(ego, buf, reg_id);
Chris@10: 	       if (slvndx == INFEASIBLE_SLVNDX)
Chris@10: 		    goto bad;
Chris@10: 	  }
Chris@10: 
Chris@10: 	  /* inter oves locum praesta */
Chris@10: 	  flags.l = l;
Chris@10: 	  flags.u = u;
Chris@10: 	  flags.timelimit_impatience = timelimit_impatience;
Chris@10: 	  flags.hash_info = BLESSING;
Chris@10: 
Chris@10: 	  CK(flags.l == l);
Chris@10: 	  CK(flags.u == u);
Chris@10: 	  CK(flags.timelimit_impatience == timelimit_impatience);
Chris@10: 
Chris@10: 	  if (!hlookup(ego, sig, &flags))
Chris@10: 	       hinsert(ego, sig, &flags, slvndx);
Chris@10:      }
Chris@10: 
Chris@10:      X(ifree0)(old.solutions);
Chris@10:      return 1;
Chris@10: 
Chris@10:  bad:
Chris@10:      /* ``The wisdom of FFTW must be above suspicion.'' */
Chris@10:      X(ifree0)(ht->solutions);
Chris@10:      *ht = old;
Chris@10:      return 0;
Chris@10: }
Chris@10: 
Chris@10: /*
Chris@10:  * create a planner
Chris@10:  */
Chris@10: planner *X(mkplanner)(void)
Chris@10: {
Chris@10:      int i;
Chris@10: 
Chris@10:      static const planner_adt padt = {
Chris@10: 	  register_solver, mkplan, forget, exprt, imprt
Chris@10:      };
Chris@10: 
Chris@10:      planner *p = (planner *) MALLOC(sizeof(planner), PLANNERS);
Chris@10: 
Chris@10:      p->adt = &padt;
Chris@10:      p->nplan = p->nprob = 0;
Chris@10:      p->pcost = p->epcost = 0.0;
Chris@10:      p->hook = 0;
Chris@10:      p->cost_hook = 0;
Chris@10:      p->wisdom_ok_hook = 0;
Chris@10:      p->nowisdom_hook = 0;
Chris@10:      p->bogosity_hook = 0;
Chris@10:      p->cur_reg_nam = 0;
Chris@10:      p->wisdom_state = WISDOM_NORMAL;
Chris@10: 
Chris@10:      p->slvdescs = 0;
Chris@10:      p->nslvdesc = p->slvdescsiz = 0;
Chris@10: 
Chris@10:      p->flags.l = 0;
Chris@10:      p->flags.u = 0;
Chris@10:      p->flags.timelimit_impatience = 0;
Chris@10:      p->flags.hash_info = 0;
Chris@10:      p->nthr = 1;
Chris@10:      p->need_timeout_check = 1;
Chris@10:      p->timelimit = -1;
Chris@10: 
Chris@10:      mkhashtab(&p->htab_blessed);
Chris@10:      mkhashtab(&p->htab_unblessed);
Chris@10: 
Chris@10:      for (i = 0; i < PROBLEM_LAST; ++i)
Chris@10: 	  p->slvdescs_for_problem_kind[i] = -1;
Chris@10: 
Chris@10:      return p;
Chris@10: }
Chris@10: 
Chris@10: void X(planner_destroy)(planner *ego)
Chris@10: {
Chris@10:      /* destroy hash table */
Chris@10:      htab_destroy(&ego->htab_blessed);
Chris@10:      htab_destroy(&ego->htab_unblessed);
Chris@10: 
Chris@10:      /* destroy solvdesc table */
Chris@10:      FORALL_SOLVERS(ego, s, sp, {
Chris@10: 	  UNUSED(sp);
Chris@10: 	  X(solver_destroy)(s);
Chris@10:      });
Chris@10: 
Chris@10:      X(ifree0)(ego->slvdescs);
Chris@10:      X(ifree)(ego); /* dona eis requiem */
Chris@10: }
Chris@10: 
Chris@10: plan *X(mkplan_d)(planner *ego, problem *p)
Chris@10: {
Chris@10:      plan *pln = ego->adt->mkplan(ego, p);
Chris@10:      X(problem_destroy)(p);
Chris@10:      return pln;
Chris@10: }
Chris@10: 
Chris@10: /* like X(mkplan_d), but sets/resets flags as well */
Chris@10: plan *X(mkplan_f_d)(planner *ego, problem *p, 
Chris@10: 		    unsigned l_set, unsigned u_set, unsigned u_reset)
Chris@10: {
Chris@10:      flags_t oflags = ego->flags;
Chris@10:      plan *pln;
Chris@10: 
Chris@10:      PLNR_U(ego) &= ~u_reset;
Chris@10:      PLNR_L(ego) &= ~u_reset;
Chris@10:      PLNR_L(ego) |= l_set;
Chris@10:      PLNR_U(ego) |= u_set | l_set;
Chris@10:      pln = X(mkplan_d)(ego, p);
Chris@10:      ego->flags = oflags;
Chris@10:      return pln;
Chris@10: }
Chris@10: 
Chris@10: /*
Chris@10:  * Debugging code:
Chris@10:  */
Chris@10: #ifdef FFTW_DEBUG
Chris@10: static void check(hashtab *ht)
Chris@10: {
Chris@10:      unsigned live = 0;
Chris@10:      unsigned i;
Chris@10: 
Chris@10:      A(ht->nelem < ht->hashsiz);
Chris@10: 
Chris@10:      for (i = 0; i < ht->hashsiz; ++i) {
Chris@10: 	  solution *l = ht->solutions + i; 
Chris@10: 	  if (LIVEP(l)) 
Chris@10: 	       ++live; 
Chris@10:      }
Chris@10: 
Chris@10:      A(ht->nelem == live);
Chris@10: 
Chris@10:      for (i = 0; i < ht->hashsiz; ++i) {
Chris@10: 	  solution *l1 = ht->solutions + i; 
Chris@10: 	  int foundit = 0;
Chris@10: 	  if (LIVEP(l1)) {
Chris@10: 	       unsigned g, h = h1(ht, l1->s), d = h2(ht, l1->s);
Chris@10: 
Chris@10: 	       g = h;
Chris@10: 	       do {
Chris@10: 		    solution *l = ht->solutions + g;
Chris@10: 		    if (VALIDP(l)) {
Chris@10: 			 if (l1 == l)
Chris@10: 			      foundit = 1;
Chris@10: 			 else if (LIVEP(l) && md5eq(l1->s, l->s)) {
Chris@10: 			      A(!subsumes(&l->flags, SLVNDX(l), &l1->flags));
Chris@10: 			      A(!subsumes(&l1->flags, SLVNDX(l1), &l->flags));
Chris@10: 			 }
Chris@10: 		    } else 
Chris@10: 			 break;
Chris@10: 		    g = addmod(g, d, ht->hashsiz);
Chris@10: 	       } while (g != h);
Chris@10: 
Chris@10: 	       A(foundit);
Chris@10: 	  }
Chris@10:      }
Chris@10: }
Chris@10: #endif