Heuristics for Symmetric Rectilinear Matrix Partitioning
Abdurrahman Yaşar Ümit V. Çatalyürek








Heuristics for Symmetric Rectilinear Matrix Partitioning
    
  Abdurrahman Yaşar
  School of Computational Science and Engineering, Georgia Institute
  of Technology, Atlanta, GA 30332
  (mailto:ayasar@gatech.eduayasar@gatech.edu,
   mailto:umit@gatech.eduumit@gatech.edu)
Ümit V. Çatalyürek[1]

    
===============================================================================================================================================================================================================================================



[R]Copyright © 20XX by SIAM

Unauthorized reproduction of this article is prohibited

 =9pt
Partitioning sparse matrices and graphs is a common and important problem
that arises in
many scientific and graph analytics applications. In this work, we are
concerned with a spatial partitioning called rectilinear partitioning
(also known as generalized block
distribution) of sparse matrices, which is necessary for tiled (or blocked) execution of sparse matrix and graph analytics kernels. More
specifically, in this work, we address the problem of symmetric rectilinear
partitioning of square matrices. By symmetric, we mean having the same
partition on rows and columns of the matrix, yielding a special tiling where
the diagonal tiles (blocks) will be squares. We propose five heuristics to
solve two different variants of this problem, and present a thorough
experimental evaluation showing the effectiveness of the proposed
algorithms.




§ INTRODUCTION



After advances in the social networks and the rise of interactions on the web,
we are witnessing an enormous growth in the volume of generated data. A large
portion of this data remains sparse and irregular. Graphs and sparse matrices
are used to store and analyze an important portion of this data. However,
analyzing data stored in that kind of irregular data structures is becoming
more and more challenging, especially for traditional architectures due to the
growing size of these irregular problems.
The sheer size of the problems necessitates parallel execution. There
have been many studies developing parallel graph and sparse matrix algorithms for
shared and distributed memory systems as well as GPUs and hybrid systems.
Effective data and computation partitioning is the
first step to propose efficient portable (parallel) algorithms <cit.>.

Two-dimensional matrix partitioning is a hard problem and has been used in
dense linear algebra <cit.> for a long time.
Checkerboard partitioning, where the partitioned matrix maps naturally onto a
2D mesh of processors, is widely used in earlier two-dimensional matrix
partitioning <cit.>. Checkerboard partitioning
is highly useful to limit the total
number of messages on distributed settings. However, these works are suited to
dense or well structured sparse matrices.

In the context of this paper, we focus on spatial, two-dimensional
checkerboard-like partitioning problem that we call Symmetric
Rectilinear Partitioning. Here, we assume that given matrix is square and
we would like to partition that matrix into p × p tiles such that
by definition diagonal blocks will be squares. This type of partitioning is
very convenient if one wants to gather information along the rows/columns
and distribute along columns/rows. Also, in the context of graphs, diagonal
tiles can be visualized as sub-graphs and any other tile represents the edges
between two sub-graphs. This type of partitioning becomes highly
useful to reason about graph algorithms.
For example, in a concurrent work, we have leveraged the
symmetric rectilinear partitioning for developing a block-based triangle counting
algorithm <cit.> that reduces data movement during, both sequential and parallel,
execution and also naturally suitable for heterogeneous architectures.

In this work, we define two variants of the symmetric rectilinear partitioning
problem and we propose refinement based and probe-based partitioning heuristics
to solve these problems. Refinement based
heuristics <cit.> apply a dimension reduction
technique to map the two-dimensional problem into one-dimension and compute partition
vector on one-dimensional data by running an optimal partitioning
algorithm <cit.>. Probe-based algorithms compute the partitioning
vector by seeking for the best cut for each point. The novelty of the proposed
approaches is to use the natural order of the matrix instead of running expensive
hypergraph models or graph partitioning algorithms to order vertices. We combine
lightweight spatial partitioning techniques with simple heuristics.

Contributions of this work are as follows:



  * We propose heuristics for symmetric rectilinear partitioning
problem that does not require row orderings.


  * We evaluate the effect of the simple, degree-based and
RCM based vertex orderings on the tile distributions.


  * We experimentally evaluate the performances of proposed
algorithms wrt. state-of-the-art algorithms in different settings.


Our experimental results show that our proposed algorithms are
very effective in finding symmetric rectilinear partitions. In all instances,
our algorithms produce similar or better load-imbalanced solutions than
Nicol's <cit.> rectilinear partitioning algorithm, which has more
freedom in choosing row and column partitions.



§ PROBLEM DEFINITION





In this paper, we are concerned with partitioning sparse matrices.
In particular, we are interested in partitioning adjacency matrix
representation of graphs. A directed graph G=(V, E), consists of a
set of vertices V and a set of edges E. A directed edge e is
referred as e = (u, v) ∈ E, where u, v ∈ V, and u is called
the source of the edge and v is called the target. The neighbor list
of a vertex u ∈ V is defined as N[u] = {v∈ V|(u, v)∈ E}.
We will use n and m for number of vertices and edges, respectively,
i.e., n=|V| and m=|E|. Let A_G be the adjacency matrix representation
of the graph G, where all edges are represented with nonzeros, and the
rest of the entries will be zero. That is, A_G is an n× n matrix,
where ∀ (u, v) ∈ E, A_G[u,v] =1, and everything else will be 0.
Without loss of generality, we will assume source vertices are represented
as rows, and target vertices represented as columns. In other words,
elements of N[u] will correspond to column indices of nonzero elements
in row u. We will also simply refer to matrix A_G as A when G is
clear in the context. <Ref> lists the notations used in
this paper.





Given adjacency matrix of G, A_G ∈𝔹^n × n and an
integer p, 1 ≤ p ≤ n. Let C be a partition vector
that consists of sequence of p+1 integers such that
0=c_0<c_1<…<c_p=n. Then C defines a
partition of [0, n] into p intervals
[c_i, c_i+1-1] for 0 ≤ i ≤ p-1.

Rectilinear Partitioning.
Given A, and two integers, p and q, a rectilinear partitioning
consists of a partition of [0, n] into p intervals
(C_c, for columns) and into q intervals (C_r, for rows) such that A is partitioned into
non-overlapping p × q contiguous tiles.


In rectilinear partitioning, a column partition vector, C_c, and a row partition vector, C_r,
together generate p × q tiles.
For i ∈ [0,p] and j∈[0,q], we denote (i,j)-th tile by T_i,j and
|T_i,j| denotes the number of nonzeros in T_i,j.
Given partition vectors, quality of partitioning can be defined using
load imbalance, λ, among the tiles, which is computed as

    λ(A, C_c, C_r) = L_max/L_avg - 1

where

    L_max = max_0≤ i,j ≤ n{ |T_i,j| }
 and

    L_avg=∑_0 ≤ i,j ≤ n|T_i,j|/p × q =
    2× m/p × q.

A solution which is perfectly balanced achieves a load imbalance,
λ, of 0. <Ref> presents a toy example for
rectilinear partitioning where C_c={0,5,8,10}, C_r={0,4,8,10}
and λ(A, C_c, C_r) = 5/3.6-1 = 0.39.


Symmetric Rectilinear Partitioning.
Given A and p, symmetric rectilinear partitioning
can be defined as partitioning [0, n] into p intervals such that
A partitions into p × p non-overlapping contiguous tiles where
diagonal tiles are squares.



In symmetric rectilinear partitioning, same partition vector, C_c=C_r,
is used for row and column partitioning.
<Ref> presents a toy example for the symmetric rectilinear
partitioning  where C_c= C_r = {0,5,8,10}
and λ(A, C_c, C_r) = 5/3.6-1 = 0.39

In the context of this work, we consider two symmetric rectilinear
partitioning problems. The first problem (minLoadImbal, or in short),
consists in finding the optimal partition vector, C, that minimizes the load
imbalance, for a given the matrix A and an integer p:


    ( A, p )    = _Cλ(A, C, C)


The second problem (minNumCuts, or in short), is the dual
of the first problem. For a given a matrix, A, and an integer, Z,
this problem consists of finding the minimum number of intervals, p,
that will partition the matrix A, where the number of nonzeros in all
tiles is bounded by Z.


    ( A, Z )    = _p∈ [1, n]λ(A, ( A, p ), ( A, p )) 
       ≤Z× p^2/ m





§ RELATED WORK



Two-dimensional matrix distributions have been widely used in dense
linear algebra. Most of these distributions are
cartesian <cit.>; where the same partitioning
vector is used for row and column partition. For sparse and irregular
problems finding a good vector that can be aligned with both dimensions
is even harder. Therefore, many non-cartesian two-dimensional matrix
partitioning methods have been
proposed <cit.> for
sparse and irregular problems.

One way to overcome the hardness of proposing one partition vector for rows and
columns is to propose a partition vector for rows and columns.
This problem is named as rectilinear partitioning <cit.>
(or generalized block distribution <cit.>).
Independently, Nicol <cit.>
and Manne and Sørevik <cit.> proposed an algorithm to solve this
problem that is based on iteratively improving a given solution by alternating
between row and column partitioning.
These algorithms
transform two dimensional (2D) rectilinear partitioning problem into
one-dimensional (1D) partitioning problem using a heuristic
and iteratively improves the solution.
The one-dimensional partitioning problem is built by
setting the load of an interval of the problem as the maximum of the load of
the interval inside each stripe of the fixed dimension. This refinement
technique is presented in Algorithm <ref>.
Here, optimal1DPartition(P) is a function that returns the optimal 1D partition
(which is also known as chains-on-chains partitioning) <cit.>.
Hence, Algorithm <ref> returns the optimal 1D row partition for
the given column partition C_c.


[t]

  *[l]Array to store max of interval sums for each vertex
  P[i] = 0, for 0 ≤ i ≤ n+1  

  *[l]for each row
  i=0 to n-1
    *[l]Array to store interval sums
    nnz[k] = 0, for 0 ≤ k ≤ p  
    each j, where A_G[i,j]=1
      k ← 0 *[l]Interval index
      j ≥ C[k+1]
        k ← k+1 *[l]Find the interval
      
      nnz[k] ← nnz[k]+1  
    
    P[i+1] ←max_k{nnz[k]} 
  

  *[l]Compute prefix sum
  i=1 to n
    P[i] ← P[i]+P[i-1]  
  

  *[l]Return the output of 1D partitioning
  optimal1DPartition(P)
  refinement(A_G, C_c, p)
  




Computing the optimal solution was shown to be NP-hard by Grigni and
Manne <cit.>. In fact, their proof shows that the
problem is NP-hard to approximate within any factor less than 2. Khanna et
al. <cit.> have shown the problem to be constant-factor
approximate.

Rectilinear partitioning may still cause high load-imbalance due to
generalization. Jagged partitions (or Semi Generalized Block
Distribution <cit.>) tries to overcome this problem by
distinguishing between the main dimension and the auxiliary dimension.
The main dimension is split into p intervals and each of these intervals
partition into
q rectangles in the auxiliary dimension.
Each rectangle of the solution must have its main dimension matching one of
these intervals. The auxiliary dimension of each rectangle is arbitrary.
Saule et al. <cit.> present multiple variants and
generalization of jagged partitioning.



§ SYMMETRIC RECTILINEAR PARTITIONING



We propose five different algorithms for two variants of the symmetric rectilinear
partitioning problem. These algorithms can be classified as refinement
based and probe-based. In this section, we explain how these algorithms
are designed.



 §.§ Heuristics for the minLoadImbal problem.


We propose three algorithms for problem. Two of these algorithms,
Pick best direction (first) (PBD) and
Pick best (in each) iteration (PBI),
leverages previously defined refinement technique (see <Ref>)
into the symmetric rectilinear partitioning problem.
Note that these two algorithms have no convergence guarantee.
The third algorithm, Probe target cut (PTC),
implements another heuristic and probes the minimum load imbalance
by moving in partition point in the diagonal of the matrix.



  §.§.§ Pick best direction (first) (PBD)



algorithm first applies row-based and column-based refinement independently
to simply find the optimal 1D row and column partitions and chooses the one
that gives the best load imbalance. Then, iteratively applies the refinement
algorithm only in this direction until it reaches the iteration limit (τ)
or partition vector does not change significantly (computed using 2-norm).
This procedure is presented in Algorithm <ref>.

The primary advantage of this algorithm is its simplicity. This algorithm can be
easily parallelizable like <cit.>.
However, choosing a direction at the beginning may cause missing information that
can be gathered from the other direction and the solution may converge to a
local optimum quickly.

[ht]

  *[l]Current (C) and previous (C') partition vectors
  C[0] = 0; C[j] = n, for 1 ≤ j ≤ p+1  
  C'[j] = n, for 0 ≤ j ≤ p+1  

  *[l]Apply 1D partitioning refinement
  C_r ← refinement(A, C, p) *[r]Row based
  C_c ← refinement(A^T, C, p) *[r]Column based

  *[l]Aligning same partition vector for rows and columns
  L_r←λ(A, C_r, C_r ) *[r]Row based imbalance
  L_c←λ(A, C_c, C_c ) *[r]Column based imbalance
   L_r<L_c 
      C ← C_r  
    
      C ← C_c  
    

  i ← 0  
   i<τ and ||C-C'||_2 > ε
    C' ← C 
    *[l]Always pick the best initial direction
     L_r<L_c 
      C ← refinement(A, C, p)  
    
      C ← refinement(A^T, C, p)  
    
    i ← i+1  
  

  C
  PBD(A, p)
  





  §.§.§ Pick best (in each) iteration (PBI)



algorithm applies refinement on both row-based and column partitions and always
chooses the best partition vector for the next iteration.
If that partition vector improves the current best solution, PBI algorithm
updates the partition vector that stores the best solution achieved.
This procedure is presented in Algorithm <ref>.

Tracking the load imbalance creates an opportunity to output a
better partitioning. However, tracking comes up with a computational cost (O(m)).
Hence, PBI algorithm is more expensive than PBD algorithm.

[ht]

  *[l]Initialize partition vectors
  C_c[0] = 0; C_c[j] = n, for 1 ≤ j ≤ p+1 *[r]Col. based
  C_r[0] = 0; C_r[j] = n, for 1 ≤ j ≤ p+1 *[r]Row based
  C_b[0] = 0; C_b[j] = n, for 1 ≤ j ≤ p+1 *[r]Best

  i ← 0  

   i<τ 
    *[l]Compute row based and column based partition vectors.
    C_r ← refinement(A, C_c, p)  
    C_c ← refinement(A^T, C_r, p)  

    *[l]Pick the best partition vector
    λ(A, C_c, C_c ) < λ(A, C_r, C_r )
      C_r ← C_c  
    
      C_c ← C_r  
    

    *[l]Update the best partition vector, if improved
    λ(A, C_c, C_c) < λ(A, C_b, C_b)
      C_b ← C_c  
    

    i ← i+1  
  

  C_b
  PBI(A, p)
  





  §.§.§ Probe target cut (PTC)



probes for the largest possible cut point in each step using a
two-dimensional probe algorithm. In refinement based algorithms (PBI and PBD)
mapping the problem from two-dimensional space to one-dimension
ends up with losing information. PTC algorithm tries to overcome this
problem by seeking better cut points in two-dimensional space. In high-level view
PTC algorithm is inspired from Nicol's <cit.> probe-based
one-dimensional partitioning algorithm however PTC algorithm operates on
two-dimensional space on diagonal direction in a greedy manner and due to
non-convex structure of the problem it doesn't guarantee an optimal solution.
Algorithm <ref> presents the
two-dimensional probe algorithm. The elements of C are found
through binary search, β, on the matrix.
In this algorithm, β(A, C, i, ℓ), searches A in the range [C[i-1],n]
to compute the largest cut point, C[i]=j such that λ(A, C, i) ≤ℓ and for
C[i]=j+1, λ(A,C,i)>ℓ. Algorithm <ref> returns true if at the end
partition vector has p intervals and the load imbalance is less than ℓ.
In each step, i, PTC algorithm seeks for the largest cut point in the range
[C[i-1],n] for which Probe Algorithm <ref> returns true. PTC stores the
load imbalances of these steps in an array, B. After p-1 steps, PTC algorithm
computes the minimum of B and constructs partition vector C using a binary
search-based approach similar to Algorithm <ref>.
This procedure is presented in Algorithm <ref>.

[ht]

  *[l]Initialize partition vector
  C[0] = 0; C[p] = n  

  i=1 to p-1
    C[i] ←β(A, C, i, ℓ) 
  

  λ(A, C, C) ≤ℓ
    true
  
    false
  

  probe(A, p, ℓ)
  



Note that, probe algorithm can return false even if there
exists a valid partitioning for a given target load imbalance ℓ. Because PTC algorithm
always searches for the largest cut point and that may cause missing the optimal solution.
Hence, this algorithm can also be stuck in a local optimum. However,
PTC algorithm considers more cases in a two-dimensional fashion. Therefore,
PTC is expected to produce a better partitioning than PBD and PBI algorithms.
The major disadvantage of this algorithm is its computational complexity.

[ht]

  *[l]Initialize temporary partition vector
  C[0] = 0; C[j] = n, for 1 ≤ j ≤ p+1  

  *[l]An array to store load imbalances
  B[j] = 0, for 0 ≤ j < p  

  i=1 to p-1
    l ← C[i-1]  
    r ← n  
    *[l]Probe in binary search fashion
    l<r
      m ← (l+r)/2  
      C[i] ← m  
      ℓ←λ(A, C, C, i)  
      probe(A, p, ℓ)
        r ← m  
        B[i-1] ←ℓ  
      
        l ← m+1  
      
    
    C[i] ← r  
  

  *[l]Find the minimum load imbalance
  B_min←min_j{B[j]} 

  
  C[0] = 0; C[p] = n  

  *[l]Construct partition vector
  i=1 to p
    C[i] ←β(A, C, i, B_min) 
  

  C
  PTC(A, p)
  






 §.§ Algorithms for the minNumCuts problem


Given a matrix, A, and an integer, Z, problem aims to output
a partition vector, C, with the minimum number of intervals, p, where the maximum
load of a tile in the corresponding partitioning is less then Z, i.e.,
max_0 ≤ i,j ≤ p{T_i,j}≤ Z.
We propose two algorithms to solve this problem. These algorithms are
variations of the ones that we propose for problem.



  §.§.§ Bound target load (BTL)



algorithm, displayed in Algorithm <ref>, searches for the minimal p
intervals using a binary search-based approach (Algorithm <ref>).
In the worst-case, the maximum loaded tile may be filled
with nonzeros. Therefore, the initial upper bound for the search
space can be defined as p∈[1,n/√(Z)].
BTL algorithm first reduces this upper bound by using cheap uniform
partitioning.
Then using this reduced search space, BTL algorithm searches for the minimal p
intervals using PBD (or PBI) algorithm in the same fashion.

[ht]

  *[l]u is the upper bound for number of intervals
  *[l]Initial upper bound: Tile is fully connected
  u ←n/√(Z)  

  *[l]Decreasing upper bound using a cheap algorithm
  u ← findUpperBound(A, Z, 1, u, UNI) 

  *[l]Using PBD (PBI) search for a lower upper bound
  u ← findUpperBound(A, Z, 1, u, PBD) 

  PBD(A, u)
  BTL(A, Z)
  



[ht]


  l<r
    p ← (l+r)/2  
    C ←  f(A, p)  
    ℓ←Z × p^2/m  
    ℓ' ←λ(A, C, C)  
    ℓ'<ℓ
      r ← p  
    
      l ← p+1  
    
  

  r
  findUpperBound(A, Z, l, r, f()̇)
  






  §.§.§ Probe target load (PTL)



algorithm slightly modifies two-dimensional probe algorithm. In each step, k, PTL
algorithm searches for the largest cut point, C[k] = l in the range [C[k-1], n]
that assures all of the tiles until that cut point have less than Z
number of nonzeros, ∀_i,j≤ k |T_i,j|≤ Z.
This procedure is described in Algorithm <ref>.

[ht]

  *[l]Initially we don't know partition vector's size
  C[0] = 0  
  i ← 1 
  C[i-1]≠ n
    C[i] ←β(A, C, i, Z) 
    i ← i+1  
  

  C

  PTL(A, Z)
  





 §.§ Complexity Analysis.



<Ref> displays the computational complexity of the algorithms we have
used in this work.
Each iteration of the iterative refinement algorithm <cit.>
(Algorithm <ref>) has a worst-case complexity of
O(q(plogn/p)^2 + p(qlog n/q)^2) <cit.> for
unsymmetric rectilinear partitioning. The Algorithm is guaranteed to converge with
at most n^2 iterations. However, as noted in these earlier work,
in our experiments we observed that algorithm converges very quickly, and hence
for the sake of fairness we have decided to use the same limit on the number of iterations, τ.
For the symmetric case, where p=q, refinement algorithm runs in
O(p^3(logn/p)^2). and this what we displayed in <Ref>.

PBD algorithm first runs Algorithm <ref> and then computes
the load imbalance. These operations can be computed in
O(p^3(logn/p)^2) and in O(m) respectively. In the worst-case,
Algorithm <ref> is called τ times. Hence, PBD algorithm runs in
O(m + τ p^3(logn/p)^2).

PBI algorithm runs Algorithm <ref> and then computes
the load imbalance in each of the τ iterations.
Hence, PBI algorithm runs in
O( τ(m+p^3(logn/p)^2) ).

Probe algorithm (Algorithm <ref>) does O(m log n) computations in
the worst-case to find a cut point. Since there are O(p) cut points
probe algorithm runs in
O(p mlog n). PTC algorithm (Algorithm <ref>) calls Probe algorithm
at most log n times for each i ∈ [1,p-1]. Hence, PTC algorithm runs
in O(m(plog n)^2).

BTL algorithm initially tries to reduce search space using cheap uniform
partitioning, however in the worst-case uniform partitioning may not be able
to reduce search space. Hence, in the worst-case BTL algorithm runs
logn/√(Z) times uniform and PBD (or PBI) algorithm.
In the worst-case p=n/√(Z)  (i.e., when there is fully dense tile).
So, BTL algorithm runs in
O( logn/√(Z) (m + τn/√(Z)^3(log√(Z))^2) )
when PBD is used as the secondary algorithm and in
O( logn/√(Z)τ ( m+ n/√(Z)^3( log√(Z))^2) )
when PBI is used as the secondary algorithm.

PTL algorithm has the same complexity with Algorithm <ref> and in the
worst-case p=n/√(Z).
Hence, PTL algorithm runs in; O( n/√(Z) mlog n ).

In addition to our proposed algorithms, we use two reference partitioning
algorithms as baselines. First,
NIC refers to Nicol's rectilinear partitioning algorithm <cit.>.
This algorithm outputs a partition vector for each dimension, hence it does not
output symmetric partitioning. The complexity of this algorithm is provided
in <cit.> and listed in <Ref>.
Second, UNI refers to uniform partitioning. This is the simplest
checker board partitioning, where each tile has an equal number of
rows and columns. UNI algorithm runs in constant time.










§ EXPERIMENTAL EVALUATION


The partitioning algorithms presented for and problems are implemented
in the C++ programming language and compiled with GCC version 7.2. The experiments
are conducted on a server that has four, 14-core Intel Xeon E7-4850 2.20GHz
processors, 2TB of memory, 1TB disk space, running Ubuntu GNU/Linux with
kernel 4.8.0.

We have performed an extensive evaluation of the proposed algorithms on 16
different real-world and 4 different synthetic (RMAT) graphs coming
from SuiteSparse Matrix Collection (formerly known as UFL) <cit.>,
SNAP [SNAP Datasets: <http://snap.stanford.edu/data>], and
DARPA Graph Challange [GraphChallenge Datasets: <https://graphchallenge.mit.edu/data-sets>].
Properties of these datasets, along with the load balances found
by different algorithms on different layouts of these matrices are presented in <Ref>.
In the experiments, we used τ=20 and ϵ=0.0001.

In our experiments, we included three different vertex ordering techniques
before giving the adjacency matrix to respective partitioner:
the natural order of the vertices (NAT), degree-based ordering of the vertices (DEG)
and Reverse Cuthill McKee (RCM) based ordering of the vertices.
<Ref> illustrates these orderings along
with the adjacency matrix representations for a toy graph.





 §.§ Load imbalance evaluation.


We first evaluate the proposed algorithms for the problem.
Table <ref> reports load imbalances of
five different algorithms; NIC, UNI, PBD, PBI, and PTC,
on three different vertex ordering techniques, NAT, DEG and RCM for each graph.
In this experiment we chose p=q=8. Hence, every graph is partitioned into
8 × 8 tiles.
In <Ref>, the best load imbalance for a graph
instance is highlighted using the green color.
As expected, PTC algorithm gives the best load imbalance in every graph instance.
The best load imbalance among with the other four algorithms; NIC, UNI, PBD,
and PBI, for each vertex ordering are highlighted using yellow color.
PBD algorithm gives the best performance among these four algorithms
in 54 of the 60 graph instances.
PBI and NIC algorithms give very similar load imbalances.

In <Ref>, the last row presents geometric means of the five
different algorithms on three different vertex orderings.
As shown in the table, the geometric means of PTC algorithm on NAT, DEG and RCM vertex
orderings are 0.7, 0.7 and 0.8, respectively.
These results show that PTC algorithm is more resistant to the vertex order,
and hence, it can produce partitioning with similar qualities.
On the other hand, we observe a significant change in load imbalance (up to 40%)
for refinement based algorithms depending on the vertex order.
As expected, uniform partitioning performs poor on DEG vertex order due
to higher density in the bottom right portion of the adjacency matrix.



 §.§ Algorithm evaluation.


We evaluate relative load imbalance performances of NIC, UNI, PBD, PBI and PTC
algorithms. The aim is to illustrate the efficiency of the proposed algorithms with
respect to NIC and UNI. In this experiment, NAT vertex ordering is used and
we choose p={2,4,8}.
<Ref> illustrates the performance profiles of the algorithms
for different p values.
In the performance profiles, we plot the number of the test
instances (y-axis) in which an algorithm obtains a load imbalance on an instance that
is no larger than x times (x-axis) the best load imbalance found by any algorithm
for that instance <cit.>. Therefore, the higher a profile at a
given x value, the better an algorithm is. We observe that in all cases
(<Ref>-<Ref>) PTC algorithm gives the best
performance in the majority of the test instances.
PBD algorithm becomes the second-best algorithm. We also observe
that NIC and PBI algorithms have almost identical profiles. Both of these
algorithms consider row-based and column-based directions in each iteration.
Therefore, probably they fall into the same local optimum and stuck.
As expected UNI algorithm performs the worst because of the irregularity of the graphs.





 §.§ Visualization of block distributions.




In <Ref>, we present density maps of NIC, UNI, PBD and PTC
based partitioning where p=8 on soc-LiveJournal1 and friendster graphs.
The same color range used in all plots and darkness of the color of a tile is
proportional to the number of nonzeros inside the tile.
In this experiment, NAT vertex ordering is used.
Percentages presented for each tile present percentage of the number of
nonzeros within a tile.
Note that in these figures tiles are mapped into a grid for a better
visualization. hence NIC algorithm outputs a rectilinear partitioning,
therefore, there can be nonzero tiles under diagonal tiles.
<Ref> visually supports our previous
findings; PTC algorithm produces better partitioning (closely colored tiles)
and UNI gives the worst (too much variance between tile colors)
and PBD performs slightly better then NIC algorithm.



 §.§ Number of cut evaluation.


In this experiment, we evaluate performances of BTL and PTL algorithms, proposed
for the problem.
<Ref> reports the number of cuts found by these algorithms
for different graphs. We choose Z=m/8 for each graph.
UNI algorithm outputs the minimum number of cuts
that can be gathered using uniform partitioning.
In <Ref> for each graph, the minimum number of cuts
is highlighted using a green color.
We observe that NAT vertex ordering has a lower number of cuts in general
due to the more uniform distribution of the nonzeros. On the other hand algorithms
output higher number of cuts using RCM vertex ordering,
because RCM creates denser regions while trying to make nonzeros appear near diagonal.
As expected UNI partitioning performs the worst.






 §.§ Execution time evaluation.



<Ref> shows the execution time behavior of six algorithms; NIC, PBI, PBD, PTC, BTL,
and PTL, on the largest five graphs in our dataset listed in <Ref>.
We observe that, for problem both BTL and PTL algorithms' execution times are similar.
For problem, for the smaller three graphs, both PBI and PBD algorithms' execution time is
about 1.5 times slower than NIC. On largest two graphs PBD is the fastest algorithm.
Given the computational complexities of the algorithms in <Ref>, as expected
PTC algorithm was the slowest of all, due to too many prefix sum lookups on the matrix.
Please note that, <Ref> assumes linear lookup for a given partition vector and
require linear space. We first developed a naive implementation with such linear lookup, which
was two to three orders of magnitude slower than than other algorithms.
To decrease the lookup complexity for PTC algorithm we used a sparse prefix sum
data structure inspired from the Fenwick tree <cit.>, that can query the matrix in
logarithmic time but requires log(n) times more space.
With that data structure PTC algorithm now 4 to 43 times slower than NIC algorithm.







§ CONCLUSION



In this work, we proposed different heuristics for symmetric
rectilinear partitioning problem and
present a thorough
experimental evaluation showing the effectiveness of the proposed
algorithms.
Even though our problem definition is more restricted, in our experiments,
we observed that our proposed algorithms give better
load-imbalance than Nicol's <cit.> state-of-the-art rectilinear
partitioning algorithm in every test instances.
PTC algorithm gives the best load imbalance in the majority of the test instances
and PBD algorithm is the second-best algorithm. PBI algorithm's performance
is almost identical to NIC.

As future work, we are working on decreasing the space requirements
of our sparse prefix sum data structure and also decreasing the number of lookups that PTC
performs. In addition, we will also investigate approximation techniques, and
parallelization of the algorithm.





§ ACKNOWLEDGE

We would like to extend our gratitude to M. Fatih Bali̇n for his valuable
comments and sharing his sparse prefix sum data structure code with us.

siam


