Atoll island hydrogeology: flow and freshwater occurrence in a tidally dominated system

Journal of Hydrology, 120 (1990) 327-340

327

Elsevier Science Publishers B.V., Amsterdam - - Printed in The Netherlands [2]

ATOLL ISLAND HYDROGEOLOGY: FLOW A N D F R E S H W A T E R OCCURRENCE IN A TIDALLY DOMINATED SYSTEM

J U N E A. OBERDORFER 1, PATRICK J. HOGAN 1 and ROBERT W. BUDDEMEIER 2.

1Department of Geology, San Jose State University, San Jose, CA 95192 (U.S.A.) ~Nuclear Chemistry Division, Lawrence Livermore National Laboratory, P.O. Box 808, Livermore, CA 94550 (U.S.A.) (Received December 1, 1989; accepted after revision J a n u a r y 25, 1990)

ABSTRACT Oberdorfer, J.A., Hogan, P.J. and Buddemeier, R.W., 1990. Atoll island hydrogeology: flow and freshwater occurrence in a tidally dominated system. J. Hydrol., 120: 327-340. A layered-aquifer model of groundwater occurrence in a n atoll island was tested with a solute-transport numerical model. The computer model used, SUTRA, incorporates density-dependent flow. This can be significant in freshwater-saltwater interactions associated with the freshwater lens of a n atoll island. Boundary conditions for the model included ocean and lagoon tidal variations. The model was calibrated to field data from Enjebi Island, Enewetak Atoll, and tested for sensitivity to a variety of parameters. This resulted in a hydraulic conductivity of 10today -1 for the surflcia] aquifer and 1000today -1 for the deeper aquifer; this combination of values gave an excellent reproduction of the tidal response data from test wells. The average salinity distribution was closely reproduced using a dispersivity of 0.02m. The computer simulation quantitatively supports the layered-aquifer model, including under conditions of density-dependent flow, and shows t h a t tidal variations are the predominant driving force for flow beneath the island. The oscillating, vertical flow produced by the tidal variations creates a n extensive mixing zone of brackish water. The layered-aquifer model with tidally driven flow is a significant improvement over the Ghyben-Herzberg-Dupuit model as it is conventionally applied to groundwater studies for many Pacific reef islands.

NOTATION C C*

Co D g / k p

Q~

Fluid solute mass fraction (mass solute per total mass) Solute concentration of fluid sources Base fluid solute concentration Apparent molecular diffusivity Dispersion tensor Gravitational acceleration Identity tensor Solid matrix permeability Fluid pressure Fluid mass source

[L2T -1] [L2T -1] [LT -2 ]

~] [ML-1T -2] [ML-ST -I]

*Present address: Kansas Geological Survey, 1930 Constant Ave., Campus West, University of Kansas, Lawrence, K S 66046-2598, U.S.A.

0022-1694/90/$03.50

© 1990 Elsevier Science Publishers B.V.

328

J.A. OBERDORFER ET AL.

v

Average fluid velocity Compressibility of solid matrix Longitudinal dispersivity Transverse dispersivity Compressibility of fluid Porosity Fluid viscosity Fluid density Base fluid density at C = Co Differential vector operator

~L ~r fl e p Po V

[LT ~] [LT2M 11 ILl [L] [LT 2M 1] [ML IT 2] [ML ~:~J [ M L :3] [L- ~1

INTRODUCTION

Reef islands provide the framework for the retention of fresh water in a sea of surrounding saline water. Movement of the fresh water and its mixing with the surrounding seawater affect the amount of potable water available and the diagenesis of the reef structure, particularly the evolution of porosity. Past hydrogeological assessments have frequently overestimated the amount of available water by neglecting the enormous amount of mixing that takes place. An improved understanding of the hydrodynamics of the flow regime would permit estimation of pore-water residence times, understanding of the development of different geochemical environments, and estimation of the rates at which diagenetic changes occur. A numerical model of the hydrogeology and solute transport of an atoll island was constructed and calibrated to field data (Hogan, 1988) from Enjebi Island,

a)

~j/-Enjebl

ol o_.L

km

Reef c r e S t ~

N

~ Lagoon

0 300 m

¢M|

elm 0o

~

!

• 18-80 m

140°E 180°140°W

Fig. 1. Location of this study: (a) Enewetak Atoll, with location of Enjebi Island indicated, and (b) Enjebi Island, with the locations of the test wells and the computer-simulation transect. After Herman et al. (1986). Open symbols indicate the wells used to calculate the average salinity profile.

ATOLL ISLAND

HYDROGEOLOGY

329

Enewetak Atoll (Fig. 1). This island was chosen because of a good set of field data (Buddemeier and Holladay, 1977; Wheatcraft and Buddemeier, 1981) and because it had been the subject of a previous modeling effort (Herman et al., 1986), which successfully simulated the tidal control of the flow patterns. As the island's groundwater is subject to significant tidal variations, the boundary conditions for both Herman's and this paper's model were chosen to reflect this effect. The model developed here was used to test hypotheses on geological control of the flow patterns in the island. In the future, the model can be used as a tool in evaluating water resource development options, in examining recovery from water-quality perturbations (such as those caused by storm-surge overwash), and in evaluating the consequences of global climatic change, particularly sea-level rise and rainfall variability, on these very low-lying islands. FRESHWATER LENS OCCURRENCE

Because of density differences, fresh groundwater in ocean islands occurs as a lens of fresh water floating atop the saline ocean water that permeates the porous geological substructure of the island. The traditional conceptual model of groundwater occurrence in an island has been the Ghyben-Herzberg lens (Fig. 2), which is generally treated mathematically in resource assessments with the Dupuit assumption of horizontal flow (Vacher, 1988; Hunt and Peterson, 1980). This conceptual model, although theoretically valid and probably appropriate for large, homogeneous islands, has traditionally been applied to small reef islands by use of a number of assumptions that appear to be inappropriate for that environment. First, outflow is. typically assumed to occur at the island margin, in or immediately below the intertidal zone. Flow is treated mathematically as one dimensional, which implicitly demands horizontal flow. Because of the relationship between freshwater and oceanwater denisities, the lens is considered to contain fresh water below the water table to a depth ~ 40 times the freshwater head (above mean sea l e v e l - - MSL). Tidal mixing is often considered negligible or small, and is generally treated as

discharge /

zone/

~:~l::::z:::::::~ ,...:.:..

:::i;:i:i:i!:ii:i!!:!:!ii~!!i:i ....

Rainfall

~

ilili!iii!!i!iiiii!!~~ iiiiiiiiiiiiiiiiiiiiiiiiiiiiiiiiiiii~

~

~

~ ~'h I

4nh

~

i*

~i!:i~======================== ........i

--

i~ilili~i~i i~i i ~i i i~ililili~il iiiiiiiii!!iiiiii

iiiiiii!iiiiiiiiiiiiiiiiiiii iiiiiiiiiiiiiiiiiiiiii iiiiiiiiiiiii!ii!iiiii!i!iii!i !!iiiiii!i iiliii!! !ii!iii!ii!i i !i!i!!ii93!i i ijiiij!iiijii3iijij;i

Fig. 2. Conceptual model of a n atoll-island Ghyben-Herzberg lens in a homogeneous medium, with Dupuit assumptions of horizontal flow.

330

J.AOBERDORFER ET AL.

significant only at the island margins. Tidal fluctuations of the water table are i expected to decrease rapidly with distance inland from the shoreline, whereas the lagtime for the tidal response correspondingly increases, Water movement within the freshwater lens is generally considered to result only from rechargeinduced head gradients, and to be independent of water movement in the underlying sea water, which is usually ignored. Field data from a number of coral-reef islands in the equatorial Pacific (Kwajalein: Hunt and Peterson, 1980; Enewetak: Wheatcraft and Buddemeier, 1981; Majuro: Hamlin and Anthony, 1987) contradict many of the above assumptions. A very broad transition zone (generally significantly greater than the potable water resource) is found to extend beneath the entire island. The tidal signal (both lag and efficiency) is found to be virtually independent of distance from shoreline, whereas it varies with depth so that stronger tidal signals occur with increasing depth. An alternative to the Ghyben-Herzberg-Dupuit conceptual model takes into account the heterogeneity of the geological materials and the recent geological history of the atolls, as well as the referenced field observations. Most atolls were subject to subaerial exposure of the carbonate reef structure during low stands of the sea that occurred repeatedly through the Pleistocene epoch. Dissolution of the exposed limestone by rainwater created a karst topography which, combined with the contructional variations evident in modern reef structures, gave the heterogeneous geological materials a very high and variably distributed permeability. Geological logging of boreholes on reefs and islands indicates that Pleistocene solution unconformities indentified with the most recent low sea levels typically occur 7-25m below sea level in most tectonically stable areas. Overlying these Pleistocene deposits are semi- to unconsolidated Holocene deposits of primarily sand and gravel-sized coral reef rubble. Where measurements have been made, this Holocene aquifer often has lower permeability (by one to two orders of magnitude) than the Pleistocene aquifer (Wheatcraft and Buddemeier, 1981; Oberdorfer and Buddemeier, 1985). Most of the fresh water is retained within these Holocene materials. In departure from the Ghyben-Herzberg-Dupuit model, the proposed alternative conceptual model (Buddemeier and Holladay, 1977; Wheatcraft and Buddemeier, 1981) consists of a two-layer, permeability-contrast system with a moderate-permeability Holocene aquifer overlying a high-permeability Pleistocene aquifer. Herman et al. (1986) previously modeled a layered-aquifer system of groundwater flow in an atoll island for Enjebi Island and found that the model produced good agreement with the tidal response data. The finiteelement simulation code FEMWATER (Yeh and Ward, 1980) was used. The model for the island's groundwater simulated flow only, however, and did not include either solute transport or the effects of density.induced flow on this variable-salinity system. The present model was constructed to include both salinity and density variations.

ATOLL ISLAND HYDROGEOLOGY

331

OBJECTIVES OF THIS STUDY The present simulation effort proposed to take modeling of atoll islands a step closer to reality than did the previous computer simulation. The effects of density variations on flow patterns would be assessed, particularly to evaluate the assumption of H e r m a n et al. (1986) that tidal driving forces dominate this system. This retesting of the layered-aquifer hypothesis with variable-density fluids would evaluate the ability of this conceptual model to reproduce hydrodynamic responses observed in the field.The creation of a calibrated model that included solute transport would provide a tool which could be used in predictive mode to track changes in water quality under different environmental conditions. M O D E L DESCRIPTION The U.S. Geological Survey computer code S U T R A (Voss, 1984) was selected for the model because it solves equations for both fluid and solute transport, including density-dependent flow. The numerical methods used to approximate these two interdependent processes are a two-dimensional, hybrid, finiteelement method and an integrated, finite-difference method. The code used includes the modifications as written by Voss for water-table conditions; these modifications permit the assignment of a value for specific storage to the uppermost nodes. Unsaturated transport was deemed to be of negligible importance in this system; therefore, that portion of the code was not used. Fluid pressure (p) is the primary variable in the flow equation, whereas the primary variable for the solute transport equation is solute concentration (C). Fluid density varies with concentration. Scientific terms are defined in the Appendix. Fluid flow is described by v =

-

[k.(e/0-1] .(vp

- pg)

(I)

Qp

(2)

fluid mass balance by

a(~p) at

=

-v

• (~#v)

+

solute mass balance by

a(epC) at

=

-v

• (epvC) q: V'[f.p(Dm/ + D)" VC] + QpC*

(3)

and equation of state for fluid density as a function of concentration by dp p(C) _-__ Po + ~-~(C -

Co)

(4)

The model has been configured (Fig. 3) to represent the geometry and geology observed at Enjebi, conceptualized as a layered aquifer system with a moderate-permeability aquifer to a depth of 12m below sea level, and an

332

J A . OBERDORFER ET AL. P nitih I IQhlnd

Bmk

Fig. 3. Conceptual model of a layered-aquifer system, indicating the boundaries used in the simulation. The consolidated aquifer consists of Pleistocene limestone, the surficial aquifer of Holocene material. underlying high-permeability aquifer to a total depth of 1277 m. The maximum elevation of the island is 3 m. Boundaries B1 and B2 are time-dependent pressure boundaries, where the fluid pressure varies with the tidal cycle. The tidal variations in sea level are represented by a sine wave with an amplitude of 1.8 m and a period of 12 h. The sine wave is expressed as pressure varying with time according to p(xB, zB, t) =

[0.9 sin (4~t) +

1277 -

zB]p~g

(5)

Boundary B3 is a no-flow boundary at the c arbonat e-basal t contact 1277 m below MSL, and B4 is a no-flow boundary located ~ 3000 m from the island so as to limit the extent of the grid. Boundary B5 is a recharge boundary; inflow at this boundary is estimated to be 0.5 m year -1, about one-third the annual precipitation of 1.5 m year-1 at Enewetak. Although rainfall at Enewetak is seasonal in nature, with most of the rain occurring from September to December, the simulated input was distributed equally over the year, to examine average conditions. This also had the effect of keeping computational times within manageable limits. The finite-element grid (Fig. 4) was set up to represent a cross-section through the island from ocean front to lagoon. A detailed description of the model configuration is given in Hogan (1988). This two-dimensional crosssection would be a good approximation for a strip island, where the length is much greater th an the width. For the triangular shape of Enjebi, however, edge effects might become significant. When this issue was addressed as part of the model study (see below), edge effects were found to be relatively unimportant. This two-dimensional simulation can therefore be taken as a close approximation of reality. The grid consists of 672 nodes and 605 elements, with greater element density in the Holocene aquifer. The geometry is very similar to t h a t used by Herman et al. (1986).

333

ATOLL ISLAND HYDROGEOLOGY

o I-

j j ~ j j ~ j j ~

M e a n sea level

)$1IIIIIllIllllIIIl111T D"IYlIIlIIIIlIIIIIIlllll] Lllll ~'JFIlIIIlllIIIllllIIlllIll~ ll'[I l l l L I I l l l l l l l l I l I I l l l l I IIIl "l I l l ' I I l ~Illllttilllllllllllllllllll [11 ii] ~lllltllttnlllllllllllllllllll dlllltttllI~tllllllllllllilllI1111TI "

5

_

.,J (n

:E

10

.9o J=

15

[ [ I I I I I I l l l l l l l l l l l l l l l l l l l l I I I I IJill

K 20

25

)lllllllllliiii'iiiiiiiiii 211110

2500

3000

3500

Distance from no-flow boundary B 4 (m)

Fig. 4. Finite-element grid in the vicinity of the simulated island. Boxed area identifies the zone contoured in Fig. 7.

Initial conditions for the simulation were a completely saltwater system with the pressures everywhere reflecting MSL. This is expressed as C(xi, zi, O) = 0.0357

(6)

and P(Xi, zi, O) =

(1277 -

(7)

zi)p,g

With a time step of 0.25 h, it required two simulated days for the pressures to reach a stable pattern of hydraulic response; 3 years of simulated time at a time step of I h were required for the salinity distribution to reach a stable configuration. Some input parameters for the model were taken from standard values in the literature; others were estimated from field data and then refined through sensitivity analysis. The input parameters determined from literature values are given in Table 1. TABLE 1 Input parameters determined from literature values Parameter

Value

Compressibility of water, fl Porosity, Fluid viscosity, Solute mass fraction, seawater, C~ Density, seawater, p~ Density, fresh water, pf

4.4 × 10 -1°m 2N -1 0.30 1.0 × 10 - a k g m - i s -1 0.0357 kg salt kg- 1 seawater 1025 kg m - 3 l O 0 0 k g m -3

334

J A O B E R D O R F E R E T AL.

TABLE 2 Input parameters determined from sensitivity analysis Parameter

Value

Upper aquifer permeability, ku Lower aquifer permeability, ki Compressibility of porous media, Longitudinal dispersivity, ~L Transverse dispersivity, ~w

1.2 x 10 t1 m~ 1.2 x 10 s m2 1.0 × 10 "~m~ N J 0.02 m 0.0

P a r a m e t e r s d e t e r m i n e d t h r o u g h s e n s i t i v i t y a n a l y s i s (Table 2) w e r e c a l i b r a t e d to the tidal r e s p o n s e o b s e r v e d in wells at Enjebi and to the o b s e r v e d salinity distribution. T h e m e a n tidal efficiency (the r a t i o of the w a t e r - l e v e l c h a n g e in the well to a c o r r e s p o n d i n g c h a n g e in the o c e a n tide) in nine shallow-pit wells was 0.09 +_ 0.04, and the m e a n tidal lag (the delay in t i m e b e t w e e n p e a k tides in the o c e a n a n d well) was 2.96 _+ 0.29 h (Buddemeier, 1981). T h e s e v a l u e s w e r e t a k e n to r e p r e s e n t the tidal r e s p o n s e at the w a t e r table. Field d a t a also s h o w e d t h a t tidal efficiency i n c r e a s e d a n d tidal lag d e c r e a s e d w i t h d e p t h b e n e a t h t h e island to w h e r e tidal efficiencies n e a r the H o l o c e n e P l e i s t o c e n e c o n t a c t were in excess of 0.40 and lags were as low as 0.25h ( W h e a t c r a f t a n d Buddemeier, 1981). One difficulty w i t h t h e s e field d a t a is t h a t the b o r e h o l e casings w e r e slotted o v e r t h e i r e n t i r e lengths. As the tidal signal was i n t e g r a t e d o v e r the t o t a l depth, the r e s u l t s are m o r e difficult to i n t e r p r e t . C o m p a r a b l e v a l u e s h a v e been found for o t h e r atoll islands (Bikini, K w a j a l e i n , M a j u r o ) w h e r e s h o r t - s c r e e n p i e z o m e t e r s w e r e used. T h e s a l i n i t y d i s t r i b u t i o n w i t h i n the island v a r i e s w i t h time, b e c a u s e f r e s h w a t e r i n p u t tends to be seasonal. To c r e a t e a g e n e r a l i z e d view of the lens c o n f i g u r a t i o n t h a t is c o n s i s t e n t w i t h o u r use of a n a v e r a g e r e c h a r g e v a l u e in modeling, all a v a i l a b l e s a l i n i t y d a t a (from s a l i n i t y profiles m e a s u r e d f r o m s u r f a c e to f u l l - s e a w a t e r s a l i n i t y in n i n e deep wells at v a r i o u s s e a s o n s o v e r a period of 2 years) w e r e a v e r a g e d . This ' a v e r a g e ' s a l i n i t y d i s t r i b u t i o n w a s t a k e n to r e p r e s e n t w h a t would be likely to o c c u r if the t o t a l r e c h a r g e were d i s t r i b u t e d e v e n l y t h r o u g h o u t the year. RESULTS A r a n g e of v a l u e s was tested for the h y d r a u l i c p a r a m e t e r s of g r e a t e s t uncertainty, to d e t e r m i n e those p a r a m e t e r s w h i c h best m a t c h e d the field observations. E a c h of the two f o r m a t i o n s was t r e a t e d i n d i v i d u a l l y as b e i n g homogen e o u s a n d isotropic. T h e h y d r a u l i c c o n d u c t i v i t y of the surficial aquifer, as d e t e r m i n e d by a q u i f e r tests, h a s been found to r a n g e from 1 to 100m day-1 ( c o r r e s p o n d i n g to a n i n t r i n s i c p e r m e a b i l i t y of 1.2 x 10--12 to 1.2 × 10 10 m2). A r a n g e of 5-50 m d a y 1 was tested in t h e s e n s i t i v i t y analysis. D i r e c t o b s e r v a t i o n of the P l e i s t o c e n e a q u i f e r c h a r a c t e r i s t i c s is v e r y limited; however, its per-

ATOLLISLANDHYDROGEOLOGY 1.0!

,

I

'

I

'

335 I

0.9 0 8r ~ , , , . ~

'

I

q

12

'

o10 o 15 20

0.6,

01

1

Upper aquifer hydraulic c o n d u ~ h ~ i t y (m/d)

•

"5

'

• 25

,

I

10

,

1 8 Depth

,

1

6 below

,

I

4

,

I

,

2

MSL (m)

Fig. 5. Simulated tidal efficiency vs. depth as a function of upper (Holocene) aquifer hydraulic conductivity. Lower (Pleistocene) acquifer hydraulic conductivity was fixed at 1000 m day -z.

meability is thought to exceed t h a t of the Holocene by one to two orders of magnitude. The sensitivity analysis represented in Fig. 5 indicates that the values of 10 m day-1 (intrinsic permeability of 1.2 × 10-1~ m 2) for the Holocene sediments and 1000m day -~ (intrinsic permeability of 1.2 × 10 -9 m 2) give the closest match to the tidal efficiency. These values give a tidal efficiency of 0.14 at the water table and of 0.70 at the Holocene-Pleistocene contact. Corresponding tidal lags were 2.75 h at the water table and 0.25 h at the Holocene-Pleistocene contact. The tidal response produced by the model using these values is considered to be in excellent agreement with the field data. The elevation of the water table above MSL, 0.26m, predicted with these hydraulic conductivity values agrees with the average water-table elevation observed in nine shallow Enjebi wells, 0.26 + 0.05m (Buddemeier, 1981). The rate of change of fluid mass is a function of the compressibility of the water (fl) and the compressibility of the porous medium (~). A value of 1 × 10 -9 m 2 N -1 for the porous medium compressibility was found to produce a reasonable response. Values of 1 × 10 -8 or less excessively dampened the tidal signal. The value used is in the midrange of literature values. The dispersion coefficient is the product of the dispersivity of the porous medium and the fluid velocity. The longitudinal dispersivity value is a fitted parameter, that which best reproduced the observed mixing of fresh water and salt water beneath the island. The value of longitudinal dispersivity obtained (0.02 m) is certainly within the range of reasonable values. The same value was used for both aquifers. The transverse dispersivity was not used for this simulation because the large variations in flow directions and magnitudes at a given point would cause the longitudinal dispersivity to mask the effects of the much smaller transverse dispersivity. No upstream weighting was used. According to Voss (1984, pp. 231-233) a dispersivity of the order of one-quarter

336

J,A. OBERDORFER ET AL. High tide

~ ~$~*$

I

~lO-5m/s

h

A

• A AAAA A

•

lO-sm/s lo-Tm/s ~

LOW tide

.

Fig. 6. Velocity v e c t o r s in t h e s i m u l a t e d region: (a) h i g h tide a n d (b) low tide. A r r o w l e n g t h s are a p p r o x i m a t e l y p r o p o r t i o n a l to log (velocity).

the mesh size is needed to obtain a stable solution without upstream weighting. Our vertical mesh spacing was 1 m, which implies a need for longitudinal dispersivity values near 0.25m. In fact, we simulated periods of>~3years without instability problems while using a longitudinal dispersivity of 0.02 m. We hypothesize that the strong vertical oscillations in flow prevented cumulative instabilities from developing, but have not tested this hypothesis. A longitudinal dispersivity of 0.05 m was tested and was found to overmix the water (i.e. surface salinities were too high to match the averaged field data). SUTRA calculates flow velocities (v) as part of its output; these were plotted using the graphics package DISPLAA and are presented in Fig. 6. The very strong effects of the tidal variations on the flow field can be seen in this figure; a complete reversal in flow direction occurs between high and low tide. Flow /--

MSL= O - -

Groundwater

-

_o '°

a 15 251111

3000

3500

Distance from B4 model boundary (m)

Fig. 7. S a l i n i t y c o n t o u r s in per c e n t s e a w a t e r for e q u i l i b r a t e d c o n s t a n t - r e c h a r g e s i m u l a t i o n .

ATOLLISLANDHYDROGEOLOGY

337

-

Simulation

\\

_

_

_0 J: m -

_,o

-

10

r-

12

14

16

18

0

20

40

60

80

100

Percent Seawater

Fig. 8. Average salinity profiles of Enjebi groundwater from field data and from SUTRA simulation.

is seen to be vertically upward beneath the island surface on the high tide and vertically downward on the low tide. The flow field is nearly uniform across the width of the island, with only small edge effects showing in the immediate vicinity of the island's margin. The computer-generated salinity distribution within the island's groundwater is given in Fig. 7. A comparison between simulation results and the average field measurements can be seen in Fig. 8. Because of strong seasonal and annual variations in the upper portions of the field curve (averages had standard deviations approaching 100%), the fit between the simulated and the measured data is satisfactory in view of the uncertainties and assumptions involved. The isochlors (percentage of sea water) can be seen to be relatively flat beneath the center of the island; this reflects the vertical nature of transport and mixing in this zone. Salinity increases much more rapidly in the horizontal direction near the margin of the island than it does in the vertical

338

•.A. OBERDORFER ET AL.

direction. The 50% isochlor occurs at ~ 9m, which is in reasonable agreement with the Ghyben-Herzberg prediction of a freshwater depth of 10.4 m based on a head of 0.26 m. DISCUSSION The model of Enjebi presented here has been simplified from the complex original in many ways: heterogeneities within each aquifer (such as the lowerpermeability reef plate on the ocean side of the surficial aquifer) have been ignored, as have the details of the mixed semi-diurnal tide signal and the water level differences between ocean and lagoon observed to exist for part of the tidal cycle (Buddemeier, 1981). Such refinements to the model must be tested in the future to evaluate the magnitude of their effect; however, the level of complexity incorporated into the current modeling attempt is close enough to reality to provide an improved understanding of the groundwater system. One important improvement to be incorporated in the future will be to examine seasonal variability by varying aquifer recharge during different periods of the year. This is needed to assess real responses to environmental changes and to establish a realistic lower limit on available groundwater resources. The ability of the layered-aquifer model to reproduce the hydraulic response in wells, the elevation of the water table, and the average distributions of salinity reinforces its accuracy in describing atoll island hydrogeology. It must be remembered that the two-layer model used here is still just an approximation for a much more complex reality, but it appears to be a much more realistic one than the model of a Ghyben-Herzberg lens with Dupuit assumptions. The sensitivity analysis has provided insight into the probable magnitude of hydraulic conductivity of the Pleistocene aquifer, an as-yet undetermined parameter. Although the hydraulic conductivity of the Holocene aquifer is known reasonably well from field tests to be of the order of tens of metres per day, the model showed that a hydraulic conductivity contrast of approximately two orders of magnitude is needed to reproduce the field data for tidal response in wells. The findings for the flow dynamics are very similar to those of Herman et al. (1986), who found t h a t tidal response of the island's water table was virtually independent of distance from the shoreline. This is because the tidal signal travels very quickly and with little attenuation through the highpermeability Pleistocene aquifer, which then acts as the source (over the relatively short and consistent vertical distance of 15 m or less) of the tidal signal to the water table through the overlying, lower-permeability Holocene aquifer. The surficial aquifer acts as a leaky, confining layer where pressure built up by the rising tide in the lower aquifer is relieved by flow upward into the overlying Holocene deposits. As hypothesized by Herman et al. (1986), the predominant driving force for flow is the ocean's tidal variation, which masks any density-driven component. The traditional Ghyben-Herzberg-Dupuit analysis treats horizontal flow as driven by recharge-generated head gradients, This mechanism can occur in the layered-aquifer model as well, but there are also other possible driving forces

ATOLL ISLAND HYDROGEOLOGY

339

for horizontal flow. Dispersion and horizontal mixing as a result of tidal pumping can result in enhanced loss of fresh water by the discharge of brackish water to the surrounding ocean, and the asymmetry of the reef and island structure may interact with the tidal oscillations to generate net gradients (e.g. compare the relative changes in velocity vectors on the lagoon and ocean sides of Fig. 6). Although neither of these mechanisms was specifically investigated in the present study, it is interesting to note t h a t the simulated water table was relatively flat over most of the central part of the island - - an observation that is consistent with field data and that suggests that recharge-generated head gradients may not completely control the horizontal flow of fresh water. Flow in the lens is not separate from flow in the underlying seawater. To the extent that the seawater experiences flow t h a t is driven by forces other than the freshwater recharge (Buddemeier, 1981; Buddemeier and Oberdorfer, 1989), the Dupuit assumptions as conventionally applied to the Ghyben-Herzberg model may not be valid. The broad transition zone t h a t exists beneath the island is produced by the oscillating flow velocities that cause mixing between the fresh water and salt water on each tidal cycle. Ignoring tidally driven flow would mean ignoring this mixing mechanism and could lead to gross overestimation of the amount of potable water ( < 3% sea water) t h a t might be available. Solute transport models of atoll islands that do not include tidally induced flow would have to use unreasonably large dispersivities to reproduce the observed broad transition zone. This could lead to unrealistic responses when the model is later used in a predictive mode. The ability of the model presented here to simulate field observations gives confidence in its use as a predictive tool. There are certainly disadvantages (i.e. long computer run times) to having to use a time step as small as 1 h, but this time step is necessary to reproduce the groundwater hydrodynamics and the resulting mixing phenomenon. The model has been used in a predictive mode to examine the effects of global climatic change on atoll islands (Oberdorfer and Buddemeier, 1989), and further simulations and model refinements are planned for the future. If the hydrodynamics modeled here are representative of those of other atoll islands, as appears likely, this model can be generalized to represent a larger number of atoll islands. CONCLUSIONS Comparison of the model simulation results with the field observations produces the following conclusions: (1) the quantitative agreement of model results with both physical and chemical field observations further supports the layered-aquifer conceptual model, including under conditions of variable salinity where density variations could induce flow; (2) tidally induced flow predominates over density-driven flow and produces vertical flow patterns and mixing of fresh water and salt water beneath the

340

J A . OBERDORFER ET AL

island surface, with the p r i m a r y p a t h w a y for tidal signal p r o p a g a t i o n being t h r o u g h the lower aquifer; (3) the use of G h y b e n - H e r z b e r g - D u p u i t a s s u m p t i o n s w i t h o u t a d e q u a t e cons i d e r a t i o n of tidally i n d u c e d m i x i n g m a y result in a significant o v e r e s t i m a t e of reliable reef island g r o u n d w a t e r r e s o u r c e s in s i t u a t i o n s where a highpermeability, o c e a n - c o n n e c t e d aquifer underlies the w a t e r - b e a r i n g zone; and (4) the layered-aquifer model c a n be applied to reef islands as a useful and r e l e v a n t predictive tool for b o t h w a t e r r e s o u r c e and g e o c h e m i c a l analyses. ACKNOWLEDGEMENTS We t h a n k Dr. H. L. V a c h e r a n d an a n o n y m o u s reviewer for c o m m e n t s t h a t significantly improved this m a n u s c r i p t . P o r t i o n s of this w o r k were performed u n d e r the auspices of the U.S. D e p a r t m e n t of E n e r g y by L a w r e n c e L i v e r m o r e N a t i o n a l L a b o r a t o r y u n d e r C o n t r a c t W-7405-Eng-48. REFERENCES

Buddemeier, R.W., 1981. The geohydrology of Enewetak Atoll islands and reefs. Proc. 4th Int. Coral Reef Syrup., 1: 339-345. Buddemeier, R.W. and Holladay, G.L., 1977. Atoll hydrology: island groundwater characteristics and their relationship to diagenesis. Proc. 3rd Int. Coral Reef Symp., 2: 167-174. Buddemeier, R.W. and Oberdorfer, J.A., 1989. Hydrology and hydrodynamics of coral reef pore waters. Proc. 6th Int. Coral Reef Symp., 2: 485-490. Hamlin, S.N. and Anthony, S.S., 1987. Groundwater resources of the Laura Area, Majuro Atoll, Marshall Islands. Water Resources Investigations Rep. 87-4047, U.S. Geological Survey, Honolulu, HI, 69 pp. Herman, M.E., Buddemeier, R.W. and Wheatcraft, S.W., 1986. A layered aquifer model of atoll island hydrology: validation of a computer simulation. J. Hydrol., 84: 303-322. Hogan, P.J.,1988. Modeling of fresh water-saltwater interaction on Enjebi Island, Enewetak Atoll, M.S. Thesis, San Jose State University, San Jose, CA, 141 pp. Hunt, C.D. and Peterson, F.L., 1980. Groundwater resources of Kwajalein Island, Marshall Islands. Water Resources Research Center Technical Rep. 126, University of Hawaii, Honolulu, 91 pp. Oberdorfer, J.A. and Buddemeier, R.W., 1985. Coral reef hydrogeology. Proc. 5th Int. Coral Reef Syrup., 3: 307-312. Oberdorfer, J.A. and Buddemeier, R.W., 1989. Climate change: effectson island resources. Proc. 6th Int. Coral Reef Syrup., 3: 523-528. Vacher, H.L., 1988. Dupuit-Ghyben-Herzberg analysis of strip-islandlenses. Geol. Soc. Am. Bull., 100: 580-591. Voss, C.I., 1984. SUTRA: A finite-element simulation model for saturated-unsaturated, fluiddensity-dependent groundwater flow with energy transport or chemically-reactive singlespecies solute transport. U.S.Geological Survey Water Resour. Invest. Rep., 84-4369:409 pp. Wheatcraft, S.W. and Buddemeier, R.W., 1981. Atoll island hydrology. Ground Water, 19: 311-320. Yeh, G.T. and Ward, D.S., 1980. F E M W A T E R : a finite element model of water flow through saturated-unsaturated porous media. Rep. ORNL-5567. Oak Ridge National Laboratory, Oak Ridge, TN, 117 pp.