基于DTW-Wavelet的混凝土坝变形特征辨识方法

李 俊，王丽蓉，李 斌，高焕焕，张 群

(中国电建集团西北勘测设计研究院有限公司，陕西 西安 710065)

0 前 言

混凝土坝变形安全监测数据是一种随时间或空间变化的具有非平稳、非线性特征的数据信号，是水压、温度、降雨、时效等因素的综合反映[1]。混凝土坝变形数据特征辨识，是大坝安全监测资料分析进而构建监控模型、拟定预警指标、进行安全综合评判等的前提和基础。

混凝土坝变形监测数据可看作由不同频率构成的数字信号，库水位、温度等所引起的分量一般呈年周期性变化，对应的频率较高；时效分量对应的信号表现为长周期低频率[2]。王超等[3]利用小波多分辨率分析将结构的时变参数在多尺度上展开，截断高频细节成分，将结构时变参数识别问题转化为时不变小波系数估计问题。

本文通过对大坝效应量监测数据执行小波多分辨率分析，辨识大坝效应量的周期性变化、趋势性变化等特征，有助于剖析大坝的真实工作性态，从而更客观地评估大坝安全状况。

1 混凝土坝变形监测数据小波多分辨分析

小波多分辨率分析的原理是构造一定的函数空间，将混凝土坝变形监测数据序列分解到此空间中进行计算处理，实现信息重构。

1.1 小波多分辨率分析

小波多分辨率的前提是在空间中L2(R)，存在一个子空间Vj⊂L2(R)，j∈Z，满足以下条件[4-8]：

(1)单调性

Vj⊂Vj-1，∀j∈Z

(1)

(2)伸缩性

f(t)∈Vj⟺f(2t)∈Vj-1

(2)

(3)平移不变性

f(t)∈V0⟺f(t+n)∈V0,∀n∈Z

(3)

(4)正交性

(4)

(5)逼近性

(5)

式中，Vj、Vj-1为尺度j、j-1上的尺度空间；Z为自然数集合；t为时间变量；m,n为时间平移量；φ(t)∈V0且积分值不为零；{φ(t-n)|n∈Z}是V0的标准正交基。

二尺度方程是小波多分辨率赋予尺度函数φ(t)和小波函数φ(t)最基本的特征，用来描述相关子空间下基函数φj-1,k(t)、φj,k(t)和φj-1,t(t)、φj,t(t)之间的相互关系，即

(6)

(7)

式中，φ(t)为尺度函数；φ(t)为小波函数。h(n)、g(n)分别为邻近二尺度空间对应基函数的相关系数，计算公式为

h(n)=〈φ,φ-1,n〉,g(n)=〈φ,φ-1,n〉

(8)

则式(6)和式(7)，对于∀j∈Z：

(9)

(10)

式中的h(n)、g(n)值不会随尺度j的变化而变化。

1.2 混凝土坝变形监测数据小波分解与重构

将二尺度方程在时间轴上按式(2)、式(3)进行伸缩、平移，则有

(11)

并进行内积处理变换，可得小波变换重构系数公式

(12)

式中，cj,k、dj,k、cj-1,m分别为小波变换前后的重构系数；k为平移量。

图1为小波多分辨率分析的分解结构。从图1可以看出，多分辨率分析首先将混凝土坝变形监测数据分解为低频部分a1和高频部分d1，再将低频部分a1进一步分解为低频部分a2和高频部分d2，以此类推。重复上述步骤，可得到任意频率下混凝土坝变形监测数据序列的低频和高频部分，即

s=an+d1+d2+…+dn

(13)

式中，s为混凝土坝变形原始监测数据序列；an为第n次分解后的低频子序列；di为第i次分解后的高频子序列，i=1，2，…，n。

图1 混凝土坝变形监测数据小波多分辨率分析结构

2 基于DTW的最优小波分解层数确定

在小波多分辨率分解过程中，不同小波分解层数所辨识的高频子序列和低频子序列有较大差异，以往多根据低频子序列的变化趋势是否符合时效分量变化特征为标准确定最终分解层数，为提高辨识精度，减小主观影响，本文提出以动态时间弯曲距离度量小波分解所得周期性变化的高频子序列与环境量序列间的差异，以此为标准确定最优分解层数。

动态时间弯曲距离(DynamicTimeWarping，DTW)最初用于识别语音的相似性[9-11]。DTW与传统距离最大的区别在于其模糊性，DTW并不是严格对同一时刻的两个测点进行距离度量，而是允许对应点间有相互错动，使得两组监测序列能够同步进行不是特别精确地匹配计算。

对于混凝土坝变形监测数据序列X={x1，x2，…，xn}和环境量监测数据序列Y={y1，y2，…，ym}，两者之间的DTW为

(14)

Dmin=min{Ddtw(X,Rest(Y)),Ddtw(Rest(X),Rest(Y)),Ddtw(Rest(X),Y)}

(15)

其中，Rest(X)={x2，x3，…，xn}，Rest(Y)={y2，y3，…，ym}。

从式(14)和式(15)可以看出，DTW是一种递归计算过程，计算较为繁琐。为提高计算效率，做如下简化：首先构造一个大小为n×m的矩阵，称为弯曲矩阵，如图2所示。每一个方格代表效应量与环境量的匹配关系，值为w(xi，yj)，称之为基距离。基距离函数有很多，通常可以表示为

(16)

式中，N为最高维度值；xin和yjn分别表示xi和yj在第n维的值；一般情况下取m=3。

图2 混凝土坝变形与环境量监测数据DTW弯曲矩阵

由图2可知，混凝土坝变形与环境量监测序列上两测值点间的对应关系转换为弯曲矩阵中从一个角点到另一个角点的弯曲路径，对应的数据点距离记为W，则W=[w1，w2，…，wk]，k为弯曲路径经过的方格数，可得如下关系

max(m,n)≤k≤m+n

(17)

定义函数fw(xi，yj)为

wl=fw(xi,yj),1≤i≤m,1≤j≤n,1≤l≤k

(18)

混凝土坝变形与环境量监测数据序列之间的DTW计算可以转换为求解弯曲矩阵中最小弯曲距离的问题，即

(19)

式中，arg min(*)表示使*最小时，其自变量取值。

混凝土坝变形与环境量之间的最小弯曲距离DTW(X,Y)可用于度量两监测序列间的差异性。

3 混凝土坝变形监测数据特征辨识实现流程和步骤

借助小波多分辨率分析，结合DTW实现混凝土坝变形监测数据特征辨识流程如图3所示，主要步骤如下：

(1)根据混凝土坝变形监测数据序列的特征，选择合适的小波基函数，对监测数据序列进行不同层数的小波分解，得到低频子序列ai～a1和高频子序列dij～di1。

(2)将不同分解层数下的子序列进行区分，筛选出周期性变化的高频子序列予以叠加作为监测数据的周期性分量Xi。

(3)对比不同分解层数下，周期性分量与环境量序列间差异性，按照式(19)计算第(2)步中所合成的各个分解层数下的周期性分量与库水位、温度监测序列的动态时间弯曲距离DTW(Xi,W)、DTW(Xi,T)。选择DTW(Xi,W)和DTW(Xi,T)均值最小的分解层数作为最优分解层数k。

(4)基于最优分解层数k，选择其中周期性变化的高频子序列相加，其结果代表由环境量引起的周期性分量；趋势性明显的低频子序列ak即为时效分量。

图3 混凝土坝变形监测数据特征辨识流程

4 工程实例

某水库大坝溢流坝段为混凝土重力坝，共分为12个坝段。在溢流坝段坝顶布置有一条激光准直测线，每个坝段布置一个测点，编号为LA-1～LA-12，用以监测各坝段坝顶水平位移情况，测点布置情况见图4。应用本文所述方法，对该溢流坝段的变形监测数据进行特征辨识分析，验证方法的可行性和有效性。

按照图3所示流程，对该大坝溢流坝段坝顶水平位移测点LA-2从2009年6月30日至2018年12月26日的监测数据进行小波多分辨率分析，辨识其特征分量。选用工程中常用的db4小波基函数对监测数据序列进行分解，分解层数n=5～15，辨识周期性分量，如图5所示。根据式(19)分别计算不同分解层数高频周期性分量与库水位、气温间的DTW(Xi,W)和DTW(Xi,T)，如图6所示。

图4 溢流坝段激光准直测点布置

图5 不同分解层数周期性分量过程线

图6 不同分解层数下周期性分量与库水位、气温偏差度

从图6中可以看出，当分解层数n=5～9时，周期性分量与库水位、气温的DTW(Xi,W)和DTW(Xi,T)均较小，且变化平稳；当n≥10时，周期性分量与库水位、气温的DTW(Xi,W)和DTW(Xi,T)增加较大，随着分解层数的增加有逐渐减小的趋势。综合考虑，选取最优分解层数k=9，所辨识的周期性分量及时效分量如图7所示。

图7 最优分解层数(k=9)下周期性分量与时效分量

5 结 论

本文所提出的变形特征辨识方法是基于小波多分辨率分解混凝土坝变形监测数据，引入动态时间弯曲距离用于度量周期性分量和环境量间的差异性，实现了混凝土坝变形监测数据更为客观、高效的特征辨识，将该方法的特点总结如下：

(1)该方法依据根据混凝土坝变形监测数据中蕴含的水压、温度、时效等不同时频特征，利用小波多分辨率特性进行辨识，具有很强的理论基础，可靠性较高。

(2)该方法引入动态时间弯曲距离用于度量周期性分量和环境量间的差异性，作为小波分解层数的判别条件，减小了人为干预，使得辨识结果更具有客观性。

(3)将该方法用于实际工程中，所辨识得到的周期性分量和时效分量符合一般规律，说明将其应用于混凝土坝变形特征辨识是可行的，有效的。