In a recent thread ( Linearize a polynomial in several variables? ) the question of linearization of a multinomial was posed.
Frank Kampas pointed out that Series has its drawbacks and that something like a Taylor-Series is missing. Frank said
In order to to a multi-variable Taylor series expansion, it’s
necessary to use the procedure Daniel describes, since Series does its
expansion sequentially in the variables. Why Mathematica doesn’t have
a TaylorSeries function is something I’ve wondered about for years.
Such a Taylor-Expansion may be written (see textbooks of calculus, x and h are vectors)
f[ x + h ] = Sum[ 1/j! ( h . \[Del] )^j f[x] , { j, 0 , n } ] + rest
A certain difficulty arises in handling the powers of the differential operator.
In the following code I try to give a procedure for building the taylor expansion for a multivariate function:
The function g1 handles products of differential operators giving back a list of powers of derivatives
( d[ k ] meaning differentiation to the k-th variable, d[ k ] ^n means the n-th derivative, g1 operates on terms like d[1], d[ j ]^m, d[ j ]^n d[ j ]^m and so forth )
g1[x_] := Module[{},
abl = Table[0, {xL}];
If[x[[0]] == d, m = x[[1]]; n = 1; abl[[m]] = n];
If[x[[0]] == Power, m = x[[1, 1]]; n = x[[2]]; abl[[m]] = n];
If[x[[0]] == Times, abl = g1 /@ (Plus @@ (x));];
abl
]
For example
In[54]:= g1[d[1] d[2]^3]
Out[54]= {1, 3, 0}
This list is later on used in Derivative.
Function g2 handles terms like ( h . [Del] )^j building the appropriate terms for Derivative, here abbreviated to der
In[55]:= g2[({h[1], h[2]}.{d[1], d[2]})^3]
Out[55]= der[3, 0, 0] h[1]^3 + 3 der[2, 1, 0] h[1]^2 h[2] + 3 der[1, 2, 0] h[1] h[2]^2 + der[0, 3, 0] h[2]^3
The main procedure uses these things to build the expression for a multivariate Taylor-Series
taylorseries[f_, x_, x0_, n_] := Module[{},
xL = Length[x];
hh = Table[h[j], {j, xL}];
dd = Table[d[j], {j, xL}];
ser = (f /. Thread[x :> x0]) + Sum[1/j! g2[(hh.dd)^j], {j, 1, n}];
ser /. (h[j_] :> (x[[j]] - x0[[j]])) /.
der[u__] :> (Derivative[u][f[[0]]][Sequence @@ x0])
]
Examples (without output)
tsa = taylorseries[fa[x, y, z], {x, y, z}, {x0, y0, z0}, 2]
tsb = taylorseries[fa[x, y, z], {x, y, z}, {0, 0, 0}, 3]
The aforementioned linearization of a polynomial is now easy
fpoly[x_, y_, z_] := gamma + 3 x + a y + z + x y + x y z
taylorseries[uu[x, y, z], {x, y, z}, {0, 0, 0}, 2];
% /. uu -> fpoly
Out[43]= gamma + 3 x + a y + x y + z
The term x y has to be dropped, but doesn’t appear if you use a 1 as 4th argument in taylorseries.
Unfortunately I couldn’t figure out why the procedure fails if fpoly is directly given to the program
In[21]:= taylorseries[fpoly[x, y, z], {x, y, z}, {0, 0, 0}, 2]
Out[21]= gamma + x + y + z
Obviously the coefficients 3 of x and a of y are lost. I am pretty sure that is due to the fact that taylorseries doesn’t extract such coefficients correctly when f is given explicitly.
Any discussion and remarks and suggestions are welcome