ベルヌーイ型微分方程式

ベルヌーイ型微分方程式 は，次の形式で表される微分方程式である．

\[ \frac{dy}{dx} + p(x)y = q(x)y^n \]

ここで、\(p(x)\) と \(q(x)\) は \(x\) の関数であり，\(n\)は定数である．見てわかるように，\(y^n\)の項が含まれており，一階非線形微分方程式である．

ベルヌーイ型微分方程式の解法

ベルヌーイ型微分方程式と名前が付いているように，この微分方程式は一階非線形微分方程式の中でも特殊な形式であり，変数変換を行うことで一階線形微分方程式として解を求めることができる．

\(n=0\) の場合

まずは\(n=0\)のときのベルヌーイ型微分方程式を考えよう．代入すると，

\[ \frac{dy}{dx} + p(x)y = q(x) \]

という形式となり，これは一階線形非同次微分方程式である．これは前回学んだ積分因子形等の方法で解くことができる．

\(n=1\) の場合

続いて，\(n=1\)のときを考える．代入すると，

\[\begin{split} \begin{aligned} \frac{dy}{dx} + p(x)y &= q(x)y \\ \frac{dy}{dx} + \left(p(x) - q(x) \right)y &= 0 \end{aligned} \end{split}\]

となり，\(p(x),q(x)\) はともに変数 \(x\) に関する項なので \(y\) でまとめることができ，結果的にこちらも単純な一階線形同次微分方程式となる．このような微分方程式は直接積分法や変数分離法などで解くことができた．

\(n\neq0,1\) の場合

ここで扱いたいのは上記のような形式ではなく，\(n\neq0,1\)のときのベルヌーイ型微分方程式である．前述したように適切に変数変換を行うと線形微分方程式に変換できる．

ここで変数 \(z\) を導入する．

\[ z = y^{-(n-1)}=\frac{1}{y^{n-1}} \]

ここで \(y,z\) は変数\(x\)に関する関数なので\(x\)について微分すると，

\[\begin{split} \begin{aligned} \frac{dz}{dx} &= \frac{d}{dx}y^{1-n} \\ &= \left(1-n \right)y^{-n} \frac{dy}{dx} \end{aligned} \end{split}\]

となる．これを使って元の微分方程式を変数変換の関数 \(z\) と \(\frac{dz}{dx}\) で書き換える．そのために両辺を \(y^n\) で割り，変数変換の結果を代入すると

\[\begin{split} \begin{aligned} \frac{dy}{dx} + p(x)y &= q(x)y^n \\ \frac{1}{y^n}\frac{dy}{dx} + p(x) \frac{1}{y^{n-1}} &= q(x) \\ \frac{1-n}{1-n}y^{-n}\frac{dy}{dx} + p(x) \frac{1}{y^{n-1}} &= q(x) \\ \frac{1}{1-n}\frac{dz}{dx} + p(x) z &= p(x) \\ \frac{dz}{dx} + (1-n) p(x) z &= (1-n) p(x) \\ \end{aligned} \end{split}\]

とシンプルな形式に書き換えることができる．これは見てわかるように，関数 \(z\) に関する線形微分方程式であり，一般解は積分因子法や定数変化法で導出できる．一般解は以下となる．

\[ \begin{align} y &= e^{-\left( 1-n \right)\int p(x) dx} \left \{ \left( 1-n \right) \int q(x)e^{\left( 1-n \right)\int p(x) dx} dx + C \right \} \end{align} \]

ベルヌーイ型微分方程式の具体例

では，次の微分方程式を解いてみよう．

\[ \frac{dy}{dx} - y = - xy^2 \]

解答はクリックで確認できる．

これは \(n=2\) のベルヌーイ型であるので，\(z=y^{-1}\) として変数変換する．

\[\begin{split} \begin{aligned} \frac{dz}{dx} &= \frac{d}{dx}y^{-1} \\ &= - y^{-2} \frac{dy}{dx} \end{aligned} \end{split}\]

これを使って元の微分方程式を変数変換の関数 \(z\) と \(\frac{dz}{dx}\) で書き換える．そのために両辺を\(y^n\)で割り，変数変換の結果を代入すると

\[\begin{split} \begin{aligned} \frac{dy}{dx} - y &= -xy^2 \\ \frac{1}{y^2}\frac{dy}{dx} - \frac{y}{y^{2}} &= -x \\ -\frac{dz}{dx} - \frac{1}{y} &= -x \\ \frac{dz}{dx} + z &= x \\ \end{aligned} \end{split}\]

となる．これは \(p(x)=1,q(x)=x\) の非同次微分方程式であるので，\(z=\frac{1}{y}\)に注意して一般解を導出すると，

\[\begin{split} \begin{align} z &= Ce^{-x} + x - 1 \\ \frac{1}{y} &= Ce^{-x} + x - 1 \\ y &= \frac{1}{Ce^{-x} + x - 1} \end{align} \end{split}\]

となる．

Pythonによる実装

では，上記の微分方程式をsympyで解く．

from sympy import symbols, Eq, Derivative, Function, dsolve

# 変数と未知関数の定義
x = symbols('x')
y = Function('y')(x)

# dy/dxの定義
dy = Derivative(y, x)

# 微分方程式の定義
eq = Eq(dy - y, -x*y**2)

# 微分方程式を解く
y_ = dsolve(eq, y)
y_

\[\displaystyle y{\left(x \right)} = \frac{e^{x}}{C_{1} + x e^{x} - e^{x}}\]

import numpy as np
from sympy import plotting
import matplotlib.pyplot as plt

N = 10

eqs = []
for c in np.linspace(-5, 5, N):
    eqs.append(y_.rhs.subs(symbols('C1'), c))
p = plotting.plot(*eqs, (x, 0, 5), ylim=(0, 20), show=False)

cm = plt.get_cmap('magma', N)
for i in range(cm.N):
    p[i].line_color = cm(i)
p.show()

グラフのプロットは \(e^{-x}\) を含むためある値で漸近していることに注意されたい．