University	  of	  Tsukuba	  	  | Center	  for	  Computational	  Sciences

Numerical	  Computation
OpenCL-­‐ready	  High	  Speed	  FPGA	  Networking	  [1]Implementation overview
I/O channel
specified access

additionally implemented

QSFP+ port
calibration is needed

FPGA board (A10PL4)

OpenCLkernelInterconnect

FPGA

Ethernet IP
Controller

Ethernet
IP CoreQSFP+

Ethernet IP
Controller

Ethernet
IP CoreQSFP+Our FPGA board(BittWare A10PL4)


QSFP+
DDR4
mem

DDR4
Controller
DDR4
Controller
PCIeController

DDR4
mem

BSP

DriverHost PC

Host Application

14• Problem space is divided into small blocks• e.g. (16, 16, 16) → 8  (8, 8, 8)• PE is assigned to each of small blocks)

MemoryReader

444

Buffer(PE Array(2x2x2)

444

Table 2: Resource usage and clock frequency of t
DDR4
Memory4424)

Memory Network
Channelbetter

MemoryWriter(16x16x16)

size
1200
(16, 16, 16)
1000(32, 32, 32)
800(64, 64, 64)
600 128, 128)
(128,1400

44

Buffermesh• ART is ray tracing based algorithm• problem space is divided
into meshes and reactions
are computed on each mesh
• Memory access patterndepends on ray direction• Not suitable for SIMD architecture

Size	  XS,	  S	  and	  M	  on	  1,	  2	  and	  4	  FPGAs(Strong	  Scaling)

DDR4
memAuthentic	  Radiation	  Transfer	  (ART)	  on	  FPGA	  [2]• Accelerated Radiative transfer on grids Oct-Tree(ARGOT) has been developer in Center for
Computational Sciences, University of Tsukuba• ART is one of algorithms used in ARGOT anddominant part (90% or more of computationtime) of ARGOT program

3D-­‐Poisson	  Equation	  Solver,	  19-­‐point	  stencil	  computation

Performance [M mesh/s]lDev.	  tooI:	  Intel	  FPGA	  for	  OpneCLlBoard	  Support	  Package	  (BSP)	  is	  a	  hardware	  component	  to	  support	  multiple	  different	  boardsØWhich	  FPGA	  chip	  is	  used	  on	  the	  board	  ØWhat	  kind	  of	  peripherals	  are	  support	  by	  the	  board	  lBasically,	  only	  minimum	  interfaces	  are	  supported	  ØTo	  perform	  inter	  FPGA	  communication,	  implementing	  network	  controller	  and	  integrating	  it	  into	  the	  BSP	  are	  required• Our implementation uses channel based approach• One of extensions to OpenCL for FPGAs by Intel400
Fig. 5: Design Outline of ART on FPGA.
• It enables inter kernel communication much faster200• No external memory (DDR) access is required0• Lower resource
utilization
thancomputes
DDR accesseach other.
Each kernel
reaction between a mesh and a ray

on its own computation space which is dedicated to each kernel. While
computing, a ray is traversed among multiple compute kernels depend• PEs are connected by channels each other
Global Memory
ing on its location. If a ray goes out from kernel’s space, its data will be(DDR4)• PE: Processing Element
Off Chip
Source kernel through a Destinationtransferred
to a neighborchannel.
• BE: Boundary Element
KernelKernelFigure 5 shows the design outline ofour implementation.
“Memory Reader”
FIFOWriteRead• Kernel of PEs and BEs are started automatically by
Channel which is seen as a global memory
reads mesh data from DDR4 memory
autorun attribute
SourceDestinationfrom OpenCL language. “Memory Writer” is a counterpart to the reader
KernelKernel• Lower control overhead and resource usage
and updates mesh data by the result of computation. It has both of readbecause of decreasing number of host controlled
and write memory access because it computes integration of gas reaction.
kernels
without channels
withtochannels“Buﬀer” is a mesh data buﬀer
improve memory access performance.
“PE Array” is an array of PEs (Processing Element). PE computes theReference	  kernel of ART method. The array is consists of multiple kernels. We show[1]	  Ryohei Kobayashi,	  Yuma	  Oobata,	  Norihisa Fujita,	  Yoshiki Yamaguchi,	  and	  Taisuke Boku,	  OpenCL-­‐ready	  High	  Speed	  FPGA	  Network	  for	  Reconfigurable	  High	  Performance	  Computing,	  HPC	  Asia	  2018,	  pp.192-­‐201,	  January	  2018
the detail of PE network in the next subsection.
[2]	  Norihisa Fujita,	  Ryohei Kobayashi,	  Yoshiki Yamaguchi,	  Yuuma Oobata,	  Taisuke Boku,	  Makito Abe,	  Kohji	  Yoshikawa,	  and	  Masayuki	  Umemura:	  Accelerating	  Space	  Radiate	  Transfer	  on	  FPGA	  using	  OpenCL	  (Accepted),	  International	  Symposium	  on	  Highly-­‐Efficient	  Accelerators	  and	  Reconfigurable	  Technologies	  (HEART	  2018)
Since our implementation is work-in-progress, it lacks some features from
Acknowledgment	  the CPU implementation. While computation in an FPGA, all mesh data
This	  research	  is	  a	  part	  of	  the	  project	  titled	  “Development	  of	  Computing-­‐Communication	  Unified	  Supercomputer	  in	  Next	  Generation”	  under	  the	  program	  of	  “Research	  and	  Development	  for	  Next-­‐Generation	  Supercomputing	  Technology”	  by	  MEXT.	  We	  thank	  Intel	  University Program	  for	  providing	  us	  both	  of	  hardware	  and	  software.
must be put into its internal BRAM (Block Random Access Memory).
The FPGA implementation does not support to replace mesh data involved by progression of its computation. Therefore, problem size which
an FPGA can solve is limited by the size of BRAM. The CPU implementation supports inter-node parallelization using MPI (Message PassingInterface), but the FPGA implementation does not support any networking functionality and uses only one FPGA.
BEBEBEPEPEBEBEPEPEBE

96bit x2
(read,write)
ChannelyBE

(%)
31%
40%
40%
40%

Registers (%)
267,828 31%
344,447 40%
CPU(14C) 40%
344,512CPU(28C) 40%
344,505
P100(x1)

BERay Data

M20K
739
796
796
796

FPGA

Table 3: Performancecomparison
between(128,128,128)
FPGA, CPU and(16,16,16)(32,32,32)(64,64,64)
GPU implementations. Themeshunitsizeis M mesh/sec.

(8x8x8)xequal	  performance# ofalmost	  PEsALMs
(2, 2, 2) 132,283(2, 2, 2) 169,882
4.9x	  faster(2, 2, 2) 169,549(2, 2, 2) 169,662

Size(16,16,16)(32,32,32)(64,64,64)(128,128,128)

CPU(14C)
112.4
158.9
175.0
95.4

CPU(28C)
77.2
183.4
227.2
165.0

P100
105.3
490.4
1041.41116.1FPGA
1282.81165.21111.01133.5weak point
through PC
run our cod
ation Intel F
than Arria
blocks and 3

9. REFE[1] K. M. G
F. K. Ha
Healpix:per link) multiple interconnection links (up to 4 channels) onand fastit. Additionally, HLS such as OpenCL programming enviAstrophyronment is provided, and there are several tyeps of research[2] K. Hill,
to involve them in FPGA computing. In [3], Kobayashi, etComparal.
show
thebasicfeature
toutilize
thehighspeedintercon4.2 Parallelization using Channel in an FPGA
image-prnection over FPGA driven by OpenCL kernels. Therefore,
We describe the structure in “PE Array” shown in Figure 5. A PE ArrayIEEE 26
Background
althoughour implementationis almost forTheaimofand BEsthisresearchisthetoperformancedevelopof numericalsoftwareis consistsof PEs(BoundaryElement) as shownprojectin Figure6.
It shows the PE Array network on the x-y dimension. We do same
not show as NVIDIA P100 GPU, the overall performance withApplicat
connections
forzdimension
tokeepthefiguresimple.
We
also
have
alarge-­‐scale eigenvalue problems
for nodes
post-­‐petascalecomputing
The fast Fourier transform (FFT) is widely used in science andProcesso
multiple computation
with FPGA to be connected
disimilar connection to x-y dimension for z dimension.
rectly caneasily overcome
the GPU implementation
which SS [3] R. Koba
environment. An eigensolver
based
on contour
integral (theengineering. Parallel FFTs on distributed-­‐memory paralleland T. B
requires host CPU control and kernel switching for intermethod) has been proposed
bycommunication.
Sakurai andSugiuraoverhead[3]. Thismethodreconfigu
computers require intensive all-­‐to-­‐all communication, which
node
Networking
on FPGAs
isProceedi
muchandlower than
one on GPUs.
To massively
improve currentparallelARThasa
hierarchicalstructureissuitablefor
affects their performance. How to overlap the computation andPerform
method implementation with such an interconnection feaAsia 201
supercomputers [2]. Moreover,behighapplicableture on theFPGASS
is method
our next stepcan
towardperformance forthe all-­‐to-­‐all communication is an issue that needs to beACM.
parallel FPGA computing.
nonlinear
eigenvalue
problem[1].
Block
Krylovmethod[4]
improves
addressed for parallel FFTs.
[4] Y. Luo,the performance of the method.
Based on these newly designed G. Mem
8. CONCLUSION
Overview
irregularInthispaper,
weoptimizedART
methodusedinARGOT
algorithms, we have developed a massively parallel software z-­‐Pares case stu
program which solves a fundamental calculation in early
We proposed an automatic tuning (AT) of computation-­‐freely-­‐available from http://zpares.cs.tsukuba.ac.jp/.
MATLAB
version Conferen
stage universe with space radiative transfer
phenomenon,
Applicatcommunication overlap for parallel 1-­‐D FFT. We used aon
anFPGAusing
IntelFPGA
SDKforOpenCL.
We
paris also available in our webpage. We have also developed CISS [5] T. Okam
allelized the algorithm using the SDK’s channel extension
computation-­‐communication overlap method that introduces aaccelerat
eigensolver in SLEPc.
in an FPGA. We achieved 4.89 times faster performance
Monthly
communication thread with OpenMP. An automatic tuning facility
than the CPU implementation using OpenMP as well as al419(4):2
mostsameperformanceas
theGPU
implementation
usingfor selecting the optimal parameters of the computation-­‐[6] S. Tanak
Hierarchical	  Parallel	  Structure CUDA. Numerical	  Example	  on	  the	  K	  Computer
K. Hasecommunication overlap, the radices, and the block size was
Although the performance of the FPGA implementation
Hardware is grouped according to is comparableApplicationforband
calculation
withreal radiation
toNVIDIA
P100GPU,ithasroomtoimimplemented.
a hierarchical structure of the prove its performance.
space density
theory(RSDFT) [2]. Publicat
Thefunctional
most important
optimization
67(4):62
isresourceoptimization.
If
wecanimplement
largernumPerformancealgorithm.
[7] H. R. Zo
ber of PEs than one of the current, we can improve perand S. M
formance. However, it is diﬃcult for us to reduce usage
To evaluate the parallel 1-­‐D FFT with AT, we compared its
Algorithm
Hardwarekernels fof ALMs and registers because we do not describe them diTop	  level
performanceagainstthoseofFFTE
6.2alpha
Proceedi
rectly in OpenCL code. Not only resource but also frequency
Γ2
Γ3
Γ1
Performis important. We suppose Arria 10 with OpenCL design can(http://www.ffte.jp/) , FFTE 6.2alpha with AT, and FFTW
Analysisrun on 200MHz or higher frequency.
z2Middle	  level
USA, 20
3.3.7.
We will implement the network functionality into the ART
z1
design to parallelize it among multiple FPGAs. We conzNThe performance results demonstrate that the proposedsider networking using FPGAs is an important feature for
Bottom	  level
parallel applications using FPGAs. Although GPUs have
implementation of a
Performance  of  parallel  1-­D  FFTs
higher computation
performance
FLOPS
and higherBand	  structure	  of	  silicon	  nanowire	  of	  9mem,924	  atoms.
（ Oakforest-­PACS， 1024  nodes）
parallel 1-­‐D FFT withory bandwidth than FPGAs, I/O including networking is a(matrix	  size	  =	  8,719,488,	  	  Number	  of	  cores	  =	  6,144)	  T(z)X=V4000jj
FFTE
6.2alpha  (nooverlap)
FFTE
6.2alpha(NDIV=4)FFTE
6.2alpha  withAT
FFTW    3.3.7

GG

16644G

1G

6M

M
64

25

M

3500
3000
2500
2000
1500
1000
500
0

16

automatic tuning isefficient for improving
the performance oncluster of Intel Xeon
Phi processors.

GFlopsAutomatic	  Tuning	  of	  Parallel	  1-­‐D	  FFT	  on	  Cluster	  of	  Intel	  Xeon	  Phi	  Processors

Length  of  transform  N

contact	  address:	  pr@ccs.tsukuba.ac.jp

Development	  of	  Parallel	  Sparse	  Eigensolver Package:	  z	  -­‐ Pares*The	  results	  are	  tentative	  since	  they	  are	  obtained	  by	  early	  access	  to	  the	  K	  computer.	  

Reference	  [1]	  J.	  Asakura,	  T.	  Sakurai,	  H.	  Tadano,	  T.	  Ikegami	  and	  K.	  Kimura,	  A	  numerical	  method	  for	  nonlinear	  eigenvalue	  problems	  using	  contour	  integrals, JSIAM	  Letters,	  1	  (2009)	  52-­‐55.	  
[2]	  Y.	  Futamura,	  T.	  Sakurai,	  S.	  Furuya and	  J.-­‐I.	  Iwata,	  Efficient	  algorithm	  for	  linear	  systems	  arising	  in	  solutions	  of	  eigenproblems and	  its	  application	  to	  electronic-­‐structure	  calculations,	  Proc.	  10th	  International	  Meeting	  on	  High-­‐Performance	  Computing	  for	  Computational	  Science	  (VECPAR	  2012),	  7851	  (2013),	  226-­‐235.	  
[3]	  T.	  Sakurai	  and	  H.	  Sugiura,	  A	  projection	  method	  for	  generalized	  eigenvalue	  problems,	  J.	  Comput.	  Appl.	  Math.,	  159	  (2003)	  119-­‐128.	  
[4]	  H.	  Tadano,	  T.	  Sakurai	  and	  Y.	  Kuramashi,	  Block	  BiCGGR:	  A	  new	  block	  Krylov subspace	  method	  for	  computing	  high	  accuracy	  solutions,	  JSIAM	  Letters,	  1	  (2009)	  44-­‐47.	  
Acknowledgment
This	  research	  was	  partially	  supported	  by	  the	  Core	  Research	  for	  Evolutional	  Science	  and	  Technology	  (CREST)	  Project	  of	  JST	  and	  the Grant-­‐in-­‐Aid	  for	  Scientific	  Research	  of	  Ministry	  of	  Education,	  Culture,	  Sports,	  Science	  and	  Technology,	  Japan,	  Grant	  number:	  23105702.	  

http://www.ccs.tsukuba.ac.jp/