李 聰，武芳芳

(沈陽工業(yè)大學 理學院，沈陽 110870)

格子Boltzmann方法[1-2]自提出之日起，便迅速發(fā)展成為流體力學、數(shù)學、物理學等相關(guān)領域內(nèi)的研究熱點.格子Boltzmann方法基于分子動力學理論而建立，與傳統(tǒng)數(shù)值方法不同，它具有清晰的物理背景，其時間、空間完全獨立，邊界易于處理，是一種完全并行化的介觀模擬算法.格子Boltzmann方法的獨特優(yōu)勢使其在求解各類非線性偏微分方程數(shù)值解[3-7]方面克服了許多傳統(tǒng)算法的不足，受到眾多專家學者的廣泛關(guān)注，發(fā)展前景日益廣闊.

期權(quán)[8-9]是金融領域中一種風險管理的核心工具，多資產(chǎn)期權(quán)作為一種奇異期權(quán)[10]，由多個風險資產(chǎn)組成.美式多資產(chǎn)看跌期權(quán)滿足下列線性互補模型[8]：

(1)

終值條件為V(S,T)=φ(S)，其中

1 有界區(qū)域上的定價模型

通過變量替換t=T-τ，xi=lnSi，u(xi,τ)=V(Si,t)(i=1,2)和懲罰法[17]將方程(1)轉(zhuǎn)化為下列等價方程

(2)

假設截斷區(qū)域是Ω=[-l1,l1]×[-l2,l2]，其中l(wèi)i(i=1,2)是未知變量，對于固定時間t∈(0,T]，固定點(x1,x2)∈R2以及給定充分小的數(shù)ζ∈(0,1)，由遠場截斷法[18-19]可知截斷長度lmax滿足以下關(guān)系

(3)

且有l(wèi)max≤li，K表示敲定價格.

2 格子Boltzmann模型

格子Boltzmann模型是基于DnQb系列速度模型(其中n表示空間維數(shù)，b表示離散速度)的單松弛格子BGK(Bhatnagar-Gross-Krook)模型[20].方程(2)的格子Boltzmann模型的局部粒子分布函數(shù)演化方程為

(4)

對方程(4)做Taylor展開，有

(5)

再實施以下多尺度Chapman-Enskog展開[22]，

(6)

其中ε是Knudsen數(shù)[23].將式(6)代入方程(5)中，對比ε的各階系數(shù)，可得

(7)

(8)

(9)

(10)

定義宏觀量u(x,t)為該點各個方向的局部粒子分布函數(shù)的和，即

(11)

由于質(zhì)量守恒，故局部平衡態(tài)分布函數(shù)滿足

(12)

因此由方程(7)和式(12)，有

(13)

為恢復宏觀方程(2)，對局部平衡態(tài)分布函數(shù)做如下限制，

(14)

同時對補償函數(shù)做如下限制，

(15)

由式(6)(15)，可得

(16)

(17)

對方程(8)關(guān)于α求和，并結(jié)合式(13)(14)(16)和(17)有

(18)

對方程(9)關(guān)于α求和，并結(jié)合式(13)(14)(16)和(17)有

(19)

對方程(10)關(guān)于α求和，并結(jié)合式(13)(14)(16)和(17)有

(20)

2.1 D2Q5模型

對方程(2)采用D2Q5速度模型，其離散速度集合為{eα,α=0,1,…,4}={(0,0),(c,0),(-c,0),(0,c),(0,-c)},由式(18)×ε+(19)×ε2可以恢復得到具有二階精度的宏觀方程

(21)

為恢復方程(2)，令

則方程(21)可變?yōu)?/p>

由式(13)(14)可得一組平衡態(tài)分布函數(shù)為

由式(15)(16)(17)可得一組補償函數(shù)為

2.2 D2Q9模型

對方程(2)采用D2Q9速度模型，其離散速度集合為{eα,α=0,1,…,8}={(0,0),(c,0),(-c,0),(0,c),(0,-c),(c,c),(-c,c),(-c,-c),(c,-c)}，由式(18)×ε+(19)×ε2+(20)×ε3可以恢復得到具有三階精度的宏觀方程

(22)

為恢復方程(2)，令

則方程(22)可變?yōu)?/p>

由式(13)(14)可得一組平衡態(tài)分布函數(shù)為

由式(15)(16)(17)可得一組補償函數(shù)為

3 數(shù)值模擬

為驗證本文算法的有效性，考慮1年期美式多資產(chǎn)期權(quán)定價模型，即T=1，取敲定價格K=1，σ1=0.2，σ2=0.3，q1=0.05，q2=0.01，在線性互補模型(1)中取終值條件φ(S1,S2)=max{K-α1S1-α2S2,0},其中αi≥0(i=1,2)表示第i個標的資產(chǎn)所占比例，且有α1+α2=1，故取α1=0.6，α2=0.4.在數(shù)值模擬過程中，初始粒子分布函數(shù)取平衡態(tài)分布函數(shù)，采用非平衡外推格式[24]處理邊界條件，同時定義全局誤差為

其中u(xi,xj,t)表示數(shù)值解，u*(xi,xj,t)表示精確解.

在式(3)中取ζ=10-5可得截斷長度lmax=2.082 5，為敘述方便取l1=l2=lmax/2，故截斷區(qū)域為Ω∈[-1,1]×[-1,1].

在數(shù)值模擬過程中，時間間隔取Δt=0.001，空間步長取Δx=2/N，N為格子數(shù)，松弛時間取τ=2.5，懲罰法系數(shù)取η=0.01.取格子數(shù)N=200.圖1和圖2分別描繪了D2Q5模型和D2Q9模型在3組無風險利率r=0.05，r=0.10和r=0.20下T=1時刻的期權(quán)價格.參照文獻[15]中有限元法得到的期權(quán)價格圖像，圖1和圖2中期權(quán)價格圖形與其形狀相吻合，證實了本文所建立格子Boltzmann模型的有效性.

(a) r=0.05 (b) r=0.10 (c) r=0.20圖1 D2Q5速度模型下的期權(quán)價格Fig.1 Option prices of LBM under D2Q5 model

(a) r=0.05 (b) r=0.10 (c) r=0.20圖2 D2Q9速度模型下的期權(quán)價格Fig.2 Option prices of LBM under D2Q9 model

由于美式多資產(chǎn)期權(quán)定價問題不存在顯式解，取D2Q9速度模型在格子數(shù)N=800以及T=1時刻下的數(shù)值解作為精確解u*(xi,xj,t)的近似.表1列出了D2Q5速度模型和D2Q9速度模型在不同無風險利率下的全局誤差.圖3描繪了在不同無風險利率下的全局誤差圖像，圖像表明在同一速度模型下格子數(shù)越大，誤差越小；在相同無風險利率下D2Q9模型比D2Q5模型誤差更小.

表1 D2Q5和D2Q9速度模型的全局誤差Tab.1 Global relative errors under D2Q5 model and D2Q9 model

(a) r=0.05 (b) r=0.10 (c) r=0.20圖3 D2Q5和D2Q9速度模型下的全局誤差圖Fig.3 Diagrams of Global relative errors under D2Q5 model and D2Q9 model

4 結(jié)論

本文實現(xiàn)了利用格子Boltzmann方法求解美式多資產(chǎn)期權(quán)定價問題.利用懲罰法和遠場截斷法將無界區(qū)域上的線性互補模型轉(zhuǎn)化為有界區(qū)域上的拋物問題，針對轉(zhuǎn)化后的模型進行多尺度Chapman-Enskog展開等一系列計算，同時確定不同速度模型下的平衡態(tài)分布函數(shù)和補償函數(shù)，最終計算得到期權(quán)價格.數(shù)值模擬結(jié)果表明，在相同的速度模型下隨著格子數(shù)的增加，誤差越小；隨著速度方向的增加，誤差越?。徊⑴c其他數(shù)值方法做對比，期權(quán)價格圖形吻合.證實了格子Boltzmann方法對于研究美式多資產(chǎn)期權(quán)定價問題是一種有效算法,該算法為求解美式多資產(chǎn)期權(quán)定價問題提供了新思路.