矩形波的傅里叶级数及代码

最新推荐文章于 2023-02-01 19:55:11 发布

November丶Chopin

最新推荐文章于 2023-02-01 19:55:11 发布

阅读量4.6k

点赞数 9

分类专栏：专栏04-数学知识文章标签： matploylib 傅里叶级数

本文链接：https://blog.csdn.net/u012762410/article/details/127021849

版权

专栏04-数学知识专栏收录该内容

2 篇文章

订阅专栏

问题：求如下矩形波的傅里叶级数

在这里插入图片描述
可见次矩形波的周期 $T = 1$ 。用傅里叶级数在区间 $(0, 1)$ 内拟合矩形波 $f (t)$ 。

傅里叶级数公式

傅里叶级数的公式为：
在这里插入图片描述
其中， $t_0$ 表示周期起点，可以设置为0， $\omega={\frac {2\pi} T}$ 。

下面求 $a_n$ 和 $b_n$ 。
设 $t_0=0$ ，
$\begin{aligned} &a_0={\frac 2 T}\int_{0}^{T}f(t)cos(n\omega t)dt\\\ &\quad=2\int_{0}^{1/2}dt \\\ &\quad=1 \end{aligned}$
$\begin{aligned} &a_n={\frac 2 T}\int_{0}^{T}f(t)cos(n\omega t)dt\\\ &\quad={\frac 2 T}\int_{0}^{T/2}cos(n\omega t)dt \qquad\qquad 因为f(x)只在这个区间有值 \\\ &\quad={\frac 2 T}{\frac 1 {n\omega}} sin(n\omega t)\Big|_{0}^{T/2} \qquad\qquad 定积分 \\\ &\quad={\frac 1 {n\pi}}sin(n\pi) \end{aligned}$
$\begin{aligned} &b_n={\frac 2 T}\int_{0}^{T}f(t)sin(n\omega t)dt\\\ &\quad={\frac 2 T}\int_{0}^{T/2}sin(n\omega t)dt \qquad\qquad 同理f(x)只在这个区间有值 \\\ &\quad={\frac 2 T}(-{\frac 1 {n\omega}}) cos(n\omega t)\Big|_{0}^{T/2} \qquad\qquad 定积分 \\\ &\quad={\frac 1 {n\pi}}(1-cos(n\pi)) \end{aligned}$
所以：
$\begin{aligned} &f(t)=a_0/2+\displaystyle\sum_{n=1}^\infin\left(\LARGE{ a_ncos(n\omega t) + b_nsin(n\omega t) }\right)\\\ &\qquad=1/2+\displaystyle\sum_{n=1}^\infin\left(\LARGE{ {\frac 1 {n\pi}}sin(n\pi)cos(n\omega t) + {\frac 1 {n\pi}}(1-cos(n\pi))sin(n\omega t) }\right)\\\ &\qquad=1/2+\displaystyle\sum_{n=1}^\infin{\frac 1 {n\pi}}\left(\LARGE{ sin(2n\pi t)+sin(n\pi(1-2t)) }\right)\\\ \end{aligned}$

绘图

根据上面求得的 $f (x)$ ，对其进行绘图，设置不同的正弦波个数n，可以获得不同的拟合程度，n越大拟合的越好。

import matplotlib.pyplot as plt
import numpy as np

def polt(n=10):
    A0=1
    N = np.arange(1,n)
    F = []
    for t in np.linspace(0,0.99,100):
        n_pi = N*np.pi
        S = 1/n_pi*(np.sin(2*n_pi*t)+np.sin(n_pi*(1-2*t)))
        f_t = A0/2+np.sum(S)
        F.append(f_t)

    t = np.linspace(-1.1,1.1,1000)

    # 方波
    rec_wav = (abs(2*t%2)<1).astype('i')

    plt.figure(figsize=(10,10))
    ax = plt.gca()
    ax.set_aspect(1) # 设置坐标轴等长
    ax.spines['right'].set_color('none')
    ax.spines['top'].set_color('none')
    ax.xaxis.set_ticks_position('bottom')
    ax.spines['bottom'].set_position(('data',0))
    ax.yaxis.set_ticks_position('left')
    ax.spines['left'].set_position(('data',0))

    ax.plot(np.linspace(0,0.99,100), F)
    ax.plot(t,rec_wav)
    plt.show()