盛慶軻 梁海琴

摘要：達(dá)芬方程作為典型的非線性振動方程，可以直觀地用來說明與線性振動的區(qū)別。分別采用有限元程序ANSYS和編程語言Python對達(dá)芬方程進(jìn)行求解。針對ANSYS界面操作效率低、數(shù)據(jù)后處理能力不足等問題。以PyQt5為平臺，封裝ANSYS和Python的求解以及數(shù)據(jù)處理過程，并制定可視化操作界面。程序包含達(dá)芬方程參數(shù)輸入、ANSYS計算、Python計算和數(shù)據(jù)后處理四個功能。以多組系統(tǒng)參數(shù)為例，使用該程序計算達(dá)芬方程，計算結(jié)果表現(xiàn)了非線性振動中初值影響、多穩(wěn)態(tài)解和諧波響應(yīng)運動特性。

關(guān)鍵詞：達(dá)芬方程;非線性振動;有限元;可視化設(shè)計

中圖分類號：TP391.9 文獻(xiàn)標(biāo)識碼：A 文章編號：1007-9416（2020）08-0076-04

0 引言

機載設(shè)備往往由多個零部件構(gòu)成，隨著科學(xué)技術(shù)的迅速發(fā)展，零件自身的復(fù)雜性和連接關(guān)系也趨于復(fù)雜。目前對于機載設(shè)備的振動強度主要采用有限元軟件進(jìn)行計算校核，然而有限元的振動計算是基于小變形線性理論。當(dāng)結(jié)構(gòu)簡單、振幅較小、連接關(guān)系變化不大時，按線性理論是計算可以得到較為滿意的結(jié)果[1]。

一切力學(xué)問題在本質(zhì)上都屬于非線性問題。對于非線性系統(tǒng)，響應(yīng)和輸出之間并不是線性關(guān)系，按照線性理論計算得到的結(jié)果精度可能不足，甚至導(dǎo)致本質(zhì)上的錯誤。非線性振動會出現(xiàn)幅頻曲線彎曲、諧波、分岔混沌共振等線性振動不存在的現(xiàn)象[2]。在非線性振動振動理論研究中，具有代表性的方程有達(dá)芬方程、Van-der-Pol方程等[3]。

達(dá)芬方程對應(yīng)的動力學(xué)表達(dá)式如式（1）所示，也稱為達(dá)芬系統(tǒng)，其中剛度包含立方項。

（1）

其中：x代表位移，m代表質(zhì)量，c代表阻尼，k1代表線性剛度，k2代表非線性剛度，f（t）代表外載荷。

ANSYS的腳本語言APDL能夠?qū)崿F(xiàn)有限元模型的參數(shù)化建模計算功能，但其可視化程度較低，一般基于其他編程語言對其進(jìn)行二次界面開發(fā)。Python是一種面向?qū)ο蟮膭討B(tài)解釋型語言，應(yīng)用領(lǐng)域廣泛[4]。Qt是由C++語言編寫的跨平臺GUI庫，PyQt5則是將Python語言與Qt融合，具有Python開發(fā)效率高的特點[5]。

本文以PyQt5為平臺，基于Python對ANSYS求解達(dá)芬方程的過程進(jìn)行封裝并設(shè)計可視化操作程序，同時加入了基于Python科學(xué)計算庫的達(dá)芬方程計算和結(jié)果顯示功能。

1 計算程序的功能設(shè)計

達(dá)芬方程計算程序主要包含四個功能，分別是參數(shù)輸入、提交ANSYS計算、提交Python計算和結(jié)果顯示模塊。由于Python語言自身的特點，許多功能的實現(xiàn)需要引用Python標(biāo)準(zhǔn)庫和第三方庫。程序的總體設(shè)計框圖如圖1所示。

（1）參數(shù)輸入包括質(zhì)量系數(shù)、阻尼系數(shù)、剛度線性系數(shù)和立方項系數(shù)、載荷幅值和頻率、初始位移和速度、求解時間和步長。（2）提交ANSYS計算前提是用戶電腦需要裝有ANSYS軟件，然后設(shè)置ANSYS.exe的啟動位置、工作目錄和APDL文件名即可，程序會在指定目錄下生成文件并且調(diào)用ANSYS主程序進(jìn)行求解。（3）提交Python計算只需前面的參數(shù)輸入完畢后，利用Python科學(xué)計算庫求解即可。（4）對于非線性微分方程，一般很難找到精確的解析解，大都是從定性和定量兩個方向來研究。相平面法[6]是常用的定性方法之一，它是描述系統(tǒng)速度與位移的對應(yīng)關(guān)系，可以直接判斷方程解的平衡性、周期性和穩(wěn)定性等。

1.1 基于Python實現(xiàn)ANSYS二次開發(fā)

利用Python對ANSYS二次開發(fā)的主要過程分為提取界面輸出參數(shù)及相關(guān)信息、生成指定參數(shù)下的APDL文件和提交ANSYS計算。

1.1.1 APDL命令流文件的生成

在APDL建模命令流之前引入相關(guān)的參數(shù)化變量，然后對該字符串進(jìn)行正確的格式化，即可得到特定輸入下特定的ADPL文件。

本文編寫的APDL腳本程序主要采用質(zhì)量單元mass21和非線性彈簧單元combin39，分析類型為瞬態(tài)非線性，根據(jù)初始條件中的位移和速度是否為0，生成不同的APDL命令流以模擬相應(yīng)的初始條件。

1.1.2 Python調(diào)用ANSYS程序

首先在計算程序設(shè)置APDL文件的存儲位置和文件名稱，程序即可在指定目錄下生成2.1中描述的APDL文件，主要代碼如下所示：

file = str #獲取存儲位置的絕對路徑

with open （file，） as f：? #打開文件

f.write（s）? #寫入2.1中的字符串

接下來需要設(shè)置啟動ANSYS.exe的路徑，依次生成執(zhí)行批處理的Python文件。在Python中采用OS標(biāo)準(zhǔn)庫完成系統(tǒng)操作，并且依據(jù)OS庫執(zhí)行的接口可以實現(xiàn)跨平臺訪問。生成批處理文件的代碼如下所示：

最后再利用os.system命令運行上述的批處理文件，即可完成Python對ANSYS的調(diào)用計算。

1.2 提交Python計算

當(dāng)達(dá)芬方程的系統(tǒng)參數(shù)和求解參數(shù)輸入完成后，即可提交Python進(jìn)行求解。整個過程基于Python的科學(xué)計算庫[7]Numpy和Scipy，求解達(dá)芬方程采用Scipy庫的odeint函數(shù)，該函數(shù)是基于龍格-庫塔算法編寫的。

1.3 結(jié)果顯示

達(dá)芬方程計算程序的結(jié)果顯示功能是基于Python科學(xué)計算庫Matplotlib完成的，通過讀取ANSYS和Python的計算結(jié)果進(jìn)行顯示。

2 可視化界面設(shè)計

達(dá)芬方程計算程序的界面設(shè)計基于Python的第三方庫—PyQt5，基于PyQt5的界面與邏輯分離方法，利用Qt Designer設(shè)計程序頁面，主程序編寫各控件的功能以及相互之間的信號傳遞函數(shù)等實現(xiàn)整個程序的開發(fā)，最終效果如圖2和圖3所示。

（1）達(dá)芬方程參數(shù)輸入。參數(shù)輸入采用Qt內(nèi)置的Qline Edit本文輸入框來接收，并添加了數(shù)據(jù)檢測功能。（2）ANSYS計算?？紤]到ANSYS計算達(dá)芬方程是一個瞬態(tài)非線性分析類型，計算步長會計算結(jié)果的收斂性和時長影響較大，程序能夠根據(jù)用戶輸入的達(dá)芬方程參數(shù)給出一個范圍供用戶選擇;為了能夠合理地利用計算機資源，程序會基于Psutil庫獲取CPU的邏輯個數(shù)，給出范圍供用戶選擇。這兩個功能都調(diào)用QSlider滑動條類完成。（3）Python計算。這部分由一個QButton按鈕類和QLabel標(biāo)簽類組成，當(dāng)提交Python計算完成時，標(biāo)簽將會顯示“計算完成”。（4）結(jié)果顯示。在PyQt5中引入Matplotlib庫實現(xiàn)對計算結(jié)果的查看。并采用QtableView表格類讀取達(dá)芬方程的時域計算結(jié)果，并提供保存功能，可對結(jié)果進(jìn)行更深入的分析與繪制圖形。

3 實例分析

本節(jié)以具體數(shù)據(jù)為例，利用本程序計算達(dá)芬方程并根據(jù)結(jié)果描述非線性振動的一些特點。

受迫振動是指系統(tǒng)在周期性的外力作用下所發(fā)生的的振動。以表1所示的三組方程參數(shù)為例，說明達(dá)芬系統(tǒng)受迫振動的一些特點，計算受迫振動時需要一定的計算時長，可根據(jù)響應(yīng)曲線大致判斷響應(yīng)的穩(wěn)定性。

對于線性振動，受迫振動響應(yīng)必然是與載荷相同頻率的簡諧運動。當(dāng)系統(tǒng)本身和載荷參數(shù)確定后，系統(tǒng)的振動響應(yīng)就已確定，響應(yīng)幅A值與載荷頻率呈現(xiàn)一一對應(yīng)的關(guān)系。

系統(tǒng)2-1的計算結(jié)果如圖4所示，可以看出當(dāng)初始位移不同時，達(dá)芬方程受迫振動的幅值可能不同。即達(dá)芬系統(tǒng)的響應(yīng)幅值與頻率不再呈現(xiàn)一一對應(yīng)關(guān)系，這是非線性振動的多穩(wěn)定解現(xiàn)象。

系統(tǒng)2-2的計算結(jié)果如圖4所示，可以看出當(dāng)初始位移為1或3時（這兩者的差異是由于初值不同造成的相位差導(dǎo)致），幅值不僅與初值為5不同，而且已經(jīng)不再是簡諧運動。

取初始位移等于1和5，進(jìn)一步查看系統(tǒng)2-2穩(wěn)態(tài)階段的相平面圖，如圖5所示。從相平面圖中可以看出初始位移為1時相軌跡已經(jīng)從橢圓發(fā)生了畸變，這是非線性振動的諧波現(xiàn)象。

4 結(jié)論

（1）以非線性振動方程—達(dá)芬方程為研究對象，基于有限元程序ANSYS和編程語言Python編寫了可視化計算程序。該程序能夠快速的完成不同參數(shù)下達(dá)芬方程的求解與結(jié)果顯示，降低了軟件使用難度并提高分析效率。（2）介紹了基于Python語言對ANSYS二次開發(fā)的主要步驟，以及該程序的主要功能和結(jié)構(gòu)。（3）以具體參數(shù)為例，說明非線性振動與線性振動的區(qū)別，包括初值影響、多穩(wěn)定解和諧波運動。

參考文獻(xiàn)

[1] 黃安基.非線性振動[M].成都：西南交通大學(xué)出版社，1993.

[2] 王海波.Duffing方程非線性振動特性的計算與分析[D].西安：西安建筑科技大學(xué)，2009.

[3] 許磊，陸明萬，曹慶杰.Van der Pol-Duffing方程的非線性動力學(xué)分叉特性研究[J].應(yīng)用力學(xué)學(xué)報，2002（4）：130-133+168.

[4] 蔡敏.Python語言的Web開發(fā)應(yīng)用分析[J].無線互聯(lián)科技，2019，16（04）：27-28.

[5] Pina-Martins F，Paulo O S.Ncbi mass sequence downloader–large dataset downloading made easy[J].SoftwareX，2016（5）：80-83.

[6] 張立軍，曾慶東，劉志平.機械工程中常見振動的相平面分析[J].機械設(shè)計與制造，2000（1）：32-33.

[7] 姚建盛，李淑梅.Python在科學(xué)計算中的應(yīng)用[J].數(shù)字技術(shù)與應(yīng)用，2016（11）：76.