石劍平,阮麗媛
(昆明理工大學(xué)理學(xué)院系統(tǒng)科學(xué)與應(yīng)用數(shù)學(xué)系,云南 昆明650500)
計(jì)算機(jī)病毒傳播是網(wǎng)絡(luò)安全中最為常見(jiàn)的問(wèn)題,其傳播速度快、范圍廣,帶來(lái)的破壞性也非常大,不加以及時(shí)控制就會(huì)造成重大的經(jīng)濟(jì)損失,嚴(yán)重干擾社會(huì)生產(chǎn)的正常進(jìn)行.因此,對(duì)計(jì)算機(jī)病毒傳播機(jī)理問(wèn)題進(jìn)行研究非常必要,其中對(duì)于病毒傳播模型的研究具有重要的理論和現(xiàn)實(shí)意義.1991年Kephart和White[1]注意到生物病毒傳播與計(jì)算機(jī)病毒傳播存在很多共性,首次將流行病學(xué)數(shù)學(xué)模型引入到計(jì)算機(jī)病毒傳播的研究中.此后,眾多學(xué)者也陸續(xù)使用生物學(xué)中的一些模型和分析方法對(duì)計(jì)算機(jī)病毒進(jìn)行深入的研究.
分?jǐn)?shù)階模型最初是在數(shù)學(xué)基礎(chǔ)理論的研究中提出的[2],因其計(jì)算的復(fù)雜性,在很長(zhǎng)一段時(shí)間沒(méi)有得到廣泛的應(yīng)用.隨著計(jì)算機(jī)技術(shù)的發(fā)展,復(fù)雜的計(jì)算過(guò)程可以通過(guò)各種軟件程序來(lái)完成,分?jǐn)?shù)階模型重新進(jìn)入研究者的視野.由于分?jǐn)?shù)階微積分可以描述帶有遺傳屬性的演化過(guò)程,其所具有的記憶特征不僅可以描述相關(guān)系統(tǒng)過(guò)去的狀態(tài),還會(huì)影響系統(tǒng)現(xiàn)在和未來(lái)的狀態(tài).近年來(lái),數(shù)學(xué)、物理學(xué)、生物學(xué)、金融學(xué)等領(lǐng)域的研究者在各自的研究工作中均發(fā)現(xiàn),某些問(wèn)題采用分?jǐn)?shù)階模型來(lái)描述更能恰當(dāng)?shù)姆从硢?wèn)題的真實(shí)情形.
“生態(tài)記憶”是生物種群模型的重要特征,它指的是“過(guò)去的狀態(tài)或經(jīng)驗(yàn)影響群落目前或未來(lái)反應(yīng)的能力”[3],這一特征與分?jǐn)?shù)階微積分所能表現(xiàn)的特性不謀而合.很多研究者將分?jǐn)?shù)階微積分方法引入生物數(shù)學(xué)的各種模型中,獲得了一批有理論價(jià)值和實(shí)際意義的結(jié)果.例如,Djordjevic等[4]采用分?jǐn)?shù)階常微分方程對(duì)一類細(xì)胞問(wèn)題建模,模型與實(shí)際數(shù)據(jù)的吻合度比整數(shù)階的模型更為準(zhǔn)確;EI-Saka等[5]研究分?jǐn)?shù)階變動(dòng)人口SIS 流行病模型,分析了平衡點(diǎn)的穩(wěn)定性及一致穩(wěn)態(tài)解的存在性條件;原三領(lǐng)等[6]則研究了一類分?jǐn)?shù)階SIR流行病模型的穩(wěn)定性問(wèn)題.基于流行病學(xué)的計(jì)算機(jī)病毒模型,因其演化過(guò)程和控制策略具有很強(qiáng)的“記憶”屬性,采用分?jǐn)?shù)階微積分理論進(jìn)行研究,可以獲得諸如混沌等新的動(dòng)力學(xué)特征[7].
本文將分?jǐn)?shù)階導(dǎo)數(shù)引入一類計(jì)算機(jī)病毒模型[8],考慮了三個(gè)因素對(duì)模型進(jìn)行改進(jìn),且基于感染過(guò)程中病毒的潛伏性,提出了一個(gè)改進(jìn)的時(shí)滯分?jǐn)?shù)階計(jì)算機(jī)病毒模型.采用文[9]總結(jié)的結(jié)論,從理論上分析模型的正平衡點(diǎn)處出現(xiàn)Hopf分支與時(shí)滯之間的聯(lián)系.為了驗(yàn)證理論分析的合理性,選擇恰當(dāng)?shù)南到y(tǒng)參數(shù)進(jìn)行數(shù)值模擬,并通過(guò)數(shù)值模擬獲得分?jǐn)?shù)階變化的情況下時(shí)滯臨界值的變化趨勢(shì).
文[8]研究了一類計(jì)算機(jī)病毒模型(SIR模型):
其中,S(t),I(t),R(t)分別表示易感染節(jié)點(diǎn)、已感染節(jié)點(diǎn)和免疫節(jié)點(diǎn)的個(gè)數(shù),0<α,β,γ,μ,θ <1是各種節(jié)點(diǎn)之間相互轉(zhuǎn)化的概率.文[8]從信息傳播機(jī)理出發(fā),以上述系統(tǒng)作為信息傳播的模型給出了信息傳播的規(guī)則,建立了相應(yīng)的平均場(chǎng)方程,并通過(guò)計(jì)算平衡點(diǎn)和基本再生數(shù)R0,從理論上證明了平衡點(diǎn)的漸近穩(wěn)定性.
為了更準(zhǔn)確地描述計(jì)算機(jī)網(wǎng)絡(luò)系統(tǒng)病毒感染演化的過(guò)程,本文將分?jǐn)?shù)階導(dǎo)數(shù)引入模型(2.1),并進(jìn)一步考慮三個(gè)因素對(duì)模型進(jìn)行改進(jìn):第一,某些感染節(jié)點(diǎn)因系統(tǒng)的差異性,被感染后無(wú)法正常工作,因而退出整個(gè)網(wǎng)絡(luò)系統(tǒng);第二,有些易感節(jié)點(diǎn)因?yàn)檐浖牟町惗哂邢忍烀庖?第三,網(wǎng)絡(luò)系統(tǒng)根據(jù)需要按一定的比例添加新的節(jié)點(diǎn).由于某些計(jì)算機(jī)病毒的發(fā)作存在潛伏期,本文在易感染節(jié)點(diǎn)被感染的過(guò)程加入了時(shí)滯項(xiàng).基于上述分析,我們建立以下時(shí)滯分?jǐn)?shù)階SIR模型:
其中,q(0 模型(2.2)的平衡點(diǎn)為: 其中 根據(jù)上述模型假設(shè),系統(tǒng)(2.2)是一個(gè)開(kāi)放模型,在系統(tǒng)演化過(guò)程中有新的節(jié)點(diǎn)加入并有節(jié)點(diǎn)退出,所有的系統(tǒng)參數(shù)均取正值.顯然平衡點(diǎn)E2的分量均不為零,此改進(jìn)模型除了零平衡點(diǎn)和正平衡點(diǎn)外不存在無(wú)病平衡點(diǎn),由于零平衡點(diǎn)的情況無(wú)實(shí)際意義,故下文僅研究正平衡點(diǎn).為敘述方便,不妨假設(shè)正平衡點(diǎn)為E?=(x?,y?,z?),在E?點(diǎn)對(duì)系統(tǒng)(2.2)做平移變換,變換公式為: 代入模型(2.2)得 研究分?jǐn)?shù)階微分方程系統(tǒng)的穩(wěn)定性和Hopf分支問(wèn)題,目前主要是基于其對(duì)應(yīng)的線性化近似系統(tǒng)在平衡點(diǎn)處的穩(wěn)定性來(lái)進(jìn)行的.文[9]給出定義,說(shuō)明系統(tǒng)的時(shí)滯τ滿足以下三個(gè)條件,則在平衡點(diǎn)處會(huì)出現(xiàn)Hopf分支,即平衡點(diǎn)的穩(wěn)定性可以由τ值來(lái)確定. 定義2.1[9]對(duì)于以下n-維分?jǐn)?shù)階時(shí)滯系統(tǒng): 0<α ≤1 和τ ≥0.如果系統(tǒng)(2.4)滿足以下三個(gè)條件: 1)當(dāng)τ= 0,(2.4)在平衡點(diǎn)x?= (x?1,x?2,...,x?n)處的線性化系統(tǒng)的系數(shù)矩陣記為A,A的所有特征根λj(j=1,2,...,n)必須都滿足 2)當(dāng)τ=τ0時(shí),(2.4)在平衡點(diǎn)x?= (x?1,x?2,...,x?n)處的線性化系統(tǒng)的特征方程存在一對(duì)純虛根s(τ0)=±iω0; 3)橫截條件 成立,其中s(τ)是(2.4)在平衡點(diǎn)x?=(x?1,x?2,...,x?n)處的線性化系統(tǒng)的特征方程的根,Re[·]指提取復(fù)數(shù)的實(shí)部. 則當(dāng)時(shí)滯τ=τ0時(shí),系統(tǒng)(2.4)在平衡點(diǎn)x?=(x?1,x?2,...,x?n)處出現(xiàn)Hopf分支. 注1定義2.1中,條件1)說(shuō)明當(dāng)時(shí)滯τ=0時(shí),分?jǐn)?shù)階微分方程系統(tǒng)(2.4)在0<α<1時(shí),分?jǐn)?shù)階數(shù)α可擴(kuò)大模型的穩(wěn)定區(qū)域,即當(dāng)特征值穿過(guò)虛軸時(shí),系統(tǒng)仍然可能穩(wěn)定[10?11]. 下面將按照上述三個(gè)條件分析系統(tǒng)(2.2)正平衡點(diǎn)的穩(wěn)定性. 將系統(tǒng)(2.3)在平衡點(diǎn)(0,0,0)處線性化,得到對(duì)應(yīng)的線性化系統(tǒng): 其中: Ⅰ 時(shí)滯為零時(shí)系數(shù)矩陣的特征根分析 當(dāng)τ=0時(shí),線性化系統(tǒng)(3.1)的系數(shù)矩陣對(duì)應(yīng)的特征方程為 若特征方程(3.2)式滿足Routh-Hurwritz條件,則所有的特征根均有負(fù)實(shí)部,因此,在給定參數(shù)條件下,只需驗(yàn)證以下三個(gè)不等式成立,即可保證定義2.1中的條件1)成立: Ⅱ 確定出現(xiàn)Hopf分支的時(shí)滯臨界值 對(duì)(3.1)式做拉普拉斯變換[12]得到以下方程組: (3.4)式的特征方程為: 分離(3.5)式的時(shí)滯項(xiàng)和非時(shí)滯項(xiàng)可得: 其中 定義2.1的條件2)指出,要使分?jǐn)?shù)階系統(tǒng)在時(shí)滯臨界值τ0處出現(xiàn)Hopf分支,其特征方程的根s必須是純虛根.為了找到τ0,不妨假設(shè)s= iω=ω(cos(π/2)+i sin(π/2)),代入(3.6)式,分離實(shí)部與虛部可得到以下方程組: 其中M1,M2,N1和N2如下所示: 易求出 聯(lián)立公式sin2(τω)+cos2(τω)=1可得到: 求解上式可得ω.設(shè)ωi(i=0,1,2,...)為ω中所有的正值,代入(3.9)式計(jì)算sin(τω)和cos(τω),可得到τ的計(jì)算公式如下所示: (I)如果W1(ωi)>0,W2(ωi)>0,則 (II)如果W1(ωi)<0,W2(ωi)>0,則 (III)如果W1(ωi)>0,W2(ωi)<0,則 (IV)如果W1(ωi)<0,W2(ωi)<0,則 按照時(shí)滯的實(shí)際意義,我們僅考慮其最小正值,取 根據(jù)定義,τ=τ?0即為系統(tǒng)(2.3)出現(xiàn)Hopf分支的時(shí)滯臨界值,對(duì)應(yīng)于τ?0的ω值記為ω?0=ωi. Ⅲ 驗(yàn)證橫截條件 在文[9]中,通過(guò)對(duì)特征方程兩邊求τ的一階偏導(dǎo)數(shù),推導(dǎo)了的計(jì)算公式,但此公式比較復(fù)雜,由此給出的判斷橫截條件的假設(shè)較為嚴(yán)格.文[13]通過(guò)構(gòu)造輔助函數(shù),提出了一種更為恰當(dāng)?shù)呐袛喾椒?下面采用文[13]提出的判別方法驗(yàn)證橫截條件. 將s=iω帶入(3.6)式,則不妨令 構(gòu)造以下的輔助函數(shù): 定義Jacobi行列式為: 綜上分析,我們得到以下結(jié)果: 定理3.1對(duì)于系統(tǒng)(2.3),在給定參數(shù)組(α,β,γ,μ,θ,a,b,c)的情況下,如果以下假設(shè)成立: (H1)(3.3)式成立; (H2)當(dāng)τ=τ?0時(shí),系統(tǒng)(2.3)在平衡點(diǎn)(0,0,0)處的線性化系統(tǒng)的特征方程(3.5)存在一對(duì)純虛根s=±iω?0; (H3)J(ω?0,τ?0)>0,J(ω,τ)為(3.11)式和(3.12)式所定義的Jacobi行列式.則當(dāng)時(shí)滯τ=τ?0時(shí),系統(tǒng)(2.3)在平衡點(diǎn)(0,0,0)處出現(xiàn)Hopf分支. 注2由定理3.1即可推出,當(dāng)時(shí)滯τ=τ?0時(shí),系統(tǒng)(2.2)在正平衡點(diǎn)E?= (x?,y?,z?)處出現(xiàn)Hopf分支. Ⅰ 平衡點(diǎn)的穩(wěn)定性模擬 為了驗(yàn)證上述理論分析結(jié)果的正確性,選擇一組恰當(dāng)?shù)膮?shù)α=0.2,β=0.2,γ=0.4,μ=0.6,θ= 0.3,a= 0.05,b= 0.06,c= 0.05代入模型(2.2),正平衡點(diǎn)E?= (0.75,1.65,1.865625),分?jǐn)?shù)階取值為q= 0.9,計(jì)算得時(shí)滯臨界值為τ?0= 0.981197,對(duì)應(yīng)的ω?0= 1.620917262,容易驗(yàn)證(3.3)式成立,且J(ω?0,τ?0)>0.為了更好的展示理論分析的結(jié)果,我們選擇兩種情況,借助于文[14]給出的預(yù)估校正法做數(shù)值模擬.一種情況為τ= 0.99>τ?0= 0.981197,模擬結(jié)果見(jiàn)圖4.1和圖4.2,此時(shí)正平衡點(diǎn)出現(xiàn)Hopf分支,是不穩(wěn)定的.另一種情況為τ= 0.93<τ?0=0.981197,模擬結(jié)果見(jiàn)圖4.3和圖4.4,顯然此時(shí)正平衡點(diǎn)是漸近穩(wěn)定的. 圖4.1 模型(2.2)狀態(tài)變量的波形圖,初值為(1.5,0.8,1.75), q =0.9, τ =0.99>0.981197,正平衡點(diǎn)不穩(wěn)定 圖4.2 初值為(1.5,0.8,1.75), q =0.9, τ =0.99>0.981197,模型(2.2)在正平衡點(diǎn)處出現(xiàn)Hopf分支 圖4.3 模型(2.2)狀態(tài)變量的波形圖,初值為(0.822,1.25,2.1), q = 0.9, τ = 0.93 < 0.981197,正平衡點(diǎn)是漸近穩(wěn)定的 圖4.4 初值為(0.822,1.25,2.1), q =0.9, τ =0.93<0.981197,模型(2.2)的正平衡點(diǎn)是漸近穩(wěn)定的 Ⅱ 分?jǐn)?shù)階對(duì)時(shí)滯臨界值的影響 通過(guò)上述理論分析和數(shù)值模擬可知,對(duì)于有病毒傳染的計(jì)算機(jī)網(wǎng)絡(luò)系統(tǒng),在模型參數(shù)給定的情況下,系統(tǒng)從不穩(wěn)定狀態(tài)進(jìn)入穩(wěn)定狀態(tài),病毒潛伏到發(fā)作的時(shí)間需要控制在一定的時(shí)間范圍τ <τ?0.反之,在知道病毒的潛伏期,即已知τ的情況下,如果需要控制系統(tǒng)的穩(wěn)定狀態(tài),就需要對(duì)系統(tǒng)進(jìn)行干擾以改變時(shí)滯臨界值τ?0,即改變穩(wěn)定域的大小.此時(shí),需要對(duì)影響τ?0的因素進(jìn)行分析. 影響時(shí)滯臨界值的因素非常多,從時(shí)滯臨界值的計(jì)算公式來(lái)看,模型參數(shù)、分?jǐn)?shù)階都是影響時(shí)滯臨界值大小,繼而影響系統(tǒng)穩(wěn)定性的直接因素.本文以分?jǐn)?shù)階為參數(shù),不改變前面模型(2.2)參數(shù)的取值,通過(guò)數(shù)值模擬獲得模型導(dǎo)數(shù)從整數(shù)階到分?jǐn)?shù)階、分?jǐn)?shù)階不斷減小的過(guò)程中,時(shí)滯臨界值的具體變化情況,詳見(jiàn)表4.1. 表4.1說(shuō)明,模型的導(dǎo)數(shù)從1到0.8的變化過(guò)程中,時(shí)滯臨界值τ?0隨著q值的減小不斷增加.在0.8至0.75之間的某一點(diǎn)處,τ?0的變化發(fā)生轉(zhuǎn)折,不再隨著q值的減小而增加,反而隨著q值的減小而不斷變小.由此可見(jiàn),通過(guò)改變分?jǐn)?shù)階的大小改變時(shí)滯臨界值,是可以影響系統(tǒng)的穩(wěn)定域的. 表4.1 分?jǐn)?shù)階的變化對(duì)ω?0 和τ?0 的影響 本文根據(jù)病毒的演變特征和潛伏特征,引入分?jǐn)?shù)階導(dǎo)數(shù)和時(shí)滯項(xiàng)對(duì)一類計(jì)算機(jī)病毒模型進(jìn)行了推廣.采用線性化方法和拉普拉斯變換對(duì)模型進(jìn)行轉(zhuǎn)換,通過(guò)分析線性化系統(tǒng)的系數(shù)矩陣特征方程和系統(tǒng)特征方程,分析了時(shí)滯和分?jǐn)?shù)階系統(tǒng)平衡點(diǎn)穩(wěn)定性之間的關(guān)系,推導(dǎo)了Hopf分支出現(xiàn)時(shí)時(shí)滯臨界值的計(jì)算公式.為了驗(yàn)證理論分析的合理性,選擇恰當(dāng)?shù)南到y(tǒng)參數(shù)進(jìn)行數(shù)值模擬,模擬結(jié)果不僅驗(yàn)證了結(jié)論的正確性,還揭示了分?jǐn)?shù)階對(duì)時(shí)滯臨界值的影響趨勢(shì).本文的研究結(jié)果一定程度上推廣了計(jì)算機(jī)病毒模型的研究范圍和方法.3.Hopf分支分析
4.數(shù)值模擬
5.結(jié)語(yǔ)