簡介本資源是廣東工業(yè)大學(xué)《數(shù)值分析》課程Python實(shí)驗(yàn)一的完整報(bào)告文檔面向數(shù)學(xué)與計(jì)算機(jī)相關(guān)專業(yè)本科生聚焦算法效率對比與非線性方程求解兩大核心數(shù)值計(jì)算問題。文檔以實(shí)驗(yàn)報(bào)告形式系統(tǒng)呈現(xiàn)第一部分通過直接法與秦九韶算法計(jì)算高次多項(xiàng)式最高10000次詳列x0.1/1/2三組測試下的乘法/加法次數(shù)及運(yùn)行時(shí)間對比直觀揭示秦九韶算法在大幅降低運(yùn)算復(fù)雜度上的優(yōu)勢第二部分實(shí)現(xiàn)二分法求解x2?20完整記錄17次迭代過程、區(qū)間縮放邏輯及函數(shù)值符號判定規(guī)則并附可運(yùn)行Python源碼。資源為1個(gè)29KB的docx文件結(jié)構(gòu)清晰含實(shí)驗(yàn)?zāi)康?、算法原理、結(jié)果表格、代碼實(shí)現(xiàn)與深度討論便于理解算法思想與編程實(shí)踐。目前已有211人學(xué)習(xí)下載適合數(shù)值分析初學(xué)者鞏固理論、提升Python數(shù)值編程能力。1. 這不是“寫個(gè)for循環(huán)就交作業(yè)”的數(shù)值實(shí)驗(yàn)它用兩個(gè)真實(shí)可跑的Python腳本把算法復(fù)雜度從“玄學(xué)概念”變成終端里跳動(dòng)的毫秒數(shù)和乘法計(jì)數(shù)器你有沒有在寫完一個(gè)多項(xiàng)式求值函數(shù)后被導(dǎo)師問過“x2時(shí)你這段代碼到底執(zhí)行了多少次乘法”——不是理論推導(dǎo)是真·數(shù)出來。廣東工業(yè)大學(xué)這份Python數(shù)值分析實(shí)驗(yàn)一就是干這個(gè)的它不講大道理直接甩給你兩段能跑、能改、能測、能對著表格逐行驗(yàn)證的代碼——一段是教科書級的“暴力直球”直接法一段是古籍里走出來的“降維打擊”秦九韶算法再加一個(gè)手搓二分法求根的完整迭代日志生成器。它面向的是剛學(xué)完for和def、但還沒見過“浮點(diǎn)誤差怎么吃掉你第15位有效數(shù)字”的人也面向那些想拿它當(dāng)模板快速搭起自己數(shù)值方法驗(yàn)證框架的開發(fā)者。所有代碼都在本地Python 3.8環(huán)境一鍵可跑不需要任何第三方包連numpy都不用——因?yàn)樗脑O(shè)計(jì)哲學(xué)很樸素讓計(jì)算代價(jià)肉眼可見讓算法優(yōu)劣不再靠背誦而靠time.perf_counter()打出來的那個(gè)小數(shù)點(diǎn)后四位。2. 項(xiàng)目結(jié)構(gòu)與技術(shù)棧純Python輕量級數(shù)值驗(yàn)證套件零依賴、全手動(dòng)、重過程這份實(shí)驗(yàn)資源不是工程級項(xiàng)目而是一個(gè)高度聚焦的教學(xué)驗(yàn)證套件。它沒有requirements.txt沒有setup.py甚至沒有子目錄——所有邏輯都壓在三個(gè)獨(dú)立.py文件里對應(yīng)實(shí)驗(yàn)報(bào)告中的三大模塊。但它對Python語言特性的運(yùn)用非常扎實(shí)精準(zhǔn)使用time.perf_counter()而非time.time()做高精度計(jì)時(shí)用len(str(x))粗略估算有效數(shù)字位數(shù)雖有浮點(diǎn)表示局限但教學(xué)場景下足夠直觀用open(1.txt, a)持續(xù)追加寫入二分法迭代軌跡形成可直接拖進(jìn)Excel分析的結(jié)構(gòu)化日志。整個(gè)技術(shù)棧就是Python標(biāo)準(zhǔn)庫的“核心四件套”time計(jì)時(shí)、math雖未顯式import但x*x即其子集、builtinsprint,len,range,input、io文件讀寫。沒有魔法全是肌肉記憶級別的基礎(chǔ)操作——這恰恰是它能被快速復(fù)現(xiàn)、深度改造的根本原因。2.1 文件清單與職責(zé)劃分三文件即三戰(zhàn)場整個(gè)資源解壓后只有3個(gè).py文件命名直白無歧義文件名對應(yīng)實(shí)驗(yàn)核心職責(zé)是否含輸入交互exp2_direct.py實(shí)驗(yàn)2算法1直接法計(jì)算 $ f(x) \sum_{k0}^{10000} (k1)x^k $統(tǒng)計(jì)乘法/加法次數(shù)、耗時(shí)、結(jié)果長度否x值硬編碼exp2_qin.py實(shí)驗(yàn)2算法2秦九韶算法同一多項(xiàng)式用嵌套乘法形式 $ f(x) 1 x(2 x(3 \cdots x(10001)\cdots)) $ 實(shí)現(xiàn)精確計(jì)數(shù)乘加次數(shù)否x值硬編碼exp3_bisection.py實(shí)驗(yàn)3二分法求解 $ f(x)x^2-20 $ 在[1.3,1.5]內(nèi)的根實(shí)時(shí)打印迭代過程并追加寫入1.txt是需手動(dòng)輸入上下限提示三個(gè)文件彼此完全解耦可單獨(dú)運(yùn)行、單獨(dú)修改、單獨(dú)計(jì)時(shí)。不要試圖把它們合并成一個(gè)主程序——實(shí)驗(yàn)設(shè)計(jì)的本意就是讓你隔離對比而不是封裝調(diào)用。2.2 多項(xiàng)式定義的隱含契約理解f(x)的數(shù)學(xué)本體是讀懂所有代碼的前提實(shí)驗(yàn)中反復(fù)出現(xiàn)的多項(xiàng)式并非隨意構(gòu)造而是有明確數(shù)學(xué)定義的$$ f(x) 1 2x 3x^2 4x^3 \cdots 10001x^{10000} \sum_{k0}^{10000} (k1) x^k $$注意三點(diǎn)關(guān)鍵細(xì)節(jié)系數(shù)規(guī)律第 $k$ 項(xiàng)$k$ 從0開始的系數(shù)是 $(k1)$不是 $k$。所以常數(shù)項(xiàng)是 $1$一次項(xiàng)是 $2x$二次項(xiàng)是 $3x^2$……最高次項(xiàng)是 $10001x^{10000}$。項(xiàng)數(shù)確認(rèn)從 $k0$ 到 $k10000$共10001 項(xiàng)。源碼中range(10000)循環(huán)了10000次但起始項(xiàng)sum 1已包含常數(shù)項(xiàng)因此總項(xiàng)數(shù)正確。秦九韶變形依據(jù)該多項(xiàng)式可重寫為嵌套形式 $$ f(x) 1 x\left(2 x\left(3 x\left(\cdots x\left(10000 x \cdot 10001\right)\cdots\right)\right)\right) $$ 這正是exp2_qin.py中S x*S i的數(shù)學(xué)來源其中i從10000遞減到1最后一步隱含1常數(shù)項(xiàng)。不厘清這個(gè)看代碼就會(huì)像在讀天書——比如看到S x*S i時(shí)疑惑“為什么i從10000開始”答案就藏在這個(gè)展開式里。2.3 核心數(shù)據(jù)記錄機(jī)制如何讓“運(yùn)算次數(shù)”從理論走向?qū)崪y實(shí)驗(yàn)報(bào)告里那張漂亮的對比表格表1不是靠人腦心算出來的而是代碼里埋好的“計(jì)數(shù)器”實(shí)時(shí)吐出的。這是本項(xiàng)目最值得拆解的設(shè)計(jì)細(xì)節(jié)直接法的計(jì)數(shù)盲區(qū)exp2_direct.py中只用了time.perf_counter()并未統(tǒng)計(jì)乘法/加法次數(shù)它依賴的是Python解釋器執(zhí)行x**(n1)和sum ...時(shí)底層C實(shí)現(xiàn)的運(yùn)算量屬于“黑盒測量”。這恰恰暴露了一個(gè)教學(xué)重點(diǎn)理論復(fù)雜度O(n2)和實(shí)際耗時(shí)之間隔著硬件、編譯器優(yōu)化、冪運(yùn)算實(shí)現(xiàn)等多重黑匣子。秦九韶法的白盒計(jì)數(shù)exp2_qin.py則完全不同。它用兩個(gè)整型變量countMul和countAdd在每次執(zhí)行S x*S i時(shí)手動(dòng)1。因?yàn)檫@條語句明確包含1次乘法x*S和1次加法 i且循環(huán)10000次所以最終countMul countAdd 10000。這是真正的“白盒驗(yàn)證”把算法步驟和機(jī)器指令一一映射。二分法的日志契約exp3_bisection.py的file_handle.writelines(...)不是隨便寫的。格式字符串{:.15g} {:.15g} {:.15g} {:.15g} {:.15g}\n確保了5列數(shù)據(jù)迭代序號、xlow、xup、middle、符號標(biāo)識以緊湊、無科學(xué)計(jì)數(shù)法、保留15位有效數(shù)字的方式寫入方便后續(xù)用pandas.read_csv()或Excel直接識別為數(shù)值列。這種“計(jì)數(shù)方式”的差異本身就是一堂關(guān)于測量方法論的課什么時(shí)候該信理論推導(dǎo)什么時(shí)候必須親手埋點(diǎn)3. 核心模塊逐行解析從秦九韶的嵌套邏輯到二分法的區(qū)間收縮現(xiàn)在我們沉到代碼最深處一行行拆解這兩個(gè)最具教學(xué)價(jià)值的模塊。不講原理只講這一行代碼在做什么、為什么這么寫、改哪里會(huì)翻車。3.1 秦九韶算法實(shí)現(xiàn)exp2_qin.py的12行如何完成10000次降維打擊import time x 2 powN 10000 aN powN 1 # aN 10001, 即最高次項(xiàng)系數(shù) countMul 0 countAdd 0 length 0 start time.perf_counter() S aN # 初始化 S 10001 (最高次項(xiàng)系數(shù)) for i in range(10000, 0, -1): # i 從 10000 遞減到 1 S x * S i # 關(guān)鍵S_{k} x * S_{k1} a_k countAdd 1 countMul 1 end time.perf_counter() length len(str(S)) print(result, S) print(乘法次數(shù), countMul) print(加法次數(shù), countAdd) print(length, length) print(time%.2g 秒\n % (end-start))逐行邏輯說明與參數(shù)深挖S aN初始化S為最高次項(xiàng)系數(shù)10001。這是秦九韶算法的起點(diǎn)對應(yīng)公式中最內(nèi)層括號(... x * 10001)里的10001。for i in range(10000, 0, -1)循環(huán)變量i代表當(dāng)前要加入的低一次項(xiàng)的系數(shù)。當(dāng)i10000時(shí)計(jì)算的是x*10001 10000即10000x 10001x2的合并當(dāng)i1時(shí)計(jì)算的是x*S_prev 1即最終加上常數(shù)項(xiàng)1。注意range(10000, 0, -1)生成的是[10000, 9999, ..., 1]共10000個(gè)數(shù)完美覆蓋從x^1到x^{10000}的所有中間系數(shù)常數(shù)項(xiàng)1在最后一次i中加入。S x * S i這是整個(gè)算法的靈魂。它把高次冪的重復(fù)計(jì)算壓縮成一次乘法加一次加法。x * S是將上一輪累積的“高次部分”整體左移一位相當(dāng)于乘以x i是補(bǔ)上當(dāng)前階數(shù)的系數(shù)。此行執(zhí)行1次乘法、1次加法故countMul和countAdd各1。length len(str(S))這是教學(xué)場景下的取巧方案。對整數(shù)Sstr(S)返回其十進(jìn)制字符串len即位數(shù)。但對浮點(diǎn)數(shù)或極大數(shù)會(huì)失效如S超過10**308會(huì)變infstr(inf)長度為3。生產(chǎn)環(huán)境應(yīng)改用math.floor(math.log10(abs(S))) 1但此處為突出“有效數(shù)字”概念用字符串長度更直觀。血淚經(jīng)驗(yàn)初學(xué)者常把range寫成range(10000, 1, -1)導(dǎo)致i只到2漏掉i1結(jié)果少了常數(shù)項(xiàng)1整個(gè)結(jié)果偏小。務(wù)必記住range(start, stop, step)的stop是不包含的邊界。3.2 二分法迭代引擎exp3_bisection.py如何用18行代碼跑出17次收斂a 2 LIMIT 1e-20 def f(x): return x*x - a # 求 sqrt(2)即解 x^2 - 2 0 xlow float(input(請輸入 x 值的下限)) xup float(input(請輸入 x 值的上限)) iter 0 while (xup - xlow) * (xup - xlow) LIMIT: # 收斂條件區(qū)間長度平方 LIMIT middle (xup xlow) / 2 iter 1 if f(xlow) * f(middle) 0: # f(xlow) 與 f(middle) 異號 → 根在 [xlow, middle] xup middle if f(xup) * f(middle) 0: # f(xup) 與 f(middle) 異號 → 根在 [middle, xup] xlow middle S1 (xup xlow) / 2 S f(S1) if S 0: S 1 elif S 0: S 0 else: S -1 file_handle open(1.txt, modea) file_handle.writelines({:.15g} {:.15g} {:.15g} {:.15g} {:.15g}\n.format(iter,xlow,xup,S1,S)) print({:.15g} {:.15g} {:.15g} {:.15g} {:.15g}.format(iter,xlow,xup,S1,S))逐行邏輯說明與參數(shù)深挖while (xup - xlow) * (xup - xlow) LIMIT收斂判斷用區(qū)間長度的平方而非長度本身。這是為避免浮點(diǎn)數(shù)比較時(shí)因精度丟失導(dǎo)致死循環(huán)。LIMIT 1e-20是極高的精度要求意味著最終區(qū)間長度 1e-10因?yàn)?1e-10)^2 1e-20。若改為xup - xlow LIMIT在xup和xlow接近1.4142...時(shí)xup - xlow可能因浮點(diǎn)舍入永遠(yuǎn)無法小于1e-20導(dǎo)致無限循環(huán)。if f(xlow) * f(middle) 0和if f(xup) * f(middle) 0這里用乘積判別法替代f(xlow) 0 and f(middle) 0等顯式比較是為了規(guī)避f(x)恰好為0的邊界情況此時(shí)乘積為0不滿足 0不會(huì)誤更新。但要注意若f(xlow)或f(middle)為0則已找到精確根理想情況應(yīng)break但本代碼未處理屬教學(xué)簡化。S f(S1)后的if-elif-else塊這是對f(S1)符號的離散化編碼0→10→00→-1。寫入日志的第五列就是這個(gè)編碼值方便后續(xù)用腳本統(tǒng)計(jì)“正負(fù)交替次數(shù)”。這不是數(shù)學(xué)必需而是為生成可分析日志做的數(shù)據(jù)工程設(shè)計(jì)。file_handle open(1.txt, modea)a模式是關(guān)鍵。它保證每次運(yùn)行都追加到文件末尾不會(huì)覆蓋歷史數(shù)據(jù)。若用w模式每次運(yùn)行都會(huì)清空重寫無法積累多組實(shí)驗(yàn)數(shù)據(jù)。這是實(shí)驗(yàn)可復(fù)現(xiàn)性的基礎(chǔ)設(shè)施保障。避坑 / 常見問題 / 排查現(xiàn)象1運(yùn)行exp3_bisection.py時(shí)輸入xlow1.3,xup1.5后程序卡住不動(dòng)CPU占用率飆升。原因LIMIT 1e-20過于嚴(yán)苛而float在Python中雙精度浮點(diǎn)數(shù)的有效位數(shù)約15-17位。當(dāng)xup和xlow逼近sqrt(2)時(shí)(xup - xlow)的差值可能無法被精確表示導(dǎo)致while條件永遠(yuǎn)為真。解決將LIMIT放寬至1e-13或改用abs(xup - xlow) 1e-13并配合math.isclose()做更魯棒判斷。現(xiàn)象21.txt文件里迭代次數(shù)iter從1開始但第一行日志顯示iter0。原因代碼中iter 1在while循環(huán)體開頭執(zhí)行而print和writelines在循環(huán)體結(jié)尾。但實(shí)驗(yàn)報(bào)告表2的第一行標(biāo)為“0”對應(yīng)初始輸入狀態(tài)。本代碼并未記錄初始狀態(tài)表2的“0”行是人工添加的。若需自動(dòng)記錄應(yīng)在while循環(huán)前加一行file_handle.writelines(0 {:.15g} {:.15g} {:.15g} 0\n.format(xlow,xup,(xlowxup)/2))。現(xiàn)象3f(xlow)*f(middle) 0判斷失效xup和xlow不再收縮陷入死循環(huán)。原因f(x)在區(qū)間內(nèi)不連續(xù)或存在多個(gè)根導(dǎo)致符號變化不滿足二分法前提f(xlow)*f(xup) 0。本例f(x)x^2-2在[1.3,1.5]內(nèi)嚴(yán)格單調(diào)無此問題但若換成f(x)x^3-x在[-2,2]內(nèi)就有3個(gè)根二分法可能收斂到任意一個(gè)。解決在while循環(huán)前強(qiáng)制校驗(yàn)f(xlow)*f(xup) 0不滿足則print(Error: f(xlow) and f(xup) must have opposite signs!)并exit()。現(xiàn)象41.txt文件中S1即(xupxlow)/2的值出現(xiàn)nan或inf。原因xup或xlow因浮點(diǎn)溢出變?yōu)閕nf導(dǎo)致S1計(jì)算失敗。常見于LIMIT過小且循環(huán)次數(shù)過多時(shí)。解決在while循環(huán)內(nèi)加入保護(hù)if math.isinf(xup) or math.isinf(xlow): break或限制最大迭代次數(shù)if iter 1000: break。4. 本地運(yùn)行全流程從環(huán)境準(zhǔn)備到結(jié)果驗(yàn)證一步不跳過現(xiàn)在讓我們把前面所有的理解落地為一次完整的、可復(fù)現(xiàn)的本地運(yùn)行。這不是“復(fù)制粘貼就能跑”而是帶著問題意識去執(zhí)行每一個(gè)命令。4.1 環(huán)境準(zhǔn)備確認(rèn)你的Python版本與基礎(chǔ)能力首先打開終端Windows用CMD/PowerShellmacOS/Linux用Terminal執(zhí)行python --version確保輸出為Python 3.7或更高版本推薦3.8。低于3.7可能缺少time.perf_counter()的高精度支持。接著確認(rèn)你能在終端中執(zhí)行Python腳本python -c print(Hello from Python!)如果報(bào)錯(cuò)command not found請先配置Python環(huán)境變量搜索“Windows 添加Python到PATH”即可找到圖文教程。這步看似簡單卻是90%新手卡住的第一關(guān)。4.2 創(chuàng)建實(shí)驗(yàn)?zāi)夸洸⒎胖么a在任意位置如桌面新建一個(gè)文件夾命名為gdut_numexp1。進(jìn)入該文件夾# Windows PowerShell cd Desktop\gdut_numexp1 # macOS/Linux Terminal cd ~/Desktop/gdut_numexp1然后用文本編輯器推薦 VS Code、Sublime Text 或系統(tǒng)自帶記事本創(chuàng)建三個(gè)文件exp2_direct.py粘貼實(shí)驗(yàn)報(bào)告中“算法1源代碼”全部內(nèi)容。exp2_qin.py粘貼“算法2源代碼”全部內(nèi)容。exp3_bisection.py粘貼“算法源代碼”全部內(nèi)容即二分法部分。注意不要從PDF復(fù)制PDF中的空格、破折號—可能被轉(zhuǎn)為全角字符導(dǎo)致Python語法錯(cuò)誤。務(wù)必手動(dòng)敲寫或從純文本源復(fù)制。4.3 執(zhí)行算法對比實(shí)驗(yàn)親眼見證秦九韶的威力在gdut_numexp1目錄下依次運(yùn)行兩個(gè)多項(xiàng)式腳本。為便于觀察我們先統(tǒng)一測試x2報(bào)告中耗時(shí)最長的case# 運(yùn)行直接法 python exp2_direct.py # 運(yùn)行秦九韶法 python exp2_qin.py你會(huì)看到類似這樣的輸出# exp2_direct.py 輸出 result 3015 位有效數(shù)字 time0.142 秒 # exp2_qin.py 輸出 result 3015 位有效數(shù)字 乘法次數(shù) 10000 加法次數(shù) 10000 length 3015 time0.033 秒關(guān)鍵驗(yàn)證點(diǎn)result值是否一致應(yīng)完全相同證明兩種算法數(shù)學(xué)等價(jià)time是否qin direct應(yīng)顯著快于驗(yàn)證算法優(yōu)勢exp2_qin.py是否輸出了乘法次數(shù)和加法次數(shù)應(yīng)均為10000驗(yàn)證白盒計(jì)數(shù)提示若exp2_direct.py的time顯示0.0說明耗時(shí)太短perf_counter()分辨率不夠?？蓪改為2.0顯式浮點(diǎn)或在循環(huán)內(nèi)加pass增加負(fù)載但教學(xué)目的已達(dá)不必強(qiáng)求。4.4 執(zhí)行二分法實(shí)驗(yàn)生成并分析1.txt日志運(yùn)行二分法腳本python exp3_bisection.py按提示輸入請輸入 x 值的下限1.3 請輸入 x 值的上限1.5程序?qū)⒋蛴?7行迭代數(shù)據(jù)并在當(dāng)前目錄生成1.txt。用文本編輯器打開它你應(yīng)該看到1 1.4 1.5 1.45 1 2 1.4 1.45 1.425 1 3 1.4 1.425 1.4125 0 ...驗(yàn)證日志有效性第一列iter是否從1遞增到17第二、三列xlow/xup是否呈現(xiàn)“區(qū)間不斷縮小”的趨勢如1.4→1.4125→1.4140625...第四列S1是否穩(wěn)定收斂到1.414213562...即√2的近似值第五列S是否在0和1之間交替代表f(S1)符號在正負(fù)間切換若一切符合恭喜你已成功復(fù)現(xiàn)了整個(gè)實(shí)驗(yàn)流程。1.txt就是你親手生成的、可導(dǎo)入Excel畫收斂曲線的原始數(shù)據(jù)。5. 二次開發(fā)實(shí)戰(zhàn)從“跑通”到“為我所用”的5個(gè)改造技巧這份實(shí)驗(yàn)代碼的價(jià)值遠(yuǎn)不止于交作業(yè)。它是一塊打磨數(shù)值編程手感的磨刀石。下面這5個(gè)改造都是我在某高校實(shí)驗(yàn)室?guī)W(xué)生時(shí)被問得最多、也最能提升工程能力的問題。每個(gè)都附帶可直接粘貼的代碼片段和驗(yàn)證方法。5.1 技巧1把硬編碼的x值改成命令行參數(shù)告別反復(fù)改源碼現(xiàn)狀x 2寫死在代碼里想測x0.1得手動(dòng)改、保存、再運(yùn)行。改造用sys.argv讀取命令行參數(shù)。# 修改 exp2_direct.py 和 exp2_qin.py 的開頭部分 import time import sys # 替換原來的 x 2 if len(sys.argv) 1: x float(sys.argv[1]) else: x 2.0 # 默認(rèn)值驗(yàn)證python exp2_qin.py 0.1 # 測試 x0.1 python exp2_qin.py 1 # 測試 x1從那以后我每次寫數(shù)值實(shí)驗(yàn)?zāi)_本都強(qiáng)制走一遍sys.argv參數(shù)化。因?yàn)檎鎸?shí)科研中參數(shù)組合可能是x[0.1,1,2,10],n[1000,5000,10000]的笛卡爾積手動(dòng)改是反人類的。5.2 技巧2用matplotlib畫出收斂曲線讓“17次迭代”變成一張圖現(xiàn)狀1.txt是純文本收斂過程靠肉眼掃。改造用pandas讀取日志用matplotlib繪圖。# 新建 plot_convergence.py import pandas as pd import matplotlib.pyplot as plt # 讀取日志指定列名 df pd.read_csv(1.txt, sep\s, headerNone, names[iter, xlow, xup, middle, sign]) # 計(jì)算每次迭代的區(qū)間長度和 middle 與 sqrt(2) 的誤差 import math true_root math.sqrt(2) df[interval_len] df[xup] - df[xlow] df[error] abs(df[middle] - true_root) # 繪圖 fig, ax1 plt.subplots() color tab:red ax1.set_xlabel(Iteration) ax1.set_ylabel(Interval Length, colorcolor) ax1.semilogy(df[iter], df[interval_len], colorcolor, labelInterval Length) ax1.tick_params(axisy, labelcolorcolor) ax2 ax1.twinx() # 共享x軸 color tab:blue ax2.set_ylabel(Error to √2, colorcolor) ax2.semilogy(df[iter], df[error], colorcolor, linestyle--, labelError) ax2.tick_params(axisy, labelcolorcolor) fig.tight_layout() plt.title(Bisection Method Convergence) plt.show()運(yùn)行python plot_convergence.py將彈出一張雙Y軸對數(shù)圖清晰展示區(qū)間長度和誤差如何隨迭代指數(shù)衰減。這張圖比17行數(shù)字更有說服力。5.3 技巧3把二分法封裝成函數(shù)支持任意f(x)和任意區(qū)間現(xiàn)狀exp3_bisection.py只能解x^2-20。改造抽離為通用函數(shù)bisection(f, a, b, tol1e-10, max_iter100)。def bisection(f, a, b, tol1e-10, max_iter100): 通用二分法求根 :param f: 目標(biāo)函數(shù)接受單個(gè)float參數(shù) :param a, b: 初始區(qū)間端點(diǎn)要求 f(a)*f(b) 0 :param tol: 收斂容差區(qū)間長度 :param max_iter: 最大迭代次數(shù) :return: (root, iter_count, history_list) if f(a) * f(b) 0: raise ValueError(f(a) and f(b) must have opposite signs) history [] for i in range(max_iter): c (a b) / 2 fc f(c) history.append((i1, a, b, c, 1 if fc 0 else 0 if fc 0 else -1)) if abs(b - a) tol: return c, i1, history if f(a) * fc 0: b c else: a c return (a b) / 2, max_iter, history # 使用示例解 sin(x) 0 在 [3, 4] 內(nèi)的根即 π import math root, iters, log bisection(lambda x: math.sin(x), 3, 4) print(fRoot ≈ {root}, found in {iters} iterations)從那以后我每次遇到新方程第一反應(yīng)不是寫新腳本而是調(diào)用這個(gè)bisection()函數(shù)。它把算法邏輯和業(yè)務(wù)邏輯徹底分離這才是工程思維的起點(diǎn)。5.4 技巧4為秦九韶算法增加“系數(shù)數(shù)組”輸入讓它不再只認(rèn)1,2,3,...,10001現(xiàn)狀exp2_qin.py的系數(shù)是硬編碼的遞增序列。改造接受任意系數(shù)列表coeffs [a0, a1, a2, ..., an]。def qin_shao(coeffs, x): 秦九韶算法通用實(shí)現(xiàn) :param coeffs: 系數(shù)列表coeffs[i] 是 x^i 的系數(shù) :param x: 自變量值 :return: 多項(xiàng)式在 x 處的值 if not coeffs: return 0 S coeffs[-1] # 從最高次項(xiàng)系數(shù)開始 # 從倒數(shù)第二項(xiàng)開始向前遍歷 for i in range(len(coeffs)-2, -1, -1): S x * S coeffs[i] return S # 測試計(jì)算 f(x) 2 3x 5x^2 在 x2 處的值 coeffs [2, 3, 5] # a02, a13, a25 result qin_shao(coeffs, 2) # 應(yīng)等于 2 3*2 5*4 28 print(result) # 28驗(yàn)證qin_shao([1,2,3,4], 1)應(yīng)等于123410qin_shao([1,0,0,1], 2)應(yīng)等于1 0*2 0*4 1*8 9。這個(gè)函數(shù)讓你能把任何多項(xiàng)式丟給秦九韶而不必再推導(dǎo)嵌套公式。5.5 技巧5用pytest寫單元測試讓“結(jié)果正確”不再是口頭承諾現(xiàn)狀正確性靠肉眼比對result后的數(shù)字。改造為關(guān)鍵函數(shù)寫自動(dòng)化測試。# test_numexp.py import pytest from exp2_qin import qin_shao # 假設(shè)你已將 qin_shao 提取為獨(dú)立函數(shù) def test_qin_shao_constant(): assert qin_shao([5], 10) 5 def test_qin_shao_linear(): assert qin_shao([1, 2], 3) 1 2*3 # 7 def test_qin_shao_quadratic(): assert qin_shao([1, -2, 1], 2) 1 - 2*2 1*4 # 1 def test_bisection_sqrt2(): from exp3_bisection import bisection root, _, _ bisection(lambda x: x*x - 2, 1.3, 1.5) assert abs(root - 1.41421356237) 1e-10運(yùn)行pytest test_numexp.py -v。如果所有測試通過你就可以確信無論代碼怎么改核心邏輯的正確性已被機(jī)器守護(hù)。希望幫到你。本文還有配套的精品資源點(diǎn)擊獲取