2017, Vol. 60
2. 海底科学与探测技术教育部重点实验室, 山东青岛 266100
2. Key Lab of Submarine Geosciences and Prospecting Techniques of Ministry of Education, Ocean University of China, Qingdao 266100, China
近年来，海洋可控源电磁(Controlled-Source ElectroMagnetic method，CSEM)作为海洋油气资源勘探的新方法取得了突飞猛进的发展(Constable and Weiss, 2006).频率域海洋CSEM方法通常使用拖曳在离海底上方几十米处的电偶极子作为发射源，并在被拖曳过程中向位于海底的电磁采集站发射低频电磁信号(0.1~10 Hz)，采集站接收到的信号主要是来自海底地层的电磁感应信号，其振幅和相位依赖于海底介质的电阻率.由于海洋CSEM法对海底高阻薄层有较强的分辨能力，海洋CSEM法已经成功应用到海洋油气勘探(Kong et al., 2002)以及天然气水合物资源勘查评价(Boswell et al., 2012; Shedd et al., 2012; Weitemeyer et al., 2006).
在海洋CSEM法勘探中，由于水平电偶极子源(HED)所产生的场在海底环境中的复杂性，定性解释和反演算法在海洋CSEM资料解释过程中显得尤为重要.借鉴地震勘探中的共中心点理论(李振春等, 2000)，由电磁勘探的发射源和接收站组成共中心域(Common Mid-Point, CMP)(Buland and KolbjØrnsen, 2012; Crepaldi et al., 2011; Mittet et al., 2008)，将原观测系统中得到的电(磁)场响应转换至不同CMP域单元中分析，该方法对高阻储层有较高的分辨能力，还可将二维地电结构的电磁响应在CMP单元内表达为一维模型的响应(Crepaldi et al., 2011).传统的定性解释是通过归一化异常或者MVO/PVO (Magnitude/Phase Versus Offset)曲线来完成(裴建新等, 2015a)，其优势在于反映衰减规律，而CMP域的响应分析更适于在异常区域的横向上直观地分辨高阻储层.在地球物理的解释方法中，海洋CSEM法一维反演已趋于成熟(Key, 2009; Ray and Key, 2012; 刘颖等, 2013; 王涛, 2015；罗鸣等，2016)，但是对于复杂的海底地电结构，传统的一维反演往往不能够满足解释的需要；然而，海洋CSEM二维(Abubakar et al., 2008)和三维(Newman and Alumbaugh, 1997; Sasaki, 2013)反演算法大多基于有限单元法(FE)或有限差分法(FD)的正演数值模拟来实现的，较于一维算法会占用相对较大的内存和消耗更多的运算时间，影响其在实际工作中现场处理、分析的效率.鉴于此，本文实现了一维频率域海洋可控源电磁资料高斯-牛顿反演算法，提出将二维海洋CSEM响应转换成一维CMP响应的方法，分析了不同天然气水合物模型CMP域的水平电场和归一化响应的特征，运用一维高斯-牛顿反演对不同CMP单元进行反演解释，成功实现了二维结构的CMP域一维快速反演算法.
2 CMP域反演 2.1 CMP域数据转换方法海洋CSEM法CMP域反演方法是借鉴地震勘探共中心点(CMP)理论来实现的.假设实际探测目标如图 1所示，海水深度为1000 m，在海底下方一定深度处存在一个厚度为100 m、宽度为4000 m (横向范围为-2000~2000 m)的天然气水合物储层.假设发射源拖曳于海底上方50 m的位置，沿测线方向从-10000~10000 m等间距地发射201次可控源电磁信号，11个海洋CSEM接收站均匀布置于测线(y轴方向)-5000~5000 m范围内.
|
图 1 天然气水合物模型CMP模拟图 Fig. 1 CMP schematic diagram of hydrate model |
将接收站(图 1中白色倒三角)及相邻接收站中点所在位置作为CMP域的分割点，11个接收站共等分为20个单元(图 1红线)，相邻两个红线之间区域为一个CMP单元.由发射源坐标ys和接收站坐标yr计算得到CMP域共中心点坐标yc和收发距yh：
|
(1) |
单一接收站接收到单一发射源的电磁信号为一个数据点，通过(1)式可计算该数据点所对应的CMP位置yc和收发距yh，频率域电磁信号由此也可转换到CMP域中.图 2为CMP域yc中的数据点与收发距yh的关系示意图，横坐标为CMP位置，纵坐标为收发距，黑点表示不同发射源与接收站的数据点.假设有一个数据(发射源坐标ys=800 m，接收站坐标yr=100 m)，其CMP位置yc=450 m，介于0 m与500 m之间，处于第11个CMP单元(yc=0~500 m，从左向右数)中，该数据点将作为第11个CMP单元(yh=700 m)的数据进行CMP单元的反演解释.在CMP域内进行单元反演时，假定单元内的电阻率在横向上是不变的，而只是在垂直方向上发生变化.由此，传统的海洋CSEM观测系统记录的电磁数据可通过等式(1)转换至CMP域中.基于2.5D海洋CSEM自适应有限元算法(Li and Key, 2007)，可将二维地电结构的海洋CSEM正演响应转化为CMP域响应.
|
图 2 CMP数据示意图 Fig. 2 Schematic diagram of CMP data |
海洋CSEM法CMP域反演是通过将二维模型响应转换为CMP域响应，再对每个CMP单元内分别进行一维反演.在此，我们采用高斯-牛顿法进行反演(罗鸣等，2016)，其目标函数为
|
(2) |
其中，m为模型参数，是由反演模型的各层的电阻率的对数值组成的向量，m=(log10ρ1, log10ρ2, …, log10ρM)T，下标M为反演模型的层数.ϕd和ϕm分别为观测数据与理论模型响应的拟合差和反演模型约束.λ为正则化因子，用于平衡数据拟合差ϕd与反演模型约束ϕm之间的权重.可将ϕd写为(3)式：
|
(3) |
其中，‖ · ‖为标准差，e=d-F(m)为残差向量；d=(d1, d2, …, dN)T为观测数据；F(m)为模型正演响应函数；W d=diag (1/σ1, 1/σ2, …, 1/σN)为数据加权矩阵(N×N的对角阵，N为观测数据的个数)，且σi为第i个数据的标准差.ϕm可表示为如下形式：
|
(4) |
其中，R为模型加权函数，本文中均使用最小支撑函数(Minimum support，MS)(Zhdanov, 2002).
2.3 目标函数极小值问题假设在第k次迭代，通过求取目标函数的极小值，得到搜索方向pk，则可得到第k+1次迭代的模型参量为
|
(5) |
其中，目标函数的梯度(g k)与海森矩阵(H k)具有如下形式：
|
(6) |
|
(7) |
式(6)、式(7)中J k为雅克比矩阵，写成元素形式为
|
(8) |
在反演迭代过程中，反演迭代模型会逐步向真实模型收敛，其目标函数拟合差ψ：
|
(9) |
也逐渐收敛于1.0，其中，δi为第i个数据的标准差.
2.4 正则化因子正则化因子λ的选择方式是能否获得合理反演结果的关键.在地球物理反演方法中，许多正则化因子的选择方式已经被提出(Constable et al., 1987; Hansen, 1992; Newman and Alumbaugh, 1997).本文使用一种衰减近似的方法选择正则化因子，我们选择矩阵乘积[(W d J)T(W d J)]中各行的总和的最大值作为正则化因子(Chen and Alumbaugh, 2009)，假设在反演第i次迭代，正则化因子λi可写作如下形式：
|
(10) |
其中，amj为反演第i次迭代时矩阵乘积[(W d J)T(W d J)]的元素，为正则化衰减因子(常设为一个常数).
3 模型算例 3.1 CMP域电磁响应为了对比分析天然气水合物储层在不同埋藏埋深及电阻率情况下的海洋CSEM响应及CMP域电磁响应的差别，在图 1所示观测系统的基础上设定图 3所示的四个天然气水合物模型：图 3a，天然气水合物储集体埋深为200 m，电阻率10 Ωm，水平位置-2000~2000 m；图 3b，天然气水合物储集体埋深为500 m，其他参数与模型(a)相同；图 3c，天然气水合物储集体电阻率5 Ωm，其他参数与模型(a)相同；图 3d，天然气水合物储集体埋深为500 m，其他参数与模型(c)相同.
|
图 3 天然气水合物模型图 (a)天然气水合物埋深200 m，电阻率10 Ωm；(b)天然气水合物埋深500 m，电阻率10 Ωm；(c)天然气水合物埋深200 m，电阻率5 Ωm；(d)天然气水合物埋深500 m，电阻率5 Ωm. Fig. 3 Hydrate model (a) Hydrate with 200 m depth and 10 Ωm; (b) Hydrate with 500 m depth and 10 Ωm; (c) Hydrate with 200 m depth and 5 Ωm; (d) Hydrate with 500 m depth and 5 Ωm. |
正演模拟时设定发射频率为5 Hz，在轴向观测模式下(inline geometry，发射源方向与测线方向相同)分别计算如图 3所示的天然气水合物模型的海洋CSEM响应.然后，通过CMP域数据的转换方法，我们将二维海洋CSEM数据转换成CMP单元内的数据进行讨论和处理，水平电场Ey分量振幅响应对比如图 4所示：蓝点为背景模型的电磁响应，红点为储层在横向上无限延伸的电磁响应(一维模型响应)，绿点为第11个CMP单元(yc=0~500 m，模型中心区域单元)的电磁响应，黑点为天然气水合物模型第6个接收站(yr=0 m，模型中心区域的接收站)的电磁响应.由图 4可以看出，所有含天然气水合物储层的水平电场响应均大于背景模型的响应.需指出，二维模型的CMP单元响应均大于相同模型的单个接收站的响应，并且趋近于一维模型的响应，表现出对目标体较高的分辨能力.
|
图 4 天然气水合物模型水平电场响应对比 Fig. 4 The amplitude of horizontal electric field for hydrate model |
在海洋CSEM勘探中，位于海底表面的接收站所接收到的电磁场信号包含有地下地电结构信息以及“空气波”(Chen and Alumbaugh, 2009)和场源的影响.为了能够突出地下地电结构所反映的信息，通常我们对接收到电磁响应做归一化处理.归一化异常和相位差的计算公式为
|
(11) |
其中，Eobs和φobs分别为含有油气储层的电场振幅和相位，E0和φ0分别为背景模型(无油气储层)的电场振幅和相位.图 5为四个天然气水合物模型(图 3)水平电场响应基于背景模型归一化的结果.由图 5可以看出，不论是对于埋深浅、电阻率大的天然气水合物模型(电磁信号更强)，还是埋深较大、电阻率较小的模型(电磁信号更微弱)，CMP域处理方式均能够在CMP单元内较好地将二维模型的响应结果表达为信号较强的一维模型电磁响应，这表明在CMP域转换能够提高分辨能力，在其基础上进一步实现海洋CSEM快速反演是可行的.
|
图 5 天然气水合物模型水平电场归一化对比 Fig. 5 The normalized field of horizontal electric field for hydrate model |
为了分析CMP域在资料处理上的优势，在整个CMP域内，将inline模式下获得的四个天然气水合物模型(图 3)响应(水平电场Ey)的振幅和相位资料转换至CMP域(裴建新等, 2015b)，得到图 6和图 7，再将图 6所示的CMP域Ey响应与背景模型响应做归一化处理得到如图 8所示的结果(图中虚线为天然气水合物储层的横向范围).由图 6-图 8(各图中四个结果分别对应于图 3中的四个模型的电磁响应)可见：
|
图 6 天然气水合物模型CMP域水平电场对比图(虚线为天然气水合物储层的横向范围) Fig. 6 The horizontal electric field of hydrate model in CMP domain |
|
图 7 天然气水合物模型CMP域水平电场相位对比图(虚线为天然气水合物储层的横向范围) Fig. 7 The phase of horizontal electric field for hydrate model in CMP domain |
|
图 8 天然气水合物模型CMP域水平电场归一化对比图(虚线为天然气水合物储层的横向范围) Fig. 8 The normalized field of horizontal electric field for hydrate model in CMP domain |
(1) CMP域水平电场Ey的振幅和相位响应均能够反映出高阻储集体的存在.异常响应明显的区域横向上正好对应于高阻储集体模型所在位置，其中异常响应最强的区域也是模型的中心区域，表明基于CMP域的海洋CSEM定性解释在横向上能够很好地分辨高阻储层；
(2)通过归一化振幅结果能够更好地体现高阻储集体的存在，并且在横向-2000~2000 m范围内能够很好地吻合储层的横向位置.对于埋深为200 m，电阻率为10 Ωm的模型中心区域的归一化异常达到800，随着埋深增大以及天然气水合物模型电阻率的减小，由于异常响应减小，归一化异常也随之减小；
需指出，收发距小于1000 m的情况下，由于场源的影响，其振幅和相位响应并不能真实反映异常的存在.
3.2 天然气水合物模型CMP域反演由于目前尚缺少海洋CSEM实测资料，在此借鉴其他地球物理方法实际资料建立海洋CSEM理论模型进行研究.南海地区关于天然气水合物的地震勘探及测井资料显示其埋藏深度较浅，储集体与围岩的电阻率差异小(苏正等, 2011)，这将造成其异常响应相对微弱，也给天然气水合物的海洋CSEM二维反演问题带来了困难.若直接将常用的一维反演方法应用于二维资料的解释，将很难真实地反映出地下实际的地电结构.鉴于此，采用折中的方案，将二维异常体海洋CSEM响应在CMP单元内表达为一维模型的响应，从而可实现模型的CMP域快速反演.
对于图 3所示的天然气水合物模型，观测系统与CMP单元剖分如图 1所示(包含20个CMP单元).建立反演初始模型为：空气-海水-围岩的三层模型，其中，空气层及海水层的电阻率和厚度已知且不参与反演，围岩电阻率的初始值为1 Ωm.由于天然气水合物通常埋藏较浅，合成数据时采用4个较高的发射频率(10, 15, 20, 30 Hz).将正演得到的水平电场Ey分量的振幅和相位加入2%的随机高斯噪声用于反演计算.图 9为CMP反演方法得到的各个水合物模型对应的反演结果，可见：当水合物真实电阻率为10 Ωm时(图 9a、9b)，其反演异常体中心位置电阻率值能够达到7~8 Ωm，当水合物真实电阻率为5 Ωm时(图 9c、9d)，其反演异常体中心位置电阻率值能够达到4~5 Ωm.以上分析可知：因为该方法是在二维模型的CMP域划分基础上实现一维计算，其实质为一维反演，其计算效率与一维反演效率相对应，较于二维反演的计算效率具有显著优势；反演结果能够明确反映出天然气水合物异常体的存在，反演得到异常区域的水平和垂直位置与真实模型(图 9中白色线框)对应较好.
|
图 9 天然气水合物模型CMP域反演剖面图 Fig. 9 The CMP inversion profile of hydrate model |
为了进一步验证反演算法的效果，我们引入相对复杂的二维模型(刘颖, 2014)(图 10).设定海水电阻率为0.3 Ωm，水深1000 m；海底下方覆盖层厚度为1000 m，电阻率为10 Ωm；其下有一厚度为3000 m、电阻率为200 Ωm的沉积层，该沉积层受低阻推覆体(电阻率10 Ωm)侵入，在分界处形成斜坡，斜坡面的水平范围为7000 m.正演模拟时，沿测线方向，在-17000~19000 m范围内等间距布置181个发射源，均位于海底上方50 m；19个海底接收站等间距地布置于-8000~10000 m范围内；将接收站所在位置及相邻接收站中点所在位置作为CMP域的分割点，19个接收站共等分为36个CMP单元.设定3个发射频率分别为0.1 Hz，0.5 Hz，0.75 Hz.
|
图 10 2D验证模型 Fig. 10 The 2D verification model |
反演时使用合成数据的水平电场Ey分量的振幅和相位(加入2%的随机高斯噪声).反演时初始模型为：空气-海水-岩层的三层模型，空气、海水的电阻率和厚度为定值，岩层电阻率的初始值为1 Ωm.最终反演剖面如图 11所示，由图可见：
|
图 11 推覆构造模型CMP域反演结果 Fig. 11 The CMP inversion profile of the 2D verification model |
(1)反演结果能基本反映出各岩体的分布形态和轮廓，对于覆盖层，CMP域反演结果的轮廓能够较好地吻合真实模型，反演电阻率也接近于真实模型电阻率；
(2)对于200 Ωm的沉积层以及低阻推覆体，反演结果能够很好地体现出倾斜界面的存在，反演电阻率与模型电阻率也吻合较好；
(3)该方法的假设条件是在各CMP单元内横向上电阻率是均匀的，各CMP单元内的反演结果为电阻率在当前深度的综合表现，当所划分的相邻CMP单元存在电阻率横向的突变时，反演结果会出现扰动，若对测站进行加密，该现象将显著改善，同时这也取决于实际勘探的方式和效率.
为了更加直观地分析CMP域反演误差及单元内数据拟合情况，定义第i个数据的拟合差Dd如下：
|
(12) |
其中，Fmod和Fobs分别为反演模型正演响应和观测数据，δ为标准差，当反演模型与真实模型完全相同时，数据拟合差Dd等于1.0.以第18个CMP单元(500~1000 m，图 11中虚线标记的CMP单元)中的反演数据为例，不同收发距的水平电场振幅(“+”线)与由反演模型的正演计算结果得到的水平电场振幅(“ 
|
图 12 CMP单元数据拟合图 Fig. 12 Data misfit in CMP inversion gather |
总体来说，CMP域快速一维反演能够有效地反映理论模型的地电结构，理论计算验证其是可行且有效的.
4 结论本文提出了将二维海洋CSEM响应转换成一维CMP域响应的方法并进行处理和快速反演，研究结果表明：
(1)由海洋CSEM资料的CMP单元响应能够直接反映出天然气水合物横向分布范围，在横向上可较好地分辨高阻储层.
(2)通过不同天然气水合物模型CMP域的水平电场振幅和相位分布，以及归一化响应特征的讨论，分析进行CMP域反演的可行性.
(3)运用一维高斯-牛顿反演对不同CMP单元响应进行反演解释，实现了二维结构的CMP域一维快速反演算法，由板状天然气水合物模型和推覆构造体模型的CMP域反演结果表明，CMP域反演能够快速合理地解释二维海洋可控源电磁资料.
| Abubakar A, Habashy T M, Druskin V L, et al. 2008. 2.5 D forward and inverse modeling for interpreting low-frequency electromagnetic measurements. Geophysics, 73(4): F165-F177. DOI:10.1190/1.2937466 | |
| Boswell R, Collett T S, Frye M, et al. 2012. Subsurface gas hydrates in the northern Gulf of Mexico. Marine & Petroleum Geology, 34(1): 4-30. | |
| Buland A, Kolbjørnsen O. 2012. Bayesian inversion of CSEM and magnetotelluric data. Geophysics, 77(1): E33-E42. DOI:10.1190/geo2010-0298.1 | |
| Chen J P, Alumbaugh D. 2009. Three methods for mitigating airwaves in shallow water marine CSEM data.//SEG Technical Program Expanded Abstracts. SEG, 785-789. | |
| Constable S, Weiss C J. 2006. Mapping thin resistors and hydrocarbons with marine EM methods:Insights from 1D modeling. Geophysics, 71(2): G43-G51. DOI:10.1190/1.2187748 | |
| Constable S C, Parker R L, Constable C G. 1987. Occam's inversion:a practical algorithm for generating smooth models from electromagnetic sounding data. Geophysics, 52(3): 289-300. DOI:10.1190/1.1442303 | |
| Crepaldi J L S, Buonora M P P, Figueiredo I. 2011. Fast marine CSEM inversion in the CMP domain using analytical derivatives. Geophysics, 76(5): F303-F313. DOI:10.1190/geo2010-0237.1 | |
| Hansen P C. 1992. Analysis of discrete ill-posed problems by means of the L-curve. SIAM Review, 34(4): 561-580. DOI:10.1137/1034115 | |
| Key K. 2009. 1D inversion of multicomponent, multifrequency marine CSEM data:Methodology and synthetic studies for resolving thin resistive layers. Geophysics, 74(2): F9-F20. DOI:10.1190/1.3058434 | |
| Kong F N, Westerdahl H, Ellingsrud S, et al. 2002. 'Seabed logging':A possible direct hydrocarbon indicator for deepsea prospects using EM energy. Oil and Gas Journal, 100(19): 30-38. | |
| Li Y G, Key K. 2007. 2D marine controlled-source electromagnetic modeling:Part 1-An adaptive finite-element algorithm. Geophysics, 72(2): WA51-WA62. DOI:10.1190/1.2432262 | |
| Li Z C, Wang H Z, Ma Z T. 2000. Migration velocity analysis of common mid point gathers. Geophysical Prospecting for Petroleum (in Chinese), 39(1): 20-26. | |
| Liu Y, Li Y G, Liu J X, et al. 2013. One-dimentional inversion of marine controlled-source electromagnetic fields. The Chinese Journal of Nonferrous Metals (in Chinese), 23(9): 2551-2556. DOI:10.1016/S1003-6326(13)62767-3 | |
| Liu Y. 2014. 2D finite element modeling and inversion for marine controlled-source electromagnetic fields[Ph. D. thesis] (in Chinese). Qingdao:Ocean University of China. | |
| Luo M, LI Y G, Li G. 2016. Frequency-domain inversion of marine CSEM data in one-dimensional vertically anisotropic structures. Chinese J. Geophys. (in Chinese), 59(11): 4349-4359. DOI:10.6038/cjg20161134 | |
| Mittet R, Brauti K, Maulana H, et al. 2008. CMP inversion and post-inversion modelling for marine CSEM data. First Break, 26(8): 59-67. | |
| Newman G A, Alumbaugh D L. 1997. Three-dimensional massively parallel electromagnetic inversion-Ⅰ.Theory.. Geophysical Journal International, 128(2): 345-354. DOI:10.1111/gji.1997.128.issue-2 | |
| Pei J X, Wang Q, Zhang X L. 2015a. Effective anomaly for gas hydrate detection with marine CSEM method. Oil Geophysical Prospection (in Chinese), 50(1): 177-183. | |
| Pei J X, Wang Q, Yuan X. 2015b. Lateral resolution for 2D horizontal Plate-Like hydrate reservoir detection using marine CSEM method. Marine Geology Frontiers (in Chinese), 31(6): 44-49. | |
| Ray A, Key K. 2012. Bayesian inversion of marine CSEM data with a trans-dimensional self parametrizing algorithm. Geophysical Journal International, 191(3): 1135-1151. | |
| Sasaki Y. 2013. 3D inversion of marine CSEM and MT data:An approach to shallow-water problem. Geophysics, 78(1): E59-E65. DOI:10.1190/geo2012-0094.1 | |
| Shedd W, Boswell R, Frye M, et al. 2012. Occurrence and nature of "bottom simulating reflectors" in the northern Gulf of Mexico. Marine & Petroleum Geology, 34(1): 31-40. | |
| Su Z, Wu N Y, Zhang K N. 2011. Assessment of gas production potential of hydrate deposits at Shenhu area on northern continental slope of south China Sea. Marine Geology Frontiers (in Chinese), 27(6): 16-23. | |
| Wang T. 2015. Research on marine CSEM inversion in CMP domain of submarine gas hydrate[Master thesis] (in Chinese). Qingdao:Ocean University of China. | |
| Weitemeyer K A, Constable S C, Key K W, et al. 2006. First results from a marine controlled-source electromagnetic survey to detect gas hydrates offshore Oregon. Geophys. Res. Lett., 33(3). | |
| Zhdanov M S. Geophysical Inverse Theory and Regularization Problems.The Netherlands: Elsevier Science, 2002. | |
| 李振春, 王华忠, 马在田. 2000. 共中心点道集偏移速度分析. 石油物探, 39(1): 20–26. | |
| 刘颖, 李予国, 柳建新, 等. 2013. 海洋可控源电磁场的一维反演. 中国有色金属学报, 23(9): 2551–2556. | |
| 刘颖. 2014.海洋可控源电磁法二维有限元正演及反演[博士论文].青岛:中国海洋大学. | |
| 罗鸣, 李予国, 李刚. 2016. 一维垂直各向异性介质频率域海洋可控源电磁资料反演方法. 地球物理学报, 59(11): 4349–4359. DOI:10.6038/cjg20161134 | |
| 裴建新, 王启, 张秀丽. 2015a. 海洋CSEM探测海底天然气水合物的有效异常研究. 石油地球物理勘探, 50(1): 177–183. | |
| 裴建新, 王启, 袁翔. 2015b. 海洋CSEM法对二维水平板状水合物储集体的横向分辨能力研究. 海洋地质前沿, 31(6): 44–49. | |
| 苏正, 吴能友, 张可霓. 2011. 南海北部陆坡神狐天然气水合物开发潜力. 海洋地质前沿, 27(6): 16–23. | |
| 王涛. 2015.海底天然气水合物的海洋CSEM法CMP域反演研究[硕士论文].青岛:中国海洋大学. | |