ref: c2c9562e3c2994d87f65ab09779190d1e7e09517
parent: e0be49d7f1dcb50048bead1b7d62633448482246
author: aiju <devnull@localhost>
date: Sat Mar 17 15:26:26 EDT 2018
add libsat
--- /dev/null
+++ b/sys/include/sat.h
@@ -1,0 +1,104 @@
+#pragma lib "libsat.a"
+
+typedef struct SATParam SATParam;
+typedef struct SATClause SATClause;
+typedef struct SATSolve SATSolve;
+typedef struct SATBlock SATBlock;
+typedef struct SATVar SATVar;
+typedef struct SATLit SATLit;
+typedef struct SATConflict SATConflict;
+
+/* user adjustable parameters */
+struct SATParam {+ void (*errfun)(char *, void *);
+ void *erraux;
+ long (*randfn)(void *);
+ void *randaux;
+
+ uint goofprob; /* probability of making a random decision, times 2**31 */
+ double varρ; /* Δactivity is multiplied by this after a conflict */
+ double clauseρ; /* Δclactivity is multiplied by this after a conflict */
+ int trivlim; /* number of extra literals we're willing to tolerate before substituting the trivial clause */
+ int purgeΔ; /* initial purge interval (number of conflicts before a purge) */
+ int purgeδ; /* increase in purge interval at purge */
+ double purgeα; /* α weight factor for purge heuristic */
+ u32int flushψ; /* agility threshold for restarts */
+};
+
+/* each block contains multiple SATClauses consecutively in its data region. each clause is 8 byte aligned and the total size is SATBLOCKSZ (64K) */
+struct SATBlock {+ SATBlock *next, *prev;
+ SATClause *last; /* last clause, ==nil for empty blocks */
+ void *end; /* first byte past the last clause */
+ uchar data[1];
+};
+
+struct SATSolve {+ SATParam;
+
+ uchar unsat; /* ==1 if unsatisfiable. don't even try to solve. */
+ uchar scratched; /* <0 if error happened, state undefined */
+
+ SATBlock bl[2]; /* two doubly linked list heads: list bl[0] contains user clauses, list bl[1] contains learned clauses */
+ SATBlock *lastbl; /* the last block we added a learned clause to */
+ SATClause *cl; /* all clauses are linked together; this is the first user clause */
+ SATClause *learncl; /* first learned clause */
+ SATClause **lastp[2]; /* this points towards the last link in each linked list */
+ int ncl; /* total number of clauses */
+ int ncl0; /* number of user clauses */
+ SATVar *var; /* all variables (array with nvar elements) */
+ SATLit *lit; /* all literals (array with 2*nvar elements) */
+ int nvar;
+ int nvaralloc; /* space allocated for variables */
+ int *trail; /* the trail contains all literals currently assumed true */
+ int ntrail;
+ int *decbd; /* decision boundaries. trail[decbd[i]] has the first literal of level i */
+ int lvl; /* current decision level */
+ SATVar **heap; /* binary heap with free variables, sorted by activity (nonfree variables are removed lazily and so may also be in it) */
+ int nheap;
+
+ uint *lvlstamp; /* used to "stamp" levels during conflict resolution and purging */
+ uint stamp; /* current "stamp" counter */
+
+ int forptr; /* trail[forptr] is the first literal we haven't explored the implications of yet */
+ int binptr; /* ditto for binary implications */
+
+ int *cflcl; /* during conflict resolution we build the learned clause in here */
+ int ncflcl;
+ int cflclalloc; /* space allocated for cflcl */
+ int cfllvl; /* the maximum level of the literals in cflcl, cflcl[0] excluded */
+ int cflctr; /* number of unresolved literals during conflict resolution */
+
+ double Δactivity; /* activity increment for variables */
+ double Δclactivity; /* activity increment for clauses */
+
+ uvlong conflicts; /* total number of conflicts so far */
+
+ uvlong nextpurge; /* purge happens when conflicts >= nextpurge */
+ uint purgeival; /* increment for nextpurge */
+ /* during a purge we do a "full run", assigning all variables and recording conflicts rather than resolving them */
+ SATConflict *fullrcfl; /* conflicts found thus */
+ int nfullrcfl;
+ int fullrlvl; /* minimum cfllvl for conflicts found during purging */
+ int *fullrlits; /* literals implied by conflicts at level fullrlvl */
+ int nfullrlits;
+ int rangecnt[256]; /* number of clauses with certain range values */
+
+ u64int nextflush; /* flush happens when conflicts >= nextflush */
+ u32int flushu, flushv, flushθ; /* variables for flush scheduling algorithm */
+ u32int agility; /* agility is a measure how quickly variables are being flipped. high agility inhibits flushes */
+
+ void *scrap; /* auxiliary memory that may need to be freed after a fatal error */
+};
+
+SATSolve *satnew(void);
+SATSolve *satadd1(SATSolve *, int *, int);
+SATSolve *sataddv(SATSolve *, ...);
+SATSolve *satrange1(SATSolve *, int *, int, int, int);
+SATSolve *satrangev(SATSolve *, int, int, ...);
+int satsolve(SATSolve *);
+int satmore(SATSolve *);
+int satval(SATSolve *, int);
+void satfree(SATSolve *);
+void satreset(SATSolve *);
+
--- /dev/null
+++ b/sys/man/2/sat
@@ -1,0 +1,231 @@
+.TH SAT 2
+.SH NAME
+satnew, satadd1, sataddv, satrange1, satrangev, satsolve, satmore, satval, satreset, satfree \- boolean satisfiability (SAT) solver
+.SH SYNOPSIS
+.de PB
+.PP
+.ft L
+.nf
+..
+.PB
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+.PB
+struct SATParam {+ void (*errfun)(char *msg, void *erraux);
+ void *erraux;
+ long (*randfn)(void *randaux);
+ void *randaux;
+ /* + finetuning parameters, see sat.h */
+};
+.PB
+struct SATSolve {+ SATParam;
+ /* + internals */
+};
+.PB
+.ta +\w'\fLSATSolve* \fP'u
+SATSolve* satnew(void);
+void satfree(SATSolve *s);
+SATSolve* satadd1(SATSolve *s, int *lit, int nlit);
+SATSolve* sataddv(SATSolve *s, ...);
+SATSolve* satrange1(SATSolve *s, int *lit, int nlit,
+ int min, int max);
+SATSolve* satrangev(SATSolve *s, int min, int max, ...);
+int satsolve(SATSolve *s);
+int satmore(SATSolve *s);
+int satval(SATSolve *s, int lit);
+void satreset(SATSolve *s);
+.SH DESCRIPTION
+.PP
+.I Libsat
+is a solver for the boolean satisfiability problem, i.e. given a boolean formula it will either find an assignment to the variables that makes it true, or report that this is impossible.
+The input formula must be in conjunctive normal form (CNF), i.e. of the form
+.IP
+.if t (x\d\s71\s10\u ∨ x\d\s72\s10\u ∨ x\d\s73\s10\u ∨ …) ∧ (y\d\s71\s10\u ∨ y\d\s72\s10\u ∨ y \d\s73\s10\u ∨ …) ∧ …,
+.if n (x1 ∨ x2 ∨ x3 ∨ ...) ∧ (y1 ∨ y2 ∨ y3 ∨ ...) ∧ ...,
+.PP
+where each
+.if t x\d\s7i\s10\u or y\d\s7i\s10\u
+.if n x_i or y_i
+can optionally be negated.
+.PP
+For example, consider
+.IP
+.if t (x\d\s71\s10\u ∨ x\d\s72\s10\u ∨ x\d\s73\s10\u) ∧ (¬x\d\s71\s10\u ∨ ¬x\d\s72\s10\u) ∧ (¬x\d\s72\s10\u ∨ ¬x\d\s73\s10\u) ∧ (¬x\d\s71\s10\u ∨ ¬x\d\s73\s10\u).
+.if n (x1 ∨ x2 ∨ x3) ∧ (¬x1 ∨ ¬x2) ∧ (¬x2 ∨ ¬x3) ∧ (¬x1 ∨ ¬x3).
+.PP
+This formula encodes the constraint that exactly one of the three variables be true. To represent this as input for
+.I libsat
+we assign positive integers to each variable.
+Negation is represented by the corresponding negative number, hence our example corresponds to the set of "clauses"
+.IP
+1, 2, 3
+.br
+-1, -2
+.br
+-1, -3
+.br
+-2, -3
+.PP
+To actually solve this problem we would create a
+.B SATSolve
+structure and add clauses one by one using
+.I satadd1
+or
+.I sataddv
+(the former takes an
+.B int
+array, the latter a variadic list terminated by 0).
+The
+.B SATSolve
+is modified inplace but returned for convenience.
+Passing
+.B nil
+as a first argument will create and return a new structure.
+Alternatively,
+.I satnew
+will create an empty structure.
+.PP
+Once clauses have been added,
+.I satsolve
+will invoke the actual solver.
+It returns 1 if it found an assignment and 0 if there is no assignment (the formula is unsatisfiable).
+If an assignment has been found, further clauses can be added to constrain it further and
+.I satsolve
+rerun.
+.I Satmore
+performs this automatically, excluding the current values of the variables.
+It is equivalent to
+.I satsolve
+if no variables have assigned values.
+.PP
+Once a solution has been found,
+.I satval
+returns the value of literal
+.I lit.
+It returns 1 for true, 0 for false and -1 for undetermined.
+If the formula is satisfiable, an undetermined variable is one where either value will satisfy the formula.
+If the formula is unsatisfiable, all variables are undetermined.
+.PP
+.I Satrange1
+and
+.I satrangev
+function like their
+.I satadd
+brethren but rather than adding a single clause they add multiple clauses corresponding to the constraint that at least
+.I min
+and at most
+.I max
+literals from the provided array be true.
+For example, the clause from above corresponds to
+.IP
+.B "satrangev(s, 1, 1, 1, 2, 3, 0);"
+.PP
+.I Satreset
+resets all solver state, deleting all learned clauses and variable assignments.
+It retains all user provided clauses.
+.I Satfree
+deletes a solver structure and frees all associated storage.
+.PP
+There are a number of user-adjustable parameters in the
+.B SATParam
+structure embedded in
+.BR SATSolve .
+.I Randfun
+is called with argument
+.I randaux
+to generate random numbers between 0 and
+.if t 2\u\s731\s10\d-1;
+.if n 2^31-1;
+it defaults to
+.IR lrand (2).
+.I Errfun
+is called on fatal errors (see DIAGNOSTICS).
+Additionally, a number of finetuning parameters are defined in
+.BR sat.h .
+By tweaking their values, the run-time for a given problem can be reduced.
+.SH EXAMPLE
+Find all solutions to the example clause from above:
+.PB
+.ta .5i 1i 1.5i
+ SATSolve *s;
+
+ s = nil;
+ s = sataddv(s, 1, 2, 3, 0);
+ s = sataddv(s, -1, -2, 0);
+ s = sataddv(s, -1, -3, 0);
+ s = sataddv(s, -2, -3, 0);
+ while(satmore(s) > 0)
+ print("x1=%d x2=%d x3=%d\\n",+ satval(s, 1),
+ satval(s, 2),
+ satval(s, 3));
+ satfree(s);
+.SH SOURCE
+.B /sys/src/libsat
+.SH "SEE ALSO"
+Donald Knuth, ``The Art of Computer Programming'', Volume 4, Fascicle 6.
+.SH DIAGNOSTICS
+.I Satnew
+returns
+.B nil
+on certain fatal error conditions (such as
+.IR malloc (2)
+failure).
+Other routines will call
+.I errfun
+with an error string and
+.IR erraux .
+If no
+.I errfun
+is provided or if it returns,
+.IR sysfatal (2)
+is called.
+It is permissible to use
+.IR setjmp (2)
+to return from an error condition.
+Call
+.I satfree
+to clean up the
+.B SATSolve
+structure in this case.
+Note that calling the
+.I satadd
+or
+.I satrange
+routines with
+.B nil
+first argument will invoke
+.I sysfatal
+on error, since no
+.I errfun
+has been defined yet.
+.SH BUGS
+Variable numbers should be consecutive numbers starting from 1, since variable data is kept in arrays internally.
+.PP
+Large clauses of several thousand literals are probably inefficient and should be split up using auxiliary variables.
+Very large clauses exceeding about 16,000 literals will not work at all.
+.PP
+There is no way to remove clauses (since it's unclear what the semantics should be).
+.PP
+The details about the tuning parameters are subject to change.
+.PP
+Calling
+.I satadd
+or
+.I satrange
+after
+.I satsolve
+or
+.I satmore
+may reset variable values.
+.PP
+.I Satmore
+will always return 1 when there are no assigned variables in the solution.
+.P
+Some debugging routines called under "shouldn't happen" conditions are non-reentrant.
+.SH HISTORY
+.I Libsat
+first appeared in 9front in March, 2018.
--- /dev/null
+++ b/sys/src/libsat/debug.c
@@ -1,0 +1,130 @@
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+#include "impl.h"
+
+static SATSolve *globsat;
+
+static int
+satclausefmt(Fmt *f)
+{+ SATClause *c;
+ char *s;
+ int i, fl;
+
+ fl = f->flags;
+ c = va_arg(f->args, SATClause *);
+ if(c == nil)
+ return fmtstrcpy(f, "Λ");
+ if(c->n == 0)
+ return fmtstrcpy(f, "ε");
+ s = "%s%d";
+ for(i = 0; i < c->n; i++){+ if((fl & FmtSign) != 0)
+ switch(globsat->lit[c->l[i]].val){+ case 1: s = "%s[%d]"; break;
+ case 0: s = "%s(%d)"; break;
+ case -1: s = "%s%d"; break;
+ default: abort();
+ }
+ fmtprint(f, s, i != 0 ? " ∨ " : "", signf(c->l[i]));
+ }
+ return 0;
+}
+
+void
+satprintstate(SATSolve *s)
+{+ int i;
+ Fmt f;
+ char buf[512];
+ SATVar *v;
+
+ fmtfdinit(&f, 1, buf, sizeof(buf));
+ fmtprint(&f, "trail:\n");
+ for(i = 0; i < s->ntrail; i++){+ v = &s->var[VAR(s->trail[i])];
+ fmtprint(&f, "%c%-8d %- 8d %-8d ", i == s->forptr ? '*' : ' ', i, signf(s->trail[i]), v->lvl);
+ if(v->isbinreason)
+ fmtprint(&f, "%d ∨ %d\n", signf(s->trail[i]), signf(v->binreason));
+ else
+ fmtprint(&f, "%+Γ\n", v->reason);
+ }
+ fmtrune(&f, '\n');
+ fmtfdflush(&f);
+}
+
+void
+satsanity(SATSolve *s)
+{+ int i, j, k, m, tl, s0, s1;
+ SATVar *v;
+ SATLit *l;
+ SATClause *c;
+
+ for(c = s->cl; c != nil; c = c->next){+ assert(c->n >= 2);
+ assert((uint)((uchar*)c->next - (uchar*)c) >= sizeof(SATClause) + (c->n - 1) * sizeof(int));
+ for(j = 0; j < c->n; j++)
+ assert((uint)c->l[j] < 2*s->nvar);
+ for(i = 0; i < 2; i++)
+ c->watch[i] = (void*)((uintptr)c->watch[i] | 1);
+ }
+ for(i = 0; i < s->nvar; i++){+ tl = -1;
+ for(j = 0; j < s->ntrail; j++)
+ if(VAR(s->trail[j]) == i){+ assert(tl == -1);
+ tl = j;
+ }
+ v = &s->var[i];
+ l = &s->lit[2*i];
+ if(l->val >= 0){+ assert(l->val <= 1);
+ assert(l[0].val + l[1].val == 1);
+ assert((uint)v->lvl <= s->lvl);
+ assert(tl != -1);
+ assert(s->trail[tl] == 2*i+l[1].val);
+ assert(tl >= s->decbd[v->lvl]);
+ assert(v->lvl == s->lvl || tl < s->decbd[v->lvl+1]);
+ }else{+ assert(l[0].val == -1 && l[1].val == -1);
+ assert(v->lvl == -1);
+ assert(v->heaploc >= 0);
+ assert(tl == -1);
+ }
+ assert(v->heaploc == -1 || (uint)v->heaploc <= s->nheap && s->heap[v->heaploc] == v);
+ for(j = 0; j < 2; j++){+ m = 2 * i + j;
+ for(c = l[j].watch; c != nil; c = c->watch[k]){+ k = c->l[1] == m;
+ assert(k || c->l[0] == m);
+ assert((uintptr)c->watch[k] & 1);
+ c->watch[k] = (void*)((uintptr)c->watch[k] & ~1);
+ }
+ }
+ }
+ for(c = s->cl; c != nil; c = c->next)
+ for(i = 0; i < 2; i++)
+ assert(((uintptr)c->watch[i] & 1) == 0);
+ if(s->forptr == s->ntrail)
+ for(c = s->cl; c != nil; c = c->next){+ s0 = s->lit[c->l[0]].val;
+ s1 = s->lit[c->l[1]].val;
+ if(s0 != 0 && s1 != 0 || s0 == 1 || s1 == 1)
+ continue;
+ for(i = 2; i < c->n; i++)
+ if(s->lit[c->l[i]].val != 0){+ satprintstate(s);
+ print("watchlist error: %+Γ\n", c);+ assert(0);
+ }
+ }
+}
+
+void
+satdebuginit(SATSolve *s)
+{+ globsat = s;
+ fmtinstall(L'Γ', satclausefmt);
+}
--- /dev/null
+++ b/sys/src/libsat/impl.h
@@ -1,0 +1,81 @@
+/* note that internally, literals use a representation different from the API.
+ * variables are numbered from 0 (not 1) and 2v and 2v+1 correspond to v
+ * and ¬v, resp. */
+#define VAR(l) ((l)>>1)
+#define NOT(l) ((l)^1)
+
+/* l[0] and l[1] are special: they are the watched literals.
+ * all clauses that have literal l on their watchlist form a linked list starting with s->lit[l].watch
+ * and watch[i] having the next clause for l[i] */
+struct SATClause {+ SATClause *next;
+ SATClause *watch[2];
+ double activity; /* activity is increased every time a clause is used to resolve a conflict (tiebreaking heuristic during purging) */
+ int n; /* >= 2 for learned clauses and > 2 for input clauses (binary input clauses are kept in the bimp tables) */
+ ushort range; /* heuristic used during purging, low range => keep clause (range 0..256) */
+ int l[1];
+};
+
+struct SATLit {+ int *bimp; /* array of literals implied by this literal through binary clauses (Binary IMPlications) */
+ SATClause *watch; /* linked list of watched clauses */
+ int nbimp;
+ char val; /* -1 = not assigned, 0 = false, 1 = true */
+};
+
+struct SATVar {+ double activity; /* activity is increased every time a variable shows up in a conflict */
+ union {+ SATClause *reason; /* nil for decision and free literals */
+ int binreason; /* used when isbinreason == 1: the reason is the clause l ∨ l->binreason */
+ };
+ int lvl; /* level at which this variable is defined, or -1 for free variables */
+ int heaploc; /* location in binary heap or -1 when not in heap */
+ uint stamp; /* "stamp" value used for conflict resolution etc. */
+ uchar flags; /* see below */
+ char isbinreason;
+};
+
+enum {+ VARPHASE = 1, /* for a free variables, VARPHASE is used as a first guess the next time it is picked */
+ VARUSER = 0x80, /* user assigned variable (unit clause in input) */
+};
+
+/* records conflicts during purging */
+struct SATConflict {+ union {+ SATClause *c;
+ uvlong b;
+ };
+ int lvl; /* bit 31 denotes binary conflict */
+};
+#define CFLLVL(c) ((c).lvl & 0x7fffffff)
+
+enum {+ SATBLOCKSZ = 65536,
+ SATVARALLOC = 64,
+ CLAUSEALIGN = 8,
+ CFLCLALLOC = 16,
+};
+
+void saterror(SATSolve *, char *, ...);
+void sataddtrail(SATSolve *, int);
+void satdebuginit(SATSolve *);
+void satprintstate(SATSolve *);
+void satsanity(SATSolve *);
+SATVar *satheaptake(SATSolve *);
+void satheapput(SATSolve *, SATVar *);
+void satreheap(SATSolve *, SATVar *);
+void satheapreset(SATSolve *);
+int satnrand(SATSolve *, int);
+void *satrealloc(SATSolve *, void *, ulong);
+SATClause *satnewclause(SATSolve *, int, int);
+SATClause *satreplclause(SATSolve *, int);
+void satcleanup(SATSolve *, int);
+void satbackjump(SATSolve *, int);
+
+#define signf(l) (((l)<<31>>31|1)*((l)+2>>1))
+#pragma varargck type "Γ" SATClause*
+
+#define ε 2.2250738585072014e-308
+#define MAXACTIVITY 1e100
--- /dev/null
+++ b/sys/src/libsat/misc.c
@@ -1,0 +1,206 @@
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+#include "impl.h"
+
+SATSolve *
+satnew(void)
+{+ SATSolve *s;
+
+ s = calloc(1, sizeof(SATSolve));
+ if(s == nil) return nil;
+ s->bl[0].next = s->bl[0].prev = &s->bl[0];
+ s->bl[1].next = s->bl[1].prev = &s->bl[1];
+ s->bl[0].end = (uchar*)&s->bl[0] + SATBLOCKSZ; /* this block is "full" */
+ s->bl[1].end = (uchar*)&s->bl[1] + SATBLOCKSZ;
+ s->lastp[0] = &s->cl;
+ s->lastp[1] = &s->learncl;
+ s->lastbl = &s->bl[1];
+ s->randfn = (long(*)(void*)) lrand;
+
+ s->goofprob = 0.02 * (1UL<<31);
+ s->varρ = 1/0.9;
+ s->clauseρ = 1/0.999;
+ s->trivlim = 10;
+ s->purgeΔ = 10000;
+ s->purgeδ = 100;
+ s->purgeα = 0.2;
+ s->flushψ = (1ULL<<32)*0.05;
+
+ s->Δactivity = 1;
+ s->Δclactivity = 1;
+
+ return s;
+}
+
+void
+satfree(SATSolve *s)
+{+ SATBlock *b, *bb;
+ int i;
+
+ if(s == nil) return;
+ for(i = 0; i < 2; i++)
+ for(b = s->bl[i].next; b != &s->bl[i]; b = bb){+ bb = b->next;
+ free(b);
+ }
+ for(i = 0; i < 2 * s->nvar; i++)
+ free(s->lit[i].bimp);
+ free(s->heap);
+ free(s->trail);
+ free(s->decbd);
+ free(s->var);
+ free(s->lit);
+ free(s->cflcl);
+ free(s->fullrcfl);
+ free(s->fullrlits);
+ free(s->scrap);
+ free(s);
+}
+
+void
+saterror(SATSolve *s, char *msg, ...)
+{+ char buf[ERRMAX];
+ va_list va;
+
+ va_start(va, msg);
+ vsnprint(buf, sizeof(buf), msg, va);
+ va_end(va);
+ s->scratched = 1;
+ if(s != nil && s->errfun != nil)
+ s->errfun(buf, s->erraux);
+ sysfatal("%s", buf);+}
+
+int
+satval(SATSolve *s, int l)
+{+ int m, v;
+
+ if(s->unsat) return -1;
+ m = l >> 31;
+ v = (l + m ^ m) - 1;
+ if(v < 0 || v >= s->nvar) return -1;
+ return s->lit[2*v+(m&1)].val;
+}
+
+int
+satnrand(SATSolve *s, int n)
+{+ long slop, v;
+
+ if(n <= 1) return 0;
+ slop = 0x7fffffff % n;
+ do
+ v = s->randfn(s->randaux);
+ while(v <= slop);
+ return v % n;
+}
+
+void *
+satrealloc(SATSolve *s, void *v, ulong n)
+{+ v = realloc(v, n);
+ if(v == nil)
+ saterror(s, "realloc: %r");
+ return v;
+}
+
+#define LEFT(x) (2*(x)+1)
+#define RIGHT(x) (2*(x)+2)
+#define UP(x) ((x)-1>>1)
+
+static SATVar *
+heapswap(SATSolve *s, int i, int j)
+{+ SATVar *r;
+
+ if(i == j) return s->heap[i];
+ r = s->heap[i];
+ s->heap[i] = s->heap[j];
+ s->heap[j] = r;
+ s->heap[i]->heaploc = i;
+ s->heap[j]->heaploc = j;
+ return r;
+}
+
+static void
+heapup(SATSolve *s, int i)
+{+ int m;
+
+ m = i;
+ for(;;){+ if(LEFT(i) < s->nheap && s->heap[LEFT(i)]->activity > s->heap[m]->activity)
+ m = LEFT(i);
+ if(RIGHT(i) < s->nheap && s->heap[RIGHT(i)]->activity > s->heap[m]->activity)
+ m = RIGHT(i);
+ if(i == m) break;
+ heapswap(s, m, i);
+ i = m;
+ }
+}
+
+static void
+heapdown(SATSolve *s, int i)
+{+ int p;
+
+ for(; i > 0 && s->heap[p = UP(i)]->activity < s->heap[i]->activity; i = p)
+ heapswap(s, i, p);
+}
+
+SATVar *
+satheaptake(SATSolve *s)
+{+ SATVar *r;
+
+ assert(s->nheap > 0);
+ r = heapswap(s, 0, --s->nheap);
+ heapup(s, 0);
+ r->heaploc = -1;
+ return r;
+}
+
+void
+satheapput(SATSolve *s, SATVar *v)
+{+ assert(s->nheap < s->nvar);
+ v->heaploc = s->nheap;
+ s->heap[s->nheap++] = v;
+ heapdown(s, s->nheap - 1);
+}
+
+void
+satreheap(SATSolve *s, SATVar *v)
+{+ int i;
+
+ i = v->heaploc;
+ if(i < 0) return;
+ heapup(s, i);
+ heapdown(s, i);
+}
+
+void
+satheapreset(SATSolve *s)
+{+ int i, n, j;
+
+ s->heap = satrealloc(s, s->heap, s->nvar * sizeof(SATVar *));
+ n = s->nvar;
+ s->nheap = n;
+ for(i = 0; i < n; i++){+ s->heap[i] = &s->var[i];
+ s->heap[i]->heaploc = i;
+ }
+ for(i = 0; i < n - 1; i++){+ j = i + satnrand(s, n - i);
+ heapswap(s, i, j);
+ heapdown(s, i);
+ }
+ heapdown(s, n - 1);
+}
--- /dev/null
+++ b/sys/src/libsat/mkfile
@@ -1,0 +1,22 @@
+</$objtype/mkfile
+
+LIB=/$objtype/lib/libsat.a
+
+OFILES=\
+ misc.$O \
+ satadd.$O \
+ satrange.$O \
+ satsolve.$O \
+ satmore.$O \
+ debug.$O \
+
+HFILES=\
+ /sys/include/sat.h\
+ impl.h\
+
+UPDATE=\
+ mkfile\
+ $HFILES\
+ ${OFILES:%.$O=%.c}\+
+</sys/src/cmd/mksyslib
--- /dev/null
+++ b/sys/src/libsat/satadd.c
@@ -1,0 +1,230 @@
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+#include "impl.h"
+
+static SATBlock *
+newblock(SATSolve *s, int learned)
+{+ SATBlock *b;
+
+ b = calloc(1, SATBLOCKSZ);
+ if(b == nil)
+ saterror(s, "malloc: %r");
+ b->prev = s->bl[learned].prev;
+ b->next = &s->bl[learned];
+ b->next->prev = b;
+ b->prev->next = b;
+ b->end = (void*) b->data;
+ return b;
+}
+
+SATClause *
+satnewclause(SATSolve *s, int n, int learned)
+{+ SATBlock *b;
+ SATClause *c;
+ int f, sz;
+
+ sz = sizeof(SATClause) + (n - 1) * sizeof(int);
+ assert(sz <= SATBLOCKSZ);
+ if(learned)
+ b = s->lastbl;
+ else
+ b = s->bl[0].prev;
+ for(;;){+ f = (uchar*)b + SATBLOCKSZ - (uchar*)b->end;
+ if(f >= sz) break;
+ b = b->next;
+ if(b == &s->bl[learned])
+ b = newblock(s, learned);
+ }
+ c = b->end;
+ memset(c, 0, sizeof(SATClause));
+ b->end = (void *)((uintptr)b->end + sz + CLAUSEALIGN - 1 & -CLAUSEALIGN);
+ b->last = c;
+ if(learned){+ if(s->lastp[1] == &s->learncl)
+ *s->lastp[0] = c;
+ s->lastbl = b;
+ }else
+ c->next = s->learncl;
+ *s->lastp[learned] = c;
+ s->lastp[learned] = &c->next;
+ s->ncl++;
+ return c;
+}
+
+/* this is currently only used to subsume clauses, i.e. n is guaranteed to be less than the last n */
+SATClause *
+satreplclause(SATSolve *s, int n)
+{+ SATBlock *b;
+ SATClause *c, **wp;
+ int f, sz, i, l;
+
+ assert(s->lastbl != nil && s->lastbl->last != nil);
+ b = s->lastbl;
+ c = b->last;
+ f = (uchar*)b + SATBLOCKSZ - (uchar*)c;
+ sz = sizeof(SATClause) + (n - 1) * sizeof(int);
+ assert(f >= sz);
+ b->end = (void *)((uintptr)c + sz + CLAUSEALIGN - 1 & -CLAUSEALIGN);
+ for(i = 0; i < 2; i++){+ l = c->l[i];
+ for(wp = &s->lit[l].watch; *wp != nil && *wp != c; wp = &(*wp)->watch[(*wp)->l[1] == l])
+ ;
+ assert(*wp != nil);
+ *wp = c->watch[i];
+ }
+ memset(c, 0, sizeof(SATClause));
+ return c;
+}
+
+static int
+litconv(SATSolve *s, int l)
+{+ int v, m, n;
+ SATVar *vp;
+ SATLit *lp;
+
+ m = l >> 31;
+ v = (l + m ^ m) - 1;
+ if(v >= s->nvaralloc){+ n = -(-(v+1) & -SATVARALLOC);
+ s->var = vp = satrealloc(s, s->var, n * sizeof(SATVar));
+ s->lit = lp = satrealloc(s, s->lit, 2 * n * sizeof(SATLit));
+ memset(vp += s->nvaralloc, 0, (n - s->nvaralloc) * sizeof(SATVar));
+ memset(lp += 2*s->nvaralloc, 0, 2 * (n - s->nvaralloc) * sizeof(SATLit));
+ for(; vp < s->var + n; vp++){+ vp->lvl = -1;
+ vp->flags = VARPHASE;
+ }
+ for(; lp < s->lit + 2 * n; lp++)
+ lp->val = -1;
+ s->nvaralloc = n;
+ }
+ if(v >= s->nvar)
+ s->nvar = v + 1;
+ return v << 1 | m & 1;
+}
+
+static void
+addbimp(SATSolve *s, int l0, int l1)
+{+ SATLit *lp;
+
+ lp = &s->lit[NOT(l0)];
+ lp->bimp = satrealloc(s, lp->bimp, (lp->nbimp + 1) * sizeof(int));
+ lp->bimp[lp->nbimp++] = l1;
+}
+
+static SATSolve *
+satadd1special(SATSolve *s, int *a, int n)
+{+ int i, l0, l1;
+
+ if(n == 0){+ s->unsat = 1;
+ return s;
+ }
+ l0 = a[0];
+ l1 = 0;
+ for(i = 1; i < n; i++)
+ if(a[i] != l0){+ l1 = a[i];
+ break;
+ }
+ if(l1 == 0){+ l0 = litconv(s, l0);
+ assert(s->lvl == 0);
+ switch(s->lit[l0].val){+ case 0:
+ s->unsat = 1;
+ return s;
+ case -1:
+ s->trail = satrealloc(s, s->trail, sizeof(int) * s->nvar);
+ memmove(&s->trail[1], s->trail, sizeof(int) * s->ntrail);
+ s->trail[0] = l0;
+ s->ntrail++;
+ s->var[VAR(l0)].flags |= VARUSER;
+ s->var[VAR(l0)].lvl = 0;
+ s->lit[l0].val = 1;
+ s->lit[NOT(l0)].val = 0;
+ }
+ return s;
+ }
+ if(l0 + l1 == 0) return s;
+ l0 = litconv(s, l0);
+ l1 = litconv(s, l1);
+ addbimp(s, l0, l1);
+ addbimp(s, l1, l0);
+ return s;
+}
+
+SATSolve *
+satadd1(SATSolve *s, int *a, int n)
+{+ SATClause *c;
+ int i, j, l, u;
+ SATVar *v;
+
+ if(s == nil){+ s = satnew();
+ if(s == nil)
+ saterror(nil, "satnew: %r");
+ }
+ if(n < 0)
+ for(n = 0; a[n] != 0; n++)
+ ;
+ for(i = 0; i < n; i++)
+ if(a[i] == 0)
+ saterror(s, "satadd1(%p, %p, %d): a[%d]==0, callerpc=%p", s, a, n, i, getcallerpc(&s));
+ satbackjump(s, 0);
+ if(n <= 2)
+ return satadd1special(s, a, n);
+ /* use stamps to detect repeated literals and tautological clauses */
+ if(s->stamp >= (uint)-6){+ for(i = 0; i < s->nvar; i++)
+ s->var[i].stamp = 0;
+ s->stamp = 1;
+ }else
+ s->stamp += 3;
+ u = 0;
+ for(i = 0; i < n; i++){+ l = litconv(s, a[i]);
+ v = &s->var[VAR(l)];
+ if(v->stamp < s->stamp) u++;
+ if(v->stamp == s->stamp + (~l & 1))
+ return s; /* tautological */
+ v->stamp = s->stamp + (l & 1);
+ }
+ if(u <= 2)
+ return satadd1special(s, a, n);
+ s->stamp += 3;
+ c = satnewclause(s, u, 0);
+ c->n = u;
+ for(i = 0, j = 0; i < n; i++){+ l = litconv(s, a[i]);
+ v = &s->var[VAR(l)];
+ if(v->stamp < s->stamp){+ c->l[j++] = l;
+ v->stamp = s->stamp;
+ }
+ }
+ assert(j == u);
+ s->ncl0++;
+ return s;
+}
+
+SATSolve *
+sataddv(SATSolve *s, ...)
+{+ va_list va;
+
+ va_start(va, s);
+ /* horrible hack */
+ s = satadd1(s, (int*)va, -1);
+ va_end(va);
+ return s;
+}
--- /dev/null
+++ b/sys/src/libsat/satmore.c
@@ -1,0 +1,24 @@
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+#include "impl.h"
+
+int
+satmore(SATSolve *s)
+{+ int *a, i, n;
+
+ if(s == nil) return 1;
+ s->scrap = a = satrealloc(s, nil, s->nvar * sizeof(int));
+ n = 0;
+ for(i = 0; i < s->nvar; i++)
+ switch(s->lit[2*i].val){+ case 0: a[n++] = i+1; break;
+ case 1: a[n++] = -(i+1); break;
+ }
+ if(n > 0)
+ satadd1(s, a, n);
+ free(a);
+ s->scrap = nil;
+ return satsolve(s);
+}
--- /dev/null
+++ b/sys/src/libsat/satrange.c
@@ -1,0 +1,68 @@
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+#include "impl.h"
+
+static SATSolve *
+satmin(SATSolve *s, int *a, int n, int *id, int *l, int m, int mul)
+{+ int i;
+
+ if(m > n) return s;
+ for(i = 0; i < m; i++)
+ id[i] = i;
+ for(;;){+ for(i = 0; i < m; i++)
+ l[i] = a[id[i]] * mul;
+ s = satadd1(s, l, m);
+ for(i = m-1; i >= 0; i--){+ if(++id[i] < n+i+1-m)
+ break;
+ if(i == 0)
+ return s;
+ }
+ while(++i < m)
+ id[i] = id[i-1]+1;
+ }
+}
+
+SATSolve *
+satrange1(SATSolve *s, int *a, int n, int min, int max)
+{+ int sz, na;
+
+ if(s == nil){+ s = satnew();
+ if(s == nil)
+ saterror(nil, "satnew: %r");
+ }
+ if(n < 0)
+ for(n = 0; a[n] != 0; n++)
+ ;
+ if(min > n || max < 0)
+ return sataddv(s, 0);
+ if(min < 0) min = 0;
+ if(max > n) max = n;
+ sz = n+1-min;
+ if(min == 0 || max != n && sz < max+1) sz = max+1;
+ if(s->cflclalloc < 2*sz){+ na = -(-2*sz & -CFLCLALLOC);
+ s->cflcl = satrealloc(s, s->cflcl, na * sizeof(int));
+ s->cflclalloc = na;
+ }
+ s = satmin(s, a, n, s->cflcl, s->cflcl+sz, max+1, -1);
+ s = satmin(s, a, n, s->cflcl, s->cflcl+sz, n+1-min, 1);
+ return s;
+}
+
+SATSolve *
+satrangev(SATSolve *s, int min, int max, ...)
+{+ va_list va;
+
+ va_start(va, max);
+ /* horrible hack */
+ s = satrange1(s, (int*)va, -1, min, max);
+ va_end(va);
+ return s;
+}
--- /dev/null
+++ b/sys/src/libsat/satsolve.c
@@ -1,0 +1,887 @@
+#include <u.h>
+#include <libc.h>
+#include <sat.h>
+#include "impl.h"
+
+/* the solver follows Algorithm C from Knuth's The Art of Computer Programming, Vol. 4, Fascicle 6 */
+
+#define verbosestate 0
+#define verboseforcing 0
+#define verboseconflict 0
+#define paranoia 0
+#define sanity(s) if(paranoia) satsanity(s)
+
+void
+sataddtrail(SATSolve *s, int l)
+{+ s->trail[s->ntrail++] = l;
+ s->lit[l].val = 1;
+ s->lit[NOT(l)].val = 0;
+ s->var[VAR(l)].lvl = s->lvl;
+ s->agility -= s->agility >> 13;
+ if(((s->var[VAR(l)].flags ^ l) & 1) != 0)
+ s->agility += 1<<19;
+ if(verbosestate) satprintstate(s);
+}
+
+/* compute watchlists from scratch */
+static void
+rewatch(SATSolve *s)
+{+ SATLit *l;
+ SATClause *c;
+ int i, j, x;
+
+ for(l = s->lit; l < s->lit + 2*s->nvar; l++)
+ l->watch = nil;
+ for(c = s->cl; c != nil; c = c->next)
+ for(i = 0; i < 2; i++){+ if(s->lit[c->l[i]].val == 0)
+ for(j = 2; j < c->n; j++)
+ if(s->lit[c->l[j]].val != 0){+ x = c->l[i], c->l[i] = c->l[j], c->l[j] = x;
+ break;
+ }
+ c->watch[i] = s->lit[c->l[i]].watch;
+ s->lit[c->l[i]].watch = c;
+ }
+}
+
+/* jump back to decision level d */
+void
+satbackjump(SATSolve *s, int d)
+{+ int l;
+ SATVar *v;
+
+ if(s->lvl == d) return;
+ while(s->ntrail > s->decbd[d + 1]){+ l = s->trail[--s->ntrail];
+ v = &s->var[VAR(l)];
+ if((v->flags & VARUSER) != 0){ /* don't delete user assignments */+ s->ntrail++;
+ break;
+ }
+ s->lit[l].val = -1;
+ s->lit[NOT(l)].val = -1;
+ v->flags = v->flags & ~1 | l & 1;
+ v->lvl = -1;
+ v->reason = nil;
+ v->isbinreason = 0;
+ if(v->heaploc < 0)
+ satheapput(s, v);
+ }
+ s->lvl = d;
+ if(s->forptr > s->ntrail) s->forptr = s->ntrail;
+ if(s->binptr > s->ntrail) s->binptr = s->ntrail;
+ if(verbosestate) satprintstate(s);
+}
+
+static void
+solvinit(SATSolve *s)
+{+ satdebuginit(s);
+ satheapreset(s);
+ s->decbd = satrealloc(s, s->decbd, s->nvar * sizeof(int));
+ s->decbd[0] = 0;
+ s->trail = satrealloc(s, s->trail, sizeof(int) * s->nvar);
+ s->fullrlits = satrealloc(s, s->fullrlits, sizeof(int) * s->nvar);
+ s->lvlstamp = satrealloc(s, s->lvlstamp, sizeof(int) * s->nvar);
+ memset(s->lvlstamp, 0, sizeof(int) * s->nvar);
+ if(s->cflclalloc == 0){+ s->cflcl = satrealloc(s, s->cflcl, CFLCLALLOC * sizeof(int));
+ s->cflclalloc = CFLCLALLOC;
+ }
+ rewatch(s);
+
+ s->conflicts = 0;
+ s->nextpurge = s->purgeΔ;
+ s->purgeival = s->purgeΔ;
+ s->nextflush = 1;
+ s->flushu = 1;
+ s->flushv = 1;
+ s->flushθ = s->flushψ;
+ s->agility = 0;
+
+ satbackjump(s, 0);
+ s->forptr = 0;
+ s->binptr = 0;
+}
+
+void
+satcleanup(SATSolve *s, int all)
+{+ SATBlock *b, *bn;
+
+ if(all){+ *s->lastp[0] = nil;
+ s->learncl = nil;
+ s->lastp[1] = &s->learncl;
+ s->ncl = s->ncl0;
+ }
+ for(b = s->bl[1].next; b != &s->bl[1]; b = bn){+ bn = b->next;
+ if(b->last != nil && !all) continue;
+ b->next->prev = b->prev;
+ b->prev->next = b->next;
+ free(b);
+ }
+ s->lastbl = s->bl[1].prev;
+ free(s->fullrlits);
+ s->fullrlits = nil;
+ free(s->lvlstamp);
+ s->lvlstamp = nil;
+ free(s->cflcl);
+ s->cflcl = nil;
+ s->cflclalloc = 0;
+}
+
+static void
+stampoverflow(SATSolve *s)
+{+ int i;
+
+ for(i = 0; i < s->nvar; i++){+ s->var[i].stamp = 0;
+ s->lvlstamp[i] = 0;
+ }
+ s->stamp = -2;
+}
+
+/* "bump" the variable, i.e. increase its activity score. reduce all score when one exceeds MAXACTIVITY (1e100) */
+static void
+varbump(SATSolve *s, SATVar *v)
+{+ v->activity += s->Δactivity;
+ satreheap(s, v);
+ if(v->activity < MAXACTIVITY) return;
+ for(v = s->var; v < s->var + s->nvar; v++)
+ if(v->activity != 0){+ v->activity /= MAXACTIVITY;
+ if(v->activity < ε)
+ v->activity = ε;
+ }
+ s->Δactivity /= MAXACTIVITY;
+}
+
+/* ditto for clauses */
+static void
+clausebump(SATSolve *s, SATClause *c)
+{+ c->activity += s->Δclactivity;
+ if(c->activity < MAXACTIVITY) return;
+ for(c = s->cl; c != nil; c = c->next)
+ if(c->activity != 0){+ c->activity /= MAXACTIVITY;
+ if(c->activity < ε)
+ c->activity = ε;
+ }
+ s->Δclactivity /= MAXACTIVITY;
+}
+
+/* pick a literal. normally we pick the variable with highest activity from the heap. sometimes we goof and pick a random one. */
+static void
+decision(SATSolve *s)
+{+ SATVar *v;
+
+ s->decbd[++s->lvl] = s->ntrail;
+ if((uint)s->randfn(s->randaux) < s->goofprob){+ v = s->heap[satnrand(s, s->nheap)];
+ if(v->lvl < 0)
+ goto gotv;
+ }
+ do
+ v = satheaptake(s);
+ while(v->lvl >= 0);
+gotv:
+ sataddtrail(s, 2 * (v - s->var) + (v->flags & VARPHASE));
+}
+
+/* go through the watchlist of a literal that just turned out false. */
+/* full == 1 records the first conflict and goes on rather than aborting immediately */
+static SATClause *
+forcing(SATSolve *s, int l, int full)
+{+ SATClause **cp, *rc, *c, *xp;
+ int v0;
+ int x, j;
+
+ cp = &s->lit[l].watch;
+ rc = nil;
+ if(verboseforcing) print("forcing literal %d\n", signf(l));+ while(c = *cp, c != nil){+ if(l == c->l[0]){+ /* this swap implies that the reason r for a literal l always has r->l[0]==l */
+ x = c->l[1], c->l[1] = c->l[0], c->l[0] = x;
+ xp = c->watch[1], c->watch[1] = c->watch[0], c->watch[0] = xp;
+ }
+ assert(c->l[1] == l);
+ v0 = s->lit[c->l[0]].val;
+ if(v0 > 0) /* the clause is true anyway */
+ goto next;
+ for(j = 2; j < c->n; j++)
+ if(s->lit[c->l[j]].val != 0){+ /* found another literal to watch for this clause */
+ if(verboseforcing) print("moving clause %+Γ onto watchlist %d\n", c, signf(c->l[j]));+ *cp = c->watch[1];
+ x = c->l[j], c->l[j] = c->l[1], c->l[1] = x;
+ c->watch[1] = s->lit[x].watch;
+ s->lit[x].watch = c;
+ goto cont;
+ }
+ if(v0 == 0){+ /* conflict */
+ if(!full) return c;
+ if(rc == nil) rc = c;
+ goto next;
+ }
+ if(verboseforcing) print("inferring %d using clause %+Γ\n", signf(c->l[0]), c);+ sataddtrail(s, c->l[0]);
+ s->var[VAR(c->l[0])].reason = c;
+ next:
+ cp = &c->watch[1];
+ cont: ;
+ }
+ return rc;
+}
+
+/* forcing() for binary implications */
+static uvlong
+binforcing(SATSolve *s, int l, int full)
+{+ SATLit *lp;
+ int i, m;
+ uvlong rc;
+
+ lp = &s->lit[l];
+ rc = 0;
+ if(verboseforcing && lp->nbimp > 0) print("forcing literal %d (binary)\n", signf(l));+ for(i = 0; i < lp->nbimp; i++){+ m = lp->bimp[i];
+ switch(s->lit[m].val){+ case -1:
+ if(verboseforcing) print("inferring %d using binary clause (%d) ∨ %d\n", signf(m), -signf(l), signf(m));+ sataddtrail(s, m);
+ s->var[VAR(m)].binreason = NOT(l);
+ s->var[VAR(m)].isbinreason = 1;
+ break;
+ case 0:
+ if(verboseforcing) print("conflict (%d) ∨ (%d)\n", -signf(l), signf(m));+ if(rc == 0) rc = (uvlong)NOT(l) << 32 | (uint)m;
+ if(!full) return rc;
+ break;
+ }
+ }
+ return rc;
+}
+
+/* check if we can discard the previously learned clause because the current one subsumes it */
+static int
+checkdiscard(SATSolve *s)
+{+ SATClause *c;
+ SATVar *v;
+ int q, j;
+
+ if(s->lastp[1] == &s->learncl) return 0;
+ c = (SATClause*) ((uchar*) s->lastp[1] - (uchar*) &((SATClause*)0)->next);
+ if(s->lit[c->l[0]].val >= 0) return 0; /* clause is a reason, hands off */
+ q = s->ncflcl;
+ for(j = c->n - 1; q > 0 && j >= q; j--){+ v = &s->var[VAR(c->l[j])];
+ /* check if literal is in the current clause */
+ if(c->l[j] == s->cflcl[0] || (uint)v->lvl <= s->cfllvl && v->stamp == s->stamp)
+ q--;
+ }
+ return q == 0;
+}
+
+/* add the clause we just learned to our collection */
+static SATClause *
+learn(SATSolve *s, int notriv)
+{+ SATClause *r;
+ int i, l, triv;
+
+ /* clauses that are too complicated are not worth it. learn the trivial clause (all decisions negated) instead */
+ if(triv = !notriv && s->ncflcl > s->lvl + s->trivlim){+ assert(s->lvl + 1 <= s->cflclalloc);
+ for(i = 1; i <= s->lvl; i++)
+ s->cflcl[i] = NOT(s->trail[s->decbd[s->lvl + 1 - i]]);
+ s->ncflcl = s->lvl + 1;
+ }
+ if(s->ncflcl == 1) /* unit clauses are handled by putting them on the trail in conflict() */
+ return nil;
+ if(!triv && checkdiscard(s))
+ r = satreplclause(s, s->ncflcl);
+ else
+ r = satnewclause(s, s->ncflcl, 1);
+ r->n = s->ncflcl;
+ memcpy(r->l, s->cflcl, s->ncflcl * sizeof(int));
+ for(i = 0; i < 2; i++){+ l = r->l[i];
+ r->watch[i] = s->lit[l].watch;
+ s->lit[l].watch = r;
+ }
+ return r;
+}
+
+/* recursive procedure to determine if a literal is redundant.
+ * to avoid repeated work, each known redundant literal is stamped with stamp+1
+ * and each known nonredundant literal is stamped with stamp+2.
+ */
+static int
+redundant(SATSolve *s, int l)
+{+ SATVar *v, *w;
+ SATClause *c;
+ int i, r;
+
+ v = &s->var[VAR(l)];
+ if(v->isbinreason){+ /* stupid special case code */
+ r = v->binreason;
+ w = &s->var[VAR(r)];
+ if(w->lvl != 0){+ if(w->stamp == s->stamp + 2)
+ return 0;
+ if(w->stamp < s->stamp && (s->lvlstamp[w->lvl] < s->stamp || !redundant(s, r))){+ w->stamp = s->stamp + 2;
+ return 0;
+ }
+ }
+ v->stamp = s->stamp + 1;
+ return 1;
+ }
+ if(v->reason == nil) return 0; /* decision literals are never redundant */
+ c = v->reason;
+ for(i = 0; i < c->n; i++){+ if(c->l[i] == NOT(l)) continue;
+ w = &s->var[VAR(c->l[i])];
+ if(w->lvl == 0)
+ continue; /* literals at level 0 are redundant */
+ if(w->stamp == s->stamp + 2)
+ return 0;
+ /* if the literal is not in the clause or known redundant, check if it is redundant */
+ /* we can skip the check if the level is not stamped: */
+ /* if there are no literals at the same level in the clause, it must be nonredundant */
+ if(w->stamp < s->stamp && (s->lvlstamp[w->lvl] < s->stamp || !redundant(s, c->l[i]))){+ w->stamp = s->stamp + 2;
+ return 0;
+ }
+ }
+ v->stamp = s->stamp + 1;
+ return 1;
+}
+
+/* "blitting" a literal means to either add it to the conflict clause
+ * (if v->lvl < s->lvl) or to increment the counter of literals to be
+ * resolved, plus some bookkeeping. */
+static void
+blit(SATSolve *s, int l)
+{+ SATVar *v;
+ int p;
+
+ v = &s->var[VAR(l)];
+ if(v->stamp == s->stamp) return;
+ v->stamp = s->stamp;
+ p = v->lvl;
+ if(p == 0) return;
+ if(verboseconflict) print("stamp %d %s (ctr=%d)\n", signf(l), p == s->lvl ? "and increment" : "and collect", s->cflctr);+ varbump(s, v);
+ if(p == s->lvl){+ s->cflctr++;
+ return;
+ }
+ if(s->ncflcl >= s->cflclalloc){+ s->cflcl = satrealloc(s, s->cflcl, (s->cflclalloc + CFLCLALLOC) * sizeof(int));
+ s->cflclalloc += CFLCLALLOC;
+ }
+ s->cflcl[s->ncflcl++] = l;
+ if(p > s->cfllvl) s->cfllvl = p;
+ /* lvlstamp[p] == stamp if there is exactly one literal and ==stamp+1 if there is more than one literal on level p */
+ if(s->lvlstamp[p] <= s->stamp)
+ s->lvlstamp[p] = s->stamp + (s->lvlstamp[p] == s->stamp);
+}
+
+/* to resolve a conflict, we start with the conflict clause and use
+ * resolution (a ∨ b and ¬a ∨ c imply b ∨ c) with the reasons for the
+ * literals to remove all but one literal at the current level. this
+ * gives a new "learned" clause with all literals false and we jump back
+ * to the second-highest level in it. at this point, the clause implies
+ * the one remaining literal and we can continue.
+ * to do this quickly, rather than explicitly apply resolution, we keep a
+ * counter of literals at the highest level (unresolved literals) and an
+ * array with all other literals (which will become the learned clause). */
+static void
+conflict(SATSolve *s, SATClause *c, uvlong bin, int full)
+{+ int i, j, l, p, *nl, found;
+ SATVar *v;
+ SATClause *r;
+
+ if(verboseconflict) satprintstate(s);
+ /* choose a new unique stamp value */
+ if(s->stamp >= (uint)-3)
+ stampoverflow(s);
+ s->stamp += 3;
+ s->ncflcl = 1;
+ s->cflctr = 0;
+ s->cfllvl = 0;
+ /* we start by blitting each literal in the conflict clause */
+ if(c != nil){+ clausebump(s, c);
+ for(i = 0; i < c->n; i++)
+ blit(s, c->l[i]);
+ /* if there is only one literal l at the current level, we should have inferred ¬l at a lower level (bug). */
+ if(s->cflctr <= 1){+ satprintstate(s);
+ print("conflict clause %+Γ\n", c);+ assert(s->cflctr > 1);
+ }
+ }else{+ blit(s, bin);
+ blit(s, bin>>32);
+ if(s->cflctr <= 1){+ satprintstate(s);
+ print("binary conflict clause %d ∨ %d\n", (int)(bin>>32), (int)bin);+ assert(s->cflctr > 1);
+ }
+ }
+ /* now we go backwards through the trail, decrementing the unresolved literal counter at each stamped literal */
+ /* and blitting the literals in their reason */
+ for(i = s->ntrail; --i >= 0; ){+ v = &s->var[VAR(s->trail[i])];
+ if(v->stamp != s->stamp) continue;
+ if(verboseconflict) print("trail literal %d\n", signf(s->trail[i]));+ if(--s->cflctr == 0) break;
+ if(v->isbinreason)
+ blit(s, v->binreason);
+ else if((r = v->reason) != nil){+ clausebump(s, r);
+ for(j = 0; j < r->n; j++)
+ blit(s, r->l[j]);
+ }
+ }
+ /* i should point to the one remaining literal at the current level */
+ assert(i >= 0);
+ nl = s->cflcl;
+ nl[0] = NOT(s->trail[i]);
+ found = 0;
+ /* delete redundant literals. note we must watch a literal at cfllvl, so put it in l[1]. */
+ for(i = 1, j = 1; i < s->ncflcl; i++){+ l = nl[i];
+ p = s->var[VAR(nl[i])].lvl;
+ /* lvlstamp[p] != s->stamp + 1 => only one literal at level p => must be nonredundant */
+ if(s->lvlstamp[p] != s->stamp + 1 || !redundant(s, l))
+ if(found || p < s->cfllvl)
+ nl[j++] = nl[i];
+ else{+ /* watch this literal */
+ l = nl[i], nl[j++] = nl[1], nl[1] = l;
+ found = 1;
+ }
+ }
+ s->ncflcl = j;
+ if(!full){+ /* normal mode: jump back and add to trail right away */
+ satbackjump(s, s->cfllvl);
+ sataddtrail(s, nl[0]);
+ }else{+ /* purging: record minimum cfllvl and literals at that level */
+ if(s->cfllvl < s->fullrlvl){+ s->fullrlvl = s->cfllvl;
+ s->nfullrlits = 0;
+ }
+ s->fullrlits[s->nfullrlits++] = nl[0];
+ }
+ r = learn(s, full);
+ if(!full && r != nil)
+ s->var[VAR(nl[0])].reason = r;
+ if(verboseconflict)
+ if(r != nil)
+ print("learned %+Γ\n", r);+ else
+ print("learned %d\n", signf(nl[0]));+ s->Δactivity *= s->varρ;
+ s->Δclactivity *= s->clauseρ;
+ s->conflicts++;
+}
+
+/* to purge, we need a fullrun that assigns values to all variables.
+ * during this we record the first conflict at each level, to be resolved
+ * later. otherwise this is just a copy of the main loop which never
+ * purges or flushes. */
+static int
+fullrun(SATSolve *s)
+{+ int l;
+ uvlong b;
+ SATClause *c;
+
+ while(s->ntrail < s->nvar){+ decision(s);
+ re:
+ while(s->binptr < s->ntrail){+ l = s->trail[s->binptr++];
+ b = binforcing(s, l, 1);
+ if(b != 0){+ if(s->lvl == 0){+ s->unsat = 1;
+ return -1;
+ }
+ if(s->nfullrcfl == 0 || s->lvl > CFLLVL(s->fullrcfl[s->nfullrcfl-1])){+ s->fullrcfl = satrealloc(s, s->fullrcfl, sizeof(SATConflict) * (s->nfullrcfl + 1));
+ s->fullrcfl[s->nfullrcfl].lvl = 1<<31 | s->lvl;
+ s->fullrcfl[s->nfullrcfl++].b = b;
+ }
+ }
+ sanity(s);
+ }
+ while(s->forptr < s->ntrail){+ l = s->trail[s->forptr++];
+ c = forcing(s, NOT(l), 1);
+ if(c != nil){+ if(s->lvl == 0){+ s->unsat = 1;
+ return -1;
+ }
+ if(s->nfullrcfl == 0 || s->lvl > CFLLVL(s->fullrcfl[s->nfullrcfl-1])){+ s->fullrcfl = satrealloc(s, s->fullrcfl, sizeof(SATConflict) * (s->nfullrcfl + 1));
+ s->fullrcfl[s->nfullrcfl].lvl = s->lvl;
+ s->fullrcfl[s->nfullrcfl++].c = c;
+ }
+ }
+ }
+ if(s->binptr < s->ntrail) goto re;
+ }
+ return 0;
+}
+
+/* assign range scores to all clauses.
+ * p == number of levels that have positive literals in the clause.
+ * r == number of levels that have literals in the clause.
+ * range == min(floor(16 * (p + α (r - p))), 255) with magic constant α. */
+static void
+ranges(SATSolve *s)
+{+ SATClause *c;
+ int p, r, k, l, v;
+ uint ci;
+
+ ci = 2;
+ memset(s->lvlstamp, 0, sizeof(int) * s->nvar);
+ memset(s->rangecnt, 0, sizeof(s->rangecnt));
+ for(c = s->learncl; c != nil; c = c->next, ci += 2){+ if(!s->var[VAR(c->l[0])].binreason && s->var[VAR(c->l[0])].reason == c){+ c->range = 0;
+ s->rangecnt[0]++;
+ continue;
+ }
+ p = 0;
+ r = 0;
+ for(k = 0; k < c->n; k++){+ l = c->l[k];
+ v = s->var[VAR(l)].lvl;
+ if(v == 0){+ if(s->lit[l].val == 1){+ c->range = 256;
+ goto next;
+ }
+ }else{+ if(s->lvlstamp[v] < ci){+ s->lvlstamp[v] = ci;
+ r++;
+ }
+ if(s->lvlstamp[v] == ci && s->lit[l].val == 1){+ s->lvlstamp[v] = ci + 1;
+ p++;
+ }
+ }
+ }
+ r = 16 * (p + s->purgeα * (r - p));
+ if(r > 255) r = 255;
+ c->range = r;
+ s->rangecnt[r]++;
+ next: ;
+ }
+}
+
+/* resolve conflicts found during fullrun() */
+static void
+fullrconflicts(SATSolve *s)
+{+ SATConflict *cfl;
+ int i;
+
+ s->fullrlvl = s->lvl;
+ s->nfullrlits = 0;
+ for(cfl = &s->fullrcfl[s->nfullrcfl - 1]; cfl >= s->fullrcfl; cfl--){+ satbackjump(s, CFLLVL(*cfl));
+ if(cfl->lvl < 0)
+ conflict(s, nil, cfl->b, 1);
+ else
+ conflict(s, cfl->c, 0, 1);
+ }
+ satbackjump(s, 0);
+ if(s->fullrlvl == 0)
+ for(i = 0; i < s->nfullrlits; i++)
+ sataddtrail(s, s->fullrlits[i]);
+ free(s->fullrcfl);
+ s->fullrcfl = nil;
+}
+
+/* note that nil > *, this simplifies the algorithm by having nil "bubble" to the top */
+static int
+actgt(SATClause *a, SATClause *b)
+{+ if(b == nil) return 0;
+ if(a == nil) return 1;
+ return a->activity > b->activity || a->activity == b->activity && a > b;
+}
+
+/* select n clauses to keep
+ * first we find the upper limit j on the range score
+ * to get the exact number, we move htot clauses from j to j+1
+ * to this end, we put them in a max-heap of size htot, sorted by activity,
+ * continually replacing the largest element if we find a less active clause.
+ * the heap starts out filled with nil and the nil are replaced during the first
+ * htot iterations. */
+#define LEFT(i) (2*(i)+1)
+#define RIGHT(i) (2*(i)+2)
+static int
+judgement(SATSolve *s, int n)
+{+ int cnt, i, j, htot, m;
+ SATClause *c, **h, *z;
+
+ cnt = 0;
+ for(j = 0; j < 256; j++){+ cnt += s->rangecnt[j];
+ if(cnt >= n) break;
+ }
+ if(j == 256) return j;
+ if(cnt > n){+ htot = cnt - n;
+ h = satrealloc(s, nil, sizeof(SATClause *) * htot);
+ memset(h, 0, sizeof(SATClause *) * htot);
+ for(c = s->learncl; c != nil; c = c->next){+ if(c->range != j || actgt(c, h[0])) continue;
+ h[0] = c;
+ m = 0;
+ for(;;){+ i = m;
+ if(LEFT(i) < htot && actgt(h[LEFT(i)], h[m])) m = LEFT(i);
+ if(RIGHT(i) < htot && actgt(h[RIGHT(i)], h[m])) m = RIGHT(i);
+ if(i == m) break;
+ z = h[i], h[i] = h[m], h[m] = z;
+ }
+ }
+ for(i = 0; i < htot; i++)
+ if(h[i] != nil)
+ h[i]->range = j + 1;
+ free(h);
+ }
+ return j;
+}
+
+/* during purging we remove permanently false literals from learned clauses.
+ * returns 1 if the clause can be deleted entirely. */
+static int
+cleanupclause(SATSolve *s, SATClause *c)
+{+ int i, k;
+
+ for(i = 0; i < c->n; i++)
+ if(s->lit[c->l[i]].val == 0)
+ break;
+ if(i == c->n) return 0;
+ for(k = i; i < c->n; i++)
+ if(s->lit[c->l[i]].val != 0)
+ c->l[k++] = c->l[i];
+ c->n = k;
+ if(k > 1) return 0;
+ if(k == 0)
+ s->unsat = 1;
+ else if(s->lit[c->l[0]].val < 0)
+ sataddtrail(s, c->l[0]);
+ return 1;
+}
+
+/* delete clauses by overwriting them. don't delete empty blocks; we're going to fill them up soon enough again. */
+static void
+execution(SATSolve *s, int j)
+{+ SATClause *c, *n, **cp, *p;
+ SATBlock *b;
+ SATVar *v0;
+ int f, sz;
+
+ b = s->bl[1].next;
+ p = (SATClause*) b->data;
+ s->ncl = s->ncl0;
+ cp = &s->learncl;
+ for(c = p; c != nil; c = n){+ n = c->next;
+ if(c->range > j || cleanupclause(s, c))
+ continue;
+ sz = sizeof(SATClause) + (c->n - 1) * sizeof(int);
+ f = (uchar*)b + SATBLOCKSZ - (uchar*)p;
+ if(f < sz){+ memset(p, 0, f);
+ b = b->next;
+ assert(b != &s->bl[1]);
+ p = (SATClause *) b->data;
+ }
+ b->last = p;
+ /* update reason field of the first variable (if applicable) */
+ v0 = &s->var[VAR(c->l[0])];
+ if(!v0->isbinreason && v0->reason == c)
+ v0->reason = p;
+ memmove(p, c, sz);
+ *cp = p;
+ cp = &p->next;
+ p = (void*)((uintptr)p + sz + CLAUSEALIGN - 1 & -CLAUSEALIGN);
+ b->end = p;
+ s->ncl++;
+ }
+ *cp = nil;
+ *s->lastp[0] = s->learncl;
+ s->lastp[1] = cp;
+ s->lastbl = b;
+ f = (uchar*)b + SATBLOCKSZ - (uchar*)p;
+ memset(p, 0, f);
+ for(b = b->next; b != &s->bl[1]; b = b->next){+ b->last = nil;
+ b->end = b->data;
+ }
+}
+
+static void
+thepurge(SATSolve *s)
+{+ int nkeep, i, j;
+ SATVar *v;
+
+ s->purgeival += s->purgeδ;
+ s->nextpurge = s->conflicts + s->purgeival;
+ nkeep = (s->ncl - s->ncl0) / 2;
+ for(i = 0; i < s->ntrail; i++){+ v = &s->var[VAR(s->trail[i])];
+ if(!v->isbinreason && v->reason != nil)
+ nkeep++;
+ }
+ if(nkeep <= 0) return; /* shouldn't happen */
+ s->nfullrcfl = 0;
+ if(fullrun(s) < 0){ /* accidentally determined UNSAT during fullrun() */+ free(s->fullrcfl);
+ s->fullrcfl = nil;
+ return;
+ }
+ ranges(s);
+ fullrconflicts(s);
+ j = judgement(s, nkeep);
+ execution(s, j);
+ rewatch(s);
+}
+
+/* to avoid getting stuck, flushing backs up the trail to remove low activity variables.
+ * don't worry about throwing out high activity ones, they'll get readded quickly. */
+static void
+theflush(SATSolve *s)
+{+ double actk;
+ int dd, l;
+
+ /* "reluctant doubling" wizardry to determine when to flush */
+ if((s->flushu & -s->flushu) == s->flushv){+ s->flushu++;
+ s->flushv = 1;
+ s->flushθ = s->flushψ;
+ }else{+ s->flushv *= 2;
+ s->flushθ += s->flushθ >> 4;
+ }
+ s->nextflush = s->conflicts + s->flushv;
+ if(s->agility > s->flushθ) return; /* don't flush when we're too busy */
+ /* clean up the heap so that a free variable is at the top */
+ while(s->nheap > 0 && s->heap[0]->lvl >= 0)
+ satheaptake(s);
+ if(s->nheap == 0) return; /* shouldn't happen */
+ actk = s->heap[0]->activity;
+ for(dd = 0; dd < s->lvl; dd++){+ l = s->trail[s->decbd[dd+1]];
+ if(s->var[VAR(l)].activity < actk)
+ break;
+ }
+ satbackjump(s, dd);
+}
+
+int
+satsolve(SATSolve *s)
+{+ int l;
+ SATClause *c;
+ uvlong b;
+
+ if(s == nil) return 1;
+ if(s->scratched) return -1;
+ if(s->nvar == 0) return 1;
+ solvinit(s);
+
+ while(!s->unsat){+ re:
+ while(s->binptr < s->ntrail){+ l = s->trail[s->binptr++];
+ b = binforcing(s, l, 0);
+ sanity(s);
+ if(b != 0){+ if(s->lvl == 0) goto unsat;
+ conflict(s, nil, b, 0);
+ sanity(s);
+ }
+ }
+ while(s->forptr < s->ntrail){+ l = s->trail[s->forptr++];
+ c = forcing(s, NOT(l), 0);
+ sanity(s);
+ if(c != nil){+ if(s->lvl == 0) goto unsat;
+ conflict(s, c, 0, 0);
+ sanity(s);
+ }
+ }
+ if(s->binptr < s->ntrail) goto re;
+ if(s->ntrail == s->nvar) goto out;
+ if(s->conflicts >= s->nextpurge)
+ thepurge(s);
+ else if(s->conflicts >= s->nextflush)
+ theflush(s);
+ else
+ decision(s);
+ }
+unsat:
+ s->unsat = 1;
+out:
+ satcleanup(s, 0);
+ return !s->unsat;
+}
+
+void
+satreset(SATSolve *s)
+{+ int i;
+
+ if(s == nil || s->decbd == nil) return;
+ satbackjump(s, -1);
+ s->lvl = 0;
+ for(i = 0; i < s->nvar; i++){+ s->var[i].activity = 0;
+ s->var[i].flags |= VARPHASE;
+ }
+ satcleanup(s, 1);
+ s->Δactivity = 1;
+ s->Δclactivity = 1;
+}
--- a/sys/src/mkfile
+++ b/sys/src/mkfile
@@ -30,6 +30,7 @@
libndb\
libplumb\
libregexp\
+ libsat\
libscribble\
libsec\
libstdio\
--
⑨