2
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?

More than 1 year has passed since last update.

Pythonで三角関数を自作してみた

Last updated at Posted at 2022-03-19

ライブラリを使うと

Pythonでは標準ライブラリで三角関数の計算ができます。

$ python3
Python 3.10.2 (main, Mar  9 2022, 20:22:18) [GCC 7.5.0] on linux
Type "help", "copyright", "credits" or "license" for more information.
>>> import math
>>> math.sin(0.1)
0.09983341664682815

これを使わずに三角関数を計算してみます。
理由は特にありません。暇つぶしです。

マクローリン展開

マクローリン展開を使うと、三角関数を多項式に置き換えることが出来ます。

f(x) = \sum_{k=0}^{\infty}f^{(k)}(0)\frac{x^k}{k!}

$\sin(x)$のマクローリン展開は以下の通りです。

\sin(x) = \sum_{k=0}^{\infty}(-1)^k\frac{x^{2k+1}}{(2k+1)!}

これを実装します。

sin.py
# -*- coding:utf-8 -*-

def fact(x):
    if x == 0: return 1
    else: return x * fact(x-1)

def sin(x):
    ans = 0.0
    for k in range(30):
        ans += (-1 if k%2 else 1) * (x**(2*k+1)) / fact(2*k+1)
    return ans

if __name__ == "__main__":
    import math
    print("自作", sin(0.1))
    print("math", math.sin(0.1))
    print("差分", abs(sin(0.1) - math.sin(0.1)))

"""
$ python3 sin.py
自作 0.09983341664682817
math 0.09983341664682815
差分 1.3877787807814457e-17
"""

それっぽい値がでました。

高階関数

Pythonでは関数も第一級オブジェクトですので、関数$g(n)=f^{(n)}(0)$を受け取って関数$f$を返すような関数を定義できます。引数または返り値に関数を含む関数を高階関数と言います。ここで、$f^{(n)}(x)$は$f(x)$のn次導関数とします。$f(x)=\sin(x)$とすると、$f^{(n)}(0)$は次のように表すことが出来ます。

f^{(n)}(0) = \left\{
\begin{array}{ll}
0 & (n \equiv 0 \pmod 4) \\
1 & (n \equiv 1 \pmod 4) \\
0 & (n \equiv 2 \pmod 4) \\
-1 & (n \equiv 3 \pmod 4) \\
\end{array}
\right.

$f(x)=\cos(x)$のときも同様に定義できます。

f^{(n)}(0) = \left\{
\begin{array}{ll}
1 & (n \equiv 0 \pmod 4) \\
0 & (n \equiv 1 \pmod 4) \\
-1 & (n \equiv 2 \pmod 4) \\
0 & (n \equiv 3 \pmod 4) \\
\end{array}
\right.

これをコードにします。

trigonometric_function.py
# -*- coding:utf-8 -*-

def fact(x):
    if x == 0: return 1
    else: return x * fact(x-1)

# マクローリン展開
def maclaurin(f):
    def g(x):
        ans = 0.0
        for k in range(30):
            ans += f(k) * (x**k) / fact(k)
        return ans
    return g

# sin(x)のn次導関数
def diff_sin(n):
    return (0, 1, 0, -1)[n % 4]

# cos(x)のn次導関数
def diff_cos(n):
    return (1, 0, -1, 0)[n % 4]

sin = maclaurin(diff_sin)
cos = maclaurin(diff_cos)

if __name__ == "__main__":
    import math
    print("---sin---")
    print("自作", sin(0.1))
    print("math", math.sin(0.1))
    print("差分", abs(sin(0.1) - math.sin(0.1)))

    print("---cos---")
    print("自作", cos(0.1))
    print("math", math.cos(0.1))
    print("差分", abs(cos(0.1) - math.cos(0.1)))

"""
$ python3 trigonometric_function.py 
---sin---
自作 0.09983341664682817
math 0.09983341664682815
差分 1.3877787807814457e-17
---cos---
自作 0.9950041652780258
math 0.9950041652780258
差分 0.0
"""

このmaclaurin関数を用いれば、n次導関数が分かっている関数を定義できます。例えば、$f(x)=e^x$のとき、任意のnについて$f^{(n)}(0)=1$なので

f = maclaurin(lambda _: 1)
print(f(1)) # 2.7182818284590455

と定義できます。

参考文献・サイト

マクローリン展開|高校数学の美しい物語

2
1
0

Register as a new user and use Qiita more conveniently

  1. You get articles that match your needs
  2. You can efficiently read back useful information
  3. You can use dark theme
What you can do with signing up
2
1

Delete article

Deleted articles cannot be recovered.

Draft of this article would be also deleted.

Are you sure you want to delete this article?