江 山，易年余，孫美玲

（1.揚州大學數(shù)學科學學院，江蘇 揚州225002；2.湘潭大學數(shù)學與計算科學學院，湖南 湘潭411105；3.南通職業(yè)大學基礎課部，江蘇 南通226007）

拋物型微分方程的多尺度有限元高效計算

江 山1*，易年余2，孫美玲1，3

（1.揚州大學數(shù)學科學學院，江蘇 揚州225002；2.湘潭大學數(shù)學與計算科學學院，湖南 湘潭411105；3.南通職業(yè)大學基礎課部，江蘇 南通226007）

提出了拋物型微分方程的高效多尺度數(shù)值計算方法.與傳統(tǒng)有限元基函數(shù)相比，多尺度有限元基函數(shù)能更好地反映問題自身的強振蕩微觀信息，結(jié)合多尺度有限元格式，可使計算結(jié)果在宏觀尺度獲得很好的數(shù)值逼近.對時間采用歐拉向后差分離散化，得到穩(wěn)定且收斂的數(shù)值結(jié)果.新方法在取得高仿真逼近的同時，節(jié)約了大量計算資源和時間，因而更具應用價值.

拋物型模型；多尺度有限元；時間離散化；歐拉向后差分；算法效率

奇異攝動微分方程廣泛應用于流體力學、熱電傳導、分子動力學、聲光學等領域［1-2］，其特點是若存在小攝動參數(shù)ε，方程的解就會因在子區(qū)域劇烈變化而產(chǎn)生邊界層現(xiàn)象［3］，使其在許多跨尺度問題中難以有效求解，故奇異攝動微分方程高效數(shù)值解法引起學者們的關(guān)注.Hou等［4］開創(chuàng)性地提出多尺度有限元法（multiscale finite element method，MsFEM），隨后的理論研究和數(shù)值應用不斷取得進展［5-8］.與時間、空間有關(guān)的多尺度拋物型微分方程是當今研究的熱點和難點，人們試圖在計算精度和計算代價之間尋求平衡，致力于探尋既能保證計算精度又可節(jié)約計算資源的新數(shù)值方法來處理復雜的時空多尺度問題.Efendiev等［9］提出了廣義多尺度有限元法，采用在脫機空間計算基函數(shù)、在聯(lián)機空間求解特征值的方法，實現(xiàn)了約化計算.Sun等［10］針對二維反應擴散方程，利用嵌入的多尺度格式并結(jié)合等級網(wǎng)格，僅在邊界層粗單元采用多尺度基，在內(nèi)部光滑粗單元采用傳統(tǒng)有限元基，從而節(jié)約了大量的空間與時間.在本文中，筆者提出時間、空間的多尺度拋物型初邊值問題：

式中u為真解，t為時間，kε（x）為強振蕩系數(shù)，依賴于尺度參數(shù)ε（ε?1），x為二維空間變量，f為右端函數(shù)，g為邊值函數(shù)，u0為初值函數(shù)，Ω為空間域及邊界?Ω.在一定條件下，該問題是適定的，但缺少普遍適用性的解析式，故其高效的數(shù)值解是本文研究的重點.

1　傳統(tǒng)有限元法與多尺度有限元法

1.1變分形式和有限元格式

問題（1）的變分形式是尋求u∈H1（Ω）滿足

將二維區(qū)域Ω分為若干子區(qū)域Ωk，且∪Ωk=Ω.令Kh是矩形網(wǎng)格剖分，在每一個單元K∈Kh定義4個節(jié)點基函數(shù)及其對應的二維坐標［xi，xi+1］×［yj，yj+1］.由有限元的雙線性標準基與等參變換可知：φ1=（1-ξ）（1-η），φ2=ξ（1-η），φ3=ξη，φ4=（1-ξ）η為4個基函數(shù)；ξ=（x-xi）/hx，η=（y-yj）/hy為等參變換；hx=xi+1-xi，hy=yj+1-yj為x，y 方向單元步長，且微分d x=hxdξ，d y=hydη，有

1.2多尺度有限元格式

采用新的多尺度有限元法處理時空拋物型問題時，對應的變分形式是尋求uh∈Uh滿足

其中uh為多尺度有限元解，而多尺度空間Uh=span｛φi：i=1，…，4，?K∈Kh｝是由多尺度基所生成的函數(shù)空間.

針對多尺度基函數(shù)，在粗網(wǎng)格單元K中求解子問題

給定線性邊界條件，即i=j時φi（xj，yj）=1，i≠j時φi（xj，yj）=0，且φi（xj，yj）在?K上呈線性變化.由于子問題（4）與原問題（1）具有相同的微分算子，故在粗網(wǎng)格K 中求得的多尺度基函數(shù)φi能自動反映宏觀單元的微觀信息，如尺度振蕩性（見圖1（a）），而不再是有限元的簡單雙線性基函數(shù)φi（見圖1（b））.

圖1　多尺度基函數(shù)（a）與有限元雙線性基函數(shù)（b）的三維直觀圖Fig.1 The 3D mesh for multiscale basis functions（a）and standard basis functions（b），respectively

傳統(tǒng)有限元法須在非常密的細網(wǎng)格上求解，二維空間全局節(jié)點數(shù)Nglobenode=（Nx×nx+1）×（Ny×ny+1），式中Nx，nx分別為x 方向上粗、細單元剖分數(shù).而多尺度有限元法僅在粗網(wǎng)格上求解，其全局粗節(jié)點數(shù)為Ncoarsenode=（Nx+1）×（Ny+1），因此計算時所需的存儲量和時間都大為縮減.

2　時間相關(guān)的全離散多尺度計算

針對原問題（1）的時間尺度離散，采用更有效的Euler向后差分隱格式

這樣可依次求得不同時間間隔的多尺度數(shù)值解un+1h.

3　數(shù)值實驗

圖2　強振蕩系數(shù)kε（x）Fig.2 The ploted oscillating coefficient kε（x）

圖3　初值u0Fig.3 Initial value u0

計算程序用Matlab2012語言編寫，傳統(tǒng)有限元法在細網(wǎng)格上所需存儲量為O（Nx×nx）2，而多尺度有限元在粗網(wǎng)格上僅需O（Nx）2的存儲量.用R存儲多尺度基函數(shù)及其對應的節(jié)點矩陣，其大小size（R）=Ncoarsenode×Nglobenode，即行、列數(shù)分別為全局粗、細節(jié)點數(shù).細網(wǎng)格的有限元解用u=A-1F求解（計算時間較長）；而用多尺度有限元計算時，利用之前得到的M，A，F(xiàn)，即粗網(wǎng)格上的總質(zhì)量矩陣Mms=R×M×R′，總剛度矩陣Ams=R×A×R′，總載荷向量Fms=R×F，則可用ums=A-1msFms很快算得多尺度有限元在粗網(wǎng)格上的解，再用umsfine=R′×ums還原得到多尺度有限元在細網(wǎng)格上的解，從而方便地比較2種數(shù)值算法的精度.

令每步時間間隔Δt=0.001 s，迭代20步的計算結(jié)果見圖4.結(jié)果表明，迭代20步后有限元解與多尺度有限元解數(shù)值模擬效果非常接近.但前者是在非常密的細網(wǎng)格上運算得到，計算存儲量和計算時間都很大；而后者僅在較粗網(wǎng)格上運算，節(jié)約了大量的計算資源并獲得很好的數(shù)值仿真結(jié)果.

圖4　迭代20步時拋物型問題的傳統(tǒng)有限元解（a）與多尺度有限元解（b）之比較Fig.4 The standard FEM solution（a）and the MsEFM solution（b）for the parabolic problem after 20 steps，respectively

L2范數(shù)的絕對誤差和相對誤差分別為

其中u表示真解，uapprox表示2種方法對應的數(shù)值解，其值的大小可作為度量的依據(jù)（見圖5）.從圖5可以看出拋物型問題的MsFEM解隨時間步增加具有數(shù)值穩(wěn)定性和收斂性的特征.數(shù)值算例很好地驗證了多尺度有限元模擬可以高效處理時空多尺度拋物型微分方程，該方法能在計算精度和計算代價之間取得理想的平衡，具備廣闊的數(shù)值應用優(yōu)勢.

圖5　拋物型多尺度有限元解隨時間變化L2范數(shù)的絕對誤差（a）與相對誤差（b）Fig.5 The L2absolute error（a）and relative error（b）of the MsFEM parabolic solutions varying with time

［1］CAI Xin，CAI Danlin，WU Ruiqian，et al.High accuracy non-equidistant method for singular perturbation reaction-diffusion problem［J］.Appl Math Mech，2009，30（2）：175-182.

［2］SU Fang，XU Zhan，CUI Junzhi，et al.Multiscale computation method for an advection-diffusion equation［J］. Appl Math Comput，2012，218（14）：7369-7374.

［3］DENG Weibing，YUN Xulai，XIE Chunhong.Convergence analysis of the multiscale method for a class of convection-diffusion equations with highly oscillating coefficients［J］.Appl Numer Math，2009，59（7）：1549-1567.

［4］HOU T Y，WU Xiaohui.A multiscale finite element method for elliptic problems in composite materials and porous media［J］.J Comput Phys，1997，134（1）：169-189.

［5］YUE Xingye，E Weinan.The local microscale problem in the multiscale modeling of strongly heterogeneous media：effect of boundary conditions and cell size［J］.J Comput Phys，2007，222（2）：556-572.

［6］JIANG Lijian，PRESHO M.A resourceful splitting technique with applications to deterministic and stochastic multiscale finite element methods［J］.Multiscale Model Simul，2012，10（3）：954-985.

［7］DENG Weibing，WU Haijun.A combined finite element and multiscale finite element method for the multiscale elliptic problems［J］.Multiscale Model Simul，2014，12（4）：1424-1457.

［8］孫美玲，江山，唐元生.多尺度有限元法在Shishkin邊界層的數(shù)值模擬［J］.揚州大學學報：自然科學版，2014，17（2）：25-28，32.

［9］EFENDIEV Y，GALVIS J，HOU T Y.Generalized multiscale finite element methods（GMsFEM）［J］.J Comput Phys，2013，251：116-135.

［10］SUN Meiling，JIANG Shan.Multiscale basis functions for singular perturbation on adaptively graded meshes［J］.Adv Appl Math Mech，2014，6（5）：604-614.

The multiscale finite element computation for
efficiently solving the parabolic differential equation

JIANG Shan1*，YI Nianyu2，SUN Meiling1，3
（1.Sch of Math Sci，Yangzhou Univ，Yangzhou 225002，China；2.Sch of Math&Comput Sci，Xiangtan Univ，Xiangtan 411105，China；3.Dept of Fundam Courses，Nantong Vocat Coll，Nantong 226007，China）

An efficient multiscale finite element computation is proposed to solve the time-space parabolic problems.By comparing the standard finite element basis functions with the multiscale basis functions，the latter has the ability to reflect the local oscillating information，and by the multiscale finite element scheme it may achieve good approximation on the macroscopical scale.For time scale to apply the Euler backward difference discretization，the author demonstrates the stability and convergence by the numerical experiment.This new method obtains the good simulation，and at the same time it saves plenty of computer resource and time，as a consequence it is available for further application values.

parabolic model；multiscale finite element；time discretization；Euler backward difference；algorithm efficiency

O 241.82

A

1007-824X（2015）02-0026-05

（責任編輯 秋 實）

2014-10-28.*聯(lián)系人，E-mail：jiangshan@yzu.edu.cn.

國家自然科學青年基金資助項目（11301462）；江蘇省高校自然科學基金資助項目（13KJB110030）；江蘇省高校研究生科研創(chuàng)新資助項目（KYLX-1332）；揚州大學博士后研究資助項目；揚州大學新世紀人才工程資助項目.

江山，易年余，孫美玲.拋物型微分方程的多尺度有限元高效計算［J］.揚州大學學報：自然科學版，2015，18（2）：26-30.