How to expand a composite function into series?

Question

I need to expand such a function $$g[y,z(x,y)]=\frac{-y (z+1)^4-z^4-4 z^3+8 z+8}{z+1},\tag{1}$$ into powers of $x$ and $y$. Among $x,y,z$ there is a constraint equation, for example $$(3 y+3) z^4+z^3 (4-8 x+12 y)+z^2 (18 y-24 x)+z (12 y-24 x)-8 x+3 y=0.\tag{2}$$ There is also the requirement, $z\rightarrow 0$ as $x,y \rightarrow 0$.

The most direct method is to solve Eq.(2) for $z$, and then substituting the $z=z(x,y)$ into Eq.(1) to obtain $f(x,y)$. At last, one can use, such as

Series[f[x,y],{x,0,2},{y,0,1}]

to expand the composite function as a series of $x$ and $y$.

The code for the two equations are:

g[y_,z_]:=(8 + 8 z - 4 z^3 - z^4 - y (1 + z)^4)/(1 + z);

and

-8 x + 3 y + (-24 x + 12 y) z + (-24 x + 18 y) z^2 + (4 - 8 x + 
 12 y) z^3 + (3 + 3 y) z^4 == 0;

I tried this direct method. It takes a very long time, but I do not get a result. If the constraint equation is not only a quartic equation, but a more complicated algebraic equation, the direct method may not work. I want to know whether there are simpler methods to deal with this kind of problems?

Additional remarks：

According to @Daniel Lichtblau's method, one can first expand the implicit function equation (the constraint equation) into a series of $x,y$. The key problem here is that we cannot first expand $z=z(x,y)$ formally as a standard Taylor series, and then substitute it into the constraint equation to solve the coefficients, because the true series of $z=z(x,y)$ may have minus power, even fractional power of $x$ or $y$. Below is an example: a simple quadratic equation,

-x + y + (4 x + y) z - (3 + 2 y) z^2==0

We certainly can solve it to obtain $$ z=\frac{\sqrt{16 x^2-12 x+9 y^2+12 y}-4 x-y}{2 (-2 y-3)}. $$ (* the first root *) Expanding the solution, Series[z, {x, 0, 2}, {y, 0, 1}], one can obtain such a result $$ \left(-\frac{\sqrt{y}}{\sqrt{3}}+\frac{y}{6}+O\left(y^{3/2}\right)\right)+x \left(\frac{1}{2 \sqrt{3} \sqrt{y}}+\frac{2}{3}-\frac{25 \sqrt{y}}{48 \sqrt{3}}-\frac{4 y}{9}+O\left(y^{3/2}\right)\right)+x^2 \left(\frac{1}{8 \sqrt{3} y^{3/2}}-\frac{19 \sqrt{3}}{64 \sqrt{y}}+\frac{333 \sqrt{3} \sqrt{y}}{1024}+O\left(y^{3/2}\right)\right)+O\left(x^3\right).$$

I also tried the methods on this post: [Calculating Taylor polynomial of an implicit function given by an equation][1] [1]: Calculating Taylor polynomial of an implicit function given by an equation, it does not work.

Answer 1

One possibility is to scale x and y and to then find the series in terms of the scaling parameter. For instance:

constraint = -8 x + 3 y + (-24 x + 12 y) z + (-24 x + 18 y) z^2 + (4 - 8 x + 12 y) z^3 + (3 + 3 y) z^4 == 0;
roots = Solve[constraint /. {x -> s x, y -> s y}, z, Quartics->False]

{{z -> Root[-8 s x + 3 s y + (-24 s x + 12 s y) #1 + (-24 s x + 18 s y) #1^2 + (4 - 8 s x + 12 s y) #1^3 + (3 + 3 s y) #1^4 &, 1]}, {z -> Root[-8 s x + 3 s y + (-24 s x + 12 s y) #1 + (-24 s x + 18 s y) #1^2 + (4 - 8 s x + 12 s y) #1^3 + (3 + 3 s y) #1^4 &, 2]}, {z -> Root[-8 s x + 3 s y + (-24 s x + 12 s y) #1 + (-24 s x + 18 s y) #1^2 + (4 - 8 s x + 12 s y) #1^3 + (3 + 3 s y) #1^4 &, 3]}, {z -> Root[-8 s x + 3 s y + (-24 s x + 12 s y) #1 + (-24 s x + 18 s y) #1^2 + (4 - 8 s x + 12 s y) #1^3 + (3 + 3 s y) #1^4 &, 4]}}

(I used the Quartics->False option because Series works much better with Root objects)

The first root is the branch where z == -4/3:

roots[[1]] /. s->0

{z -> -(4/3)}

So, you are interested in the other 3 roots:

roots[[2 ;; 4]] /. s->0

{{z -> 0}, {z -> 0}, {z -> 0}}

Mathematica can find the series expansion of Root objects. For example:

zser = Series[z /. roots[[2]], {s, 0, 1}];
zser //TeXForm

Root::sbr: Because of branch cuts, the series may represent a different root of -8 s x+3 s y+(-24 s x+12 s y) #1+(-24 s x+18 s y) #1^2+(4-8 s x+12 s y) #1^3+(3+3 s y) #1^4& for some values of {x,y,s}.

$\frac{\sqrt[3]{s} \sqrt[3]{8 x-3 y}}{2^{2/3}}+\frac{s^{2/3} (24 x-13 y)}{8 \sqrt[3]{2} \sqrt[3]{8 x-3 y}}+\frac{1}{64} s (56 x-57 y)+O\left(s^{4/3}\right)$

Plugging this into your g expression:

(8 + 8 z - 4 z^3 - z^4 - y (1 + z)^4)/(1 + z) /. {x -> s x, y -> s y, z -> zser} //Simplify //TeXForm

$8+s (2 y-8 x)+O\left(s^{4/3}\right)$

Answer 2

Basically the same approach as in this recent MSE thread.

The derivatives can be found iteratively using the constraint.

The set-up:

g[y_, z_] := (8 + 8 z - 4 z^3 - z^4 - y (1 + z)^4)/(1 + z);
poly[x_, y_, z_] := -8 x + 
   3 y + (-24 x + 12 y) z + (-24 x + 18 y) z^2 + (4 - 8 x + 
      12 y) z^3 + (3 + 3 y) z^4;
zconstraint = poly[x, y, z[x, y]];

First solve for z at x=0,y=0. We use the constraint only (and similar when solving for derivatives).

zsolns = 
 Solve[(zconstraint /. {x -> 0, y -> 0}) == 0, z[0, 0]]

(* Out[659]= {{z[0, 0] -> -(4/3)}, {z[0, 0] -> 0}, {z[0, 0] -> 
   0}, {z[0, 0] -> 0}} *)

We'll work with the branch hitting -4/3.

zderivsolns = 
 Solve[{D[zconstraint, x] == 0, 
     D[zconstraint, y] == 0}, {D[z[x, y], x], D[z[x, y], y]}] /. {x ->
      0, y -> 0} /. zsolns[[1]]

(* Out[660]= {{Derivative[1, 0][z][0, 0] -> 1/24, 
  Derivative[0, 1][z][0, 0] -> 1/192}} *)

Iterate to get the higher order derivatives needed for the series.

zderiv2solns = 
 Solve[{D[zconstraint, {x, 2}] == 0, 
      D[zconstraint, {x, 1}, {y, 1}] == 0, 
      D[zconstraint, {x, 2}, {y, 1}] == 0}, {D[z[x, y], {x, 2}], 
      D[z[x, y], {x, 1}, {y, 1}], 
      D[z[x, y], {x, 2}, {y, 1}]}] /. {x -> 0, y -> 0} /. 
   zsolns[[1]] /. zderivsolns

(* Out[686]= {{{Derivative[2, 0][z][0, 0] -> -(3/128), 
   Derivative[1, 1][z][0, 0] -> 
         -(11/3072), Derivative[2, 1][z][0, 0] -> 135/32768}}} *)

Now take a series for g, convert to a Taylor polynomial, and plug in values for z and the needed partial derivatives.

ser = Normal[Series[g[y, z[x, y]], {x, 0, 2}, {y, 0, 1}]];
ser /. zsolns[[1]] /. zderivsolns /. zderivsolns /. zderiv2solns

(* Out[695]= {{{{-(296/27) + x^2 (1/9 - (13 y)/1152) + 
     x (-(8/9) + y/36) - (2 y)/27}}}} *)