On Time Optimal Supernode Shape.
AbstractWith the objective of minimizing the total execution time of a parallel program on a distributed memory parallel computer, this paper discusses the selection of an optimal supernode shape of a supernode transformation (also known as tiling). We identify three parameters of a supernode transformation: supernode size, relative side lengths, and cutting hyperplane directions. For supernode transformations on algorithms with perfectly nested loops and uniform dependencies, we prove the optimality of a constant linear schedule vector and give a necessary and sufficient condition for optimal relative side lengths. We also prove that the total running time is minimized by a cutting hyperplane direction matrix from a particular subset of all valid directions and we discuss the cases where this subset is unique. The results are derived in continuous space and should be considered approximate. Our model does not include cache effects and assumes an unbounded number of available processors, the communication cost approximated by a constant, uniform dependences, and loop bounds known at compile time. A comprehensive example is discussed with an application of the results to the Jacobi algorithm.
Introduction
Supernode partitioning is a transformation technique
that groups a number of iterations in a nested
loop in order to reduce the communication startup
cost. This paper addresses the problem of selecting
optimal cutting hyperplane directions and optimal
supernode relative side lengths with the objective of
minimizing the total running time, the sum of communication
time and computation time, assuming a large
number of available processors which execute multiple
supernodes.
A problem in distributed memory parallel systems
is the communication startup cost, the time it takes
a message to reach transmission media from the mo-
This author was supported by AT&T Labs.
y This research was supported in part by the National Science
Foundation under Grant CCR-9502889 and by the Clare Boothe
Luce Professorship from Henry Luce Foundation.
ment of its initiation. The communication startup
cost is usually orders of magnitude greater than the
time to transmit a message across transmission media
or to compute data in a message. Supernode transformation
was proposed in [14], and has been studied
in [1, 2, 3, 4, 15, 18, 20, 25, 26, 27] and others to
reduce the communication startup cost. Informally,
in a supernode transformation, several iterations of
a loop are grouped into one supernode and this supernode
is assigned to a processor as a unit for ex-
ecution. The data of the iterations in the same su-
pernode, that need to be sent to another processor,
are grouped as a single message such that the number
of communication startups is reduced from the number
of iterations in a supernode to one. A supernode
transformation is characterized by the supernode size,
the relative lengths of the sides of a supernode, and
the directions of hyperplanes which slice the iteration
index space of the given algorithm into supernodes.
All the three factors affect the total running time. A
larger supernode size reduces communication startup
cost, but may delay the computation of other processors
waiting for the message and therefore, result in
a longer total running time. Also, a square supernode
may not be as good as a rectangular supernode
with the same supernode size. In this paper, selection
of optimal cutting hyperplane directions and optimal
relative side lengths is addressed.
The rest of the paper is organized as follows. Section
presents necessary definitions, assumptions,
and terminology. Section 3 discusses our results in
detail. Section 4 briefly describes related work and
the contribution of this work compared to previous
work. Section 5 concludes this paper. A bibliography
of related work is included at the end.
Basic definitions, models and assumptions
The architecture under consideration is a parallel
computer with distributed memory. Each processor
has access only to its local memory and is capable of
communicating with other processors by passing mes-
sages. In our model, the cost of sending a message
is represented by t s , the message startup time. The
computation speed of a single processor is characterized
by the time it takes to compute a single iteration
of a nested loop. This parameter is denoted by t c .
Algorithms under consideration consist of a single
nested loop with uniform dependencies [22]. Such algorithms
can be described by a pair (J; D), where J
is an iteration index space and D is an n \Theta m dependence
matrix. Each column in the dependence matrix
represents a dependence vector. The cone generated
by dependence vectors is called dependence cone. The
cone generated by the vectors orthogonal to the facets
of the dependence cone is called tiling cone. We assume
that m - n, matrix D has full rank (which is
equal to the number of loop nests n), and all elements
on the main diagonal of the Smith normal form of D
are equal to one. As discussed in [23], if the above assumptions
are not satisfied, then the iteration index
space J contains independent components and can be
partitioned into several independent sub-algorithms
with the above assumptions satisfied.
In a supernode transformation, the iteration space
is sliced by n independent families of parallel equidistant
hyperplanes. The hyperplanes partition iteration
index space into n-dimensional parallelepiped supernodes
(or tiles). Hyperplanes of one family can be
specified by a normal vector orthogonal to the hy-
perplanes. The square matrix consisting of n normal
vectors as rows is denoted by H. H is of full rank
because the n hyperplanes are assumed to be inde-
pendent. These parallelepiped supernodes can also
be described by the n linearly independent vectors
which are supernode sides. As described in [20], column
vectors of matrix are the n side vec-
tors. A supernode template T is defined as one of
the full supernodes translated to the origin 0, i.e.,
1g.
Supernode index space, J s , obtained by the supernode
transformation H is:
Supernode dependence matrix, D s 1 , resulting from
supernode transformation H consists of elements of
the set:
We use D to denote either a matrix or a set consisting of
the column vectors of matrix D. Whether it is a matrix or a
set should be clear from the context.
As discussed in [20, 26] 2 , partitioning hyperplanes defined
by matrix H have to satisfy HD - 0, i.e. each
entry in the product matrix is greater than or equal
to zero, in order to have (J s ; D s ) computable. This
implies that the cone formed by the column vectors in
E has to contain all dependence vectors in D. There-
fore, components of vectors in D s defined above are
nonnegative numbers. In our analysis, throughout the
paper, we further assume that all dependence vectors
of the original algorithm are properly contained in the
supernode template. Consequently, components of D s
are only 0 or 1, and I ' D s . This is a reasonable assumption
for real world problems [5, 26].
To present our analysis, the following additional
notations are introduced. The column vector l =
called supernode side length vector. Let
L be an n \Theta n diagonal matrix with vector l on its
diagonal and E u be a matrix with unit determinant
and column vectors in the directions of corresponding
column vectors of matrix E. Then,
components of vector l are supernode side lengths in
units of the corresponding columns of E u . We define
the cutting hyperplane direction matrix as:
The supernode size, or supernode volume, denoted by
g, is defined as the number of iterations in one su-
pernode. The supernode volume g, matrix of extreme
vectors E and the supernode side length vector l are
related as
. The relative supernode
side length vector,
sg
\Theta l;
and clearly
lengths of supernodes relative to the supernode size.
For example, if H I, the identity matrix,
and then the supernode is a square. How-
ever, if
with the same H u and n, then
the supernode is a rectangle with the same size as
the square supernode but the ratio of the two sides
being 2:1. We also use R to denote diagonal n \Theta n
matrix with vector r on its diagonal, and
gE u R and
transformation is completely specified by H u , r, and g,
and therefore, denoted by (H g). The advantage
of factoring matrix H this way is that it allows us to
study the three supernode transformation parameters
separately.
Implication 2 of Corollary 1 of [26]
For an algorithm A = (J; D), a linear schedule
[22] is defined as oe - : J ! N , such that oe -
is a linear schedule vector,
a row vector with n rational components, minf-d
Jg. A linear
schedule assigns each node j 2 J an execution step
with dependence relations respected. We approximate
the length of a linear schedule
Note that j 1 and j 2 for which oe - (j 1
are always extreme points in the iteration index
space.
The execution of an algorithm (J; D) is as follows.
We apply a supernode transformation (H
obtain (J s ; D s ). A time optimal linear schedule - can
be found for (J s ; D s ). The execution based on the
linear schedule alternates between computation and
communication phases. That is, in step i, we assign
with the same oe - available
processors. After each processor finishes all the computations
of a supernode, processors communicate by
passing messages in order to exchange partial results.
After the communication is done, we go to step i + 1.
Hence, the total running time of an algorithm depends
on all of the following: (J; D), H u , g, r, -, t c , t s .
The total running time is a sum of the total computation
time and the total communication time which
are multiples of the number of phases in the execu-
tion. The linear schedule length corresponds to the
number of communication phases in the execution.
We approximate the number of computation phases
and the number of communication phases by the linear
schedule length (2). The total running time is then
the sum of the computation time T comp and communication
time Tcomm in one phase, multiplied by the
number of phases P . Computation time is the number
of iterations in one supernode multiplied by the
time it takes to compute one iteration, T
The cost of communicating data computed in one supernode
to other dependent supernodes is denoted by
Tcomm . If c is the number of processors to which the
data needs to be sent, then ct s . This model
of communication greatly simplifies analysis, and is
acceptable when the message transmission time can
be overlapped with other operations such as computations
or communication startup of next message, or
when the communication startup time dominates the
communication operation. Thus, the total runningv1v2 denotes vector dot product of vectors v1 and v2
time is:
ct s
3 Optimal Supernode Shape
In this section, we present the results pertaining
to the time optimal supernode shape, i.e. supernode
relative side length matrix, R, and cutting hyperplane
direction matrix, H u , derived in the model and under
the assumptions set in the previous section.
In the model with constant communication cost,
only linear schedule vector and linear schedule length
in the expression (3) depend on the supernode shape.
Therefore, in order to minimize total running time, we
need to choose supernode shape that minimizes linear
schedule length of the transformed algorithm.
The problem is a non-linear programming problem:
diagonal matrix
det(H
H
where the scalar n
1=g is a constant that can be computed
independent of H u and Q, and without loss of
generality, we can exclude it from the objective func-
tion. We studied selection of supernode size in [9].
The floor operator of (1) has been droped in the objective
function to simplify the model. It can be shown
that the error in the linear schedule length is bounded
by
which is insignificant for components of -
close to 1 and large iteration index spaces.
Theorem 1 gives a closed form of the optimal linear
schedule vector for the transformed algorithm.
Theorem 1 An optimal supernode transformation,
with I ' D s , has an optimal linear schedule -
Proof: As defined in section 2, min- s d
D s . Since I ' D s , and in order to have feasible linear
schedule, i.e. - s D s ? 0, we must have - s - 1.
If all extreme projection vectors have non-negative
components, then their linear schedule length is minimized
with the linear schedule vector with smallest
components, i.e. -
If there are extreme projection vectors with negative
components, then an initial optimal linear schedule
vector may be different from 1. We still must
have in order to satisfy the definition of linear
schedule vector, i.e. min- s d 1. Let the i-th
component of - s be greater than 1. Then, we
can set - 0
modify the supernode shape
by setting is a diagonal matrix
with
Linear schedule of all vectors in the transformed algorithm
remains the same, i.e. linear schedule of vector
and we can divide - 0
s with min- 0
s which
will shorten linear schedule of all points. Therefore,
we got a shorter linear schedule for the algorithm and
got one more component in the linear schedule vector
equal 1. Continuing the process, we eventually get to
the linear schedule with all ones. 2
Theorem 2 gives a necessary and sufficient condition
for an optimal relative side length matrix R, and
consequently its inverse matrix Q, assuming the optimal
linear schedule vector 1.
Theorem 2 Let g and H u be fixed, let the linear
schedule vector be 1 and let M be the set of maximal
projection vectors in the transformed space. Relative
side lengths vector r is optimal if and only if vector
with equal components, v, belongs to the cone generated
by maximal projection vectors, of the transformed
algorithm:
Proof: Let m be linear schedule length, and without
loss of generality, let v be a vector with equal
components such that
Sufficient condition. Let vector v be included in
cone(M) and let the corresponding relative supernode
side length matrix be R. Consider another supernode
transformation close to the original with
slightly different from R.
Suppose the image of v of transformation with R 0
is and the schedule length of v 0 is 1v 0 .
Then, based on relation between geometric arithmetic
. The new
maximal projection vectors' linear schedule length can
only be greater than or equal to 1v 0 , which is greater
than Therefore, supernode relative side
length matrix R is optimal.
Necessary condition: We prove by contradiction.
Let R be optimal and assume v is not in cone(M ).
Then there exists a separating hyperplane Z:
for all x 2 Z, such that za ! 0, for all a 2 M , and
1. We can select the normal
vector z arbitrarily close to being orthogonal to vector
1 and of arbitrary length in order to ensure s
4 for convenience, we abuse notation and write P (v) to mean
-sv, where the choice of - is clear from the context0000000000000000000000000000000000000000000000001111111111111111111111111111111111111111111111111111
Z
z
s
cone(M)

Figure

1: Construction of vector s. Illustration for
the proof of Theorem 2.
and
1. The former is ensured by selecting
sufficiently small length of vector z. The latter is ensured
by selecting an appropriate angle between z and
1 such that s sinks on the curve
Based
on relation between arithmetic and geometric mean,
we must have s1 ?
The latter is the case
by the construction of Z, and thus construction of s
is feasible. Figure 1 illustrates construction of vector
s in two dimensional space. Then, by further scaling
supernode index space by diag(s), i.e. by choosing
improve linear schedule length
of vectors in M :
for all a 2 M . By choosing vector z arbitrarily close
to being orthogonal to vector 1, we can ensure that
no extreme projection vector becomes a new maximal
projection vector. By improving R to R 0 we contradict
the hypothesis that R is optimal. 2
Theorem 2 also implies the relation between Q and
H u , and enables analysis of H u independent of Q, as
stated by the following corollary.
Corollary 1 Given H u and the vector u in the convex
hull of the original iteration index space which
maps into vector with equal components and with maximal
linear schedule length in the convex hull of the
supernode index space, optimal Q has components:
With optimal selection of Q, objective function of (4):
1QH u u (6)
reduces to the expression:
Y
Proof: The relation (5) is readily derived from
const1. The expression (7)
follows by substituting expression (5) for Q in (6). 2
In the special case of a single maximal projection
vector, optimal Q can be easily computed based on
(5). For example, if the iteration index space is a
hyperrectangle and H optimal Q is the one
which turns supernode index space into a hypercube,
and makes supernode similar to the original iteration
index space.
Based on (7), our objective function for optimal H u
is:
Y
where u is the vector in the convex hull of the original
iteration index space as defined by Corollary 1.
The following shows that by positively combining
rows of a matrix, i.e. by taking any matrix with rows
inside the cone generated by rows of the original ma-
trix, we obtain no better cutting hyperplane direction
matrix.
u1 be a cutting hyperplane direction
matrix. Let U be a square matrix such that
and U - 0. Then, matrix H u 2
does not
give a cutting hyperplane matrix with shorter schedule
length.
Proof: It is enough to show that F (H u 1
with the vector u in (8) being the vector as defined
by Corollary 1, corresponding to H u 1
Let W be a square non-negative matrix of unit determinant
such that
for H u1 and u, and let 1Q 1 H u 1
Y
Y
Y
Y
Y
Y
where we used inequality between sum of non-negative
numbers and root of the sum of their squares, and
Hadamard inequality in the deductive sequence above.Based on Lemma 1, we can state the following regarding
the optimality of the choice of H u in general.
Theorem 3 Optimal hyperplane direction matrix H u
assumes row vectors from the surface of the tiling
cone.
Proof: If any of the row vectors of H u are from
the interior of the tiling cone, then there is another
hyperplane direction matrix H 0
u with row vectors from
the surface such that H
i.e. the cone generated by the rows of H u is included
in the cone generated by H 0
which takes row vectors
from the surface of the tiling cone. Based on Lemma
u is a better choice than H u . 2
In the case of algorithm with exactly n extreme
dependence vectors, we can state a stronger result,
provided by the following theorem.
Theorem 4 Optimal hyperplane direction matrix H u
for an algorithm with n extreme dependence directions
is uniquely defined, up to uniform rescaling, by the n
extreme directions of the tiling cone.
Proof: The dependence cone having exactly n
extreme directions, implies that there are exactly n
extreme directions of the corresponding tiling cone.
Then each row of any hyperplane direction matrix is
a non-negative linear combination of the n extreme
directions of the tiling cone, and we are sure that the
hyperplane direction matrix with extreme directions
of the tiling cone as its rows is the best choice, based
on Lemma 1. 2
An equivalent statement for two dimensional algorithms
gives an even stronger result.
Theorem 5 Optimal extreme direction matrix E u for
two dimensional algorithms has column vectors in the
directions of the two extreme dependence vectors.
Proof: Two dimensional algorithms always have
exactly two extreme dependence vectors. Since
Theorem 4 applies always. 2
The following example confirms that the optimal
hyperplane directions matrix H u can take row vectors
from all of the surface of the tiling cone, in the general
case, and not only from the set of extreme directions
of the tiling cone.
Example 1 In this example, we apply supernode
transformation to Jacobi algorithm [19]. We select
optimal Q, and discuss the selection of H u , assuming
(the selection of g is discussed in [9]).
The core of the Jacobi algorithm, with constant
number of iterations, can be written as follows:
do
do all
Iteration index space is a three dimensional rectan-
gle. There are eight extreme points in the iteration
index space, shown as column vectors of the matrix
X:
500 500 500 500C A :
There are four dependencies, caused by access to
elements of array a in the code above. They are represented
in matrix D as column vectors:
The tiling cone generating vectors, i.e. extreme vectors
of the tiling cone, for this example, include the
following four row vectors:
We can construct four different hyperplane direction
matrices, H ui , 4, from the four row vectors
of matrix B. However, none of those four matrices,
constructed from extreme directions of the tiling cone,
is necessarily optimal as we will see through the rest
of this example.
From the set of extreme points, X, we can construct
a set of extreme projection vectors. Those are 28 vector
differences between all different pairs of column
vectors of X. Applying an iterative non-linear optimization
procedure we obtain optimal Q i for each
H ui above. Applying the supernode transformations
obtained in this way to the set of extreme projection
vectors and finding the maximum, we get the same
linear schedule length of 606:2445 for each transfor-
mation. Thus, all H ui 's are equally good. This is due
to regularity of the iteration index space and the dependence
vectors.
Let us consider another hyperplane direction matrix

0:5 \Gamma0:5 \Gamma0:5
constructed from two row vectors of B and a sum of
the other two row vectors of B, and then properly
normalized. Corresponding optimal supernode relative
side lengths are given by:
Corresponding tiling matrix is:
0:125 \Gamma0:125 0:125
0:125 \Gamma0:125 \Gamma0:125
Matrix
5 gives the shape of the supernode:
Applying H 5 to the set of extreme index points, we get
the extreme index points in the supernode index space:
62:5 \Gamma312:5 \Gamma187:5
Similarly, we get the transformed extreme projection
vectors, but we don't show all 28 of them. We show
only the maximal projection vectors, there are nine of
them:
The nine maximal projection vectors sink on a single
two dimensional plane, and the vector
with equal components and the same linear schedule
length, belongs to the cone formed by the maximal projection
vectors and sinks on the same plane, ensuring
that we have selected the optimal relative side lengths
based on Theorem 2. This can be easily verified by
computing normal vectors to the faces of the cone generated
by maximal projection vectors and showing that
vector v's projection onto these normals indicates that
the vector v is inside the cone.
4 Related work
Irigoin and Triolet [14] proposed the supernode partitioning
technique for multiprocessors in 1988, as
a new restructuring method. Ramanujam and Sadayappan
[20] studied tiling multidimensional iteration
spaces for multiprocessors. They showed equivalence
between the problem of finding a partitioning
hyperplane matrix H, and the problem of finding a
containing cone for a given set of dependence vectors.
In [4], the choice of supernode shape is discussed with
the goal of minimizing a new objective function. The
key feature of the new objective function is its scal-
ability. It is defined in such a way that it is independent
of the supernode size. In [24], the authors
discussed the problem of finding the optimum wave-front
(optimal linear schedule vector, in our terminol-
ogy) that minimizes the total execution time for two
dimensional (data) arrays executed on one or two dimensional
processor arrays. In [18], the optimal tile
size is studied under different model and assumptions.
In [26], an extended definition of supernode transformation
is given. It is an extension of the definition
originally given in [14]. In [27] the choice of matrix
H with the criterion of minimizing the communication
volume is studied. Similar to [4], the optimization criterion
does not include iteration index space, and thus
the model does not include the linear schedule effects
on the execution time. In [3], the choice of optimal tile
size that minimizes total running time is studied. The
authors' approach is in two steps. They first formulate
an abstract optimization problem, which does not
include architectural or program characteristics, and
partially solve the optimization problem. In the second
step they include the architectural and program
details into the model and then solve the problem for
the optimal tile size yielding a closed form solution.
Recently, an international seminar [6] was held with
the topic of tiling where twenty five lectures were pre-
sented. The lectures covered many issues related to
the tiling transformation.
Selection of optimal supernode size was studied in
our previous work within a similar model in [9]. We
studied the choice of cutting hyperplane directions in
two dimensional algorithms in [12], selection of supernode
shape for the case of dependence cone with
extreme directions in [13] and the results presented
in this paper are an extenssion to the results of [13].
Compared to the related work, our optimization criterion
is to minimize the total running time, rather
than communication volume or ratio between communication
and computation volume. In addition, we
use a different approach where we specify a supernode
transformation by the supernode size, the relative
side length vector r, and the cutting hyperplane direction
such that the three variables become
independent, and can be studied separately.
5 Conclusion
We build a model of the total running time based
on algorithm characteristics, architecture parameters
and the parameters of supernode transformation. The
supernode transformation is specified by three independent
parameters: supernode size, supernode relative
side lengths, and the cutting hyperplane direc-
tions. The independence of the parameters allows
us to study their selection separately. In this paper,
two of the three parameters are studied. We give a
necessary and sufficient condition for optimal relative
side lengths, show that the optimal cutting hyper-plane
directions from the surface of the tiling cone,
and show that linear schedule - is an optimal
linear schedule in the transformed algorithm. If the
final supernode transformation violates the assumption
that I ' D s , then the results do not hold. The
results are derived in continuous space and should for
that reason be considered approximate.



--R

"Scanning Polyhedra with Do Loops,"
"Optimal Orthogonal Tiling,"
"Optimal Orthogonal Tiling of 2-D Iterations,"
"(Pen)- Ultimate Tiling,"
"Practical Dependence Testing,"
"Tiling for Optimal Resource Utiliza- tion,"
"Linear Scheduling is Nearly Optimal,"
"Eval- uating Compiler Optimizations for Fortran D,"
"On Supernode Transformations with Minimized Total Running Time,"


"On Supernode Partitioning Hyperplanes for Two Dimensional Algorithms,"
"Time Optimal Supernode Shape for Algorithms with n Extreme Dependence Di- rections,"
"Supernode Partitioning,"
"The Cache Performance and Optimizations of Blocked Al- gorithms,"
New Jersey
"Modeling Optimal Granularity When Adapting Systolic Algorithms to Transputer Based Supercomput- ers,"
"Op- timal Tile Size Adjustment in Compiling General DOACROSS Loop Nests,"
Parallel computing: theory and practice
"Tiling Multidimensional Iteration Spaces for Multicomputers,"
"Automatic Blocking of Nested Loops"
"Time Optimal Linear Schedules for Algorithms With Uniform Dependen- cies,"
"Independent Partitioning of Algorithms with Uniform Dependencies,"
"Finding Optimum Wavefront of Parallel Computation,"
"More Iteration Space Tiling,"
"On Tiling as a Loop Transformation,"
"Communication-Minimal Tiling of Uniform Dependence Loops,"
--TR

--CTR
Georgios Goumas , Nikolaos Drosinos , Maria Athanasaki , Nectarios Koziris, Automatic parallel code generation for tiled nested loops, Proceedings of the 2004 ACM symposium on Applied computing, March 14-17, 2004, Nicosia, Cyprus
Sriram Krishnamoorthy , Muthu Baskaran , Uday Bondhugula , J. Ramanujam , Atanas Rountev , P Sadayappan, Effective automatic parallelization of stencil computations, ACM SIGPLAN Notices, v.42 n.6, June 2007
Georgios Goumas , Nikolaos Drosinos , Maria Athanasaki , Nectarios Koziris, Message-passing code generation for non-rectangular tiling transformations, Parallel Computing, v.32 n.10, p.711-732, November, 2006
Maria Athanasaki , Aristidis Sotiropoulos , Georgios Tsoukalas , Nectarios Koziris , Panayiotis Tsanakas, Hyperplane Grouping and Pipelined Schedules: How to Execute Tiled Loops Fast on Clusters of SMPs, The Journal of Supercomputing, v.33 n.3, p.197-226, September 2005
Saeed Parsa , Shahriar Lotfi, A New Genetic Algorithm for Loop Tiling, The Journal of Supercomputing, v.37 n.3, p.249-269, September 2006
distributed memory multicomputer;minimizing running time;parallelizing compilers;algorithm partitioning;tiling;supernode transformation