Comments 11
Шучу. А если серьезно, так оно и есть, и в силу удобства и мощи Rcpp крайне рекомендуется для всех пользователей R. В сложных проектах часто даже «векторные» эрные функции переписывал на C++, и работали они потом на порядок быстрее.
Ребята, вот же красивое решение безо всяких циклов. Строит ганкелеву матрицу размера LxK, состоящую из номеров побочных диагоналей, в три строчки безо всяких циклов, Rcpp, Vectorize и прочего:
K <- N - L + 1 i <- rep(1:L, K) j <- rep(1:K, each = L) matrix(i + j - 1, nrow = L)
Естественно, если последнюю строчку перед взятием matrix поменять на что-то типа:
matrix(min(1:N, N:1)[i + j - 1], nrow = L)
то получим решение вашей задачи.
Ой, переносы строк поехали. И pmin() вместо min(), конечно же.
Спасибо, хороший вариант. Задача, как я упоминал, отчасти учебная — на ней удобно показывать спектр различных приёмов. А вариантов её решить миллион, и всегда найдётся простор для ещё большей оптимизации. Так, например, заменить первый rep
на rep.int
— будет ещё быстрее.
Быстрый и безболезненный способ повысить производительность циклов — использовать пакет compiler
. Он присутствует в R из коробки, с версии R 3.4 обещают включать компиляцию по умолчанию. На моей машине строчка library(compiler);enableJIT(3)
перед всеми расчетами сразу ускоряет build_naive
раза в четыре — после включения компиляции она вычисляется даже чуточку быстрее скомпилированного же build_recursive
.
Спасибо за статью!
Еще немного оптимизации — можно воспользоваться тем, что матрица является персимметричной и работать с побочными диагоналями. В таком случае сокращается число итераций в цикле и не надо использовать min
IntegerMatrix build_cpp_PerSym(const int n) {
IntegerMatrix x(n, n);
for (int i=0; i<n; i++)
for (int j=0; j<=i; j++)
x(j, i-j) = x(n-1-j, n-1-(i-j)) = i+1;
return x;
}
У меня microbenchmark показывает, что этот метод немного быстрее, чем build_cpp
.
В R также можно использовать из коробки замечательные cpp библиотеки для работы с линейной алгеброй — Armadillo и Eigen. В Armadillo требуемая матрица генерируется, можно сказать, в одну строчку
umat flipToeplitz(const int n) {
return fliplr(toeplitz(linspace<uvec>(n, 1, n)));
}
Здесь столбцы матрицы Тёплица от вектора n:1 записываются в обратном порядке.
Этим методом матрицы не очень большого размера (~ 100) получается генерировать примерно за то же время, что и методом build_cpp
.
Спасибо за идеи — абсолютно согласен. У меня есть опыт работы с обеими библиотеками, и он весьма положительный, так что присоединяюсь к рекомендациям.
Для полноты картины прямые ссылки на пакеты, реализующие соответствующие интерфейсы: RcppArmadillo, RcppEigen.
Низкоуровневая оптимизация и измерение производительности кода на R