用Python绘制傅里叶级数和泰勒级数逼近已知函数的动态过程

news2024/9/20 1:24:18

文章目录

- Taylor级数
- Fourier级数

本文代码： Fourier级数和Taylor级数对原函数的逼近动画

Taylor级数

级数是对已知函数的一种逼近，比较容易理解的是Taylor级数，通过多项式来逼近有限区间内的函数，其一般形式为

$f(x)=\sum_{n=0}^N a_nx^n$

其中最著名的应该是自然指数，根据其导数不变的特点，我们可以很容易得到其表达式

$e^x=\sum_{n=0}^N \frac{x^n}{n!}$

随着N的不断增加，其逼近过程如图所示

在这里插入图片描述
其中，Taylor级数的实现方法如下，除了exp函数之外，还包括sin和cos函数的Taylor级数。

def Taylor(x,funcType='exp',n=0):
    func = {
        'exp' : lambda x,n : x**n/fac(n),
        'sin' : lambda x,n : (-1)**n*x**(2*n+1)/fac(2*n+1),
        'cos' : lambda x,n : (-1)**n*x**(2*n)/fac(2*n)
    }
    return func[funcType](x,n)

绘图代码如下

def approxGif(num=30, funcType='exp'):
    funcType = func if type(func)==str else 'auto'
    if type(func)==str:
        func = funcDict[func]
    
    x = np.linspace(0,10,1000)
    Y =func if type(func)==type(x) else func(x)
    
    if method in ['Taylor','taylor']:
        y = Taylor(x,funcType,0)
    elif method in ['Fourier','fourier']:
        y = Fourier(x,funcType,0)
    num = range(num)

    #画图初始化
    fig = plt.figure()
    ax = fig.add_subplot(111,autoscale_on=False,
        xlim=(0,10),ylim=(min(Y)-0.5,max(Y)+0.5))
    ax.plot(x,Y,color='g',lw=0.2)
    ax.grid()

    line, = ax.plot([],[],lw=0.5)
    time_text = ax.text(0.1,0.9,'',transform=ax.transAxes)

    # 动画初始化
    def init():
        line.set_data([],[])
        time_text.set_text('level:'+str(0))
        return line, time_text

    # 动画迭代
    def animate(n):
        nonlocal y
        y = Taylor(x,funcType,n) if n==0 else y+Taylor(x,funcType,n)
        line.set_data(x,y)
        time_text.set_text('level:'+str(n))
        print(n)
        return line, time_text
    
    ani = animation.FuncAnimation(fig,animate,
        num,interval=200,blit=False,init_func=init)
    #ani.save(funcType+'.gif',writer='pillow')
    plt.show()

Fourier级数

Fourier级数也是本着相同的思维，只不过采用了不同频率的三角函数作为其空间中的基底。

对于以 $2\pi$ 为周期的方波信号

$f(x)=\left\{\begin{aligned} 1\quad &x \in[0,\pi)\\ -1\quad &x\in[-\pi,0) \end{aligned}\right.$

其Fourier级数为

$f(x)=\frac{4}{\pi}\sum_{n=0}^N\frac{\sin x}{2n+1}$

实现为

def Fourier(x,funcType='square',n=0):
    func = {
        'square' : lambda x,n : 4/np.pi*np.sin((2*n+1)*x)/(2*n+1),
        'tri' : lambda x,n: np.pi/2 if n == 0 \
            else -4/np.pi*np.cos((2*n-1)*x)/(2*n-1)**2,
        'oblique': lambda x,n : 2*np.sin((n+1)*x)/(n+1)*(-1)**n
    }
    return func[funcType](x,n)

绘图代码为

def square(x):
    x = np.mod(x,np.pi*2)
    x[x>np.pi] = -1
    x[x!=-1] = 1
    return x 

def tri(x):
    return np.pi-np.abs(np.mod(x,2*np.pi)-np.pi)

funcDict = {
    'exp':np.exp,
    'sin':np.sin,
    'cos':np.cos,
    'square':square,
    'tri':tri,
}

# func支持三种输入模式，即字符串，函数以及numpy数组
def approxGif(func='square',method='fourier',num=30):
    funcType = func if type(func)==str else 'auto'
    if type(func)==str:
        func = funcDict[func]
    
    x = np.linspace(0,10,1000)
    Y =func if type(func)==type(x) else func(x)
    
    if method in ['Taylor','taylor']:
        y = Taylor(x,funcType,0)
    elif method in ['Fourier','fourier']:
        y = Fourier(x,funcType,0)
    num = range(num)

    #画图初始化
    fig = plt.figure()
    ax = fig.add_subplot(111,autoscale_on=False,
        xlim=(0,10),ylim=(min(Y)-0.5,max(Y)+0.5))
    ax.plot(x,Y,color='g',lw=0.2)
    ax.grid()

    line, = ax.plot([],[],lw=0.5)
    time_text = ax.text(0.1,0.9,'',transform=ax.transAxes)

    # 动画初始化
    def init():
        line.set_data([],[])
        time_text.set_text('level:'+str(0))
        return line, time_text

    # 动画迭代
    def animate(n):
        nonlocal y
        if method in ['taylor','Taylor']:
            y = Taylor(x,funcType,n) if n==0 \
                else y+Taylor(x,funcType,n)
        elif method in ['fourier','Fourier']:
            y = Fourier(x,funcType,n) if n==0 \
                else y+Fourier(x,funcType,n)
            pass
        line.set_data(x,y)
        time_text.set_text('level:'+str(n))
        print(n)
        return line, time_text
    
    ani = animation.FuncAnimation(fig,animate,
        num,interval=200,blit=False,init_func=init)
    #ani.save(funcType+'.gif',writer='pillow')
    plt.show()