錐體海洋平臺(tái)結(jié)構(gòu)冰荷載的離散單元分析

李 曄

（中國(guó)船級(jí)社 海洋工程技術(shù)中心，天津300457）

錐體海洋平臺(tái)結(jié)構(gòu)冰荷載的離散單元分析

李 曄

（中國(guó)船級(jí)社 海洋工程技術(shù)中心，天津300457）

在冰區(qū)油氣開(kāi)發(fā)中，錐體結(jié)構(gòu)可以有效降低冰力，避免強(qiáng)烈的冰激振動(dòng)，是目前渤海油氣平臺(tái)的主要結(jié)構(gòu)形式。為研究海冰與錐體結(jié)構(gòu)的相互作用過(guò)程，文章建立了適用于模擬海冰破碎特性的離散單元模型。該模型將海冰離散為若干個(gè)具有粘接-破碎功能的顆粒單元，并通過(guò)海冰彎曲試驗(yàn)確定了單元間的粘接強(qiáng)度；然后對(duì)海冰與錐體結(jié)構(gòu)的作用過(guò)程進(jìn)行了數(shù)值計(jì)算，獲得了相應(yīng)的動(dòng)冰荷載及冰振響應(yīng)；在此基礎(chǔ)上討論了不同錐角影響下冰荷載及結(jié)構(gòu)振動(dòng)響應(yīng)的變化規(guī)律。結(jié)果表明，水平方向冰荷載及結(jié)構(gòu)冰振響應(yīng)隨錐角的增加明顯增加，而豎直方向冰荷載則顯著降低。該離散單元模型還可進(jìn)一步應(yīng)用于不同類型抗冰結(jié)構(gòu)的冰荷載分析，有助于解決冰區(qū)結(jié)構(gòu)物的抗冰結(jié)構(gòu)設(shè)計(jì)和冰致疲勞分析。

海冰；錐體結(jié)構(gòu)；冰荷載；冰激振動(dòng)；離散單元模型

0 引 言

海冰與海洋結(jié)構(gòu)相互作用時(shí)產(chǎn)生的冰荷載主要受海冰性質(zhì)與結(jié)構(gòu)形式的影響。在冰區(qū)油氣開(kāi)發(fā)的導(dǎo)管架式海洋平臺(tái)結(jié)構(gòu)設(shè)計(jì)中，人們最初采用的是直立結(jié)構(gòu)形式，但直立結(jié)構(gòu)存在冰力偏大和冰激振動(dòng)嚴(yán)重等問(wèn)題[1-3]。由于海冰是一種典型的脆性材料，其彎曲強(qiáng)度要明顯低于壓縮強(qiáng)度。為此，人們提出了在直立結(jié)構(gòu)的水面處安裝破冰錐體以降低冰荷載。海冰與斜面或錐體結(jié)構(gòu)相互作用過(guò)程的研究一直是人們所關(guān)注的問(wèn)題[4-6]。特別是近年來(lái)，渤海錐體結(jié)構(gòu)的廣泛采用，錐體結(jié)構(gòu)冰荷載及相應(yīng)結(jié)構(gòu)冰振響應(yīng)引起了人們更加廣泛的關(guān)注[7-10]。

近年來(lái)，人們通過(guò)理論分析、試驗(yàn)測(cè)試和數(shù)值方法研究了海冰與錐體結(jié)構(gòu)的相互作用過(guò)程[10-12]。在海冰與錐體結(jié)構(gòu)相互作用過(guò)程中，海冰以彎曲破壞為主，其中海冰彎曲強(qiáng)度是決定冰力的重要參數(shù)；此外，錐體結(jié)構(gòu)冰荷載還與冰厚和冰速等海冰參數(shù)、錐角和錐徑等結(jié)構(gòu)形式等有著密切的關(guān)系。由于海冰與錐體結(jié)構(gòu)作用時(shí)的破壞過(guò)程非常復(fù)雜，影響因素較多，其冰荷載特征還需要進(jìn)一步深入研究。

海冰荷載的數(shù)值模擬大多采用有限元方法[11,13]。該方法能夠很好地分析海冰與海洋結(jié)構(gòu)相互作用中冰內(nèi)應(yīng)力分布以及海冰破碎特性，但尚不能合理地模擬海冰連續(xù)破碎過(guò)程中的動(dòng)冰力。在有限元方法的基礎(chǔ)上，人們發(fā)展了一種粘接單元模型用以模擬海冰的破碎特征，并由此確定相應(yīng)的冰荷載[14]。最近，基于離散元模型的海冰與結(jié)構(gòu)物相互作用研究取得了很大的進(jìn)展[15-18]。它能夠合理地模擬海冰與結(jié)構(gòu)作用中的破碎過(guò)程，并獲得相應(yīng)的冰力時(shí)程。離散元模型在處理海冰與錐體結(jié)構(gòu)作用時(shí)，將冰排離散為若干個(gè)具有粘結(jié)—破碎功能的顆粒單元，通過(guò)顆粒材料性質(zhì)及單元間粘接強(qiáng)度模擬海冰的破碎模式[16]。海洋平臺(tái)可簡(jiǎn)化為具有等效的質(zhì)量—?jiǎng)偠取枘岬挠?jì)算模型。由于顆粒單元模型具有物理意義明確、計(jì)算模型簡(jiǎn)便、接觸判斷快速等優(yōu)點(diǎn)，可較好地模擬海冰與海洋結(jié)構(gòu)相互作用中的破碎特性，并由此確定冰荷載的動(dòng)力特性以及海冰作用下的結(jié)構(gòu)振動(dòng)。

為此，本文采用粘接顆粒離散單元模型構(gòu)造海冰材料，通過(guò)對(duì)海冰彎曲試驗(yàn)的離散元分析確定單元間的粘結(jié)強(qiáng)度；在此基礎(chǔ)上對(duì)海冰與錐體之間的相互作用進(jìn)行數(shù)值分析，確定相應(yīng)的冰荷載及結(jié)構(gòu)冰振響應(yīng)；最后對(duì)錐角影響下的冰荷載和冰振響應(yīng)進(jìn)行討論。

1 海冰的離散單元模型及彎曲強(qiáng)度數(shù)值模擬

海冰在與錐體結(jié)構(gòu)的相互作用過(guò)程中主要呈現(xiàn)彎曲破碎特性。這里首先建立海冰的離散單元模型，并通過(guò)海冰的彎曲試驗(yàn)進(jìn)行模型的檢驗(yàn)和計(jì)算參數(shù)的確定。

1.1 海冰離散單元模型

在采用離散單元模型計(jì)算海冰的動(dòng)力作用及破碎過(guò)程時(shí)，將海冰離散為若干具有一定質(zhì)量和大小的顆粒單元，單元間具有相應(yīng)的粘接作用。顆粒之間可視為在法向上由彈簧與阻尼器串聯(lián)連接，在切向上由彈簧、阻尼器和滑動(dòng)摩擦器連接。在顆粒相互作用過(guò)程中，一般考慮單元間因相對(duì)速度和彈性變形而引起的粘彈性作用力，并采用Mohr-Coulomb摩擦定律計(jì)算剪切力。

海冰彎曲強(qiáng)度的試驗(yàn)表明，海冰彎曲強(qiáng)度是海冰鹵水體積（溫度、鹽度）、加載速率等因素的函數(shù)[20-21]。由于海冰的彎曲強(qiáng)度主要由顆粒間的粘接強(qiáng)度決定，而粘結(jié)強(qiáng)度主要受鹵水體積的影響，可將海冰單元間的粘接強(qiáng)度用單元間最大粘接強(qiáng)度σmbax表示，即

式中：β（ vb）為鹵水體積vb影響下的海冰強(qiáng)度折減系數(shù)。海冰單元間的最大粘接強(qiáng)度σmbax可通過(guò)海冰彎曲強(qiáng)度的敏度分析進(jìn)行確定?？紤]海冰的壓縮和彎曲強(qiáng)度與鹵水體積有相似的對(duì)應(yīng)關(guān)系[21]，則有

式中：vb可設(shè)為海冰溫度和鹽度的函數(shù)[22]，即

式中：T為海冰溫度（℃），S為海冰鹽度（‰）。

1.2 海冰彎曲強(qiáng)度的離散單元模擬

采用以上海冰離散元模型對(duì)海冰三點(diǎn)彎曲強(qiáng)度試驗(yàn)進(jìn)行數(shù)值模擬，其主要計(jì)算參數(shù)列于表1中。海冰單元按一定排開(kāi)方式組合后粘接為海冰彎曲試樣，如圖1（a）所示。試樣下側(cè)的兩個(gè)支撐點(diǎn)固定，其上側(cè)加載點(diǎn)以勻速0.01 m/s向下運(yùn)動(dòng)。計(jì)算得到的海冰彎曲斷裂過(guò)程及所對(duì)應(yīng)的冰力時(shí)程曲線分別如圖1和圖2所示。圖中顆粒顏色代表其平均接觸力，由此可發(fā)現(xiàn)海冰受彎時(shí)中性層的位置及其上下兩側(cè)的對(duì)稱受力情況。從圖2中可以發(fā)現(xiàn)，應(yīng)力隨著時(shí)間大體呈線性增長(zhǎng)，并存在一定的波動(dòng)性。當(dāng)t=0.033 s時(shí)，應(yīng)力達(dá)到最大值，此時(shí)海冰試樣斷裂，其對(duì)應(yīng)的海冰彎曲強(qiáng)度σf為1.57 MPa。

表1 海冰彎曲強(qiáng)度的離散單元模擬主要計(jì)算參數(shù)Tab.1 Main computational para meters in numerical simulation of sea ice flexural strength with DEM

圖1 離散元數(shù)值模擬的海冰試樣彎曲破碎過(guò)程Fig.1 Bending failure process of sea ice sample simulated with DEM

圖2 離散元模擬的海冰彎曲過(guò)程中的應(yīng)力變化時(shí)程Fig.2 The stress versus time in sea ice bending process simulated with DEM

圖3 不同粘接強(qiáng)度下離散元計(jì)算的海冰彎曲強(qiáng)度Fig.3 Sea ice flexural strengths simulated with DEM under various bonding strengths

圖4 不同鹵水體積下計(jì)算的海冰彎曲強(qiáng)度Fig.4 Sea ice flexural strengths simulated under various brine volumes

2 錐體結(jié)構(gòu)冰荷載的離散元模擬

采用離散單元模型對(duì)海冰與錐體結(jié)構(gòu)的相互作用過(guò)程進(jìn)行數(shù)值分析，由此確定冰力時(shí)程及平臺(tái)結(jié)構(gòu)的冰振響應(yīng)，并在此基礎(chǔ)上進(jìn)一步討論錐角對(duì)冰力的影響。

2.1 海冰與錐體結(jié)構(gòu)的離散單元模擬

在海冰與錐體結(jié)構(gòu)的相互作用過(guò)程中，海冰主要以彎曲破碎為主，顆粒間發(fā)生拉伸和剪切破壞。遼東灣現(xiàn)場(chǎng)監(jiān)測(cè)的海冰與錐體結(jié)構(gòu)的作用狀況如圖5所示。下面采用離散元模型對(duì)海冰與錐體結(jié)構(gòu)作用時(shí)的破碎過(guò)程、冰荷載時(shí)程和冰激結(jié)構(gòu)振動(dòng)響應(yīng)進(jìn)行數(shù)值計(jì)算，主要計(jì)算參數(shù)列于表2中。

圖5 遼東灣JZ20-2油氣田海冰與錐體相互作用Fig.5 Conical offshore structure in the JZ20-2 oil field of Liaodong Bay and the interaction between ice cover and cone

采用離散元方法模擬的海冰與錐體結(jié)構(gòu)相互作用過(guò)程如圖6所示，圖中顆粒顏色代表其平均作用力，也反應(yīng)了冰內(nèi)應(yīng)力場(chǎng)分布的大致規(guī)律。從海冰破碎后的冰塊尺度來(lái)看，采用平行粘接模型可模擬出海冰與錐體結(jié)構(gòu)作用時(shí)的破碎長(zhǎng)度。計(jì)算得到的x向樁腿冰力時(shí)程如圖6（a）所示，其最大和平均冰力分別為62.68 kN和18.41 kN。該冰力時(shí)程呈現(xiàn)出明顯的周期性特征，這與海冰和錐體結(jié)構(gòu)相互作

用時(shí)的彎曲破碎長(zhǎng)度和冰速密切相關(guān)[7,10]。計(jì)算得到的y向冰力如圖6（b）所示，其最大和平均值分別為18.51 kN和0.67 kN。這表明海冰對(duì)錐體側(cè)面也會(huì)產(chǎn)生一定的沖擊力，但其作用力是基本對(duì)稱的。z方向的冰力時(shí)程如圖6（c）所示，其最大值和均值分別為49.09 kN和16.20 kN，其峰值也顯示出一定的周期性。錐體結(jié)構(gòu)的冰振加速度時(shí)程如圖6（d）所示，其中最大加速度為9.25 gal（1.0 gal=1.0 cm/s2）。以上結(jié)果表明，該離散元模型可以較好地模擬海冰與錐體結(jié)構(gòu)的相互作用過(guò)程，確定動(dòng)冰荷載及錐體結(jié)構(gòu)的冰振響應(yīng)。

表2 海冰與錐體結(jié)構(gòu)作用離散元模擬的主要計(jì)算參數(shù)Tab.2 Main computational parameters of DEM simulation of dynamic ice load on conical structures

圖6 離散元模擬的海冰與錐體結(jié)構(gòu)作用過(guò)程Fig.6 The interaction between sea ice and the conical structure simulated with DEM

圖7 離散元模擬的錐體海洋結(jié)構(gòu)冰荷載及冰速方向冰振加速度Fig.7 Dynamic ice loads and the ice-induced acceleration of conical offshore structure simulated with DEM

2.2 錐角對(duì)錐體結(jié)構(gòu)冰荷載的影響

在海冰與錐體結(jié)構(gòu)相互作用中，海冰以彎曲破壞為主，并由此產(chǎn)生相應(yīng)的動(dòng)冰力。海冰發(fā)生彎曲破壞的模式及其相應(yīng)的破碎長(zhǎng)度在很大程度上受錐體角度的影響。錐體角度越大，海冰受到的錐體作用力的水平分量也就越大，豎向分量則越小。為進(jìn)一步分析錐角對(duì)冰荷載的影響，這里將錐角分別設(shè)為65.8°、55.0°、46.5°和39.0°，其他計(jì)算參數(shù)采用表1和表2中的數(shù)值。當(dāng)海冰與錐體的作用時(shí)間為6.5 s時(shí)，計(jì)算得到的不同錐角下海冰破碎現(xiàn)象以及冰力時(shí)程分別如圖8和圖9所示?？梢钥闯觯１扑榈某叽缭诓煌F角下有明顯差異。x方向的冰力時(shí)程更清晰地表明，錐角越大，冰力周期也就越小，即海冰的破碎尺寸也就越小。

圖8 不同錐角下海冰與錐體的相互作用Fig.8 Interaction between sea ice and the conical structure with different conical angles

圖9 離散元模擬的不同錐角下的水平和豎向冰力Fig.9 Ice loads in horizontal and vertical directions simulated with DEM under various conical angles

圖10 不同錐角下的冰力及結(jié)構(gòu)冰振加速度Fig.10 Mean ice loads and accelerations of ice-induced structure vibration simulated under various conical angles

在不同錐角下，錐體結(jié)構(gòu)的水平和豎向冰力均值、最大值如圖10（a）和（b）所示?？梢园l(fā)現(xiàn)水平方向冰力的均值和最大值隨錐角的增加而增大，而豎向冰力則均呈減小趨勢(shì)。這主要是由于在錐角由小到大變化過(guò)程中，海冰破碎由彎曲模式向擠壓模式轉(zhuǎn)變，并由此導(dǎo)致冰力水平分量的增加和豎向分量的降低。不同錐角下平臺(tái)結(jié)構(gòu)的冰振加速度均值和最大值如圖10（c）所示，其變化規(guī)律與水平方向冰力是一致的。這主要是因?yàn)槠脚_(tái)結(jié)構(gòu)的冰振響應(yīng)主要受控于水平方向冰力。

3 結(jié) 語(yǔ)

為研究海冰與錐體海洋平臺(tái)結(jié)構(gòu)的作用過(guò)程，本文采用離散單元模型對(duì)海冰的彎曲強(qiáng)度、錐體結(jié)構(gòu)動(dòng)冰力及冰振響應(yīng)進(jìn)行了數(shù)值計(jì)算。通過(guò)海冰彎曲強(qiáng)度的離散元數(shù)值模擬，分析了鹵水體積對(duì)海冰彎曲強(qiáng)度的影響，并由此確定了海冰單元間的粘結(jié)強(qiáng)度，驗(yàn)證了離散元模擬的可行性；在此基礎(chǔ)上對(duì)海冰與錐體結(jié)構(gòu)的作用過(guò)程進(jìn)行了離散元計(jì)算，得到了動(dòng)冰力及結(jié)構(gòu)冰振響應(yīng)；通過(guò)對(duì)不同錐角下結(jié)構(gòu)冰力和冰振響應(yīng)的對(duì)比，發(fā)現(xiàn)錐角越大水平冰力和冰振響應(yīng)就越大，而豎向冰力則相應(yīng)減少。

在海冰與錐體結(jié)構(gòu)相互作用中，海冰的斷裂長(zhǎng)度和正倒錐結(jié)合處的海冰破碎模式是兩個(gè)重要研究?jī)?nèi)容。這需要進(jìn)一步對(duì)海冰離散元模型進(jìn)行改進(jìn)并確定合理的計(jì)算參數(shù)以更加全面地揭示錐體結(jié)構(gòu)的冰力特性。此外，對(duì)海洋平臺(tái)結(jié)構(gòu)冰荷載及冰振響應(yīng)的離散元模擬還需要通過(guò)室內(nèi)模型試驗(yàn)和現(xiàn)場(chǎng)監(jiān)測(cè)結(jié)果進(jìn)行驗(yàn)證，以提高其計(jì)算精度和可靠性。該海冰離散元模型對(duì)深入研究海冰與海洋結(jié)構(gòu)的作用特性具有很好的應(yīng)用前景，可對(duì)不同類型海洋平臺(tái)的冰荷載、冰激結(jié)構(gòu)振動(dòng)機(jī)理和多樁腿結(jié)構(gòu)的冰力掩蔽效應(yīng)進(jìn)行分析。

致謝：感謝大連理工大學(xué)工業(yè)裝備結(jié)構(gòu)分析國(guó)家重點(diǎn)實(shí)驗(yàn)室季順迎教授的有益討論。

[1]Karna T,Qu Y,Kuhnlein W,Yue Q,Bi X.A spectral model of ice forces due to ice crushing[J].Journal of Offshore Mechanics&Arctic Engineering,2007,129:138-145.

[2]Huang Y,Shi Q,Song A.Model test study of the interaction between ice and compliant vertical narrow structure[J].Cold Regions Science and Technology,2007,49:151-160.

[3]Sodhi D.Crushing failure during ice-structure interaction[J].Engineering Fracture Mechanics,2001,68:1889-1921.

[4]季順迎,王安良,車嘯飛,岳前進(jìn).錐體導(dǎo)管架海洋平臺(tái)冰激結(jié)構(gòu)振動(dòng)響應(yīng)分析[J].海洋工程,2011,29(2):32-39.Ji S Y,Wang A L,Che X F,et al.Analysis of ice induced structure vibration of offshore jacket platform with ice breaking cone[J].Ocean Engineering,2011,29(2):32-39.

[5]Huang Y,Li X.Dynamic ice loads on conical structures[J].Theoretical&Applied Mechanics Letters,2011,1:022007.

[6]Timco G.Ice loads on conical-shaped structure[J].Cold Regions Science and Technology,2009,55:1-2.

[7]岳前進(jìn),畢祥軍,于 曉,等.錐體結(jié)構(gòu)的冰激振動(dòng)與冰力函數(shù)[J].土木工程學(xué)報(bào),2001,36(2):16-19.Yue Q J,Bi X J,Yu X,et al.Ice-induced vibration and ice force function of conical structure[J].China Civil Engineering Journal,2001,36(2):16-19.

[8]史慶增,黃 焱,宋 安,等.錐體冰力的實(shí)驗(yàn)研究[J].海洋工程,2004,22(1):88-92.Shi Q Z,Huang Y,Song A,et al.Experimental study of ice force on conical structures[J].Ocean Engineering,2004,22(1):88-92.

[9]李 峰,岳前進(jìn).不同破壞模式下錐體結(jié)構(gòu)冰力定性分析[J].大連理工大學(xué)學(xué)報(bào),2005,45(6):785-789.Li F,Yue Q J.Qualitative analysis of ice forces on conical structures under different failure modes[J].Journal of Dalian University of Technology,2005,45(6):785-789.

[10]Xü N,Yue Q,Bi X,et al.Experimental study of dynamic conical ice force[J].Cold Regions Science and Technology,2015,120:21-29.

[11]王 剛,武文華,岳前進(jìn).錐體接觸寬度對(duì)冰排彎曲破壞模式影響的有限元分析[J].工程力學(xué),2008,25(1):235-240.Wang G,Wu W H,Yue Q J.FEM analysis on ice-bending failure mode with width effect of ice-cone interaction[J].Engineering Mechanics,2008,25(1):235-240

[12]Huang Yan.Model test study of the nonsimultaneous failure of ice before wide conical structures[J].Cold Regions Science and Technology,2010,63:87-96.

[13]Konuk I,Gurtner A,Yu S.Study of dynamic ice and cylindrical structure interaction by the cohesive element method[C]//Proceedings of the 20th International Conference on Port and Ocean Engineering under Arctic Conditions.Lulea,Sweden,2009.

[14]Gagnon R E.A numerical model of ice crushing using a foam analogue[J].Cold Regions Science and Technology,2011,65:335-350.

[15]Lau M,Lawrence K P,Rothenburg L.Discrete element analysis of ice loads on ships and structures[J].Ships and Offshore Structures,2011,6(3):211-221.

[16]狄少丞,季順迎.海冰與自升式海洋平臺(tái)相互作用GPU離散元模擬[J].力學(xué)學(xué)報(bào),2014,46(4):561-571.Di S C,Ji S Y.GPU-based discrete element modelling of interaction between sea ice and jack-up platform structure[J].Chinese Journal of Theoretical and Applied Mechanics,2014,46(4):561-571.

[17]Poloj?rvi A,Tuhkuri J,Pustogvar A.DEM simulations of direct shear box experiments of ice rubble:Force chains and peak loads[J].Cold Regions Science&Technology,2015,116:12-23.

[18]季順迎,李春花,劉 煜.海冰離散元模型的研究回顧及展望[J].極地研究,2012,24(4):315-329.Ji S Y,Li C H,Liu Y.A review of advances in sea ice discrete element models[J].Chinese Journal of Polar Research,2012,24(4):315-329.

[19]Ji S,Shen H H.Internal parameters and regime map for soft polydispersed granular materials[J].Journal of Rheology,2008,52(1):87-103.

[20]Timco G W,O’Brien S.Flexural strength equation for sea ice[J].Cold Regions Science and Technology,1994,22:285-298.

[21]王安良,許 寧,畢祥軍,等.鹵水體積和應(yīng)力速率影響下海冰強(qiáng)度的統(tǒng)一表征[J].海洋學(xué)報(bào),2016,38(9):126-133.Wang A L,Xü N,Bi X J,et al.Unified representation of sea ice strengths under influences of brine volume and stress rate[J].Acta Oceanolo Gica Sinica,2016,38(9):126-133.

[22]Frankenstein G,Garner R.Equations for determining the brine volume sea ice from-0.5℃ to-22.0℃[J].Journal of Glaciology,1967,6(48):943-944.

Analysis of ice load on conical structure with discrete element model

LI Ye

(Offshore Technology Center of China Classification Society,Tianjin 300457,China)

In oil/gas exploitations of ice covered cold regions,the conical structure can reduce ice force and avoid intense structure vibration.It is applied most widely in the Bohai Sea recently.In order to investigate the interaction between ice cover and conical structure,a discrete element method is established to model the sea ice breaking characteristics.The ice cover is dispersed into a series of bonded spherical elements.The bonding strength among ice particles is determined with the flexural test of sea ice.The interaction between the ice cover and the conical structure is simulated with the discrete element method to determine the dynamic ice load and structure vibration.Moreover,the influence of conical angle is discussed based on the simulated results.The ice load and the ice-induced structure response in horizontal direction increase obviously with the increase of cone angle,while the ice load in vertical direction decreases significantly.This discrete element method can also be applied to simulate ice loads of other type structures,and to be aided for the design of ice-resistance structure and the analysis of ice induced structure fatigue.

sea ice;conical structure;ice load;ice-induced vibration;discrete element method

P731.15 P751

A

10.3969/j.issn.1007-7294.2017.10.009

1007-7294（2017）10-1254-09

2017-03-27

國(guó)家重點(diǎn)研發(fā)計(jì)劃（2106YCF1401505）和國(guó)家自然科學(xué)基金項(xiàng)目（41576179，51639004）

李 曄（1975-），男，博士，高級(jí)工程師，E-mail:liye@ccs.org.cn。