Try to improve massiv benchmark performance

bench/Main.hs

 {-# LANGUAGE TypeApplications  #-}
 {-# LANGUAGE OverloadedStrings #-}
+{-# LANGUAGE NumericUnderscores #-}
+{-# LANGUAGE BangPatterns #-}
+{-# LANGUAGE ScopedTypeVariables #-}
 
 module Main where
 
@@ -10,14 +13,17 @@ import Gargantext.Core.Viz.Phylo.API.Tools (readPhylo)
 import Gargantext.Core.Viz.Phylo.PhyloMaker (toPhylo)
 import Gargantext.Core.Viz.Phylo.PhyloTools
 import Gargantext.Prelude.Crypto.Auth (createPasswordHash)
-import qualified Gargantext.Core.LinearAlgebra as LA
-import qualified Gargantext.Core.Methods.Similarities.Accelerate.Distributional as Accelerate
-
-import Test.Tasty.Bench
 import Paths_gargantext
 import qualified Data.Array.Accelerate as A
-import qualified Data.Massiv.Array as Massiv
 import qualified Data.Array.Accelerate.Interpreter as Naive
+import qualified Data.List.Split as Split
+import qualified Data.Massiv.Array as Massiv
+import qualified Gargantext.Core.LinearAlgebra as LA
+import qualified Gargantext.Core.Methods.Similarities.Accelerate.Distributional as Accelerate
+import qualified Gargantext.Core.Methods.Matrix.Accelerate.Utils as Accelerate
+import Test.Tasty.Bench
+import qualified Data.Array.Accelerate.Interpreter as LLVM
+import qualified Data.Array.Accelerate as Accelerate
 
 phyloConfig :: PhyloConfig
 phyloConfig = PhyloConfig {
@@ -43,32 +49,29 @@ phyloConfig = PhyloConfig {
   , exportFilter = [ByBranchSize {_branch_size = 3.0}]
   }
 
+matrixValues :: [Int]
+matrixValues = [ 1 .. 10_000 ]
+
+matrixDim :: Int
+matrixDim = 100
+
 testMatrix :: A.Matrix Int
-testMatrix = A.fromList (A.Z A.:. 14 A.:. 14) $
-  [  30,  36, -36, -16,   0,   7,  34,  -7,   5,  -4,   0,  21,   6, -35,
-      0, -31,  20, -15, -22,  -7, -22, -37, -29, -29,  23, -31, -29, -23,
-    -24, -29,  19,  -6,  16,   7,  15, -27, -27, -30,  -9, -33,  18, -23,
-      7, -36,  12,  26, -17,  -3,  -2, -15,  -4,  26,  24,   9,  -4,   4,
-     32,  28,  -2, -10,  34,  -3,  20,  -9, -22,  20, -26,  34,  18, -21,
-      7, -12,  12,  -2,  36,  10,  34, -37,  13,  -9, -28,  34,  33, -18,
-     -4, -32,  -1,  29,  29, -28,  24,  28,  35,  19,   8, -18,  25, -35,
-    -14,  -4, -24,  -1,   7,  34, -37, -28, -12, -32,  -5, -23,  27,  33,
-    -36, -28,  21, -29,  -2, -26,  -4, -31, -26, -21,  33, -11, -33,  20,
-     25,  14,   5,  -7,   5,  24,  37,   1,  -3,  23,  25, -16,  17,   5,
-    -35,  36,  -2,  -2,   1, -14,  34, -30, -10,  12,  25,  21,   0,  34,
-     17,  -1,  20, -19,  15,  20,  -5, -30, -35, -13,   5,  17, -10, -19,
-    -34, -11, -18,  26, -29, -28,   0,   3,  23,  -6,  36,   4,  16,  28,
-     13, -37, -16,   2,   7, -13,  21, -10, -33, -33, -26, -19,  -1,  29]
+testMatrix = A.fromList (A.Z A.:. matrixDim A.:. matrixDim) $ matrixValues
 {-# INLINE testMatrix #-}
 
-testMassivMatrix :: Massiv.Matrix Massiv.D Int
-testMassivMatrix = LA.accelerate2MassivMatrix testMatrix
+testMassivMatrix :: Massiv.Matrix Massiv.U Int
+testMassivMatrix = Massiv.fromLists' Massiv.Par $ Split.chunksOf matrixDim $ matrixValues
 {-# INLINE testMassivMatrix #-}
 
 main :: IO ()
 main = do
   _issue290Phylo     <- force . setConfig phyloConfig <$> (readPhylo =<< getDataFileName "bench-data/phylo/issue-290.json")
   issue290PhyloSmall <- force . setConfig phyloConfig <$> (readPhylo =<< getDataFileName "bench-data/phylo/issue-290-small.json")
+  let !accInput       = force testMatrix
+  let !massivInput    = force testMassivMatrix
+  let !(accDoubleInput :: Accelerate.Matrix Double) = force $ Naive.run $ Accelerate.map Accelerate.fromIntegral (Accelerate.use testMatrix)
+  let !massivInput    = force testMassivMatrix
+  let !(massivDoubleInput :: Massiv.Matrix Massiv.U Double) = force $ Massiv.computeP $ Massiv.map fromIntegral testMassivMatrix
   defaultMain
     [ bgroup "Benchmarks"
       [ bgroup "User creation" [
@@ -79,10 +82,22 @@ main = do
       , bgroup "Phylo" [
           bench "toPhylo (small)" $ nf toPhylo issue290PhyloSmall
       ]
+      , bgroup "diag" [
+          bench "Accelerate (Naive)" $ nf (Naive.run . Accelerate.diag . Accelerate.use) accInput
+      ,   bench "Accelerate (LLVM)"  $ nf (LLVM.run  . Accelerate.diag . Accelerate.use) accInput
+      ,   bench "Massiv " $ nf (LA.diag @_) massivInput
+      ]
+      , bgroup "termDivNan" [
+          bench "Accelerate (Naive)" $
+            nf (\m -> Naive.run $ Accelerate.termDivNan (Accelerate.use m) (Accelerate.use m)) accDoubleInput
+      ,   bench "Accelerate (LLVM)" $
+            nf (\m -> LLVM.run $ Accelerate.termDivNan (Accelerate.use m) (Accelerate.use m)) accDoubleInput
+      ,   bench "Massiv " $ nf (\m -> LA.termDivNan @Massiv.U m m) massivDoubleInput
+      ]
       , bgroup "distributional" [
-          bench "Accelerate (Naive)" $ nf (Accelerate.distributionalWith @Double Naive.run) testMatrix
-      ,   bench "Accelerate (LLVM)"  $ nf Accelerate.distributional testMatrix
-      ,   bench "Massiv " $ nf (Massiv.computeAs Massiv.U . LA.distributional @Double) testMassivMatrix
+          bench "Accelerate (Naive)" $ nf (Accelerate.distributionalWith @Double Naive.run) accInput
+      ,   bench "Accelerate (LLVM)"  $ nf Accelerate.distributional accInput
+      ,   bench "Massiv " $ nf (Massiv.computeP @Massiv.U . LA.distributional @_ @Double) massivInput
       ]
       ]
     ]

gargantext.cabal

@@ -309,6 +309,7 @@ library
       Gargantext.Utils.SpacyNLP.Types
       Gargantext.Utils.Tuple
       Gargantext.Utils.Zip
+      Paths_gargantext
   other-modules:
       Gargantext.API.Admin.Auth
       Gargantext.API.Admin.FrontEnd
@@ -477,7 +478,6 @@ library
       Gargantext.Utils.Aeson
       Gargantext.Utils.Servant
       Gargantext.Utils.UTCTime
-      Paths_gargantext
   ghc-options: -Wall -Wincomplete-uni-patterns -Wincomplete-record-updates -Wmissing-signatures -Wunused-binds -Wunused-imports -Wunused-packages -Werror -freduction-depth=300 -fprint-potential-instances
   hs-source-dirs:
       src
@@ -908,8 +908,9 @@ benchmark garg-bench
                , deepseq
                , gargantext
                , gargantext-prelude
+               , split
                , tasty-bench
-  ghc-options:   "-with-rtsopts=-T -A32m"
+  ghc-options:   -threaded "-with-rtsopts=-N -T -A32m"
   if impl(ghc >= 8.6)
-    ghc-options: "-with-rtsopts=--nonmoving-gc"
+    ghc-options: -threaded "-with-rtsopts=-N --nonmoving-gc"
 

src/Gargantext/Core/LinearAlgebra.hs

@@ -58,10 +58,10 @@ createIndices = set2indices . map2set
     set2indices s = foldr (uncurry Bimap.insert) Bimap.empty (zip [0..] $ S.toList s)
 
 -- | Converts an accelerate matrix into a Massiv matrix.
-accelerate2MassivMatrix :: (A.Unbox a, Acc.Elt a) => Acc.Matrix a -> Matrix D a
+accelerate2MassivMatrix :: (A.Unbox a, Acc.Elt a) => Acc.Matrix a -> Matrix A.U a
 accelerate2MassivMatrix m =
   let (Acc.Z Acc.:. _r Acc.:. c) = Acc.arrayShape m
-  in A.delay $ A.fromLists' @A.U A.Par $ Split.chunksOf c (Acc.toList m)
+  in A.fromLists' @A.U A.Par $ Split.chunksOf c (Acc.toList m)
 
 -- | Converts a massiv matrix into an accelerate matrix.
 massiv2AccelerateMatrix :: (Acc.Elt a, A.Source r a) => Matrix r a -> Acc.Matrix a
@@ -80,13 +80,11 @@ massiv2AccelerateVector m =
 
 
 -- | Computes the diagnonal matrix of the input one.
-diag :: Num e => Matrix D e -> Vector D e
+diag :: (A.Unbox e, A.Manifest r e, A.Source r e, Num e) => Matrix r e -> Vector A.U e
 diag matrix =
   let (A.Sz2 rows _cols) = A.size matrix
       newSize = A.Sz1 rows
-  in A.backpermute' newSize (\j -> j A.:. j) matrix
-
-
+  in A.makeArrayR A.U A.Seq newSize $ (\(A.Ix1 i) -> matrix A.! (A.Ix2 i i))
 
 -- | `distributional m` returns the distributional distance between terms each
 -- pair of terms as a matrix.  The argument m is the matrix $[n_{ij}]_{i,j}$
@@ -124,30 +122,31 @@ diag matrix =
 --
 -- /IMPORTANT/: As this function computes the diagonal matrix in order to carry on the computation
 -- the input has to be a square matrix, or this function will fail at runtime.
-distributional :: forall e. (Ord e, Fractional e, Num e)
-               => Matrix D Int
-               -> Matrix D e
+distributional :: forall r e. (A.Manifest r e, A.Unbox e, A.Source r Int, A.Size r, Ord e, Fractional e, Num e)
+               => Matrix r Int
+               -> Matrix r e
 distributional m' = result
  where
-    m :: Matrix D e
-    m = A.map fromIntegral m'
+    m :: Matrix A.U e
+    m = A.compute$ A.map fromIntegral m'
     n = dim m'
 
-    diag_m :: Vector D e
-    diag_m = diag m
+    diag_m :: Vector A.U e
+    diag_m = A.computeP $ diag m
+
+    diag_m_size :: Int
+    (A.Sz1 diag_m_size) = A.size diag_m
 
     -- replicate (constant (Z :. n :. All)) diag_m
     d_1 :: Matrix D e
-    d_1 = let (A.Sz1 r) = A.size diag_m
-          in A.backpermute' (A.Sz2 n r) (\(_ A.:. i) -> i) diag_m
+    d_1 = A.backpermute' (A.Sz2 n diag_m_size) (\(_ A.:. i) -> i) diag_m
 
     -- replicate (constant (Z :. All :. n)) diag_m
     d_2 :: Matrix D e
-    d_2 = let (A.Sz1 r) = A.size diag_m
-          in A.backpermute' (A.Sz2 r n) (\(i A.:. _) -> i) diag_m
+    d_2 = A.backpermute' (A.Sz2 diag_m_size n) (\(i A.:. _) -> i) diag_m
 
-    mi :: Matrix D e
-    mi = (.*) (termDivNan m d_1) (termDivNan m d_2)
+    mi :: Matrix A.U e
+    mi = (.*) (termDivNan @A.U m d_1) (termDivNan @A.U m d_2)
 
     -- The matrix permutations is taken care of below by directly replicating
     -- the matrix mi, making the matrix w unneccessary and saving one step.
@@ -167,33 +166,42 @@ distributional m' = result
     -- The matrix ii = [r_{i,j,k}]_{i,j,k} has r_(i,j,k) = 0 if k = i OR k = j 
     -- and r_(i,j,k) = 1 otherwise (i.e. k /= i AND k /= j).
     -- generate (constant (Z :. n :. n :. n)) (lift1 (\( i A.:. j A.:. k) -> cond ((&&) ((/=) k i) ((/=) k j)) 1 0))
-    ii :: Array D Ix3 e
-    ii = A.makeArrayR A.D A.Par (A.Sz3 n n n) $ \(i A.:> j A.:. k) -> if k /= i && k /= j then 1 else 0
+    ii :: Array A.U Ix3 e
+    ii = A.makeArrayR A.U A.Par (A.Sz3 n n n) $ \(i A.:> j A.:. k) -> if k /= i && k /= j then 1 else 0
 
-    z_1 :: Matrix D e
+    z_1 :: Matrix A.U e
     z_1 = sumRows ((.*) w' ii)
 
-    z_2 :: Matrix D e
+    z_2 :: Matrix A.U e
     z_2 = sumRows ((.*) w_1 ii)
 
     result = termDivNan z_1 z_2
 
 -- | Term by term division where divisions by 0 produce 0 rather than NaN.
-termDivNan :: (Eq a, Fractional a) => Matrix D a -> Matrix D a -> Matrix D a
-termDivNan = A.zipWith (\i j -> if j == 0 then 0 else i / j)
-
-sumRows :: Num e => Array D A.Ix3 e -> Array D A.Ix2 e
+termDivNan :: (A.Manifest r3 a, A.Source r1 a, A.Source r2 a, Eq a, Fractional a)
+           => Matrix r1 a
+           -> Matrix r2 a
+           -> Matrix r3 a
+termDivNan m1 m2 = A.compute $ A.zipWith (\i j -> if j == 0 then 0 else i / j) m1 m2
+
+sumRows :: (A.Load r A.Ix2 e
+           , A.Source r e
+           , A.Strategy r
+           , A.Size r
+           , Num e
+           ) => Array r A.Ix3 e
+             -> Array r A.Ix2 e
 sumRows matrix =
   let A.Sz3 rows cols z = A.size matrix
   in A.makeArray (A.getComp matrix) (A.Sz2 rows cols) $ \(i A.:. j) ->
        A.sum (A.backpermute' (A.Sz1 z) (\c -> i A.:> j A.:. c) matrix)
 
 -- | Matrix cell by cell multiplication
-(.*) :: (A.Index ix, Num a)
-     => Array D ix a
-     -> Array D ix a
-     -> Array D ix a
-(.*) = A.zipWith (*)
+(.*) :: (A.Manifest r3 a, A.Source r1 a, A.Source r2 a, A.Index ix, Num a)
+     => Array r1 ix a
+     -> Array r2 ix a
+     -> Array r3 ix a
+(.*) m1 m2 = A.compute $ A.zipWith (*) m1 m2
 
 -- | Get the dimensions of a /square/ matrix.
 dim :: A.Size r => Matrix r a -> Int

src/Gargantext/Core/Methods/Matrix/Accelerate/Utils.hs

@@ -62,16 +62,14 @@ import qualified Gargantext.Prelude as P
 (./) = zipWith (/)
 
 -- | Term by term division where divisions by 0 produce 0 rather than NaN. 
-termDivNan :: ( Shape ix
-        , Slice ix
-        , Elt a
-        , Eq a
-        , P.Num (Exp a)
-        , P.Fractional (Exp a)
-        )
-     => Acc (Array ((ix :. Int) :. Int) a)
-     -> Acc (Array ((ix :. Int) :. Int) a)
-     -> Acc (Array ((ix :. Int) :. Int) a)
+termDivNan :: ( Elt a
+              , Eq a
+              , P.Num (Exp a)
+              , P.Fractional (Exp a)
+              )
+           => Acc (Matrix a)
+           -> Acc (Matrix a)
+           -> Acc (Matrix a)
 termDivNan = zipWith (\i j -> cond ((==) j 0) 0 ((/) i j))
 
 (.-) :: ( Shape ix

test/Test/Core/LinearAlgebra.hs

@@ -19,10 +19,8 @@ import Gargantext.Core.Methods.Matrix.Accelerate.Utils qualified as Legacy
 import Gargantext.Core.Methods.Similarities.Accelerate.Distributional qualified as Legacy
 import Gargantext.Core.Viz.Graph.Index qualified as Legacy
 import Gargantext.Orphans.Accelerate (sliceArray)
-import Gargantext.Prelude (headMay)
 import Prelude hiding ((^))
 import qualified Data.Array.Accelerate as A
-import qualified Data.List.Split.Internals as Split
 import Test.Tasty
 import Test.Tasty.QuickCheck
 import Data.Proxy
@@ -124,21 +122,22 @@ compareTermDivNan :: (Array TermDivNanShape Double)
                   -> (Array TermDivNanShape Double)
                   -> Property
 compareTermDivNan i1 i2
-  = let massiv     = LA.termDivNan (accelerate2MassivMatrix i1) (accelerate2MassivMatrix i2)
+  = let massiv     = LA.termDivNan @Massiv.U (LA.accelerate2MassivMatrix i1) (LA.accelerate2MassivMatrix i2)
         accelerate = Naive.run (Legacy.termDivNan (use i1) (use i2))
-    in accelerate === massiv2AccelerateMatrix massiv
+    in accelerate === LA.massiv2AccelerateMatrix massiv
 
 compareDiag :: SquareMatrix Int -> Property
 compareDiag (SquareMatrix i1)
-  = let massiv     = Massiv.computeAs Massiv.U $ LA.diag (accelerate2MassivMatrix i1)
+  = let massiv     = LA.diag (LA.accelerate2MassivMatrix i1)
         accelerate =  Naive.run (Legacy.diag (use i1))
-    in accelerate === massiv2AccelerateVector massiv
+    in accelerate === LA.massiv2AccelerateVector massiv
 
 compareDistributional :: forall e.
                       ( Eq e
                       , Show e
                       , FromIntegral Int e
                       , Prelude.RealFloat e
+                      , Massiv.Unbox e
                       , A.Ord e
                       , Ord e
                       , Prelude.Fractional (Exp e)
@@ -147,7 +146,7 @@ compareDistributional :: forall e.
                       -> SquareMatrix Int
                       -> Property
 compareDistributional Proxy (SquareMatrix i1)
-  = let massiv     = Massiv.computeAs Massiv.B $ LA.distributional @e (accelerate2MassivMatrix i1)
+  = let massiv     = Massiv.computeAs Massiv.B $ LA.distributional @_ @e (LA.accelerate2MassivMatrix i1)
         accelerate = Legacy.distributionalWith Naive.run i1
         expected   = map conv (A.toList accelerate)
         actual     = map conv (mconcat (Massiv.toLists2 massiv))