Passive control and damage assessment with a pendulum tuned mass damper using mean square equations and second-order eigen-perturbation techniques

Abstract

This paper proposes a novel framework for damage detection and retuning of a Pendulum Tuned Mass Damper (PTMD) addressing effects of stochasticity via mean-square equations and establishing rapid change detection. The framework was made possible via a Stochastic Differential Equation (SDE) formulation of the equations of motion considering a two degree of freedom base and structure system equipped with a PTMD. The system was assumed to be excited at the base with Gaussian white noise and the solutions were computed using a stochastic integration scheme. A Second Order Eigen-Perturbation (SOPS) estimate is used to carry out detection of sudden stiffness reduction in the base floor of the structure. Once damage is detected, mean square equations on monitored displacements were used to retune the PTMD by adapting the length of the pendulum. It was found that damage could be detected rapidly and successive retunings were able to minimize the mean square displacement of the damaged structure and stabilise it over time.

1. Introduction

Tuning is a challenge in passive or semi-active vibration control (Wagg (2021)), especially for Tuned Mass Dampers (TMDs) Yang et al. (2022), impacting design of control and understanding natural frequency shifts of the host structure due to environmental conditions Ubertini et al. (2017) or damage Buckley et al. (2023). Retuning of TMDs over time is necessary for resilient control. Here, Pendulum Tuned Mass Dampers (PTMDs) are important as they conveniently retune by changing the pendulum length and (Roffel et al., 2013) for improved performance.

Retuning can occur in multiple ways. In passive approaches, they are often based on adjusting the pendulum frequency through a change in pendulum length. In Roffel et al. (2011), a frame is used to adjust the pendulum length as it moves along mass-supporting steel cables. Passive tuning is also employed in Wang et al. (2019) where a motor modified the pendulum length. Other solutions like adjustable stiffness devices (Nagarajaiah and Varadarajan, 2005) exist.

Solutions exists around forward and backward techniques. The forward techniques were amongst the first to originate from mathematical formulations of physical models (Hartog, 1956), allowing to derive control variables for TMD tuning via, amongst several others, frequency response functions (Gerges and Vickery, 2005a; Ghosh and Basu, 2007). Recently, a more generalized optimal parameter results for a range of systems and excitations were presented (Soltani and Deraemaeker, 2021). A thorough review is found in Elias and Matsagar (2017).

Optimal parameters for PTMDs led to control frameworks allowing for their adaptive retuning. For example, Roffel et al. (2011) investigated simulated detuning as change in system properties and Bertollucci Colherinhas et al. (2021) investigated optimal design for vibration reduction of offshore wind turbines. These studies investigated frequency domain objective functions assuming detuning, independent to damage identification.

An inverse or backward approach can also address damage detection via system identification (Farrar and Worden, 2007). Several data-driven methods exist for estimating system damage, including Empirical Mode Decomposition, analysis (Li et al., 2007; Xu and Chen, 2004) and monitoring nonlinear dynamics (Worden et al., 2008). Time domain analysis is capable of studying complex nonlinear dynamics like a Coriolis-damped pendulum (Viet, 2015). Recent interest in stochastic calculus show developments in design and simulations under stochastic excitation (Tiwari et al., 2021; Tripura et al., 2020a). Applying such methods to PTMDs indicate that mean-square equations can be determined a priori for operation as long as some assumptions are valid. Mean-square equations can also serve as a tuning objective function with satisfactory accuracy when compared to Monte Carlo-based tuning Mucchielli et al. (2022). Despite the advantages, the possibility of retuning PTMDs over time after damage not established, which we address here. There are some studies considering both damage and control, including improving seismic performance Shahraki et al. (2023) via the Park-Ang damage index (Khatibinia et al. (2024)). Combined damage and control for TMD was also observed for moving loads on beams Pakrashi et al. (2010); Chawla and Pakrashi (2022), via retuning with empirical wavelet transform Yuan et al. (2018) and in designs of Azizi et al. (2021) semi-active adjustable device for retuning Karami et al. (2019). However, they do not consider our current approach combining retuning, optimal control and also rapid damage assessment considering stochasticity. Conseuqently, this paper also allows existing devices in this sector to expand their abilities.

Passive control of a PTMD-equipped structure requires some insights to be effective, integrating with previously developed simulation tools, damage detection scenarios and required tuning. First, an integrated damage detection and control framework via this stochastic differential equation approach relies on damage detection metrics to assess relative severity and reasonably accurate estimate of when it occurs. Secondly, the framework should also be able to retune for optimal minimization of its mean-square displacement for the damage severity range considered. Finally, the framework must satisfy the requirement of rapid to close to real-time analysis.

This paper presents a hybrid forward-inverse framework of damage detection and retuning for a 2-story building model equipped with a PTMD on second floor. The system is simulated using the Ito-Taylor 1.5 stochastic integration scheme, with displacement and acceleration responses obtained for stochastic excitation. Damage detection abilities via a Second-Order Eigen-Perturbation (SOP) approach is demonstrated with associated Damage Sensitive Features (DSFs). The lower computational cost of this technique allows for rapid processing of acceleration responses and damage detection. Accuracy and computational costs of SOP are detailed in Mucchielli et al. (2020). While damage detection delivers essential information on the nature of damage, a PTMD becomes useless once a thereshold of damage is exceeded, requiring retuning.

To address this, Ito’s lemma is used to compute a priori mean square equations for each state (floor displacement and velocities, pendulum angle and angular velocity) for PTMD tuning. It is demonstrated, how minimization of second floor’s mean square displacement as a function of pendulum length and damping leads to optimal PTMD performance, allowing for retuning to be aligned to damage or degradation.

A key contribution of this paper is thus in introducing the stochastic framework for a PTMD, combining damage detection and retuning. Another contribution is the use of mean square estimates of displacements to understand PTMD retuning conditions, expanding its operational boundaries and performance. Finally, the paper demonstrates the efficacy of a combined detection, control and retuning approach through a two degree of freedom system with base excitation, which is a canonical form many real systems undergoing random vibration and requiring vibration control, including multistory buildings (Djedoui and Ounis (2022); Love et al. (2023); Taha (2021)), wind response mitigation (Stanikzai et al. (2022), pounding mitigation ( Djerouni et al. (2023)) and related design (Challagulla et al. (2022); Karakostas et al. (2024)) of passive and semi-active control devices Arrigan et al. (2011); Cruze et al. (2021); Suthar et al. (2022); Jangid (2022); Love et al. (2020)), even on curved spaces ( Panda et al. (2023)).

2. Methodology

The proposed framework has several modular parts. The word simulation indicates that the Ito-Taylor 1.5 numerical integration scheme was employed to obtain it (Roy and Rao, 2017):

• Initial design: Appropriate PTMD parameters are determined via derivation of the mean-square displacement equations of an Ito SDE, containing the

 – Derivation of the equations of motion of the PTMD system.

 – Derivation of mean-square equations for an equivalent Ito SDE.

 – Computation (through optimization of the mean-square) of the optimal length and auxiliary damping at pristine state.

 – Computation of a realistic mean-square threshold.

• Damage detection: A combination of forward (using mean-square equations) and inverse (using eigen-perturbation) solutions detects damage by

 – Forward problem (Mean-square module):

 * Mean-square displacement is monitored from the simulation. Several operations are performed to check whether it is significantly away from its statistically most probable region.

 – Inverse problem (Eigen-perturbation module):

 * Eigen-perturbation update is performed at each sample and damage sensitive features are monitored for damage detection.

2.1. System description

We consider a 2 degree of freedom (DOF) building based on Ramallo et al. (2004). The bottom DOF, referred to as the base supports a heavier structure equipped with a PTMD as presented in this paper. Both top and bottom floors are not rigid bodies. A schematic rendering presented in Figure 1. The base is isolated by a spring-dashpot system with spring constant k _b and damping c _b , the structural stiffness and damping are noted k _s and c _s and the pendulum mass and length are m and l, respectively. The constant h, determining the point of attachment of the auxiliary damper (with damping c), was fixed to 80 % of length l. Both floors have a damping ratio of 2 % of critical viscous damping, with a period of 2.5 s and 0.3 s for the base and structure, respectively. The displacement of the base and structure in the u-direction are x _b and x _s , respectively.

Figure 1.

Schematic of the system under study.

2.2. Equations of motion

The equations of motion for this system are derived from its Lagrangian. The pendulum mass position has independent coordinates,

$\begin{array}{ll}\hfill p_x& =x_s+l\sin \theta \hfill \end{array}$

(1)

$\begin{array}{ll}\hfill p_y& =-l\cos \theta \hfill \end{array}$

(2)

This leads to the following equations for the kinetic, potential and dissipative energies, respectively

$\begin{array}{ll}\hfill V& =\frac{1}{2}{m}_b{\dot{x}}_b^2+\frac{1}{2}{m}_s{\dot{x}}_s^2+\frac{1}{2}m[{\dot{x}}_s^2+2{\dot{x}}_{s}l\cos \theta \dot{\theta }+l^2{\dot{\theta }}^2]\hfill \\ \hfill E& =\frac{1}{2}{k}_b{x}_b^2+\frac{1}{2}{k}_s{(x_s-x_b)}^2-m{g}_{a}l\cos \theta \hfill \\ \hfill R& =\frac{1}{2}{c}_b{\dot{x}}_b^2+\frac{1}{2}{c}_s{({\dot{x}}_s-{\dot{x}}_b)}^2+\frac{1}{2}c{h}^2{\cos }^2\theta {\dot{\theta }}^2\hfill \end{array}$

(3)

where g _a is the acceleration due to gravity. The equations of motion can then be derived using the Lagrangian L = V − E and, considering inclusion of dissipative terms after the minimisation

$\frac{d}{dt}\frac{\partial L}{\partial \dot{\mathbf{q}}}-\frac{\partial L}{\partial \mathbf{q}}+\frac{\partial R}{\partial \dot{\mathbf{q}}}=\mathbf{Q}$

(4)

Where q = (x _b , x _s , θ) and Q is the sum of all forces acting on the system. This leads to,

$\begin{array}{l}\hfill \underset{\mathbf{M}}{\underbrace{\left(\begin{array}{ccc}\hfill m_b\hfill & \hfill 0\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill m_s+m\hfill & \hfill ml\cos \theta \hfill \\ \hfill 0\hfill & \hfill ml\cos \theta \hfill & \hfill m{l}^2\hfill \end{array}\right)}}\left(\begin{array}{c}\hfill {\ddot{x}}_b\hfill \\ \hfill {\ddot{x}}_s\hfill \\ \hfill \ddot{\theta }\hfill \end{array}\right)\\ \hfill +\underset{\mathbf{C}}{\underbrace{\left(\begin{array}{ccc}\hfill c_b+c_s\hfill & \hfill -c_s\hfill & \hfill 0\hfill \\ \hfill -c_s\hfill & \hfill c_s\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill c{h}^2{\cos }^2\theta \hfill \end{array}\right)}}\left(\begin{array}{c}\hfill {\dot{x}}_b\hfill \\ \hfill {\dot{x}}_s\hfill \\ \hfill \dot{\theta }\hfill \end{array}\right)\\ \hfill +\underset{\mathbf{K}}{\underbrace{\left(\begin{array}{ccc}\hfill k_b+k_s\hfill & \hfill -k_s\hfill & \hfill 0\hfill \\ \hfill -k_s\hfill & \hfill k_s\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 0\hfill \end{array}\right)}}\left(\begin{array}{c}\hfill x_b\hfill \\ \hfill x_s\hfill \\ \hfill \theta \hfill \end{array}\right)\\ \hfill =\underset{\mathbf{f}}{\underbrace{\left(\begin{array}{c}\hfill 0\hfill \\ \hfill 0\hfill \\ \hfill -m{g}_{a}l\sin \theta \hfill \end{array}\right)}}+\underset{{\mathbf{B}}_{\mathit{init}}\dot{\mathbf{w}}}{\underbrace{\left(\begin{array}{c}\hfill -m_b{b}_1\hfill \\ \hfill -m_s{b}_2\hfill \\ \hfill 0\hfill \end{array}\right)\left(\begin{array}{c}\hfill {\dot{w}}_1(t)\hfill \end{array}\right)}}\end{array}$

(5)

The noise term ${\dot{w}}_1(t)$ represents white noise, the generalized derivative of Wiener process w₁(t). Here white noise is modeled as a stationary process in a wide sense where the expected value is zero and autocorrelation function is σ²δ(t) where σ is a scalar. Note, b₁ and b₂ represent the excitation amplitude. The system of equations were cast into an SDE format analogous to that presented below,

$\begin{array}{ll}\hfill \mathrm{d}\mathbf{y}& ={\widehat{\mathbf{M}}}^{-1}(-\widehat{\mathbf{C}}\mathbf{y}-\widehat{\mathbf{K}}\mathbf{y}+\widehat{\mathbf{f}}+{\widehat{\mathbf{B}}}_{\mathit{init}}\dot{\mathbf{w}})\mathrm{d}t\hfill \\ \hfill \mathrm{d}\mathbf{y}& ={\widehat{\mathbf{M}}}^{-1}(-\widehat{\mathbf{C}}\mathbf{y}-\widehat{\mathbf{K}}\mathbf{y}+\widehat{\mathbf{f}})\mathrm{d}t+{\widehat{\mathbf{M}}}^{-1}{\widehat{\mathbf{B}}}_{\mathit{init}}\mathrm{d}\mathbf{w}\hfill \end{array}$

(6)

where $d\mathbf{w}=\dot{\mathbf{w}}dt$ , $\mathbf{y}=(x_b,{\dot{x}}_b,x_s,{\dot{x}}_s,\theta ,\dot{\theta })$ and $\widehat{\mathbf{M}},\widehat{\mathbf{C}},\widehat{\mathbf{K}},\widehat{\mathbf{f}},{\widehat{\mathbf{B}}}_{\mathit{init}}$ are the augmented versions of their namesakes established in equation (5). Gathering time-independent terms in equation (6) as vector a and matrix B results in,

$\mathrm{d}\mathbf{y}=\mathbf{a}\mathrm{d}t+\mathbf{B}\mathrm{d}\mathbf{w}$

(7)

Note that this system is base excited. The excitation amplitudes were chosen to be 2.5 % of the floor mass and the mass ratio between the structure and the pendulum was set to μ = m/m _s = 5%.

2.3. The forward problem solution

2.3.1. Mean square equations derivation

To derive mean square equations, linearisation of the equations of motion was done by reasonably assuming that the stroke angle is of limited amplitude. This small-angle approximation results in, up to first order,

$\begin{array}{l}\hfill \left(\begin{array}{ccc}\hfill m_b\hfill & \hfill 0\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill m_s+m\hfill & \hfill ml\hfill \\ \hfill 0\hfill & \hfill ml\hfill & \hfill m{l}^2\hfill \end{array}\right)\left(\begin{array}{c}\hfill {\ddot{x}}_b\hfill \\ \hfill {\ddot{x}}_s\hfill \\ \hfill \ddot{\theta }\hfill \end{array}\right)\\ \hfill +\left(\begin{array}{ccc}\hfill c_b+c_s\hfill & \hfill -c_s\hfill & \hfill 0\hfill \\ \hfill -c_s\hfill & \hfill c_s\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill c{h}^2\hfill \end{array}\right)\left(\begin{array}{c}\hfill {\dot{x}}_b\hfill \\ \hfill {\dot{x}}_s\hfill \\ \hfill \dot{\theta }\hfill \end{array}\right)\\ \hfill +\left(\begin{array}{ccc}\hfill k_b+k_s\hfill & \hfill -k_s\hfill & \hfill 0\hfill \\ \hfill -k_s\hfill & \hfill k_s\hfill & \hfill 0\hfill \\ \hfill 0\hfill & \hfill 0\hfill & \hfill 0\hfill \end{array}\right)\left(\begin{array}{c}\hfill x_b\hfill \\ \hfill x_s\hfill \\ \hfill \theta \hfill \end{array}\right)\\ \hfill =\left(\begin{array}{c}\hfill 0\hfill \\ \hfill 0\hfill \\ \hfill -m{g}_{a}l\theta \hfill \end{array}\right)+\left(\begin{array}{c}\hfill -m_b{b}_1\hfill \\ \hfill -m_s{b}_2\hfill \\ \hfill 0\hfill \end{array}\right)\left(\begin{array}{c}\hfill {\dot{w}}_1(t)\hfill \end{array}\right)\end{array}$

(8)

For algebraic ease, we define,

${\omega }_s=\sqrt{\frac{k_s}{m_s}}{\omega }_p=\sqrt{\frac{g_a}{l}}$

(9)

From equation (8), the derivation of mean square equations was based on Socha (2008) and Roberts and Spanos (2003), making use of Ito’s lemma to derive mean square derivatives for stationary stochastic excitation. For a scalar-valued function g of t and the stochastic process X _t , that is mathematically, g: = g(t, X _t ),

$\begin{array}{l}\hfill \mathrm{d}g=\\ \hfill \left(\frac{\partial g}{\partial t}+\sum _{i=1}^p{a}_i\frac{\partial g}{\partial y_i}+\frac{1}{2}\sum _{i=1}^n\sum _{j=1}^n\sum _{k=1}^m{B}_{ik}{B}_{jk}\frac{{\partial }^2}{\partial y_i\partial y_j}\right)\mathrm{d}t\\ \hfill +\sum _{i=1}^p\sum _{j=1}^n{B}_{ij}\frac{\partial g}{\partial y_i}\mathrm{d}{w}_j\end{array}$

(10)

where p and n are respectively the number of states in vector y and the number of Wiener processes involved in the system formulation. The value of g is taken for each of the unique combinations of y _i y _j , i, j = 1, …, p. Applying the expected value operator leads to a mean square derivative formulaton dE[y _i y _j ]/dt leading to 21 unique derivative equations. These were solved considering wide sense stationarity of each state as

$\begin{array}{c}\hfill \frac{\mathrm{d}E[y_i]}{\mathrm{d}t}=0,\frac{\mathrm{d}E[y_i{y}_j]}{\mathrm{d}t}=0\hfill \end{array}$

(11)

using symbolic computations in MATLAB^® Symbolic Toolbox.

2.3.2. Initial tuning

Initial tuning was performed as an optimization of the main structure’s mean square ${\sigma }_s^2$ as a function of the frequency ratio ξ and auxiliary damping ratio ζ through a performance index J as

$\begin{array}{cc}\hfill {\min }_{\xi ,\zeta }& J\hfill \\ \hfill s.t.& 0\le \xi \le 0.3\hfill \\ \hfill & 0\le \zeta \le 1\hfill \end{array}$

(12)

where,

$\begin{array}{ll}\hfill J& =\frac{{\sigma }_s^2}{{\sigma }_{\mathit{uncontrolled}}^2}\hfill \end{array}$

(13)

$\begin{array}{ll}\hfill \xi & =\frac{{\omega }_p}{{\omega }_s}\hfill \end{array}$

(14)

$\begin{array}{ll}\hfill \zeta & =\frac{c{h}^2}{2m{l}^2{\omega }_p}\hfill \end{array}$

(15)

and ${\sigma }_{\mathit{uncontrolled}}^2$ is the mean square of the uncontrolled structure.

Here the MATLAB^® fmincon function with the Interior Point algorithm was used, yielding an optimal point at ξ = 0.1211 and ζ = 0.0996. The performance index shown in Figure 2 reveals a best reduction in the mean-square response to J _opt = 30.22% compared to an uncontrolled system at this optimal point.

Figure 2.

Mean-square displacement along the x-axis as a function of frequency ratio and auxiliary damping ratio. The data tip shows the optimal point at which this system response is evaluated at the beginning of the analysis.

2.3.3. Determination of mean square thresholds

Damage was introduced as a sudden reduction in base stiffness, characterised by a damage index κ where κ = 0.25 is equivalent to 75 % stiffness reduction. Good estimation of large damages are often useful for a post-event assessment of a structure.

Handling damage detection using mean-square equations and subsequent retuning requires establishing detection thresholds considering natural variations of mean-square responses with stochastic excitation considered.

Several damage indices were mapped to their effect on the mean-square responses (Table 1), demonstrating a nonlinear relationship between both quantities. For example, a reduction of 5% in base stiffness would cause a 9% increase in the mean-square value. Several thresholds were tested but the threshold of 10% (of the current target mean-square value) was retained as the best compromise between sensitivity and robustness against false positives. Damage at or below 5% could still be detected due to small-scale variations in the mean-square response. Other relevant damage indices should also cause a mean-square increase of sufficient amplitude for damage detection.

Table 1.

Table of damage indices and their associated multiplication factor effect on a pristine mean-square response.

Damage index	0.25	0.30	0.40	0.50	0.60	0.70	0.80	0.90	0.95
Mean-square multiplication factor	12.34	10.11	7.18	5.22	3.73	2.57	1.73	1.23	1.09

2.4. The inverse problem solution

Damage detection via an inverse approach uses SOP to identify significant shift in data covariance of the eigenspace of the matrix at every time sample and subsequently compute damage sensitive features.

2.4.1. Second Order Eigen-Perturbation

Given a covariance matrix C _k at step k, one can update it with a vector of new samples X _k as,

${\mathbf{C}}_k=\frac{k-1}{k}{\mathbf{C}}_{k-1}+\frac{1}{k}{X}_k{X}_k^T$

(16)

Both covariance matrices can be expressed with their eigenvalue (Ω_(.)) and eigenvector matrices (E_(.)),

${\mathbf{E}}_{k}k{\mathbf{\Omega }}_k{\mathbf{E}}_k^T={\mathbf{E}}_{k-1}((k-1){\mathbf{\Omega }}_{k-1}+\delta \mathbf{\Omega }){\mathbf{E}}_{k-1}^T$

(17)

where $\delta \mathbf{\Omega }={\mathbf{E}}_{k-1}^T{X}_k{X}_k^T{\mathbf{E}}_{k-1}$ . Matrix ((k − 1)Ω_k−1 + δΩ) will increasingly approach and achieve diagonal dominance as k increases. The perturbation update, detailed in Mucchielli et al. (2020), allows to easily compute the eigenvectors and eigenvalues of this matrix (V _k and Λ _k respectively) and to link them to the update of Ω _k and E _k ,

${\mathbf{E}}_{k}k{\mathbf{\Omega }}_k{\mathbf{E}}_k^T={\mathbf{E}}_{k-1}{\mathbf{V}}_k{\mathbf{\Lambda }}_k{\mathbf{V}}_k^T{\mathbf{E}}_{k-1}^T$

(18)

${\mathbf{E}}_k={\mathbf{E}}_{k-1}{\mathbf{V}}_k$

(19a)

${\mathbf{\Omega }}_k=\frac{{\mathbf{\Lambda }}_k}{k}$

(19b)

This update rule makes E _k and Ω _k available at every k as the data is supplied. A sudden change in structural properties will cause a sudden eigenvalue shift or through rotation of the orthogonal basis described by the eigenvector matrix. Such modal indicators have been used before ((Farrar and Worden, 2007; Feeny, 2002)).

2.4.2. Damage detection features

An array of damage sensitive features were selected from their established performance in previous studies (Bhowmik et al., 2019a; Hot et al., 2012; Krishnan et al., 2018). They include

• Recursive Residual Error - 1 $({\chi }_{r{r}_1})$

$\begin{array}{l}\hfill {\chi }_{r{r}_1}(k)=\Vert X_k^{*}-{\mathbf{V}}_k^{1^T}{X}_k{\Vert }^2\\ \hfill X_k^{*}={\mathbf{V}}_{k-1}^1{\mathbf{V}}_k^{1^T}{X}_k\end{array}$

(20)

where ${\mathbf{V}}_{k-1}^1$ , ${\mathbf{V}}_k^1$ are the collection of eigenvectors whose corresponding eigenvalues represent more than 90 % of the variance at sample k.

• • Recursive Residual Error - 2 $({\chi }_{r{r}_2})$

${\chi }_{r{r}_2}(k)=\Vert X_k^{*}-X_k{\Vert }^2$

(21)

${\chi }_{r{r}_2}$ is the difference between ${\mathbf{V}}_{k-1}^1$ and ${\mathbf{V}}_k^1$ through the difference between the reconstructed sample and original data vector. A higher ${\chi }_{r{r}_2}$ denotes a shift of underlying orthogonal basis, and therefore the covariance of accelerations.

• • Recursive Residual Error - 3 $({\chi }_{r{r}_3})$

$\begin{array}{ll}\hfill {\chi }_{r{r}_3}(k)& =\Vert {\overline{X}}_k-{\widehat{X}}_k\Vert \hfill \end{array}$

(22)

$\begin{array}{ll}\hfill {\overline{X}}_k& ={\mathbf{V}}_k^1{X}_k\hfill \end{array}$

(23)

$\begin{array}{ll}\hfill {\widehat{X}}_k& ={\mathbf{V}}_{k-1}^1{X}_k\hfill \end{array}$

(24)

${\chi }_{r{r}_3}$ quantifies the difference between ${\mathbf{V}}_{k-1}^1$ and ${\mathbf{V}}_k^1$ . A higher ${\chi }_{r{r}_3}$ denotes shift in underlying orthogonal basis.

• • Recursive Residual Error - 4 $({\chi }_{r{r}_4})$

${\chi }_{r{r}_4}(k)=1-\Vert {\mathbf{V}}_{k-1}^1{\mathbf{V}}_k^{1^T}\Vert$

(25)

${\chi }_{r{r}_4}$ quantifies the difference between ${\mathbf{V}}_{k-1}^1$ and ${\mathbf{V}}_k^1$ through the rotation in the orthogonal basis they represent. A lower ${\chi }_{r{r}_4}$ denotes a shift in the basis.

• • Eigenvector change (EVC)

$\begin{array}{l}\hfill EVC(k)\\ \hfill =\mathrm{m}\mathrm{a}\mathrm{x}\left(\sqrt{\mathrm{d}\mathrm{i}\mathrm{a}\mathrm{g}\left(({\mathbf{V}}_k^1-{\mathbf{V}}_{k-1}^1){({\mathbf{V}}_k^1-{\mathbf{V}}_{k-1}^1)}^T\right)}\right)\end{array}$

(26)

• • Mahalanobis distance (m _d )

$m_d(k)=\sqrt{(X_k-{\mathit{\mu }}_k){\mathbf{C}}_k^{-1}(X_k-{\mathit{\mu }}_k)}$

(27)

Where μ _k and C _k are the mean vector and covariance matrix respectively of the acceleration signal for kth sample, quantifying the position of each data sample within the distribution of the acceleration signal. Statistical outliers will have a greater Mahalanobis distance than a sample consistent with the signal’s distribution.

• • Maximum eigenvalue (λ _max )

${\lambda }_{\mathit{max}}(k)=\mathrm{m}\mathrm{a}\mathrm{x}({\mathbf{\Sigma }}_k)$

(28)

which will shift from damage incidence.

• • First eigenvalue (λ₁)

${\lambda }_1(k)={\mathbf{\Sigma }(k)}_1$

(29)

where Σ is the eigenvalue matrix corresponding to the covariance matrix of the collected acceleration responses. This eigenvalue corresponds to the frequency of oscillation of the base floor and can change significantly due to a stiffness change. The sensitivity of these features to 20% stiffness change (185.6 kN/m here) are shown in Figure 3. Each feature was further processed to highlight the change more significantly. For a generalized feature q this led to

$q_n={\mathrm{l}\mathrm{o}\mathrm{g}}_{10}\left(\frac{1}{500}\sum _{i=1}^{500}{\left.\left(q_{n-i}-q_{n-(i+1)}\right)\right)}^2\right)$

(30)

where the difference equation gives information about the rate of variation of q. The averaging suppresses high frequency variations from stochasticity of the signal. The logarithm enables a clearer representation of the quantities over time.

Figure 3.

Plots of the discussed damage detection features for a simulation of the undamaged and damaged system under study over 100 s.

It is observed that ${\chi }_{r{r}_2}$ , ${\chi }_{r{r}_3}$ , the EVC and λ₁ are sensitive to the stiffness reduction at 60s. Eigenvalue features are prone to false positives, and consequently one eigenvalue and one eigenvector feature with the best performance are used ( ${\chi }_{r{r}_2}$ and λ₁) and after the processing through equation (30), named ν and λ respectively.

2.5. Proposed framework for damage detection

The damage detection framework can be described as:

1. For each time step, the simulation module produces a new sample of displacement, velocity and acceleration for each of the three DOFs of this system using design parameters obtained during tuning phase. Exactly after 500s of simulation, the damage index κ was modified to reflect damage severity.

2. Eigen-perturbation module:

 (a) In damage detection module, nothing occurs as long as the number of initial acceleration samples is not reached. These are used to calculate the initial covariance matrix for SOP. Once the initial number of 100 samples is reached (n = 100.), initial covariance matrix and eigenspace are computed.

 (b). If the n > 100, each new sample is directly fed into the SOP to update the eigenspace. Damage sensitive features λ and ν are computed using the updated eigenspace.

 (c) If within an interval of 10 samples (equivalent to 0.1s) both features increase more than 5%, damage is flagged with its time of occurrence. If not, the module waits for the next acceleration data from simulation

• 3. Mean square analysis module:

• (b) Mean-square response of displacement is computed using a sample-wise update (Pébay et al., 2016).

• (c) To eliminate high frequency changes in mean-square response, a moving-average filter is applied with a window of 1000 samples (10s).

• (d) This averaged mean-square is used for graphical interface display and fed to a derivative filter.

• (e) If average mean-square exceeds the determined threshold (here 10%) a damage flag is raised and the damage instant is recorded. If not, the module waits for the next displacement sample.

• (f) If a damage event is flagged, the framework next inspects the averaged mean-square to determine whether it reached a plateau and steady enough to correspond to the asymptotic, stationary value of the system’s response with the new/reduced stiffness. Othereise the module waits for the next displacement sample.

• (g) If average mean-square has reached a plateau, the mean-square equations from equation (8) are used in a backward way to determine the new base stiffness based on the mean-square value at time step n assuming that no other system parameter has been modified. Provided a plateau is reached, the number of retuning attempts is limited to within at least 50,000 samples (500 s ≈ 8 min) of each other, to preserve computing resources.

• (h) When the new stiffness is determined, the system is re-tuned to establish new optimal pendulum length, as auxiliary damping is unchanged throughout the study.

• (i) The new pendulum length is included in simulations to reflect passive control adaptation applied in this framework.

• (j) Next, the mean-square value is monitored to detect when the mean-square value approaches its prospected asymptotic value as computed through mean square equations. When the last 250 s of data is close enough to that value for a certain precision, the system is deemed to have recovered to a pristine state and the framework is reset to its initial state, until the new threshold is detected.

A sampling frequency of 100 Hz was applied for 5000s simulation. The fidelity of the model is thus linked to the performance of this framework along with the modal analysis accuracy in experiments (Figure 4). With more sensors, damage location identification may be possible while a separate module is required to characterise the type of damage, including discriminating and assessing different nonlinearities present.

Figure 4.

Detection and re-tuning flow chart.

3. Results

Figure 5 represents an example output. In the top horizontal axes, the sample mean-square displacement over time (solid thick black) is represented with the instant of damage introduction (vertical dashed thick magenta line) and damage instant detected by the established mean-square threshold (red dotted vertical line). The vertical thick dashed-dotted blue lines show the attempts at system retuning. The horizontal black dashed-dotted line shows the expected value of the mean-square and the two horizontal blue dashed lines show the mean-square threshold at the 10% limit. The solid thick black line represents the sample mean-square displacement over time when retuning is performed. The solid red line shows the mean-square when no retuning is performed. In the bottom left axes, a contour plot displays the expected value of the mean-square displacement as a function of pendulum length and auxiliary damping with tuning parameters at the end of simulation (t = 5000s). Dimensionless variables (ξ, ζ and J) used for the initial tuning were not used as both would vary with length changes, making interpretation and tracking of tuning over time difficult. Instead, length l, auxiliary damping c and mean square value are shown. In the bottom right corner, a plot of the ν (red) and λ (thick blue) over time is presented. The purple vertical dashed line shows the instant of damage detected by these metrics. After damage introduction, two events occur. First, the mean-square displacement crosses the threshold at ’Damage detected’ line.

Figure 5.

Final state of the interface used during this study of the end of the recursion for 20% stiffness damage.

Second, in the bottom right axes, the DSFs spike at t = 500s, with damage detection marked by the vertical purple line. Subsequently, the first retuning attempt occurs at the first blue vertical dashed line (1300s). Here, the retuned (solid thick black) and not retuned (solid red) lines start deviating, since the pendulum length is changed. Over time and several other retuning attempts, the gap between controlled and uncontrolled systems increases as the mean-square is optimized.

Figure 6 shows a timeline of occurrence of key events within the same run displayed in Figure 6. The damage introduction instant and its detection instant (by the two proposed modalities: mean-square, DSFs) are shown with orange milestones. Subsequent retuning attempts are then matched with time instants and newly computed stiffness, error in its estimate and the corresponding optimal pendulum length are shown in teal color.

Figure 6.

Timeline of detection and retuning for a run with 20 % damage.

Table 2 shows, the damage detection instant for both mean-square threshold and DSFs. The DSFs are sensitive to damage for tested severities and consistent with their use in other works ( Bhowmik et al. (2019a); Tripura et al. (2020a)). The mean-square detection indicates that detection instant is guided by damage. Higher severity leads to a lower delay in detection and is consistent with the fact that higher damage means a steeper increase in mean-square displacement, leading to a quicker crossing of the threshold. Damage severity and the magnitudes of DSFs amplitude show consistent, monotonic relationship. Stiffness estimate errors during each retuning attempt, as well as the new optimal length estimates indicate that the highest errors occur at the first retuning attempt before consistently decreasing with each new retuning attempt.

Table 2.

Table presenting the detection and retuning results at corresponding damage severities.

Damage severity		10 %		20 %		30 %
		1123.00		544.66		542.85
Detection instant (s)	Mean-square DSFs	500.05		500.04		500.00
Increase at damage instant (%)	λ	5.9		7.4		14.0
	ν	12.7		16.6		29.1
		Error (%)	Length (m)	Error (%)	Length (m)	Error (%)	Length (m)
Retuning attempt	1	4.06	1.77	6.13	1.95	15.36	2.05
	2	0.35	1.83	1.71	2.10	0.99	2.33
	3	0.94	1.86	0.18	2.06	4.08	2.46
	4	0.81	1.85	0.15	2.07	4.04	2.46
	5	0.50	1.85	0.48	2.06	3.07	2.43
	6	0.12	1.84	0.66	2.08	1.92	2.40
	7	-	-	0.04	2.06	-	-

Figure 7 shows stiffness updating over time for 20% damage. Figure 7(a) represents the interface at the first time step, upon the start of the recursion. The optimal parameters shown in the contour plot are equivalent to the results from the optimization procedure detailed in equation (12). Figure 7(b) shows the interface after 500s, at the moment of damage introduction where the DSFs detect damage instantaneously. Figure 7(c) at the instant of detection of damage using mean-square threshold, occurring at 544.66s. Figures 7(d)–9(j) show successive stiffness estimation and retuning attempts The optimal parameters at the bottom left axes change for each figure as a new optimization with new estimated stiffness is carried out. From the first retuning step (Figure 7), the mean-square corresponding to the system that is not retuned (red solid line), starts deviating from the retuned system’s mean-square. The redefinition of the mean-square threshold after each retuning is shown through change in y-value of the horizontal blue dashed lines.

The impact of mass ratio on the mean-square value was investigated next (Table 3). The 30.22% mean square value reported before could be improved to a 16.95% of the uncontrolled mean-square. This would although require a significant amount of mass and therefore auxiliary damping will be required, both of which can be unrealistic. Alternatively, a 1% mass ratio, closer to many real-world examples, led to an optimal reduction of 51.94 % of the uncontrolled mean-square, in line with Roffel et al. (2013). For mass ratios considered here, for root mean square responses of earthquakes such values are typically 10% (Gerges and Vickery (2005b)) while Bhowmik et al. (2019b) reports 33% of peak response reduction of transfer function using adaptive TMD. Bridge vehicle interactions with damage and TMD give 10% peak displacement reduction maximum (Pakrashi et al. (2010)), including double beam systems with multiple TMDs (Chawla and Pakrashi (2022)). For tuned liquid column dampers, peak response reduction is around what was obtained in this paper for wind turbine type structures (Buckley et al. (2018)) and 20% for floating platforms (Jaksic et al. (2015)). So, the framework obtains control performance in line with what is typically expected from passive systems, but with the advantage of rapid assessment of damage, as well as handling stochasticity properly.

Table 3.

Table of mass ratios and their corresponding mean-square as a percentage of the uncontrolled mean-square.

Mass ratio	1%	5%	10%	15%	20%
Mean-square value	51.94 %	30.22 %	22.88 %	19.25 %	16.95 %

4. Discussions

Detection via DSFs were at most delayed by tens od milliseconds from damage introduction and were sensitive to small damages (10%). While, mean-square alone is capable of detecting damage, it comes with a greater delay - sometimes >10 min. The delay is inversely related to damage severity and damages <10% might not be flagged for a long time. The use of equation (30) may cause some issues in detecting successive damages. For two damages spaced at 500s, DSFs do not lose accuracy. However, the nature of moving average operation prevents discriminating damage events spaced less than the window length. This can be resolved by shortening the moving average window with a tradeoff between removal of high frequency components against time separation of successive damage events.

The features λ and ν. were similar in damage sensitivity, and computing just one in SOP can work. However, ν provides an extra layer of detection and helps corroborate detection results with eigenvalue feature which is more prone to false positives. As stiffness is reduced significantly, the damping of base may increase, introducing nonlinearities and eventually compromising the performance SOP (Mucchielli et al., 2020). For example, at the end of the update for a 30% damage case, the damping ratio was 10.5 %, slightly over the established performance of SOP.

The observed maximum stiffness estimation error for first retuning (Table 2) is consistent with expectations since plateau detection triggers retuning when the derivative is overall constant over several samples. However, early retuning trigger provides the retuning module a mean sqaure value far from its stationary value, causing higher initial error.

For Figure 5, the first two retunings therefore cause significant decrease in mean-square without targeting the most ideal mean-square value due to high estimated stiffness error. Subsequent retunings fine-tune the pendulum length causing the mean-square to approach more closely and hover about the ideal value.

Higher damage severity corresponded with higher stiffness estimation error for all retuning and a greater increase in mean-square. So the number of retuning attempts to estimate stiffness to sufficient accuracy differ for different damage severities. Also, from Figure 9 in appendix, the absolute rate of change of the mean-square with 30% damage is observed to be much higher than in Figure 5. This implies that variations in the mean-square provided to the retuning module may incorporate more variation at higher severities, causing this increase in error in stiffness estimation.

These observations can be made in reverse via Figure 8, in appendix. Delay in damage detection via mean-square is 623 s (≈10 min), due to low damage severity (κ = 0.90, 10%). Compared to higher damage, retuning attempts perform fine tuning right from the first occurrence. This is because stiffness was better estimated due to the smaller rate of change of the mean-square, also evidenced by low initial error in stiffness estimates in Table 2.

Rapid detection is possible for many scenarios as the average total time of processing was 1000s, equivalent to 0.2s/sample. Removing graphics, this number decreases to 0.05 s/sample and can be further reduced via parallel computation for mean-square/damage detection modules. Less time intensive computing facilities (e.g. C/C++,Python) can also make the framework faster. Note, response covariance matrix can be sensitive to excitation changes Noel et al. (2017); Li et al. (2019). Here, measurement rates are much higher than changes of interest which come from slower time scales. Also, ambient vibrations are typically wide banded and not far away from a white noise. The use of covariance matrix involving white noise excitation Law et al. (2012) works well even with external noise trying to mask damage features. Yan et al. (2003) indicates that establishing a correlation between the measured features and corresponding environmental conditions (Buckley et al. (2021); Cross et al. (2011) can be useful. Another variant, is the covariance of covariance of acceleration Li and Law (2010) which may be investigated in future.

Some limitations are mentioned. Damage is introduced as instantaneous base stiffness reduction, but other complex and nonlinear forms (Ramallo et al., 2004; Tripura et al., 2020) or slow degradation (Tripura et al., 2020b) exist. Excitation noise is Gaussian but other forms should be considered. We assume sampling frequency of synchronous acceleration and displacement inputs equal to the simulation period. More work is needed to determine the influence of sampling frequency on performance. Very low stiffness change (typically <10%) is difficult to detect. Varying random excitation spectrum may result eigenvalue variation based on correlation functions and these should be calibrated separately. Strong non-stationary situations will reduce performance due to stationarity assumption of the analytical derivatives. Finally, apart from displacement there are markers of control linked to acceleration and force and these markers should be investigated individually in future, including combined acceleration conditions, system uncertainties, interstory/inter-structure coupling and damage classes.

5. Conclusions

A hybrid forward-inverse framework incorporating stochastic simulation, damage detection and retuning of a 2-DOF base-isolated PTMD-equipped structure was presented. The performance of the framework was analysed for 10%, 20% and 30% instantaneous stiffness loss damage in base floor. Two damage sensitive features computed from the eigenvalue and eigenvector outputs from second order eigen-perturbation were used as damage detection metrics. Also, the mean square equations of the system were derived and initial PTMD tuning via optimization of first floor mean-square displacement was performed. Monitoring of sample mean-square displacement from the system enabled damage detection through by establishing a statistical threshold. Following damage detection, the mean-square displacement was minimized over time through a change in pendulum length whose optimal value was computed via optimisation of the mean square equation of the first floor. Results show

• near instantaneous damage detection is possible by selecting the proposed damage sensitive features

• a novel conjunction of the forward and inverse approaches provides an advantage for combined control and damage detection

• a stochastic differential equation approach is useful for damage detection and control

• a detection via mean-square displacement threshold is a function of and delayed by damage severity

• retunings of PTMD can be efficiently carried out by the proposed method leading to consistently decreased error in estimating the new reduced stiffness as well as providing the appropriate optimal length of pendulum for adequate mean-square minimization.