本發(fā)明涉及含區(qū)間參數(shù)結(jié)構(gòu)響應(yīng)的不確定性傳播領(lǐng)域，特別涉及一種基于超體積迭代策略的含區(qū)間參數(shù)結(jié)構(gòu)響應(yīng)區(qū)間的分析方法，該方法是基于超體積迭代策略全局優(yōu)化算法而提出的。

背景技術(shù)：

現(xiàn)代化工程結(jié)構(gòu)復(fù)雜而精細(xì)，但是由于制造工藝的限制、認(rèn)知的匱乏、測(cè)量的誤差以及環(huán)境因素的變化，不確定性問題廣泛存在于工程問題中，因此在分析和設(shè)計(jì)時(shí)需要更加注重先進(jìn)數(shù)值算法和設(shè)計(jì)技術(shù)的研究。對(duì)復(fù)雜系統(tǒng)而言，不確定因素會(huì)因?yàn)橄到y(tǒng)的非線性而使響應(yīng)的區(qū)間上/下界較難確定，現(xiàn)在解決該問題的方法主要有區(qū)間參數(shù)攝動(dòng)法、頂點(diǎn)組合法等。前者往往基于一階Taylor級(jí)數(shù)展開理論，將問題在很小的范圍內(nèi)線性化以此達(dá)到簡(jiǎn)化問題的目的，但是其不足也顯而易見，即只能針對(duì)非線性不太強(qiáng)的問題在較小不確定度的情況，對(duì)于機(jī)翼等中復(fù)合材料具有較大分散性的區(qū)間響應(yīng)分析問題就會(huì)有較大誤差。頂點(diǎn)組合法針對(duì)單調(diào)問題可以精確地得到響應(yīng)區(qū)間上/下界，如果不是單調(diào)的問題，則該方法將不再適用。以上都是基于區(qū)間數(shù)學(xué)的理論方法，蒙特卡洛法也可以解決該問題并得到較為精確的區(qū)間上/下界，甚至可以可作為對(duì)比其他方法精確度的一個(gè)參考基準(zhǔn)，該方法結(jié)果越精確意味著隨機(jī)抽樣次數(shù)越多，對(duì)于復(fù)雜結(jié)構(gòu)將帶來較大的計(jì)算負(fù)擔(dān)。因此針對(duì)非線性較強(qiáng)且具有較大不確定度的問題，研究一種在保證一定計(jì)算效率的條件下可以精確地確定響應(yīng)量區(qū)間上下界的方法，十分有意義。

本發(fā)明基于優(yōu)化方法的思路將響應(yīng)量作為優(yōu)化的目標(biāo)函數(shù)，通過求解其最大/小值即可得到區(qū)間的上/下界。這樣將一個(gè)區(qū)間響應(yīng)的不確定性傳播轉(zhuǎn)化為一個(gè)優(yōu)化問題，不僅使問題簡(jiǎn)化，無需考慮區(qū)間數(shù)學(xué)的理論方法，同時(shí)也克服上述諸多問題。在諸多優(yōu)化方法中，牛頓法由于其較快的收斂速度廣泛應(yīng)用于各優(yōu)化問題中，但是牛頓法針對(duì)的是局部?jī)?yōu)化問題，當(dāng)針對(duì)全局優(yōu)化問題時(shí)，牛頓法以及其他基于梯度信息的算法將陷入局部最優(yōu)解區(qū)域無法跳出，從而無法得到全局最優(yōu)解。啟發(fā)式算法是現(xiàn)在發(fā)展較好的一類全局優(yōu)化算法，包括模擬退火法、遺傳算法、粒子群算法等，該類方法通過模擬或揭示某些自然現(xiàn)象或過程而得到發(fā)展，其思想和內(nèi)容涉及數(shù)學(xué)、物理學(xué)、生物進(jìn)化、人工智能、統(tǒng)計(jì)力學(xué)等方面，為解決復(fù)雜問題提供新的思路和手段，但是該類方法屬于隨機(jī)算法，沒有明確的收斂準(zhǔn)則。本發(fā)明基于超體積迭代策略全局優(yōu)化算法，該方法以目標(biāo)函數(shù)在可行域的超體積作為迭代的收斂準(zhǔn)則，可以快速準(zhǔn)確地收斂到全局最優(yōu)解。

技術(shù)實(shí)現(xiàn)要素：

本發(fā)明要解決的技術(shù)問題是：克服現(xiàn)有方法和技術(shù)的不足，提供一種基于超體積迭代策略的含區(qū)間參數(shù)結(jié)構(gòu)響應(yīng)區(qū)間的分析方法，本發(fā)明充分考慮了不確定性因素對(duì)結(jié)構(gòu)響應(yīng)帶來的影響，建立了一種基于優(yōu)化思路的求解不確定性傳播的響應(yīng)區(qū)間上/下界的方法，相較于傳統(tǒng)的方法拓展了不確定性傳播分析的方法，并可以快速準(zhǔn)確地確定非線性較強(qiáng)、具有較大不確定度系統(tǒng)中結(jié)構(gòu)響應(yīng)區(qū)間的上/下界。

本發(fā)明采用的技術(shù)方案為：一種基于超體積迭代策略的含區(qū)間參數(shù)結(jié)構(gòu)響應(yīng)區(qū)間的分析方法，具體流程見圖1，其實(shí)現(xiàn)步驟如下：

步驟(1)基于MSC.Patran構(gòu)建結(jié)構(gòu)的有限元模型，包括建立幾何模型，賦予材料屬性，施加荷載和邊界條件。

步驟(2)確定由步驟(1)構(gòu)建的模型中的不確定參數(shù)x及其不確定區(qū)間。由于制作工藝、測(cè)量誤差等因素導(dǎo)致結(jié)構(gòu)中的材料和荷載具有一定的不確定性，結(jié)構(gòu)中材料的彈性模量作用在結(jié)構(gòu)上的荷載這里E和P分別表示彈性模量和荷載的下界，和分別表示彈性模量和荷載的上界，然后在MSC.Patran中修改相應(yīng)不確定參數(shù)的數(shù)值，并且該修改過程需要用到MSC.Patran的宏錄制功能最后形成記錄該修改過程的ses文件。其中區(qū)間不確定性參數(shù)向量x可以表示為：

其中，x^c＝(E^c,P^c,…)為區(qū)間中點(diǎn)，x^r＝(E^r,P^r,…)為區(qū)間半徑，e∈Ξⁿ，Ξⁿ定義為所有元素包含在[-1,1]內(nèi)的n維向量集合，符號(hào)“×”定義為兩個(gè)向量各對(duì)應(yīng)元素乘積仍為n維向量。

步驟(3)修改結(jié)構(gòu)有限元模型中對(duì)應(yīng)的不確定參數(shù)。在步驟(2)中形成的ses文件中找到修改模型中不確定參數(shù)的關(guān)鍵性命令，利用MATLAB的讀寫功能修改ses文件中修改模型的命令，使得可以在MATLAB中賦值來修改有限元模型中對(duì)應(yīng)的不確定參數(shù)x(如彈性模量、荷載)，完成有限元模型可被MATLAB調(diào)用并被修改的功能。

需要說明的是：第一次執(zhí)行該步驟時(shí)，首先取不確定參數(shù)區(qū)間上/下界及區(qū)間中點(diǎn)修改結(jié)構(gòu)有限元模型，第二次執(zhí)行及之后將按照步驟(9)中增加的積分點(diǎn)來修改不確定參數(shù)。

步驟(4)生成新的有限元模型。通過步驟(3)生成新的結(jié)構(gòu)有限元模型(db文件)，新的有限元模型除了上述修改的部分改變，其它都保持原來屬性。

步驟(5)利用MATLAB調(diào)用有限元求解器MSC.Nasran讀取步驟(4)中db文件，計(jì)算步驟(3)中結(jié)構(gòu)的響應(yīng)，形成xdb文件。

步驟(6)利用MATLAB調(diào)用MSC.Patran讀取步驟(5)中的xdb文件，對(duì)計(jì)算結(jié)果進(jìn)行后處理生成結(jié)果報(bào)告(rpt文件)，以MATLAB讀取相關(guān)的響應(yīng)量。

步驟(7)刪除計(jì)算的過程文件，因?yàn)榻Y(jié)算過程中產(chǎn)生的過程文件對(duì)下一次計(jì)算沒有意義，甚至可能會(huì)和下次循環(huán)中的文件難以區(qū)分，而影響下次循環(huán)中對(duì)相關(guān)文件的調(diào)用、讀取或修改，同時(shí)為了節(jié)約計(jì)算機(jī)儲(chǔ)存空間同時(shí)不影響后續(xù)計(jì)算，刪除與后續(xù)計(jì)算無關(guān)的過程文件，如步驟(5)中的xdb文件和步驟(6)中的rpt文件。

步驟(8)以結(jié)構(gòu)的響應(yīng)為被積函數(shù)，在不確定參數(shù)x的區(qū)間范圍內(nèi)進(jìn)行自適應(yīng)辛普森積分，所謂“自適應(yīng)”即積分點(diǎn)的分布和被積函數(shù)在可行域內(nèi)的變化趨勢(shì)相關(guān)，如果在某區(qū)間內(nèi)被積函數(shù)變化平緩，則積分點(diǎn)分布較為稀疏；如果某段區(qū)間被積函數(shù)變化劇烈，則該區(qū)間內(nèi)積分點(diǎn)則較密，如圖2所示，以一維問題來作說明，首先在積分域[a,b]的兩端和中點(diǎn)c布置3個(gè)積分點(diǎn)，對(duì)函數(shù)f(·)進(jìn)行辛普森積分，如下式所示：

為了驗(yàn)證3個(gè)積分點(diǎn)計(jì)算結(jié)果是否滿足精度，分別在[a,c]和[c,b]域內(nèi)中點(diǎn)增加積分點(diǎn)d和e，分別在這兩個(gè)區(qū)間采用上式積分并求和為Q_s，如果Q_s和Q的相對(duì)誤差δ滿足容許誤差ε，即3個(gè)積分點(diǎn)計(jì)算結(jié)果滿足計(jì)算精度；如果Q_s和Q的相對(duì)誤差δ不滿足容許誤差ε，則在在[a,c]和[c,b]域內(nèi)分別按照上述方法驗(yàn)證各自區(qū)間內(nèi)是否滿足計(jì)算精度，如此循環(huán)，直至滿足精度為止。這樣就可以保證被積函數(shù)的最大/小值兩側(cè)的積分點(diǎn)之間的函數(shù)為凸函數(shù)，以便進(jìn)一步進(jìn)行牛頓法迭代得到最優(yōu)解，如圖3所示。

步驟(9)判斷步驟(8)中積分的結(jié)果是否滿足精度，判斷是否滿足精度的方法是判斷增加積分點(diǎn)前后計(jì)算得到的積分值的相對(duì)誤差δ是否滿足容許誤差ε要求。如果積分結(jié)果精度滿足容許誤差ε進(jìn)入下一步；如果不滿足計(jì)算精度則返回到步驟(3)進(jìn)行新一輪的循環(huán)。其中，容許誤差ε一般設(shè)置為10^-3，而相對(duì)誤差δ是迭代過程中前后兩次計(jì)算的積分值只差與前一次積分值的比值，即為：

其中，Q和Q_s分別是迭代前后的積分值。

步驟(10)選取現(xiàn)有的所有積分點(diǎn)集合中所對(duì)應(yīng)目標(biāo)函數(shù)的最大/小值及對(duì)應(yīng)的積分點(diǎn)，如圖3所示，圖3(a)為全局的積分點(diǎn)分布情況，選取這些積分點(diǎn)中目標(biāo)函數(shù)的最大/小值的積分點(diǎn)，如圖3(b)所示。

步驟(11)以(11)選取的積分點(diǎn)作為起點(diǎn)，由于該點(diǎn)到最優(yōu)解處之間為凸函數(shù)甚至是單調(diào)函數(shù)，故可以采用采用收斂速度較快的牛頓法進(jìn)行局部尋優(yōu)，得到的最優(yōu)解即為結(jié)構(gòu)響應(yīng)的上/下界。

步驟(12)結(jié)束。

本發(fā)明與現(xiàn)有技術(shù)相比的有點(diǎn)在于：

本發(fā)明以優(yōu)化思路為切入點(diǎn)，提出了基于超體積迭代策略全局優(yōu)化算法的結(jié)構(gòu)區(qū)間響應(yīng)的不確定性傳播分析方法，拓展了解決結(jié)構(gòu)區(qū)間響應(yīng)不確定性傳播問題的思路?；谠撍悸房梢院芎玫亟鉀Q非線性較強(qiáng)系統(tǒng)在不確定度較大時(shí)的區(qū)間響應(yīng)傳播問題，可以較為精確地確定區(qū)間響應(yīng)的上下界，相較于區(qū)間參數(shù)攝動(dòng)法精度更高，相較于蒙特卡洛法計(jì)算效率更高，相較于遺傳算法等優(yōu)化算法具有明確的收斂準(zhǔn)則。

附圖說明

圖1是本發(fā)明一種基于超體積迭代策略的含區(qū)間參數(shù)結(jié)構(gòu)響應(yīng)區(qū)間的分析方法流程圖；

圖2是一維自適應(yīng)辛普森積分示意圖；

圖3是一維問題自適應(yīng)積分點(diǎn)配置結(jié)果，其中，圖3(a)為全局的積分點(diǎn)分布情況，圖3(b)為選取這些積分點(diǎn)中目標(biāo)函數(shù)的最大/小值的積分點(diǎn)；

圖4是桁架結(jié)構(gòu)及相關(guān)物理參數(shù)。

具體實(shí)施方式

下面結(jié)合附圖以及具體實(shí)施方式進(jìn)一步說明本發(fā)明。

如圖3所示，本發(fā)明提出了一種基于超體積迭代策略的含區(qū)間參數(shù)結(jié)構(gòu)響應(yīng)區(qū)間的分析方法，包括以下步驟：

步驟(1)基于MSC.Patran構(gòu)建結(jié)構(gòu)的有限元模型，包括建立幾何模型，賦予材料屬性，施加荷載和邊界條件。

步驟(2)確定由步驟(1)構(gòu)建的模型中的不確定參數(shù)x及其不確定區(qū)間。由于制作工藝、測(cè)量誤差等因素導(dǎo)致結(jié)構(gòu)中的材料和荷載具有一定的不確定性，結(jié)構(gòu)中材料的彈性模量作用在結(jié)構(gòu)上的荷載這里E和P分別表示彈性模量和荷載的下界，和分別表示彈性模量和荷載的上界，然后在MSC.Patran中修改相應(yīng)不確定參數(shù)的數(shù)值，并且該修改過程需要用到MSC.Patran的宏錄制功能最后形成記錄該修改過程的ses文件。其中區(qū)間不確定性參數(shù)向量x可以表示為：

其中，x^c＝(E^c,P^c,…)為區(qū)間中點(diǎn)，x^r＝(E^r,P^r,…)為區(qū)間半徑，e∈Ξⁿ，Ξⁿ定義為所有元素包含在[-1,1]內(nèi)的n維向量集合，符號(hào)“×”定義為兩個(gè)向量各對(duì)應(yīng)元素乘積仍為n維向量。

步驟(3)修改結(jié)構(gòu)有限元模型中對(duì)應(yīng)的不確定參數(shù)。在步驟(2)中形成的ses文件中找到修改模型中不確定參數(shù)的關(guān)鍵性命令，利用MATLAB的讀寫功能修改ses文件中修改模型的命令，使得可以在MATLAB中賦值來修改有限元模型中對(duì)應(yīng)的不確定參數(shù)x(如彈性模量、荷載)，完成有限元模型可被MATLAB調(diào)用并被修改的功能。

需要說明的是：第一次執(zhí)行該步驟時(shí)，首先取不確定參數(shù)區(qū)間上/下界及區(qū)間中點(diǎn)修改結(jié)構(gòu)有限元模型，第二次執(zhí)行及之后將按照步驟(9)中增加的積分點(diǎn)來修改不確定參數(shù)。

步驟(4)生成新的有限元模型。通過步驟(3)生成新的結(jié)構(gòu)有限元模型(db文件)，新的有限元模型除了上述修改的部分改變，其它都保持原來屬性。

步驟(5)利用MATLAB調(diào)用有限元求解器MSC.Nasran讀取步驟(4)中db文件，計(jì)算步驟(3)中結(jié)構(gòu)的響應(yīng)，形成xdb文件。

步驟(6)利用MATLAB調(diào)用MSC.Patran讀取步驟(5)中的xdb文件，對(duì)計(jì)算結(jié)果進(jìn)行后處理生成結(jié)果報(bào)告(rpt文件)，以MATLAB讀取相關(guān)的響應(yīng)量。

步驟(7)刪除計(jì)算的過程文件，因?yàn)榻Y(jié)算過程中產(chǎn)生的過程文件對(duì)下一次計(jì)算沒有意義，甚至可能會(huì)和下次循環(huán)中的文件難以區(qū)分，而影響下次循環(huán)中對(duì)相關(guān)文件的調(diào)用、讀取或修改，同時(shí)為了節(jié)約計(jì)算機(jī)儲(chǔ)存空間同時(shí)不影響后續(xù)計(jì)算，刪除與后續(xù)計(jì)算無關(guān)的過程文件，如步驟(5)中的xdb文件和步驟(6)中的rpt文件。

步驟(8)以結(jié)構(gòu)的響應(yīng)為被積函數(shù)，在不確定參數(shù)x的區(qū)間范圍內(nèi)進(jìn)行自適應(yīng)辛普森積分，所謂“自適應(yīng)”即積分點(diǎn)的分布和被積函數(shù)在可行域內(nèi)的變化趨勢(shì)相關(guān)，如果在某區(qū)間內(nèi)被積函數(shù)變化平緩，則積分點(diǎn)分布較為稀疏；如果某段區(qū)間被積函數(shù)變化劇烈，則該區(qū)間內(nèi)積分點(diǎn)則較密，如圖2所示，以一維問題來作說明，首先在積分域[a,b]的兩端和中點(diǎn)c布置3個(gè)積分點(diǎn)，對(duì)函數(shù)f(·)進(jìn)行辛普森積分，如下式所示：

為了驗(yàn)證3個(gè)積分點(diǎn)計(jì)算結(jié)果是否滿足精度，分別在[a,c]和[c,b]域內(nèi)中點(diǎn)增加積分點(diǎn)d和e，分別在這兩個(gè)區(qū)間采用上式積分并求和為Q_s，如果Q_s和Q的相對(duì)誤差δ滿足容許誤差ε，即3個(gè)積分點(diǎn)計(jì)算結(jié)果滿足計(jì)算精度；如果Q_s和Q的相對(duì)誤差δ不滿足容許誤差ε，則在在[a,c]和[c,b]域內(nèi)分別按照上述方法驗(yàn)證各自區(qū)間內(nèi)是否滿足計(jì)算精度，如此循環(huán)，直至滿足精度為止。這樣就可以保證被積函數(shù)的最大/小值兩側(cè)的積分點(diǎn)之間的函數(shù)為凸函數(shù)，以便進(jìn)一步進(jìn)行牛頓法迭代得到最優(yōu)解，如圖3所示。

步驟(9)判斷步驟(8)中積分的結(jié)果是否滿足精度，判斷是否滿足精度的方法是判斷增加積分點(diǎn)前后計(jì)算得到的積分值的相對(duì)誤差δ是否滿足容許誤差ε要求。如果積分結(jié)果精度滿足容許誤差ε進(jìn)入下一步；如果不滿足計(jì)算精度則返回到步驟(3)進(jìn)行新一輪的循環(huán)。其中，容許誤差ε一般設(shè)置為10^-3，而相對(duì)誤差δ是迭代過程中前后兩次計(jì)算的積分值只差與前一次積分值的比值，即為：

其中，Q和Q_s分別是迭代前后的積分值。

步驟(10)選取現(xiàn)有的所有積分點(diǎn)集合中所對(duì)應(yīng)目標(biāo)函數(shù)的最大/小值及對(duì)應(yīng)的積分點(diǎn)，如圖3所示，圖3(a)為全局的積分點(diǎn)分布情況，選取這些積分點(diǎn)中目標(biāo)函數(shù)的最大/小值的積分點(diǎn)，如圖3(b)所示。

步驟(11)以步驟(11)選取的積分點(diǎn)作為起點(diǎn)，由于該點(diǎn)到最優(yōu)解處之間為凸函數(shù)甚至是單調(diào)函數(shù)，故可以采用采用收斂速度較快的牛頓法進(jìn)行局部尋優(yōu)，得到的最優(yōu)解即為結(jié)構(gòu)響應(yīng)的上/下界。

步驟(12)結(jié)束。

實(shí)施例：

為了充分地了解本發(fā)明的特點(diǎn)及其對(duì)工程實(shí)際的適用性，本發(fā)明建立如圖所示的包含彈性模量和荷載不確定性的桁架結(jié)構(gòu)，并對(duì)結(jié)構(gòu)區(qū)間位移傳播進(jìn)行分析。該桁架結(jié)構(gòu)及相關(guān)物理參數(shù)如圖4所示，其中，桿①、⑤和⑥的彈性模量的不確定度為恒荷載的不確定度為活荷載的不確定度為本發(fā)明的方法和蒙特卡洛法、區(qū)間參數(shù)攝動(dòng)法的結(jié)果作對(duì)比，計(jì)算結(jié)果如表1所示。通過對(duì)比發(fā)現(xiàn)三者的計(jì)算結(jié)果相差無幾，并且從表3看出，區(qū)間參數(shù)攝動(dòng)法具有一定優(yōu)勢(shì)。但是這是針對(duì)較小的不確定度的工況，系統(tǒng)的非線性較弱。當(dāng)結(jié)構(gòu)中不確定度較大時(shí)，就會(huì)使系統(tǒng)的非線性增強(qiáng)，區(qū)間參數(shù)攝動(dòng)法對(duì)小參數(shù)的強(qiáng)烈依賴性就體現(xiàn)出來了。當(dāng)桿①、⑤和⑥的彈性模量的不確定度為恒荷載的不確定度為活荷載的不確定度為時(shí)，計(jì)算的結(jié)果如表2所示。本發(fā)明方法較蒙特卡洛法調(diào)用計(jì)算有限元計(jì)算次數(shù)更少，計(jì)算量較少具有很大優(yōu)勢(shì)。

表1小不確定度工況下桁架3、4節(jié)點(diǎn)位移(單位：×10^-5m)

表2大不確定度工況下桁架3、4節(jié)點(diǎn)位移(單位：×10^-5m)

表3各方法調(diào)用有限元程序次數(shù)對(duì)比

本發(fā)明未詳細(xì)闡述部分屬于本領(lǐng)域技術(shù)人員的公知技術(shù)。

以上僅是本發(fā)明的具體步驟，對(duì)本發(fā)明的保護(hù)范圍不構(gòu)成任何限制；其可擴(kuò)展應(yīng)用于非線性系統(tǒng)區(qū)間相應(yīng)的不確定性傳播分析中，凡采用等同變換或者等效替換而形成的技術(shù)方案，均落在本發(fā)明權(quán)利保護(hù)范圍之內(nèi)。