Modeling of Groundwater Flow and Water Use for San Luis Potosí Valley Aquifer System

Land use changes are currently one of the indisputable factors in the alteration of processes and cycles of the aquifer system in the San Luis Potosí Valley. Due to its importance, is considered indispensable to investigate this detrimental factor of aquifers. The aim of this research is to use a numerical flow model to analyze the impact that land use changes have had on the aquifer. A finite differences numerical model was adapted to the size and hydrological properties of the aquifer system. It consisted of a regular grid with 30 columns and 34 rows with constant spacing of 1000 meters. It has two layers; the first includes the shallow aquifer and the second, the deep aquifer. The initial hydraulic head of the model corresponds to 1986 and was verified for 1995 and 2007. The model shows the development of a drawdown cone (central valley) extending toward the industrial area (southern valley). Piezometric water levels revealed a decrease of 0.6 to 1.6 meters annually during a period from 1977 to 2007. This work demonstrates that it is the consequence of land use changes and of the incessant overall decline in groundwater reserves. Based on the flow model, population growth projections and water use change, the calculated predictions indicate that by 2021, the total established volume of 136 Mm/year for consume will be reached. The flow model of the San Luis Potosí Valley aquifer system shows a clear effect of the risks associated with aquifer mining.


Introduction
The San Luis Potosí Valley (SLPV) aquifer is located in the southwest region of the state of San Luis Potosí (SLP), with an approximate area of 1,980 km 2 . In terms of water management, it is part of "El Salado" Watershed Number 37 (Figure 1), and represents a zone in the SLPV with favorable hydrogeological conditions that have enabled the intensive exploitation of groundwater. The study area has a temperate climate with a warm-arid summer. The mean annual rainfall is 402.6 mm and the mean annual temperature in the region is 17.5 °C; mean annual potential evaporation is 2038.7 mm.
In the early 19th century, wells supplied water to the population of SLP (Sheridan-Prieto, 2001). By 1960, 59% of domestic water was supplied by the shallow aquifer and 41% by the deep aquifer. The supply sources have changed dramatically in recent years so that 92% of the water is obtained from the subsoil and only 8% from surface water (COTAS, 2005). The combined municipalities of SLP and Soledad de Graciano Sánchez (SGS), which embraces the urban area of the valley and they require the main demand for water, grew 30% between 195030% between and 197030% between . Between 197030% between and 1980, a staggering population increase of more than 50% was recorded ( Figure 2). In 2000, the amount of inhabitants in the two municipalities was equivalent to 38% of the total state population, and between 2000 and 2005 (Figure 2), the SLP metropolitan area increased 32%, exceeding the growth of 29% between 199029% between and 200029% between (INEGI, 198129% between , 199129% between , 200129% between , 2008. This rapid urban growth strongly impacts land use, potentially affecting recharge and groundwater extraction (Noyola et al., 2009).
Flow models have been used to carry out prognostic studies and to develop management schemes for the effects of groundwater pumping (Anderson & Wang, 1982;Sanford et al., 2004;Milzow et al., 2009;Yustres et al., 2013). Early works to develop mathematical modeling in Mexico were performed by Cruikshank (1982) and Herrera (1982 and1989) for the Basin of Mexico. Over the last decade, 18 mathematical model applications were achieved to manage groundwater in the state of Guanajuato (Chávez et al., 2006). In SLP, similar studies were performed (SARH-CNA-UNAM, 1992;Difurt et al., 1995 andCNA, 1996-a) adapting coarse mathematical models to the SLPV aquifer system. Kohn-Ledesma (2010) developed a flow model for the SLP aquifer to propose extraction alternatives that would increase the balance of the aquifer.
Groundwater models were also applied to study the dynamics response of the ground water systems (Bredehoeft, 2002). The purpose of this article is to analyze an important portion of the groundwater flow system in the SLPV to explore the impact of water and land use changes on the deep aquifer and on the groundwater resource management. For this, a transient groundwater flow model was adapted and calibrated in the MODFLOW computer code (MacDonald & Harbaugh, 1988) to quantify changes in the flow system. Although this kind of studies have been conducted in other parts of the world, there are not research works on the effects of water use changes on the aquifer system of SLPV and on how this factor may affect the groundwater.

Economy and Population
The SLPV aquifer is one of the most important in the state of San Luis Potosí because it supplies the capital city, which is responsible for 80% of the state's production. It also partially includes the municipalities of SLP, Soledad de Graciano Sanchez (SGS), Mexquitic de Carmona, Cerro de San Pedro and Zaragoza. The total population of these municipalities is 1, 122.502 (INEGI, 2010), which depends on groundwater to supply 95% of (COTAS, 2005) the various uses and consumes. Urban and geographic charts reveal that in the regional development in the last decades, social and economic considerations are among the most important drivers of landscape change (Turner et al., 1996).

Water and Land Use
Water and land use changes can substantially alter hydrologic ecosystem services (Kepner et al., 2012). In recent decades, land use change in the SLPV has become one of the indisputable factors in the alteration of the processes and cycles that have existed since the early years of the last century. Figure 1 shows the distribution of the three main land uses in the study area.

Urban Land Use
One of the first and principal factors for changes in the SLPV is the urban land use. Nowadays it remains growing and changing over time. Population growth coupled with rapid urban expansion and increasing economic and productive (industrial and agricultural) activities have put severe pressure on natural resources in the SLP Valley and on groundwater consume, in particular.
Therefore, the greatest urban demand has been concentrated in the central valley (urban zone), greatly stressing and impacting the aquifer system and causing a large drawdown cone.
Some of the conditions created by the urban area include domestic and industrial wastewater discharges which travel through open channels. These wastewaters are used in the agriculture, impacting directly the perched aquifer. In fact, it is wrong handling of wastewater, since their use occurs primarily in the agriculture with the argument that the water quality is not suitable for human consumption. The accelerated urban and population growth in recent years and has led invariably to a significant change that greatly impacts the aquifer ( Figure 2). The most significant change in urban land occurred between 1970and 1993(Noyola et al., 2009, when the area of the urbanized zone quadrupled and the suburbs in SLP experienced a 15-fold increase in size from 1959 to 2005 ( Figure 2). According to Peña (2006), the distribution of the use of extracted water was 90.82 Mm 3 /year for domestic use and 3 Mm 3 /year for other uses.

Agricultural Land Use
Historically, agricultural development in the valley has been directed towards increasing the productivity and exploitation of natural resources while ignoring the complex interactions and implications of these processes. As a result, agriculture has represented for a long time the second largest land use and has had a high impact on the shallow aquifer. Artificial recharge induced by the use of wastewater for agricultural irrigation has resulted in hydrogeochemical anomalies in the shallow aquifer, which may result from the impact of the agricultural zone on the aquifer. The volume water for agriculture use was 8.57 Mm 3 /year.

Industrial Land Use
Although mining was the dominant industrial activity in SLP for centuries, other industries have developed in recent decades, including food, automotive, chemical, textile, paper, steel and metal mechanics. The increasing industrial activity and its corresponding displacement of the peri-urban, agriculture have intensified the demand for water, since agricultural demand for the resource is cyclical. Today, the industrial and urban changes in land use require a constant volume and permanent supply of water (SEMARNAT, 2008). The consequence of this land use change is the development of a large drawdown cone in the urban area, extending to the industrial zone. Industrial demand increase in nine years (1996)(1997)(1998)(1999)(2000)(2001)(2002)(2003)(2004)(2005), from 7.6 to 13.24 Mm 3 /year, while in the same period agriculture decreased from 19.55 to 8.57 Mm 3 /year (Table 1). SEMARNAT (2008) recognized that the growing industry currently occupies spaces that were intended for agriculture, claiming 14% of the water extracted from the aquifer. The distribution of the use of extracted water was 13.24 Mm 3 /year for industrial use.

Geology
The SLP tectonic valley was filled with alluvial deposits that originated in the Sierra de San Miguelito (SSM), which resulted in a graded distribution, with preferential northeast direction. It is bounded to the west by staggered faults on the edges of the SSM and the Sierra de Alvarez (SA) (Figures 3 and 4).
The stratigraphic column in the region includes sedimentary and volcanic rocks ranging in age from Cretaceous to Recent (Figure 3). In the valley, the Cretaceous rocks are mainly thin-layered limestone of Indidura and Cuesta del Cura formations, which have low hydraulic conductivity because these units are from the marine www.ccsenet.org/jgg Journal of Geography and Geology Vol. 6, No. 3; basin. Unconformable covering these units are Cretaceous Tertiary volcanic rocks consisting of basalt spills, latites and ignimbrites. These fractured volcanic units reach a thicknesses up to 500 m and are covered by semi-consolidated granular material with greatly variable thicknesses (from 80 m to 350 m and, in some cases, over 500 m thick). Observed in some areas are interbedded basalts with Quaternary sediments, rhyolitic ignimbrites and tuff sands from the Middle and Upper Oligocene (Labarthe et al., 1982).
Groundwater in the study area is present in Quaternary aquifers with porous alluvial deposits, which consist of laterally discontinuous layers of gravel, sand and clay, and Tertiary volcanic rocks (Figure 3). The aquifer system can be divided according to shallow and deep units. The shallow unit or granular aquifer consists of sand and gravel and extends to a depth of approximately 15 to 200 m below the ground surface. This unit has good hydraulic conductivity and water production. Above the shallow aquifer, there is a perched aquifer which is poorly distributed in the valley. Its depth is variable with a maximum of 40 m (Figure 4). The perched aquifer is not hydraulically connected with the deep unit.
The thickness varies in the deep unit; its bottom is between 300 and 800 m below the ground surface and it is mainly composed of fractured volcanic rocks with good hydraulic conductivity. Because the shallow aquifer unit has been exploited and almost exhausted, pumping wells have been extracting groundwater from the deep unit. Therefore, for the purpose of this study, the deep unit is considered a flow system for a confined aquifer.

Conceptual Model
Large aquifers are relevant due to their considerable size and as potential water sources. Regional groundwater modeling of this kind of aquifers remains a challenge because the geologic complexity and data scarcity in space and time (Rodriguez et al., 2012). The model includes an area that extends from the airport in the northern region of the city of San Luis Potosí to the town of La Pila in the southern part of the valley (groundwater divide), and from the SSM border in the west to the SA border in the east ( Figure 4). The topography was obtained from digital terrain models (INEGI, 2007). The conceptual model (aquifer system) consists of two main units; a shallow aquifer which consists of a granular media (Layer 1). The second unit, the deep aquifer, consists of fractured volcanic rocks (Layer 2). The hydrogeological basis consists of argillaceous limestone, laminated with low hydraulic conductivity ( Figure 3).

Hydrogeological Evolution and Hydraulic Heads
The initial hydraulic head corresponds to the year 1986. Water table (  Based on the 2003 piezometric network, the drawdown cone is observed to be predominant in the city of SLP. The hydraulic gradients converge towards the regional drawdown cone and the depth of the cone increased 60 m from 1971 to 1995 (CNA, 1996-b). The evolution of the water table for the period 1995-2001, its depth increased as much as 25 m in the center of the drawdown cone and developed towards the north, while in certain areas there were no water level drawdowns (COTAS, 2005). In 2007, a new direction in the expansion of the drawdown cone was observed, which goes from the center of the SLP urban area toward the southern region of the city.
For calibration, historical data were used from water tables from 21 wells in the monitoring network of the deep aquifer. The depth of the wells varies from 200 to 726 m.

Aquifer Properties
The values used for the initial estimates of hydraulic conductivities, the storage coefficient, transmissivity and specific storage for the domain model were based on those reported in previous studies (Cardona-Benavides, 1990, 2007Carrillo-Rivera et al., 1992, 2002, as well as values commonly accepted in the literature for materials similar to those in the study area (Figures 3, 4, 6 and Table 1).
Hydraulic conductivities in layer 1 are related to the distribution of materials in the alluvial valley ( Figure 6). They have a range of conductivities from 0.1 to 50 m/d. The lowest values are located in the northeast part of the valley, while the highest are distributed in areas that coincide with the presence of geological faults (Figures 3, 4, and 6).
In layer 2, hydraulic conductivities have a range of conductivities from 2 to 10 m/d. The highest are distributed in areas that coincide with geological faults too (Figures 3, 4 and 6).

Numerical Model
The numerical flow model was constructed using the MODFLOW (McDonald & Harbaugh, 1988) algorithm based on block-centered finite differences and a three-dimensional graphic interface. The flow model was calibrated in stationary conditions.

Equations
A general solution for the response of multiple aquifer systems to pumping stress requires solving the three-dimensional flow equations (Bredehoeft et al., 1970). The three-dimensional flow of constant-density groundwater through a porous medium can be described by the partial differential equation: (1) Where h is the head hydraulic (L), t is time (T), Kxx, and KZZ Kyy are values of the hydraulic conductivity along the x, y and z coordinate axes, which are assumed to be parallel to the principal axes of the hydraulic conductivity (LT-1), W is the volumetric flow rate per unit volume and represents water sources and/or sinks (T-1) and Ss is the specific storage of porous material (L-1).
The MODFLOW computer code developed by the USGS (McDonald & Harbaugh, 1996) was used to solve the flow model and the pre-and post-processor Visual MODFLOW (VM) Pro version 4.1 (a product owned by Waterloo Hydrogeologic Inc., 2005) was used to plot the input data for analysis and presentation of the output data.

Spatial and Temporal Discretization
The discretization of the finite-difference flow model is shown in Figure 6. The grid contains 30 columns and 34 rows with a constant spacing of 1000 m. Vertically, two layers were defined in the model to simulate the two-dimensional regional flow system. The first layer comprises alluvial deposits and granular medium (shallow aquifer) and the second represents fractured volcanic material (deep aquifer). The elevations of the roof and the bottom of the layers were interpolated from a geoelectric model of the study area. In the present study, we simulated the groundwater flow in SLPV for the time period from 1986 to 2007.

Boundary Conditions
The physical and hydraulic boundaries in the study area were defined from available geological and hydrological information ( Figure 6). On the east and west edges of the interface model is located the interface between the plains and the mountainous terrain of the SSM and SA, respectively. The mountain-front recharge enters the study area laterally, while the east and west boundaries of the model were defined as no-flow boundaries. Dirchlet (constant head) boundary conditions were imposed on the northeast and southwest borders. The values of the constant heads for 1986 in the northeast and southwest borders were 1753 and 1744 m, respectively (Martínez & Aguirre, 1987).

Recharge and Lateral Input
The main sources of groundwater recharge for the shallow aquifer in the study area are lateral flow from the surrounding mountains and irrigation return flow. The fact that evapotranspiration and precipitation rates are similar limits the deep aquifer recharge.
The rivers that run through the valley are intermittent and therefore were not considered as elements of aquifer recharge. Nevertheless, open channels are located in the agricultural zone of SGS, through which wastewater flows constantly. Although farming generates artificial recharge due to irrigation returns, it goes into the perched aquifer, but only in areas where the aquifer behaves freely, since the rest of the area is urbanized, covered by a thin layer of clay or layer caliche which prevents infiltration into the perched aquifer.
Three distinct recharge areas were determined. The most influential was for the agricultural area (15 mm/year), the second the urban area (10 mm/year) and the third were for the rest of valley (7 mm/year) (Figure 7).
Lateral recharges from the SSM and the SA, which constitute the most important orographic areas in the SLPV, contribute little due to low permeability, making infiltration into the subsoil very difficult. Regional recharge is mainly generated in the northwest portion of the study area. The main discharge in the study area is from pumping for domestic and industrial uses in the city of SLPV and pumping for irrigation in neighboring municipalities. The maximum volume allowed by the authorities to be extracted from the aquifer system is 136 Mm 3 /year (INTERAPAS, 2004). Individual point wells were used to represent pumping in the city of SLPV.
An analysis of the total volume extracted from the aquifer system by pumping for different uses was conducted by INTERAPAS (2004) for the period 1995-2005. When considering the volume of water demand according to type of use, the majority of the volume extracted was determined to be for human use and consumption, followed by industrial and agricultural uses.
Since 1950, the population has grown at an annual average rate of 3% with a demand of 200 liters per capita daily. Thus by the year 2005 the urban use represented 78.2% of the total extraction from the deep aquifer. Industry demanded 11.4%, with an annual growth rate of 0.63% through 2005. Meanwhile, the agricultural zone has experienced an annual average decrease of 1.22% in 2005, with a demand of 7.4% of the total extraction from the deep aquifer (Table 1).

Calibration
The numerical model was calibrated to reproduce the measured groundwater levels. This calibration includes two sequential steps. First, the steady-state model was calibrated and developed by adjusting the input parameters, such as hydraulic conductivity and boundary conditions that were able to reproduce behaviors whose characteristics represented the developed aquifer system, until obtaining the optimal condition for the modeling system. The criteria used for the best calibration was mean squared error (

Results
During the calibration tests, several experiments were performed with meshes of 250×250 and 500×500 m. However the mean square error (MSE) was high and the results were unrealistic. The best results in the calibration process were obtained with a grid spacing of 1000×1000 m. This mesh is suitable for the available information and for the dimensions of this basin. Smaller grid spacing does not lead necessarily to better results. In other studies for similar basin dimensions, the same grid spacing of 1000 m was also applied (Sanford et al., 2004;Milzow et al., 2009;Yustres et al., 2013). Numerical modeling has allowed simulating an MSE of 0.26 for the transient-state, for the initial conditions of 1986 (Figure 8a).
MSE of 1.44 and 1.70 were obtained for the transient states in 1995 and 2007, respectively (Figures 8b and 8c).
Using the flow model, it was possible to determine a preferential direction from north to south during initial conditions, which were very different from nowadays conditions. The flow directions have been changing as the aquifer has been evolving due to extraction so that the contour of water tables now shows the direction of flow converging towards the urban area, resulting in a drawdown cone (Figure 9).
During the simulation period of 1986, the drawdown cone was observed to be located mainly in the urban areas (central valley) and another in the northeastern portion of the valley, while for the year 2007, the cone in the urban zone began to develop towards the industrial area (southern portion of the valley) as a result of increased extraction due to a change in land use by the agricultural to industry.
Based on the flow model and considering population growth projections and the analysis of land use changes, predictive simulations were performed to define groundwater management scenarios, which found that the amount allowed by CONAGUA would be reached by 2021 (Figure 9).
Considering the recent developments in the urban and peri-urban valley since 1950, predictive simulations were performed in which urban and industrial changes in land use were inversely proportional to the change in agricultural use (Figure 9). It impacted directly the recharge rate in the agricultural zone. Another way to see the exchange of water used by agriculture with that used in the industry is to consider the percentage change between 1996 and 2005.  Figure 10 shows that human consumption impacts more on the total extraction system of SLPV, it has increased steadily since the end of last century. In the period under review, there has been an increase in volume for industrial use. This volume was compensated by a decrease in the extracted volume for agriculture.
These changes are proportional, so that they do not represent a major impact on the trend of the total extraction. But in time, this proportionality does not occur. There will be an impact on the trend of the total extraction ( Figure 10).
industrial land use. In the former case, areas considered to be natural recharge (SSM) for the shallow aquifer have been altered by the urban use of land with the consequence of little or no recharge in the shallow aquifer, causing its depletion and promoting greater extraction of groundwater from the deep aquifer to supply the population.
Similarly, areas that were traditionally of agricultural land use have gradually become urban areas.
Additional projections with the 50-year model calculation show that the drawdown cone develops throughout the valley (Figure 10).
Looking for alternative sources of water supply for the different uses by the population is a priority, given this scenario and the continued increase in the population and industrial activity.
Based on the model results, it is concluded that the distribution of extraction wells, located in the neighborhood of the hydrological boundaries without flow contributions, causes a large drawdown cone and land subsidence in the urban area (López-Álvarez et al., 2013).
The local water sources in the basin are almost depleted. Inside the VSLP, it is necessary to take steps to reduce the extraction in the system and to reuse wastewater in the industry and agriculture. Since extraction for human use impacts predominantly the system, it would be desirable that the water authorities promote campaigns of water culture (water use efficiency and appropriate use). It is also desirable to establish the endowment reduction of the volume of water for each person as well as political and water managements. There should be also considered the water extraction of new areas.

Conclusions
The SLPV has undergone changes in land use over the past 20 years which impacted greatly the aquifer system in the study area.
The main land use changes took place from agricultural to industrial and urban land use.
With the application of the flow model, it was possible to show that the induced recharge has a greater effect www.ccsenet.org/jgg Journal of Geography and Geology Vol. 6, No. 3; than the natural recharge, especially on the shallow aquifer, while for the deep aquifer both types of recharge have minimal or no effect.
Numerical simulations show that land use changes have had the greatest impact on the growth of the drawdown cone towards the industrial area.
For the predictive phase, it was identified that the volume established by CONAGUA will be reached in 2021.
The aquifer of VSLP had a drawdown of 199 meters in the period 1997-2007. If this scenario continues, the drawdown will continue in the valley, causing an overall decline in groundwater levels.
The flow model for the SLPV aquifer system shows a clear effect of the risks associated with mining aquifer. Therefore, groundwater management policies appropriate to the precarious conditions of the system need to be implemented.