A329369 на PARI/GP: три цикла против матриц

[kthxbye]

Пусть задана относительно простая последовательность (A329369)

$$a(2^m(2k+1))=\sum\limits_{j=0}^{m}\binom{m+1}{j}a(2^jk), a(0)=1$$

Чтобы вычислять ее без рекурсии с использованием небольшого числа относительно маленьких векторов, я (интуитивно) придумал следующий алгоритм (функция a1):

a(n) = {if(n == 0, 1,
my(A = valuation(n, 2), B = n\2^(A+1));
sum(j=0, A, binomial(A+1, j)*a(2^j*B)))}
a1(n) = {my(n = 2*n, A = 1, v1, v2, v3); v1 = [];
for(i=0, logint(n, 2), 
if(bittest(n, i) == bittest(n, i+1), A++,
v1 = concat(v1, A); A = 1));
v1 = Vecrev(v1);
for(i=2, #v1/2, v1[2*i] += v1[2*(i-1)]);
v2 = vector(v1[#v1], i, 1); v3 = v2;
v4 = vector(#v2+1, i,
vector(i, j, j == 1 || j == i));
for(i=3, #v4, for(j=2, i-1, 
v4[i][j] = v4[i-1][j] + v4[i-1][j-1]));
for(i=1, #v1/2, 
A = if(i == 1, 0, v1[2*(i-1)]);
for(j=1, v1[2*i-1],
for(k=A+1, #v2,
v3[k] = sum(q=1, k-A, v4[k - A + 1][q]*v2[q+A]));
v2 = v3));
v2[#v2]}
test(n) = a1(n) == a(n)

Потом я задал вопрос на сайте MathOverflow (вот он). Там некий Peter Taylor подтвердил (пусть и без деталей), что мой алгоритм эквивалентен другому алгоритму с использованием матриц. Этот другой алгоритм он вроде бы подробно описал в своем ответе, но к сожалению я ничего не понял.

В ходе дальнейшей дискуссии Питер (насколько я понял) объявил мой алгоритм крайне сложным в плане вычислений, а вот матричный метод по его мнению довольно хорош. Прав ли он в этом вопросе? Будет ли алгоритм с использованием матриц работать быстрее? Если да, то как написать его на PARI/GP для сравнения?