30

在 Repa 中，我想d在我的数组的最里面的维度上并行应用某个维线性变换，即在所有“列”向量上。

一般来说，这样的变换可以表示为一个矩阵M，而 的每个条目M*v只是 与 的适当行的点M积v。所以我可以只使用traverse一个计算适当点积的函数。这成本d^2工作。

然而，myM很特别：它采用线性工作顺序算法。例如，M可能是一个下三角矩阵，1s 贯穿整个下三角。然后M*v只是v（又名“扫描”）的部分和的向量。这些总和可以以一种明显的方式顺序计算，但是需要(i-1)结果的第 st 项才能i有效地计算第 th 项。（我有几个这样M的，所有这些都可以在线性顺序时间内以一种或另一种方式计算。）

我没有看到任何明显的方式来使用traverse（或任何其他 Repa 函数）来利用M. 可以做到吗？d^2当有如此快速的线性工作算法可用时，使用 -work 算法（即使具有丰富的并行性）将是非常浪费的。

（我看过一些旧的 SO 帖子（例如，here）提出了类似的问题，但没有什么与我的情况完全相符。）

更新

根据要求，这里有一些说明性代码，用于M计算部分和（如上所述）。d正如我所料，运行时（工作）在数组范围的第二个参数（）中超线性增长ext。这是因为mulM'仅指定如何计算i输出的第 th 项，与所有其他项无关。即使在数组的总大小中有一个线性时间算法，但我不知道如何在 Repa 中表达它。

有趣的是，如果我从 中删除定义清单的行array'，main那么运行时只会在数组的总大小中线性缩放！因此，当数组“一直向下”延迟时，融合/优化必须以某种方式提取线性工作算法，但没有我的任何明确帮助。这太棒了，但对我来说也不是很有用，因为实际上，我需要调用mulM清单数组。

{-# LANGUAGE TypeOperators, ScopedTypeVariables, FlexibleContexts #-}

module Main where

import Data.Array.Repa as R

-- multiplication by M across innermost dimension
mulM arr = traverse arr id mulM'
    where mulM' _ idx@(i' :. i) =
              sumAllS $ extract (Z:.0) (Z:.(i+1)) $ slice arr (i' :. All)

ext = Z :. (1000000::Int) :. (10::Int) -- super-linear runtime in 2nd arg
--ext = Z :. (10::Int) :. (1000000::Int) -- takes forever

array = fromFunction ext (\(Z:.j:.i) -> j+i)

main :: IO ()
main = do
  -- apply mulM to a manifest array
  array' :: Array U DIM2 Int <- computeP $ array
  ans :: Array U DIM2 Int <- computeP $ mulM array'
  print "done"
4

0 回答 0