日本ソフトウェア科学会第 31 回大会 (2014 年度) 講演論文集

SGL Hierarchical Bridging Modelfor Parallel Programming
LI Chong§†

HAINS Gaétan†‡

Parallel programming and data-parallel algorithms have been the main techniques supporting highperformance computing for many decades. Like all non-functional properties of software, the conversionof computing resources into scalable and predictable performance involves a delicate balance of abstraction
and automation with semantic precision. A major conceptual step was taken by L. Valiant who introduced
the Bulk-Synchronous Parallel (BSP) model. Parallel algorithms on BSP can be designed and measured bytaking into account not only the classical balance between time and parallel space but also communicationand synchronization. However, the flat view of a parallel machine as a set of communicating sequential machines remains true but is more and more incomplete. Inspired by BSP, we introduce here the SGL bridging
model, aim to improve the simplicity of parallel program development, the portability of parallel programcode, and the precision of parallel algorithm performance prediction on both classical parallel machinesand novel hierarchical machines. SGL can be used later as a part of model transformation for generating
automatically eﬃcient parallel program with its cost model.

and scalable performances as predicted by the

1 Introduction

model. BSP is thus a bridging model relating par-Parallel programming and data-parallel algo-allel algorithms to hardware architectures.

rithms have been the main techniques supporting

From the mid-1990’s, it became clear that BSP

high-performance computing for many decades. A

is the model of choice for implementing algorith-major conceptual step was taken by L. Valiant who

mic skeletons: its view of the parallel system in-

introduced the Bulk-Synchronous Parallel (BSP)

cluded explicit processes and added a small set

model [23]. Inspired by the complexity theory ofof network performance parameters to allow pre-PRAM model of parallel computers [6], Valiant

dictable performance. Hains et al. designed BS-

proposed that parallel algorithms can be designed

lambda [18] as a minimal model of computationand measured by taking into account not only thewith BSP operations. BS-lambda became the basis

classical balance between time and parallel spacefor Bulk-Synchronous ML (BSML) [17] a variant(hence the number of processors) but also com-of CAML under development by Loulergue et al.

munication and synchronization.

The BSP per-since 2000. BSML has a simplified programming in-formance model is both realistic and tractable soterface of only four operations: mkpar to constructthat researchers like McColl et al. [19] were able toparallel vectors indexed by processors, proj to map

define BSP versions of all important PRAM algo-them back to lists/arrays, apply to generate asyn-

rithms, implement them and verify their portablechronous parallel computation and put to generate§ , 国立情報学研究所, National Institute of Informatics.

communication and global synchronization from a† , バリ大学（東）, Université Paris-Est, France.

two-dimensional processor-processor pairing. As a

‡ , 国立応用科学院サントル校, INSA Centre, France.

result, parallel performance mathematically follows

from program semantics and the BSP parametersof the host architecture.

The communication on SGL is conceptually centralized. SGL does not express “horizontal” com-

While BSML was evolving and practical experi-munication patterns, even this defect aﬀects a mi-ence with BSP algorithms was accumulating, one

nority of algorithms and can be compensated byof its basic assumptions about parallel hardware

automated compilation/interpretation. Therefore,

was changing. The flat view of a parallel machine

the GPS theorem was proposed in [16] that can

as a set of communicating sequential machines re-

be implemented later in a compiler to optimize themains true but is more and more incomplete. Re-

SGL’s “horizontal” all-to-all communication.

cent supercomputers like IBM Blue Gene [21] fea-

The programming model of SGL replaces thetures multi-processors on one card, multi-core pro-

BSML programming primitives with scatter, gather

cessors on one chip, multiple-rack clusters etc. The

and pardo. The SGL cost model improves the clar-

Cell/B.E. [9], Cell-based supercomputer RoadRun-

ity of algorithms performance analysis. At the samener [2] and GPU-equipped hybrid machines featuretime, it allows benchmarking machine performancea CPU with Master-Worker architecture. More-and algorithm scalability. SGL provides the pos-over, [11] observes that heterogeneous chip mul-sibility to divide parallel programming into three

tiprocessors present unique opportunities for im-

diﬀerent dimensions in order to reduce the devel-

proving system throughput and reducing processor

opment complexity. Then the tasks on each di-consumption. The trend towards Green Comput-

mension can be conquered independently by either

ing puts even more pressure on the optimal use of

system engineer, parallelism specialist, or end-use

architectures that are not only highly scalable but

domain expert.

hierarchical and non-homogeneous.

The rest of this paper is organized as follow: Sec-Li and Hains [14] [15] proposed a hierarchicaltion 2 gives an overview of the SGL model; Sec-

bridging model named SGL for modelling hetero-tion 3 presents skeletal parallelism on SGL; Sec-

geneous parallel architectures and providing a pro-tion 4 introduces the GPS theorem for horizontal

gramming platform for parallel algorithms devel-

communication on SGL; and Section 5 is the con-opment. SGL attempts to improve the simplicity

clusion.

of parallel program development, the portability ofparallel program code, and the precision of parallel

2 The SGL Modelalgorithm performance prediction on both classical2. 1 The SGL Computer

parallel machines and novel hierarchical machines.

The PRAM model [6] of parallel computation be-

Several skeletons such as parallel scan, parallelgins with a set of p sequential Von-Neumann ma-sort, distributed homomorphism, etc. have beenchines (processors). It then connects them with

implemented on SGL [15] [13] to attempt to mo-a shared memory and an implicit global controller

tivate and support the view that BSP’s advan-(as in a SIMD architecture). The BSP model [23]tages for parallel software can be enhanced by therelaxes the global control and direct memory ac-recursive hierarchical and heterogeneous machinecess hypotheses. It also begins with a set of p se-

structure of SGL, while simplifying the program-quential machines but assumes asynchronous ex-

ming interface even further by replacing point-to-

ecution except for global synchronisation barriers

point messages with logically centralised communi-

and point-to-point communications that complete

cations.

between two successive barriers. It is common toassume that the PRAM/BSP processors are identical, but it is relatively easy to adapt algorithms toa set of heterogeneous processors with varied processing speeds. In both cases the set of processors
has no structure other than its numbering from 0
to p−1. This absence of structure on the set of processors can be traced back to the failure of a trendof popular research in the 1980s: algorithms for
specific interconnect topologies such as hypercube

図1

A multi-level SGL computer

algorithms, systolic algorithms, rectangular grid algorithms, etc. We trace this failure to the excessive

ters themselves or leaf-workers. The number ofvariety of algorithms and architectures without achildren is unbounded so that the BSP/PRAM

model to bridge the portability gap between algo-concept of a flat p-vector of processors is easily

rithms and parallel machines.

simulated in SGL. The tree structure oﬀers a

The SGL model [15] we proposed does introduce

vertical scalability to SGLC that the BSP con-a degree of topology on the set of processors. But

cept can scale only horizontally the number ofthis topology is purely hierarchical and is not in-processors.

tended to become an explicit factor of algorithm
diversity. Just as PRAM/BSP algorithms use the p

2. 2 SGL Execution Modelparameter to adapt to systems of all sizes, our SGL

An SGL program execution is a sequence of super-algorithms use the machine parameters to adapt tosteps. The initial computing data and final result

all possible systems.

can be either distributed in workers or centralised

The abstract machine or an SGL Computer(SGLC) is defined as the combination of three attributes (図 1):

in the root-master. Each superstep is composed offour phases:
1. A scatter communication phase initiated

1. A set of workers, which are sequential pro-by the master. The master scatters data to

cessors composed of a computation element

its children if the data is not yet distributed,(“computing core”) and local memory unit.

then it engages its children to start the second

The workers provide the primary computing

phase.

power.

2. An asynchronous computation phase per-

2. A set of masters, which are also sequential

formed by the children. The children executeprocessors composed of a coordination element

the task initialized by their master in the first(“coordinating core”) and central memory unit

phase. A child can be either a master or awhose data is accessible by its own children

worker. The workers accomplish its master ’s

via the coordinating core, i.e. the coordinat-

task; and the child-master can start a sequenceing core can scatter (resp. gather ) message toof supersteps with its own children nested in(resp. from) its children. The masters limit

this child computation phase so-called sub-the communication cost.

supersteps.

3. A tree structure that the root is a master3. A gather communication phase centred on

and the master ’s children may be either mas-

the master. The master synchronizes its chil-dren’s task of the second phase, then gath-ter and its children. As always, local computation

ers the computation results from its children

cost is the number of instructions executed divided

if necessary.

by processing speed.

4. A local computation phase on the master.

The machine-dependent parameters can be madeThe master post-processes the gathered data

as abstract or concrete as desired.

of the third phase and terminates the super-words, theoretical SGL algorithms can be inves-

step.

tigated with respect to their asymptotic complex-

In otherWith these four phases, we can usually write an

ity classes, while portable concrete SGL algorithms

SGL program recursively. The choice of using (orcan be analysed for more precise cost formulae us-not) children in the recursive program depends on

ing bytecode-like instruction counts and normalized

performance parameters that combine (known) pa-

communication parameters, and finally SGL pro-rameters of: communication costs, synchronization

grams can be compiled and measured for actual

costs, computation costs, and load balancing.

performance on a given architecture.
Concrete cost estimations and measurements use

2. 3 SGL Cost Modelthe following machine and algorithm parameters.

Like all bridging models, SGL’s main goal of ex-

Each level of a hierarchy can have diﬀerent param-pressing parallel algorithms requires a precise notion of the execution time for a program. We give

eter values.
• Machine parameters:

here the mathematical form of this cost model with

– p : the number of children processors thatthe understanding that a full definition should bea master has. We use P for the total num-

based on the operational semantics and will be de-ber of leaf-workers of a machine.

fined in a further document. The cost equations– c : computation speed of processors, c0

below are nevertheless suﬃcient for informal algo-denotes the time interval for performing arithm design/understanding and comparison withunit of work on the master and ci (i=1..p)other parallel models.

denotes the time interval for performing

The cost of an SGL algorithm (i.e. of its execu-

a unit of work on children processor. Pa-tion on given input data) is the sum of all its super-rameter c without an index refers to a local

steps’ cost. The cost of an individual superstep isquantity.

composed of computation cost and communication– g : the gap, g↓ defined as the minimum

cost. Computation cost is the sum of local com-time interval for transmitting one word

putation time in the children (parallel, hence com-from master to its children, and g↑ for chil-bined with max) and of the local computation indren to its master. We use a single g in

the master. Addition of both terms realizes the hy-

case of symmetric communication cost.

pothesis that the master’s local work may not over-– l : the latency to perform a gather com-

lap with the children’s work. Communication costmunication synchronization. i.e. the timeis split into the time for performing the scatter andto execute a 1-bit gather. Since a scat-

the gather operation. Communication costs are es-

ter communication should be initialized by

timated by a linear term similar to BSP’s g × h + L,

the master, synchronisation is not neces-

based on machine parameters for (the inverse of)sary for scatter phase. We use L for the

bandwidth and synchronization between the mas-

total cost of all-level synchronisations of

one global superstep.

tion.
−
→
• V ranges over vectorial locations VecLoc,• Algorithm parameters:– w : work or number of local operations
performed by processors, w0 denotes themaster ’s work and wi (i=1..p) denote thework of a child i. A parameter w without
an index refers to a local quantity.
– k : number of words to transmit, k↓ denotes number of words that master scatteri.e. names of memory elements to store ar−→rays. Here V i=pid denotes master /children lo−→
cations; V without index denotes master location.
f ranges over vectorial vectorial locations• W
VVecLoc, i.e. names of memory elements tostore arrays of arrays.

to its children, and k↑ denotes number ofExpressions are relatively standard with thewords that master gathers from its chil-convenience of scalar-to-vector (sequential) oper-dren.

ations.

Thus, the cost formula is:[
Costnode =w×c+⊙ ranges over binary arithmetic oper-a ranges over scalar arithmetic expres−→sions::= n | X | a ⊙ a | V [a] b ranges over scalarations.

max {Costchildi }

i=1..p]+k↓ × g↓ + k↑ × g↑ + l
which clearly covers the possibility of a heteroge-

and w ranges over vectorial vectorial expressions::=f.
< v1 , v2 , ..vℓ > | W

neous architecture.

The language’s commands include classical se-

2. 4 SGL Programming Model
Bougé advocated in his paper [3], that an abstract computing model needs an execution model,
but also a programming model.

boolean expressions. v ranges over vectorial expres−→f [a]sions ::= < a1 , a2 , . . . aℓ > | V | v ⊙ a | v ⊙ v | W

We enrich

Winskel’s basic imperative language IMP [24] toyield our deterministic parallel programming language – SGL.

quential constructs with SGL’s three primitives:
scatter, pardo and gather. Their exact meaningsof parallel statements are defined in the semantic
rules below.
• c ranges over primitive commands.
Com ::=−→
f := w | c ; cskip | X := a | V := v | W

Values are integers, booleans and arrays (vectors)| if b then c else cbuilt from them. Vectors of vectors are necessary| for X from a to a do c−→
| scatter w to V | scatter v to X−→→
f | gather X to −| gather V to W
V

for building blocks of work to be scattered among
workers. The scatter operation takes a vector ofvectors in the master and distributes it to workers/children. The gather operation inverts this process.
Imperative variables are abstractions of memory
positions and are called locations. They are manysorted like the language’s values.
• X ranges over scalar locations NatLoc, i.e.
names of memory elements to store numbers. Here Xi=pid denotes master /children locations; X without index denotes master loca-| pardo c | if master then c else c• Auxiliary commands are also needed.
−→f
Aux ::= numChd | len V | len W
States (or environments) are maps from imperative variables (locations) to values of the corresponding sort. Like values they are many-sorted
and we use the following notations and definitions
for them. The functions in States Σ are defined asfollow:• σ : N atLoc → N at, thus σ(X) ∈ N at−→
• σ : V ecLoc → V ec, thus σ( V ) ∈ V ecf) ∈• σ : V V ecLoc → V ecV ec, thus σ(WV ecV ec
Here P os ∈ N at is what we call the (relative)• Gather−→⟨ V :=< X1 , X2 , . . . XnumChd >, σ⟩ → σ ′−→⟨gather X to V , σ⟩ → σ ′position:
P os = 0 denotes master position (same as above),

The values in location X on the children are

and P os = i ∈ {1..p} denotes position in isent to the master, then they are stored as a−→vector in the location V on the master. A syn-thchild.

It is the recursive analog of BSP’s (or MPI’s) pids.
• σ : N atLoc → P os → N atchronisation is presented here to ensure the re-• σ : V ecLoc → P os → V ec

ception of all values. The communication cost

The semantics of vector expressions is standardof this statement is (p words × g↑ + l).

and deserves no special explanations. The evaluation of parallel primitive commands are defined asfollow:→ −→−→f :=< −
⟨W
V 1 , V 2 , . . . V numChd >, σ⟩ → σ ′−→f , σ⟩ → σ ′⟨gather V to W−→
Same as the statement (gather X to V ) but• Scatterfor gathering vectorial vectorial expression.
⟨v, σ⟩→<n1 ,n2 ,...np > ∀i=1..numChd ⟨Xi :=ni , σ⟩→σi′⟨scatter v to X, σ⟩→σ ′The communication cost of this statement is∑( number of words of v × g↑ + l).

The vector expression v is evaluated in themaster ’s local environment which is a part ofthe initial global environment σ to get thescalar value of each element of v; and the length• Parallel∀i=1..numChd ⟨c, σi ⟩ → σi′⟨pardo c, σ⟩ → σ ′of v shall be the same as the number of chilAfter

Each child of the master evaluates indepen-

that, the scalar value of each element of v isdently the command c in its own local environ-sent to a child ’s environment that the position

ment σi where the position i = 1 . . . p and p isof the child is the same as the indexation of

the number of children. This statement costs

the sending element of v. We thus have a new

maxi=1..p {Cost(c)i }.

dren processors that the master has.

′

global environment σ composed of all new local environments of the master and its children.

⟨numChd = 0, σ⟩ → true ⟨c2 , σ⟩ → σ ′⟨if master then c1 else c2 , σ⟩ → σ ′The communication cost of this statement is (pwords × g↓ ).
⟨w, σ⟩→<v1 ,v2 ,...vp >−→∀i=1..numChd ⟨ V i :=vi , σ⟩→σi′⟨numChd = 0, σ⟩ → false ⟨c1 , σ⟩ → σ ′⟨if master then c1 else c2 , σ⟩ → σ ′−→⟨scatter w to V , σ⟩→σ ′If the local environment is on a master, i.e.

Same as the statement (scatter v to X) but

the number of children is zero in the local en-

for scattering vectorial vectorial expression.

vironment, the command c1 will be evaluated;

The communication cost of this statement is∑( number of words of v × g↓ ).

otherwise, the command c2 will be evaluated.

repl x n=

[x, . . . , x]

map f [x1 , . . . , xn ]=[(f x1 ), . . . , (f xn )]mapidx g [x1 , . . . , xn ]=[(g 1 x1 ), . . . , (g n xn )]zip ⊕ [x1 , . . . , xn ] [y1 , . . . , yn ]=[x1 ⊕ y1 , . . . , xn ⊕ yn ]reduce ⊕ [x1 , . . . , xn ]=x1 ⊕ · · · ⊕ xn

scan ⊕ [x1 , . . . , xn ]=[x1 , (x1 ⊕ x2 ), · · · , ((x1 ⊕ x2 ) · · · ⊕ xn )]表1

Simple data-parallel skeletons

him and use a SGL programming style when new

3 Skeletal Parallelism

patterns are needed.

3. 1 Skeleton paradigm
There exist two kinds of algorithmic skele-

3. 2 Data-parallel skeletons

tons [10]: tasks and data-parallel ones. The former

表 1 defines the functional semantics of a set of

can capture parallelism that originates from exe-

data-parallel skeletons [4] [1]. It can also be seen

cuting several tasks, i.e. diﬀerent function calls, in

as a naive sequential implementation using lists.

parallel. They mainly describe various patterns for

The skeletons work as follow. Skeleton repl cre-organizing parallelism, including pipelining, farm-ates a new list containing n times element x. Hereing, client-server, etc. The latter parallelize compu-

we speak of lists for the specification but paral-

tation on a data structure by partitioning it amonglel implementations would use more eﬃcient data-processors and performing computation simultane-

structures as arrays (in this article) or a stream (inously on diﬀerent parts of it.

a client/server or grid environment) since the size

A well-know disadvantage of skeleton languagesof the lists remain constant.

is that the only admitted parallelism is usually that

The map, mapidx and zip skeletons are equiv-of skeletons, while many parallel applications arealent to the classical Single-Program-Multiple-Data

not easily expressible as instances of known skele-(SPMD) style of parallel programming, where a sin-

tons. Skeleton languages must be constructed as togle program f is applied to diﬀerent data in parallel.

allow the integration of skeletal and ad-hoc paral-

Parallel execution is obtained by assigning a sharelelism in a well defined way [4]. In this light, hav-of the input list to each available processor.

ing skeletons in SGL would combine the expressive

reduce is an elementary data-parallel skeleton:

power of collective communication patterns withthe function reduce ⊕ e l computes the “sum” ofthe clarity of the skeleton approach.

all elements in a list l, using the associative binary

In this work we consider the implementation ofoperator ⊕ and a neutral element e. Reduction has

well-known data-parallel skeletons because they are

traditionally been very popular in parallel program-

simpler to use than task-parallel ones for coarse-

ming and is provided as the collective operation

grained models and also because they encode many

MPI Reduce in the MPI standard. Note that thescientific computation problems and scale natu-

binary operator ⊕ may itself be time-consuming.

rally. Even if the SGL’s implementation is certainly

To parallelize the reduce skeleton, the input list is

less eﬃcient compared to a dedicated skeleton lan-divided into sub-lists that are assigned to each pro-guage (using MPI send/receive [5]), the program-cessor. The processors compute the ⊕-reductions

mer can compose skeletons when it is natural forof their elements locally in parallel, and the local (F F T u) ⊕results are then combined either on a single processor or using a tree-like pattern of computation andcommunication, making use of associativity in thebinary operator.

i(F F T x)i =i,n (F F T v)i

 (F F T u)j ⊗j,n (F F T v)jif i < n2otherwise

where u = [x0 , x2 , . . . , xn−2 ], v = [x1 , x3 , . . . , xn−1 ],j = i − n2 , and a ⊕i,n b = a + ωni b and a ⊗j,n b =The scan skeleton is similar to reduce (and isa − ωnj b. This formulation is close to the dh skele-

provided as the collective operation MPI Scan), but

ton except ⊕i,n and ⊗j,n being parametrized withrather than the single “sum” produced by reduce,
scan computes the partial (prefix) sums for all list
elements. Parallel implementation is done as for
reduce.

i and n.
These operators repeatedly compute the rootsof unity.

Instead of computing them for ev-

ery call, they can be computed once a priori

reduce and scan are basic skeletons that arelike the MPI’s collective operations. A more complex data-parallel skeleton, the Distributable Homomorphism dh presented in [1], is used to expressnand stored in a list Ω = [ωn1 , . . . , ωn2 ] accessible by both operators.

a scan. FFT can thus be expressed as follow:(F F T l) = let Ω = scan + 1 (repl (ω n) n2 )in map π1 (dh ⊕ ⊗ (mapidx triple l))a special class of divide-and-conquer algorithms.
dh ⊕ ⊗ l transforms a list l = [x1 , · · · , xn ] of size n =2m into a result list r = [y1 , · · · , yn ] of the same size,
whose elements
 are recursively computed as follows: u ⊕v
if i ≤ n2iiyi = ui− n ⊗ vi− n otherwise22

where u = dh ⊕ × [x1 , . . . , x n2 ], i.e.

dh ap-

plied to the left half of the input list l and v =dh ⊕ × [x n2 +1 , . . . , xn ], i.e. dh applied to the righthalf of l. The dh skeleton provides the well-known
butterfly pattern of computation which can be used
to implement many computations.

3. 3 Program Example
We give here an example of classical parallel numerical computation that can be performed using
the skeletons presented above.
The Fast Fourier Transformation (FFT) of a list
x = [x0 , . . . , xn−1 ] of length n = 2m yields a list
whose ith element is defined as:n−1∑(F F T x)i =xk ωnkik=0
where ωn denotes the nth complex root of unity√e2π −1/n .
The FFT can be expressed in a divide-andconquer form:

For this, we first use

where:






x2
x1
x1 ⊕in11 x2  i2  i1
=
 i1  ⊕ 
 n2 
n1
2n1
t2
(⊗ is defined similarly). The first element of each
triple contains the input value, the second one its
position and the last one the current list length. In
each dh step, these operators are applied elementwise to two lists of length n1 = n2 , resulting
in a list of length 2n1 .

x1 ⊕in11 x2 (resp.

for⊗) is defined as x1 + x2 × (th (n × ni12 ) Ω) whereth n [x1 , . . . , xn , . . .] = xn .
3. 4 Implementation of DH on SGL
The implementation of basic skeletons on SGL
can be found in [12]. Here we show the implementation of dh presented in [13].
表 2 shows the algorithmic of dh. All data is

distributed in p workers and the algorithm performas recursively as follow: first of all, each worker
performs a sequential dh with its own local data;
then, the master gathers the computed data, permutes them according the position, and scatters
the permuted data to the workers; after that, eachworker performs either ⊕ operation or ⊗ operationaccording its position. After log(p) times aboveachievements, we obtain the final result.

lel, and lines 22 - 31, the no-children case, represent
a local sequential loop. The cost of the super-stepDH(IN ⊕, IN ⊗, INOUT data)
1: if master thenis below:
CostM aster=log(p) × (2n × (g↑ + g↓ ) + 2l

par do2:DH(⊕, ⊗, data);

3:
4:

end par5:for n from 1 to log2(numChd) do

maxi=1..p (DHChildi )++2n × max(c⊕ , c⊗ ))CostW orker=

2n × n ×c⊕ +c⊗
2

We can see from the pseudo code, dh imposes

6:

gather data to tmp;

7:

for i from 1 to (len(numChd)) do

8:if ((i-1) % exp2(n))/2 = 0 thenrience in teaching Master course “Parallel Algo-9:

Swap(tmp[i], tmp[i+exp2(n)])

rithm” at Université Paris-Est Créteil [16] and the10:

end ifabove example confirm our claims that scatter and11:

end forgather can be used to program a large subset of all

12:

scatter tmp to list;

parallel functions. However, an optimal execution

13:

par doof all-to-all communication is still desirable. In theif ((PID-1) % exp2(n))/2 = 0 then14:

data := data ⊕ list;

15:

Our expe-

following section will discuss how to express such
communication pattern on SGL.

else

16:

date := data ⊗ list;

17:
18:

end if19:

end par20:

many “horizontal” communications.

end for

4 The GPS Theorem
4. 1 Gather-Scatter Communication[14] shows that many basic algorithms can be ef-

ficiently and cleanly programmed with the model,
while satisfying a simple performance model. But

21: else
22:

for n from 1 to log2(len(data)) do

the case of the distributed homomorphism algo-23:

for i from 1 to (len(data)) dorithm poses a fundamental problem for SGL. The

if ((i-1) % exp2(n))/2 = 0 then24:

tmp := data[i] ⊕ data[i+exp2(n)];

25:

else

26:

tmp := data[i-exp2(n)] ⊗ data[i];27:tree topology architecture suggests a limited bandwidth around the master. Either this is experimentally true, in which case sample-sort cannot be
generally eﬃcient (for which SGL’s simplicity could28:

end if

not be blamed) or the topology abstracts the actual29:

end forcommunication bandwidth assuming that all links

30:

data := tmp;

in the architecture have the same “width”. This31:

end for32: end if
表2

Distributable Homomorphism (dh)may hold if:(a) indeed the network is a fat-tree; or(b) the architecture’s level below the master allows for horizontal communications that SGL
does not express.

In the pseudo code, line 3 is a recursive call to

Case (a) is already covered by our SGL cost model

the algorithms, lines 14 - 18 are executed in paral-

but it was not observed in our experiments: themaster ’s bandwidth can be overwhelmed when an
algorithm such as sorting moves a large subset of
all data across the master. Case (b) is not only
likely, but is the norm for flat architectures.

−→
f and S ≡ scatter W
f
Let G ≡ gather V to W−→
to V in a system with one master and p workers.
−→
Assume also that values for V are all vectors off are equivalent tolength p. As a result values for W

The pending issue with SGL is thus: can we pro-p×p matrices of scalars. SGL code for reorganizing

gram general communication patterns such as thosesuch a parallel matrix of values is a sequence G; P; S

involved in parallel sample-sort. The if master

where P is a sequential program in the master thatconditional selects an instruction branch depend-

realizes a permutation of the matrix.

ing on the local node’s number of children: zeroor more. Also, SGL’s communications do not provide direct “horizontal” communications because of

GPS-theorem
Let G; P; S be as above, P a sequential pro→f, −gram whose non-local variables are W
V , andcan be extracted from SGL semantics either stati-

π a permutation of {1, . . . , p} such that (π):
f )i,j = σ ′ (Wf(π(i,j)) ) whenever ⟨P, σ ′ ⟩ →∀i, j. σ ′′ (W−→−→σ ′′ . Then σ ′′′ ( V )(i,j) = σ( V )π(i,j) whenevercally or dynamically.

⟨(G; P; S), σ⟩ → σ ′′′ .

their hierarchical logical structure. But if the hardware supports them, inter-worker communications

Proof.

The sub-programs must be evaluated

4. 2 The GPS Theorem

through steps: (g) ⟨G, σ⟩ → σ ′ ; (p) ⟨P, σ ′ ⟩ → σ ′′Solutions to SGL’s algorithmic incompletenessand (s) ⟨S, σ ′′ ⟩ → σ ′′′ .

Recall that environ-are proposed: a compromise between simple cen-ments σ are maps from identifiers and machine

tralised communications and more general hori-

positions (master, child 1, child 2, child of childzontal communications that are more diﬃcult toi ...)

program and error-prone. We start from the no-dices. The semantics translates step (g) into (g’):→
f /σ(−
σ ′ = σ[WV )i | i = 1, . . . , p] and step (s) into
→
′′′
′′ −
f )i |i = 1, . . . , p]. We thus(s’): σ = σ [ V /σ ′′ (Wtion that we may program communications thatare logically centralised but physically “horizontal”. To support this claim we first present a theoretical result showing how successive supersteps
can be partially compressed into a single horizontal communication-synchronization phase.
A first type of solution to SGL’s lack of horizontal
communications can come in the form of compiler
optimisations. Sub-programs that concentrate data
(gather), then locally process it on the master, thenredistribute it (scatter) can be compiled to horizontal communications without the need for concentration. This possibility was not implemented so far
for lack of a complete SGL language with syntax,
compiler and code generation, a rather long subproject to develop. But we can prove that this kindof optimisation is possible in the form of a theorem
based on our operational semantics.

to values.

The former are written as in-

have:
−
→
−
→
σ ′′′ ( V )(i,j) = (σ ′′′ ( V )i )jf )i )j= (σ ′′ (W′′ f= σ (W )(i,j)
f )π(i,j)= σ ′ (W−→= σ( V )π(i,j)

(s′ )

(π)
(g ′ ).

□
If the proposition’s hypotheses are satisfied thenpermutation π can be applied locally to the subset
of matrix data available on one worker node, andthen given as local argument to a collective communication operation, thus combining two verticalcommunications into a single horizontal one. The
local interpretation of sub-program P into π is beyond the scope of this chapter and would requirea more complex language than what is covered byour current semantics for SGL.

4. 3 BSML Programminga given distant processor.

BSML developed at Université d’Orléans and

The execution of mkpar is purely local and so

Université Paris XII (now called Université Paris-is the execution of the apply primitive. The ex-

Est Créteil ) is a library for OCaml implementing

ecution of proj uses an all-to-all communication

partially the Bulk Synchronous Parallel ML lan-and the execution of put is a general BSP commu-

guage [7]. There is in BSML an abstract polymor-

nication (any processor-processor relation can bephic type α par which represents the type of p-wide

implemented with it). Experience with BSML forparallel vectors of values of type α, one per process.

more than a decade has shown that proj is much

It is very diﬀerent from usual SPMD programmingeasier to use than put , that proj can be used towhere messages and processes are explicit, and pro-

program a large subset of all parallel functions, but

grams may be non-deterministic or may containthat algorithms such as sample-sort cannot be im-

deadlocks. In fact a large subset of BSML par-plemented without put .

allel programs are purely functional. The newest
version (0.5) of core BSML library is based on the

4. 4 Simplifying BSML’s PUT

following primitives:

We now propose one solution to SGL’s lack of

mkpar : (int → α) → α par

horizontal communication – the GPS function – forproj : α par → (int → α)parallel programming to simplify BSML’s general

apply : (α → β) par → α par → β parcommunication put function. It is based on theput : (int → α) par → (int → α) par

following simulation of SGL by BSML for flat BSP

The semantics of BSML primitives is described by
the use of parallel values. Parallel value <x0 , x1 ,. . .,xp−1 > represents a set of local values of a given
type, such that xi is stored on processor i and p is
the number of processors.
In BSML, mkpar is the parallel constructor:
mkpar f computes the value <f 0, f 1,. . ., f (p1)>. proj is the parallel destructor: proj <x0 ,machines.
There is a natural correspondence between a twolevel SGL machine and a BSML program:• the master corresponds to all non-par types
in the BSML program i.e. all values that are
replicated on each processor;• the workers correspond to all local elements ofpar types in the BSML program;

x1 ,. . ., xp−1 > computes a function f such that (f i)• scatter corresponds to mkpar (seq-to-par );= xi . apply is the asynchronous parallel trans-• gather corresponds to proj (par-to-seq); andformer:• pardo corresponds to apply .

apply <f0 , f1 ,. . ., fp−1 > <x0 , x1 ,. . .,xp−1 > computes <f0 x0 , f1 x1 ,. . ., fp−1 xp−1 >. Fi-

An algorithm such as parallel sample-sort cannot

nally, put is the synchronous (communicating) par-be programmed with pure SGL primitives i.e. withallel transformer: put <g0 , g1 ,. . ., gp−1 > computes

only mkpar, apply, proj and without the generala parallel vector of functions that contain the trans-

communication primitive put. Our proposed solu-

ported messages that were specified by the gi . Thetion is a simplified form of put that leads to the

input local functions are used to specify the outgo-

same parallel performance but resembles a G;P;S

ing messages thus: gi j is the value that processor

program as in the GPS theorem presented before.

i wishes to send to processor j. The result of ap-

The general BSML communication primitive isplying put is a parallel vector of functions dual toput and it has the following type:the gi : they specify which value was received fromput: (int -> ’a) par -> (int -> ’a) parwith the input parallel vector containing adestination-value map at each processor, and theoutput parallel vector containing a sender-value(* === BSML l i b r a r i e s === *)open Bsml
open Base
open Tools

map at each processor.
The simplified version that we propose is called
sgl gps and has the following type:(int -> ’a -> ’b list)-> (int -> ’b list -> ’c)-> ’a par -> ’c par
The first argument relates to data input for thecommunication phase. It specifies how each processor (rank an integer) splits a value into a listof p values, one for every destination. The second
argument relates to data reception after communication. It specifies how each processor (rank an
integer) aggregates a list of p received values into a
local value. The third argument is the input parallel vector and the function produces its output
parallel vector by applying a single communication-(* === A u x i l i a r y f u n c t i o n s === *)let rec ( from_to : int -> int -> int list
) = fun debut fin ->if debut > fin then [] else debut ::(from_to ( debut +1) fin ) ;;let ( procs_list : int list ) = from_to 0 (bsp_p -1) ;;(* === SGL g -p - s f u n c t i o n === *)let ( sgl_gps : ( int -> ’a -> ’b list ) ->( int -> ’b list -> ’c ) -> ’a par ->’c par ) = fun split assemble indata->
let splitted = apply ( mkpar split )indata inlet exchange = put ( parfun List . nthsplitted ) in
let permuted = parfun ( fun l -> List .
map l procs_list ) exchange inapply ( mkpar assemble ) permuted;;

synchronisation superstep, like BSML’s put.
The new function is considered to be simpler to

表3

Code of GPS function over BSML

be used than put because: (a) it separates metadata (first two arguments) from the actual data to

We implemented Tiskin-McColl parallel sample-

be communicated (third argument); and (b) thesort [22] with this new communication primitive inmeta-data is sequential.

BSML, then experimented our GPS implementation on the French Tier-0 GENCI/TGCC Curiesu-

4. 5 Experiments

percomputer’s Thin nodes: B510 bullx. Each node

表 3 shows the BSML implementation for the newis composed of two 8-core 2.7 GHz Inter Sandy

communication function. Here parfun is an abbre-

Bridge EP (E5-2680) CPUs, and the nodes are

viation forinterconnected via an InfiniBand QDR Full Fat

fun f x -> apply (replicate f) x

Tree network.

We fixed the total size of input

provided by BSML’s module Bsmlbase, wheredata at 12 800 000 integer numbers, and varied

replicate is an abbreviation for

the number of processors from 2 to 512. The ex-

fun f -> mkpar (fun i -> f)from the same module. This function is similar to

ecution time was measured by OCaml’s function
Unix.gettimeofday.

The formula of speedup we used is:
ExeT imeSeq
Speedup =parfun f < x0 , . . . , x(p−1) > = < f x0 , . . . , f x(p−1) >ExeT imeP ar
where ExeT imeSeq = 47.024 sec is based on
And bsp p is a machine parameter accessor of
List.map for parallel vectors:

BSML giving the number p of processors in the

OCaml’s sequential function List.sort processing

parallel machine.

12 800 000 random integer numbers generated by

the function Random.int.
Since the bound of sequential sorting algorithm# of proc

2

4

8

16

32

64

is O(n log n), the formula for the eﬃciency we used

Time (s)

27.5

14.3

7.5

3.7

1.9

1.1

is slightly diﬀerent from the standard one which is

Speedup

1.7

3.3

6.2

12.7

25.3

44.2

simply speedup × p1 . Here it is:Eﬃciency

0.85

0.82

0.78

0.79

0.79

0.69

n log nm log mwhere n is the size of input data per processor and

表5

Ef f iciency = Speedup ×Speedup and eﬃciency of parallel
sample-sort using BSML put

m is the total size of data: 12 800 000. This formula is coherent with Shi & Schaeﬀer’s method pre-keeps the same speedup as the standard BSML one,sented in the paper [20] introducing Parallel Sort-

because the sorting algorithms are the same foring by Regular Sampling (PSRS). Tiskin-McColl’sboth implementations. These results confirm ourclaim that GPS can simplify programming whilealgorithm is a BSP adaptation of PSRS.
表 4 shows the results of our experiment using 2

keeping a good performance.

to 64 processors. The eﬃciency is about 0.80 until
the communication cost overtakes the computation
64

cost when we use a large parallel system to compute

56
48

cies with 128, 256 and 512 processors are decreased40

to 0.34, 0.05 and 0.009.

Speedupa small amount of number. That’s why the eﬃcien-Linear
GPS
put

32
24
16
8

# of proc

2

4

8

16

32

64

Time (s)

27.3

14.3

7.5

3.7

1.9

1.1

Speedup

1.7

3.3

6.2

12.8

25.2

44.4

0.86

0.80

0.78

0.80

0.79

0.69

Eﬃciency
表4

Speedup and eﬃciency of parallel

0
2

4

8

16

32

64

Num of Proc

図2

Speedup of GPS and put parallel
sample-sort implementations

sample-sort in GPS implementation

The main diﬀerence between these two impleWe also tested on the same machine the samementations is that the sgl gps function provides

sorting algorithm implemented in standard BSMLa sequential view for coding. Developers need only(with put), which can be found in L. Gesbert’s

to focus on how to split the local data and how to

PhD thesis [8]. We used the same method as be-

assemble the received data on a local node. The all-fore for measuring the execution time, speedup andto-all communication and synchronisation are per-eﬃciency of this put implementation. The results

formed implicitly by sgl gps. The cost of commu-

are shown in 表 5.

nication may be optimised during compilation (if図 2 shows the speedup of the GPS imple-

the compiler supports natively the sgl gps func-

mentation and the put implementation of paral-tion). Contrariwise, the put function focusses onlel sample-sort. Our measurement shows that thethe communication. A destination-value map for

GPS implementation is as fast as the put imple-

each processor before put and a sender-value map

mentation. Furthermore, the GPS implementationfor each processor after put must be designed by

the developer using mkpar, apply, etc. Moreover,the function’s design makes it easier to use thandevelopers must work on a parallel environment

BSML’s put function.

even before applying put itself.

Future work will experiment multi-level and heterogeneous versions of this extension to SGL. A

5 Conclusion

study on automatic detection of machine’s commu-

The definition of SGL is only a first step towards

nication throughput with a formal network perfor-

its safe and eﬃcient application to high-level par-mance formula is also desired.

allel programming. Adding a general communications primitive to SGL may be a radical solution.

謝辞

But this would have hurt the practical and theo-

The first author thanks EXQIM company for

retical simplicity of SGL. Experience and researchfunding his research for the SGL model during 2009

with BSML shows that the semantics and practicalto 2013 with French National Association of Re-

use of such a primitive (put) is more complex andsearch and Technology. The second author thankserror-prone than pure SGL. This is obvious because

W.F. McColl for introducing him to the BSP model

the semantics of put is based on a communicationin the early 1990s. SGL is inspired by the work ofmatrix, while scatter and gather are described bycolleagues who are too numerous to cite here, but

1-dimensional vectors of messages.

special thanks are due to Prof. Frédéric Louler-

Further analysis has led us to propose two ele-gue et al. for their work on BSML. Both authors

ments of a general solution to this dilemma. The

thank Partnership for Advanced Computing in Eu-

first one is the GPS theorem: a semantic equiva-

rope (PRACE) Research Infrastructure for access-lence between a gather; P; scatter sub-programs

ing the Curie supercomputer for our experiments(where P is a local sequential program in the mas-through the MHBSP project (RA1634). We thankter ) and horizontal put-like operations. This re-sincerely Prof. Zhenjiang Hu for supporting thesult is the basis of future compiler optimisations

presentation of this paper on JSSST 2014.

whereby a family of clean but ineﬃcient SGL programs can be compiled to lower-level but more ef-

参 考 文 献ficient programs using horizontal communications.

[ 1 ] Alt, M.: Using Algorithmic Skeletons for Efficient Grid Computing with Predictable Performance, PhD Thesis, Universität Münster, 2007.
[ 2 ] Barker, K. J., Davis, K., Hoisie, A., Kerbyson,
D. J., Lang, M., Pakin, S., and Sancho, J. C.: Entering the petaflop era: the architecture and performance of Roadrunner, Proceedings of the 2008
ACM/IEEE conference on Supercomputing, SC ’08,
Piscataway, NJ, USA, IEEE Press, 2008, pp. 1:1–
1:11.
[ 3 ] Bougé, L.: The data parallel programming
model: A semantic perspective, The Data Parallel Programming Model, Perrin, G.-R. and Darte,A.(eds.), Lecture Notes in Computer Science,Vol. 1132, Springer Berlin Heidelberg, 1996, pp. 4–26.
[ 4 ] Cole, M.:
Bringing Skeletons out of the
Closet: A Pragmatic Manifesto for Skeletal Parallel Programming, Parallel Computing, Vol. 30,
No. 3(2004), pp. 389–406.

The second proposed solution is a simplified formof put that we have designed and experimented in
BSML. With it, the Tiskin-McColl BSP samplesort algorithm is programmed without having toencode a general communication matrix but only
its sender-side 1-to-n and receiver-side n-to-1 communication relations.
As a result our overall position w.r.t. the communication primitives is to promote SGL-style simplified programming for most programs, while allowing to “escape” in the more complex use of horizontal communication with the GPS function. The
GPS theorem reduces the need for the latter and[ 5 ] Falcou, J., Serot, J., Chateau, T., and Lapreste,
J. T.: QUAFF : Eﬃcient C++ Design for Parallel Skeletons, Parallel Computing, Vol. 32, No. 78(2006), pp. 604–615.
[ 6 ] Fortune, S. and Wyllie, J.: Parallelism in random access machines, Proceedings of the tenth annual ACM symposium on Theory of computing,
STOC ’78, New York, NY, USA, ACM, 1978,
pp. 114–118.
[ 7 ] Gava, F.:
Approches fonctionnelles de laprogrammation parallèle et des méta-ordinateurs.
Sémantiques, implantations et certification, PhD
Thesis, Université Paris XII-Val de Marne, LACL,
2005.
[ 8 ] Gesbert, L.: Développement systématique etsûreté d’exécution en programmation parallèlestructurée, PhD Thesis, Université Paris-Est, 2009.
[ 9 ] Johns, C. R. and Brokenshire, D. A.: Introduction to the cell broadband engine architecture, IBM
J. Res. Dev., Vol. 51(2007), pp. 503–519.
[10] Kuchen, H. and Cole, M.: The Integration of
Task and Data Parallel Skeletons, Parallel Processing Letters, Vol. 12, No. 2(2002), pp. 141–155.
[11] Kumar, R., Tullsen, D. M., Jouppi, N. P., andRanganathan, P.: Heterogeneous Chip Multiprocessors, Computer, Vol. 38(2005), pp. 32–38.
[12] Li, C.: A Software-Hardware Bridging Modelfor Simplifying Parallel Programming, PhD Thesis,
Université Paris-Est, 2013.
[13] Li, C., Gava, F., and Hains, G.: Implementationof Data-Parallel Skeletons: A Case Study Using a
Coarse-Grained Hierarchical Model, Proceedings ofthe 2012 11th International Symposium on Parallel and Distributed Computing, ISPDC ’12, Washington, DC, USA, IEEE Computer Society, 2012,
pp. 26–33.
[14] Li, C. and Hains, G.: A simple bridging
model for high-performance computing, High Performance Computing and Simulation, 2011 International Conference on, HPCS 2011, Washington,

DC, USA, IEEE Computer Society, 2011, pp. 249–
256.
[15] Li, C. and Hains, G.: SGL: towards a bridging model for heterogeneous hierarchical platforms,
Int. J. High Perform. Comput. Netw., Vol. 7,
No. 2(2012), pp. 139–151.
[16] Li, C. and Hains, G.: GPS: Towards Simplified Communication on SGL Model, Proceedings ofIEEE International Symposium on Parallel & Distributed Processing, Workshops, IPDPS ’14, Washington, DC, USA, IEEE Computer Society, 2014.
[17] Loulergue, F.: Conception de langages fonctionnels pour la programmation massivement parallèle,thèse de doctorat, Université d’Orléans, LIFO, 4 rue
Léonard de Vinci, BP 6759, F-45067 Orléans Cedex
2, France, January 2000.
[18] Loulergue, F., Hains, G., and Foisy, C.: A
Calculus of Functional BSP Programs, Science ofComputer Programming, Vol. 37, No. 1-3(2000),
pp. 253–277.
[19] McColl, W. F. and Walker, D.: Theory andAlgorithms for Parallel Computation, Euro-Par,
Pritchard, D. J. and Reeve, J.(eds.), Lecture Notes
in Computer Science, Vol. 1470, Springer, 1998,
pp. 863–864.
[20] Shi, H. and Schaeﬀer, J.: Parallel sorting by regular sampling, J. Parallel Distrib. Comput., Vol. 14,
No. 4(1992), pp. 361–372.
[21] staﬀ, I. j. o. R. and Development: Overview of
the IBM Blue Gene/P project, IBM J. Res. Dev.,
Vol. 52(2008), pp. 199–220.
[22] Tiskin, A.: The Design and Analysis of BulkSynchronous Parallel Algorithms, PhD Thesis, University of Oxford, 1999.
[23] Valiant, L. G.: A bridging model for parallel computation, Commun. ACM, Vol. 33(1990),
pp. 103–111.
[24] Winskel, G.: The formal semantics of programming languages: an introduction, MIT Press, Cambridge, MA, USA, 1993.