# How to get series coefficients from functional equations?

$p(z),t(z)$ are two mutually defined functional equations, while $\widehat{G}(z)$ is the exponential generation function of A182173 (maybe, I am not sure...lol):

$$\begin{cases} p(z)=e^{t(z)}-t(z)+2 z-1\\ t(z)=2\left(e^{p(z)}-e^{p(z)/2}+z\right)-p(z)\\ \widehat{G}(z)=p(z)+t(z)-2 z \end{cases}$$

Is there any way to get series coefficients effectively?

p[z]->2z+E^t[z]-1-t[z]
t[z]->2*( z+ E^p[z] -E^(p[z]/2) ) - p[z]
Ge=p[z]+t[z]-2z
A182173[n_]:=n!SeriesCoefficient[Ge,{z,0,n}]

• For completeness, do you have a reference for the expressions for $p$, $t$, and $\widehat{G}$? – J. M.'s ennui Mar 4 '18 at 14:51
• @J.M. 四则运算表达式计数问题, with some methods from 《Analytic Combinatorics》. – GalAster Mar 5 '18 at 1:28

Find Ge[z]:

r1 = p[z] -> 2 z + E^t[z] - 1 - t[z]

r2 = t[z] -> 2*(z + E^p[z] - E^(p[z]/2)) - p[z]

Ge[z_] = p[z] + t[z] - 2 z /. r1 /. r2 // FullSimplify

(*   -1 + E^(2 (-E^((p[z]/2)) + E^p[z] + z) - p[z])   *)


SeriesCoefficient produces derivatives of p at z == 0. Therefore we have to solve for the p[0], p'[0], p''[0] ...

der1[i_] := Derivative[i][p][0] == Derivative[i][2 # + E^t[#] - 1 - t[#] &][0]

der2[i_] :=
Derivative[i][t][0] == Derivative[i][2*(# + E^p[#] - E^(p[#]/2)) - p[#] &][0]

sol[0] = First @
Solve[{der1[0], der2[0]}, {Derivative[0][p][0],
Derivative[0][t][0]}, Reals]

(*   {p[0] -> 0, t[0] -> 0}   *)


Lower derivatives are needed to solve for higher ones:

sol[j_ /; j > 0] := sol[j]=
First@Solve[{der1[j], der2[j]} //.
Flatten[Table[sol[k], {k, 0, j - 1}]], {Derivative[j][p][0],
Derivative[j][t][0]}, Reals]


Inserting these found derivatives of p into the series coefficients gives the desired A numbers.

A[n_] := n! SeriesCoefficient[Ge[z], {z, 0, n}] //.
Flatten[Table[sol[i], {i, 0, n}]]

Table[A[i], {i, 1, 10}]

(*   {2, 10, 94, 1466, 31814, 887650, 30259198, 1218864842, 56644903958,
2983300619410}   *)


As Akku14 notes in his answer, one can determine from the conditions that both functions $p,t$ have a zero constant term. With that additional piece of information, we can use SolveAlways[] to determine the rest of the required coefficients:

n = 10; (* desired expansion order *)
p = Sum[pc[k] z^k, {k, 0, n}] + O[z]^(n + 1);
t = Sum[tc[k] z^k, {k, 0, n}] + O[z]^(n + 1);

vals = SolveAlways[{p == Exp[t] - t + 2 z - 1, t == 2 (Exp[p] - Exp[p/2] + z) - p} /.
{pc[0] -> 0, tc[0] -> 0}, z];

Rest[CoefficientList[Sum[(pc[k] + tc[k]) z^k, {k, 0, n}] - 2 z + O[z]^(n + 1) /.
Join[{pc[0] -> 0, tc[0] -> 0}, First[vals]], z]] Range[n]!
{2, 10, 94, 1466, 31814, 887650, 30259198, 1218864842, 56644903958,
2983300619410}