Revision of archive recovery tests using analytical and numerical methods on thermal water wells in sandstone and fractured carbonate aquifers in the vicinity of Budapest, Hungary

This study offers a reinterpretation of archive aquifer tests, predominantly on the basis of recovery data, from an original datasheet of thermal water wells located in carbonate and sandstone aquifer units in the vicinity of Budapest, Hungary. The study compares the hydraulic conductivity (K) and specific storage (Ss) values derived in the first instance from an aquifer test evaluation. This included an initial application of the classical analytical Cooper and Jacob method. Subsequently, the visual two-zone (VTZ) numerical method was applied, then third, a more complex model, namely, WT software. It was found that the simple analytical solution is not able to represent the field conditions accurately, while in the course of the application of the VTZ model, it proved possible to alter the various hydraulic parameters within reasonable limits to fit the field data. In the case of the VTZ model, the researcher is required to calculate the accuracy of the fitted model separately, while with the WT model, this is automatic, the software seeks out the best fit. In addition to VTZ parameters, the WT model can efficiently incorporate data on up to 500 model layers, water level, and pressure. The optimization of the parameters may be achieved by automatic calibration, improving the accuracy of the numerical results. Recovery tests for 12 wells were numerically simulated to obtain values for vertical and horizontal hydraulic conductivity and specific storage for Triassic and Eocene fractured carbonate and the Upper-Miocene-Pliocene granular sandstone aquifer units. When an analytical solution is applied, only average values could be obtained. The conclusion reached was that the results of the analytical solution can be improved by the use of numerical methods. These methods are able to incorporate basic information on well design, aquifer material and the hydrogeological environment in the course of the evaluation. The revision of the archive recovery data using numerical methods may assist in the quest for better data for numerical flow and transport simulations without the need to perform new tests. In addition, the methods employed here can explain cases in which the original analytical interpretations proved unable to yield reliable data and predictions.


Introduction and goals
An important step in numerical simulations for hydrogeological purposes is the provision of adequate input of parameter values for any given geological unit of reference. Depending on the objectives, such values are usually taken from one of two sources: first, progressive drill-and-test borehole (slug test) characterization made in the course of their construction, providing information on discrete hydraulic properties plotted against depth (Spane and Newcomer 2009); or second, bulk values from specific literature on the geological unit, or from surface geophysical and/or logging measurements, as well as parameters based on laboratory studies. The characteristic values of in situ formations are more commonly obtained from an evaluation of an aquifer test, carried out in wells in the region of interest. During their evaluation, the best estimated values must be detected by integrating the hydrogeological functioning of the setting. A further constraint is that during the evaluation the fact that the parameters are scale-dependent, that is to say, the results are dependent on the duration of the aquifer test and consequently the involved rock volume, has to be taken into account.
With regard to aquifer tests, it should be noted that the basic drawdown/time (s-t) response for a granular aquifer has been solved analytically by various researchers, (Theis 1935;Boulton 1954;Hantush 1956), for confined, water-table, and leaky semi-confined conditions, respectively. A contrasting challenge is presented for fractured test interpretation at the aquifer scale (Francis 2010). There has been a certain level of application of classical analytical methods made for granular porous media to fractured media. On the other hand, most of the equations to be solved for fractured media require a knowledge of the special features of the fractures; additionally, for blocks formed by fractures, data on their distribution and porosity are required. This, in turn, implies that these data are usually not available on most wells constructed for groundwater abstraction, and especially in the case of archive datasets. However, there are various advantages to the available methodology for double porosity media, such as the possibility of computing the anisotropy (K v and K h ) and the storage coefficient of the matrix (S m ) of blocks, as well as that of the fissures (S f ). A solution to flow in vertical or horizontal fractures, or where well-loss and heterogeneity are included, may be found in Gringarten (1982). Moreover, the real producing thickness of a lithological unit needs to be identified if a more realistic K value for the fractured media is sought. There are several methods to obtain K, η and S from step drawdown tests (Birsoy and Summers 1980;Clark 1977;Eden and Hazel 1973). Important contributions have been published on how to cope with the limitations of the Theis (1935) method: a well of finite diameter in which the stored water present in the well is incorporated in the test (Papadopulos and Cooper 1967); well losses, as initially included by Clarke (1977); a case in which a well might be subject to variable discharge yield, as solved by Birsoy and Summers (1980); the solution to low hydraulic and constant head boundary found by Ferris et al. (1962); the question of the presence of vertical and horizontal variation of the hydraulic conductivity and storage coefficient, as solved by Neuman (1975); leakage input to the drawdown cone and the partial penetration case successfully examined by Hantush (1956Hantush ( , 1961; the delay yield response solved by Boulton (1954), Boulton and Streltsova (1978), which presented a solution for an aquifer test carried out in a fractured medium. Certainly, neither the identification of the various flow components, nor the data related to the physical media affecting the flow response in an abstraction well can be involved in an aquifer test data analysis carried out using analytical solutions, and especially in the case of a double porosity medium.
The main goal of the application of numerical modelling codes such as VTZ (Rushton and Redshaw 1979;Rathod and Rushton 1991) and WT (Székely 1992(Székely , 2006, respectively, to archive datasets previously handled using analytical methods, is the more accurate representation of the field conditions and flow in double porosity media. The archive data series of aquifer tests, and especially recovery tests, were carried out in the course of the drilling of the wells bored in carbonate and sandstone aquifers in the testing area. In the vicinity of Budapest, recovery data on deep thermal wells were gathered for this analysis. The aim of the analysis was to use different methods and compare results, thus hopefully gaining better understanding of the differences between analytical and numerical methods in the definition of aquifer properties. Such an evaluation implied a search for reliable values for hydraulic conductivity (K) and specific storage (S s ) for sandstone and carbonate aquifers in the study area. Over the timespan of the study, the question of how the heterogeneity and anisotropy of the fractured carbonate formations can be considered during aquifer tests was also addressed. In the case of the wells screening sandstone formation, it was desirable for the test to answer questions raised by the model concerning the hydraulic response to the numerous strata. The further significance of the present study lies in the demonstration of the process of obtaining best-fit values using numerical solutions on data provided by an analytical evaluation. This should be the result of an understanding of groundwater flow into the abstraction well. It should be noted that this appears to be key to the interpretation involved in gathering and reviewing further available aquifer test data. The study area and the dataset

Hydrogeological conditions
The study area is located in the vicinity of Budapest, central Hungary. Geologically, the Pre-Neogene basement of the study area is considered part of the Transdanubian range (TR) in the ALCAPA (Alps, Carpathians, and Pannonian Basin) mega-unit of Hungary (Csontos and Vörös 2004). The aquifer system of the area is characterized by slightly metamorphosed Variscian formations overlain by nonmetamorphosed Alpine sequences dating from the Middle and Late Triassic to Early Jurassic (Wein 1977;Haas et al. 2000). The study area is built up Triassic carbonates, Middle Triassic cherty dolomite and limestone and mostly from Upper Triassic dolomite and limestone, to a total thickness of several thousand metres. These formations are overlain by Upper Eocene limestones and an Upper Eocene/Early Oligocene calcareous marl, the Buda Marl. The sequence has been affected by normal and reverse faults and behaves as a fractured and slightly karstified aquifer system (Wein 1977;Haas 1988;Fodor et al. 1994). During the Middle Eocene-Early Oligocene, a generally deep, underfilled basin evolved here (Tari et al. 1993). In the Late Eocene, terrestrial siliciclastic rocks, shallow-water limestone, and deep-water marl, were deposited here. These are overlain by Oligocene anoxic deep marine shales (Báldi and Báldi-Beke 1985). In the Late Oligocene to Early Miocene, the basin filled up with shallow marine-to-continental siliciclastics. From the Early Miocene, the formation of the Pannonian basin commenced in an area characterized by marine siliciclastic formations, shallow marine limestones, and volcanics and volcanoclastics. The basin was structurally reactivated beginning in the Late Miocene-Pliocene (Horváth and Royden 1981;Horváth and Cloetingh 1996;Bada et al. 2007). In this area, the isolated Lake Pannon evolved, and was gradually filled up with sediments derived from the uplifting Alps and Carpathians. The deep-to-shallow lacustrine-to-alluvial sedimentary succession is represented by deep water marls, turbidite, sandstones, and shales overlain by an alluvial sequence (Juhász et al. 2007;Sztanó et al. 2013Sztanó et al. , 2016. The Quaternary sediments in the area are represented by alluvial gravel, sand, clay, and aeolian sand and loess (Kőrössy 2004). The inversion of the Pannonian basin caused the uplift of the Transdanubian Range and the gradual evolution of topography-driven fluid flow systems and related karstification (Wein 1977;Kele et al. 2011;Leél-Őssy 1995;Erőss 2010;Havril et al. 2016). The whole area is characterized by elevated heat flux (~ 100 mW/m 2 ), as compared to the surrounding regions in the Pannonian basin (Lenkey et al. 2002).
The study area is divided by the River Danube into western and eastern parts (Fig. 1a, b). From the hydrogeological point of view, the Danube represents the boundary between the semi-unconfined (Buda Hills), and the confined carbonate aquifer systems (Gödöllő Hills and Pest side) . The dominantly Triassic carbonate formations partially outcrop and are partially covered by an Oligocene aquitard and Early-Middle Miocene aquifer-aquitard formations to the west of the Danube. The Triassic formations, however, are completely covered by up to 3000 m thick of these latter formations to the east of the river.
The covering formations are east of the Danube the Upper-Miocene-Pliocene and Quaternary formations behave as aquifer units (Fig. 1b). The Triassic and Eocene carbonate formations of the area, made up of dolomite, limestone and calcareous marl, are characterized by double porosity. In contrast, the siliciclastic Upper Miocene-Pliocene sandstone formations behave as a granular aquifer and display primary porosity.
The natural discharge of the system takes place near the Danube in the form of lukewarm and thermal springs (Papp 1942;Alföldi et al. 1968;Lorberer 2002;. However, the many thermal wells operating in the area since the mid-nineteenth century have influenced the natural discharge (Déri-Takács et al. 2015). The total dissolved solids (TDS) of water in the examined carbonate and sandstone aquifers is less than 10,000 mg/l .

Database
Aquifer test data from the archives of 46 tested wells, all of which yield thermal water in the research area were examined. All tested wells reported, exclusively, recovery test data series collected after a relatively short abstraction period, based on regular practice during the preformance of tests for new wells. These wells were screened to different aquifers described as Upper Miocene-Pliocene sandstone, as well as in Middle Triassic dolomite, Upper Triassic dolomite, limestone and Upper Eocene-Lower Oligocene calcareous marl, respectively. The tested wells were unevenly distributed throughout the research area (Fig. 1a).
Most of the data series could not be evaluated due to the lack of appropriate documentation of static water level or hydraulic head, production yield, production time, etc. Moreover, with regard to wells sunk into carbonate rocks, it was a common occurrence that when the drilling of the well reached an open fracture, its water level recovery was almost immediate. Therefore, after disregarding inadequate data series, 12 wells were identified as suitable for processing ( Fig. 1a; Table 1). The density contrast is negligible in thermal wells due to the low TDS content of the examined aquifers. Consequently, no density contrast was taken into consideration during the carrying out of the aquifer test.

Aquifer test methodology
On the basis of the geological setting, the analysis of available aquifer tests from wells constructed in these two types of formations present in the study area requires different evaluation methods and distinct approaches. For this project, the classical analytical method developed by Cooper and Jacob (1946) was used for recovery test evaluation and its results were compared to those deriving from the application of two numerical models, the Visual Two Zone (Rathod and Rushton 1991) interface designed by Noel Laloth (Faculty of Engineering, UNAM), and the WT model (Székely 1992(Székely , 2006(Székely , 2015. Cooper and Jacob (1946) proposed a straight-line solution to the estimation of the hydraulic properties of a confined aquifer unit yielding water to a well through the analysis of time-drawdown data collected in an observation well during a constant-rate abstraction test.

Analytical, Cooper and Jacob
where d is thichness of the layer [L], K is hydraulic conductivity [L/T], Q is the abstraction rate [L 3 /T] and Δs is the slope of the fitted line (change in drawdown per log cycle time).
This application provides a simplified method for Theis (1935) solution with its basic assumptions. These assumptions are difficult to meet in the actual field conditions. Therefore, individual corrections have been devised by different authors to make field data to comply with the requirements for the application of the analytical solution.
In the application of the Theis solution by Cooper and Jacob (1946), the Aquitest Software Package by Waterloo Hydrogeologic Inc. (Roehrich 1996), was applied as a "recovery test" evaluation procedure. Available aquifer test data were initially reviewed through the application of this method. On the one hand, there was a lack of data from observation wells. Furthermore, their construction characteristics (i.e., depth, screen interval, well diameter, etc.) and well losses could not be integrated into the calculations. Additionally, limitations were found on defining how the different nature of the lithology crossed by the well (i.e., fractured, porous, etc.) could assist in defining the groundwater flow conditions in the tested well. Originally, this method was developed to investigate the hydraulic properties of a one-layered aquifer unit bounded above and below by confining units. Finally, this method avoids the question of the nature of the groundwater flow (fissured, granular) and the handling of multilayered flow to the abstraction well.

Numerical visual two zone (VTZ)
The Two Zone numerical model for aquifer test data interpretation was developed for two-layered aquifer units (Rushton and Redshaw, 1979;Rathod and Rushton, 1991). This program was further developed by N. Hernandez Laloth in 2008, via the creation of a piece of Windows software with a user-friendly graphical interface. This numerical model solves the flow equation radially centred on the abstraction well, creating a logarithmic mesh of cells whose size increases with distance from the abstraction well. The calculated water-flow to the well is based on the extended version of the Thiem equation (Thiem 1906) for radial flow.
The equations of the numerical model can be found in the original paper by Rathod and Rushton (1991). The standard model consists of 4 layers: two aquifer units (upper and lower), an aquitard layer in between, and a top layer, which may technically be described as a leaky aquitard. The structure of the basic model and the parameters of the different layers are shown in Fig. 2. One of the most important improvements in this theoretical approach is that it considers the horizontal flow in two aquifer units, and also the resulting vertical components of flow due to abstraction, as well as leakage from the top layer. These are significant improvements because in most cases it is found that in the analysed hydrogeological setting, the vertical components of flow are not negligible (Rathod and Rushton 1991). Furthermore, the model can be used to understand and inspect water quality changes in those cases when groundwater abstraction causes saline-water upwelling (Carrillo-Rivera et al. 1996).
The VTZ model can only handle horizontal strata. For systems with double porosity, the model can be applied to a limited degree, since it was developed for systems made up of granular aquifers. However, it has been successfully used in fractured rocks, on the basis of an approach that represents the average flow through the system as the sum of the water released from fractures and matrix, respectively (Rathod and Rushton 1991). The average hydraulic conductivity (K) of a fractured double porosity aquifer cell unit can be estimated as indicated in Fig. 3, where Q is the flow, i is the hydraulic gradient, A is the cross-section of flow, and K is the hydraulic conductivity.

Analytical and numerical, WT Software
The WT software is designed to simulate the hydraulic response caused by discharge/recharge/recovery, constant head, intermittent free flow, and slug as well as packer tests considering linear cross flow between the model layers or diffusive cross flow through aquitards (Székely 2006(Székely , 2015. The description of the mathematical model and the details of the original software can be found in Székely (1992Székely ( , 2015. The WT is able to use an analytical technique based on the numerical Laplace inversion (Hemker 1999) and also a numerical technique which is a finite difference method improved by Székely and Galsa (2006). The wellbore simulator of the WT software includes the evaluation of the following effects: laminar and/or turbulent well-loss with depth and/or time variant parameters of the well screen; turbulent axial friction loss; a non-uniform static level in the screened model layers; induced flow-controlled drawdown in observation wells; and wellbore storage. Simulation options of a thin skin of infinitesimal extent, a ring-shaped skin zone and a radially two-zone formation model are also available in the software (Székely 2015).
The WT software allows for the incorporation of geological features, permitting the design of a multi-layered hydrostratigraphic model and the selection of the most appropriate solution option. Automatic parameter calibration may be carried out by fitting the simulated drawdown, recovery, flow metering and overflow data to the measured ones. The software has been successfully applied to the evaluation of field tests conducted in sedimentary aquifers (Székely 1992(Székely , 2012Mukhopadhyay et al. 1994), fissured granite (Székely 2013;Székely and Galsa 2006) and in volcanic aquifers , as well. Table 2 shows the comparison of the VTZ and the WT models.

Results of the study
As may be seen, the screened intervals of the aquifers for the examined thermal wells of the region may be classified into carbonate and sandstone groups (Table 1). In the present paper, of the evaluations of aquifer tests performed on these groups, one example from each group will be presented.

VTZ model
The depth of the well is 20 m, it was deepened in the Main Dolomite (T3 Main Dolomite Fm.) aquifer (Figs. 1a, b, 4).
The upper part of the dolomite (5.8 m) is intact which followed by 5.7 m fractured section. After this point, intact dolomite is found down to the bottom of the well. There is an absence of a low hydraulic conductivity layer covering this aquifer unit. The static water level was at a depth of 4.5 m-consequently it is an unconfined aquifer. Since the Main Dolomite aquifer below the water level is open throughout to the well, a double-layered model structure was created that is capable of analysed using the VTZ model ( Fig. 4; Table 3).
The upper-7.0 m thick-aquifer unit was defined from the static water level down to the clogged section of the well, so the aquifer unit partially includes intact and fractured dolomite, too. The lower aquifer unit was defined as a 50 m thick intact dolomite. The initial parameters of the numerical simulations derived from archive recovery test data can be seen in Table 3. The horizontal hydraulic conductivity of the upper aquifer unit was three orders of magnitude higher than that of the lower unit. However, in the fractured section, anisotropy only occurred to a smaller degree (K h =2K v ). Contrarily, in the lower aquifer unit it was found that K v = 3K h , which as expected, suggests the presence of vertical fracturing in the dolomite block.
The value of specific yield (S y ) for the upper aquifer came to be 0.1, which is one order of magnitude larger than that obtained for typical dolomite, but a value of this size is still admissible. The lower aquifer, S s , had an acceptable value of 6 × 10 -6 m −1 . These specific yield values must be treated with a great deal of caution. They depend on the length of the test, and there were no continuous water level  5).

WT model
In the WT model, the upper aquifer was represented by two layers, as opposed to the one layer in the VTZ, and the parameters were set accordingly. The uppermost 1.3 m thick layer and the lowest, third layer (50 m thick), were considered as intact dolomite without fractures, and with similar parameters of hydraulic conductivity values (Table 3).
Through an independent analysis of data, the modelled parameters in this model were found to coincide with those obtained in the VTZ model. A value of 0.1 was obtained for S y of the upper layer using the WT model. In pursuit of a more refined result, zones characterized by well-loss were also defined for this well. The hydraulic conductivity was found to increase in the surroundings of the well; the best match of data was obtained when the skin zone radius was 2.0 m, in which case K was 2 orders of magnitude greater close to the well (Table 3). When the WT model was run with these parameters, the adjustment of the model was excellent; the modelled values showed an average variance of approximately 1 cm from the measured values (Fig. 5). In the abstraction period, the shape of the drawdown curve differs in the two models.

VTZ model
The thermal water wells screened in the Upper Miocene-Pliocene sandstone aquifer produce water from more than two strata. Since the VTZ software can only deal with a maximum of two aquifers units at once, numerous screened  (Figs. 1a, 6). In the surroundings of the screened section, fine and medium sand, silt, silt loam and silty clay alternate. In Well 8, three screens can be found. The upper screens are open to fine and medium sand, while the lower two screens to medium sand. Above and between the screened sections, primarily silty clay is reported in the wellbore log.
An archive data series was available for Well 8, Pécel, which was observed after a 3-day abstraction. In the course of the aquifer test, wellbore flow measurements were made and the log revealed that the uppermost screen supplied 42%, while the one below supplied 58% of the obtained water during the test. The lowermost screen supplied nothing at all. It can be estimated that in the course of the 3-day test, the 15 m thick silty clay could have behaved as a hydraulically active layer between the two screened sand layers.
Accordingly, a confined 4-strata model was defined, a screened section was included in both the upper and lower aquifers (Fig. 6). Between the two, a layer with low conductivity was placed. The estimated well-loss value had a factor of 4. The calibrated model parameters are presented in Table 4. The hydraulic conductivity of the upper aquifer had a value of ~ 10 -6 m/s, while for the lower it was 10 -5 m/s. These results are consistent with the lithology, i.e. with the upper layer screening finer sediments than the lower. The model response to these parameters showed a high degree of agreement with measured values (Fig. 7).

WT model
The wellbore flow measurements were taken during the test of Well 8, Pécel. These data were used to construct a model based also on the capabilities of the WT software. Segments of the strata, comprising a pair of aquifer and aquitard layers, were introduced into the model. An adequate mapping of the surrounding screens was defined, including 5 segments for the well. Figure 6 shows the segments within the aquifer and aquitard sections, and also which layers matched the drilled strata. The first and fifth segments are truncated. The aquitard part of the first segment and the aquifer part of the fifth segment was 0 m thick. The aquifer section of the second and third segments was screened. Based on the wellbore flow measurements, the third screen in the fourth segment was not considered. In the segment-based model, the aquifer part had only horizontal hydraulic conductivity (K haf ), while the aquitard part had only vertical (K vat ), in accordance with the preferred flow orientation. The values used to define the hydraulic parameters can be seen in Table 4. A similar value for K haf , 10 -5 m/s, was chosen for the first two segments screening a sandstone formation. The K haf value in the third segment was smaller, to represent the effect of the thin silt layer in the screened strata. Since there are no available data on the aquitard formations, a K v value of 5 × 10 -10 m/s and a S s value of 5 × 10 -7 1/m were chosen after the calibration of the model using the wellbore flow measurement results.
The laminar skin zone was set on the two screened sections with the additional use of the above results. Table 4 shows the first screen with a poor value, while the second screen implies an improved value in the well surroundings. The running of the model with these parameters displayed a strong correlation with the measured values (Fig. 7). The mean absolute deviation is under 19 cm, which, compared to the maximum drawdown value (30 m), may be considered acceptable.

Interpretation and discussion
In this study, a completed revision of previously evaluated historical recovery tests in carbonate and sandstone thermal aquifer units has been carried out. The goal of the revision was to compare the application of the analytical solution of the Theis equation using the Cooper and Jacob method with further numerical solutions. In addition, the results may represent the potential for further archive data re-evaluation based on straightforward modern analytical and numerical solutions.
The first analytical method, the Cooper and Jacob, produces an average K value which is meant to represent the screened strata. It is a subjective straight-line solution into which the parameters of well performance and those of the related dimensions of the aquifer units are not incorporated.
On the other hand, the parameters of the model built on each tested well by the VTZ were adjusted to well structure and lithology data. Consequently, the modelled results match the measured ones. One disadvantage here is that the accuracy of the solution using the VTZ is not measured by the model. In this regard, the WT model does provide numerical information on the correctness of the results. As for providing an appropriate initial estimate, the calibration is automatic and seeks the best match during the calibration. Moreover, with the VTZ and WT models, the K h and K v of the productive layers may be obtained, as well as other factors affecting the aquifer test. Since the archive data set was based on short-term aquifer tests, the K h factor was adoptable with a greater degree of certainty compared to K v from among the tested parameters. Therefore, such values may be viewed as the primary basis of comparison. In Table 5, these columns have been highlighted in bold.
The use of VTZ for the interpretation of drawdown data for carbonate aquifers seems to imply a response similar to a Theis solution, in the absence of additional drawdown data. However, in the recovery test both models provide comparably accurate representations of reality. Values obtained for hydraulic conductivity are characteristic of the Main Dolomite aquifer for both the fractured and solid material as an average value for hydraulic conductivity. In the cases of Well 3, Budapest, Well 6, Budaörs, and Well 12, Nagykáta, an underestimation of up to one order of magnitude of the original results (Theis solution) above the VTZ model can be observed. In the first two cases, this deviation presumably arises from the peculiarities of the double porosity of the carbonate rocks. In the wells screening the Upper Miocene Pliocene sandstone, although it could be argued that numerous aspects were simplified in the VTZ model, the conditions of the Theis model were only met to a limited degree. For those wells in which a WT model was used, e.g. Well 8 Pécel, a VTZ model had been constructed for the purposes of comparison, and in most cases the values were close to be identical. That is to say, the better options of the WT software could not offer more reliable results. Table 6 shows the mean vertical and horizontal hydraulic conductivity values for each well.
At the same time, the results summarized in Table 6 can only be evaluated exclusively in terms of local hydraulic conductivity, since they are derived from a short-term evaluation, maximum 120-min aquifer tests (except Well 8 and 9). It must be noted that the parameters obtained from the short-term well tests may be, on average, 1 to 2 orders of magnitude smaller or greater than the parameters of the same formation on the regional scale and in a long-term test. In the case of sandstone formations, the basin-scale hydraulic conductivity values are smaller than those for local scale hydraulic conductivity (Garven 1989). In the case of double porosity formations (Main Dolomite (T 3 ), Middle Triassic-Carnian dolomite (T 2-3 ), Buda Marl (E 3 -Ol 1 ), on the basis of work by Király (1975) basin-scale hydraulic Table 5 Comparison of hydraulic parameters obtained from original analytical data and from the numerical re-evaluation of the aquifer tests on wells of the study area using the VTZ and WT

Summary and conclusions
This study represents a re-interpretation of recovery aquifer tests from archive historical data of thermal water wells belonging to carbonate and sandstone aquifer units in the vicinity of Budapest. The resulting K and S s obtained values were derived using one analytical and two independent numerical methods. The aquifer test evaluation began with an initial application of the classical analytical Cooper and Jacob (1946) method, followed by the Visual Two-Zone (VTZ) numerical method (Rathod and Rushton 1991) and a more complex model, the WT software (Székely 1992(Székely , 2015. In the course of the application of the VTZ model, the hydraulic parameters for the well were changed within reasonable limits to achieve the best fit with the field data s-t reports. In the case of the VTZ model, the researcher is required to calculate the accuracy of the fitted model separately, while with the WT model, this takes place automatically, the software seeks out the best fit. The later model can also incorporate some 500 model layers, water level, and pressure data. The optimization of the parameters could be achieved by automatic calibration improving the exactness of the numerical results. Finally, 12 recovery tests were numerically simulated to obtain hydraulic conductivity and specific storage values for the Main Dolomite, as well as for the Buda Marl, and the Upper-Miocene-Pliocene sandstone aquifer units. When applying the Cooper and Jacob (1946) model its original assumptions were difficult to meet. The application of the analytical model yields only average hydraulic conductivity values, while from the other two numerical methods horizontal and vertical hydraulic conductivity values could also be derived. Due to the short length of the aquifer tests K h values may well be more reliable than K v values; the former may be compared with existing data from archives. Values of K h obtained from numerical simulations in this study display a similar range to that which might be expected in analytical solutions. However, in Wells 3, 6 and 12, a range difference was revealed in the derived values between the Cooper-Jacob and the VTZ models. For the first two wells this can be the consequence of the double-porosity of the carbonate aquifer. However, in case of the Upper Miocene/Pliocene sandstone aquifer, many simplifications were made to use the VTZ model, and these may have caused the observed difference. On the whole, the comparison of the WT and VTZ models provided K h , K v and S s values within more or less the same range. Therefore, both models may be used to reach a better understanding of the examined groundwater movement between the abstraction well and the geological media. Due to the short-term nature of the evaluated aquifer tests the better options of the WT software were not able to offer more reliable results, because in most cases during the test some useful parameters (wellbore temperature, salinity and flow) were not measured. Derived data reflect local hydraulic parameters limited by the cone of depression in the immediate vicinity of the wells tested. Given that studies carried out on a regional scale (Király 1975;Hartmann et al. 2014) have proposed that such parameters for carbonate rocks could be two orders of magnitude larger than those obtained from aquifer test analysis, a set of longer abstraction aquifer tests would be desirable in the interests of addressing this question. funding from the European Union's Horizon 2020 research and innovation program under grant agreement No 810980.
Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit http://creat iveco mmons .org/licen ses/by/4.0/.