腐蚀管道可靠性评价方法研究

帅义，帅健*，刘朝阳

1 中国石油大学(北京) 机械与储运工程学院，北京 102249 2 中国石油管道科技研究中心，廊坊 065000

腐蚀管道可靠性评价方法研究

帅义1，帅健1*，刘朝阳2

1 中国石油大学(北京) 机械与储运工程学院，北京 102249 2 中国石油管道科技研究中心，廊坊 065000

针对常用含腐蚀缺陷管道可靠性评价方法较多且适用性和保守程度不一的问题，对一种新的腐蚀管道失效压力预测公式进行误差分析，证实该模型比常用ASME B31G-2009、DNV RP-F101、PCORRC、CSA Z662-07、RAM PIPE REQUAL、NETTO ET AL等方法具有更高的预测精度。基于该预测公式，建立腐蚀缺陷管道失效极限状态方程，采用蒙特卡洛模拟方法计算了管道失效概率，并与其它常用工业标准与方法的计算结果进行对比分析。结果表明：基于该模型计算的失效概率与PCORRC, ASME B31G- 2009以及DNV RP-F101 许用应力法计算结果相近，进行可靠性分析具有较好的适用性。通过参数敏感度分析得出结论：管道壁厚、直径、内压、径向腐蚀速率等随机变量的变异系数的增大对管道失效概率的影响具有双向性，可能造成管道失效概率的增加，也可能导致管道失效概率的减小。指出径向腐蚀速率和管道内压是影响腐蚀管道可靠性的两个关键因素。

腐蚀管道；内检测；失效概率；极限状态；敏感性

0 引言

腐蚀缺陷是油气管道常见的一种缺陷，对管道的结构完整性构成严重威胁，容易造成管道的破裂和泄漏，给社会经济及环境造成重大损失。基于可靠性的分析方法能够综合考虑各种因素的不确定性对管道失效概率的影响，例如缺陷尺寸、材料属性、操作工况的不确定性等，这些不确定性对腐蚀管道的剩余寿命产生重要影响[1]。因此，采用可靠性分析方法对管道上的腐蚀缺陷进行管理成为国内外普遍采用的管理模式，对管道的安全运营具有重要意义[2-4]。

目前，工程上有许多预测含腐蚀缺陷管道失效压力的公式。例如ASME B31G-2009[5]、DNV RP-F101 MOP[6]、RSTRENG[7]、PCORRC[8]、CSA Z662-07[9]、RAM PIPE REQUAL[10]、NETTO ET AL[11]。国内外学者基于不同的预测公式针对腐蚀管道的可靠性进行了大量研究[12-14]。但是这些研究都只是基于其中某一种预测公式进行可靠性分析。由于不同预测公式的保守程度不一样，基于这些预测模型的可靠性分析结果进行维修和检测决策，就可能会造成检测、维修费用浪费或者管道安全问题。因此，有必要对不同的预测模型进行可靠性分析的适用性进行分析。尽量采用精度高的预测公式进行腐蚀管道可靠性分析。

针对上述问题，笔者基于一种新的预测精度较高的爆破压力预测公式，建立了腐蚀管道失效极限状态方程，对两次管道内检测腐蚀数据进行匹配，获取了管道腐蚀速率，采用蒙特卡洛模拟方法计算了腐蚀管道的失效概率。并与其它常用模型计算的失效概率进行对比，分析各模型的保守性。最后分析了该新的模型中各随机变量对管道失效概率的敏感性。

1 含腐蚀缺陷管道失效压力预测模型

有限元方法是对含缺陷管道剩余强度进行精确评定的有效方法，它与全尺寸管道爆破实验相结合，是改进含缺陷管道剩余强度评定方法的有效手段。帅健在文献[15]中采用20节点的六面体等参单元，建立含腐蚀缺陷管道有限元数值模型。认为当腐蚀区域的应力达到材料后屈服终点时，即腐蚀区域最小等效应力达到材料的抗拉强度UTS时，管道发生失效。采用该有限元模型计算了多组爆破试验腐蚀管道的失效压力。结果表明，有限元计算得到的管道失效压力非常接近于爆破试验结果，其计算最小相对误差为0.49%，最大相对误差为10.78%，验证了有限元模型的准确性。在此基础上计算了39组不同情况下的腐蚀管道爆破压力，并进行拟合得到如下公式(1)，以下讨论简称“CUP”模型。

式中，D为管道直径，t为管道壁厚，L为腐蚀缺陷长度，d为腐蚀缺陷深度，w为腐蚀缺陷宽度，UTS为管材拉伸强度。

为了验证公式(1)的准确性与可靠性。采用公式(1)以及常见ASME B31G -2009、DNV RP-F101许用应力法、PCORRC、RSTRENG、CSA Z662-07、RAM PIPE、NETTO ET AL等方法分别计算所收集14组爆破试验中腐蚀缺陷管道的失效压力，并与实验爆破压力值[16-21](见表1)进行对比，误差分析结果如图1所示。可以看出，该预测公式比ASME B31G-2009、DNV RP-F101、PCORRC、CSA Z662-07、Netto ETAL以及RAM PIPE等公式具有更高的预测精度。

2 腐蚀管道可靠性分析原理

2.1 腐蚀速率模型

腐蚀缺陷的长度、宽度、深度取决于缺陷的腐蚀速率。随着时间的增长，腐蚀缺陷尺寸逐渐增长，现场经验表明，管道最初的腐蚀增长速率较大，随后逐渐下降，最终有可能达到稳定。由于最初增长阶段时间较短，且缺陷较小，对管道的完整性威胁不大，所以主要考虑的就是稳态的增长速率，假设该状态下增长速率为线性，则：

表1 含腐蚀缺陷管道全尺寸爆破试验数据Table 1 Full scale burst experiment datas of corroded pipelines

图1 几种方法的计算误差对比Fig. 1 The comparison of calculation errors using fi ve method

2.2 极限状态方程

失效概率为

依照上面的理论，含腐蚀缺陷管道极限状态函数可表示为

其中，LP为腐蚀缺陷管道的失效压力；0P是管道的工作压力；当时，管道失效；管道不发生失效；若腐蚀缺陷管道失效压力LP采用公式(1)进行计算，结合式(6)～(8)，则含腐蚀缺陷管道极限状态方程为(12)所示。同理可以得到基于其他预测公式的腐蚀管道极限状态方程。

2.3 蒙特卡洛模

工程中直接求解式(12)来获取腐蚀管道的失效概率比较困难，采用蒙特卡洛模拟方法可以有效的解决这个复杂的概率求解问题。蒙特卡洛模拟方法是通过随机变量的统计试验，求得数学物理和工程技术问题的近似解的数值方法，也称为统计试验法或随机模拟法。它只需要已知随机变量的概率分布即可，而不会因为状态函数的非线性、状态函数中随机变量是否为非正态以及随机变量之间是否相关而发生求解困难。根据随机变量的分布函数，随机生成109样本数据，代入极限状态方程，利用计算机模拟计算可以得到腐蚀管道的失效概率。采用蒙特卡洛模拟方法计算腐蚀管道失效概率具体步骤流程如图2所示。

3 实例应用与分析

3.1 参数概率分布统计

对国内某段(库尔勒-马兰)原油管道进行可靠性评估。该段管线长153 km，管材为API 5L X65，管径610 mm，壁厚7.1 mm，工作压力8 MPa。为了解该管线腐蚀状况，分别于2011年12月1-4号和2013年5月3-6号采用漏磁内检测器(两次内检测均为同一服务商，采用的内检测器型号、性能一致)对该段管线进行了两次内检测，检测时间间隔519天，获取了该段管道腐蚀缺陷数据，其中2011年检测到腐蚀缺陷823处，2013年检测到腐蚀缺陷1 358处。通过漏磁检测信号识别得到了每个腐蚀缺陷点的位置里程、深度、宽度、长度的具体尺寸及时钟方位等信息。对两次内检测腐蚀缺陷数据进行对比、匹配通常比较困难，这是由于检测时间间隔内会有新的腐蚀缺陷产生，并且腐蚀深度较小时，漏磁信号比较微弱，造成匹配困难。因此设置阀值，对腐蚀缺陷深度阀值5%wt以上的腐蚀缺陷进行匹配，获取了两次检测相同的缺陷，剔除未能匹配的缺陷，共匹配成功腐蚀缺陷305处，基本可以满足该段管线的腐蚀速率统计需求。通过2.1节的腐蚀速率计算模型就可以对匹配成功的腐蚀缺陷点进行腐蚀速率计算。对现场采集的该管线压力波动数据以及内检测到的腐蚀缺陷深度、长度、宽度数据以及匹配成功的腐蚀缺陷点的腐蚀速率进行统计分析，得到了各随机变量的分布类型和统计参数如表2所示。

图2 蒙特卡洛模拟计算方法流程图Fig. 2 The fl ow chart depicts the calculation procedure of Monte Carlo method

3.2 失效概率计算

采用上述理论计算该段管线的失效概率，为保证计算精度，蒙特卡洛模拟次数取109。各模型的失效概率计算结果如图3、4所示。从图中可以看出，各种模型计算出的该管道失效概率随着时间增长，逐渐增大，但是各模型计算结果存在差异。可将各模型计算结果大致分为三个等级：RAM PIPE REQUAL和CSA Z662-2007计算的失效概率相对其他模型明显偏大；然而Netto计算的失效概率明显偏小；其他模型ASME B31G-2009、RSATRENG、CUP、PCORRC和DNV RP-F101许用应力法计算结果为一个等级。其中，基于 CUP模型所计算的失效概率与ASME B31G-2009以及PCORRC计算结果相近。

为分析各个模型预测的失效压力与其计算的失效概率之间的关系，采用剩余强度定量分析方法分析不同腐蚀深度下各模型所计算失效压力的变化规律，如图5、6所示。可以看出，随着腐蚀缺陷深度的增长，各个模型的保守程度发生变化。同一腐蚀缺陷深度下，RAM PIPE REQUAL 和 CSA Z662-2007预测的失效压力低，所以其计算的失效概率较大。Netto预测的失效压力比较大，导致其计算的失效概率比其他模型低。从图6可以看出，CUP、 DNV RP-F101许用应力法、RSTRENG 以及 PCORRC方法所预测的失效压力较为相近，因此，其计算得到的失效概率也较为相近。其中，PCORRC 和 DNV RP-F101许用应力法计算的失效压力更为接近。

3.3 参数敏感性分析

参数敏感性分析的目的是研究随机变量的分散性和随机性对结构失效概率的影响，找出影响结构安全的关键变量，并在工程中尽可能减小关键变量的分散性和随机性，即减小关键参量的变异系数(如尽可能的减小压力波动等)，以提高结构的可靠性。参数敏感性分析可以采用变异系数来表征，它是随机变量不确定性程度、分散性程度的度量，其定义表达式为[23]

式中，u为平均值，s为标准差。固定平均值，改变标准偏差，就可以计算不同变异系数下腐蚀管道的失效概率。基于此，以下对基于CUP预测公式的可靠性模型中各随机变量对管道失效概率的敏感性进行分析，计算结果如图7-16。

表2 变量分布类型及参数Table 2 Input probability distributions and parameter

图3 不同计算模型的失效概率计算结果(20年)Fig. 3 Failure probability calculated from different models (20 years)

图4 不同计算模型的失效概率计算结果(40年)Fig. 4 Failure probability calculated from different models (40 years)

图7-8分别是随机变量D，t的变异系数cov( D)、cov( t)变化对管道逐年失效概率的影响，从图中可以看出，变异系数cov( D)、cov( t)在0.01～0.05之间变化对管道失效概率的影响较小。从最近一次检测时间30年内，失效概率随着cov( D)、cov( t)的增大而增大，超过30年以后管道失效概率随着cov( D)、cov( t)的增大反而减小。同时可以看出，壁厚比管径对腐蚀管道失效概率的影响更大。

图5 各模型计算的失效压力随缺陷深度变化规律Fig. 5 The law of burst pressures to defect depth calculated by different models

图6 CUP, DNV, RSTRENG, PCORRC计算的失效压力随缺陷深度变化规律Fig. 6 Comparison of burst pressure of model CUP, DNV, RSTRENG, PCORRC

图7 管径D变异系数变化对管道失效概率的影响Fig. 7 Effect of the cov(D) on the failure probability

图8 壁厚t变异系数变化对管道失效概率的影响Fig. 8 Effect of the cov(t) on the failure probability

图9是管道内压变异系数cov( P0)变化对腐蚀管道失效概率的影响。当时间间隔小于33年时，失效概率随着cov( P0)增大明显增大，当超过33年之后开始呈负增长趋势。从图中可以明显看出，操作压力的波动对腐蚀管道的失效概率影响十分显著。

图10是管材拉伸强度变异系数cov( UTS)变化对腐蚀管道逐年失效概率的影响。管道失效概率随着cov( UTS)的增大而增大，这种趋势一直持续到第40年。

图11～13为第一次内检测腐蚀缺陷深度、长度、宽度即初始缺陷尺寸的变异系数cov( d0)、cov( L0)、cov( w0)对管道失效概率的影响。可以看到，cov( d0)、cov( L0)、cov( w0)的变化几乎对管道的失效概率没有影响，这是因为该段管线初次(2011年)内检测缺陷深度普遍较浅(平均深度0.559 mm)。

图14～16显示了腐蚀速率的变异系数cov( Vd)、cov( Vl)、cov( Vw)对腐蚀管道失效概率的影响。从图14中可以看到，在管道的失效概率小于0.5的范围内，概率随着cov( Vd)的增大而增大。当管道失效概率为0.2时，若腐蚀深度速率变异系数从0.1增加到0.4，则管道的服役寿命将会减少4年。在管道的失效概率大于0.5的范围内，概率随着cov( Vd)的增大反而减小。这是因为当失效概率小于0.5时，随机变量有较小的可能性使得极限状态函数g( x)的值落在失效区域(g( x)＜0)内。增加cov( Vd)的值，却增加了g( x)＜0的可能性，从而提高了腐蚀管道的失效概率。然而，当腐蚀管道失效概率大于0.5时，管道失效概率随cov( Vd)的变化呈现与上面相反的变化规律。这是因为失效概率大于0.5，意味着随机变量有较大可能性导致管道失效，增大cov( Vd)的值，反而减小了g( x)＜0的可能性，从而减小了管道的失效概率。从图15中可以看到，距离最近一次检测11年之内，管道的失效概率对轴向腐蚀速率不敏感。11年之后，管道的失效概率随着cov( Vl)增大而增大。图16显示，cov( Vw)对管道失效概率几乎没有影响，这是由于腐蚀缺陷宽度对管道承压能力几乎没有影响的原因。综合可以得知，三个(径向、环向、轴向)腐蚀速率中，径向腐蚀速率对腐蚀管道失效概率的影响最大。

图9 操作压力0P变异系数变化对管道失效概率的影响Fig. 9 Effect of the cov(P0) on the failure probability

图10 拉伸强度UTS变异系数变化对管道失效概率的影响Fig. 10 Effect of the cov(UTS) on the failure probability

图11 初始腐蚀深度d0变异系数变化对管道失效概率的影响Fig. 11 Effect of the cov(d0) on the failure probability

图12 初始腐蚀长度L0变异系数变化对管道失效概率的影响Fig. 12 Effect of the cov(L0) on the failure probability

上面的分析是通过分别改变变量d0，L0，w0，D，t，P0，UTS，Vd，Vl和Vw的中的某个变量标准偏差，固定其他变量来进行分析的。通过分析可以得知，腐蚀管道的失效概率对变量t，D，U TS，P0，Vd和 Vl比较敏感。因此，准确确定各随机变量分布参数是对腐蚀管道进行可靠性分析的重要前提。

图13 初始腐蚀宽度w0变异系数变化对管道失效概率的影响Fig. 13 Effect of the cov(w0) on the failure probability

图14 腐蚀深度速率变异系数变化对管道失效概率的影响Fig. 14 Effect of the cov(Vd) on the failure probability

图15 腐蚀长度速率变异系数变化对管道失效概率的影响Fig. 15 Effect of the cov(Vl) on the failure probability

图16 腐蚀宽度速率变异系数变化对管道失效概率的影响Fig. 16 Effect of the cov(Vw) on the failure probability

4 结论

(1) 误差分析显示CUP方法预测腐蚀管道失效压力具有较好的预测精度。基于此预测公式，建立了腐蚀管道失效的极限状态方程，计算了腐蚀管道的失效概率，并与常用的预测模型进行对比，对各种方法的适用性和保守性进行了论证和分析，将保守程度分为了三个等级，分析了保守程度差异性的原因，认为该新的预测公式用于腐蚀管道可靠性分析具有较好的适用性。可以为腐蚀管道可靠性分析提供参考和选择。

(2) cov( D)、cov( t)、cov( P0)、cov( Vd)的变化对管道失效概率的影响具有双向性，变异系数的增大可能造失效概率的增大或减小。管道失效概率对cov( P0)和cov( Vd)最为敏感，因此，在管道尺寸及材料不可变的情况下，可采取措施减小管道径向腐蚀速率和内压波动幅度以提高管道的安全性。另外，相对而言，cov( t)比cov( D)对腐蚀管道失效概率的影响更大一些。cov( Vl)、cov( Vw)对腐蚀管道的失效概率的影响较小。其中，腐蚀缺陷宽度的影响几乎可以忽略，这与常用的、多数学者开发腐蚀管道失效压力公式未考虑缺陷宽度影响是相一致的。

(3) 本文针对腐蚀管道的可靠性分析方法是基于各随机变量独立性假设。实际上，有些随机变量之间是存在一定联系的，例如缺陷深度、长度、壁厚与管材强度，壁厚与直径之间等。因此，综合分析各随机变量之间的相关性对腐蚀管道失效概率的影响是下一步研究工作的重点。

[1] CALEYO F, GONZAL J L, HALLEN J M. A study on the reliability assessment methodology for pipelines with active corrosion defects [C]. International Journal of Pressure Vessels and Piping, 2002, 79(1)： 77-86.

[2] PAPADRAKAKIS M, LAGAROS N D. Reliability-based structural optimization using neural networks and Monte Carlo simulation [J]. Computer Methods in Applied Mechanics and Engineering. 2002, 191(32)： 349-357.

[3] STEPHENS M J, AND VAN ROODSELAAR A. 2008. Developments in Reliability-based Corrosion Management and the Signif i cance of In-line Inspection Uncertainties [C]. Proceedings of IPC 2008, 7th International Pipeline Conference, Paper No. IPC2008-64384, Calgary, Alberta, September 29-October 3, 2008.

[4] 罗金恒, 赵新伟, 白真权, 路民旭. 输油管道腐蚀剩余寿命预测[J]. 石油机械, 2000. 28(2)： 29-32. [LUO J H, ZHAO X W, BAI Z Q, LU M X. Residual Life Prediction Method of Corroded Buried Oil Pipeline [J]. China Petroleum Machinery, 2000. 28(2)： 29-32.]

[5] ASME B31G-2009. Manual for determining the remaining strength of corroded pipelines[S]. American Society of Mechanical Engineers, 2009.

[6] Recommended Practice DNV RP-F101. Corroded pipelines[S], 2004.

[7] KIEFNER J F, VIETH P H. A Modif i ed criterion for evaluation the remaining strength of corroded pipe [C]. Final report on project PR 3-805 to the Pipeline Research Committee of the American Gas Association, 1989.

[8] DENNY R, STEPHENS B N. Development of alternative criterion for residual strength of corrosion defects in moderate-to-high-toughness pipe [C]. Battelle Energy Products Division, 505 King Ave Columbus, Ohio 43201-2693. In： 2000 international pipeline conference, vol. 2. ASME 2000.

[9] CSA Z662-07. Limit state equation for burst of large leaks and rupture for corrosion defect [S]. Oil and Gas Pipeline Systems, Canadian Standards Association, 2007, 554-555.

[10] BEA R, XU T. Corrosion effects on burst pressures RAM pipe requal[R]. University of California at Berkeley, 1999. Pipeline Requalif ication Guidelines Project Report 1, 103-104.

[11] NETTO T A, FERRAZ U S, ESTEFAN S F. The effect of corrosion defect on the burst pressure of the pipeline[J]. Journal of Constructional Steel Research 2005(61)： 185-204.

[12] SIKDER H, FAISAL K, SHAWN K. Probability assessment of burst limit state due to internal corrosion[J]. International Journal of Pressure Vessels and Piping, 2012, 89(2)： 48-58.

[13] NESSIM M A, ZHOU W. Reliability based design and assessment for location-specif i c failure threats[C]. Proceedings of IPC 2006, International Pipeline Conference, Paper No. IPC06-10097, Calgary, October.

[14] MCCALLUM M, FORD G, FRANCIS A. Extend the in-line inspection interval for a gas pipeline using direct assessment[C]. 7th International Pipeline Conference, Canada, 2008.

[15] 帅健, 张春娥, 陈福来. 非线性有限元法用于腐蚀管道失效压力预测[J]. 石油学报, 2008(06)： 933-937. [SHUAI J, ZHANG C E, CHENG F L. Prediction of failure pressure in corroded pipelines based on non-linear fi nite element analysis [J]. Acta Petrolei Sinica, 2008(06)： 933-937.]

[16] CHOUCHAOUI B A, PICK R J, YOST D B. Burst pressure predictions of line pipe containing single corrosion pits using fi nite element method [J]．Pipeline-technology．ASME OMAE 1992, Vol. V-A： 203-210.

[17] CHOUCHAOUI B A, PICK R J. Behavior of circumferentially aligned corrosion pits [J]．Journal of pressure vessel and piping. 1994(57)： 187-200.

[18] BJORNOY O H, OLE R, STEIN F. Residual strength of dented Pipelines[C], DNV Test Results. Proceedings of the 10th (2000) International Offshore and Polar Engineering Conference Seattle, USA, May 28-June 2, 2000.

[19] ADILSON C, BENJAMIN R D, VIEIRA J L F, et al. Burst tests on pipeline with long external corrosion [C]. International pipeline conference, 2000.

[20] DARUO B N J, ASILSON C B. Finite element models for the prediction of the failure pressure of pipeline with long corrosion defects [C]. Proceedings of IPC’ 02, 4th International Pipeline Conference. Sep. 29-Oct. 3, 2002, Canada.

[21] MOK D R B, PICK R. J, GLOVER A G. Bursting of line pipe with long external corrosion[J]. Journal of pressure vessel and piping, 1991 (46)： 195-216.

[22] 赵新伟, 罗金恒. 油气管道完整性评价技术[M]. 陕西出版集团, 2010. [ZHAO X W, LUO J H. Integrity assessment technology of oil and gas pipelines[M]. Shanxi Publishing Group, 2010.]

[23] ZHOU W. System reliability of corroding pipelines[J]. International Journal of Pressure Vessels and Piping, 2010： 587-595.

[24] QIAN G, NIFFENEGGER M, LI S X. Probabilistic analysis of pipelines with corrosion defects by using FITNET FFS procedure[J]. Journal of corrosion Science, 2011(53)： 855-861.

Research on the reliability methods of corroded pipeline

SHUAI Yi1, SHUAI Jian1, Liu Chaoyang2
1 College of Mechanical and Transportation Engineering, China University of Petroleum-Beijing, Beijing 102249, China 2 Petrochina pipeline R&D center, Langfang 065000, China

There are many reliability evaluation methods for the corroded pipelines. But there is no unif i ed standard and the versatility is poor. In order to solve the problem, this paper studied a prediction formula mentioned in one literature. The error analysis shows that it has higher precision than the models of ASME B31G-2009, DNV RP-F101, PCORRC, CSA Z662-07, RAM PIPE REQUAL and NETTO ET AL. Based on the structural reliability theory, a limit-state function is def i ned. Then the probabilistic damage model is formulated and numerically calculated by using Monte Carlo Simulation (MCS). The failure probability of a corroded pipeline was calculated, and compared with the candidate standards. The results show that the failure probability calculated by this model is similar to that calculated by the allowable stress methods of PCORRC, ASME B31G-2009 and DNV RP-F101. This indicated that the new prediction formula can provide a choice for the reliability analysis of corroded pipelines. Through the parameter sensitivity analysis, it is found that an increase in the variation coeff i cient of random variables as wall thickness, diameter, internal pressure and radial corrosion rate may lead to an increase or a decrease in the failure probability. It is pointed out that the radial corrosion rate and internal pressure are the two key factors affecting the reliability of corroded pipelines.

corroded pipeline; in-line inspection; failure probability; limit state; sensitivity

10.3969/j.issn.2096-1693.2017.02.027

(编辑 马桂霞)

帅义, 帅健, 刘朝阳. 腐蚀管道可靠性评价方法研究. 石油科学通报, 2017, 02： 288-297

SHUAI Yi, SHUAI Jian, Liu Chaoyang. Research on the reliability methods of corroded pipeline. Petroleum Science Bulletin, 2017, 02：288-297. doi： 10.3969/j.issn.2096-1693.2017.02.027

*通信作者, shuaij@cup.edu.cn

2016-11-06

“十二五”国家科技支撑计划课题(2015BAK16B00)资助