Hydrologic control on natural land subsidence in the shallow coastal aquifer of the Ravenna coast, Italy

Multiple processes contributing to natural land subsidence in a shallow coastal aquifer near Ravenna (Italy) were identified by analysing the relationships among different data set time series (water table level, rainfall, land reclamation drainage, sea level, etc.) and establishing the correlations with vertical ground motion observed at a high-resolution settlement gauge. Our study highlights the presence of three deformation components related to different processes controlling land subsidence: elastic, delayed-elastic, and irreversible (plastic) components. The elastic and delayed-elastic components are closely related to water table fluctuations that change the effective stress in two portions of the coastal aquifer at a daily (in the sandy unconfined portion) and seasonal time scales (in the layered clay-rich semi-confined prodelta portion), respectively. The irreversible component represents the trend in the land subsidence time series and is due to primary consolidation (pore pressure dissipation) of the fine-grained prodelta levels above where the settlement gauge is located. The amplitudes of the elastic component can be up to 0.2–0.3 mm whereas the amplitude of the delayed-elastic component reaches 0.89 mm. The primary consolidation rate of deformation is 0.9 mm yr⁻¹ and constrains the likely age of prodelta sediments deposition to 1300–2800 years before present. The delayed-elastic subsidence rate has similar magnitude to that due to primary consolidation and is connected to poroelastic effects in the prodelta sequence following seasonal variations in water table. Our findings are important for planning land subsidence management and monitoring strategies especially where the surface aquifer structure is heterogeneous due to different depositional settings. The natural land subsidence rate in the Holocene sediments of the shallow coastal aquifer of Ravenna (North eastern Italy) that we measured in this study accounts for 10 %–20 % of the total current land subsidence rate observed in this portion of Ravenna coastal area (10–20 mm yr⁻¹).

1 Introduction

By using data from a high-resolution settlement gauge, our present work aims to define the processes governing shallow ground settlement (magnitude and development over time since the installation of the instrument) and verify land subsidence and water table fluctuations interactions in the shallow coastal aquifer of Ravenna. Our work highlights the contribution of natural processes such as primary consolidation (compaction of sediments under their own weight via expulsion of interstitial pore water) and water table fluctuations in the Holocene shallow coastal aquifer of Ravenna to the cumulative land subsidence rate observed in the area (Bertoni et al., 1995; Teatini et al., 2005; Baldi et al., 2009). The processes in the shallow coastal aquifer are uncoupled from other processes contributing to the cumulative land subsidence rate and that are primary and secondary consolidation of deep aquifers (Teatini et al., 2011), fluids extraction from reservoirs (Teatini et al., 2006), tectonism, and isostasy (Carminati et al., 2003). The time series analysis and correlation of hydrologic parameters with settlement from high-resolution gages are powerful tools to investigate the coupling of different processes within a shallow aquifer and to constrain the time of its deposition (Antonellini et al., 2019).

Study area

The area addressed by the present study includes the coastal area of the Ravenna city, in the Emilia-Romagna coastland, south of the Po River Delta (Northeastern Italy, Fig. 1). It is a lowland coastal area not exceeding 2 m a.s.l. (above sea level), with a large portion below mean sea level, because of the combined effects of natural and anthropogenic land subsidence, land reclamation, and sea level rise.

Figure 1 (a) Study area and location of the monitoring stations. The green dots represent the locations of data acquisition: MF1_3 piezometer, “Settlement gauge + SCPT (Standard Cone Penetration Test)” see Fig. 2; (b) Subsidence rate and isokinetic isolines (2006–2011) of the study area (source: ARPAe database: https://arpae.it/cartografia/, last access: 15 July 2018; from Antonellini et al., 2019).

In most studies performed in the area there is evidence of a close relationship between land subsidence, deep groundwater pumping, and deep gas field development (Teatini et al., 2005, 2006).

The shallow coastal aquifer, where the instruments of this research are installed, consists of Holocene sediments; there are three fundamental units (Amorosi et al., 1999a, b; Campo et al., 2017), as represented in the stratigraphic core log and Standard Cone Penetration Test (SCPT) in Fig. 2. Two sandy units, one at the surface, made up by medium- to fine-grained sand and with a thickness of about 8 m, and a deeper one (from −20 to −25 m a.s.l.), consisting of fine sand with very thin clay inter-layers. The intermediate unit is a package of alternating clay-silt and sandy-silt strata (prodelta sediments). The basement of the unconfined coastal aquifer consists of Pleistocene impermeable clay at a depth of 26–30 m a.s.l. In this paper, when we refer to the “shallow coastal aquifer of Ravenna”, we mean the whole sediment package described above (depth 0–30 m). It is apparent, however, that the upper part (0–8 m depth) is unconfined (phreatic), whereas below 8 m depth the aquifer is semi-confined or consists of aquitards.

Figure 2 Stratigraphic core log of the borehole where the settlement gauge is installed and Standard Cone Penetration Test (SCPT) carried out next to the log. Depositional environments (dep. env.): floodplain (F), levee and crevasse (L/C), marsh (M), dune backshore and foreshore (D/B/F), upper shoreface (US), lower shoreface (LS), offshore and prodelta (OP) and lagoon (L). Lithology: P = peat, C = clay, S = silt, S = sand f = fine, m = medium, c = coarse. CPT parameters: qt = cone resistance, fs = lateral friction, U = pore pressure, Vs = shear wave velocity, z = depth, H_a = depth of the settlement gauge, H_b = depth of the bottom of the first sand package, and H_dr = equivalent thickness of clay level (from Antonellini et al., 2019).

2 Methods

We used a classical decomposition by moving averages based on an additive model (Kendall et al., 1983) of data time series (subsidence-settlement, water table level, precipitation, drainage, sea level) into the trend, seasonality and noise components and find correlation coefficients between each analysed component, especially between settlement and changes in water levels measured in a nearby piezometer (Fig. 1). We also used the effective stress principle (Geertsma, 1957) and Terzaghi (1923) equation in the closed-form solution given by Singh and Swamee (2008) to model the irrecoverable and elastic components of land subsidence.

A TS 377 AGISCO settlement gauge is anchored 30 m below ground with the aim to evaluate the contribution of the Holocene compaction to the land subsidence (Fig. 2). The gauge is equipped with a 15 mm inductive transducer and the recorded time series of settlement data consists of 4 measurements per day. In addition, the dominant parameters controlling water level in the piezometer (one at a depth of 3 m open at the bottom used for the correlations and one in the upper sandy portion of the surface aquifer used for control) were analyzed: the stage level (m) of the Bevano river; the atmospheric pressure (Pa), and rainfall (mm) at the “ARPAe” Ravenna weather station; the daily discharge at the “V Basin” pumping station; the water table elevation (m), groundwater temperature (^∘C), and salinity (g L⁻¹) at the MF1–3 piezometer; the sea level (m) at the Porto Garibaldi tide gauge station of the Emilia-Romagna.

Given that the frequency of data is constant over time, we decomposed the daily time series using an additive model into seasonal, trend and irregular components using moving averages (Kendall et al., 1983).

Correlation tables among all data series were elaborated for daily and monthly time series (Table 1); lag times among correlated parameters were also computed.

Table 1 Monthly correlation coefficients (from Antonellini et al., 2019). Bold values are for a very good correlation and NS means not significant.

3 Results

Monthly correlation coefficients calculated for time series dataset are listed in Table 1. The time series data for the most relevant (and correlated) parameters recorded – ground movement (mm), water table level (m), drainage (mm), sea level (m) – are shown in Fig. 3a.

Figure 3 (a) Time series of recorded data. (b) Trend component of the settlement time series data (solid black line) compared with the 1-D best fit (R²=0.99) consolidation model (dashed red line line). (c) Seasonal and noise (d) components of the settlement, land reclamation drainage in mm equivalent of water per unit area, sea level and water table time series data (from Antonellini et al., 2019).

The results of the decomposition of the settlement daily time series are presented in Fig. 3b–d. The seasonal component (Fig. 3c) of the settlement has an absolute amplitude of about 0.89 mm (from −0.41 to 0.48 mm) with a maximum positive value (0.48 mm) at the end of May and a secondary maximum (−0.05 mm) in mid-November, whereas the principal minimum (−0.41 mm) is at the beginning of September and a secondary one (−0.38 mm) is at the end of December. The seasonal component of the water table fluctuation has an amplitude of about 0.60 m (from −0.19 to 0.41 m) with a maximum positive value (0.41 m) at the beginning of March and a secondary maximum (0.3 m) at the end of October, whereas the principal minimum (−0.19 m) is at the end of August and a secondary one (−0.08 m) is at the end of December. The noise components of settlement and water table are correlated (Fig. 3d).

The trend component of settlement (Fig. 3b) fits (R²=0.99) the 1-D consolidation model based on Terzaghi (1923) equation when a consolidation coefficient $c_v=\mathrm{7.1}\times {\mathrm{10}}^{-\mathrm{3}}$ [m² d⁻¹] is used. The trend settlement rate is −0.91 [mm yr⁻¹].

Table 1 shows the correlation at monthly level between settlement and other hydrologic parameters measured, which helps in smoothing lags due to the daily observations. The settlement has a positive correlation with water table level (0.44), Bevano River stage (0.49), drainage (0.55), rainfall (0.46), and water salinity (0.49).

4 Discussion and conclusions

The data we present are important because they allow to separate and quantify the effects of natural processes (primary consolidation, water table fluctuations, etc.) in the heterogeneous Holocene aquifer sediments from the cumulative land subsidence rate observed in the area, which includes other natural (secondary consolidation in deep aquifers, tectonism, isostasy, etc.) and anthropogenic effects (gas and fluids extraction from deep aquifers and reservoirs).

One of the important findings of our research is that the daily ground movement variations are related to the water table fluctuations (Fig. 3a and Table 1). This is shown by the similarity of the two curves in Fig. 3a and by the correlation coefficients in Table 1. The positive correlations observed in Table 1 and described in the results section are expected since the daily water table fluctuations modify the effective stress, which induces an almost immediate elastic response in the upper sandy portion (8 m-thick) of the coastal aquifer (upward-directed and therefore in the positive direction). The other correlations follow from the direct links among precipitation, drainage, river stage, and water table level; the water salinity correlation is probably related to dissolution of salts from the vadose zone as the water table rises. To the best of our knowledge, the elastic effect described in the aquifer was not documented before.

Another important result from time series analysis is the identification of a seasonal recoverable (time-dependent) component in vertical displacement with an amplitude of 0.89 mm (Fig. 3c), which is similar in magnitude to the yearly irrecoverable consolidation rate (0.9 mm yr⁻¹) (Fig. 3b). We propose that this displacement component is delayed-elastic given its time-dependency and its delayed coupling to effective stress variations. This poroelastic effect occurs in the prodelta aquitard following seasonal water table changes. The elastic seasonal component of vertical displacement has a positive (uplift) value that initiates with the seasonal water table rise (Fig. 3c) and then persists after the water table has fallen (for about 3 months). This behaviour points out an increase in pore pressure within the prodelta section of the aquifer, starting in the period when the water table is high, followed by a slow seepage and outflow after the water table has decayed. The relationships between daily water table fluctuations and daily land subsidence are almost instantaneous as shown in Fig. 3a; the best r value achieved after a lag of 5 d is only slightly higher than that obtained with a 1 d lag. This is an indication of elastic deformation in the high-hydraulic-conductivity (k>5 m d⁻¹) sandy part of the aquifer (upper 8 m).

The trend component of land subsidence (Fig. 3b) allows us to compute the total irrecoverable settlement, which occurred in the shallow coastal aquifer (from 0 to 30 m below the ground surface) during the monitoring interval time (May 2014–November 2016) and that is 2.29 mm. The consolidation rate of the Holocene succession, therefore, is 0.9 mm yr⁻¹ and is comparable to the lowest rates observed in the Holocene sequence of the Po Delta (Teatini et al., 2011). The consolidation rate that we measured is likely connected to residual primary consolidation of the prodelta clays (Teatini et al., 2011). By using consolidation coefficients in the range varying from 1.58 to 3.15 m² yr⁻¹, which is what proposed by Teatini et al. (2011), we obtain an average degree of consolidation U_m varying from 0.8 to 0.99. This average degree of consolidation is compatible with values reported in the literature (Teatini et al., 2011), with the geotechnical characteristics of the prodelta clays (pore index), and the irrecoverable land subsidence rate that we identified with time series analysis. Therefore, we can constrain the age of the prodelta clays to an interval of 1300–2800 years before present, which is compatible with the range given in the literature (Amorosi et al., 1999b; Campo et al., 2017).

Our data, however, highlight the contribution of the shallow Holocene coastal aquifer to the overall land subsidence rate, which is not negligible. It represents about 10 %–20 % of the total land subsidence rate in the area (10–20 mm yr⁻¹; Teatini et al., 2005). The elastic daily response of land subsidence to water table fluctuations might have a maximum amplitude of 0.2 mm and the seasonal delayed-elastic component of land subsidence has an amplitude (0.89 mm) that is similar to the irrecoverable yearly rate (0.9 mm yr⁻¹).

Recognition of the recoverable elastic and delayed-elastic vertical displacement components at different time periods (daily to seasonal) is complementary to and may help constraining uncertainties associated to data collected with areal monitoring land subsidence techniques (GPS, InSAR) at sparse time intervals (seasonal to multiple years). The methodology and the findings that we present here show the advantages of using settlement gauges to characterize the effect of aquifer heterogeneity on land subsidence and can be used in other similar settings around the world where multiple processes contribute to the cumulative land subsidence rate.