[心得] Graeffe's method
http://en.wikipedia.org/wiki/Dandelin%E2%80%93Gr%C3%A4ffe_method
最近在拿數值方法教科書玩樂。以Graeffe法找n次多項式的根的精神是
「不斷製造新的多項式,使得新的根=舊的根平方,使絕對值最大的根脫穎而出。」
有點像讓一群"根"一起跑電泳...有何陰謀呢?是為了讓根與係數關係中,忽略掉
較小根的影響,讓形式變得很簡單。
這回因為要模擬浮點精度學到了SetPrecision/SetAccuracy的用法是為收穫
= = = = = = = =
(*方便產生多項式,Mathematica除了較複雜的FromCoefficientRules,竟然別無方式
由List->Polynomial 真奇怪*)
ToPolynomial[k_List, v_] := v^Range[0, Length@k - 1].k
(*iter是根據維基p_odd,p_even法寫出來的*)
iter[pol_] :=
Module[{n, po, pe},
n = Length@pol; {pe, po} =
ToPolynomial[#,
x] & /@ (CoefficientList[pol, x][[Range[#, n, 2]]] & /@ {1, 2});
Expand[pe^2 - x po^2] (-1)^n]
Rand = RandomReal[ChiDistribution[1.5], 4] (*分布是亂配的可隨意*)
(*模擬精度有限,根依照絕對值大小排序*)
root = SetPrecision[Sort[Rand, Greater], 6]
(*產生多項式f*)
f = Product[x - i, {i, root}]
f = Expand@f
(*max iteration,可改*)m = 6;
(*輸出*)
iterlist = NestList[iter, f, m]; Column[iterlist]
Table[Reverse@
Abs[1/Ratios[(#[[;; -2]]/#[[-1]]) &[
CoefficientList[iterlist[[#]], x]]~Join~{1}]]^(1/
2^(# - 1)) &[i], {i, 2, m}] // Column
root
= = = = = = = = = = = = = = =
Graeffe法最被詬病的就是精度損失很大,SetPrecision設成誇張小就會發現
輸出全變成0,Indeterminate,或無限大。相反的優點是全局性的找根,精度的問題
就用牛頓法補救......有點有一好沒兩好之意。
話說虛根或重根不知道會發生什麼事
--
※ 發信站: 批踢踢實業坊(ptt.cc)
◆ From: 140.112.213.88
Mathematica 近期熱門文章
PTT數位生活區 即時熱門文章