PROCEDURE GenerateGraph (coins: Random.T): Graph =

  VAR 
    (* Adjacency matrix: *)
    g: Graph;
    (* Invariant: 
         /\ g[u,v] = TRUE iff there is an edge from vertex u to vertex v
         /\ g[u,u] = FALSE
         /\ g[u,v] = g[v,u]
    *)

      (* Degree table: *)
      deg: ARRAY Vertex OF CARDINAL;     
      (* Invariant: deg[v] is the degree of vertex v in g. *)
      
      (* Induced subgraph tables: *)
      ne: ARRAY Subgraph OF CARDINAL;   (* Number of edges in subgraph s *)
      nv: ARRAY Subgraph OF CARDINAL;   (* Number of vertices in subgraph s *)
      maxe: ARRAY Subgraph OF CARDINAL; (* Max num of edges allowed in subgraph s *) 
      (* 
        "Subgraph s" means the subgraph of g induced by the set of vertices
        whose indices are the bit positions of the "1" bits in the binary 
        representation of integer s (with bit 0 in the units place).  Thus,
        for example, subgraph 0 is the empty graph, and subgraph TN-1 is
        g itself.
        
        Note that nv[s] is simply the number of "1" bits in the binary representation
        of s.
        
        Invariant:
          s IN [0..TN-1] ==> ne[s] < D * (nv[s] DIV 2)
      *)
        
      (* Candidate edge table: *)
      ncands: CARDINAL;
      candK: ARRAY Vertex OF ARRAY Vertex OF INTEGER; 
      candE: ARRAY [0..M-1] OF Edge;
      (* Invariant: 
          /\ candK[u,v] = candK[v,u]
          /\ candK[u,u] = -1
          /\ (0 < k < ncands /\ candE[k] = (u,v)) ==> candK[u,v] = k
          /\ candK[u,v] # -1 ==> 
              /\ candE[candK[u,v]] = (u,v)
              /\ g[u,v] = FALSE
              /\ candK[u,v] < ncands
              /\ degree[u] < D
              /\ degree[v] < D
      *)

  PROCEDURE InitializeGraph() =
    (* Initializes g  deg, ncands, candK, candE, ne, nv: *)
    BEGIN
      FOR u := 0 TO N-1 DO
        g[u,u] := FALSE;
        candK[u,u] := -1;
        deg[u] := 0;
        FOR v := u+1 TO N-1 DO
          g[u,v] := FALSE; g[v,u] := FALSE;
          candK[u,v] := ncands; candK[v,u] := ncands;
          candE[ncands] := Edge{u,v};
          INC(ncands);
        END
      END;
      <* ASSERT ncands = M *>
      nv[0] := 0;  ne[0] := 0; maxe[0] := 0;
      FOR s := 1 TO TN-1 DO 
        nv[s] := nv[s DIV 2] + (s MOD 2);
        ne[s] := 0;
        maxe[s] := D * (nv[s] DIV 2);
      END;
    END InitializeGraph;

  PROCEDURE PickCandEdge(): Edge =
    (* Picks a candidate edge at random: *)
    BEGIN
      <* ASSERT ncands > 0 *>
      RETURN candE[Random.Subrange(coins, 0, ncands-1)];
    END PickCandEdge;
    
  PROCEDURE EliminateCand(u, v: Vertex) =
    (* 
      Removes a candidate edge (u,v) from the candidate pool, by
      clearing its entries in candK and candE. 
      A no-op if (u,v) is not in the pool.
    *)
    BEGIN
      IF u > v THEN VAR t := u; BEGIN u := v; v := t END END;
      <* ASSERT u <= v AND v < N *>
      WITH 
        kuv = candK[u,v],
        kvu = candK[v,u]
      DO
        <* ASSERT kuv = kvu *>
        IF kuv # -1 THEN
          <* ASSERT ncands > 0 *>
          Msg("      eliminating edge " & FE(u, v));
          WITH
            ek = candE[kuv]
          DO
            <* ASSERT ek.org = u AND ek.dst = v *>
            IF kuv # ncands-1 THEN
              (* Fill hole candE[kuv] with candE[ncands-1]:*)
              WITH
                en = candE[ncands-1],
                un = en.org,
                vn = en.dst,
                kunvn = candK[un, vn],
                kvnun = candK[vn, un]
              DO
                ek := en;
                <* ASSERT kunvn = ncands-1 *>
                <* ASSERT kunvn = kvnun *>
                kunvn := kuv;
                kvnun := kvu;
              END
            END;
          END;
          kuv := -1;
          kvu := -1;
          DEC(ncands)
        END;
      END
    END EliminateCand;
    
  PROCEDURE BumpDegree(u: Vertex) =
    (* 
      Increments the degree of vertex u.  
      If it becomes D, all edges incident TO u are disqualified.
    *)
    BEGIN
      <* ASSERT deg[u] < D *>
      INC(deg[u]);
      IF deg[u] = D THEN
        Msg("    vertex " & FV(u) & " saturated");
        (* Eliminate all edges incident on u from the candidate pool: *)
        FOR v := 0 TO N-1 DO EliminateCand(u, v) END
      END
    END BumpDegree;
    
  PROCEDURE BumpEdgesInSubgraph(s: Subgraph) =
    (* 
      Updates ne[s] to account for the new edge u, v, assuming both
      are in s.  If s becomes full (one less edge than overfull),
      then eliminates all pairs of vertices in s from the candidate 
      pool. 
    *)
    BEGIN
      <* ASSERT ne[s] < maxe[s] *>
      INC(ne[s]);
      IF ne[s] = maxe[s] THEN
        Msg("    subgraph " & FS(s) & " full");
        FOR i := 0 TO N-2 DO
          WITH mi = Word.Shift(1, i) DO
            IF Word.And(mi, s) = mi THEN
              FOR j := i+1 TO N-1 DO
                WITH mj = Word.Shift(1, j) DO
                  IF Word.And(mj, s) = mj THEN
                    EliminateCand(i, j)
                  END
                END
              END
            END
          END
        END
      END
    END BumpEdgesInSubgraph;
  
  PROCEDURE AddEdge(u, v: Vertex) =
    BEGIN
      <* ASSERT u < v AND v < N *>
      <* ASSERT (NOT g[u, v]) AND (NOT g[v,u]) *>
      g[u,v] := TRUE; g[v,u] := TRUE;
      Msg("  added edge " & FE(u, v));
      EliminateCand(u, v);
      BumpDegree(u);
      BumpDegree(v);
      WITH 
        muv = Word.Or(Word.Shift(1, u), Word.Shift(1, v))
      DO
        FOR s := 0 TO TN-1 DO 
          IF Word.And(muv, s) = muv THEN BumpEdgesInSubgraph(s) END
        END
      END
    END AddEdge;

  PROCEDURE MinDegree(): CARDINAL =
    VAR d := N-1;
    BEGIN
      FOR u := 0 TO N-1 DO IF deg[u] < d THEN d := u END END;
      RETURN d
    END MinDegree;
    
  BEGIN
    REPEAT
      InitializeGraph();
      WHILE ncands > 0 DO
        WITH
          e = PickCandEdge()
        DO
          AddEdge(e.org, e.dst);
        END
      END;
    UNTIL MinDegree() = D AND NOT IsClassOne(g);    
    RETURN g
  END GenerateGraph;