Modelling Water Transferences from the Almádena Odeáxere Aquifer System to Vale do Barão and Bensafrim Creeks (Algarve – South Portugal) J. MARTINS & J. P. MONTEIRO Algarve University Geo Systems Centre UALG/ CVRM, Marine and Environmental Sciences Faculty, Campus de Gambelas 8005-139 Faro, Portugal e-mail: [email protected]

Abstract The Almádena-Odeáxere aquifer system (Algarve, Portugal) is crossed by temporary watercourses which are mostly influent. In order to understand existent groundwater-surface water interactions, a numeric finiteelement groundwater flow model was used to analyse transferences at the short effluent reaches of the two most important creeks, where surface water and groundwater are hydraulically connected. Hydrographs simulating these transferences, considering homogeneous transmissivity and a heterogeneous distribution of that parameter, obtained through inverse calibration, were compared and analysed. Overall results allow the discussion of the influence of the heterogeneity in the balance of transferences toward the two discharge areas of the aquifer. Key words groundwater-surface water interactions; inverse modelling; Algarve

INTRODUCTION Tourism and agriculture activities in Algarve (the southernmost region of Portugal) rely heavily on the use of the available water resources, particularly during dry periods, posing a threat to the public water supply. The multimunicipal public water supply system became recently, entirely dependent on surface water provided by large reservoirs. During the occurrence of an extreme drought in 2004-2005, serious consequences of this single-source strategy were felt when several reservoirs reached their exploitation limit and could not satisfy the water demand. In this period an important part of the urban supply was supported again by groundwater use. In order to guarantee public water supply in the long-term, a strategy based on an integrated water resources management is therefore needed. However, for an adequate integration of groundwater resources into the water supply system, the existent interactions between groundwater and surface water occurring at Algarve aquifers must be better understood, since scenarios of overexploitation of the subterranean resources may pose additional risk on natural systems which depend on groundwater reserves. With the aim to understand these processes, a numeric finite-element groundwater flow model was applied to the study of groundwater-surface water interactions occurring at the short effluent reaches of the two most important creeks on the Almádena-Odeáxere aquifer system (AO) area: “ribeira de Bensafrim” and “ribeira de Vale do Barão”, Algarve region (Fig. 1). On these reaches, surface water and groundwater are hydraulically connected and an effort was made to understand the influence of the distribution of transmissivity (T) on the water balance of transferences,

268

by comparing simulated hydrographs with gauging data. Simulations were based on a calibrated version of the AO model, where T was estimated by inverse modelling for zones where the behaviour of piezometers allows a reasonable fitting of field data using a single value of T. A similar approach has already been used by Monteiro et al. (2007b) on the improvement of data collection and definition of boundary conditions towards a better conceptualisation of river-aquifer relations. Similar numeric models were also used by Monteiro et al. (2006) and Monteiro et al. (2007a) to study Algarve aquifer systems with different objectives.

CONCEPTUAL FLOW MODEL AND INVERSE CALIBRATION OF THE NUMERIC FLOW MODEL The AO has an area of 63,5 km2 across Odeáxere (East) and Almádena (West). It develops in carbonate Lias-Dogger lithologies, limestone, dolomitic limestone and dolomite rock, which show, in some places, a well developed karst with a thickness in the order of 750 meters (Reis, 1993; Almeida et al., 2000). As can be seen from Fig. 1, regional water flow occurs predominantly from NE (recharge area) to SW (discharge occurs at “ribeira de Vale do Barão” creek). An N (recharge area) to S (discharge occurs at “ribeira de Bensafrim” creek) flow can also be identified on the right section of AO. A hydrometric monitoring network collects data throughout the aquifer area and is composed of 7 gauging stations, located at both effluent (stations 6, 2, 1 and 3) and influent (stations 7, 5 and 4) stream reaches. On the last two decades, hydraulic parameters have been obtained for AO from

4 Bensafrim 7 Odeáxere Sargaçal B.S. João 5

Mte. Judeu

2 1 Chinicato

3 Portelas

Meia Praia

B.S. Miguel Almádena

Espiche Lagos

6 V. Boi Luz Burgau

0

20000

40000

60000

80000 m

Fig. 1 – Stream network (dotted grey lines) and effluent reaches of the two most important creeks “ribeira de Bensafrim and “ribeira de Vale do Barão” (thick black dots) on the AO. The 7 gauges composing the hydrometric network are represented as black squares. Dashed lines represent the equipotential surface obtained from measurements at boreholes of the hydraulic head monitoring network (black crosses).

269

pumping tests in individual boreholes (Reis, 1993; Almeida et al., 2000); however these methods cannot provide the necessary data to carry out realistic representations of regional aquifers. In order to overcome these limitations, simulations were performed, using a synthetic bi-dimensional numerical representation of the AO, where T was estimated by inverse modelling. These results allowed a very important improvement on the reliability of the simulation of the regional flow pattern observed at the AO, namely in a best characterization of the spatial distribution of hydraulic head. The physical principles on the basis of the simulation of the hydraulic behaviour of the aquifer system are expressed by the equation (1):

S

uuuuur ∂h + div −[T ]gradh = Q ∂t

(

)

(1)

Where: T [L 2 T -1] is transmissivity, h is the hydraulic head [L], Q is a volumetric flux per volume unit [L 3T -1L -3], representing sources and/or sinks and S is the storage coefficient [-]. In steady state conditions the variables are time-independent. In this case the equation (1) is reduced to equation (2): uuuuur div −[T ]gradh = Q

(

)

(2)

The direct solution of equations (1) and (2) was implemented using a standard finite element model based on the Galerkin method of weighted residuals. The model has a mesh with 7494 nodes and 14533 triangular finite elements for which the defined conceptual flow model was translated. The use of finite element models is currently a standard approach for solving problems in hydrogeology and is described in textbooks such as Huyakorn & Pinder (1983), Kinzelbach (1986), Wang & Anderson (1982) and Bear & Verruijt (1987). The calibration of T was performed by inverse modelling, using the GaussMarquardt-Levenberg method, implemented in the nonlinear parameter estimation software PEST (Doherty, 2002). This calibration approach has already been thoroughly described and implemented by authors such as Carrera et al. (2005) and Poeter & Hill (1997) and its application has widespread across academic circles. The present work is therefore also a contribution to put into effect a model-supported water resources management using this methodology as common practice for calibration. Individual values of T were defined for 16 zones where the behaviour of piezometers allows a reasonable fitting of field data using a single value for this parameter. The optimisation of the results was based in numerous model runs performed by inverse modelling, whilst several variants of the model were tested to search the best possible reproduction of the aquifer equipotential surface, accommodating values of water balance. The process of calibration of the AO numeric flow model is thoroughly described in Martins (2007) and Martins & Monteiro (Accepted). Calibration results (Fig. 2) ranged between 86 m2 day-1 and 8158 m2 day-1, throughout the defined zones. The calibration revealed a good fit between simulated and measured head values, the correlation coefficient, R, value was 0,9967, well above the minimum acceptable limit, 0,9, and the sum-of-squared weighted residuals between model outcomes and corresponding field data (i.e. objective function, Φ) was 4,56 m.

270

271

Simulated hydraulic head, in meters above sea level

593/5 23 Fict22 22 Fict20 21 Fict21 20 19 18 Fict18 17 Fict19 Fict16 16 15 14 Fict17 13 Fict14 12 11 Fict15 10 9 AO-02 Fict12 8 Fict11 602/242 AO-06 7 603/38 AO-08 Fict10 Fict13 6 Fict9 602/187 5 AO-10 AO-14 AO-16 4 Fict7 AO-01 Fict6 Fict8 3 Fict4 Fict5 2 Fict1Fict3 1 0 Fict2 0 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23

Observation points Fict2 Fict3 Fict5 Fict1 Fict4 Fict6 Fict7 Fict8 AO-16 Fict9 Fict10 AO-08 AO-06 Fict11 AO-02 602/242 Fict12 Fict14 Fict16 AO-14 Fict18 AO-01 Fict20 602/187 AO-10 Fict17 593/5 Fict15 Fict21 Fict22 Fict19 Fict13 603/38

Hydraulic head computed from measurements, in meters above sea level

Measured (m) 0.50 1.50 2.50 0.50 1.50 2.50 3.50 3.50 4.42 4.40 4.70 5.14 4.69 5.10 7.23 6.98 8.50 11.50 14.50 5.27 18.00 3.79 21.50 5.33 3.95 14.50 22.70 11.50 21.50 23.00 18.00 8.50 5.00

Calculated (m) 0.66 1.37 2.53 0.97 2.03 3.02 3.84 3.48 4.16 4.96 5.29 5.64 6.01 6.26 7.64 6.82 7.57 12.04 15.47 5.09 17.52 4.11 21.66 5.43 4.42 14.36 22.99 10.81 21.47 22.91 17.07 6.28 5.83

| Residuals | (m) 0.16 0.13 0.03 0.47 0.53 0.52 0.34 0.02 0.26 0.56 0.59 0.50 1.32 1.16 0.41 0.16 0.93 0.54 0.97 0.18 0.48 0.32 0.16 0.10 0.47 0.14 0.29 0.69 0.03 0.09 0.93 2.22 0.83

Fig. 2 – Hydraulic heads computed from measurements versus simulated hydraulic heads (left and above) and residuals (absolute value of the difference between measured hydraulic head and model-generated hydraulic head) associated to each observation point (table). Spatial distribution of T values for each of the 16 predefined zones obtained using inverse calibration (left and below).

An interpolation of simulated hydraulic head values, obtained using the calibrated version of the model (T values were distributed according to the results shown on Fig. 2) is presented on Fig. 3. A comparison is made between this version

0

20000

40000

60000

80000 m

0

20000

40000

60000

80000 m

Fig. 3 Simulated hydraulic head contours for the non-calibrated (homogeneous T) version of the model (left) and simulated hydraulic head contours for the calibrated version (right) represented as solid lines, versus measured head values (dashed lines).

272

(Fig. 3, right) and a previous one (Fig. 3, left), where the homogeneous value of T = 3.92×102 m2 s-1 (an average of values obtained on pumping tests of individual boreholes), was used. It can be observed that the calibrated version is able to supply a more realistic spatial distribution of hydraulic head throughout the aquifer, since there is a fairly better fit between measured and simulated values in this case.

TRANSIENT SIMULATIONS

In order to perform transient simulations, recharge was calculated using average daily precipitation values, from January 2005 until November 2006, obtained at the Lagos (31E/01UC) meteorological station. These values consist of abrupt point inputs and were converted into smoother continuous recharge inputs before being introduced in the AO finite-element model in an attempt to trim down a potential source of numeric instability into the model (Fig. 4). Given an episode of precipitation (consisting of one or more days of significant precipitation), equation (3) is used to calculate the highest value of recharge, hrec, attained during the episode: hrec =

2× R × P n

(3)

Where: R is recharge rate, n is the number of days of the precipitation episode and P is the total precipitation occurring at the episode.

45

45 40

35

35

30

30

25

25

20

20

15

15

10

10

5

5

0

Recharge (mm)

Precipitation (mm)

Precipitation Recharge

Recharge Peak

40

0 0

2

4

6

Time (days)

0.4

Recharge (l.h-1.m2)

0.35 0.3 0.25 0.2 0.15 0.1

0

Jan-2005 Feb-2005 Mar-2005 Apr-2005 May-2005 Jun-2005 Jul-2005 Aug-2005 Sep-2005 Oct-2005 Nov-2005 Dec-2005 Jan-2006 Feb-2006 Mar-2006 Apr-2006 May-2006 Jun-2006 Jul-2006 Aug-2006 Sep-2006 Oct-2006 Nov-2006 Dec-2006

0.05

Fig. 4 Illustrated example of the transient recharge input method based on measured precipitation values of a typical precipitation episode (above) and calculated recharge values used on transient simulations (below).

273

The recharge rate, R = 0.41, was obtained from the conceptual recharge model, proposed by Vieira & Monteiro (2003) and Monteiro et al. (2003). It corresponds to a weighted average having into account the presence of sub-areas where recharge classes range between 5 % and 50 %, depending on the outcropping geologic structures. Values of transient recharge were distributed linearly increasing from the beginning of the episode until hrec and decreasing afterwards until reaching the last day of the episode. Several simulations were performed with varying storage coefficient, S, values, ranging between 0.05 and 0.15. This sensitivity analysis determined which value would be more suitable for model simulations to reproduce the varying pattern of the hydrographs obtained at the gauges positioned close to the effluent reaches of the two selected creeks (Fig. 6). The used S was 0.1 due to the need to maintain simulation parameters inside realistic physical boundaries.

ANALYSIS OF GROUNDWATER-SURFACE WATER INTERACTIONS

-0.5

-0.5

-0.4

-0.4

Q (m3.s-1)

-0.3

-0.2

-0.3

-0.2

-0.1

Nov-2006

Jul-2006

Jan-2006

Nov-2005

Jul-2005

Sep-2005

May-2005

Jan-2005

0

Mar-2005

Nov-2006

Jul-2006

Sep-2006

May-2006

Jan-2006

Mar-2006

Nov-2005

Jul-2005

Sep-2005

May-2005

Jan-2005

Mar-2005

0

Sep-2006

Bensafrim_hetero VBarão_hetero

Bensafrim_homo VBarão_homo

Mar-2006

-0.1

May-2006

Q (m3.s-1)

Balanced water fluxes, Q, ocurring from the aquifer to the “ribeira de Bensafrim” and “ribeira de Vale do Barão” creeks were simulated considering both homogeneous T and a heterogeneous distribution of that parameter, obtained through inverse calibration, similarly to the permanent simulations. Streams were represented as specified-head boundaries; model nodes where streams are located were simulated with an unchanging head. Head values were set at the stage of the stream. It implies that there is no head loss between the stream and the groundwater system and that the flow of groundwater into or from the stream will not affect the stage of the stream. The amount of water flowing between the stream and the aquifer depends upon the groundwater heads in the nodes that surround the specifiedhead boundary representing the stream. This representation was found appropriate since the selected reaches are expected to be well connected to the groundwater system

Fig. 5 Hydrographs of balanced water fluxes of the effluent reaches of “ribeira de Bensafrim” (black line) and “ribeira de Vale do Barão” (cyan line) creeks, between January 2005 and November 2006. On the left, a homogeneous distribution of transmissivity. On the right, a heterogeneous distribution, obtained through inverse calibration.

and the stream stage is not expected to change (Reilly, 2001). The resulting hydrographs can be observed on Fig. 5 (fluxes out of AO have negative values).

274

1.2 0.9

Gauge 6 Gauge 1

Q (m3 . s-1)

0.8 0.7 0.6 0.5 0.4 0.3 0.2

Dec-2006

Oct-2006

Aug-2006

Jun-2006

Apr-2006

Feb-2006

Dec-2005

Feb-2005

0

Dec-2004

0.1

Fig. 6 Hydrographs of the “ribeira de Bensafrim” (Gauge 1, cyan line) and “ribeira de Vale do Barão” (Gauge 6, black line) creeks. Hydrographs are truncated on the “x“ and “y” axis to keep the figure legible (Gauge 6 has no entries before October 2005 and Gauge one has only “0” values registered between March 2005 and October 2005.

Through the observation of measured values of streamflow discharge, Q, on Fig. 6., and considering that groundwater is the most important contribution for flow occurring on these creeks, it can be stated that simulated transferences (Fig. 5) from the aquifer to the stream at “ribeira de Bensafrim” are expected to be higher than at “ribeira de Vale do Barão”, given the relative importance of measured values at the corresponding gauges 1 and 6. However, calibrated simulation results suggested the opposite (Fig. 5, right). In fact, simulation results achieved considering homogeneous T on the flow domain (Fig. 5, left) are more coherent with the balance of proportions of Q, existing between these creeks. The objective of simulating a realistic spatial distribution of hidraulic head has been entirely accomplished through the use of inverse calibration. Nevertheless, the subsequent transient simulations have revealed that this result does not contribute necessarily for a more realistic description, in terms of the distribution of water outputs from the aquifer on different discharge areas, than a model whose flow domain is characterized as homogeneous, in terms of T. It must also be considered that both simulated balanced water fluxes have higher values than streamflow discharge values measured at gauges 6 and 1. Simulated Q quantifies groundwater transferences occurring on the full extent of the selected creek reaches, i.e., the sum of all the contributions for base flow at a linear section of the creeks, whereas data collected at gauging stations merely quantifies streamflow discharge at a given point on the creek’s course. In order to fully understand and quantify the groundwater component of flow occurring at creeks, the use of surface models coupled with groundwater models is recommended for future simulations.

CONCLUSIONS

The developed work was based on a sequence of operations started by the implementation of a bi-dimensional numerical representation of AO. The conductive parameter T was estimated by inverse modelling for 16 zones, using the GaussMarquardt-Levenberg method, implemented in a non-linear parameter estimation software, coupled with a standard finite-element model based on the Galerkin method of weighted residuals.

275

Calibration results revealed a good fit between simulated and measured hydraulic head values, however, transient simulations of stream-aquifer hydraulic transferences allowed verifying that a calibration of T providing a reliable representation of the spatial distribution of hydraulic head, does not imply necessarily a better distribution of observed discharges on different aquifer output areas, in comparison to what is simulated considering homogeneous T. In order to improve the process of calibration of the AO numeric model, simulations coupling surface models, targeting the achievement of more realistic hydrograph results, should be performed.

Acknowledgements The authors gratefully acknowledge the Portuguese Science and Technology Foundation (FCT) for the financial support of the POCTI/AMB/57432/2004 Project “Groundwater flow modelling and optimisation of groundwater monitoring networks at the regional scale in coastal aquifers – The Algarve case study” (CVRM Geo-Systems Centre, Algarve University).

REFERENCES Almeida, C., Mendonça, J. L., Jesus, M. R. & Gomes, A. J. (2000) Sistemas Aquíferos de Portugal Continental. [Aquifer Systems of Mainland Portugal]. INAG. Relatório Técnico [Technical Report], CD-ROM, 661pp. Bear, J. & Verruijt, A. (1987) Modeling Groundwater Flow and Pollution. Dordrech. The Netherlands. D. Reidel Publ. Co., 414pp. Carrera, J., Alcolea, A., Medina, A., Hidalgo, J. & Slooten, L. J. (2005) Inverse problem in hydrogeology. Hydrogeology Journal 13(1), 206-222. Doherty, J. (2002) PEST, Model-Independent Parameter Estimation. 4th Edition, Watermark Numerical Computing, Australia, 279 pp. Huyakorn, P. S. & Pinder, G. F. (1983) Computational Methods in Subsurface Flow. New York. U.S.A.. Academic. Press., 473pp. Kinzelbach, W. (1986) Groundwater modelling. An introduction with sample programs in Basic. Developments in Water Science nº25. New York. U.S.A.. Elsevier, 333pp. Martins, J. (2007) Inverse Calibration of a Groundwater Flow Model for the Almádena-Odeáxere Aquifer System (Algarve – Portugal). MSc. Thesis, Faculdade de Ciências do Mar e do Ambiente, Universidade do Algarve. Martins, J. & Monteiro, J. P (Accepted) Coupling Monitoring Networks and Regional Scale Flow Models for the Management of Groundwater Resources. The Almádena-Odeáxere Aquifer Case Study (Algarve-Portugal). In: HydroPredict2008. Prague. Czech Republic. Monteiro, J. P., Nunes, L., Vieira, J., Stigter, T., Martins, R., Santos, I. J. & Reis, M. E. (2003) Sintese Bidimensional dos Modelos Conceptuais de Funcionamento Hidráulico de Seis Sistemas Aquíferos do Algarve. [Bidimensional Sinthesis of Hydraulic Behaviour Conceptual Models of Six Algarve Aquifer Systems]. In: Proceedings of Jornadas Luso-Espanholas Sobre as Águas Subterrâneas no Sul da Península Ibérica. Universidade do Algarve, Faro. Portugal, 10pp. Monteiro, J. P., Vieira, J., Nunes, L. & Younes, F. (2006) Inverse Calibration of a Regional Flow Model for the QuerençaSilves Aquifer System (Algarve-Portugal). In: Integrated Water Resources Management and Challenges of the Sustainable Development. Marrakech. Morocco, CD-ROM, 6pp. Monteiro, J. P.; Oliveira, M. M. & Costa, J. P. (2007a) Impact of the Replacement of Groundwater by Dam Waters in the Albufeira-Ribeira de Quarteira and Quarteira Coastal Aquifers. In: XXXV IAH Congress. Groundwater and Ecosystems. Lisbon. Portugal, pp 489-480, CD-ROM, 10pp. Monteiro, J. P., Ribeiro, L., Reis, E., Martins, J. & Silva, J. M. (2007b) Modelling Stream-Groundwater Interactions in the Querença-Silves Aquifer System. In: XXXV AIH Congress. Groundwater and Ecosystems. Lisbon. Portugal, pp 4142, CD-ROM, 10pp. Poeter, E. P. & Hill, M. C. (1997) Inverse models: a necessary next step in ground-water modelling. Ground Water 35(2), 250-260. Reilly T. E. (2001) System and Boundary Conceptualization in Ground-Water Flow Simulation. Techniques of waterResources Investigations of the United States Geological Survey. Book 3, Applications of Hydraulics. Reston. Virginia. U.S.A., 30pp. Reis, E. (1993) Estudo Hidrogeológico das Formações do Lias-Dogger situadas a ocidente do Rio Arade (Algarve), [Hydrogeologic Study of the Lias-Dogger Formations West of River Arade (Algarve)]. MSc. Thesis, Departamento de Geologia, FCUL. Vieira, J. & Monteiro, J. P. (2003) Atribuição de Propriedades a Redes Não Estruturadas de Elementos Finitos Triangulares (Aplicação ao Cálculo da Recarga de Sistemas Aquíferos do Algarve). [Attribution of Properties to Non-Structured Triangular Finite Element Networks (Application to the Calculation of Recharge of Algarve Aquifer Systems)]. In: Proceedings of Jornadas Luso-Espanholas Sobre as Águas Subterrâneas no Sul da Península Ibérica.

276

Faro. Portugal, 10pp. Wang, H. F. & Anderson M. P., (1982) Introduction to groundwater modelling. Finite difference and finite element Methods. Freeman. New York, 237pp.

277