1. 過剰決定系
最初の例として、1 つの従属変数 y(x) を持つ 2 種類の不等式の系を過剰決定します。
>
|
![sys1 := [x*(diff(y(x), x))^2*(diff(y(x), x, x))^2-2*x*(diff(y(x), x))*(diff(y(x), x, x))*y(x)*(diff(y(x), x, x, x))+x*y(x)^2*(diff(y(x), x, x, x))^2-y(x)*(diff(y(x), x, x))+(diff(y(x), x))^2 = 0, -(diff(y(x), x))*(diff(y(x), x, x))+y(x)*(diff(y(x), x, x, x))+2*y(x)^2*(diff(y(x), x, x))^2-4*y(x)*(diff(y(x), x, x))*(diff(y(x), x))^2+2*(diff(y(x), x))^4 = 0]](/support/helpjp/helpview.aspx?si=1448/file00448/math112.png)
|
![sys1 := [x*(diff(y(x), x))^2*(diff(y(x), `$`(x, 2)))^2-2*x*(diff(y(x), x))*(diff(y(x), `$`(x, 2)))*y(x)*(diff(y(x), `$`(x, 3)))+x*y(x)^2*(diff(y(x), `$`(x, 3)))^2-y(x)*(diff(y(x), `$`(x, 2)))+(diff(y(x), x))^2 = 0, -(diff(y(x), x))*(diff(y(x), `$`(x, 2)))+y(x)*(diff(y(x), `$`(x, 3)))+2*y(x)^2*(diff(y(x), `$`(x, 2)))^2-4*y(x)*(diff(y(x), `$`(x, 2)))*(diff(y(x), x))^2+2*(diff(y(x), x))^4 = 0]](/support/helpjp/helpview.aspx?si=1448/file00448/math115.png)
| (4.1) |
ここで次のように rifsimp を呼び出します。
>
|
|
>
|
|
| (4.2) |
ここで系がより単純な方程式に簡素化され (y(x) の場合は完全にゼロになるわけではありません)、Maple の dsolve で処理できるようになります。
>
|
|
| (4.3) |
厳密解を得るために dsolve を支援するだけでなく、この簡素化された形式は、initialdata 関数と組み合わせて使用することで、式を解くのに必要な初期データを得ることができます。
>
|
|
| (4.4) |
この場合は、2 階 ODE で想定された結果が得られましたが、さらに複雑な PDE 系のための必須初期データを計算することもできます (詳細は initialdata を参照してください)。これで、initialdata で計算した型の初期条件を持つ縮小された系に対し、数値計算方法を正しく適用できるようになりました。
>
|
|
| (4.5) |
最終的には、この問題の局所解の形式べき級数 (Taylor 級数) を得ることもできます (詳細については rtaylor を参照してください)。
>
|
|
![tay_ser := [y(x) = y(x[0])+(D(y))(x[0])*(x-x[0])+(1/2)*(D(y))(x[0])^2*(x-x[0])^2/y(x[0])+(1/6)*(D(y))(x[0])^3*(x-x[0])^3/y(x[0])^2+(1/24)*(D(y))(x[0])^4*(x-x[0])^4/y(x[0])^3]](/support/helpjp/helpview.aspx?si=1448/file00448/math183.png)
| (4.6) |
前述の初期データを与えると、上記のように級数が完全に決定されます:
>
|
|
![tay_ser := [y(x) = _C1+_C2*(x-x[0])+(1/2)*_C2^2*(x-x[0])^2/_C1+(1/6)*_C2^3*(x-x[0])^3/_C1^2+(1/24)*_C2^4*(x-x[0])^4/_C1^3]](/support/helpjp/helpview.aspx?si=1448/file00448/math192.png)
| (4.7) |
2. 矛盾する系
rifsimp は、系が矛盾していないか判別するのによく使用されます。
>
|
|
>
|
|
| (4.8) |
3. 拘束機械系
この例では、rifsimp を拘束機械系 (代数方程式、つまり DAE 系) のプリプロセッサとして使用します。Lagrange 法の公式を使用すると、以下の重力下における摩擦なし形状 Phi(x,y) のワイヤのビーズ運動の非ゼロ質量 m を求めることができます:
>
|
)(x(t), y(t)), m*(diff(y(t), t, t)) = -lambda(t)*(D[2](Phi))(x(t), y(t))-m*g]](/support/helpjp/helpview.aspx?si=1448/file00448/math224.png)
|
)(x(t), y(t)), m*(diff(y(t), `$`(t, 2))) = -lambda(t)*(D[2](Phi))(x(t), y(t))-m*g]](/support/helpjp/helpview.aspx?si=1448/file00448/math227.png)
| (4.9) |
空気抵抗を考慮せずに、重力下で降下する質量は次のように表されます:
>
|
|
| (4.10) |
ここでは振子を例に挙げ、ビーズは円を描くようにワイヤ上を移動します。
>
|
|
| (4.11) |
>
|
|
| (4.12) |
このような拘束系を数値ソルバで表すことは非常に困難です。もちろん、極座標を使用した拘束を排除することはできますが、通常、複雑な拘束を受ける機械系にとってこれは不可能です。たとえば、拘束を別の関数に置き換えることでワイヤの形状を変えてみます。
>
|
|
| (4.13) |
(ページ Rif の参考文献およびパッケージ情報 (Reidら著: 1996 年) を参照してください)。振子はこのような系における代表的な例で、アプリケーションでの重要性が高いことから近年の中心的な研究対象です。
>
|
|
![SimpSys := TABLE([Case = [[x(t) <> 0, diff(x(t), t)], [(y(t)-1)*(y(t)+1) <> 0, diff(lambda(t), t)]], Solved = [diff(y(t), `$`(t, 2)) = -(2*lambda(t)*y(t)+m*g)/m, diff(lambda(t), t) = -(3/2)*(diff(y(t), t))*m*g, diff(x(t), t) = -y(t)*(diff(y(t), t))/x(t)], Pivots = [m <> 0, x(t) <> 0, y(t)-1 <> 0, y(t)+1 <> 0], Constraint = [m*(diff(y(t), t))^2+2*y(t)^2*lambda(t)-y(t)*m*g+y(t)^3*m*g-2*lambda(t) = 0, x(t)^2+y(t)^2-1 = 0]])](/support/helpjp/helpview.aspx?si=1448/file00448/math276.png)
| (4.14) |
>
|
|
![TABLE([Finite = [lambda(t[0]) = _C1, x(t[0]) = _C2, y(t[0]) = _C3, (D(y))(t[0]) = _C4, g = _C5, m = _C6], Pivots = [_C6 <> 0, _C2 <> 0, _C3-1 <> 0, _C3+1 <> 0], Infinite = [], Constraint = [_C6*_C4^2+2*_C3^2*_C1-_C3*_C6*_C5+_C3^3*_C6*_C5-2*_C1 = 0, _C2^2+_C3^2-1 = 0]])](/support/helpjp/helpview.aspx?si=1448/file00448/math283.png)
| (4.15) |
前述した例ほどの簡潔さはないものの、rifsimp によって初期条件で満たす必要があるさらに別の拘束 (出力 initialdata における Constraint および Pivots) が見つかり、時間導関数の方程式 lambda(t) を得ることができます。ここで Maple の dsolve[numeric] を使用すると前述した例のような数値解が得られます。
4. Lie 対称性決定
この例では、rifsimp で ODE の Lie 点の対称性を決定する方法を示します:
>
|
|
| (4.16) |
以下のようにして、Lie 点の対称性を求める PDE 系を得ます:
>
|
|

| (4.17) |
rifsimp を使用してこの系を次のように大幅に簡約します:
>
|
|
![rifsys := TABLE([Solved = [diff(eta(x, y), `$`(x, 2)) = (-eta(x, y)*x+xi(x, y)*y)/x^3, diff(xi(x, y), x) = eta(x, y)/y, diff(eta(x, y), y) = eta(x, y)/y, diff(xi(x, y), y) = 0]])](/support/helpjp/helpview.aspx?si=1448/file00448/math332.png)
| (4.18) |
この時点で pdsolve を使用すると以下のように Lie 点の対称性を得ることができます:
>
|
|

| (4.19) |