Analysis of Seismic Site Effects in Plio-Quaternary Intermontane Basin (L'Aquila, Central Italy)

1. Introduction

In Italy, Level 3 analysis involves a detailed representation of the local seismic response using a strictly quantitative technique to assess local seismic impacts. This technique relies on laboratory and site investigations, as well as the processing of analytical and/or numerical data. It is now well-known that, at a local site scale, the seismostratigraphic and topographic effect can induce the seismic wave amplification [Borcherdt 1970; Kramer 1996]. These phenomena are known as local seismic effects. It is easy to see how these effects have a significant impact on infrastructures and buildings [Bisch et al. 2012¸ CS.LL.PP. 2018]. Therefore, they need to be assessed for the best practice of seismic risk mitigation [SM Working Group 2008; Ansal et al. 2010; Aversa and Crespellani 2016]. Their evaluation is typically obtained through numerical simulations of seismic wave propagation and the estimation of the amplification factors (AFs) [Pergalani et al. 2020]. For a given period interval (e.g., 0.1 – 0.5 s), AF is the ratio between the average spectral acceleration (averaged over the period interval itself) measured at the ground surface and that measured on outcropping seismic bedrock, given the same seismic input. Modeling requires, as input, the geometry of the seismo-layers and their geotechnical and geophysical properties [Ciancimino et al. 2019]. These primarily include the depth of the seismic bedrock, the dip of the seismostratigraphic boundaries, Vs, and the dynamic behavior of each seismo-layer.

Investigating these local effects becomes especially challenging when the subsoil model and the dynamic properties of the seismo-layers are complex and deviate from typical patterns, such as not in the simplest case of flat, parallel layers where Vs increases with depth [Bard and Bouchon 1985; Lanzo and Silvestri 1999; Smerzini et al. 2011; Giallini et al. 2024].

The intermontane plains of central Italy, subject to a high seismicity, as evidenced by the 2009 L'Aquila earthquake (Mw 6.3) and 2016-17 central Italy seismic sequence (Mw 6.0 and 6.5), exhibit particularly complex geometries and geophysical and geotechnical characteristics [MS-AQ Working Group 2010; Evangelista et al. 2011; Lanzo et al. 2011; Monaco and Amoroso 2016¸ Alleanza et al. 2024]. These can be schematized as follows: (i) thin shallow seismo-layers with strong seismic impedance contrasts; (ii) sub-vertical or strongly inclined geological boundaries between different lithologies such as faults or erosion surfaces; (iii) inversions of seismo-layer velocities with depth; (iv) at least two or more seismic impedance contrasts at different depths; (v) buried paleomorphology of deep and shallow paleo-valleys; (vi) positive topographies such as ridges and hills; (vii) channelling of seismic waves within fault rocks; (viii) regolith of outcropping rock with Vs < 800 m/s; (ix) infilling basin deposits that instead show Vs > 800 m/s.

Therefore, for proper prevention and management of urban areas of central Italy, the study of local seismic effects plays a significant role, to the extent that numerous pilot projects on seismic microzonation and the study of local seismic effects have been carried out.

One of these was carried out in the areas of L'Aquila Municipality, following the near-source event of April 6th, 2009 [Chiarabba et al. 2009]. A portion of the results of this seismic microzonation project is presented in this work. It is focused on two western (Preturo-Sassa) and eastern (Bazzano-Monticchio) pilot areas of L'Aquila Municipality (Figure 1 and Figure 2).

This seismic microzonation project, aimed at best practice on the management of urban and land territory for seismic risk mitigation, was funded by the Abruzzo Region's Department of Government of the Territory and Environmental Policies - Civil Protection Risk Prevention Service.

The pilot areas of about 7 km² are geologically situated within the L'Aquila intermontane basin, characterized by thick accumulations of coarse- and fine-grained sediments deposited during the Plio-Quaternary primarily in lacustrine, slope, and fluvial settings [Nocentini et al. 2017].

In compliance with the Italian Civil Protection Department's guidelines [SM Working Group 2008], the seismic microzonation project concluded with the production of numerical maps in which the AFs are reported. The AFs were calculated via 1D and 2D numerical simulations of local seismic response on fifteen 3500-5000 m long geological sections representative of the pilot areas. They were specifically elaborated within the seismic microzonation project via many borehole data. Furthermore, the geophysical and geotechnical properties of the seismo-layers represented in them were assigned based on many in-situ and lab investigations. Three-dimensional local seismic response simulations were not performed due to the complexity and inherent uncertainties associated with reconstructing the 3D subsurface structure. Additionally, the substantial computational demands required for modeling such large areas made this approach infeasible. The 2D modelling and the comparison between AFs and the related spectra obtained with 1D and 2D modelling permitted also to recognize several of the above-mentioned local seismic effects. In some cases, these could be reasonably predicted (e.g., the presence of thin shallow seismo-layers with strong seismic impedance contrasts), in other cases not (e.g., buried paleomorphology of deep paleo-valleys). Among the aims of the seismic microzonation project, we find proposals for warnings for the best practices of numerical simulations to estimate AFs used in seismic microzonation maps.

The results presented herein offer twofold utility, serving both speculative and practical purposes. On the one hand, the studied area provides an excellent natural laboratory exhibiting a complex geological and seismostratigraphic setting. Extensive field surveys and geophysical investigations, coupled with numerical simulations, have enabled a detailed reconstruction of the subsurface characteristics and local seismic response, so providing interesting insights into the factors controlling site-specific seismic amplification. On the other hand, this study provides a detailed reconstruction of ground stratigraphic and geotechnical features, as well as of the local seismic response, over several sections in the area, serving as valuable tools for seismic risk mitigation and future studies in this area.

2. Geological Background

The Plio-Quaternary extensional tectonic regime in central Italy, characterized by SW-dipping normal faults, has shaped the landscape, forming graben structures and intermontane basins such as the L'Aquila Basin [Faure Walker et al. 2021]. This latter has been filled with various detrital deposits since the Upper Pliocene (3.6 million years) to the present [Cavinato and De Celles 1999; Nocentini et al. 2017]. The ongoing activity of these faults poses a significant seismic hazard to the region, with the potential for moderate to strong earthquakes (awaited maximum magnitudes of up to 6.5-7) [Devoti et al. 2010; Pondrelli et al. 2006; Galadini and Galli 2000].

The L'Aquila Basin sedimentary record begins in the upper Piacenzian (3.0–2.6 Myr) with the deposition of breccias and conglomerates: Colle Cantaro–Cave Formation in the Preturo–Sassa area and Valle Orsa Formation and Valle Valiano Formation in Bazzano–Monticchio area (all1). These units were followed, via an unconformity boundary above all in the Preturo–Sassa area, by Calabrian (1.8–0.8 Myr) alluvial floodplain silts and sands (Madonna della Strada Synthem, all2). A significant hiatus, marked by an unconformity, separates all2 from the overlying Middle Pleistocene (0.8–0.12 Myr) alluvial fan and plain deposits of the Fosso Genzano Synthem composed by gravels and sands (at1).

The L'Aquila historic center is situated on a hill composed of late Middle Pleistocene (0.3–0.12 Myr) calcareous breccias (Colle Macchione-L’Aquila Synthem, dbf). These breccias overlie the older units all2 and at1 and the Meso–Cenozoic substratum [Sciortino et al. 2024a; 2024b]. The surrounding L'Aquila Basin landscape is lastly characterized by recent slope (fal) and colluvial (col) deposits, as well as the alluvial deposits of the Aterno River and its tributaries such as the Raio Stream (all3).

In the cross-sections reported in Figure 3 the geological relationships of the main units are shown. The Plio-Quaternary L’Aquila basin-cover units are ant, fra, col, at1, all3 S, all3 G, all2, and all1. The units UAP and UAM refer to the terrigenous bedrock; the unit CCD to the carbonate bedrock.

The above-mentioned units of the Preturo-Sassa and Bazzano-Monticchio areas are summarized in Table 1 [Tallini et al., 2024]. The unit codes adopted in the modeling sections are based on the Abruzzo Seismic Microzonation Working Group (2012).

3. Materials and Methods

3.1. Geotechnical and Geophysical Database

A database was set up for storing and managing borehole, geotechnical and geophysical data, and seismic microzonation numerical maps via (i) QGIS (version 3.30.3-‘s-Hertogenbosch) and (ii) the plugin “mzs-tools” by CNR-IGAG (https://github.com/CNR-IGAG/mzs-tools, accessed on 15 November 2024) [Cosentino et al. 2022]. The stored investigations were carried out for the anti-seismic urban reconstruction activities following the April 6, L’Aquila earthquake. Furthermore, ad hoc investigations (above all 2D microtremor arrays, boreholes, static cone penetration tests, ERTs) were performed within the activities of the seismic microzonation project. They were placed in significative sites of the pilot areas for investigating the subsoil model and the coseismic instabilities [Tallini et al. 2024].

Several borehole logs, geotechnical and geophysical surveys, were stored in the QGIS database (Table 2). As an example, in Figure 4 and Figure 5 the investigations and the resonance frequency map obtained through the QGIS database are shown. The maps highilights the high density of data used within the seismic microzonation project. In the Sassa area, crossed by section A-A', numerous microtremor measurements have highlighted resonance frequency values between 2 and 4 Hz.

3.2. Geological Sections

Seven and eight 3500-5000 m-long geological sections were elaborated at fine scale (1:5000) in the Preturo-Sassa and Bazzano-Monticchio pilot areas, respectively and 1D and 2D simulations were conducted on them. The sections were designed to intersect the significant geological features within the pilot areas, such as fault lines, terrace scarps, landslides, and anthropic deposits and most microzones. The microzone is defined as a geologically homogeneous zone susceptible to local amplification or/and to coseimic instabilities at fine - 1:5000 scale [SM Working Group 2008]). Additionally, the sections were placed near the geophysical investigation sites, primarily microtremor measurements and down-hole tests, as well as boreholes, to enhance the accuracy of the subsoil model. The sections start and end at least 400 m outside their extremes located into the seismic bedrock for minimizing in the simulations the edge effects caused by the lateral dispersion of seismic energy [Tallini et al. 2024]. The seismic bedrock is defined by the Italian technical standards for construction [CS.LL.PP. 2018] as having the shear wave velocity Vs > 800 m/s. Furthermore, to validate simulations of possible 2D phenomena arising from seismic directionality, orthogonal sections have been elaborated and 2D simulated [Lanzo and Silvestri 1999].

In this paper we discuss in detail five of these sections (sections A-A’, B-B’, C-C’ in Figure 1, and D-D’, E-E’ in Figure 2), as they provide compelling examples of local seismic effects related to 2D subsurface conditions. Sections A-A’, B-B’, C-C’, and D-D’ exhibit both conditions well approximated by one-dimensional stratigraphy and scenarios where the influence of two-dimensional subsurface geometry (2D effects) is significant. They also point out how distinct sites with identical stratigraphy but varying seismic layer thicknesses can exhibit differing amplification factors. Section E-E’, conversely, presents a noteworthy case where seismic amplification, albeit moderate, is observed on outcropping bedrock.

3.3. Geotechnical and Geophysical Characteristics of the Units

The geophysical and geotechnical parameters of the units listed in Table 3 and used for numerical simulations represent average values derived from extensive on-site investigations conducted across the L'Aquila Basin area [e.g., MS-AQ Working Group 2010; Monaco et al. 2012; de Magistris et al. 2013] or stored in the QGIS database (Table 3).

Examples of results of cross-hole prospections conducted by previous studies [Del Monaco et al. 2013; Chiaradonna et al. 2024] in representative locations of the study areas (cross-hole A and B in Figure 1 and Figure 2, respectively) are illustrated in Figure 6, in terms of shear wave velocities Vs measured for the main seismo-layers. Finally, in the two studied areas no cases of inversion of seismo-layer velocities with depth were found.

3.4. Seismic Input

The Italian technical standards for construction, [CS.LL.PP. 2018] recommends using seven appropriately selected natural accelerograms as input signals for numerical simulations. This involved running each 1D and 2D simulation seven times, once for each accelerogram, and subsequently averaging the resulting seven spectra to obtain the final output. In accordance with the mentioned technical standards, seven natural accelerograms were selected from the online REXELite database (https://itaca.mi.ingv.it/ItacaNet_40/#/rexel, accessed on November 15, 2024) [Sgobba et al. 2019] using a spectrum matching approach. REXELite is a simplified version of the Rexel database [Iervolino et al. 2009].

The seven seismic events were selected based on the following criteria: (i) magnitude (Mw or Ml) between 5.5 and 7; (ii) extensional fault focal mechanism; (iii) epicentral distance of 0-30 km, and (iv) recording site on seismic bedrock.

The seven inputs are acceptable for the simulations if the average value, calculated for each frequency or period, is compatible with a target spectrum. Following the Italian regulation [CS.LL.PP. 2018], the upper and lower tolerance between the average value and the target spectrum is equal to 30% and 10%, respectively, in the 0.1-1.1 second period range. The target response spectrum used for L’Aquila (latitude: 42.377; longitude: 13.310) is associated with a reference earthquake having a 475-year return period on seismic bedrock (Vs > 800 m/s), assuming horizontal topography and 5% damping. This spectrum is automatically calculated by REXELite when querying the database.

Table 4 shows the set of the seven accelerograms, which comply with the criteria above outlined and were employed as inputs for 1D and 2D numerical simulations for the estimation of AFs in the study of the local seismic response. Figure 7 illustrates the spectra of those accelerograms.

3.5. Approach to Simulation of Real Soil Behavior

Soils subjected to cyclical loading exhibit nonlinear elastic dissipative behavior, manifested as hysteresis loops in the shear stress (τ) versus shear strain (γ) plane. This nonlinear elastic behavior cannot be neglected in earthquake related strain conditions and can be effectively represented using a Kelvin-Voigt viscoelastic model [Kramer 1996; Yoshida 2015]. Numerical analyses addressing soil nonlinearity typically employ one of two approaches [Kramer 1996; Yoshida 2015]: (1) the equivalent linear method or (2) the incremental nonlinear method. The equivalent linear method involves a sequence of linear analyses, where soil nonlinearity is approximated using the secant shear modulus (G), defined as the ratio τ/γ for an opportunely chosen γ value, and the damping ratio (D), which is proportional to the area enclosed by the hysteresis loop (i.e., to the energy loss at each loop).

Conversely, the incremental nonlinear method entails a step-by-step time-domain integration of the motion equations based on a nonlinear Kelvin-Voigt soil constitutive model. This approach, however, presents certain limitations. For instance, established linear analysis techniques, such as the superposition principle, transfer functions, and Fourier analysis, become inapplicable. Moreover, solving nonlinear differential equations is inherently more complex than solving linear ones. Numerical solutions of these nonlinear equations can frequently encounter convergence issues, necessitating numerous iterations. Consequently, a fully nonlinear treatment introduces greater algorithmic complexity and often a substantial computational demand, making codes implementing this approach more difficult to develop, validate, and maintain, resulting in increased costs. To mitigate these challenges, many simulation codes adopt the equivalent linear approach. This method involves integrating the dynamic motion equations (in the time or frequency domain), driven by an input accelerogram representing the target earthquake, using initial estimates for the elastic moduli (usually, the maximum value, denoted by G₀, deduced by in situ V_s measurements) and damping ratios [Kramer 1996]. This integration yields the maximum strain (γ_Max) at each point within the model. Subsequently, effective strain values (γ_eff = α γ_Max), where α is a magnitude-dependent coefficient typically ranging from 0.6 to 0.7) are used to update G and D based on experimental curves describing the nonlinear relationship between shear modulus G, damping ratio D, and strain γ. These updated G and D values are then compared with the previous iteration's estimates. The iterative process continues until the relative differences between successive G and D values fall below a predefined tolerance, at which point the solution is considered converged.

3.6. Numerical Simulation: Code Utilized, Discretization and Boundary Conditions

To balance analysis accuracy and computational efficiency, while considering nonlinear soil behavior, local seismic effects, and the complex geology of L'Aquila Basin, a 2D equivalent linear model was adopted.

While the equivalent linear approach uses simplifications in the calculation, it significantly reduces the complexity of modeling wave propagation in heterogeneous, nonlinear media, enabling efficient and cost-effective software solutions.

Therefore, for the AFs evaluation, the LSR 2D code was utilized to assess the amplification factors (https://www.stacec.com/lsr-2d_pp92.aspx, accessed on 15 November 2024) [STACEC, 2021]. LSR 2D performs 2D equivalent linear modeling using the finite element method in the time domain, total stress, and Kelvin-Voigt model [Kuhlemeyer and Lysmer 1973].

Furthermore, in the conducted 2D analysis using an equivalent linear and concentrated masses approach, the subsoil model was discretized into a mesh composed of triangular or preferably quadrangular elements. In this case study, an adaptive mesh generation approach was used to optimize computational resources by focusing on control points relevant to the desired output results. The mesh step would rise from higher values at the seismic bedrock (about 4 meters) to smaller values (about 1 meter) closer to the control points. An example of the developed mesh is illustrated in Figure 8.

LSR 2D demands as input, for each seismo-layer, the volume weight, the shear modulus, the damping at low strain, the Poisson’s ratio, the G/G₀ vs. γ and D vs. γ curves and lastly the constant α for the calculation of the shear deformation value starting from the maximum one of γ (t) (usually equal to 0.65).

LSR 2D furnishes as output the maximum accelerations at all nodes, the maximum tangential stresses and strains in each element the acceleration time histories at selected nodes for the vertical and horizontal components.

The finite element model thus assembled is subjected to boundary conditions consisting of the application of the input seismic signal (Sect. 3.4) to the base nodes and free-field conditions to the lateral boundary nodes. At the base, the input seismic motion is applied simultaneously to all base nodes as vertically propagating shear waves (in-plane motion). To account for energy dissipation through radiation damping at the model base, linear viscous dampers are incorporated.

Free-field conditions are implemented on the lateral boundaries by coupling viscous dampers between the model's lateral boundary nodes and nodes of appropriate one-dimensional soil columns (free-field columns) capable of representing free-field motion. The coordinates of the model's lateral boundary nodes and the free-field column nodes can coincide. For these boundary conditions to be applicable, the model's lateral boundaries must be vertical. The velocity components of the lateral column nodes at each time step are obtained through a one-dimensional lumped-mass numerical solution of the wave propagation equation. This solution is computed in parallel with the solution of the main model. In other words, for each time integration step, the velocity components of the free-field columns are first calculated via the 1D solution and subsequently converted into loads applied to the main 2D model.

4. Results

The subsoil model for the pilot areas was elaborated by integrating detailed geological information and borehole data (among which several reach the seismic bedrock) with single-station microtremor measurements [Nocentini et al. 2017] and 2D passive microtremor array measurements [Tallini et al. 2024]. This enabled us to estimate the shear wave velocity (Vs) profile for this area and to conjecture also the seismic bedrock depth.

Figure 9, Figure 10, Figure 11, Figure 12 and Figure 13 show the outcomes of 1D and 2D numerical simulations carried out on the sections whose locations are reported in Figure 1 and Figure 2.

The AFs were computed using the LSR 2D code along the sections, spaced every 30-50 meters. For each section, three calculation points were defined within a zone ranging from approximately 70 to 120 meters in width.

Along the sections, the class values of the AFs for the period range of 0.1 to 0.5 seconds were reported taking into account this is evaluated the most critical of the three required intervals (0.1-0.5 s, 0.04-0.8 s, and 0.7-1.1 s) as defined by Italian construction standards [CS.LL.PP. 2018] (Figure 9, Figure 10, Figure 11, Figure 12 and Figure 13).

The red, green, and blue curves refer to AFs in the period ranges of 0.1–0.5 s, 0.04–0.8 s, and 0.7–1.1 s, respectively. The diamond symbol refers to AFs computed using 1D simulation for comparative analysis with the results of 2D simulations.

The formations denoted by CCD, CCF, SCZ, SCZb, and CBZ are carbonatic units, while UAP, UAM1b and UAM3 are terrigenous units. All together pertain to the seismic rocky bedrock (Vs > 800 m/s). The other codes refer to the Plio-Quaternary basin fill detrital deposits (Table 1 and Table 3).

5. Discussion

The investigated area, showing a complex underground structure and stratigraphy, represents a valuable case study which provides interesting insights into the factors controlling site-specific seismic amplification. Below we discuss some sections, together with seismic response simulation results, which offer illustrative examples (Figure 9, Figure 10, Figure 11, Figure 12 and Figure 13).

While amplification factors (AFs) provide a concise summary of local site response, they do not fully capture the complexities of the associated seismic spectra. For instance, two spectra may exhibit identical AFs within a specific period range (e.g., 0.1-0.5 s) but have significantly different peak values, as illustrated in Figure 14, which shows a detail of spectra related to progressive distances (PD) of 1040 and 1263 m in the section A’-A (Figure 9). This is because the AF represents an average amplification over the considered frequency range. Consequently, a comprehensive site-specific seismic response analysis is essential for engineering design, and reliance solely on seismic microzonation data is insufficient. Nevertheless, a combined analysis of AFs across different period ranges can provide a valuable initial assessment of the amplification characteristics, indicating the frequency distribution of amplified ground motion. Based on this assessment, a more detailed site-specific seismic response analysis can be planned.

In section A’-A (Figure 9), a comparison of AFs at PD 311 m and PD 1494 m reveals the influence of all2 layer thickness on the spectral shape. The all2 layer is located in between shallow all3 and the deep located all1, and it shows different thickness. The greater all2 layer thickness, at PD 1494 m with respect to PD 311 m, results in a broader frequency response with lower amplifications at high and mid frequencies (red and green lines) together with higher amplification at low frequencies (blue line). Additionally, at both PDs, the highest AFs are observed at mid frequencies, which can be attributed to impedance contrasts at depths greater than 30 m. Furthermore, 2D effects, induced by both lateral surface, deep-seated variations and buried morphology, are observed at several points within this section. Figure 15a compares the 1D and 2D response spectra at PD 162 m in section A'-A, where a lateral variation in surface lithology within the flat plain is present. The transition from an outcropping all2 layer (Vs ≈ 450 m/s) to an outcropping all3 layer (Vs ≈ 250 m/s) induces 2D effects, as evidenced by the differences between the 1D and 2D spectra and their corresponding AFs. Conversely, at PDs 481 m and 1145 m within the same section (Figure 15b,c), the observed differences in the 1D and 2D spectra are attributable to lateral variations in deep-seated structures, such as variable thickness of the all2 layer and the presence of a steeply inclined bedrock and/or buried morphology.

Figure 10 points out an inverse relationship between high-frequency AF values (red line) and the depth of the shallowest impedance contrast. As we move towards the Raio Stream thalweg, the decreasing thickness of the upper layer due to the presence of an alluvial terrace edge (ATE) leads to an increase in AF values. This is evident in the transition from lower AFs (1.7) in the left section side to higher AFs (2.3) in the right side.

The seismic response estimated for the shallow valley in Figure 11, is predominantly 1D, with minimal 2D effects. The close agreement between the 1D and 2D AFs indicates that 1D modeling adequately captures the seismic response, except for regions near the valley edges where subtle 2D effects may arise due to the gradual transition to the rocky substrate. The 1D AF values at the valley center exhibit greater amplification at middle (green line) and lower frequencies (blue line) and are primarily influenced by the impedance contrast between all2 and all1.

Figure 12 illustrates a pronounced 2D basin edge effect near the southwestern edge of the valley. At PD 856 m of section D-D’, the 2D AF values significantly exceed the 1D AF values and the related 1D and 2D spectra appear substantially different (Figure 16a), indicating a strong influence of the valley boundary. For comparison, Figure 16b,c illustrate, in the same section, two cases in which 1D and 2D spectra exhibit modest differences (PD 635 m) or are identical (PD 20 m). This trend is observed in various sections within the Bazzano-Monticchio area, where 2D effects are more prominent in the Aterno River plain.

Figure 13 suggests a potential topographic effect on a slope where the terrigenous rocky substratum outcrops. The higher 2D AF value of 1.3 (which although being modest is not negligible) compared to the expected value of 1.0 indicates local amplification due to topographic irregularities.

To validate our numerical simulations of local site response, we compared the results with microtremor analysis using the Horizontal to Vertical Spectral Ratio (HVSR) technique [Nakamura 1989, 2019]. The HVSR method provides a proxy for the dominant resonant frequencies of a site by analyzing the Fourier spectra of the horizontal and vertical components of ground motion associated with microtremors or seismic events. Figure 17 shows a comparison between the transfer function calculated from a 1D simulation at PD 1494 m in section A'-A and the HVSR plot for a nearby site (066049-P183 in Figure 5). The transfer function exhibits a primary peak at approximately 1.5 Hz, followed by a secondary peak at around 3.5 Hz. The HVSR plot clearly shows a primary peak at 2 Hz, with a more complex pattern above this frequency, possibly indicating overlapping secondary peaks. The overall shape of the HVSR plot is similar to the transfer function, with a slight rightward shift of the primary peak. This shift can be attributed to the lower deformation levels associated with HVSR waves compared to those assumed in numerical simulations of local site response. The significant deformations hypothesized in these simulations lead to reductions in elastic moduli and shear wave velocities, consequently lowering resonant frequencies [Lanzo and Silvestri 1999]. Since HVSR analysis provides a direct, albeit approximate, measure of site resonant frequencies, the compatibility between HVSR results and numerical simulation results confirms the reliability of the subsurface model used in the simulations.

Several key considerations emerge for best practices in local seismic response analysis and seismic microzonation studies in Italy. AFs provide a valuable tool for urban planning by assigning an AF value to each microzone, so pointing out for them the extent of amplification [SM Working Group 2008]. Conversely, site-specific engineering design requires detailed seismic response analyses at the site scale, by using response spectra, because sites characterized by the same AF may show different spectrum shape (Figure 14). Furthermore, for microzones showing high AF values, numerical simulation is preferable to simplified method allowed by the Italian technical standards for construction [CS.LL.PP. 2018]. Given the prevalence of 2D basin edge effects in Plio-Quaternary geologically complex intermontane basins in central Italy, as demonstrated in the L'Aquila Basin case study (Figure 9, Figure 12, Figure 15 and Figure 16), if the stratigraphy of the site under investigation is not confidently known to be horizontally layered, then 2D simulation models are recommended. Furthermore, our analysis has revealed topographic effects, even in areas with moderate rocky slopes (Figure 13). In such cases, 2D simulations are also advisable. Finally, to validate the numerical models and their associated spectra and AFs, it is recommended to compare transfer functions with HVSR microtremor measurements at control points along the simulated section. Figure 17 illustrates a comparison between the transfer function obtained for the Progressive Distance 1494 m of section A-A’ and the HVSR plot. Here we use the results from 1D simulation, as this site is in the middle of the valley and subsoil exhibits a horizontally layered structure.

6. Conclusions

We illustrate a detailed analysis of site effects in the highly seismic region of L'Aquila (central Italy) characterized by recent and significant seismicity, such as the 2009 L'Aquila earthquake (Mw 6.3) and the 2016 Central Italy seismic sequence (Mw 6.0 and 6.5). An accurate evaluation of local seismic amplification was assessed by integrating detailed geophysical and geotechnical data with numerical modelling. Amplification factors for different period ranges were estimated through 2D numerical simulations using the LSR 2D code on many (n. 15) 3500-5000 m long geological sections. Five of these sections have been illustrated and discussed in this paper.

Several key considerations emerge for best practices in local seismic response analysis and seismic microzonation studies in Italy. While AFs provide valuable tools for urban planning, site-specific engineering design always requires detailed seismic local response analyses, as sites characterized by the same AF may show substantially different spectrum shape. Furthermore, given the prevalence of 2D effects in Plio-Quaternary geologically complex intermontane basins in central Italy, as evidenced by the L'Aquila Basin case study, 2D simulation models are generally preferable against 1D ones. In several cases, we found that the local seismic response is significantly affected by structures located at depths greater than 30 m. Finally, to validate the utilized numerical models and their associated spectra and AFs, it is suggested to compare transfer functions with HVSR microtremor measurements at control points along the studied sections.