Chihuahua-Sacramento and Tabalaopa-Aldama aquifers of Chihuahua state, Mexico: linkage and importance of geostatistical and hydrogeochemical analysis

Over the evolutionary period of life on the planet, and especially since the appearance of man, water resources were exploited in a balanced and stable way in regard to the ecosystem. The gradual evolution of human intelligence did not produce enough awareness for the rational use of resources. In other words, resources are being used in accordance to population growth and its needs. In the State of Chihuahua, the availability of the liquid in its vast arid and semi-arid regions is limited in quantity and quality, and it determines the State's economic activity and life. This problem has worsened over the years due to prolonged periods of drought in large areas of the region (Herrera Peraza et al. 2016).

Groundwater is widely used in many countries. It is often the primary source of drinking water (supplying half of the world's population) and it contributes significantly to irrigation; thus, it provides food security in arid and semiarid regions (Margat & Van der Gun 2013).

Groundwater systems have relatively large volumes of water in storage providing a unique buffer capacity, the groundwater systems' major strength. It allows the periodic and seasonal conditions for survival in semi-arid and arid regions and it helps to reduce the risk of temporary water shortages (enhancing the value of groundwater), specially groundwater, which represents 97% of the available freshwater in the world and is the main source of water for one third of the world's population (Hunter et al. 2016).

Renewable groundwater resources in the arid and semi-arid zone are less widespread and variable but with a generally limited magnitude (CGMWW/UNESCO Geological map of the world 2000). In the countries with highest withdrawal of non-renewable groundwater, this resource is significant and it is usually the predominant source of water. It covers a major part of the total water demand (predominated by irrigation water demand). Even though groundwater provides nearly 50% of all drinking water (Blanco & Donoso 2021), there is much uncertainty and inconsistency in available data on non-renewable groundwater extraction; however, in Figure 1 the distribution of mean annual groundwater recharge in mm/year (1961–1990) is shown (Foster & Loucks 2006)

According to Figure 1, we can see that the study region has an estimated recharge of up to 20 mm/year; it is estimated that 60% of the total volume of water used in the region is extracted from groundwater to be used in all socioeconomic activities (Plan Estatal Hídrico del Estado de Chihuahua PEH 2040a, 2040b).

More than 90% of the drinking water supply for the population of the state of Chihuahua is obtained from groundwater, the most important cities being Ciudad Juarez and Chihuahua.

Derived from the analysis carried out by the international association of hydrogeologists on water utilities and groundwater, it is necessary to raise these concerns:

What are the risks to groundwater resource sustainability that need to be addressed by water utilities? For waterwell use to be sustainable (and available at times of severe water stress), groundwater resources need to be pro-actively managed to avoid:

• Permanent depletion due to overexploitation

• Pollution from uncontrolled on-site sanitation, agricultural land-use and industrial activities. (IAH 2021)

Water utilities, being the major stakeholder in potable water supplies, should formally recognize and embrace their co-responsibility for groundwater resource management and protection to avoid the aforementioned impacts. That is why the Municipal Board of Chihuahua supports this work in search of elements that define the relationships between the factors that intervene in the use of groundwater (quantity and quality).

In the state of Chihuahua there are 61 identified aquifers; of those, overexploited aquifers increased from 11 in 2011 to 40 aquifers in 2020 (SINA-CONAGUA 2020). Figure 2 shows (in red) the overexploited aquifers, according to the latest publication of the regulatory entity in Mexico (CONAGUA).

In terms of aquifer balance, five of them: Chihuahua-Sacramento, Laguna de Mexicanos, Cuauhtémoc, Villa Ahumada-Flores Magón and Jiménez-Camargo, report a severe condition.

The water that supplies the City of Chihuahua is extracted from two different aquifers:

Chihuahua-Sacramento (ACHS) or 0830 aquifer according to the classification of the Groundwater Management of the National Water Commission (CONAGUA), which currently has 69 active wells and is the city's main and oldest drinking water supplier. ACHS contributes the largest volume of water: 49% of the total annual volume extracted for drinking water.

The second aquifer, Tabalaopa-Aldama (ATA) or 0835 aquifer, has 42 active wells for drinking water and provides 30.14% of the total annual volume for supply of the city.

The objectives of this study are to investigate the available links of information to determine the combination of factor analysis and cluster analysis as a tool to find:

• (a)

  The variables that control the regional and intermediate flow systems.

• (b)

  The way to differentiate the flows in a heterogeneous aquifer of great thickness.

These will create, from the information obtained, a solid base that will contribute to the Groundwater Protection Program in the State of Chihuahua.

The Environmental Protection Agency defines the hydrogeologic protection zone as ‘the surface or subsoil area surrounding a water catchment or well field that serves as a drinking water supply to a community’ (Seur Spencer & Drake 1987).

According to the update of the annual groundwater average, published from 2015 to 2020 by CONAGUA, there are several technical studies and their results conclude that a modification in the availability of groundwater in the analyzed regions has occurred due to the change in the natural recharge regime, concessioned volume and committed natural recharge, which has modified the value of the average annual water availability.

The elements that intervene in the quantification of the condition of availability or overexploitation of aquifers are identified below:

KEY: This is ID granted by the National Water Commission (CONAGUA)

ACUIFER: The acronym of the aquifer

R: Represents the volume of the average annual natural recharge (Mm³/year)

DNCOM: The volume of the committed annual discharge (Mm³/year)

VCAS: The volume granted or registered with CONAGUA for the exploitation of groundwater (Mm³/year)

VEXTET: Magnitude of the volume of groundwater extraction identified in technical studies (Mm³/year)

DAS: This means the average annual availability of groundwater (Mm³/year)

DEFICIT: This refers to the real imbalance for the availability of water consumption (Mm³/year)

Table 1 shows the values that were published by CONAGUA to determine the exploitation conditions of the aquifers.

All definitions of the above terms are contained in the Mexican Official Standard NOM-011-CONAGUA-2015. They are also based on the Secretary of Environment and Natural Resources Decree (2020), which endorses the National Water Program 2020–2024 and on the aquifer maps reported by CONAGUA in the Manual of Average Annual Groundwater Availability (DOF 17/09/2020).

The results of the physicochemical analysis of 59 groundwater samples obtained during the study performed from 2003 to 2016 by the Municipal Water and Sanitation Board of Chihuahua (JMAS) in wells for drinking water supply were considered. The parameters considered for determining the relationships between water quality and flow patterns are: iron, manganese, fluorides, chlorides, temperature, electrical conductivity, pH, total hardness and total dissolved solids.

One of the most important parameters to determine the directions or patterns of flow is the concentration of total dissolved solids. In general, concentrations come from the north and from the foothills of the mountains towards the center of the aquifers. This shows us that the recharge zones are in the water part of the mountains and their discharge zones occur in the outcrops of the rivers that cross the aquifers. It can be observed that the concentration of Total Dissolved Solids, in most cases, is less than 600 ppm, less than the 1,000 ppm established by the Official Mexican Standard NOM-127-SSA1–1994 for water intended to be consumed by humans.

The aquifers ACHS and ATA have a design of a hydrogeological monitoring network, carried out in the year 2011, and it can be consulted in the study called ‘Evaluation of the sources current supply to the city of Chihuahua, feasibility study of alternative sources and preliminary draft of the necessary hydraulic infrastructure’, (JCAS/CONAGUA 2012(a), 2012(b)).

In the period from 2007 to 2009 and from 2011 to 2020, with depth data at groundwater level in meters (m), elevation and evolution of the water level, the depth settings were measured and built. These configurations show the local and regional behavior of the aquifer system in the present study.

The following figures show the static level configurations and the preferential flow directions for the year 2020, Figure 3 for ACHS and Figure 4 for ATA. In both cases, it is observed that the parallelism is maintained with the surface hydrology to the rivers that cross the aquifers, in the case of ACHS from north-west to south-east and in ATA north-west and from south-west to east.

Both Figures 3 and 4 show the outline of the city of Chihuahua (in yellow) and the limits of the aquifers. In ACHS there are also some small recoveries, located in the south central portion of the aquifer, precisely in the urban area from the city of Chihuahua, with values from 1 to 4.4 m. In general, the aquifer presents abatement of around −25.58 m in the period of 10 years, which would represent an average annual drop of −2.55 m. In ATA it is observed that the minor abatements occur in the central portion of the aquifer, along the entire length of the Chuvíscar River and its surroundings and it increases as it moves away from that area towards the mountains, both Name of God as Saint Ignatius. In general, the aquifer shows an abatement of around −27.7 m in a period of 9 years, it represents a low annual average of −3.1 m.

In the period from 2007 to 2009, the evolution of the Chihuahua-Sacramento and Tabalaopa-Aldama aquifers was mapped with piezometric information. However, in the case of the ATA, there is not much piezometric information to evaluate in detail the configuration of the evolution of the static level. But it is known, and it can be observed in Figure 3, that the incipient extraction that has been taking place has not apparently caused any alteration to the original conditions of the subsurface flow regime.

Figure 3 shows the boundaries of the Chihuahua-Sacramento aquifer. In the aforementioned figure, the hydrogeodynamics of the underground water flow are represented with small black arrows. The hydrogeochemical and hydrogeodynamic aspects for both aquifers are the same and some of them are described in a few paragraphs below. The border with the Tabalaopa-Aldama aquifer can also be seen. Figure 4 shows the static level elevation (masl) updated to 2014, according to the data reported by the JMAS of Chihuahua and taken up by the National Water Commission in the study to update the availability of both aquifers (CONAGUA 2015). Also, in this Figure, the similarity of the cuts of the geological sections A-B and C-D is shown.

Figure 5 shows the C-D geologic section located in the foothills of La Haciendita and Nombre de Dios in the tectonic basin, with a probable maximum thickness of no more than 600 m in the Chihuahua City area.

According to the results of the physicochemical analyses, it can be observed that the concentration of Total Dissolved Solids, in most cases, is less than 600 ppm, less than the 1,000 ppm established by the Official Mexican Standard NOM-127-SSA1-1994 for water intended to be consumed by humans. In general, concentrations increase from the north and from the foothills of the mountains towards the center of the aquifer.

The presence of its distribution and behavior of groundwater, initially, is dependent on the hydrological cycle of the area, and later on the management and use that occurs in the region. As to the hydrological cycle, we can highlight the aridity index. The level of aridity is defined by the Aridity Index (AI) which is the amount of precipitation divided by the potential evapotranspiration (Poeter & Ying 2021).

Arid areas in Chihuahua have an AI between 0.05 and 0.2, they lose a large fraction of precipitation to evapotranspiration and only a small amount reaches the water table (from near 0% up to about 4% of precipitation); temperature plays an important factor in the magnitude of recharge.

The present work uses multivariate analysis as a tool to find the relationships between the factors (water quality, climate and hydrogeology), which objectively make visible where (flow order or system of groundwater flows) the exploitation of groundwater occurs.

The present work uses multivariate analysis as a tool to find the relationships between the factors (water quality, climate and hydrogeology) and it objectively shows where (flow order) the exploitation of groundwater occurs. This indicates the degree of exploitation and problems at region.

Remembering that in most flow networks and in most field areas, it is possible to differentiate between local groundwater flow systems, intermediate groundwater flow systems and regional groundwater flow systems, as illustrated schematically (Freeze & Cherry 1979),

Analytical or mathematical simulation methods are normally used to determine groundwater flow systems; however, they are not very assertive in the differentiation of this order of groundwater flow systems. In these analysis tools, a strong weight is the expertise of the modeler, whereas, in a multivariate analysis it is independent of the human component to determine the relationships between factors.

In several published works on systems of groundwater flows and in contrast to the geological environment, as in the valley of San Luis Potosi (Carrillo-Rivera et al. 1996), (Carrillo-Rivera et al. 2007) and (Carrillo-Rivera & Cardona 2008). The statistical methodology to determine the components of the intermediate flows (IF), as well as the mixed or intermediate (MF) and regional (RF) flows; concluded that for the identification of groundwater flow systems, the factors that determine a change between flow systems are fluoride and water temperature. Based on these previous works, the process determined the variables or factors that have an association, relationship or impact on the expected response or ‘groundwater flow systems’.

The statistical elements used are: Correlation Matrix (CM), Principal Component Analysis (PCA) and Factor Component Analysis (FCA) with factorial rotation (VARIMAX). FCA is a data dimensionality reduction technique whose purpose is to find the minimum number of factors or variables that explain the information contained in the sampled data. In Factor Component Analysis, all variables have the same role, since they are independent and, a priori, there is no conceptual dependence between them.

For the statistical analysis, the time series from 2003 to 2016 was used. The following parameters are related from both aquifers:

Temp, Temperature (°C)

DSL, Static Level depth (m)

ESL, Evolution Static Level (m)

F⁻, Fluoride (ppm)

TDS, Total dissolved solids (mg/L or ppm)

SO4-, Sulfates (mg/L or ppm)

EC, Electric Conductivity (μS/cm.)

Cl-, Chloride (mg/L or ppm)

CaCO3, Calcium Carbonate (total hardness) (mg/L or ppm)

Mn, Manganese (mg/L or ppm)

Fe, Iron (mg/L or ppm)

NO_2, Nitrites (mg/L or ppm)

NO₃ Nitrates (mg/L or ppm)

The results of the factorial analysis can be seen graphically, where the position and magnitude of the vectors of each parameter define the interaction in the groundwater flow system. The goodness of this analysis is that the parameters that have a relationship between them can be observed grouped, as shown in Table 1 and reported graphically in Figure 6 by the remarkable paper of Tóth (1963).

Figure 7 shows clearly in the graphical result of the multivariate factor loading, the dependence of water temperature on fluoride contents. It also gives an idea of the quality of the extraction water in ATA aquifer.

Figure 7 shows clearly in the graphical result of the multivariate factor loading, the dependence of water temperature on fluoride contents, due to their similarity in eigenvalues and eigenvectors, in addition to forming an acute angle between them. This graphical result is similar to that of the two aquifers studied and it was obtained from the water quality reports generated from 2003 to 2016.

The data were obtained from 60 and 40 wells in ACHS and ATA respectively, of the mentioned water quality parameters, as well as the depth levels of the groundwater, the temperature was recorded in the CONAGUA weather station. These data were processed annually and reviewed for their behavior in the observation period using MINITAB version 17 (MINITAB 2013).

In Table 2, the colors follow the strict order of statistically determined flow types: Yellow represents local flow (LF), orange, medium flow (MF) and in red, regional flow (RF).

All the analysis that follows is carried out through the calculations in Table 2 are strictly related to the results shown in Figure 7 through the Factorial Component Analysis

For more information on the multivariate analysis processes, see the book Multivariate Analysis of Hair et al. 1999, 5th Edition.

In the table above it is shown in summary how the parameters of water quality, temperature and hydrogeology were grouped (DSL and ESL). In the observation period in the monitored wells it can be seen that those with the greatest change in magnitude (percentage) are fluoride and temperature, in addition to chlorides and sulfates. In the particular case of ATA, the iron also presents significant variability. This grouping of water quality parameters are related to the depth of the static level and the evolution of the static level that occur annually. In the case of LF, the DSL monitored in the wells, in both aquifers, is less than or equal to 80 m and ESL on average is equal to or less than 1 m, with a water temperature of 24 °C.

In the case of MF, the depths of the DSL in the wells monitored in both aquifers are between 80 to 100 m and ESL, on average, is between 1 and 2 m, with a water temperature of 27 to 28 °C both ACHS and ATA.

In the case of RF, the depths of the DSL in the wells monitored in both aquifers are between 100 to 120 m and ESL on average is between 2 and more than 3 m, with a water temperature of 32 °C.

It is essential to exercise great caution in the analysis and interpretation of the results in order to obtain acceptable conclusions that reflect a clear picture of the performance of the flow systems also both ACHS and ATA respectively, because changing a single value in the individual steps of the analysis procedure may lead to different results and conclusions, especially with respect to the chemical and physical components of the flow systems.

After the statistical analysis, it is important to follow up by monitoring the hydrogeochemical characteristics of the water obtained in order to observe and analyze, in detail, the evaluation of the main regional flow indicators. It is good to note that all these indicators have been detected and proposed for study in Plan Estatal Hídrico del Estado de Chihuahua PEH 2040 (State Water Plan of the State of Chihuahua PEH 2040).

The procedures are shown using the hydrogeochemical parameters as witnesses in the hydraulic flow simulation; that is, considering only the groups of the wells that were identified by VISUAL MODFLOW (Waterloo Hydrogeologic 2012).

By following this methodology, the layout of the flow systems of the ACHS was estimated. Figure 8 shows the estimated rock types and flows.

After the statistical analysis, it is important to follow up by monitoring the hydrogeochemical characteristics of the water obtained in order to observe and analyze, in detail, the evaluation of the main indicators of regional flow. By following this methodology, the layout of the flow systems of the ACHS was estimated. Figure 8 shows the estimated rock types and flows.

Figure 9 shows the static levels in detail, based on the sampling carried out with support from the JMAS of Chihuahua and also confirming the flow behavior based on the geostatistical and hydrogeochemical analysis performed in this study.

During the study process, a horizontal (x, y) grid model was developed. This grid completely covers the aquifer surface, with square cells of 500 m per side, 120 columns and 155 rows. With this spacing, the following is achieved: good demarcation of the system geometry, fine regular parameterization and a good resolution of the spatial distribution of the simulated hydraulic loads. With this model, the interpolation error towards the observation points is reduced, as shown in Figure 10 for the ACHS aquifer.

The groundwater levels considered as the initial level in both ACHS and ATA aquifers (with key 0830 and 0835 respectively) are those reported in the 1971 study. These are complemented with hydrogeological criteria, such as topography support. An interpolation of these data by the Kriging method was also carried out so that the initial load values were spread and vary throughout the model. Colors represent piezometric elevation, as seen in Figures 11 for ACHS, and Figure 13 for ATA.

Since the most recent reports made by the Technical Sub-directorate General of Groundwater Management of CONAGUA have not presented new inputs since 2016 (Update of the Average Annual Water Availability in the Tabalaopa-Aldama Aquifer (0835), State of Chihuahua, Mexico, December 2020).

Figures 10–14 were interpolated to cover the entire aquifers zone.

ACHS aquifer comprises an estimated area of 1,889 km² (SINA-CONAGUA, 2020), and is considered a free to semi-confined aquifer type; it is adjacent to the east with the ATA aquifer.

It is located within the parallels 28° 26 ‘and 28° 56′ of North latitude, and 105° 58 ‘and 106° 32′ of West longitude, which positions it at the center of the state of Chihuahua.

The Tabalaopa-Aldama aquifer, defined with key 0835 in CONAGUA's Geographic Information System for Groundwater Management (SIGMAS), is located in the central portion of the State of Chihuahua, between parallels 28°32′30″ and 29°0′0″ north latitude and meridians 106° 10′ 0″ and 105° 52′ 0″ west longitude, covering an area of 728 km².

This aquifer is part of the Sierras y Llanuras del Norte physiographic province and the Bolsón de Mapimí subprovince. It presents its vertexes to the northwest, reaching the end of the vertex at Cerro El Colorado, to the northwest. From the north it crosses a mountainous area to the south crossing the Sierra El Cuervo-Peña Blanca, bordering the Sierra San Ignacio in the same direction from north to south until reaching the towns of Santa Eulalia and San Guillermo.

The ATA descends on the northwest side to the south, enters Cerro Chilicote and continues through the Nombre de Dios mountain range and Cerro el Coronel in a north-south direction. It passes through Cerro Grande and ends its apex in the valley, reaching the highest elevation on the north side and descending to the lower elevation south side.

This is verified by the fact that it runs with the El Mimbre stream and connects with two other streams: Los Nogales and Santa Eulalia, crossing the valley from north to south.

The proposed model better simulates the interaction among surface water and groundwater, and reduces the calculation error caused by the concentration of pumping in the center of the cells. In addition, the model was divided into four layers in vertical orientation: the first three in granular media and the last one in low conductivity fractured media. The three layers follow the subsurface geometry observed in the geological sections of Figure 13.

Figure 12 shows the cross section of the piezometric elevation levels obtained by means of the VISUAL MODFLOW program.

The groundwater levels considered as the initial level in the Tabalaopa-Aldama aquifer (ATA) are those reported in the 1971 study. These are complemented with hydrogeological criteria, such as topography support. An interpolation of these data by the Kriging method was also carried out so that the initial load values were spread and vary throughout the model. Colors represent piezometric elevation, as seen in Figure 14.

To determine the value and distribution of the input parameters to the model, information from the hydrogeological units, surface geology, pumping tests, land use and vegetation, census and groundwater balance were verified.

The Tabalaopa-Aldama aquifer is considered a no-flow aquifer, since there are no groundwater inflows or outflows at the aquifer boundaries, except for the border located in the northeastern portion of the aquifer corresponding to the Chuviscar River. The Chuviscar River flows mainly from southwest to northwest within the aquifer.

A total of 589 pumping wells were assigned to the model, totaling an extraction of 78.5 hm³ per year⁻¹, including 72 observation wells.

The groundwater levels considered as the initial level in the Tabalaopa-Aldama aquifer are those reported in the 1971 study. These are complemented with hydrogeological criteria, such as topography support. An interpolation of these data by the Kriging method was also carried out so that the initial load values are spread and vary throughout the model.

In the present study, 44 ACHS and ATA effort periods were assigned to the model. The calibration period spans from 1971 to 2011.

However, in both models, the simulations were run for the period from 2003 to 2016, considering in each period the simulations by grouping of the wells that were identified by the flow systems in the statistical analysis.

As a result, in both aquifers, the wells that were considered the LF had a better fit than in the other two (MF and RF). The simulation considering only the wells grouped in the RF presents a very large error (greater than 50%).

Table 3 shows a summary of the stress periods and extraction volumes between August 1971 and December 2011.

Through the study, it was possible to estimate the flow behavior of the Chihuahua-Sacramento and Tabalaopa-Aldama aquifers. The fluorine content values reported in the three average flow levels, mixed and regional, serve as a temporal and spatial guide to monitor the behavior of water quality and its temperature. This can tell us how to handle extractions, affecting the population's consumption as little as possible.

All relevant data are included in the paper or its Supplementary Information.