Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018jh170049-ISJ

Theory and Practice of Vector Processing
for Data and Memory Centric Applications
Hiroaki Kobayashi（Graduate School of Information Sciences, Tohoku University）
Abstract
Since vector data processing is one of the key technologies for the next-generation
supercomputing system, it is strongly required to clarify the capability of vector dataprocessing. This project conducts the research on vector data processing potentialthat can be used both for development of new HPC systems and for the developmentof algorithms for extremely large-scale numerical problems. To this end, theresearch group of Lomonosov Moscow State University (MSU) that has experience
with graph problems and the research group of Tohoku University that has
experience with program development for vector supercomputers work together toclarify the potential of vector data processing on graph problems. The MSU group
has developed various graph algorithms for an ultra-high scale that contains morethan 109 vertices and edges. The research group of Tohoku University has knowhow
to exploit the potential of vector supercomputers. By the international collaborationof the two groups, the theory and practice of vector data processing for ultra-highscale graph problems aims to be clarified. In this report, by the performance analysis
of four graph problems on SX-ACE, the potential of the vector data processing to thegraph problems is shown.

1.

Basic Information

statistics)(1) Collaborating JHPCN CentersoCyberscience Center, Tohoku University

analysis)(2) Research Areaso

Very large-scale numerical computationoVery large-scale data processingoVery large capacity network technology■

Very large-scale information systems

Alexander Antonov (Informational structureoAlexey Teplov (Scalability and performance
analysis)oIlyaAfanasyev(Codedesignandcode

optimization)
Tohoku University (Japan)(3) Roles of Project MembersoThis collaborative work is undertaken by

Hiroaki Kobayashi (Code optimization andperformance analysis)

Lomonosov Moscow State University (MSU)oand Tohoku University. The researchers of

Hiroyuki Takizawa (Code optimization andperformance analysis)

MSU mainly focus on algorithm study andodevelopment of graph problems that can

Akihiro

Musa

(Code

designandcode

optimization)effectively utilize vector data processing. Theoresearchers of Tohoku University conduct

Ryusuke Egawa (Code optimization andperformance analysis)

performance analysis and optimizations foroobtaining aspects for development of next

Kazuhiko Komatsu (Code design and code
optimization)

generation vector computer architectures.

2.

Lomonosov Moscow State University (Russia)

Purpose and Significance of the ResearchoVladimir Voevodin (Algorithm Development)

The overall goal of this project is theoDmitry Nikitenko (Performance analysis anddeep research on vector processing potential
1Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

on large-scale numerical methods. Thetheobtained knowledge of this project can beexperimental data is open to be madeuseful not only for developments of newcommonly and widely available. These

vector supercomputing systems, but also forresults are also positively utilized fordevelopments

development of other scientific applicationsofextremely

large-scale

numerical solutions.

as

To this end, this project focuses on a

Russianwellasgroup,

betterthevaluableunderstanding

therequirements for the next generation vector

graph problem that is one of the important

processor and memory architecture.

extremely large-scale numerical solutions. A

Theeffectiveuseofvectordata

graph problem is well known as the basis

processing becomes more essential for thefor solution of the wide range of applied andfuture design of applications as vector data

scientificasprocessing is widely used in x86 processors,communicational and transport network

IBM Power processors as well as vector

optimization, social network analysis, andprocessors. Therefore, the deep study on

web data analysis. To satisfy the demands ofvector data processing is important for boththese challenges, a graph problem withsupercomputer

more than 109 vertices and edges should beapplications.

challenges

suchdesignandlarge-scalesolved. The Russian research group tries to3.

solve such huge-scale graph problems byeffectively

utilizingprocessing

on

the

vectorSignificance as a JHPCN Joint Research
Project

data

supercomputing

From the scientific aspects, this project cansystems in this project. The Russian researchshow the effectiveness of vector datagroup considers the efficient algorithm thatprocessing in large-scale graph algorithms.

is suitable for the vector data processing and

The development in the vector-friendly

implements

graph algorithm is important. Thus, thethevectoralgorithm

onvector

JHPCN

supercomputing systems.

Joint

Research

Tohoku university research group that has ainterdisciplinary

long

Russian and Japan teams.

experiencedevelopment

of

with

operation

vectorsupercomputing

Project

isanexcellent opportunity to conduct our joint

Furthermore, by collaborating with the
and

From

the

betweencomputational

aspects,

Russiais further optimized and tuned so as to fullysystems from Japan. Regarding vector

exploitsupercomputing systems, currently, there

potentialofvectordifferentthesystems, the vector-friendly graph algorithm
the

has

researchsupercomputing

are no vector systems in Russia. Thus, thesupercomputing systems.
Through this collaborative work of both

JHPCN project helps to investigate andgroups about theory and practice, thesedevelop a graph algorithm that is suitable

research efforts become public to sharefor vector data processing on vector

valuableresults

andknowledge.

supercomputing systems.

For

4.

example, as a part of the AlgoWiki project

Outline of the Research Achievements up
to FY2016that is mainly developed and maintained in2

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

As this year is the first year, there are nosome pre-selected pairs of vertices.

research achievements up to FY2016.

In this project, the SSSP problem is
examined, since it is the simplest and most5.

Details of FY2017 Research Achievements

basic for other problems: for example, theIn the period of FY2017, we have conducted

APSP problem for large-scale graphs can be

performance

evaluation

of

four

graph

solved by repeated calls of the SSSP

an

Intel

Xeon

Gold

operations for each source vertex, since

Intel

Xeon

Phi

problems

using

processor,anKnights

traditionalalgorithms,

such

as

Landing processor. Then, we have been

Floyd-Warshal, cannot be applied because ofworking

the high memory costs.

on

deep

performanceandbottleneck analysis of graph algorithms on avector processor, especially NEC SX-ACE.

5.1.2

5.1 Performance evaluation of Single Source
Shortest Path (SSSP) problem
5.1.1

traditional

Problem description of SSSP
An

undirectedverticesgraphalgorithms:

Dijkstra,Bellman-Ford, and Delta-stepping.

𝐺 = 𝑉, 𝐸

𝑉 = (𝑣! , 𝑣! , … , 𝑣! )Algorithm description
The SSSP problem can be solved with three

and

with

•

Dijkstra's algorithm is designed tosolve the problem in graphs with edges,

edges

𝐸 = (𝑒! , 𝑒! , … , 𝑒! ) is given. Each edge 𝑒 ∈ 𝐸havinghas a weight value 𝑤 𝑒 . The path between

Variation

vertices 𝑢 and 𝑣 is defined as edges

implemented with Fibonacci heap has

sequence 𝜋!,! = (𝑒! , … , 𝑒! ) , beginning inthevertex 𝑢 and ending in vertex 𝑣, so each

𝑂 |E| + |V|log|V| . The computations ofedge follows another one. The path length

thiscan be defined as 𝑤 𝜋!,! =!
!!! 𝑤 𝑒! .

only

non-negativeoftheasymptoticallyalgorithm

weights.
algorithm,

fastestincludetime

sequentialA

traversal of vertices, starting from the∗path 𝜋!,!
, which has a minimal possible

source vertex, while putting adjacent

length between vertices 𝑢 and 𝑣, is calledvertices to the stack (heap) to betheshortestpath:∗
𝑑 𝑢, 𝑣 = 𝑤 𝜋!,!
=min 𝑤 𝜋!,! .
Depending on the choice of a vertex pair

processed

later.

Asaalgorithm

canbeexecuted

result,theonlysequentially.

between which the search is performed, the•The

Bellman-Fordalgorithm

isshortest paths problem can be formulated indesigned to solve the problem inthree different ways:

graphs, including edges with negative

•

•

•

SSSP (single source shortest paths)

weights. The computational kernel ofcomputes the shortest paths from athesingle selected source vertex.

iterations, which require traverses of allalgorithm

consists

of

a

few

paths)

graph edges; the array of distances iscomputes the shortest paths between

computed from the data of each edge.

all pairs of graph vertices.

The computations continue until there

APSP

SPSP

(all

(some

pairs

pairs

shortest

shortest

are no changes in the distance array.

paths)

computes the shortest paths between

The3

algorithm

hassequential

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure

5.1

Informational

graph

Figure 5.2 Graph edges reordering example

offor graph with 5 vertices and 16 edges

Bellman-Ford algorithm
complexity equal to 𝑂 𝑝|E| , where 𝑝computations in 𝑂

is the maximum possible length of the
5.1.3

shortest path from the source vertex to.

Implementation details of SSSP

any other. As a consequence, theFirst of all, it is very important to select theworst-casetomost suitable algorithm for the target NEC

many

SX architecture. Since Dijkstra’s algorithm isreal-world graphs, the algorithm isnaturally sequential and requires complex

finished in a much smaller number ofdata structures as well as the Delta-stepping

steps. Moreover, the algorithm has aalgorithm,

significant parallel potential: its parallel

implemented on NEC SX architecture. On

𝑂 |V||E|complexity,ishowever,

equal
for

complexity is equal to 𝑂 𝑝

|!|
!

The

delta-stepping

algorithm

thosealgorithms

are

not

the other hand, the Bellman-Ford algorithm, wheredoes not require complex data structure
processing. Moreover, each iteration of thisN is the number of processors used.
•(!"#|!|)^!
(!"# !"#|!|)

algorithm is a parallel traversal of all graphisdesigned to solve the single sourceedges,

shortest paths problem using novel

vectorization can be applied.

whereparallelprocessing

anddata-structures called buckets. It can be

During the selection of the most suitableviewed as a generalization of thealgorithm for particular parallel architecture,

Bellman-Fordtoparallel algorithm properties have to be

improve its complexity in a averagestudied. For this purpose, informational

sense. For random directed graphs withgraphs, introduced in [1], can be used.

!

Figure|!|algorithm,

aimingand uniformly distributed edge

weights,

thealgorithm

5.1

demonstrates

informationalgraphs of the Bellman-Ford algorithm.

completes

4Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Table 5.1.

System configurations.

Model
name

Intel Xeon
Gold 6126

Intel Xeon
Phi 7230

NEC
SX-ACE

Socket1

1

1

Cores persocket12

64

4

CPU MHz

2601

1297

1000

L1d cache

32K

32K

1M

L2 cache

1024K

1024K

-

L3 cache:19712K--

Figure 5.3 Performance (in MTEPS) of multicore

Parallel operations as shown in named 3
updates of graph edges on each

iteration,

the data with indirect memory accesses areparallel initialization of distances array. Bothoperations

can

be

placed more locally and stored in caches for

successfully

the longer period of time. The reordering isparallelized and vectorized.
Beforeimplementing

theimplemented in the following way: an arraychosenof distances is divided into segments (red,algorithm, it is important to determine thegreen and blue colors on Figure 5.2), whose

storage data format for input graphs. For thesizes are equal to the size of the lowest level

Bellman-Ford algorithm, the most suitablecache on the target architecture. Then, theformat is an edges list, where each edge isedges are placed into the array in thestored as a triple {vertex-start, vertex-end,following way: from the beginning of theedge’s weight}; all edges are stored in aarray, an edge, whose source vertices belong

single array in any order. Moreover, thisto the first segment of distances array, isformat allows simpler vectorization andsorted. Then, the second and the third edges

possible sorting-optimization, which aims toare sorted in order. Edges with the same

improve data locality. Thus, this format issegment number are sorted with the similar

going to be used for all implementations instrategy,

this research project.

accesses,

updates(indirected

parallel

is

graph

edges

reordering.

theirdestinationof #pragma cdir nodep directive and sxc++

oncompiler.

informational graph), the main optimization
used

to

Vectorization itself is achieved with help

memoryoperations

appliedvertices.

To achieve data locality improvement on
per-edge

the Bellman-Ford

to achieve much higher performance, since

while operations in named 2 correspond to
those

implementations ofalgorithm.

in Figure 5.1 correspond to independent

This

5.1.4

optimization significantly improves memory

Evaluation

of

the

Bellman-Fordalgorithm on NEC SX, Intel KNL, andaccess pattern, which allows the algorithm

Skylake systems
5

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure 5.4 Performance (in MTEPS) of single core implementations of the Bellman-Ford algorithm.

This section describes the comparisons ofisperformances of the Bellman-Ford algorithm

processes a fixed graph with a total amountwidely

used.

When

the

algorithm

among three different architectures: NEC

of |E| edges in T seconds, the performance

SX-ACE, Intel Knight Landing (KNL), andcan be calculated in this way.

Intel Skylake.
All results presented in this report have

𝑝𝑒𝑟𝑓𝑜𝑟𝑚𝑎𝑛𝑐𝑒(𝑖𝑛𝑇𝐸𝑃𝑆) =been obtained on two supercomputers: NEC
SX-ACEin

TohokuUniversity,

and=Lomonosov-2 in Lomonosov Moscow State
University.
partitions

Lomonosov-2
withdifferent

has

|𝐸|
𝑇

𝑒𝑑𝑔𝑒𝑠𝑐𝑜𝑢𝑛𝑡𝑒𝑥𝑒𝑐𝑢𝑡𝑖𝑜𝑛𝑡𝑖𝑚𝑒

severalconfigurations,Thus, the performance metric indicates

including nodes with Intel Xeon Gold 6126

the amount of Traversed Edges Per Second

and Intel Xeon Phi 7230. The specifications(TEPS) of a fixed graph during algorithmof those systems are presented in Table 5.1.

execution time. This metric can be calculated

In general, the KNL has the highest number

forof lightweight cores among the reviewed

implementations; Figure 5.3 demonstrates

architectures, in the meantime with theperformance comparison between multicore

smallest

amountimplementations

Moreover,inofcurrent

cachememory.

configuration

on

both

single-core

for

and

three

multicoredifferentarchitectures.

Lomonosov-2, KNL has high-bandwidth

Since

all

three

architectures

have

MCDRAM memory, which can be used as adifferent numbers of cores, it is interesting tolarge last-level cache on a variety ofcompare performance values obtained on a

problems.

single core of each target processor. For

In order to compare different graphalgorithms

implementations

these purposes we use 2 metrics: performance

betweenper core, which is a performance in TEPS

different architectures, a performance metric

metric, obtained on a single core of each6

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure 5.5 Performance (in MTEPS) of single core implementations of the Bellman-Ford algorithm.

target processor, and performance per core

important metrics for the vector architecture,normalized, which is usual performance inare99.62%

and255.9,TEPS obtained on maximum number ofrespectively. Therefore, NEC

SX

could

cores on a target processor, divided by theachieve such performance.

very

high

of

number of cores used during computations

In order to further evaluate parallel(4 for NEC SX, 12 for Skylake, 64 for Intelefficiency of the developed implementations,

KNL). In general, a performance per core

we use another 2 widely used metrics:

metric represents how well an algorithm

parallel efficiency and parallel acceleration,utilizes

which can be calculated using the following

singlecore

computational

capabilities of the target system, while aformulas:performance per core normalized metric takes
into an account how well all the cores work

𝑝𝑎𝑟𝑎𝑙𝑙𝑒𝑙 𝑎𝑐𝑐𝑒𝑙𝑒𝑟𝑎𝑡𝑖𝑜𝑛 =together on a multicore system.
Figure 5.4 demonstrates performance
per

core

implementations

comparison
for

three

=

between𝑠𝑒𝑞𝑢𝑒𝑛𝑡𝑖𝑎𝑙 𝑡𝑖𝑚𝑒𝑝𝑎𝑟𝑎𝑙𝑙𝑒𝑙 𝑡𝑖𝑚𝑒𝑝𝑎𝑟𝑎𝑙𝑙𝑒𝑙 𝑝𝑒𝑟𝑓𝑜𝑟𝑚𝑎𝑛𝑐𝑒
𝑠𝑒𝑞𝑢𝑒𝑛𝑡𝑖𝑎𝑙 𝑝𝑒𝑟𝑓𝑜𝑟𝑚𝑎𝑛𝑐𝑒different

reviewed architectures. These figures clearlyshow that the performance of NEC SX is𝑠𝑒𝑞𝑢𝑒𝑛𝑡𝑖𝑎𝑙 𝑡𝑖𝑚𝑒𝑝𝑎𝑟𝑎𝑙𝑙𝑒𝑙 𝑡𝑖𝑚𝑒
∗ 100
𝑐𝑜𝑟𝑒𝑠 𝑛𝑢𝑚𝑏𝑒𝑟
𝑝𝑎𝑟𝑎𝑙𝑙𝑒𝑙 𝑝𝑒𝑟𝑓𝑜𝑟𝑚𝑎𝑛𝑐𝑒
𝑠𝑒𝑞𝑢𝑒𝑛𝑡𝑖𝑎𝑙 𝑝𝑒𝑟𝑓𝑜𝑟𝑚𝑎𝑛𝑐𝑒=𝑐𝑜𝑟𝑒𝑠 𝑛𝑢𝑚𝑏𝑒𝑟

comparable to or higher than those of𝑝𝑎𝑟𝑎𝑙𝑙𝑒𝑙 𝑒𝑓𝑓𝑖𝑐𝑖𝑒𝑛𝑐𝑦 =Skylake and KNL in the case of the single
core performance. This is because thealgorithm implemented is suitable for thevector architecture. Thus, the average vector
length and vector operation ratio, which are

Figure 5.5 demonstrates the parallel
7

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018acceleration and parallel efficiency derived

algorithm

forsolvingtheBFS

by these equations. The parallel efficiencies

problemonprocessors.

Theof NEC SX is higher than or comparable toalgorithm is based around storing

those of Skylake. More investigation isvertices, reached on the current steprequired when the graph scale is more thanin a separate data-structure called

24 in the case of NEC SX.

frontier. Also, there have been many
implementations

5.2 Performance evaluation of BFS problem
5.2.1

Problem description of BFS
An

undirectedverticesgraph

𝑉 = (𝑣! , 𝑣! , … , 𝑣! )𝐸 = (𝑒! , 𝑒! , … , 𝑒! )•𝐺 = 𝑉, 𝐸andformulti-coreCPUs and GPU.
Direction-optimizing BFS is a novel

withapproachedges

nowadays. This approach uses two

to

BFS,

widelyusedpathdifferent bottom-up and top-down

between vertices 𝑢 and 𝑣 is defined as astages; the top-down stage is usually

set arcs𝜋!,! = (𝑒! , … , 𝑒! ), beginning in vertex

corresponding𝑢 and ending in vertex 𝑣 , so each edge

frontier-based approach, while thefollows another one. A length 𝑑 𝑢, 𝑣

bottom-up

isgiven.

The

of

stage

to

the
is

a

basic
reversepath 𝜋!,! between 𝑢 and 𝑣 vertices is the

direction search. Combining thesenumber of its arcs (edges). The breadth-firsttwo approaches allows one to checksearch finds the shortest 𝑑 𝑢, 𝑣

from asignificantly fewer edges, which

specified vertex 𝑢 to all the other nodes. (1)results in significant accelerations

There are several possible modifications

obtained [2][3].

of the BFS problem:•

it is required only to check which
vertices are reachable from specified
source vertex; (2)•it is required to find out «parent»
vertex for each one. (3)In this report, we have investigated andimplemented the basic version (1) of the
problems, such as SCC and transitive
closure. The BFS employs the version (2).
Figure 5.6 Performance (in MTEPS) of multi

5.2.2

core

Algorithms for BFS

implementations

of

the

BFS is a fundamental problem of graphBellman-Ford-based BFS algorithm for NEC

processing, so there are a lot of algorithms

SX, Intel KNL, and Skylake architectures.

that have been designed so far.
•

Frontier-based BFS is a classical8

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018•onarchitectures. It is shown that thissparse-matrix dense-vector (SpMV)

approach demonstrates comparable

multiplications, to which the BFS

results

problem

implementations.

SpMV-based

can

BFS

be

is

based

reduced.

A

implementation of this method for•toBellman-ford

thestate-of-art

based
the

BFS

has

vector architectures (Intel KNL) is

approximately

same

described in [4], where authors claim

computational complexity as thethat even though this approach has

SpMV-based BFS, but has different

great computational complexity, it isdata-structures used. The idea of thisstill more suitable for modern vector

approach is very similar to the

Figure 5.7 Performance (in MTEPS) of single core implementations of the Bellman-Ford-based BFS
algorithm for NEC SX, Intel KNL, and Skylake architectures.

Figure 5.8 Parallel acceleration and parallel efficiency comparison for NEC SX, Intel KNL, andSkylake architectures.

9

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Bellman-Ford

paths

parallel acceleration and parallel efficiency.

algorithm, where weight-distances

shortest

From these figures, the parallel efficiencies

updateswithof NEC SX are higher than those of the

distancesothers, and the high efficiency can bearereachabilityorreplaced
integerarrays.

obtained even when the graph scale is small.
A high vector operation ratio of 99.62% and

5.2.3

Implementation details of BFS

a high average vector length of 255.9

In this report, the Bellman-Ford based BFS iscontribute to the high efficiency in the case

implemented. It has very similar potentialof the single core performance.

for vectorization and concurrency. The#pragma

cdir

vprefetch

directive,the

5.3 Performance evaluation of the SCC problem#pragmacdir

vovertakedirective,

the

5.3.1

Problem description of SCC#pragma cdir vob directive, and manual loop

A directed graph 𝐺 = 𝑉, 𝐸 with vertices

unrolling are used to effectively issue the𝑉 = (𝑣! , 𝑣! , … , 𝑣! )vector𝐸 = (𝑒! , 𝑒! , … , 𝑒! ) is given. Edges may notgather/scatter

instructions

byovertaking other vector instructions.

and

edges

have any data assigned (so graphs withoutedges weights are discussed in current

5.2.4

Evaluation of the Bellman-Ford based

section). A strongly connected component

BFS on NEC SX, Intel KNL, and Skylake(SCC) of a directed graph G is a strongly

systems

connected subgraph, which is maximal

Similar to the evaluation of SSSP, Figure 5.6

within the following property: no additional

demonstrates the performance comparisonedges or vertices from G can be included in

betweenthe subgraph without breaking its property

multicore

implementations

forthree different architectures. Figure 5.7
demonstrates

performanceperof being strongly connected.

core

comparison. Figure 5.8 demonstrates the

Figure 5.9 Informational graph of Forward-Backward-Trim algorithm: top-level (left), trim-step(right).

10

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

5.3.2

Algorithm description of SCC

of the DCSC algorithms. Different variations

SCCs in the BFS can be found with one ofcan be more or less effective for different

the following algorithms.

types of graphs; the paper [7] concluded that•

Tarjan's algorithm is based on a

the forward-backward-trim algorithm is the

single depth first search (DFS) andmost efficient way to process graphs.

uses 𝑂(|E|) operations. Due to thefact that the algorithm is based on
the

•

DFS,

only

a

5.3.3

sequential

Implementation details of SCC
During the selection of the most suitable

implementation is possible.

algorithm

The DCSC algorithm (Divide andarchitectures, parallel algorithm properties

Conquer Strong Components), or FB

have to be studied. For this purpose,(Forward-backward) is based on theinformational graphs, introduced in [1], can

BFS and requires 𝑂(|E| ∗ log(|E|))

be used.

operations. This algorithm is initiallydesignedforforparticular

Figure 5.9 demonstrates informationalofthe

implementations: in each step, it

algorithm.

Asfinds a single strongly connected

informational

component

Forward-Backward-Trimandparallelallocates

parallel

graph

the

Forward-Backward-Trim
one

can

see

from

the

graph,

thealgorithm

iscomponent up to three subgraph,designed in the following way: on the firsteach of which may contain otherstep (trim step), the removal of the strongly

strongly connected components, and,

connected

as a result, this algorithm can be

performed. After that, in each step, theprocessed in parallel.

algorithm finds one nontrivial strongly

Variations of the DCSC algorithm

connected

are proposed such as Coloring andallocates

FB with step-trim. These modified

subgraphs, each of which contains otherversions of the DCSC algorithm arecomponents, and, more importantly, thisdescribed in detail in papers [5][6].

algorithm can be processed in parallel. ThisFor obvious reasons, Tarjan's algorithm

step heavily relies on breadth-first search tois not suitable for solving a problem onfind all vertices, which can be reached fromparallel architectures, since it is based on thethe selected pivot vertex, and all vertices,

depth first search as well as complex datafrom which pivot vertex can be reached.

structures (stack and queue) processing,Thus, this algorithm has two levels of

which cannot be implemented efficiently onparallelism: BFS level and parallel subgraphs

GPUs.

handling level, which are a big advantage for•

A large number of papers such as [7]components

component
thecomponent

parallel target architectures.

have already investigated different variation

11

of

size

1

is

(FB-step)

and

up

three

to

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

It is important to describe that both the
Trim and FB steps can be successfully
vectorized, since they both have a lot of
independent operations on bottom level asshown in boxes 3, 4, 5, 7, 8, and inside BFS in
Figure 5.9.
Asaresult,Forward-backward-Trimtheparallelalgorithm

isimplemented in the following way. First, in
order

to

avoid

buildingaFigure 5.10 Performance (in MTEPS) of multi

reverse(transposed) graph, which is necessary for

core

efficient backward search implementation, aalgorithm for NEC SX, Intel KNL and Skylake

graph is converted to the storage format ofarchitectures.

an edge list. To further increase theimplementationsof

Forward-Backwardof SCC on NEC SX, as the tendency of theperformance, the input edges list can bebehavior is similar to the BFS, the similar

pre-sorted as described in Section 5.1.3. The
approach greatly improves data locality andcache usage. Moreover, this optimization

optimizationsareappliedtotheBellman-Ford

based

BFS

such

as

the

directives

andmanualinsertions

also significantly improves the trim stepofunrolling in order to enhance the overtake ofefficiency, since it has the similar memory

vector gather/scatter instructions.

access pattern.
By the first analysis of the performance

Figure 5.11 Performance (in MTEPS) of single core implementations of the Forward-Backward
algorithm for NEC SX, Intel KNL, and Skylake architectures.

12

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure 5.12 Parallel acceleration and parallel efficiency comparison for NEC SX, Intel KNL, andSkylake architectures.

5.3.4

Evaluation of the Forward-Backward

Forward-Backward algorithm, especially inalgorithm on NEC, KNL, and Skylake

the case of the larger graph scale.

systems
Figure 5.10 demonstrates the performancecomparison

between

implementations
architectures.

for

Figurethree
5.11

5.4 Performance

multicoreevaluationofTransitiveofTransitive

Closure problem

different

5.4.1

demonstrates

Problem

description

Closureperformance per core comparison. Figure

A directed graph 𝐺 = 𝑉, 𝐸 with vertices

5.12 demonstrates the parallel acceleration

𝑉 = (𝑣! , 𝑣! , … , 𝑣! )and parallel efficiency. These figures show

𝐸 = (𝑒! , 𝑒! , … , 𝑒! )that the performance of NEC SX is good

between vertices 𝑢 and 𝑣 is defined aswhen the graph scale is small such as 18, 19,edges sequence 𝜋!,! = (𝑒! , … , 𝑒! ), beginningand 20. However, as the graph scale increase,in vertex 𝑢 and ending in vertex 𝑣, so eachthe performance and efficiency decrease. Inedge follows another one. Vertex 𝑣 isthe case of KNL, the performance andreachable from vertex 𝑢, if at least a single

efficiency are stable and in the case ofpath 𝑃 𝑢, 𝑣 between vertices u and v exists

Skylake, the performance and efficiency(every vertex is considered reachable fromincrease. More analysis and optimizations

itself).

arenecessary

for

NEC

SX

for

the

13

and
is

given.

edges
The

path

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure 5.13 Informational graph of Purdom's algorithm.

Computing the transitive closure of

5.4.2

Algorithm of Transitive Closure

graph 𝐺 = 𝑉, 𝐸 means obtaining graph

The transitive closure computation problem

𝐺 ! = 𝑉, 𝐸 ! , where an edge 𝐸 𝑣, 𝑤 fromin directed graph G can be solved using

𝐺 belongs to 𝐸 ! if and only if vertex 𝑤 isthreereachable form vertex 𝑣 in graph 𝐺. As a

described below.

result,thetransitive

closureproblem1.

different

traditional

approaches,

The transitive closure computation can be

solution requires the |V|! storage space.

reduced to the shortest paths computation

Thus, the memory capacity of a node of

in a corresponding graph with identicalmodern computers is not enough even forweights. Consequently, it can be solved

medium-sized graphs (starting with around

with

2!"

introduced in [3] and [4]. This algorithm

vertices).

generalization

For
of

the

this

reason,transitive

theclosure

has

the
the

Floyd-WarshallO(|V|3)algorithm,

computational

problem is used: only the specified pairs ofcomplexity, and historically is the firstvertices (𝑢! , 𝑣! ) , (𝑢! , 𝑣! ) , …, (𝑢! , 𝑣! ) aredeveloped algorithm for the transitive

checked to belong to the transitive closure.

closureThe number of the pairs required to check 𝑛important property of the Floyd-Warshallbecomes an additional flexible algorithm

algorithm

parameter: varying it may greatly affect therequirement for computations, which

overall algorithm performance.

immediately reduces its applicability only

problemissolution
O(|V|2)to small-scale graphs.

14

[5].

An

memoryJoint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Table 5.2 The comparison of sizes between original input graphs and corresponding intermediate
representation graphs.

2.

The transitive closure can be obtained

Tartan’s algorithm (O(|E|)) is used forusing several multiple BFS, executed fromstrongly

each vertex of the graph. The BFS

computation; if another algorithm, suchalgorithm has been first described in [6];as the DCSC is used, computational

thecomplexity can be different.

algorithm

allowschecking

connected

components

reachability between the selected source
vertex and other vertices of the graph
with

only

required3.

O(|E|+|V|)(usinga

5.4.3

operations

Implementation details of Transitive
Closurequeue-based

Based

on

computational

complexity

implementation approach). As a result,estimates, the most suitable approach tothe full transitive closure computation

solve the transitive closure problem isrequires O(|V| ∗ (|V| + |E|)) operations

Purdom’s algorithm. However, during theAmong the reviewed approaches, theselection of the most suitable algorithm formost optimal computational complexity

particular

has Purdom’s algorithm, introduced in [7].

properties of the parallel algorithm have to

Purdom’s algorithm is based on the

be studied. Figure 5.13 demonstrates anfollowing idea: the transitive closure

informational graph of Purdom’s algorithm.

computation for graph G can be reduced

The presented graph is rather complicatedto the transitive closure computation forand therefore includes two subgraphs, eachgraph G−, obtained from graph G byone

collapsingG’sstronglyalgorithm building-block: the BFS and SCC

components

into

G−

connected

vertices.

based

on

an

importantdivided into 4 separate stages (steps):
1. Detecting

u is the number of edges in graph G−. The
is

to

the

Thus, the developed algorithm can beu|V|) computational complexity, where
estimate

correspondingarchitectures,computations.

The

described approach provides O(|E|+providedparallelthestrongly

connected

components in input graph G,

assumption that asymptotically optimal

2. Creating intermediate representation
15

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure 5.14 Performance (in MTEPS) of multi core implementations of the Purdom’s algorithm
for NEC SX, Intel KNL, and Skylake architectures.

Figure 5.15 Performance (in MTEPS) of single core implementations of the Purdom’s algorithm for
NEC SX, Intel KNL, and Skylake architectures.

no repeated SCC computation.

graph G−,
3. Computing an answer for all vertex

The edges sorting is not required forpairs, related to the same stronglyintermediate representation graph G−, sinceconnected component,

usually these graphs have a significantly

4. Computing an answer for the rest

smaller number of vertices compared to thepairs, using parallel BFS in theoriginal input graphs (while it is still used forintermediate representation graphthe original graph G in order to improve the

G −.

SCC

If

required,anintermediate

step

performance).

corresponding

reachabilityAsaresult,arraysofrepresentation graph can be saved to hardintermediate representation graphs usually

drive after stage 1. Then, the transitive

fit into caches. Table 5.2 demonstrates

closure can be found more efficiently due tocomparison between sizes of original input
16

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Figure 5.16 Parallel acceleration and parallel efficiency comparison for NEC SX, Intel KNL, andSkylake architectures.

graphsandgraphsofcorresponding

intermediate representations.

[2] Beamer, Scott, Krste Asanovi , andIn the loops, the container functions are
utilized

andvectorization.

these
By

BHV, St. Petersburg (2002). (in Russian)functionsavoiding

prevent

using

the

David

Patterson.

"Direction-optimizing

breadth-first search." Scientific Programming
21.3-4 (2013).

functions, the compiler could vectorize the[3] Zhong, Jianlong, and Bingsheng He.

code.

"Medusa: Simplified graph processing on
GPUs." IEEE Transactions on Parallel and

5.4.4

Evaluation of Puron’s algorithm on NEC,Distributed Systems 25.6 (2014).

KNL, and Skylake systems[4]

Besta,

Maciej,"SlimSell:

A

Representation

al.

for

Vectorizable

comparison

betweenmulticore

Breadth-First Search." Parallel and Distributed

implementationsfordifferent

Processing Symposium (IPDPS), 2017 IEEEarchitectures.

Figurethree
5.15

demonstrates

Graph

et

Figure 5.14 demonstrates the performance

International. IEEE, 2017.

performance per core comparison. Figure[5] Fleischer, Lisa K, Bruce Hendrickson,

5.16 demonstrates the parallel acceleration

and Ali Pınar. «On Identifying Strongly

and parallel efficiency. These figures show

Connected Components in Parallel». In

that the performances of NEC SX and KNL

Lecture Notes in Computer Science, Volume

are not high compared to Skylake. In this

1800, Springer, 2000, pp. 505–11.

algorithm, the deep analysis and further[6] Hong, Sungpack, Nicole C Rodia, andoptimization are mandatory to exploit the

Kunle Olukotun. «On Fast Parallel Detection

potential of NEC SX-ACE.

of Strongly Connected Components (SCC) in
Small-World

Graphs»,ProceeedingsofReferences

SC'13, 1–11, New York, New York, USA:

[1] Voevodin, V.V. Parallel Computing. 608p.
17

ACM Press, 2013.

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018[7] Jiˇr ı Barnat, Petr Bauch, Lubosˇ Brim,

Fujimoto, Hiroyuki Takizawa, Akihiro

and Milan Cˇesˇka. «Computing Strongly

Musa

Connected“Performance and Power Analysis of

Components

inParallelon

CUDA». Faculty of Informatics, Masaryk
University,

Botanicka

68a,

60200

and

SX-ACE

Brno,

Czech Republic.

Hiroakiusing

Kobayashi,

HP-X

Benchmark

Programs,”

Proceedings

of

IEEE

International

Conferenceon

Cluster

Computing 2017, pp. 693-700, Sep. 2017.
3. Hiroyuki Takizawa, Thorsten Reimann,

6.

Progress of FY2017 and Future Prospects

Kazuhiko

As we planned, four graph algorithms have

Ryusuke Egawa, Akihiro Musa andbeen implemented on NEC SX-ACE andHiroaki Kobayashi, “Vectorization-aware

conducted

performance

Loop Optimization with User-definedconsideringtheanalysischaracteristics

byofTakashi

Soga,

Code Transformations”, Proceedings ofeachIEEE International Conference on Cluster

algorithms. Through the initial analysis ofComputing 2017, pp. 685 – 692, Sep. 2017.

the first versions of implementation, the

4. Hiroyuki Takizawa, Kenta Yamaguchi,

optimization candidates are clarified. Then,

TakashiSoga,

Thorsten

Reimann,

optimizations for NEC SX-ACE have been

Kazuhiko Komatsu, Ryusuke Egawa,

applied to these four graph algorithms. By

Akihiro Musa and Hiroaki Kobayashi,

appropriate optimizations, the potential of"Migrating an old vector code to modern

the vector data processing to the graphvectormachines,"

acceptedforpresentation at The 30th International

problems is clarified.

Conference on Parallel Computational

For future work, further optimizations

Fluid Dynamics, 2018.

for Transitive Closure, which is one of the

5. Hiroyuki Takizawa, Thorsten Reimann,graph algorithms, need to be performed.

Kazuhiko

Komatsu,

Takashi

Soga,

Furthermore, as there are other graph

Ryusuke Egawa, Akihiro Musa, andalgorithms that are not implemented on a

Hiroakivector architecture, these graph algorithms

differences

are the target to find more vector-friendly

between Intel Knights Landing and NEC

Kobayashi,

SX-ACE

algorithms.

7.

Komatsu,in"Expressing

code

processors,"

theoptimizations
acceptedforpresentation at the 13th World Congress
on Computational Mechanics / 2nd Pan

List of Publications and Presentations

American Congress on Computational(1) Journal Papers

Mechanics, 2018.

(2) Conference Papers

6. Afanasyev

1. Ilya Afanasyev, Vladimir Voevodin, “Thecomparisonoflarge-scale

Ilya,

“An

Efficient

Implementation of the Transitive Closure

graph

Problem on Intel KNL Architecture,” theprocessing algorithms implementation

3rd

methods for Intel KNL and NVIDIA

Distributed, and Cloud Computing forGPU”,

large-scale

graph

processing

7. Afanasyev

2. Ryusuke Egawa, Kazuhiko Komatsu,
Isobe,

Toshihiro

Kato,

WorkshopYoung Scientists, 2017,

Workshop, Sep. 2017.
Yoko

Ural

Comparison

Souya

I.,

on
. 10-19.

Voevodinof

Parallel,

Vl,

Large-Scale

“The
Graph

Processing Algorithms Implementation
18

Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

Methods for Intel KNL and NVIDIA,”

efficiency,” Invited talk at the 27th

Supercomputing. RuSCDays 2017, series"Communications

inComputerWorkshopandp. 80-94.
(3) Oral Presentations
1. IlyaAfanasyev,“Efficient

Implementation of the Transitive Closure
Problem on Intel KNL Architecture”, 3rd
Ural Workshop on Parallel, Distributed,
and

Cloud

Computing

for

Young

Scientists (Ural-PDC 2017), Oct. 2017.
2. Hiroaki

Kobayashi,“Two-Year

Experiences with Vector Supercomputer
SX-ACE and Design Space Exploration of
the Next Generation Vector System”,
KeynoteSpeechin

RussianSupercomputing Days 2017, Sep. 2017.
3. Alexander

Antonov,

Jack

Dongarra,

DmitryNikitenko,

Pavel

Shvets,VladimirVoevodin,Architectures

RequireChangesin

Changes

in

Implementation, not in

Algorithms,”

Russian

Days,

Supercomputing

Sep.

2017.
4. Vladimir Voevodin , The AlgoWiki
Project

andChallenges

of

the

Well-Known Area,” Invited talk at the
3rd

InternationalHigh

Performance

Computing Forum, Guanzhou, China,
Sept. 2017.
5. Vladimir

Voevodin.,“The

AlgoWiki

Project: An Algorithmic Pillar of ExascaleComputing,”

Invited

talk

at

ATIP

workshop "International exascale andnext-generation computing programs",
Nov. 2017.
6. Vladimir

Voevodin,“Algorithms,

Computing Platforms and Unlimited
Freedom of Comparison,” Invited talk at
the

27th

WorkshopononSustained

Performance, March 2018.
(4) Others

Information Science", Springer, Vol. 793,

Sustained

Simulation Performance, March 2018
7. Vadim Voevodin, “HPC software formass analysis of parallel applications
19

Simulation