常寧，梁習(xí)鋒，黃尊地, ，尹志春，黃俞銘，鄭炯杰，莫廣杏

?

基于橫風(fēng)非定常的CRH6橫截面優(yōu)化研究

常寧1，梁習(xí)鋒2, 3，黃尊地1, 2, 3，尹志春4，黃俞銘1，鄭炯杰1，莫廣杏1

(1. 五邑大學(xué) 軌道交通學(xué)院，廣東 江門(mén) 529020； 2. 中南大學(xué) 交通運(yùn)輸工程學(xué)院，湖南 長(zhǎng)沙 410075；3. 中南大學(xué) 軌道交通安全教育部重點(diǎn)實(shí)驗(yàn)室，湖南 長(zhǎng)沙 410075；4. 中車廣東軌道交通車輛有限公司，廣東 江門(mén) 529100)

結(jié)合Lisp語(yǔ)言參數(shù)化建模和列車空氣動(dòng)力學(xué)橫風(fēng)非定常理論，對(duì)高架上運(yùn)行的CRH6型城際動(dòng)車組橫斷面進(jìn)行優(yōu)化研究。在對(duì)正交試驗(yàn)的參數(shù)化模型劃分網(wǎng)格的基礎(chǔ)上，進(jìn)行充分的橫風(fēng)定常計(jì)算，進(jìn)行非定常計(jì)算記錄監(jiān)測(cè)數(shù)據(jù)，研究列車周圍的非定常流場(chǎng)特性，對(duì)數(shù)據(jù)進(jìn)行時(shí)域和頻域特性的分析，得到優(yōu)化方案。研究結(jié)果表明：定常橫風(fēng)作用下，列車氣動(dòng)力及周圍流場(chǎng)都呈現(xiàn)明顯的非定常特性，并有周期性；非定常主振頻率與來(lái)流橫風(fēng)風(fēng)速呈線性遞增關(guān)系，在0~3 Hz范圍內(nèi)；綜合考慮傾覆力矩的最大值和幅值對(duì)列車運(yùn)行安全的影響，最終確定CRH6橫斷面優(yōu)化結(jié)果。

橫風(fēng)；非定常；橫斷面；優(yōu)化研究；傾覆力矩

和諧號(hào)CRH6型城際動(dòng)車組是為滿足中國(guó)區(qū)域經(jīng)濟(jì)快速發(fā)展和城市群崛起對(duì)城際軌道交通的需求而研制的一種新型運(yùn)輸工具。在珠三角臺(tái)風(fēng)等橫風(fēng)環(huán)境下運(yùn)行，列車的橫向運(yùn)行穩(wěn)定性至關(guān)重要[1]。開(kāi)展列車橫斷面優(yōu)化設(shè)計(jì)，可有效縮短新車體的設(shè)計(jì)周期，提高車體設(shè)計(jì)的效率，在城際動(dòng)車組的研發(fā)制造中占有重要地位。在列車的外形優(yōu)化中，SUN等[2?3]提出了結(jié)合遺傳算法或混合遺傳算法進(jìn)行優(yōu)化的方法，在一定程度上減少流場(chǎng)計(jì)算次數(shù)，縮短設(shè)計(jì)周期；李明等[4]提出了基于Isight軟件建立參數(shù)化驅(qū)動(dòng)的高速列車頭型氣動(dòng)性能自動(dòng)計(jì)算的優(yōu)化分析方法，提高了優(yōu)化計(jì)算效率?；趨?shù)化建模、正交試驗(yàn)及遺傳算法等較好的優(yōu)化組合也可以很好應(yīng)用在列車橫斷面的優(yōu)化上。張潔等[5]對(duì)強(qiáng)側(cè)風(fēng)作用下客車的氣動(dòng)外形進(jìn)行了優(yōu)化，利用CFD技術(shù)和類比法對(duì)列車橫斷面氣動(dòng)外形進(jìn)行了詳細(xì)的研究。但越來(lái)越多的研究表明，在強(qiáng)橫風(fēng)作用下，車體的氣動(dòng)力及傾覆力矩等呈現(xiàn)出明顯的非定常特性[6?12]，所以基于橫風(fēng)作用非定常特性下，對(duì)城際動(dòng)車CRH6橫斷面的參數(shù)化優(yōu)化設(shè)計(jì)具有現(xiàn)實(shí)意義。

1 計(jì)算模型

1.1 參數(shù)化建模

lisp語(yǔ)言最早為麻省理工學(xué)院為研究人工智能而開(kāi)發(fā)的，這種編程語(yǔ)言主要為括號(hào)加命令構(gòu)成的列表，在CAD軟件二次開(kāi)發(fā)中得到很好的應(yīng)用，所以既可以為設(shè)計(jì)工程師轉(zhuǎn)化為圖紙，又可以輸出通用IGES格式進(jìn)行自動(dòng)網(wǎng)格劃分。

橫截面外形的分析建立在原型橫截面的基礎(chǔ)上，由圓弧和線段組成。采用lisp語(yǔ)言能夠保證圓弧與圓弧之間的平滑連接，還可以通過(guò)多段圓弧逼近法來(lái)達(dá)到圓弧與線段之間的平滑過(guò)渡，程序編寫(xiě)完成后封裝，得到的CRH6橫截面參數(shù)化設(shè)計(jì)封裝界面如圖1所示。橫斷面中每一個(gè)參數(shù)在允許范圍內(nèi)均可改變，并保證過(guò)渡平滑符合制造工藝要求；車輛斷面在一線還是二線可以選擇；并以通用格式保存的指定位置，供網(wǎng)格劃分軟件調(diào)用。為了控制后續(xù)工作的計(jì)算量，暫取3個(gè)參數(shù)變量，即5，2和6，表示的含義分別為列車最寬處距軌面高度、車頂側(cè)墻的過(guò)渡弧半徑和底架側(cè)墻的過(guò)渡弧半徑。

圖1 橫斷面參數(shù)化設(shè)計(jì)封裝

1.2 網(wǎng)格模型

CRH6型城際動(dòng)車組是8車編組，近軌面運(yùn)行，設(shè)備艙距離軌面0.2 m，應(yīng)充分考慮道床和鋼軌對(duì)動(dòng)車組氣動(dòng)性能的影響；同時(shí)，珠三角城際線路基本為高架形式，應(yīng)建立車輛、鋼軌、道床和高架橋的耦合模型，軌面距離地面12.75 m。車輛運(yùn)行工況對(duì)橫斷面的影響是類似的，選擇二維模型可以更好的分析橫斷面的氣動(dòng)特性并進(jìn)行優(yōu)化。

圖2 橫斷面及高架橋網(wǎng)格剖分圖

將參數(shù)化模型自動(dòng)導(dǎo)入網(wǎng)格離散軟件中，生成附面層網(wǎng)格，其中第一層附面層網(wǎng)格厚度為5 mm，同時(shí)對(duì)動(dòng)車組橫斷面、鋼軌、道床、高架橋及空氣域進(jìn)行網(wǎng)格劃分，最后網(wǎng)格總數(shù)為60萬(wàn)左右，局部網(wǎng)格如圖2所示。

2 橫風(fēng)非定常算法

2.1 計(jì)算方法

最大風(fēng)力設(shè)置為珠三角地區(qū)有記錄的最大臺(tái)風(fēng)風(fēng)速，為17級(jí)風(fēng)，即來(lái)流風(fēng)速為60 m/s，經(jīng)計(jì)算得到馬赫數(shù)為0.176，小于0.3，按不可壓流處理；最小風(fēng)速取6級(jí)風(fēng)，風(fēng)速為12 m/s，計(jì)算得到雷諾數(shù)約為3*106，大于105，按湍流處理。壓力與速度耦合處理方式選用SIMPLE算法，方程離散格式采用二階迎風(fēng)格式，時(shí)間差分格式為二級(jí)隱式，穩(wěn)態(tài)計(jì)算時(shí)選用Standard-雙方程湍流模型，非定常計(jì)算時(shí)采用分離渦湍流模型。

為了保證車體橫斷面附近流場(chǎng)充分湍流，首先進(jìn)行穩(wěn)態(tài)計(jì)算，迭代次數(shù)為30 000次，從殘差的監(jiān)控圖中可以明顯的看出，前面幾千次迭代時(shí)逐漸趨于穩(wěn)定，后面會(huì)發(fā)生周期性波動(dòng)的非定常規(guī)律；穩(wěn)態(tài)的湍流計(jì)算作為非定常的初始解，經(jīng)反復(fù)試算，最終確定非定常時(shí)間步長(zhǎng)為0.01 s，每一個(gè)時(shí)間步長(zhǎng)內(nèi)進(jìn)行30次子迭代，通過(guò)監(jiān)測(cè)列車傾覆力矩氣動(dòng)參數(shù)等保證每一個(gè)時(shí)間步長(zhǎng)內(nèi)計(jì)算的收斂。先計(jì)算500個(gè)時(shí)間步長(zhǎng)，讓非定常湍流再次充分發(fā)展；然后再計(jì)算500個(gè)時(shí)間步長(zhǎng)，輸出橫向力、升力和傾覆力矩等氣動(dòng)力參數(shù)變化及附近監(jiān)測(cè)點(diǎn)的壓力波動(dòng)情況。

(a) 橫向力時(shí)程曲線；(b) 橫向力頻譜分布情況；(c) 升力時(shí)程曲線；(d) 升力頻譜分布情況；(e) 傾覆力矩時(shí)程曲線；(f) 傾覆力矩頻譜分布情況

2.2 列車非定常氣動(dòng)特性分析

圖3為CRH6原型車的各項(xiàng)氣動(dòng)力參數(shù)的時(shí)域和頻域的變化情況。從時(shí)域變化圖中可以看出橫向力、升力和傾覆力矩都呈現(xiàn)明顯的非定常特性，而且有固定的周期性；從頻域分布圖中，可以清醒的看到列車在橫風(fēng)作用下所受到的非定常氣動(dòng)力的特征頻率主要集中在低頻范圍內(nèi)，主要峰值為1.923 Hz。

表1中列出了列車非定常氣動(dòng)特性與風(fēng)速的關(guān)系，隨著來(lái)流風(fēng)速的增加，列車側(cè)向力、升力和傾覆力矩的最大值、最小值和幅值均呈現(xiàn)遞增趨勢(shì)；列車的激振頻率也隨之線性增加，范圍在0~3 Hz，與文獻(xiàn)[6]等有很好的對(duì)應(yīng)關(guān)系，且與列車系統(tǒng)本身存在的一些固有振動(dòng)頻率接近，存在強(qiáng)橫風(fēng)下列車運(yùn)行安全及引發(fā)列車脫軌傾覆的可能性。

斯特勞哈爾數(shù)(Strouhal number)是流體力學(xué)中討論物理相似與模擬時(shí)引入的相似準(zhǔn)則，考慮具有特征頻率的圓周運(yùn)動(dòng)時(shí)要使用斯特勞哈爾數(shù)進(jìn)行評(píng)判，其具體表達(dá)式為：

其中：為漩渦脫落頻率，Hz；為特征長(zhǎng)度，本文中為車體高度，m；為來(lái)流速度，即橫風(fēng)風(fēng) 速，m/s。

不同橫風(fēng)風(fēng)速下，斯特勞哈爾數(shù)的計(jì)算值如圖4(d)所示，其數(shù)值在0.1~0.2之間，這與文獻(xiàn)[9]等有很好的對(duì)應(yīng)關(guān)系。

表1 非定常氣動(dòng)特性與風(fēng)速的關(guān)系

(a) 氣動(dòng)力最值變化曲線；(b) 氣動(dòng)力幅值變化曲線；(c) 主振頻率變化曲線；(d) 斯特勞哈爾數(shù)變化曲線

綜上所示，本文得到的氣動(dòng)力數(shù)據(jù)、激振頻率及計(jì)算得到的斯特勞哈爾數(shù)等跟國(guó)內(nèi)外其他文獻(xiàn)有很好的驗(yàn)證，同時(shí)具有更好的周期性。

2.3 列車周圍非定常流場(chǎng)特性分析

列車各項(xiàng)氣動(dòng)力參數(shù)的非定常周期性波動(dòng)，必定與橫風(fēng)作用下列車周圍的流場(chǎng)變化有關(guān)系。圖5列出了CRH6橫斷面周圍流場(chǎng)的監(jiān)測(cè)點(diǎn)以及不同時(shí)刻下流場(chǎng)的變化情況，主要反映其渦量云圖的變化。可以看出，在定常橫風(fēng)作用下，列車背風(fēng)側(cè)流動(dòng)非常復(fù)雜，產(chǎn)生了許多不同空間尺度的渦。從記錄有效數(shù)據(jù)開(kāi)始，在0.5 s內(nèi)，CRH6橫斷面周圍的流場(chǎng)也發(fā)生周期性的變化，脫落渦不斷的從車體、鋼軌、道床和高架橋前流場(chǎng)分離處產(chǎn)生，并不斷向后運(yùn)動(dòng)，同時(shí)伴隨著合并和脫落。渦的運(yùn)動(dòng)規(guī)律具有隨機(jī)性，但也有周期性，這種非定常的流動(dòng)規(guī)律也是列車非定常氣動(dòng)力產(chǎn)生的原因。

(a) 監(jiān)測(cè)點(diǎn)布置圖；(b) 0 s；(c) 0.2 s；(d) 0. 4 s

從記錄有效數(shù)據(jù)開(kāi)始，對(duì)各測(cè)點(diǎn)的靜壓力隨時(shí)間的變化歷程進(jìn)行了記錄和存儲(chǔ)。受本文篇幅限制，僅列出離列車壁面較近的1和2測(cè)點(diǎn)，以及在二線外側(cè)的3和4測(cè)點(diǎn)。圖6顯示了4個(gè)測(cè)點(diǎn)的靜壓時(shí)程曲線，可以看出，4個(gè)測(cè)點(diǎn)的靜壓力不是靜止不變的，而是呈現(xiàn)了典型非定常特性，具有周期性的波動(dòng)規(guī)律。同時(shí)對(duì)其進(jìn)行傅里葉頻譜分析，得到其振動(dòng)主頻約為2 Hz，與圖3中列車所受的氣動(dòng)力參數(shù)的主頻具有很好的對(duì)應(yīng)關(guān)系。

3 優(yōu)化結(jié)果分析

本文3個(gè)參數(shù)變量，即5，2和6，表示的含義分別為列車最寬處距軌面高度、車頂側(cè)墻的過(guò)渡弧半徑和底架側(cè)墻的過(guò)渡弧半徑。根據(jù)初步的設(shè)計(jì)需要，每個(gè)因素選擇3個(gè)水平，1個(gè)三因素三水平的實(shí)驗(yàn)，按全面實(shí)驗(yàn)要求，須進(jìn)行3^3=27種組合的實(shí)驗(yàn)；若按L9(3^4)正交表安排實(shí)驗(yàn)，只需作9次。三因素三水平構(gòu)成表2的前9個(gè)試驗(yàn)組，都是正交分布的，沒(méi)有2個(gè)水平是交叉的；但是構(gòu)建三元二次回歸方程求解10個(gè)未知參數(shù)，需要10個(gè)封閉的方程組，所以表2中又加上了原型車的數(shù)據(jù)。完成10個(gè)完整參數(shù)組合的流場(chǎng)計(jì)算后，得到表2中氣動(dòng)力數(shù)據(jù)。

(a) 監(jiān)測(cè)點(diǎn)1靜壓力時(shí)程曲線；(b) 監(jiān)測(cè)點(diǎn)2靜壓力時(shí)程曲線；(c) 監(jiān)測(cè)點(diǎn)3靜壓力時(shí)程曲線；(d) 監(jiān)測(cè)點(diǎn)4靜壓力時(shí)程曲線

表2 計(jì)算結(jié)果統(tǒng)計(jì)表

側(cè)向力和升力對(duì)最后的傾覆力矩都有貢獻(xiàn)，所以在后期的尋優(yōu)中，主要分析傾覆力矩的數(shù)值。分別將傾覆力矩最大值和幅值作為目標(biāo)函數(shù)求極值，得到的優(yōu)化結(jié)果見(jiàn)表3。

傾覆力矩絕對(duì)值最大值的擬合函數(shù)為：

傾覆力矩絕對(duì)值幅值的擬合函數(shù)為：

表3 CRH6橫斷面優(yōu)化結(jié)果

從擬合公式可知，傾覆力矩的最大值與5和6的關(guān)聯(lián)系數(shù)的絕對(duì)值分別為19.428 8和11.671 1，同時(shí)幅值與5和6的關(guān)聯(lián)系數(shù)的絕對(duì)值分別為2.573 6和3.972，相比較2，5和6對(duì)列車傾覆力矩的貢獻(xiàn)值較高，即列車最寬處距軌面高度和底架側(cè)墻的過(guò)渡弧半徑對(duì)列車傾覆穩(wěn)定性的影響因子較大。經(jīng)尋優(yōu)后得到的結(jié)果中，傾覆力矩最大值最小和幅值最小的結(jié)果分別為方案1和2所示，在傾覆力矩最大值優(yōu)化結(jié)果同在5%以內(nèi)的情況下，2中傾覆力矩的波動(dòng)幅值優(yōu)化減小了21.4%，所以最終CRH6橫斷面優(yōu)化結(jié)果5，2和6分別為1 425，400和800 mm。

5 結(jié)論

1) 定常橫風(fēng)作用下，高架上運(yùn)行的列車橫向力、升力、傾覆力矩以及周圍流場(chǎng)都呈現(xiàn)明顯的非定常特性，有周期性；非定常氣動(dòng)力及周圍流場(chǎng)的特征頻率主要集中在低頻范圍內(nèi)，主振頻率在0~3 Hz范圍內(nèi)。

2) 列車最寬處距軌面高度、底架與側(cè)墻的過(guò)渡弧半徑對(duì)列車傾覆穩(wěn)定性的影響較大。

3) 相對(duì)比定常分析中只考慮傾覆力矩最大值，非定常分析要同時(shí)考慮傾覆力矩最大值和幅值對(duì)列車運(yùn)行安全的影響。在傾覆力矩最大值優(yōu)化結(jié)果同在5%以內(nèi)的情況下，方案2中傾覆力矩的波動(dòng)幅值優(yōu)化減小了21.4%，所以最終CRH6橫斷面優(yōu)化結(jié)果5，2和6分別為1425，400和800 mm。

[1] 田紅旗. 列車空氣動(dòng)力學(xué)[M]. 北京: 中國(guó)鐵道出版社, 2007: 160?161. TIAN Hongqi. Train aerodynamics[M]. Beijing: China Railway Press, 2007: 160?161.

[2] SUN Zhenxu, SONG Jingjing, AN Yiran. Optimization of the head shape of the CRH3 high speed train[J]. China Technological Sciences, 2010, 53(12): 3356?3364.

[3] 姚拴寶, 郭迪龍, 楊國(guó)偉. 基于混合遺傳算法的高速列車截面變化率優(yōu)化設(shè)計(jì)[J]. 科學(xué)技術(shù)與工程, 2013, 13(28): 8349?8355. YAO Shuanbao, GUO Dilong, YANG Guowei. Optimal cross-sectional area distribution of a high-speed train nose based on a hybrid genetic algorithm[J]. Science Technology and Engineering, 2013, 13(28):8349?8355.

[4] 李明, 李明高, 李國(guó)清, 等. 參數(shù)化驅(qū)動(dòng)的高速列車頭型氣動(dòng)外形優(yōu)化設(shè)計(jì)[J]. 鐵道學(xué)報(bào), 2013, 35(11): 14? 20.LI Ming, LI Minggao, LI Guoqing, et al. Optimized design of parametric-driven aerodynamic shape of high-speed EMU head-type[J]. Journal of the China Railway Society, 2013, 35(11): 14?20.

[5] 張潔, 梁習(xí)鋒, 劉堂紅, 等. 強(qiáng)側(cè)風(fēng)作用下客車車體氣動(dòng)外形優(yōu)化[J]. 中南大學(xué)學(xué)報(bào)(自然科學(xué)版), 2011, 42(11): 3578?3584. ZHANG Jie, LIANG Xifeng, LIU Tanghong, et al. Optimization research on aerodynamic shape of passenger car body with strong crosswind[J]. Journal of Central South University (Science and Technology), 2011, 42(11): 3578?3584.

[6] 劉加利, 張繼業(yè), 張衛(wèi)華. 基于大渦模擬的橫風(fēng)下高速列車非定常氣動(dòng)載荷特性[J]. 鐵道學(xué)報(bào), 2013, 35(6): 13?21. LIU Jiali, ZHANG Jiye, ZHANG Weihua. Study on characteristics of unsteady aerodynamic loads of a high-speed train under crosswinds by large eddy simulation[J]. Journal of the China Railway Society, 2013, 35(6): 13?21.

[7] 郗艷紅, 毛軍, 高亮, 等. 基于分離渦方法的高速列車橫風(fēng)非定常氣動(dòng)特性[J]. 中南大學(xué)學(xué)報(bào)(自然科學(xué)版), 2015, 46(3): 1129?1139.XI Yanhong, MAO Jun, GAO Liang, et al. Aerodynamic force/moment for high-speed train in crosswind field based on DES[J]. Journal of Central South University (Science and Technology), 2015, 46(3): 1129?1139.

[8] 苗秀娟, 高廣軍. 基于DES的車輛橫風(fēng)氣動(dòng)性能模擬[J]. 中南大學(xué)學(xué)報(bào)(自然科學(xué)版), 2012, 43(7): 2855? 2860. MIAO Xiujuan, GAO Guangjun. Aerodynamic performance of train under cross-wind based on DES[J]. Journal of Central South University (Science and Technology), 2012, 43(7): 2855?2860.

[9] Hassan Hemida, Sinisa Krajnovic, Lars Davidson. Large-eddy simulation of the flow around a simplified high speed train under the influence of a cross-wind[C]// AIAA. 17thAIAA Computational Fluid Dynamics Conferece 2005. Toronto, Ontario, Canada: AIAA, 2005: 5354?5370.

[10] Hassan Hemida, Sinisa Krajnovic. Numerical study of the unsteady flow structures around train-shaped body subjected to side winds[DB/OL]. ECCOMAS CFD 2006. http://proceedings.fyper.Com/eccomascfd2006/documents/218.pdf.

[11] Ben Diedrichs. Aerodynamic calculations of crosswind stability of a high-speed train using control columes of arbitrary polyhedral shape[C]// BBAA VI International Colloquium on: Bluff Bodies Aerodynamic & Applications Milano, Italy, 2008: 20?24.

[12] Asress Mulugeta Biadgo, Aleksandar Simonovic, Jelena Svorcan,et al. Aerodynamic characteristics of high speed train under turbulent cross winds: A numerical investigation using Unsteady-RANS method[J]. Article in FME Transactions, 2014, 42(1): 10?18.

Optimization research on cross-section of CRH6by unsteady aerodynamic characteristics under crosswind

CHANG Ning1, LIANG Xifeng2, 3, HUANG Zundi1, 2, 3, YIN Zhichun4, HUANG Yuming1, ZHENG Jiongjie1, MO Guangxing1

(1. School of Rail Transportation, Wuyi University, Jiangmen 529020, China; 2. School of Traffic & Transportation Engineering, Central South University, Changsha 410075, China; 3. Key Laboratory of Traffic Safety on Track, Ministry of Education, Central South University, Changsha 410075, China; 4. CRRC Guangdong Co., Ltd, Jiangmen 529100, China)

Combined parametric modeling based on Lisp language and unsteady aerodynamic theory of train under crosswinds, cross section of CRH6intercity EMU which is common running on an elevated bridge is optimized. On the basis of orthogonal test, parametric models are meshed; firstly steady calculation under crosswind is fully simulated, and then unsteady calculation is carried out; monitoring date is recorded to research the unsteady flow field characteristics of the train around; finally, optimization results are got by the data of time domain and frequency domain characteristics. The results show that: under steady cross-wind effect, train aerodynamics’ force and flow field around are all showed significantly unsteady characteristics, and cyclical; the relationship between unsteady main vibration frequency and flow cross-wind wind speed is linearly increasing, whose range is between 0 to 3 Hz; considering the operating safety effect of the train by the maximum and the magnitude of overturning moment, the final cross-section optimization results are determined.

crosswinds; unsteady; cross section; optimization research; overturning moment

10.19713/j.cnki.43?1423/u.2018.11.002

U270.2

A

1672 ? 7029(2018)11 ? 2730 ? 08

2017?09?30

高速鐵路基礎(chǔ)研究聯(lián)合基金重點(diǎn)資助項(xiàng)目(U1334205)；廣東省青年創(chuàng)新人才類資助項(xiàng)目(2016KQNCX172)；廣東大學(xué)生科技創(chuàng)新培育專項(xiàng)資金(“攀登計(jì)劃”專項(xiàng)資金)立項(xiàng)項(xiàng)目(pdjh2017b0520)

黃尊地(1987?)，男，山東嘉祥人，講師，博士，從事軌道交通空氣動(dòng)力學(xué)研究；E?mail：wyuhzd@163.com

(編輯 蔣學(xué)東)