python フ―リエ級数展開した式をグラフに出力する。

Question

元の関数ｆ（ｔ）＝１（ーπ<t<=0）かつ　０（０<t<=π）のフ―リエ級数展開した式をグラフに出力したいです。
フ―リエ級数展開は1/2+(Σ（k=1~∞）1-(-1)のk乗*sin(nπx)/nπ)
def f(t,k)の所を何度も変えても波は重なりません。ただ、最初の1.0を無くすと重なりますが、グラフの形が違う上に、フ―リエ級数展開は1.0もシグマの中に含んでいるので、絶対この位置だと思います。

# -*- coding: utf-8 -*-
import numpy as np
import matplotlib.pyplot as plt

def f(t,k):
    return -(1.0-np.power(-1,k)*np.sin(k*np.pi*t))/k*np.pi
    #return -np.sin((2.0*k-1.0)*t)/(2.0*k-1.0)*np.pi

dt =  0.01
start   = -2.0*np.pi
end     =  2.0*np.pi

# 横軸 t の作成．-6.28 から 6.28 を0.01 おきに．t は横ベクトル．
t = np.arange(start, end, dt)  

K=1
s1 = (1.0/2.0)*np.ones(t.size)
for k in range(K):
    s1 = s1 + f(t,k+1)
K=3
s3 = (1.0/2.0)*np.ones(t.size)
for k in range(K):
    s3 = s3 + f(t,k+1)

K=10
s10 = (1.0/2.0)*np.ones(t.size)
for k in range(K):
    s10 = s10 + f(t,k+1) # k+1 にしているのは，ranng(5) は 0 はじまりのため．

K=100
s100 = (1.0/2.0)*np.ones(t.size)
for k in range(K):
    s100 = s100 + f(t,k+1)

K=1000
s1000 = (1.0/2.0)*np.ones(t.size)
for k in range(K):
    s1000 = s1000 + ft,k+1)


fig=plt.figure(0)
plt.plot (t, s3, linewidth=1.0, color="r",linestyle="solid",label="$ K=3 $")
plt.plot (t, s10, linewidth=1.0, color="g",linestyle="solid",label="$ K=10 $")
plt.plot (t, s100, linewidth=1.0, color="b",linestyle="solid",label="$ K=100 $")
plt.plot (t, s1000, linewidth=1.0, color="m",linestyle="solid",label="$ K=1000 $")
plt.grid()
plt.legend()
plt.show()
fig.savefig('f102.pdf')

Answer 1

元の関数で表される矩形波を表示したいってことですね。
下記の出力が期待結果であってますか？

# 修正前：return -(1.0-np.power(-1,k)*np.sin(k*np.pi*t))/k*np.pi

#修正後
return (-1) * (1.0 - np.power(-1,k)) * (np.sin(k * t) / (k * np.pi))

あと、下記が"("漏れでした。

s1000 = s1000 + ft,k+1)
                ^^^^^^^ f(t,k+1)

修正後の出力

元の関数ｆ（ｔ）＝１（ーπ<t<=0）かつ　０（０<t<=π）
フ―リエ級数展開は1/2+(Σ（k=1~∞）1-(-1)のk乗*sin(nπx)/nπ)

あまり数学は得意ではありませんが、フーリエ級数展開は、1/2 - (Σ（k=1~∞）(1-(-1)^k) * sin(kt) / (k * pi)) ですかね？

Comments

[tiroha]

回答ありがとうございます！
こんなようなグラフを出力したかったんです！ありがとうございます！
フ―リエ級数展開は惜しいです。
正しくは1/2 - (Σ（k=1~∞）(1-(-1)^k) * sin(nπt) / (n * pi))
sin(nt)のところにπが抜けています。
なので、return (-1) * (1.0 - np.power(-1,k)) * (np.sin(k * t) / (k * np.pi))
のnp.sin(k*np.pi*t)にしたら正しくグラフが出力されました！

Answer 2

かっこの位置が間違っているようです。

質問者様のコード

    return -(1.0-np.power(-1,k)*np.sin(k*np.pi*t))/k*np.pi

正しそうなコード

    return -(1.0-np.power(-1,k))*np.sin(k*np.pi*t)/k*np.pi

Comments

[tiroha]

回答ありがとうございます！
（）の位置でこんなに変わるんですね！