There is a nice method to compute the Narayana polynomials. With Maple we can write
P := n -> simplify(hypergeom([-n,-n+1], [2], 1/x));
seq(expand(x^k*P(k)), k=0..5);
and get the answer
1, x, x^2+x, x^3+3*x^2+x, x^4+6*x^3+6*x^2+x, x^5+10*x^4+20*x^3+10*x^2+x.
Trying the same with Sage
P = lambda n: simplify(hypergeometric([-n,-n+1],[2], 1/x))
[expand(x^k*P(k)) for k in (0..5)]
gives the answer
[1, x, ..., x^n*hypergeometric((-n, -n-1),(2),1/x)]
This is really disappointing. Is there a workaround?
