问题标签 [repa]

在 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清单数组。

在 repa-3.2.3.3 中是否有一种紧凑的方法可以通过常数来缩放矩阵？现在我正在使用这个：

我很困惑，因为这样的函数通常包含在数组库中。

我编写了一些代码，使用 repa计算距离矩阵：

它似乎工作。但后来我添加了一个快速检查：

换句话说，被距离 D 分开的两个点应该产生一个看起来像 [0,D,D,0] 的距离矩阵。在我的临时手动测试中，确实如此。但 QuickCheck 很快发现 5.0e-324 的距离会产生 [0,0,0,0] 的距离矩阵

这仅仅是因为双打的精度吗？我是否需要限制 QuickCheck 将发送的可能值？或者这是一个真正的错误？

背景

我repa更多地用作“管理”工具。我在:中传递reactive-bananas 。AddHandlersArrayArray D DIM2 (AddHandler Bool)

目前我正在使用这个kludge：

所以我可以做这样的事情：

问题

有没有理由不包括在内repa？

如何使用haskell对图像进行离散傅里叶变换。我相信这两个库repa-devil和repa-fftw可能会有所帮助，但我不知道如何整合它们。我之所以喜欢这个，是因为我可以尝试混合图像（我有一个绝妙的想法，即我可以简单地根据频率进行加权平均，而不是单独过滤图像然后对其进行平均。）

注意：此外，一旦完成此操作，人们将如何将其重新解释为图像。

注意：我实际上不确定如何准确地应用 dft 图像处理。我已经看过方程式并且知道如何评估它们，并且知道一维的“含义”和二维的粗略想法。我也知道如何使用它以某种方式过滤掉图像的高或低“频率”。我只是不确定使用 dft 生成的复杂频率图像到底有多复杂。

编译ghc -O2 -threaded test.hs

但是，当我执行可执行文件时./test +RTS -N1

结果是1.6777216e7

同时./test +RTS -N2或具有更多核心，

结果是正确的：2.0e7

我的代码有什么问题？

通常，在位图中绘制矩形的一种好方法是循环遍历 2 个边界尺寸并设置单个像素。例如，在伪代码上：

Haskell 的 REPA 等价物是什么？

我编写了算法的两种实现来计算对称矩阵的所有特征值和特征向量。一种实现使用 REPA 库

https://github.com/felipeZ/Haskell-abinitio/blob/master/Science/QuantumChemistry/NumericalTools/JacobiMethod.hs

而其他实现使用可变向量和 ST monad

https://github.com/felipeZ/Haskell-abinitio/blob/master/Science/QuantumChemistry/NumericalTools/JacobiMethodST.hs

该算法是 Jacobi 方法，其描述可在 http://en.wikipedia.org/wiki/Jacobi_eigenvalue_algorithm找到

我已经使用 100 个尺寸为 100 x 100 的矩阵测试了这两种实现，按顺序运行代码，我发现了以下时间：

Jacobi 算法需要对一些矩阵条目进行迭代更新，这意味着大多数矩阵在迭代之间保持不变。因此，我会（错误地）猜到，在 REPA 实现中为每次迭代复制一个新矩阵的成本将大于使用 monad ST 改变矩阵的成本，因为据我所知，REPA 不会改变数组，但复制它。

是 REPA 融合了所有操作并避免在每次迭代中复制新数组吗？还是别的什么？

有人可以评论这个结果吗？

我在理解 repa 的Slice模块时遇到了一些麻烦，特别是它的类型和系列如何与该slice功能一起使用。我相信以下问题的答案会为我解决问题。

假设x有形状shapeOfList [m_0, m_1, …, m_{d-1}]。给定整数i < d，j < mi我如何构建具有合适类型作为slice第二个参数的东西，从冻结索引号定义的形状x中提取切片？shapeOfList [m_0, …, m_{i-1}, m_{i+1}, …, m_{d-1}]ij

案例中的例子d=2是：

• i=0: 提取行j
• i=1: 提取列j

我有一个 3D 坐标 (x,y,z) 列表，它在每个位置指定一个标量值 (rho)，我已将其加载到 repa 数组中：

即每个数组元素是一个4元组，代表(x,y,z,rho)。此外，我不依赖于使用 repa。

我想将其转换为 DIM3 数组，每个位置都包含适当的 rho 值（或默认值，或者可能是插值）。

我想不出一种方法来做到这一点traverse（至少效率不高），因为traverse使用来自新数组的索引，而不是源数组。repa 是否有可能，还是我需要使用不同的数组库（甚至是可变的）？