簡(jiǎn)介這是一份面向數(shù)學(xué)、物理及工程領(lǐng)域?qū)W習(xí)者的偏微分方程數(shù)值解入門(mén)資料內(nèi)容從橢圓型、拋物型、雙曲型三類(lèi)典型方程出發(fā)系統(tǒng)講解泊松方程、熱傳導(dǎo)方程和波動(dòng)方程的定解問(wèn)題并重點(diǎn)展示如何在MATLAB中通過(guò)有限差分法、pdepe等工具進(jìn)行數(shù)值求解適合需要結(jié)合代碼理解偏微分方程數(shù)值方法的本科高年級(jí)學(xué)生或科研人員。壓縮包內(nèi)包含1個(gè)PDF文件體積僅1.67MB便于下載后隨時(shí)查閱。目前已有222人學(xué)習(xí)瀏覽內(nèi)容結(jié)構(gòu)清晰完整覆蓋網(wǎng)格剖分、差分格式構(gòu)造、邊值問(wèn)題處理、解的收斂性與誤差分析等關(guān)鍵環(huán)節(jié)同時(shí)給出差分方法的實(shí)現(xiàn)思路與MATLAB函數(shù)使用要點(diǎn)對(duì)初值和邊界條件的分類(lèi)也做了細(xì)致梳理。整體而言這是一份兼具理論梳理和實(shí)操參考價(jià)值的緊湊型學(xué)習(xí)手冊(cè)。 手邊這本從網(wǎng)上下載的《偏微分方程-matlab.pdf》封面已經(jīng)有點(diǎn)舊了但我還是愿意把它放在“常翻”文件夾里。原因很簡(jiǎn)單偏微分方程的數(shù)值求解MATLAB 仍然是我見(jiàn)過(guò)最省心的環(huán)境。不管你是做課程作業(yè)、畢業(yè)論文還是工程項(xiàng)目里要快速驗(yàn)證某個(gè)擴(kuò)散、波動(dòng)或傳熱過(guò)程用 MATLAB 搭數(shù)值模型基本都能在很短時(shí)間里得到能直接分析的結(jié)果。這篇文章不打算堆數(shù)學(xué)公式而是把我在使用這份 PDF 資料和實(shí)際項(xiàng)目過(guò)程中積累的求解方法、代碼模板、以及踩過(guò)的坑整理出來(lái)給同樣被偏微分方程折磨過(guò)的人一份可以直接參考的實(shí)操記錄。1. 拿到《偏微分方程-matlab.pdf》后我重新梳理的求解路線很多人把偏微分方程當(dāng)成一門(mén)純數(shù)學(xué)課來(lái)學(xué)動(dòng)輒翻公式翻到懷疑人生。但如果你手里有 MATLAB其實(shí)可以換一種思路先把它當(dāng)成工具跑通一個(gè)結(jié)果再回頭補(bǔ)理論。我最早看這本 PDF 的時(shí)候就是被里面“用數(shù)值方法求解 PDE本質(zhì)上是把無(wú)限維問(wèn)題變成有限維問(wèn)題”這句話點(diǎn)醒的。從那以后我遇到實(shí)際的偏微分方程問(wèn)題基本都按同一條流程走先判斷方程類(lèi)型再確認(rèn)定解條件然后選求解路徑最后寫(xiě)代碼、調(diào)精度、驗(yàn)證結(jié)果。1.1 為什么 MATLAB 值得用來(lái)啃偏微分方程MATLAB 在這類(lèi)問(wèn)題上有三個(gè)天然優(yōu)勢(shì)。第一矩陣運(yùn)算是它的看家本領(lǐng)而有限差分、有限元、譜方法最終都落在矩陣運(yùn)算上第二內(nèi)置可視化工具非常成熟求解完直接就能畫(huà)溫度場(chǎng)、波場(chǎng)、濃度場(chǎng)不需要額外裝庫(kù)第三工具箱和內(nèi)置函數(shù)覆蓋了大量常用場(chǎng)景比如 pdepe、PDE Toolbox以及更底層的 ode 求解器用來(lái)處理半離散化的 PDE 問(wèn)題也足夠順手。我在項(xiàng)目里經(jīng)常要快速驗(yàn)證一個(gè)傳熱模型在 Python 里要用 numpy、scipy、matplotlib 搭一堆環(huán)境在 MATLAB 里只需要一個(gè)腳本文件就能從定義方程到出圖全部完成。當(dāng)然Python 也很優(yōu)秀但如果你手上的資料就是《偏微分方程-matlab.pdf》這種教材按它的思路走是最省力的。1.2 從連續(xù)到離散貫穿所有 PDE 求解的一條主線偏微分方程描述的是連續(xù)變量之間的關(guān)系但計(jì)算機(jī)只能處理有限個(gè)數(shù)字。所以數(shù)值求解的核心思路就是把連續(xù)問(wèn)題離散成代數(shù)方程組。比如一維熱傳導(dǎo)方程[ \frac{\partial u}{\partial t} \kappa \frac{\partial^2 u}{\partial x^2} ]計(jì)算機(jī)并不知道這個(gè)偏導(dǎo)數(shù)怎么算我們得把空間 (x) 和時(shí)間 (t) 切成很多小段在網(wǎng)格點(diǎn)上用差分公式近似導(dǎo)數(shù)。這個(gè)過(guò)程叫離散化。離散化之后原來(lái)的 PDE 就變成了一系列關(guān)于網(wǎng)格點(diǎn)上的值 (u_{i}^{n}) 的代數(shù)關(guān)系式。MATLAB 做的無(wú)非就是幫你高效地解這些代數(shù)關(guān)系式。理解這條主線比記住某一類(lèi)方程的特定解法重要得多。因?yàn)椴还苁?pdepe、PDE Toolbox還是自己寫(xiě)有限差分本質(zhì)上都在做這件事。1.3 拿到方程先做類(lèi)型判斷再?zèng)Q定用哪個(gè)求解器不是所有偏微分方程都能用同一個(gè)函數(shù)解。拿到一個(gè) PDE我一般先問(wèn)三個(gè)問(wèn)題方程是什么類(lèi)型橢圓、拋物還是雙曲自變量是哪些維度有多高邊界條件和初始條件給的是完整這三個(gè)問(wèn)題直接決定選型。我把常用場(chǎng)景整理成了下面的表格方便對(duì)照方程類(lèi)型經(jīng)典例子數(shù)學(xué)特征MATLAB 常用路徑橢圓方程拉普拉斯方程、泊松方程穩(wěn)態(tài)問(wèn)題沒(méi)有時(shí)間項(xiàng)PDE Toolbox或自寫(xiě)迭代法拋物方程熱傳導(dǎo)方程、擴(kuò)散方程含一階時(shí)間導(dǎo)數(shù)擴(kuò)散特性pdepe或自寫(xiě)隱式差分雙曲方程波動(dòng)方程、對(duì)流方程含二階時(shí)間導(dǎo)數(shù)或一階時(shí)間導(dǎo)數(shù)強(qiáng)對(duì)流項(xiàng)pdepe 部分適用復(fù)雜情況自寫(xiě)差分/FVM我見(jiàn)的比較多的坑是有人拿橢圓方程的求解器去解雙曲方程結(jié)果波峰處產(chǎn)生劇烈振蕩還以為是參數(shù)設(shè)置問(wèn)題。其實(shí)第一步方程類(lèi)型沒(méi)判斷對(duì)后面全是白費(fèi)。判斷類(lèi)型的方法不復(fù)雜就看最高階導(dǎo)數(shù)全是空間二階導(dǎo)、沒(méi)有時(shí)間項(xiàng)的基本是橢圓含時(shí)間一階導(dǎo)和空間二階導(dǎo)的基本是拋物含時(shí)間二階導(dǎo)或者有強(qiáng)一階對(duì)流項(xiàng)的基本是雙曲。2. 動(dòng)手前必須想清楚的事定解條件、離散化與邊界處理我見(jiàn)過(guò)不少同學(xué)方程拿到手就急著寫(xiě)代碼結(jié)果算出來(lái)一張漂亮的圖拿去和解析解一比完全對(duì)不上。原因往往不是差分格式不對(duì)而是邊界條件或者初始條件沒(méi)理清楚。定解條件就像題目的“額外規(guī)則”少了它方程的解不唯一錯(cuò)了它整個(gè)數(shù)值解就是一場(chǎng)大型幻覺(jué)。2.1 三類(lèi)邊界條件必須熟到條件反射偏微分方程的邊界條件常見(jiàn)有三種MATLAB 里 pdepe 的寫(xiě)法也基于這三類(lèi)第一類(lèi)邊界條件Dirichlet直接給定邊界上的函數(shù)值。比如桿的一端溫度固定為 100 度就是 (u(0,t)100)。第二類(lèi)邊界條件Neumann給定邊界上的導(dǎo)數(shù)或通量。比如桿的末端絕熱即傳熱通量為 0也就是 (\frac{\partial u}{\partial x}0)。第三類(lèi)邊界條件Robin函數(shù)值和導(dǎo)數(shù)的線性組合常見(jiàn)于對(duì)流換熱。比如邊界上有熱交換可寫(xiě)成 (\alpha u \beta \frac{\partial u}{\partial x} \gamma)。在 pdepe 里面邊界條件統(tǒng)一寫(xiě)成 (p q \cdot f 0) 的形式其中 (p) 和 (q) 是我們需要返回的函數(shù)值(f) 是 PDE 方程里的通量項(xiàng)。一開(kāi)始接觸會(huì)有點(diǎn)繞但多寫(xiě)幾次就順了。2.2 離散化基礎(chǔ)前向、后向、中心差分離散化是所有數(shù)值解法的地基。以一維網(wǎng)格點(diǎn) (x_i) 上某個(gè)函數(shù) (u(x_i)) 為例常見(jiàn)的差分近似有前向差分(\frac{du}{dx} \approx \frac{u_{i1} - u_i}{\Delta x})后向差分(\frac{du}{dx} \approx \frac{u_i - u_{i-1}}{\Delta x})中心差分(\frac{du}{dx} \approx \frac{u_{i1} - u_{i-1}}{2\Delta x})二階導(dǎo) (\frac{d^2u}{dx^2} \approx \frac{u_{i1} - 2u_i u_{i-1}}{\Delta x^2})其中中心差分的精度通常更高而前向或后向在顯式時(shí)間推進(jìn)時(shí)實(shí)現(xiàn)更簡(jiǎn)單。我個(gè)人習(xí)慣空間二階導(dǎo)一律用中心差分時(shí)間一階導(dǎo)根據(jù)情況選擇顯式歐拉或隱式格式。2.3 在 MATLAB 里“定義微分方程”的標(biāo)準(zhǔn)姿勢(shì)不少人搜過(guò)“matlab 中定義微分方程”這類(lèi)詞其實(shí)在 MATLAB 里定義微分方程的關(guān)鍵不是像寫(xiě)數(shù)學(xué)公式一樣寫(xiě)出一整行 (dt/dx)而是把方程拆成標(biāo)準(zhǔn)形式再寫(xiě)成函數(shù)。以 pdepe 為例它要求方程寫(xiě)成[ c(x,t,u,\frac{\partial u}{\partial x}) \cdot \frac{\partial u}{\partial t} \frac{\partial}{\partial x} \left( f(x,t,u,\frac{\partial u}{\partial x}) \right) s(x,t,u,\frac{\partial u}{\partial x}) ]這里三個(gè)函數(shù)分別叫 (c)、(f)、(s)然后在代碼里對(duì)應(yīng)定義三個(gè)子函數(shù)。一開(kāi)始我覺(jué)得這個(gè)格式很死板后來(lái)發(fā)現(xiàn)它把“擴(kuò)散項(xiàng)、源項(xiàng)、時(shí)間系數(shù)”拆得清清楚楚反而有利于排查錯(cuò)誤。2.4 步長(zhǎng)和網(wǎng)格設(shè)計(jì)的黃金法則網(wǎng)格步長(zhǎng)不是越小越好因?yàn)椴介L(zhǎng)太小會(huì)顯著增加計(jì)算量而步長(zhǎng)太大則會(huì)導(dǎo)致精度不足甚至數(shù)值發(fā)散。對(duì)一個(gè)拋物型方程如果使用顯式格式推進(jìn)時(shí)間時(shí)間步長(zhǎng) (\Delta t) 和空間步長(zhǎng) (\Delta x) 之間必須滿足穩(wěn)定性限制典型條件形如[ \frac{\kappa \Delta t}{\Delta x^2} \le \frac{1}{2} ]這個(gè)條件叫 CFL 條件或擴(kuò)散穩(wěn)定性條件。我第一次寫(xiě)顯式差分程序時(shí)完全沒(méi)在意這個(gè)限制結(jié)果溫度場(chǎng)越算越奇怪甚至出現(xiàn)了負(fù)溫度最后查了很久才發(fā)現(xiàn)是步長(zhǎng)匹配的問(wèn)題。所以網(wǎng)格設(shè)計(jì)要綜合考慮穩(wěn)定性、精度和計(jì)算資源不能盲目加密。3. 手把手跑通一個(gè)模型熱傳導(dǎo)方程的 pdepe 實(shí)現(xiàn)看完前面那些基礎(chǔ)我們來(lái)點(diǎn)能直接上手的。以經(jīng)典的一維熱傳導(dǎo)方程為例目標(biāo)區(qū)間是 (x \in [0,1])時(shí)間從 0 到 0.5擴(kuò)散系數(shù) (\kappa 1)初始溫度分布為 (u(x,0) \sin(\pi x))兩端邊界溫度固定為 0。這是一個(gè)教科書(shū)式的拋物型方程非常適合用 pdepe 跑通。3.1 三個(gè)子函數(shù)把數(shù)學(xué)方程翻譯成 MATLAB 代碼先給整段代碼后面我再逐行解釋關(guān)鍵點(diǎn)。新建一個(gè)腳本命名為heat_conduction_pdepe.m內(nèi)容如下function heat_conduction_pdepe m 0; x linspace(0, 1, 101); t linspace(0, 0.5, 101); sol pdepe(m, pdefun, icfun, bcfun, x, t); u sol(:, :, 1); surf(x, t, u) xlabel(x) ylabel(t) zlabel(u) title(一維熱傳導(dǎo)方程的 pdepe 求解結(jié)果) end function [c, f, s] pdefun(x, t, u, DuDx) kappa 1; c 1; f kappa * DuDx; s 0; end function u0 icfun(x) u0 sin(pi * x); end function [pl, ql, pr, qr] bcfun(xl, ul, xr, ur, t) pl ul; ql 0; pr ur; qr 0; end3.2 逐段解釋pdepe 到底在解什么pdepe的第一個(gè)參數(shù)m代表坐標(biāo)系類(lèi)型0 代表直角坐標(biāo)1 代表柱坐標(biāo)2 代表球坐標(biāo)。我們這個(gè)問(wèn)題是一維直角坐標(biāo)所以填 0。接著傳進(jìn)去的四個(gè)函數(shù)句柄分別對(duì)應(yīng) PDE 系數(shù)、初始條件、邊界條件。最后的x和t是網(wǎng)格點(diǎn)返回的sol是一個(gè)三維數(shù)組第一維是時(shí)間點(diǎn)第二維是空間點(diǎn)第三維是方程解的個(gè)數(shù)。pdefun里的寫(xiě)法非常直觀方程的標(biāo)準(zhǔn)形式 (c \cdot u_t \partial_x(f) s)這里 (c1)(f \kappa u_x)(s0)。icfun就是初始條件返回一個(gè)值sin(pi*x)。bcfun里返回的四項(xiàng)分別代表左右邊界的 (p) 和 (q)因?yàn)檫吔鐥l件是 (p q \cdot f 0)對(duì)于固定溫度邊界 (u0)可以直接令 (pu)(q0)這樣邊界條件就變成了 (u0)。3.3 結(jié)果怎么驗(yàn)證解析解是最快的照妖鏡這個(gè)熱傳導(dǎo)方程的解析解很容易求[ u(x,t) \sin(\pi x) \cdot e^{-\pi^2 \kappa t} ]我通常會(huì)在代碼里加一個(gè)對(duì)照部分把解析解算出來(lái)然后算最大絕對(duì)誤差u_analytical sin(pi * x) * exp(-pi^2 * t); error max(max(abs(u - u_analytical))); disp([最大絕對(duì)誤差: , num2str(error)])實(shí)測(cè)下來(lái)pdepe 在默認(rèn)容差下精度很不錯(cuò)最大誤差通常在 (10^{-4}) 量級(jí)。這一步非常關(guān)鍵一旦代碼有小錯(cuò)誤對(duì)比解析解能第一時(shí)間發(fā)現(xiàn)而不是等畫(huà)出來(lái)的圖變得稀奇古怪才去排查。4. 踩過(guò)的坑與調(diào)試筆記從數(shù)值振蕩到矩陣運(yùn)算代碼跑通只是第一步真正花時(shí)間的是調(diào)試。這里我把遇到過(guò)的幾個(gè)典型問(wèn)題展開(kāi)講講全是實(shí)際操作中容易踩的地方。4.1 顯式差分振蕩與穩(wěn)定性條件用 pdepe 很少遇到顯式差分那種“炸場(chǎng)”但自己寫(xiě)顯式格式時(shí)穩(wěn)定性問(wèn)題幾乎是繞不開(kāi)的。我之前寫(xiě)過(guò)一段求解擴(kuò)散方程的顯式差分代碼核心迭代就三行u(2:end-1) u(2:end-1) ... kappa * dt / dx^2 * (u(3:end) - 2*u(2:end-1) u(1:end-2));如果 ( \kappa dt / dx^2 0.5 )數(shù)值解就會(huì)開(kāi)始出現(xiàn)高頻振蕩幅度越來(lái)越大最后全是 NaN。這不是代碼邏輯錯(cuò)誤而是格式本身就發(fā)散。解決辦法有兩條把時(shí)間步長(zhǎng)壓縮到穩(wěn)定范圍內(nèi)或者改用隱式格式。隱式格式?jīng)]有這個(gè)穩(wěn)定性限制是實(shí)際工程中更穩(wěn)妥的選擇。4.2 點(diǎn)乘和直接乘的區(qū)別一個(gè)不起眼但致命的錯(cuò)誤很多剛接觸 MATLAB 的人會(huì)搜“matlab 中點(diǎn)乘和直接乘有啥區(qū)別”。在偏微分方程數(shù)值解中這個(gè)問(wèn)題特別容易跳出來(lái)咬你一口。比如在有限差分迭代里如果某個(gè)系數(shù)向量要和網(wǎng)格上的值逐點(diǎn)相乘必須用.*如果誤寫(xiě)成*MATLAB 就會(huì)按照矩陣乘法去計(jì)算輕則報(bào)維度錯(cuò)誤重則矩陣形狀碰巧匹配但結(jié)果完全錯(cuò)了。我印象最深的一次是在求解帶變系數(shù)的擴(kuò)散方程時(shí)把kappa .* (u(i1)-u(i))寫(xiě)成了kappa * (u(i1)-u(i))而kappa恰好是一個(gè)向量另一個(gè)是向量MATLAB 沒(méi)報(bào)錯(cuò)因?yàn)閮烧叩木S度正好匹配但計(jì)算結(jié)果完全不對(duì)整個(gè)濃度場(chǎng)都扭曲了。從那次以后我在涉及逐點(diǎn)運(yùn)算的代碼里都會(huì)特意檢查乘號(hào)前面有沒(méi)有點(diǎn)。4.3 pdepe 報(bào)錯(cuò)和邊界條件函數(shù)的分工pdepe 報(bào)錯(cuò)常見(jiàn)的有幾類(lèi)。一類(lèi)是“Spatial discretization has failed”大概率是邊界條件函數(shù)返回的pl、ql、pr、qr某些項(xiàng)寫(xiě)錯(cuò)了導(dǎo)致離散后的 Jacobian 矩陣奇異。還有一類(lèi)是求解到中途出現(xiàn)NaN一般是系數(shù)函數(shù)pdefun里傳入了u和DuDx結(jié)果你在里面做了類(lèi)似1/u的操作而某些時(shí)刻u會(huì)正好等于 0。調(diào)試邊界條件時(shí)我習(xí)慣先把左右邊界條件都簡(jiǎn)化成 Dirichlet 形式跑通之后再逐步加復(fù)雜條件。這樣能把問(wèn)題范圍縮小不會(huì)一上來(lái)就是一團(tuán)亂麻。4.4 驗(yàn)證數(shù)值解的幾板斧寫(xiě)偏微分方程求解代碼我基本都會(huì)做這幾步驗(yàn)證找解析解對(duì)比沒(méi)有解析解就想辦法構(gòu)造一個(gè)帶有已知解的“人造源項(xiàng)”。檢查守恒量。比如熱傳導(dǎo)方程無(wú)源無(wú)熱交換時(shí)若邊界是絕熱的總熱量應(yīng)該不隨時(shí)間變化。做網(wǎng)格收斂性測(cè)試步長(zhǎng)減半看誤差是否按預(yù)期比例下降。這三個(gè)方法在絕大多數(shù)情況下都能把隱蔽的 bug 挖出來(lái)。我自己的體會(huì)是越早做驗(yàn)證后面返工越少。5. 進(jìn)階用 PDE Toolbox 和自寫(xiě)有限差分應(yīng)對(duì)更復(fù)雜的問(wèn)題pdepe 雖然好用但它主要解決一維問(wèn)題。遇到二維、三維或者不規(guī)則幾何區(qū)域時(shí)就要升級(jí)工具了。5.1 PDE Toolbox 適合的場(chǎng)景PDE Toolbox 是 MATLAB 專(zhuān)門(mén)用來(lái)求解 PDE 的圖形化工具箱支持有限元方法能處理二維和三維區(qū)域上的橢圓方程、拋物方程、雙曲方程也支持多個(gè)方程聯(lián)立。對(duì)于傳熱、結(jié)構(gòu)力學(xué)、電磁場(chǎng)這類(lèi)問(wèn)題我通常直接用工具箱里預(yù)定義的方程模塊來(lái)建模型設(shè)置好幾何區(qū)域、材料參數(shù)、邊界條件和初始條件然后生成網(wǎng)格并求解。尤其當(dāng)幾何區(qū)域復(fù)雜時(shí)手寫(xiě)網(wǎng)格生成和組裝有限元矩陣是非常費(fèi)時(shí)的事情PDE Toolbox 能把這部分時(shí)間省下來(lái)。5.2 什么時(shí)候值得自寫(xiě)有限差分自寫(xiě)差分最大的優(yōu)勢(shì)是透明、可控適合教學(xué)演示也適合處理一些工具箱覆蓋不到的非標(biāo)準(zhǔn)方程。比如帶上強(qiáng)非線性源項(xiàng)或者需要人為控制離散格式細(xì)節(jié)的時(shí)候?qū)懽约旱牡a往往比強(qiáng)行套工具箱更靈活。對(duì)于一維問(wèn)題我習(xí)慣先寫(xiě)一個(gè)顯式模板做快速試探確認(rèn)方程形態(tài)沒(méi)問(wèn)題后再改寫(xiě)成隱式格式用于正式計(jì)算。5.3 一維擴(kuò)散方程的隱式差分模板這里給一個(gè)最常用的隱式差分模板使用 Crank-Nicolson 格式兼顧精度和穩(wěn)定性function u implicit_diffusion_1d(kappa, L, T, Nx, Nt) dx L / (Nx - 1); dt T / Nt; u zeros(Nt, Nx); u(1, :) sin(pi * linspace(0, L, Nx)); r kappa * dt / dx^2; A diag((1 r) * ones(Nx, 1)) ... diag(-r/2 * ones(Nx-1, 1), 1) ... diag(-r/2 * ones(Nx-1, 1), -1); A(1, :) 0; A(1, 1) 1; A(end, :) 0; A(end, end) 1; for n 1:Nt-1 b u(n, :).; b(2:end-1) b(2:end-1) ... r/2 * (u(n, 3:end) - 2*u(n, 2:end-1) u(n, 1:end-2)).; b(1) 0; b(end) 0; u(n1, :) (A \ b).; end end這個(gè)模板雖然不復(fù)雜但已經(jīng)可以處理很多實(shí)際擴(kuò)散問(wèn)題了。如果你要處理非線性項(xiàng)只需要在組裝A或計(jì)算右端項(xiàng)時(shí)加入對(duì)應(yīng)的線性化處理。5.4 后續(xù)可以繼續(xù)擴(kuò)展的方向偏微分方程本身不是終點(diǎn)算出來(lái)的結(jié)果經(jīng)常還要對(duì)接后續(xù)分析。比如用lsqcurvefit根據(jù)實(shí)驗(yàn)觀測(cè)數(shù)據(jù)反演方程中的未知擴(kuò)散系數(shù)或者把 PDE 模型的輸出和卷積、系統(tǒng)辨識(shí)等方法結(jié)合用于更復(fù)雜的動(dòng)態(tài)系統(tǒng)建模。這些方向需要單獨(dú)寫(xiě)文章來(lái)講但核心還是先把 PDE 數(shù)值解做穩(wěn)后面才有底氣做反向題。我自己這些年用 MATLAB 解 PDE最大的體會(huì)是“先跑通再優(yōu)化”。第一版代碼往往不是最高效的但它能讓你快速知道方程行為、參數(shù)范圍和潛在問(wèn)題。等一切都在掌控之中再去考慮更精細(xì)的格式、更快的算法或者更復(fù)雜的模型完全來(lái)得及。本文還有配套的精品資源點(diǎn)擊獲取