Experimental investigations and numerical modelling: a fruitful interaction for the nonlinear dynamical analysis of masonry structures

Abstract

This paper describes the experiments carried out on a mediaeval masonry tower in the historic centre of Lucca and some finite element numerical simulations of the tower’s experimental response. The Guinigi Tower, one of the most iconic monuments in Lucca, has been continuously monitored by high-sensitivity seismic stations that recorded the structure’s response to the dynamic actions of the surrounding environment. The monitoring campaign results have been analysed to show the effectiveness of dynamic monitoring as a valuable source of information on the structural properties of the tower. The dynamic analyses of the tower and the surrounding palace subjected to some seismic events recorded during the experiments have highlighted the capabilities of experiment-based finite element modelling. The calibration of the finite element model and the numerical analysis have been carried out by resorting to procedures developed at ISTI-CNR and able to consider the nonlinear behaviour of masonry materials.

1 Introduction

Towers and bell towers spread worldwide are historical symbols of social value whose safeguarding requires a comprehensive approach encompassing the definition of programmes and initiatives on behalf of national and local authorities, long-term monitoring campaigns and effective numerical modelling and prediction tools.

Long-term vibration monitoring is an efficient, scarcely invasive technique to investigate the dynamic behaviour and check the health status of historic structures. Data recorded by sensor networks are processed using suitable numerical procedures to determine the structure’s dynamic properties, such as frequencies, damping ratios and mode shapes. This approach, known as operational modal analysis (OMA) [3], allows tracking of the dynamic properties and makes damage detection possible. Monitoring the variation of frequencies over time, assessing environmental effects and analysing changes and anomalies are essential ingredients of damage detection. This operation can be accomplished via regression and output-only methods [1, 10, 11, 25, 26] and machine learning techniques [4, 19].

The interaction between ambient vibration tests and numerical modelling constitutes a fruitful approach to structural health monitoring. Using the experimental modal properties of a structure makes it possible to calibrate its finite element (FE) model by adopting model updating procedures [9, 12, 22]. The FE model suitably calibrated can be used to compare the experimental and numerical behaviours of the structure subjected to accelerations recorded during a seismic event or to predict its dynamic response to any time-dependent loads.

Several papers on the structural behaviour of masonry towers are available, addressing different issues. Recent contributions focus on the static and dynamic behaviour of masonry towers, from experimental and numerical points of view. A review of the principal methods adopted for analysing slender structures subjected to seismic actions and an extensive collection of references are provided in Masciotta and Lourenco [16]. The soil–structure interaction is considered in Romero et al. [24] to model the settlements of the Giralda tower in Sevilla (Spain). The results of a 2-month dynamic monitoring survey carried out after a seismic event on a historic masonry tower in France are presented in Azzara et al. [7], where the dependence of modal properties on environmental conditions is also investigated. In Monchetti et al. [18], the authors propose a Bayesian model updating method for confined masonry towers, and in Ponsi et al. [21], a comparison between Bayesian-based and deterministic approaches to the calibration of masonry towers is presented. In Mehrotra et al. [17], the earthquake-induced collapse of masonry towers is investigated via an integrated strategy relying on FE analyses and rocking dynamics implemented in an open-source platform.

This paper is focused on analysing the effects of different vibration sources (anthropic activities, wind and earthquakes) on the dynamic behaviour of an ancient masonry tower in an urban setting. The path from acquiring experimental data to numerical simulation is described in detail. We make use of original numerical procedures for FE modelling and model updating. In particular, the low tensile strength that characterizes masonry materials is taken into account via the constitutive equation of masonry-like materials, implemented in the NOSA-ITACA code (www.nosaitaca.it), a FE software developed at ISTI−CNR and devoted to the structural analysis of ancient masonry constructions. The paper is organized as follows.

In Sect. 2 we describe the experiments carried out on the Guinigi Tower, a mediaeval masonry structure in the historic centre of Lucca, whose ambient vibrations were continuously monitored from June 2021 to October 2022 [2]. Experimental data are then used in Sect. 3 to calibrate FE models of the tower via model updating techniques. Finally, FE linear and nonlinear dynamic analyses are performed by assigning to the tower’s model the signal of the Viareggio earthquake recorded at the base of the structure in February 2022. Experimental and numerical results are then compared.

The main goals of the paper are to show the effectiveness of the ambient vibration tests as a valuable source of information on the structural properties of the tower and highlight the capabilities of experiment-based FE modelling. The experimental results are useful not only for measuring the structure’s dynamic behaviour over time and under earthquakes and exceptional loadings but also for building a realistic FE model via model updating and validating the numerical response via a direct comparison between numerical simulation and recorded signals.

2 Dynamic monitoring of the Guinigi Tower

The results of a long-term dynamic monitoring campaign conducted on the Guinigi Tower in Lucca are presented in this Section. The Guinigi Tower is a mediaeval structure located in the restricted vehicular traffic area of the historic centre (Fig. 1) and managed by the Municipality of Lucca. The tower, one of the most famous buildings in Lucca, is about 44 m high and dates to the XIV century. It is regularly open to the public and accessible via a metallic staircase from the level of 23 m to the roof terrace, at the height of about 43 m. From the ground floor to the level of 23 m, the tower is incorporated into the surrounding Guinigi Palace (Fig. 2) and is entirely free from 23 to 40 m, without inner diaphragms. The terrace is accessed from a massive masonry vault at 40 m, and the thickness of the tower’s walls is about 1 m, constant along the height.

The tower was monitored from June 2021 to October 2022 using seismic stations produced by SARA Electronic Instruments. Four triaxial velocimeters (three SS45s with eigenfrequency 4.5 Hz and one SS20, with eigenfrequency 2 Hz), each coupled with a 24-bit digitizer (SL06), were installed on the tower according to the sensors’ layout shown in Fig. 3. One station was placed on the underground floor (SS20 2045), one at 17.9 m (SS45 2542), and the remaining two stations at 39.88 m (SS45 2896 and SS45 2898, Fig. 4a). The sampling frequency was 100 Hz. Two thermo-hygrometers with a sampling time of 60 s completed the monitoring system. A Virtual Private Network allows sending the data recorded from the instruments to a server hosted at ISTI-CNR for storage and processing (Fig. 4b).

Before installing the permanent monitoring system, three preliminary tests, each lasting about 40 min, were conducted in January 2021, adopting three different sensor layouts to assess the tower’s natural frequencies and mode shapes reported in Table 1.

Velocities recorded on the tower were processed via the Covariance-driven Stochastic Subspace Identification method (SSI/Cov), an OMA technique in the time domain implemented in the MACEC code [23] for calculating the modal parameters of a structure using vibration data measured under operational conditions. The structure, subjected to unknown input, is modelled in the time domain as a discrete linear time-invariant system whose dynamic behaviour is governed by the following state-space model:

$$\begin{aligned}{} & {} \mathrm {\textbf{x}}_{k+1}=\mathrm {\textbf{A}}\mathrm {\textbf{x}}_{k}+\mathrm {\textbf{w}}_{k}, \end{aligned}$$

(1)

$$\begin{aligned}{} & {} \mathrm {\textbf{y}}_{k}=\mathrm {\textbf{C}}\mathrm {\textbf{x}}_{k}+\mathrm {\textbf{v}}_{k}, \end{aligned}$$

(2)

where $\mathrm {\textbf{x}}_{k}\in {\mathbb {R}}^{n}$ is the state of the system at the $k^{th}$time, $\mathrm {\textbf{y}}_{k}\in {\mathbb {R}}^{m}$ is the measured output vector, and $\mathrm {\textbf{w}}_{k}\in {\mathbb {R}}^{n}$ and $\mathrm {\textbf{v}}_{k}\in {\mathbb {R}}^{m}$ are the process and output noises, modelled as white noise random processes. Integers n and m are the system’s order and the number of channels of the velocimeters installed on the structure. $\mathrm {\textbf{C}}\in {\mathbb {R}}^{m\mathrm {\times }n}$ is the output matrix and $\mathrm {\textbf{A}}\in {\mathbb {R}}^{n\mathrm {\times }n}$ is the state matrix whose eigenvalues characterize the systems’ modal properties.

Fig. 1

Panoramic view of Lucca and the Guinigi Tower (courtesy of the Municipality of Lucca)

Fig. 2

The Guinigi Palace and its tower (courtesy of the Municipality of Lucca)

Fig. 3

Positions of the seismic stations placed on the Guinigi Tower for the long-term monitoring

Fig. 4

a Seismic station installed on the tower at 40 m during the long-term monitoring. b The server hosted at ISTI-CNR

A detailed analysis of the velocities recorded on the tower and the experimental frequencies calculated via the SSI/Cov algorithm was presented in Azzara et al. [2], where the dependence of frequencies on temperature and humidity, as well as the effect of the visitors on the tower’s dynamic behaviour, was investigated. Here we limit ourselves to recalling the main results of such analysis. Table 1 reports the first five frequencies and their corresponding damping ratios calculated in the preliminary test in January 2021. The first four frequencies correspond to flexural mode shapes, and the fifth relates to a torsional mode shape.

The main frequencies of the towers exhibit a variation in the order of 5–6% from August 2021 to April 2022. Frequencies and temperature are scarcely correlated. The presence of tourists strongly influences the dynamics of the building and does not allow evaluation of the effects of the environmental parameters on the tower’s modal properties [2].

Table 1 Experimental frequencies and damping ratios calculated in the preliminary test, January 2021

The average velocities measured on the tower, which is located in a restricted traffic area, are very low, on the order of $0.3 \times \,10^{\mathrm {-3}}$ m/s. However, considerable variations in the velocity levels were observed during the study period, with peaks of up to 5 mm/s, corresponding to the windiest days, the hours with high numbers of visitors and the seismic events. The tower’s velocities and the number of tickets sold from 2 to 6 July 2021 (opening local time from 9:00 to 19:00) are plotted in Fig. 5 versus time.

Figure 6 shows the spectrograms of the signals recorded by SS45 2896 from 13 September (Monday, Santa Croce feast) to 19 September (Sunday) 2021 in the x (upper) and y (lower) directions. The tower’s first two main frequencies (corresponding to bending mode shapes) are clearly identifiable in the figure: $f_{1} = 1.23$ Hz in the x direction and $f_{\textrm{2}} = 1.34$ Hz in the y direction. Every day from 9:00 to 19:00, there is an evident increase in the power spectral energy (characterized by a − 200 dB/Hz peak) of the signal recorded on the tower. This increase is due to the people visiting the tower during the opening period. An energy increase is also detectable in correspondence to the Santa Croce feast and the weekends. Interestingly, the effect of the visitors inside the tower is mainly visible on the two first frequencies, thus indicating an increase in the oscillation velocities of the top of the building. In contrast, the crowd moving in the historical centre during the feast days affects the whole frequency spectrum of the structure.

Fig. 5

Daily trend of velocities in the x direction recorded by the sensor SS 2542 ($+$ 17.90 m) and the number of presences inside the tower from 2 to 6 July 2021 (right); zoom of the first day (left)

Between July 2021 and May 2022, about 1300 earthquakes have been located by the Italian National Seismological Network of INGV. In this period, the seismic monitoring stations recorded 17 earthquakes of magnitude between 3.0 and 4.2 in a circular area centred in the Guinigi Tower with a radius of 108 km. In May 2022, a seismic sequence occurred near Florence, about 60 km from Lucca, producing more than 200 earthquakes. Figure 7 shows the recorded events and the map of their epicentres, along with their magnitude and distance from the Guinigi Tower.

Fig. 6

Spectrogram of the signal recorded by SS45 2896 in the x (up) and y (bottom) directions from 13 (Monday, Santa Croce feast) to 19 (Sunday) September 2021 (UTC time)

Fig. 7

Map of the seismic events occurred from July 2021 to May 2022 in a circular area centred on Lucca. For each event, magnitude and distance from the Guinigi Tower are listed

The Viareggio earthquake (event no. 5 in Fig. 7), with magnitude 3.7, located 20 km from Lucca, induced the maximum acceleration on the top of the tower. The maximum acceleration recorded in the x direction at Level 5 (see Fig. 3) is just over 8 mg, while for the earthquake of maximum magnitude 4.2, located just over 100 km from the tower, and for one of the strongest events of the Impruneta sequence (magnitude 3.7 at about 60 km from Lucca) the maximum acceleration does not exceed 2 mg.

Figure 8 reports the x, y and z velocities and accelerations recorded by stations SS20 2045 (green), SS45 2542 (red) and SS45 2896 (blue) during the Viareggio earthquake on 6 February 2022. The maximum ground acceleration and velocity recorded at the tower’s base are $3.6 \times 10 ^{\mathrm {-2}}$ m/s$^{\textrm{2}}$ and $1.8 \times 10^{\mathrm {-3}}$ m/s, while the maximum values recorded at the top are $8.1 \times 10^{\mathrm {-2}}$ m/s$^{\textrm{2}}$ and $6.3 \times 10^{\mathrm {-3}}$ m/s. Strong amplification of the signal along the tower’s height can be observed, particularly in the horizontal directions, along which the velocity at the top of the tower is more than 3.5 times that recorded at the base.

Fig. 8

Seismograms of the Viareggio (LU) earthquake (magnitude 3.7, on 6 February 2022, at about 20 km from the tower), velocities $v_{\textrm{x}}$, $v_{\textrm{y}}$ and $v_{\textrm{z}}$ (m/s) (left) and accelerations $a_{\textrm{x}}$, $a_{\textrm{y}}$ and $a_{\textrm{z}}$ (m/s$^{\textrm{2}}$) (right) versus time t (s) recorded by stations SS20 2045 (green), SS45 2542 (red) and SS45 2896 (blue) (colour figure online)

Figure 9 shows the Standard Spectral Ratio (SSR), i.e. the ratio between the Fast Fourier Transform (FFT) of the signals recorded by station SS45 2896 located on the top of the tower (Level 5, $+$ 40 m) and the station at the base, before (blue), during (black) and after (red) the Viareggio earthquake for the x, y and z directions. The curves are normalized with respect to their maximum value. It is worth noting that the blue and red curves almost coincide, thus indicating a substantially unchanged behaviour of the tower before and after the event in terms of frequencies. The plots reported in Fig. 10 refer to the signals recorded by station SS45 2542 ($+$ 17.90 m). The effect of the earthquake on the tower’s frequencies is here less evident, likely because station 2542 is located in the portion of the tower surrounded by the Guinigi Palace.

Figure 11 shows the spectrograms of the signals recorded by SS45 2896 from 1 to 7 February 2022 in the x (upper) and y (lower) directions. As in Fig. 6, every day there is an evident increase in the signal’s power spectral energy due to the people visiting the tower. In addition, the occurrence of the Viareggio earthquake is visible in the spectrogram, which exhibits an energy peak of about − 100 dB/Hz on 6 February at 1:36. It is worth noting that the energy increase due to the visiting people is comparable to that induced by the Viareggio earthquake.

Fig. 9

Normalized spectral ratio of the Viareggio earthquake (black line) and the signals recorded before (blue line) and after (red line) the earthquake versus the frequency (Hz), station SS45 2896 ($+$ 40 m) (colour figure online)

Fig. 10

Normalized spectral ratio of the Viareggio earthquake (black line) and the signals recorded before (blue line) and after (red line) the earthquake versus the frequency (Hz), station SS45 2542 ($+$ 17.90 m) (colour figure online)

Fig. 11

Spectrogram of the signal recorded by SS45 2896 in the x (up) and y (bottom) directions from 1 to 7 February 2022. The Viareggio earthquake occurred on 6 February at 1:36 (UTC) is highlighted in the black box

3 FE model updating and dynamic analysis

The availability of the experimental frequencies makes it possible to calibrate a FE model of the Guinigi Tower via model updating procedures. The goal is to estimate some unknown parameters (Young’s moduli and mass densities of the materials constituting the structure) and thus develop a reliable model to conduct FE analyses and simulate the tower’s structural response.

The model updating of the numerical model was conducted via the global optimization algorithm [12] implemented in the FE code NOSA-ITACA developed by ISTI-CNR for the analysis and calibration of masonry structures (www.nosaitaca.it/software).

The algorithm minimizes the distance $\upphi$(x)

$$\begin{aligned} {\upphi }\left( \mathrm {\textbf{x}} \right) =\sum \limits _{k=1}^q w_{k}^{2} \left( \textrm{f}_\textrm{exp}^{k}-\textrm{f}_\textrm{num}^{k}\left( \mathrm {\textbf{x}} \right) \right) ^{2}, \mathrm {\textbf{x}}\in {\Omega } \end{aligned}$$

(3)

between the first q structure’s natural frequencies $f^{k}_{\textrm{exp}}$ evaluated experimentally and their numerical counterparts $f^{k}_{\textrm{num}}$(x) ($k = 1$, ..., q) evaluated by modal analysis, depending nonlinearly on the vector x of the p model’s unknown parameters (with $p \le q)$ belonging to a feasible set ${\Omega }$.

The scalar $w_{k}$ is the weight that should be given to frequency $f^{k}_{\textrm{exp}}$ in the optimization scheme. If the goal is to minimize the distance between the vectors of the measured and computed frequencies in the usual Euclidean norm, $w_{k}=1$, should be chosen. If, instead, relative accuracy on the frequencies is desired, $w_{k}={\mathrm {1/f}}_\textrm{exp}^{k}$ is a natural choice.

The algorithm is based on a recursive procedure for constructing local parametric reduced-order models; the minimization problem is solved via a trust-region scheme. Given a feasible set ${\Omega }$ into which the parameters vector x is allowed to range, the procedure provides a set of local minimum points, including the global one. Generally, experimental frequencies may not be accurate, since they are derived by manipulating measured data that may be contaminated by noise. Thus, when minimizing objective function $\upphi$(x), one must ensure that the optimal parameters are well-defined and robust to perturbations in the data $f^{k}_{\textrm{exp}}$. To this purpose, for x* a local minimum point of function $\upphi$(x), we recall the quantities $\upzeta _{i}$ and $\upeta _{i}$ ($i = 1$, ..., p) introduced in Girardi et al. [12], which involve the Jacobian of the numerical frequencies $f^{k}_{\textrm{num}}$(x*) and measure how trustworthy the single parameter $x_{i}$* is. We have $\upeta _{i} \le \upzeta _{i}$ and the following situations can occur.

(i) $\upeta _{i} \le \upzeta _{i}<< 1$: parameter $x_{i}$ cannot be reliably determined, as no information on it is encoded in the optimization problem.

(ii) $0<< \upeta _{i} \le \upzeta _{i}$: parameter $x_{i}$ can be reliably determined from the data, even if it is subject to noise.

(iii) $\upeta _{i}<< 1$, but $\upzeta _{i}>> 0$: there is some information on parameter $x_{i}$ encoded in the problem, but the result will not be free of noise.

Moreover, quantities $\upzeta _{i}^{\mathrm {-1}}$ and $\upeta _{i}^{\mathrm {-1}}$ estimate the minimum and maximum percentage error in assessing the parameter’s optimal value $x_{i}$* under the hypothesis of a 1% error in identifying the experimental frequencies.

The information provided by the quantities $\upzeta _{i}$ and $\upeta _{i}$ helps the user to assess the parameters’ reliability avoiding preliminary sensitivity analyses [20]. Such analyses are usually performed before calibrating the numerical model to choose the number of updating parameters and to exclude some uncertain parameters from the model updating process. The computational cost of sensitivity analysis is very high (of the order of hundreds of modal analyses runs) with respect to the cost of the minimization procedure implemented in NOSA-ITACA, which provides the global minimum point and an assessment of its reliability in few iterations [12].

A FE model of the complex constituted by the tower and the palace was created with NOSA-ITACA. The mesh of the structure, assumed to be perfectly clamped at the base (no information related to the foundations and the soil is available, so the soil–structure interaction is neglected), is shown in Fig. 12 and consists of 23,161 iso-parametric thick shell and beam elements (element no. 10 and no. 9 of the NOSA-ITACA library) with 22,665 nodes, for a total of 135,990 degrees of freedom.

The finite elements of the tower and the palace have a mean dimension of 0.6 m and 1.2 m, respectively. The thickness of the tower is 0.9 m for the northern and southern walls and 1.1 m for the western and eastern walls. The palace’s walls are 1.0 m thick, and the slabs have a mean thickness of 0.35 m.

Fig. 12

FE discretization of the Guinigi Palace and Tower (left), materials considered in the mesh (right). Different colours correspond to different materials (colour figure online)

The model updating of the structure is carried out by considering the first five experimental frequencies $f^{k}_{\textrm{exp}}$ determined in the preliminary tests conducted in January 2021 and reported in Table 1 and setting $w_{k}={\mathrm {1/f}}_\textrm{exp}^{k}$ in (3). Bearing in mind that the number of parameters to be optimized is expected to be no greater than the number of frequencies to match, we start by considering a first FE model (named M$^{\mathrm {5-5}})$ in which the structure is divided into three parts: the palace (in blue in Fig. 12), the portion of the tower contained within the palace (in red in Fig. 12) and the remaining free portion of the tower (in brown in Fig. 12). In the absence of information on the materials constituting the structure, this choice seems natural and reflects the actual characteristics of the complex. The unknown parameters are the Young’s modulus E of the palace, Young’s moduli $E_{\textrm{1}}$ and $E_{\textrm{2}}$ and mass densities $\rho _{\textrm{1}}$ and $\rho _{\textrm{2}}$ of the portions of the tower outside the palace and contained in it, respectively (Fig. 12). The elastic moduli $E_{\textrm{1}}$ and $E_{\textrm{2}}$ are allowed to vary within the interval [1.0, 6.0] GPa, while E ranges in [0.5, 6.0] GPa. The mass densities $\rho _{\textrm{1}}$ and $\rho _{\textrm{2}}$ vary in the interval [1200, 2000] kg/m$^{\textrm{3}}$, and the Poisson’s ratio of all materials is 0.2. Table 2 summarizes the structure’s numerical frequencies ($f^{\textrm{1}}_{\textrm{num}}$, $f^{\textrm{2}}_{\textrm{num}}$, $f^{\textrm{3}}_{\textrm{num}}$, $f^{\textrm{4}}_{\textrm{num}}$, $f^{\textrm{5}}_{\textrm{num}})$ corresponding to the optimal value of the parameter vector x = ($E, E_{\textrm{1}}$, $\rho _{\textrm{1}}$, $E_{\textrm{2}}$, $\rho _{\textrm{2}})$ and their relative error |$\Delta f^{\textrm{k}}$| with respect to the experimental counterparts. The maximum value of the relative error for the first frequency is 0.81%. Table 2 also reports the optimal values of the parameters recovered by the NOSA-ITACA code. It shows, for the optimal value $x_{i}$*, the quantities $\upzeta _{i}$ and $\upeta _{i}$ and the percentage error interval [$\upzeta _{i}^{\mathrm {-1}}$, $\upeta _{i}^{\mathrm {-1}}$]. The table shows that, in the worst-case scenario, at most, the Young’s modulus E of the palace will be affected by a 2.14% error and the Young’s moduli of the tower by 10.60%. The mass density $\rho _{\textrm{1}}$ will be affected by a 2.10% error; the mass density $\rho _{\textrm{2}}$ of the material constituting the portion of the tower contained in the palace exhibits a very high percentage error (175.44%), likely because $\rho _{\textrm{2}}$ has a low influence on the evaluation of the frequencies and cannot be reliably estimated by the model updating process.

Table 2 Model M$^{\mathrm {5-5}}$: experimental and numerical frequencies, optimal values of the parameters calculated by NOSA-ITACA and their percentage error intervals

To highlight the effectiveness of the minimization procedure implemented in NOSA-ITACA, and because of the high degree of uncertainty that affects the calculation of $\rho _{\textrm{2}}$ in model M$^{\mathrm {5-5}}$, we considered a second model of the structure (named M$^{\mathrm {5-4}})$ obtained from M$^{\mathrm {5-5}}$ removing $\rho _{\textrm{2}}$ (set equal to the mass density of the palace, 1800 kg/m$^{\textrm{3}})$ from the unknown parameters. Table 3 summarizes the results of the minimization process for model M$^{\mathrm {5-4}}$.

Table 3 Model M$^{\mathrm {5-4}}$: experimental and numerical frequencies, optimal values of the parameters calculated by NOSA-ITACA and their percentage error intervals

The numerical frequencies calculated using the optimal values of the parameters coincide for both models. For the sake of comparison, Table 4 reports the optimal values of the parameters for models M$^{\mathrm {5-5}}$ and M$^{\mathrm {5-4}}$.

Table 4 Comparison between models M$^{\mathrm {5-5}}$ and M$^{\mathrm {5-4}}$

The table shows that the optimal values of the unknown parameters for models M$^{\mathrm {5-5}}$ and M$^{\mathrm {5-4}}$ are very close, and their relative differences are at most 5%. This comparison confirms the robustness of the proposed model updating procedure to the parameters’ choice and its capability to highlight unreliable results. Figure 13 reports the first five mode shapes of the structure calculated by NOSA-ITACA for the optimal values of the parameters summarized in Table 2. The MAC values calculated between the numerical mode shapes and their experimental counterparts are greater than 0.9.

Fig. 13

Mode shapes calculated by NOSA-ITACA using the optimal values of the parameters

The proposed calibrated model is an acceptable trade-off between the relative errors on the frequencies and the reliability of the parameters. It can be used to conduct dynamic analyses by assigning the signals of seismic events recorded at the tower’s base to the model and comparing experimental and numerical velocities.

The nonlinear analyses are conducted adopting the constitutive equation of masonry-like materials [8], which models masonry as an isotropic nonlinear elastic material with zero tensile strength and infinite compressive strength. Assumptions underlying the model are that the infinitesimal strain tensor E is the sum of an elastic part E$^{e}$ and a fracture part E$^{f}$ and that the stress tensor T, negative semidefinite, depends linearly on the elastic part, T $=$ C[E$^{e}$], with C the isotropic fourth-order elasticity tensor containing the mechanical properties of the material (the Young’s modulus and the Poisson’ ratio). The fracture strain E$^{f}$ is positive semidefinite and orthogonal to T.

Stress T can be also characterized as the projection of C[E] onto the cone of negative semidefinite symmetric tensors with respect to the inner product defined by the inverse of C [15]. It is possible to define the stress function ${\hat{\textbf{T}}}$ from the linear space of symmetric tensors to the cone of negative semidefinite symmetric tensors defined by ${\hat{\textbf{T}}}\left( \mathrm {\textbf{E}} \right) \mathrm {=\textbf{T}}$, for any symmetric tensor E. The function ${\hat{\textbf{T}}}$ is nonlinear, homogeneous of degree 1, monotone, Lipschitz continuous and Fréchet differentiable on an open dense subset of the space of symmetric tensors.

The explicit expression of ${\hat{\textbf{T}}}(\mathrm {\textbf{E}})$ can be found in Lucchesi et al. [15], together with its derivative $D_{\mathrm {\textbf{E}}}{\hat{\textbf{T}}}(\mathrm {\textbf{E}})$ with respect to E.

The constitutive equation of masonry-like materials can realistically describe the mechanical behaviour of masonry, in particular, its inability to withstand tensile stresses. Generalizations of the equation aimed at considering the limited compressive strength of masonry and the presence of thermal dilatations due to temperature variations are described in Lucchesi et al. [15]. The masonry-like model has been implemented in NOSA-ITACA and adopted to simulate the static behaviour of several historical buildings. Limit analysis of arches and vaults has been conducted to investigate their behaviour near the collapse and calculate collapse loads and associated mechanisms.

It is worth noting that the equation of masonry-like materials, inspired by the work of Heyman on masonry arches [13], is not able to predict the collapse of a structure due to shear stresses and the lack of bounds on shear stress may affect the response of masonry structures to horizontal loads; this is not, in general, the case of towers (like the Guinigi Tower) that behave as slender structures. In any case, the shear stresses calculated by the FE analysis can be used to assess the safety conditions according to seismic codes and regulations. Furthermore, considering the masonry’s limited compressive strength imposes indirect bounds on the shear stress.

In NOSA-ITACA the nonlinear dynamic problem is solved by integrating the system of ordinary differential equations obtained by discretizing the structure into finite elements. The presence of damping in the equations of motion is modelled by the damping matrix, which depends on the elastic stiffness and mass matrices according to the Rayleigh assumption [6, 14]. The NOSA-ITACA code adopts the Newmark method for solving the system of differential equations and the Newton–Raphson scheme for solving the nonlinear algebraic system obtained at each time step. The tangent stiffness matrix is calculated using the explicit expression of $D_{\mathrm {\textbf{E}}}{\hat{\textbf{T}}}(\mathrm {\textbf{E}})$.

In the FE analysis, the palace was constituted by a linear elastic material with Young’s modulus E given in Table 2; instead, the tower was made of two masonry-like materials with properties ($E_{\textrm{1}}$, $\rho _{\textrm{1}})$ and ($E_{\textrm{2}}$, $\rho _{\textrm{2}})$ provided in Table 2. The analysis of the complex palace tower was conducted using the accelerogram recorded at the tower’s base (SS20 2045) on 6 February 2022, corresponding to the 1:36 a.m. (UTC) Viareggio seismic event. After the dead loads were assigned, the seismic signal was applied to the model, whose numerical response to the dynamic excitation was compared to that recorded on the tower. This analysis provided a check of the numerical method implemented in NOSA−ITACA. The damping matrix has been calculated using the experimental damping ratios for the first two mode shapes estimated after the tests conducted in January 2021 (Table 1). The duration of the quaking was 60 s, and the time step for numerical integration was set at 0.01 s. The results of a linear dynamic analysis conducted by modelling the masonry as a linear elastic material are also reported for comparison.

Figure 14 shows the maximum eigenvalue of the fracture strain E$^{f}$ in the inner and outer surfaces of the building subjected to its weight, calculated in the preliminary static analysis. The overall structure of the tower does not exhibit evident cracks: this is confirmed by the small values of the maximum eigenvalue of E$^{f}$. The local crack patterns visible in the tower’s wall and in the nearby palace close to the Level 0 (Fig. 3) area not evidenced by the model, whose aim is catching the overall dynamic behaviour of the complex. Applying the seismic load does not significantly modify the fracture strain distribution and, as expected given the features of the seismic signal recorded by the velocimeter 2045 during the Viareggio earthquake, the values of the maximum eigenvalue of E$^{f}$ during the dynamic analysis do not increase significantly.

Fig. 14

Maximum eigenvalue of the fracture strain tensor E$^{f}$ in the intrados (left) and extrados (right) of the building subjected to its weight

Figure 15 shows the response at Level 0 (Fig. 3) in the x and y directions versus time recorded by the instrument SS45 2542 ($+$ 17.9 m, red line), together with that calculated by NOSA−ITACA in the nonlinear (black) and linear (blue) cases. Figure 16 shows the same quantities at Level 5 for the instrument SS45 2896 ($+$ 40 m). Figure 17 reports the FFT of the velocity (m/s/Hz) in the x and y directions versus frequency (Hz) recorded by SS45 2542 (red), calculated by modelling the masonry as a nonlinear elastic material (black) and a linear elastic material (blue). Figure 18, showing the same quantities as Fig. 17, refers to the sensor SS45 2896.

Fig. 15

Viareggio earthquake. Left: velocities $v_{\textrm{x}}$ and $v_{\textrm{y}}$ (m/s) in the x and y directions versus time t (s) recorded by SS45 2542 ($+$ 17.9 m, red) and calculated by NOSA-ITACA modelling the masonry as a nonlinear elastic material (black). Right: velocities $v_{\textrm{x}}$ and $v_{\textrm{y}}$ (m/s) in the x and y directions versus time t (s) recorded by SS45 2542 ($+$ 17.9 m, red) and calculated by NOSA-ITACA modelling the masonry as a linear elastic material (blue) (colour figure online)

Fig. 16

Viareggio earthquake. Left: velocities $v_{\textrm{x}}$ and $v_{\textrm{y}}$ (m/s) in the x and y directions versus time t (s) recorded by SS45 2896 ($+$ 40 m, red) and calculated by NOSA-ITACA modelling the masonry as a nonlinear elastic material (black). Right: velocities $v_{\textrm{x}}$ and $v_{\textrm{y}}$ (m/s) in the x and y directions versus time t (s) recorded by SS45 2896 ($+$ 40 m, red) and calculated by NOSA-ITACA modelling the masonry as a linear elastic material (blue) (colour figure online)

Fig. 17

Viareggio earthquake. Left: FFT of the velocities (m/s/Hz) in the x and y directions versus the frequency (Hz) recorded by SS45 2542 ($+$ 17.9 m, red) and calculated by NOSA-ITACA modelling the masonry as a nonlinear elastic material (black). Right: FFT of the velocities (m/s/Hz) in the x and y directions versus the frequency (Hz) recorded by SS45 2542 (red) and calculated by NOSA-ITACA modelling the masonry as a linear elastic material (blue) (colour figure online)

Fig. 18

Viareggio earthquake. Left: FFT of the velocities (m/s/Hz) in the x and y directions versus frequency (Hz) recorded by SS45 2896 ($+$ 40 m, red) and calculated by NOSA-ITACA modelling the masonry as a nonlinear elastic material (black). Right: FFT of the velocities (m/s/Hz) in the x and y directions versus the frequency (Hz) recorded by SS45 2896 (red) and calculated by NOSA-ITACA modelling the masonry as a linear elastic material (blue) (colour figure online)

For comparison, Table 5 reports the peak component particle velocities (PCPV) (PCPV denotes the maximum absolute value of the experimental or numerical velocities components) calculated in the nonlinear and linear case along with their experimental counterparts. As shown in Figs. 15 and 16 and highlighted by Table 5, the numerical analysis can estimate the peak velocity values reached by the structure at the monitored points. Except for the x component related to the SS 2896 sensor, the nonlinear dynamic analyses can better approximate the structural response than the simulations made with a linear elastic material. Furthermore, both nonlinear and linear analyses underestimate the velocity value in the x direction; this fact is probably due to the modelling of the wooden floors that laterally restrain the tower. In fact, instead of considering the actual beam grid, the floors were modelled by shell elements with equivalent thickness and mechanical properties; such modelling led to overestimating the flexural and axial stiffness and constituted a more rigid constraint for the tower in the x-direction.

Table 5 Comparison between experimental and numerical PCPV

Figures 19 and 20 show the experimental displacements recorded by the stations SS45 2542 and SS45 2896 (red) during the Viareggio earthquake plotted along with the numerical displacements calculated by NOSA-ITACA in the nonlinear case (black). Experimental displacements have been obtained by integrating the velocities after applying a 0.1 Hz high-pass filter.

Fig. 19

Viareggio earthquake. Displacements $u_{\textrm{x}}$ and $u_{\textrm{y}}$ (m) in the x and y directions versus time t (s) recorded by SS45 2542 ($+$ 17.9 m, red) and calculated by NOSA-ITACA in the nonlinear case (black) (colour figure online)

Fig. 20

Viareggio earthquake. Displacements $u_{\textrm{x}}$ and $u_{\textrm{y}}$ (m) in the x and y directions versus time t (s) recorded by SS45 2896 ($+$ 40 m, red) and calculated by NOSA-ITACA in the nonlinear case (black) (colour figure online)

Figure 21 shows the velocities in the y direction recorded by the stations SS20 2045 (underground), SS45 2542 ($+$ 17.9 m) and SS45 2896 ($+$ 40 m) during the Viareggio earthquake.

The velocity spectra demonstrate that the first two modes dominate the tower’s velocity responses, with minor contributions from the other modes. An initial upward travelling impulse is observed at $t = 8.6$ s (Fig. 20). The first crest of the impulse is recorded on the top at $t = 8.84$ s, after 0.24 s, about 1/3 of the building’s period $T_{\textrm{2}}=1/f_{\textrm{2}}=0.746$ s [5]. As shown also in Table 5, the peak velocity reaches 6 mm/s at the top of the tower.

Fig. 21

Numerical (black) and experimental (red) velocities $v_{\textrm{y}}$ (m/s) recorded on the tower in the y direction versus time t (s) during the Viareggio earthquake (left) at the underground floor (bottom), at $+$ 17.9 m (centre), at $+$ 40 m (top) and corresponding FFT (m/s/Hz) versus the frequency (Hz) (right) (colour figure online)

4 Conclusions

This chapter presents the results of long-term vibration monitoring campaigns conducted on the Guinigi Tower, the most iconic tower in the historic centre of Lucca, whose dynamic monitoring was conducted from June 2021 to October 2022. A FE model of the tower is presented and fine-tuned through a model updating procedure implemented in the NOSA-ITACA code. This procedure allows calibrating the unknown parameters with good accuracy and evaluating the reliability of the optimal solution without resorting to expensive sensitivity analyses. Finally, a nonlinear dynamic analysis is conducted by assigning to the tower’s FE model the signals recorded by the instrument installed at the structure’s base during the Viareggio earthquake (magnitude 3.7, on 6 February 2022). Experimental and numerical results are directly compared in the cases of linear elastic and masonry-like solutions.

The paper describes the whole pathway from acquiring and analysing the experimental data to numerical simulation. The lessons learned during this process can be summarized as follows:

• Continuous long-term dynamic monitoring of ancient towers and monuments allows to collect valuable information regarding the measurement of the vibration levels and sources, the trend in time of the dynamic properties of the structure, and the assessment of the structural response to earthquakes and exceptional loadings.

• The model updating procedure implemented in NOSA-ITACA allows to estimate the materials’ mechanical properties and provide information on the solution’s reliability with low computational effort.

• Direct comparison between numerical simulation and experimental response provides interesting results. The dynamic analysis conducted with NOSA-ITACA shows a good agreement between numerical and experimental results and proves that the tower behaves nonlinearly for low-amplitude earthquakes. We point out that the effect of the Viareggio earthquake did not induce any visible damage on the tower, as also shown by comparing the tower’s spectra before and after the seismic event.

The good agreement between numerical and experimental results represents a validation of the algorithms implemented in the NOSA-ITACA code for dynamic analysis. It shows that the constitutive equation of masonry-like materials can catch the tower’s dynamic behaviour for low to moderate seismic events.