Python で微分積分などの記号計算をする場合は,以下のように必要なパッケージを import します。(弘大 JupyterHub では必要なパッケージはインストール済みです。)
セルをクリックして,上のメニューから「▶︎ Run」をクリックするか,Shift キーを押しながら Enter キー(Return キー)を押すと実行します。
from sympy.abc import *
from sympy import *
# SymPy Plotting Backends (SPB): グラフを描く際に利用
from spb import *
# グラフを SVG で Notebook にインライン表示
%config InlineBackend.figure_formats = ['svg']
導関数の定義
$f(x) = x^2$ を微分する。
Python では関数の定義は以下のように書きます。
# Python での関数 f(x) の定義
def f(x):
return x**2
まずは以下の導関数の定義に基づいて極限を計算する。
$$ \frac{df}{dx} \equiv \lim_{h \rightarrow 0} \frac{f(x+h) – f(x)}{h}$$
limit((f(x+h)-f(x))/h, h, 0)
ともかく答えを求めたい場合はこれで OK です。極限の計算をする前の表示をしたい場合は Limit()
を使います。計算を実行させるには .doit()
をつけます。
Eq(Limit((f(x+h)-f(x))/h, h, 0),
Limit((f(x+h)-f(x))/h, h, 0).doit())
微分を計算するのにいちいち極限をとっているわけにはいきません。次に微分を行う関数 diff()
を使って直接 f(x)
を微分してみます。
diff(f(x), x)
以下のように微分を計算する前の式を表示させたい場合は,Derivative()
を使います。
Eq(Derivative(f(x), x),
Derivative(f(x), x).doit())
微分法の公式
微分法の公式 にある公式を確認します。
$\displaystyle 1.\ \left\{ f(x) \pm g(x) \right\}’ = f'(x) \pm g'(x)$
# from sympy.abc import * していれば不要ですが
# 念のため
x = Symbol('x')
# 「未定義」な関数の定義
f = Function('f')
g = Function('g')
# これらの関数を使うときの呼び出し方
f(x)
# 「未定義」関数はこのようにも定義できる
F = Function('F')(x)
# これを呼び出すときには...
F
Eq(Derivative(f(x) + g(x), x),
Derivative(f(x) + g(x), x).doit())
$\displaystyle 2. \ \left\{ c f(x) \right\}’ = c f'(x) $
Eq(Derivative(c*f(x), x),
Derivative(c*f(x), x).doit())
上の例では左辺に $\displaystyle \frac{\partial}{\partial x}$ という微分記号が見えていますが,これは偏微分記号です。別の授業で習います。
$\displaystyle 3. \ \left\{ f(x)\, g(x) \right\}’ = f'(x)\,g(x) + f(x)\,g'(x) $
ライプニッツルールを知っているか?
Eq(Derivative(f(x)*g(x), x),
Derivative(f(x)*g(x), x).doit())
$\displaystyle 4. \ \left\{ \frac{1}{g(x)} \right\}’ = – \frac{g'(x)}{\left\{g(x)\right\}^2} $
Eq(Derivative(1/g(x), x),
Derivative(1/g(x), x).doit())
$\displaystyle 5. \ \left\{ \frac{f(x)}{g(x)} \right\}’ = \frac{f'(x)\,g(x) – f(x)\, g'(x)}{\left\{g(x)\right\}^2} $
Eq(Derivative(f(x)/g(x), x),
Derivative(f(x)/g(x), x).doit())
simplify()
は通分などを行なって「簡単化」してくれます。
simplify(diff(f(x)/g(x), x))
合成関数の微分法
$y = y(u), \ u = u(x)$,つまり $y$ が $u$ を通して $x$ の関数であるとき,
$$\frac{dy}{dx} = \frac{dy}{du} \frac{du}{dx}$$
ということを理解しているか確認。
y = Function('y')
u = Function('u')
Eq(Derivative(y(u(x)), x),
Derivative(y(u(x)), x).doit())
例: $ y = (a x^2 + b x + c)^5 $ の微分。
u = a*x**2 + b*x + c
y = u**5
Eq(Derivative(y, x),
Derivative(y, x).doit())
べき関数の微分
べき乗(累乗)は **
で表します。$x^p$ は,x**p
。べき関数 x**p
を微分してみます。
_
は直前の計算結果をあらわします。
diff(x**p, x)
simplify(_)
例として,べき関数の合成関数である,以下のような関数の微分。
$$\frac{d}{dx} \frac{1}{\sqrt{1 + x^2}}$$
$\sqrt{x} =$ sqrt(x)
と書きます。
f = 1/sqrt(1 + x**2)
Eq(Derivative(f, x),
Derivative(f, x).doit())
参考:べき関数のグラフ
$y = x^2$ を $ -3 \leq x 3$ の範囲でグラフを描く例。2次元グラフは plot()
関数で描きます。弘大 JupyterHub では SymPy Plotting Backends もインストールされているので以下のように import していることを前提に説明します。
# SymPy Plotting Backends (SPB): グラフを描く際に利用
from spb import *
plot(x**2, (x, -3, 3));
$y = x, x^2, x^3$ を $ -3 \leq x 3$ の範囲でグラフを描く例。
- $y$ 軸(縦軸)の表示範囲を設定するときは
ylim = (ymin, ymax)
などのようにします。
plot(x, x**2, x**3, (x, -3, 3),
ylim = (-3, 7));
$y = x, x^{1/2}, x^{1/3}$ のグラフを描く例です。
(細かなことですが,Python では 1/3 = 0.3333333333333333
なので分数のままにするために Rational(1, 3)
などとしています。)
plot(x, x**(Rational(1, 2)), x**(Rational(1,3)), (x, 0, 3));
$y = x^{-2}, \ x^{-3}$ のグラフを描く例。
plot(1/x**2, 1/x**3, (x, -3, 3), ylim = (-5, 5));
$y = x^{-3}$ のグラフが変です。$y \rightarrow -\infty$ と $y \rightarrow +\infty$ をつなぐ縦の線が見えています。不連続点の処理を適切に行うために,detect_poles = True
オプションをつけてみます。
plot(1/x**2, 1/x**3, (x, -3, 3), ylim = (-5, 5),
detect_poles = True);
指数関数の微分
念のため,SymPy が $\displaystyle \lim_{h \rightarrow 0}(1 + h)^{\frac{1}{h}}$ がわかるか確認します。
Eq(Limit((1 + h)**(1/h), h, 0),
Limit((1 + h)**(1/h), h, 0).doit())
ネイピア数 $e$ を底とした指数関数 $e^x$ は exp(x)
で表します。
$$\frac{d}{dx} e^x = \cdots$$
Eq(Derivative(exp(x), x),
Derivative(exp(x), x).doit())
微分した答えだけを知りたいなら以下のようにするだけ。
diff(exp(x), x)
参考:SymPy の基本的な定数として,以下のようなものがあります。
円周率 $\pi$ | ネイピア数 $e$ | 虚数単位 $i$ | 無限大 $\infty$ |
---|---|---|---|
pi | E | I | oo |
$\exp(1) = e^1 =$ E
であることを確認します。
exp(1) - E
float()
関数でネイピア数の値を小数点表示してみます。
float(E)
参考:指数関数のグラフ
$y = e^x, e^{-x}$ のグラフ。
plot(exp(x), exp(-x), (x, -3, 3));
対数関数
定義
自然対数は log(x)
です。自然対数は,ネイピア数 $e$ を底とした指数関数 $e^x$ の逆関数であるから,
$$e^{\log x} = x$$
が成り立ちます。確認してみます。
exp(log(x))
微分
$$\frac{d}{dx} \log x = \cdots$$
diff(log(x), x)
参考:対数法則
- $\log(x y) = \log x + \log y$
- $\log(x^y) = y \log x$
これらを確認するには expand_log()
関数に force=True
オプションが必要。
# x, y = symbols('x y')
var('x y')
expand_log(log(x*y), force=True)
expand_log(log(x**y), force=True)
さらに,$\log x$ は $e^x$ の逆関数であることから
$$\log e^x = x$$
も成り立ちます。確認するには,force=True
オプションをつけて…
expand_log(log(exp(x)), force=True)
参考:対数関数のグラフ
$y = \log x$ のグラフを $0 < x < 5$ の範囲で描く。$y$ の表示範囲を $ -3 < y < 3$ くらいにして。
plot(log(x), (x, 0.01, 5), ylim = (-3, 3));
三角関数
性質
三角関数は sin(x), cos(x), tan(x)
のように書きます。
$\sin^2 x + \cos^2 x = 1$ を確認します。
simplify(sin(x)**2 + cos(x)**2)
極限公式の確認
三角関数の微分を求めるためには,以下の極限
$$ \lim_{x \rightarrow 0} \frac{\sin x}{x} = \cdots$$
の値が必要でした。この極限ができるか確認します。
Eq(Limit(sin(x)/x, x, 0),
limit(sin(x)/x, x, 0))
$\displaystyle \lim_{x \rightarrow 0} \frac{1 – \cos x}{x} = 0$ であることも確認します。
Eq(Limit((1-cos(x))/x, x, 0),
limit((1-cos(x))/x, x, 0))
加法定理の確認
$\sin (x + h) = \sin x \, \cos h + \cos x \, \sin h$ および
$\cos (x + h) = \cos x \, \cos h – \sin x \, \sin h$ の確認。
三角関数を「展開」するときは,expand()
関数に trig = True
オプションをつけます。
expand(sin(x + h), trig = True)
expand(cos(x + h), trig = True)
微分
まずは導関数の定義に従って $\sin x$ の微分をしてみます。
$\displaystyle \frac{d}{dx} \sin x = \lim_{h \rightarrow 0} \frac{\sin(x+h) – \sin (x)}{h} = \cdots$
Eq(Limit((sin(x+h)-sin(x))/h, h, 0),
limit((sin(x+h)-sin(x))/h, h, 0))
diff()
関数で三角関数を微分しても同じ答えが出ることを確認します。
Eq(Derivative(sin(x), x),
diff(sin(x), x))
$\cos x, \tan x$ の微分についても diff()
関数で計算します。
Eq(Derivative(cos(x), x),
diff(cos(x), x))
Eq(Derivative(tan(x), x),
diff(tan(x), x))
右辺はこれでもよいですが,simplify()
関数で「簡単化」してみます。
simplify(diff(tan(x), x))
参考:三角関数のグラフ
$y = \sin x, \cos x, \tan x$ のグラフを描きます。$ -2 \pi \leq x \leq 2 \pi$ の範囲で。$\tan x$ は非常に大きい値になる場合があるので,$y$ の表示範囲を適宜設定します。
また,$y = \tan x $ の不連続点を適切に処理するために detect_poles = True
オプションをつけておきます。
plot(sin(x), cos(x), tan(x), (x, -2*pi, 2*pi), ylim = (-3, 3),
detect_poles = True);
参考:横軸の目盛を $\pi/2$ ごとに
$x$ 軸(横軸)の目盛を $\pi/2$ ごとにしてみる。また SPB/Matplotlib のデフォルトではグリッド線が実線でちと目障りだったり,$x$ 軸,$y$ 軸が表示されなかったりなので,その辺も設定してみる。
p = plot(sin(x), cos(x), tan(x), (x, -2*pi, 2*pi),
xlim = (-2*pi, 2*pi), ylim = (-5, 5),
detect_poles = True, size = (6.28, 6.4),
show = False)
ax = p.ax
# x 軸の目盛は 0.5π ごとに
Pi = float(pi)
ax.set_xticks(
[(0.5*i*Pi) for i in range(-4, 5)])
# y 軸の目盛は 1 ごとに
ax.set_yticks(
[i for i in range(-5, 6)])
# x 軸の目盛のラベル。既約分数表記で。
ax.set_xticklabels(
[(str(Rational(i,2))+" $\pi$") for i in range(-4, 5)])
# グリッドを灰色点線に
ax.grid(True, which="major", color="lightgray",
dashes=(3, 5), linewidth=0.5)
# x軸 y軸。表示オプションはないのでこれで代用。
ax.axhline(0, color='black', dashes=(3, 5), linewidth=0.5)
ax.axvline(0, color='black', dashes=(3, 5), linewidth=0.5);
逆三角関数
定義
逆三角関数は三角関数の逆関数として定義され,それぞれ
$\sin^{-1} x = \arcsin x = $ asin(x)
$\cos^{-1} x = \arccos x = $ acos(x)
$\tan^{-1} x = \arctan x = $ atan(x)
と書きます。逆関数ですから以下の関係が成り立つはずです。
$\sin\left(\sin^{-1} x\right) = x$
$\cos\left(\cos^{-1} x\right) = x$
$\tan\left(\tan^{-1} x\right) = x$
確認してみます。
sin(asin(x)), cos(acos(x)), tan(atan(x))
しかし,$\sin^{-1} \left(\sin x\right) = x$ などとは一般にはなりません。
反例は $x = \pi$ とすれば…
asin(sin(pi))
微分
Eq(Derivative(asin(x), x),
diff(asin(x), x))
Eq(Derivative(acos(x), x),
diff(acos(x), x))
Eq(Derivative(atan(x), x),
diff(atan(x), x))
参考:逆三角関数のグラフ
$\sin^{-1} x$ と $\cos^{-1} x$ の定義域は $ -1 < x < 1$,$\tan^{-1} x$ の定義機は $-\infty < x < \infty$ であることに留意して,逆三角関数のグラフを描く。
$y$ 軸(縦軸)の目盛を $\pi/4$ ごとにしてみる。
p = plot((asin(x), "arcsin $x$"), (acos(x), "arccos $x$"),
(x, -1, 1), xlim = (-1, 1), show = False)
ax = p.ax
# y 軸の目盛
Pi = float(pi)
ax.set_yticks(
[(0.25*i*Pi) for i in range(-2, 5)])
# y 軸のラベル
ax.set_yticklabels(
["%.2f $\pi$" % (0.25*i) for i in range(-2, 5)]);
p = plot((atan(x), "arctan $x$"), (x, -10, 10), legend = True,
xlim = (-10, 10), ylim = (-pi/2, pi/2), show = False)
ax = p.ax
# y 軸の目盛
Pi = float(pi)
ax.set_yticks(
[(0.25*i*Pi) for i in range(-2, 3)])
# y 軸のラベル
ax.set_yticklabels(
["%.2f $\pi$" % (0.25*i) for i in range(-2, 3)]);
双曲線関数
定義
双曲線関数は指数関数を使って以下のように定義されます。
$\displaystyle \sinh x \equiv \frac{e^{x} – e^{-x}}{2} =$ sinh(x)
$\displaystyle \cosh x \equiv \frac{e^{x} + e^{-x}}{2} =$ cosh(x)
$\displaystyle \tanh x \equiv \frac{\sinh x}{\cosh x} = \frac{e^{x} – e^{-x}}{e^{x} + e^{-x}} =$ tanh(x)
双曲線関数を指数関数で書き直すには rewrite()
を使います。
sinh(x).rewrite(exp)
cosh(x).rewrite(exp)
tanh(x).rewrite(exp)
性質
$\cosh^2 x – \sinh^2 x = 1$ を確認してみます。
simplify(cosh(x)**2 - sinh(x)**2)
加法定理の確認
双曲線関数にも,三角関数の場合と似たような加法定理があります。
Eq(expand(cosh(x + y)),
expand(cosh(x + y), trig = True))
Eq(expand(sinh(x + y)),
expand(sinh(x + y), trig = True))
Eq(expand(tanh(x + y)),
expand(tanh(x + y), trig = True))
上式の右辺を通分してみます。
expand(tanh(x + y), trig = True).factor()
微分
Eq(Derivative(sinh(x), x),
diff(sinh(x), x))
Eq(Derivative(cosh(x), x),
diff(cosh(x), x))
Eq(Derivative(tanh(x), x),
diff(tanh(x), x))
上式の右辺を「簡単化」してみます。
diff(tanh(x), x).simplify()
参考:双曲線関数のグラフ
$y = \sinh x, \ \cosh x$ のグラフ。
plot((sinh(x), "sinh $x$"), (cosh(x), "cosh $x$"), (x, -3, 3),
xlim = (-3, 3));
$ y = \tanh x$ のグラフ。
plot((tanh(x), "tanh $x$"), (x, -5, 5), legend = True,
xlim = (-5, 5));
逆双曲線関数
定義
逆双曲線関数は双曲線関数の逆関数として定義され,それぞれ
$\sinh^{-1} x = \mbox{arsinh}\ x = $ asinh(x)
$\cosh^{-1} x = \mbox{arcosh}\ x = $ acosh(x)
$\tanh^{-1} x = \mbox{artanh}\ x = $ atanh(x)
と書きます。(ar であり,arc ではない。)逆関数ですから以下の関係が成り立つはずです。
$\sinh\left(\sinh^{-1} x\right) = x$
$\cosh\left(\cosh^{-1} x\right) = x$
$\tanh\left(\tanh^{-1} x\right) = x$
確認してみます。
sinh(asinh(x)), cosh(acosh(x)), tanh(atanh(x))
双曲線関数は指数関数を使って表すことができましたが,逆双曲線関数は対数関数を使って表すことができます。
asinh(x).rewrite(log)
acosh(x).rewrite(log)
atanh(x).rewrite(log)
微分
Eq(Derivative(asinh(x), x),
Derivative(asinh(x), x).doit())
Eq(Derivative(acosh(x), x),
Derivative(acosh(x), x).doit())
Eq(Derivative(atanh(x), x),
Derivative(atanh(x), x).doit())
参考:逆双曲線関数のグラフ
$y = \sinh^{-1} x$ の定義域は $-\infty < x < \infty$
plot(asinh(x), "sinh$^{-1}\, x$", (x, -10, 10), legend = True);
$y = \cosh^{-1} x$ の定義域は $x \geq 1$
plot(acosh(x), "cosh$^{-1}\, x$", (x, 1, 10), legend = True,
xlim = (1, 10));
$y = \tanh^{-1} x$ の定義域は $-1 < x < 1$
plot(atanh(x), "tanh$^{-1}\, x$", (x, -0.999, 0.999), legend = True,
xlim = (-1, 1));