# A System Identification-Based Damage-Detection Method for Gravity Dams

Authors: Masoud Mirtaheri; Mojtaba Salkhordeh; Masoud Mohammadgholiha
Journal: Shock and Vibration
 (2021)

---

## Abstract

Dams are essential infrastructures as they provide a range of economic, environmental, and social benefits to the local populations. Damage in the body of these structures may lead to an irreparable disaster. This paper presents a cost-effective vibration-based framework to identify the dynamic properties and damage of the dams. To this end, four commonly occurred damage scenarios, including (1) damage in the neck of the dam, (2) damage in the toe of the structure, (3) simultaneous damage in the neck and the toe of the dam, and (4) damage in the lifting joints of the dam, are considered. The proposed method is based on processing the acceleration response of a gravity dam under ambient excitations. First, the random decrement technique (RDT) is applied to determine the free-vibration of the structure using the structural response. Then, a combined method based on Hilbert–Huang Transform (HHT) and Wavelet Transform (WT) is presented to obtain the dynamic properties of the structure. Next, the cubic-spline technique is used to make the mode shapes differentiable. Finally, Continuous Wavelet Transform (CWT) is applied to the residual values of mode shape curvatures between intact and damaged structures to estimate the damage location. In order to evaluate the efficiency of the proposed method in field condition, 10% noise is added to the structural response. Results show promising accuracy in estimating the location of damage even when the structure is subjected to simultaneous damage in different locations.


---

## Body

## 1. Introduction

Civil engineering infrastructures provide the essential instruments for a community to have an immediate operation. These infrastructures include buildings, vehicular and pedestrian bridges, factories, transmission towers, and dams. Concrete dams provide economic benefits by supplying water for irrigation, energy for hydroelectricity generation, and obstacle for flood control. The performance of such infrastructures under ambient vibrations as well as the operational loads may reduce over their life-time. Failure of the dams leads to irreparable losses for beneficiary communities near these structures [1, 2]. The most important catastrophes that occurred in the recent years include the failure in 1864 of a 30 m embankment dam which claimed 254 people in Sheffield, UK [3], Canyon lake dam that claimed the lives of 236 people in Texas, USA [4], and Austin dam disaster that claimed the lives of 78 people [5]. Hence, a continuous Structural Health Monitoring (SHM) framework is needed to prevent such sudden failures in the dams.The process of identifying and tracking the structural damage in a civil, mechanical, or electrical infrastructure is known as structural health monitoring [6]. The existence of damage in a structural system changes its dynamic properties such as damping, natural frequencies, and mode shapes. Therefore, the damage state of a structure can be determined by identifying and evaluating these dynamic parameters [7]. In the following section, a comprehensive literature review of damage detection in structural systems as well as recent advances in SHM of dams is presented.


### 1.1. Structural Damage Detection

Damage detection methods could be divided into two general categories, model-based and data-driven methods. In the former approach, a system identification algorithm should be employed for identifying the dynamic properties of the structure. Then, the structural health condition will be determined by capturing the variation of these dynamic properties. On the other hand, data-driven methods apply features extracted from the vibration response of the structure through a statistical learning algorithm to predict the damage [8]. Cui et al. [9] proposed a damage diagnosis procedure based on the natural excitation technique to identify the structural damage under ambient conditions. The authors combined the strain response with the Eigensystem realization algorithm to identify the strain modal parameters of a structure under ambient vibration. Montazer and Seyedpoor [10] presented a damage index based on flexibility index to detect the location of damage in truss systems. They implemented their proposed method on two numerical models, including a planar truss having 31 elements and a 47-bar planar power line tower. Lee and Eun [11] provided a static damage detection procedure by comparing the variation of stress distribution before and after damage. They showed the validity of their method utilizing two truss structures. Ding et al. [12] developed a method based on Modal Strain Energy (MSE) to identify the damage in the grid structures. They employed an updated mode shape expansion method to ensure that the modal shape obtained from the reference baseline model is reliable. Liu et al. [13] proposed a method based on Stationary Wavelet Transform (SWT) for the purpose of localization and severity detection of cracks in the cantilever beams. Montejo [14] evaluated different vibration-based damage diagnosis methods for the purpose of detecting the damage in the structures under ambient excitation. They concluded that the method based on Continuous Wavelet Transform (CWT) is much more efficient than uncovering the high-frequency spikes of the structural response obtained by Hilbert–Huang transform, high-pass filtering, or Discrete Wavelet Transform (DWT). Pu et al. [15] presented the use of Frequency Response Functions (FRFs) along with the model updating theory for identifying the damage occurred in concrete beams. They formulated an optimization algorithm to set the analytical FRFs from a benchmark finite element model with those obtained through the experimental response. Worden et al. [16] implemented a statistical approach based on Mahalanobis squared distances for healthy and damaged structures. This method was actually a novelty detection problem where a damaged structure should be diagnosed from the healthy cases. Ditommaso et al. [17] presented a method for damage detection of the framed-structures subjected to earthquake motions by evaluating the variation of modal curvature in the natural frequencies. They validated their method by using numerical models and implementing some experiments. Zhang et al. [18] proposed an innovative damage index according to the macrostrain modal shapes for the purpose of detecting the potential damage in steel bridges. Yazdanpanah et al. [19] presented a damage indicator based on modal information for damage detection in beam-like structures. Their method was founded on three factors, namely, the mode shape, the slope of the mode shape, and the curvature of the mode shape. Figueiredo et al. [20] compared four different machine learning techniques by applying autoassociative neural network, factor analysis, Mahalanobis distance, and singular value decomposition to identify the structural damage. They have implemented their study on an experimental base-excited three-story frame structure. Ghiasi et al. [21] proposed a least square support vector machine based on the thin-plate spline Littlewood–Paley wavelet kernel function to identify the structural damage. Ghannadi and Kourehli [22] explored the efficiency of the moth flame optimization algorithm to identify the location and severity of structural damage. They used the Eigenvalues obtained from the equation of motion along with Modal Assurance Criteria (MAC) flexibility as damage indices. Nguyen et al. [23] presented a novel framework based on the transmissibility function fed to ANNs to identify the structural damage in bridges. They verified their method using data measured from an existing bridge in Taiwan.


### 1.2. Structural Health Monitoring of Dams

Wang and He [24] presented a methodology based on the reduction of natural frequencies to identify the crack location in arch dams. They proposed a statistical neural network to detect the crack by measuring the reductions of natural frequencies. Turker et al. [25] explored a method based on model updating to identify the global damage that occurred in the arch dams. They have validated their method on an experimental model under ambient excitation. Pirboudaghi et al. [26] proposed a multistage method based on the extended finite element method (XFEM) and continuous wavelet transform (CWT) to detect the seismic cracks of concrete dams. Sevim et al. [27] implemented an ambient vibration test to identify the dynamic properties of a real concrete arch dam, including modal frequencies, damping, and mode shapes. Bianchi and Bremen [28] proposed various procedures to capture the long-term performance of Ferden and Roggiasca dams in Switzerland. In this regard, they captured variations of the dam behavior before the snowfall and after the break period. Mata [29] conducted a comparative study on two statistical models named multiple linear regression (MLR) and neural network (NN) for monitoring of Alto Rabagão dam performance in the environmental condition. Seyed-Kolbadi et al. [30] evaluated the stability of the Boostan earth dam by processing its long-term performance and interpreting the measured data.All the previous research assumed that the modal parameters of the dam are known. In fact, they directly used the modal parameter obtained from Finite Element (FE) analysis to identify the structural damage. It seems that there is a research gap in identifying the damage that occurred in dams especially when the modal parameters are calculated through the system identification process. In the light of previous research, in this paper, a combined method based on wavelet and Hilbert–Huang transform is proposed to identify the structural dynamics of the dam. In this regard, four commonly occurred damage scenarios, including damage in the neck, damage in the toe, simultaneous damage in the neck and the toe, and damage in the lifting joints of the dam, are considered to verify the capability of the proposed method. An optimization technique based on minimizing the nondiagonal entries of the MAC matrix is implemented to optimize the location of sensors needed to install at the body of the dam. Next, by adding 10% noise to the structural response under the ambient excitation, the dynamic properties of the structure are obtained for predefined locations, and the cubic-spline technique is used to make the mode shapes differentiable. In the following, damage location is estimated through employing CWT on the difference between the curvature of the mode shapes of the intact and damaged structure.


## 1.1. Structural Damage Detection

Damage detection methods could be divided into two general categories, model-based and data-driven methods. In the former approach, a system identification algorithm should be employed for identifying the dynamic properties of the structure. Then, the structural health condition will be determined by capturing the variation of these dynamic properties. On the other hand, data-driven methods apply features extracted from the vibration response of the structure through a statistical learning algorithm to predict the damage [8]. Cui et al. [9] proposed a damage diagnosis procedure based on the natural excitation technique to identify the structural damage under ambient conditions. The authors combined the strain response with the Eigensystem realization algorithm to identify the strain modal parameters of a structure under ambient vibration. Montazer and Seyedpoor [10] presented a damage index based on flexibility index to detect the location of damage in truss systems. They implemented their proposed method on two numerical models, including a planar truss having 31 elements and a 47-bar planar power line tower. Lee and Eun [11] provided a static damage detection procedure by comparing the variation of stress distribution before and after damage. They showed the validity of their method utilizing two truss structures. Ding et al. [12] developed a method based on Modal Strain Energy (MSE) to identify the damage in the grid structures. They employed an updated mode shape expansion method to ensure that the modal shape obtained from the reference baseline model is reliable. Liu et al. [13] proposed a method based on Stationary Wavelet Transform (SWT) for the purpose of localization and severity detection of cracks in the cantilever beams. Montejo [14] evaluated different vibration-based damage diagnosis methods for the purpose of detecting the damage in the structures under ambient excitation. They concluded that the method based on Continuous Wavelet Transform (CWT) is much more efficient than uncovering the high-frequency spikes of the structural response obtained by Hilbert–Huang transform, high-pass filtering, or Discrete Wavelet Transform (DWT). Pu et al. [15] presented the use of Frequency Response Functions (FRFs) along with the model updating theory for identifying the damage occurred in concrete beams. They formulated an optimization algorithm to set the analytical FRFs from a benchmark finite element model with those obtained through the experimental response. Worden et al. [16] implemented a statistical approach based on Mahalanobis squared distances for healthy and damaged structures. This method was actually a novelty detection problem where a damaged structure should be diagnosed from the healthy cases. Ditommaso et al. [17] presented a method for damage detection of the framed-structures subjected to earthquake motions by evaluating the variation of modal curvature in the natural frequencies. They validated their method by using numerical models and implementing some experiments. Zhang et al. [18] proposed an innovative damage index according to the macrostrain modal shapes for the purpose of detecting the potential damage in steel bridges. Yazdanpanah et al. [19] presented a damage indicator based on modal information for damage detection in beam-like structures. Their method was founded on three factors, namely, the mode shape, the slope of the mode shape, and the curvature of the mode shape. Figueiredo et al. [20] compared four different machine learning techniques by applying autoassociative neural network, factor analysis, Mahalanobis distance, and singular value decomposition to identify the structural damage. They have implemented their study on an experimental base-excited three-story frame structure. Ghiasi et al. [21] proposed a least square support vector machine based on the thin-plate spline Littlewood–Paley wavelet kernel function to identify the structural damage. Ghannadi and Kourehli [22] explored the efficiency of the moth flame optimization algorithm to identify the location and severity of structural damage. They used the Eigenvalues obtained from the equation of motion along with Modal Assurance Criteria (MAC) flexibility as damage indices. Nguyen et al. [23] presented a novel framework based on the transmissibility function fed to ANNs to identify the structural damage in bridges. They verified their method using data measured from an existing bridge in Taiwan.


## 1.2. Structural Health Monitoring of Dams

Wang and He [24] presented a methodology based on the reduction of natural frequencies to identify the crack location in arch dams. They proposed a statistical neural network to detect the crack by measuring the reductions of natural frequencies. Turker et al. [25] explored a method based on model updating to identify the global damage that occurred in the arch dams. They have validated their method on an experimental model under ambient excitation. Pirboudaghi et al. [26] proposed a multistage method based on the extended finite element method (XFEM) and continuous wavelet transform (CWT) to detect the seismic cracks of concrete dams. Sevim et al. [27] implemented an ambient vibration test to identify the dynamic properties of a real concrete arch dam, including modal frequencies, damping, and mode shapes. Bianchi and Bremen [28] proposed various procedures to capture the long-term performance of Ferden and Roggiasca dams in Switzerland. In this regard, they captured variations of the dam behavior before the snowfall and after the break period. Mata [29] conducted a comparative study on two statistical models named multiple linear regression (MLR) and neural network (NN) for monitoring of Alto Rabagão dam performance in the environmental condition. Seyed-Kolbadi et al. [30] evaluated the stability of the Boostan earth dam by processing its long-term performance and interpreting the measured data.All the previous research assumed that the modal parameters of the dam are known. In fact, they directly used the modal parameter obtained from Finite Element (FE) analysis to identify the structural damage. It seems that there is a research gap in identifying the damage that occurred in dams especially when the modal parameters are calculated through the system identification process. In the light of previous research, in this paper, a combined method based on wavelet and Hilbert–Huang transform is proposed to identify the structural dynamics of the dam. In this regard, four commonly occurred damage scenarios, including damage in the neck, damage in the toe, simultaneous damage in the neck and the toe, and damage in the lifting joints of the dam, are considered to verify the capability of the proposed method. An optimization technique based on minimizing the nondiagonal entries of the MAC matrix is implemented to optimize the location of sensors needed to install at the body of the dam. Next, by adding 10% noise to the structural response under the ambient excitation, the dynamic properties of the structure are obtained for predefined locations, and the cubic-spline technique is used to make the mode shapes differentiable. In the following, damage location is estimated through employing CWT on the difference between the curvature of the mode shapes of the intact and damaged structure.


## 2. Models and Materials

### 2.1. Model Description

In this study, a well-known reservoir dam named “Koyna” was numerically modeled in Ansys finite element software, and further investigations were implemented on the response obtained through the numerical analysis. The following sections introduce the structural behavior and modeling procedure of the understudy dam. Fluid–Structure Interaction (FSI) is an important issue in both linear and nonlinear analyses of hydraulic structures. There are three main solutions to resolve this problem, including added mass method, Eulerian method, and Lagrangian method. The added mass method is related to additional mass considered in the system because of the existence of an accelerating or decelerating body that must move some volume of surrounding fluid as it moves through it. The Lagrangian method is related to individual particles, which calculates the trajectory of each particle independently. The Eulerian method is based on the concentration of particles and computes the general convection of the particles [31]. The coupled fluid–structure equation of motion can be easily solved by applying some proper boundary conditions for the reservoir object. Readers are encouraged to study the paper presented by Park et al. [32] for more information about fluid–structure interaction modeling.There are various types of foundation modeling procedures, such as rigid foundation and massless and massed foundation. The rigid foundation model disregards the interaction between the dam and its foundation. Actually, the stiffness and mass of the foundation are not considered in the global coupled equation of motion. Besides, the massless foundation model considers the flexibility of foundation in the FE model. In this problem, the foundation model must stretch to a large length so that its effects on stresses and strains of the dam approach to an enough small value [33]. Generally, in FE models, the meshing system includes solid elements to model the dam body and reservoir. In order to model the reservoir domain, pressure-based fluid elements were applied. The reservoir of the dam is modeled three-fifths of the dam’s height. The water is defined as a linearly compressible and small-amplitude material. The Helmholtz equation is implemented as governing equation in the reservoir medium, as shown in Figure 1. In addition, the boundary conditions considered for the modeling procedure are shown in this figure [34]. The coupled equation of the dam-foundation-reservoir system can be written as follows:(1)Ms0ρwQFSITMFU¨P¨+Cs00CFU˙P˙+Ks−QFSI0KFUP=FstrFacc,where K, C, and M are stiffness, damping, and mass matrices, Fstr is the summation of external forces and ground excitation on the solid domain, QFSI is fluid–structure coupling matrix, Facc is the summation of the forces that arise from ground motions on solid–fluid boundaries and the total acceleration acts on the other boundaries, P and U are the hydrodynamic pressures and displacements, and the subscripts “F” and “S” refer to fluid and structure, respectively.Figure 1
Numerical model of Koyna dam. (a) Schematic view of 3D model. (b) Mesh distribution at considered middle section. (c) Overall view of coupled dam–reservoir foundation system.
(a)(b)(c)


### 2.2. Case Study

Koyna gravity dam is founded in India with a crest length equal to 807 m and a height of 103 m. Its thickness at the base and the crest is 70.2 and 12.1 m, respectively. The 1967 Koynanagar event imposed intensive damages to the dam, including horizontal cracks on the upstream and downstream faces, as a result of which the slope of the downstream face changes significantly [35]. Table 1 presents the mechanical properties of the concrete and the acoustic properties used to model the water dynamics.Table 1
Mechanical properties of the dam body and reservoir water.
Model partCharacteristicUnitStaticDynamicDam bodyModulus of elasticityGPa31.035.7Poisson’s ratio—0.020.14DensityKg/m326432643Reservoir waterSonic velocitym/s——DensityKg/m310001000According to the large dimensions in the cross-stream direction, the assumption of plane strain is applicable for gravity dams. Therefore, these structures are usually modeled and analyzed as 2D structures [36–38].


### 2.3. Modal Analysis

By discarding the damping term in the coupled equation of motion, one can derive the Eigenvalues of the undamped FE model as defined by Hariri-Ardebili and Saouma [33]. Since the mass matrix and stiffness matrix are unsymmetric, the Lanczos formulation is implemented to solve this equation. This technique uses a secondary variable and the undamped coupled equation is written as [39](2)KS−QFSI0−QFSIT−MFKF0KFT0−MS0000000KFUPΘ=000,where Θ is equal to 1/ω2P. Readers are encouraged to study the paper presented by Hariri-Ardebili and Saouma [33] for more information about FSI and modal analysis of the dam.


## 2.1. Model Description

In this study, a well-known reservoir dam named “Koyna” was numerically modeled in Ansys finite element software, and further investigations were implemented on the response obtained through the numerical analysis. The following sections introduce the structural behavior and modeling procedure of the understudy dam. Fluid–Structure Interaction (FSI) is an important issue in both linear and nonlinear analyses of hydraulic structures. There are three main solutions to resolve this problem, including added mass method, Eulerian method, and Lagrangian method. The added mass method is related to additional mass considered in the system because of the existence of an accelerating or decelerating body that must move some volume of surrounding fluid as it moves through it. The Lagrangian method is related to individual particles, which calculates the trajectory of each particle independently. The Eulerian method is based on the concentration of particles and computes the general convection of the particles [31]. The coupled fluid–structure equation of motion can be easily solved by applying some proper boundary conditions for the reservoir object. Readers are encouraged to study the paper presented by Park et al. [32] for more information about fluid–structure interaction modeling.There are various types of foundation modeling procedures, such as rigid foundation and massless and massed foundation. The rigid foundation model disregards the interaction between the dam and its foundation. Actually, the stiffness and mass of the foundation are not considered in the global coupled equation of motion. Besides, the massless foundation model considers the flexibility of foundation in the FE model. In this problem, the foundation model must stretch to a large length so that its effects on stresses and strains of the dam approach to an enough small value [33]. Generally, in FE models, the meshing system includes solid elements to model the dam body and reservoir. In order to model the reservoir domain, pressure-based fluid elements were applied. The reservoir of the dam is modeled three-fifths of the dam’s height. The water is defined as a linearly compressible and small-amplitude material. The Helmholtz equation is implemented as governing equation in the reservoir medium, as shown in Figure 1. In addition, the boundary conditions considered for the modeling procedure are shown in this figure [34]. The coupled equation of the dam-foundation-reservoir system can be written as follows:(1)Ms0ρwQFSITMFU¨P¨+Cs00CFU˙P˙+Ks−QFSI0KFUP=FstrFacc,where K, C, and M are stiffness, damping, and mass matrices, Fstr is the summation of external forces and ground excitation on the solid domain, QFSI is fluid–structure coupling matrix, Facc is the summation of the forces that arise from ground motions on solid–fluid boundaries and the total acceleration acts on the other boundaries, P and U are the hydrodynamic pressures and displacements, and the subscripts “F” and “S” refer to fluid and structure, respectively.Figure 1
Numerical model of Koyna dam. (a) Schematic view of 3D model. (b) Mesh distribution at considered middle section. (c) Overall view of coupled dam–reservoir foundation system.
(a)(b)(c)


## 2.2. Case Study

Koyna gravity dam is founded in India with a crest length equal to 807 m and a height of 103 m. Its thickness at the base and the crest is 70.2 and 12.1 m, respectively. The 1967 Koynanagar event imposed intensive damages to the dam, including horizontal cracks on the upstream and downstream faces, as a result of which the slope of the downstream face changes significantly [35]. Table 1 presents the mechanical properties of the concrete and the acoustic properties used to model the water dynamics.Table 1
Mechanical properties of the dam body and reservoir water.
Model partCharacteristicUnitStaticDynamicDam bodyModulus of elasticityGPa31.035.7Poisson’s ratio—0.020.14DensityKg/m326432643Reservoir waterSonic velocitym/s——DensityKg/m310001000According to the large dimensions in the cross-stream direction, the assumption of plane strain is applicable for gravity dams. Therefore, these structures are usually modeled and analyzed as 2D structures [36–38].


## 2.3. Modal Analysis

By discarding the damping term in the coupled equation of motion, one can derive the Eigenvalues of the undamped FE model as defined by Hariri-Ardebili and Saouma [33]. Since the mass matrix and stiffness matrix are unsymmetric, the Lanczos formulation is implemented to solve this equation. This technique uses a secondary variable and the undamped coupled equation is written as [39](2)KS−QFSI0−QFSIT−MFKF0KFT0−MS0000000KFUPΘ=000,where Θ is equal to 1/ω2P. Readers are encouraged to study the paper presented by Hariri-Ardebili and Saouma [33] for more information about FSI and modal analysis of the dam.


## 3. Methodology

### 3.1. Proposed Framework

This study presents a cost-effective method for the purpose of identifying dynamic properties and damage location of the dams. The optimum locations for sensor placements were determined through minimizing the nondiagonal entries of the MAC matrix. In consequence, the dynamic properties of the dam were determined by processing the output signals captured from the structure under ambient vibration excitation. The dynamic properties of a structure will change by causing a damage to the body of the dam. The proposed method is based on the difference between curvature of the mode shapes in the healthy and defective structures. Eventually, by applying the continuous wavelet transform (CWT) on the residual of the curvature of the mode shapes, the location of the damage is estimated. Figure2 shows the general view of the proposed framework.Figure 2
Main algorithm of the proposed method.


### 3.2. Theoretical Background of Continuous Wavelet Transform

Since the foundation of the proposed method is based on continuous wavelet transform, it is essential to present a brief theoretical background for the readers who are less familiar with the topic. A wavelet is a functionψt∈L2ℜ with zero mean and defined as follows [40, 41]: (3)ψa,bt=a−0.5ψt−ba,a>0,where L2ℜ is the set of integrable functions, a and b are scale and translation parameters, respectively, and a−0.5 is used to ensure that ψa,b is equal to unity. Also, ψ is the wavelet function (the mother wavelet), which must satisfy the admissibility condition(4)Cψ=∫−∞+∞ψ^ω2ωdω<∞,where ω is the frequency, and ψ^ω is the Fourier transform of ψ. The continuous wavelet transform of xt is the inner product of the family of wavelets ψa,bt with the signal xt as follows:(5)FWa,b=<xt,ψa,bt>=∫−∞+∞xtaψ∗t−badt,where ψ∗ is the complex conjugate of ψ. This condition guarantees the workability of the signal to be used as the mother wavelet in continuous wavelet transform. In this study, two different mother wavelets were used, i.e., Morlet and Daubechies. Complex Morlet wavelet has a close form to the following equation and applied to identify the modal properties of the dam:(6)ψt=1πe−t2eiω0t,where t is the time, and ω0 is the central frequency of the wavelet. By applying a as dilation and b as translation, a son wavelet is written as follows:(7)ψa,bt=1aπe−t−b/a2eiω0t−b/a,On the other hand, Daubechies wavelets are defined by computing the running average and difference of scalar products with wavelets and scaling signals. As the Daubechies wavelets apply overlapped windows, the high-frequency coefficient spectrum presents all high-frequency variations. It can be used as an appropriate wavelet to identify the high-frequency jumps resulting from damage in the response signal of the structures [42]. Therefore, the Daubechies mother wavelet was implemented to identify the structural damage that occurred in the dam.


### 3.3. Mathematics of Cubic Spline

The cubic-spline interpolation was presented to the best curve fitting between some discrete values. The cubic spline includes a set of weights connected to a plane surface where the points are connected. A flexible strip is curved through each of these weights, which forms a straight and smooth curve. The foundation of this interpolation method is based on curve fitting of piecewise-defined function with the below form:(8)hx=h1x,ifx1≤x<x2,h2x,ifx2≤x<x3,…hnx,ifxn−1≤x<xn,where hi is a third polynomial function. The cubic spline must have the below conditions:(1)
The piecewise-defined function should interpolate all the points.(2)
hx,h′x, and h″x should be continuous in x1,xn; the cubic spline is widely applied to determine the slope (or cumulative variation) of a finite number of variables over an interval [43].


## 3.1. Proposed Framework

This study presents a cost-effective method for the purpose of identifying dynamic properties and damage location of the dams. The optimum locations for sensor placements were determined through minimizing the nondiagonal entries of the MAC matrix. In consequence, the dynamic properties of the dam were determined by processing the output signals captured from the structure under ambient vibration excitation. The dynamic properties of a structure will change by causing a damage to the body of the dam. The proposed method is based on the difference between curvature of the mode shapes in the healthy and defective structures. Eventually, by applying the continuous wavelet transform (CWT) on the residual of the curvature of the mode shapes, the location of the damage is estimated. Figure2 shows the general view of the proposed framework.Figure 2
Main algorithm of the proposed method.


## 3.2. Theoretical Background of Continuous Wavelet Transform

Since the foundation of the proposed method is based on continuous wavelet transform, it is essential to present a brief theoretical background for the readers who are less familiar with the topic. A wavelet is a functionψt∈L2ℜ with zero mean and defined as follows [40, 41]: (3)ψa,bt=a−0.5ψt−ba,a>0,where L2ℜ is the set of integrable functions, a and b are scale and translation parameters, respectively, and a−0.5 is used to ensure that ψa,b is equal to unity. Also, ψ is the wavelet function (the mother wavelet), which must satisfy the admissibility condition(4)Cψ=∫−∞+∞ψ^ω2ωdω<∞,where ω is the frequency, and ψ^ω is the Fourier transform of ψ. The continuous wavelet transform of xt is the inner product of the family of wavelets ψa,bt with the signal xt as follows:(5)FWa,b=<xt,ψa,bt>=∫−∞+∞xtaψ∗t−badt,where ψ∗ is the complex conjugate of ψ. This condition guarantees the workability of the signal to be used as the mother wavelet in continuous wavelet transform. In this study, two different mother wavelets were used, i.e., Morlet and Daubechies. Complex Morlet wavelet has a close form to the following equation and applied to identify the modal properties of the dam:(6)ψt=1πe−t2eiω0t,where t is the time, and ω0 is the central frequency of the wavelet. By applying a as dilation and b as translation, a son wavelet is written as follows:(7)ψa,bt=1aπe−t−b/a2eiω0t−b/a,On the other hand, Daubechies wavelets are defined by computing the running average and difference of scalar products with wavelets and scaling signals. As the Daubechies wavelets apply overlapped windows, the high-frequency coefficient spectrum presents all high-frequency variations. It can be used as an appropriate wavelet to identify the high-frequency jumps resulting from damage in the response signal of the structures [42]. Therefore, the Daubechies mother wavelet was implemented to identify the structural damage that occurred in the dam.


## 3.3. Mathematics of Cubic Spline

The cubic-spline interpolation was presented to the best curve fitting between some discrete values. The cubic spline includes a set of weights connected to a plane surface where the points are connected. A flexible strip is curved through each of these weights, which forms a straight and smooth curve. The foundation of this interpolation method is based on curve fitting of piecewise-defined function with the below form:(8)hx=h1x,ifx1≤x<x2,h2x,ifx2≤x<x3,…hnx,ifxn−1≤x<xn,where hi is a third polynomial function. The cubic spline must have the below conditions:(1)
The piecewise-defined function should interpolate all the points.(2)
hx,h′x, and h″x should be continuous in x1,xn; the cubic spline is widely applied to determine the slope (or cumulative variation) of a finite number of variables over an interval [43].


## 4. System Identification

### 4.1. Input Excitation

The input excitation considered to simulate the ambient vibration excitation to the dam was a Gaussian white noise signal. White noise is a random signal enjoying the same intensity at different frequencies, which provides a constant power spectral density. Therandn function of MATLAB [44] was utilized to generate a white noise signal. The total time of the input signal was 1200 s, with the time step of 0.01 s. It should be noted that the input signal was applied to the whole structure of the dam for simulating the real condition. The dynamic characterization of the dam was performed by implementing ambient vibration analysis. It should be noted that the key reason for using such a long-duration signal was to reach a stationary condition, which was essential for proper ambient vibration analysis. Figure 3 illustrates the input excitation subjected to the dam.Figure 3
Input excitation signal.


### 4.2. Sensor Location Optimization

One of the main objectives of the present study was to optimize the location of sensors. Therefore, it was necessary to know how to place these instruments in such a way that the maximum agreement between the numerical modal shapes and those obtained through the finite number of sensors be obtained. In order to achieve an economical plan, as well as effective coverage of the whole dam body in terms of detecting the probable damage, the team decided to place only six sensors in this study. The following steps were implemented to determine the optimal locations:(i)
First, the height of the dam was divided into 60 locations, where it is possible to place the sensors.(ii)
Second, the mode shapes of the dam were determined for all the possible nodes by applying modal analysis. Figure4 illustrates the first three mode shapes of the dam for 60 possible nodes.(iii)
Third, six intended sensors should be placed at the possible nodes to achieve the maximum linear independency. In fact, the instrumentation should be conducted so that the nondiagonal entries of the MAC matrix converge to their minimum values [45]. To do so, the first sensor was placed at the summit of the dam. The next sensor was placed at one of the 59 remaining locations so that the nondiagonal entries of the MAC matrix became minimum.(iv)
This operation was repeated for other remained sensors until the best agreement was reached. Finally, the optimal location of the sensors was determined, as shown in Figure5.Figure 4
Mode shapes of the dam according to 60 possible nodes for sensor placement.Figure 5
Optimal sensor location.


### 4.3. Free Vibration

In order to extract information about the dynamic properties of a structure from acceleration sensors, the vibration response of those structures must be measured. To this end, in this study, the random decrement technique (RDT) [46] was used. The concept behind the RDT is that by extracting a large number of time segments from the generated signal and averaging such signals, the random part of the response will continue to disappear from the system and only the response to the initial conditions will remain. For example, Figures 6(a) and 6(b) show the response of sensor no. 5 under ambient vibration excitation and corresponding free vibration obtained through RDT.Figure 6
(a) Response recorded by sensor no. 5. (b) Free vibration obtained by RDT.
(a)(b)


### 4.4. Mode Decomposition

The Continuous Wavelet Transform (CWT) of the free-vibration response gives time and scales coefficients of the wavelet coefficient contours. Considering this, a complex Morlet wavelet with a central frequency of 4.57 was used. As an example, Figure7 depicts an instance of wavelet coefficients obtained for sensors no. 1 and no. 6. In the next step, the free vibration associated with each vibrating mode can be achieved by a plotting summation of the wavelet coefficient at a point where the wavelet coefficients were maximum (maximum of contours). In other words, one should find the local maximum of the wavelet contours presented in Figure 7. Then, by applying the complex Morlet wavelet to the signal obtained from these peak values and plotting the imaginary parts of the wavelet coefficient, free vibration of each mode will be obtained according to each peak value.Figure 7
Wavelet coefficient contours. (a) Sensor no. 1; scale: 1,200. (b) Sensor no. 5; scale: 1,250.
(a)(b)


### 4.5. Modal Damping and Natural Frequencies

The Hilbert transform of the free vibration was then calculated to determine the natural frequency and modal damping ratio of each mode. Figure8 displays an algorithm for computing the dynamic properties. The steps in this algorithm are as follows:(1)
Hilbert transform was applied to the free vibration of each mode(2)
The phase and amplitude of the signal obtained in the previous step was computed(3)
The rate of change of the amplitude and the phase signal are considered as−ζiωi andωDi, respectively(4)
Finally, the fundamental concepts from dynamics of structures were used to findωi and ζiFigure 8
Algorithm for computingζi and ωi.For example, in an attempt to calculate the first and second natural frequencies, as well as modal damping ratios, it is required to employ the scaling of 60 in wavelet coefficient for the first mode and the scale of 145 in wavelet coefficient for the second mode as shown in Figure7. According to Figure 8, the slope of L21 in the linear region is ωD1, which is equal to ωD1=24.4rad (see Figure 9(a)). In addition, the slope of L11 in the linear region is −ζ1ω1, and as shown in Figure 9(b), the natural frequency and damping ratio of the first vibrating mode are obtained by implementing the following simple calculations:ωD1=244.0/10=24.4Hz and –ζ1ω1=−1.05⟶ω1=24.4/1−1.05/ω12⟶ω1=24.42rad⟶f1=3.88Hz and ζ1=1.05/24.42=0.043Figure 9
Curve fitting to extract the dynamic properties. First row for the first mode (i = 1) and second row for the second mode (i = 2).
(a)(b)(c)(d)It is worth noting that using the same procedure, the second mode properties can be calculated. As illustrated in Figure8, the slope of L22 in the linear region is ωD2, which is obtained as ωD2=45.9rad (see Figure 9(c)). Furthermore, the slope of L12 in the linear region is −ζ2ω2, and according to Figure 9(d), the natural frequency and damping ratio of the second mode can be computed by the following calculations:ωD2=459.0/10=45.9Hz and –ζ2ω2=−1.125⟶ω2=45.9/1−1.125/ω12⟶ω2=45.92rad⟶f2=7.31Hz and ζ2=1.125/45.92=0.025The first three natural frequencies and modal damping of the dam are calculated as shown in Table2.Table 2
Estimated modal parameters.
ModeFinite element modelWithout noiseWith 10% noiseFrequency (Hz)Damping ratio (%)Frequency (Hz)Damping ratio (%)Frequency (Hz)Damping ratio (%)13.674.03.884.33.924.227.954.07.312.57.453.1311.804.011.914.111.824.2


### 4.6. Noise Effect

In general, noise is an inevitable part of any measurement. Hence, there is a concept known as the Signal-to-Noise Ratio (SNR), which is defined as the ratio of signal power to the noise power, and commonly stated in decibels (dB). A ratio higher than 1 denotes a higher amount of signal than noise [41]:(9)SNR=σsignal2σnoise2,where σsignal2 is the variance of signal and σnoise2 represents the variance of the noise. For a more realistic simulation, an appropriate level of noise needs to be applied to the acceleration time history. In the present study, the structural response is polluted by 10% noise to evaluate the practical application of the proposed method. It is proved that this value is the upper bound value of noise used in the previous studies [47]. Table 2 also presents the effect of noise on the obtained modal parameter of the dam. The results show that the proposed method had a promising ability to identify the modal parameters of the dam in the presence of noise.


### 4.7. Mode Shape Identification

In an attempt to identify the structural mode shapes, a method presented by Yang et al. [43] is utilized. From the previous sections, it can be concluded that while the damping ratio and natural frequencies can be estimated by only one sensor, the acceleration time histories at all DOFs are required and should be measured to identify mode shapes. The absolute values of Eigenvectors (mode shapes) can be obtained using the below relation:(10)φpiφqi=expApi′t0−Aqi′t0,where φpi is the ith mode shape at the pth degree of freedom. Similar definition is valid for φqi. The parameters Api′t0 and Aqi′t0 are the values calculated from curve fitting of L2i for pth and qth degrees of freedom at time t0, respectively. It is notable that t0 is the time at the median value of L2i. McKinley and Levine [48] presented the below equation to calculate the sign of the mode shape at each degree of freedom:(11)φpi,q=θpi′t0−θqi′t0,where φpi,q denotes the difference between the phase value of two signals in pth and qth DOF and ith mode. In addition, θpi′t0 and θqi′t0 are the values obtained through curve fitting to L1i for pth and qth degree of freedom at time t0, respectively. As the previous equation defines only the absolute value of mode shapes, in order to determine the sign of the mode shape entries, the following relations are presented:(12)ifφpi,q=∓2mπ⟶φpiφqi>0and  ifφpi,q=∓2m+1π⟶φpiφqi<0.By applying the abovementioned relations, the mode shape entries can be determined on the basis of free-vibration response of each mode with significant modal participating mass ratio. Comparison between finite element analysis mode shapes and corresponding signal processing (SP) values is shown in Figure10. The results indicate that the mode shapes are identified with a similar trend. However, it is clear that the third mode shape was estimated with more errors as compared to the first two modes. In fact, variations in the mode values due to damage were less than the errors observed in identification process. This indeed makes this mode shape useless for identifying the structural damage.Figure 10
Difference between the mode shapes obtained by Finite Element (FE) analysis and System Identification (SI).In order to evaluate the efficiency of the cubic-spline technique, mode shapes obtained through the system identification (SI) process are compared with those modified by cubic spline (CS). Figure11 illustrates the effectiveness of cubic spline in smoothing and interpolating the mode shapes.Figure 11
Comparison between identified model shapes and interpolated ones.


## 4.1. Input Excitation

The input excitation considered to simulate the ambient vibration excitation to the dam was a Gaussian white noise signal. White noise is a random signal enjoying the same intensity at different frequencies, which provides a constant power spectral density. Therandn function of MATLAB [44] was utilized to generate a white noise signal. The total time of the input signal was 1200 s, with the time step of 0.01 s. It should be noted that the input signal was applied to the whole structure of the dam for simulating the real condition. The dynamic characterization of the dam was performed by implementing ambient vibration analysis. It should be noted that the key reason for using such a long-duration signal was to reach a stationary condition, which was essential for proper ambient vibration analysis. Figure 3 illustrates the input excitation subjected to the dam.Figure 3
Input excitation signal.


## 4.2. Sensor Location Optimization

One of the main objectives of the present study was to optimize the location of sensors. Therefore, it was necessary to know how to place these instruments in such a way that the maximum agreement between the numerical modal shapes and those obtained through the finite number of sensors be obtained. In order to achieve an economical plan, as well as effective coverage of the whole dam body in terms of detecting the probable damage, the team decided to place only six sensors in this study. The following steps were implemented to determine the optimal locations:(i)
First, the height of the dam was divided into 60 locations, where it is possible to place the sensors.(ii)
Second, the mode shapes of the dam were determined for all the possible nodes by applying modal analysis. Figure4 illustrates the first three mode shapes of the dam for 60 possible nodes.(iii)
Third, six intended sensors should be placed at the possible nodes to achieve the maximum linear independency. In fact, the instrumentation should be conducted so that the nondiagonal entries of the MAC matrix converge to their minimum values [45]. To do so, the first sensor was placed at the summit of the dam. The next sensor was placed at one of the 59 remaining locations so that the nondiagonal entries of the MAC matrix became minimum.(iv)
This operation was repeated for other remained sensors until the best agreement was reached. Finally, the optimal location of the sensors was determined, as shown in Figure5.Figure 4
Mode shapes of the dam according to 60 possible nodes for sensor placement.Figure 5
Optimal sensor location.


## 4.3. Free Vibration

In order to extract information about the dynamic properties of a structure from acceleration sensors, the vibration response of those structures must be measured. To this end, in this study, the random decrement technique (RDT) [46] was used. The concept behind the RDT is that by extracting a large number of time segments from the generated signal and averaging such signals, the random part of the response will continue to disappear from the system and only the response to the initial conditions will remain. For example, Figures 6(a) and 6(b) show the response of sensor no. 5 under ambient vibration excitation and corresponding free vibration obtained through RDT.Figure 6
(a) Response recorded by sensor no. 5. (b) Free vibration obtained by RDT.
(a)(b)


## 4.4. Mode Decomposition

The Continuous Wavelet Transform (CWT) of the free-vibration response gives time and scales coefficients of the wavelet coefficient contours. Considering this, a complex Morlet wavelet with a central frequency of 4.57 was used. As an example, Figure7 depicts an instance of wavelet coefficients obtained for sensors no. 1 and no. 6. In the next step, the free vibration associated with each vibrating mode can be achieved by a plotting summation of the wavelet coefficient at a point where the wavelet coefficients were maximum (maximum of contours). In other words, one should find the local maximum of the wavelet contours presented in Figure 7. Then, by applying the complex Morlet wavelet to the signal obtained from these peak values and plotting the imaginary parts of the wavelet coefficient, free vibration of each mode will be obtained according to each peak value.Figure 7
Wavelet coefficient contours. (a) Sensor no. 1; scale: 1,200. (b) Sensor no. 5; scale: 1,250.
(a)(b)


## 4.5. Modal Damping and Natural Frequencies

The Hilbert transform of the free vibration was then calculated to determine the natural frequency and modal damping ratio of each mode. Figure8 displays an algorithm for computing the dynamic properties. The steps in this algorithm are as follows:(1)
Hilbert transform was applied to the free vibration of each mode(2)
The phase and amplitude of the signal obtained in the previous step was computed(3)
The rate of change of the amplitude and the phase signal are considered as−ζiωi andωDi, respectively(4)
Finally, the fundamental concepts from dynamics of structures were used to findωi and ζiFigure 8
Algorithm for computingζi and ωi.For example, in an attempt to calculate the first and second natural frequencies, as well as modal damping ratios, it is required to employ the scaling of 60 in wavelet coefficient for the first mode and the scale of 145 in wavelet coefficient for the second mode as shown in Figure7. According to Figure 8, the slope of L21 in the linear region is ωD1, which is equal to ωD1=24.4rad (see Figure 9(a)). In addition, the slope of L11 in the linear region is −ζ1ω1, and as shown in Figure 9(b), the natural frequency and damping ratio of the first vibrating mode are obtained by implementing the following simple calculations:ωD1=244.0/10=24.4Hz and –ζ1ω1=−1.05⟶ω1=24.4/1−1.05/ω12⟶ω1=24.42rad⟶f1=3.88Hz and ζ1=1.05/24.42=0.043Figure 9
Curve fitting to extract the dynamic properties. First row for the first mode (i = 1) and second row for the second mode (i = 2).
(a)(b)(c)(d)It is worth noting that using the same procedure, the second mode properties can be calculated. As illustrated in Figure8, the slope of L22 in the linear region is ωD2, which is obtained as ωD2=45.9rad (see Figure 9(c)). Furthermore, the slope of L12 in the linear region is −ζ2ω2, and according to Figure 9(d), the natural frequency and damping ratio of the second mode can be computed by the following calculations:ωD2=459.0/10=45.9Hz and –ζ2ω2=−1.125⟶ω2=45.9/1−1.125/ω12⟶ω2=45.92rad⟶f2=7.31Hz and ζ2=1.125/45.92=0.025The first three natural frequencies and modal damping of the dam are calculated as shown in Table2.Table 2
Estimated modal parameters.
ModeFinite element modelWithout noiseWith 10% noiseFrequency (Hz)Damping ratio (%)Frequency (Hz)Damping ratio (%)Frequency (Hz)Damping ratio (%)13.674.03.884.33.924.227.954.07.312.57.453.1311.804.011.914.111.824.2


## 4.6. Noise Effect

In general, noise is an inevitable part of any measurement. Hence, there is a concept known as the Signal-to-Noise Ratio (SNR), which is defined as the ratio of signal power to the noise power, and commonly stated in decibels (dB). A ratio higher than 1 denotes a higher amount of signal than noise [41]:(9)SNR=σsignal2σnoise2,where σsignal2 is the variance of signal and σnoise2 represents the variance of the noise. For a more realistic simulation, an appropriate level of noise needs to be applied to the acceleration time history. In the present study, the structural response is polluted by 10% noise to evaluate the practical application of the proposed method. It is proved that this value is the upper bound value of noise used in the previous studies [47]. Table 2 also presents the effect of noise on the obtained modal parameter of the dam. The results show that the proposed method had a promising ability to identify the modal parameters of the dam in the presence of noise.


## 4.7. Mode Shape Identification

In an attempt to identify the structural mode shapes, a method presented by Yang et al. [43] is utilized. From the previous sections, it can be concluded that while the damping ratio and natural frequencies can be estimated by only one sensor, the acceleration time histories at all DOFs are required and should be measured to identify mode shapes. The absolute values of Eigenvectors (mode shapes) can be obtained using the below relation:(10)φpiφqi=expApi′t0−Aqi′t0,where φpi is the ith mode shape at the pth degree of freedom. Similar definition is valid for φqi. The parameters Api′t0 and Aqi′t0 are the values calculated from curve fitting of L2i for pth and qth degrees of freedom at time t0, respectively. It is notable that t0 is the time at the median value of L2i. McKinley and Levine [48] presented the below equation to calculate the sign of the mode shape at each degree of freedom:(11)φpi,q=θpi′t0−θqi′t0,where φpi,q denotes the difference between the phase value of two signals in pth and qth DOF and ith mode. In addition, θpi′t0 and θqi′t0 are the values obtained through curve fitting to L1i for pth and qth degree of freedom at time t0, respectively. As the previous equation defines only the absolute value of mode shapes, in order to determine the sign of the mode shape entries, the following relations are presented:(12)ifφpi,q=∓2mπ⟶φpiφqi>0and  ifφpi,q=∓2m+1π⟶φpiφqi<0.By applying the abovementioned relations, the mode shape entries can be determined on the basis of free-vibration response of each mode with significant modal participating mass ratio. Comparison between finite element analysis mode shapes and corresponding signal processing (SP) values is shown in Figure10. The results indicate that the mode shapes are identified with a similar trend. However, it is clear that the third mode shape was estimated with more errors as compared to the first two modes. In fact, variations in the mode values due to damage were less than the errors observed in identification process. This indeed makes this mode shape useless for identifying the structural damage.Figure 10
Difference between the mode shapes obtained by Finite Element (FE) analysis and System Identification (SI).In order to evaluate the efficiency of the cubic-spline technique, mode shapes obtained through the system identification (SI) process are compared with those modified by cubic spline (CS). Figure11 illustrates the effectiveness of cubic spline in smoothing and interpolating the mode shapes.Figure 11
Comparison between identified model shapes and interpolated ones.


## 5. Damage Detection

### 5.1. Proposed Method for Damage Detection

The presented method is based on applying CWT on the residual curvature of the mode shapes. Therefore, the curvature of the mode shapes was determined by using the following relation [49]: (13)ϕi″=ϕi+1−2ϕi+ϕi−1dh,where ϕi is the mode shape value at the ith height increment and dh is the corresponding height step. In general, a large enough number of data are needed to calculate the slope and the curvature of the mode shapes. Nevertheless, the present study aims to identify the structural properties as well as the damage by implementing only six sensors. This leads to only six data points for each mode shape, which makes it impossible to calculate the curvature of the signal. To address this issue, cubic spline technique was implemented in MATLAB [33] to create enough data for computing the rate of change as well as the curvature of the mode shapes. Next, the difference between the curvature of the interpolated mode shape in the healthy and damaged structures was obtained. Eventually, by applying CWTon the residual vector of the mode shape’s curvature, the location of damage was determined. As mentioned before, Daubechies wavelets were used to capture the changes occurring in the mode shape’s curvature because of the structural damage.


### 5.2. Damage Scenarios

The general potential failure mechanisms in concrete gravity dams are overstressing, sliding along the weak surface. The intense earthquake may create tensile cracks at the base or near the downstream slope change discontinuity in concrete gravity dams. In certain earthquake motions, depending on whether the earthquake is near-field or far-field, simultaneous tensile cracks are observed in the neck and toe of the dam. In addition, rigid sliding at the lifting joints of dam was reported as a potential damage scenario for near-field motions [33]. To validate the presented method, four damage scenarios were considered as follows:(i)
Scenario 1. Neck: stiffness of a small region in the neck of the dam is reduced to 10% of its value in the healthy condition, as shown in Figure 12(a). To understand the extent of damage, the ratio of the defective elements to the total number of elements of the dam can be used. For example, in this scenario, 0.25% of the dam’s elements have been injured (light damage on the neck). It is notable that the stiffness of elements is decreased by reducing the Young modulus of them.(ii)
Scenario 2. Toe: stiffness of a small region (0.35% damage) in the toe of the dam is decreased to 10% of its initial value in the intact condition, as shown in Figure 12(b) (light damage in the toe).(iii)
Scenario 3. Toe and neck together: stiffness of a small number of elements of the toe as well as the neck of the dam (0.6% damage) is reduced 90%, as illustrated in Figure 12(c) (light damage in the toe and neck).(iv)
Scenario 4. Lift joints: stiffness of a small number of elements in the lift joints (0.2% damage) are decreased to 10% of their initial value in the intact condition as shown in Figure 12(d).Figure 12
Damage scenarios. (a) Damage in the neck. (b) Damage in the toe. (c) Simultaneous damage. (d) Damage in lift joint elements.
(a)(b)(c)(d)


### 5.3. Results

Figure13 illustrates the residual mode shape curvature as well as wavelet coefficients corresponding to different scales vs. the dam height for the first scenario. As mentioned before, since the third mode shape was obtained with a bias with respect to its real vector, this mode shape was not efficient in damage detection scheme. Therefore, results are shown only for the first and second mode shapes.Figure 13
Residual value of mode-shape curvature and corresponding scalogram of CWT for the first scenario (damage in the neck). (a) First mode. (b) Second mode.
(a)(b)Results of CWT on the residual vector of curvature between the healthy and damaged dam for the first and second mode shape are shown in Figures13(a) and 13(b), respectively. It is clear that at about 70% of the dam’s total height, a mutation in the wavelet coefficients is revealed. In fact, this jump in the scalogram of wavelet coefficients occurs because of the existing damage in this location.Figure14 illustrates the result of damage detection for the second scenario. As it can be seen, the absolute value of wavelet coefficients is higher at the base of the dam as compared to other locations. Therefore, these spikes in the CWT at the base of the dam represent the location of the damage. Also, both mode shapes are applicable for identifying the location of the damage in this scenario.Figure 14
Residual value of mode-shape curvature and corresponding scalogram for the second scenario (damage in the toe). (a) First mode. (b) Second mode.
(a)(b)Figure15 presents the results of damage detection when the damage occurs in the toe and the neck of the dam simultaneously. Although the result of CWT on the first mode shape had not significant effect in detecting the damage in the toe of the dam, a jump in scalogram is observed in the height corresponds to the dam’s neck. Moreover, for the second mode, by applying the CWT on the residual mode shape curvature of the healthy and defective structure, there are two jumps in the absolute values of CWT where the damage has occurred. It can be concluded that the second mode shape is much sensitive to the damage scenario where multiple instances of damage have occurred in the body of the dam.Figure 15
Residual value of mode-shape curvature and corresponding scalogram for the third scenario (damage in the toe and the neck). (a) First mode. (b) Second mode.
(a)(b)Finally, Figure16 shows the result of damage detection for the fourth damage scenario where slight damage was created in the elements located in the lift joints. As it can be seen, a sharp jump is revealed at the height where the damage has occurred. Moreover, it is evident that both mode shapes (first and second modes) were efficient for detecting the location of the damage.Figure 16
Residual value of mode-shape curvature and corresponding scalogram for the fourth scenario (damage in the lift nodes elements). (a) First mode. (b) Second mode.
(a)(b)


## 5.1. Proposed Method for Damage Detection

The presented method is based on applying CWT on the residual curvature of the mode shapes. Therefore, the curvature of the mode shapes was determined by using the following relation [49]: (13)ϕi″=ϕi+1−2ϕi+ϕi−1dh,where ϕi is the mode shape value at the ith height increment and dh is the corresponding height step. In general, a large enough number of data are needed to calculate the slope and the curvature of the mode shapes. Nevertheless, the present study aims to identify the structural properties as well as the damage by implementing only six sensors. This leads to only six data points for each mode shape, which makes it impossible to calculate the curvature of the signal. To address this issue, cubic spline technique was implemented in MATLAB [33] to create enough data for computing the rate of change as well as the curvature of the mode shapes. Next, the difference between the curvature of the interpolated mode shape in the healthy and damaged structures was obtained. Eventually, by applying CWTon the residual vector of the mode shape’s curvature, the location of damage was determined. As mentioned before, Daubechies wavelets were used to capture the changes occurring in the mode shape’s curvature because of the structural damage.


## 5.2. Damage Scenarios

The general potential failure mechanisms in concrete gravity dams are overstressing, sliding along the weak surface. The intense earthquake may create tensile cracks at the base or near the downstream slope change discontinuity in concrete gravity dams. In certain earthquake motions, depending on whether the earthquake is near-field or far-field, simultaneous tensile cracks are observed in the neck and toe of the dam. In addition, rigid sliding at the lifting joints of dam was reported as a potential damage scenario for near-field motions [33]. To validate the presented method, four damage scenarios were considered as follows:(i)
Scenario 1. Neck: stiffness of a small region in the neck of the dam is reduced to 10% of its value in the healthy condition, as shown in Figure 12(a). To understand the extent of damage, the ratio of the defective elements to the total number of elements of the dam can be used. For example, in this scenario, 0.25% of the dam’s elements have been injured (light damage on the neck). It is notable that the stiffness of elements is decreased by reducing the Young modulus of them.(ii)
Scenario 2. Toe: stiffness of a small region (0.35% damage) in the toe of the dam is decreased to 10% of its initial value in the intact condition, as shown in Figure 12(b) (light damage in the toe).(iii)
Scenario 3. Toe and neck together: stiffness of a small number of elements of the toe as well as the neck of the dam (0.6% damage) is reduced 90%, as illustrated in Figure 12(c) (light damage in the toe and neck).(iv)
Scenario 4. Lift joints: stiffness of a small number of elements in the lift joints (0.2% damage) are decreased to 10% of their initial value in the intact condition as shown in Figure 12(d).Figure 12
Damage scenarios. (a) Damage in the neck. (b) Damage in the toe. (c) Simultaneous damage. (d) Damage in lift joint elements.
(a)(b)(c)(d)


## 5.3. Results

Figure13 illustrates the residual mode shape curvature as well as wavelet coefficients corresponding to different scales vs. the dam height for the first scenario. As mentioned before, since the third mode shape was obtained with a bias with respect to its real vector, this mode shape was not efficient in damage detection scheme. Therefore, results are shown only for the first and second mode shapes.Figure 13
Residual value of mode-shape curvature and corresponding scalogram of CWT for the first scenario (damage in the neck). (a) First mode. (b) Second mode.
(a)(b)Results of CWT on the residual vector of curvature between the healthy and damaged dam for the first and second mode shape are shown in Figures13(a) and 13(b), respectively. It is clear that at about 70% of the dam’s total height, a mutation in the wavelet coefficients is revealed. In fact, this jump in the scalogram of wavelet coefficients occurs because of the existing damage in this location.Figure14 illustrates the result of damage detection for the second scenario. As it can be seen, the absolute value of wavelet coefficients is higher at the base of the dam as compared to other locations. Therefore, these spikes in the CWT at the base of the dam represent the location of the damage. Also, both mode shapes are applicable for identifying the location of the damage in this scenario.Figure 14
Residual value of mode-shape curvature and corresponding scalogram for the second scenario (damage in the toe). (a) First mode. (b) Second mode.
(a)(b)Figure15 presents the results of damage detection when the damage occurs in the toe and the neck of the dam simultaneously. Although the result of CWT on the first mode shape had not significant effect in detecting the damage in the toe of the dam, a jump in scalogram is observed in the height corresponds to the dam’s neck. Moreover, for the second mode, by applying the CWT on the residual mode shape curvature of the healthy and defective structure, there are two jumps in the absolute values of CWT where the damage has occurred. It can be concluded that the second mode shape is much sensitive to the damage scenario where multiple instances of damage have occurred in the body of the dam.Figure 15
Residual value of mode-shape curvature and corresponding scalogram for the third scenario (damage in the toe and the neck). (a) First mode. (b) Second mode.
(a)(b)Finally, Figure16 shows the result of damage detection for the fourth damage scenario where slight damage was created in the elements located in the lift joints. As it can be seen, a sharp jump is revealed at the height where the damage has occurred. Moreover, it is evident that both mode shapes (first and second modes) were efficient for detecting the location of the damage.Figure 16
Residual value of mode-shape curvature and corresponding scalogram for the fourth scenario (damage in the lift nodes elements). (a) First mode. (b) Second mode.
(a)(b)


## 6. Conclusions

This study presents a system identification-based damage detection procedure for the purpose of identifying the damage that occurred in the body of dam structures. An optimization algorithm based on minimizing the nondiagonal entries of the MAC matrix is applied to estimate the best location of sensors, which leads to the best possible precision of damage detection procedure. The proposed method requires to implement only six sensors along the height of dam. The conducted method present in this study is founded on identifying the modal properties of the structure under ambient excitation. To simulate the practical condition (effect of arbitrary noise), the structural response obtained through FE analysis was polluted with 10% noise. By applying the random decrement technique, the free vibration of the dam was identified at each sensor location. Then, the dynamic characteristics of the dam are detected using the Hilbert–Huang method. The mode shapes of the structure were calculated based on the method proposed by Yang et al. [43]. Four potential damage scenarios were created, and mode shapes are obtained for both healthy and defective structures.Eventually, by implementing CWT to the residual vector of the mode shape’s curvature of the fully operational and damaged structures, the wavelet coefficients in terms of the scale parameter were calculated along the height of the dam. This research illustrated that the first, second, and fourth damage scenarios were easily detectable from the first two mode shapes. In addition, the damage scenario no. 3 was detected by only the second mode shape, and the first mode shape did not reveal the damage in the toe of the dam. Although the third mode shape of the structure was also identified, it did not contribute in detecting the location of the damage due to the high level of discrepancy compared to the ideal mode shape. The results show a promising and accurate framework in detecting the damage location even when the damage occurred in multiple locations. Also, regarding the high level of importance of such structures, it is recommended that the procedure presented in this study be accompanied by a visual inspection of experienced structural engineers. Finally, the authors propose additional research conducted by applying further experimental tests for dams with various size characteristics and element configurations.



---

*Source: 6653254-2021-02-26.xml*