分裂的Euler-Maruyama誤差修正法

殷政偉, 王天軍(河南科技大學 數學與統計學院, 河南 洛陽 471003)

?

分裂的Euler-Maruyama誤差修正法

殷政偉, 王天軍
(河南科技大學 數學與統計學院, 河南 洛陽 471003)

在Euler-Maruyama方法基礎上提出一種新的顯式的分裂法,該方法可以用來求解隨機常微分方程的數值解.與Euler-Maruyama方法相比,該方法具有較大的均方穩定域,因此該方法可以用來求解一些剛性隨機微分方程.最后,通過數值算例說明該方法的有效性及可實現性.

隨機微分方程數值解; 誤差修正法; 分裂法; Euler-Maruyama方法

隨機微分方程由于考慮了隨機噪聲對系統變化的影響,因此,較確定性微分方程更能準確地刻畫、描述現實生活中的各種現象.隨機微分方程模型在眾多應用領域,如生物、化學、微電子學、經濟、金融等領域都有廣泛的應用[1-5].1951年,日本數學家K. It發表了劃時代的著作[6],給出了隨機微分方程(SDE)的嚴格數學描述,建立了隨機積分,創立了隨機微分方程理論.此后,隨機微分方程理論得到快速發展.然而對于隨機微分方程而言,它的解析解是非常少見的,因此對于它的數值解的研究具有非常重要的現實意義和應用價值.

隨機微分方程的數值方法[7]主要分為2大類,一類是顯式方法,例如經典的Euler-Maruyama方法[8]、Milstein方法[9]和Runge-Kutta方法[10].顯式方法的特點是計算量較小,計算方便,容易實現,只要給定初值后便可以直接計算.但顯式方法的穩定性較差,在遇到剛性問題時,對步長的要求極為嚴格.相比較而言,隱式方法具有較好的穩定性,如平衡法[11]、強1階平衡Milstein方法[12]、分步平衡Milstein方法[13]、隱式Taylor方法[14]、全隱的1階和1.5階的Runge-Kutta型方法[15]等.隱式方法雖然在求解剛性隨機微分方程時有著較好的穩定性,但是由于是隱格式,每一步計算時都要求解一個大型的非線性方程組,這對于一些較復雜的問題,計算量無疑要增加很多,甚至不一定能夠實現,因此構造一些具有較好穩定性的顯式數值方法就顯得尤為重要.在這方面,近年來,也有一些成果,可參考文獻[16-17].

本文主要研究自治的隨機微分方程,其形式為

(1)

其中,T是一個正常數,f:Rd→Rd是漂移項系數,g=(g1,g2,…,gm):Rd→Rd×m是擴散項系數.w(t)是定義在完備概率空間(Ω,F,{Ft}t≥0,P)上的一個m維的維納過程,概率空間的濾子{Ft}t≥0滿足通常的條件,即{Ft}t≥0是遞增右連續的,且F0包含所有的零測集.初值x0與維納過程w(t)是相互獨立的,且滿足E|X0|2<∞.

受到文獻[18]構造求解常微分方程的誤差校正法思想的啟發,提出一種分裂的Euler-Maruyama誤差修正法,該方法是一個顯式方法,但卻有較好的穩定性.

1 Euler誤差修正法

Euler誤差修正法是求解常微分方程初值問題

(2)

的一種顯式方法,收斂階較高,且穩定域較大,所以可以用來求解剛性微分方程.Euler誤差修正法的原理是利用切比雪夫配點法求解Euler多項式與常微分方程初值問題(2)的準確解的差所滿足的微分方程,然后將所求得的解,即修正項加到Euler法上即可.下面簡要介紹Euler誤差修正法的推導過程.

在區間[tn,tn+1]上,設y(t)為Euler多項式

y(t)=yn+(t-tn)f(tn,yn).

(3)

定義φ(t)為方程(2)的準確解x(t)與Euler多項式y(t)的差,即

φ(t)=x(t)-y(t),

則φ(t)滿足如下的線性微分方程

φ′(t)=G(t)φ(t)+F(t),

t∈[tn,tn+1],

(4)

F(ts)),s∈[-1,1],

(5)

ξ.

令C=max|f″xx(t,χ)|,可以得出

).

因此(5)式可以轉化為

sj=cosm-jmπ,j=0,1,…,m

(7)

d=[β1,β2,…,βm]T,

F=[F(ts1),F(ts2),…,F(tsm)]T,

yn+1=yn+hf(tn,yn)+βm.

(8)

2 分裂的Euler-Maruyama誤差修正法及其穩定性

將第一節中所提到的Euler誤差修正法稍做變動應用到隨機微分方程(1),得出分裂的Euler-Maruyama誤差修正法.在離散的小區間[tn,tn+1]上,該方法的形式如下

(9)

其中,yn是方程(1)的解析解x(t)在tn時刻的近似值,h是步長,βm是Euler法求解常微分方程初值問題時誤差的修正項.

同常微分方程一樣,隨機微分方程數值解法的穩定性也是通過一線性測試方程進行分析.考慮如下形式的線性標量測試方程

(10)

其中λ,μ∈R是常數,x0依概率1不等于零,方程(10)的準確解為

).

用單步方法求解測試方程(10)得到隨機差分方程

(11)

定義 1 對于固定的λ,μ,h,一種數值方法稱為均方穩定的,如果

引理 1[18]將方法(8)應用到Dahlquist測試方程x′(t)=λx(t)可以得到差分方程

yn+1=Sm(p)yn,n≥0,

p=λh,其中穩定函數Sm(p)可以表示為

注 1 當m=1,2,3,4時,穩定函數Sm(p)可分別表示為

定理 1 對于固定的λ、μ、h,分裂的Euler-Maruyama誤差修正法是均方穩定的,如果有

Sm(p)+q2<1,

證明 由引理1可知,將分裂的Euler-Maruyama誤差修正法應用到測試方程(10)可得

yn+1=Sm(p)yn+μyn△wn,

為了更直觀的看到分裂的Euler-Maruyama誤差修正法的均方穩定域,利用Matlab軟件分別畫出2種方法的均方穩定域.圖1中圓內部是Euler-Maruyama方法的均方穩定區域,另一條曲線內部是m=2時分裂的Euler-Maruyama誤差修正法的均方穩定區域,其中.從圖1中可以明顯看出,經過誤差修正后的分裂Euler-Maruyama方法的均方穩定區域比Euler-Maruyama方法的均方穩定域增大了很多.

3 數值實驗

取m=2,用分裂的Euler-Maruyama誤差修正法求解測試方程(10),其中λ=2,μ=1,x0=1.圖2中曲線是測試方程(10)的精確解,帶*曲線是用分裂的Euler-Maruyama誤差修正法求得的數值解,從圖2可以看出,數值解與真解吻合得很好.

4 結語

本文在切比雪夫配點法的基礎上,提出一種新的求解隨機常微分方程數值解的顯式方法.該方法在求解隨機常微分方程時,算法簡單,計算量小,且與同階的Euler-Maruyama方法相比,該方法的均方穩定域增大很多,因此也可以考慮用該方法求解一些具有剛性的隨機常微分方程.

[1] Oksendal B. Stochastic Differential Equations:an Introduction with Applications[M]. 5th ed. Berlin:Springer-Verlag,1998.

[2] Mao X. Stochastic Differential Equations and Applications[M]. Chichester:Horwood,1997.

[3] Sobczyk K. Stochastic Differential Equations with Application to Physics and Engineering[M]. Dordrecht:Kluwer Academic,1991.

[4] Hobson D, Rogers L. Complete models with stochastic volatility[J]. Math Finance,1998,8:27-48.

[5] 楊鵬. 基于確定繳費型養老金最優投資的隨機微分博弈[J]. 四川師范大學學報:自然科學版,2015,38(2):194-200.

[7] Kloeden P E, Platen E. Numerical Solution of Stochastic Differential Equations[M]. Berlin:Springer-Verlag,1992.

[8] Maruyama G. Continuous Markov processes and stochastic equations[J]. Rend Circ Mat Palermo,1955,4:48-90.

[9] Milstein G N. Approximate integration of stochastic differential equations[J]. Theory Prob Appl,1974,19(3):557-562.

[10] Rumelin W. Numerical treatment of stochastic differential equations[J]. SIAM J Numer Anal,1982,19(3):604-613.

[11] Milstein G N, Platen E, Schurz H. Balanced implicit methods for stiff stochastic systems[J]. SIAM J Numer Anal,1998,35(3):1010-1019.

[12] Kahl C, Schurz H. Balanced Milstein methods for ordinary SDEs[J]. Monte Carlo Meth Appl,2006,12(2):143-170.

[13] Wang P, Liu Z. Split-step backward balanced Milstein methods for stiff stochastic systems[J]. Appl Numer Math,2009,59:1198-1213.

[14] Tian T, Burrage K. Implicit Taylor methods for stiff stochastic differential equations[J]. Appl Numer Math,2001,38:167-185.

[15] Burrage K, Tian T. Implicit stochastic Runge-Kutta methods for stochastic differential equations[J]. Bit Numerical Math,2004,44:21-39.

[16] Abdulle A, Cirilli S. Stabilized methods for stiff stochastic systems[J]. C R Acad Sci,2007,345(10):593-598.

[17] Abdulle A, Cirilli S. S-ROCK methods for stiff stochastic problems[J]. SIAM J Sci Comput,2008,30(2):997-1014.

[18] Kim P, Piao X, Kim S D. An error corrected Euler method for solving stiff problems based on Chebyshev collocation[J]. SIAM J Numer Anal,2011,49(6):2211-2230.

2010 MSC：34A05; 65L05

(編輯 鄭月蓉)

An Error Corrected Split Euler-Maruyama Method

YIN Zhengwei, WANG Tianjun
(SchoolofMathematicsandStatistics,HenanUniversityofScienceandTechnology,Luoyang471003,Henan)

In this paper, a new explicit split method is proposed based on Euler-Maruyama method. This method can be applied to solve the stochastic differential equations. Compared to Euler-Maruyama method, this method has larger mean-square stability region. Therefore, the method can be used to solve stiff stochastic differential equations. Some numerical experiments illustrate the effectiveness and operability.

numerical solution of stochastic differential equations; error corrected method; split method; Euler-Maruyama method

2014-07-08

國家自然科學基金(11371123)和教育部博士點基金(20120162120096)

殷政偉(1980—),女,博士生,主要從事隨機微分方程數值解的研究,E-mail:yzw_0379@163.com

O241.8

A

1001-8395(2015)06-0830-04

10.3969/j.issn.1001-8395.2015.06.007