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

High-performance Randomized Matrix Computations
for Big Data Analytics and Applications
Takahiro Katagiri (Nagoya University)

Abstract
In this project, we evaluate performance of randomized matrix algorithm for big data
analysis with several kinds of supercomputers with providing JHPCN. Target of thealgorithm is singular value decomposition (SVD) and related eigenproblem. In thisyear, we focus on evaluating of prototyping for SVD. The integrated SVD (iSVD),
which is a new algorithm by members of this project, is evaluated with provided
supercomputer resources.

1.

JHPCN joint research project are two-folded.

Basic Information

First, the project aims at dealing with very

(1) Collaborating JHPCN Centers
This research is planned to use multiple

large-scale numerical computations and data

supercomputer resources in the centers ofprocessing. Supercomputers offered by the

JHPCN. The consideration for division of joint

JHPCN are critical and extremely necessary. On

research are thus covered. Necessity of multiple

with the accesses of these supercomputers, wesupercomputer in the proposal is summarized ascan (i) deploy the proposed random sketching

follows.

algorithms

GPUs for TSUBAME2.5

Implementationalarge-scaleparallelenvironments, (ii) develop and study autoon

tuning schemes of the computer codes for theselarge-scale cluster for SVD (singularadvanced supercomputers to achieve highvalue decomposition) and LS (Linear

computing capabilities, and (iii) handle the very

System)large-scale data sets which are generated from

1-CPU-1-GPU; multiple-CPU-Multiple-

important and timely scientific, engineering,GPU;social network applications.

Multiple-node:

CPU

of

parallelism

cluster(MPI+

Second, the project can be conducted only if

OpenMP);

inresearchers

withvariousexpertisework

Advanced CPUs for the CX400, thetogether. A strong international team has been

FX10 and the FX100.

formed to conduct this project. The team

Performance comparison with the FX10 andmembers include experts in high-performancethe FX100 with respect to hardware

computations,improvements, such as bandwidths forstatistical analysis and computations. The

memory accesses.

support from JHPCN is essential to establish

Supercomputer adaptation for SVD and LS.

the collaborations of these researches to achieve

Performance evaluation and profiling withthe goals of the project.

numerical

linearalgebra,Fujitsu performance profiler.


Evaluating auto-tuning effects.

(2) Research Areas

The necessity for implementing the project as a1Very large-scale numerical computation

☑

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

Matrix computation is one of the most

Very large-scale data processingVery large capacity network technology

important computational kernels in manyVery large-scale information systems

applications of numerical computations anddataprocessing.

Becausematrixrepresentations are almost everywhere in(3) Roles of Project Members
Takahiro Katagiri researches large scalescientific simulations, engineering innovations,

implementation, and adaptation of auto-tuning.

data

analysis,statisticalinferences,andWeichung Wang researches parallel algorithm

knowledge extractions, just to name a few, thesedevelopment, surrogate-assisted turning, andtimely applications lead to strong needs ofbig data applications.

efficient parallel numerical solvers based on

Su-Yun Huang researches random sketchingalgorithm

development,

mathematicalmatrices. However, few numerical solvers,

andespecially randomized algorithms, are designed

statistical analysis.

to tackle very large-scale matrix computations

Kengo Nakajima provides knowledge of highperformance

implementationofon the latest supercomputers. This project plans

iterative

to fill the gap by focusing on algorithmic andmethods.

software aspects.

Osni Marques provides knowledge of parallel

In the algorithmic side, we intend to develop

eigenvalue algorithms. He provides knowledge

efficientof SVD and implementation of numericalapproximate SVD and LS solutions of large-

libraries.

scale matrices. The main idea is to sketch thesketchingschemes

tocomputeFeng-Nan Hwang provides knowledge ofmatrices by randomized algorithms to reduceparallel eigenvalue algorithms. He provides

the computational dimensions and then suitably

knowledge

integrate the sketches to improve the accuracy

of

SVD

and

its

parallelimplementation.
ToshioEndoand to lower the computational costs. Suchprovidesknowledgeoftechniques will benefit many scenarios that low-

optimizations for hierarchical memory andrank SVD approximation or approximate LS

adaptation of its auto-tuning. He provides

solutions with quick sketches are sufficient. Our

knowledge of hierarchical memory optimization.

schemes will also benefit the applications that it
is not possible to compute full SVD or highly-

2.

Purpose and Significance of the Research

accurate LS solutions provided the matrices areThis project aims at developing random

very-large or the matrices are provided in thesketching algorithms with high-performance

form of streaming.

implementations on supercomputers to compute

In

the

softwareaspect,we

intendtoalgorithmson

SVD and LS solutions of very large-scale

implement

matrices.

resultingsupercomputers. Consequently, the softwarealgorithms and packages can be used to solve

packages can be used to solve very large-scale

important and demanding problems arising inreal world problems and advance the frontier ofvery large-scale numerical computations andsciencedataessential component of this project is to develop

Consequently,

processing

nowadays.

that

the

remainchallenges
2

and

the

proposedtechnology

innovations.

Oneeffective automatic software auto-tuning (AT)Joint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018technologies, so that the package can fully takesoftware can fully take advantage of theadvantage of the computational capabilities ofcomputational

the target supercomputers that include CPUsupercomputers.

homogeneous and CPU-GPU heterogeneous

capabilities

of

the

target(3) Implications of applications with large-

parallel computers.

scale data sets: We expect the resulting
algorithms and software packages will impact

Significance as a JHPCN Joint ResearchseveralProject

Singular value decomposition of matrices and

The significance is summarized as follows:

solutions of linear systems are two essential(1) Novel and efficient randomized sketchingcomponents of big data analytics and various

3.

important

big

data

applications.

algorithms: We anticipate the international

practical

interdisciplinary joint efforts in this project will

include unsupervised dimension reduction (e.g.

lead to highly efficient randomized sketching

matrixalgorithms for computing very large-scale

analysis, spectral clustering) and supervised

singular value decompositions and for solving

dimension reduction (e.g. linear discriminantlarge-scale linear systems. We have observed

analysis, canonical correlation analysis, inversethat the size of matrices increases rapidly

regression),

nowadays. These large matrices arising in, forsimulations in geophysics, energy, human

example, finer meshes in discrete differential

genetic variation, and many others.

applications.
factorization,

machineThese
principal

applications
component

learning,

numerical

equations, larger number of sensors, andcollections of activities on the Internet. To tackle

4.

these large-scale problems, random sketching is

Outline of the Research Achievements up
to FY2016

one of the most powerful approach. Our

Randomized Algorithm (A Prototyping)
Fig. 1 shows the algorithm of Integrated SVD

integration type algorithms are expected towith multiple sketches (iSVD)

improve the accuracy and accelerate theWe have parallelized iSVD with MPI. The line

computation on parallel computers.
software3 in Fig. 1 can be parallelized easily. We need to

packages with auto-tuning mechanisms: The

distribute Ω[i] to each process in the line 1 byproposedbeusing MPI_Alltoall. We also need to gather

implemented on the JHPCN supercomputers.

distributed Q[i] by using MPI_gather. In thisWhile the data size keeps increasing rapidly,

version, matrix A is not distributed.

(2)

High-performancealgorithms

arescalable
intendedtofortunately, the number of computing nodes also
keeps increasing and more and more recent
parallel computers are equipped with many-core
co-processors. It is thus essential to develop
scalable algorithms to take advantage of suchlatest computer architectures for tackling verylarge-scale problems. On the other hand, thanks

Fig.1 Algorithm of Integrated SVD withto the auto-tuning technologies to be studied andmultiple sketches (iSVD)

deployed within the software, we expect the3

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

Accuracy evaluation of iSVD



analysis. After permutation and partition by H-

In interim report, we reported parallel

matrices, which is one of other low rank

performance of iSVD. In this report, we evaluate

approximations for matrices, we obtain a target

accuracy of iSVD in viewpoint of low rank

sub-matrix. We use the sub-matrix with N = 300

approximation to show ability of performance forx 300 for the test matrix.
Error of the sub-matrix can be evaluatedprocessing big data.

with the following three viewpoints:

Our computer environment in this experience
is as follows:

The CX400 (ITC, Nagoya U.)

A − UΛ V T

Node: 2 sockets of Xeon E5-2697v3(Haswell) (2.6GHz).

Interconnection: Infiniti Band, full

/ A F,TF

i =1Intel composer_xe_2015.5.223, MKL



T

i =1

Intel MPI version 5.0.3.049.

/ A F )i / T ,(2)

Square root (Sqrt):∑ ( A − UΛ V

Intel compiler: version 15.0.4.

(1)

Mean:
T

with same version is used.
F∑ ( A − UΛ V

connection.


Frobenius Norm (FB):/ A F − mean )i ,
2

T
F

(3)

where T is the number of test trials in iSVD.

The test matrices used in here are as follows:
Matrix 1: A Randomized Matrix

In this matrix 1, matrices U and V are made



Results

by using normal random generator [0,1]. A is

Matrix 1

made by using the generated U and V with A =In this chapter, we evaluate errors of low rank

U S V^T, where, singular value (true answer) isapproximation by iSVD. Fig.2 shows errors by

set to i = 1 ~ k : 1.0 / i, i = k + 1 ~ m : 1e-2 / (i -Matrix 1 with respect to parallel execution.

k). Hence, true answer of U and V is the U and

Fig.2

V.



the

accuracy

improves

according to the number of cores. Hence theWe can check the answer with:

indicatesiSVD has scalability of accuracy. Note that theCompute mean value of diag(S’’), (U U’^T) =execution time is almost same, but the number

U’’ S’’ V’’^T, where, A = U’ S’ V’^T.

of samplings for low rank approximations are

Error(max, mean, min): sqrt (diag(S’’) -increasing to the number of cores.

mean).


Matrix2:

H-Matrix

(Low

Rank

Approximation)
This matrix 2 is an example to show abilityof low rank approximation by using a wellapproximated problem. The projection withrandomized

algorithm,

hence,

canbeestablished in theoretically.
We use a problem from static electric field

Fig.3 Accuracy of the Matrix 1.
4

Matrix 2:

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

Fig. 3 shows errors of low rank approximation
by equation (1)-(3). The number of trials T is set
to 100. The k is dimension number of low rank
approximation.
Fig. 3 indicates that accuracy depends
dimensions

of

row

rank

approximation.

According to Fig. 3 (a), if we accept accuracy
around 2e-02, we can reduce matrix dimension(a) Case of k = 3x3

from 300x300 to 3x3. This is 1/100 scale. Hence
this means 1/1003 reduction of computation
complexity, since computation complexity ofdense eigenproblem is O ( n3 ).
As summary of this experiment, we conclude
that iSVD provides reasonable accuracy inviewpoint of low rank approximation. Thismeans that iSVD is one of candidates to process
big data analysis for the eigenvalue problem.
Several publications for the other topics are(b) Case of k = 10x10

now under preparation.

5.

Details of FY2017 Research Achievements

The original proposal is planned with thefollowing three years.
Year 1 (FY2016): Algorithm development
and testing environments deployment. (Aprototyping)



Year

2

(FY2017):

Large-scale(c) Case of k = 30x30

implementation and software integrations.
Year 3: Auto-tuning of large-scale codes andtests of applications.

(d) Case of k = 300x300
Fig.3 Accuracy of the Matrix 2.
5

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

Mathematical and statistical algorithms:

Parallel Implementation of iSVD

In the Year 1, we have evaluated several

In

this

year,

a

prototypingforparallelcontents as follows: (1) SVD and LS based on

implementation has been finished for Prof.

random projection and random sampling; (2)Wang’sAlgorithms integrating multiple randomized results

implementation for iSVD. Main contribution isbased on optimization and numerical linear

to parallelize input matrix A with row-block

algebra techniques; (3) Algorithm analysis from

distribution with parallel reduction for MPI to

viewpoints of geometry, statistics, and complexity;

reduce communication time. Fig. 4 shows a(4) Preliminary numerical experiments for the

typical parallel performance for iSVD with row-

purpose of proof-of-concept.

block distribution.

team

to

establish

large

scale

Parallelism on large-scale cluster:In the Year 2 (this year), we focus on thefollowings: (1) Algorithm parallelism including task
separation, task distribution, data structure, data
communication; (2) Parallel implementations: (2a) Single-node: 1-CPU-multi-cores (OpenMP);
Multiple-CPU (OpenMP); 1-CPU-1-GPU; multipleCPU-Multiple-GPU; CPU with Xeon Phi; (2-b)
Multiple-node: CPU cluster (MPI + OpenMP);
CPU-GPU cluster (MPI + OpenMP);

Fig. 4 A Typical Parallel Performance for

Auto-tuning and applications:

iSVD with row-block distribution.
According to Fig. 4, communication time forIn the end of the Year 2 and the Year 3, we
focus on the followings: (1) Adaptation of general

Row-Block distribution can be dramatically

performance model for AT, in particular, surrogate

reduced according to number of processors,models; (2) Development of auto-tuning

while communication time in conventional

methodologies for computation kernels and MPI(naive)communications of SVD and LS, in particular,

because, parallel reduction operation based onhierarchical memory optimizations; (3) Extension

row-block distribution is crucial to conventionalof AT functions for ppOpen-AT with respect toimplementation.

implementation

increases.

This

is

nature of processes for SVD and LS. (4) Human
I/O evaluation

Genetic Variation; Underground water; Airpollution; Images analysis;

To treat with big data, I/O time cannot beignored in general. Even in iSVD, I/O time

Progress of FY2017 and Future Prospects

inside processes will be a serious problem if we

We have several results for eigenvalue

do not implement efficient I/O. To show evidence

program and sparse linear equations solvers inof this, breakdown among each procedure andthis project. Due to page limitation, we only

I/O time for previous version of iSVD is shownshow result of iSVD.

in Fig. 5.

6.

In Fig. 5, 50.1% to total time is occupied with
6

I/O execution in this environment. It should beJoint Usage/Research Center for Interdisciplinary Large-scale Information Infrastructures
Final Report for JHPCN Joint Research of FY 2017, May 2018

parallelized for the I/O process. To show effect ofiSVD. We apply iSVD algorithm to Facebook

the parallelization, a test for parallel I/O in the

100k data. This application requires 108,585 x

example in Fig. 5 is shown in Fig. 6.

108,585 matrix 𝑨𝑨. For the iSVD, we summarize
settings as follows.


Row-block

Gaussianprojection

sketching (CPU/GPU),

Row-block Gramian orthogonalization,Row-block Wen-Yin integration, andRow-block Gramian former.

We use C++ for iSVD. Target machine is theReedbush-H at ITC, the University of Tokyo.
Using number of nodes is 4 nodes, 144 cores,
includingFig.5 A Breakdown of Execution Time of iSVD

8

GPUs

(Pascal).

Numerical

parameters of iSVD are set to 𝑁𝑁 = 256, 𝑘𝑘 = 20, 𝑝𝑝

in one node of the Oakforest-PACS.

= 12, 𝑃𝑃 = 8.

The size of I/O is 24 GB. The size of input

For MATLAB execution, we use a PC cluster

matrix is 5,000 x 200,000 , and size of sampled(WLabmatrix is 500 x 500.

Cluster)

with

1

Node,

24

cores.

Numerical parameters of 𝑁𝑁 = 256, 𝑘𝑘 = 20, 𝑝𝑝 = 12.

Fig. 7 shows total execution time in seconds

between

MATLAB

execution

and

iSVD

execution.

1x
10x

Fig.6 Parallel I/O Effect with 4 MPI processes.

Fig. 7 Total Runtime of Facebook Data.

According to Fig. 6, we can establish 4x

Time in seconds.

speedup of I/O with 4 MPI processes. This is

According to Fig.7, our developed parallel

perfect scaling. Please note that the result is not

iSVD has a merit from 10x to 19x. Hence, we canusing fast I/O technologies, such as burst buffer.

conclude that our developed iSVD code is very

Hence, we have enough room to speed up the I/O

useful for actual application.

implementation.


19x



Application Adaptation

Future Prospects

There are several topics to progress the research

We have adapted several applications, but due

in this financial year. The research topics areto page limit, we show an adaptation example of7summarized as follows.



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

Detailed

evaluation

supercomputers,

and

on

severalimprovement

Conjugate Gradient Eigensolvers”, DPMAT2017,of(Sep. 2017)[5] Takahiro Katagiri, Mu Yang, Weichungimplementations for MPI.


Adaptation of auto-tuning for selection of

Wang, Akihiro Ida, “Performance Evaluation of

implementations that are inside processes

Integration of Multiple Randomized Low-Rank

in iSVD. However, this is main research

Singular Value Decompositions”, 2017 年並列／

topic for the Year 3.

分散／協調処理に関する『秋田』サマー・ワーク

Additional adaptations and evaluations toショップ (SWoPP2017), 日本応用数理学会 「行many application software.

列・固有値問題の解法とその応用」研究部会(MEPA) (Jul. 2017)



[6] Weichung Wang, “Parallel Singular Value

Related Workshop

In this year, we held in an international

Decomposition via Multiple Random Sketches forconference related to the project, named Second

Large Matrices”, International Conference on

International

Econometrics and Statistics (EcoSta 2017) (JuneWorkshoponDeepening

Performance Models for Automatic Tuning

2017)

Nagoya[7] Weichung Wang, “Parallel SVD by Random

University. Many members in this project

Sketches and Big Data Intelligent Analytics Forumpresented their work in DPMAT2017. See theon Artificial Intelligence for Data Science”,

attached list of publications.

National Sun Yat-sen University, Kaohsiung,(DPMAT)in

September

2017

at

Taiwan (Dec. 2017)

7.

List of Publications and Presentations[8](1) Journal Papers

Decompositions for Photonic Device Simulations

None.

Weichung

Computational[1] Weichung Wang: “Parallel Singular Valueof Tokyo (Dec. 2017)

Decomposition for Large Matrices by Multiple
Random Sketches”, DPMAT2017, (Sep. 2017)
Takahiro

Katagiri:“Auto-tuning(4) Others
to

None

Scientific Applications –Traditional Approachby Code Transformations and New Approach by
AI –”, DPMAT2017, (Sep. 2017)
Feng-Nan

Preconditionerfor

Hwang,“Nonlinear

Full-space

Lagrange-

Newton-Krylov Algorithms with applications in
Large-scale

PDE-constrained

Optimization

Problems”, DPMAT2017, (Sep. 2017)
Osni

Marques,

Science(RECS),

The

Computational Science Alliance of the University(3) Oral Presentations[4]

MatrixSymposium on Research and Education of

None.

[3]

“Efficient& Big Data Analytics”, Second International(2) Conference Papers

[2]

Wang,

“Unconstrained

Functionals for Efficient Parallel Scaling of8