Portal de Periódicos da CAPES

Machine learning‐based automatic operational modal analysis: A structural health monitoring application to masonry arch bridges

2022; Wiley; Volume: 29; Issue: 10 Linguagem: Inglês

10.1002/stc.3028

1545-2263

Marco Civera, Vezio Mugnaini, Luca Zanotti Fragonara,

Infrastructure Maintenance and Monitoring

Structural Control and Health MonitoringVolume 29, Issue 10 e3028 RESEARCH ARTICLEOpen Access Machine learning-based automatic operational modal analysis: A structural health monitoring application to masonry arch bridges Marco Civera, Corresponding Author Marco Civera marco.civera@polito.it orcid.org/0000-0003-0414-7440 Department of Structural, Building and Geotechnical Engineering, Politecnico di Torino, Turin, Italy Correspondence Marco Civera, Department of Structural, Building and Geotechnical Engineering, Politecnico di Torino, Turin 10129, Italy. Email: marco.civera@polito.itSearch for more papers by this authorVezio Mugnaini, Vezio Mugnaini Department of Structural, Building and Geotechnical Engineering, Politecnico di Torino, Turin, ItalySearch for more papers by this authorLuca Zanotti Fragonara, Luca Zanotti Fragonara orcid.org/0000-0001-6269-5280 School of Aerospace, Transport and Manufacturing, Cranfield University, Cranfield, UKSearch for more papers by this author Marco Civera, Corresponding Author Marco Civera marco.civera@polito.it orcid.org/0000-0003-0414-7440 Department of Structural, Building and Geotechnical Engineering, Politecnico di Torino, Turin, Italy Correspondence Marco Civera, Department of Structural, Building and Geotechnical Engineering, Politecnico di Torino, Turin 10129, Italy. Email: marco.civera@polito.itSearch for more papers by this authorVezio Mugnaini, Vezio Mugnaini Department of Structural, Building and Geotechnical Engineering, Politecnico di Torino, Turin, ItalySearch for more papers by this authorLuca Zanotti Fragonara, Luca Zanotti Fragonara orcid.org/0000-0001-6269-5280 School of Aerospace, Transport and Manufacturing, Cranfield University, Cranfield, UKSearch for more papers by this author First published: 19 June 2022 https://doi.org/10.1002/stc.3028Citations: 1AboutSectionsPDF ToolsRequest permissionExport citationAdd to favoritesTrack citation ShareShare Give accessShare full text accessShare full-text accessPlease review our Terms and Conditions of Use and check box below to share full-text version of article.I have read and accept the Wiley Online Library Terms and Conditions of UseShareable LinkUse the link below to share a full-text version of this article with your friends and colleagues. Learn more.Copy URL Share a linkShare onFacebookTwitterLinked InRedditWechat Summary Structural health monitoring (SHM) is one of the main research topics in civil, mechanical and aerospace engineering. In this regard, modal parameters and their trends over time can be used as features and indicators of damage occurrence and growth. However, for practical reasons, output-only techniques are particularly suitable for the system identification (SI) of large civil structures and infrastructures, as they do not require a controlled source of input force. In this context, these approaches are typically referred to as operational modal analysis (OMA) techniques. However, the interpretation of the OMA identifications is a labour-intensive task, which could be better automated with artificial intelligence and machine learning (ML) techniques. In particular, clustering and cluster analysis can be used to group unlabelled datasets and interpret them. In this study, a novel multi-stage clustering algorithm for automatic OMA (AOMA) is tested and validated for SHM applications—specifically, for damage detection and severity assessment—to a masonry arch bridge. The experimental case study involves a 1:2 scaled model, progressively damaged to simulate foundation scouring at the central pier. 1 INTRODUCTION Structural health monitoring (SHM) is one of the fastest-growing topics among structural engineers, especially when applied in a machine learning (ML) framework.1 This is also thanks to public opinion, which is becoming more and more sensitive to the importance of the maintenance of strategic infrastructures, the preservation of cultural heritage buildings, and the safeguarding of human life. For instance, people interact with ageing infrastructures more than ever as modernisation cannot keep up with the increasing demand for transport flows. The constant monitoring of structures and infrastructure can give some important information about the possible damages which may arise during their service life. Moreover, an assessment of the residual life can be estimated, and maintenance strategies may be planned in time to avoid risks to human life and structural failure which may generate damage to other buildings. To extract useful features that characterise the behaviour of a structure or the changing of these with time, different methods of investigation can be pursued for system identification (SI).2 The choice of the most appropriate method and the planning of study strategies depend on the typology of the object of study, the budget and many other factors.1 A common method is the use of vibration-based SHM, aimed at extracting useful features from the dynamic response of a structure,3 especially in the form of modal parameters (i.e., natural frequencies, modal damping ratios and mode shapes). This is the basic concept of modal analysis, which is arguably the most common approach for SI and damage detection. When the structure can be excited with a known input, experimental modal analysis (EMA) can be applied to this aim. Several EMA approaches exist; a recent method was applied as well for masonry bridges in recent studies.4, 5 However, for massive structures and infrastructures and/or for continuous monitoring over a prolonged time, output-only methods should be preferred. When ambient vibrations (AVs) are utilised for output-only MA, the term operational modal analysis (OMA) is used.6 Some examples of OMA applications for large civil infrastructures can be found in the works of several authors.7-11 The main problem for large and complex systems is that the number of relevant modes is not known a priori. Therefore, it is not possible to guess the model order for SI. The common practice is thus to re-run the SI procedure repeatedly for an arbitrary range of model orders. The result is called a stabilisation diagram, where both physical and numerical poles (due to the overfitting of the higher order models over noisy data) are included. Thus, the actual vibrational modes need to be selected among all the identified poles, discarding the rest. This is generally performed by an expert user. Consequently, the main limitation to the extensive use of OMA for continuous monitoring is the requirement of constant supervision by a skilled professional. This is incompatible with the need for uninterrupted structural surveillance. Therefore, the automation of the whole modal extraction procedure is a necessary step for the future generation of civil buildings; this is generally known as automatic (or automated) OMA, i.e., AOMA for short. In this regard, two extensive reviews of OMA and AOMA approaches can be found in Peeters and de Roeck and Reynders.12, 13 A recent example of application in the ambit of continuous SHM for civil engineering is the AOMA system deployed at the Sanctuary of Regina Montis Regalis in Vicoforte, Italy.14 Similar applications were also deployed for the famous case studies of the Z24 bridge in Switzerland15 and the Tamar suspension bridge in Southwest England.16 Some interesting aspects can be found in Rainieri and Fabbrocino.17 Very recently, Cheema et al.18 applied a Dirichlet process Gaussian mixture model (DP-GMM) clustering approach to discern true physical modes from the mathematically spurious modes on a cable-stayed bridge in New South Wales, Australia. He et al.19 applied a modified version of fuzzy C-means (FCM) clustering to the Fourth Nanjing Yangtze River Bridge. Zeng and Hoon Kim20 tested a self-adaptive clustering approach with a weighted multi-term distance on the Z24 and Downing Hall benchmarks. The classic hierarchical clustering-based AOMA was applied by Anastasopoulos et al.21 to almost one year of acquisitions from a steel single-span tied arch railway bridge. Indeed, the Structural Mechanics Section of KU Leuven developed a MatLab® Toolbox,22 which, similarly to the algorithm applied here for this study, includes (among some other options) an SSI algorithm.23 The method was tested on footbridges24 and other steel, concrete or even composite infrastructures (e.g., Reynders et al. and Liu et al.25, 26). Among all the large civil constructions, bridges are arguably the ones of major interest, for several reasons. The main factor derives from their strategic importance as the commercial arteries of any country. They also pose a potentially serious threat to human lives due to the intense traffic which they sustain daily. Some recent events, for instance, are well-descriptive of the risks of bridge structural failures. The collapse of the Morandi Bridge in Genoa on 14 August 2018 caused tens of victims, hundreds of displaced inhabitants and millions of euros of direct economic losses, even without considering the long-term socio-economic consequences of the disruptions. On a smaller scale, these were the same aftermaths faced after the 19th April 2017 viaduct collapse in Fossano (Piedmont), the 9th March 2017 viaduct collapse in Camerano (Marche), the 23rd January 2017 event of Fiumara Allaro (Calabria), the 28th October 2016 event of Annone and the very recent collapse of the bridge of Albiano Magra (Tuscany) happened the 8th April 2020, just to cite some of the most recent bridge and viaduct failures in Italy. This is even more accentuated for historical bridges, as they do not only have an economic and practical value but also an architectural and cultural one. More specifically, masonry structures are notoriously difficult to handle, model27 and investigate.28 Their dynamic response is also very sensitive to damage.29 All these aspects are mainly due to the uncertainties and complex interactions which are characteristic of their building materials.30 The natural ageing of these materials (with most of the cultural heritage being centuries-old) can only aggravate the situation. In this regard, the non-destructive testing (NDT) and the continuous monitoring of existing masonry arch bridges are a complex yet necessary requirement. An investigation about a recent and interesting case of scour-induced collapse for a multi-span masonry arch bridge can be found in Scozzes et al.31 Two examples concerning railway masonry bridges are reported in Costa et al.32, 33 All the previous studies reported in the scientific literature point out the effectiveness and efficiency of OMA and AOMA for the SHM of bridges of any size and use,34-36 ranging from long-span cable-stayed bridges37 to iron arch ones38 and pedestrian walkways.39 However, to the best of the authors' knowledge, there is still a lack of AOMA applications to masonry arch bridges in the scientific literature to this day. The procedure described here is characterised by a novel approach for discerning the physical modes from the computational ones. This overcomes the current issue of spurious modes, which can affect the results. In the scenario investigated in Politecnico di Torino's laboratories, the scouring of the central pile was emulated by means of imposed foundation settlements and rotations; further detail will be provided in the next sections. The rest of this paper is organised as follows. Section 2 briefly discusses the AOMA procedure. Section 3 describes the case study, i.e., the 1:2 scaled model of the masonry two-span arch bridge. Section 4 reports the results, and the paper ends with the Conclusions. 2 THE AOMA PROCEDURE The algorithm applied and tested here for SHM purposes consists of the preliminary application of the stochastic subspace identification (SSI) algorithm, followed by four consecutive steps: application of SSI for an arbitrary range of model orders; a preliminary partition by hard validation criteria (HVC); a further selection by soft validation criteria (SVC); the cluster-based modal identification; and the modal parameters identification. After the initial identification of the poles, the four-stepped concept follows closely the classic procedure defined by Reynders and colleagues40; this is intended for the cleansing of the stabilisation diagram. The main aspects are here recalled for completeness. The HVC are physically based and enforce commonsense rules. The SVC, on the other hand, are purely data driven and do not rely on any further assumptions. In this aspect, clustering and cluster analysis, which are important ML techniques, can be used. The concept is to binary classify the remaining poles as stable or unstable and discard the latter ones, while the elements of the first group are passed to the next phases. The same ML approaches can be applied for the following steps as well. In the clustering phase, the poles that passed both HVC and SVC selections are grouped together according to their distances. These clusters are then further sifted: The ones classified as ‘most probably physical’ are considered as representing a distinct vibrational mode, while the others (assumed as certainly numerical) are discarded. Finally, in the last step, the modal parameters which best characterise any remaining cluster (i.e., its natural frequency, damping ratio, and mode shape) are statistically defined, based on the poles included in it. The four steps can be further detailed as reported in the following subsections. 2.1 Hard validation criteria Three rules are implemented for the HVC. These reflect the suggestions made by Reynders et al.13, 40 The poles with negative or excessively high damping ratios are immediately discarded as unrealistic. The higher boundary is set to 20% according to Cabboi et al. and Demarie and Sabia.38, 41 Real-valued poles are rejected as well since a zeroed imaginary part would mean a non-vibrational free decay, which is certainly spurious for an oscillatory response. 2.2 Soft validation criteria Five parameters are considered for the definition of the SVC: the eigenvalues λ , the natural frequencies f , the damping ratios ξ , the complementary of the modal assurance criterion ν = 1 − MAC , and the mean phase deviation MPD . These will be described in more detail later on. For each one of these parameters, the absolute differences between the poles for different model orders are computed. These are then utilised to define if the pole is stable or not. The rationale for preferring absolute differences over the more traditional approach (relative differences) was investigated and described in previous works.42 In this discussion, the model order is defined accordingly to the number of corresponding poles. Since poles are complex conjugated, the model order ranges between the (arbitrarily defined) n min and n max with increments of two ( n min , n min + 2 , n min + 4 , …, n max ). l n n = n − n min indicates the number of poles belonging to the nth model order. Therefore, the comparison parameters are defined for any i = 1 , 2 , … , l n pole belonging to the model order n when coupled to any m = 1 , 2 , … , l n + 2 pole belonging to the next model order n + 2. For a generic parameter p , that means ∆ p n i : ∆ p n i , m = p n i − p n + 2 m , (1)where ∆ p n i denotes the array of absolute differences for the ith pole of the nth model order, while ∆ p n i , m is its mth term. Equation (1) can be rewritten to accommodate p = λ , f, ξ , 1 − MAC , or MPD , considering that ∆ λ n i ∈ ℂ 1 × l n + 2 for the complex-valued eigenvalues, while the arrays of the other parameters are subsets of ℝ 1 × l n + 2 . The modal assurance criterion (MAC) is computed between the mode shapes ϕ i n and ϕ m n + 2 according to its classic definition43 MAC ϕ i n ϕ m n + 2 = ϕ i n * ϕ m n + 2 2 ϕ i n * ϕ i n · ϕ m n + 2 * ϕ m n + 2 , (2)while MPD is defined as in Philips and Allemang.44 In conclusion, these comparison parameters form a five-dimensional (5D) space, where the poles can be clustered. To do so, once these five arrays have been defined, the distance between poles of subsequent model orders is evaluated as d i n = Δf * + ν * , d i n ∈ ℝ 1 × l n + 2 , (3)with the first term computed as in Equation (1) (for ∆ p = ∆ f ) and then min-max normalised, that is to say, Δf * = ∆ f n i − min ∆ f n i max ∆ f n i − min ∆ f n i . (4) Such that ∆ f * is defined between 0 and 1 and thus comparable with ν * ≡ ν i n (in both cases, zero corresponds to no correlation, and maximum similarity tends to 1). Please note that Equation (3) is similar yet different to the classic procedure described in Reynders et al.,13 where the (min-max normalised) eigenvalue distance in the complex plane ∆ λ is utilised in lieu of the (min-max normalised) difference of eigenfrequencies ∆ f . The rationale for this design choice implemented here was to not include the dependence on the damping estimates, which are well-known to be less reliable and subject to greater variability than frequency and mode shape estimates.45 At this point, for the ith pole of the nth model order, the ‘neighbouring pole’ is defined as the single pole belonging to model order n + 2 corresponding to the minimum value of d i n , i.e. ∆ p i n = ∆ p i n arg min d i n . (5) Therefore, five scalars are extracted from the five arrays ∆ p n i . This process is reiterated for all the poles belonging to a model order and then to all model orders. The resulting values are then concatenated (for any of the five parameters) into a unique vector ∆ p , defined as ∆ p = ∆ p n min ∆ p n min + 2 ∆ p n min + 4 … ∆ p n max − 2 T , (6)where for any model order ∆ p n = ∆ p 1 n ∆ p 2 n … ∆ p l n n n T , (7)with ∆ p 1 n , ∆ p 2 n , … ∆ p l n n n defined as for Equation (5). Remembering that ∆ p is intended here as a shorthand for Δ λ , Δ f , Δ ξ , ν, or Δ MPD , the stable and unstable poles will be defined through k-means clustering (with k = 2 ) in the 5D space defined by the comparison parameters, according to their distance. However, before doing so, another intermediate step is required. The clustering algorithm assumes a Gaussian distribution of these parameters.46 Since this may be not satisfied by the raw parameters (as highlighted by Neu et al.47), a Box–Cox transformation is applied to each one of the five arrays Δ p = Δ λ , Δ f , Δ ξ , ν , Δ MPD , as h γ = Δp − 1 / γ ln Δp , , γ ≠ 0 γ = 0 (8)where the transformation parameter γ is defined through profile log-likelihood maximisation as described in the original paper of Box and Cox.48 At this point, the z scores of the obtained h arrays are computed, separately for the five cases, as z i n = h i n − μ p σ p , (9)where μ p is the global mean and σ p the global standard deviation, respectively computed as μ p = 1 L ∑ n ∈ N ∑ i = 1 l n h i n , (10)and σ p = 1 L ∑ n ∈ N ∑ i = 1 l n h i n − μ p 2 , (11)where N = n min n min + 2 n min + 4 … n max is the array of the considered model orders. Note that the length of the resulting z array is equal to L = ∑ n ∈ N l n n . At this point, is finally possible to define the two clusters—stable and unstable poles—in the 5D space defined by the z scores of the Box–Cox transforms of the five Δ p arrays. 2.3 Cluster-based modal identification In this third step, the poles are firstly clustered according to their distances and by means of hierarchical clustering. Then, the obtained clusters are in turn divided between physical and non-physical (computational) modes. Finally, for any remaining cluster, the poles flagged as outliers are removed. 2.3.1 Hierarchical clustering At this point, the vector d i n described in Equation (3) is extended to include all model orders, resulting in the L-by-L distance matrix D, the generic term of which can be defined (for the uth row and vth column) as D u , v = ∆ f u , v − min ∆ f u , v max ∆ f u , v − min ∆ f u , v + ν u , v , (12) considering ∆ f u , v = f u − f v , ν u , v = 1 − MAC ϕ u , ϕ v . (13) The symmetry of the square matrix D is exploited by considering only its lower triangular part. Moreover, the distances among the poles belonging to the same model order are not considered. This reflects the fact that, by definition, two poles from the same identification cannot represent the same mode. This is intended to manage pole splitting at high model orders, as evidenced in Mugnaini et al.,42 Section 2.6. Thus, D u , v = NaN was forcefully imposed in these cases. The terms in D are then used to hierarchically cluster the poles according to their neighbourhood. The agglomerative (bottom-up) process is stopped when the inter-cluster distance, defined with the average linkage method49 as ℓ q , t = 1 n q n t ∑ i = 1 n q ∑ j = 1 n t D u = x qi v = x tj (14)reaches a threshold value d ~ , defined according to the algorithm detailed in Neu et al.47 It is worth mentioning that Reynders et al.40 proposed a stricter threshold definition; however, this alternative was tested here as well and found to be almost inconsequential. In Equation (14), x qi and x tj indicate the ith and jth elements of the two clusters q and t, where n q and n t are the number of elements included in these two groups. The algorithm of the hierarchical clustering described in Algorithm 6.4.1 from Aggarwal50 was only slightly modified to account for the presence of the NaN terms. The complete implementation is reported in detail in Mugnaini et al.42 (Algorithm 1) and is here omitted for brevity. 2.3.2 Selection of physical clusters This step represents the main novelty (even if not the only one) with respect to the similar, already-existing AOMA procedures such as the ones of Reynders et al.40 and Neu et al.47 As in these two well-known procedures, the clusters are labelled as ‘possibly physical’ (PP) or ‘certainly mathematical’ (CM) according to the number of poles included. The rationale is that clusters with very few poles are almost surely computational. However, in all classic implementations, this is performed (again) through 2-means clustering.40 This approach is understandable since it is intended to distinguish between only two categories—CM and PP. However, the use of the adjective possibly physical in Reynders et al.40 (and probably physical in Neu et al.47) implicitly states that the use of 2-means clustering is largely conservative, in the sense that it tends to preserve non-physical modes rather than risking to discard actual ones. On the one hand, this unbalance makes practical sense since the risk of missing a weakly excited mode is far more consequential than the outcomes of retaining a few spurious identifications. However, on the other hand, this was found to be excessively unbalanced, often leaving too many mathematical modes to be discarded in the next subsequent step or by user manual selection. This was found as well in previous numerical and experimental applications,42 which motivated this research for a stricter sifting strategy. To bypass this problem, a tuneable method is here proposed. This can be seen as a variation of the previous technique. The procedure is as follows: To improve the stability of the process, a number w of empty clusters is artificially added (with w corresponding to 20% of the number of poles included in the largest cluster, as suggested by Reynders et al.40) The k-means clustering is performed, with k = 2, to obtain the two CM and PP centroids. For any potential mode, the distances from the CM and PP centroids are calculated as D I CM and D I PP . The total index is computed as D I TOT = D I CM max D I CM + 1 − D I PP max D I PP / 2 . For any candidate mode, the percentage D I % is obtained as D I % = D I TOT max D I TOT · 100 . At this point, D I % can be seen as a probability of physical meaningfulness. The upper limit (100%) is set by the point(s) closer to the PP centroid and farther from the CM one and thus certainly physical. The user can then manually set the threshold of probability as desired. It derives from the classic implementation of the k-means clustering that for k = 2, the threshold is implicitly set to 50%. The proposed variant allows selecting a different value to require a stricter selection of physical modes, discarding the other (less probably physical) ones. Regarding Step (1), the addition of the empty clusters did not produce any noticeable improvement, at least for the case study here inspected. The resulting modes did not change significantly even by adding empty clusters as suggested by Reynders et al.40 However, this can be a specific aspect of the case analysed, and no detrimental effects were caused by it. Therefore, the 20% value suggested in the literature was kept in use. After these, only the physically meaningful clusters are passed to the next step. 2.3.3 Outlier removal An outlier detection procedure is run on the remaining clusters. Specifically, the mean absolute deviation (MAD)51 is calculated among the poles belonging to the same group. The poles that have at least one among ∆ f , ∆ ξ , or ν higher than 3 scaled MAD from the cluster median are removed. The remaining poles are used to characterise the modal parameters of their corresponding cluster in the next step. 2.4 Identification of the modal parameters Four statistical measures were investigated to better define the cluster-wide parameters. These four options were found to be equivalent with no major differences in the results. These methods can be briefly described as follows: Method 1: Assigning to each cluster the mean of frequency, damping ratio, and (normalised) mode shape ( f ¯ , ξ ¯ , and ϕ ¯ , in this order). Method 2: Assigning to each cluster the modal parameters of the pole which minimise the distance between the mean of the damping ratios ξ ¯ and its own ξ . Method 3: Assigning to each cluster the modal parameters of the pole which minimise the sum of the three distances f ¯ − f , ξ ¯ − ξ , and ϕ ¯ − ϕ . Method 4: Assigning to each cluster the modal parameters of the pole which minimise ϕ ¯ − ϕ . Having defined the whole AOMA procedure, the shifts in frequency and the other modal parameters can be directly linked to the developing damage, as this induced a stiffness reduction in the monitored system. This is here validated on the 1:2 scale model of the masonry arch bridge. 3 THE MASONRY TWO-SPAN ARCH BRIDGE The test structure is not an exact reproduction of an existing bridge but was purposely designed taking into account common features, geometric proportion, and historical design codes. The intent was to recreate the main characteristics of Italian masonry arch bridges. 3.1 The experimental setup The model (portrayed in Figure 1a) was built with clay bricks sized according to the adapted scale law (130 × 65 × 30 mm). The geometry of the two-span arch bridge is described in Table 1; its mechanical properties (Young's modulus E, Poisson ratio ν, and the S and P wave velocities Vs and Vp) are reported in Table 2. Construction materials (bricks and mortar) with poor mechanical resistances were purposely chosen to be similar to the aged materials of historical buildings.52 The arch barrels were stabilised with a backfill made up of a mix of gravel, sand, and debris, topped by a levelled 10 cm-thick concrete layer. The whole structure rested on two reinforced concrete slabs, fixed to the floor. Further details can be found in Ruocci et al.52, 53 FIGURE 1Open in figure viewerPowerPoint (a) Photo portrait of the 1:2 scaled bridge model. (b) The settlement application system, front view (adapted from Ruocci et al.53, 54) TABLE 1. Geometry of the bridge model Parameter Value Measurement unit Total width 1.55 m Total height 1.73 m Total length 5.87 m Arch span 2.00 m Arch radius 2.00 m Arch angular opening 30 ° Note: Values were retrieved from Ruocci.53 TABLE 2. Construction materials Material E (MPa) ν (−) ρ (kg/m3) Vp (m/s) Vs (m/s) Masonry 1500 0.20 1900 937 574 Backfill material 50 0.10 2000 160 107 Concrete 5000 0.15 2200 1549 994 Reinforced concrete 30,000 0.15 2400 3633 2331 Note: Values were retrieved from Ruocci.53 A settlement application system was inserted under the central pile to investigate the effect of the displacements and rotations imposed at the masonry pier. These were intended to emulate riverbank erosion and calibrated accordingly to ad hoc hydraulic flume tests.53 A thin layer of polystyrene, set in between the top steel plate of the settlement application system and the bottom of the mid-length pier, simulated the riverbed sediment. This was progressively taken off to simulate the progressive erosion. Special attention was dedicate