\datethis
\def\adj{\mathrel{\!\mathrel-\mkern-8mu\mathrel-\mkern-8mu\mathrel-\!}}
\def\dadj{\mathrel{\!\mathrel-\mkern-8mu\mathrel-\mkern-12mu\to\!}}
@s mod and
\let\Xmod=\bmod % this is CWEB magic for using "mod" instead of "%"
@i gb_types.w
@*Intro. This program finds all the Hamiltonian cycles of a given graph.
It uses an interesting algorithm that chooses the edges of subpaths
without knowing where those edges will appear in the final cycle
until they are all linked together. (That important idea was introduced in
Geoffrey Selby's thesis (1970), and extended by Silvano Martello in 1983
to incorporate the {\mc MRV} branching heuristic. I rediscovered Selby's
approach independently in 2001, in a slightly more symmetrical form;
in my variant, all subpaths have equal status, so that we needn't branch only
on the ways of extending one end of a principal subpath.) My first
implementation, {\mc HAMDANCE}, used dancing links;
here I'm using sparse-set data structures instead.
As in other programs such as {\mc SSXCC}, this one reports the running time
in ``mems.'' One mem is counted whenever we read or write a 64-bit word of
memory, but not when we access data that's already in a register.
(The number of mems reported does not include the work that we do when
inputting the graph or printing the results.)
@d o mems++ /* count one mem */
@d oo mems+=2 /* count two mems */
@d ooo mems+=3 /* count three mems */
@d O "%" /* used for percent signs in format strings */
@d mod % /* used for percent signs denoting remainder in \CEE/ */
@d maxn 1000 /* at most this many vertices */
@d infty maxn /* larger than any vertex number */
@c
#include "gb_graph.h" /* use the Stanford GraphBase conventions */
#include "gb_save.h" /* and its routine for inputting graphs */
#include "gb_flip.h" /* and its random number generator */
@h
typedef unsigned long long ullng; /* a convenient abbreviation */
@<type definitions@>;
@<global variables@>@;
Graph *g; /* the given graph */
@<subroutines@>@;
int main(int argc,char *argv[])
{
 register int i,j,k,d,t,u,v,w,mm;
 @<process the command line, inputting the graph@>;
 @<prepare the graph for backtracking@>;
 imems=mems, mems=0;
 @<backtrack through all solutions@>;
done: @<print the results@>;
 exit(0);
}
@ The command line names the graph, which is supplied in a file
such as &#124;"foo.gb"&#124; in Stanford GraphBase format. Other options
may follow this file name, in order to cause printing of some or
all of the solutions, or to provide diagnostic information.
Here's a list of the available options:
\smallskip
\item{$\bullet$}
`\.v$\langle,円$integer$,円\rangle$' enables or disables various kinds of verbose
 output on &#124;stderr&#124;, specified as a sum of binary codes such as &#124;show_choices&#124;;
\item{$\bullet$}
`\.m$\langle,円$integer$,円\rangle$' causes every $m$th solution
to be output (the default is \.{m0}, which merely counts them);
\item{$\bullet$}
`\.s$\langle,円$integer$,円\rangle$' causes the algorithm to randomize
the input graph data (thus providing some variety, although
the solutions are by no means uniformly random);
\item{$\bullet$}
`\.d$\langle,円$integer$,円\rangle$' sets &#124;delta&#124;, which causes periodic
state reports on &#124;stderr&#124; after the algorithm has performed approximately
&#124;delta&#124; mems since the previous report (default 10000000000);
\item{$\bullet$}
`\.t$\langle,円$positive integer$,円\rangle$' causes the program to
stop after this many solutions have been found;
\item{$\bullet$}
`\.T$\langle,円$integer$,円\rangle$' sets &#124;timeout&#124; (which causes abrupt
termination if &#124;mems>timeout&#124; at the beginning of a level);
\item{$\bullet$}
`\.n$\langle,円$integer$,円\rangle$' ignores all but the first &#124;n&#124; vertices
of the input graph (default &#124;maxn&#124;);
@d show_basics 1 /* &#124;vbose&#124; code for basic stats; this is the default */
@d show_choices 2 /* &#124;vbose&#124; code for backtrack logging */
@d show_details 4 /* &#124;vbose&#124; code for further commentary */
@d show_raw_sols 64 /* &#124;vbose&#124; code to show solutions in order of arcs added */
@d show_profile 128 /* &#124;vbose&#124; code to show the search tree profile */
@d show_full_state 256 /* &#124;vbose&#124; code for complete state reports */
@<glob...@>=
int random_seed=0; /* seed for the random words of &#124;gb_rand&#124; */
int randomizing; /* has option `\.s' been specified? */
int vbose=show_basics; /* level of verbosity */
int spacing; /* solution $k$ is output if $k$ is a multiple of &#124;spacing&#124; */
int maxl; /* maximum level actually reached */
ullng count; /* solutions found so far */
ullng imems,rmems,mems; /* mem counts */
ullng delta=10000000000; /* report every &#124;delta&#124; or so mems */
ullng thresh=10000000000; /* report when &#124;mems&#124; exceeds this, if &#124;delta!=0&#124; */
ullng maxcount=0xffffffffffffffff; /* stop after finding this many solutions */
ullng timeout=0x1fffffffffffffff; /* give up after this many mems */
ullng nodes; /* total size of search tree */
ullng profile[maxn]; /* number of nodes at each level of the search tree */
int nn; /* number of vertices in the given graph */
int nmax=maxn; /* cutoff for the largest vertex number */
int mind,maxd; /* smallest and largest degree in the given graph */
@ If an option appears more than once on the command line, the first
appearance takes precedence.
@<process the command line...@>=
for (j=argc-1,k=0;j>1;j--) switch (argv[j][0]) {
case 'v': k&#124;=(sscanf(argv[j]+1,""O"d",&vbose)-1);@+break;
case 'm': k&#124;=(sscanf(argv[j]+1,""O"d",&spacing)-1);@+break;
case 's': k&#124;=(sscanf(argv[j]+1,""O"d",&random_seed)-1),randomizing=1;@+break;
case 'd': k&#124;=(sscanf(argv[j]+1,""O"lld",&delta)-1),thresh=delta;@+break;
case 't': k&#124;=(sscanf(argv[j]+1,""O"lld",&maxcount)-1);@+break;
case 'T': k&#124;=(sscanf(argv[j]+1,""O"lld",&timeout)-1);@+break;
case 'n': k&#124;=(sscanf(argv[j]+1,""O"d",&nmax)-1);@+break;
default: k=1; /* unrecognized command-line option */
}
if (argc<2) k=1; if (k==0) { g=restore_graph(argv[1]); if (!g) { fprintf(stderr,"I couldn't reconstruct graph "O"s!\n",argv[1]); k=1; }@+else { nn=g->n;
 if (nn>nmax) nn=nmax;
 if (nn>maxn) {
 fprintf(stderr,"Sorry, graph "O"s has too many vertices ("O"d>"O"d)!\n",
 argv[1],nn,maxn);
 exit(-2);
 }
 }
}
if (k) {
 fprintf(stderr,
 "Usage: "O"s foo.gb [v<n>] [m<n>] [s<n>] [d<n>] [t<n>] [T<n>] [n<n>]\n",
 argv[0]);
 exit(-1);
}
if (randomizing) gb_init_rand(random_seed);
@ @<print the results@>=
if (vbose&show_profile) @<print the profile@>;
if (vbose&show_basics) {
 fprintf(stderr,"Altogether "O"llu solution"O"s, "O"llu nodes,",@&#124;
 count,count==1?"":"s",nodes);
 fprintf(stderr," "O"llu+"O"llu mems.\n",imems,mems);
}
@ To help detect faulty reasoning,
we provide a routine that we hope is never invoked.
@<sub...@>=
void confusion(char *m) {
 fprintf(stderr,"This can't happen: "O"s!\n",m);
 exit(666);
}
@*Data structures. This program can be regarded as an algorithm that
finds Hamiltonian cycles by starting with
a graph~&#124;g&#124; and removing edges until only a cycle is left. 
We use a sparse-set representation for &#124;g&#124;, because such structures are
an especially attractive way to maintain
the current status of a graph that is continually getting smaller.
The idea is to have two arrays, &#124;nbr&#124; and &#124;adj&#124;, with one row for
each vertex~&#124;v&#124;. If &#124;v&#124; has &#124;d&#124; neighbors in~&#124;g&#124;, they're listed (in any order)
in the first &#124;d&#124; columns of &#124;nbr[v]&#124;. And if &#124;nbr[v][k]=u&#124;, where 0ドル\le k<d,ドル we have &#124;adj[v][u]=k&#124;; in other words, there's an important invariant relation, $$\hbox{&#124;nbr[v][adj[v][u]]=u&#124;.}$$ Neighbors can be deleted by moving them to the right and decreasing~&#124;d&#124;; neighbors can be undeleted by simply increasing~&#124;d&#124;. Furthermore, if &#124;u&#124; is not a neighbor of~&#124;v&#124;, &#124;adj[v][u]&#124; has the impossible value &#124;infty&#124;; thus the &#124;adj&#124; matrix functions also as an adjacency matrix. @<glob...@>=
int nbr[maxn][maxn], adj[maxn][maxn]; /* sparse-set representation of &#124;g&#124; */
int degree[maxn]; /* vertex degree in the input graph (for diagnostics only) */
@ The edges of &#124;g&#124; are considered to be pairs of arcs that run in
opposite directions. (In other words, the edge $u\adj v$ is actually
treated as two arcs, $u\dadj v$ and $v\dadj u$.) When an edge is deleted,
we often need to delete only one of those arcs, because our algorithm
doesn't always depend on both of them.
The algorithm proceeds not only by removing unwanted edges but also by
choosing edges that will {\it not\/} be removed. Those edges appear
in an array &#124;e&#124; of \&{edge} structs, each of which has two fields
&#124;u&#124; and~&#124;v&#124;. If the $k$th chosen edge is $u\adj v,ドル we have &#124;e[k].u=u&#124; and
&#124;e[k].v=v&#124;.
@<type...@>=
typedef struct edge_struct {
 int u,v; /* the vertices joined by this edge */
} edge;
@ @<glob...@>=
edge e[maxn]; /* the edges chosen so far */
int eptr; /* we've currently chosen this many edges */
@ @<sub...@>=
void print_edges(void) {
 register int k;
 for (k=0;k<eptr;k++) printf(""O"s--"O"s\n",name(e[k].u),name(e[k].v)); } @ The chosen edges form one or more subpaths of the final cycle. If $$v_0\adj v_1\adj \cdots\adj v_k$$ is a chosen subpath, where $v_0$ and $v_k$ participate in only one of the edges chosen so far, we say that $v_0$ and $v_k$ are ``outer'' vertices, while $\{v_1,\ldots,v_{k-1}\}$ are ``inner.'' A vertex that's neither outer nor inner is called ``bare.'' Every vertex begins bare, and is eventually clothed. At the end all vertices will have become inner, except for the two vertices of the last-chosen edge; and the chosen edges will be a Hamiltonian cycle. As the algorithm proceeds, two crucial integer values are associated with every vertex~&#124;v&#124;, namely &#124;mate(v)&#124; and &#124;deg(v)&#124;. In the chosen subpath above, we have $&#124;mate&#124;(v_0)=v_k$ and $&#124;mate&#124;(v_k)=v_0$; that rule defines &#124;mate(v)&#124; for all outer vertices~&#124;v&#124;. We also define &#124;mate(v)=-1&#124; if and only if &#124;v&#124; is bare. The value of &#124;mate(v)&#124; is undefined when &#124;v&#124; is an inner vertex, except for the fact that it's nonnegative. If $v$ is an outer vertex or a bare vertex, the value of &#124;deg(v)&#124; is the number of unchosen edges touching~$v$ that haven't yet been ruled out for the final path. (Again, &#124;deg(v)&#124; is undefined if &#124;v&#124; is inner; an inner vertex is essentially invisible to the algorithm.) The current values of &#124;mate(v)&#124; and &#124;deg(v)&#124; are maintained in a \&{vert} struct, so that we can conveniently access both of them at once. @d mate(v) vrt[v].m @d deg(v) vrt[v].d @<type...@>=
typedef struct vert_struct {
 int m,d; /* the &#124;mate&#124; and &#124;deg&#124; of this vertex */
} vert;
@ @<glob...@>=
vert vrt[maxn];
@ @<sub...@>=
void print_vert(int v) {
 register int k;
 printf(""O"s:",name(v));
 for (k=0;;k++) {
 if (k==deg(v)) printf("&#124;");@+else printf(" ");
 if (k==degree[v]) break;
 printf(""O"s",name(nbr[v][k]));
 }
 if (mate(v)<0) printf(" bare\n"); else if (ivis[v]>=visible) printf(" mate "O"s, inner\n",name(mate(v)));
 else printf(" mate "O"s\n",name(mate(v)));
}
@#
void print_verts(void) {
 register int v;
 for (v=0;v<nn;v++) { printf(""O"d,",v); print_vert(v); } } @ As mentioned above, an inner vertex is essentially invisible to the algorithm. An array &#124;vis&#124; lists the visible vertices --- those that are either bare or outer. It's a sparse-set representation, containing a permutation of the vertices, with the invisible ones at the end. The inverse permutation appears in &#124;ivis&#124;, a companion array, so that we have $$\hbox{&#124;vis[k]=v&#124; \qquad $\Leftrightarrow$ \qquad &#124;ivis[v]=k&#124;.}$$ Vertex &#124;v&#124; is visible if and only if &#124;ivis[v]<visible&#124;; thus &#124;v&#124; is inner if and only if &#124;ivis[v]>=visible&#124;.
@d makeinner(v) {@+register int vv,k;
 o,vv=vis[--visible];
 o,k=ivis[v];
 oo,vis[visible]=v,ivis[v]=visible;
 oo,vis[k]=vv,ivis[vv]=k;
 }
@ @<glob...@>=
int vis[maxn],ivis[maxn]; /* sparse-set representation of visibility */
int visible; /* this many vertices are currently visible */
@ @<initialize@>=
for (k=0;k<nn;k++) oo,vis[k]=ivis[k]=k; visible=nn; @ Here's how we remove an existing arc from &#124;u&#124; to &#124;v&#124;. We assume that &#124;u&#124; is visible, and that &#124;v&#124; is currently a neighbor of &#124;u&#124;, namely that &#124;adj[u][v]<deg(u)&#124;. In a production version of this program, the &#124;remove_arc&#124; subroutine can be declared \&{inline}. Thus we don't charge any extra mems for invoking it. The test for $k=d$ in this case saves six mems, at the cost of possibly fouling up branch prediction. So it may not be wise. @<sub...@>=
void remove_arc(int u,int v) {
 register int d,k,w;
 o,d=deg(u)-1;
 oo,k=adj[u][v]; /* we assume that $k\le d$ */
 if (k!=d) {
 oo,w=nbr[u][d];
 o,nbr[u][d]=v;
 o,nbr[u][k]=w;
 o,adj[u][v]=d;
 o,adj[u][w]=k;
 }
 o,deg(u)=d;
} 
@ We maintain a doubly linked list of all the outer vertices in
the current partial solution.
Each entry of this list
is a \&{pair} struct, containing two pointers &#124;llink&#124; and &#124;rlink&#124;.
The list head is a \&{pair} struct called &#124;head&#124;.
Vertex &#124;v&#124; is an outer vertex if and only if the pair &#124;act[v]&#124; is
currently in the list reachable from &#124;head&#124;. Putting it into this
list is called ``activating'' &#124;v&#124;; taking it out is ``deactivating'' it.
@<type...@>=
typedef struct pair_struct {
 int llink,rlink; /* links to left and right in a doubly linked list */
} pair;
@ @d head maxn /* address of the list header in the &#124;act&#124; array */
@d activate(v) {@+register int l=(o,act[head].llink);
 oo,act[l].rlink=act[head].llink=v;
 o,act[v].llink=l,act[v].rlink=head;@+}
@d deactivate(v) {@+register int l=(o,act[v].llink),r=act[v].rlink;
 oo,act[l].rlink=r,act[r].llink=l;
 makeinner(v);@+}
@<glob...@>=
pair act[maxn+1];
@ @<initialize@>=
o,act[head].llink=act[head].rlink=head; /* active list starts empty */
@ @<sub...@>=
void print_actives(void) {
 register int v;
 for (v=act[head].rlink;v!=head;v=act[v].rlink)
 printf(" "O"s",name(v));
 printf("\n");
}
@ One of the command-line options listed above allows randomization of
the input graph. Vertex &#124;k&#124; of our graph corresponds to vertex &#124;perm[k]&#124;
of that one, where &#124;perm[0]&#124;, \dots,~&#124;perm[nn-1]&#124; is a random permutation.
@d name(v) (g->vertices+iperm[v])->name
@<prepare the graph for backtracking@>=
@<initialize@>;
if (randomizing) {
 for (j=0;j<nn;j++) { mems+=4,k=gb_unif_rand(j+1); ooo,perm[j]=perm[k],perm[k]=j; } for (j=0;j<nn;j++) iperm[perm[j]]=j; }@+else@+for (j=0;j<nn;j++) perm[j]=iperm[j]=j; @<set up the &#124;nbr&#124; and &#124;adj&#124; arrays@>;
for (mind=infty,maxd=u=0;u<nn;u++) { if (o,deg(u)<mind) mind=deg(u),curv=u; if (o,deg(u)>maxd) maxd=deg(u);
 if (deg(u)==2) o,trigger[trigptr++]=u;
 for (v=0;v<nn;v++) if (u!=v) { if (adj[u][v]!=infty && adj[v][u]==infty) { fprintf(stderr,"graph "O"s is directed at "O"s and "O"s!\n", argv[1],name(u),name(v)); exit(-5); } } } if (mind<2) { printf("There are no Hamiltonian cycles, because "O"s has degree "O"d!\n", name(curv),mind); exit(0); } fprintf(stderr,"OK, I've got a graph with "O"d vertices, "O"d edges,\n", nn,mm/2); fprintf(stderr," minimum degree "O"d, and maximum degree "O"d.\n",mind,maxd); @ @<glob...@>=
int perm[maxn]; /* vertex mapping between this program and the input graph */
int iperm[maxn]; /* the inverse mapping */
@ @<set up the &#124;nbr&#124; and &#124;adj&#124; arrays@>=
for (i=0;i<nn;i++) for (o,j=0;j<nn;j++) o,adj[i][j]=infty; for (mm=v=0;v<nn;v++) { register int up,vp; register Arc *a; rmems++,vp=perm[v]; oo; /* mems to fetch &#124;nbr[vp]&#124; and &#124;adj[vp]&#124;, needed in the following loop */ for (d=0,o,a=(g->vertices+v)->arcs;a;o,a=a->next) {
 o,u=a->tip-g->vertices;
 if (u>=nmax) continue;
 if (u==v) {
 fprintf(stderr,"graph "O"s has a self loop "O"s--"O"s!\n",
 argv[1],(g->vertices+v)->name,(g->vertices+u)->name);
 exit(-44);
 }
 rmems++,up=perm[u];
 if (adj[vp][up]!=infty) {
 fprintf(stderr,"graph "O"s has a repeated edge "O"s--"O"s!\n",
 argv[1],(g->vertices+v)->name,(g->vertices+u)->name);
 exit(-4);
 }
 oo,nbr[vp][d]=up,adj[vp][up]=d++,mm++;
 }
 o,mate(vp)=-1,deg(vp)=degree[vp]=d;
 if (randomizing) { /* permute the list of neighbors */
 for (j=1;j<d;j++) { mems+=4,k=gb_unif_rand(j+1); oo,u=nbr[vp][j],w=nbr[vp][k]; oo,nbr[vp][j]=w,nbr[vp][k]=u; oo,adj[vp][w]=j,adj[vp][u]=k; } } } if (randomizing) mems+=rmems; /* &#124;rmems&#124; are ignored if &#124;perm&#124; is the identity */ @ Here's a subroutine that painstakingly doublechecks the integrity of the data structures in their current state. It does not, however, attempt to be bulletproof. For example, it assumes that links in the &#124;act&#124; array aren't out of bounds. It doesn't even bother to check that &#124;vis&#124; and &#124;ivis&#124; are inverse permutations. @d sanity_checking 0 /* set this to 1 if you suspect a bug */ @<sub...@>=
void sanity(void) {
 register int u,v,pv,k;
 for (pv=head,v=act[pv].rlink;v!=head;pv=v,v=act[pv].rlink) {
 if (act[v].llink!=pv) fprintf(stderr,"llink of "O"s is bad!\n",name(v));
 if (ivis[v]>=visible) fprintf(stderr,"active "O"s is invisible!\n",name(v));
 u=mate(v);
 if (u<0) fprintf(stderr,"active "O"s has no mate!\n",name(v)); else if (u>=nn) fprintf(stderr,"active "O"s has bad mate!\n",name(v));
 else if (mate(u)!=v)
 fprintf(stderr,"mate(mate("O"s))!="O"s!\n",name(v),name(v));
 else if (adj[v][u]<deg(v)) fprintf(stderr,"there's an arc from "O"s to its mate!\n",name(v)); } if (act[head].llink!=pv) fprintf(stderr,"llink of head is bad!\n"); for (v=0;v<nn;v++) { for (k=0;k<degree[v];k++) if (adj[v][nbr[v][k]]!=k) fprintf(stderr,"Bad nbr["O"s]["O"d]!\n",name(v),k); for (u=0;u<nn;u++) if (adj[v][u]!=infty && nbr[v][adj[v][u]]!=u) fprintf(stderr,"Bad adj["O"s]["O"s]!\n",name(v),name(u)); if (ivis[v]<visible && eptr<nn) { /* &#124;v&#124; is outer or bare */ for (k=0;k<deg(v);k++) { u=nbr[v][k]; if (ivis[u]>=visible)
 fprintf(stderr,"inner "O"s is touched by "O"s!\n",name(u),name(v));
 else if (adj[u][v]>=deg(u))
 fprintf(stderr,"arc "O"s to "O"s is missing!\n",name(u),name(v));
 }
 }
 }
}
@*Nodes, stacks, and the trigger list.
Our backtrack process corresponds to traversing the nodes of a search tree,
and we control that traversal by maintaining status information in
an array of \&{node} structs. On the current level, that info is in
&#124;nd[level]&#124;; and we'll eventually be resuming what we were doing in
&#124;nd[level-1]&#124;, \dots, &#124;nd[0]&#124;.
Thus &#124;nd&#124; is a stack that helps to control this algorithm.
@<type...@>=
typedef struct node_struct {
 int v; /* the active vertex &#124;curv&#124; on which we're branching */
 int m; /* the number of edges chosen so far */
 int i; /* the index &#124;curi&#124; of &#124;curv&#124;'s current neighbor &#124;curu&#124; */
 int d; /* the total number &#124;deg(curv)&#124; of possibilities for &#124;curi&#124; */
 int s; /* the number of visible vertices */
 int t; /* base position in the trigger list (see below) */
 int a; /* base position in the active stack (see below) */
} node;
@ Two more stacks act in parallel with &#124;nd&#124;, but grow at different
rates, namely &#124;savestack&#124; (which records the mates and degrees of vertices)
and &#124;actstack&#124; (which records which vertices were active).
The &#124;savestack&#124; grows by exactly &#124;nn&#124; entries at each level.
@<glob...@>=
int level; /* the depth of branching */
node nd[maxn]; /* nodes between current level and the search tree root */
int trigger[maxn*maxn]; /* vertices whose degree became 2 while bare */
int trigptr; /* the number of vertices in the &#124;trigger&#124; lists */
vert savestack[maxn*maxn]; /* data for the visible vertices at each level */
int saveptr; /* number of entries on &#124;savestack&#124; */
int actstack[maxn*maxn]; /* lists of active vertices at each level */
int actptr; /* number of entries on &#124;actstack&#124; */
@ If a bare vertex &#124;v&#124; has degree~2 in the current graph, every Hamiltonian
cycle must contain the two edges that touch~&#124;v&#124;. We put &#124;v&#124; into a list
called &#124;trigger&#124;, because we want it to force those edges to be chosen
as soon as we have a chance to do so.
@ We call &#124;removex&#124; instead of &#124;remove_arc&#124; when &#124;u&#124; might be bare, because
&#124;removex&#124; will make &#124;u&#124; a trigger at the appropriate time.
(Note: Sometimes &#124;remove_arc&#124; is called when &#124;u&#124; is bare but will soon become
outer. It's more efficient to do that than to use &#124;removex&#124; in all cases.)
@<sub...@>=
void removex(int u,int v) {
 register int d,k,w;
 o,d=deg(u)-1;
 if (mate(u)<0 && d==2) o,trigger[trigptr++]=u; oo,k=adj[u][v]; /* we assume that $k\le d$ */ if (k!=d) { oo,w=nbr[u][d]; o,nbr[u][d]=v; o,nbr[u][k]=w; o,adj[u][v]=d; o,adj[u][w]=k; } o,deg(u)=d; } @*Marching forward. Here we follow the usual pattern of a backtrack process (and I follow my usual practice of &#124;goto&#124;-ing). In this particular case it's a bit tricky to get the whole process started, so I'm deferring that bootstrap calculation until the program for nonroot levels is in place and understood. @<backtrack through all solutions@>=
@<bootstrap the backtrack process@>;
advance:@<clothe everything on the trigger list, or &#124;goto try_again&#124;@>;
if (sanity_checking) sanity();
nodes++,level++;
if (level>maxl) maxl=level;
if (vbose&show_profile) profile[level]++;
if (vbose&show_details) fprintf(stderr,"Entering level "O"d:\n",level);
if (eptr>=nn-1) @<check for solution and &#124;goto backup&#124;@>;
@<do special things if enough mems have accumulated@>;
@<set &#124;curv&#124; to an outer vertex of minimum degree &#124;d&#124;@>;
if (d==0) goto backup;
e[eptr].u=curv; /* no mem charged for the &#124;e&#124; array */
o,trigptr=nd[level-1].t;
@<promote &#124;curv&#124; from outer to inner@>;
if (sanity_checking) sanity();
curi=0;
@<record the current status, for backtracking later@>;
try_move:@<choose the edge from &#124;curv&#124; to &#124;nbr[curv][curi]&#124;@>;
goto advance;
backup:@+if (--level<0) goto done; if (vbose&show_details) fprintf(stderr,"Back to level "O"d\n",level); try_again:@<restore &#124;d&#124; and &#124;curi&#124;, increasing &#124;curi&#124;@>;
if (curi>=d) goto backup;
@<undo the other changes made at the current level@>;
if (level) {
 if (sanity_checking) sanity();
 goto try_move;
}
@<advance at root level@>;
@ @<glob...@>=
int curt,curu,curv,curw; /* current vertices of interest */
int curi; /* index of the neighbor currently chosen */
@ Here's where we force edges to be in the cycle, because some bare
vertex of degree~2 had entered the trigger list.
As we work through that list, the situation might have changed,
because the formerly bare vertex may have become active.
Indeed, giving clothes to one bare vertex might have a ripple effect, causing
other bare vertices to enter the trigger list. The value of &#124;trigptr&#124;
in the following loop might therefore be a moving target.
One case needs to handled with special care: If the two neighbors of &#124;v&#124;
are mates of each other, we are forced to complete a cycle. This is
legitimate only if the cycle includes all vertices.
When this loop has finished, every remaining bare vertex will have
degree 3 or more.
@d vprint() if (vbose&show_choices) fprintf(stderr,
 " "O"s--"O"s\n",name(e[eptr-1].u),name(e[eptr-1].v));
@<clothe everything...@>=
for (o,j=(level?nd[level-1].t:0);j<div class="naked_tr">=0) continue; /* &#124;v&#124; is no longer bare */
 if (deg(v)<2) { if (vbose&show_details) fprintf(stderr,"oops, low degree at "O"s\n",name(v)); goto try_again; } ooo,u=nbr[v][0],w=nbr[v][1]; if ((o,mate(u)==w) && eptr!=nn-2) { if (vbose&show_details) fprintf(stderr,"oops, short cycle "O"s--"O"s--"O"s--..."O"s forced\n", name(u),name(v),name(w),name(u)); goto try_again; } /* now &#124;v&#124; is bare and connected only to &#124;u&#124; and &#124;w&#124;, which aren't mates */ e[eptr].u=u,e[eptr++].v=v;@+vprint(); /* the &#124;e&#124; array is memfree */ e[eptr].u=v,e[eptr++].v=w;@+vprint(); o,mate(v)=v;@+makeinner(v); if (o,mate(u)<0) { if (o,mate(w)<0) @<promote BBB to OIO@>@;
 else @<promote BBO to OII@>;
 }@+else if (o,mate(w)<0) @<promote OBB to IIO@>@;
 else @<promote to OBO to III@>;
}
 
@ A subtle point ought to be explained here: Suppose the input graph
is simply the complete graph $K_3,ドル so that its vertices are $\{0,1,2\}$
and the edges are 0ドル\adj 1\adj 2\adj 0$. Then all vertices go immediately
onto the trigger list. The first promotion, &#124;trigger[0]=0&#124;,
will generate two forced edges 1ドル\adj0$
and 0ドル\adj2$. Then 1 and~2 will no longer be bare; so they won't actually
trigger anything. Moreover, they will be mates; and the edge between
them will have been removed (by &#124;makemates&#124;), leaving them with
degree~0! But that won't be a problem, because the algorithm
never branches after &#124;nn-1&#124; edges have been chosen.
@<promote BBB to OIO@>=
{
 activate(u);
 activate(w);
 remove_arc(u,v);
 remove_arc(w,v);
 makemates(u,w);
}
@ @<sub...@>=
void makemates(int u,int w) {
 if (ooo,adj[w][u]<deg(w)) { /* &#124;u&#124; is a neighbor of &#124;w&#124; */ remove_arc(w,u); remove_arc(u,w); } oo,mate(u)=w,mate(w)=u; } @ @<promote OBB to IIO@>=
{
 activate(w);
 remove_arc(w,v);
 for (oo,k=deg(u)-1;k>=0;k--) {
 o,t=nbr[u][k]; /* the mem for fetching &#124;nbr[u]&#124; was charged above */
 if (t!=v) removex(t,u);
 }
 o,makemates(mate(u),w);
 deactivate(u);
}
@ @<promote BBO to OII@>=
{
 activate(u);
 remove_arc(u,v);
 for (oo,k=deg(w)-1;k>=0;k--) {
 o,t=nbr[w][k]; /* the mem for fetching &#124;nbr[w]&#124; was charged above */
 if (t!=v) removex(t,w);
 }
 o,makemates(mate(w),u);
 deactivate(w);
}
@ In the final case, &#124;u&#124; and &#124;w&#124; are outer vertices. We have cleverly arranged
things so that they are mates if and only if &#124;eptr=nn&#124;, in which case
all edges of a Hamiltonian cycle have already been chosen.
(Consider, for example, the case where the input graph is the cyclic
graph $C_4,ドル with edges 0ドル\adj1\adj2\adj3\adj0$. Then &#124;trigger[0]=0&#124;
will choose edges 1ドル\adj0$ and 0ドル\adj3$. And &#124;trigger[1]=1&#124;
will do nothing, because 1 is no longer bare. Then &#124;trigger[2]=2&#124; will
choose edges 1ドル\adj2$ and 2ドル\adj3,ドル as desired, in spite of the
fact that 1 and~3 are mates.
The subsequent &#124;trigger[3]=3&#124; will again do nothing.)
@<promote to OBO to III@>=
if (eptr!=nn) {
 for (oo,k=deg(u)-1;k>=0;k--) {
 o,t=nbr[u][k]; /* the mem for fetching &#124;nbr[u]&#124; was charged above */
 if (t!=v) removex(t,u);
 }
 for (oo,k=deg(w)-1;k>=0;k--) {
 o,t=nbr[w][k]; /* the mem for fetching &#124;nbr[w]&#124; was charged above */
 if (t!=v) removex(t,w);
 }
 oo,makemates(mate(u),mate(w));
 deactivate(u);
 deactivate(w);
}
@ @<set &#124;curv&#124; to an outer vertex of minimum degree &#124;d&#124;@>=
for (oo,curv=k=act[head].rlink,d=deg(curv);k!=head;o,k=act[k].rlink) {
 if (vbose&show_details) fprintf(stderr," "O"s("O"d)",name(k),deg(k));
 if (o,deg(k)<d) curv=k,d=deg(k); } if (vbose&show_details) fprintf(stderr,", branching on "O"s("O"d)\n",name(curv),d); @ The &#124;d&#124; neighbors of &#124;curv&#124; will remain in &#124;curv&#124;'s list. But &#124;curv&#124; will be removed from {\it their\/} lists. @<promote &#124;curv&#124; from outer to inner@>=
for (o,k=0;k<d;k++) { o,u=nbr[curv][k]; /* the mem for fetching &#124;nbr[curv]&#124; was charged above */ removex(u,curv); } deactivate(curv); @ @<choose the edge from &#124;curv&#124; to &#124;nbr[curv][curi]&#124;@>=
o,curu=nbr[curv][curi];
o,curw=mate(curv);
e[eptr++].v=curu;
if (vbose&show_choices) fprintf(stderr,""O"3d: "O"s--"O"s ("O"d of "O"d)\n",
 level,name(e[eptr-1].u),name(e[eptr-1].v),curi+1,d);
o,curt=mate(curu);
if (curt<0) { /* &#124;curu&#124; is bare */ makemates(curu,curw); activate(curu); }@+else { /* &#124;curu&#124; is outer */ makemates(curt,curw); for (oo,k=deg(curu)-1;k>=0;k--) {
 o,u=nbr[curu][k]; /* the mem for fetching &#124;nbr[curu]&#124; was charged above */
 removex(u,curu);
 }
 deactivate(curu);
}
@*Backtracking. As we explore the search tree, we often want to go back and
investigate the branches not yet taken.
Only one mem is needed to access both &#124;nd[level].v&#124; and &#124;nd[level].m&#124;
simultaneously, because those 32-bit \&{int}s occupy the same 64-bit word.
A similar remark applies to other pairs of fields.
@<record the current status, for backtracking later@>=
{
 o,nd[level].d=d,nd[level].i=curi;
 o,nd[level].s=visible,nd[level].t=trigptr;
 if (d>1) {
 o,nd[level].v=curv,nd[level].m=eptr;
 saveptr=level*nn;
 for (k=0;k<visible;k++) { o,u=vis[k]; oo,savestack[saveptr+u]=vrt[u]; } for (o,u=act[head].rlink;u!=head;o,u=act[u].rlink) actstack[actptr++]=u; } o,nd[level].a=actptr; } @ Here, as suggested by Peter Weigel, we restore only the two most crucial state variables --- because they might tell us that we needn't bother to restore any more. @<restore &#124;d&#124; and &#124;curi&#124;, ...@>=
oo,d=nd[level].d, curi=++nd[level].i;
@ @<undo the other changes made at the current level@>=
for (o,actptr=nd[level].a,v=head,k=(level?o,nd[level-1].a:0);k<actptr;k++) { o,u=actstack[k]; oo,act[v].r v=u; } oo,act[v].r o,visible=nd[level].s,trigptr=nd[level].t; saveptr=level*nn; for (k=0;k<visible;k++) { o,u=vis[k]; oo,vrt[u]=savestack[saveptr+u]; } o,curv=nd[level].v,eptr=nd[level].m; @*Reaping the rewards. Once all vertices have been connected up, no more decisions need to be made. In most such cases, we'll have found a valid Hamiltonian cycle, although its last link usually still needs to be filled in. At this point, exactly two vertices should be active. @ @<check for solution and &#124;goto backup&#124;@>=
{
 if (eptr<nn) { @<if the two &#124;outer&#124; vertices aren't adjacent, &#124;goto backup&#124;@>;
 e[eptr].u=act[head].llink,e[eptr++].v=act[head].rlink;@+vprint();
 }
 count++;
 if (spacing && count mod spacing==0) {
 nd[level].i=0,nd[level].d=1;
 nd[level].m=eptr;
 if (vbose&show_raw_sols) {
 printf("\n"O"llu:\n",count);@+print_state(stdout);
 }@+else @<unscramble and print the current solution@>;
 fflush(stdout);
 }
 if (count>=maxcount) goto done;
 goto backup;
}
@ @<sub...@>=
void print_state(FILE *stream)
{
 register int i,j,l;
 for (j=l=0;l<=level;j++,l++) { while (j<nd[l].m) { fprintf(stream," "O"s--"O"s\n",name(e[j].u),name(e[j].v)); j++; } if (l) { if (j<nn) fprintf(stream," "O"3d: "O"s--"O"s ("O"d of "O"d)\n", l,name(e[j].u),name(e[j].v), nd[l].d==1? 1: nd[l].i+1,nd[l].d); }@+else @<print the state line for the root level@>;
 }
}
@ At this point we've formed a Hamiltonian path, which will be a
Hamiltonian cycle if and only if its two &#124;outer&#124; vertices are
neighbors.
@<if the two &#124;outer&#124; vertices aren't adjacent, &#124;goto backup&#124;@>=
{
 o,u=act[head].llink,v=act[head].rlink;
 if (oo,adj[u][v]==infty) goto backup;
}
@ @d index(v) ((v)-g->vertices)
@<unscramble and print the current solution@>=
{
 register int i,j,k;
 for (k=0;k<nn;k++) v1[k]=-1; for (k=0;k<nn;k++) { i=e[k].u,j=e[k].v; if (v1[i]<0) v1[i]=j; else v2[i]=j; if (v1[j]<0) v1[j]=i; else v2[j]=i; } path[0]=0, path[1]=v1[0]; for (k=2;;k++) { if (v1[path[k-1]]==path[k-2]) path[k]=v2[path[k-1]]; else path[k]=v1[path[k-1]]; if (path[k]==0) break; } for (k=0;k<=nn;k++) printf(""O"s ",name(path[k])); printf("#"O"llu\n",count); } @ @<glob...@>=
int v1[maxn],v2[maxn]; /* the neighbors of a given vertex */
int path[maxn+1]; /* the Hamiltonian cycle, in order */
@*Getting started. Our program is almost complete, but we still need to
figure out how to get the ball rolling by setting things up properly
at backtrack level~0.
There's no problem if the graph has at least one vertex of degree 2,
because the &#124;trigger&#124; list will provide us with at least two active vertices
in such a case. But if all vertices have degree 3 or more, we've got to
have some outer vertices as seeds for the rest of the computation.
In the former (easy) case, we set &#124;curv&#124; to $-1$. In the latter case,
we take a vertex &#124;curv&#124; of minimum degree &#124;d&#124;; we set &#124;nd[0].v=curv&#124;,
and try each neighbor of &#124;curv&#124; in turn. (More precisely, after we've found
all Hamiltonian cycles that contain an edge from &#124;curv&#124; to some other vertex,
&#124;u&#124;, we'll remove that edge permanently from the graph, and repeat
the process until &#124;curv&#124; or some other vertex has only two neighbors left.)
@<bootstrap...@>=
level=0;
d=mind-1;
if (d==1) curv=-1;
else {
 curi=0;
force: oo,curu=nbr[curv][d-curi];
 e[0].u=curv,e[0].v=curu,eptr=1;
if (vbose&show_choices) fprintf(stderr," 0: "O"s--"O"s ("O"d of "O"d)\n",
 name(e[0].u),name(e[0].v),curi+1,d);
 activate(curu);
 activate(curv);
 makemates(curu,curv);
}
@<record the current status...@>;
@ Back at root level, all vertices are again bare.
Since the edge that was previously tried at root level is now no longer
present, one or both of its vertices might now have degree~2;
and in such a case the trigger list will provide a way to finish
the final round.
@<advance at root level@>=
o,curu=mate(curv); /* the previous edge &#124;curv&#124; to &#124;curu&#124; is now gone */
o,act[head].llink=act[head].rlink=head,actptr=0; /* nothing is active */
oo,mate(curu)=mate(curv)=-1,visible=nn; /* everything is bare */
if (deg(curu)==2) trigger[0]=curu,trigptr=1;@+else trigptr=0;
if (deg(curv)==2) trigger[trigptr++]=curv;
if (trigptr==0) goto force;
oo,nd[0].v=-1,nd[0].a=eptr=0;
if (vbose&show_choices) fprintf(stderr," 0: ("O"d of "O"d)\n",
 curi+1,d);
goto advance;
@ If &#124;nd[0].v&#124; is negative, the root level began with its first edges supplied
by the trigger list, so there was no ``chosen'' edge.
@<print the state line for the root level@>=
{
 if (nd[0].v>=0 &#124;&#124; nd[0].d>1)
 fprintf(stream," 0: ("O"d of "O"d)\n",nd[0].i+1,nd[0].d);
 j--; /* compensate for &#124;j++&#124; in the &#124;for&#124; loop */
}
@*Progress reports. It's interesting to watch this algorithm in operation,
and we provide several ways for a user to do that.
@<do special things if enough mems have accumulated@>=
if (delta && (mems>=thresh)) {
 thresh+=delta;
 if (vbose&show_full_state) print_state(stderr);
 else print_progress();
}
if (mems>=timeout) {
 fprintf(stderr,"TIMEOUT!\n");@+goto done;
}
@ During a long run, it's helpful to have some way to measure progress.
The following routine prints a string that indicates roughly where we
are in the search tree. The string consists of character pairs, separated
by blanks, where each character pair represents a branch of the search
tree. When a node has $d$ descendants and we are working on the $k$th,
the two characters respectively represent $k$ and~$d$ in a simple code;
namely, the values 0, 1, \dots, 61 are denoted by
$$\.0,\ \.1,\ \dots,\ \.9,\ \.a,\ \.b,\ \dots,\ \.z,\ \.A,\ \.B,\ \dots,\.Z.$$
All values greater than 61 are shown as `\.*'. Notice that as computation
proceeds, this string will increase lexicographically.
Following that string, a fractional estimate of total progress is computed,
based on the na{\"\i}ve assumption that the search tree has a uniform
branching structure. If the tree consists
of a single node, this estimate is~.5; otherwise, if the first choice
is `$k$ of~$d$', the estimate is $(k-1)/d$ plus 1ドル/d$ times the
recursively evaluated estimate for the $k$th subtree. (This estimate
might obviously be very misleading, in some cases, but at least it
tends to grow monotonically.)
@<sub...@>=
void print_progress(void) {
 register int l,k,d,p;
 register double f,fd;
 fprintf(stderr," after "O"lld mems: "O"lld sols,",mems,count);
 for (f=0.0,fd=1.0,l=0;l<level;l++) { d=nd[l].d; k=(d==1? 1: nd[l].i+1); fd*=d,f+=(k-1)/fd; /* choice &#124;l&#124; is &#124;k&#124; of &#124;d&#124; */ fprintf(stderr," "O"c"O"c", k<10? '0'+k: k<36? 'a'+k-10: k<62? 'A'+k-36: '*', d<10? '0'+d: d<36? 'a'+d-10: d<62? 'A'+d-36: '*'); } fprintf(stderr," "O".5f\n",f+0.5/fd); } @ @<print the profile@>=
{
 fprintf(stderr,"Profile:\n");
 for (level=1;level<=maxl;level++)
 fprintf(stderr,""O"3d: "O"lld\n",
 level,profile[level]);
}
@*Index.
</div><div class="naked_ctrl">
<form action="/index.cgi/speech" method="get" name="gate">
<p><a href="http://altstyle.alfasado.net">AltStyle</a> によって変換されたページ <a href="https://cs.stanford.edu/~knuth/programs/ssham.w">(-&gt;オリジナル)</a>
/ <label>アドレス: <input type="text" name="naked_post_url" value="https://cs.stanford.edu/~knuth/programs/ssham.w" size="22" /></label> <label>モード: <select name="naked_post_mode">
<option value="default">デフォルト</option>
<option value="speech" selected="selected">音声ブラウザ</option>
<option value="ruby">ルビ付き</option>
<option value="contrast">配色反転</option>
<option value="larger-text">文字拡大</option>
<option value="mobile">モバイル</option>
</select>
<input type="submit" value="表示" />
</p>
</form>
</div>