EVALUATION AND OPERATION POLICY FOR AQUIFERS IN FOLDED LIME STONES

Instituto de Ingeniería de la Universidad Nacional Autónoma de México

Ciudad Universitaria, México D.F., 04510, MéXICO

Tel. (52) 56 223326; Fax (52) 56 162798; E-mail ccv@pumas.iingen.unam.mx

Abstract: A model for aquifers in folded limestone in highlands based on simple linear reservoir concepts is developed; this is possible due to the high permeability of the formations where the piezometric levels evolve similarly even for large changes. The model was applied to several structures and explains, in the mean, water level evolution in response to rain and pumping. From its analysis a relation between pumping depth and sustainable discharge is obtained, thus giving a realistic view of water cost.

 

Keywords: aquifer evaluation, folded lime stones, karst, linear reservoir

The chain of mountains in the east of the Mexican Republic called “Sierra Madre Oriental” is the product of a large and very thick sequence of intercalated maritime and continental sediment deposits formed in the Mesozoic era, at the end of the Jurassic and during the Cretassic periods. It is clear that the intercalation took place in alternating sinking and emerging phases of the continent. In the early tertiary age these sediments suffered strong compression stresses which folded and raised them to elevations of thousands of meters over sea level. Later these great rock masses suffered tension stresses with the resulting faults, horsts and graven, the erosion of the soft and elevated sediments with their further deposition at lower levels.

Maritime sediments deposited in shallow waters formed mostly permeable limestone’s with dissolution channels (karst) originating in the fossil nodules while deposition in deep waters formed mostly compact and sometimes clayey limestone of low permeability. The alternation of permeable and impermeable strata and its subsequent folding and faulting facilitated the existence of large water deposits in permeable zones contained by impermeable ones in syncline structures (Fig. 1).

The northeast part of Mexico is an arid region with very variable extreme temperatures (registered –13℃ in winter and 39℃ in summer). Mean annual rainfall goes from 200 mm in low continental valleys to 1000 mm in highlands nearer to the Gulf of Mexico. Most populated centres are in the lowlands while highlands are dedicated to weather dependent agriculture with some irrigation by wells and small impoundment. Water demand for human consumption is high in the cities and the resource scarce; this is why nearly 40 years ago the perforation of very deep wells began, searching to exploit the referred permeable limestone. Results were in many cases favourable with good quality water and high transmissivity of the formation. Nevertheless, piezometric levels were very sensitive to prolonged exploitation and to heavy rainy seasons, which indicated low storage capacity.

In this study it was possible to access registers of monthly observations of piezometric levels and extracted volumes of several wells whose water is used to provide potable water to the city of Saltillo with 700,000 inhabitants. The largest record has 17 years with some missing intervals. The objective was to determine the sustainable rate of water extraction for a well field located in a specific structure (bottom or side of a syncline).

The evolution of registered water levels was first examined. It was noted that wells belonging to the same structure have very similar, parallel evolutions even when they are apart as much as 10 km; instead, when the wells are in different synclines, even being very close to each other, their water level evolutions could not be correlated. This suggests that each structure conforms a reservoir which is fairly independent from adjacent ones since they are separated by impermeable layers. These reservoirs must receive their recharge from rain and surface flow over the surface of permeable layers and have their natural discharge as flow to other formations and also as springs. The nature of their recharge must be very variable according to the rain regime while their discharge shall change fairly slowly. In natural conditions and in the long run discharge should equal recharge. Pumping from the reservoir introduces a perturbation whose nature will be analysed in what fallows.The hypothesis was made that these reservoirs could be treated as linear ones, that is to say that their natural discharge is proportional to their storage.

A quantitative expression for a linear reservoir is given by the continuity equation which for a given time interval is:

(V_i - V_i-1)/D t = I_i – a (V_i + V_i-1)/2 – B_i                                   (1)

Where:

V_i - Volume stored at the end of the interval i

I_i - Volume of water infiltrated from time i-1 to time i (Infitration rate)

a - Constant that defines the natural discharge as proportional to storage.[T^-1]

B_i - Volume of water pumped from time i-1 to time i. (Pumping rate)

D t - Time interval

If one considers the subsurface reservoir having a constant plan area and a constant storage coefficient the stored water is expressed as:

V = H A S                                                 (2)

Where :

A - Plan area of the reservoir

S - Storage coefficient of the permeable formation

H - Piezometric elevation over a base level

On the other hand, infiltration should be related to the amount of rain falling upon the surface of the permeable formation but also to the antecedent humidity of that surface. It has been shown in previous studies [1] that surface runoff and infiltration are a non linear function of precipitation that may be expressed as:

I_i = (b P_i + c P_i-1 P_i ) A_c                                           (3)

Where:

b - Linear infiltration factor for the rain from time i-1 to time i

P_i - Rain height from time i-1 to time i [L T^-1]

c - Factor that multiplies previous rain to form the antecedent humidity index. [L^-1 T].

A_c - Capturing area of the permeable surface of the formation. [L²]

Substituting equations (2) and (3) in (1) one obtains:

H_i A S (1/D t+a/2) = H_i-1 A S (1/D t-a/2) + A_c b P_i + A_c c P_i-1 P_i – B_i                 (4)

which divided by A S (1/D t +a/2) gives:

H_i = d₁ H_i-1 + d₂ P_i + d₃ P_i-1 P_i – d₄ B_i                                   (5)

Where:

d₁ = (1/D t -a/2) / (1/D t +a/2) dimensionless;

d₂ = A_c b / (A S (1/D t +a/2)) [T];

d₃ = A_c c / (A S (1/D t +a/2)) [T² L^-1];

d₄ = 1 / (A S (1/D t +a/2)) [T L^-2 ];

To obtain the parameters d_i one depends highly on having trustworthy observations of pumped volumes and of water levels over the largest time interval possible that includes rainy and dry years. For the case studied these variables are measured in general monthly, so that the chosen time interval in the above equations was one month. It should be beard in mind that parameter “a” is linked to the time interval, its value being actually 2(1 – d₁ ) / ( 1 + d₁ ) D t .

If an error term is added to equation 5, there is a procedure to determine coefficients d_i in that equation; that is the dynamic least squares fitting by the Kalman filter which also allows to define tendencies of the parameters if they are not constant [1], [2]. This procedure was applied to several syncline structures of the region under study obtaining very similar results in all of them: d₁ is smaller but very close to unity and fairly constant with variations of (1-d₁) of ± 10%; d₂ is variable but always around zero with insignificant values (10^-5); d₃ is very variable around a mean value with fluctuations of ± 300%; d₄ is constant with fluctuations of only ± 0.3%.

Two conclusions may be drawn from the described behaviour of parameters d_i: one is that the absence of tendencies (specially the very low variation of d₄) during the observation time with water level variations of around 50 meters, validates the hypothesis of a constant product AS; the second one is that the great variability of the parameter that multiplies the non linear precipitation term d₃, avoids the precise determination of its value with the procedure employed. The variability can be blamed to the fact that precipitation is not measured directly over the formation and probably also to the influence of the type of rain that occurs during a month, its distribution in time, its intensity, etc.. This has to be further investigated. With these facts in mind an alternative approach was used to obtain the d_I parameters.

As was already pointed out, in the absence of pumping and in the long run, inflow to the subsurface reservoir must be in equilibrium with its outflow; the latter is determined by aV = aASH; if a mean initial value of H is taken, the one before exploitation begins, and is named H_in, the integral of the inflow volume over a long time span is calculated as NaASH_in where N is the number of time intervals considered (in the present case, months). On the other side, the inflow integrated over the same time span according to the model may be obtained by cA_c S P_i-1P_i; but cA_c from the definition of d₃ is equal to d₃AS(1/D t+a/2); substituting this value and equating integrated outflow and inflow one gets:

d₃ = a N H_in / [(1/D t +a/2) S P_i-1 P_i ]                          (6)

When pumping, one can also integrate equation (1) for the whole known time interval; the total storage change is:

AS (H_fin – H_in) = cA_{c D} t S P_i-1P_i – a D t AS S H_i - D t S B_i

Where:

H_fin - Piezometric level in month N at the end of the observed pumping time.

Again, substituting here the recharge in the long run by the initial natural discharge one obtains:

AS = S B_i / [ a N H_in – a S H_i – ( H_fin – H_in )/ D t ]                            (7)

From equations (6) and (7) d₃ and AS depend only of observations and of parameter “a”, so that another condition has to be imposed to the solution and that is to have the minimum square error between observed and predicted water elevations. A direct formulation of this condition was obtained but it did not resulted in a reliable value for “a” (probably due to approximation errors). So, a search method was preferred; the parameter was given different values to obtain d₃ and AS for simulations until the minimum square error was obtained.

The model itself indicates that no sustainable pumping is possible unless the reservoir natural discharge diminishes and that is precisely accomplished by pumping as the piezometric level is lowered. Thus, the extraction policy is very clear: one has to over exploit the reservoir until a desirable pumping level H_des is reached and then diminish extraction to a rate that equals the rate of discharge that has been avoided to happen naturally; the sustainable pumping rate is:

B_sus = a AS ( H_in – H_des)                                    (8)

The described evaluation method was applied to five well fields in structures not related to one another, even when some were close together.

The best well field is located on the northern slope of a mountain ridge named Zapalinamé where 12 wells were perforated in a length of 10 km. Their water levels show a definite parallel behaviour in spite of their distance and of having experienced a draw down of about 85 meters. Another good well field is that of Jagüey defile with only three wells which receives flow contributions from the two adjacent structures.

The other three fields are located at the nose of ridges that deepen into valley sediments. Their evolution although similar is not of the same magnitude, indicating a less uniform structure than the former.

The main characteristics of each well field are given in table 1. There, the last two columns refer to the relative mean quadratic error in fitting calculated to observed water level evolutions (MQE) and to the head draw dawn necessary to obtain 100 l/s of sustainable pumping. This last value shows the great difference among aquifers.

In figure 2 a graph is drawn to show the evolution of water levels at a well in Zapalinamé structure, observed an calculated. At the end of the observation time interval, pumping was diminished to the sustainable rate and it may be seen that predicted water levels fluctuate around a mean value when the model is fed with previous rains.

The developed model for karst aquifers in folded lime stone allows the evaluation of its recharge, response to pumping and infiltration and is the basis for its operation policy. Based on the concept of a linear reservoir and the continuity equation the model has parameters with physical significance. They can be determined by integrating the continuity equation for natural conditions and for the pumping time intervals. The latter requires the observation of water levels evolution and of extracted volumes of water. The sustainable exploitation volume equals the amount of natural discharge that is avoided through lowering piezometric levels. This will affect naturally flows to other formations and springs.

The philosophy that may be drawn from this experience is the same that applies to other more common aquifers like those in unconsolidated sediments. In any case they are natural reservoirs that receive an inflow and have naturally an outflow to other water bodies. Pumping for local use will always affect resources downstream.

It is recognised the support and the facilities to access their information to the Gerencia Regional Río Bravo of the Comisión Nacional del Agua, in the consecution of this study.

[1] Cruickshank C. and Verde C. (1996), “Some Examples of Simple Rainfall-Runoff Models with Dynamic Parameter Identification”, XI Intern. Conf. on Comput. Meths. in Water Resour., V2, pp 445-452, Cancún, México.

[2] Bras R.L. and Rodríguez-Iturbe I. (1986) Random functions in hydrology. Academic Press

               Table 1  Main characteristics of the analysed structures

Fig. 1  Example of the folding in the region

Fig.2  Water level evolution at a well in Zapalinamé field