Distributed Graph Layout for Scalable Small-world Network Analysis

The in-memory graph layout or organization has a considerable impact on the time and energy efficiency of distributed memory graph computations. It affects memory locality, inter-task load balance, communication time, and overall memory utilization. Graph layout could refer to partitioning or replication of vertex and edge arrays, selective replication of data structures that hold meta-data, and reordering vertex and edge identifiers. In this work, we present DGL, a fast, parallel, and memory-efficient distributed graph layout strategy that is specifically designed for small-world networks (low-diameter graphs with skewed vertex degree distributions). Label propagation-based partitioning and a scalable BFS-based ordering are the main steps in the layout strategy. We show that the DGL layout can significantly improve end-to-end performance of five challenging graph analytics workloads: PageRank, a parallel subgraph enumeration program, tuned implementations of breadth-first search and single-source shortest paths, and RDF3X-MPI, a distributed SPARQL query processing engine. Using these benchmarks, we additionally offer a comprehensive analysis on how graph layout affects the performance of graph analytics with variable computation and communication characteristics.

INTRODUCTION
Layouts of graphs and sparse matrices in distributed memory and shared memory have been well-studied for regular graphs that arise in the scientific computing domain. "Layout" in this instance refers to how vertices and edges are partitioned in distributed-memory and how the vertex identifiers are ordered in sharedmemory. Recently, several new open-source distributedmemory graph processing frameworks have emerged into mainstream usage. These include GraphLab [35] and its derivatives PowerGraph [27] and PowerLyra [17], Giraph [19], Trinity [47], and PEGASUS [31], among others. The primary goal of these frameworks is to analyze real-world graphs such as web crawls and social networks, which tend to be low-diameter graphs with skewed vertex degree distributions. Most of these frameworks assume an initial topology-agnostic vertex and edge partitioning and ordering. With these frameworks • G. M. Slota  and small-world irregular graphs in mind, this paper attempts to answer the following questions: 1) Will the layout of the graphs impact the performance of irregular, data-analytic algorithms and frameworks ? 2) Can such a layout be computed in a scalable and efficient fashion to be applicable in graph analytics ? 3) What kind of graph computations will be impacted by the graph layouts and how ? As has been observed, the impact of partitioning and ordering on irregular graph computations can be considerable [10], [21], [26]. However, using traditional layout strategies based on graph/hypergraph partitioners and orderings for data layout of highly irregular small-world graphs may not be appropriate for the following reasons: 1) Traditional partitioners and even some ordering methods, for example nested dissection, are heavyweight tools that are expensive both in terms of memory usage and time. They are appropriate when followed by more expensive linear solvers or when the partitioning results can be used for multiple solves. In contrast, graph analytic workloads are constantly evolving and a typical analytic operation is typically cheaper than a linear solver. 2) Previous ordering algorithms are designed for metrics appropriate for linear solvers such as minimizing a bandwidth [22] or minimizing the fill-in in a LU factorization [2], [33]. In contrast, ordering methods that improve the layouts in a shared memory context for small-world graphs are needed.
3) The performance of distributed-memory graph algorithms can be dependent on both local and global graph topology. Global topology affects the number of parallel phases and synchronization overhead, while local topological changes impact per-phase load balance. Optimizing for aggregate measures such as conductance or edge cut would ignore local topology changes and may not account for dynamic variations in per-phase execution. Graph computations on highly irregular graphs require a layout that depends on parallel partitioners and ordering methods that are highly scalable for very large graphs. Label propagation-based partitioners are shown to be useful for partitioning small-world graphs [38].
We utilize such partitioning algorithms (PULP [50]) to compute the distributed memory layout. Label propagation exploits the community structure inherent in many real small-world graphs to quickly partition even multibillion edge networks. Label propagation also allows for optimization of various objectives under multiple constraints, which enables us to explore the impact of these objectives and constraints on total execution and communication times for our various test analytics. In addition, we also introduce a breadth-first searchbased ordering that is more scalable than other ordering schemes and suitable for small-world graphs in the shared-memory layout. In case of distributed graph processing, we consider various partitioning-ordering possibilities, a simultaneous global partitioning and ordering of all vertices, and a local ordering of vertices after the partitioning phase.
In short, we propose a "distributed-memory graph layout" based on vertex partitioning using label propagation and a BFS-based parallel ordering strategy. The proposed DGL (Distributed Graph Layout) is a fast, memory-efficient, and scalable graph layout strategy. We demonstrate the new DGL layout scheme is about 10-12× times faster to compute than METIS partitioning [33], and about 2.3× faster to compute than Reverse Cuthill-McKee (RCM)-based orderings.
We demonstrate the impact of DGL and present detailed analysis on the end-to-end performance of distinct graph analytic workloads. The graph analysis routines include subgraph counting, breadth-first search (BFS), single-source shortest paths (SSSP), resource description framework (RDF) queries, and PageRank. The five algorithms were chosen to be representative of the diversity in modern graph analytics. We chose a recent algorithm for subgraph counting [49] which is a randomized parallel algorithm to generate approximate counts of treestructured subgraphs. Although recent related work [15], [16] primarily looks at strong scaling of BFS and related computations on massive synthetic Graph 500 networks, our work examines the subgraph counting algorithm, an analytic that is computationally very different from BFS. However, we also do an in-depth evaluation of BFS and SSSP performance. The fourth benchmark evaluates a distributed-memory implementation of the popular RDF store RDF-3X [42]. Our final included algorithm is a highly scalable implementation of PageRank [51], which is a popular and more computationally-intensive implementation than BFS for benchmarking performance of frameworks and systems.
We use the end-to-end graph analysis times for partitioning-ordering-workload in both single-threaded (MPI) and multi-threaded (MPI+OpenMP) distributed programming models. We also consider computation and communication times of the analytic separately, in order to better isolate the effects of partitioning and ordering on performance. We primarily consider realworld rather than synthetic graphs in our study. We use tuned implementations, all developed by us, in order to ensure consistency. We also analyze trade-offs between partitioning quality on computational load balance and communication overhead for several large real-world networks. The following is a summary of the key observations and findings from this workload analysis. We finally mention that DGL is not limited to the MPI processing models considered in this work, and can therefore be utilized as a preprocessing step while running under other graph engines and parallel execution environments.

DISTRIBUTED GRAPH LAYOUT
In this section, we discuss the distributed graph layout using label propagation-based partitioning and BFSbased ordering methods. We define a distributed graph layout as the pair of partitioning×ordering. The partitioning part of the layout affects the number of parallel phases and synchronization overhead in a graph computation. It is important to balance the computation in different parallel phases as well as minimize the communication overhead. We explore trade-offs in work and memory balance and communication minimization between tasks with different partitioning strategies. Work performed and memory utilization per-task roughly correlates with the number of vertices and adjacent edges stored on each task. The communication requirements roughly correlates with the number of inter-task edges, or edge cut resulting from partitioning. The ordering part of the layout affects the per-phase computation time in graph computations. We ideally want to increase intranode memory access locality to reduce cache misses and improve execution times. In order to be practical the partitioning×ordering pair must be computed in parallel, scalable fashion.

Partitioning
We utilize three partitioners in this work. We use a random partitioning to establish a baseline for benchmarking, which randomly assigns part assignments to each vertex. We use the well-known METIS [33] partitioner as a representation of the state-of-the-art.
We also utilize the PULP partitioner, which is specifically optimized to partition the small-world graphs we are considering in this work. We consider both single constraint and multi-constraint partitioning scenarios, where we either balance partitions for vertices or for both vertices and edges. We attempt to minimize total edge cut for both PULP and METIS. Additionally, for PULP, we also attempt to balance communication among parts by minimizing the maximal number of cut edges coming out of any single part.
The PULP partitioner is based off of the community detection label propagation algorithm [45]. Label propagation methods are attractive as they have low computational overhead, low memory utilization, are easy to parallelize, and demonstrate scaling to graphs with billions of vertices.
An overview of the basic label propagation algorithm is as follows: Initially, each vertex in a graph is initialized to having a unique label. Iteratively, each vertex then examines all of its neighbors' labels then updates to its label the label that appears most frequently among its neighbors, with ties broken randomly. The loop over all vertices can be parallelized without any explicit synchronizations or locking with minimal effect on solution quality [50]. We continue to loop over all vertices until no labels are updated, or, more commonly, after some number of iterations of the outermost While loop (usually 10 or fewer iterations is sufficient).

Algorithm 1 PULP Multi-Constraint Multi-Objective Algorithm Overview
Initialize p parts Execute degree-weighted label propagation. for k 1 iterations do Balance parts for vertex constraint. Refine parts to minimize edge cut. for k 2 iterations do Balance parts to satisfy edge constraint and minimize max per-part cut. Refine parts to minimize edge cut.
PULP's subroutines essentially use variants label propagation that limit the number of possible labels to the number of desired parts and impose additional weighting criteria to create balanced partitions. This weighted form of label propagation is utilized in two separate stages during execution of PULP. Algorithm 1 gives a very broad overview of the PULP multiconstraint (vertices and edge per part) multi-objective (minimize total edge cut and maximum edge cut per part) algorithm that demonstrates these two stages. After initialization, we first utilize weighted label propagation in k 1 alternating stages to balance the initial parts for our vertex constraint and then refining to minimize the total edge cut. Next, we perform k 2 alternating stages of balancing for our edge constraint while minimizing the secondary objective of max per-part cut and then again refining to minimize the total edge cut. In prior work [50], we describe the algorithm in considerably greater detail and demonstrate the approach's effectiveness in terms of cut quality and runtime with respect to other traditional partitioners. However, it is critical to show that such label propagation-based partitionings are not only easy to compute, but that they also improve the end-to-end runtimes of graph analytic applications. With DGL, we are able to utilize such a partitioner in the layout strategy and demonstrate its applicability for the first time.

Ordering
For a distributed graph computation, a good graph partitioning will reduce inter-node communication cost. The goal of on-node vertex ordering is to increase locality of intra-node memory references, and thereby reduce intra-node computation time. This is done by relabeling vertex identifiers so that consecutive accesses of pervertex specific information occur with greater spatial and temporal locality. As many graph computations access per-vertex data based on adjacencies, and per-vertex data is commonly stored in a flat array, minimizing the numeric difference between the vertex identifiers of adjacent vertex pairs can greatly improve access locality and therefore cache utilization.
Reverse Cuthill-McKee (RCM) is a commonly-used vertex ordering strategy in sparse matrix and graph applications. We propose a BFS-based ordering (see Algorithm 2) which can be considered an approximation to RCM. It avoids the costly sorting step used in RCM where it tries to order the nodes with the same parent in terms of the degree. Recently, a similar ordering was proposed for improving the matrix-vector multiply time and bandwidth reduction [32]. The primary focus of that approach was to arrive at parallel orderings to improve the linear solver time. Our focus is to improve the graph computations' end-to-end time.
We demonstrate our approach in Algorithm 2. We randomly choose a minimal-degree vertex as the root and perform a standard BFS routine, tracking visitation status with visited and the current level with level. We add vertices to level sets L when they are visited, as with RCM. We avoid explicit sorting by assuming that each L 0···M ax level , where M ax level is the maximum BFS level, is mostly sorted in the order of decreasing vertex degree, as there is a higher likelihood of encountering high-degree vertices sooner in any given level for most real world graphs. We assign new labels using an incrementing value of n by starting with the vertices in the highest level and working backwards to those in the lowest level. As we will show in the next section, this approach performs better than both random and RCM orderings in applications that have a high number of irregular memory accesses. As with RCM in [32], Algorithm 2 can be straightforwardly parallelized.
We compare the performance of our ordering method and RCM to a random ordering, where vertices are shuffled into a random order. Another common ordering strategy is using the "natural" ordering, which is using the original vertex assignments as given in the graph data file. We avoid comparison to the natural ordering (and associated vertex block partitioning) for a couple of reasons. The quality of the natural layout for graphs retrieved from a database can be highly variable. E.g. a natural layout for a web graph that's based on crawling methodology might give considerably better performance than the layout for a social graph that's based off of user ID (we observed a spread of both 2× speedup and slowdown for a natural ordering vs. random ordering on our PageRank test; similarly, we observed consistent variance in vertex block vs. random partitioning). Additionally, a "natural" layout for certain graphs is not necessarily well-defined. For example, the user IDs for a social network might be considerably larger than the number of vertices in the graph. To create and efficiently store a graph in-memory, these values need to be mapped to vertex identifiers 0 · · · |V | − 1). This mapping can be done by compressing the user IDs and retaining their original order, or through a more efficient first-come-first-served mapping where IDs are mapped to vertex identifiers as they are encountered in the data file. The method used would be applicationspecific. Because the extensive number of tests we performed limited the number of graphs we could use in our experiments, we wished to eliminate any potential impacts due to the above graph creation methodologies for the sake of uniformity in our comparative results. Hence why we view a random ordering to be the best baseline for relative comparison.
With the five partitioning methods (random, METIS {single constraint, single objective and multiple constraint, single objective} and PULP {multiple constraint, single objective and multiple constraint, multiple objective}) and three ordering methods (random, RCM, and DGL) we evaluate all the combinations of partitioning×ordering pairs and demonstrate that the DGL layout with PULP partitioner and DGL-based ordering performs the best in irregular graph computations.

PARALLEL GRAPH COMPUTATIONS
In this section, we will give an overview of the five distributed graph analytics used during our experimental analysis of the impact of partitioning and ordering on analytic performance. In an attempt to best understand the general effects of varying partitioning and ordering on the performance, the graph analytics were selected as to represent a wide range of execution characteristics. The test suite includes an implementation which is relatively computation-heavy, PageRank, algorithms which are relatively more communication-heavy, breadth-first search and single source shortest paths, an algorithm which is both very computation and communication intensive, color-coding subgraph counting, as well as an algorithm whose performance is dependent on the sizes of the n-hop neighborhoods of each partition, distributed query processing of Resource Description Framework stores.

Distributed PageRank
Our distributed PageRank uses an MPI+OpenMP approach and an |V | p partitioning, with each of p MPI tasks calculating the counts for an equivalent portion of the |V | vertices in the graph G. With one MPI task per node, we then use thread parallelism while updating the counts of owned vertices. With the exception of the single MPI communication call on each iteration, all per-task work can be done in parallel. Updates are passed among neighbors using an MPI all-to-all exchange. In practice, this specific implementation has been observed to be very efficient and scalable, giving per-iteration costs of less than a few seconds for networks of over 100 billion edges while running on 256 compute nodes. The specific technical details of the implementation are omitted, but please see [51] for a more in-depth discussion.

Subgraph Counting
Subgraph counting is a computationally challenging task, with the naïve approach scaling as O(n k ), where n is the number of vertices in a graph and k the number of vertices in the subgraph being counted. The best known exact algorithm [25] improves the exponent by a factor of α 3 , where α is the exponent for fast matrix multiplication. Because of these extremely high execution time bounds, recent work has focused on approximation algorithms. One such approach for counting tree-structured subgraphs utilizes the color-coding technique of Alon et al. [1].
In prior work, we developed a fast parallel implementations of color-coding subgraph counting in both shared-memory and distributed-memory environments [49]. The distributed version of our approach uses several optimizations, including fully partitioning and compressing the memory-intensive dynamic programming table to decrease memory requirements across all tasks, further compressing the table during communication to reduce the total transfer volume, and using all-to-all exchanges in lieu of broadcasts to reduce communication times. These optimizations demonstrate good scaling and enable us to count subgraphs of 10 and 11 vertices on billion-edge networks in minutes on a modest number of 16 nodes. For space consideration, we omit a detailed description of our implementation. Instead, please refer to [49] for an in-depth discussion of the stages and execution of the algorithm.

SSSP and BFS
We also assess the performance impact of layout on tuned implementations for parallel breadth-first search (BFS) and single-source shortest paths (SSSP) computation in this paper. Our parallel BFS approach can take advantage of both 1D and 2D graph distributions [10], [11], [12]. We use a 1D distribution in this work, as it is easier to correlate communication time with edge cut after partitioning with a 1D distribution. Recent BFS and SSSP implementations use a 1D partitioning and direction-optimizing search [4] for work-efficient and highly scalable execution on Graph 500 test instances. For an overview of the current state-of-the-art in performance optimizations for these routines, we refer the reader to [15], [16].
We use an optimized parallel implementation [44] of the ∆-stepping algorithm [37] for parallel SSSP in this paper. Each BFS iteration and ∆-stepping phase is comprised of three main steps: local discovery, all-to-all exchange, and local update. To aid adjacency queries, we use a distributed compressed sparse row representation for a graph. The distance array is also partitioned and distributed along with the distributed vertices (for ∆stepping). In the local discovery step, both algorithms expand their frontiers by listing all corresponding adjacencies and their proposed distance based on vertices in a queue of recently-visited vertices (for BFS) or in a current bucket (for ∆-stepping). Note that BFS visits each reachable vertex only once while ∆-stepping may visit each reachable vertex multiple times before it is settled.
Once all vertices in the queue are processed or the current bucket is empty (with no more vertex reinsertions), all p tasks exchange vertices in these generated lists to make them local to the owner tasks. This step is the same for both BFS and ∆-stepping, and uses an all-to-all collective communication routine. At the end of each BFS iteration and ∆-stepping phase, each task locally updates the distance of its own vertices using the exchanged information. The update in BFS is only on unvisited vertices, while ∆-stepping updates all vertices whose distances can be decreased. Thus, the ∆-stepping algorithm performs more computation and has a higher communication complexity.
Since our goal is to analyze and evaluate the effect of graph partitioning and vertex reordering, we have not yet implemented all the optimizations in [15], [16]. However, our approach has three new optimizations: (i) A semi-sort of vertex adjacencies based on weights is used prior to execution of the algorithm. (ii) Memoryoptimized queues are used to represent the bucket data structure. This decreases the algorithm memory requirement, while slightly increasing the running time. (iii) An array of all local unique adjacencies is created and locally used to track tentative distance of adjacencies. This array improves efficiency by filtering out unnecessary requests to be added in the new frontiers.

Distributed RDF Stores and SPARQL Query Processing
Resource Description Framework (RDF) [46] is a popular data format for storing web data sets. Informally, the RDF format specifies typed relationships between entities, and the basic record in an RDF data set is a triple. There are a growing number of publicly-available RDF data sets that contain billions of triples. Thus, database methodologies for storing these RDF data sets, also called triple stores, are becoming popular. We have developed a distributed MPI-based implementation of an open-source triple store called RDF-3X [42]. Our distributed RDF store is called RDF3X-MPI [20].
An alternate approach to viewing an RDF data set is as a directed graph with edge types. RDF data sets can be queried using a language called SPARQL. We extend the distributed RDF store methodology of RDF-3X to the SPARQL querying phase as well. Thus our RDF3X-MPI tool has two phases, a load phase and a query phase. In the load phase, the given triple data set is partitioned into several independent files, one per task, and each task then constructs an index for helping answering SPARQL queries. It is possible to parallelize some query evaluation in a purely data-parallel manner (i.e., with no communication between tasks), provided there is sufficient replication of triples among partitions. Formally, if the triple partitions satisfy an n-hop guarantee, then SPARQL queries in which all pairs of join variables are at distance of less than n hops from each other can be solved without any inter-task communication [30]. So the role of graph partitioning in this application is to reorder vertices such that the number of triples that are replicated between tasks after applying an n-hop guarantee are minimized. If the number of triples that are replicated is reduced, then the database indexes are smaller, making them potentially faster to query. For this application, we study the impact of partitioning on the number of replicated triples. A smaller value of replication is desired, and further, smaller index sizes should translate to faster query times.

EXPERIMENTAL SETUP
We evaluate performance of our new partitioning and ordering strategy DGL and the graph analytics workload on a collection of nine large-scale low diameter graphs, listed in Table 1. LiveJournal, Orkut, and Twitter (follower network) are crawls of online social networks obtained from the SNAP Database and the Max Planck Institute for Software Systems [14], [52]. uk-2005 and sk-2005 are crawls of the United Kingdom (.uk) and Slovakian (.sk) domains performed in 2005 using UbiCrawler and downloaded from the University of Florida Sparse Matrix Collection [6], [8], [23]. WebBase is similarly a crawl obtained in 2001 by the Stanford WebBase crawler. We created the BSBM and LUBM graphs from RDF data sets generated using the Berlin SPARQL benchmark [5] and Lehigh University Benchmark [28] generators. DBpedia was created from RDF triples extracted from Wikipedia [40].
The Orkut graph is undirected and the remaining graphs are directed. For the web and social graphs, we preprocessed the graphs before executing PageRank, BFS, SSSP, and subgraph counting. Specifically, we removed all degree-0 vertices, multi-edges, and extracted the largest (weakly) connected component. Further, edge directivity was ignored when partitioning the graphs using PULP and METIS and reordering with RCM and DGL. Table 1 lists the sizes of these nine graphs after preprocessing. The scalability studies for subgraph counting, BFS, SSSP, and RDF query processing were done primarily on Blue Waters, a large petascale supercomputer at the National Center for Supercomputing Applications (NCSA). Each XE compute node of Blue Waters is a dual-socket system with 64 GB main memory and AMD 6276 Interlagos processors at 2.3 GHz. The system uses a Cray Gemini 3D torus interconnect. We built our programs with the GNU C++ compiler (version 4.8.2), using OpenMP for multithreading and the -O3 optimization parameter during compilation. For the pre-processing phases of DGL (partitioning and reordering) and some scalability runs, we utilized Compton, a testbed cluster. Compton has a dual socket setup with Intel Xeon E5-2670 (Sandy Bridge) CPUs at 2.60 GHz and 64 GB main memory. Due to the large memory requirements of partitioning with METIS, we also had to use the large memory nodes on Carver at NERSC for partitioning the larger networks (Twitter, uk-2005, Webbase, and sk-2005). Carver's large memory nodes have 1024 GB main memory and four Intel Xeon X7550 ("Nehalem-EX") CPUs at 2.00 GHz. We performed k-way partitioning with METIS using version 5.1.0.

DGL Performance Evaluation
We first evaluate our DGL label propagation-based partitioning methodology, PULP, against METIS partitioning by examining total running time for generating 16 and 64 partitions. We consider two versions of both PULP and METIS. For PULP, we have an implementation that has both maximal vertex and edge balance constraints and minimizes both total edge cut and maximal per-part edge cut. We consider this our baseline implementation, and label it in figures as PULP-MM (PULP multi-objective multi-constraint). We also have a dual constraint version that only attempts to minimize the total edge cut, which we call PULP-M. Similarly for METIS, the dual constraint single objective version is termed METIS-M, while the single constraint (vertex balance) and single-objective version is termed simply as METIS. METIS-M and PULP-M are solving the same problem. For our constraints, we fix the maximal vertex imbalance ratio at 1.10 and the edge imbalance ratio at 1.50. The results will show that the multi-constraint, multi-objective mode of PULP-MM can be important for irregular graph computations. Table 2 shows the partitioning time of PULP-MM along with METIS-M running on Compton. Due to METIS's large memory requirements (close to 500GB for Twitter), only LiveJournal, Orkut, and the RDF graphs could be partitioned on Compton. The larger web graphs and Twitter were all partitioned on a large memory node of Carver. We also report the relative speedup of PULP to METIS. We omit time comparison to ParMETIS, as the only graphs it was able to successfully partition on any system were LiveJournal and Orkut. Further, ParMETIS's speedups relative to METIS for those two instances were minimal (less than 2× with 16-way parallelism). From Table 2 we observe considerable speedup for PULP, with a geometric mean speedup of 12.4× for 16 parts and 10.1× for 64 parts.  The partitioning quality in terms of both vertex and edge balance constraints and edge cut and maximal perpart edge cut objectives for the different partitioners is shown in Table 3 as geometric averages. We also note that aggregate measures don't fully capture the wide spread of results among different tests, so include min and max improvements for edge cut and max per-part cut as well. E.g., the improvement METIS has relative to random partitioning varies from 1.5× for 64-way partitioning of Twitter to 107× improvement for 16-way partitioning of uk-2005.
In terms of the total edge cut (EC), the singleconstraint, single-objective METIS does the best, but it performs poorly in the maximum per-part edge cut (EC max ) and edge balance (E max ). PULP-MM also performs better than all the methods in the EC max metric without sacrificing a lot in EC and still respecting the vertex balance and edge balance constraints. Also note the much larger E max of single constraint METIS. As we will demonstrate, this can have a considerably impact of execution time for the applications in our benchmarks. We bold the best values for each column for edge cut and max per-part cut, and note that METIS performs consistently best overall in the edge cut metric and PULP-MM performs best overall in the max per-part metric by a wide margin.
We additionally compare our DGL vertex ordering strategy to RCM. Table 4 gives the average running  [32]. However, their timings are insignificant in the end-to-end performance of complex analytics such as our subgraph counting benchmark. We demonstrate that the DGL and PULP strategies are also able to efficiently compute the layout for larger numbers of parts beyond 16 and 64. Table 5 gives the execution times of DGL ordering and PULP partitioning when computing the layout of the various test graphs with 256, 512, and 1024 parts. We observe flat scaling of DGL ordering with increasing part counts due to its intrinsic work efficiency and O(n + m) expected execution time. We observe an increase in PULP times for higher part counts. This is due to PULP having a per-iteration workload of O(np + m) [50], where p is the number of parts being computing. However, we still observe sub-linear scaling relative to p for almost all test cases. Overall, we observe no intrinsic scalability bottlenecks of these methods at this higher scale.
We include one more table to demonstrate how our DGL ordering strategy might improve cache performance of executing codes. To improve the performance of linear solvers, a common ordering metric to optimize for is graph bandwidth, which is the maximum integer distance between vertex identifiers for vertices that share a single neighbor. RCM is an effective means of bandwidth reduction for regular matrices. However, for small-world graphs, the bandwidth is usually going to be large, on the order of d max , where d max is the maximal degree of any vertex in the graph. Comparing bandwidth measures between different orderings therefore won't show any global improvements in compaction for rows of much lesser degree vertices. As such, we look at other metrics to give an indication of the possible cache efficiency in practice. Across the entire adjacency array, we measure how often edges listed in order also have identifiers within a single integer value of each other. This indicates that these edges would be neighboring nonzeros in the same row of an adjacency matrix. Co-located edges improve cache utilization of per-vertex information accesses, such as checking visitation status for BFS or PageRank value lookups. To quantify how many co-located vertex identifiers for the edges are in the adjacency list, we report two values. First, we report a ratio of how co-located all edges are, where a value of zero indicates that no edges are co-located and a value of one indicates that all edges are co-located. Second, we report a running "cost" as the sum of the distances, or gaps, between vertex identifiers in the adjacency list. We scale the distances by their log, as a distance of one or close to it indicates that vertices are closely co-located and would have minimal cache cost for their subsequent accesses, and the cost difference between large and very large distances is minimal, since it's likely a new cache line would need to be loaded in both instances. Finding an ordering that minimizes this sum is referred to in the literature as the Minimum Logarithmic Gap Arrangement problem, which is NP-hard [18]. In an attempt to give a graph-independent ratio, we further scale the log sum by a worst-case possible value of m log n. The true dependence of cache utilization on distance would be architecture-specific, but the approximation of this cost gives enough insight for comparative purposes when examining ordering quality. Table 6 gives both the co-location ratio (Co-loc. Ratio) as well as the log sum of gap distances ratios (Gap Sum Ratio) for all ordering combinations across all graphs for 16 and 64 parts. We report the geometric mean values across all five partitioning strategies (Random, METIS, METIS-M, PULP-M, PULP-MM). For co-location ratio, higher indicates better locality, while for the log gap sum ratio, lower indicates better locality. We omit reporting the co-location ratio for Random ordering in Table 6, as all values are close to zero, a few orders of magnitudes less than RCM and DGL. We observe nearly that DGL ordering results in the best co-location ratio and lowest log gap sum ratio across almost all instances. The computational timings results we'll report next in our benchmarks will demonstrate that these measurements translate into real performance benefits across a wide range of graph analytics. For our first set of experimental benchmarks results, we examine the effect of partitioning and ordering on a distributed PageRank implementation. We will first show the effect that different partitionings have on communication times, and then we will show the effect that orderings have on computation times. For these experiments we use the three social network graphs (LiveJournal, Orkut, and Twitter) as well as the three web crawls (uk-2005, WebBase, sk-2005). Figure 1 (top) gives the speedups relative to random partitioning for METIS single and multiple constraint partitionings and PULP multiple constraint with single and multiple objective partitionings. Figure 1 (bottom) gives the speedups relative to random ordering for DGL and RCM. Table 7 gives the explicit speedup values and overall geometric means across the six test graphs. These value are for 20 iterations of PageRank executing on 16 nodes of Blue Waters.

PageRank Performance
We observe that all partitionings offer considerable speedups relative to random. In general, the web crawls show even greater speedups than the social networks. This is due to the web crawls being greater in diameter and more separable than social networks, resulting in a decrease in the number of cut edges and subsequently greater performance improvements relative to random. Averaged across all six test graphs, PULP multiple constraint and multiple objective partitioning offers the greatest speedup. The performance benefit is due to the implementation's use of an iterative bulk synchronous model and moderately low required communication, so the improved communication balance resulting from PULP-MM's decrease in max per-part cut becomes apparent in the timings.
Additionally, we note that both RCM and DGL offer considerable speedups for total computation times relative to random ordering. On the uk-2005 and WebBase graphs, the speedups for DGL are about 5×. Again we observe that the social networks generally show    less performance benefit relative to random, and this is again due to their lower diameter and small-world characteristics, which makes effective ordering more difficult. However, we still observe a consistent 20%-40% speedups with the improved orderings. Overall, DGL gives a greater performance speedup over RCM by about 10%, a result we expected based on our measurement of potential locality and cache performance as demonstrated in Table 6.

Subgraph Counting Performance
We next compare the impact of various partitioning and ordering strategies with regards to the running times of our subgraph counting implementation. We run on 16 node of Blue Waters. We compare communication times resulting from each of the 5 partitioners with a fixed random ordering. We also compare the computation times resulting from the 3 ordering strategies with fixed PULP-MM partitioning. The speedups for each strategy on the 6 test graphs are given in Figure 2 and Table 8. We also look at total end-to-end execution time for the five partitioning strategies with random ordering in terms of total time spent in the communication, computation, and partitioning steps. Note that the results with single constraint METIS for the uk-2005 and WebBase graphs are absent. This is due to execution times taking longer than 24 hours for these instances. Several trends can be observed in Figure 2. The top subfigure gives the speedup of the communication phase of subgraph counting for each of the partitioning strate-  gies relative to random partitioning. We again note considerable speedup for all partitioners. We note that the PULP methods give the best improvement for five out of the six tested graphs. Overall PULP-M gives the highest speedup overall. This implementation doesn't benefit as highly from the more communication-balanced PULP-MM partitioning due to the overall higher communication requirements (the Twitter graph requires compression and transfer of several terabytes of data in total for the Count table exchanges between tasks) and lower overall synchronization cost relative to PageRank, so total edge cut is observed to have a greater effect in practice. This emphasizes the fact that a one-size-fits all solution is not optimal in practice, and implementation knowledge is required to extract the best performance for any given running application when utilizing a layout strategy.
The middle subfigure of Figure 2 plots the speedup relative to random ordering for the DGL and RCM reordering strategies with PULP-MM partitioning. Overall, we note about a 6% improvement for DGL and 5% improvement for RCM ordering relative to random. These improvements are much lower than PageRank's improvement due to considerably more information stored per-vertex in the stored counts table, so greater cache locality has less of an effect in preventing re-accesses to main memory; however, we note even a modest 5%-6% consistent improvement can be noteworthy in this instance. On processors with larger cache, this relative improvement would be expected to increase.
Finally, the bottom subfigure of Figure 2 shows the total end-to-end execution time in seconds for initial partitioning plus running of the subgraph counting application. We further split subgraph counting into the sum of time spent in each of its computation and communication phases. We observe that our partitioning and ordering strategies result in the fastest end-to-end running times for all test instances. The time spent for partitioning is considerable relative to execution time for METIS, as is the extra communication costs that result with random partitioning. The additional partitioning time cost for METIS might be amortized in practice by re-using the same partitions for subsequent analysis, but we note that PULP partitioning shows an immediate decrease in total end-to-end time after a single analytic run.

Execution Timelines
To offer visual explanation of the performance of balanced constraint partitioning on total execution time, we give execution timelines in Figure 3 of a single run of counting a 10 vertex template on the LiveJournal graph and 10 iterations of PageRank on the Webbase graph. We used the Compton system for these tests and random, single-constraint METIS, multi-constraint METIS, and PULP-MM partitioning (from left to right, respectively) with random ordering. On looking at subgraph counting (top), we note first the two extreme cases. Random shows the lowest total computation times at a high cost of communication, while single objective METIS results in low communication times but high total times during the execution stages. This is due to unbalanced work among each task, which is directly proportional to the edge balance among each part. We observe that balanced multi-objective PULP partitioning gives the best tradeoff in terms of work balance and communication requirements. For PageRank, we note a large performance gap between Random partitioning and the other strategies. This is due to the implementation's computational and communication requirements for each task being dependent on the one hop neighborhood and per-part cut. These values are much higher with Random partitioning. We observe that PULP gives the best performance, due to the fact that the multiple objectives are explicitly optimizing for these metrics while keeping work balance very consistent. Overall, we notice about a 5-10% total execution time improvement for PULP versus the METIS variants by using multi-objective partitioning. As noted, considering total end-to-end execution time with partitioning costs, this speedup would be even more dramatic.

SSSP and BFS Performance
In this section, we analyze the performance of our SSSP and BFS implementation when using the different par-  titioning and ordering layouts. These benchmarks were run on 64 nodes of Blue Waters. While the running time of distributed subgraph counting is dominated by largescale data transfers during the communication phases, SSSP's performance is more dependent on intra-task computation, similar to PageRank, but has considerably less communication. BFS has the overall lowest communication and computation requirements out of all of the benchmarks thus far. Figure 4 and Table 9 show the speedups for communication and computation for SSSP performance with            64 MPI tasks. The top subfigure of Figure 4 shows the communication speedups relative to random partitioning for the other partitioning strategies. Due to the relatively lower communication requirements for this SSSP implementation, we correspondingly observe lower speedups relative to what was observed in Figure 2 with subgraph counting. We note that METIS gives the highest communication speedup, due to the lower overall communication load and total synchronization costs which emphasize a lower total workload than explicit balance. We observe speedups for computation times on all graphs, and especially the web crawls, with both ordering strategies. DGL ordering gives around 30% speedup overall.
The two subfigures of Figure 5 and Table 10 give the speedup in communication time with different partitioners and speedups in computation time with different orderings for the BFS implementation. We notice similar trends to SSSP in these plots. Overall, the lower total computation and communication workload of BFS contributes to lower speedups when using better ordering and partitioning strategies compared to random.

SPARQL Query Processing
In this final section, we study the impact of partitioning and ordering on the performance of RDF stores and SPARQL querying, a benchmark algorithm that is very different than the previous ones. In Table 11, we report replication ratios observed when an undirected 2-hop guarantee is enforced. Our RDF3X-MPI implementation uses a 2-hop guarantee to partition the graph, and a lower replication ratio indicates a smaller index size, which should translate to faster query times in practice. Table 11 compares PULP-MM with METIS-M and random partitioning for 16 and 64 parts on the 3 RDF graphs. Out of the 6 total graph-part count scenarios, the PULP-MM approach shows the lowest replication ratio for half of them. Note that none of these partitioners are explicitly optimizing for this metric, so the performance of PULP-MM in this instance is indirect. In Table 12, we report sum of query times of RDF3X-MPI averaged over the BSBM, LUBM, and DBpedia data sets. We use a selection of queries from the Berlin SPARQL Benchmark. We use the 16 part partitions for this test and additionally look at the performance affects of the three ordering strategies. PULP-MM partitioning with random ordering yields the best performance, while PULP-MM further demonstrates the highest performance when using the other two ordering strategies as well. This corresponds to PULP-MM having the lowest replication ratios. We note that since PULP-MM is faster and much more memory-efficient than METIS, this is a promising result, and future work can attempt to optimize PULP for the one and two hops replication ratio metrics for further improvements. The effect of ordering strategy on query times is interesting, in that the higher-locality orderings demonstrate correspondingly worse performance. To store the RDF data, RDF3X-MPI converts the input RDF graph structure into multiple indexes, which are created by sorting the RDF data, creating B+ trees, and then performing compression. We note that the worsened performance with localityoptimized ordering is most likely an artifact of this preprocessing stage. This further indicates that knowledge of a graph analytic's algorithmic details is important when determining an optimal graph layout, as unexpected and counter-intuitive performance impacts are a real possibility.

RELATED WORK
We selected METIS and RCM for comparison to the partitioning and ordering aspects of DGL, as they represent the most popular and current state-of-the-art approaches for these problems in terms of both speed of computation and quality produced. There are various other partitioning algorithms and methods, including multi-level partitioners similar to METIS [3], [29], [38], coordinate and geometry-based partitioners [9], and hypergraph partitioners [13]. Hypergraph partitioners can often calculate higher quality partitions than graph partitioners for regular matrices, but at a considerably higher cost to compute. Other graph partitioners have utilized label propagation in single or multilevel approaches [39], [54], [55], [56], demonstrating improved algorithm execution times with these partitions versus naïve methods. However, while some of these partitioners produce very high partition quality with good computational efficiency [39], they only consider single constraint partitioning scenarios. As we've demonstrated, using multi-constraint partitioning is important for optimal algorithm performance with high computation loads. Though we acknowledge that the tradeoff of pre-processing time for runtime performance with communication-bound algorithms is going to be application-specific. Other recent work [24] has correspondingly demonstrated that complex partitioning scenarios beyond single objective partitioning optimizing for edge cut and/or communication volume is a necessary consideration for optimal performance in other distributed computations.
Recently, other partitioning methods have been developed with goals similar to that of PULP, in terms of striking the balance between both high scalability and high quality of part computation. We compare to the most notable of these, FENNEL [53], using their published results. On a large Twitter dataset, FENNEL reports an edge cut relative to METIS of 0.56×, 1.19×, and 1.33× (lower is better) for computing 2, 4, and 8 parts, respectively. However, the relative imbalance for FENNEL is higher than the balance constraint imposed on METIS, which makes a lower edge cut considerably easier to achieve. Using the same graph, we impose equivalent balance constraints for PULP and METIS and compute relative edge cuts of 1.19×, 1.05×, and 1.13× again for 2, 4, and 8 parts. FENNEL reports a time of 40 minutes to partition the Twitter graph. On this graph, which is smaller than the Twitter graph we used in our primary results, we compute these partitions with PULP in about 5 minutes. For 16 parts on LiveJournal, use of FENNEL speeds up PageRank execution relative to random (hash) partitioning by about 1.18×. On the same graph and number of parts, we reported an improvement of 2.83× with PULP relative to random partitioning.
In addition to RCM ordering, Cuthill-McKee (CM) [22], nested dissection [33], and Approximate Minimum Degree (AMD) [2] are a few examples of sparse matrix reordering strategies used in the past. Some techniques, such as space-filling curves [48] or spectral bisection and orderings based on calculated eigenvectors [43] have been utilized for both partitioning and ordering of sparse matrices. Ordering methods on irregular networks such as social and Internet graphs has been studied for the purposes of visualization [41] and compression [7], [18]. Although the authors know of no performance analysis of the effects of applying these ordering techniques to distributed graph computation in literature, promising future work might involve utilizing and optimizes these ordering for such purposes.
We omit direct comparison using the distributedmemory graph processing frameworks mentioned previous [17], [19], [27], [31], [35], [47], as prior work has demonstrated a several orders-of-magnitude performance difference between them and optimized code [36], [51]. However, we acknowledge the goal of these frameworks is often programmer efficiency rather than pure performance. Implementing our methods within such frameworks would make for interesting future work.

CONCLUSIONS
In this paper, we present a methodology for distributed graph layout (partitioning, vertex ordering). The partitioning method for our layout, PULP, is based on the scalable label propagation community detection method. The partitions produced are comparable in quality to the k-way multilevel partitioning scheme in METIS, but only take a fraction of the execution time. Our DGL vertex ordering strategy can also improve computational performance of graph computations that consist of a high proportion of irregular accesses.
Additionally, we give a comprehensive performance analysis by examining the effect of graphs layouts on a graph analytics workload. In general, we note that graph analytics which have a high relative computation cost can benefit greatly from a locality-optimizing vertex ordering strategy. Graphs analytics that have a relatively high communication volume but few synchronizations might benefit most from a partitioning that optimizes explicitly for edge cut, while computations which consist of numerous synchronizations would benefit from a more balanced partitioning in terms of per-task communication loads. In general, the higher the computation to communication ratio for an analytic, the greater the impact of partitioning and ordering.


extfieldsofstudy:
Computer Science
Mathematics

text:
An Exact Spectrum Formula for the Maximum Size of Finite Length Block Codes

An exact information spectrum-type formula for the maximum size of finite length block codes subject to a minimum pairwise distance constraint is presented. This formula can be applied to codes for a broad class of distance measures. As revealed by the formula, the largest code size is fully characterized by the information spectrum of the distance between two independent and identically distributed (i.i.d.) random codewords drawn from an optimal distribution. A new family of lower bounds to the maximal code size is thus established, and the well-known Gilbert-Varshamov (GV) lower bound is a special case of this family. By duality, an explicit expression for the largest minimum distance of finite length block codes of a fixed code size is also obtained. Under an arbitrary uniformly bounded symmetric distance measure, the asymptotic largest code rate (in the block length $n$) attainable for a sequence of $(n,M,n\delta)$-codes is given exactly by the maximum large deviation rate function of the normalized distance between two i.i.d. random codewords. The exact information spectrum-type formula also yields bounds on the second-order terms in the asymptotic expansion of the optimum finite length rate for block codes with a fixed normalized minimum distance.

I. INTRODUCTION
The determination of the maximal size M * n (d) of a finite length block code with pairwise distance d and block length n has been a long-standing problem in information and coding theory. In this paper, we establish an exact formula for this fundamental quantity for all block lengths n ≥ 1. More specifically, we show that M * n (d) is equal to the maximal reciprocal of the cumulative distribution function Pr[µ(X n , X n ) < d] over all i.i.d. pairs of random vectorsX n and X n . This tight characterization gives a new way to construct good codes, and allows us to derive improved bounds on the optimal code size for a given distance.
Before describing our results, we introduce some of the previously known solutions and bounds, both in the asymptotic regime and for finite block length.

A. Related Works on Asymptotic Bounds
A common approach to understanding the maximum size of a block code is to allow the block length n to grow, and to apply asymptotic analysis of the largest code rate (1/n) log M * n (d), subject to a normalized distance constraint d/n ≥ δ. Many bounds on this asymptotic quantity have been derived, for which interested readers are referred to [1]- [10] and the references therein and thereafter. The best-known upper and lower bounds are perhaps the linear programming upper bound [1] and the Gilbert-Varshamov (GV) lower bound [7].
Both of these bounds were stated originally for binary alphabets and the Hamming distance, but they can be generalized [11]. The GV lower bound remained the best lower bound for decades until the invention of Goppa codes [12], based on which a better lower bound, called the algebraic geometry (AG) bound, was established for certain range of δ when the code alphabet size is an even power of a prime and no less than 49 [13], [14]. Later, Zinoviev and Litsyn proved that a better lower bound than the GV lower bound is actually possible for any code alphabet size larger than 45 [15]. Nonetheless, a gap remains between these asymptotic bounds.
In 2000, Chen, Lee and Han [16] used ideas from information spectrum analysis [17], [18] to establish an exact formula for the asymptotic largest minimum distance (in the block length n) of deterministic block codes under generalized distance measures. Dually, the largest minimum distance considered in [16] is formulated as lim n→∞ (1/n)d * n (e nR ) with d * n (M) max C⊂X n : |C|≥M min x n ,x n ∈C andx n =x n µ(x n , x n ), where | · | denotes the size of a set and µ(x n , x n ) represents the distance betweenx n and x n .
A major advantage of the exact formula in [16], as with most information spectrum results [18], is that its validity is not limited to the usual Hamming distance but can be applied to arbitrary distance measures. Furthermore, the derivation in [16] can be applied to arbitrary code alphabets, including continuous ones. This exact formula shows that the asymptotic behavior of the largest minimum distance of block codes is completely determined by the statistical properties of the normalized distance function (1/n)µ(X n , X n ) evaluated under optimally chosen i.i.d.X n and X n . This result complements what has been obtained for channel capacity [17] and its optimistic version [19], where the information density is shown to be fundamental in the determination of the general formulae for channel capacities.

B. Related Works on Finite Length Bounds
Although the results stated above are for block lengths tending to infinity, in practice, encoding and decoding complexities often grow with the block length. Hence, it may not be practical to use long block lengths in order to approach asymptotic performance guarantees.
Most of the asymptotic bounds mentioned at the beginning of this section, such as the linear programming upper bound [1] and the GV lower bound [7], also admit finite-length expressions, and there are many other well-known upper bounds, including the Singleton, Plotkin, and Elias bounds [20]. In general, however, these bounds are not tight.
In recent years, there has been renewed interest in understanding the performance of block codes with a fixed, finite block length (see, e.g., [21]- [27]). As a particular example, the authors in [21] established bounds on the maximal achievable transmission rate subject to a tolerable error probability and a finite block length. As with previous bounds, the results cited are not exact, and gaps generally remain between the upper and lower bounds.

C. Main Contributions
The main contributions of this work are as follows: 1) We provide an exact information spectrum formula for the maximal code size M * n (d), given in Theorem 1. This simple formula completely characterizes M * n (d) by an information spectrum expression. Our formula stands in stark contrast to existing finite length studies (e.g., [21]) in which the non-asymptotic bounds (e.g, dependence testing, meta-converse) obtained are, in general, not tight. An immediate implication of our result is a new family of lower bounds to M * (d), which includes the well-known finite length GV bound as a special case.
2) The Huffman-coding-like codeword-combining procedure used in the proof of Theorem 1 suggests an iterative algorithm for the construction of a code subject to a minimum pairwise distance constraint. Numerical experiments show that our procedure yields codes that outperform the finite length GV lower bound.
3) We extend our results to the determination of the largest asymptotic code rate and for the largest asymptotic minimum relative distance. These provide alternative, and arguably simpler, expressions compared to the asymptotic exact spectrum formulas in [16]. The limiting rate and distance are shown to be characterized by the large deviation rate function for (1/n)µ(X n , X n ) [28], [29]. Novel bounds on the second-order terms in the asymptotic expansion of the optimum finite length rate for block codes with a fixed normalized minimum distance are obtained.

D. Paper Organization
The rest of the paper is organized as follows. The exact spectrum formula for M * n (d) is stated and proven in Section II. A generalization of the finite length GV lower bound is proposed in Section III. The development of an iterative algorithm for the construction of a block code satisfying a given minimum distance is also presented in this section. The spectrum expression for the largest minimum distance of finite length block codes is addressed in Section IV. Extensions to the asymptotic regime, including bounds on the second-order terms, are studied in Section V. Open problems are discussed in Section VI.

II. LARGEST CODE SIZE ATTAINABLE UNDER A FIXED MINIMUM PAIRWISE DISTANCE
We introduce the notation in this paper. Let X be a finite set. 1 . . x m,n ) belongs to X n [30]. An (n, M, d)-code denotes an (n, M)-code with the minimum pairwise distance among codewords equal to d, i.e.,

d = min
x n i ,x n j ∈C and x n i =x n j µ(x n i , x n j ).
The maximal code size M * n (d) subject to a fixed pairwise distance lower bound d is where N is the set of positive integers. The distance measure µ(·, ·) is assumed to be symmetric, non-negative, and zero between a point and itself. In other words, for any two elementsx n and x n in X n , µ(x n , x n ) = µ(x n ,x n ), and µ(x n , x n ) ≥ 0 with equality holding if and only ifx n = x n . In order to simplify notation, we will sometimes use x and X, instead of the conventional notations x n and X n in the information spectrum literature [18], to respectively denote a single element and a random variable taking values in X n .
A. An Exact Spectrum Formula for M * n (d) Definition 1: Define the distance spectrum corresponding to a distance measure µ(·, ·) and a distribution P X on X n as where in this definition (and also throughout the paper),X and X are used to denote two independent random variables over X n having common distribution P X .
If we consider attempting to construct a code of distance d by drawing i.i.d. samples from X, the quantity F X (d) yields the probability that any two of the selected codewords violate the distance constraint. It turns out, perhaps surprisingly, that this quantity provides an exact information spectrum characterization of M * n (d). This key result is summarized in the next theorem.
Theorem 1: For all n ≥ 1 and d > 0, Before giving the proof of Theorem 1, which appears in Section II-B, we provide some remarks on the theorem itself, including some of its implications.
The determination of M * n (d)-that is, finding an "information" expression in terms of d and µ(·, ·) that is exactly equal to the operational quantity M * n (d)-has been a long-standing problem in information and coding theory. Up until now, only upper and lower bounds on this quantity were known. The above theorem beautifully shows that M * n (d) can be fully characterized by a simple distance spectrum formula and confirms definitively that µ(X n , X n ) is a "sufficient statistic" for the determination of M * n (d). It is also interesting to point out that the reciprocal of F X (d) is often not an integer for a general distribution P X . However, after we have optimized 1/F X (d) over all distributions, the optimal value is indeed an integer for all n ≥ 1.
An immediate consequence of Theorem 1 is that a family of lower bounds to M * n (d) can be obtained in a straightforward manner by evaluating 1/F X (d) for any X. We summarize this observation in the following corollary.

Corollary 1 (Distance spectrum (DS) lower bound):
For all n ≥ 1 and d > 0, for any distribution P X over X n , .
Note that if there does not exist a code of size two, in which the pairwise distance is no less than d, then M * n (d) = 1. This observation coincides with (2) and (3) as F X (d) = 1 for any X if d exceeds maxx ,x∈X n µ(x, x).
The above corollary implies that a good lower bound can be obtained by employing a good distribution P X over X n . If the optimizer to (2) is used to evaluate 1/F X (d), the associated bound is indeed tight! A question that naturally follows is what we can possibly say about the set of optimizing distributions. Despite the challenge of determining optimizers to the optimization problem in Theorem 1, extensive numerical studies for binary codes of small block lengths n leads to the following observations. First, the optimizer of the optimization problem in Theorem 1 may not be unique. Secondly, we have observed in numerical experiments that the set of optimizers seems to always include , of which all one-dimensional marginal distributions are uniform over X , i.e., P X * i (x) = 1/|X | for every x ∈ X and every 1 ≤ i ≤ n. For example, an optimal distribution for codes of block length n = 3 subject to a minimum It can be verified that , for x ∈ {001, 010, 100, 110}; 0, otherwise, which again yields Third, the two examples just introduced indicate that the optimizer may vary with d and is uniformly distributed over a subset of X n . These observations motivate the main idea behind the proof of Theorem 1.

B. Proof of Theorem 1
The theorem will be proven in two steps: For an (n, , let Z be the distribution which places probability mass 1/M on each codeword of C. Then, .
This immediately implies that Suppose that the support S(Z) of a distribution P Z consists of ℓ distinct elements, i.e., , there exists another distribution P W , whose support contains only M * n (d) elements, satisfying As a result of this claim, inf X F X (d) can be simplified to We thus continue from (5) to show that for any distribution P X with |S(X)| = ℓ ≤ M * n (d), where (6) follows from the Cauchy-Schwarz inequality, 2 and (7) is due to the assumption . 2 The Cauchy-Schwarz inequality states that which recovers (6) by setting ri = 1 for 1 ≤ i ≤ ℓ.
It remains to substantiate the claim in (4), which can be proven through a Huffman-codinglike [31] codeword-combining procedure as follows. Among the support S(Z) = {x 1 , x 2 , . . ., x ℓ } of size ℓ > M * n (d), there must exist two distinct elements with distance smaller than d. Without loss of generality, we let this pair be x 1 and x 2 , and denote for 1 ≤ k ≤ 2, where 1{·} is the indicator function. Assume, without loss of generality, that q 1 ≥ q 2 .
Construct P W (1) as Note that the support of P W (1) consists of only (ℓ − 1) elements, i.e., where (10) holds because equalities (11) and (12) follow from the definition ofp i in (9), and the inequality in (13) follows from the assumption that q 1 ≥ q 2 .
If |S(W (1) )| = ℓ − 1 is still larger than M * n (d), we can similarly construct P W (2) satisfying By repeating such a construction (ℓ − M * n (d)) times, we obtain the desired distribution This completes the proof of Theorem 1.
Remarks. We conclude this section by pointing out that the validity of the proof does not require the assumption that the alphabet X is finite but only requires the finiteness of M * n (d). As long as M * n (d) < ∞, the claim in (4) holds. In addition, the non-negativity assumption on the distance measure µ(·, ·) is not necessary. We only require the property that As a result, Theorem 1 can be extended to arbitrary code alphabets and arbitrary symmetric distance measure satisfying (14), provided that d > µ 0 and M * n (d) is finite.
In this section, further explorations based on the theoretical result in the previous section are conducted. Specifically, Section III-A verifies that the finite length GV lower bound is a special case of the DS lower bound, Section III-B explores the achievability (to M * n (d)) of a subclass of the DS lower bound under uniform X, and Section III-C presents an algorithmic construction of block codes with code sizes exceeding the finite length GV bound.

A. Special Cases of the DS Lower Bound for
In the literature, the most well-known lower bound to M * n (d) is perhaps the GV lower bound [7]. Under a finite code alphabet with |X | = Q and the Hamming distance measure, the bound states that When Q is a prime power and the code is linear with respect to modulo-Q addition and modulo-Q multiplication, G n (d) has been improved [7] to The next two examples show that both G n (d) and G lin n (d) are special cases of the DS lower bound in Corollary 1.
Example 1: Consider a finite code alphabet X with |X | = Q and the Hamming distance measure. Let the components of X = (X 1 X 2 . . . X n ) be i.i.d. random variables with a common distribution P X . This choice yields according to the Cauchy-Schwarz inequality. Since (17) is monotonically increasing in β, the best DS bound L X (d) is to take the smallest β, i.e., β = 1/Q. Hence, X is uniformly distributed over X n and which is exactly G n (d) in (15). ✷ Example 2: Continue from Example 1 and consider a binary code alphabet X = {0, 1}.
Instead of taking X to be uniformly distributed over X n as was done in the previous example, we now let it be uniformly distributed over a linear space spanned by k linearly independent Hence, X can be written as where w H (x) is the Hamming weight of vector x, and U i B i ⊕ B i . Accordingly, from which we can see that the distance spectrum is reduced to the probability that the Hamming weights of linear codewords are less than d. Accordingly, As anticipated, L X (d) = 2 k .
As a specific example, when taking {z 1 , z 2 } = {0011, 1100} for n = 4, we have d = 2 and This yields L X (2) = G lin 4 (2) = 4, which is an improvement over G 4 (2) = 16/5 = 3.2 from (18). ✷ The two examples indicate that the finite length GV lower bounds could be regarded as special cases of the DS lower bound by adopting different uniform X. While (15) is based on uniform X over the entire X n , the improved finite length GV lower bound in (16) can be obtained by adopting uniform X over a linear subspace of X n . The former observation in (15) has actually been stated by Kolesnik and Krachkovsky in [5, pp. 1446]. As we will see in the next subsection, this may help determining distributions P X which yield tighter bounds on M * n (d).

B. Uniform Distribution (UD) Lower Bounds
Theorem 1 implies the sufficiency of employing uniform X (over an appropriate subset of The examples in previous subsection suggest that adopting uniform X may result in a good lower bound. Thus, L X (d) based on uniform X forms an important subclass of the DS lower bounds, which is referred to as the uniform distribution (UD) lower bound for convenience. The examples below for specific distance measures confirm the significance of the UD lower bound.
Example 3: For a finite code alphabet X with |X | = Q and a so-called "probability-oferror" distance measure [16], defined as which can be validated by where X is the uniform distribution on X n . ✷ Example 4: Let the distance measure be given by wherex and x are in {0, 1} n , and κ n (x) x n 2 n−1 + x n−1 2 n−2 + ... + x 2 2 1 + x 1 is the binary representation of x = (x 1 x 2 . . . x n ). By its definition, µ(x, x) is the absolute difference between two decimal numbers κ n (x) and κ n (x), and is a separable distance measure [16, Def. 1].
Since κ n (x) is an integer in {0, 1, 2, . . . , 2 n − 1}, it can be easily seen that for d > 0, where ⌈·⌉ is the ceiling function. One of those uniform X's that equate L X (d) with M * n (d) has support {0, ⌈d⌉, 2⌈d⌉, . . . , (M * n (d) − 1)⌈d⌉}, and there are exactly ⌈d⌉ optimizers for M * n (d) = sup X (1/F X (d)). We then recall that Example 1 has illustrated that G n (d) can be regarded as a special case of the DS lower bound under uniform X over X n . Along this perspective, the "GV lower bound" for this separable distance measure is given by As a consequence, G n (d) is strictly less than M * n (d) except when ⌈d⌉ = 1 and ⌈d⌉ ≥ 2 n . This result confirms that the finite length GV lower bound is not tight in general.
In this special example, an upper bound U n (d) for M * n (d) can also be provided based on Theorem 1. If there exists U n (d) such that U n (d) ≥ 1/F X (d) for all X's, then Now setting j = j(n, d) 2 n /⌈d⌉, we derive where (19) and (20) hold becauseX and X are i.i.d., and (21) follows from the Cauchy-Schwarz inequality (similar to (6)). This immediately gives which is precisely M * n (d). ✷

C. Algorithmic Construction of Block Codes Subject to a Given Minimum Distance Criterion
The direct proof (i.e., step 2) of Theorem 1 suggests an iterative procedure to remove elements from an initial support of size larger than M * n (d). Between any two elements selected from the current support, the one with a larger q-value (cf. (8)) should be removed, as implied from (13). Through such an iterative procedure, it can be anticipated that and the DS lower bound is iteratively improved until either M * n (d) is reached or no further improvement can be obtained. An experiment on this algorithmic construction is thus conducted via the procedures outlined below. We observe that a number of block codes of sizes exceeding the GV lower bound are obtained.

Algorithm 1 (Algorithmic lower bound
Step 1. Set i = 1. Let P X (i) be the uniform distribution over (finite) support X n . Initialize for every z and x in X n , η(z, x) = 1{µ(z, x) < d}, which indicates whether or not one of z and x needs to be removed. Record the number of such pairs as Step 2. If Γ (i) = 0, then output A n (d) = |S(X (i) )| and stop the algorithm.
Step 3. Select arbitrarily two distinct elements z (i) and x (i) from S(X (i) ), satisfying and q x (i) x∈X n P X (i) (x) · η(x (i) , x).
Step 4. Construct P X (i+1) as and update Step 5. Set i = i + 1 and go to Step 2.
For clarity, an example is provided to demonstrate the algorithmic establishment of A n (d).
Example 5: Let µ(·, ·) be the Hamming distance measure, and initialize P X (1) as the uniform distribution with support {0, 1} 3 . Then, a lower bound A 3 (2) to M * 3 (2) can be established through the proposed algorithmic procedure as follows.
Because q 000 and q 001 are equal, we just choose to keep 001 with its probability being updated as 1/8 + 1/8 = 1/4, and remove 000 by letting the corresponding probability mass be zero.
So, we keep 111, update P X (5) (111) as 1/4, and remove 110. By leveraging duality and Theorem 1, this quantity is characterized in the following theorem.
Theorem 2: For all n ≥ 1 and M ≥ 2, where ℜ + is the set of positive real numbers.
Similar to the remarks made at the end of Section II, if µ(·, ·) is not a non-negative distance measure, the range of a in (22) should be changed to (µ 0 , ∞) (where µ 0 is defined in (14)).
An immediate consequence of Theorem 2 is that any distribution P X provides a lower bound to d * n (M). We summarize this observation in the next corollary.

Corollary 2:
for any n-letter distribution P X .
Based on this corollary, lower bounds to d * n (M) for binary block codes have been computed for n = 6, 7, 8, 11 and M = 2 k for any integer 1 ≤ k < n; this was presented in the conference version of the present work [32, Fig. 1]. In contrast to the algorithmic construction in Section III-C, the distributions employed therein are uniformly distributed over a k-dimensional linear subspace with k properly chosen length-n vectors as their bases, which satisfy the condition that all pairwise distances among these k basis vectors are no less than d (as was similarly done in Example 2). The results confirm, again, the lower bound obtained via Corollary 2 can considerably outperform the finite length GV lower bound.

A. The Maximal Code Rate and the Largest Minimum Relative Distance
In this section, we extend the results in Theorems 1 and 2 to a form that is amenable to asymptotic analyses. Distance spectrum-based formulae for the relative distance δ = d/n and the code rate R = log(M)/n are provided in the next theorem. We adopt the natural logarithm in the following. In order to reflect clearly its dependence on block length n, X n rather than X will be used to denote an n-tuple of random variables taking values in X n .
The formula in (24) provides a different form compared to the general asymptotic formula of δ * n (R) from [16], which we recapitulate here to facilitate subsequent comparison: 3 where X {X n } ∞ n=1 and ℜ denotes the set of real numbers. When compared to (25) and (26), the new way to characterize the largest minimum relative distance in (24) exhibits the following advantages.
1) The cumbersome exponent of e nR in (25) and (26) no longer exists.
2) The task of taking the supremum over all distributions is now done before the functional optimization of parameter a over the real line. This reduces the set of all distributions to be optimized over from being of infinite dimension for X to being of finite dimension for X n , and hence may help in evaluating R * n (δ) via the numerical characterization of the optimal distribution P * X n , particularly when n is small. 3) With the exact formula for every n, we exclude the necessity of sandwiching the quantities of (25) and (26), and remove the mathematical peculiarity, such as equality holding for all rates R, except for (at most) countably many points.

B. Bounds for the largest code rate attainable for an (n, M, nδ)-code
The formula of R * n (δ) in Theorem 3 provides a new quantitative characterization of the largest code rate attainable for an (n, M, nδ)-code. As a result, R * n (δ) can be lower-bounded by the large deviation rate function of (1/n)µ(X n , X n ).

Theorem 4:
The largest code rate R * n (δ) and the largest minimum relative distance δ * n (R), are respectively lower bounded as and δ * n (R) ≥ max a ∈ ℜ + : sup where J X n (δ) inf a≤δ I X n (a), and ϕ X n (θ) 1 n log E e θµ(X n ,X n ) .
Proof: First, we note from basic properties 4 of the large deviation rate function [29] that Therefore, it suffices to prove (27) under the condition that Let Y nδ − µ(X n , X n ) and note that E[Y ] < 0. It can then be derived from Markov's inequality that for θ > 0, By using the fact that and the convexity of M Y (θ) over θ ∈ ℜ, we obtain This completes the proof of (27). The bound in (28) can be easily obtained by duality.
When the distance measure µ(·, ·) is additive and the components of X n are i.i.d. with generic distribution P X , I X n (a) exhibits a single-letter expression for all block lengths n as: where ϕ X (θ) = log E e θµ(X,X) .
This reduces (27) to Further assuming that X = {α 1 , α 2 , . . . , α Q } is finite and µ(·, ·) is the Hamming distance measure, we obtain and is the binary Kullback-Leibler divergence [33]. This leads to The same argument used to verify (6) implies that 0 ≤ b X ≤ (Q − 1)/Q. Consequently, for and for 0 < R < log Q, As a result, when we take X n = (X 1 X 2 . . . X n ) to have i.i.d. components, (27) is reduced to (31) under the Hamming distance measure, which yields exactly the asymptotic GV lower bound [7], [34]. This suggests that the asymptotic GV lower bound may be improved with more general X n , particularly X n with dependent components.
We then have the following result, which strengthens (31) and (32) by the addition of a logarithmic term.
Before we prove this result, we remark that Jiang and Vardy [35,Thm. 1] proved an asymptotic improvement to the GV lower bound. They showed, by using a graph-theoretic framework, that the achievable second-order term in (33) is at least (log n)/n. This is slightly stronger than our result in (33) because we showed that the achievable secondorder term in (33) is at least (log n)/(2n) but our methods, based primarily on information spectrum analysis [18], are completely different from Jiang and Vardy's [35] and provide some additional insight into the suboptimality of choosing X n with i.i.d. components since our evaluation of the relevant distance spectrum is asymptotically tight. See (35) to follow.
Proof: Per the proof of Theorem 4, to lower bound R * n (δ), it suffices to provide an upper bound on Pr[µ(X n , X n ) < nδ]. However, instead of using Markov's inequality in (30), we will evaluate this using exact asymptotics. Choose X n (and alsoX n ) to also be i.i.d. with generic distribution Bernoulli(1/2). Then it is clear that K i µ(X i , X i ) is also a sequence of i.i.d. Bernoulli(1/2) random variables. Since the span 5 of K i is 1, by the Bahadur-Rao theorem [36] (see also [28,Thm. 3.7.4]) applied to lattice random variables, we know that where f n ∼ g n means that lim n→∞ f n /g n = 1. See [37, Eqn. (5.41)] for a detailed derivation of (35). Uniting (23) and (35), we obtain R * n (δ) ≥ D δ which is (33).
For (34), we simply have to "invert" the former result carefully. From (24) and (35), we note that for n sufficiently large, where is some continuous function of a. Let {a * n } n∈N be the sequence of maximizers in (36). Then because 0 < R < log 2, the sequence {a * n } n∈N converges to some limit a = δ GV (R) ∈ (0, 1/2); so, g(a) ∈ (0, ∞). Due to the additional term −(1/2) log n in the exponent in (36), we have The final result in (34) follows by Taylor expanding δ GV (·).
In parallel to Theorems 4 and 5, we also provide an upper bound to R * n (δ) with the aid of the so-called "twisted distributions" technique [38]. In the following, we again do not impose any structure (e.g., Hamming or additive) on the distance measure µ(·, ·).
Given that θ * is the maximizer to achieve I X n (δ) = sup θ∈ℜ {aθ − ϕ X n (θ)} with P X n being the optimizer of sup X n J X n (δ), we have that under (29), provided that µ(X n , X n ) has no point mass at nδ, i.e., Pr[µ(X n , X n ) = nδ] = 0.
We conclude from all the above derivations that In fact, the preceding theorems give us much more, assuming the conditions for them are satisfied. Theorem 5 refines the lower bound to R * n (δ) in Theorem 4 yielding R * n (δ) ≥ max sup X n J X n (δ), D δ 1 2 + log n 2n + Θ 1 n .

VI. CONCLUSION AND OPEN PROBLEMS
This paper provides an exact formula for maximal size of a code with pairwise distance d. As a result, we have made progress towards resolving a long-standing open problem in information and coding theory, albeit the formula may be difficult to compute for large block lengths. Future work includes: 1) Section III-C presents an algorithmic construction of block codes subject to a given minimum distance. It would be fruitful to conduct systematic theoretical studies to examine properties of the resultant code produced by this algorithm.
2) In [40], the covering analogue of the minimum distance (packing) problem was considered from an information spectrum perspective. Naturally, one would wonder whether it is possible to prove a similar finite length formula for the covering analogue.