Power series expansion

[Julie]

Hi all,

I am attempting to obtain coefficients of a generating function to
obtain probabilites, but in order to obtain the coefficients, I first
need to expand a power series, which is necessary for my paricular
function.
Is there a simple way to expand such a series in sage e.g. 2^x?

(For my exact problem, the generating function contains 2 variables (p
and y), and when expanded up to terms where p=1 for example, I have
the formula:
(0.030*0.248244^y)y+0.05721*(0.248244^y)p +0.08838*(0.248244^y)

Thus, before being able to extract the coefficients of p^0*y^0, p1,
y1,p*y etc, I need to expand 0.248244^y as a power series - will the
same programming also hold for this problem?)

Many thanks,
Julie

[Carl Eberhart]

Hi Julie,

I'm new here, but if I understood your question, something like works to give the power series expansion of a nice function f

about x=a 

def taylor(f,a,n):

    tay=f(a)

    g = diff(f(x),x)

    for i in range(1,n):

        tay=tay+g(a)/factorial(i)*(x-a)^i    

    return tay

#now get the first 4 terms of the expansion of .248244^x about a=0   

y=var('y')

taylor((.248244)^x,0,4)(x=y)

Then in order to get the coefficients you want, you'll need to multiply out and collect like terms

Carl

Julie

--
To post to this group, send email to sage-s...@googlegroups.com
To unsubscribe from this group, send email to sage-support...@googlegroups.com
For more options, visit this group at http://groups.google.com/group/sage-support
URL: http://www.sagemath.org

[achrzesz]

I'm not familiar with generating functions but to extract coefficients
you can try:
sage: var('y p')
(y, p)
sage: s=taylor(0.248244^y,y,0,6)
sage: w=expand(0.030*s*y+0.05721*s*p +0.08838*s)
sage: w.polynomial(RR).coefficient([0,0])
0.0883800000000000
sage: w.polynomial(RR).coefficient([1,1])
-0.0797131613559558

Andrzej Chrzeszczyk

[Carl Eberhart]

Oops! I left out a line.  Sorry

def taylor(f,a,n):

    tay=f(a)

    g = diff(f(x),x)

    for i in range(1,n):

        tay=tay+g(a)/factorial(i)*(x-a)^i

        g = diff(g(x),x)    

    return tay

y=var('y')

taylor((.324)^x,0,4)(x=y)

Carl

On Fri, Dec 2, 2011 at 7:24 AM, Julie <juliewil...@googlemail.com> wrote:

Julie

[Julie]

Hi,

Thanks for such quick responses!
Unfortunately, having the Tayor series approach out, don't think it's
really appropriate for my problem afterall, as what I esentially need
to do is find the coefficientsof p^0*y^0, p, y, p^2*y etc in the
formula
(0.030*0.248244^y)y+0.05721*(0.248244^y)p +0.08838*(0.248244^y)
and was hoping that there may be a simple way of evaluating 0.248244^y
as 0.248244^0 + 0.248244^1*y^1+ 0.248244^2*y^2...as this would have
subsequently allowed me to look up the coefficients in the next step
Yet, the correct application of the Taylor series, as you pointed out,
is 0.248244^y = 0.248244^0 + ln(0.248244)*0.248244y + .... and so I'm
still in need of guidance for where to go next!

Thinking back to finding the coefficients, relating to Andrzej's
advice, if we call the formula:
i(p,y) = (0.030*0.248244^y)y+0.05721*(0.248244^y)p
+0.08838*(0.248244^y)
By using the function
i(p,y).coefficients()
sage outputs these coefficients as
[[0.03*0.248244^y*y + 0.08838*0.248244^y, 0], [0.05721*0.248244^y, 1]]
With more basic two-variable problems, I have used the command similar
to Andrzej's:
i(p,y).polynimial(SR).coefficient([0,0]) which outputs the specific
coefficient for p^0*y^0, and we coud then maybe substitute y=0 into
the coefficient after it is pulled out
But as there are still powers of y involved with the coefficients at
the primary stage, the i(p,y).polynimial(SR).coefficient([0,0])
command gives an error message.

With ony the i(p,y).coefficients() command, I can see by inspection
that to find the coefficients of p^0*y^0, I need first pull out the
1st bracket, as above, and evaluate 0.248244^0 to find the coefficient
And similarly to find the coefficients of p^0*y^1, I need first pull
out the 1st bracket, as above, and evaluate 0.248244^1 to find the
coefficient

Are you aware if there is a way to do this in sage?

Thank you,
Julie

[David Joyner]

On Fri, Dec 2, 2011 at 8:24 AM, Julie <juliewil...@googlemail.com> wrote:
> Hi all,
>
> I am attempting to obtain coefficients of a generating function to
> obtain probabilites, but in order to obtain the coefficients, I first
> need to expand a power series, which is necessary for my paricular
> function.
> Is there a simple way to expand such a series in sage e.g. 2^x?

sage: f = 2^x
sage: f.taylor(x,0,5)
1/120*x^5*log(2)^5 + 1/24*x^4*log(2)^4 + 1/6*x^3*log(2)^3 +
1/2*x^2*log(2)^2 + x*log(2) + 1

>
> (For my exact problem, the generating function contains 2 variables (p
> and y), and when expanded up to terms where p=1 for example, I have
> the formula:
> (0.030*0.248244^y)y+0.05721*(0.248244^y)p +0.08838*(0.248244^y)
>
> Thus, before being able to extract the coefficients of p^0*y^0, p1,
> y1,p*y etc, I need to expand 0.248244^y as a power series - will the
> same programming also hold for this problem?)

sage: x,y = var("x,y")
sage: f = 2^x*sin(y)
sage: f.taylor((x,0), (y,0) ,5)
1/24*x^4*y*log(2)^4 + 1/6*x^3*y*log(2)^3 - 1/12*x^2*y^3*log(2)^2 +
1/2*x^2*y*log(2)^2 - 1/6*x*y^3*log(2) + 1/120*y^5 + x*y*log(2) -
1/6*y^3 + y

Hope that helps.

[Anton Sherwood]

On 2011-12-02 08:17, Julie wrote:
> Unfortunately, having the Tayor series approach out, don't think it's
> really appropriate for my problem afterall, as what I esentially need
> to do is find the coefficientsof p^0*y^0, p, y, p^2*y etc in the
> formula
> (0.030*0.248244^y)y+0.05721*(0.248244^y)p +0.08838*(0.248244^y)

> [...]

Since this is not a polynomial in p and y, what does it mean to obtain
the coefficients of p^j y^k *if not* those of something very similar to
a Taylor series?

--
Anton Sherwood *\\* www.bendwavy.org *\\* www.zazzle.com/tamfang

[Julie]

Thank you very much to everyone for all your help.
I've now solved the issue I was having trouble with - the reason
finding the coefficients of the y terms didn't give me the required
results was because the generating function was really in terms of one
variable (p), not two, and required values to be substituted in for
y's rather than finding the various coefficients.

Now I have completed this work, wondered if you knew the best way to
export the outputted coefficients by the internet version of sage to
another package such as Excel? Or indeed, if this is even possible?

Many thanks,
Julie

On Dec 3, 8:57 pm, Anton Sherwood <bro...@pobox.com> wrote:
> On 2011-12-02 08:17, Julie wrote:> Unfortunately, having the Tayorseriesapproach out, don't think it's

> > really appropriate for my problem afterall, as what I esentially need
> > to do is find the coefficientsof p^0*y^0, p, y, p^2*y etc in the
> > formula
> > (0.030*0.248244^y)y+0.05721*(0.248244^y)p +0.08838*(0.248244^y)
>
>  > [...]
>
> Since this is not a polynomial in p and y, what does it mean to obtain

> thecoefficientsof p^j y^k *if not* those of something very similar to

[Robert Bradshaw]

If you create an actual power series element, you can easily write the
coefficients to a file:

sage: f = taylor(sin(x), x, 0, 10); f
1/362880*x^9 - 1/5040*x^7 + 1/120*x^5 - 1/6*x^3 + x
sage: power_series = RR[['x']](f); power_series
0.000000000000000 + 1.00000000000000*x + 0.000000000000000*x^2 -
0.166666666666667*x^3 + 0.000000000000000*x^4 +
0.00833333333333333*x^5 + 0.000000000000000*x^6 -
0.000198412698412698*x^7 + 0.000000000000000*x^8 +
2.75573192239859e-6*x^9
sage: list(power_series)
[0.000000000000000, 1.00000000000000, 0.000000000000000,
-0.166666666666667, 0.000000000000000, 0.00833333333333333,
0.000000000000000, -0.000198412698412698, 0.000000000000000,
2.75573192239859e-6]
sage: ",".join(str(a) for a in list(power_series))
'0.000000000000000,1.00000000000000,0.000000000000000,-0.166666666666667,0.000000000000000,0.00833333333333333,0.000000000000000,-0.000198412698412698,0.000000000000000,2.75573192239859e-6'
sage: open("file.txt", "w").write("\n".join(str(a) for a in list(power_series)))
sage: open("file.txt", "w").write("\n".join(str(a) for a in list(power_series)))
sage: !cat file.txt
0.000000000000000
1.00000000000000
0.000000000000000
-0.166666666666667
0.000000000000000
0.00833333333333333
0.000000000000000
-0.000198412698412698
0.000000000000000
2.75573192239859e-6

Note that ZZ[[x]], QQ[[x]], RR[[x]] are *much* faster to compute with as well

sage: R.<t> = PowerSeriesRing(ZZ, 't', 1000)
sage: f = 1/(t+1)
sage: f + O(t^10) # this is really O(t^1000), don't want to print it all
1 - t + t^2 - t^3 + t^4 - t^5 + t^6 - t^7 + t^8 - t^9 + O(t^10)
sage: timeit("f^100")
25 loops, best of 3: 28.2 ms per loop

- Robert