Перманент матрицы

[maxal]

wrest в сообщении #1600055 писал(а):

Скажите пож-ста, 40 секунд Sage и 8 секунд PARI/GP -- это один и тот же комп? В коде Sage как я понимаю включена мемоизация (хотя там мне кажется и хитов нету), так что рекурсия против прямого цикла вроде не проблема. Или разные компы?

Комп один и тот же. Просто Sage проигрывает в производительности PARI/GP, ничего удивительного.

[Null]

Geen в сообщении #1600042 писал(а):

На самом деле меньше

Длинна чисел $O(n\log(n))$

[vicvolf]

maxal в сообщении #1600061 писал(а):

Комп один и тот же. Просто Sage проигрывает в производительности PARI/GP, ничего удивительного.

Конечно время счета одного и того же алгоритма зависит от программы и производительности компа. Сложность алгоритма зависит от количества длинных и коротких операций или от средней длины числа. В данном случае наилучший показатель:

-- 06.07.2023, 08:09 --

Null в сообщении #1600064 писал(а):

$O(n\log(n))$

[kthxbye]

Ну что ж, наверное стоит добавить результаты в OEIS. Предлагаю следующие варианты программ на PARI/GP:

A204262

upto(n)=my(v1); v1=vector(n+1, i, i--; i!*x^i); for(i=1, n, for(j=i, n, my(A=intformal((j-i)^2*v1[j])); v1[j+1] = A + subst(v1[j+1] - A, x, i))); v1

A204264

a(n)=if(n==0, 1, my(v1); v1=vector(n+1, i, i--; i!*x^i); for(i=1, n, for(j=i, n, my(A=intformal((j-i)^2*v1[j])); v1[j+1] = A + subst(v1[j+1] - A, x, n - i + 1))); v1[n+1])

Добавил бы сам, но у меня все 3 слота правок заняты. Желательно также добавить ссылку на эту тему. Если кто-нибудь откроет правку к A204262, то я смогу добавить туда свою формулу (если успею до аппрува). Результат участника Null тоже желательно оформить формулой (при этом я не уверен как это правильно сделать).

[wrest]

kthxbye
Сделайте так, чтобы код не зависел от того, определён ранее x или нет.
Иначе может быть так:

? x=6;
? upto(3)
  ***   at top-level: upto(3)
  ***                 ^-------
  ***   in function upto: ...1[j+1]=A+subst(v1[j+1]-A,x,i)));v1
  ***                                                 ^---------
  ***   incorrect type in evaluator [variable name expected] (t_INT).
  ***   Break loop: type 'break' to go back to GP prompt
break>

И это в лучшем случае. В худшем - может неверно считать.

Например, вместо upto(n)=my(v1); сделайте upto(n)=my(v1, x='x);

Да, и по правилам хорошего тона, вместо A=intformal((j-i)^2*v1[j]) лучше писать A=intformal((j-i)^2*v1[j],x)

P.S. Но в целом, хотел вам сказать что я весьма впечатлён. Хороший код, ничего лишнего. Прогресс налицо!
P.P.S. Меня смущает intformal() но я не знаю, ускорит ли тут что-то "прямой" код.

[wrest]

Ещё процентов на 20 побыстрее:

{upto(n)=
my(v,x='x,A);
v1=vector(n+1,i,i--;i!*x^i);
for(i=1,n,
  for(j=i,n,
    if(i==j,A=0,A=intformal((j-i)^2*v[j],x));
    if(j>(2*i),v[j+1]=A,v1[j+1]=A+subst(v[j+1]-A,x,i))
  )
 );
return(v)}

[Null]

wrest в сообщении #1600404 писал(а):

if(j>(2*i),v[j+1]=A,v1[j+1]=A+subst(v[j+1]-A,x,i))

Да. Матрица с большим квадратом из 0 имеет нулевой перманент.

[wrest]

Нельзя ли снизить сложность до $O(n^2)$?

[kthxbye]

Задал вопрос на MathOverflow и получил ответ от Теренса Тао. Оказывается, моя формула непосредственно связана с формулой участника Null.