海水涡动黏性与湍流结构密切相关, 在上层海洋运动中起重要作用^[1]。在海水混合、海洋环流、中尺度涡等动力过程的数值模型研究中, 垂向涡动黏性系数(vertical eddy viscosity coefficient, VEVC)的参数化具有重要意义^[2-4]。但VEVC难以被直接测量, 通常需要根据湍流模型计算, 比如依赖于Richardson数的模型、k-ε模型, M-Y模型, Prandtl混合长模型、KPP模型等^[5-9]。由于海洋环境往往会随空间和时间发生变化, 再加上观测数据的不足, 各个湍流模型的应用容易因普适性问题而受到限制。不同模型计算出的垂向涡动黏性系数差别较大, 没有某一个模型能精准适用于全球所有海域。

本文基于数学物理反问题的思想, 采用数据同化方法, 由观测资料反演VEVC。数据同化方法是在严格的数学理论基础上发展起来的, 比如反问题理论和最优控制理论^[10]。基于最优控制方法和摄动理论的四维变分同化是海洋数据同化未来的发展方向之一^[11]。伴随同化方法作为四维变分同化方法的一种, 具有提取不同种观测数据有效信息, 同时与模型保持动力学和物理一致性的优势, 是参数反演的有力工具^[12-13]。

Ekman理论在物理海洋学的研究中有悠久的历史, 奠定了风生海流的理论基础。Ekman模型具有不少优点, 例如: (1)不需要显式地表达浮力强迫和层化效应, 而是通过黏性和边界层参数的取值反应其影响; (2)容易求解, 便于探究具体参数的影响^[14]。Ekman层的VEVC垂向结构对Ekman螺旋结构和Ekman输运都有重要影响^[15-17]。近年来, 伴随同化方法被广泛应用于VEVC的估计^[18-22]。在过去的研究中, VEVC往往被设为常数或随深度变化的量, 而与时间无关。然而海洋是一个时变系统, 海洋动力学参数也应该是随时间变化的。例如, 在不同风况下, 海面拖曳系数会随时间和风速变化^[23-24]。底摩擦系数等其他参数在时间域上也会发生变化^[25-27]。海表与海底摩擦、升降流、大小潮等因素导致的热量传递过程和垂向混合过程, 都与海水垂向涡动黏性密切相关^[28-30]。显然, 如果能在海洋模型中实现时变的VEVC, 模拟的准确性可以得到提高。因此, 本文旨在利用伴随法, 将观测资料同化到Ekman模型中, 对时变VEVC的伴随同化反演进行研究。

1 模型描述 1.1 控制方程(正向模型)

在本文中, Ekman模型的控制方程如公式(1)~公式(8)^[18]:

$ \left\{ \begin{array}{l} \frac{{\partial u}}{{\partial t}} - fv = \frac{\partial }{{\partial z}}\left( {A\frac{{\partial u}}{{\partial z}}} \right)\\ \frac{{\partial v}}{{\partial t}} + fu = \frac{\partial }{{\partial z}}\left( {A\frac{{\partial v}}{{\partial z}}} \right) \end{array}, \right. $

(1)

相应的边界条件为

$ \begin{array}{l} {\left. {\frac{{\partial u}}{{\partial z}}} \right|_{z = 0}} = \frac{{{C_{\rm{D}}}}}{A}\frac{{{\rho _{\rm{a}}}}}{{{\rho _w}}}\sqrt {u_{\rm{a}}^{\rm{2}} + v_{\rm{a}}^{\rm{2}}} {u_{\rm{a}}}\;, \;\;{\left. {\frac{{\partial v}}{{\partial z}}} \right|_{z = 0}} = \\ \;\;\;\;\;\;\;\;\;\;\;\;\frac{{{C_{\rm{D}}}}}{A}\frac{{{\rho _{\rm{a}}}}}{{{\rho _{\rm{w}}}}}\sqrt {u_{\rm{a}}^2 + v_{\rm{a}}^{\rm{2}}} {v_{\rm{a}}}, \; \end{array} $

(2)

$ A{\left. {\frac{{\partial u}}{{\partial z}}} \right|_{z = - {H_0}}} = 0\;, \;\;A{\left. {\frac{{\partial v}}{{\partial z}}} \right|_{z = - {H_0}}} = 0, $

(3)

相应的初始条件为

${\left. u \right|_{t = 0}} = {u_0}\;, \;\;{\left. v \right|_{t = 0}} = {v_0},$

(4)

其中, z轴向上为正且海表面处z=0, u(z)(向东为正)和v(z)(向北为正)为各深度水平方向的流速分量; f为科氏参数(常值), C_D为风应力拖曳系数, A为涡动黏性系数; ρ_a和ρ_w分别为空气和水的密度, u_a和v_a是风速的分量, H₀为Ekman层的深度, u₀和v₀是初始流速。

在本文中, 模型参数设置如下:

f=10^–4 s^–1(中纬度海域);

u₀(j)=0.01cos[π/4–(j–0.5)π/20]e^{–(j–0.5)π/20} m/s, v₀(j)= 0.01sin[π/4–(j–0.5)–π/20]e^{–(j–0.5)π/20} m/s(生成厘米级初始流速, 其中j为层数);

u_a=10sin(2πiΔt/T₀) m/s(生成周期性变化风速, 其中i为时间步), v_a=0(假设只有东西向风速);

T₀=10 h(风速变化周期), T=10 d(总积分时间), Δt=0.5 h(时间步长);

C_D=1.2×10^–3; A(i)=0.02×[1.2+0.5sin(2πiΔt/24T₀)+ 0.5cos(2πiΔt/4T₀)] m²/s, (其中i为时间步);

ρ_a=1.2 kg/m³, ρ_w=1.025×10³ kg/m³;

H₀=100 m, Δz=5 m。

由于控制方程中变量具有不同的单位和量级, 为了使方程具有更好的通用性, 我们参考前人工作^[18], 进行了变量替换。由此得到方程:

$\left\{ \begin{array}{l} \frac{{\partial u}}{{\partial t}} - v = \frac{\partial }{{\partial z}}\left( {A\frac{{\partial u}}{{\partial z}}} \right)\\ \frac{{\partial v}}{{\partial t}} + u = \frac{\partial }{{\partial z}}\left( {A\frac{{\partial v}}{{\partial z}}} \right) \end{array}, \right.$

(5)

相应的边界条件为

${\left. {\frac{{\partial u}}{{\partial z}}} \right|_{z = 0}} = {C_D}\sqrt {u_{\rm{a}}^2 + v_{\rm{a}}^2} {u_{\rm{a}}}\;, \;\;{\left. {\frac{{\partial v}}{{\partial z}}} \right|_{z = 0}} = {C_D}\sqrt {u_{\rm{a}}^2 + v_{\rm{a}}^2} {v_{\rm{a}}},$

(6)

$ A{\left. {\frac{{\partial u}}{{\partial z}}} \right|_{z = - {H_0}}} = 0\;, \;\;A{\left. {\frac{{\partial v}}{{\partial z}}} \right|_{z = - {H_0}}} = 0, $

(7)

相应的初始条件为

${\left. u \right|_{t = 0}} = {u_0}\;, \;\;{\left. v \right|_{t = 0}} = {v_0}.$

(8)

1.2 伴随方程(反向模型)

构造代价函数

$J\left( {u, v, A} \right) = \frac{1}{2}{K_m}\int_t {\int_z {\left( {{{\left( {u - \hat u} \right)}^2} + {{\left( {v - \hat v} \right)}^2}} \right){\rm{d}}z} {\rm{d}}t}, $

(9)

其中, $\hat u$和$\hat v$为流速观测值; K_m为权重矩阵, 在本文中设为单位矩阵^[31]。

构造拉格朗日函数

$\begin{array}{l} L\left( {u, v, A, \lambda , \mu } \right) = J + \\ \int_t {\int_z {\left[ {\lambda \left( {\frac{{\partial u}}{{\partial t}} - v - \frac{\partial }{{\partial z}}\left( {A\frac{{\partial u}}{{\partial z}}} \right)} \right) + \mu \left( {\frac{{\partial v}}{{\partial t}} + u - \left( {A\frac{{\partial v}}{{\partial z}}} \right)} \right)} \right]{\rm{d}}z} {\rm{d}}t}, \end{array}$

(10)

其中, λ和μ分别是u和v的拉格朗日乘子(又称为伴随变量), 代价函数在控制方程约束下的最小化问题, 变成寻找拉格朗日函数梯度消失时的u, v, A, λ和μ的驻点问题。由此得到以下方程:

$\frac{{\partial L\left( {u, v, A, \lambda , \mu } \right)}}{{\partial \lambda }} = 0\;, \;\;\frac{{\partial L\left( {u, v, A, \lambda , \mu } \right)}}{{\partial \mu }} = 0,$

(11)

$\frac{{\partial L\left( {u, v, A, \lambda , \mu } \right)}}{{\partial u}} = 0\;, \;\;\frac{{\partial L\left( {u, v, A, \lambda , \mu } \right)}}{{\partial v}} = 0,$

(12)

$\frac{{\partial L\left( {u, v, A, \lambda , \mu } \right)}}{{\partial A}} = 0.$

(13)

方程组(11)对应原始的控制方程, 方程组(12)得到伴随方程:

$\left\{ \begin{array}{l} \frac{{\partial \lambda }}{{\partial t}} - \mu + \frac{\partial }{{\partial z}}\left( {A\frac{{\partial \lambda }}{{\partial z}}} \right) = {K_m}\left( {u - \hat u} \right)\\ \frac{{\partial \mu }}{{\partial t}} + \lambda + \frac{\partial }{{\partial z}}\left( {A\frac{{\partial \mu }}{{\partial z}}} \right) = {K_m}\left( {v - \hat v} \right) \end{array}, \right.$

(14)

相应的边界条件为

${\left. {\frac{{\partial \lambda }}{{\partial z}}} \right|_{z = 0}} = 0\;, \;\;{\left. {\frac{{\partial \mu }}{{\partial z}}} \right|_{z = 0}} = 0,$

(15)

${\left. {\frac{{\partial \lambda }}{{\partial z}}} \right|_{z = - {H_0}}} = 0\;, \;\;{\left. {\frac{{\partial \mu }}{{\partial z}}} \right|_{z = - {H_0}}} = 0,$

(16)

相应的初始条件为

${\left. \lambda \right|_{t = T}} = 0\;, \;\;{\left. \mu \right|_{t = T}} = 0.$

(17)

对控制方程进行正向积分, 得到模拟流速值; 对伴随方程进行反向积分, 得到伴随变量值。一般来说, 反向积分时间等于正向积分时间^[32]。

涡动黏性系数A(t)的值, 可由方程(13)结合优化算法进行反演得到。

1.3 优化算法

利用方程(13), 可求解代价函数关于涡动黏性系数的梯度。

若将VEVC视为恒定值:

$ \frac{{\partial J}}{{\partial A}} = - \int_t {\int_z {\left( {\frac{{\partial \lambda }}{{\partial z}}\frac{{\partial u}}{{\partial z}} + \frac{{\partial \mu }}{{\partial z}}\frac{{\partial v}}{{\partial z}}} \right){\rm{d}}z} {\rm{d}}t} . $

(18)

若考虑VEVC的深度变化:

$ \frac{{\partial J}}{{\partial A\left( z \right)}}{\rm{ = }} - \int_t {\left( {\frac{{\partial \lambda }}{{\partial z}}\frac{{\partial u}}{{\partial z}} + \frac{{\partial \mu }}{{\partial z}}\frac{{\partial v}}{{\partial z}}} \right){\rm{d}}t} . $

(19)

在本文中, 考虑VEVC的时间变化:

$ \frac{{\partial J}}{{\partial A\left( t \right)}}{\rm{ = }} - \int_z {\left( {\frac{{\partial \lambda }}{{\partial z}}\frac{{\partial u}}{{\partial z}} + \frac{{\partial \mu }}{{\partial z}}\frac{{\partial v}}{{\partial z}}} \right){\rm{d}}z} . $

(20)

结合梯度与优化算法, 即可反演涡动黏性系数:

${A_{n + 1}} = {A_n} + {\alpha _n}{d_n},$

(21)

其中, A_n和A_n+1为第n次迭代中涡动黏性系数的当前值和调整值, α_n和d_n分别代表迭代步长和搜索方向。经过反复调试^[33], 在本文中α_n设为4.0×10^–4。确定d_n有多种可行的优化算法^[34], 本文应用了其中三种, 分别是梯度下降法(gradient descent, GD)、共轭梯度法(quasi-Newton conjugate gradient, CG)和有限记忆BFGS(limited-memory Broyden-Fletcher-Goldfarb-Shanno, L-BFGS)法。

在梯度下降法中, 搜索方向为:

${d_n} = - {g_n} = - \left( {\frac{{\partial J}}{{\partial A}}} \right)_n^0 = - \frac{{{{\left( {\partial J/\partial A} \right)}_n}}}{{\left\| {{{\left( {\partial J/\partial A} \right)}_n}} \right\|}},$

(22)

其中, (∂J/∂A)_n是第n次迭代中代价函数关于涡动黏性系数的梯度, 而||(∂J/∂A)_n||是其L₂范数。

在共轭梯度法中, 搜索方向被定义为共轭方向:

$ {d_n} = - {g_n} + \frac{{{{\left\| {{g_n}} \right\|}^2}}}{{{{\left\| {{g_{n - 1}}} \right\|}^2}}}{d_{n - 1}}. $

(23)

L_BFGS法是有限记忆的BFGS法, 它要求搜索方向满足^[35]:

$ {d_n} = - {H_n}{g_n}, $

(24)

其中,

$\begin{array}{l} {H_n} = \left( {V_{n - 1}^T \cdots V_{n - m}^T} \right){H_0}\left( {{V_{n - m}} \cdots {V_{n - 1}}} \right)\\ \;\;\;\;\;\;\; + {\rho _{n - m}}\left( {V_{n - 1}^T \cdots V_{n - m + 1}^T} \right){s_{n - m}}s_{n - m}^T\left( {{V_{n - m + 1}} \cdots {V_{n - 1}}} \right)\\ \;\;\;\;\;\;\; + {\rho _{n - m + 1}}\left( {V_{n - 1}^T \cdots V_{n - m + 2}^T} \right){s_{n - m + 1}}s_{n - m + 1}^T\left( {{V_{n - m + 2}} \cdots {V_{n - 1}}} \right)\\ \;\;\;\;\;\;\; + \cdots \\ \;\;\;\;\;\;\; + {\rho _{n - 1}}{s_{n - 1}}s_{n - 1}^T.\; \end{array}$

(25)

在方程(25)中, ρ_n=1/s_ny_n^T, V_n=1–ρ_ny_ns_n^T, s_n=A_n+1–A_n=α_nd_n, y_n=g_n+1–g_n; m是校正对的数量, 通常取在3~7, 本文中设为5^[36]。

2 理想试验

为了验证模型的合理性和可行性, 设计若干组理想孪生实验, 分别探讨优化算法、初始猜测、观测深度等因素对实验结果的影响。具体步骤如下:

(1) 以涡动黏性系数的给定值运行正向模型, 将模拟流速视为观测值;

(2) 以涡动黏性系数的猜测值(初始猜测值A₁=0.002 m²/s)运行正向模型, 得到流速模拟值;

(3) 以流速观测值和模拟值的差驱动反向模型, 对猜测值进行优化;

(4) 重复步骤(2)和(3), 直至结果满足收敛性判别标准, 得到涡动黏性系数的最终反演值。

2.1 优化算法的影响

GD、CG和L-BFGS是三种常用的优化算法。一些研究表明, L-BFGS法相对于其他算法具有明显优越性^[37-39]; 不过, 也有研究认为GD法反演的结果更好^[40-41]。为此, 我们设计了一组实验探讨三种算法的优劣。

在参数优化的过程中, 迭代了足够多步, 选取其中最好的结果进行展示分析。本组垂向涡动黏性系数的给定值和反演值如图 1a所示, 代价函数下降情况如图 1b所示。同化后的代价函数和均方根误差见表 1。

在本组实验中, 三个算例的代价函数都下降了4个数量级, 均方根误差都下降了1个数量级, 均得到了不错的结果。其中2号和3号的收敛速度高于1号, 1号和3号的反演结果明显好于2号。综合考虑, 和另外两种算法相比, GD法处理本问题的效果更好。因此, 在以下的实验中, 我们选取GD法来进行反演。

2.2 初始猜测的影响

合理的初猜测值可以提高反演结果, 加快收敛速度^[42]。因此, 我们设计了本组实验以探究不同初始猜测值对本模型反演结果的影响。本组垂向涡动黏性系数的给定值和反演值如图 2a所示, 代价函数下降情况如图 2b所示。同化后的代价函数和均方根误差见表 2。

本组实验中, 三个算例的反演结果都接近一开始的给定值。在同化后, 代价函数下降了3~4个数量级, 均方根误差下降了1个数量级。1号实验的初始猜测值最小, 优化效果最显著。4号实验的初始猜测值与给定分布的平均值最接近, 所需迭代步数最少, 反演误差也最小。实验表明, 与实际值接近的初始猜测值能提高反演结果。处理实测数据时, 可以参考其他文献中的涡动黏性估计结果^{[29, 43-44]}, 给定合适的初始猜测值。

在以下的实验中, 为了展示明显的同化前后结果对比, 仍采用较小的初始猜测值。

2.3 观测深度的影响

在以上几组实验中, 我们用到了所有层的流速数据。然而在实际操作中, 由于观测仪器盲区等原因, 在某些深度上很可能没有观测数据。吕咸青^[45]曾仅凭对海洋表层和次表层的“观测”, 反演出海洋涡动黏性剖面。在本组实验中, 我们也进行了类似的尝试, 即通过部分层的数据来进行反演, 探究数据缺失对结果的影响。

本组垂向涡动黏性系数的给定值和反演值如图 3a所示, 代价函数下降情况如图 3b所示。同化后的代价函数和均方根误差见表 3。

从代价函数下降情况来看, 使用的数据越少, 收敛速度越快。就反演结果而言, 6号优于7号, 但与使用全部层的流速相比仍略逊一筹。这说明反演结果对表层和次表层的流速比较敏感, 而较深层流速对反演结果的贡献较小, 与Ekman理论相符合。不过整体而言, 本组实验都成功地反演了VEVC的时间变化, 说明本模型能在一定程度上适应观测比较缺乏的海洋环境。

3 实际实验 3.1 数据描述

本文研究的区域为百慕大群岛邻近海域(如图 4a所示)。这片区域极端事件频发, 1851~2005年期间受到188次热带风暴和飓风袭击^[46], 是研究上层海洋响应风应力变化的理想区域。

本文的研究数据来自百慕大试验站锚泊系统(Bermuda testbed mooring, BTM)第24次部署(如图 4b所示)。BTM是位于百慕大东南80 km处的深水平台, 自1994年开始部署。它能在传统采样无法进行的恶劣天气保持观测, 收集高时间分辨率的长时间序列数据, 可用于开发测试新仪器、收集跨学科数据, 以及为卫星遥感传感器提供标定和验证数据^[47-48]。

风速数据由海表面的浮标塔测量, 时间分辨率为5 min。风速及风应力时间序列如图 5a所示。流速数据由200 m深处的声学多普勒流速剖面仪(acoustic Doppler current profiler, ADCP)测量, 时间分辨率为15 min, 垂向分辨率为3 m, 采样深度为26.52~ 191.52 m。26.52 m处流速时间序列如图 5b所示。

我们选取2006年8月9日0点~2006年8月19日0点时段的数据, 拟从中提取Ekman流分量, 并将流速插值到25~100 m深度用于本文的同化模型。

3.2 数据处理 3.2.1 滤除潮流

利用带通滤波可以从观测流速中提取近惯性流分量^{[22, 49]}。对流速时间序列进行频谱分析(如图 6a所示), 近惯性频段与全日潮P₁、K₁两个分潮频率比较接近。对观测流速进行希尔伯特变换(Hilbert transform)带通滤波处理, 提取频段范围取为[0.96 1.2]f_l, 其中f_l为当地惯性频率。选择该范围既能够避免全日潮信号污染, 又能够最大限度提取近惯性流。过滤前后的流速对比如图 6b所示。

3.2.2 分离地转流

$ {U_{{\rm{ni}}}} = {U_{{\rm{ek}}}} + {U_{{\rm{geo}}}}. $

(26)

如公式(28)所示, 近惯性流速U_ni由Ekman流速(风生流速)U_ek和地转流速U_geo两大组分构成。在以往的研究中, 常选取一个地转参考零面, 将该深度以下的流速(视为恒定值)作为地转参考速度^[50-51]。然而, 由于地转剪切的存在, 恒定地转流速的假设很可能导致对一部分地转流速与Ekman流速发生混叠, 造成对Ekman流速的错误估计^[52-53]。为此, 我们参考Roach等的做法, 将地转剪切也考虑进来, 把地转流速拆分为两部分^[54]:

$ {U_{{\rm{ni}}}} = {U_{{\rm{ek}}}} + {U_{{\rm{deep}}}} + \frac{{{\rm{d}}{U_{{\rm{geo}}}}}}{{{\rm{d}}z}}z, $

(27)

其中, U_deep是恒定的地转参考速度, dU_geo/dz是恒定的地转剪切。

提取出的近惯性流的典型流速剖面如图 7a所示, 140 m以下流速变化很小, 因此U_deep取为140 m以下的近惯性流速平均值; 考虑到模型深度设为100 m, 因此dU_geo/dz取为100~140 m近惯性流速切变平均值。提取出的Ekman流的典型流速剖面如图 7b所示。

3.3 优化算法

对抛物化Navier-Stokes方程模型的研究发现, 有无黏性可能对优化算法的反演结果产生影响^[3]。而基于内潮模型反演开边界条件的研究发现, 当同化的观测包括未知边界附近区域时, L-BFGS法表现最好, 而当同化的观测仅包括远离未知边界的区域时, 只有GD法才能取得令人满意的结果^[55]。这说明理想实验和实际实验的结果可能略有不同。

为了验证在实际实验中, 三种优化算法的表现是否跟理想实验一致, 我们设置了本组实验。本组VEVC的反演值如图 8a所示, 代价函数下降情况如图 8b所示。同化后的代价函数和VEVC平均值见表 4。

同化一开始, 三个算例的代价函数下降速度差不多。但9号和10号的下降很快陷入停滞, 只有8号代价函数下降最为显著, 为同化前的67%。结果表明, GD法在实际实验中的表现最好, 这与理想实验结果一致。因此, GD法虽然比较简单, 但因其可控且易于实施的特性, 在实际应用中仍具有强大竞争力。在本文接下来的实际实验中, 我们同样采用了GD法。

3.4 反演策略

对风暴潮的数值模拟研究发现, 与将风应力拖曳系数设为常数的参数化方案相比, 将其设为随风速变化的方案模拟结果误差更小^[56]。二维潮波模型的底摩擦系数反演研究说明, 反演策略假定的参数分布情况与实际分布一致时, 取得的反演结果最好^[57]。以往反演Ekman模型VEVC的研究, 一般把VEVC设为随深度变化而忽视其时间变化^{[18, 42, 46]}。但大量研究表明, VEVC同样有很强的时间变化特征^[28-30]。

为了对比不同方案的优劣, 我们在本组实验中采取了三种反演策略: (1)恒定; (2)深度变化; (3)时间变化。本组时间平均的流速剖面如图 9a所示, 25 m处u向流速时间序列如图 9b所示。同化后的代价函数和VEVC平均值见表 5。

由表 5可知, 11号同化后的代价函数为同化前的95%, 下降很小, 说明恒定型VEVC的反演策略优化效果很弱。12号的这一数值为78%, 比11号小, 但不及8号的67%。这说明深度变化型的反演结果相比于恒定型有一定进步, 却还是比不上时间变化型。从图 9可以看出, 采用时间变化型VEVC反演策略的流速模拟结果, 最接近观测值。由此可以得出结论, 时间变化型反演策略相较于其他两种策略具备一定优越性。相对而言, 时间变化结合深度变化更为合理, 但独立变量的增加将会加剧反问题不适定性的影响, 如何克服这一限制, 将有待未来研究。

根据实验结果, 研究区域的VEVC在0.01 m²/s左右, 且具有较大的波动。不同研究中得到的VEVC大小如表 6所示, 量级一般为O(10^–1~10^–3 m²/s)。对比可知, 本文得到的VEVC量级与它们大体范围一致。不少研究发现, 按流速随深度的大小衰减和按流速随深度的相位变化计算得到的VEVC相差8~10倍, 这被认为是层化现象导致的^{[51-52, 58-59]}。而在南大洋的研究证明, 两种方式计算结果的差异可以用地转剪切来解释^{[53, 55]}。本文在提取Ekman流时, 已经考虑了地转剪切的影响, 然而得到的VEVC值与他们的结果相比仍偏小。可能的原因是本文使用的观测数据来自中低纬度的Bermuda试验站锚泊系统, 层化现象比南大洋要强。

4 结论

本文提出了基于伴随同化技术反演Ekman模型VEVC时间变化的方法。通过一系列理想实验和实际实验中对优化算法、初始猜测值、观测深度、反演策略等因素进行了探究, 验证了方法的有效性。主要结论包括:

(1) L-BFGS法和CG法代价函数下降速度快于GD法, L-BFGS法和GD法的准确性高于CG法, 综合来看优化效果最好的是GD法; 而在实际应用中, 可控且易于实施的GD法相较于L-BFGS法和CG法仍然具有一定优势。(2)合理的初始猜测值能改善实验结果。初始猜测值与实际值接近时, 反演所需的步数减少, 最终误差变小。(3)使用不同观测深度的流速反演时, 结果有所不同。本模型对表层和次表层的流速更为敏感。(4)将VEVC设为时间变化型的反演策略, 比恒定型和深度变化型更具优越性。

此外, 本文利用从BTM数据中提取的Ekman流成功反演了该地区的VEVC, 并将其与南大洋等地区的数值进行了对比。总的来说, 本文提出的方法能适应各种不同的实验情况, 并且在百慕大地区得到了初步应用。这项工作增进了我们对Ekman模型的理解, 并提供了一种有效的方法来确定VEVC的时间变化。