Geomechanical modelling and stability analysis of the shallow underground water reservoir ‘Palombaro Lungo’ (Matera-Italy)

The city of Matera, as several other ancient historical cities and towns, is recovering a large part of its ancient parts developing a policy of urban regeneration. In the early 90s, during the restoration works of the main square, Piazza Vittorio Veneto, it was re-discovered a huge underground space used as water reservoir named Palombaro Lungo. It was built linking ancient pre-existing smaller caves and completed at the end of the nineteenth century making it impermeable with the cocciopesto (opus signinum) technique. The hypogeum is dug in a Plio-Pleistocene rock called Gravina Calcarenite. It is a calcareous sandstone weak rock with good mechanical properties, low permeability and easy to be dug. The internal geometry of the hypogeum and the rock thicknesses above the caves was reconstructed by the integration of topographic and GPR surveys. The paper presents a geo-mechanical model of Palombaro Lungo underground reservoir and evaluation about the stability of this structure. The stability analysis were developed using the numerical code UDEC for several transect of Palombaro Lungo assuming the calcareous sandstone blocks as deformable, trying to address the evolution of the stress–strain conditions. The case study can be considered as an example of general interest for the study of rupestrian underground reservoir of a longed form, where one planimetric dimension strictly prevails over the other. Numerical simulations showed a stress–strain state compatible with the fractures detected in situ and confirmed the absence of instability problems in the groundwater reservoir.


Introduction
Several ancient historical cities and towns are developing a policy of urban regeneration with the recovery and the restoration of their underground spaces. The use of dug underground space for several purpose was a really common purpose in Mediterranean country where the rupestrian civilization develop (Margottini et al. 2017). Stability studies for restoration work of rupestrian settlement in urban area have studied in several other Mediterranean areas (Deniz and Topal 2021;Hatzor et al. 2006Hatzor et al. , 2010 that have geological conditions similar to that of Matera. These restoration work have been really relevant for the city of Matera (Southern Italy) whose ancient district called Sassi, can be considered as a sort of capital of rupestrian settlement. This ancient part of the city is completely dug in a Plio-Pleistocene calcareous sandstone, a weak rock called Gravina Calcarenite. It is a weak rock with good mechanical properties with low permeability easy to be dug.
In this paper is presented a study of the stability of a huge underground space used as water reservoir named Palombaro Lungo as an example of the use and reliability of the 2D approach for the study of the stability of rupestrian underground reservoir of a longed form, where one planimetric dimension strictly prevails over the other. So, the result here presented are of general interest for the study of the stability of underground reservoir not only in Matera, but also for other sites where there are similar artefacts.

Historical notes and geometrical features of the Palombaro Lungo
The hypogeum of the Palombaro Lungo is a large hydrodrinking reservoir built in the outcropping calcareous sandstone of Piazza Vittorio Veneto starting from the sixteenth century, whose last documented extension dates back to 1870 (Fig. 1).
Its current planimetric and volumetric articulation arises from the progressive coalescence of several underground cavities and pre-existing cisterns, and it is only one of the many underground environments dug into the calcareous sandstone rock mass present in the subsoil of the square. The tank walls are waterproofed with a 2-3 cm thick crushed earthenware cocciopesto lining and, at the deepest point, it is about fifteen meters high. It is estimated that, at the maximum filling level, it can contain about 5000 m 3 of water (Fig. 2).

Geological-structural features of the Matera area
From a geological point of view, Matera's area represents a singular crossing point between the domain of Foredeep and Foreland. The stratigraphic sequence (Fig. 3) is characterized by a substratum of Cretaceous limestones of the Apulian carbonate platform, covered by the sequence of the Bradanic sedimentary cycle and Quaternary deposits (Cotecchia and Grassi 1975;Ciaranfi et al. 1988;Pieri et al. 1996;Simeone et al. 2019).
Cretaceous limestones of Altamura Limestone Formation is the geological substratum. These are compact dolomitic limestone and limestone significantly fractured that give evidence of a more intense and significant tectonic history than imaginable from the current seismicity of the area. The lithotype has layers that on average dip towards SW with inclinations of about 10° and variable thicknesses from a few centimetres up to 2 m. Plio-Quaternary deposits of the Bradanic Trough sedimentary cycle rest on the limestone substratum. Calcareous sandstone are the lower Plio-Pleistocene sedimentary term. They are locally called Fig. 1 Overview of Piazza Vittorio Veneto with the entrance to the underground reservoir of the Palombaro Lungo and a detail of the area during the work of arranging the pavement of the square "tuffs", ascribable to the geological formation of Gravina Calcarenite (Auct.). These calcareous sandstone are weak rock with a variable color from yellowish to gray-whitish. Generally they are massive, but significantly fractured, probably both as a result of the tectonic uplift and the complex recent structural tectonic events that have led to the current geological configuration of the area (Doglioni et al. 2020). The calcareous sandstone are almost impermeable when not fractured and are self-supporting to excavation. They allow easy and precise excavation and thanks to these characteristics they have allowed the development in the rock mass of the typical rock settlements of the city, including the historic district Sassi.
On the calcarenites, in the westernmost portions of the city, lay, in stratigraphic continuity, the Sub-Apennine Clays (Auct.) which form the base of the hills and then the coarser lithotypes that closes the sedimentary cycle of the Bradanica Foredeep.
The Matera's area has evident and important tectonic dislocations resulting from a polyphasic tectonics that deeply disturbed both the Mesozoic substrate and the Plio-Quaternary sedimentary coverings. Tectonic fractures are clearly evident in stone lithotypes such as the Altamura limestone and the Gravina Calcarenite but are hardly or not at all distinguishable in soils with plastic behavior such as Sub-Apennine clays.
The tectonic structures have also significantly affected the morphological and hydrographic evolution of the Matera's area, as in the case of the "Gravina di Matera" watercourse, a deep morphological incision profoundly conditioned by tectonics (Doglioni et al. 2020).

Geo-morphological layout of the hypogeum area
The area of Piazza Vittorio Veneto, below which the hypogeum of the Palombaro Lungo develops, is characterized in its entirety by outcrops of the Gravina Calcarenite which, from the western edge of the square and proceeding towards the west, pass to the overhanging ones Sub-Apennine clays whose thickness tends to grow moving towards WNW. The area of the current Piazza Vittorio Veneto was characterized by a very articulated morphology; as evidenced by the two views left after the renovation of the square in the early 90s of the twentieth century with the presence of numerous underground rooms, which are still present today below the floor of the square. The area was crossed by some outflow routes, which conveyed both the runoff waters and those that flowed into the surface coulters towards the area of the current square and then into the Grabiglione of the Sasso Barisano. In fact, in this area there were activities that required an abundant availability of water, such as tanneries and the iron warehouse.
The area of Piazza Vittorio Veneto occupies the upper portion of the morphological incision where the complex system of hypogea of the Sasso Barisano develops. Currently it is sub-flat, with a slight inclination towards the SSW, but originally it had a somewhat articulated morphology as mentioned in the introduction. The strong anthropic reworking of the subsoil was significant, through the construction of a succession of hypogea and the neighborhoods. In the square, all around there are other hypogea ( Fig. 1): San Domenico, Spirito Santo and Fondaco di Mezzo (Bruno et al. 2016), dug in the same calcareous sandstone of the subsoil of the square.
The calcarenitic rock mass in this area is mainly massive and has a maximum thickness of about 40 m. The geological-structural arrangement of the calcarenitic rock masses shows a generalized dip of about 10° towards the West and WSW; the rock mass is also characterized by a significant state of fracturing whose intensity is locally variable. The carbonate content of the rock, according to the most recent findings in scientific literature, varies from 80 to 99% and richly fossiliferous horizons are present locally. Overall, the rock has a highly variable degree of diagenesis and, consequently, of mechanical resistance also in relation to its natural water content.

Geometrical survey
For the reconstruction of the geometry of the Palombaro Lungo and the adjacent hypogea it was possible to refer to the plans of previous topographical surveys, the last of which developed for the interventions of arrangement of the hypogea of Piazza Vittorio Veneto as well as the geophysical surveys carried out for the study of the hypogea of San Domenico (Bruno et al. 2016). However, to improve the level of reliability of the survey, given the complex geometry of the hypogea in the area in question, a revision was carried out through precise measurements and direct observations in situ of the real shapes and extensions of the underground environments.
From the measurements and observations on site it was possible to define the plan of the hypogeum and four sections A-A′ (longitudinal) and B-B′, C-C′, D-D′ (transversals) used for the construction of the geomechanical models of the stability verification sections of the Palombaro Lungo (Fig. 4). The precise values of the thicknesses of the calcarenitic roof covering the hypogeum in correspondence with the area of the former fountain (section D-D′), installed in 1993 and removed in 2018, were deduced from a mesh of GPR strips carried out in correspondence with this area (Cardone and Cristallo 2018).

Geomechanical characterization of the rock masses in the Palombaro Lungo area
The stress-strain behavior of the rock mass of the geomechanical model, used for stability checks, is strongly conditioned both by the mechanical strength of the intact rock, by the presence of the less cemented material zones and by the discontinuities that cross it. A further complication is determined by the numerous voids in the hypogea adjacent to the one under study, which cause alterations in the stress field that create local concentrations of stresses. The widespread presence in the area of underground cavities and discontinuities favors, among other things, the infiltration of water into the subsoil and this gives rise to a degradation of the mechanical characteristics of the rock mass, especially the calcarenitic one.
For the geomechanical characterization of calcarenite lithotypes, reference was made to the extensive literature on Matera calcarenites and similar lithotypes from Southern Italy (Cotecchia and Grassi 1975;Baldassarre 1990;Cherubini et al. 1996Cherubini et al. , 2007Bruno and Cherubini 2005;Bruno and Rotolo 2018;Simeone et al. 2019), which highlighted how the physical-mechanical characteristics of Gravina Calcarenites in the Matera area, are highly variable depending on the content of calcareous cement, the granulometric assortment and the degree of saturation (Table 1). In particular, the value of the uniaxial compressive strength of these calcarenites allows them to be classified as a soft rock (Geological Society of London, Engineering Group Working Party 1970).
According to the fracturing state of the calcareous sandstone (calcarenite) rock mass, where Palombaro Lungo was dug, it was modeled as a discontinuous medium using the mechanical parameters of both the intact rock and the discontinuities. For this purpose an in situ geomechanical survey was performed, with the generalized scanline method   (Bruno 2012), which allowed to define the geometric characteristics of the discontinuities and the parameters necessary to classify the rock masses of the calculation model. In particular, through on-site surveys, the numerous discontinuities present in the hypogeal environments around and within the Palombaro Lungo have been identified and characterized, which can be classified into 5 families K 1 , K 2 , K 3 , K 4 , K 5 (Fig. 5), including one oriented NW-SE, affects the Palombaro Lungo for almost its entire length ( Fig. 2 red arrow). The geomechanical classification index used is the basic RMR '79-'89 (Bieniawski 1989), i.e. without considering the discontinuity orientation parameter, as required by the specific calculation code used for the numerical stability analyses of the case study. From the classification index RMR '79-'89 (Tables 2), through the correlations of the literature between the geomechanical classifications, the GSI index was obtained and by linearizing the Hoek-Brown resistance criterion the mechanical parameters of the rock masses (Table 3).
Geomechanical parameters of both substratum limestone and dug calcareous sandstone were evaluated using literature data (Bruno and Carucci 2019;Lupo and Gallipoli 2012) and present study survey data (Tables 1 and 2). These data were interpreted using formulas suggested in the UDEC code user manual to obtain reliable, conservative values (Tables 4,5). The values assumed for the dynamic parameters of the lithotypes are shown in Table 6.

Geomechanical models for distinct elements analysis
To evaluate the degree of stability of the underground cavities that make up the hypogeum of the Palombaro Lungo, three of the four sections made were considered, one longitudinal A-A′ and two transversal B-B′ and D-D′, which are extend laterally beyond the limits of the hypogeum (Fig. 5). This is to ensure that the mechanical response of the calculation models was not influenced by the edge effects and, at the same time, to take into account the influence of the adjacent hypogea which represent a relevant buried morphological variation of the rock mass both for checks static and, above all, for the dynamic ones and the consequent effects of local seismic response. For the purposes of stability analyses, section D-D′ is considered particularly important, because it is located at the point of greatest depth of the Palombaro Lungo and, at the same time, where the coverage of the underground environments is minimal. The section A-A′, although significant for the description of the underground environment, is certainly much less significant for the description of the stability conditions because it does not take into account the three-dimensional effects of the distribution of stresses in the rock mass that, for elongated structures such as the one under study are of great importance (Simeone et al. 2019). Observing the geomechanical models of the three verification sections (Fig. 6), it is noted that the calcarenitic rock mass has fewer fractures in the section A-A′ with respect to sections B-B′ and D-D′. This is not due to a change in the state of fracturing of the rock mass but to the subparallel orientation that section A-A′ presents with respect to the discontinuity families K 1 and K 4 . Another peculiarity of the different environments that make up the Palombaro Lungo is that in sections A-A′ and D-D′, in correspondence with the SE environment, the thickness of the roof is reduced to a minimum value of a few tens of centimeters (Figs. 4 and 6). The variability and the extreme reduction of the thickness of the hypogeum roof together with the disproportionate longitudinal dimension of the cavity compared to the transversal one mean that the section A-A′ is little or not at all suitable for 2D analysis of the stability of the underground reservoir of the Palombaro Lungo. Section A-A′, however, was equally investigated in order to clearly verify the feasibility and limits of a 2D analysis of such a geomechanical model.

Stability analysis and discussion
The assessment of the stability and geomechanical behavior of the Palombaro Lungo underground water reservoir was developed with the two-dimensional calculation code UDEC 7.0 "Universal Distinct Element Code" (Itasca 2019). The analyses were conducted considering the deformable blocks and adopting a constitutive model for intact rock, respectively, of the strain-softening type for the calcarenitic mass and linear elastic-plastic for the limestone mass. As regards the discontinuities of the rock mass, however, the Mohr-Coulomb "contact area" type was adopted with the residual strength option exclusively for those involving the calcarenitic rock mass (Fig. 6).
The static checks of the analyzed sections were performed in three steps to reproduce the evolution of the stress-strain state suffered by the rock masses of the geomechanical   In particular, in the first step the stress state existing in the rock mass before excavation was simulated, but taking into account the numerous mechanical discontinuities detected, in the second step the stress state resulting from the construction of the hypogea was calculated, and finally in the third step the effect of the application of anthropogenic loads was evaluated.
Step 1-Condition before the realization of the hypogea (natural tensional state).
The first step of the stability analysis consisted in bringing the models of the verification sections into the stress-strain state prior to the excavation of the currently existing underground cavities. The models of the three sections, without cavities and subject only to the lithostatic pressure of the rocks, undergo slight deformations and reach equilibrium as can be seen from the zeroing of the unbalanced forces and maximum speeds (Fig. 7 on the right). The presence of the discontinuities in the system does not generate significant shifts in the blocks delimited by them; in fact, the maximum values obtained for the total displacements are of the order of 10 -3 ÷ 10 -4 m in all three sections (Fig. 7 on the left).
Step 2-Excavation of the hypogea and analysis of the induced tensional state.
The second step simulated the excavation of the hypogea. The excavation of the rock was simulated as performed all at the same time, even if the creation of this underground environment is the result of the union of numerous pre-existing underground environments, joined to form a single volume. The excavation therefore developed over several centuries and then unified only at the end of the nineteenth century, however, it is believed that the simulation of a single excavation does not affect the assessments made.
It has been decided to develop only a 2D stability analysis according to the potentialities of available calculation software and also because elongated shape of the dug underground space made this sort of evaluation really reliable. The assessments of the stress and displacement state resulting from the excavation of the hypogeum were, as could be expected, really similar to the cross sections B-B' and D-D' (Fig. 8) that can be considered representative for a 2D analysis of the stability of the dug underground space. Moreover, the longitudinal-section A-A′ has been analyzed to show how this section cannot be considered representative in a 2D analysis.
In the two cross sections, the stress relief gives rise to a generalized phenomenon of convergence of the walls and base of the hypogea, with maximum upward displacements of the order of 10 -3 ÷ 10 -4 m ( Fig. 8 on the left). After the displacements, the two sections reach equilibrium as can be seen from the zeroing of unbalanced forces and maximum speeds (Fig. 8 on the right).
In the longitudinal section A-A′, however, the displacements are concentrated only on the roof of the Palombaro Lungo with maximum downward displacements of the order of 10 -1 m (Fig. 8 on the left). The extent of the displacements in this section is clearly incompatible with stability conditions as evidenced by the trend of the unbalanced forces and the maximum speeds that do not tend to zero (Fig. 8 on the right). Furthermore, the excavation of the cavities determines a redistribution of the stress acting in the rock mass which induce states of localized plasticization exclusively to the hypogeum roof of the longitudinal section A-A′ (Fig. 9). In numerical models fracture inclination angles were corrected as a function of the angle between the section trace and the fracture strike according to the approach of Bruno (2012).
The state of instability that the simulation returns, for the roof of the longitudinal section A-A′ of the Palombaro Lungo is not, however, reflected in reality. The deformations related to the decompression of the overall excavation, despite having determined the onset of a longitudinal fracture in the vault of the hypogeum (Fig. 2), did not induce an instability and/or collapse of the vault. Moreover, even over time it has not undergone any changes, as can be easily seen from the state of the crushed earthenware cocciopesto lining present, despite the extensive use of the spaces above.
The verification of the stability of section A-A′ was also performed by simulating the presence of a 3 cm thick crushed earthenware cocciopesto lining. The results obtained are displayed in terms of displacement (Fig. 10  on the left), unbalanced forces, state of plasticization and maximum velocities (Fig. 10 on the right). The result give evidence of the not representativeness of the cross section A-A′ in a 2D stability analysis. It show that representative section are only that ones normal to the elongation direction.
In consideration of what has been said, in the third step, therefore, it was decided to omit the longitudinal section A-A′ of the Palombaro Lungo considering, for further stability analyses, only the cross sections B-B′ and D-D′.  6 Verification section locations and relative geomechanical models adopted for intact rock blocks and discontinuities Step 3-Application of permanent and accidental overloads to the calculation model.
In the third step, the loads that currently affect the environments under study were introduced for the two sections considered B-B′ and D-D′ (Fig. 11).
• permanent loads, due to the pavement of Piazza Vittorio Veneto and the underlying reinforced concrete plate; • accidental overloads, generated by the passage of vehicles and pedestrians; • hydrostatic permanent loads, determined by the water present inside the Palombaro Lungo and whose level, for safety reasons, is kept at a constant level by means of a drainage system and reaches a depth of 1.64 m at the deepest point.
The stability checks with these load conditions have returned trends of the maximum displacements which are slightly different for the two sections but, in both cases, largely compatible with hypogeum stability conditions (Fig. 12 on the left). Even the zeroing of the unbalanced forces and the maximum speeds (Fig. 12 on the right) Fig. 7 Step 1: trend of total displacements (left) and maximum speeds (right) in the rock mass before the quarrying of the hypogea indicate that the models reach stability conditions at the end of the simulation.
Furthermore, the application of loads in sections B-B′ and D-D′ does not determine redistribution of stresses compatible with the onset of yielding states (Fig. 13) and the values of the safety factors, calculated considering the Hoek-Brown resistance criterion, are largely precautionary (Fig. 14).

Conclusions
The paper illustrates, based on results from the Palombaro Lungo case study in Matera, how the stability of this excavated cave can be reliably studied using a 2D approach, according to cross sections perpendicular to the main length direction. This type of approach does not take into account the contribution to stability of the three-dimensional effects Fig. 8 Step 2: trend of total displacements (left) and maximum speeds (right) in the rock mass after the quarrying of the hypogea Fig. 9 Step 2: yield points in the rock mass after the excavation of the cavities Fig. 10 Step 2: trend of total displacements (left) and maximum speeds (right) in the rock mass after excavating the hypogea and applying cocciopesto coating of the vault type but, in the authors' opinion, is conservative enough to be considered a reliable approach for this type of problem.
The stress-strain state analysis of the hypogeum was based on two main representative cross sections perpendicular to the main direction of the hypogeum (B-B′ and D-D′). The stability analyses were based on a detailed numerical model analysis of the site and provided stress-strain results in good agreement with the really stability conditions, even when anthropogenic loads due to plaza use and of water inside the reservoir are applied. The maximum displacements were obtained in the vault covering the reservoir and are on the order of a tenth of a millimeter in section B-B' and of the millimeter in section D-D′, where the thickness of the vault cover is equal to the minimum evaluated value of 0.41 m. It has been also developed an analysis according to the longitudinal section of the hypogeum A-A′, that which to highlight conditions of instability. This section is not relevant from a static point of view and the results confirms that a longitudinal two-dimensional analysis is absolutely not representative of the stability conditions of a such an underground environment.
The in situ geo-mechanical surveys of the Palombaro Lungo made it possible to detect a fracture on the vault of the hypogeum which develops approximately parallel to the longitudinal section A-A′, without however, evidence of Fig. 11 Step 3: diagram of the loads imposed on the model of the cross sections B-B′ and D-D': permanent (square paving) and accidental (people and vehicles) overloads acting on Piazza Vittorio Veneto (in brown); permanent (water) inside the Palombaro Lungo (in light blue)

Fig. 12
Step 3: total displacements (left) and maximum speeds (right) in the rock mass after the application of the water load inside the Palombaro Lungo and the loads acting on the square above Page 15 of 16 302 Fig. 13 Step 3: yielding points in the rock mass after the application of loads Fig. 14 Step 3: safety factors in the rock mass, according Hoek-Brown strength criterion