IIR フィルターを使って sin、cos を再帰的に求める / Python で再現する

参考記事: Interface(インターフェース) 2017年10月号, p.159
　
計算手順はわかるが理窟は全然わからない。
　
手順は下のとおりである。

1. a = 2 cos(Δθ)、c = -cos(Δθ)、s = sin(Δθ) の 3 つだけを組み込み函数か何かでちゃんと計算しておく。
2. u[n] = a * u[n-1] - u[n-2] という漸化式を考える。n は 0, 1, 2, ......。初期値は u[0] = 0、u[1] = 1 とする。
3. sin、cos は、この漸化式 u[n] からそれぞれ sin(n Δθ) = s * u[n]、cos(n Δθ) = u[n+1] + c * u[n] で求まる。n は 0, 1, 2, ......

import math
import matplotlib.pyplot as plt

def SinCos(Δθ):
    A = 2 * math.cos(Δθ)
    C = -math.cos(Δθ)
    S = math.sin(Δθ)
    u = [0, 1]
    def sin_cos():
        sin = S * u[0]
        cos = u[1] + C * u[0]
        u[0], u[1] = u[1], A * u[1] - u[0]
        return sin, cos

    return sin_cos


########
# test #
########
Δθ      = math.pi/2**3
sincos  = SinCos(Δθ)
θList   = []
sinList = []
cosList = []

for i in range(50):
    sin, cos = sincos()

    θList.append(i * Δθ)
    sinList.append(sin)
    cosList.append(cos)
    #print("θ: %.2f, sin: %.6f, cos: %.6f, mag: %.16f" % (i*Δθ, sin, cos, np.sqrt(sin**2 + cos**2)))
    #print(i*Δθ, sin, cos, sin**2 + cos**2, sep="\t")

plt.plot(θList,sinList,
         θList,cosList)
plt.xlabel("θ")
plt.grid()
plt.show()

実行結果: