數值積分范文

導語：如何才能寫好一篇數值積分，這就需要搜集整理更多的資料和文獻，歡迎閱讀由公務員之家整理的十篇范文，供你借鑒。

篇1

關鍵詞：高斯數值積分方法；積分變量；積分區間；非圓弧拱；多拱梁法

1問題的提出

高斯數值積分方法是一種節點很少、精確度很高的方法，它的特點是節點不等距，計算精度很高，一般利用正交多項式的有關關系式來確定其節點位置和系數.當節點為n時，其代數準確度可達2n-1次.如節點數為3,則求積公式對于任意5(=2×3-1)次多項式都是準確的，這樣的精度完全可應用于拱壩程序中的拱段及梁段的計算.筆者曾在圓弧拱多拱梁法程序中，廣泛采用3節點高斯數值積分，取得成功[1].但將這一方法推廣應用于非圓弧拱多拱梁法程序，卻有一定的難度.

由計算數學可知[2]，一般定積分式與高斯積分式的變換形式為：

式中:l=(b-a)/2，為積分限變換系數；n為高斯節點數；ξi為高斯積分節點坐標；gi為與ξi對應的高斯積分系數.

如所周知，拱圈形、載常數計算公式為：

(1)

(2)

式中:S為弧長.

當采用高斯積分公式且為等截面圓拱時，上式相應改為：

(3)

(4)

式中:φ為與弧長S對應的中心角；r為中心半徑；E為壩體彈模；Ii,MLi為與高斯節點i對應的截面慣性矩及靜定力矩.從上式可見，由于等截面圓拱r為定值，存在簡單的dS=rdφ,S=rφ的關系，所以積分變量由S變為φ，積分區間由弧段0～S變為中心角0～φ.在等截面圓拱計算中，許多參數(如坐標及靜定力系等)可直接由中心角用顯式求得.積分變量由S改為φ后，可大大簡化計算.

然而，對于非圓弧拱，問題要復雜得多.因為非圓弧拱的曲率半徑和曲率中心處處都在變化，往往不能用簡單的顯式來表示某一拱段中心角與弧長的關系.并且某些曲線計算弧長也很麻煩，甚至不易用顯式求得.如何采用適宜的積分變量和積分區間，是迫切需要解決的問題.

2不同曲線積分變量區間的選取

2.1基本資料5種非圓弧拱示意見圖1.幾種非圓弧拱曲線方程及有關公式見表1.

圖15種非圓弧拱平面示意

表1幾種非圓弧拱曲線公式

類型拋物線對數螺旋線橢圓雙曲線

曲線方程y=x2/2Rρ=ρ0eaψ

a=cosβ

曲率半徑r=(R2+x2)3/2/R2r=Reaψ

r=a2b2x

[(a-y)2/a4+x2/b4]3/2r=a2b2x

[(a+y)2/a4+x2/b4]3/2

弧長S=1/2{xA/cosψA+RLn[(sinψA+1)/cosψA]}S=R/a(eaψA-1)

dSReaψdψ

備注R為拱冠曲率半徑β為切線角,R為拱冠曲率半徑a為長軸之半,b為短軸之半a為實軸之半,b為虛軸之半

2.2弧的微分公式表1所列弧的微分(dS)算式，除了對數螺旋線稍簡單外，其余3種曲線的算式都比較復雜.若從曲線的一般性質來看，當曲線方程y=f(x)時，則弧的微分為：

式中:y′=tgφ，為函數y=f(x)在點x的導數，φ為過點(x,y)的切線與x軸的交角，不難看出，此角與該點的中心角相等.將y′=tgφ代入上式：

所以dS=dx/cosψ(5)

將式(1)作為拋物線、橢圓、雙曲線3種曲線弧的微分一般公式，比表1中所列dS算式要簡捷得多.式(1)也適用于其他一階導數存在的任何曲線.

2.3積分變量區間的選取根據各類曲線的性質，選取3種積分變量區間，分述于后.

(1)對于拋物線、橢圓、雙曲線3種曲線，采用式(1)所列弧的微分一般公式.此時積分變量為x，積分區間為與S相應的x變化區間.高斯積分時各項均應乘以1/cosφi,φi為與節點i對應的中心角.如A1算式改為:

(6)

(2)對于對數螺旋線，采用表1所列算式：dS=Reaφdφ.此時積分變量為φ，積分區間為與S相應的φ的變化區間，各高斯積分項i均應乘以Reaφi，如A1算式為：

(7)

上式除了高斯積分項乘以eaφi外，其余與圓弧拱相似.

圖2五心拱平面示意

(3)對于五心拱，其中弧段為圓弧，算法與等截面圓拱相同，當其邊弧段為變截面時，則中心拱弧線為非圓弧曲線，且不便用顯式表達.仔細考察該段曲線，發現它與以R3=(RM+RD)/2為半徑的圓弧很相近(RM,RD分別為邊弧外、內半徑)，此圓弧中心O3在邊弧起始截面外弧中心O1與內弧中心O2聯線中點處，如圖2所示.

邊弧段計算時，積分變量為φ，積分區間為與邊弧S相應的φ的變化區間φ3，如A1算式為：

(8)

上式形式上與圓弧拱一樣.

3算例

以上述3種積分變量區間的選取方式，計算各類曲線的半拱弧長，舉例于下.

設采用3節點高斯數值積分方法，節點坐標及高斯積分系數列于表2.如以中心角φ為積分變量，以Δφ為積分區間，則與節點i相對應的φi=Δφ(1+ξi)/2.又如以水平坐標x為積分變量，以Δx為積分區間，則與節點i相對應的xi=Δx(1+ξi)/2.

表2節點坐標及高斯積分系數

節點號i節點坐標ξi高斯積分系數gi

1-0.7745966910.555555582

200.888888896

30.745966910.555555582

3.1求拋物線、橢圓、雙曲線拱半拱弧長例1：設拋物線拱拱冠曲率半徑R=140m，拱端中心角φA=45.32°，拱端坐標xA=141.5726m，求半拱弧長S，可以采用兩種方法.一種是按表1所列S的公式直接計算，這是理論積分后的公式，是精確的.另一種是采用高斯數值積分方法，以x為積分變量，xA為積分區間.

(1)按理論公式

S=1/2{xA/cosψA+RLn[(sinψ\-A+1)/coxψA]}

算得S=162.921371m.

(2)按高斯數值積分方法，應有

算得gi/cosφi=2.301610646,S=(141.5726×2.301610646)/2=162.922485m.

該數值積分值與理論計算值162.921371m相比，僅相差0.001114m，相對誤差僅為6.8×10-6.

例2：設橢圓拱拱冠曲率半徑R=164.51m，拱端中心角φA=51°，坐標xA=141.5843m，長軸之半236.9m，短軸之半197.42m，求半拱弧長.例3：設雙曲線拱拱冠曲率半徑R=152.71m，拱端中心角φA=41°，坐標xA=141.58421m，實軸之半954.41m，虛軸之半381.77m，求半拱弧長.

上兩例因兩種曲線無理論積分公式直接用顯式計算弧長S，只能與表1中dS算式的高斯數值積分值相比.用dS算式直接數值積分時，積分變量、積分區間與前述方法一樣，仍為x及xA，但每一高斯積分項不用除以cosφi，而是乘以dx前的算式等，可見后一算法較繁.兩例成果S及比較見表3.

表3橢圓拱、雙曲線拱計算成果比較

類別例2橢圓拱例3雙曲線拱

原dS算式成果164.143661158.610046

dS=dx/cosψ成果164.143737158.610031

兩種算法差值/m0.0000760.000015

相對誤差4.63×10-79.457×10-8

由表3可見兩種算法成果非常接近.

3.2求對數螺旋線拱半拱弧長例4：設對數螺旋線拱拱冠曲率半徑R=154m，拱端中心角φA=46°，切線角β=52°，a=ctgβ=0.781285626,以中心角φ為積分變量，積分區間為拱端中心角φA，求半拱弧長S.

(1)按理論公式計算

S=R/a(eaψA-1)=171.9726625m

(2)按高斯數值積分計算

兩種算法成果差值為0.0000061m，相對誤差僅為3.547×10-8.

3.3用近似方法求五心拱邊例5：設五心拱半拱邊弧夾角10°，外半徑RM=290m，內半徑RD=193.514m，邊弧起點拱厚6.466m，拱端厚度8.466m，求半拱邊弧長S.

(1)用較精確的計算公式S=φ3×(RA+2×R3)/3.式中3=(RM+RD)/2;φ3為以O3為近似中心的邊弧夾角，以弧度計；RA為邊弧拱端點與O3聯線的長度，見圖2.算得RA=241.56728m，φ3=0.206893551弧度，R3=241.757m，從而算出S=50.00488034m.

(2)用近似計算公式S＝φ3×R3=50.01796429m.

兩種算法所得的差值為0.013m，相對誤差為2.6×10-4.

由以上成果可知，以R3為半徑，以φ3為中心角所得的圓弧與邊弧的中心弧很近似，因此在高斯數值積分計算時，可以采用較簡單的積分變量φ及積分區間φ3.

4結語

綜上可知，高斯數值積分應用于非圓弧拱時，需針對非圓弧拱曲線的性質和特點，選取適宜的積分變量與相應的積分區間，這樣常可收到事半功倍之效.如對于拋物線拱、橢圓拱、雙曲線拱，采用dS=dx/cosφ的通式，既避免了各種曲線的繁復計算，也便于程序規格化.對于對數螺旋線拱，則利用對數函數微分、積分都簡單以及極坐標方程弧的微分的特點，積分變量選取φ而不選取x，既大大簡化了計算，也很容易利用圓弧拱的算法稍加變換.對于五心拱，在控制誤差足夠小的前提下，邊弧線采用近似圓弧，大大簡化了計算.

將3節點高斯數值積分方法推廣應用于各種非圓弧拱壩多拱梁法程序，對于提高拱壩程序的計算速度和精度，有很大價值；對于將來推廣應用高斯數值積分于各類復雜結構的分析、計算，也有一定的啟發和借鑒作用.

參考文獻

篇2

關鍵詞：曲面面積 數值計算 數值微分 積分

中圖分類號：O172 文獻標識碼：A 文章編號：1672-3791（2014）03（c）-0238-02

二重積分的數值計算方法有很多，但是在實際應用中，曲面面積的很重要，而曲面面積計算的數值方法卻不多，目前還沒有找到一種高效、精確的計算其表面積方法。文獻[1][2]模型的建立是基于多網格化下小區域內曲面積近似等于平面面積，因此計算結果存在一定誤差，且計算精度不易分析。為了減小誤差，提高精度，我們建立利用積分中值定理和數值微分公式，建立一個新的計算表面積的數值計算公式―― “四點”插值算法。

1 單元構造和數值計算公式

已知曲面函數為，則考慮曲面在矩形區域內的表面積。對區域進行分割，首先考慮如圖1網格單元區域：

利用積分中值定理[3]

則

若，充分小，則由偏導數的連續性有：

，

于是

由三點數值微分公式[4]

，

于是

2 誤差估計

其中

由三點數值微分公式[4]

，其中

由二階泰勒展開：

，其中

于是

其中

同理

所以

3 復化公式

計算矩形區域內函數的表面積，在格網化區域計算表面積。首先對區域進行劃分，把目標區域劃分成個方格，則有：

，

取如圖2的方格，則在每個方格上應用表面積的近似計算公式，只需計算4個信息點。

4 算例分析

例：曲面函數在矩形區域內的表面積。

其表面積計算的精確值為：

在相同的分割網格下：

“四點”插值算法節點數：

三角形法需要的節點數：

數值計算結果如表1。

5 結語

通過實驗的matlab仿真，可知基于本文的方法求解曲面面積的算法誤差和傳統的“三角形法”誤差雖然都是，但本文方法的誤差是“三角形法”的，計算時間是“三角形法”的二十分之一。由此可以看出本算法需要信息點少，精度較好，運算速度快，具有較大的實用價值。

參考文獻

[1] 陳吉龍，武偉，劉洪斌.DEM在林地表面積計算中的應用研究[J].西南農業學報，2008，21（5）.

[2] 魏東，張秀程.基于遞歸算法的三維地形面積計算方法研究[J].沈陽：沈陽工業大學信息科學與工程學院，2007（3）.

篇3

關鍵詞 地層熱阻 熱響應測試地埋管換熱器

中圖分類號：TK172 文獻標識碼： A

Thermal resistance calculation on vertical ground heat exchanger

Li Jinghui1Wang Jiankui2 Lu Lin2 Fang Xugen1

1 Zhejiang Construction Division Building Energy Technology Co., Ltd.

2 Zhejiang Academy of Building Research & Design.ltd

Abstract: Using the data of grock-soil thermal response test, andaccordingto , this paper presents the numericalcalculationofthe vertical ground heat exchanger's thermal resistance calculation formula,andalso analyzeshow heat exchanger resistance affect the performance of heat exchanger. The calculation and analysis have a certain reference valueinguidingthe design and construction of ground heat exchanger.

Key word: ground heat resistance，heat responsetest，ground heat exchanger

0 引言

地層熱阻是決定土壤源熱泵合理設計和科學應用的核心因素，是影響地下埋管換熱器傳熱性能及土壤源熱泵系統節能與經濟性的重要原因。換熱器的熱阻工程實際計算及熱阻的影響對于設計人員設計及方案的優化有著重要的影響。

1 地埋管換熱器的熱阻計算

地埋管地源熱泵系統地埋管換熱器的設計計算主要有采用單位井深換熱量、專業軟件和依據《地源熱泵系統工程技術規范》[1]（以下簡稱《規范》）作為目前指導地源熱泵工程設計和施工最重要的規范，在其中第4.3.5條明確指出“豎直地埋管換熱器的設計也可按本規范附錄B的方法進行計算。

但是《規范》附錄B中對地層熱阻及短期連續脈沖負荷引起的附加熱阻的計算公式的表述存在一定的歧義，也不完整，在實際應用中會產生較大的誤差甚至錯誤[2]。本文以供熱工況下的換熱器的熱阻為例各熱阻值對進行計算分析。

1.1巖土熱響應測試數據

表1 巖土熱響應測試基礎數據表

項目 測試孔 項目 測試孔

鉆孔深度（m） 100 鉆孔直徑（mm） 135

埋管形式 單U型 埋管材質 PE管

埋管內徑（mm） 26 埋管外徑（mm） 32

鉆孔回填材料 原漿 細沙 主要地質結構 粉質粘土

巖土熱響應測試測試報告結果：埋管區域的平均綜合導熱系數為1.667 W/ m·℃，平均容積比熱為2.034×106 J/m3·℃。巖土體初始溫度20.1℃。

1.2值計算

公式1

式中：

——傳熱介質與U形管內壁的對流換熱熱阻，（m·℃/W）;

——U形管的內徑（m）；

——傳熱介質與U形管內壁的對流換熱系數[W/ m2·℃]。

換熱器溫差＜20℃，根據熱物性測試報告定性溫度按6℃，查表得=56.3×10-2W/ m·℃，ν=1.553×10-6m2/s ，=11.60 =1547×10-6N·s /m2

d：管內徑=26mm管內媒質流速按0.9 m/s計算；=0.9×0.026/1.553×10-6=15067＞104；=0.023×15067×11.6=4020；K==4020×0.563÷0.026=87045 W/ m2·℃；=1÷(3.14×0.026×87045)=1.41×10-4 m·℃/W

1.3 值計算

公式2

公式3

式中：

——U形管的管壁熱阻，（m·℃/W）;

——U形管導熱系數[w/(m·℃)]；

——U形管的外徑(m)；

——U形管的當量直徑(m)；對單u形管，n=2；對雙U形管，n=4。

=0.42~0.45 w/(m·℃)；=0.032 m； ——0.0452m

=0.0539 m·℃/W

1.4計算

公式4

式中：

——鉆孔灌漿回填材料的熱阻(m·℃/w)；

——灌漿材料導熱系數[w/(m·℃)]；

——鉆孔的直徑(m)。

=1.47w/(m·℃)， =0.135m，= 0.0809m·℃/w

1.5計算

對應于單個鉆孔：

公式5

公式6

式中：

——地層熱阻(m·℃/w)；

——指數積分公式；

——巖土體的平均導熱系數[w/(m·℃)]；

——巖土體的熱擴散率m2/s；

——鉆孔的半徑（m）；

——運行時間（S）；

根據文獻[2-4]對進行修正簡化為下式：

公式7

= 公式8

公式9

=1.667w/(m·℃)，=/cρ，cρ=2.034×106J/m3℃，則=0.8196×10-6m2/s，鉆孔的半徑=0.0675（m）,運行時間=68×24×3600=5.8752×106（S）；z=0.0006

==0.3712 m·℃/w

1.6計算

公式10

式中：

——短期連續脈沖負荷引起的附加熱阻(m·℃/W)；

——短期脈沖負荷連續運行的時間（s）。

根據文獻[2-4]修正簡化為下式：

公式11

=10×3600s

=0.1297 m·℃/w

2 熱阻影響分析

根據《規范》中公式B.0.2-4中供熱運行份額Fn取0.2。根據以上計算得各熱阻值及比例見下表，豎直埋管換熱器分項熱阻計算分析。

表2 豎直埋管換熱器分項熱阻計算分析表

項目 ×Fn ×(1-Fn)

熱阻值

（m·℃/w） 0.001 0.0539 0.0742 0.1038 0.1184

份額（%） 0.04 15.39 21.18 29.60 33.79

根據以上數據可知，單孔U型管回填材料的熱阻、地層熱阻和短期連續脈沖負荷引起的附加熱阻為構成換熱器熱阻的主要部分，U型管內壁的對流換熱熱阻很小，對換熱器熱阻影響幾乎可以忽略。

3 結論

3.1對豎直地埋管換熱器各熱阻值進行了詳細的計算，結合相關文獻對地層熱阻及短期連續脈沖負荷引起的附加熱阻進行了簡化計算。

3.2 計算結果表明，單孔U型管回填材料的熱阻、地層熱阻和短期連續脈沖負荷引起的附加熱阻為構成換熱器熱阻的主要部分，對于優化換熱器設計方案具有一定的參考價值。

參考文獻

[1] 中國建筑科學研究院.GB 50366-2005 地源熱泵系統工程技術規范 [S].北京：中國建筑工業出版社,2009版

[2] 雷建平,於仲義. 關于《地源熱泵系統工程技術規范》地層熱阻計算式的商榷及應用分析[J], 暖通空調,2009,39(6):27-30

[3] 余延順. 土壤耦合熱泵系統地下埋管換熱器傳熱模型的研究[J].暖通空調,2005,35(1):26-31

篇4

在計算機上運用數值分析解決實際問題的過程為：實際問題→數學模型→數值計算方法→程序設計→上級計算求出結果[1]。數值實驗的設計應當充分體現這一過程，同時也應當充分體現數學建模思想. 而目前地方高校在數值分析課程的教學與實驗中，普遍存在重理論、輕實踐、純數學十足的問題，部分即使避免這些問題，做到利用計算機進行可視化教學和算法編程實踐，但學生也缺乏分析問題、解決問題的能力，因為學生的算法與程序多是百度而來，缺乏思考. 因此，數值分析課程的數值實驗因當因時而異，因專業而異，緊跟時展，做到充分吸引學生的注意力，激發學生的興趣，調動學生積極性，讓學生積極主動的去想辦法解決問題，而非被動的接受。這樣的數值實驗才能有效培養學生應用能力。

數值分析的內容廣泛，包含插值擬合、數值積分與數值微分、數值代數、微分方程數值解法、非線性方程與方程組的數值解法，最優化等等。下面以插值擬合、微分方程數值解法、非線性方程數值解法、最優化的相關實際問題為例，研究數值實驗的設計。

2 結束語

通過以上案例發現：如果知道解決實際問題的數值分析方法，利用MATLAB能很方便的求出實際問題的結果。但是對于多數學生而言，由實際問題無法得到數學模型，更無法知道相應的數值分析方法。因此在學生有一定數值分析基礎后，才能引入這些實際問題的數值實驗；同時選取的數值實驗必須為學生深入研究預留了空間，因為在講解相關的知識背景，詳細分析問題，建立模型后，對學生進行的是分層次指導解決這些問題（即或利用MATLAB函數命令簡單編程計算，或利用MATLAB的專用工具箱計算，或設計算法流程利用MATLAB編程計算）；這樣每一個學生都能參與到實驗中來，每一個學生在實驗中都有收獲；最后撰寫實驗報告，闡明實驗的目的、要求、過程、收獲等.如此，通過這些數值實驗，學生既自己動手解決實際問題，培養了其綜合應用能力，又讓其充分體會到數值分析的魅力，為大學生數學建模競賽和后續專業學習打下堅實基礎.當然，隨著技術的進步，不同專業的數值分析應用實驗需要教師不停的探索和研究。

參考文獻：

[1] 徐翠薇，孫繩武.計算方法引論[M]. 3版.北京：高等教育出版社，2007.

[2] 呂國英.算法設計與分析[M]. 2版.北京：清華大學出版社，2009.

[3] 張德豐等. MATLAB數值分析[M]. 2版.北京：機械工業出版社，2012.

[4] 杜廷松.關于《數值分析》課程教學改革研究的綜述和思考[J]. 大學數學，2007，23（2）：8-15.

[5] 張濤，. 數值實驗在《數值分析》教學中的重要性分析[J]. 長江大學學報：自然科學版，2009，6（1）：371-372.

篇5

關鍵詞：質子交換膜燃料電池；數值模擬；兩相流；水傳遞

中圖分類號：U473.4文獻標文獻標志碼：A文獻標DOI：10.3969/j.issn.2095-1469.2012.06.03

質子交換膜燃料電池（Proton Exchange Membrane Fuel Cell, PEMFC）可以看作是一個正極不斷生成水的發電裝置。如果正極生成的水不能有效排出，就會造成正極水淹。一方面造成催化層氧氣供應不足，此時原本在高電位下進行的氧氣得電子產生還原反應，就會被低電位下進行的H+得電子還原反應取代，電池端電壓由正轉負，形成所謂“負差效應”；另一方面催化層活性區域減小，導致催化劑利用率降低。目前大部分工作放在系統設計和工程方面，這些策略常不可避免地會帶來嚴重的寄生功率損耗。對小功率便攜式PEMFC而言，由于所容許的空間和能量有限，不可能采用主動型排水裝置。對于大功率的燃料電池系統，雖然允許采用增加氣體流速，強制排除水分措施[1-5]，但是，這些方法只能排除雙極板流場中的水分，無法有效解決夾藏在多孔電極孔隙中的水，而PEMFC生成的第1滴水恰恰在氣體多孔電極的催化層中。

多孔電極內液態水量過多，液態水會淹沒電極減少電極有效面積，同時堵塞傳質通道嚴重阻礙氧氣傳輸，使陰極氧氣供應不足，濃差極化增大，從而導致電池性能大幅度下降[6-9]。正極水淹造成催化層氧氣供應不足的核心問題是：非極性的O2分子在極性水分子中的溶解度太低，在PEMFC工作溫度（60～80?℃）下，由于水的蒸汽壓增大，O2分壓減小，O2在水中的溶解度還會進一步減小。我們基于“相似相溶原理”，在電極微孔中預先加入對O2有很高溶解度的非極性硅油（DMS），代替孔隙作為氣體不可侵占的傳輸通道，確保氧氣的供應不會因為水的堵塞而中斷，而未被硅油占據的孔隙則作為反應生成水的固定排除通道，從而實現了氣體傳輸與水排出通道的有序分離[10]。圖1為傳統氣體電極（Conventional Gas Pt/C Electrode, CGE）和抗水淹電極（Anti-Flooding Pt/C Electrode, AFE）示意圖，圖2展示了CGE催化層內氧氣的通道為空氣和液態水，AFE催化層內氧氣的通道為空氣、DMS和液態水。AFE的優勢主要體現在以下幾方面：（1）氧氣在DMS內的溶解度是水中的10倍，DMS提供了輸送氧氣的專用通道。（2）DMS主要填充在孔徑在20?～70?nm之間的孔隙中，因而解決了最難克服的水淹問題，該尺度內的傳質是傳統的對流的方法難以解決的。（3）氣體多孔電極催化層微孔而非雙極板流場的水淹問題才是影響燃料電池性能的首要問題，而AFE直接解決了發生在多孔電極內的水淹問題。

本文在文獻[11]和[12]模型基礎上，針對抗水淹電極內DMS在多孔電極內的作用，發展一個三維兩相流的數學模型，用于研究抗水淹電極單電池，然后討論多孔電極內液態水飽和度分布和電流密度，分析抗水淹電極內DMS對燃料電池輸出特性的影響。

1 幾何模型與數學方程

圖3中試驗與數值計算都是基于圖3（a）所示的蛇形流道。圖3（b）的單電池由陽極極板、陽極流道、陽極擴散層、陽極催化層、質子交換膜、陰極催化層、陰極擴散層、陰極流道和陰極極板依次組合而成，圖3（c）為單電池的計算域，電極活性面積為20?mm×21?mm。三維數值計算域是根據流場板的流道外形來建立流道的實體幾何域作為流場計算的流體域。由于流場板為固體導電區域，只有電傳導和熱傳導的計算，不涉及多組分的流體傳輸和各種源項的計算，因此流場板的網格尺寸相對大些，而催化層是電化學反應發生的區域，涉及多種源項的計算和多組分流體的傳輸，所以對催化層網格進行了局部加密。

本文的數值模擬工作基于一個穩態、多組分、非等溫和三維的兩相流傳遞模型，模型考慮水的相變，模型簡化了液態水的生成和傳輸過程，將液態水飽和度分布描述為水蒸汽過飽和度的簡單函數，這樣就可以通過氣體擴散傳質方程得到的水蒸汽分壓力直接計算出液態水飽和度分布。液態水飽和度分布的結果將被用來修正反應氣體擴散系數和電化學反應速率方程，這樣模型就能夠模擬液態水生成堵塞傳質通道和淹沒電極有效面積的影響。

流場邊界條件：陰陽極入口設為質量入口邊界，先設定進口邊界上總的質量流量，然后再指定不同組分的摩爾分數；陰、陽極出口設為壓力出口，設定相對壓力為0?Pa。流場板和流道的交界面以及流場板和擴散層的交界面設為壁面邊界條件；流道和擴散層的交界面、擴散層和催化層的交界面以及催化層和質子交換膜的交界面均設為內部連續性邊界，其余邊界為壁面邊界條件。

溫度場邊界條件：模型外部邊界面及流道入口和出口均設定為燃料電池的工作溫度；內部交界面設為連續性邊界條件。

電勢場邊界條件：陰極流場板電位邊界設為電池輸出電壓，陽極流場板電位邊界設為0?V，內部各交界面設為連續性邊界條件，其余邊界設為絕緣邊界條件。

數學模型應用了以下假設：（1）氣體混合物為理想氣體。（2）電池內流體的流動為層流。（3）多孔介質為各向同性，即在同一層中具有相同的特征參數，如滲透率、組分擴散系數和黏性系數等。（4）流道

足夠長，即流動在出口處達到充分發展狀況。（5）催化層內催化劑和DMS含量比例恒定且均勻。

質量守恒方程：

，

x、y和z 3個方向的動量守恒方程：

，

組分守恒方程：

，

溫度控制方程：

，

式中：為孔隙率；ρ為密度；為速度矢量；Sm是由于電化學反應發生而導致的混合氣體質量源項；ui為x、y和z 3個方向的速度；為動力學粘度；p為流體壓力；Su為動量源項；Ci為組分體積濃度； Dm,i為組分擴散系數；Si是由于電化學反應發生而導致的各組分質量源項；Cp為等壓比熱容；T為熱力學溫度；keff為有效導熱系數；ST為能量源項。各方程具體源項和參數計算方法參見文獻[11]和[12]。

以下是針對抗水淹電極建立的局部有效擴散系數修正方程。

傳統多孔電極里組分氣體的有效擴散系數由孔隙率和水的飽和度s共同決定，水的飽和度定義為多孔電極內局部液態水的體積與孔隙體積的比值。

，

，

式中：Dm,i為組分氣體的擴散系數；VWL為局部液態水的體積；Vpore為孔隙體積。

有效介質理論[13]用于估計有效擴散系數，對于重疊的球狀顆粒

．

多孔電極液態水的飽和度，通過減小擴散面積和增加擴散路徑的彎曲阻礙了氣體的擴散

g(s)=(1-s)x ．

如果孔隙的結構不隨著水飽和程度改變，根據有效介質理論，x=1.5。Mezedur等[14]針對他們模型中數值計算得到x=0.71。Martys[15]用數值方法研究了多孔介質中液體以球狀形式存在的影響問題，通過網格計算得到液態水的分布，結果表明x=2。但對于接觸角為90°的情況，x遠遠大于2。因此多數理論分析中x一般大于2[16]。所以x的值是與具體的多孔介質的結構有關，是不同燃料電池多孔電極的物性參數，本文通過數值仿真結合試驗結果推算出x值為2。

針對抗水淹多孔電極，多孔電極里組分氣體的有效擴散系數由孔隙率、水的飽和度S和DMS在多孔電極里的體積分數Soil共同決定。假設DMS均勻分布在催化層內，在式（8）基礎上發展為式（9）。

，

，

，

式中：為氧氣在DMS內的擴散系數；koil為擴散系數的比值；Voil為局部DMS的體積。

氧氣在DMS內的擴散系數由Wilke-Chang公式（12）[17]計算

，

式中：為溶質A在溶劑B里的擴散系數；MB和分別為溶劑B的摩爾質量和粘度mPa?s；T為絕對溫度 K；VbA是溶質A在沸點時的摩爾體積 （O2: 25.6?cm3/ mol）；是與溶劑B相關的因子，DMS是1。在溫度333?K時，氧氣在DMS里的擴散系數是1.0×10-5?cm2 /s。

氧氣在水蒸汽里的擴散系數由式（13）計算

，

式中：Dij為雙組元氣體擴散系數。當 p0=1.01×105?Pa，T0=352.3?K時，DO2-H2O=0.352?cm2/s；p0=1.82×105?pa，T0=333?K時，氧氣在水蒸氣里的擴散系數DO2-H2O=

0.182?cm2/s。則由式（10）可得koil=5.5×10-5，所以在式（9）中koil Soil可以被忽略掉。

根據試驗測量的DMS填充的面密度值（如2.5?mg/cm2）乘以填充的面積4.2?cm2計算出DMS的填充總質量，然后將填充總質量除以DMS的密度就可以得到DMS的體積Voil，使用式（11）就可以計算得到DMS在多孔電極里的體積分數Soil。模擬計算DMS填充的面密度值2.5?mg/cm2對應的Soil值為0.333。

將發展后的氣體有效擴散系數修正式（9）繪制成曲線圖，如圖4所示。圖中展示了電極分別是抗水淹電極（其中DMS填充的面密度值2.5?mg/cm2）和傳統電極（其中DMS填充的面密度值0?mg/cm2）時液態水飽和度對組分氣體擴散系數的影響，單個曲線體現的是隨著液態水在多孔介質里的增加，因此阻礙了氣體的擴散，所以有效擴散系數減少。曲線間的不同體現的是，多孔電極填充了DMS后，由于DMS占據了部分氣體通道同樣也會阻礙氣體的擴散，所以氣體的有效擴散系數整體降低了。

電極反應速率：電極反應速率采用Butler-Volmer動力學方程。

，

，

式中：Ra和Rc分別為陽極和陰極電極反應速率，下標a、c分別代表陽極、陰極；Aa和Ac為催化層比表面積；i0,a和i0,c為交換電流密度；和分別為氫氣和氧氣的當地摩爾濃度；和分別為氫氣和氧氣的參考摩爾濃度；和為濃度指數，對于陽極=0.5，對于陰極=1；和為電化學反應能量傳遞系數。

在陰極催化層內反應生成的液態水不僅占據了傳質通道而且減少了催化劑的有效表面積，傳統多孔電極通過式（16）修正，考慮液態水對催化劑有效表面積的影響，如圖5所示。

k(s)=1-s．

為了描述液態水和DMS對電化學反應速率的影響，水的飽和度S和DMS在多孔電極里的體積分數Soil共同修正陰極的電化學反應方程，在式（16）基礎上發展為式（17）。

k(s)=[Soil+(1-Soil)(1-s)]y ，

式中：電化學反應速率修正函數指數y為3.5。

邊界條件和數值計算方法：在流道入口處解動量平衡方程時采用質量流量邊界條件，入口質量流量、組分濃度和溫度已知；在流道出口處解動量平衡方程時采用壓力出口邊界條件，在該面上的所有變量（壓力除外），梯度都為0；在流道進出口之間的固體壁面（除燃料電池各部分交界面）采用組分流量為0，無滑移壁面邊界條件，并假定極板外壁面溫度恒定為電池工作溫度，通過給定兩極板處的電勢求電流。模型仿真求解使用FLUENT12.0，應用求解壓力耦合方程組的半隱式方法SIMPLE(Semi-implicit Method for Pressure-linked Equations)[18] 對所有的守恒方程進行迭代求解，直到連續兩次迭代的各個物理量相對誤差為10-5。

2 抗水淹電極模型驗證

為了通過數值模擬計算來研究抗水淹電極PEMFC，使用上述考慮了DMS影響得到的式（9）和式（17）來分別模擬傳統電極PEMFC和抗水淹電極PEMFC。模擬計算中傳統電極PEMFC的Soil值為0，抗水淹電極PEMFC中的DMS含量為2.5?mg/cm2。

本文根據單電池試驗的實際情況，設置了仿真數學模型的物理參數和運行工況列于表1，其它參數主要從文獻[11]得到。

圖6對比了PEMFC試驗值和模擬計算值。圖中顯示了不論是傳統電極還是抗水淹電極的電池極化特性曲線，模擬計算和試驗值都吻合得較好，模擬計算和試驗結果一樣都顯示出了活化極化區、歐姆極化區和濃差極化區；還顯示出傳統電極PEMFC和抗水淹電極PEMFC性能的差異，試驗中單電流密度值大于1.32 A/cm2后陰極為抗水淹電極的單電池性能要好于傳統電極單電池性能，模擬計算結果中當電流密度值大于1.3 A/cm2后陰極為抗水淹電極的單電池性能要好于傳統電極單電池性能；在輸出電壓為0.35 V時，試驗中抗水淹電極較傳統電極PEMFC的電流密度值提高了0.87 A/cm2，模擬計算中抗水淹電極較傳統電極PEMFC的電流密度值提高了0.74 A/cm2。

圖7和圖8展示的是當輸出電壓為0.35 V時陰極催化層中間截面的電化學反應速率云圖。由圖7可見傳統電極電化學反應速率在圖的右邊低左邊高，也就是說產生的電荷分布右端低左端高，整體分布不均勻；由圖8可見抗水淹電極電化學反應速率在局部上成高低的交錯變化，但是在整體上產生的電荷分布較均勻，沒有出現圖 7中產生的電荷電池的一端低而另一端高的現象。這說明陰極電極填充了DMS后，能夠使電荷分布均勻，從而改善電池輸出性能。

3 抗水淹電極性能預測

為了進一步研究DMS填充量對電池輸出特性的影響，應用上述的抗水淹電極模型對不同的DMS含量的燃料電池極化特性進行模擬計算。模擬計算中設置了7種DMS含量，從0?mg/cm2(即電極不含DMS的傳統電極)到6?mg/cm2。將不同的DMS面密度值轉換成DMS的在孔隙里體積分數后，分別代入到式（9）和式（17），得到不同DMS含量時所對應的擴散系數修正方程和電化學反應速率修正方程，見表2。

圖9展示的是將不同DMS含量的擴散系修正函數的曲線圖，由圖可以看出隨著DMS的含量由0?mg/cm2增大到6?mg/cm2的過程中，有效擴散系數的值出現整體下降，這反映的是DMS對氣體擴散的阻礙；在每種DMS含量下，隨著液態水飽和度的增大，氣體有效擴散系數不斷下降，這反映的是液態水對氣體擴散的阻礙，所以計算模型同時反映了DMS和液態水共同對氣體擴散的阻礙效果。

圖10展示的是不同DMS含量的電化學反應速率修正方程的曲線圖。由圖可以看出隨著DMS的含量由0?mg/cm2增大到6?mg/cm2的過程中，催化劑的有效面積整體提高了，這反映的是DMS起到增大催化劑有效面積的作用；在每種DMS含量下，隨著液態水飽和度的增大，催化劑有效面積不斷減少，這反映的是液態水使催化劑有效面積減少，所以計算模型同時反映出了DMS和液態水共同對催化劑有效面積的影響。

綜合圖9和圖10的分析可知，增大DMS含量，一方面增大了催化劑的有效面積，另一方面也增大了對氣體的阻礙作用，前者有利于電池輸出特性的提高，后者不利于電池輸出特性的提高。因此，有必要對DMS的含量進行選擇，采用表2中不同DMS含量時的擴散系數修正方程和電化學反應速率修正方程，對質子交換膜燃料電池進行了模擬計算。

圖11展示的是不同DMS含量時PEMFC的極化特性。對比DMS含量為0?mg/cm2、1?mg/cm2、2?mg/cm2、3 mg/cm2和4?mg/cm2各極化特性曲線發現，隨著DMS含量的增大，在濃差極化段的輸出電壓不斷增大，即隨著DMS含量的增大PEMFC的輸出特性得到了提高。但是當DMS含量為5?mg/cm2而輸出電壓為0.45?V時，極化特性曲線出現了劇烈的下降，并且電流密度值有倒退的趨勢，這說明增大DMS含量對燃料電池輸出特性的負面影響開始占主導地位，即由于氣體擴散能力的降低，使多孔電極更早地出現了傳質極限，因此輸出電壓急劇下降，當DMS含量為6?mg/cm2時輸出電壓為0.5?V時就出現了電壓的急劇下降，可見當DMS含量超過4?mg/cm2后在高電流輸出時電池性能開始下降。為了更清楚地對比不同DMS含量所對應的電流密度值，取輸出電壓為0.35?V時不同DMS含量所對應的電流密度值繪制成圖12，圖中顯示開始電流密度值隨著DMS含量的增加而增加，當DMS含量超過4?mg/cm2后電流密度值開始下降。這說明抗水淹電極內DMS的含量不是越大越好，DMS的含量對于提高PEMFC性能存在著最優值，在本次試驗電池所對應的仿真模擬計算中DMS含量的最優值為4?mg/cm2。

文中采用DMS作為電極多孔介質中氧氣額外的固定傳輸通道，是在電極微孔中預先加入進去的，它不參與化學反應，不會被消耗，因此其生命周期與PEM同步，具有耐久性。

4 結論

篇6

關鍵詞: 圖像去噪；全變分模型；數字TV濾波器；擴散系數

中圖分類號:TP391.41

文獻標識碼:A文章編號:1672-8513(2010)03-0211-05

A New Numerical Realization Algorithm Based on the Total Variational Model

WANG Yiyan

(Department of Physics and Engineering Technology，Sichuan University of Arts and Science，Dazhou 635000，China)

Abstract:

Based on the classical total variational restoration model，this research proposes a new numerical realization algorithm of digital TV filter for image denoising from the viewpoint of local diffusion. The proposed algorithm can select the diffusion coefficients adaptively according to the local gradient of the image，so as to reduce the diffusion in edge regions while enhancing the diffusion in flat regions. It is used as an anisotropic filter for the noise image. So the proposed algorithm can not only suppress noises very effectively but also preserve image edge details very well. Experimental results show that this algorithm has better performance and higher computation efficiency than the conventional digital TV filter.

Key words: image denoising; total variational model; digital TV filter; diffusion coefficient

圖像在獲取或傳輸過程中不可避免地會遭受噪聲的干擾，從而導致圖像質量下降，對圖像的后續處理產生嚴重的影響.為了消除這種影響，就必須在后續處理前，對圖像進行去噪預處理.因此，圖像去噪一直是圖像處理領域的研究熱點.其難點在于：去噪的同時要求保留圖像的邊緣等細節信息.目前，基于偏微分方程的變分模型方法是圖像去噪的一種流行方法，其高質量的處理效果已引起了國內外廣大學者的重點關注［1-6］.Rudin等［1-2］提出的Total Variation(TV),模型是其中的典型代表.后來，Chan等［4］在文獻［1］的基礎上提出了數字TV濾波器，并證明該濾波器是基于TV模型的最小能量泛函在半點格式下［7］的等價數值求解格式.由于數字TV濾波器是建立在TV模型的基礎上，因此它具有良好的邊緣保持特性.同時，數字TV濾波器的設計特點還使之具有特殊的信息反饋機制.國內肖亮等［5］也對數字TV濾波器的設計及性能進行了深入的分析與研究，驗證了數字TV濾波器的上述性能.然而，傳統數字TV濾波器本質上是一種倒數梯度加權濾波器［5］，對圖像去噪后的像素值是通過原來噪聲圖像中對應點的周圍像素的梯度值平均加權得到，因此，它是一種線性濾波器，在去噪時容易導致圖像模糊.本文根據全變分最小的能量泛函，從局部擴散的角度出發，設計了一種新的應用于圖像去噪的數值算法.該算法可以根據圖像的局部梯度信息自適應地選取擴散系數，使邊緣處擴散較弱，在平坦區域擴散較強.實現了對噪聲圖像的各向異性濾波處理，在去除噪聲的同時能夠較好地保護圖像的細節信息.實驗結果表明，該算法去噪后的圖像比傳統數字TV濾波器具有更好的客觀評價指標和主觀視覺效果，同時運算效率也更高.

1 全變分模型及其數字濾波器

1.1 全變分模型

令u為原始清晰的圖像，u0為被噪聲污染的圖像，即u0=u+n.式中n是均值為0，方差為σ2的高斯噪聲.即：E(n)=0，E2(n)=σ2.（1）

全變分模型是由Rudin等［1-2］提出的，并且現在是圖像恢復中最成功的方法之一.Rudin等認為，有噪聲圖像的的全變分比無噪聲圖像的全變分明顯大，最小化全變分可以達到去除噪聲的目的.全變分模型是一種變換，可用下式表示：

TV(u)=∫ΩΔudxdy.（2）

其中Δ u=uxuy，Δu=ux2+uy2是圖像的梯度場，Ω為圖像的定義域，(x，y)∈Ω為像素點的二維坐標.

并假設所有可觀測的圖像具有有界變分，有界變分空間定義為：

艷：基于全變分模型的新型數值實現算法

BV(Ω)={u:u∈L1(Ω)，并且TV(u)

BV(Ω)是在BV范數uBV=uL1+TV(u)下的Banach空間.

由（1）式可知，TV范數最小化應滿足以下2個條件：

∫Ωudxdy=∫Ωu0dxdy，1Ω∫Ω(u－u0)2dxdy=σ2.（4）

其中，Ω為圖像區域Ω的面積.對于常數c，有TV［u+c］=TV［u］.因此，第1個條件通?？偸浅闪⒌模?］.于是，從u0中恢復u的問題就演變為在條件1Ω∫Ω(u－u0)2dxdy=σ2下，求TV(u)=∫ΩΔudxdy的最小化問題.通過引入參數λ，可以定義全變分最小的能量泛函：

minuJ(u)=∫ΩΔudxdy+λ2∫Ωu－u02dxdy.（5）

上式中第1項為圖像的正則化項，第2項為數據保真項.λ為Lagrange乘子，是平衡正則化項和保真項之間的權重參數.

由于目標函數J(u)是u的凸函數，極值存在的一階充分必要條件為：Δ J(u)=0.然而由于當Δu=0，一階條件不能滿足.為了克服這種退化問題，一般對圖像正則化項引入小擾動進行修正：

Δuε=ux2+uy2+ε2.（6）

只要ε保持足夠小，便不會影響噪聲濾波的性能.所以式（5）改寫為：

minuJ(u)=∫ΩΔuεdxdy+λ2∫Ωu－u02dxdy.（7）

這樣，就建立了圖像恢復的全變分模型.

1.2 數字TV濾波器

由一階最優條件，式（7）可歸結為具有黎曼邊界條件的歐拉-拉格朗日方程：

－Δ•ΔuΔuε+λ(u－u0)=0.（8）

記Ω為Ω的邊界，n為邊界Ω的法向量，黎曼邊界條件為：

un=0，x∈Ω . （9）

即式（8）成立時J(u)最小.由梯度下降法得：

ut=Δ•ΔuΔuε+λ(u0－u).（10）

數字TV濾波器［4］實質上就是對上述極值問題的數值逼近.對于數字圖像，設α是待處理的中心像素點，像素值為uα，β等是α的十字鄰域，如圖1所示.

圖1中hαα(u)和hαβ(u)等為加權系數，分別定義為：

hαα(u)=λλ+∑β~αωαβ(u)，hαβ(u)=ωαβ(u)λ+∑β~αωαβ(u),（11）

其中：ωαβ(u)=1Δαuε+1Δβuε.（12）

顯然，hαα(u)+∑β~αhαβ(u)=1.從式（11）和（12）可以看出，加權系數與中心像素和它鄰域像素的正則化梯度的倒數和有關.定義中心像素α和它的十字鄰域的正則化梯度分別為：

Δαuε=∑β~α(uβ－uα)2+ε2，

Δβuε=(uβ－uα)2+ε2.（13）

結合式（11）~（13），則可得到TV數字濾波器，記為Fε，λα(u)：

Fε，λα(u)=∑β~αhαβ(u)uβ+hαα(u)u0α.（14）

由式（14）可知，TV數字濾波器包含2個參數：

1） 正參數ε稱為規整化參數，為一個很小的正數，目的是避免圖像平坦區域梯度為0.ε取值很小，當ε≤10－4時，對濾波器影響很小.

2） 正參數λ稱為匹配參數，即Lagrange乘子，決定了濾波器如何在邊緣保護和平滑噪聲之間取得平衡，其選取依賴于噪聲水平σ.噪聲越大，λ越小，平滑作用越明顯，邊緣損失越大；反之則λ取值越大，邊緣保護越明顯.

TV數字濾波器采用迭代策略對圖像中的每個像素點進行處理，整個濾波過程如下：

初始化u(0)，通常取u(0)=(u(0)α1，u(0)α2，…，u(0)αN)

For k=1，2，….

For α1α2α3…αN

根據式（11）~（13）分別計算以α為中心點的Δαuε及其相鄰Δβuε等4點、ωαβ(u)、hαβ(u)、hαα(u)，然后利用式（14）計算u(k)αj=Fε，λαj(u(k－1)).

特別地，由于圖像在迭代過程中不斷變化，因此，需要適時更新λ值［4］：

λ=1σ21Ω∑α∈Ω∑β~αωαβuβ－uαuα－u0α,（15）

其中，σ2是圖像的噪聲方差，需要估計，Ω是圖像的面積.通常λ被初始化為1σ2.

由上述分析可知，數字TV濾波器實質上是一種倒數梯度加權濾波器［5］，對圖像去噪后的像素值是通過原來噪聲圖像中對應點的周圍像素的梯度值平均加權得到.因此，它是一種線性濾波器，在去噪時容易導致圖像模糊.

2 本文算法

針對數字TV濾波器的不足，我們對全變分模型給出了一種新的結合梯度加權最速下降法和半點格式的數值迭代算法.

首先由變分法［8］可推導出式（5）對應的歐拉-拉格朗日方程為：

－Δ•ΔuΔu+λ(u－u0)=0.（16）

由梯度下降法，引入人工時間變量t，則可得如下方程：

ut =Δ• ΔuΔu + λ(u0 - u),

unΩ=0,u(x,y,0)=u0 (x,y) . （17）

令η為圖像的梯度方向，ξ為與η垂直的切線方向（即邊緣方向），分別表示如下：

η=ΔuΔu=(ux，uy)u2x+u2y，ξ=ΔuΔu=(－uy，ux)u2x+u2y.（18）

定義uξξ和uηη分別為u在ξ和η方向上的二階導數，表達式如下：

uηη=u2xuxx+2uxuyuxy+u2yuyyu2x+u2y,uξξ=u2yuxx－2uxuyuxy+u2xuyyu2x+u2y . （19）

則式（17）可化為η－ξ方向上的偏微分方程(PDE)：

ut=1Δu•uξξ+λ(u0－u).（20）

由式（20）可看出，邊緣方向對應的擴散系數為1Δu，它依賴于圖像的局部梯度信息.在不同區域，擴散程度不一樣.在邊緣處，Δu較大，擴散系數較小，因此沿邊緣方向的擴散較弱，從而有利于保護邊緣；在平坦區域，Δu較小，擴散系數較大，因此在平坦區域的擴散能力較強，從而可以有效去除噪聲，整個過程相當于對噪聲圖像進行各向異性濾波處理.因此，該方法能在去噪的同時保護圖像的邊緣細節.

采用半點格式的中心差分來離散化式（20）中的PDE，為了利用像素周圍更多的信息，本文選擇8鄰域系統，如圖2.

由一階和二階導數的差分格式可以得到Δu，uξξ的離散形式，進一步對式（20）利用時間步進法求解，其迭代形式為：

un+1i，j=uni，j+δt•1Δu•uξξni，j+λu0－uni，j.（22）

在實際計算過程中，設圖像大小為M×N，每個像素位置點(i，j)∈1，…，M×1，…，N，采取對稱邊界條件為：u0，j=u1，j，uM+1，j=uM，j，ui，0=ui，1，ui，N+1=ui，N.并且為了避免上述方程中分母處的梯度為0，將梯度提升為Δuε=Δu2+ε2，其中，ε>0.Lagrange乘子λ可通過對式（16）兩邊同時乘以(u－u0)，然后在整個圖像區域Ω上積分求得：

λ=1σ2Ω∫ΩΔ•ΔuΔu•(u－u0)dxdy.（23）

3 實驗結果與分析

為了驗證本文算法的有效性，我們進行了大量的仿真實驗.采用平均絕對誤差(MAE)和峰值信噪比(PSNR)作為去噪性能的客觀評價尺度，MAE和PSNR分別定義如下：

MAE=1MN∑i，jyi，j－xi，j，（24）

PSNR=10 lg(M×N)max(xi，j)2∑i，j(yi，j－xi，j)2，（25）

式中，x是大小為M×N的原始圖像，y是去噪后的圖像.MAE值越小越好，而PSNR值越大越好.

以Woman灰度圖像為測試對象，對其加入均值為0，不同噪聲方差（依次為5，10，15，20）的高斯噪聲，將本文算法與數字TV濾波器［4］進行了去噪比較.圖3給出了方差為20時，Woman圖像去噪后的視覺效果；表1為不同噪聲條件下，2種方法去噪前后的客觀評價指標對比；表2為2種算法在去噪過程中所需時間,運算在Matlab 6.5平臺上進行.

從實驗結果來看，本文算法去噪效果、客觀評價指標和運算效率方面都比傳統數字TV濾波器具有明顯的優勢.

4 結語

本文對全變分模型提出了一種新的數值實現方法，與傳統數字TV濾波器不同，該算法從局部擴散的角度，充分利用了圖像想素點的梯度信息，自適應的選取擴散系數，對噪聲圖像進行各向異性濾波處理.實驗結果表明，本文算法在去除噪聲的同時能較好的保護圖像細節，與傳統數字TV濾波器相比，本文算法在去噪效果、客觀評價指標和運算效率等方面均具有明顯的優勢.

參考文獻:

［1］RUDIN L，OSHER S，FATEMI E.Nonlinear total variation based noise removal algorithms［J］.Physica D，1992，60(1-4): 259-268.

［2］RUDIN L，OSHER S.Total variation based image restoration with free local constraints［C］//Proc 1st IEEE Int Conf Image Process，Austin，TX: IEEE Press，1994，1: 31-35.

［3］COBSON D，VOGEL C.Convergence of an iterative method for total variation denoising［J］.SIAM Journal of Numerical Analysis，1997，34(5): 1779-1791.

［4］CHAN T F，OSHER S， SHEN J.The digital TV filter and nonlinear denoising［J］.IEEE trans on Image Processing，2001，10(2): 231-241.

［5］肖亮，吳慧中，韋志輝，等.基于總體變差模型的數字濾波器設計及其性能研究［J］.信號處理，2003，19(3): 247-251.

［6］莊紅林，王文斌，段鵬，等.整體變分最小化方法在圖像處理中的應用［J］.云南民族大學學報：自然科學版，2008，17(2): 176-180.

篇7

【關鍵詞】應用型人才；模塊化；創新能力；教學改革

0 引言

應用型人才培養模式的具體內涵是隨著高等教育的發展而不斷發展的“應用型人才培養模式是以能力為中心，以培養技術應用型專門人才為目標的”[1]。本科應用型是本科層次教育，既有著普通本科教育的共性，又有別于普通本科的自身特點，它更加注重的是實踐性、應用性和技術性。

隨著現代科學技術的發展和計算機的廣泛應用，科學計算已經成為平行于理論分析和科學實驗的第三種科學手段。而《數值分析》這門課程是科學計算的基礎和核心，逐漸成為本科生和研究生的必修課。是理論與實踐結合緊密的一門科學。那么在我們講授《數值分析》課程的過程中，將理論分析和實驗教學相結合，不斷地改進我們的教學方法，吸取其他數值分析教育工作者的改革措施和經驗[2-4]激發學生學習的興趣和熱情， 提高學生對算法的實際應用能力。

《數值分析》是一門介紹科學計算的核心理論和基本方法的數學課程[5]，分析《數值分析》課程的特點，根據在該課程的具體教學實踐中的體會，對該課程的內容進行分類整理，論述以MatLab 為教學平臺，并對以數學建模為實驗內容的教學方式的改革作初步探索。

1 《數值分析》課程的特點及傳統教學中的不足

《數值分析》課程有不同于其他數學專業課的特點， 有其自身的理論體系。數值分析要面向計算機。數值分析是為了解決實際問題而產生的科學，每個算法除了在理論上可行外， 還要通過數值試驗證明是可行的。這對學生的編程能力以及數學軟件的應用能力提出了較高的要求。

基于《數值分析》課程以上的這些特點可以看出《數值分析》是一門與計算機使用密切結合的，實用性和實踐性都很強的課程。但在傳統的教學中，只注重講授數值方法的原理以及誤差估計等理論。另外，在教學中由于缺乏對數值實驗的重視以及學生數學軟件運用能力和編程能力的限制， 使學生不能把這些算法在計算機上實現并應用到解決實際問題中去。這幾年我們一直在不斷地對教學過程進行思考和探索，不斷地總結教學經驗，克服教學過程中的不足，不斷地對教學方法進行改革。

2 對《數值分析》課程教學改革的幾點體會

2.1 數值分析的模塊化教學內容改革

“數值分析”是研究各種數學問題求解的數值計算方法及其理論的一門課程[5]。在理、工科數學類科目教學體系中，“數值分析”起著承上啟下的作用，是高等理工科院校的重要基礎課程，同時也具有培養學生創新思維、創新能力的特點。數值分析的教學目的不是讓學生僅僅會利用己有的算法去解決某些問題，而應有理有據地、系統地傳授其基本思想、基本方法和基本原理。數值分析的宗旨應是提高學生數學素養和使學生掌握實用算法并重。這就是說要使學生能夠知其然，更得知其所以然，并能讓學生在遇到新問題時有能力依靠數值方法去解決，注重培養學生舉一反三的能力。因此，理論方法和應用都不可忽視，這就要求對教學內容、環節等重新優化設計。

經過我們的分析發現《數值分析》課程主要內容可以總結為三個方面：①離散問題的連續化， 即根據離散數據構造一個簡單且易于計算的函數代替原有的復雜函數或數據， 例如函數插值、函數逼近、數據擬合等內容；②連續問題的離散化，即在一系列離散的點上求解數值解，例如非線性方程求根、數值積分、微分方程數值求解等內容；③對于方程組的求解以及矩陣的特征值計算等內容。在教學過程中，我們要從宏觀上把課程內容分析給學生，讓學生沿著這種思路把課本中的內容系統化， 把各個分散的內容都建立聯系進行模塊化教學。模塊化的教學結構能夠更好地適應和啟發引導學生積極思維，大膽創新。

2.2 引入MatLab 軟件為主的計算機平臺，利用現代化的教學手段

在教學過程中，根據本課程的特點，不能一味地只做理論分析，而是要打破傳統的教學體系，采用多媒體輔助教學設備，引入采用數值解法、解析解法和圖形可視化相結合的方法。

在實際教學中， 我們在利用傳統的教學手段的同時，引入現代化教學手段―――電子教學，將一些推導過程、實例應用利用數值計算軟件，以可視化的方式來表現數值計算的過程，通過分步實現、誤差數據、幾何圖形直觀地展示給學生， 從而使學生更容易理解抽象的理論， 也能讓學生學習到如何將算法在計算機上實現。我們將計算機軟件MatLab 引入數值分析的各個教學內容之中，在教學中充分利用MatLab 軟件為主的計算機平臺。MatLab 是最為普遍和功能強大的科學計算工具，我們針對具體的數值算法，編寫MatLab 數值算法程序和圖形可視化程序， 進行數值數據和圖形的比較。引導學生理解所學的內容，激勵學生的興趣，培養他們自己動手解決問題的能力和創新能力。

2.3 從實際問題出發，將數學建模案例和工程應用實例引入課堂教學

《數值分析》是一門應用很強的課程，以數學建模案例和工程應用實例引導學生運用數學方法解決實際問題，因此選擇與授課核心內容相對應的數學建模案例和工程實例進行教學具有必要性和可能性.例如用2011年全國大學生數學建模競賽的C題“山東省職工工資預測問題”作為實例講解曲線擬合.職工工資增長有兩個基本模型：一個是指數增長模型，一個是阻滯增長模型.按照兩種不同的模型，使用非線性最小二乘法分別擬合自1978年至2010年的工資數據，預測從2015年至2035年的工資數據，并將預測數據與實際數據比較，可以讓學生形象直觀地用預測結果分析兩個模型的優劣.又如網絡流量、核輻射量的測量與計算與數據擬合理論，根據山區地形采樣點數據來繪制地形圖與插值理論，導彈追蹤問題與微分方程數值解，核廢料的妥善處理、油罐刻度設計問題與方程組的求根問題等，都是把數值分析與數學建模和工程實例相結合的典型例子，會給《數值分析》課程的教學帶來蓬勃

3 結語

應用型本科院?；谀K化的數值分析教學，增強了課堂教學的直觀性，優化了課堂教學與實驗教學，能真正實現教與學的良性互動?；诂F代化的教學手段的數值分析教學，培養學生應用數值計算方法解決實際工程問題的能力，開拓學生思維，培養具有創新能力與素質的應用型現代化建設人才。

【參考文獻】

[1]莊華潔，周金其.本科應用型人才培養模式的研究與實踐[J].高等教育研究，2004（6）.

[2]李大美，等.基于創新能力培養的計算方法課程改革[J].中國大學教學，2007.

[3]張光輝，任敏.MATLAB平臺上《數值分析》課程教學的幾點思考[J].甘肅聯合大學學報（自然科學版），2012，26（5）：103-105.

篇8

關鍵詞：建筑照明；智能控制；ZigBee技術

中圖分類號：TM923 文獻標識碼：A

隨著國民經濟的快速發展，能源消耗成為國家和社會關注的熱點問題，同時隨之而來的節能減排工作成為研究的重要內容之一。據不完全統計在現代建筑中照明系統對能源的消耗占的比例較大[1-2]。因此，建筑電氣照明節能技術的理論及應用研究成為了建筑電氣綜合節能研究內容之一。本文以學校智能照明控制系統設計為例，研究并設計出一套基于Zigbee技術智能照明控制系統，該系統詳細的介紹了控制系統的總體設計方案、硬件電路設備選型、Zigbee無線網絡的組建以及軟件系統設計等內容。本文所提方案和做法對于建筑電氣照明節能起到一定的參考和借鑒作用。

1 智能照明控制系統總體設計方案

針對學校的具體情況本文提出一種基于ZigBee 技術的智能照明控制系統。該系統主要由調光控制電路、ZigBee 節點模塊、ZigBee 無線通信網絡及遠程控制上位機四個部分組成。智能照明控制系統結構圖如圖 1 所示。

從圖中可以看出，整個控制過程如下：信息采集電路將所采集的信息經單片機處理后傳遞給ZigBee 收發節點模塊將其轉換為無線電信號發送到 ZigBee 無線網絡中，經過ZigBee 收發協調器模塊傳輸給上位機，在上位機中將獲得相關的現場采集數據。

圖1 智能照明控制系統結構圖

2 智能照明硬件系統設計

本文的智能照明控制系統主要以單片機為核心，單片機采用AT89S52芯片，系統同時還設有參數修改、數據采集與處理以及聯動控制等功能，其結構如圖 2 所示。從圖2可以看出，該控制電路主要由光照傳感器模塊、紅外傳感器模塊、ZigBee節點模塊、LED 燈管驅動電路模塊、保護控制電路、輸入/輸出等相關部分構成。當系統檢測到室內光照度后，首先將其數字化處理，然后送給單片機。單片機根據所檢照度值，將所產生的脈沖寬度調制信號送到驅動電路中，從而控制流過 LED 燈管上的電流，以調整室內亮度使之與其相匹配。其硬件選型如表1所示。

圖 2 智能照明控制電路結構示意圖

表1 智能照明控制硬件選型表

序號 硬件名稱 硬件選型

1 單片機 AT89S52

2 通信網絡節點 CC2430

3 光照傳感器 TSL2561 光傳感器芯片

4 紅外傳感器 BIS0001傳感器芯片

3 ZigBee無線網絡的組建與智能照明控制系統軟件設計

3.1ZigBee無線網絡的組建

由于本文的智能照明系統設計對象主要指學校，而且在系統中主要以傳感器作為主要接入對象，其數據傳輸量相對而言不大，采用一般的網絡拓撲結構進行數據通信即可。由于學校的直徑一般在幾公里范圍內，然而ZigBee 節點數據傳輸距離僅有 100m 左右，因此為了增加系統數據的傳輸距離，在整個網絡布置時需增加ZigBee 路由器。正是由于基于以上考慮，本文所設計的無線網絡拓撲結構采用簇樹型網絡拓撲結構。該網絡結構具有成本低、實現容易、可擴充、易操作等優點。本文所設計的智能照明控制網絡拓撲結構如圖 3 所示。從圖3中可以看出，學校網絡協調器和各樓宇子網路由器是其中非常重要的兩個設備，網絡協調器完成網絡建立、協調整個網絡運行。相關各個教學樓路由器主要負責數據傳輸作用，目的是為保證將現場所采集的數據經過長距離傳送后到達協調器。教室中的各個終端控制節點模塊，主要負責現場的數據采集和初步處理以及對各個LED控制管理。

圖3 智能照明控制系統網絡拓撲結構

3.2智能照明控制系統軟件設計

智能控制系統程序流程圖如圖 4 所示，從圖4看出，系統首先進行初始化，之后可以通過參數的錄入設定系統所需要的光照值，現場采集信息系統含有兩個部分，一部分采集光照強度，另一部分通過人體紅外傳感器采集是否有人員的信息。上述這些信息進入到單片機中執行相應的算法程序，本文所采用的自動控制算法為比例-積分-微分（PID）算法，通過該算法的運算輸出合適的控制量以減少偏差可以實現恒定光照控制，通過幾次循環反饋系統進入較穩定的狀態，教室的光照強度基本上穩定在預先設定的數值左右。本文所采用的調光方式為 PWM方式，由于單片機系統不帶PWM模塊，因此需要模擬一個PWM信號，本文利用數字電路中的定時器來模擬產生一個 PWM 信號。

圖4 智能照明控制系統程序流程圖

4結論

本文將ZigBee 技術應用于學校智能照明控制系統中，對智能照明控制系統的軟、硬件進行設計，其中重點對系統方案、智能照明控制節點以及ZigBee無線網絡的組建進行了詳細地闡述，本系統可以根據室內光照調節照度，從而使學校教室內具有很好的照度。同時本系統采用ZigBee方式，減少了布線鋪設，使得工程的施工更加容易，達到國家所提倡的節能減排，該系統的設計方案為其它智能照明控制系統提供了參考與借鑒。

參考文獻

[1]任紅.綠色照明與建筑照明節能設計[J].電氣應用，2009（7）：21-25.

[2]陳歲生，盧建剛 . 基于智能檢測技術的室內照明系統 [J]. 機電工程技術， 2008，37（5）：53-55.

[3]朱明碩，陳校.淺談照明節能技術的應用[J].電氣應用，2010（3）：54-58.

[4]樂鈺唯.照明控制技術的發展及應用現狀[J].照明工程學報，2010（3）：14-18.

篇9

[關鍵詞] 旅游收入 增加值 旅游收入乘數

一、旅游乘數的研究回顧

1882年,經濟學家巴奇霍特(Bagehot)分析了緊縮產業對經濟中其他產業引起的負面影響后,不少經濟學家便開始對乘數理論展開了研究。1931年,美國經濟學家卡恩(Kahn)首次把乘數概念引入經濟理論,他的乘數計量的是總就業量的增量和與之相關的投資品行業的初級就業量的增量之間的比。1936年,凱恩斯(Keynes)對其進行補充,形成了乘數模型。之后現代乘數理論大多沿著凱恩斯的乘數模型和里昂惕夫投入-產出模型發展。1982年馬西森(Mathieson)和沃爾(Wall)提出了旅游乘數的概念,即“旅游乘數是這樣一個數值,最初旅游消費和它相乘后能在一定時期內產生總收入效應”。

1998年,張凌云在《南開經濟研究》上發表的《旅游業乘數效應的幾個問題》是我國旅游業乘數研究的早期成果。西方學者弗萊徹(Fletcher)和斯尼(Snee)鑒定了六種旅游乘數,包括產出乘數、銷售乘數、收入乘數、就業乘數、政府收入和進口乘數。左冰(2002)、石培華(2003)、吳國新(2003)在各自的論文中研究了我國的產出乘數和就業乘數。但是,更多的學者偏向研究乘數作用對區域旅游的影響,李志清(2001)對上海、李興緒(2004)對云南、譚佩琴(2006)對澳門、王瑜(2006)對福建、張小利(2007)對西部等。乘數計算方面,我國學者田里(2002)、李天元(2002)、羅明義(2004)等研究旅游乘數的時候,在各自出版的旅游經濟學書籍里面均提出了不盡相同的公式?？锪?1996)著重介紹了旅游業乘數的類型和計算,同時介紹了經濟分析模型。對于旅游收入以李江帆(1999)為先,認為計算增加值代替總收入有其合理性,其后有魏小安(2000)、宋子千(2001)對此繼續深入討論。智瑞芝(2003)、鄒煒(2008)在其論文中分別對黑龍江和云南旅游產業進行了增加值的計算研究。

二、旅游收入乘數的測算

1.增加值法的引入

由國民經濟核算理論上看,增加值是生產單位從總產出的價值中間除去其中所含貨物、服務消耗的價值之后的余額,它代表的是該生產單位會及各種生產要素在生產過程中新創造的價值(高敏雪,2006)。眾所周知,國民生產總值GDP的核算就是從增加值的角度出發,衡量一個國家或地區在一定時期內新創造的價值。在宏觀上,增加值算法可以避免重復計算,加總結果可以表現整個經濟生產的真實成果。旅游收入是旅游經營者從旅客們在旅游過程中購買服務、貨物來滿足他們旅行中的游覽、住宿、飲食、娛樂、購物等方面獲得的貨幣收入。但是旅游總收入中旅游產品所消耗的中間產品,并不是旅游業所生產的,而是與旅游相關的其他部門生產的。因此,嚴格上來說,若直接用旅游收入計算在國民經濟中所占份額的話,可能會扭曲旅游消費對國民經濟的實際貢獻。所以,應該從旅游總收入中扣除實物產品和中間消耗的價值,來得到旅游增加值。在實際問題中,由于統計的難度和問題,我國并沒有相關單獨的統計數據。旅游活動發生在乘坐交通工具、購物、用食、娛樂活動之中,數據上很難區分是屬于日常生活還是旅游,這樣由旅游產生的經濟價值被統計作為相關的交通、商業、餐飲和其他服務業的價值中。難以獨立統計,這給相關研究帶來了難度,不過,運用旅游增加值剝離測算方法,可以將相關產業中屬于旅游創造的價值提取出來。

按照李江帆(1999)對旅游增加值的界定,認為旅游主要涉及交通、郵電、商業、飲食業和社會服務業這五個第三產業部門。這些部門的產品只有一部分投入旅游消費,其增加值需要按一定的比例剝離出,計入旅游增加值中。這些含有游客消費的服務行業提供的增加值中旅游者消費部分所占比重稱為該行業的旅游消費剝離系數。利用旅游消費剝離系數分別計算出各個行業的旅游增加值,然后將這些加總就得到該年的旅游總收入。計算步驟如圖1所示。在計算的時候基于了兩個假設,第一,假定行業內部的增值率相同,也就是比如在餐飲消費中不論是花在旅游上還是在其他方面帶來的銷售額增長都相同。第二,假定旅游消費剝離系數不隨時間的變化而變化。該假設一方面由于統計局限,旅游收入中的行業構成無法得到各年的數據,另一方面也是為了化簡計算。這樣我們可以認為交通、郵電、商業銷售、飲食業和社會服務業這五個行業的旅游消費剝離系數為24.8%、7.2%、3.0%、19.5%和54.1%。再根據統計年鑒上每年相關數據計算各個行業的增加值,乘上剝離系數,得到各個行業的旅游增加值,最后加總獲得該年旅游總收入,即:∑(i行業旅游消費剝離系數×i行業增加值)=某年全國旅游增加值。具體計算結果在下表2中。

2.模型的建立和求解

根據何佳梅(2005)關于旅游收入乘數測定的方法,按照柯布-道格拉斯生產函數Y=A?Ka?Lβ來構造關于旅游收入的生產函數GDP=A?Ia,其中GDP為國內生產總值,I為增加值算法下的旅游收入,A為方程的常系數,而a則實際上表示的是旅游收入產出彈性。

由經濟學的知識可以知道,旅游收入產出彈性指的就是旅游收入I的變化引起GDP產出的變化,即:

那么對旅游收入的生產函數GDP=A?Ia,根據微積分可知,a表示的就是旅游收入產出彈性,即:Ei=a

給旅游收入生產函數GDP=A?Ia兩邊取自然對數則可得到:

Ln(GDP)=Ln(A)+a*Ln(I)

利用最小二乘法(OLS),就可以估計出a的值,也就是旅游收入產出彈性。在估算的時候,首先,不直接取GDP的值,而將其指數化,以便保證數據的可比性。G是以1978年為基年計算的GDP指數,替代GDP進行估算。第二,在公式的右邊加入時間趨勢變量t,用來剔除時間趨勢。最后估算公式調整為:

Ln(G)=Ln(A)+a*Ln(I)+b*t

Ln(G)=c+a*Ln(I)+b*t

利用Eviews5.0軟件,運用OLS進行參數估計:

Ln(G)=6.92-0.095*Ln(I)+0.10*t

(17.61)(-1.82)(15.60)

R2=0.9966,調整R2=0.9958,F-statistic=1320.488,Prob(F-statistic)=0.000

旅游收入產出彈性為負,這與經濟事實不相符合,因而懷疑方程具有多重共線性。根據Klein判別法進行分析,給出解釋變量間的簡單相關系數矩陣,發現其中一個簡單相關系數大于R2,說明模型中存在多重共線性。因此,必須對模型進行修正。利用修正Frish法進行修正得到:

由對比可見,修正后旅游收入產出彈性都得到了比修正前更好的擬合

數據來源:相關年度的《中國統計年鑒》

3.測算結果小結

從模型測算出來乘數的結果上來看,呈現自1996年逐年下降,但是2003年突然增加,再到2004年下降后到2007年逐步上升的走勢。自1998年我國出于刺激消費擴大內需的愿望,中央經濟工作會議將旅游業列為我國國民經濟發展新的增長點后,各省市把旅游業列入支柱產業或主導產業的名單中。在這樣的宏觀背景下,旅游收入增長率在隨后幾年超過了國內生產總值的增長率,造成乘數的逐漸變小。但是2003年我國遇上非典的非常時期,旅游業在一定程度上遭到打擊,總收入相較于上年銳減了20%左右,但是國內生產總值的增長卻仍然保持了一貫的高水平,這樣計算出的乘數自然會突增。到2004年旅游業基本恢復到2003年以前的水平,2004年開始我國的國內生產總值繼續保持高增長,這也就造成了旅游收入乘數的回升。

三、對旅游收入乘數的認識

測算乘數可以從宏觀上反映某行業部門對整體經濟的影響程度,體現各部門相互聯系,任何部門最終需求的變動都會引起整個經濟中產出、收入、就業等變動,即現代經濟的特點。一方面,對于乘數較大的行業可以選擇著重或優先發展,因為其可較大程度的拉動整體經濟,另一方面要看清乘數作用的發揮有一定的局限性。即如旅游業是個產值貢獻大、關聯性強,俗稱“興一業,旺百業”的行業,其發揮聯動增長作用的前提之一是社會上有一定數量可利用的閑置資源、存貨、非充分就業等,這樣才可以保證有擴張的后勁。對于乘數具體數量上的看待同樣要帶有理性的眼光,不論多精密的模型所討論的總是基于若干假設,簡化不可量化或者很難度量的條件,可能造成與實際上的偏差。如閆敏(1999)在驗證旅游產業和經濟發展水平之間關系的文章中,就涉及到了23各部門,這種復雜的產業交織的關系增加了問題的復雜性,要完全的從相關產業中剝離出旅游業的增加值幾乎不可能。再如2003年的旅游收入乘數,數值上均比前后兩年的大,但我們不可以直接得出當年旅游業發展突出的好的結論,還要看當年的具體宏觀背景以及相應發生的不可抗拒的事件。

參考文獻:

[1]葉秀霜:旅游經濟學[M].北京:北京大學出版社,2005:152

[2]朱沁夫:旅游經濟學[M].長沙:湖南大學出版社,2005:154～155

[3]師守祥:旅游業乘數研究辯證[J].旅游學刊,2007,(10):30～33

[4]匡 林:關于旅游乘數理論的幾個問題[J].華僑大學學報(社會科學版),2001,(1):9～11

[5]李江帆 李美云:旅游產業與旅游增加值的算法[J].旅游學刊,1999,(5):16～19

篇10

關鍵詞：告警數據 Hadoop Spark

1 引言

隨著電信網絡的不斷演進，全省數據網、交換網、接入網設備單月產生告警原始日志近億條。以上告警通過網元網管、專業綜合網管、智能網管系統[1]三層收斂，監控人員每月需處理影響業務或網絡質量的告警事件為20萬條，但一些對網絡可能造成隱患的告警信息被過濾掉。如何從海量告警數據中獲取與網絡性能指標、運維效率相關的有價值的數據，對于傳統的關系型數據庫架構而言，似乎是一個不可能完成的任務。

在一般告警量情況下，ORACLE數據處理能力基本可以滿足分析需求，但當告警分析量上升到億級，如果采用傳統的數據存儲和計算方式，一方面數據量過大，表的管理、維護開銷過大，要做到每個字段建索引，存儲浪費巨大；另一方面計算分析過程耗時過長，無法滿足實時和準實時分析需求。因此必須采用新的技術架構來分析處理海量告警信息，支撐主動維護工作顯得尤為必要，為此我們引入了大數據技術。

2 分析目標

（1）數據源：電信運營商網絡設備告警日志數據，每天50 G。

（2）數據分析目標：完成高頻翻轉類（瞬斷）告警分析；完成自定義網元、自定義告警等可定制告警分析；完成被過濾掉的告警分析、TOPN告警分析；核心設備和重要業務監控。

（3）分析平臺硬件配置：云計算平臺分配8臺虛擬機，每臺虛機配置CPU16核；內存32 G；硬盤2 T。

3 制定方案

進入大數據時代，行業內涌現了大量的數據挖掘技術，數據處理和分析更高效、更有價值。Google、Facebook等公司提供可行的思路是通過類似Hadoop[2]的分布式計算、MapReduce[3]、Spark[4]算法等構造而成的新型架構，挖掘有價值信息。

Hadoop是Apache基金會用JAVA語言開發的分布式框架，通過利用計算機集群對大規模數據進行分布式計算分析。Hadoop框架最重要的兩個核心是HDFS和MapReduce，HDFS用于分布式存儲，MapReduce則實現分布式任務計算。

一個HDFS集群包含元數據節點（NameNode）、若干數據節點（DataNode）和客戶端（Client）。NameNode管理HDFS的文件系統，DataNode存儲數據塊文件。HDFS將一個文件劃分成若干個數據塊，這些數據塊存儲DataNode節點上。

MapReduce是Google公司提出的針對大數據的編程模型。核心思想是將計算過程分解成Map（映射）和Reduce（歸約）兩個過程，也就是將一個大的計算任務拆分為多個小任務，MapReduce框架化繁為簡，輕松地解決了數據分布式存儲的計算問題，讓不熟悉并行編程的程序員也能輕松寫出分布式計算程序。MapReduce最大的不足則在于Map和Reduce都是以進程為單位調度、運行、結束的，磁盤I/O開銷大、效率低，無法滿足實時計算需求。

Spark是由加州伯克利大學AMP實驗室開發的類Hadoop MapReduce的分布式并行計算框架，主要特點是彈性分布式數據集RDD[5]，中間輸出結果可以保存在內存中，節省了大量的磁盤I/O操作。Spark除擁有Hadoop MapReduce所具有的優點外，還支持多次迭代計算，特別適合流計算和圖計算。

基于成本、效率、復雜性等因素，我們選擇了HDFS+Spark實現對告警數據的挖掘分析。

4 分析平臺設計

4.1 Hadoop集群搭建

基于CentOS-6.5系統環境搭建Hadoop集群，配置如表1所示。

4.2 Spark參數設置[6]

Spark參數設置如表2所示。

4.3 數據采集層

數據采集：由于需采集的告警設備種類繁多，故采取分布式的告警采集，數據網設備、交換網設備、接入網設備分別通過IP綜合網管、天元綜合網管、PON綜合網管進行采集，采集周期5分鐘一次。采集機先將采集到的告警日志文件，通過FTP接口上傳到智能網管系統文件服務器上，再對文件進行校驗，通過Sqoop推送到Hadoop集群上。

4.4 邏輯處理層

（1）建立高頻翻轉告警監控工作流程

先將海量告警進行初步刪選，通過數量、位置和時間三個維度的分析，得出高頻翻轉類告警清單列表，最后由專業工程師甄別確認，對某類告警進行重點關注和監控。

（2）差異化定制方案

按組網架構細分，針對核心重要節點的所有告警均納入實時監控方案；

按業務網絡細分，針對不同業務網絡設計個性化的監控方案；

按客戶業務細分，針對客戶數字出租電路設計個性化的監控方案。

4.5 數據分析層

Spark讀取Hive[7]表的告警數據，然后在Spark引擎中進行SQL統計分析。Spark SQL模K在進行分析時，將外部告警數據源轉化為DataFrame[8]，并像操作RDD或者將其注冊為臨時表的方式處理和分析這些數據。一旦將DataFrame注冊成臨時表，就可以使用類SQL的方式操作查詢分析告警數據。表3是利用Spark SQL對告警工單做的一個簡單分析：

5 平臺實踐應用

探索運維數據分析的新方法，利用大數據分析技術，分析可能影響業務/設備整體性能的設備告警，結合網絡性能數據，找到網絡隱患，實現主動維護的工作目標。

5.1 高頻翻轉類告警監控

首先制定了高頻翻轉類告警分析規則，將連續7天每天原始告警發生24次以上定義為高頻翻轉類告警，并基于大數據平臺開發了相應的分析腳本，目前已實現全專業所有告警類型的分析。表4是全省高頻翻轉類TOP10排名。

5.2 核心設備和重要業務監控

目前以設備廠商或專家經驗評定告警監控級別往往會與實際形成偏差，主要表現在以下幾個方面：監控級別的差異化設定基于已知的告警類型，一旦網絡重大故障上報未知的告警類型就無法在第一時間有效監控到；同一類型的故障告警出現在不同網絡層面可能影響業務的程度是完全不同的；不同保障級別的客戶對故障告警監控的實時性要求也是不同的。

通過大數據分析平臺對差異化監控提供了靈活的定制手段，可根據告警關鍵字，分專業、地市、網管、機房、告警頻次等維度自主定制需要的告警數據，實現日、周、月、某個時間區等統計分析。

應用案例：省NOC通過大數據分析出一條編號為CTVPN80113的中國平安大客戶電路在一段時間內頻繁產生線路劣化告警，但用戶未申告，省NOC隨即預警給政企支撐工程師，政支工程師與用戶溝通后，派維護人員至現場處理，發現線路接頭松動，緊急處理后告警消除、業務恢復。

5.3 被過濾告警分析

全省每天網絡告警數據300萬條～500萬條，其中99%都會根據告警過濾規則進行過濾篩選，把過濾后的告警呈現給網絡監控人員。過濾規則的準確性直接影響告警數據的質量。一般來說告警過濾規則可以從具有豐富運維經驗的網絡維護人員獲得，但是這個過程非常繁瑣，而且通過人工途徑獲得的告警過濾規則在不同的應用環境可能存在差異，無法滿足網絡維護的整體需要。采用大數據技術對被過濾的告警進行分析可以很好地完善過濾規則，讓真正急迫需要處理的告警優先呈現給維護人員及時處理，真正做到先于客戶發現故障。表5是動環專業被過濾的告警情況分布。

5.4 動環深放電分析

動環網管通過C接口采集蓄電池電壓數據，在停電告警產生之后，電壓數據首次下降到45 V，表示該局站電池出現深放電現象，通過計算這一放電過程的持續時間，記為深放電時長，該時長可以初步反映電池的放電性能。一個局站每天產生幾十萬條電壓等動環實時數據。

在告警數據分析的基礎上，實現對蓄電池電壓變化數據的分析，提醒分公司關注那些深放電次數過多和放電時長過短的局站，核查蓄電池、油機配置、發電安排等，并進行整治。利用Spark SQL統計了一個月內撫州、贛州、吉安三分公司幾十億條動環數據，分析了其中深放電的情況如表6所示。

6 結論

本文利用HDFS+Spark技術，實驗性地解決告警數據存儲和分析等相關問題：一是通過數據分析，從海量告警數據中發現潛在的網絡隱患；二是結合資源信息和不同專業的告警，最終為用戶提供綜合預警；三是轉變網絡監控思路和方式，通過數據匯聚、數據相關性分析、數據可視化展示，提高了網絡監控效率；最后還擴展到對動環實時數據、信令數據進行分析。

從實際運行效果來看，HDFS和Spark完全可以取代傳統的數據存儲和計算方式，滿足電信運營商主動運維的需求。

參考文獻：

[1] 中國電信股份有限公司. 中國電信智能網管技術規范-總體分冊[Z]. 2015.

[2] Tom white. Hadoop權威指南[M]. 4版. 南京： 東南大學出版社， 2015.

[3] RP Raji. MapReduce： Simplified Data Processing on Large Clusters[Z]. 2004.

[4] Spark. Apache Spark?[EB/OL]. [2016-11-27]. http：///.

[5] Matei Zaharia， Mosharaf Chowdhury， Tathagata Das， et al. Resilient Distributed Datasets： A Fault-Tolerant Abstraction for In-Memory Cluster Computing[J]. Usenix Conference on Networked Systems Design & Implementation， 2012，70（2）： 141-146.

[6] S鵬. Apache Spark源碼剖析[M]. 北京： 電子工業出版社， 2015.

[7] Hive. Apache HiveTM[EB/OL]. [2016-11-27]. http：///.

[8] Holden Karau， Andy Konwinski， Patrick Wendell， et al. Learning Spark： Lightning-Fast Big Data Analysis[M]. Oreilly & Associates Inc， 2015.

[9] 員建廈. 基于動態存儲策略的數據管理系統[J]. 無線電工程， 2014，44（11）： 52-54.