Shared Memory Parallelization of MTTKRP for Dense Tensors BLIS - - PowerPoint PPT Presentation

Shared Memory Parallelization of MTTKRP for Dense Tensors

Koby Hayashi, Grey Ballard, Yujie Jiang, Michael Tobia hayakb13@,ballard@,jiany14@,tobiamj@wfu.edu

BLIS Retreat 2017, September 18th

Neuroimaging Application

ensor: Time by Subjects by

• xel Correlation Matrix

est: Rest à Activity à ecovery Subjects: Control, MDD, SAD, COMO Time

Quick Introduction to Tensors

Multidimensional arrays, an N-dimensional nsors is said to be N-way or order-N.

ay 2-way 5-way 3-way 4-way

CP Decomposition

nonical Polyadic composition (CP): Decomposes a tensor into a sum of rank 1 tensors

𝒴≈∑𝑑=0↑𝐷−1▒​𝑣↓𝑗𝑑 ∘​𝑤↓𝑘𝑑 ∘​𝑥↓𝑙𝑑 𝒴≈⟦𝑉,𝑊,𝑋⟧

CP via Alternating Least Squares

Hadamard Product

Element wise matrix product denoted *

𝐷 = 𝐵∗𝐶 ​𝐷↓𝑗𝑘 =​𝐵↓𝑗𝑘 ∗​𝐶↓𝑗𝑘 ​𝐷↓𝑗𝑘 ​𝐵↓𝑗𝑘 ​𝐶↓𝑗𝑘

*

= 𝐽 𝐾 𝐾 𝐾 ↑𝑈 ​𝑉↓0 )∗…∗​(𝑉↓𝑜−1↑𝑈 ​𝑉↓𝑜−1 )∗​(𝑉↓𝑜+1↑𝑈 ​𝑉↓𝑜+1 )∗…∗​(𝑉↓𝑂−1↑𝑈 ​𝑉↓𝑂−1 )

Khatri Rao Product

Khatri Rao Product (KRP):

𝐿 = 𝐵⊙𝐶

• lumn-wise Kronecker Product

𝐿(:,𝑗)= 𝐵(:,𝑗)⊙𝐶(:,𝑗)

r Hadamard Product of Rows

𝐿(​𝑠↓𝐶 +​𝑠↓𝐵 ​𝐽↓𝐶 ,:)= 𝐵(​𝑠↓𝐵 ,:)∗𝐶(​𝑠↓𝐶 ,:) 𝐷 𝐿(​𝑠↓𝐶 +​𝑠↓𝐵 ​𝐽↓𝐶 ,:) =

A(​𝑠↓𝐵 ,:) B(​𝑠↓𝐶

⊙ 𝐷 𝐷 ​𝐽↓𝐵 ∙​𝐽↓𝐶 =​𝒀↓ 𝒀↓(𝒐) (​𝑽↓ 𝑽↓𝑶−𝟐 ⨀…⨀​𝑽↓ 𝑽↓𝒐+𝟐 ⨀​𝑽↓ 𝑽↓𝒐−𝟐 ⨀…⨀​𝑽↓ 𝑽↓𝟏 )

Tensor Fibers

​𝑜=0, 𝒴↓(:𝑘𝑙) ​𝑜=1, 𝒴↓(𝑗:𝑙) ​𝑜=2, 𝒴↓(𝑗𝑘:)

Unfolding Tensors

• The nth mode matricization of a N-

way tensor 𝒴 that is ​𝐽↓0 ×​𝐽↓1 ×…×​𝐽↓𝑂−1 is denoted ​𝑌↓(𝑜) and is ​𝐽↓𝑜 ×​𝐽↓≠𝑜

• ​𝐽↓≠𝑜 =∏𝑜≠𝑙∈[𝑂]↑▒​𝐽↓𝑙
• ​𝑌↓(𝑛:𝑜) denotes a matricization

where {𝑛,𝑛+1,…,𝑜} are the row modes

​𝐽↓≠𝑜 ​𝐽↓𝑜 ∏𝑙={𝑛,𝑛+1,…,𝑜}↑▒​𝐽↓𝑙 ​𝑌↓(𝑜) ​𝑌↓(𝑛:𝑜) 𝐍=​𝒀↓ 𝒀↓(𝒐) (​𝑽↓ 𝑽↓𝑶−𝟐 ⨀…⨀​𝑽↓ 𝑽↓𝒐+𝟐 ⨀​𝑽↓ 𝑽↓𝒐−𝟐 ⨀…⨀​𝑽↓ 𝑽↓𝟏 )

Matricized Tensor Times Khatri Rao Product

𝑁=​𝑌↓(𝑜) (​𝑉↓0 ⨀…⨀​𝑉↓𝑜−1 ⨀​𝑉↓𝑜+1 ⨀…⨀​𝑉↓𝑂−1 ) Naïve algorithm 1. Permute 𝒴 to ​𝑌↓(𝑜) 2. Form K=​(𝑉↓0 ⨀…⨀​𝑉↓𝑜−1 ⨀​𝑉↓𝑜+1 ⨀…⨀​𝑉↓𝑂−1 ) 3. Call DGEMM 1-Step and 2-Step MTTKRP 1. Avoid permuting 𝒴 2. Efficiently form the KRP § 1Step

• (​𝑉↓𝑂−1 ⨀…⨀​𝑉↓𝑜+1 ⨀​𝑉↓𝑜−1 ⨀…⨀​𝑉↓0 )

§ 2Step

• ​𝐿↓𝑀 = ​(𝑉↓0 ⨀…⨀​𝑉↓𝑜−1 )
• ​𝐿↓𝑆 =​(𝑉↓𝑜+1 ⨀…⨀​𝑉↓𝑂 )

3. Utilize BLAS = ​𝐽↓𝑜 𝐷 ​𝐽↓≠𝑜 𝐷 ​𝐽↓𝑜 ​𝑌↓(𝑜) K

Computing the KRP

Consider 𝐿=𝐵⨀𝐶⨀𝐷

• 𝐿(𝑘,:)=𝐵(𝑏,:)∗𝐶(𝑐,:)∗𝐷(𝑑,:)

𝐵 𝐶 𝐷 ⨀ ⨀ = ​𝐽↓𝐵 ​𝐽↓𝐶 ​𝐽↓𝑑 𝐵(0,:)∗𝐶(0,:)⨀𝐷 𝐿

Timings for KRPs of naïve and reuse algorithms.

1-Step MTTKRP

void permuting tensor entries ast computation as matmul y observation: the nth mode tricization of a tensor can be tained by chunking the tensor contiguous submatrices of ual size.

!"#$% !#$% !"& !& !'

) blocks

!' 1 !'

(

2(4) 2(6) 2(7$8) 9

!'

) blocks

Parallel 1-Step MTTKRP

Form ​𝐿↓𝑀 Form ​𝐿↓𝑆 (𝑘,:) Form 𝐿(𝑘,:) MatMul Reduce

2-Step MTTKRP

• First Compute a Partial MTTKRP
• 1. Compute ​𝐿↓𝑀 and ​𝐿↓𝑆
• 2. ℒ​← 𝑌↓(0:𝑜−1)↑𝑈 ∙​𝐿↓𝑀
• ℒ is ​𝐽↓𝑜 ×…×​𝐽↓𝑂−1 ×𝐷
• Second Compute a Series of ___?___ operations.
• a. Tensor Times Vector (TTVs)
• b. Tensor Times Matrix (TTMs)
• c. Quasi-Tensor Times Matrix (q-TTMs)

2-Step MTTKRP: ℒ

• First Compute a Partial MTTKRP

!"

#

!" $ !"

%

& !"

#

'(

!" $ !"

%

((*:,-.-/)

&

1(*:.-/)

2

$ =

2-Step MTTKRP: ℒ

• Second Compute a series of TTVs

! blocks )*

+(-)[0] 2(: , 0)

! )* !

56(: , 0)

)*

7

)*

7

8 =

Parallel 2-Step MTTKRP

Call Parallel BLAS WOW!!!

60×60×60×60×60

Per iteration time of a CP decomposition via ALS. Matlab used the Tensor Toolbox cp_als function, version 2.6. [1]

Findings

wo interesting networks

• Positive affect
• Negative affect

Tobia M., Hayashi K., Ballard G., Gotlib I. Dynamic Functional Connectivity and Individual Differences in Emotions During Social Stress - to appear in uman Brain Mapping

References

Tamara G. Kolda and Bre8 W. Bader. 2009. Tensor DecomposiAons and ApplicaAons. SIAM Rev. 51, 3 (Septembe 2009), 455–500. h8ps://doi.org/10.1137/ 07070111X Jiajia Li, Jee Choi, Ioakeim Perros, Jimeng Sun, and Richard Vuduc. 2017. Model Driven Sparse CP DecomposiAon for Higher-Order Tensors. In IEEE InternaAonal Parallel and Distributed Processing Symposium (IPDPS). 1048–10 h8ps://doi.org/10.1109/IPDPS.2017.80 Shaden Smith, Niranjay Ravindran, Nicholas D. Sidiropoulos, and George Karypis. 2015. SPLATT: Efficient and Parallel Sparse Tensor-Matrix MulAplicaAon. In Proceedings of the 2015 IEEE InternaAonal Parallel and Distribute Processing Symposium (IPDPS ’15). IEEE Computer Society, Washington, DC, USA, 61–70. h8ps://doi.org/10.1109/ IPDPS.2015.27 D.C. Van Essen, K. Ugurbil, E. Auerbach, D. Barch, T.E.J. Behrens, R. Bucholz, A. Chang, L. Chen, M. Corbe8a, S.W. CurAss, S. Della Penna, D. Feinberg, M.F. Glasser, N. Harel, A.C. Heath, L. Larson-Prior, D. Marcus, G. Michalareas

• S. Moeller, R. Oostenveld, S.E. Petersen, F. Prior, B.L. Schlaggar, S.M. Smith, A.Z. Snyder, J. Xu, and E. Yacoub. 20

The Human Connectome Project: a data acquisiAon perspecAve. Neuroimage 62, 4 (2012), 2222–2231. h8ps:// doi.org/10. 1016/j.neuroimage.2012.02.018 Anh-Huy Phan, Petr Tichavsky, and Andrzej Cichocki. 2013. Fast AlternaAng LS Algorithms for High Order CANDECOMP/PARAFAC Tensor FactorizaAons. IEEE TransacAons on Signal Processing 61, 19 (Oct 2013), 4834–

• 4846. h8ps://doi. org/10.1109/TSP.2013.2269903

End

Thanks for listening