\datethis
@i gb_types.w
@*Introduction. This program inputs a directed graph.
It outputs a not-necessarily-reduced binary decision diagram
for the family of all simple oriented cycles in that graph.
The format of the output is described in another program,
{\mc SIMPATH-REDUCE}. Let me just say here that it is intended
only for computational convenience, not for human readability.
I've tried to make this program simple, whenever I had to
choose between simplicity and efficiency. But I haven't gone
out of my way to be inefficient.
(Notes, 30 November 2015: My original version of this program,
written in August 2008, was hacked from {\mc SIMPATH}. I~don't
think I used it much at that time, if at all, because
I made a change in February 2010 to make it compile without
errors. Today I'm making two fundamental changes:
(i) Each ``frontier'' in {\mc SIMPATH} was required
to be an interval of vertices, according to the vertex numbering.
Now the elements of each frontier are listed explicitly; so
I needn't waste space by including elements that don't really
participate in frontier activities. (ii)~I do {\it not\/}
renumber the vertices.
The main advantage of these two changes is
that I can put a dummy vertex at the end, with arcs to and from
every other vertex; then we get all the simple {\it paths\/} instead
of all the simple {\it cycles}, while the frontiers stay the same size
except for the dummy element. And we can modify this program to get all
the oriented {\it Hamiltonian\/} paths as well.)
@d maxn 90 /* maximum number of vertices; at most 126 */
@d maxm 2000 /* maximum number of arcs */
@d logmemsize 27
@d memsize (1<
#include
#include
#include "gb_graph.h"
#include "gb_save.h"
char mem[memsize]; /* the big workspace */
unsigned long long tail,boundary,head; /* queue pointers */
unsigned int htable[htsize]; /* hash table */
unsigned int htid; /* ``time stamp'' for hash entries */
int htcount; /* number of entries in the hash table */
int wrap=1; /* wraparound counter for hash table clearing */
Vertex *vert[maxn+1];
int f[maxn+2],ff[maxn+2]; /* elements of the current and the next frontier */
int s,ss; /* the sizes of |f| and |ff| */
int curfront[maxn+1],nextfront[maxn+1]; /* inverse frontier map */
int arcto[maxm]; /* destination number of each arc */
int firstarc[maxn+2]; /* where arcs from a vertex start in |arcto| */
char mate[maxn+3]; /* encoded state */
int serial,newserial; /* state numbers */
@@;
@#
main(int argc, char* argv[]) {
register int i,j,jj,jm,k,km,l,ll,m,n,p,t,hash,sign;
register Graph *g;
register Arc *a,*b;
register Vertex *u,*v;
@;
@;
@;
}
@ @=
if (argc!=2) {
fprintf(stderr,"Usage: %s foo.gb\n",argv[0]);
exit(-1);
}
g=restore_graph(argv[1]);
if (!g) {
fprintf(stderr,"I can't input the graph %s (panic code %ld)!\n",
argv[1],panic_code);
exit(-2);
}
n=g->n;
if (n>maxn) {
fprintf(stderr,"Sorry, that graph has %d vertices; ",n);
fprintf(stderr,"I can't handle more than %d!\n",maxn);
exit(-3);
}
if (g->m>maxm) {
fprintf(stderr,"Sorry, that graph has %ld arcs; ",(g->m+1)/2);
fprintf(stderr,"I can't handle more than %d!\n",maxm);
exit(-3);
}
@ The arcs will be either $j\to k$ or $j\gets k$ between vertex number~$j$
and vertex number~$k$, when $j=
@;
for (m=0,k=1;k<=n;k++) {
firstarc[k]=m;
v=vert[k];
printf("%d(%s)\n",k,v->name);
for (a=v->arcs;a;a=a->next) {
u=a->tip;
if (u>v) {
arcto[m++]=u-g->vertices+1;
if (a->len==1) printf(" -> %ld(%s) #%d\n",u-g->vertices+1,u->name,m);
else printf(" -> %ld(%s,%ld) #%d\n",u-g->vertices+1,u->name,a->len,m);
}
}
for (a=v->invarcs;a;a=a->next) {
u=a->tip;
if (u>v) {
arcto[m++]=-(u-g->vertices+1);
if (a->len==1) printf(" <- %ld(%s) #%d\n",u-g->vertices+1,u->name,m);
else printf(" <- %ld(%s,%ld) #%d\n",u-g->vertices+1,u->name,a->len,m);
}
}
}
firstarc[k]=m;
@ To aid in the desired sorting, we first create an inverse-arc
list for each vertex~|v|, namely a list of vertices that point to~|v|.
@d invarcs y.A
@=
for (v=g->vertices;vvertices+n;v++) v->invarcs=NULL;
for (v=g->vertices;vvertices+n;v++) {
vert[v-g->vertices+1]=v;
for (a=v->arcs;a;a=a->next) {
register Arc *b=gb_virgin_arc();
u=a->tip;
b->tip=v;
b->len=a->len;
b->next=u->invarcs;
u->invarcs=b;
}
}
@*The algorithm.
Now comes the fun part. We systematically construct a binary decision
diagram for all simple paths by working top-down, considering the
arcs in |arcto|, one by one.
When we're dealing with arc |i|, we've already constructed a table of
all possible states that might arise when each of the previous arcs has
been chosen-or-not, except for states that obviously cannot be
part of a simple path.
Arc |i| runs from vertex |j| to vertex |k=arcto[i]|,
or from |k=-arcto[i]| to~|j|.
Let $F_i=\{v_1,\ldots,v_s\}$ be the {\it frontier\/} at arc~|i|,
namely the set of vertex numbers |>=j| that appear in arcs~|=
for (t=1;t<=n;t++) mate[t]=t;
@;
for (i=0;i;
}
printf("%x:0,0\n",
serial);
@ Each state for a particular arc gets a distinguishing number, where
its ZDD instructions begin.
Two states are special: 0 means the losing state, when a simple path
is impossible; 1 means the winning state, when a simple path has been
completed. The other states are 2 or more.
Initially |i| will be zero, and the queue is empty. We'll want
|jj| to be the the |j| vertex of arc |i+1|, and |ss| to be the
size of~$F_{i+1}$. Also |serial| is the identifying number for
arc~|i+1|.
@=
jj=1,ss=0;
while (firstarc[jj+1]==0) jj++; /* unnecessary unless vertex 1 is isolated */
tail=head=0;
serial=2;
@ The output format on |stdout| simply shows the identifying numbers of a state
and its two successors, in hexadecimal.
@d trunc(addr) ((addr)&(memsize-1))
@=
if (ss==0) head++; /* put a dummy byte into the queue */
boundary=head,htcount=0,htid=(i+wrap)<0?sign:-sign),s=ss;
for (p=0;p;
while (tail;
@;
printf(",");
@;
printf("\n");
}
@ Here we set |nextfront[t]| to |i+1| whenever $t\in F_{i+1}$.
And we also set |curfront[t]| to |i+1| wheneer $t\in F_i$;
I~use |i+1|, not~|i|, because the |curfront| array is initially zero.
@=
while (jj<=n && firstarc[jj+1]==i+1) jj++;
for (p=ss=0;p=jj) {
nextfront[t]=i+1;
ff[ss++]=t;
}
}
if (j==jj && nextfront[j]!=i+1) nextfront[j]=i+1,ff[ss++]=j;
if (k>=jj && nextfront[k]!=i+1) nextfront[k]=i+1,ff[ss++]=k;
@ This step sets |mate[t]| for all $t\in F_i\cup\{j,k\}$, based on a
queued state, while taking |s| bytes out of the queue.
@=
if (s==0) tail++;
else {
for (p=0;p=
if (sign>0) {
jm=mate[j],km=mate[k];
if (jm==j) jm=-j;
if (jm>=0 || km<=0) printf("0"); /* we mustn't touch a saturated vertex */
else if (jm==-k)
@@;
else {
mate[j]=0,mate[k]=0;
mate[-jm]=km,mate[km]=jm;
printstate(j,jj,i,k);
}
}@+else {
jm=mate[j],km=mate[k];
if (km==k) km=-k;
if (jm<=0 || km>=0) printf("0"); /* we mustn't touch a saturated vertex */
else if (km==-j)
@@;
else {
mate[j]=0,mate[k]=0;
mate[jm]=km,mate[-km]=jm;
printstate(j,jj,i,k);
}
}
@ @=
printstate(j,jj,i,k);
@ See the note below regarding a change that will restrict consideration
to Hamiltonian paths. A similar change is needed here.
@=
{
for (p=0;p=
void printstate(int j,int jj,int i,int k) {
register int h,hh,p,t,tt,hash;
for (p=0;pmemsize) {
fprintf(stderr,"Oops, I'm out of memory: memsize=%d, serial=%d!\n",
memsize,serial);
exit(-69);
}
@;
@;
h=trunc(hh-boundary)/ss;
printf("%x",
newserial+h);
}
}
@ @=
for (p=0,h=trunc(head),hash=0;p=
for (hash=hash&(htsize-1);;hash=(hash+1)&(htsize-1)) {
hh=htable[hash];
if ((hh^htid)>=memsize) @;
hh=trunc(hh);
for (t=hh,h=trunc(head),tt=trunc(t+ss-1);;t=trunc(t+1),h=trunc(h+1)) {
if (mem[t]!=mem[h]) break;
if (t==tt) goto found;
}
}
found:
@ @=
{
if (++htcount>(htsize>>1)) {
fprintf(stderr,"Sorry, the hash table is full (htsize=%d, serial=%d)!\n",
htsize,serial);
exit(-96);
}
hh=trunc(head);
htable[hash]=htid+hh;
head+=ss;
goto found;
}
@*Index.