miles_span.w

模拟器提供了一个简单易用的平台

are competing for minimum-mem prizes. They must be ready to
submit their programs to inspection by impartial judges. A good
algorithm will not need to abuse the spirit of realistic mem-counting.

Mems can be analyzed theoretically as well as empirically.
This means we can attach constants to estimates of running time, instead of
always resorting to $O$~notation.

@*Kruskal's algorithm.
The first algorithm we shall implement and instrument is the simplest:
It considers the edges one by one in order of nondecreasing length,
selecting each edge that does not form a cycle with previously
selected edges.

We know that the edge lengths are less than $2^{12}$, so we can sort them
into order with two passes of a $2^6$-bucket radix sort.
We will arrange to have them appear in the buckets as linked lists
of |Arc| records; the two utility fields of an |Arc| will be called
|from| and |klink|, respectively.

@d from a.V /* an edge goes from vertex |a->from| to vertex |a->tip| */
@d klink b.A /* the next longer edge after |a| will be |a->klink| */

@<Put all the edges into |bucket[0]| through |bucket[63]|@>=
o,n=g->n;
for (l=0;l<64;l++) oo,aucket[l]=bucket[l]=NULL;
for (o,v=g->vertices;v<g->vertices+n;v++)
  for (o,a=v->arcs;a&&(o,a->tip>v);o,a=a->next) {
    o,a->from=v;
    o,l=a->len&0x3f; /* length mod 64 */
    oo,a->klink=aucket[l];
    o,aucket[l]=a;
  }
for (l=63;l>=0;l--)
  for (o,a=aucket[l];a;) {@+register long ll;
    register Arc *aa=a;
    o,a=a->klink;
    o,ll=aa->len>>6; /* length divided by 64 */
    oo,aa->klink=bucket[ll];
    o,bucket[ll]=aa;
  }

@ @<Glob...@>=
Arc *aucket[64], *bucket[64]; /* heads of linked lists of arcs */

@ Kruskal's algorithm now takes the following form.

@<Sub...@>=
unsigned long krusk(g)
  Graph *g;
{@+@<Local variables for |krusk|@>@;@#
  mems=0;
  @<Put all the edges...@>;
  if (verbose) printf("   [%ld mems to sort the edges into buckets]\n",mems);
  @<Put all the vertices into components by themselves@>;
  for (l=0;l<64;l++)
    for (o,a=bucket[l];a;o,a=a->klink) {
      o,u=a->from;
      o,v=a->tip;
      @<If |u| and |v| are already in the same component, |continue|@>;
      if (verbose) report(a->from,a->tip,a->len);
      o,tot_len+=a->len;
      if (--components==1) return tot_len;
      @<Merge the components containing |u| and |v|@>;
    }
  return INFINITY; /* the graph wasn't connected */
}

@ Lest we forget, we'd better declare all the local variables we've
been using.

@<Local variables for |krusk|@>=
register Arc *a; /* current edge of interest */
register long l; /* current bucket of interest */
register Vertex *u,*v,*w; /* current vertices of interest */
unsigned long tot_len=0; /* total length of edges already chosen */
long n; /* the number of vertices */
long components;

@ The remaining things that |krusk| needs to do are easily recognizable
as an application of ``equivalence algorithms'' or ``union/find''
data structures. We will use a simple approach whose average running
time on random graphs was shown to be linear by Knuth and Sch\"onhage
@^Knuth, Donald Ervin@>
@^Sch\"onhage, Arnold@>
in {\sl Theoretical Computer Science\/ \bf 6} (1978), 281--315.

The vertices of each component (that is, of each connected fragment defined by
the edges selected so far) will be linked circularly by |clink| pointers.
Each vertex also has a |comp| field that points to a unique vertex
representing its component. Each component representative also has
a |csize| field that tells how many vertices are in the component.

@d clink z.V /* pointer to another vertex in the same component */
@d comp y.V /* pointer to component representative */
@d csize x.I /* size of the component (maintained only for representatives) */

@<If |u| and |v| are already in the same component, |continue|@>=
if (oo,u->comp==v->comp) continue;

@ We don't need to charge any mems for fetching |g->vertices|, because
|krusk| has already referred to it.
@^discussion of \\{mems}@>

@<Put all the vertices...@>=
for (v=g->vertices;v<g->vertices+n;v++) {
  oo,v->clink=v->comp=v;
  o,v->csize=1;
}
components=n;

@ The operation of merging two components together requires us to
change two |clink| pointers, one |csize| field, and the |comp|
fields in each vertex of the smaller component.

Here we charge two mems for the first |if| test, since |u->csize| and
|v->csize| are being fetched from memory. Then we charge only one mem
when |u->csize| is being updated, since the values being added together
have already been fetched. True, the compiler has to be smart to
realize that it's safe to add the fetched values |u->csize+v->csize|
even though |u| and |v| might have been swapped in the meantime;
but we are assuming that the compiler is extremely clever. (Otherwise we
would have to clutter up our program every time we don't trust the compiler.
After all, programs that count mems are intended primarily to be read.
They aren't intended for production jobs.) % Prim-arily?
@^discussion of \\{mems}@>

@<Merge the components containing |u| and |v|@>=
u=u->comp; /* |u->comp| has already been fetched from memory */
v=v->comp; /* ditto for |v->comp| */
if (oo,u->csize<v->csize) {
  w=u;@+u=v;@+v=w;
} /* now |v|'s component is smaller than |u|'s (or equally small) */
o,u->csize+=v->csize;
o,w=v->clink;
oo,v->clink=u->clink;
o,u->clink=w;
for (;;o,w=w->clink) {
  o,w->comp=u;
  if (w==v) break;
}
  
@* Jarn{\'\i}k and Prim's algorithm.
A second approach to minimum spanning trees is also pretty simple,
except for one technicality: We want to write it in a sufficiently
general manner that different priority queue algorithms can be plugged in.
The basic idea is to choose an arbitrary vertex $v_0$ and connect it to its
nearest neighbor~$v_1$, then to connect that fragment to its nearest
neighbor~$v_2$, and so on. A priority queue holds all vertices that
are adjacent to but not already in the current fragment; the key value
stored with each vertex is its distance to the current fragment.

We want the priority queue data structure to support the four
operations |init_queue(d)|, |enqueue(v,d)|, |requeue(v,d)|, and
|del_min()|, described in the {\sc GB\_\,DIJK} module. Dijkstra's
algorithm for shortest paths, described there, is remarkably similar
to Jarn{\'\i}k and Prim's algorithm for minimum spanning trees; in
fact, Dijkstra discovered the latter algorithm independently, at the
@^Dijkstra, Edsger Wijbe@>
same time as he came up with his procedure for shortest paths.

As in {\sc GB\_\,DIJK}, we define pointers to priority queue subroutines
so that the queueing mechanism can be varied.

@d dist z.I /* this is the key field for vertices in the priority queue */
@d backlink y.V /* this vertex is the stated |dist| away */

@<Glob...@>=
void @[@] (*init_queue)(); /* create an empty priority queue */
void @[@] (*enqueue)(); /* insert a new element in the priority queue */
void @[@] (*requeue)(); /* decrease the key of an element in the queue */
Vertex *(*del_min)(); /* remove an element with smallest key */

@ The vertices in this algorithm are initially ``unseen''; they become
``seen'' when they enter the priority queue, and finally ``known''
when they leave it and enter the current fragment.
We will put a special constant in the |backlink| field
of known vertices. A vertex will be unseen if and only if its
|backlink| is~|NULL|.

@d KNOWN (Vertex*)1 /* special |backlink| to mark known vertices */

@<Sub...@>=
unsigned long jar_pr(g)
  Graph *g;
{@+register Vertex *t; /* vertex that is just becoming known */
  long fragment_size; /* number of vertices in the tree so far */
  unsigned long tot_len=0; /* sum of edge lengths in the tree so far */
  mems=0;
  @<Make |t=g->vertices| the only vertex seen; also make it known@>;
  while (fragment_size<g->n) {
    @<Put all unseen vertices adjacent to |t| into the queue,
      and update the distances of the other vertices adjacent to~|t|@>;
    t=(*del_min)();
    if (t==NULL) return INFINITY; /* the graph is disconnected */
    if (verbose) report(t->backlink,t,t->dist);
    o,tot_len+=t->dist;
    o,t->backlink=KNOWN;
    fragment_size++;
  }
  return tot_len;
}

@ Notice that we don't charge any mems for the subroutine call
to |init_queue|, except for mems counted in the subroutine itself.
What should we charge in general for subroutine linkage when we are
counting mems? The parameters to subroutines generally go into
registers, and registers are ``free''; also, a compiler can often
choose to implement a procedure in line, thereby reducing the
overhead to zero. Hence, the recommended method for charging mems
with respect to subroutines is: Charge nothing if the subroutine
is not recursive; otherwise charge twice the number of things that need
to be saved on a runtime stack. (The return address is one of the
things that needs to be saved.)
@^discussion of \\{mems}@>

@<Make |t=g->vertices| the only vertex seen; also make it known@>=
for (oo,t=g->vertices+g->n-1;t>g->vertices;t--) o,t->backlink=NULL;
o,t->backlink=KNOWN;
fragment_size=1;
(*init_queue)(0L); /* make the priority queue empty */

@ @<Put all unseen vertices adjacent to |t| into the queue,
      and update the distances of the other vertices adjacent to~|t|@>=
{@+register Arc *a; /* an arc leading from |t| */
  for (o,a=t->arcs; a; o,a=a->next) {
    register Vertex *v; /* a vertex adjacent to |t| */
    o,v=a->tip;
    if (o,v->backlink) { /* |v| has already been seen */
      if (v->backlink>KNOWN) {
        if (oo,a->len<v->dist) {
          o,v->backlink=t;
          (*requeue)(v,a->len); /* we found a better way to get there */
        }
      }
    }@+else { /* |v| hasn't been seen before */
      o,v->backlink=t;
      o,(*enqueue)(v,a->len);
    }
  }
}

@*Binary heaps.
To complete the |jar_pr| routine, we need to fill in the four
priority queue functions. Jarn{\'\i}k wrote his original paper before
computers were known; Prim and Dijkstra wrote theirs before efficient priority
queue algorithms were known. Their original algorithms therefore
took $\Theta(n^2)$ steps. 
Kerschenbaum and Van Slyke pointed out in 1972 that binary heaps could
@^Kerschenbaum, A.@>
@^Van Slyke, Richard Maurice@>
do better. A simplified version of binary heaps (invented by Williams
@^Williams, John William Joseph@>
in 1964) is presented here.

A binary heap is an array of $n$ elements, and we need space for it.
Fortunately the space is already there; we can use utility field
|u| in each of the vertex records of the graph. Moreover, if
|heap_elt(i)| points to vertex~|v|, we will arrange things so that
|v->heap_index=i|.

@d heap_elt(i) (gv+i)->u.V /* the |i|th vertex of the heap; |gv=g->vertices| */
@d heap_index v.I
 /* the |v| utility field says where a vertex is in the heap */

@<Glob...@>=
Vertex *gv; /* |g->vertices|, the base of the heap array */
long hsize; /* the number of elements currently in the heap */

@ To initialize the heap, we need only initialize two ``registers'' to
known values, so we don't have to charge any mems at all. (In a production
implementation, this code would appear in-line as part of the
spanning tree algorithm.)
@^discussion of \\{mems}@>

Important Note: This routine refers to the global variable |g|, which is
set in |main| (not in |jar_pr|). Suitable changes need to be made
if these binary heap routines are used in other programs.

@<Priority queue subroutines@>=
void init_heap(d) /* makes the heap empty */
  long d;
{
  gv=g->vertices;
  hsize=0;
}

@ The key invariant property that makes heaps work is
$$\hbox{|heap_elt(k/2)->dist<=heap_elt(k)->dist|, \qquad for |1<k<=hsize|.}$$
(A reader who has not seen heap ordering before should stop at this
point and study the beautiful consequences of this innocuously simple
set of inequalities.) The enqueueing operation turns out to be quite simple:

@<Priority queue subroutines@>=
void enq_heap(v,d)
  Vertex *v; /* vertex that is entering the queue */
  long d; /* its key (aka |dist|) */
{@+register unsigned long k; /* position of a ``hole'' in the heap */
  register unsigned long j; /* the parent of that position */
  register Vertex *u; /* |heap_elt(j)| */
  o,v->dist=d;
  k=++hsize;
  j=k>>1; /* |k/2| */
  while (j>0 && (oo,(u=heap_elt(j))->dist>d)) {
    o,heap_elt(k)=u; /* the hole moves to parent position */
    o,u->heap_index=k;
    k=j;
    j=k>>1;
  }
  o,heap_elt(k)=v;
  o,v->heap_index=k;
}

@ And in fact, the general requeueing operation is almost identical to
enqueueing.  This operation is popularly called ``siftup,'' because
the vertex whose key is being reduced may displace its ancestors
higher in the heap. We could have implemented enqueueing by first
placing the new element at the end of the heap, then requeueing it;
that would have cost at most a couple mems more.

@<Priority queue subroutines@>=
void req_heap(v,d)
  Vertex *v; /* vertex whose key is being reduced */
  long d; /* its new |dist| */
{@+register unsigned long k; /* position of a ``hole'' in the heap */
  register unsigned long j; /* the parent of that position */
  register Vertex *u; /* |heap_elt(j)| */
  o,v->dist=d;
  o,k=v->heap_index; /* now |heap_elt(k)=v| */