Rmpfr包中outer等矩阵乘法函数运行缓慢的原因咨询

我最近在用R的Rmpfr包处理高精度数值时，发现一个很明显的性能差异：outer()、tcrossprod()和%*%这些常规矩阵乘法函数，居然比先把向量重复再做逐元素乘法慢超多！

先给大家看一下我做的测试代码和结果：

library(Rmpfr)

m <- 100
n <- 50
x <- mpfr(runif(m), 128)
y <- mpfr(runif(n), 128)

bench::mark(
  outer = outer(x, y),
  tcrossprod = tcrossprod(x, y),
  "%*%" = x%*%t(y),
  repmult = matrix(rep(x, n)*rep(y, each = m), m, n),
  relative = TRUE
)

运行后的性能对比结果：

A tibble: 4 × 6 #

expression min median itr/sec mem_alloc gc/sec
bch:expr
1 outer 89.5 67.7 1.51 525. 1
2 tcrossprod 83.7 63.3 1.61 513. 1.18
3 %*% 135. 102. 1 517. 1.36
4 repmult 1 1 101. 1 2.00

从结果能看到，repmult的速度比最慢的%*%快了100倍左右，比outer也快了近60倍！而且我发现这个性能差距还会随着最终矩阵的规模增大而大致线性增长，我用的是默认的R BLAS。

我想请教一下，为啥会出现这么大的性能差异呢？我猜测可能是这些矩阵操作函数对于Rmpfr的高精度数值类型没有做针对性优化？比如outer()这类函数内部可能是用循环逐个处理元素，而重复向量后做逐元素乘是利用了R的向量化优势？或者是默认BLAS对mpfr类型的矩阵乘法没有优化支持，导致这些函数只能走低效的实现路径？

备注：内容来源于stack exchange，提问作者jblood94