葛新廣 盧嘉康 張梨榮 羅臻

摘 要：基于頻域法研究幾類隨機激勵下工程結(jié)構(gòu)的動力響應(yīng)，需要求解4種含符號函數(shù)在無窮區(qū)間的積分解，而該4種積分式目前存在計算效率和精度的問題。首先，根據(jù)積分運算法則將4種積分式的計算轉(zhuǎn)換為一種積分式的計算；其次，基于留數(shù)定理和高斯-雅克比積分提出了4種積分式的高效簡明數(shù)值解法；最后，通過算例研究了傳統(tǒng)方法受積分上限和積分區(qū)間的影響、本文方法受高斯-雅克比積分節(jié)點數(shù)的影響。結(jié)果表明：傳統(tǒng)方法的計算結(jié)果受積分區(qū)間和積分上限的影響較大，而本文方法取25個高斯-雅克比積分節(jié)點數(shù)即可獲得較高的精度，且具有較高的計算效率。

關(guān)鍵詞：無窮積分；符號函數(shù)；高斯-雅克比積分；高效數(shù)值解

中圖分類號：TU311.3 DOI：10.16375/j.cnki.cn45-1395/t.2024.03.007

0 引言

隨機振動現(xiàn)象廣泛存在于工程領(lǐng)域，如航空器在飛行過程中因風阻導(dǎo)致的振動[1]、土木工程領(lǐng)域的結(jié)構(gòu)地震動[2]和風振動[3]、車輛在路面運行過程中因路面不平順導(dǎo)致的振動[4-5]等。頻域法是分析各類隨機激勵下結(jié)構(gòu)構(gòu)件動力響應(yīng)的重要方法之一[6-7]。林家浩等[8]基于頻域法提出了虛擬激勵法，該方法廣泛應(yīng)用于各種結(jié)構(gòu)的隨機激勵下的響應(yīng)分析[9]，是研究隨機地震動激勵下土木工程結(jié)構(gòu)的主要方法之一。實際地震動過程具有典型的非平穩(wěn)特性，土木工程結(jié)構(gòu)基于非平穩(wěn)地震動激勵下的動力響應(yīng)分析[10]是評估結(jié)構(gòu)地震過程安全的重要內(nèi)容。大跨度土木工程結(jié)構(gòu)在地震動激勵下具有典型的行波效應(yīng)[11]，考慮行波效應(yīng)比不考慮行波效應(yīng)的此類結(jié)構(gòu)地震響應(yīng)更加顯著?；疖嚮谲壍啦黄巾樇頪12]和汽車基于路面不平順激勵下動力響應(yīng)[13]分析需要考慮相干函數(shù)，考慮相干函數(shù)的影響會比不考慮行波效應(yīng)的交通工具更加顯著。上述工程領(lǐng)域的分析需要計算一類含有符號函數(shù)[sgn（ω）egiωtkiω-q-n]（sgn為符號函數(shù)，[ω]為積分變量，g=±1，i=[-1]，t為振動時長，k=±1，q為實部為負值的復(fù)數(shù)或者負實數(shù)，n為正整數(shù)）的無窮積分[14-16]，而目前針對這種積分所采用的方法較為復(fù)雜[15]，或者采用數(shù)值解的計算精度不易控制。張志超等[17]基于虛擬激勵法-精細積分法研究了車橋耦合系統(tǒng)的非平穩(wěn)激勵下的車輛系統(tǒng)響應(yīng)，所得解為數(shù)值解，其分析精度受時間步長和頻率上限的影響較大，需要試算才能確定精度。Barbato等[14]研究了非平穩(wěn)隨機地震動激勵下結(jié)構(gòu)時變可靠度的非幾何譜矩的計算，所得封閉解含有積分表達式，而該積分表達式是針對實部較小的復(fù)特征值，容易出現(xiàn)偏差較大的結(jié)果。葛新廣[16]在Barbato[14]所提方法的基礎(chǔ)上提出了上述4個表達式的改進計算方法，并將其應(yīng)用于分數(shù)導(dǎo)數(shù)阻尼器耗能結(jié)構(gòu)基于Conte-Peng完全非平穩(wěn)激勵下的響應(yīng)分析，但所獲得的表達式依然復(fù)雜，且存在對q值的實部較小時計算失真的情況。因此，有必要針對此類含有符號函數(shù)的積分提出更高效的計算方法。

高斯-雅克比積分法具有分析精度高和效率高的特點，廣泛應(yīng)用于各種復(fù)雜表達式的積分計算[16，18]。本文根據(jù)[sgn（ω）egiωtkiω-q-n]式中g(shù)、k的取值，首先研究了4種式子之間的關(guān)系式，并將4種積分式的計算轉(zhuǎn)換為一種積分式的計算；其次，基于留數(shù)定理和高斯-雅克比積分提出了一種積分式的表達式的解，同時將其推廣至其他3種表達式的計算；最后通過算例驗證了所提方法的正確性和高效性。

1 4個無窮積分表達式的關(guān)系

根據(jù)g、k的取值，可將[sgn（ω）egiωtkiω-q-n]展開為4個無窮積分表達式

[A1q=-∞∞sgn（ω）eiωtiω-qndω ，] （1）

[A2q=-∞∞sgn（ω）e-iωtiω-qndω]， （2）

[A3q=-∞∞sgn（ω）eiωt-iω-qndω]， （3）

[A4q=-∞∞sgn（ω）e-iωt-iω-qndω]. （4）

對式（2）積分進行變換得

[A2q=-∞∞sgn（-ω）eiωt-iω-qnd-ω]. （5）

調(diào)換積分上下限及利用符號函數(shù)性質(zhì)，式（5）改寫為

[A2q=--∞∞sgn（ω）eiωt-iω-qndω]. （6）

為獲得A2（q）與A1（q）的關(guān)系，對式（6）進一步改寫為

[A2q=（-1）n-1-∞∞sgn（ω）eiωtiω+qndω]. （7）

對比式（7）與式（1）的關(guān)系，則存在

[A2q=（-1）n-1A1-q]. （8）

式（8）給出了A2（q）與A1（q）的關(guān)系式。

對式（3）進行積分變換可得

[A3q=-∞∞sgn（-ω）e-iωtiω-qnd-ω]. （9）

調(diào)整積分上下限并利用符號函數(shù)性質(zhì)，對式（9）改寫為

[A3q=--∞∞sgn（ω）e-iωtiω-qndω]. （10）

比較式（2）與式（10），則

[A3q=-A2q]. （11）

利用式（8），則建立A1（q）與A3（q）的關(guān)系式

[A3q=（-1）nA1-q]. （12）

對式（4）進行積分變換可得

[A4q=-∞∞sgn（-ω）eiωtiω-qnd-ω]. （13）

調(diào)整積分上下限并利用符號函數(shù)性質(zhì)，對式（13）改寫為

[A4q=--∞∞sgn（ω）eiωtiω-qndω]. （14）

比較式（14）與式（1），則

[A4q=-A1q]. （15）

式（15）給出了A1（q）與A4（q）的關(guān)系式。

式（8）、式（12）、式（15）、A2（q）—A4（q）均與A1（q）建立了關(guān)系，為此若獲得式（1）—式（4）的積分解，只需要獲得式（1）的積分解。

2 4個無窮積分的近似解

2.1 A1（q）的近似解

對式（1）進行改寫得

[A1q=-in-∞∞sgn（ω）eiωtω+iqndω=-inA11+A12，]? （16）

式中：

[A11=0∞eiωtω+iqndω]， （17）

[A12=-∞0eiωtω+iqndω]. （18）

利用留數(shù)定理并結(jié)合式（17）及式（18）的特點，A11和A12的計算如圖1所示。

[C5][C1][C3][C4][C2][實軸][虛軸] [O]

圖1 留數(shù)定理應(yīng)用示意圖

[A11=0∞eiωtω+iqndω=C1eiωtω+iqndω，] （19）

[A12=-∞0eiωtω+iqndω=C5eiωtω+iqndω]. （20）

由留數(shù)定理，可構(gòu)造積分函數(shù)為

[（21）]

[（22）]

[[式中：Rey=max0，ReyRey；????????? ???????Imy=max0，ImyImy；] ]

Re（y）、Im（y）分別表示求y的實部和虛部。

由式（21）及式（22），則A11與A12可表示為

[（23）]

[（24）]

對式（23）及式（24），C2—C4積分為

[C2eiztz+iqndz=limr→∞0π2eireiθtrieiθreiθ+iqndθ=0]，?? （25）

[（26）]

[C4eiztz+iqndz=limr→∞π2πeireiθtrieiθreiθ+iqndθ=0]，??? （27）

式中：r為留數(shù)定理中的半徑。

把式（25）—式（27）代入式（23）和式（24），則A11和A12表示為

[A11=Re-iqIm-iq2πintn-1eqtn-1！+-in-1D]，?? （28）

[A12=-2πintn-1eqtn-1！ReiqIm-iq+-in-1D]，?? （29）

式中：

[D=0∞e-yty+qndy]. （30）

對式（30）利用高斯-雅克比積分進行高精度近似計算，首先對積分區(qū)間進行變換，令[x=1-y1+y]，則[y=1-x1+x]，代入式（30），則

[D=-11e-1-x1+xt1-x1+x+qn-21+x2dx]. （31）

對式（31）進一步簡化可得

[D=2-11（1+x）n-2e-1-x1+xt1+q+q-1xndx]. （32）

利用高斯-雅克比積分對式（30）進行改寫，即

[D≈2j=1N（1+xj）n-2Aje-1-xj1+xjt1+q+q-1xjn ，] （33）

式中：[xj、Aj]為權(quán)函數(shù)[ρ=1-x01+x0]的高斯-雅克比積分高斯點及系數(shù)；N為積分節(jié)點數(shù)。

把式（28）、式（29）及式（33）代入式（16），則

[A1q≈Im-iqRe-iq-Reiq2πtn-1eqtn-1！+4-1nij=1N（1+xj）n-2Aje-1-xj1+xjt1+q+q-1xjn.] （34）

由于工程上的q為負實數(shù)或者實部為負數(shù)的復(fù)數(shù)，為此可對式（34）進一步簡化為

[A1q≈2πts-1eqtμs-1！+???????????????????? 4-1sij=1N（1+xj）n-2Aje-1-xj1+xjt1+q+q-1xjn，] （35）

式中：[μ=Imq/Imq]。

由式（35）可知，本文所提方法將無窮積分轉(zhuǎn)化為高斯-雅克比的積分，相對于文獻[14]和文獻[16]具有簡潔性。

2.2 A2（q）、A3（q）和A4（q）的近似解

由式（8）及式（34），則A2（q）表示為

[A2q≈（-1）n-1ImiqReiq-Re-iq×2πtn-1e-qtn-1！-4ij=1N（1+xj）n-2Aje-1-xj1+xjt1-q-q+1xjn.] （36）

由于工程上的q為負實數(shù)或者實部為負數(shù)的復(fù)數(shù)，為此可將式（36）進一步簡化為

[A2q=-4ij=1N（1+xj）n-2Aje-1-xj1+xjt1-q-q+1xjn]. （37）

由式（11）及式（37），則當q為負實數(shù)或者實部為負數(shù)的復(fù)數(shù)，則A3（q）表示為

[A3q=4ij=1N（1+xj）n-2Aje-1-xj1+xjt1-q-q+1xjn]. （38）

由式（15）及式（35），則當q為負實數(shù)或者實部為負數(shù)的復(fù)數(shù)，則A4（q）表示為

[A4q≈-2πtn-1eqtμn-1！+4-1n-1×???????????????????? ij=1N（1+xj）n-2Aje-1-xj1+xjt1+q+q-1xjn.] （39）

至此，獲得了當q為負實數(shù)或者實部為負數(shù)的復(fù)數(shù)，則得到A1（q）—A4（q）的近似解。

3 算例

3.1 計算方法精度及效率驗證

為驗證本文所提方法的精度和分析效率，利用傳統(tǒng)的數(shù)值積分進行了對比驗證。數(shù)值積分計算A1（q）—A4（q）的計算表達式為

[A1q=12k=0ωu/Δω-1ei×k×Δωti×k×Δω-qn-e-i×k×Δωt-i×k×Δω-qn+ei×k+1×Δωti×k+1×Δω-qn-e-i×k+1×Δωt-i×k+1×Δω-qnΔω]， （40）

[A2q=12k=0ωu/Δω-1e-i×k×Δωti×k×Δω-qn-ei×k×Δωt-i×k×Δω-qn+e-i×k+1×Δωti×k+1×Δω-qn-ei×k+1×Δωt-i×k+1×Δω-qnΔω]， （41）

[A3q=12k=0ωu/Δω-1ei×k×Δωt-i×k×Δω-qn-e-i×k×Δωti×k×Δω-qn+ei×k+1×Δωt-i×k+1×Δω-qn-e-i×k+1×Δωti×k+1×Δω-qnΔω]， （42）

[A4q=12k=0ωu/Δω-1e-i×k×Δωt-i×k×Δω-qn-ei×k×Δωti×k×Δω-qn+e-i×k+1×Δωt-i×k+1×Δω-qn-ei×k+1×Δωti×k+1×Δω-qnΔω]， （43）

式中：[ωu、Δω]分別為積分上限和積分步長。

從式（40）—式（43）可知，采用數(shù)值方法需要給定積分上限和積分間距。為驗證所提方法的正確性，A1（q）—A4（q）中的參數(shù)為n、q及t，故本算例的參數(shù)取值為：1）q=-50±100i，t=0.1 s， n=9和2） q=-300±400i， t=10.0 s， n=9。分別從積分上限和積分步長來研究本文方法的正確性，本文方法取30個高斯積分節(jié)點。

3.1.1 基于積分上限的驗證

經(jīng)試算積分步長[Δω]取10-4 rad/s，3種工況：工況1—工況3積分上限[ωu]取102、103、104 rad/s，計算結(jié)果如表1—表4。工況3的耗時為619.640 s，本文方法的耗時為0.032 s，本文方法具有較高的計算效率。

從表1—表4可知，當積分步長[Δω]=10-4 rad/s、q=-50±100i時，計算A1（q）—A4（q）時積分上限在3種工況下均一致，且與本文方法誤差很小，說明[ωu]=102 rad/s時計算結(jié)果即可獲得精確解；當q=-300±400i時，計算A1（q）—A4（q）時3種工況積分上限均不同，但當[ωu]=104 rad/s時最接近本文方法。

綜上所述，本文方法在計算A1（q）—A4（q）時不受積分上限的影響，具有精度高的特點，也說明了本文方法的正確性。

3.1.2 基于積分步長的驗證

通過3.1.1可知，積分上限[ωu]=104 rad/s時，計算q=-50±100i和q=-300±400i可獲得高精度解。為此，本算例的計算基于積分步長的驗證[ωu]取104 rad/s。3種工況積分步長[Δω]?。?0-2、10-3、10-4 rad/s，計算結(jié)果如表5—表8所示。

由表5—表8可知，在積分上限[ωu]取104 rad/s、q=-50±100i時，計算A1（q）—A4（q）積分上限[Δω=10-3]和[10-4 rad/s]的計算結(jié)果一致且與本文方法比較接近，說明獲得了精確解；在q=-300±400i時，計算A1（q）—A4（q）時3種工況積分上限均不同，但當[Δω]=10-4 rad/s時最接近本文方法。綜上所述，本文方法在計算A1（q）—A4（q）時不受積步長的影響，具有精度高的特點，也說明了本文方法的正確性。

3.2 計算精度與高斯-雅克比積分節(jié)點數(shù)的關(guān)系

高斯系列積分的精度受積分節(jié)點個數(shù)的影響，積分節(jié)點個數(shù)越大，精度越高；為此本文研究了高斯積分節(jié)點個數(shù)對計算精度的影響。根據(jù)3.1節(jié)的結(jié)論，本文方法在計算參數(shù)分別為1） q=-50+100i、t=0.1 s， n=9和2）q=-300+400i、t=10.0 s、n=10時計算精度高，故下面的研究取上述參數(shù)進行。高斯積分節(jié)點數(shù)取5～30個。圖2給出了q=-50+100i、t=0.1 s、n=9的積分節(jié)點數(shù)對計算結(jié)果的誤差分析。圖3給出了q=-300+400i、t=10.0 s、n=9的積分節(jié)點數(shù)對計算結(jié)果的誤差分析。

由圖2可知，在q=-50+100i、t=0.1 s、n=9時，在計算A2（q）和A3（q）時，積分節(jié)點數(shù)較少時精度較差，且實部的影響比虛部更大；而計算A1（q）和A4（q）時計算精度受高斯積分節(jié)點數(shù)影響不大。由圖3可知，q=-300+400i、t=10.0 s、n=10時，計算所有的積分式時精度受高斯積分節(jié)點個數(shù)影響不大。對比圖2和圖3，高斯積分節(jié)點的個數(shù)受q的影響會有所不同，q的實部或虛部絕對值越小受高斯積分節(jié)點個數(shù)的影響越大；高斯積分節(jié)點的個數(shù)取25個時可具有較好的精度。

4 結(jié)論

本文針對振動領(lǐng)域里4個含有符號函數(shù)[sgn（ω）egiωtkiω-q-n]的無窮積分的表達式求解方法存在的問題，提出了一種高效算法，通過算例分析獲得如下結(jié)論：

1）利用本文方法分析4個積分式的解時，無需考慮積分間距和積分上限的影響。而使用傳統(tǒng)數(shù)值方法計算時，則受積分間距和積分上限的影響，且不同q值的積分解的精度受積分間距和積分上限的影響也會不同。因此，采用傳統(tǒng)數(shù)值積分時可能會出現(xiàn)偏差很大的解，這一點需要引起注意。

2）本文所提解法的精度受高斯-雅克比積分節(jié)點個數(shù)的影響較大，q的實部或虛部絕對值越小受高斯積分節(jié)點個數(shù)的影響越大。從論文所獲得結(jié)果來看，當取25個高斯積分節(jié)點，積分精度可不受q值實部或虛部的影響。為此，建議采用本文方法計算此類含有符號函數(shù)積分式的解時取25個高斯積分節(jié)點。

參考文獻

[1] 衡楠，鄭鋒華，李金鵬．直升機機載掃描反射鏡組件隨機振動仿真分析[J]．機械設(shè)計與制造，2020（6）：290-293，297.

[2] 李創(chuàng)第，楊雪峰，李宇翔，等．組合體結(jié)構(gòu)隨機地震動響應(yīng)的簡明解法[J]．地震工程與工程振動，2023，43（1）：24-33.

[3] 葛新廣，張夢丹，龔景海，等．頻響函數(shù)二次正交法在Davenport風速譜下結(jié)構(gòu)系列響應(yīng)簡明封閉解的應(yīng)用研究[J]．振動與沖擊，2021，40（21）：207-214.

[4] 謝娟娟，李晉，田震，等．考慮路面不平順隨機性的汽車過橋動力響應(yīng)分析[J]．振動與沖擊，2021，40（14）：299-306.

[5] 李創(chuàng)第，李宇翔，葛新廣．路面不平順激勵下車輛系統(tǒng)振動響應(yīng)譜矩的簡明封閉解[J]．廣西科技大學(xué)學(xué)報，2021，32（3）：13-20.

[6] 李創(chuàng)第，賀王濤，葛新廣.卷積型非粘滯阻尼結(jié)構(gòu)隨機地震動系列響應(yīng)求解的虛擬激勵法[J].廣西科技大學(xué)學(xué)報，2021，32（1）：78-84.

[7] 趙中偉，張永高.基于本征-虛擬激勵法的大跨鋼結(jié)構(gòu)風振響應(yīng)分析[J].空間結(jié)構(gòu)，2020，26（1）：15-23.

[8] 林家浩，張亞輝.隨機振動的虛擬激勵法[M].北京：科學(xué)出版社，2004.

[9] 林家浩，張亞輝，趙巖.虛擬激勵法在國內(nèi)外工程界的應(yīng)用回顧與展望[J].應(yīng)用數(shù)學(xué)和力學(xué)，2017，38（1）：1-32.

[10] ALDERUCCI T，MUSCOLINO G．Fully nonstationary analysis of linear structural systems subjected to multicorrelated stochastic excitations[J]．ASCE-ASME Journal of Risk and Uncertainty in Engineering Systems，Part A：Civil Engineering，2016，2（2）：C4015007.

[11] 張文首，林家浩.大跨度結(jié)構(gòu)考慮行波效應(yīng)時平穩(wěn)隨機地震響應(yīng)的閉合解[J].固體力學(xué)學(xué)報，2004，25（4）：446-450.

[12] 何彬彬，馮磊，馮玉林，等.不同軌道不平順譜激勵下高鐵列車走行性分析[J].鐵道工程學(xué)報，2022，39（7）：36-43.

[13] 桂水榮，陳水生，萬水.基于路面一致激勵車橋耦合非平穩(wěn)隨機振動分析[J].振動.測試與診斷，2018，38（5）：908-915，1077.

[14] BARBATO? M，CONTE J P.Time-variant reliability analysis of linear elastic systems subjected to fully nonstationary stochastic excitations[J].Journal of Engineering Mechanics，2015，141（6）：04014173.

[15] BARBATO? M，VASTA? M．Closed-form solutions for the time-variant spectral characteristics of non-stationary random processes[J].Probabilistic Engineering Mechanics，2010，25（1）：9-17.

[16] 葛新廣.分數(shù)導(dǎo)數(shù)型阻尼耗能結(jié)構(gòu)地震動響應(yīng)的分析方法研究[D].上海：上海交通大學(xué)，2022.

[17] 張志超，張亞輝，林家浩.車橋耦合系統(tǒng)非平穩(wěn)隨機振動分析的虛擬激勵-精細積分法[J].工程力學(xué)，2008，25（11）：197-204.

[18] DIETHELM K．An improvement of a nonclassical numerical method for the computation of fractional derivatives[J]．Journal of Vibration and Acoustics，2009，131（1）：014502.

Efficient numerical solutions for a class of infinite

integrals with signed functions

GE Xinguang1， LU Jiakang1， ZHANG Lirong1， LUO Zhen*2

（1. School of Civil and Architecture Engineering， Liuzhou Institute of Technology，

Liuzhou 545616， China; 2. School of Civil and Architecture Engineering， Guangxi University of

Science and Technology， Liuzhou 545006， China）

Abstract： To study the dynamic response of engineering structures under some random excitation based on frequency domain method， it is necessary to solve four infinite integrals with signed functions for which， however， there is computational efficiency and accuracy problem. Firstly， the calculation of four integral equations was converted into that of one integral equation according to the integral arithmetic. Secondly， an efficient and concise numerical solution for the above integral equations was presented based on the residue theorem and Gaussian-Jacobian integration. Finally， through numerical examples， the effect of the integration upper limit and integration interval on the traditional method and that of the number of Gaussian-Jacobian integration nodes on the method in this paper were studied. The results show that the calculation results of the traditional method are greatly affected by the integration interval and integration upper limit. However， the method in this paper can achieve high accuracy and computational efficiency by taking 25 Gaussian-Jacobian integration nodes.

Keywords： infinite integral; signed functions; Gaussian-Jacobian integral; efficient numerical solution

（責任編輯：羅小芬）

收稿日期：2023-10-11；修回日期：2023-10-25

基金項目：國家自然科學(xué)基金項目（51368005）；柳州工學(xué)院高層次人才項目（202201）資助

第一作者：葛新廣，博士，副教授，研究方向：土木工程結(jié)構(gòu)抗震、抗風研究

*通信作者：羅臻，碩士，高級工程師，研究方向：土木工程結(jié)構(gòu)抗震、抗風研究，E-mail：20507257@qq.com