Applying Downscaled Global Climate Model Data to a Groundwater Model of the Suwannee River Basin , Florida , USA

The application of Global Climate Model (GCM) output to a hydrologic model allows for comparisons between simulated recent and future conditions and provides insight into the dynamics of hydrology as it may be affected by climate change. A previously developed numerical model of the Suwannee River Basin, Florida, USA, was modified and calibrated to represent transient conditions. A simulation of recent conditions was developed for the 372-month period 1970-2000 and was compared with a simulation of future conditions for a similar-length period 2039-2069, which uses downscaled GCM data. The MODFLOW groundwater-simulation code was used in both of these simulations, and two different MODFLOW boundary condition “packages” (River and Streamflow-Routing Packages) were used to represent interactions between surface-water and groundwater features. The hydrologic fluxes between the atmosphere and landscape for the simulation of future conditions were developed from dynamically downscaled precipitation and evapotranspiration (ET) data generated by the Community Climate System Model (CCSM). The downscaled precipitation data were interpolated for the Suwannee River model grid, and the downscaled ET data were used to develop potential ET and were interpolated to the grid. The future period has higher simulated rainfall (10.8 percent) and ET (4.5 percent) than the recent period. The higher future rainfall causes simulated groundwater levels to rise in areas where they are deep and have little ET in either the recent or future case. However, in areas where groundwater levels were originally near the surface, the greater future ET causes groundwater levels to become lower despite the higher projected rainfall. The general implication is that unsaturated zone depth could be more spatially uniform in the future and vegetation that requires a range of conditions (substantially wetter or drier than average) could be detrimentally affected. This vegetation would include wetland species, especially in areas inland from the coast. How to cite this paper: Swain, E. and Davis, J.H. (2016) Applying Downscaled Global Climate Model Data to a Groundwater Model of the Suwannee River Basin, Florida, USA. American Journal of Climate Change, 5, 526-557. http://dx.doi.org/10.4236/ajcc.2016.54037 Received: February 23, 2016 Accepted: November 15, 2016 Published: November 18, 2016 Copyright © 2016 by authors and Scientific Research Publishing Inc. This work is licensed under the Creative Commons Attribution International License (CC BY 4.0). http://creativecommons.org/licenses/by/4.0/ Open Access E. Swain, J. H. Davis 527


Introduction
Global Climate Models (GCMs) are important tools for simulating historical climate and projecting future climate, including precipitation [1].The precipitation predicted by some GCMs has been "downscaled"; that is, converted to a finer resolution.The downscaling method can be statistical or dynamic, and both methods have varying degrees of uncertainty.Comparisons of several models over the conterminous United States indicate that root-mean square errors of precipitation predictions differ by less than 0.1 mm/day between statistical and dynamic methods [2].Analyses of uncertainty in statistical downscaling methods indicate significant variations between stochastic and regression-based techniques [3].Dynamic downscaling involves embedding a smaller-scale regional climate model within the GCM [4].This approach resolves atmospheric processes on a smaller scale and with physically consistent processes, but is computationally intensive and sensitive to uncertainties in the GCM-derived boundary conditions.
The Florida State University Center for Ocean-Atmospheric Prediction Studies (COAPS) has dynamically downscaled GCM-simulated precipitation for the southeastern United States [5].The downscaled precipitation rates are bias-corrected by using the quantile-matching approach [6], which assigns corrections by percentile of the precipitation's cumulative distribution function.COAPS has downscaled other GCM variables including evapotranspiration (ET).Precipitation and ET are both essential to hydrologic simulations, and a predictive model can use the downscaled GCM values to represent future conditions [7] [8] [9].The downscaled data can be used to compare recent and future simulations to estimate changes in hydrology attributable to climate change.
The Suwannee River Basin, occupying about 9950 mi 2 in north-central Florida and southern Georgia (Figure 1), is an area where the effects of future precipitation changes are of concern.Forested and agricultural lands account for much of the current land use.Parts of the basin are subject to future population growth and increases in groundwater withdrawals [10].These factors increase concerns about the effects of precipitation changes and water availability in the future.The Suwannee River, originating in the Okefenokee Swamp in southeastern Georgia (Figure 1) and terminating at the Gulf of Mexico, 12 mi northwest of Cedar Key, is a dominant surface-water feature in the basin (Figure 1).The Alapaha and northern Withlacoochee Rivers are major tributaries to the Suwannee River.The other major tributary is the Santa Fe River, which flows westward from its headwaters to join the Suwannee near Branford, Florida (Figure 1).The lower reaches of these tributaries are incised into the highly transmissive Upper Floridan aquifer, as is most of the Suwannee River (the reach downstream from White Springs), and substantial interactions occur between the aquifer and these river reaches (and the numerous springs occurring along these reaches).To support analyses of the groundwater system in the Suwannee River Basin, a numerical model was constructed using the MODFLOW code that simulates groundwater flow in the Upper Floridan aquifer and the interaction of the river system with the aquifer [11].The Suwannee River Basin Model (SRBM) simulates a single aquifer layer under steady-state conditions approximate to September 1990 [11].Data for calibration came from September 1990 measurements of groundwater levels, discharge at gaging stations for the Suwannee, Alapaha, northern Withlacoochee, Santa Fe, Fenholloway, Aucilla, Econfina, and Steinhatchee Rivers, and spring flows at seven first-magnitude springs [11] (Figure 2).It should be noted that the study area includes two distinct Withlacoochee Rivers that are described herein as "northern Withlacoochee," referring to the major tributary with headwaters in Georgia, and "southern Withlacoochee," referring to the river discharging to the Gulf of Mexico approximately 30 mi. southeast of the Suwannee River mouth (Figure 2).

Study Area
The Suwannee, Santa Fe, Alapaha, and northern Withlacoochee Rivers are interconnected and drain the Suwannee River Basin.Outside of the Suwannee River Basin are the Aucilla, Econfina, Fenholloway, and Steinhatchee Rivers to the west and the Waccasassa and southern Withlacoochee Rivers to the south, all of which provide substantial groundwater drainage to the Gulf of Mexico (Figure 2).Large areas of the Suwannee River Basin are devoid of channelized or surface drainage, and most of the drainage directly infiltrates the karst topography of the area, which is generally flat and contains numerous sinkholes.The highly permeable aquifer allows numerous springs, including several major springs and spring groups to augment streamflow (Figure 2).The Upper Floridan aquifer is unconfined in much of the study area, but is generally covered by

Gulf of Mexico
Model Area Spring surficial sediments that confine or partially confine the Upper Floridan aquifer in the northern and eastern parts of the study area.
The climate of the Suwannee River Basin ranges from temperate to humid subtropical.Temperatures typically range from 39˚F to 50˚F in the winter and from 77˚F to 95˚F in the summer.Average annual precipitation in the study area ranges from about 51˚F to 59˚F inches per year (in/yr); about half of the annual precipitation occurs from June through September.Summer precipitation is generally associated with localized thunderstorms that can produce small-scale, intense rainfall.Winter precipitation is generally associated with cold fronts and is more evenly distributed geographically.For the purposes of this study, precipitation is equivalent to rainfall.

Methods
The hydrologic modeling for recent and future conditions in the Suwannee River Basin is implemented with MODFLOW-2000 [12] and is based on the original steady-state model by Planert [11].

Groundwater Model Parameterization
The groundwater model was discretized with 1 layer, 163 rows, and 148 columns and has a uniform grid-spacing of 5000 ft (Figure 3).The Suwannee, Santa Fe, Alapaha, northern Withlacoochee, Fenholloway, Waccasassa, and the southern Withlacoochee Rivers are represented by the MODFLOW River Package RIV1 as in the original steady-state simulation [11].The upper Fenholloway, Aucilla, Ecofina, and Steinhatchee Rivers and Spring Warrior Creek are represented with the MODFLOW Drain Package (Figure 2 and Figure 3).A total of 696 groundwater withdrawal wells are simulated in the study area.Details of the original steady-state model parameterization can be found in Planert [11].

Groundwater and Surface-Water Boundaries
The groundwater head boundaries follow the configuration in Planert [11].The northern boundary represents the effects of a potentiometric high near Valdosta, Georgia (Figure 4) and is represented by the specified-head package [12].This representation allows the groundwater head to vary linearly over each stress period of the simulation.Stage in the Suwannee, Santa Fe, Alapaha, and northern Withlacoochee Rivers was interpolated from eight stream gages (Figure 4).

Gainesville
were estimated for the lower Fenholloway, Waccasassa, and the southern Withlacoochee Rivers by using topographic maps.The difference in bottom and shoreline topography was used to estimate river stage.In order to predict river flows for both the recent and future simulations, a method was developed (described below) that used measured rainfall to generate boundary flows at the Withlacoochee River near Pinetta, Georgia; the Alapaha River at Statenville, Georgia; the Suwannee River at White Springs, Florida; and the Santa Fe River near High Springs, Florida (Figure 4).The net recharge or discharge to the aquifer for the river reaches downstream from the gages on these reaches was added or subtracted to these boundary values to determine river discharge from simulation results.

Transient Simulations
The SRBM was modified from the original steady-state version [11] to become a transient simulation tool, necessary for the simulation of different rainfall time series.For the recent rainfall simulations, the January 1970-December 2000 (372 months) period is used, and for the future rainfall simulation, the January 2039-December 2069 (372 months) period is used.A 1-month time step was chosen for the calibration of the transient model for the recent period.For comparisons between recent and future periods, the time step was reduced to 1 day.As well as providing better temporal resolution, a shorter time step was chosen for stability concerns when ET is transiently computed in the future simulations, as described below.
For the purposes of determining river stages for the future simulation period, the MODFLOW Streamflow-Routing (SFR2) Package [13] was applied to allow discharge in the rivers to be simulated and stages to be calculated (rather than being specified as input to the model).Only the Suwannee, Santa Fe, Alapaha, and northern Withlacoo-

Precipitation and Evapotranspiration in the Recent Period
Rainfall in the recent period was based on data from 15 rainfall stations (Figure 5).For each station, monthly average rainfall and temperature were computed from the field data.These monthly values were used to compute an initial estimate of net recharge (the difference between rainfall and ET) for each station based on the Thornthwaite method [14] as follows.The potential ET is calculated first using the Hamon equation [15]: where PET is potential evapotranspiration in millimeters per month (mm/month); d is the number of days in the month, D is the mean monthly hours of daylight in units of 12 hrs, and W t is a saturated water vapor density term:

4.95e 100
where T is the monthly mean temperature.In order to develop a more realistic actual ET, at times when the precipitation P is less than computed PET, the actual ET is set to the precipitation plus the available soil moisture (STW): ( ) where ST i−1 is the soil-moisture storage for the previous month and STC is the soilmoisture storage capacity.
Because the Thornthwaite method net-recharge values do not take into account loca-  4) and groundwater levels at monitoring wells (Figure 6) were used as criteria in a trial and error method to estimate multipliers for the net recharge at each rainfall station zone (Figure 5).A multiplier greater than 1.0 could indicate that ET is limited due to unsaturated zone depth, whereas a multiplier less than 1.0 could indicate losses due to runoff.A reasonably good fit was found with the following recharge multipliers:  The net recharge is calculated for the rainfall zones surrounding each rainfall station (Figure 5) by multiplying the initial estimate of recharge for the zone by its respective net-recharge multiplier.As expected, the lowest net recharge values occur along the coast (Perry, Steinhatchee, and Usher Tower; Figure 5) where coastal runoff losses would predominate and reduce the water available for net recharge.

Parameter Estimation
Adjusting the model-input hydraulic conductivity values is functionally identical to adjusting hydraulic transmissivity because the model aquifer thickness is defined as constant.Therefore, hydraulic conductivity is the MODFLOW input for calibration, but field values of hydraulic transmissivity are used for initializing and comparison, and the results are discussed in terms of transmissivity.The SRBM transient groundwater simulations were calibrated to reasonable values of hydraulic conductivity (or transmissivity) and storage properties (specific yield and specific storage) based on the limited aquifer testing available and properties derived from the previous steady-state modeling, using the parameter-estimation techniques implemented in the PEST software suite [16].PEST accesses MODFLOW input and output to evaluate the sensitivity of computed values to perturbations in model-input parameters [17].Planert [11] calibrated the steady-state version of the SRBM through trial-and-error comparisons of measured and computed groundwater wells, springs, and river discharge.Through experimentation, and guided by the known hydrostratigraphy, Planert [11] divided the aquifer into zones of differing transmissivity and specific storage.For the purpose of further calibrating the SRBM transient groundwater simulation, the parameter-estimation application is based on this same zonation, with seven zones chosen for adjustment by PEST (Figure 6).An additional four zones of varying size were not adjusted in the parameter-estimation application and retain the same aquifer parameters as the steadystate model [11].
The targets for the parameter estimation were measured aquifer heads at the groundwater wells shown in Figure 6.In order to reasonably limit the number of PEST parameters and accounting for previous calibration efforts, river discharges and spring flows were not used as targets for PEST.Upper and lower bounds for estimating transmissivity were defined differently for each zone (Table 1), but the specific storage estimation is given bounds of 0.000005 to 0.000200 for all zones (Table 2).The estimated transmissivity values range from 22,000 to 10,000,000 ft 2 /d (Table 1), whereas transmissivity values in the nonparameterized zones are as low as 1000 ft 2 /d (Figure 7).The range of transmissivity field values from pumping tests [11] is the same order of magnitude as the estimated values (Table 1) with the exception of zone 6 where the estimated value is 257,000 ft 2 /d and the single aquifer test in the zone yielded 25,000 ft 2 /d.
The estimated specific storage values, which were not needed for the earlier steady-state model, range from 0.000030 to 0.000200 (Table 2).

Gulf of Mexico
rivers (S021516001, S051311001, S061629001) and the poorest match is for a well near the coast (S111117007) where the simulated values are too high (Figure 8).The actual coastal groundwater exchange to offshore may be larger than the simulated value, assuming that the specified values of net recharge and the computed values of direct ET from groundwater are reasonably accurate.

Incorporation of Downscaled Global Climate Model Data and Simulation of Future Conditions
To gain insight into possible effects of changes in future rainfall, dynamically downscaled GCM rainfall data were scaled appropriately for input to the SRBM.The GCM used for this application is the Community Climate System Model (CCSM) [18]  increasing-greenhouse scenario for terrestrial carbon emissions was also used [19].

Future Rainfall and Evapotranspiration
In order to represent ET for future conditions, the preprocessed values representing present conditions were not considered useful.However, the COAPS effort also downscaled latent heat flux values that can be converted directly to potential ET (PET) by dividing by the latent heat of vaporization of water (2450 kJ/kg at 20 °C).This was done for the SRBM area, and the resulting monthly PET values are used in the MODFLOW Evapotranspiration Package.The ET rate is set to the PET value for water levels 2 ft or less below land surface.Typical ET extinction depths, to which the ET rates linearly decrease to zero, for grass and forest on loamy soils can be 20 ft [20], which was used due to the substantial percentage of forested land.This value would normally not be as deep in the less arid conditions of this study area, but the high connectivity of the karst aquifer warrants a deep extinction depth.
The rainfall and PET values were defined at the dynamic downscaled points shown in Figure 9.These monthly values were bilinearly interpolated to the Suwannee River model cell values though the scheme: where cellval is the value of rainfall or PET at the model cell, val i is the value of rainfall or PET at location i out of 4 downscaled points around the model cell, and dist i is the distance from the downscaled point i to the model cell.
The Suwannee River model incorporating the ET computation previously discussed had substantial stability problems when running on a monthly time step.The ET rate is a strong function of the groundwater depth, and the iterative scheme tends to oscillate: a high water level produces a high ET rate, causing a lower computed water level for the next iteration, resulting in a low ET rate, which causes the cycle to repeat and the model to not converge.Large changes in ET volumes are possible over the monthly time step; therefore the time step was reduced to 1 day (while maintaining a uniform recharge rate for all time steps within a given month), allowing for consistent solution convergence.The forcing functions (stress periods) are defined in monthly increments in the future simulation of conditions for the years 2039-2069, just as in the 1970-2000 simulation of recent historical conditions.

Future Surface-Water Inflows
Measured values for the surface-water boundary flows at the Withlacoochee River near Pinetta, Georgia, the Alapaha River at Statenville, Georgia, the Suwannee River at White Springs, Florida, and the Santa Fe River near High Springs, Florida, stream gages (Figure 4) are available for the recent-period simulation , but can also be used to develop relations between rainfall and the measured flows to generate boundary inflows for future periods.A catchment basin is associated with each rainfall station, with one or more rainfall stations with catchments representing inflows to each of the river boundary locations (Figure 10 and Table 3).The monthly rainfall measured at the station is multiplied by the catchment area and a runoff coefficient to produce the boundary flows.The runoff coefficients are adjusted to calibrate the computed flows to

Gulf of Mexico
To predict future (2039-2069) runoff flow in these catchments, the time series of average GCM-generated future rainfall in the SRBM area was used.Attempting to use a more localized predicted rainfall time series to drive these river boundaries was considered unnecessary because of the inherent uncertainty in the GCM predictions.The computed river boundary inflows for recent and future simulation periods are shown in Figure 11.

Future Surface-Water Stages
The primary objective for further developing the SRBM under transient conditions was to quantify the difference in streamflow during relatively recent times with those predicted and resulting from future climatic conditions.To contrast with the RIV1 Package used in the original SRBM, which computes groundwater/surface-water exchanges based on the difference between user-defined surface-water levels and model-simulated aquifer heads, the SFR2 Package was applied to help define differences in river stage between recent and future conditions.The SFR2 Package does not require user-defined surface-water levels (it computes them as part of the model solution), but is more complex and prone to stability problems than RIV1.Accordingly, the SRBM is applied for a shorter simulation period (113 months), and the resulting mean differences in stage between recent and future conditions are used to develop the user-defined stages that are required for RIV1 in the full-length future-conditions simulation.The SFR2 Package computes discharge with the steady uniform equations of flow and mass continuity [13], which can improve understanding of river discharge.The flow at the end of a stream reach is equal to the inflow plus or minus any groundwater leakage.Depth at the beginning and end of reaches is computed with Manning's equation in the form: The model application with SFR2 was compared with field measurements of stage and discharge at the Ellaville, Branford, and Wilcox stream gages on the Suwannee River and the Fort White stream gage on the Santa Fe River (Figure 4).The simulated variations in stage exceeded the upper range of the defined cross-sectional geometry; therefore, for the Suwanee River south of Branford (Figure 4), the cross-sectional  1970 1971 1972 1973 1974 1975 1976 1977 1978 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998 1999 1970 1971 1972 1973 1974 1975 1976 1977 1978 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998 1999 2000 1970 1971 1972 1973 1974 1975 1976 1977 1978 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998 1999 2000 1970 1971 1972 1973 1974 1975 1976 1977 1978 1979 1980 1981 1982 1983 1984 1985 1986 1987 1988 1989 1990 1991 1992 1993 1994 1995 1996 1997 1998 1999  geometry data were extended vertically by 10 ft, which expands the top width in the overbank from 620 ft to 720 ft.However, the cross-sectional geometry in the river remained the same as in RIV1.A better fit was obtained at all stream gages, especially for discharge, although the stage at the Wilcox stream gage on the Suwannee River and the Fort White stream gage on the Santa Fe River remained high (Figure 12).These errors may be due to nonuniform flow dynamics that the SFR2 cannot represent.
The SFR2 Package did not produce particularly accurate fluctuations in simulated stages nor was the full future 2039-2069 simulation period included, so the SFR2 simulated differences in stages between the recent and future simulation were used as a guide to modify the recent-period measured stages to estimate future-period stages for RIV1.The SFR2 simulation was applied to both the recent (1981-1990) and future (2050-2059) periods, and the mean difference in stage between the two simulations was determined at four locations (Figure 4): −0.89 ft at Fort White; 0.64 ft at Ellaville; 0.21 ft at Branford; and −0.29 ft at Wilcox.Sea-level differences are not considered, so zero change in stage at the coastline is defined for the interpolation across the Suwannee, Withlacoochee, and Santa Fe Rivers, producing the change values shown in Figure 13.
The surface-water depths in other rivers were not modified from the recent simulation.

Results
The recent  simulation period and the future (2039-2069) simulation period model results were compared to evaluate the predicted effects of future rainfall and ET changes.No efforts were made to represent other potential differences between recent and future conditions, such as groundwater withdrawal, sea level, or land cover changes.

Comparison of Recent and Future Rainfall and Evapotranspiration Rates
Actual ET rates for the future simulation were calculated during the simulation and can be tabulated from the model output.However, the recent simulation used ET computed from the Thornthwaite method summed with rainfall to be input as net recharge.Determining the actual rainfall and ET values in the recent simulation is complicated by the multipliers applied to the net recharge to correct for uncertainties, as discussed previously.When determining the rainfall and actual ET for the recent simulation, it was assumed that the multipliers account for uncertainty in the actual ET and the rainfall does not change, so the actual ET is determined by: ( )   (Table 4).Simulated average annual ET between these periods increases from 36.2 to 37.9 inches; this 4.5-percent increase is due to higher mean temperatures and latent heat predicted by the CCSM.Although both the actual ET and the rainfall are higher in the future period, the ET as a percentage of rainfall decreased slightly from 65 percent to 62 percent with a 4.4-inch increase in net annual recharge (Table 4).A higher spatial variation in the average annual net recharge is seen in the future period when compared to that in the recent period (Figure 14(a) and Figure 14(b)), but this can be reasonably attributed to the methods used in each simulation.The net recharge is preprocessed in the recent simulation, but in the future simulation the actual ET is computed during the simulation and must be combined with the preprocessed rainfall to determine net recharge.The recent simulation does not have any negative net recharge values as are seen in the future simulation along the coast and in the far  west (Figure 14(a) and Figure 14(b)).Because the average net recharge is higher in the future simulation (Table 4), the locations with the highest net future recharge areas have values substantially higher than in the recent simulation, mostly in the north and northeast portions of the model domain (Figure 14(c)).This is due to very little ET in these areas even though the average future ET is higher than the recent ET.The effect on groundwater levels is discussed in the next section.

Comparison of Simulated Recent and Future Groundwater Levels
With no changes in sea level represented, groundwater levels simulated in the future period are lower than in the recent period at locations nearer to the coast and substantially higher further inland and eastward (Figure 15).These areas of highest future groundwater increase correspond largely to areas of highest simulated future net recharge (Figure 14 no-flow condition representing a groundwater divide (Figure 4).Groundwater contours over the larger region [21] indicate that this groundwater divide could only be maintained as groundwater heads increased if the potentiometric high points outside the model area were also increased.The model boundary condition assumptions therefore depend on the hydrology outside the model area, which is not examined in this study.Certainly the groundwater heads in the northeastern area of the model domain have the highest uncertainty.
Smaller differences in groundwater levels between recent and future simulations are seen nearer the river system.In both simulations, the groundwater levels generally slope towards the river systems, indicating that the rivers drain the groundwater system.Water levels near most springs (Figure 2) do not differ appreciably between the simulations.The exception is with the westernmost spring near the Aucilla River (Figure 2 and Figure 15), which experiences lower groundwater levels; simulated mean spring flow decreased from 392 ft 3 /s to 188 ft 3 /s.
An important trend is observed when the simulated groundwater levels are expressed as a depth below land surface (Figure 16).These groundwater levels are consistently shallower in the future period in areas where they were deepest in the recent period.
Conversely, groundwater levels that were shallowest in the recent period become deeper in the future simulation period.Simulation results suggest future groundwater levels will have a greater conformity to land surface (Figure 16); the standard deviation (σ) of water-table depths changed from σ = 36.2ft in the recent simulation to σ = 28.8ft in the future simulation.This trend toward conformity is from a combined effect of higher rates of rainfall and ET in the future simulation (Table 4).Areas where the groundwater levels are relatively lower in the recent simulation received more rainfall recharge in the future simulation, and ET rates are always lower relative to rainfall at these greater depths to water.Areas where groundwater levels are relatively higher also get more rainfall recharge, but the higher ET rate near the surface offsets the higher rainfall recharge.Although the average increase in ET for the entire model domain is less than the average increase in rainfall, ET is relatively greater in wetland areas near the coast where water levels are much closer to surface (Figure 14(b)).The resulting future net recharge rates are more negative in coastal areas, which are consequently noticeably drier in the future simulation (Figure 16).Conversely, the easternmost areas of the model domain become substantially wetter in the future compared to those in the recent simulations.

Comparison of Discharge in Recent and Future Simulation Periods
River levels for the entire 372-month future simulation period were estimated based on testing of the SFR2 Package with the 113-month future period.These future changes in river levels were applied, and the total leakage exchanges for all the river reaches were summed to compute net flows at points along the rivers.A comparison of monthly average flow rates for the recent and future periods at the Ellaville, Branford, and Wilcox stream gages on the Suwannee River and the Fort White stream gage on the Santa Fe River (Figure 4) indicates differing patterns of discharge with similar peak magnitudes (Figure 17).The ratio of the standard deviations (σ) between the future and recent periods (σ future /σ recent ) was determined at each station to quantify discharge variability between the recent and future simulations: Ellaville = 0.924, Branford = 0.934, Wilcox = 0.981, and Fort White = 1.040.Note that the three sites on the Suwannee River all show a reduction in the simulated discharge variability with more reduction in variability Figure 17.Discharge for stream gages on the Suwannee and Santa Fe Rivers.
Figure 18.Cumulative discharge for stream gages on the Suwannee and Santa Fe Rivers.
surface-water inflows are based on rainfall-runoff relations and do not take into account the complex interactions that were simulated within the study area.
5) Simulation of groundwater levels was simplified as one layer with no vertical gradients in hydraulic head.The change in aquifer flow area with changing head over time is neglected because it is assumed that the saturated thickness of the aquifer is constant.
Boundary conditions assume that the potentiometric high near Valdosta, Georgia, be- recharge assure some confidence that a groundwater level increase would be expected.
The predicted groundwater differences near the coast have fewer limitations that induce uncertainty, so ultimately the future trend of higher rainfall and higher ET caus-ing the groundwater table depth to become more uniform seems reasonable.Even with the limitations described herein, the model simulations presented in this study can provide useful insight to Suwannee River Basin hydrology under potential future conditions of climate change.

Summary
A numerical model of the Suwannee River Basin area was adapted to compare recent hydrology to future hydrologic conditions resulting from predicted changes in rainfall and evapotranspiration (ET).The steady-state model of Planert [11] was modified to represent transient conditions and was calibrated with the PEST parameter-estimation code [16].The surface-water system was represented using the MODFLOW River Package, which required user-defined river levels.Some variables for simulation of the future period were developed from dynamically downscaled precipitation and potential ET produced by the Community Climate System Model (CCSM).The precipitation and potential ET data were interpolated for the Suwannee model grid from this downscaled dataset.The conditions simulated in the recent period  represent ET estimated with the Thornthwaite method based on measured data.Downscaled potential ET was applied with a standard function to decrease actual ET at lower groundwater levels under future conditions.This method produced undesirable oscillations at the month-long simulation time step originally used, and consequently the time step was reduced to 1 day.
To help guide the estimation of future river levels, the MODFLOW Streamflow-Routing (SFR2) Package was applied in a shorter recent simulation period (1981)(1982)(1983)(1984)(1985)(1986)(1987)(1988)(1989)(1990) and future simulation period (2050-2059), and average differences in river stage were determined.These differences in mean river stage were then used to create a time series of river stages for the River Package in the 2039-2069 simulation.The recent and future flows were calculated at river locations by adding or subtracting the simulated leakages for all upstream river segments to or from the initial inflows.
Comparisons of simulation results for the recent and future periods provide useful insight into the potential effects of future climate changes.Limitations from the representation of no-flow model boundaries and the spatial distribution of net recharge must be considered when interpreting model results.Even with these limitations, the simulations indicate that increasing rainfall and ET together in the future period caused groundwater levels to rise in areas where they were low in the recent period, and groundwater levels declined in areas where they were higher in the recent period.The general implication is that unsaturated zone depth would be more spatially uniform in the future, and vegetation that requires a range of conditions (substantially wetter or drier than average) would be detrimentally affected.This vegetation would include wetland species, especially in areas near the coast.

Figure 1 .
Figure 1.Location of the Suwannee River Basin study area.
This paper documents the modification and refinement of the SRBM to create a transient version of the model that is suitable for implementing downscaled GCM rainfall and ET data.The calibration of the transient simulations is discussed, along with the adaptation of the downscaled future precipitation to the SRBM, and the representation of future ET.Because the objective for the modeling was to compare the effects of future climate conditions to relatively recent simulated hydrologic conditions, no attempt was made to predict and simulate future pumping scenarios.Analysis of river and spring flow contributions to the Suwannee Estuary can also be used to predict the effects of changes in precipitation and ET on total freshwater input and salinity.

Figure 2 .
Figure 2. Suwannee River Basin Model area and hydrologic features.
hlaco oche e Riv er W it h la c o o c h e e R iv e r W a c c a s a s s a R iv e r S te in h a tc h e e R iv e r Spr ing War rior Cr eek F en ho llo w a y R iv e r A u c il la R iv e r E c o n f in a R iv e r Su wa nn ee R ive r S u w a n n e e R iv e r A la p a h a R iv e r Modifications to the original model include the addition of parameters for transient simulations and future climatic conditions.Parameter estimation techniques are used to define the spatial and temporal distribution of parameter values required for transient simulation, and potential ET is defined by several methods for the recent and future simulations.The datasets and model files and code are available at Swain [22].

Figure 4 .
Figure 4. Suwannee River Basin Model groundwater and surface-water boundaries.

chee
Rivers are represented by SFR2 for two 113-month periods, August 1981-December 1990 and August 2050-December 2059.The differences in computed stages are used to estimate river stages in the RIV1 Package for the future January 2039-December 2069 period.

Figure 5 .
Figure 5. Location of rainfall stations and rainfall zones in the Suwannee River Basin Model area.

Figure 6 .
Figure 6.Groundwater monitoring wells and zones for parameter estimation of hydraulic parameters.

Figure 6 Figure 7 .
Figure 6 indicate that the closest matches are for wells in the central areas near the Program), and a sea ice model (Community Sea Ice Model).The A2

Figure 9 .
Figure 9. Location of dynamically downscaled rainfall and PET from CCCSM.
depth, Q = flow rate, n = Manning's friction factor, C = Manning's constant (1.0 in SI units or 1.49 in English units), w = channel width, and S o is slope of the stream channel.The periods of simulation for the SFR2 Package are 113 months each-August 1981 through December 1990 for recent conditions and August 2050 through December 2059 for future conditions.This recent period has sufficient measured river flow data to allow a reasonable simulated-to-measured comparison.The SFR2 Package was applied only to the Suwannee, Santa Fe, Withlacoochee, and Alapaha Rivers (Figure1), and the RIV1 Package was still used in the Fenholloway River, the Waccasassa River, and the lower Withlacoochee River.River-bottom elevations and widths, and surface-water/ groundwater leakage coefficients are the same as used in the RIV1-only simulations.

Figure 12 .
Figure 12.Stream discharge and stage determined by the MODFLOW Streamflow-outing Package.
(c) and Figure15(c)).Although in this northeast area the model predicts markedly large recent to future increases in groundwater head, as much as 40 ft (Figure15(c)), it can be noted that slightly further east, at station S021516001 (Figure6and Figure8), observed groundwater heads do fluctuate as much as 35 ft over the recent time period.Therefore, the mean annual net recharge gain in this northeast area could support a large localized increase in groundwater heads if the different methods of computing ET in the recent and future simulations, and their effect on net recharge, are accepted as appropriate for comparison.Additional factors to be considered are the model boundary conditions near this northeast area, including a

6 ) 7 )
haves the same in the recent and future simulation periods, and the eastern boundary continues to coincide with the flow lines and acts as a no-flow boundary.This can be a factor in the large increase in simulated future groundwater heads in the northeastern model area.Similarly, the boundary at the southeast corner of the study area representing drainage to the Ocklawaha River is the same in the recent and future simulation periods.Sea-level change in offshore areas is not accounted for in the future simulation period.Changes in land use and water use due to human activities are not accounted for in the future-conditions simulation.Increased groundwater pumping, for example, can be expected as population increases, but no changes to groundwater withdrawals were simulated in the future-conditions simulation.Changing land use would induce changes in surface runoff and net-recharge rates, which have been shown to have an important effect on groundwater levels.The scope of the present study included a single GCM, the CCSM, to predict future rainfall and ET.Distinct differences in predicted rainfall and ET exist among GCM results because the assumptions needed to predict future climate have inherent uncertainty.Nevertheless, the SRBM simulation procedure described herein can be repeated with downscaled data from other GCMs.Our study was not designed to evaluate the variations between climate models but rather to demonstrate the hydrologic implications of one possible future climate scenario.This study demonstrates how downscaled GCM data can be combined with a groundwater flow model to gain insight into potential effects of climate change.8) The accuracy of the model results depends on the sensitivity of the model output to uncertainty in the model input parameters.The results of the parameter estimation and the scenario simulations suggest that groundwater heads are notably sensitive to net recharge values; at least for the range of values experienced in the scenario simulations.This indicates that uncertainty in aquifer transmissivity is a smaller source of error compared to factors affecting net recharge, such as the different methods of computing net recharge mentioned in limitation 3 above.The limitations listed above help guide the interpretation of model results and the relative confidence of the model results.The recharge and boundary variability issues listed in limitations 3 and 5, respectively, indicate that the high simulated groundwater heads in the northeastern model domain are quite uncertain, although the trends in net lized topographic dynamics such as losses due to rapid surficial runoff or limitations in ET due to unsaturated zone depth, the net-recharge values were further adjusted through model experimentation.The discharge values at the Ellaville, Branford, Fort White, and Wilcox stream gages (Figure

Table 2 .
Specific storage in parameter-estimation zones.
, which was developed by the University Corporation for Atmospheric Research (UCAR).The coupled components of the CCSM include an atmospheric model (Community At- mosphere Model), a land-surface model (Community Land Model), an ocean model Figure 8.Comparison of measured and simulated groundwater levels.

Table 4 .
Average rainfall and evapotranspiration rates from model simulations.