翁 順, 左 越, 朱宏平, 陳 波, 趙會(huì)賢, 田 煒, 顏永逸
(1.華中科技大學(xué) 土木工程與力學(xué)學(xué)院,控制結(jié)構(gòu)湖北省重點(diǎn)實(shí)驗(yàn)室,武漢 430074; 2.武漢理工大學(xué) 土木工程與建筑學(xué)院,武漢 430070)
基于子結(jié)構(gòu)的有限元模型修正方法
翁 順1, 左 越1, 朱宏平1, 陳 波2, 趙會(huì)賢1, 田 煒1, 顏永逸1
(1.華中科技大學(xué) 土木工程與力學(xué)學(xué)院,控制結(jié)構(gòu)湖北省重點(diǎn)實(shí)驗(yàn)室,武漢 430074; 2.武漢理工大學(xué) 土木工程與建筑學(xué)院,武漢 430070)
提出了一種基于子結(jié)構(gòu)的有限元模型修正方法。該方法將整體結(jié)構(gòu)有限元模型劃分為多個(gè)子結(jié)構(gòu)模型,求解獨(dú)立子結(jié)構(gòu)的主模態(tài)特征解和特征靈敏度;通過(guò)位移協(xié)調(diào)條件和能量方程,約束相鄰獨(dú)立子結(jié)構(gòu),得到整體結(jié)構(gòu)的特征解和特征靈敏度;并以整體結(jié)構(gòu)模態(tài)和結(jié)構(gòu)試驗(yàn)?zāi)B(tài)的殘差為目標(biāo)函數(shù),通過(guò)調(diào)整子結(jié)構(gòu)單元參數(shù),完成有限元模型修正。當(dāng)結(jié)構(gòu)局部參數(shù)發(fā)生變化,通過(guò)分析某一個(gè)或幾個(gè)子結(jié)構(gòu)即可求解整體結(jié)構(gòu)特征解和特征靈敏度,而不需要分析其他未發(fā)生變化的子結(jié)構(gòu)。由于子結(jié)構(gòu)模型尺寸遠(yuǎn)小于整體結(jié)構(gòu),該方法能夠極大地提高有限元模型修正方法的精度和效率。
有限元模型修正; 損傷識(shí)別; 子結(jié)構(gòu);特征解; 特征靈敏度
精確的有限元模型是結(jié)構(gòu)健康評(píng)估、動(dòng)態(tài)分析以及優(yōu)化設(shè)計(jì)的基礎(chǔ)。建模過(guò)程中參數(shù)的不確定性、邊界條件的假定等,導(dǎo)致有限元模型的動(dòng)態(tài)響應(yīng)和實(shí)測(cè)的試驗(yàn)數(shù)據(jù)存在誤差。設(shè)計(jì)規(guī)范規(guī)定, 有限元模型必須通過(guò)振動(dòng)模態(tài)試驗(yàn)或者地面共振試驗(yàn)來(lái)檢驗(yàn)。因此,近年來(lái), 有限元模型修正技術(shù)得到了長(zhǎng)足的發(fā)展[1-8]。根據(jù)修正對(duì)象的不同可將修正方法分為矩陣型方法和設(shè)計(jì)參數(shù)型方法。前者直接重建分析模型的剛度和質(zhì)量矩陣。后者重復(fù)的修正有限元模型的物理參數(shù)來(lái)使模型的模態(tài)特性(頻率和振型)和現(xiàn)場(chǎng)試驗(yàn)?zāi)B(tài)的殘差最小化。設(shè)計(jì)參數(shù)型方法物理意義明確,并且修正后的系統(tǒng)矩陣在對(duì)稱性、正定性和稀疏性等方面與修正前的矩陣保持一致。其中,基于靈敏度分析的設(shè)計(jì)參數(shù)型修正方法為模型修正過(guò)程提供有效地優(yōu)化搜索方向,效率高,被廣泛的應(yīng)用于實(shí)際工程結(jié)構(gòu)中[9-13]。
大多數(shù)設(shè)計(jì)參數(shù)型修正方法采用優(yōu)化技術(shù),通過(guò)重復(fù)計(jì)算模型的特征解和特征靈敏度矩陣,尋找最接近實(shí)際工程結(jié)構(gòu)的物理參數(shù),實(shí)現(xiàn)模型修正過(guò)程。但是,土木工程結(jié)構(gòu)體積龐大,為了有效地模擬實(shí)際工程結(jié)構(gòu),其有限元模型通常由大量的單元組成,并且包含許多需要被修正的結(jié)構(gòu)設(shè)計(jì)參數(shù),優(yōu)化迭代過(guò)程通常需要花費(fèi)很多時(shí)間,特別是從龐大的有限元模型的系統(tǒng)矩陣中重復(fù)計(jì)算特征解和特征靈敏度[14-16]。
子結(jié)構(gòu)方法將龐大的整體結(jié)構(gòu)分解為若干個(gè)較小的獨(dú)立子結(jié)構(gòu),通過(guò)位移協(xié)調(diào)條件和力的協(xié)調(diào)條件建立整體結(jié)構(gòu)和子結(jié)構(gòu)之間的關(guān)系。由于獨(dú)立子結(jié)構(gòu)系統(tǒng)矩陣較小,未知參數(shù)少,在大型結(jié)構(gòu)的模型修正及其相關(guān)應(yīng)用中具有以下優(yōu)勢(shì)[17]: ①子結(jié)構(gòu)系統(tǒng)矩陣較小,分析精度和效率高;②各子結(jié)構(gòu)相互獨(dú)立,可以對(duì)獨(dú)立子結(jié)構(gòu)進(jìn)行修正和重新分析,不需要分析整體結(jié)構(gòu);③針對(duì)當(dāng)前的多核計(jì)算機(jī)硬件,進(jìn)行并行計(jì)算,提高計(jì)算效率。
本文研究基于子結(jié)構(gòu)的有限元模型修正方法。該方法通過(guò)組集子結(jié)構(gòu)主模態(tài)的特征解和特征靈敏度,求解整體結(jié)構(gòu)的特征解和特征靈敏度,用于有限元模型修正。在計(jì)算整體結(jié)構(gòu)對(duì)結(jié)構(gòu)單元參數(shù)的特征靈敏度時(shí),只需要計(jì)算包含該單元參數(shù)的一個(gè)子結(jié)構(gòu)的特征靈敏度,其他子結(jié)構(gòu)的特征靈敏度為0。由于子結(jié)構(gòu)的尺寸遠(yuǎn)遠(yuǎn)小于整體結(jié)構(gòu),該子結(jié)構(gòu)方法極大的提高結(jié)構(gòu)特征值和特征靈敏度的計(jì)算效率,從而有效地提高了模型修正過(guò)程的效率。基于正向子結(jié)構(gòu)的有限元模型修正方法應(yīng)用于一個(gè)實(shí)際工程橋梁的有限元模型修正中,驗(yàn)證其精度和效率。
將N個(gè)自由度的整體結(jié)構(gòu)分解為NS個(gè)子結(jié)構(gòu)。以第j(j=1, 2, …,NS)個(gè)子結(jié)構(gòu)為例,假如第j個(gè)子結(jié)構(gòu)有n(j)個(gè)自由度,剛度矩陣為K(j),質(zhì)量矩陣為M(j),其n(j)對(duì)特征值和特征向量為
[Φ(j)]TK(j)Φ(j)=Λ(j)
[Φ(j)]TM(j)Φ(j)=I(j)
(j=1, 2, …,NS)
(1)
KRON[18]的子結(jié)構(gòu)方法基于虛功原理和幾何相容性原則,通過(guò)在子結(jié)構(gòu)相鄰界面處施加約束,得到整體結(jié)構(gòu)的特征方程為
(2)
其中,
Γ=[CΦp]T
Λp=Diag[Λ(1),Λ(2),…,Λ(NS)]
Φp=Diag[Φ(1),Φ(2),…,Φ(NS)]
(3)

從能量守恒的觀點(diǎn)出發(fā),全部子結(jié)構(gòu)的所有階模態(tài)都會(huì)對(duì)整體結(jié)構(gòu)的特征模態(tài)有貢獻(xiàn),也就是說(shuō),要組裝得到Λp和Φp需要計(jì)算全部子結(jié)構(gòu)的所有階特征解。實(shí)際工程結(jié)構(gòu)往往只計(jì)算部分低階模態(tài),求解子結(jié)構(gòu)的全部模態(tài)效率很低。為了克服這個(gè)困難,本文通過(guò)引入模態(tài)截?cái)喾椒▉?lái)提高KRON子結(jié)構(gòu)方法的效率。選取每個(gè)子結(jié)構(gòu)的低階模態(tài)作為“主模態(tài)”,剩下的高階模態(tài)作為“從模態(tài)”。在組裝整體結(jié)構(gòu)的特征方程時(shí)只需要計(jì)算主模態(tài),從模態(tài)的能量貢獻(xiàn)用剩余柔度來(lái)補(bǔ)充。


(j=1,2,…,NS)

(4)
整體結(jié)構(gòu)特征方程式(2)可以根據(jù)主模態(tài)和從模態(tài)重新分解為
(5)
根據(jù)式(5)的第二行,從模態(tài)的模態(tài)參與系數(shù)可以表示為
(6)
將式(6)代入式(5)中得到
(7)
(8)
(9)
式(9)的第二行將τ表示為zm,并且代入式(9)的第一行,可得到
(10)

(11)
(12)
簡(jiǎn)化后的特征方程式(10)的尺寸等于NPm×NPm,遠(yuǎn)小于原始特征方程式(2)的尺寸NP×NP,因此將極大地提高子結(jié)構(gòu)方法求解特征解的效率。
本節(jié)以第R個(gè)子結(jié)構(gòu)的第r單元參數(shù)為例,推導(dǎo)整體結(jié)構(gòu)第i階模態(tài)對(duì)參數(shù)r的特征靈敏度矩陣。式(10)用第i階模態(tài)表示為
(13)

(14)

(15)
其中,
(16)



(18)
將式(18)對(duì)參數(shù)r求偏導(dǎo),可以得到第i階模態(tài)的特征向量靈敏度為
(19)

(20)
式中,ci為第i階模態(tài){zi}的參與系數(shù)。
將式(20)代入式(14)得
(21)
求解式(21),得到常向量{vi}。
求解式(13)的特征解,特征向量{zi}滿足正交條件
{zi}T{zi}=1
(22)
式(22)對(duì)r求偏導(dǎo)后得
(23)
將式(20)代入式(23)中,得到第i階模態(tài)的參與系數(shù)ci為
(24)
求解到常向量{vi}和第i階模態(tài)的參與系數(shù)ci后即可得到
(25)

在基于靈敏度分析的模型修正過(guò)程中,將有限元模型的計(jì)算模態(tài)同結(jié)構(gòu)試驗(yàn)?zāi)B(tài)的殘差作為目標(biāo)函數(shù)
(26)

基于靈敏度分析的設(shè)計(jì)參數(shù)型修正方法,求解目標(biāo)函數(shù)關(guān)于結(jié)構(gòu)設(shè)計(jì)參數(shù)的靈敏度矩陣,提供最優(yōu)的優(yōu)化搜索方向。目標(biāo)函數(shù)對(duì)結(jié)構(gòu)設(shè)計(jì)參數(shù)的靈敏度矩陣為
(27)
其中特征值和特征向量對(duì)于參數(shù)r的靈敏度矩陣表示為

(28)
本研究采用子結(jié)構(gòu)方法求解所有單元參數(shù)的一階偏導(dǎo)。在有限元模型修正過(guò)程中,基于子結(jié)構(gòu)方法通過(guò)式(26)構(gòu)建目標(biāo)函數(shù),通過(guò)式(29)建立靈敏度矩陣為優(yōu)化算法提供搜索方向,采用常用的基于Trust-region的優(yōu)化算法,完成有限元模型修正過(guò)程。
為了驗(yàn)證基于正向子結(jié)構(gòu)的有限元模型修正方法在實(shí)際結(jié)構(gòu)中的可行性和計(jì)算效率,將其應(yīng)用于西澳大利亞的巴拉巴拉河大橋的有限元模型修正。根據(jù)設(shè)計(jì)圖紙可以建立如圖1所示有限元模型。該有限用模型的橋包含907個(gè)單元、947個(gè)節(jié)點(diǎn)。每個(gè)節(jié)點(diǎn)有6個(gè)自由度,共5 420個(gè)自由度。在現(xiàn)場(chǎng)的振動(dòng)測(cè)試中,加速度計(jì)放置在7個(gè)縱梁對(duì)應(yīng)的7列上,每列有19個(gè)測(cè)量點(diǎn),共133個(gè)測(cè)量點(diǎn)。通過(guò)現(xiàn)場(chǎng)模態(tài)試驗(yàn)提取前10階自振頻率和振型(見(jiàn)圖2)。

圖1 Balla Balla Bridge有限元模型和子結(jié)構(gòu)劃分方法Fig.1 FE model of Balla Balla Bridge and the substructures
分別采用傳統(tǒng)的基于整體結(jié)構(gòu)的有限元模型修正方法和本文所提出的基于子結(jié)構(gòu)的有限元模型修正方法對(duì)有限元模型進(jìn)行修正。一共選取了包括膈、梁、板的楊氏模量(E),剪力連接件的拉壓剛度(EA)和抗彎剛度(EIxx,EIxy)在內(nèi)的1 289個(gè)物理參數(shù)作為修正參數(shù)。目標(biāo)函數(shù)根據(jù)式(26)定義為有限元模型的頻率振型與試驗(yàn)?zāi)B(tài)的頻率振型之間的殘差。由于振型測(cè)量的誤差大于頻率,將振型的權(quán)重設(shè)定為0.1而頻率的權(quán)重設(shè)定為1.0。
首先應(yīng)用傳統(tǒng)的整體結(jié)構(gòu)有限元模型修正方法完成有限元模型修正。采用Lanczos方法從整體結(jié)構(gòu)系統(tǒng)矩陣中計(jì)算結(jié)構(gòu)特征解構(gòu)建目標(biāo)函數(shù),用NELSON方法計(jì)算特征靈敏度用于提供優(yōu)化搜索方向。由于試驗(yàn)?zāi)B(tài)只測(cè)到部分模態(tài),在每次優(yōu)化迭代中,計(jì)算有限元模型前30階模態(tài)并與實(shí)測(cè)的前10階模態(tài)進(jìn)行匹配。通過(guò)模態(tài)置信準(zhǔn)則MAC(Modal Assurance Creterion)
(29)
從有限元模型計(jì)算的30階模態(tài)中找出與試驗(yàn)?zāi)B(tài)匹配的10階模態(tài)。
模型修正過(guò)程迭代69次后收斂到預(yù)先設(shè)定的閾值,在普通個(gè)人電腦上每一步迭代大約花費(fèi)1.26 h,模型修正過(guò)程共花費(fèi)86.16 h,其收斂過(guò)程如圖3所示。
然后,應(yīng)用所提出的基于子結(jié)構(gòu)的方法有限元模型修正方法完成有限元模型修正,并采用和上述整體結(jié)構(gòu)方法相同的優(yōu)化算法、修正參數(shù)、收斂準(zhǔn)則等前提條件,并且在相同的個(gè)人電腦上運(yùn)行程序。用本文所提出的子結(jié)構(gòu)方法計(jì)算整體結(jié)構(gòu)的特征值和特征靈敏度。將整體結(jié)構(gòu)沿縱向分成11個(gè)子結(jié)構(gòu),如表1所示。每一步優(yōu)化迭代過(guò)程中,通過(guò)組集獨(dú)立子結(jié)構(gòu)的特征值求解整體結(jié)構(gòu)的特征值。計(jì)算整體結(jié)構(gòu)對(duì)某一單元參數(shù)的特征靈敏度,只需要求解包含這個(gè)單元參數(shù)的一個(gè)子結(jié)構(gòu)的主模態(tài)特征靈敏度,而其他的子結(jié)構(gòu)的特征靈敏度設(shè)定為0。

圖2 Balla Balla Bridge的測(cè)量頻率和振型Fig.2 Measured frequencies and mode shapes

圖3 有限元模型修正過(guò)程Fig.3 The finite element model udating process
比較有限元模型修正前后結(jié)構(gòu)的頻率和振型,如表2所示。該基于子結(jié)構(gòu)的有限元模型修正方法可以得到與整體結(jié)構(gòu)非常相近的結(jié)果。修正后模型頻率和實(shí)驗(yàn)頻率之間的平均差異<1%。修正后振型的相關(guān)系數(shù)MAC從0.85提高到0.93。對(duì)于這個(gè)中型結(jié)構(gòu),在達(dá)到相同計(jì)算精度的條件下,基于子結(jié)構(gòu)的有限元模型修正過(guò)程為傳統(tǒng)方法的大約一半,極大地提高了有限元模型修正的效率。子結(jié)構(gòu)方法應(yīng)用于大型結(jié)構(gòu)對(duì)效率的提高將更明顯。
由于子結(jié)構(gòu)中保留的主模態(tài)的數(shù)量會(huì)影響求解特征值和特征靈敏度的精度和效率,本文在模型修正過(guò)程中動(dòng)態(tài)選取主模態(tài)。在模型修正初始階段,選取每個(gè)子結(jié)構(gòu)的前40個(gè)模態(tài)作為主模態(tài)來(lái)計(jì)算整體結(jié)構(gòu)的前30個(gè)特征值和特征靈敏度。然后,隨著參數(shù)逐漸接近最優(yōu)解,逐步增加子結(jié)構(gòu)的主模態(tài)的數(shù)量。直到在最后幾步中,在每個(gè)子結(jié)構(gòu)中保留90個(gè)主模態(tài)來(lái)提高特征值和特征向量計(jì)算的精度。該基于子結(jié)構(gòu)的有限元模型修正過(guò)程迭代76次收斂到相同的閾值,如圖3所示。基于子結(jié)構(gòu)的有限元模型修正方法所花費(fèi)的計(jì)算時(shí)間如表3所示,共花費(fèi)48.07 h,所需時(shí)間大約是傳統(tǒng)的整體結(jié)構(gòu)模型修正方法的56%。

表1 整體結(jié)構(gòu)劃分為11個(gè)子結(jié)構(gòu)后各子結(jié)構(gòu)的信息

表2 修正前后橋的頻率和振型

表3 基于整體結(jié)構(gòu)有限元模型修正方法和基于子結(jié)構(gòu)的有限元模型修正方法計(jì)算時(shí)間對(duì)比
本文提出了一種基于子結(jié)構(gòu)的有限元模型修正方法。將整體結(jié)構(gòu)有限元模型劃分為多個(gè)獨(dú)立子結(jié)構(gòu)有限元模型,通過(guò)求解一個(gè)或幾個(gè)發(fā)生變化的獨(dú)立子結(jié)構(gòu)特征解,即可求解整體結(jié)構(gòu)特征解,并用于有限元模型修正和損傷識(shí)別。并且,采用子結(jié)構(gòu)方法只需要計(jì)算某一個(gè)子結(jié)構(gòu)特征解靈敏度矩陣即可完成對(duì)整體結(jié)構(gòu)特征靈敏度的求解。當(dāng)結(jié)構(gòu)局部發(fā)生損傷,只需要重復(fù)某一個(gè)或幾個(gè)子結(jié)構(gòu)模型,避免對(duì)整體結(jié)構(gòu)模型重復(fù)分析,從而有效地提高大型結(jié)構(gòu)有限元模型修正的精度和效率。基于子結(jié)構(gòu)的有限元模型修正方法可以有效地用于對(duì)大型結(jié)構(gòu)的動(dòng)力分析和健康監(jiān)測(cè)中。
[ 1 ] MOTTERSHEAD J E, FRISWELL M I. Model updating in structural dynamics: a survey[J]. Journal of Sound and Vibration,1993,167(2): 347-375.
[ 2 ] 朱宏平, 黃民水. 基于環(huán)境激勵(lì)的橋梁結(jié)構(gòu)動(dòng)力有限元模型修正研究[J]. 華中科技大學(xué)學(xué)報(bào), 2009,26(1):1-11. ZHU Hongping, HUANG Minshui. Study on dynamic finite element model updating of bridge structures based on ambient excitation[J]. Journal of Huazhong University of Science and Technology, 2009, 26(1):1-11.
[ 3 ] 方圣恩, 基于有限元模型修正的結(jié)構(gòu)損傷識(shí)別方法研究[D]. 長(zhǎng)沙:中南大學(xué),2010.
[ 4 ] 姜東, 費(fèi)慶國(guó), 吳邵慶. 基于攝動(dòng)法的不確定性有限元模型修正方法研究[J].計(jì)算力學(xué)學(xué)報(bào), 2014, 31(4): 431-437. JIANG Dong,FEI Qingguo, WU Shaoqing. A study on stochastic finite element model updating based on perturbation approach[J]. Chinese Journal of Computational Mechanics, 2014,31(4): 431-437.
[ 5 ] 任偉新, 陳華斌. 基于響應(yīng)面的橋梁有限元模型修正[J]. 土木工程學(xué)報(bào), 2008,41(12): 73-78. REN Weixin, CHEN Huabin. Response-surface based on finite element model updating of bridge structures[J].China Civil Engineering Journal, 2008, 41(12): 73-78.
[ 6 ] 郭勤濤, 張令彌, 費(fèi)慶國(guó). 結(jié)構(gòu)動(dòng)力學(xué)有限元模型修正的發(fā)展-模型確認(rèn)[J].力學(xué)進(jìn)展, 2006, 36(1): 36-42. GUO Qintao, ZHANG Lingmi, FEI Qingguo. From FE model updating to model validation: advances in modeling of dynamic structures[J]. Advances in Mechanics, 2006, 36(1): 36-42.
[ 7 ] BROWNJOHN J M W, MOYO P, OMENZETTER P, et al. Assessment of highway bridge upgrading by dynamic testing and finite-element model updating[J]. Journal of Bridge Engineering,2003,8(3): 162-172.
[ 8 ] JAISHI B, REN W X. Structural finite element model updating using ambient vibration test results[J]. Journal of Structural Engineering,2005, 131(4): 617-628.
[ 9 ] ZHU H P, MAO L, WENG S. A sensitivity-based structural damage identification method with unknown input excitation using transmissibility concept[J]. Journal of Sound and Vibration, 2014, 333(26): 7135-7150.
[10] WENG S, ZHU A Z, ZHU H P, et al. Dynamic condensation approach to the calculation of eigensensitivity[J]. Computers and Structures, 2014, 132(1): 55-64.
[11] WENG S, ZHU H P, XIA Y, et al. Substructuring approach to the calculation of higher-order eigensensitivity[J]. Computers and Structures, 2013, 117(2): 23-33.
[12] WENG S, ZHU H P, XIA Y, et al. Damage detection using the eigenparameter decomposition of substructural flexibility matrix[J]. Mechanical Systems and Signal Processing, 2013, 34(1/2): 19-38.
[13] WENG S, XIA Y, XU Y L, et al. Improved substructuring method for eigensolutions of large-scale structures[J]. Journal of Sound and Vibration, 2009, 323(3/4/5): 718-736.
[14] CHOI K K, KIM N H. Structural sensitivity analysis and optimization 1: linear systems[M].New York: Springer Science and Business Media, 2005.
[15] FOX R L, KAPOOR M P. Rate of change of eigenvalues and eigenvectors[J]. AIAA Journal,1968, 6(12): 2426-2429.
[16] NELSON R B. Simplified calculation of eigenvector derivatives[J].AIAA Journal,1976,24(9):823-832.
[17] CRAIG R R. Coupling of substructures for dynamic analysis:an overview[C]//Atlanta: Proceedings of 41st AIAA/ASME/ASCE/ASC Structures, Structural Dynamics, and Materials Conference, 2000:171-179.
[18] KRON G. Diakoptics[M].London:Macdonald and Co., 1963: 46-87.
Model updating based on a substructuring method
WENGShun1,ZUOYue1,ZHUHongping1,CHENBo2,ZHAOHuixian1,TIANWei1,YANYongyi1
(1. School of Civil Engineering and Mechanics,Huazhong University of Science and Technology,Wuhan 430074, China; 2. School of Civil Engineering and Architecture,Wuhan University of Technology,Wuhan 430070, China)
This paper proposes a substructure-based model updating method. The global structure was divided into independent manageable substructures, which were analyzed to calculate the substructural eigensolutions and eigensensitivity. Afterwards, the independent substructures were constrained to calculate the eigensolutions and eigensensitivity of the global structure. Finally, the eigensolutions were used for the objective function and the eigensensitivity was used to indicate the searching direction to achieve the model updating. When the local area of a structure was changed, only one or more substructures were analyzed whereas the other substructures were untouched. Since the substructures were much smaller than the global structure, the proposed method significantly improved the accuracy and efficiency of the model updating process.
model updating; substructure; damage identification; substructure; eigensolution; eigensensitivity
國(guó)家自然科學(xué)基金(51328802;51108205);中央高校基本科研業(yè)務(wù)費(fèi)專項(xiàng)資金(2014TS130;2015MS064);武漢城建委科研項(xiàng)目(201511)
2015-11-23 修改稿收到日期:2016-01-31
翁順 女,博士,副教授,1982年生
TH212;TH213.3
A
10.13465/j.cnki.jvs.2017.04.016