複雑な方程式の解を得るにはどうすればよいですか?
はじめに
Update6/2014
もともと、私は三つの解決策を見つけました。 方程式を合理化すると次数24の多項式が得られ、NSolve
は24個の根を見つけ、元の方程式のゼロではない根を取り除くので、3つあることを知っていたと思 それは方程式が数学の観点から面倒だったことが判明し、私は三つを逃しました。
その他の改造:
-
この場合、係数を
SetPrecision
よりも正確な数値に変換するためにRationalize
を使用する方が自然なようです。 -
同様に、私は
First @ eq
をSubtract @@ eq
に置き換えました。
方程式を合理化する
一つは、上のドキュメントによると、ことを覚えておく必要がありますNSolve
NSolve
主に線形方程式と多項式を扱います。
OPの方程式は、多項式に変換することができる代数方程式であり、NSolve
は簡単に解くことができます。 NSolve
はこれを単独で行うわけではないので、”手で合理化する必要があります。「合理化は無関係な解決策を導入する傾向があるので、NSolve
によって返された答えを確認する必要があります。 eq
のやや厄介な係数は、eq
を合理化する際に、丸め誤差が平方根の一部を必要に応じてキャンセルしないようにするため、難しくします。 これを修正するには、eq
の精度をInfinity
に設定します。
式expr
で平方根を見つけることができます
DeleteDuplicates @ Cases, Infinity]
次に、各平方根について、expr
に平方根を負に置き換えて式で乗算し、単純化することができます。
以下は、eq
のEqual
をSubtract
でApply
(@@
)に置き換えた関数にこの方法を適用した結果です。 方程式を正確に解くことができます。
eqExact = Rationalize;rationalized = Simplify @ Fold &, eqExact, DeleteDuplicates@Cases, Infinity]];rootsRatExact = Solve;rootsRat = NSolve(* 24 roots {{ζ -> -24559.6 + 24559.6 I}, <<22>>, {ζ -> 24559.6 - 24559.7 I}}*)
二つの解の間にはいくつかの大きな違いがあります:
(ζ /. rootsRat) - (ζ /. rootsRatExact) // Abs // Max(* 3.15226*10^-10*)
与えられた方程式の根を選択する
更新:ここで私はいくつかのゼロを逃した場所です。
私は、元の方程式eq
が小さい値を持つ根を選択するためにPick
を使用しました。10^-1
と言う小さなしきい値よりも小さい。 それぞれが根であることを後で確認することができますが、他の根は確認しませんでした。 元のルーツは次のとおりでした:
rootsEq = Pick < 10^-1 /. rootsRat](* {{ζ -> -3.78042 - 5.50655 I}, {ζ -> -3.20562 + 5.39914 I}, {ζ -> 6.98478 + 0.493405 I}}*)
これらは、根7、9、および20に対応します。rootsRat
:
Position < 10^-1 &)](* {{7}, {9}, {20}}*)
私がPossibleZeroQ
でゼロの正確な方程式をチェックすると、より多くの根が得られます:
Position(* {{1}, {7}, {9}, {14}, {20}, {21}}*)
これらの潜在的な根に永遠にかかり、または私が待って喜んでいたよりも少なくとも長くなります。]
興味深いことに、これらの追加の根は、NSolve
によって返された結果とNSolve
によって返された結果との間の差に接続されています。Solve
:
Position, _?Positive](* {{1}, {2}, {3}, {14}, {21}, {22}, {23}, {24}}*)
私たちの新しいルーツを次のようにしましょう、そして私たちはそれらを以下でチェックすることができます:
rootsEqNew = Pick];rootsEqNew // N(* {{ζ -> -24559.6 + 24559.6 I}, {ζ -> -3.78042 - 5.50655 I}, {ζ -> -3.20562 + 5.39914 I}, {ζ -> -0.0832786 - 722.827 I}, {ζ -> 6.98478 + 0.493405 I}, {ζ -> 722.642 - 0.100823 I}}*)
チェック
最初は、新しい追加のためにあまりにも良い見ていません:
eqExact /. rootsEqNew // N(* { 0. + 9.7614*10^12 I, (* {1} *) -1.49012*10^-8 + 3.91155*10^-8 I, (* {7} *) -1.39698*10^-8 - 3.63216*10^-8 I, (* {9} *) -2. + 1536. I, (* {14} *) 1.49012*10^-8 + 1.11759*10^-8 I, (* {20} *) 0. + 1024. I} (* {21} *)*)
問題は、機械精度で方程式を評価することです。
N(* N::meprec warning about $MaxExtraPrecision being reached *)(* {0``69.62973568978663 + 0``69.69302870899077 I, 0``90.5174054423328 + 0``90.55817837498498 I, 0``90.50250822096415 + 0``90.54414468499085 I, 0``80.1824915073549 + 0``79.76650578675965 I, 0``90.47483650216002 + 0``90.49782363232914 I, 0``80.17292602755023 + 0``79.76710897249409 I}*)
警告は重要ではないようです。 $MaxExtraPrecision
を2000
にバンプすると、2020
から2040
の周りの精度で警告とゼロが生成されます。 それらがゼロの場合、は
Precision
(MachinePrecision
以外)を持つことができないため、N
はおそらく常に警告を生成します。