Fastest way to compute the coefficients of this generating function?

res0001
55 ●3 ●5 ●10

I need the first N coefficients of a certain generating function, with N around 30000. Using the following code, about 2.8 hrs were needed for N=10000, and 20 hrs for N=20000; unfortunately, it looks like N=30000 would require days:

N = 30000
R.<z> = PowerSeriesRing(ZZ, default_prec=N+1)
gf = (1-z)^2 * sum( z^n/(1 - 2*z + z^(n+1))^2 for n in range(N))
coeffs = gf.coefficients()

Question: Is there a better (faster!) way to do this in SageMath? E.g., some alternative to PowerSeriesRing, or some improvement to how I'm using it?

(Aside: The coefficients of this gf are theoretically expected to exhibit a certain phenomenon that doesn't show up until N is around 30000, and I'm just interested in confirming it.)

Comments

It is possible to derive a formula for a particular coefficient of gf. However, getting them all unavoidably requires ~ $N^2$ operations, which would be time-consuming.

Max Alekseyev ( 1 year ago )

add a comment

answered 1 year ago

Max Alekseyev

6932 ●9 ●52 ●162

updated 1 year ago

Let $c_k$ be the coefficient of $z^k$ in sum( z^n/(1 - 2*z + z^(n+1))^2 for n in range(N)), then for $0 < k < N$ we have $c_k = [z^k] \sum_{n\geq 1} z^{n} (1-2z+z^{n+1})^{-2} = [z^k] \sum_{n\geq1} \sum_{m\geq0} z^{n+m} \binom{-2}{m} (-2+z^n)^m$ $=\sum_{n=1}^k \sum_{m=0\atop m\equiv k\pmod{n}}^{k-n} \binom{-2}{m} \binom{m}{\frac{k-n-m}n}(-2)^{m-\frac{k-n-m}n}.$ Here is a corresponding code:

def c(k):
    return 1 if k==0 else sum(binomial(-2,m) * binomial(m,(t:=(k-m-n)//n)) * (-2)^(m-t) for n in range(1,k+1) for m in range(k%n,k-n+1,n))

Then, the coefficient of $z^k$ in (1-z)^2 * sum( z^n/(1 - 2*z + z^(n+1))^2 for n in range(N)) can be obtained simply as c(k) - 2*c(k-1) + c(k-2).

link

Comments

Thanks for this answer! I hadn't realized one might find such a formula. Using your c(k) in a program to compute the final coefficients for a selected range of k-values, I find that it takes about 7 minutes to compute the ten of them for k=30000..30009. (I'll have to study the method you used to derive that formula. BTW, I'm sure you intended to write $[z^k]$ .)

res0001 ( 1 year ago )

Yes, it was a typo. Now fixed.

Max Alekseyev ( 1 year ago )

add a comment

Your Answer

Please start posting anonymously - your entry will be published after you log in or create a new account.

Fastest way to compute the coefficients of this generating function?

Comments

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

Fastest way to compute the coefficients of this generating function? savecancel

Comments

1 Answer

Comments

Your Answer

Question Tools

Stats

Related questions

Fastest way to compute the coefficients of this generating function?