沈敏捷，曾振柄，林 望，楊爭(zhēng)峰

(1.上海市高可信計(jì)算重點(diǎn)實(shí)驗(yàn)室(華東師范大學(xué))，上海200062; 2.上海大學(xué) 數(shù)學(xué)系，上海200444;3.溫州大學(xué)數(shù)學(xué)與信息科學(xué)學(xué)院，浙江溫州325035)

(*通信作者電子郵箱zbzeng@shu.edu.cn)

0 引言

嵌入式系統(tǒng)是一種嵌入機(jī)械或電氣系統(tǒng)內(nèi)部，具有專一功能和實(shí)時(shí)計(jì)算性能的計(jì)算機(jī)系統(tǒng)，被廣泛應(yīng)用于電信、交通、商業(yè)、工業(yè)、醫(yī)療等重要領(lǐng)域。區(qū)別于通用計(jì)算機(jī)，嵌入式系統(tǒng)是為某些特定任務(wù)設(shè)計(jì)的。其中一些系統(tǒng)具有低功耗低成本的要求，因此對(duì)其軟硬件資源要求較為嚴(yán)格。另一些系統(tǒng)對(duì)實(shí)時(shí)性有很高的要求，如數(shù)控機(jī)床、列車控制系統(tǒng)、飛機(jī)導(dǎo)航系統(tǒng)等。一旦系統(tǒng)的響應(yīng)時(shí)間未能達(dá)到要求，可能會(huì)造成嚴(yán)重的后果。

在嵌入式實(shí)時(shí)系統(tǒng)中，廣泛存在一類既包含連續(xù)動(dòng)態(tài)行為又包含離散事件行為，且兩者互相作用的系統(tǒng)。其動(dòng)力學(xué)行為既隨時(shí)間而連續(xù)變化，又受事件而離散驅(qū)動(dòng)。以列車自動(dòng)控制(Automatic Train Control，ATC)系統(tǒng)［1－2］為例，ATC 系統(tǒng)是一個(gè)完整的列車速度監(jiān)督系統(tǒng)，它可以控制列車行進(jìn)的速度。當(dāng)列車行駛速度超過容許的速度時(shí)，剎車設(shè)備應(yīng)強(qiáng)制減慢其速度，以保證行車安全。在這樣的系統(tǒng)中，離散變化的邏輯控制和實(shí)時(shí)連續(xù)行為結(jié)合在一起，共同構(gòu)成一種復(fù)雜的動(dòng)力系統(tǒng)，稱為混成系統(tǒng)［3］?；斐上到y(tǒng)通常采用微分方程來描述其連續(xù)動(dòng)態(tài)系統(tǒng)，采用狀態(tài)自動(dòng)機(jī)來描述系統(tǒng)模式之間的變遷?；斐上到y(tǒng)的安全性驗(yàn)證是指，從給定初始狀態(tài)區(qū)域出發(fā)的任意軌線不會(huì)到達(dá)指定的不安全區(qū)域。由于混成系統(tǒng)在安全攸關(guān)系統(tǒng)中廣泛存在，對(duì)混成系統(tǒng)安全性的研究已經(jīng)成為一個(gè)非常重要的課題［4］。

目前大部分關(guān)于混成系統(tǒng)的工作都集中在確定性模型上。然而，由于未知輸入、外界干擾、內(nèi)部錯(cuò)誤等因素的存在，嵌入式系統(tǒng)中不可避免存在不確定性因素。如何對(duì)相關(guān)不確定或隨機(jī)因素建模，是混成系統(tǒng)中的一個(gè)重要的研究方向。不確定性混成系統(tǒng)可以對(duì)離散變遷或連續(xù)演化進(jìn)行一些非確定性的選擇。一般而言，其局限在于在很多情況下結(jié)論過于粗糙，無法給出令人滿意的回答。相比之下，隨機(jī)混成系統(tǒng)［5］可以對(duì)安全性給出一些概率性的保證，且在實(shí)際當(dāng)中應(yīng)用范圍更廣。通過為初始狀態(tài)、離散變遷與連續(xù)演變添加隨機(jī)因素，近年來已經(jīng)提出一系列隨機(jī)或不確定性混成系統(tǒng)模型。文獻(xiàn)［6－8］考慮了具有隨機(jī)初始狀態(tài)的隨機(jī)混成系統(tǒng)，文獻(xiàn)［9－10］討論了具有隨機(jī)性離散變遷的混成系統(tǒng)。在文獻(xiàn)［11］中，使用隨機(jī)微分方程來刻畫系統(tǒng)的連續(xù)動(dòng)態(tài)。通過將幾種隨機(jī)因素結(jié)合起來，可以得到較為一般的隨機(jī)混成系統(tǒng)［12－13］。除此以外，還有各種其他形式的混成系統(tǒng)，例如文獻(xiàn)［14］在仿射混成自動(dòng)機(jī)中加入不確定性的向量場(chǎng)。文獻(xiàn)［15］在實(shí)時(shí)自動(dòng)機(jī)中引入離散概率分布，稱為概率實(shí)時(shí)自動(dòng)機(jī)。

對(duì)于隨機(jī)混成系統(tǒng)來說，安全性驗(yàn)證問題是指對(duì)計(jì)算系統(tǒng)狀態(tài)到達(dá)不安全區(qū)域的概率。文獻(xiàn)［16］考慮離散時(shí)間的隨機(jī)混成系統(tǒng)，首先將其近似為有限狀態(tài)馬爾可夫鏈，然后計(jì)算安全性概率。文獻(xiàn)［7］利用逐段決定馬爾可夫過程來研究具有隨機(jī)初始狀態(tài)的混成系統(tǒng)。文獻(xiàn)［17］使用 SSMT(Stochastic Satisfiability Modulo Theory)來計(jì)算隨機(jī)或不確定混成自動(dòng)機(jī)安全性概率的上下界。對(duì)于具有隨機(jī)初始條件的混成系統(tǒng)，ProbReach［6］將初始區(qū)域分割成足夠小的區(qū)間然后分段進(jìn)行驗(yàn)證，其不足在于計(jì)算負(fù)擔(dān)較重，而相比之下SReach［18］雖然利用統(tǒng)計(jì)檢驗(yàn)可以處理大規(guī)模問題，但只能得到統(tǒng)計(jì)意義下的結(jié)果。此外，以上兩種工具都只能驗(yàn)證在有界時(shí)間有限步內(nèi)的安全性。

障礙驗(yàn)證［19］是混成系統(tǒng)安全性驗(yàn)證中的一種重要的方法，其主要思想是尋找一個(gè)障礙函數(shù)，將可達(dá)集與不安全區(qū)域分開，從而確保系統(tǒng)的安全性，與其他方法相比，障礙驗(yàn)證具有計(jì)算簡(jiǎn)單，易于驗(yàn)證的特點(diǎn)。文獻(xiàn)［11］將其推廣到連續(xù)動(dòng)態(tài)由隨機(jī)微分方程刻畫的混成系統(tǒng)中。對(duì)于同時(shí)具有隨機(jī)連續(xù)動(dòng)態(tài)和離散動(dòng)態(tài)的隨機(jī)混成系統(tǒng)［12］，文獻(xiàn)［20］進(jìn)一步給出了障礙驗(yàn)證條件。對(duì)于具有隨機(jī)初始狀態(tài)的隨機(jī)連續(xù)系統(tǒng)，文獻(xiàn)［8］給出了一種基于初始集選擇的障礙驗(yàn)證方法。然而，對(duì)于同時(shí)在初始狀態(tài)、離散變遷與連續(xù)演變方面具有隨機(jī)性的一般隨機(jī)混成系統(tǒng)，目前尚未有文獻(xiàn)給出基于障礙驗(yàn)證的方法。為了對(duì)這類系統(tǒng)進(jìn)行研究，一種自然的想法是對(duì)單個(gè)隨機(jī)連續(xù)系統(tǒng)進(jìn)行研究，然后結(jié)合這些連續(xù)系統(tǒng)的性質(zhì)和離散跳轉(zhuǎn)的條件得出隨機(jī)混成系統(tǒng)的性質(zhì)。

受文獻(xiàn)［8］啟發(fā)，本文基于初始集選擇的方法提出一種隨機(jī)連續(xù)系統(tǒng)障礙驗(yàn)證方法，從而驗(yàn)證隨機(jī)連續(xù)系統(tǒng)的安全性。本文的主要工作為:給出隨機(jī)連續(xù)系統(tǒng)安全性驗(yàn)證的充分條件，并對(duì)半代數(shù)隨機(jī)連續(xù)系統(tǒng)給出隨機(jī)障礙函數(shù)求解的方法。

1 預(yù)備知識(shí)

本章主要介紹本領(lǐng)域的基礎(chǔ)知識(shí)與本文中需要用到的基本概念，包括連續(xù)動(dòng)力系統(tǒng)、隨機(jī)連續(xù)系統(tǒng)、障礙驗(yàn)證、平方和分解等;然后闡明本文希望解決的問題。在后文中，R代表實(shí)數(shù)集，x=(x1，x2，…，xn)T表示變量， x表示變量x對(duì)時(shí)間t的導(dǎo)數(shù)。除非特別說明，以下提到的X、I以及XU均默認(rèn)為有界閉集，f和g為連續(xù)函數(shù)。

1.1 連續(xù)系統(tǒng)

自治的連續(xù)系統(tǒng)通常由以下常微分方程組刻畫:

其中f:Rn→Rn稱為向量場(chǎng)，或微分規(guī)則。通常要求f滿足Lipschitz連續(xù)性，即存在常數(shù)K，使得:

連續(xù)系統(tǒng)的定義如下:

定義1 連續(xù)系統(tǒng)。連續(xù)系統(tǒng)S由多元組〈X，D，I，XU〉構(gòu)成，其中:

1)X Rn稱為位置不變集;

2)D定義了X上的微分規(guī)則 x=f(x)，其中f滿足X上的Lipschitz連續(xù)性;

3)I Rn表示初始狀態(tài)集合;

4)XURn為不安全區(qū)域。

對(duì)于連續(xù)系統(tǒng)S，從初始點(diǎn)x0∈I出發(fā)的軌線指的是方程限制在X中的解，所有從初始狀態(tài)集出發(fā)的軌線的并集稱為連續(xù)系統(tǒng)S的可達(dá)集。

在連續(xù)系統(tǒng)中，常常關(guān)心如下的安全性驗(yàn)證問題:對(duì)于給定的連續(xù)系統(tǒng)S，從初始區(qū)域I出發(fā)的軌線是否會(huì)進(jìn)入不安全區(qū)域XU。如果任意從初始區(qū)域出發(fā)的軌線不會(huì)進(jìn)入不安全區(qū)域，則稱該連續(xù)系統(tǒng)是安全的。

1.2 障礙驗(yàn)證

對(duì)于連續(xù)系統(tǒng)的安全性驗(yàn)證問題，一種簡(jiǎn)單的想法是尋找一條曲線，將系統(tǒng)的可達(dá)集與不安全區(qū)域分開，從而證明系統(tǒng)的軌線不可能進(jìn)入不安全區(qū)域，如圖1所示。為此，只需要尋找一個(gè)實(shí)函數(shù)B(x)，使其將可達(dá)集中的所有狀態(tài)映成負(fù)數(shù)，將不安全區(qū)域映成正數(shù)，從而證明可達(dá)集與不安全區(qū)域交集為空。在障礙驗(yàn)證中，障礙函數(shù)B(x)可利用系統(tǒng)的一些性質(zhì)求出，從而避免了直接計(jì)算系統(tǒng)可達(dá)集的困難。文獻(xiàn)［19］給出了連續(xù)系統(tǒng)障礙驗(yàn)證的一個(gè)充分條件，現(xiàn)陳述如下。

定理1 對(duì)于連續(xù)系統(tǒng)S=〈X，D，I，XU〉，如果存在連續(xù)可微函數(shù)B:Rn→R滿足:

一般來說，障礙函數(shù)的構(gòu)造與驗(yàn)證并非一件容易的事，但是當(dāng)連續(xù)系統(tǒng)為半代數(shù)系統(tǒng)，即微分規(guī)則D由多項(xiàng)式給出，并且位置不變集X，初始區(qū)域I與不安全區(qū)域XU由多項(xiàng)式的等式或不等式給出時(shí)，可以利用平方和規(guī)劃求解出多項(xiàng)式形式的障礙函數(shù)。

圖1 連續(xù)系統(tǒng)障礙驗(yàn)證Fig.1 Barrier certificates of continuous system

1.3 平方和

多項(xiàng)式的非負(fù)性判定是很多領(lǐng)域中的基本問題，例如動(dòng)力系統(tǒng)的穩(wěn)定性分析問題可歸結(jié)為驗(yàn)證某些條件非負(fù)。然而，多項(xiàng)式的非負(fù)性判定問題是NP-hard問題［21］，這促使人們使用一些易于判斷的充分條件。平方和(Sum of Squares，SOS)分解就是這樣一種可在多項(xiàng)式時(shí)間內(nèi)驗(yàn)證的條件。

首先介紹平方和多項(xiàng)式的基本概念。對(duì)于給定的多項(xiàng)式p(x)，如果存在一組多項(xiàng)式 p1(x)，p2(x)，…，pm(x)，滿足p(x)的平方和分解。顯然，如果p(x)為平方和多項(xiàng)式，那么p(x)≥0，反之則不一定成立，例如考慮Motzkin多項(xiàng)式［22］:

該多項(xiàng)式是非負(fù)多項(xiàng)式，但是無法表示成平方和的形式。事實(shí)上，只有在以下三種情況下非負(fù)多項(xiàng)式才一定能寫成平方和多項(xiàng)式的形式:一元多項(xiàng)式，2次多項(xiàng)式，或者次數(shù)為4的二元多項(xiàng)式，在其他情況下非負(fù)多項(xiàng)式不一定存在平方和分解。

對(duì)于給定的多項(xiàng)式是否可以存在平方和分解的問題，有以下結(jié)論［23］:設(shè) p(x) ∈ R［x］，次數(shù)為 deg p=2d，則判斷p(x)是否為平方和多項(xiàng)式等價(jià)于判斷p(x)能否寫成二次型ZTQZ的形式，其中Z是由次數(shù)不大于d的單項(xiàng)式構(gòu)成的向量，Q為半正定矩陣，稱為p(x)的Gram矩陣。半正定Gram矩陣的存在性可由半定規(guī)劃(Semidefinite Programming，SDP)的方法求解。

例1［24］試將多項(xiàng)式 p(x1，x2) =+－－－表示平方和多項(xiàng)式的形式。

于是p(x1，x2)=ZTQZ轉(zhuǎn)化為以下半定規(guī)劃可行性問題:q11=5，q22=2，2q12+q33= － 1，q13= － 1，q23= － 1，Q 0。利用SeDuMi［25］等軟件包求解上述半定規(guī)劃問題，得到Q的一個(gè)可行解:當(dāng)q12=－1且q33=1時(shí)Q為半正定矩陣，且Q=LTL，其中L。于是有以下平方和分解:

根據(jù)上面的討論可知，多項(xiàng)式的非負(fù)性判斷非常困難，而判定多項(xiàng)式是否為平方和多項(xiàng)式可以轉(zhuǎn)化為半定規(guī)劃問題，從而可以在多項(xiàng)式時(shí)間內(nèi)完成驗(yàn)證。

以下基于平方和松弛方法給出半代數(shù)連續(xù)系統(tǒng)障礙驗(yàn)證的充分條件。

定理2 設(shè)半代數(shù)連續(xù)系統(tǒng)S的位置不變集X、初始區(qū)域I與不安全區(qū)域XU分別由下式給出:

其中θX(x)、θI(x)、θU(x)表示多項(xiàng)式向量。如果存在多項(xiàng)式B(x)∈R［x］，常數(shù)ε＞0以及平方和多項(xiàng)式向量σX(x)、σI(x)、σU(x)，滿足:

其中，S表示平方和多項(xiàng)式全體構(gòu)成的集合，那么連續(xù)系統(tǒng)S是安全的，且B(x)為障礙函數(shù)。

上述定理的正確性依賴于以下事實(shí):如果 p(x)=σ1(x)+σ2(x)q(x)，σ1(x)，σ2(x)∈S，那么由q(x) ≥0可推出 p(x)≥ 0。類似地，如果 p(x) = ε+σ1(x) +σ2(x)q(x)，ε ＞ 0，σ1(x)，σ2(x) ∈ S，那么 q(x) ≥ 0 p(x) ＞0。

為了將定理2應(yīng)用到實(shí)際的半代數(shù)連續(xù)系統(tǒng)的安全性驗(yàn)證中，首先要確定所要尋找的障礙函數(shù)次數(shù)的上界(可以隨意決定，但是次數(shù)的提高會(huì)導(dǎo)致計(jì)算復(fù)雜度的提高)，假設(shè)要尋找次數(shù)為2d的障礙函數(shù)，則可以將障礙函數(shù)寫成B(x)=的形式，其中，bi(x)為次數(shù)不超過d

的單項(xiàng)式，ci∈R是未知量。接下來，人為取定乘子σX(x)、σI(x)與σU(x)的次數(shù)，并將乘子表示為單項(xiàng)式的線性組合，單項(xiàng)式的系數(shù)為未知數(shù)。將障礙函數(shù)B(x)以及乘子σX、σI、σU代入式(3)并將問題寫成二次型的形式，得到關(guān)于未知系數(shù)的線性方程組以及關(guān)于平方和多項(xiàng)式Gram矩陣正定的條件，于是可以通過半定規(guī)劃求解得到障礙函數(shù)B(x)。由于半定規(guī)劃問題具有多項(xiàng)式時(shí)間的解法，故上述的步驟可在多項(xiàng)式時(shí)間內(nèi)完成。事實(shí)上，可以通過專門的平方和規(guī)劃軟件，如SOSTOOLS［26－27］等，直接對(duì)式(3) 進(jìn)行求解。只需給定障礙函數(shù)與平方和多項(xiàng)式的次數(shù)上界，軟件會(huì)自動(dòng)將問題轉(zhuǎn)化為半定規(guī)劃的形式，并使用內(nèi)點(diǎn)法求解。一旦計(jì)算得到滿足式(3)的B(x)，則系統(tǒng)的安全性得證。如果沒有找到滿足條件的B(x)，則可以嘗試提高障礙函數(shù)次數(shù)的上界以擴(kuò)大搜索范圍。同時(shí)也應(yīng)考慮到，由于平方和多項(xiàng)式與非負(fù)多項(xiàng)式之間存在天然的不等價(jià)性，對(duì)于半代數(shù)連續(xù)系統(tǒng)，可能存在多項(xiàng)式障礙函數(shù)B(x)滿足定理1的條件卻不滿足定理2的條件的情況，此時(shí)，通過本文中描述的算法無法證明系統(tǒng)的安全性。

1.4 隨機(jī)連續(xù)系統(tǒng)

在混成系統(tǒng)安全性驗(yàn)證的研究中，對(duì)單個(gè)連續(xù)系統(tǒng)的研究是極為重要的。為了研究一般性隨機(jī)混成系統(tǒng)的安全性驗(yàn)證問題，必須考慮一個(gè)具有一定隨機(jī)性的連續(xù)系統(tǒng)。例如，可以考慮帶有不確定參數(shù)的連續(xù)系統(tǒng)，即系統(tǒng)的初始狀態(tài)由隨機(jī)變量定義，并且系統(tǒng)的連續(xù)動(dòng)態(tài)行為由確定性微分方程刻畫。對(duì)于這類隨機(jī)連續(xù)系統(tǒng)，已有相關(guān)文獻(xiàn)對(duì)其進(jìn)行研究。

本文中所考慮的隨機(jī)連續(xù)系統(tǒng)是一類具有隨機(jī)初始值，并且其連續(xù)動(dòng)態(tài)由隨機(jī)微分方程刻畫的連續(xù)系統(tǒng)。具體定義如下。

定義2 隨機(jī)連續(xù)系統(tǒng)。隨機(jī)連續(xù)系統(tǒng)S由多元組〈X，D，Ω，XU〉構(gòu)成，其中:

1)X Rn稱為位置不變集;

2)D=〈f，g〉定義了 X上的隨機(jī)微分方程 dx=f(x)dt+g(x)dw(t)，其中w(t)為n維布朗運(yùn)動(dòng);

3)Ω是系統(tǒng)初始值所服從的分布;

4)XU是給定的不安全區(qū)域。

通常要求定義中的f和g滿足某些條件以確保解的存在唯一性［28］，例如 Lipschitz條件:

以及線性增長(zhǎng)條件:

其中K1和K2為某常數(shù)。

設(shè)x0～Ω是給定的初始隨機(jī)變量，隨機(jī)連續(xù)系統(tǒng)S由x0出發(fā)的軌線指的是如下隨機(jī)微分方程限制在X上的解:

正式地說，設(shè)x(t)是X上滿足式(4)的隨機(jī)過程，停時(shí)τ為x(t)第一次離開X的內(nèi)部的時(shí)間，則隨機(jī)連續(xù)系統(tǒng)的軌線是滿足如下條件的隨機(jī)過程:

1.5 問題描述

本文考慮半代數(shù)隨機(jī)連續(xù)系統(tǒng)，其約束條件由多項(xiàng)式給出。

定義3 半代數(shù)隨機(jī)連續(xù)系統(tǒng)。半代數(shù)隨機(jī)連續(xù)系統(tǒng)是一個(gè)隨機(jī)連續(xù)系統(tǒng)〈X，D，Ω，XU〉，滿足:

1)隨機(jī)微分方程的漂移項(xiàng)f和擴(kuò)散項(xiàng)g都是多項(xiàng)式;

2)位置不變集X和不安全區(qū)域XU都是半代數(shù)集(由多項(xiàng)式等式與不等式給出)。

對(duì)于給定的半代數(shù)隨機(jī)連續(xù)系統(tǒng)S，本文所研究的安全性驗(yàn)證問題指的是從x0～Ω出發(fā)的軌線不進(jìn)入不安全區(qū)域XU的概率的下界。換句話說，即尋找p∈［0，1］使得:

2 隨機(jī)連續(xù)系統(tǒng)安全性分析

由于本文所討論的隨機(jī)連續(xù)系統(tǒng)同時(shí)在初始狀態(tài)和連續(xù)動(dòng)態(tài)方面具有隨機(jī)性，故對(duì)其進(jìn)行安全性驗(yàn)證時(shí)，需要同時(shí)考慮這兩方面的隨機(jī)性，較為復(fù)雜。為了解決概率安全性驗(yàn)證問題，本文首先將問題轉(zhuǎn)化為相應(yīng)的具有確定性初始狀態(tài)區(qū)域的隨機(jī)連續(xù)系統(tǒng)，其連續(xù)動(dòng)態(tài)由隨機(jī)微分方程刻畫。然后利用隨機(jī)障礙函數(shù)驗(yàn)證的方法求解該隨機(jī)連續(xù)系統(tǒng)并使用平方和松弛的方法計(jì)算出安全性概率下界。具體而言包含以下幾個(gè)步驟:

1)根據(jù)隨機(jī)連續(xù)系統(tǒng)S的初始分布Ω確定帶參數(shù)的初始狀態(tài)集I(α);

2)計(jì)算具有確定性初始狀態(tài)集I(α)的隨機(jī)連續(xù)系統(tǒng)的安全性概率的下界;

3)通過調(diào)整α的取值，可以找到安全性概率盡可能好的下界。

2.1 初始狀態(tài)集

為了克服隨機(jī)連續(xù)系統(tǒng)同時(shí)具有隨機(jī)性的初始狀態(tài)和隨機(jī)微分方程的難點(diǎn)，首先將系統(tǒng)轉(zhuǎn)化為具有確定性初始狀態(tài)集的隨機(jī)連續(xù)系統(tǒng)，即考慮軌線從初始狀態(tài)集I中的某個(gè)狀態(tài)出發(fā)的隨機(jī)連續(xù)系統(tǒng)，其中I Rn是有界閉集。與本文中定義的隨機(jī)連續(xù)系統(tǒng)S相比，該系統(tǒng)僅僅在刻畫其連續(xù)動(dòng)態(tài)時(shí)引入了隨機(jī)微分方程。以下將初始狀態(tài)集為I的隨機(jī)連續(xù)系統(tǒng)記為SI。對(duì)于給定的初始分布Ω，本文考慮如下形式的初始區(qū)域:

其中:α是待定的參數(shù);I(α)是有界閉集。

考慮到不同的概率分布具有不同的密度函數(shù)，一種自然的想法是根據(jù)其密度函數(shù)來選擇系統(tǒng)的初始區(qū)域。以下通過四種不同的初始概率分布為例來闡述其思想。

2.1.1 均勻分布

均勻分布，或稱為連續(xù)型均勻分布，是一種簡(jiǎn)單的連續(xù)型概率分布，其密度函數(shù)在支集上處處相等。均勻分布在概率統(tǒng)計(jì)中有許多應(yīng)用，例如根據(jù)最大熵原理，如果對(duì)某個(gè)隨機(jī)變量的性質(zhì)一無所知，則可假設(shè)其服從均勻分布。

設(shè)初始狀態(tài) x0= ［x1，x2，…，xn］T為隨機(jī)向量，其中 x1，x2，…，xn為獨(dú)立隨機(jī)變量。對(duì)于每個(gè) i∈ {1，2，…，n}，xi服從［li，ui］上的均勻分布，即xi～ U(li，ui)，則x0的密度函數(shù)為:

隨機(jī)向量x0在超立方體{x∈Rn|li≤xi≤ui，1≤i≤n}中等可能性地取值。因此，對(duì)于均勻分布，選擇初始狀態(tài)集:

其中ci=(ui+li)/2，ri=(ui－ri)/2，α為0到1之間的參數(shù)。顯然，初始狀態(tài)落在IU(α)中的概率為:

2.1.2 指數(shù)分布

在概率論中，指數(shù)分布常用來表示獨(dú)立事件發(fā)生的間隔，如旅客進(jìn)機(jī)場(chǎng)的時(shí)間間隔等。此外，許多電子設(shè)備的壽命也服從指數(shù)分布，有些系統(tǒng)的壽命也近似服從指數(shù)分布。指數(shù)分布的特點(diǎn)在于，隨著間隔時(shí)間變長(zhǎng)，事件發(fā)生的概率急劇下降。

設(shè) x1，x2，…，xn為獨(dú)立隨機(jī)變量，xi－ μi服從參數(shù)為 λi的指數(shù)分布，即xi－ μi～ E(λi)，其中i=1，2，…，n。則xi的密度函數(shù)可以表示為:

其中λi為大于0的常數(shù)。初始狀態(tài)x0的聯(lián)合密度函數(shù)為:

為了反映指數(shù)型分布的特性，考慮如下形式的初始區(qū)域:

其中β為給定的閾值。通過簡(jiǎn)單的計(jì)算可知，ρ(x)≥β等價(jià)于:

其中α≥0。

根據(jù)文獻(xiàn)［8］，初始狀態(tài)x0落在IE(α)中的概率可由式(11)計(jì)算:

2.1.3 正態(tài)分布

在概率論中，正態(tài)分布是一種極為重要的連續(xù)型概率分布，其密度函數(shù)曲線常被稱為鐘形曲線。正態(tài)分布具有許多良好的性質(zhì)，例如，兩個(gè)正態(tài)隨機(jī)變量的和差服從正態(tài)分布。根據(jù)中心極限定理，在某些條件下，大量相互獨(dú)立的隨機(jī)變量，其均值的分布以正態(tài)分布為極限。因此，在自然科學(xué)與社會(huì)科學(xué)的領(lǐng)域中，有時(shí)常用正態(tài)分布來近似某些未知的分布。另外，由于正態(tài)分布的熵在所有的已知均值及方差的分布中最大，這使得它作為一種均值以及方差已知的分布的自然選擇。

以下考慮x0的每個(gè)分量是獨(dú)立的正態(tài)隨機(jī)變量的情況。設(shè)x1，x2，…，xn為獨(dú)立的隨機(jī)變量，分別服從均值為μi、方差為σ2i的正態(tài)分布，其密度函數(shù)為:

其中1≤i≤n。于是x0的聯(lián)合密度函數(shù)可表示為:

考慮ρ(x)的β上水平集:

上式等價(jià)于:

即:

其中α為大于0的參數(shù)。

2.1.4 聯(lián)合正態(tài)分布

聯(lián)合正態(tài)分布可以看成是正態(tài)分布的一般情況。若隨機(jī)向量珘g是k維標(biāo)準(zhǔn)高斯隨機(jī)向量，即其元素為獨(dú)立的隨機(jī)變量，服從均值為0、方差為1的標(biāo)準(zhǔn)正態(tài)分布，如果存在矩陣A ∈Rn×k及向量 μ ∈Rn，使n元隨機(jī)向量x滿足x=A珘g+μ，則稱x為高斯向量，或稱其滿足聯(lián)合正態(tài)分布。向量μ為x的均值，協(xié)方差矩陣為Σ =AAT。

現(xiàn)假設(shè)初始狀態(tài)x0服從n元聯(lián)合正態(tài)分布，其均值為μ∈Rn，協(xié)方差矩陣為Σ∈Rn×n。根據(jù)定義，如果Σ為正定矩陣，則x0的密度函數(shù)為:

其中，|Σ|表示Σ的行列式。

考慮ρ(x)的β上水平集:

等價(jià)于:

不難看出，2.1.3節(jié)定義的IN(α)是IJN(α)的特殊情況。事實(shí)上，如果獨(dú)立的隨機(jī)變量x1，x2，…，xn分別服從均值為μ1，μ2，…，μn和方差為，，…，的正態(tài)分布，則 x0=［x1，x2，…，xn］T服從均值為μ = ［μ1，μ2，…，μn］T及協(xié)方差矩陣為 Σ =diag(，σ22，… ，σ2n)的聯(lián)合正態(tài)分布。此時(shí)IN(α)與IJN(α)定義一致。

根據(jù)文獻(xiàn)［29］，初始狀態(tài)x0落在IJN(α)的概率為:

其中Γ為伽瑪函數(shù)。注意到式(14)中的概率與μ和Σ無關(guān)，故對(duì)于x0的分量為獨(dú)立的正態(tài)隨機(jī)變量的情況，同樣有:

2.2 安全性驗(yàn)證

2.2.1 隨機(jī)障礙驗(yàn)證

對(duì)于隨機(jī)連續(xù)系統(tǒng)，也可以借鑒障礙驗(yàn)證的思想，求出滿足特定條件的隨機(jī)障礙函數(shù)以驗(yàn)證系統(tǒng)的安全性。然而，由于此時(shí)系統(tǒng)的軌線是一個(gè)隨機(jī)過程，因此使用障礙函數(shù)分離初始區(qū)域和不安全區(qū)域的思想在此不適用。類比于確定性連續(xù)系統(tǒng)S障礙驗(yàn)證的條件中要求B(x(t))隨著時(shí)間t增長(zhǎng)而遞減的性質(zhì)，在隨機(jī)連續(xù)系統(tǒng)障礙驗(yàn)證中要求B(x(t))是上鞅(相關(guān)內(nèi)容可參考文獻(xiàn)［28］)，即對(duì)于任意s≤t，B(x(t))關(guān)于s之前的所有時(shí)刻的期望小于B(x(s))。直觀地說，B(x(t))在期望的意義下關(guān)于時(shí)間t遞減或保持不變。

以下陳述對(duì)于具有確定性初始狀態(tài)的隨機(jī)連續(xù)系統(tǒng)障礙驗(yàn)證的結(jié)論。

定理 3［11］對(duì)于隨機(jī)連續(xù)系統(tǒng) SI= 〈X，D，I，XU〉，其中D=〈f，g〉表示隨機(jī)微分方程 dx=f(x)dt+g(x)dw(t)，且軌線x(t)從x0∈I出發(fā)。如果存在二次連續(xù)可微函數(shù)B:Rn→R滿足:

則系統(tǒng)不安全的概率滿足:

由式(15)可知，任意從初始集I出發(fā)的軌線滿足0≤B(x0)≤γ，且B(x)不安全區(qū)域U上取值至少為1。此外，隨著時(shí)間演變，隨機(jī)過程B(x(t))始終非負(fù)。另外事實(shí)上由式(16)可以推出B(x(t))是上鞅，因此從直觀上說，如果γ越小，那么隨著時(shí)間t的演變，隨機(jī)過程B(x(t))的值從γ增長(zhǎng)到1的概率越小，從而軌線x(t)進(jìn)入不安全區(qū)域XU的概率越小。

利用定理3，可以得到關(guān)于隨機(jī)連續(xù)系統(tǒng)S安全性驗(yàn)證的結(jié)論。

定理4 設(shè)S=〈X，D，Ω，XU〉是給定的隨機(jī)連續(xù)系統(tǒng)，其中 D =〈f，g〉表示隨機(jī)微分方程:dx = f(x)dt+g(x)dw(t)，區(qū)域I X是根據(jù)初始分布Ω所選擇的初始狀態(tài)區(qū)域。如果存在二次連續(xù)可微函數(shù)B:Rn→R滿足式(15)與(16)，則系統(tǒng)的安全性概率可由如式(19)計(jì)算:

證明 由于I是系統(tǒng)的初始狀態(tài)區(qū)域，由定理3知，任意從I出發(fā)的軌線滿足:

P{x(t) XU， t≥0|x(0)=x0}≥1－γ從而有:

由I X可知:

因此安全性概率下界得證。

2.1 節(jié)中分別對(duì)隨機(jī)連續(xù)系統(tǒng)S的幾種常見的初始分布進(jìn)行討論。如果初始狀態(tài)服從獨(dú)立的均勻分布，則選取初始狀態(tài)集I(α)為式(8)。如果初始狀態(tài)服從其他三種分布，則考慮其密度函數(shù)的上水平集，分別選取初始狀態(tài)集為式(10)、(12)、(13)。初始狀態(tài)落在選定的初始狀態(tài)集的概率可由相應(yīng)公式計(jì)算。安全性概率的下界為:

其中γα表示初始集為I(α)時(shí)由定理3得出的不安全的上界。若記:

則根據(jù)式(19)，系統(tǒng)安全性概率的下界計(jì)算式如下:

2.2.2 求解隨機(jī)障礙函數(shù)

由于式(16)和(17)是障礙驗(yàn)證的充分條件，這意味著可以通過尋找不同的B(x)與γ使得Psafe的下界盡可能大，因此問題可以表述成如下形式:

其中:I(α) Rn是根據(jù)分布Ω選定的初始集;A R是參數(shù)α的取值范圍，例如對(duì)于均勻分布有A=(0，1］，而對(duì)于其他三種分布有A=(0，+∞)。

由于式(21)中的隨機(jī)障礙函數(shù)B(x)∈C2(X)，故其本質(zhì)上是個(gè)無窮維的優(yōu)化問題。為了便于求解，將問題轉(zhuǎn)化為給定次數(shù)上界的平方和規(guī)劃問題，從而可以使用SOSTOOLS等求解平方和規(guī)劃的軟件對(duì)問題進(jìn)行求解。

設(shè)半代數(shù)集X、XU和I(α)分別定義如下:

則可以對(duì)式(21)進(jìn)行平方和松弛:

其中S表示平方和多項(xiàng)式全體構(gòu)成的集合。不難看出，如果p和分別是問題(21)與(22)的解，那么安全性概率滿足:

根據(jù)式(22)及(23)，設(shè)計(jì)求解隨機(jī)連續(xù)系統(tǒng)S安全性概率Psafe下界的算法1如下。

算法1 計(jì)算安全性概率的下界。

輸入 半代數(shù)隨機(jī)連續(xù)系統(tǒng)S，多項(xiàng)式隨機(jī)障礙函數(shù)B(x)的次數(shù) k，平方和多項(xiàng)式 σX(x)，σU(x)，σI(x)(x) 的次數(shù) d;

輸出 安全性概率的下界^p。

1) 根據(jù)初始分布Ω給出初始狀態(tài)集的模板I(α)并確定α的取值范圍A。

2)根據(jù)d與k分別設(shè)定σX(x)，σU(x)，σI(x)，珘σX(x)與B(x)的模板。

3) 將A離散化并在其中取若干α值作循環(huán):

b) 利用SOSTOOLS求解以下問題:

c) 對(duì)第b)步中得到的結(jié)果進(jìn)行平方和檢驗(yàn):如果σX(x)，σU(x)，σI(x)，珘σX(x)可以進(jìn)行平方和分解，則記^p←(1－γ)Pα。

4) 將第3)步中選取的最大的^p值作為安全性概率的下界。

需要特別說明的是，算法中的第c)步的平方和檢驗(yàn)是為了避免在b)步中由于數(shù)值誤差可能得到不正確的結(jié)果。這種情況經(jīng)常在deg(B(x))較高、計(jì)算誤差較大的情況下出現(xiàn)。對(duì)于次數(shù)較低的情況，省略平方和檢驗(yàn)的步驟對(duì)結(jié)果的正確性影響不大。

2.2.3 算法復(fù)雜度分析

算法1可以分成兩個(gè)部分:確定初始集I(α)并計(jì)算初始狀態(tài)的概率Pα，以及基于平方和松弛的方法求解初始集為I(α)的隨機(jī)連續(xù)系統(tǒng)。

關(guān)于前者，如果選定的初始集I(α) X，那么可以直接利用2.1 節(jié)中給出的式(9)、(11)、(14)、(15) 計(jì)算 Pα，否則需要對(duì)Pα另行計(jì)算，例如可以通過蒙特卡羅法求概率的近似值，復(fù)雜度為O(1)。

后者的關(guān)鍵在于求解一個(gè)平方和問題。通過將問題中的多項(xiàng)式用二次型的形式表示，判斷多項(xiàng)式是否為平方和多項(xiàng)式可以轉(zhuǎn)化為判斷對(duì)應(yīng)的Gram矩陣是否半正定的問題，據(jù)此可以將第b)步中的平方和問題轉(zhuǎn)化為半定規(guī)劃可行性問題，運(yùn)用內(nèi)點(diǎn)法可以對(duì)其進(jìn)行求解，復(fù)雜度為O(+)，其中N代表半定規(guī)劃問題中優(yōu)化變vars量個(gè)數(shù)，Ncons代表限制條件的個(gè)數(shù)。為簡(jiǎn)便起見，以下假設(shè)θX、θU、θI的維數(shù)至多為K。不難看出，問題的變量個(gè)數(shù)為:

約束條件個(gè)數(shù)為:

其中:

因此有 Nvars=O(nk+Knd)以及 Ncons=O(n「d'/2?+Knd/2)。如果位置不變集X和不安全區(qū)域XU由超立方體或者橢球給出，則對(duì)于本文中給出的四種初始狀態(tài)集均有K≤n，此時(shí)有 Nvars～ nk，Ncons～ n「d'/2?。設(shè)s=max{k，deg f，2deg g －1}，此時(shí)第b)步平方和規(guī)劃的復(fù)雜度可以表示為O(n4s)。類似可得，第c)步的平方和檢驗(yàn)復(fù)雜度為O(n3.5d+1)。如果算法在第3)步中取了m個(gè)α的值，那么整個(gè)算法的復(fù)雜度可表示為O(mn4s)。

3 數(shù)值實(shí)驗(yàn)與分析

本章中選取多個(gè)隨機(jī)連續(xù)系統(tǒng)的例子求解其安全性下界。實(shí)驗(yàn)中使用的計(jì)算機(jī)的處理器為2.5 GHz Intel Core i7，內(nèi)存為16 GB，開發(fā)環(huán)境為 Matlab R2015b，SOSTOOLS的版本為 3.01。

例2［11］考慮以下隨機(jī)連續(xù)系統(tǒng):

其中w(t)是1維Wiener過程。設(shè)不變式為:

初始狀態(tài)為:

不安全區(qū)域?yàn)?

取系統(tǒng)的初始集為:

設(shè)定隨機(jī)障礙函數(shù)的次數(shù)，使用算法1進(jìn)行安全性驗(yàn)證，結(jié)果如表1所示，“*”號(hào)表示SOSTOOLS得到的解未通過平方和檢驗(yàn)，這代表結(jié)果可能存在一定的誤差。

表1 例2的安全性驗(yàn)證結(jié)果Tab.1 Safety verification results of Eg.2

在表1中，當(dāng)系統(tǒng)式(24)的擴(kuò)散項(xiàng)越小時(shí)，安全性概率的下界也就越高。而對(duì)于固定的σ來說，提高B(x)的次數(shù)有助于尋找到更高的安全性下界。與此同時(shí)，deg B的提高會(huì)導(dǎo)致較高的計(jì)算復(fù)雜度與較大的數(shù)值誤差，因此在進(jìn)行數(shù)值實(shí)驗(yàn)時(shí)并未使用較高的次數(shù)。

例3［30］考慮隨機(jī)連續(xù)系統(tǒng):

不變式為 X=［－1，1］2，不安全區(qū)域?yàn)?

初始分布分別為:

表2給出了對(duì)于四種分布的安全性驗(yàn)證結(jié)果。

表2 例3的安全性驗(yàn)證結(jié)果Tab.2 Safety verification results of Eg.3

為了進(jìn)一步測(cè)試算法的性能，選擇多個(gè)隨機(jī)連續(xù)系統(tǒng)的例子，算法的結(jié)果與耗時(shí)如表3所示。

表3 隨機(jī)連續(xù)系統(tǒng)安全性驗(yàn)證結(jié)果Tab.3 Safety verification results of stochastic continuous system

表3中:n代表變量的個(gè)數(shù);m代表實(shí)驗(yàn)中取的α的個(gè)數(shù);^p代表隨機(jī)連續(xù)系統(tǒng)安全性的下界;T代表實(shí)驗(yàn)的耗時(shí)。表3中前四行的結(jié)果來自同一個(gè)連續(xù)系統(tǒng)，而后四行的結(jié)果來自不同的連續(xù)系統(tǒng)，從表3中可以看出，初始分布對(duì)算法的性能影響不大，而隨機(jī)障礙函數(shù)次數(shù)的提高或是變量的增加都會(huì)顯著增加耗時(shí)。

4 結(jié)語(yǔ)

在以往混成系統(tǒng)或連續(xù)系統(tǒng)障礙驗(yàn)證的工作中，只有對(duì)具有隨機(jī)初始狀態(tài)或隨機(jī)微分方程的系統(tǒng)的研究，然而對(duì)同時(shí)具有多種隨機(jī)因素的一般性系統(tǒng)進(jìn)行驗(yàn)證也是有必要的。本文針對(duì)具有隨機(jī)初始狀態(tài)以及隨機(jī)微分方程的隨機(jī)連續(xù)系統(tǒng)進(jìn)行安全性驗(yàn)證，基于初始集選擇的方法將問題轉(zhuǎn)化為具有確定性初始狀態(tài)的隨機(jī)連續(xù)系統(tǒng)，然后基于障礙驗(yàn)證以及平方和松弛的方法得出系統(tǒng)安全性的理論下界，基于此可以對(duì)連續(xù)系統(tǒng)安全性進(jìn)行分析。下一步，我們將研究如何利用本文所提方法對(duì)隨機(jī)混成系統(tǒng)進(jìn)行安全性驗(yàn)證。