Search code examples
pythonsympytaylor-series

Trouble with taylor series in python (sympy)


I'm trying to get the Taylor series for this function

Which should be similar to this, considering that d is centered or around rs

However when I try to take the example of @Saullo for my problem, enter image description here

As you can see the result is eliminating "d" from the series of Taylor, which should not be my goal.

Another additional info about the function in fact is: enter image description here

I'm doing something wrong ??, is there a way to get my result without deleting "d" ??

Any help is appreciated

The code

Thank you for your response and interest in helping me, here is my code until nowdays @asmeurer

import sympy as sy
#import numpy as np

from sympy import init_printing
init_printing(use_latex=True)

# Define the variable and the function to approximate
z, d, r_s, N_e, r_t, r_s, r_b = sy.symbols('z  d r_s N_e r_t r_s r_b')

# Define W_model
def W_model(r_t=r_t, r_b=r_b, r_s=r_s, z=z):

    s_model = sy.sqrt(pow(r_t, 2) - pow(r_s*sy.sin(z), 2)) - sy.sqrt(pow(r_b, 2) - pow(r_s*sy.sin(z), 2))
    d_model = r_t - r_b

    STEC_approx = N_e * s_model
    VTEC_approx = N_e * d_model

    return STEC_approx/VTEC_approx

f = W_model() 
# printing Standard model
f

# Some considerations for modify Standard model
rb = r_s - d/2
rt = r_s + d/2

f = W_model(r_b=rb, r_t=rt, r_s=r_s, z=z) 
# printing My model
f

## Finding taylor series aproximmation for W_model
num_of_terms = 2
# creates a generator
taylor_series = f.series(x=d, n=None)

# takes the number of terms desired for your generator
taylor_series = sum([next(taylor_series) for i in range(num_of_terms)])
taylor_series

Solution

  • The issue is that your expression is complicated enough that series doesn't know that the odd order terms are zero (you get complicated expressions for them, but if you call simplify() on them, they go to 0). Consider

    In [62]: s = f.series(d, n=None)
    
    In [63]: a1 = next(s)
    
    In [64]: a2 = next(s)
    
    In [65]: simplify(a0)
    Out[65]:
           rₛ
    ────────────────
       _____________
      ╱   2    2
    ╲╱  rₛ ⋅cos (z)
    
    In [66]: simplify(a1)
    Out[66]: 0
    

    If you print a0 and a1 they are both complicated expressions. In fact, you need to get several terms (up to a3) before series gets a term that doesn't cancel to 0:

    In [73]: simplify(a3)
    Out[73]:
          _____________
     2   ╱   2    2        2
    d ⋅╲╱  rₛ ⋅cos (z) ⋅sin (z)
    ───────────────────────────
               3    6
           8⋅rₛ ⋅cos (z)
    

    If you do f.series(d, n=3), it gives the expansion up to d**2 (n=3 means + O(d**3)). You can simplify the expression quite a bit using

    collect(expr.removeO(), d, simplify)
    

    Internally, when you give series an explicit n, it uses the term-by-term generator to get as many terms as it needs to give the proper O(d**n) expansion. If you use the generator yourself (n=None), you need to do this manually.

    In general, the iterator is not guaranteed to give you the next order term. If you want guarantees that you have all the terms, you need to provide an explicit n. The O term returned by series is always correct (meaning all the lower order terms are complete).