朱 博, 湯更生

(中國空氣動力研究與發(fā)展中心 設(shè)備設(shè)計(jì)及測試技術(shù)研究所, 四川 綿陽 621000)

聲學(xué)風(fēng)洞流場低湍流度及頻譜測量研究

朱 博*, 湯更生

(中國空氣動力研究與發(fā)展中心 設(shè)備設(shè)計(jì)及測試技術(shù)研究所, 四川 綿陽 621000)

采用熱線風(fēng)速儀對5.5m×4m低湍流航空聲學(xué)風(fēng)洞閉口試驗(yàn)段低湍流度流場進(jìn)行測量，根據(jù)對干擾信號的分析，提出了高通慣性衰減濾波方法，并與一般數(shù)據(jù)處理方法進(jìn)行了比較，給出了流場測量方案、方法和結(jié)果。采用功率譜方法和斯特羅哈數(shù)方法分析脈動速度信號中的干擾噪聲，發(fā)現(xiàn)40Hz～10kHz頻譜范圍內(nèi)同時存在電磁干擾噪聲和支架干擾噪聲。比較分析了0.5Hz～5kHz帶通濾波方法、電磁噪聲解耦方法和0.5Hz高通慣性衰減濾波方法對干擾信號的濾除效果，采用0.5Hz高通慣性衰減濾波方法獲得了流場低湍流度數(shù)據(jù)，流場速度30～100m/s的湍流度平均值小于0.05%。實(shí)驗(yàn)結(jié)果表明，高通慣性衰減濾波方法可以有效控制干擾信號對測量結(jié)果的影響程度，為低湍流度流場信號處理提供了一種方法。

低速風(fēng)洞；聲學(xué)風(fēng)洞；低湍流度；信號處理；熱線風(fēng)速儀

0 引 言

5.5m×4m低湍流航空聲學(xué)風(fēng)洞(簡稱聲學(xué)風(fēng)洞)屬于我國重大科技基礎(chǔ)設(shè)施建設(shè)項(xiàng)目，是一座新建的連續(xù)式單回流低速風(fēng)洞，主要用于聲學(xué)試驗(yàn)、低湍流度試驗(yàn)和低速試驗(yàn)。湍流度對飛行器氣動力特性影響較大，因此，湍流度是風(fēng)洞流場校測的重要指標(biāo)，在飛行器模型的風(fēng)洞試驗(yàn)中，需要對風(fēng)洞流場湍流度進(jìn)行準(zhǔn)確測量，并對其影響進(jìn)行修正，以提高風(fēng)洞試驗(yàn)數(shù)據(jù)精準(zhǔn)度。如果湍流度測量不準(zhǔn)確，風(fēng)洞試驗(yàn)得到的模型氣動系數(shù)也不精確，就會在很大程度上影響飛行器氣動性能的預(yù)測精度[1-3]。目前，測量流場湍流度的最理想儀器是熱線風(fēng)速儀(簡稱熱線)，但是由于測量目標(biāo)為動態(tài)信號，導(dǎo)致測量結(jié)果容易受到干擾，特別對于低湍流度測量而言，更容易受到電磁干擾[3-4]和支架干擾[5]。因此，在聲學(xué)風(fēng)洞低湍流度測量中，需要對測量數(shù)據(jù)進(jìn)行信號處理。

國內(nèi)外文獻(xiàn)對低速流場低湍流度測量數(shù)據(jù)進(jìn)行信號處理的方法主要有帶通濾波[6-9]和干擾信號解耦[10-11]的方法，這些風(fēng)洞的尺寸、阻尼網(wǎng)層數(shù)、收縮比、測量位置、試驗(yàn)風(fēng)速、湍流度和采用的信號處理方法比較如表1所示。

可見，表1中各風(fēng)洞采用的帶通濾波參數(shù)不盡相同。其中，高通濾波截止頻率取值范圍為0.1Hz～10Hz，低通濾波截止頻率取值范圍為1kHz～10kHz。濾波參數(shù)如何選取，不同的濾波參數(shù)和數(shù)據(jù)處理方法對低湍流度測量結(jié)果的影響程度如何，是測量低湍流度需要分析的問題。為此，本文在測量聲學(xué)風(fēng)洞低湍流度時，對熱線信號頻譜中的干擾信號進(jìn)行初步分析，探討識別電磁干擾和支架干擾的方法，并對幾種濾波參數(shù)的數(shù)據(jù)處理結(jié)果和文獻(xiàn)[10]的電磁噪聲解耦處理結(jié)果進(jìn)行比較，根據(jù)脈動速度頻譜能量的指數(shù)衰減特征，提出基于高通慣性衰減的濾波方法。

表1 低速風(fēng)洞低湍流度信號處理方法和效果Table 1 Low turbulence intensity signal process methods and results in low speed wind tunnel

1 設(shè)備條件及熱線儀

1.1 設(shè)備條件

本文主要試驗(yàn)在聲學(xué)風(fēng)洞的閉口試驗(yàn)段和開口試驗(yàn)段上進(jìn)行，風(fēng)洞穩(wěn)定段設(shè)置6層阻尼網(wǎng)，部分關(guān)鍵部段設(shè)置了消聲部件，開口試驗(yàn)段設(shè)置了消聲室。閉口試驗(yàn)段流場湍流度設(shè)計(jì)指標(biāo)為0.05%，開口試驗(yàn)段流場湍流度設(shè)計(jì)指標(biāo)為0.2%。

1.2 熱線儀

熱線儀采用丹麥丹迪公司的STREAMLINE系統(tǒng)，測量探頭主要有55P11一維探頭和55P61二維探頭，測點(diǎn)布置于風(fēng)洞試驗(yàn)段的中心，系統(tǒng)配置的采集卡分辨率為16位。數(shù)據(jù)采樣頻率為20kHz，采樣時間為10s。

2 低速低湍流度流場的干擾信號分析

判斷流場非低湍流度測量結(jié)果是否受到干擾，可觀察流場的速度頻譜是否有干擾信號，如圖1所示(速度40m/s)。由于充分發(fā)展的湍流能譜按照一定的標(biāo)度率呈現(xiàn)指數(shù)衰減[12-13]，因此，無明顯干擾的速度幅值隨著頻率的增加呈現(xiàn)按指數(shù)衰減的趨勢，而且其頻譜無明顯的尖峰信號。但是，對于低湍流度測量而言，由于測量目標(biāo)的量值較小(湍流度0.05%)，已經(jīng)接近測量系統(tǒng)的自生噪聲量值(丹迪熱線儀16位數(shù)采系統(tǒng)對流場速度40m/s的熱線信號僅有0.04%的分辨率，而且常規(guī)數(shù)采有1～3個最小有效位的電噪聲)，這時儀器的本底電噪聲和測量支架振動干擾對測量結(jié)果的影響是不容忽視的，因此需要進(jìn)一步識別和濾除這類干擾量。

圖1 低通濾波10kHz，湍流度為0.053%的流場速度頻譜

Fig.1 Flow velocity spectrum with low pass of 10kHz, turbulence intensity of 0.053%

圖2 低通濾波10kHz，湍流度為0.053%的流場脈動速度頻譜

Fig.2 Flow velocity fluctuation spectrum with low pass of 10kHz, turbulence intensity of 0.053%

圖3 低通濾波10kHz，40m/s的流場脈動速度平均頻譜

Fig.3 Flow velocity fluctuation mean spectrum at flow speed of 40m/s with low pass of 10kHz

從湍流度值的計(jì)算方法來看：

(1)

圖4 不同流場速度的脈動速度功率譜

Fig.4 Flow velocity fluctuation power spectrum at different flow speed

圖5 頻帶40Hz～1kHz的流場脈動速度功率譜

Fig.5 Flow velocity fluctuation power spectrum at frequencies of 40Hz～1kHz

圖6 頻帶500Hz～10kHz的流場脈動速度功率譜

Fig.6 Flow velocity fluctuation power spectrum at frequencies of 500Hz～10kHz

根據(jù)斯特羅哈數(shù)常數(shù)的定義：

(2)

在相同風(fēng)速v的條件下，改變漩渦發(fā)生體的尺寸l就可以產(chǎn)生漩渦頻率f和幅值的變化。據(jù)此，調(diào)整熱線探針支桿尺寸縮短三分之一，測量獲得圖7，可見相同風(fēng)速下信號尖峰頻率稍微變大，90～100m/s風(fēng)速的尖峰幅值有明顯減小。當(dāng)流場風(fēng)速達(dá)到110m/s時，采用長支桿測量獲得的脈動速度，比采用短支桿測量獲得的脈動速度在整個頻譜范圍內(nèi)的幅值都明顯偏大(見圖8)，可見特定頻率的干擾信號可能對整個頻域的湍流能譜標(biāo)度率都有影響。由此，可判斷700Hz～2kHz頻帶內(nèi)隨流速變化的尖峰信號為支架干擾，不是風(fēng)洞試驗(yàn)段流場固有的脈動速度信號，應(yīng)該予以濾除。

圖7 短支桿測量的流場脈動速度功率譜

Fig.7 Flow velocity fluctuation power spectrum with short probe support

圖8 長支桿和短支桿測量的110m/s流場脈動速度頻譜

Fig.8 Flow velocity fluctuation spectrum with short or long probe support at speed of 110m/s

3 低速低湍流度流場的信號處理方法

3.1 帶通濾波方法

帶通濾波方法是低速低湍流度流場的干擾信號處理的常用方法[6-9]，其對低頻干擾信號采用高通濾波，同時對高頻干擾信號采用低通濾波。對于大型風(fēng)洞而言，低頻壓力脈動對湍流度的貢獻(xiàn)是不容忽視的，但是，低頻脈動對轉(zhuǎn)捩的影響還有待研究[6]，因此，測量湍流度時通常對低頻信號進(jìn)行了高通濾波。低通截止頻率通常根據(jù)測量現(xiàn)場的高頻干擾信號頻率確定。

對脈動信號U′(n)進(jìn)行頻域帶通濾波[14]可以表示為：

(3)

(4)

式中：X(k)為U′(n)的頻譜幅值，Δf為頻譜的頻率分辨率，fl為通帶的最低頻率，fh為通帶的最高頻率。

3.2 電磁噪聲解耦方法

電磁噪聲解耦方法[10]是根據(jù)電磁噪聲信號與流場速度無關(guān)的假設(shè)，對干擾信號進(jìn)行解耦的方法，其解析表達(dá)式為：

(5)

3.3 高通慣性衰減濾波方法

若對圖2進(jìn)行0.5Hz～5kHz帶通濾波和尖峰濾波，效果如圖9所示，濾波后，在5kHz～10kHz的范圍內(nèi)不僅干擾信號濾除了，而且有效信號的衰減也很大；在40Hz～5kHz的范圍內(nèi)干擾信號不能有效濾除。因此，干擾信號的影響程度并未得到準(zhǔn)確控制。

圖9 帶通濾波前后的脈動速度譜比較

Fig.9 Flow velocity fluctuation spectrum with or without band pass filtering

文獻(xiàn)[3]指出，為準(zhǔn)確控制干擾信號的影響程度，濾波頻率的選擇不能太大也不能太小，應(yīng)該是最小湍流脈動對應(yīng)的Kolmogorov頻率。實(shí)際觀察聲學(xué)引導(dǎo)風(fēng)洞[15]蜂窩器后較大湍流度的流場頻譜(見圖10，湍流度10%，平均速度17m/s，頻譜平均數(shù)25)，也可見在干擾信號不影響有效信號頻譜的情況下，有效信號的頻譜能量是呈指數(shù)衰減的。據(jù)此分析，圖7所示信號的40Hz～5kHz頻譜在無干擾的情況下也應(yīng)該呈指數(shù)衰減。因此，對干擾信號進(jìn)行基于0.5Hz高通慣性衰減的帶通濾波(效果見圖11)，可以控制干擾信號對測量結(jié)果的影響程度。

圖10 聲學(xué)引導(dǎo)風(fēng)洞蜂窩器后湍流度為10%的流場脈動速度頻譜

Fig.10 Flow velocity fluctuation spectrum with turbulence intensity of 10% behind honeycomb in pilot aeroacoustic wind tunnel

Fig.11 Flow velocity fluctuation spectrum with or without 0.5Hz high pass inertia attenuation filtering

高通慣性衰減濾波方法，指對低頻干擾信號采用高通濾波，同時對高頻干擾信號采用阻帶慣性衰減和尖峰中值濾波的方法，其阻帶衰減斜率為無干擾尖峰的脈動速度信號指數(shù)衰減斜率。

尖峰中值濾波在頻域可以表示為：

(6)

(7)

式中：X(k)為干擾信號頻率ft的幅值，m?Δf為干擾信號頻率的阻尼半寬。

阻帶慣性衰減可采用巴特沃斯低通濾波器進(jìn)行設(shè)計(jì)，使得圖2的40Hz～10kHz頻譜沿1Hz～40Hz的斜率進(jìn)行慣性衰減。本文的巴特沃斯低通濾波器設(shè)計(jì)參數(shù)為通帶角頻率1，阻帶角頻率0.31，通帶最大衰減1dB，阻帶最小衰減50dB。

4 測量結(jié)果和分析

4.1 不同數(shù)據(jù)處理方法的比較

圖12是采用不同方法對閉口試驗(yàn)段的同一組數(shù)據(jù)進(jìn)行處理獲得的湍流度結(jié)果比較。由于存在測量支架振動，而且隨著風(fēng)速的增加，振動對流場脈動速度的影響加大(見圖6和7)，因此，電磁噪聲解耦方法獲得的湍流度結(jié)果隨著風(fēng)速的增加逐漸接近原始值。帶通濾波方法雖然對低于0.5Hz和高于5kHz的信號進(jìn)行了濾波，但是對高于5kHz的高頻干擾信號未予有效濾除，因此，湍流度測量結(jié)果較大。高通慣性衰減濾波方法對低于0.5Hz的信號進(jìn)行了濾波，并根據(jù)前文的電磁干擾分析和支架干擾分析的結(jié)果，對高頻干擾信號進(jìn)行了慣性衰減濾波，獲得的測量結(jié)果比其他兩種方法的測量結(jié)果小。

圖12 閉口試驗(yàn)段湍流度不同數(shù)據(jù)處理方法的比較

Fig.12 A comparison between results with different kinds of data process method in closed test section

圖13是采用不同方法對開口試驗(yàn)段的同一組數(shù)據(jù)進(jìn)行處理獲得的湍流度結(jié)果比較。可見，由于開口試驗(yàn)段低頻壓力脈動較大，0.5Hz高通濾波對試驗(yàn)數(shù)據(jù)的影響較大，高頻干擾信號對測量結(jié)果的影響相對較小。因此，電磁噪聲解耦方法獲得的湍流度結(jié)果與原始值接近，0.5Hz～5kHz帶通濾波與0.5Hz高通慣性衰減濾波結(jié)果接近。

Fig.13 Comparison between results with different kinds of data process method in open test section

4.2 聲學(xué)風(fēng)洞閉口試驗(yàn)段流場湍流度分析

圖14是采用0.5Hz高通慣性衰減濾波方法獲得的聲學(xué)風(fēng)洞閉口試驗(yàn)段流場湍流度分布，圖15是湍流度平均值。測量結(jié)果表明，閉口試驗(yàn)段流場湍流度值隨著流場速度的增加而呈上升趨勢(文獻(xiàn)[9]的數(shù)據(jù)也為上升趨勢)，湍流度值在流場速度100m/s以下優(yōu)于設(shè)計(jì)指標(biāo)0.05%。當(dāng)流場速度達(dá)到110m/s時，由于支架干擾信號已經(jīng)對整個頻域的有效信號的頻譜幅值有較大影響(見圖8)，這時采用高通慣性衰減的濾波方法也不能有效去除干擾信號，因此采用現(xiàn)有測量方法在110m/s流場速度獲得的湍流度值不宜作為流場考核指標(biāo)。

圖14 閉口試驗(yàn)段流場湍流度分布

Fig.14 Turbulence intensity distribution in flow field of closed test section

圖15 閉口試驗(yàn)段流場湍流度平均值Fig.15 Turbulence mean intensity in flow field of closed test section

5 結(jié) 論

本文給出了5.5m×4m低湍流航空聲學(xué)風(fēng)洞流場湍流度測量方法、分析方法和測量結(jié)果。

采用功率譜方法和斯特羅哈數(shù)方法分析了脈動速度信號中的干擾噪聲，發(fā)現(xiàn)40Hz～10kHz頻譜范圍內(nèi)同時存在電磁干擾噪聲和支架干擾噪聲。

提出了高通慣性衰減濾波方法，比較分析了0.5Hz～5kHz帶通濾波方法、電磁噪聲解耦方法和0.5Hz高通慣性衰減濾波方法對干擾信號的濾除效果。比較情況表明，高通慣性衰減濾波方法在一定風(fēng)速范圍內(nèi)，可以有效控制干擾信號對測量結(jié)果的影響程度。

采用0.5Hz高通慣性衰減濾波方法獲得了5.5m×4m低湍流航空聲學(xué)風(fēng)洞閉口試驗(yàn)段流場低湍流度數(shù)據(jù)，流場速度30～100m/s的湍流度平均值優(yōu)于設(shè)計(jì)指標(biāo)0.05%。

當(dāng)流場速度達(dá)到110m/s時，由于支架干擾信號已經(jīng)對整個頻域的有效信號的頻譜幅值有較大影響，這時采用0.5Hz高通慣性衰減的濾波方法也不能有效去除干擾信號，因此需要在測量支架設(shè)計(jì)、系統(tǒng)電磁兼容設(shè)計(jì)中進(jìn)一步做好抗干擾措施，以提高流場湍流度的測量精度。

[1] 何克敏, 白存儒, 郭渠渝, 等. 較低湍流度范圍湍流度對風(fēng)洞實(shí)驗(yàn)結(jié)果的影響[J]. 流體力學(xué)實(shí)驗(yàn)與測量, 1997, 11(3): 11-17.

He Kemin, Bai Cunru, Guo Quyu, et al. The effect of turbulence on wind tunnel results in the range of low turbulence[J]. Experiments and Measurements In Fluid Mechanics, 1997, 11(3): 11-17.

[2] 李峰, 白存儒, 周偉, 等. 湍流度對飛行器模型大迎角氣動特性影響的初步研究[J]. 實(shí)驗(yàn)流體力學(xué), 2006, 20(3): 45-52.

Li Feng, Bai Cunru, Zhou Wei, et al. Primal research of the effect of flow turbulence on aerodynamic characteristics of a aircraft model at high angles of attack[J]. Journal of Experiments in Fluid Mechanics, 2006, 20(3): 45-52.

[3] 米建春, 馮寶平. 平面射流沿軸線的特征尺度及其對測量信號過濾程度的依賴[J]. 物理學(xué)報(bào), 2010, 59(7): 4748-4755.

Mi Jianchun, Feng Baoping. Centerline characteristic scales of a turbulent plane jet and their dependence on filtration of measured signals[J]. Acta Physica Sinica, 2010, 59(7): 4748-4755.

[4] 侯志勇, 王連澤, 周建和, 等. 低(變) 湍流度風(fēng)洞設(shè)計(jì)再探討[J]. 實(shí)驗(yàn)流體力學(xué), 2011, 25(1): 92-96.

Hou Zhiyong, Wang Lianze, Zhou Jianhe, et al. The further research on the design of low (varying) turbulence wind tunnel [J]. Journal of Experiments in Fluid Mechanics, 2011, 25(1): 92-96.

[5] Harald Quix, Jürgen Quest, Christian Brzek. Hot-wire measurements in cryogenic environment[R]. AIAA-2011-880.

[6] Saric W S, Takagi S, Mousseux M. The ASU unsteady wind tunnel and fundamental requirements for freestream turbulence measurements[R]. AIAA-1988-0053.

[7] Mueller T J, Scharpf D F, Batill S M, et al. The design of a subsonic low-noise, low-turbulence wind tunnel for acoustic measurements[R]. AIAA-1992-3883.

[8] Neuhart Dan H, McGinley Catherine B. Free-stream turbulence intensity in the Langley 14- by 22-foot subsonic tunnel[R]. NASA/TP-2004-213247.

[9]Premi Amandeep, Maughmer Mark D, Brophy Christopher. Flow-quality measurements and qualification of the Pennsylvania State University low-speed, low-turbulence wind tunnel[R]. AIAA-2012-1214.

[10] Stainback P Calvin, Owen F Kevin. Dynamic flow quality measurements in the Langley low turbulence pressure tunnel [R]. AIAA-1984-0621.

[11] Kohama Y, Kobayashi R, Ito H. Tohoku university low-turbulence wind tunnel[R]. AIAA-1992-3913.

[12] 李士心, 姜楠, 舒瑋. 自由湍射流多尺度湍渦結(jié)構(gòu)標(biāo)度律的實(shí)驗(yàn)研究[J]. 實(shí)驗(yàn)力學(xué), 1999, 14(4): 409-413.

Li Shixin, Jiang Nan, Shu Wei. Experimental Research on the scaling law for the multi-scale eddy structure in free turbulent jet[J]. Journal of Experimental Mechanics, 1999, 14(4): 409-413.

[13] 畢衛(wèi)濤, 陳凱, 魏慶鼎. 射流和圓柱尾流湍流層次結(jié)構(gòu)的實(shí)驗(yàn)研究[J]. 流體力學(xué)實(shí)驗(yàn)與測量, 2001, 15(3): 13-18.

Bi Weitao, Chen Kai, Wei Qingding. Experimental study on the turbulent hierarchical structures of jet and cylindrical wake flow[J]. Experiments and Measurements in Fluid Mechanics, 2001, 15(3): 13-18.

[14] 朱博, 王元興, 余永生. 風(fēng)洞模型-支撐系統(tǒng)渦激振動測量與分析[J]. 實(shí)驗(yàn)流體力學(xué), 2014, 28(5): 59-64.

Zhu Bo, Wang Yuanxing, Yu Yongsheng. Vortex-induced vibration measurement and analysis of model-sting system in wind tunnel[J]. Journal of Experiments in Fluid Mechanics, 2014, 28(5): 59-64.

[15] 朱博, 湯更生, 顧光武, 等. 低速航空聲學(xué)風(fēng)洞背景噪聲測試技術(shù)研究[J]. 實(shí)驗(yàn)流體力學(xué), 2009, 23(4): 56-64.

Zhu Bo, Tang Gengsheng, Gu Guangwu, et al. An investigation of background noise measurement technology in aeroacoustic wind tunnel[J]. Journal of Experiments in Fluid Mechanics, 2009, 23(4): 56-64.

(編輯：楊 娟)

Low turbulence intensity and spectrum measurement research in aeroacoustic wind tunnel

Zhu Bo*, Tang Gengsheng

(Facility Design and Instrumentation Institute, China Aerodynamics Research and Development Center, Mianyang Sichuan 621000, China)

The 5.5m×4m low turbulence aeroacoustic wind tunnel is a low speed single-return wind tunnel, capable of performing aeroacoustic, low turbulence and general low speed tests. The flow of this tunnel presents the important quality of low turbulence intensity, which needs to be measured accurately. The band pass filter and signal noise separation are the two main methods to process the low turbulence intensity data. However, the parameters of the two methods are different, and thus there is no direct comparison between their results. In this paper the low turbulence intensity of the flow field is measured in the closed test section of the 5.5m×4m low turbulence aeroacoustic wind tunnel with the hot wire anemometry. A 0.5Hz high pass inertia attenuation filter is adopted to process the data according to the noise analysis, and the results are compared with other two methods, that the 0.5Hz～5kHz band pass filter and the signal noise separation. Noises in the flow field fluctuation are analyzed based on the power spectrum and the Strouhal numbers. The electromagnetic noise and support disturbing noise are found in spectrum domain of 40Hz～10kHz. Using a 0.5Hz high pass inertia attenuation filter, the flow field turbulence intensity is measured and found to be less than 0.05% at the speed of 30～100m/s. Experimental results show that, the high pass inertia attenuation filter method can be used to control the disturbance to the results induced by noises, which is a practical method for processing the low turbulence intensity digital signal.

low speed wind tunnel;aeroacoustic wind tunnel;low turbulence intensity;signal process;hot wire anemometer

1672-9897(2015)04-0058-07

10.11729/syltlx20140118

2014-10-19;

2015-01-09

ZhuB,TangGS.Lowturbulenceintensityandspectrummeasurementresearchinaeroacousticwindtunnel.JournalofExperimentsinFluidMechanics, 2015, 29(4): 58-64. 朱 博, 湯更生. 聲學(xué)風(fēng)洞流場低湍流度及頻譜測量研究. 實(shí)驗(yàn)流體力學(xué), 2015, 29(4): 58-64.

V211.7

A

朱 博(1973-)，男，廣西百色人，碩士，高級工程師。研究方向：風(fēng)洞測試、聲學(xué)測試和熱線技術(shù)。通信地址：四川省綿陽市二環(huán)路南段6號14信箱404分箱 (621000)。E-mail：bobjou@139.com

*通信作者 E-mail: bobjou@139.com