maximaの初歩的な操作19･･･定数係数２階常微分方程式

== maximaの初歩的な操作19･･･定数係数２階常微分方程式 == ○Xmaxima，wxMaxima のインストール方法，基本操作については［この頁］参照
○wxmaximaを起動すると，初めは入力行も何もないので，とりあえず「何かキーを押します：例えばEnterキーを押す」と入力行が現れます．
○wxmaximaでは，入力行の末尾でShift+Enterを押します．

○wxMaximaで１階または２階の常微分方程式の一般解を求めるには，

(1)　ode2( )関数を使う方法
…ordinary differential equation of 1st and 2nd order
(2)　ラプラス変換を使ってdesolve( )関数で解く方法

があります．この頁では(1)のode2( )関数を使う方法を扱います．

１階微分方程式では積分定数は%cで表され，
２階微分方程式では積分定数は%k1 , %k2で表されます．

○wxmaximaで，

１階微分方程式で初期条件に対応する特殊解を求めるにはic1( )関数が，
…initial value problems for 1st order differential equations
２階微分方程式で初期条件に対応する特殊解を求めるにはic2( )関数が
…initial value problems for 2nd order differential equations
使えます．

※微分方程式の厳密解は求め方が知られていない場合や，解けても特殊関数の複雑な組合せとなる場合があります．
wxmaximaのヘルプ情報（英文）によれば，
１階常微分方程式については，

線形，変数分離形，同次形，完全微分形，ベルヌーイ形などが試され，解けなければfalseが返される

とされています．
２階常微分方程式については，

定数係数，定数係数に変換できる同次線形，オイラーの微分方程式，定数変化法によって解けるもの，２つの１階微分方程式に分けられるものなどが試されます．解けなければfalseが返される

　級数を用いた解については，このような制限はあまりないがこの頁では扱いません．

≪頭の片隅に！≫
　厳密解･･･有限回の代数的変形や積分操作によって得られる解･･･は，コンピュータでも求められない場合がある．（変数変換のアイデア次第で解ける問題は，人の方が強い場合がある）

　解き方の分かっている問題に当てめる場合や，理論上無限回（実際上は巨大な回数）の操作を用いる級数解や数値解を求める場合コンピュータの方が強い．

■定数係数２階常微分方程式

≪例１≫ …［定数係数］特性方程式が異なる２つの実数解をもつ場合
【微分方程式】
$\frac{d^2y}{dx^2}-3\frac{dy}{dx}+ 2y=2x+ 3$
(1) $x=0$ のとき $y=1,\frac{dy}{dx}=2$
(2) $x=0$ のとき $y=1$ ， $x=1$ のとき $y=3$
【wxMaximaでの入力】
siki1 : ode2( 'diff(y,x,2)-3*'diff(y,x)+2*y=2*x+3 , y , x);
(1)ic2(siki1 , x=0 , y=1 , diff(y,x)=3);
(2)bc2(siki1 , x=0 , y=1 , x=1 , y=2);

【出力される結果】
$y=k_1e^{2x}+ k_2e^x+ x + 3$ …(*)
(1) $y=4e^{2x}-6e^x+ x + 3$
(2) $y=2e^{2x-1}-(2e+ 2)e^x-1+ x+ 3$

【補足説明】
○関数ode2( )の引数は
ode2( 微分方程式 , 従属変数 , 独立変数);
の順に指定します．

メニューから対話型メニューで行うには，「方程式」→「微分方程式を解く」と進んで，上記の引数３個を書き込みます．

○その場合において，注意すべきこととして導関数を表す記号diff(y , x)にはシングルクォートを付けて'diff(y , x)と書きます．（単にdiff(y , x)と書くと，この文脈の中では文字yを文字xで（偏）微分した0に変わります．ここでは，この式自体を生かしたいので，名詞化処理処理としてシングルクォートを付けるようです．）
　次のように，xの関数であることを明示して使う場合には，シングルクォートを付けても付けなくても同じ結果が得られます．
ode2( diff(f(x), x , 2) ... diff(f(x) , x) ... );

○微分方程式の部分で等号を省略して単に式を指定した場合は，=0が省略されたものとみなされます．（右辺＝0）

○diff( )の第３引数の1を省略すれば，第１次導関数を表しますが，第２次導関数を表す場合はdiff(y , x , 2 )のように第３引数を2にします．

○２階微分方程式の一般解は２個の積分定数を含みますが，Maximaの出力においてこれらは%k1, %k2で表されます．　入力においても，自然対数の底は%e，円周率は%piという形の記号定数を使います．

○２階微分方程式の一般解(*)に対して，条件を追加して関数を確定するには，次の方法があります．

(1)　【初期条件を指定する方法】

x=aのときy=pとなりy’=qとなるという形で指定する．すなわちf(a)とf ’(a)を指定する方法

(2)　【境界条件を指定する方法】

x=aのときy=pとなり，x=bのときy=qとなるという形で指定する．すなわちf(a)とf(b)を指定する方法

○Maximaにおいて上記の(1)の初期条件を満たす特殊解を求めるには，
ic2( 一般解 , x=a , y=p , diff(y,x)=q);
の形で指定します．（この関数ic2( )の命名についてHELP情報に明示されている訳ではないが，initial conditions 2と考えれば覚えやすい）
　このように，先に一般解を求めてから，続いてある初期条件を満たす特殊解を求めるには，一般解の式に名前を付けておいて後でそれを引用するとよいので，上記の一般解を求めるコマンドを
siki1 : ode2('diff(y,x,2)-3*'diff(y,x)+2*y=2*x+3 , y , x);
などとしておいてから
ic2( siki1 , x=0 , y=1 , diff(y,x)=3);
とすればよいことになります．

このとき
　siki1には名詞化処理（シングルクォート）はしない．
　第２，３，４引数は各々単なる値ではなく，x=... , y= ... , diff(y,x)= ...のように（どの変数）＝（値が幾ら）であると書かなければならない．
　第４引数は単にdiff(y,x)= ...と書けばよく，名詞化処理'diff(y,x)= ...はあってもなくてもよい．

　Maximaでは単に%と書けば直前の出力式を表すことになっているので，一般解を求めた直後に初期条件を加える場合には，式にな名前がなくても
ic2( % , x=0 , y=1 , diff(y,x)=3);
でできます．

メニューから対話型メニューで行うには，「方程式」→「特殊解を求める（2階微分方程式）」と進んで，上記の引数４個を書き込みます．

○Maximaにおいて上記の(2)の境界条件を満たす特殊解を求めるには，
bc2( 一般解 , x=a , y=p , x=b , y=q);
の形で指定します．（この関数bc2( )の命名についてHELP情報に明示されている訳ではないが，boundary conditions 2と考えれば覚えやすい）

メニューから対話型メニューで行うには，「方程式」→「境界条件から特殊解を求める」と進んで，上記の引数５個を書き込みます．

【補足説明２】･･･筆算で確かめるには
○微分方程式
$\frac{d^2y}{dx^2}+ a\frac{dy}{dx}+ by=P(x)$ …(*1)
の一般解は，同次方程式
$\frac{d^2y}{dx^2}+ a\frac{dy}{dx}+ by=0$ …(*2)
の一般解と非同次方程式
$\frac{d^2y}{dx^2}+ a\frac{dy}{dx}+ by=P(x)$ …(*3)
の１つの特殊解の和で表すことができます．

(*2)は解が $y=e^{rx}$ の形をしているものを探すとすると
$\frac{d^2y}{dx^2}=r^2e^{rx},\frac{dy}{dx}=re^{rx},y=e^{rx}$
より
$e^{rx}(r^2+ ar + b)=0$
ここで $e^{rx}$ $\neq$ $0$ だから特性方程式 $r^2+ ar + b=0$ が異なる２つの実数解 $r=\alpha,\beta$ をもつとき，各々
$y_1=e^{\alpha x},y_2=e^{\beta x}$
となりますが，(*2)の一般解は
$y=k_1 y_1+ k_2 y_2=k_1e^{\alpha x}+ k_2e^{\beta x}$
となります．
　もとの問題では， $r^2- 3r + 2=0$ の解は $r=2,1$ だからの一般解は
$y=k_1e^{2x}+ k_2e^{x}$
となります．

(*3)の１つの特殊解は，与えられた関数から連想すればよく， $P(x)$ がｎ次式の場合 $y''$ はｎ-２次， $y'$ はｎ-１次， $y$ はｎ次なので， $P(x)$ がｎ次式と想定して係数を比較すると求められることが多い．
もとの問題では， $P(x)=2x+ 3$ なので，１つの特殊解を
$y_0=Ax+ B$
と予想して， $A,B$ を求めてみると
$y_0''=0,y_0'=A,y_0=Ax+ B$
により
$0-3A+2(Ax+ B)=2x+ 3$
より
$A=1,B=3$
１つの特殊解は
$y_0=x+ 3$

結局，もとの方程式(*1)の一般解は
$y=k_1e^{2x}+ k_2e^{x}+ x+ 3$
となります．

≪例２≫ …［定数係数］特性方程式が実数の重解をもつ場合
【微分方程式】
$\frac{d^2y}{dx^2}-4\frac{dy}{dx}+ 4y=\sin 2x$
【wxMaximaでの入力】
ode2( 'diff(y,x,2)-4*'diff(y,x)+4*y=sin(2*x) , y , x);

【出力される結果】
$y=\frac{\cos 2x}{8}+ (k_2 x+ k_1)e^{2x}$

【補足説明】･･･筆算で確かめるには
　特性方程式 $r^2-4x+4=0$ が重解 $r=2$ をもつ．
　このような場合には， $y_1=e^{2x}$ と $y_2=xe^{2x}$ を用いて，同次方程式の一般解
$y=(k_2 x+ k_1)e^{2x}$
が得られます．

　１つの特殊解は，右辺の形から連想して（微分すると $\sin 2x, \cos 2x$ が入れ替わるので，両方とも準備する）
$y_0=A\sin 2x + B\cos 2x$
とおいて係数比較するとよい．
$y'=2A\cos 2x-2B\sin 2x$
$y_0''=-4A\sin 2x -4B\cos 2x$
だから
$(-4A+ 8B J)\sin 2x +(-4b-8A+ 4B)\cos 2x=\sin 2x$
より
$A=0,B=\frac{1}{8}$
したがって，１つの特殊解は
$y_0=\frac{1}{8}\cos 2x$

結局，もとの方程式の一般解は
$y=(k_2 x+ k_1)e^{2x}+ \frac{1}{8}\cos 2x$
になります

≪例３≫ …［定数係数］特性方程式が虚数解をもつ場合
【微分方程式】
$\frac{d^2y}{dx^2}-2\frac{dy}{dx}+ 2y=e^x$
【wxMaximaでの入力】
ode2( 'diff(y,x,2)-2*'diff(y,x)+2*y=%e^x , y , x);

【出力される結果】
$y=e^x(k_1\sin x+ k_2 \cos x)+ e^x$

【補足説明】･･･筆算で確かめるには
　特性方程式 $r^2-2x+ 2=0$ が虚数解 $r=1\pm i$ をもつ．
　このような場合には， $y_1=e^{(1+ i)x}$ と $y_2=e^{(1-i)x}$ を用いて，同次方程式の一般解
$y=k_1 e^{(1-i)x}+ k_2e^{(1+ i)x}$
が得られます．
　複素関数まで習っている場合はこのままでも平気ですが，もとが実変数の微分方程式なのに複素関数で答えるのは気味が悪いと思うときは，オイラーの公式
$e^{ix}=\cos x + i\sin x$
を用いて，実数の関数に直すことができます．
$e^{x+ iy}=e^xe^{iy}=e^x(\cos y + i\sin y)$
だから
$y_1=e^{(1+ i)x}=e^x(\cos x + i\sin x)$
$y_2=e^{(1- i)x}=e^x(\cos x - i\sin x)$

ここで新たに
$y_3=\frac{y_1+ y_2}{2}=e^x\cos x$
$y_4=\frac{y_1- y_2}{2i}=e^x\sin x$
を２つの独立な解とすれば
$y=k_2y_3+ k_1y_4=k_2 e^x\cos x+ k_1 e^x \sin x$
$=e^x(k_2 \cos x+ k_1 \sin x)$
が同次方程式の一般解になります．

　１つの特殊解は，右辺の形から連想して
$y_0=Ae^x$
とおいて係数比較するとよい．
$y'=Ae^x,y''=Ae^x$
だから
$Ae^x-2Ae^x+ 2Ae^x=Ae^x$
より
$A$ は任意の数 $\rightarrow A=1$ とおく
１つの特殊解は
$y_0=e^x$

結局，もとの方程式の一般解は
$y=e^x(k_2 \cos x+ k_1 \sin x)+ e^x$
になります

...メニューに戻る