Exponential Generating Function
Peter Luschny
# It is often useful to describe a sequence a_n by a generating function.
# Now given an *exponential* generating function A(x) for a_n, how can we
# recover the sequence efficiently?
# Here is my solution, but there are certainly better ones, which I would
# like to learn!
egf2seq := proc(f,n) local i,h,l,u;
u := series('f(x)',x,n);
l := convert(u,list);
h := [];
for i from 1 by 2 to nops(l)-2 do
h := [op(h),l[i]*l[i+1]!] od;
h end:
# Two test cases:
egf2seq(exp, 12);
egf2seq(x->ln(x+1), 12);
Regards Peter
Joe Riel
> # It is often useful to describe a sequence a_n by a generating function.
> # Now given an *exponential* generating function A(x) for a_n, how can we
> # recover the sequence efficiently?
> # Here is my solution, but there are certainly better ones, which I would
> # like to learn!
>
> egf2seq := proc(f,n) local i,h,l,u;
> u := series('f(x)',x,n);
> l := convert(u,list);
> h := [];
> for i from 1 by 2 to nops(l)-2 do
> h := [op(h),l[i]*l[i+1]!] od;
> h end:
>
Essentially the same, but more efficient, is
egf2seq := proc(f,n::posint)
local i,l,x;
l := convert(series(f(x),x,n),'list');
[seq](l[i]*l[i+1]!,i=1..nops(l)-2,2);
end proc:
--
Joe Riel
Peter Luschny
>> Peter Luschny <spam...@luschny.de> writes:
>> # It is often useful to describe a sequence a_n by a generating function.
>> # Now given an *exponential* generating function A(x) for a_n, how can we
>> # recover the sequence efficiently?
>> # Here is my solution, but there are certainly better ones, which I would
>> # like to learn!
> Essentially the same, but more efficient, is
> egf2seq := proc(f,n::posint)
> local i,l,x;
> l := convert(series(f(x),x,n),'list');
> [seq](l[i]*l[i+1]!,i=1..nops(l)-2,2);
> end proc:
Thank you Joe! But seq(..,..,2) is not available for
Maple V Release 5. So I have to use (for ... by 2 ...).
And my first version did not work in all cases as it should.
egf2seq(x->exp(x^2),12);
egf2seq(x->exp(x^2)*x,12);
egf2seq(x->exp(x^3),12);
To fix this, I changed it to:
egf2seq := proc(f,n::posint) local i,h,l,u,j;
u := series(f(x),x,n); l := convert(u,list);
h := []; j:= 0;
if l[2] = 1 then h := [op(h),0]; j := j + 1 fi;
for i from 1 by 2 to nops(l)-2 do
while j < l[i+1] do h := [op(h),0]; j := j + 1 od;
h := [op(h),l[i]*l[i+1]!];
j := j + 1;
od;
h end:
Regards Peter
Joe Riel
> But seq(..,..,2) is not available for
> Maple V Release 5. So I have to use (for ... by 2 ...).
I didn't realize you were using release 5.
> And my first version did not work in all cases as it should.
>
> egf2seq(x->exp(x^2),12);
> egf2seq(x->exp(x^2)*x,12);
> egf2seq(x->exp(x^3),12);
>
> To fix this, I changed it to:
>
> egf2seq := proc(f,n::posint) local i,h,l,u,j;
> u := series(f(x),x,n); l := convert(u,list);
> h := []; j:= 0;
> if l[2] = 1 then h := [op(h),0]; j := j + 1 fi;
> for i from 1 by 2 to nops(l)-2 do
> while j < l[i+1] do h := [op(h),0]; j := j + 1 od;
> h := [op(h),l[i]*l[i+1]!];
> j := j + 1;
> od;
> h end:
While it won't help you, this can be done in later versions of Maple with
egf2seq4 := proc(f,n::posint)
local l,x,i;
uses PolynomialTools;
l := CoefficientList(convert(series(f(x),x,n),'polynom'),x);
[seq](l[i]*(i-1)!, i=1..nops(l));
end proc:
PolynomialTools:-CoefficientList does pretty much what it sounds like:
returns a list of all the coefficients of a given polynomial, starting
with the constant term.
--
Joe Riel
Vladimir Bondarenko
PL> It is often useful to describe a sequence a_n by a
PL> generating function. Now given an *exponential*
PL> generating function A(x) for a_n, how can we recover
PL> the sequence efficiently?
PL> egf2seq := proc(f,n::posint) local i,h,l,u,j;
PL> u := series(f(x),x,n); l := convert(u,list);
PL> h := []; j:= 0;
PL> if l[2] = 1 then h := [op(h),0]; j := j + 1 fi;
PL> for i from 1 by 2 to nops(l)-2 do
PL> while j < l[i+1] do h := [op(h),0]; j := j + 1 od;
PL> h := [op(h),l[i]*l[i+1]!];
PL> j := j + 1;
PL> od;
PL> h end:
Just for fun, some QA ideas... ;)
> egf2seq(x->sqrt(x)*exp(x),12);
Error, (in egf2seq) cannot determine if this expression is true or
false: 0 < x^(3/2)
> egf2seq(x->0,12);
Error, (in egf2seq) invalid subscript selector
> egf2seq(x->ln(x/(1+x)),12);
[ln(x), -1, 1, -2, 6, -24, 120, -720, 5040, -40320, 362880, -3628800]
> egf2seq(x->x^x,12);
[1, ln(x), ln(x)^2, ln(x)^3, ln(x)^4, ln(x)^5, ln(x)^6, ln(x)^7,
ln(x)^8, ln(x)^9, ln(x)^10, ln(x)^11]
> egf2seq(x->product(n^x, n= 1..infinity),12);
[1, infinity, infinity, infinity, infinity, infinity, infinity,
infinity, infinity, infinity, infinity, infinity]
> egf2seq(x->product(2*n^x, n= 1..infinity),12);
[infinity, infinity, infinity, infinity, infinity, infinity,
infinity, infinity, infinity, infinity, infinity, infinity]
> length(egf2seq(x->product(1+n^x, n= 1..infinity),12));
> length(egf2seq(x->sqrt(x+2^x), 12));
> length(egf2seq(x->tanh(f(x)),12));
> length(egf2seq(x->tan(1+x)^csc(1+x),12));
> length(egf2seq(x->f(f(f(f(x)))),12));
7742
110151
215718
437064
1156941
Gruss,
Vladimir
P.S.
I was wondering, is there a way to drop you a
private message?
(I've send you some stuff at spamgr...@luschny.de
but 'fraid, grube means pit :)
My regular email address is
v b @ c y b e r t e s t e r . c o m
where the blanks are to bypass Google.
Peter Luschny
> Peter Luschny <spam...@luschny.de> writes:
>> egf2seq := proc(f,n::posint) local i,h,l,u,j;
>> u := series(f(x),x,n); l := convert(u,list);
>> h := []; j:= 0;
>> if l[2] = 1 then h := [op(h),0]; j := j + 1 fi;
>> for i from 1 by 2 to nops(l)-2 do
>> while j < l[i+1] do h := [op(h),0]; j := j + 1 od;
>> h := [op(h),l[i]*l[i+1]!];
>> j := j + 1;
>> od;
>> h end:
> While it won't help you, this can be done in later versions of Maple with
>
> egf2seq4 := proc(f,n::posint)
> local l,x,i;
> uses PolynomialTools;
> l := CoefficientList(convert(series(f(x),x,n),'polynom'),x);
> [seq](l[i]*(i-1)!, i=1..nops(l));
> end proc:
>
> PolynomialTools:-CoefficientList does pretty much what it sounds like:
> returns a list of all the coefficients of a given polynomial, starting
> with the constant term.
>
Thank you Joe! As you can see, it did help:
egf2seqJ := proc(f,n::posint) local l,p,x,i;
p := convert(series(f(x),x,n),'polynom');
l := [seq(coeff(p,x,i),i=0..n)];
[seq(l[i]*(i-1)!, i=1..nops(l))];
end:
Regards Peter
Peter Luschny
> Just for fun, some QA ideas... ;)
Great!
>> egf2seq(x->sqrt(x)*exp(x),12);
egf2seq is supposed to work for 'generating functions'. We have
to define this term first. The exponents should be nonnegative
integers because they are regarded as the 'indices' of the sequence.
If you agree to this, this is just another case of garbage in, garbage
out, which I would not fight. (Yes, error handling would be fine...)
>> egf2seq(x->0,12);
> Error, (in egf2seq) invalid subscript selector
Fixed with egf2seqJ.
>> egf2seq(x->ln(x/(1+x)),12);
> [ln(x), -1, 1, -2, 6, -24, 120, -720, 5040, -40320, 362880, -3628800]
>> egf2seq(x->x^x,12);
> [1, ln(x), ln(x)^2, ln(x)^3, ln(x)^4, ln(x)^5, ln(x)^6, ln(x)^7,..]
Again the question is, what is the function supposed to compute?
You cannot blame a fish for not being able to fly. Implied by my
definition is a sequence of constants (not depending on x).
However, note that my second version *does* give the above answers.
>> length(egf2seq(x->product(1+n^x, n= 1..infinity),12));
>> length(egf2seq(x->sqrt(x+2^x), 12));
>> length(egf2seq(x->tanh(f(x)),12));
>> length(egf2seq(x->tan(1+x)^csc(1+x),12));
>> length(egf2seq(x->f(f(f(f(x)))),12));
> 7742
> 110151
> 215718
> 437064
> 1156941
Nice :) But an unexpected reply by Maple has always to be expected.
What is really missing here is something like that: If the output
is more then, say, 20 lines, Maple should ask first: "Do you really
want to see all my symbols?" :))
Gruss Peter
Joe Riel
> egf2seq is supposed to work for 'generating functions'. We have
> to define this term first. The exponents should be nonnegative
> integers because they are regarded as the 'indices' of the sequence.
> If you agree to this, this is just another case of garbage in, garbage
> out, which I would not fight. (Yes, error handling would be fine...)
Seems reasonable to me. Here's a simple approach to error handling,
you'll want to modify it accordingly for Maple V R5.
egf2seq4 := proc(f,n::posint)
local l,x,i,p;
uses PolynomialTools;
p := convert(series(f(x),x,n),'polynom');
if not type(p, 'polynom')) then
error "cannot convert to a polynomial";
end if;
l := CoefficientList(p,x);
return [seq](l[i]*(i-1)!,i=1..nops(l));
end:
Joe
Peter Luschny
> Seems reasonable to me. Here's a simple approach to error handling,
> you'll want to modify it accordingly for Maple V R5.
> egf2seq4 := proc(f,n::posint)
> local l,x,i,p;
> uses PolynomialTools;
> p := convert(series(f(x),x,n),'polynom');
> if not type(p, 'polynom')) then
> error "cannot convert to a polynomial";
> end if;
> l := CoefficientList(p,x);
> return [seq](l[i]*(i-1)!,i=1..nops(l));
> end:
With Maple V R5:
egf2seqJE := proc(f,n::posint) local l,p,x,i;
p := convert(series(f(x),x,n),'polynom');
if not type(p, 'polynom') then
ERROR("Cannot convert ",p," to a polynomial!") fi;
l := [seq(coeff(p,x,i),i=0..n)];
[seq(l[i]*(i-1)!, i=1..nops(l))];
end:
Still not found a reason to upgrade :)
Regards Peter
Thomas Richard
> [...]
> What is really missing here is something like that: If the output
> is more then, say, 20 lines, Maple should ask first: "Do you really
> want to see all my symbols?" :))
A similar feature is available as "term elision" in the Standard
Worksheet Interface of Maple 9.5 and 10.
You can supply a threshold, a number of leading terms, and a number of
trailing terms.
Furthermore, there is a LargeExpressions package for hiding
subexpressions of more complicated expressions.
--
Thomas Richard
Maple Support
Scientific Computers GmbH
http://www.scientific.de