Hur kan jag få lösningen av komplicerad ekvation?

introduktion

Uppdatering 6/2014

ursprungligen hittade jag tre lösningar; nu är det sex. Jag trodde att jag visste att det fanns tre eftersom rationalisering av ekvationen resulterar i ett polynom av grad 24 och NSolve hittar 24 rötter och eliminerar sedan rötter som inte var nollor i den ursprungliga ekvationen. Det visar sig att ekvationen var besvärlig ur numerikens synvinkel och jag missade tre.

andra ändringar:

  • det verkar mer naturligt att använda Rationalize för att konvertera koefficienterna till exakta siffror i detta fall än SetPrecision.

  • på samma sätt ersatte jag First @ eq med Subtract @@ eq.

rationalisera ekvationen

man bör komma ihåg att enligt dokumentationen om NSolve

NSolve behandlar främst linjära och polynomekvationer.

OP: s ekvation är en algebraisk ekvation som kan omvandlas till en polynomekvation, som NSolve sedan enkelt kan lösa. Eftersom NSolve inte gör detta på egen hand, måste vi rationalisera “för hand.”Rationalisering tenderar att införa externa lösningar, så vi måste kontrollera svaren som returneras av NSolve. De något otäcka koefficienterna på eq gör det svårare,eftersom rationalisering av eq leder till att några av kvadratrötterna inte avbryter som de borde. Vi kan fixa det genom att ställa in precisionen eqtill Infinity.

vi kan hitta kvadratrötterna i ett uttryck expr med

DeleteDuplicates @ Cases, Infinity]

sedan för varje kvadratrot kan vi multiplicera expr med uttrycket med kvadratroten ersatt av dess negativa och förenkla.

nedan är resultatet av att metoden tillämpas på funktionen som erhålls genom att ersätta Equal i eq med Subtract med Apply(@@). Vi kan lösa ekvationen exakt.

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}}*)

det finns några signifikanta skillnader mellan de två lösningarna:

(ζ /. rootsRat) - (ζ /. rootsRatExact) // Abs // Max(* 3.15226*10^-10*)

välja rötterna för den givna ekvationen

Uppdatering: här missade jag några nollor.

jag använde Pick för att välja rötter där den ursprungliga ekvationen eqhar ett litet värde, mindre än något litet tröskelvärde, säg 10^-1. Vi kan kontrollera senare att var och en är en rot, men jag kollade inte efter andra rötter. De ursprungliga rötterna var:

rootsEq = Pick < 10^-1 /. rootsRat](* {{ζ -> -3.78042 - 5.50655 I}, {ζ -> -3.20562 + 5.39914 I}, {ζ -> 6.98478 + 0.493405 I}}*)

dessa motsvarar rötterna 7, 9 och 20 i rootsRat:

Position < 10^-1 &)](* {{7}, {9}, {20}}*)

om jag kontrollerar den exakta ekvationen för nollor med PossibleZeroQ får jag fler rötter:

Position(* {{1}, {7}, {9}, {14}, {20}, {21}}*)

tar evigt på dessa potentiella rötter, eller åtminstone längre än jag var villig att vänta.]

intressant är att dessa ytterligare rötter är kopplade till skillnaderna mellan resultaten som returneras av NSolve och Solve:

Position, _?Positive](* {{1}, {2}, {3}, {14}, {21}, {22}, {23}, {24}}*)

Låt oss låta våra nya rötter vara följande och vi kan kontrollera dem nedan:

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}}*)

kontroller

först, at ser inte för bra ut för de nya tilläggen:

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} *)*)

problemet är att utvärdera ekvationen med maskinprecision.

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}*)

varningen verkar inte vara signifikant. Bumping up $MaxExtraPrecision till 2000 producerar fortfarande varningen och nollorna med noggrannhet runt 2020till 2040. Om de är nollor kommer N förmodligen alltid att producera varningen, eftersom inte kan ha en Precision (annat än MachinePrecision).

Lämna ett svar

Din e-postadress kommer inte publiceras.