Groups | Search | Server Info | Login | Register

Groups > comp.soft-sys.math.maple > #1375

Re: Bug in series expansion

Show all headers | View raw

On Wednesday, April 19, 2023 at 7:42:05 AM UTC-4, peter....@gmail.com wrote:
> > Your use of `expand` is inadequate in that case.
> Dear acer, 
> could you please explain this in more detail?

The result series(egf1, z, 9) is of type `polynom(anything,z)`, but its coefficients are not of type `polynom(anything,x)`. They are rational polynomials with `x+1` as denominator. 

egf1 := (1 + x*exp(x*z)*exp(z))/(x + 1):
series(egf1, z, 9);
series(1+(x^2+x)/(x+1)*z+(1/2*x+x^2+1/2*x^3)/(x+1)*z^2+(1/6*x+1/2*x^2+1/2*x^3+1
/6*x^4)/(x+1)*z^3+(1/24*x+1/6*x^2+1/4*x^3+1/6*x^4+1/24*x^5)/(x+1)*z^4+(1/120*x+
1/24*x^2+1/12*x^3+1/12*x^4+1/24*x^5+1/120*x^6)/(x+1)*z^5+(1/720*x+1/120*x^2+1/
48*x^3+1/36*x^4+1/48*x^5+1/120*x^6+1/720*x^7)/(x+1)*z^6+(1/5040*x+1/720*x^2+1/
240*x^3+1/144*x^4+1/144*x^5+1/240*x^6+1/720*x^7+1/5040*x^8)/(x+1)*z^7+(1/40320*
x+1/5040*x^2+1/1440*x^3+1/720*x^4+1/576*x^5+1/720*x^6+1/1440*x^7+1/5040*x^8+1/
40320*x^9)/(x+1)*z^8+O(z^9),z,9)

However, if `normal` is applied to them then the results are of type `polynom(anything,x)` (due to factorization). The calls `coeff(c(n), x, k)` will error out if those coefficients have not been thusly simplified/normalized, because those `coeff` calls expect a form that is of type `polynom(anything,x)`.

normal(series(egf1, z, 9));
series(1+x*z+1/2*x*(x+1)*z^2+1/6*(x^2+2*x+1)*x*z^3+1/24*x*(x^3+3*x^2+3*x+1)*z^4
+1/120*x*(x^4+4*x^3+6*x^2+4*x+1)*z^5+1/720*x*(x^5+5*x^4+10*x^3+10*x^2+5*x+1)*z^\
6+1/5040*x*(x^6+6*x^5+15*x^4+20*x^3+15*x^2+6*x+1)*z^7+1/40320*x*(x^7+7*x^6+21*x
^5+35*x^4+35*x^3+21*x^2+7*x+1)*z^8+O(z^9),z,9)

The same thing occurs if you take the coefficients one at a time. For example,

b := n -> n!*coeff(series(egf1, z, 9), z, n):
b(2);
2*(1/2*x+x^2+1/2*x^3)/(x+1)

type(b(2), polynom(anything,x));
                             false

coeff(b(2), x, 1);
Error, unable to compute coeff

Now, let's repeat that, but with `normal`. (Using `simplify` or `factor` would also work here.)

normal(b(2));
x*(x+1)

type(normal(b(2)), polynom(anything,x));
                              true

coeff(normal(b(2)), x, 1);
                               1

Or, as an adjustment to your original `c`,

c := n -> normal(n!*coeff(series(egf1, z, 9), z, n)):
c(2):lprint(%);
x*(x+1)

type(c(2), polynom(anything,x));
                              true

coeff(c(2), x, 1);
                               1

Your original `c` tried it using `expand` instead of `normal`. That is not sufficient, ie.

expand(b(2));
1/(x+1)*x+2/(x+1)*x^2+1/(x+1)*x^3

It happens that the coefficients of `series(egf, z, 9)` happen to come out in a form that is of type `polynom(anything,x)`.

egf := (1 + x*exp(x*z+z))/(x + 1):
series(egf, z, 9);
series(1+x*z+1/2*x*(x+1)*z^2+1/6*x*(x+1)^2*z^3+1/24*x*(x+1)^3*z^4+1/120*x*(x+1)
^4*z^5+1/720*x*(x+1)^5*z^6+1/5040*x*(x+1)^6*z^7+1/40320*x*(x+1)^7*z^8+O(z^9),z,\
9)

An alternative is to use `MultiSeries:-series`, producing the desired form of the coefficients for both `egf` and `egf1`.

MultiSeries:-series(egf1, z, 9);
series(1+x*z+1/2*x*(x+1)*z^2+1/6*(x^2+2*x+1)*x*z^3+1/24*x*(x^3+3*x^2+3*x+1)*z^4
+1/120*x*(x^4+4*x^3+6*x^2+4*x+1)*z^5+1/720*x*(x^5+5*x^4+10*x^3+10*x^2+5*x+1)*z^\
6+1/5040*x*(x^6+6*x^5+15*x^4+20*x^3+15*x^2+6*x+1)*z^7+1/40320*x*(x^7+7*x^6+21*x
^5+35*x^4+35*x^3+21*x^2+7*x+1)*z^8+O(z^9),z,9)

for n from 0 to 6 do seq(coeff(d(n), x, k), k = 0..n) od;
                               1

                              0, 1

                            0, 1, 1

                           0, 1, 2, 1

                         0, 1, 3, 3, 1

                        0, 1, 4, 6, 4, 1

                     0, 1, 5, 10, 10, 5, 1

Back to comp.soft-sys.math.maple | Previous | Next — Previous in thread | Next in thread | Find similar

Thread

Bug in series expansion "peter....@gmail.com" <peter.luschny@gmail.com> - 2023-04-18 03:31 -0700
  Re: Bug in series expansion acer <maple@rogers.com> - 2023-04-18 13:27 -0700
    Re: Bug in series expansion "peter....@gmail.com" <peter.luschny@gmail.com> - 2023-04-19 04:42 -0700
      Re: Bug in series expansion acer <maple@rogers.com> - 2023-04-19 07:03 -0700
        Re: Bug in series expansion "peter....@gmail.com" <peter.luschny@gmail.com> - 2023-04-19 12:03 -0700
        Re: Bug in series expansion "peter....@gmail.com" <peter.luschny@gmail.com> - 2023-05-01 00:56 -0700

csiph-web