申请个人网站,做seo必须有自己网站吗,做服装网站要那些照片,东莞搜索引擎网站推广优于立方复杂度的 Rust 中矩阵乘法 迈克克维特 跟随 发表于 更好的编程 6 分钟阅读 7月 143 中途#xff1a;三次矩阵乘法 一、说明 几年前#xff0c;我在 C 年编写了 Strassen 矩阵乘法算法的实现#xff0c;最近在 Rust 中重新实现了它#xff0c;因为我继续… 优于立方复杂度的 Rust 中矩阵乘法 迈克·克维特 · 跟随 发表于 更好的编程 · 6 分钟阅读 · 7月 143 中途三次矩阵乘法 一、说明 几年前我在 C 年编写了 Strassen 矩阵乘法算法的实现最近在 Rust 中重新实现了它因为我继续学习该语言。这是学习 Rust 性能特征和优化技术的有用练习因为尽管 Strassen 的算法复杂性优于朴素方法但它在算法结构中的分配和递归开销中具有很高的常数因子。 通用算法换位以获得更好的性能次立方斯特拉森算法的工作原理排比标杆分析和性能优化 二、通用算法 一般朴素矩阵乘法算法是每个人在他们的第一堂线性代数课上学习的三个嵌套循环方法大多数人会将其识别为 On³ pub fn
mult_naive (a: Matrix, b: Matrix) - Matrix {if a.rows b.cols {let m a.rows;let n a.cols;// preallocatelet mut c: Vecf64 Vec::with_capacity(m * m);for i in 0..m {for j in 0..m {let mut sum: f64 0.0;for k in 0..n {sum a.at(i, k) * b.at(k, j);}c.push(sum);}}return Matrix::with_vector(c, m, m);} else {panic!(Matrix sizes do not match);}
} 这种算法很慢不仅因为三个嵌套循环还因为按列通过而不是按行的内部循环遍历对于 CPU 缓存命中率来说是可怕的。Bb.at(k, j) 三、换位以获得更好的性能 转置朴素方法允许 B 上的乘法迭代在行而不是列上运行将矩阵 B 的乘法步幅重新组织为更有利于缓存的格式。从而变成A x BA x B^t 它涉及一个新的矩阵分配无论如何在这个实现中和一个完整的矩阵迭代一个 On² 操作更准确地说这种方法是 On³ On²——我将进一步展示它的性能有多好。它如下所示 fn multiply_transpose (A: Matrix, B: Matrix):C new Matrix(A.num_rows, B.num_cols)// Construct transpose; requires allocation and iteration through BB’ B.transpose()for i in 0 to A.num_rows:for j in 0 to B.num_rows:sum 0;for k in 0 to A.num_cols:// Sequential access of B[j, k] is much faster than B[k, j]sum A[i, k] * B[j, k]C[i, j] sumreturn C 四、次立方斯特拉森算法的工作原理 要了解 Strassen 算法的工作原理此处为 Rust 代码首先考虑矩阵如何用象限表示。要概念化它的外观 在朴素算法中使用此象限模型结果矩阵 C 的四个象限中的每一个都是两个子矩阵乘积的总和总共产生 8 次乘法。 考虑到这八个乘法每个乘法都在一个块矩阵上运行其行和列跨度约为 A 和 B 大小的一半复杂性相同 斯特拉森算法定义了由这些象限组成的七个中间块矩阵 仅通过 7 次乘法而不是 8 次乘法计算。这些乘法可以是递归斯特拉森乘法并可用于组成最终矩阵 由此产生的亚立方复杂度 五、排比 中间矩阵 M1 的计算 ...M7 是一个令人尴尬的并行问题因此也很容易检测算法的并发变体一旦你开始理解 Rust 关于闭包的规则。 /*** Execute a recursive strassen multiplication of the given vectors, * from a thread contained within the provided thread pool.*/
fn
_par_run_strassen (a: Vecf64, b: Vecf64, m: usize, pool: ThreadPool) - ArcMutexOptionMatrix {let m1: ArcMutexOptionMatrix Arc::new(Mutex::new(None));let m1_clone Arc::clone(m1);pool.execute(move|| { // Recurse with non-parallel algorithm once were // in a working threadlet result mult_strassen(mut Matrix::with_vector(a, m, m),mut Matrix::with_vector(b, m, m));*m1_clone.lock().unwrap() Some(result);});return m1;
} 六、标杆 我编写了一些快速的基准测试代码该代码在不断增加的矩阵维度范围内运行四种算法中的每一种进行几次试验并报告每种算法的平均时间。 ~/code/strassen ~ ./strassen --lower 75 --upper 100 --factor 50 --trials 2running 50 groups of 2 trials with bounds between [75-3750, 100-5000]x y nxn naive transpose strassen par_strassen
75 100 7500 0.00ms 0.00ms 1.00ms 0.00ms
150 200 30000 6.50ms 4.00ms 4.00ms 1.00ms
225 300 67500 12.50ms 9.00ms 8.50ms 2.50ms
300 400 120000 26.50ms 22.00ms 18.00ms 5.50ms
[...]
3600 4800 17280000 131445.00ms 53683.50ms 21210.50ms 5660.00ms
3675 4900 18007500 141419.00ms 58530.00ms 28291.50ms 6811.00ms
3750 5000 18750000 154941.00ms 60990.00ms 26132.00ms 6613.00ms 然后我通过以下方式可视化结果pyplot 此图显示了矩阵从 7.5k 元素 到大约 19 万 的乘法时间。你可以看到朴素算法在计算上变得不切实际的速度有多快在高端需要两分半钟。N x M 75 x 100N x M 3750 x 5000 相比之下Strassen 算法的扩展更平滑并行算法计算两个 19M 个元素的矩阵的结果而朴素算法只处理 3.6M 个元素所花费的时间。 对我来说最有趣的是算法的性能。如前所述缓存性能的改进以牺牲完整矩阵副本为代价在这些结果中得到了清楚地证明 - 即使使用与该方法渐近等效的算法也是如此。transposenaive 七、分析和性能优化 这个文档是理解 Rust 性能基础知识的绝佳资源。在 Mac OS 上启动并运行仪器进行分析是微不足道的这要归功于货运仪器的 Rust 指南。这是调查分配行为、CPU 热点和其他事情的绝佳工具。 在此过程中发生了一些变化 Strassen 代码通过分而治之策略递归调用自己但是一旦矩阵达到足够小的大小其高常数因子使其比一般矩阵算法慢。我发现这个点是大约 64 的行宽或列宽;通过提高吞吐量提高几个因素来增加此阈值2斯特拉森算法要求矩阵填充到最接近的指数 2;减少这种情况以懒惰地确保矩阵只有偶数行和列 通过减少昂贵的大分配将吞吐量提高了大约两倍将小矩阵回退算法从 更改为 导致大约 20% 的改进naivetranspose添加和添加到 Cargo.toml 发布构建标志大约提高了 5%。有趣的是性能持续恶化codegen-units 1lto thinlto “true”一丝不苟地删除所有可能的副本大约提高了~10%Vec提供一些提示并删除随机访问查找中的向量边界检查又提高了大约 20%#[inline] /*** Returns the element at (i, j). Unsafe.*/#[inline]pub fn at (self, i: usize, j: usize) - f64 {unsafe {return *self.elements.get_unchecked(i * self.cols j);}} 参考资料 迈克·克维特 线性代数 算法 锈
