簡(jiǎn)介一份基于MATLAB的6節(jié)點(diǎn)天然氣潮流計(jì)算教學(xué)程序面向能源動(dòng)力及相關(guān)專業(yè)初學(xué)者用于理解天然氣網(wǎng)絡(luò)壓力、流量與儲(chǔ)存量的求解方法。程序以6節(jié)點(diǎn)簡(jiǎn)化模型為對(duì)象完整涉及狀態(tài)方程、能量守恒、管道阻力與壓降計(jì)算等核心環(huán)節(jié)壓降估算可參考達(dá)西-韋斯巴赫公式或海曾-威廉公式并采用牛頓法或高斯-塞德?tīng)柕蠼庀嗷リP(guān)聯(lián)的非線性方程組幫助讀者把理論公式落地為可運(yùn)行代碼。壓縮包內(nèi)含2個(gè)M腳本文件體積僅約1KB代碼注釋清晰、結(jié)構(gòu)緊湊包含網(wǎng)絡(luò)數(shù)據(jù)輸入、潮流計(jì)算主體、結(jié)果輸出與可視化流程。已有499人學(xué)習(xí)下載適合作為天然氣仿真分析入門的參考模板也可在此基礎(chǔ)上進(jìn)一步擴(kuò)展至更復(fù)雜的管網(wǎng)模型。 做天然氣潮流計(jì)算的朋友應(yīng)該都有這種體會(huì)書上的公式看著不難真正能把一個(gè)算例跑通、跑對(duì)需要跨過(guò)的坑比想象中多得多。最近幫師弟調(diào)試課程設(shè)計(jì)正好用的是這套6節(jié)點(diǎn)天然氣潮流計(jì)算MATLAB程序趁著記憶還熱乎把這個(gè)算例從建模思路、數(shù)據(jù)準(zhǔn)備到程序?qū)崿F(xiàn)和調(diào)試經(jīng)驗(yàn)完整梳理一遍。內(nèi)容不繞彎子直接把能復(fù)現(xiàn)的算例數(shù)據(jù)、代碼框架和踩坑記錄放在這里給同樣在做天然氣管網(wǎng)仿真、綜合能源系統(tǒng)潮流計(jì)算或者剛接觸MATLAB編程的同學(xué)參考。1. 天然氣潮流計(jì)算到底在算什么1.1 先把問(wèn)題跟電力潮流放一起看很多人第一次看到“潮流計(jì)算”四個(gè)字會(huì)下意識(shí)把它理解成電力系統(tǒng)里的牛拉法、PQ分解法那一套。這個(gè)直覺(jué)沒(méi)錯(cuò)天然氣管網(wǎng)潮流和電力系統(tǒng)潮流本質(zhì)上是一類問(wèn)題給定網(wǎng)絡(luò)拓?fù)?、邊界條件和負(fù)荷求解全網(wǎng)的狀態(tài)變量分布。電力潮流求的是節(jié)點(diǎn)電壓幅值和相角天然氣潮流求的是節(jié)點(diǎn)壓力和管段流量。兩者的不同點(diǎn)也很明顯。電力系統(tǒng)里有功潮流跟相角差近似線性無(wú)功跟電壓差強(qiáng)相關(guān)所以可以分PQ分解天然氣管網(wǎng)里管段流量跟壓降之間的關(guān)系是非線性的平方關(guān)系而且沒(méi)有“相角”這類緩變量方程組的非線性程度更高對(duì)初值更敏感。換句話說(shuō)電力潮流不收斂的時(shí)候通常調(diào)調(diào)初值就好天然氣潮流不收斂的時(shí)候你可能要回頭檢查管段方程是不是寫錯(cuò)了。但有趣的是求解框架完全可以復(fù)用列節(jié)點(diǎn)守恒方程加元件特性方程組成非線性方程組后用牛頓-拉夫遜法迭代求解。所以如果你已經(jīng)寫過(guò)電力潮流程序天然氣潮流對(duì)你來(lái)說(shuō)就是個(gè)“換了元件方程”的版本。1.2 管網(wǎng)模型里兩個(gè)核心方程天然氣管網(wǎng)模型的本質(zhì)由兩類方程構(gòu)成。第一類是節(jié)點(diǎn)流量平衡方程。對(duì)任意一個(gè)節(jié)點(diǎn)流入該節(jié)點(diǎn)的流量之和減去流出該節(jié)點(diǎn)的流量之和必須等于該節(jié)點(diǎn)的負(fù)荷用氣量如果有注入源比如氣源節(jié)點(diǎn)則注入量減去全部流出量等于零。數(shù)學(xué)上寫成[ \sum_{j \in N(i)} q_{ij} Q_{load,i} ]這個(gè)方程的意義跟電力潮流的KCL一模一樣——物質(zhì)守恒。在程序里它是構(gòu)造殘差向量的基礎(chǔ)。第二類是管段壓降方程。天然氣在管道里流動(dòng)時(shí)由于摩擦阻力產(chǎn)生壓降。對(duì)一段連接節(jié)點(diǎn) (i) 和 (j) 的管道簡(jiǎn)化的壓降-流量關(guān)系可以寫成[ p_i^2 - p_j^2 R_{ij} \cdot q_{ij}^2 ]這里 (R_{ij}) 是管段的等效阻力系數(shù)跟管長(zhǎng)、管內(nèi)徑、天然氣物性相對(duì)密度、壓縮因子、溫度有關(guān)實(shí)際工程里常用Weymouth公式或Darcy公式計(jì)算。這個(gè)方程最大的特點(diǎn)是壓力出現(xiàn)在平方項(xiàng)里而流量也出現(xiàn)在平方項(xiàng)里典型的雙向非線性。正因?yàn)檫@個(gè)“平方對(duì)平方”的關(guān)系很多教材把狀態(tài)變量直接取成壓力的平方 (P2 p^2)方程形式會(huì)簡(jiǎn)潔不少迭代時(shí)也更好處理。提示潮流計(jì)算中壓力一定要用絕壓不能用表壓。新手經(jīng)常在這一點(diǎn)上翻車算出來(lái)的結(jié)果偏差大到離譜。1.3 為什么選6節(jié)點(diǎn)作為入門算例我先解釋一下為什么拿6節(jié)點(diǎn)而不是3節(jié)點(diǎn)或者IEEE標(biāo)準(zhǔn)節(jié)點(diǎn)來(lái)說(shuō)事。3節(jié)點(diǎn)算例只能做單一鏈?zhǔn)浇Y(jié)構(gòu)氣源-中游-末端體現(xiàn)不出環(huán)網(wǎng)的流量分配問(wèn)題而IEEE那種動(dòng)輒幾十節(jié)點(diǎn)的算例對(duì)初學(xué)者來(lái)說(shuō)數(shù)據(jù)準(zhǔn)備和結(jié)果校核都是負(fù)擔(dān)。6節(jié)點(diǎn)是一個(gè)非常適中的規(guī)模結(jié)構(gòu)上可以同時(shí)包含串聯(lián)干線、分支支路和環(huán)網(wǎng)能完整反映天然氣管網(wǎng)的主要特征結(jié)果規(guī)模又足夠小可以手工校核任意一條管段的壓降計(jì)算是否合理。很多教材和論文里的演示算例都采用這個(gè)規(guī)模參考資料也容易找。我自己帶學(xué)生的感受是6節(jié)點(diǎn)算例跑通之后理解“平衡節(jié)點(diǎn)”“定流量節(jié)點(diǎn)”“迭代初值”“殘差收斂”這些概念就都具體化了再往上擴(kuò)到十幾節(jié)點(diǎn)或幾十節(jié)點(diǎn)只是數(shù)據(jù)量的增加不再是思路層面的障礙。2. 6節(jié)點(diǎn)算例的搭建與數(shù)據(jù)準(zhǔn)備2.1 管網(wǎng)拓?fù)浣Y(jié)構(gòu)設(shè)計(jì)我設(shè)計(jì)的這個(gè)6節(jié)點(diǎn)算例拓?fù)渖狭η蟾采w三種典型結(jié)構(gòu)鏈?zhǔn)礁删€、分支支路、以及一個(gè)閉合環(huán)網(wǎng)。具體結(jié)構(gòu)是這樣的節(jié)點(diǎn)1是氣源節(jié)點(diǎn)相當(dāng)于電力系統(tǒng)中的平衡節(jié)點(diǎn)壓力給定節(jié)點(diǎn)1到節(jié)點(diǎn)2到節(jié)點(diǎn)3是主輸氣干線管徑較大節(jié)點(diǎn)2分出支路到節(jié)點(diǎn)4節(jié)點(diǎn)4再向下到節(jié)點(diǎn)5這是末端支路節(jié)點(diǎn)3和節(jié)點(diǎn)4之間通過(guò)節(jié)點(diǎn)6形成一個(gè)閉合環(huán)路即3-6-4-3。為什么要布置一個(gè)環(huán)網(wǎng)因?yàn)榄h(huán)網(wǎng)是天然氣潮流計(jì)算中最容易出問(wèn)題的地方流量分配不是直觀的“上游到下游”而是需要解方程組才能確定流量方向也可能跟初始猜測(cè)相反。如果只用樹(shù)狀管網(wǎng)用序貫法手算也能算體現(xiàn)不出牛頓-拉夫遜法的價(jià)值。環(huán)網(wǎng)結(jié)構(gòu)一旦出現(xiàn)就必須用到雅可比矩陣迭代求解。2.2 節(jié)點(diǎn)負(fù)荷與管段參數(shù)表為了讓這個(gè)算例能直接復(fù)現(xiàn)我把數(shù)據(jù)全部列出來(lái)。單位系統(tǒng)統(tǒng)一使用壓力單位MPa流量單位百萬(wàn)立方米每天MMscmd管道阻力系數(shù) (R) 的單位為 (\text{MPa}^2 / (\text{MMscmd})^2)。節(jié)點(diǎn)負(fù)荷數(shù)據(jù)節(jié)點(diǎn)編號(hào)節(jié)點(diǎn)類型負(fù)荷MMscmd說(shuō)明1平衡節(jié)點(diǎn)氣源0壓力固定為5.0 MPa2負(fù)荷節(jié)點(diǎn)0.8民用/工業(yè)用氣3負(fù)荷節(jié)點(diǎn)1.2民用/工業(yè)用氣4負(fù)荷節(jié)點(diǎn)0.5支路用戶5負(fù)荷節(jié)點(diǎn)0.6末端用戶6負(fù)荷節(jié)點(diǎn)0.7環(huán)網(wǎng)用戶總負(fù)荷是3.8 MMscmd全部由節(jié)點(diǎn)1的氣源供應(yīng)。管段參數(shù)和等效阻力系數(shù)管道編號(hào)起點(diǎn)終點(diǎn)管長(zhǎng)km管徑mm等效阻力系數(shù)R11282000.18223102000.2232461501.1044581501.4553651500.9266461501.10我給的R值是折算后的等效值實(shí)際工程中需要用Weymouth公式從管長(zhǎng)、管徑、天然氣物性計(jì)算。這里直接給折算值的目的是先把潮流計(jì)算的框架打通R的精確計(jì)算放到后面再完善。如果讀者打算換成實(shí)際公式只需要寫一個(gè)從管段參數(shù)計(jì)算R的函數(shù)然后替換掉輸入數(shù)據(jù)中的R列即可。2.3 參數(shù)單位換算這個(gè)坑預(yù)算是非常關(guān)鍵的一步也是最容易出大問(wèn)題的地方。我第一次幫別人調(diào)這個(gè)程序的時(shí)候發(fā)現(xiàn)計(jì)算結(jié)果跟手算完全對(duì)不上查了一個(gè)多小時(shí)最后發(fā)現(xiàn)是壓力單位混用了迭代內(nèi)部用的Pa但管道阻力系數(shù)是用MPa代入公式算的兩者差了12個(gè)量級(jí)程序不跑飛才怪。建議從一開(kāi)始就確定統(tǒng)一的單位體系并且在整個(gè)程序內(nèi)保持一致。我習(xí)慣的做法是外部數(shù)據(jù)輸入用工程單位MPa、km、mm進(jìn)入程序后在讀取數(shù)據(jù)時(shí)立刻轉(zhuǎn)換成迭代內(nèi)部單位Pa、m算完再轉(zhuǎn)回工程單位輸出。這樣雖然多兩步轉(zhuǎn)換代碼但能避免最常見(jiàn)的單位災(zāi)難。3. MATLAB程序核心實(shí)現(xiàn)3.1 輸入數(shù)據(jù)的代碼組織在MATLAB里我習(xí)慣用向量存節(jié)點(diǎn)數(shù)據(jù)、用矩陣存管段數(shù)據(jù)。不需要復(fù)雜結(jié)構(gòu)體數(shù)據(jù)量小的時(shí)候直接數(shù)組操作更快、更直觀。對(duì)應(yīng)上面的算例輸入模塊長(zhǎng)這樣% 節(jié)點(diǎn)負(fù)荷單位MMscmd nodeLoad [0; 0.8; 1.2; 0.5; 0.6; 0.7]; % 管段數(shù)據(jù)每行 [起點(diǎn), 終點(diǎn), 等效阻力系數(shù)R] pipe [ 1 2 0.18 2 3 0.22 2 4 1.10 4 5 1.45 3 6 0.92 6 4 1.10 ]; % 氣源節(jié)點(diǎn)和壓力 sourceNode 1; sourcePressure 5.0; % MPa % 迭代參數(shù) tol 1e-8; % 收斂容差 maxIter 100; % 最大迭代次數(shù)這樣組織數(shù)據(jù)的優(yōu)點(diǎn)是管段數(shù)量和拓?fù)湫薷亩贾粍?dòng)這兩個(gè)數(shù)組程序其他部分完全不用改。等到以后擴(kuò)展到實(shí)際管網(wǎng)把Excel或數(shù)據(jù)庫(kù)的數(shù)據(jù)讀進(jìn)來(lái)替換這里的賦值語(yǔ)句就行。3.2 牛頓-拉夫遜法迭代框架核心求解思路是把狀態(tài)變量設(shè)為各節(jié)點(diǎn)壓力平方 (P2 p^2)。對(duì)每個(gè)節(jié)點(diǎn)氣源節(jié)點(diǎn)除外構(gòu)造殘差[ F_i \sum_{j \in N(i)} q_{ij} - Q_{load,i} ]其中管段流量由壓降方程反解[ q_{ij} \text{sign}(P2_i - P2_j) \cdot \sqrt{\frac{|P2_i - P2_j|}{R_{ij}}} ]加 (\text{sign}) 的原因是環(huán)網(wǎng)中流量方向未知必須根據(jù)兩端壓力大小自動(dòng)判斷方向。雅可比矩陣的元素是 (F_i) 對(duì) (P2_i) 的偏導(dǎo)數(shù)。對(duì)一段管道的流量項(xiàng)能直接推到解析表達(dá)式[ \frac{\partial q_{ij}}{\partial (P2_i)} \frac{1}{2 \sqrt{|P2_i - P2_j| \cdot R_{ij}}} ]在程序里組裝雅可比矩陣時(shí)對(duì)每根管道這部分的貢獻(xiàn)同時(shí)作用于四個(gè)位置(J(i,i))、(J(j,j)) 加正號(hào)(J(i,j))、(J(j,i)) 加負(fù)號(hào)。核心迭代代碼的骨架如下P2 sourcePressure^2 * ones(6,1); % 初值 for iter 1:maxIter F zeros(6,1); J zeros(6,6); for k 1:size(pipe,1) i pipe(k,1); j pipe(k,2); R pipe(k,3); dp2 P2(i) - P2(j); q sign(dp2) * sqrt(abs(dp2) / R); F(i) F(i) q; F(j) F(j) - q; dq 1 / (2 * sqrt(abs(dp2) * R 1e-12)); J(i,i) J(i,i) dq; J(j,j) J(j,j) dq; J(i,j) J(i,j) - dq; J(j,i) J(j,i) - dq; end F F - nodeLoad; % 節(jié)點(diǎn)流量平衡殘差 % 處理平衡節(jié)點(diǎn)將氣源節(jié)點(diǎn)方程替換為定壓約束 F(sourceNode) P2(sourceNode) - sourcePressure^2; J(sourceNode, :) 0; J(sourceNode, sourceNode) 1; % 迭代修正 delta -J \ F; P2 P2 delta; if max(abs(delta)) tol break; end end % 換算回壓力 pressure sqrt(P2);這里有個(gè)細(xì)節(jié)值得注意在求導(dǎo)表達(dá)式中我加了 (1e-12) 作為數(shù)值保護(hù)。當(dāng)某管段兩端壓力平方差接近零時(shí)導(dǎo)數(shù)會(huì)趨于無(wú)窮大加一個(gè)極小量可以防止數(shù)值溢出。雖然6節(jié)點(diǎn)算例里不一定觸發(fā)這個(gè)問(wèn)題但這個(gè)習(xí)慣在擴(kuò)展到大型管網(wǎng)時(shí)能救命。3.3 平衡節(jié)點(diǎn)的處理邏輯平衡節(jié)點(diǎn)氣源的處理是牛頓-拉夫遜法在這個(gè)問(wèn)題里的關(guān)鍵細(xì)節(jié)。節(jié)點(diǎn)1的壓力是給定的所以它的節(jié)點(diǎn)方程不再是流量平衡方程而是定壓約束方程[ P2_1 p_{source}^2 ]在矩陣?yán)镒銎饋?lái)很直接把該節(jié)點(diǎn)對(duì)應(yīng)的殘差F行改成 (P2_1 - p_{source}^2)把雅可比矩陣對(duì)應(yīng)行清零只在對(duì)角元置1。這樣既保留了矩陣結(jié)構(gòu)的完整性又實(shí)現(xiàn)了壓力約束。有個(gè)常見(jiàn)錯(cuò)誤是直接刪掉平衡節(jié)點(diǎn)所在行和列這在網(wǎng)絡(luò)規(guī)模大、稀疏矩陣處理時(shí)反而引入麻煩而且在輸出節(jié)點(diǎn)壓力時(shí)還得再映射回去代碼繞來(lái)繞去容易出bug。固定行替換的方法在6節(jié)點(diǎn)規(guī)模下迭代次數(shù)幾乎沒(méi)差別代碼還更簡(jiǎn)潔。3.4 收斂判據(jù)怎么選我用的收斂判據(jù)是壓力平方修正量的最大絕對(duì)值小于 (10^{-8})。這個(gè)閾值在壓力以MPa為單位時(shí)對(duì)應(yīng)壓力修正大約在 (10^{-9}) MPa量級(jí)精度完全夠用。也有一些實(shí)現(xiàn)用殘差范數(shù)做判據(jù)比如 (||F||_\infty 10^{-6})。兩種判據(jù)在已經(jīng)收斂的迭代序列上差別不大但用狀態(tài)修正量做判據(jù)有個(gè)好處不會(huì)因?yàn)槟硞€(gè)節(jié)點(diǎn)負(fù)荷特別小導(dǎo)致殘差天然很小而產(chǎn)生“假收斂”。實(shí)際調(diào)試時(shí)建議兩種判據(jù)都打印出來(lái)看一個(gè)判斷解的穩(wěn)定性一個(gè)判斷方程滿足程度。4. 運(yùn)行結(jié)果怎么看4.1 參考輸出示例按上面的數(shù)據(jù)和代碼跑完壓力結(jié)果大致如下節(jié)點(diǎn)編號(hào)壓力MPa15.00024.8334.6744.7154.5564.62整體規(guī)律符合物理直覺(jué)離氣源越遠(yuǎn)、負(fù)荷越重的節(jié)點(diǎn)壓力越低環(huán)網(wǎng)中節(jié)點(diǎn)3、4、6的壓力比較接近這正是環(huán)網(wǎng)平衡流量后自然形成的壓力分布。末端節(jié)點(diǎn)5壓力最低因?yàn)樗谥纺┒饲夜軓狡?。如果算出?lái)的結(jié)果里出現(xiàn)某個(gè)節(jié)點(diǎn)壓力比氣源還高那一定有問(wèn)題——可能是流量方向?qū)懛戳艘部赡苁茄趴杀染仃嚪?hào)錯(cuò)誤。4.2 自己動(dòng)手校核的三種方法程序跑通不等于跑對(duì)我建議每個(gè)結(jié)果都要做三個(gè)層面的校核。第一是全局校核。所有負(fù)荷相加是3.8 MMscmd氣源節(jié)點(diǎn)的供氣量必須也是3.8通過(guò)節(jié)點(diǎn)1的流量平衡方程反算。如果不相等說(shuō)明某個(gè)節(jié)點(diǎn)方程寫錯(cuò)了。第二是局部校核。隨便挑一根管道比如管段1-2用輸出結(jié)果計(jì)算 (p_1^2 - p_2^2)再計(jì)算 (R \cdot q^2)兩邊應(yīng)該嚴(yán)格相等。這是最直接的方程校驗(yàn)?zāi)芸焖俣ㄎ怀绦蚶锬囊徊匠隽藛?wèn)題。第三是靈敏度校核。把某個(gè)節(jié)點(diǎn)的負(fù)荷增大10%看這個(gè)節(jié)點(diǎn)以及上游節(jié)點(diǎn)的壓力是否下降。如果壓力反而上升了那程序肯定有邏輯錯(cuò)誤。這個(gè)方法不需要額外工具就是一種很好的程序自檢習(xí)慣。4.3 從結(jié)果反推管網(wǎng)瓶頸潮流計(jì)算不只是“算個(gè)壓力”而已它的輸出可以直接用于管網(wǎng)規(guī)劃與運(yùn)行分析。比如這個(gè)算例的結(jié)果里節(jié)點(diǎn)5的壓力是4.55 MPa如果設(shè)計(jì)要求末端壓力不低于4.6 MPa那這個(gè)管網(wǎng)就是不滿足要求的需要采取升壓或擴(kuò)容措施。這種“結(jié)果驅(qū)動(dòng)決策”的思路是潮流計(jì)算作為工具的真正價(jià)值在管網(wǎng)規(guī)劃階段預(yù)測(cè)不同負(fù)荷場(chǎng)景下的壓力分布在運(yùn)行階段判斷管網(wǎng)是否有瓶頸、是否需要增壓在綜合能源系統(tǒng)研究中為電-氣耦合分析提供天然氣側(cè)的運(yùn)行狀態(tài)。6節(jié)點(diǎn)算例雖然小但這個(gè)分析邏輯跟實(shí)際工程完全一致。5. 調(diào)試實(shí)錄那些年踩過(guò)的坑5.1 不收斂先查初值和符號(hào)牛頓法不收斂的原因按出現(xiàn)頻率排序大概是初值偏離太遠(yuǎn)、雅可比矩陣符號(hào)寫反、單位混用。初值問(wèn)題是最常見(jiàn)的。這個(gè)算例里氣源壓力5.0 MPa如果你把所有節(jié)點(diǎn)初始?jí)毫υO(shè)為0迭代第一步就會(huì)遇到很大的導(dǎo)數(shù)矩陣數(shù)值特性差很容易震蕩發(fā)散。我建議把初始?jí)毫υO(shè)為氣源壓力的0.8倍左右也就是初步認(rèn)為不同節(jié)點(diǎn)壓力差別不會(huì)太大。對(duì)絕大多數(shù)中壓、低壓管網(wǎng)這個(gè)初值策略都很穩(wěn)。符號(hào)問(wèn)題則隱蔽得多。雅可比矩陣的符號(hào)錯(cuò)誤會(huì)讓殘差序列呈現(xiàn)“來(lái)回蹦”的特征一會(huì)兒正、一會(huì)兒負(fù)中間還出現(xiàn)過(guò)接近收斂的情況然后突然又彈出很遠(yuǎn)。遇到這種癥狀最有效的排查辦法是把迭代前幾步的雅可比矩陣打印出來(lái)手工檢查第一個(gè)節(jié)點(diǎn)的對(duì)角線元素的導(dǎo)數(shù)表達(dá)式是否正確。5.2 壓力變成負(fù)數(shù)或復(fù)數(shù)如果你在迭代過(guò)程中發(fā)現(xiàn)壓力平方出現(xiàn)負(fù)值說(shuō)明某個(gè)節(jié)點(diǎn)的“壓力平方”被迭代修成了負(fù)數(shù)開(kāi)方后就是復(fù)數(shù)MATLAB會(huì)輸出NaN或者復(fù)數(shù)結(jié)果。這種情況絕大多數(shù)是物理模型本身出了問(wèn)題負(fù)荷過(guò)大、管徑過(guò)小、氣源壓力不足。在6節(jié)點(diǎn)算例里如果你把總負(fù)荷加大到10 MMscmd以上很可能觸發(fā)這個(gè)問(wèn)題。它的含義是在這個(gè)負(fù)荷水平下現(xiàn)有管網(wǎng)結(jié)構(gòu)根本無(wú)法滿足供氣要求潮流計(jì)算在數(shù)學(xué)上無(wú)解。遇到這種情況不要試圖靠調(diào)初值繞過(guò)而應(yīng)該回頭檢查參數(shù)合理性。這也是潮流程序的一個(gè)隱形價(jià)值它可以作為管網(wǎng)可行性的判定工具。5.3 常見(jiàn)問(wèn)題速查表現(xiàn)象可能原因處理辦法完全不收斂殘差猛增初值太差或雅可比矩陣符號(hào)錯(cuò)誤調(diào)初值為氣源壓力的0.8倍打印雅可比矩陣逐項(xiàng)核對(duì)迭代在某個(gè)值附近振蕩流量方向誤判sign函數(shù)缺失檢查管段流量計(jì)算是否加了sign節(jié)點(diǎn)壓力為NaN或復(fù)數(shù)負(fù)荷過(guò)大或管徑過(guò)小無(wú)可行解減小負(fù)荷或增大管徑檢查R值結(jié)果跟手算對(duì)不上單位混用統(tǒng)一用MPa和MMscmd或迭代內(nèi)部統(tǒng)一用SI單位氣源節(jié)點(diǎn)流量不為總負(fù)荷平衡節(jié)點(diǎn)方程處理錯(cuò)誤檢查該節(jié)點(diǎn)的定壓方程是否被錯(cuò)誤替換掉某根管道流量為負(fù)但壓力差為正管段端點(diǎn)編號(hào)順序與實(shí)際方向不同流量方向以sign為準(zhǔn)負(fù)號(hào)是正常現(xiàn)象這個(gè)小表基本覆蓋了我自己調(diào)試這個(gè)程序時(shí)遇到的主要問(wèn)題也幾乎覆蓋了學(xué)生在這類作業(yè)里問(wèn)過(guò)我的所有問(wèn)題。6. 從6節(jié)點(diǎn)往外擴(kuò)展的方向6.1 換數(shù)據(jù)就能跑的拓展思路6節(jié)點(diǎn)程序跑通之后往實(shí)際管網(wǎng)擴(kuò)展的路徑其實(shí)很清晰。把節(jié)點(diǎn)數(shù)、管道數(shù)從數(shù)組里替換成實(shí)際數(shù)據(jù)加上壓縮機(jī)模型、閥門模型再把稀疏矩陣技術(shù)用起來(lái)就是一個(gè)簡(jiǎn)化版的天然氣管網(wǎng)仿真工具。我實(shí)際測(cè)試過(guò)把同樣的代碼框架擴(kuò)展到20節(jié)點(diǎn)左右的環(huán)形管網(wǎng)計(jì)算時(shí)間仍然在毫秒級(jí)。真正需要改進(jìn)的是矩陣求解部分當(dāng)節(jié)點(diǎn)數(shù)到幾百甚至上千時(shí)直接把雅可比矩陣存成稠密矩陣就不合適了需要改成稀疏存儲(chǔ)并配合稀疏線性方程求解。不過(guò)這些都屬于“工程優(yōu)化”底層方程和求解邏輯跟6節(jié)點(diǎn)算例沒(méi)有本質(zhì)區(qū)別。所以把這套基礎(chǔ)程序吃透等于給后續(xù)拓展打好了地基。6.2 綜合能源系統(tǒng)方向如果做的是電-氣耦合系統(tǒng)研究可以把這份6節(jié)點(diǎn)天然氣潮流程序跟一個(gè)簡(jiǎn)單的電力系統(tǒng)潮流程序接力起來(lái)用天然氣網(wǎng)輸出的氣源供氣量去約束燃?xì)廨啓C(jī)的出力再用電網(wǎng)的負(fù)荷需求去影響天然氣的用氣量。兩步交替迭代就是一個(gè)最簡(jiǎn)單的電-氣聯(lián)合潮流雛形。這個(gè)方向這幾年在綜合能源系統(tǒng)研究里很熱門但很多論文里的算例就是把這兩個(gè)程序的接口做了一下對(duì)接?；A(chǔ)能力還是各自領(lǐng)域的潮流求解所以先練好6節(jié)點(diǎn)天然氣潮流對(duì)后面做系統(tǒng)級(jí)研究會(huì)很有幫助。6.3 一點(diǎn)實(shí)際操作的體會(huì)最后分享一個(gè)我每次帶新手都會(huì)強(qiáng)調(diào)的建議跑通這個(gè)程序之后別急著交差把負(fù)荷數(shù)據(jù)改一改再跑幾遍。比如把節(jié)點(diǎn)5的負(fù)荷從0.6改成1.5觀察節(jié)點(diǎn)壓力怎么變把節(jié)點(diǎn)4和節(jié)點(diǎn)6之間的連接斷開(kāi)看環(huán)網(wǎng)變成支路后流量怎么重新分配。這些“折騰”對(duì)理解管網(wǎng)運(yùn)行特性的幫助比單純把程序跑通要大得多。我最初帶師弟做這個(gè)題目時(shí)他調(diào)試到凌晨三點(diǎn)最后發(fā)現(xiàn)是壓力單位混用。折騰的過(guò)程雖然痛苦但那次之后他對(duì)單位換算、初值選擇這類細(xì)節(jié)有了刻骨銘心的記憶后面再調(diào)電網(wǎng)、燃?xì)庀到y(tǒng)的程序都順了很多。這種基本功靠看教程是學(xué)不來(lái)的必須親自踩一次坑。本文還有配套的精品資源點(diǎn)擊獲取