Hydrogeology and Simulated Groundwater Availability in Reaches 3 and 4 of the Washita River Aquifer, Southern Oklahoma, 1980–2017

Scientific Investigations Report 2023-5072
Prepared in cooperation with the Oklahoma Water Resources Board
By: , and 

Links

Acknowledgments

The project documented in this report was conducted in cooperation with the Oklahoma Water Resources Board (OWRB) as part of the U.S. Geological Survey (USGS) Water Availability and Use Science Program. The authors value the contributions of many OWRB and USGS staff that led to the successful completion of the project. The authors thank the OWRB for support on this project, especially Division Chief (Water Rights Administration Division) Christopher Neel and Technical Studies Manager (Water Rights Administration Division) Derrick Wagner, who provided hydrogeologic data and helped with defining study objectives and deliverables. The authors also thank Alan LePera and Byron Waltman, who assisted in collecting synoptic water-table-altitude measurements.

The authors express gratitude to USGS employees who performed data-collection activities in the field. Stephen Bradford, Kyle Cothren, Aaron Moyer, Kevin Smith, Steve Smith, and Adam Trevisan measured synoptic base flows during 2018. Leland Fuhrig and Kyle Rennell measured bedrock depths with the Geoprobe hydraulic profiling tool and collected synoptic water-table-altitude measurements. Derek Ryter contoured preliminary aquifer-base maps. In addition, Waylon Marler facilitated data entry to the USGS National Water Information System database. The authors also thank USGS employees Moussa Guira, Cory Russell, Chris Braun, and Martha Watt, who performed detailed technical reviews of this report and the accompanying model archive data release. The authors acknowledge and appreciate the professionalism, experience, and dedication of these helpful and resourceful colleagues.

Abstract

The 1973 Oklahoma Groundwater Law (Oklahoma Statutes §82–1020.5) requires that the Oklahoma Water Resources Board conduct hydrologic investigations of the State’s aquifers to determine the maximum annual yield for each groundwater basin. Because more than 20 years have elapsed since the final order was issued, the U.S. Geological Survey, in cooperation with the Oklahoma Water Resources Board, conducted an updated hydrologic investigation and evaluated the effects of potential groundwater withdrawals on groundwater flow and availability in reaches 3 and 4 of the Washita River aquifer in southern Oklahoma for a study period spanning 1980–2017. A hydrogeologic framework and conceptual model were developed to guide the construction and calibration of a numerical model of the Washita River aquifer. The numerical model was calibrated to water-table-altitude observations at selected wells, base-flow observations at selected U.S. Geological Survey streamgages, and the conceptual-model recharge.

Three types of groundwater-availability scenarios were run using the calibrated numerical model. These scenarios were used to (1) estimate equal-proportionate-share pumping rates, (2) quantify the potential effects of projected well withdrawals on groundwater storage over a 50-year period, and (3) simulate the potential effects of a hypothetical 10-year drought. With Washita River main-stem inflows, the 20-, 40-, and 50-year equal-proportionate-share pumping rates under normal recharge conditions were about 3.08 acre-feet per acre per year for reach 3 and about 3.80 acre-feet per acre per year for reach 4. Projected 50-year pumping scenarios were used to simulate the effects of modified well withdrawal rates. Because well withdrawals were less than 1 percent of the calibrated numerical-model water budget, changes to the well pumping rates had little effect on Washita River base flows and groundwater storage in the Washita River aquifer. A hypothetical 10-year drought scenario was used to simulate the potential effects of a prolonged period of reduced recharge on groundwater storage. Groundwater storage at the end of the drought period was 4.6 percent less than the groundwater storage of the calibrated numerical model at the end of the drought period.

Introduction

The 1973 Oklahoma Groundwater Law (Oklahoma Statutes §82–1020.5 [Oklahoma State Legislature, 2021b]) requires that the Oklahoma Water Resources Board (OWRB) conduct hydrologic investigations of the State’s aquifers (called groundwater basins in the statutes) to determine the maximum annual yield (MAY) for each groundwater basin. The MAY is defined as the total amount of fresh groundwater that can be annually withdrawn while ensuring a minimum 20-year life of that groundwater basin (OWRB, 2020). For alluvium and terrace groundwater basins, the life requirement is satisfied if, after 20 years of MAY withdrawals, 50 percent of the groundwater basin (hereinafter referred to as “aquifer”) retains a saturated thickness of at least 5 feet (ft) (OWRB, 2014). Although 20 years is the minimum period required by law, the OWRB can and often does consider management scenarios with longer periods. Once a MAY has been established, the amount of land owned or leased by a groundwater-use permit applicant determines the annual volume of water allocated to that groundwater-use permit applicant. The annual volume of groundwater allocated per acre of land is known as the equal-proportionate-share (EPS) pumping rate.

The OWRB issued a final order on November 13, 1990, that established the MAY (81,840 and 46,935 acre-feet per year [acre-ft/yr]) and EPS pumping rate (1.5 and 1.0 acre-foot per acre per year [(acre-ft/acre)/yr]) for reaches 3 and 4, respectively, of the Washita River aquifer in southern Oklahoma (OWRB, 2020). The MAY and EPS pumping rate were based on hydrologic investigations by Kent and others (1984) and Patterson (1984) that used a numerical groundwater-flow model (Trescott and others, 1976) to evaluate the effects of potential groundwater withdrawals on groundwater availability in the Washita River aquifer in southern Oklahoma. Every 20 years, the OWRB is statutorily required to update the hydrologic investigation on which the MAY and EPS were based. Because more than 20 years have elapsed since the final order was issued, the U.S. Geological Survey (USGS), in cooperation with the OWRB, conducted an updated hydrologic investigation and evaluated the effects of potential groundwater withdrawals on groundwater flow and availability in the Washita River aquifer of southern Oklahoma for a study period spanning 1980–2017.

Purpose and Scope

The purpose of this report is to describe a hydrologic investigation of the Washita River aquifer in southern Oklahoma. This description includes (1) an updated summary of the hydrogeology with a definition of the hydrogeologic framework (including an updated geographic extent) of the aquifer, (2) development of conceptual and calibrated numerical groundwater-flow models for the aquifer representing the 1980–2017 study period, and (3) results of simulations of groundwater-availability scenarios. The groundwater-availability scenarios use the calibrated numerical groundwater-flow model to (1) estimate the EPS pumping rate that ensures a minimum 20-, 40-, and 50-year life of the aquifer, (2) quantify the potential effects of projected well withdrawals on groundwater storage over a 50-year period, and (3) simulate the potential effects of a hypothetical (10-year) drought on groundwater storage. The calibrated numerical groundwater-flow model and groundwater-availability scenarios were archived and released in a USGS data release (Rogers and others, 2023).

The entire Washita River aquifer in Oklahoma consists of four administrative sections, or reaches, which are managed independently by the OWRB. Reach 1 extends from the Texas border to downstream from Foss Reservoir in western Oklahoma (fig. 1; Ellis and others, 2020), and reach 2 is hydraulically connected to and managed as part of the Rush Springs aquifer in western Oklahoma (fig. 1; Ellis, 2018a); both reach 1 and reach 2 are outside the scope of this report. The modeled area of this report is a 279.6-square-mile (mi2) (178,938-acre) subarea of the 411.5-mi2 extent of alluvium and terrace deposits (Heran and others, 2003; Miller and Stanley, 2004; Chang and Stanley, 2010; Stanley and Chang, 2012) of reaches 3 and 4 of the Washita River aquifer in southern Oklahoma (fig. 1). Reach 3 extends from near Anadarko, Okla., to Alex, Okla., and includes 113.6 mi2 (72,697 acres, or 40.6 percent) of the modeled area. Reach 4 extends from near Alex to south of Davis, Okla., and includes 166.0 mi2 (106,241 acres, or 59.4 percent) of the modeled area. The aquifer becomes narrow and thin in smaller tributary valleys, so this investigation focused on the main body of the aquifer along the Washita River main stem and selected large tributaries (modeled aquifer area, fig. 1). Though sometimes referred to as the “Washita River alluvial aquifer” (Schipper, 1983; OWRB, 2019b), the aquifer is referred to as the “Washita River aquifer” in this report because it consists not only of alluvium deposits but also of terrace and dune deposits. Selected sections in this report are modified from Ellis (2018a, b), Ellis and others (2020), and Smith and others (2021). Though the study areas are different, the organization and wording of this report are largely based on Smith and others (2021).

Figure 1. Maps of the Washita River aquifer study area including streamgage and observation
                        well locations in southwestern Oklahoma.
Figure 1.

Selected data-collection stations in the study area and the extent of A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma.

Description of Study Area

The Washita River aquifer is a long, narrow, unconfined aquifer that includes alluvium and terrace deposits of the Washita River and tributaries in Caddo, Grady, Garvin, McClain, and Murray Counties, Oklahoma (fig. 1). The Washita River is perennial in the study area except during the most extreme droughts as documented by historical USGS streamgaging records for the stream obtained from the USGS National Water Information System (NWIS) database (USGS, 2020; table 1). The Washita River generally flows from northwest to southeast for about 207 miles (mi) (Horizon Systems Corporation, 2015) in the study area. The largest tributaries in the study area are Little Washita River and Sugar, Rush, and Wildhorse Creeks (fig. 1). Many of the smaller tributaries to the Washita River in the study area are dammed by hundreds of Natural Resources Conservation Service (NRCS) floodwater-retarding structures (fig. 2), which were mostly built in the 1950s through 1970s (U.S. Army Corps of Engineers, 2021). Floodwater-retarding structures were designed to impound runoff from headwater catchments and slowly release it downstream with the purpose of reducing the frequency and magnitude of damaging flood flows on the Washita River main stem (Smith and Esralew, 2010). The Washita River main stem, however, is not dammed in the study area.

Table 1.    

Selected data-collection stations in and near the Washita River aquifer study area, southern Oklahoma.

[U.S. Geological Survey (USGS; 2020) data can be accessed using the 8- or 15-digit station number or other identifier. M/D/Y, month/day/year; NAVD 88, North American Vertical Datum of 1988; --, unknown or not applicable; SFR2, Streamflow-Routing package; SWB, Soil-Water Balance; WTF, water-table fluctuation method; OWRB, Oklahoma Water Resources Board]

Station number or identifier (figs. 1 and 2) Short name for station or other identifier Station name Latitude, in decimal degrees Longitude, in decimal degrees County Period of record (may contain gaps) (M/D/Y) Land-surface altitude, in feet above NAVD 88 Well or hole depth, in feet Contributing drainage area, in square miles Use in numerical groundwater-flow model
Begin End
Continuous-record streamgages (USGS, 2020)
07326500 Anadarko gage Washita River at Anadarko, Okla. 35.0842 −98.2434 Caddo 10/1/1963 Present (2022) -- -- 3,628 SFR2 inflow
07327000 Gracemont gage Sugar Creek near Gracemont, Okla. 35.1751 −98.2559 Caddo 10/1/1955 9/29/1974 -- -- 208 SFR2 inflow
07327550 Ninnekah gage Little Washita River east of Ninnekah, Okla. 34.9634 −97.8995 Grady 2/24/1992 9/29/2005 -- -- 232 SFR2 inflow
07328070 -- Winter Creek near Alex, Okla. 34.9931 −97.7614 Grady 10/1/1964 5/14/1987 -- -- 33.0 SFR2 inflow
07328100 Alex gage Washita River at Alex, Okla. 34.9259 −97.7739 Grady 10/1/1964 Present (2022) -- -- 4,756 Calibration
07328180 Criner gage North Criner Creek near Criner, Okla. 34.9715 −97.5848 McClain 10/15/1993 Present (2022) -- -- 7.19 SFR2 inflow
07328500 Pauls Valley gage Washita River near Pauls Valley, Okla. 34.7548 −97.2514 Garvin 10/1/1961 12/8/2020 -- -- 5,294 Calibration
07328550 -- Washington Creek near Pauls Valley, Okla. 34.8259 −97.2022 Garvin 7/3/1991 3/10/1994 -- -- 7.56 SFR2 inflow
07329500 Maysville gage Rush Creek near Maysville, Okla. 34.7434 −97.4053 Garvin 10/1/1954 9/29/1976 -- -- 206 SFR2 inflow
07329700 Hoover gage Wildhorse Creek near Hoover, Okla. 34.5415 −97.2472 Garvin 10/1/1969 6/30/2002 -- -- 600 SFR2 inflow
07329780 Davis gage Honey Creek below Turner Falls near Davis, Okla. 34.4318 −97.1472 Murray 10/22/2004 Present (2022) -- -- 16.4 SFR2 inflow
07331000 Dickson gage Washita River near Dickson, Okla. 34.2334 −96.9758 Carter 10/1/1961 Present (2022) -- -- 7,160 Calibration
Climate stations (Oklahoma Mesonet, 2019)
CHIC CHIC Chickasha 35.0324 −97.9145 Grady 1/1/1994 Present (2022) 1,076 -- -- Recharge (SWB, WTF)
PAUL PAUL Pauls Valley 34.7155 −97.2292 Garvin 1/1/1994 Present (2022) 955 -- -- Recharge (SWB, WTF)
Climate stations (National Centers for Environmental Information, 2022)
ANA USC00340224 ANADARKO, OK US 35.0619 −98.1988 Caddo 12/31/1892 11/29/2016 1,168 -- -- Recharge (SWB)
ANA2 US1OKCD0011 ANADARKO 3.8 E, OK US 35.0719 −98.1769 Caddo 12/3/2010 3/4/2017 1,160 Recharge (SWB)
ARD USC00340292 ARDMORE, OK US 34.1773 −97.1617 Carter 12/31/1900 Present (2022) 841 -- -- Recharge (SWB)
CES USC00341750 CHICKASHA EXPERIMENTAL STATION, OK US 35.0488 −97.9158 Grady 5/31/1953 Present (2022) 1,085 -- -- Recharge (SWB)
CNRA USC00341745 CHICKASAW NRA, OK US 34.4967 −96.9705 Murray 10/31/1978 Present (2022) 1,047 -- -- Recharge (SWB)
COX USC00342196 COX CITY 2 NE, OK US 34.7423 −97.7038 Grady 9/30/1980 6/29/2012 1,232 -- -- Recharge (SWB)
ELM USC00342872 ELMORE CITY 3 SW, OK US 34.6100 −97.4222 Garvin 7/9/1947 Present (2022) 1,020 -- -- Recharge (SWB)
HEA USC00344001 HEALDTON 3 E, OK US 34.2332 −97.4202 Carter 12/31/1893 Present (2022) 902 -- -- Recharge (SWB)
HEN USC00344052 HENNEPIN 5 N, OK US 34.5797 −97.3510 Garvin 8/31/1993 Present (2022) 966 -- -- Recharge (SWB)
LIN USC00345216 LINDSAY 2 W, OK US 34.8261 −97.6386 Garvin 3/31/1938 3/30/2010 980 -- -- Recharge (SWB)
MAD USC00345468 MADILL, OK US 34.0919 −96.7708 Marshall 11/30/1936 Present (2022) 770 -- -- Recharge (SWB)
PAUX USC00346931 PAULS VALLEY 1 SSW MESONET, OK US 34.7155 −97.2292 Garvin 1/31/2009 Present (2022) 954 -- -- Recharge (SWB)
PV4 USC00346926 PAULS VALLEY 4 WSW, OK US 34.7253 −97.2814 Garvin 12/31/1892 2/14/2011 940 -- -- Recharge (SWB)
TIS USC00348884 TISHOMINGO NATIONAL WLR, OK US 34.2070 −96.6424 Johnston 11/30/1902 Present (2022) 664 -- -- Recharge (SWB)
Continuous water-level recorder wells (USGS, 2020)
W01 344916097360001 04N-04W-15-ADA 1 W34001 34.8210 −97.5999 Garvin 3/28/2018 6/12/2019 972.68 51.5 -- --
W02 350514097593201 07N-08W-13 ADA 1 W34002 35.0874 −97.9923 Grady 11/23/2016 6/12/2019 1,099.23 36.8 -- Recharge (WTF)
W03 345824097484201 06N-06W-22 DDA 1 W34003 34.9733 −97.8116 Grady 11/23/2016 4/15/2019 1,044.56 62.8 -- --
W04 345137097414501 05N-05W-35 DBB 1 W34004 34.8602 −97.6957 Grady 7/2/2017 11/30/2017 999.84 -- -- --
W05 344704097183601 04N-01W-28 CDA 1 W34005 34.7846 −97.3100 Garvin 11/22/2016 6/12/2019 905.43 32.6 -- Recharge (WTF)
W06 343341097102901 01N-01E-14 BBD 1 W34006 34.5613 −97.1748 Garvin 4/6/2017 6/12/2019 823.34 -- -- --
Water-table-altitude observation wells (OWRB, 2019a)
25848 345017097205701 04N-01W-07 BBD 1 34.8381 −97.3495 Garvin 2/20/2001 1/26/2010 930 -- -- Calibration
25851 345102097252601 04N-02W-05 ADA 1 34.8506 −97.4242 Garvin 2/20/2001 1/26/2010 942 -- -- Calibration
9428 345304097421301 05N-05W-26 BBB 1 34.8845 −97.7039 Grady 3/10/1983 2/3/1993 1,010 -- -- Calibration
Table 1.    Selected data-collection stations in and near the Washita River aquifer study area, southern Oklahoma.
Figure 2. Map of Oklahoma water-management basins, flood-retarding structures, and
                        selected data-collection stations in the study area.
Figure 2.

Selected data-collection stations and the major geographic and surface-water features in and near the Washita River aquifer study area, southern Oklahoma.

Land-cover data for the Washita River aquifer study area were obtained from the CropScape database (National Agricultural Statistics Service, 2019) (fig. 3), which includes land-cover characteristics at 30-meter (m) resolution for land overlying the Washita River aquifer for the 2010–17 period. During this period, land-cover type was primarily cropland (44.1 percent), grass or pasture (36.3 percent), forest and shrubland (9.3 percent), or developed (7.7 percent). Winter wheat was the major crop-cover type in the study area, accounting for 57.6 percent of cropland by area. Alfalfa or hay (22.0 percent) was the next largest crop-cover type by area. Corn and soybeans accounted for 7.5 and 6.6 percent of cropland by area, respectively. Cotton (1.3 percent) and fallow or idle cropland (1.2 percent) were the only other crop-cover types that accounted for more than 1 percent of cropland by area (fig. 3). Crop-cover types typically change in response to economic conditions and hydrologic factors, but the percentages of total cropland cover and individual crop types did not change substantially during the 2010–17 period.

Figure 3. Maps and pie charts showing the type of land and crop cover with associated
                        quantities over the Washita River aquifer.
Figure 3.

Land- and crop-cover types for land overlying A, reach 3, B, reach 4, and C, the study area of the Washita River aquifer, southern Oklahoma, 2010–17.

Climate Characteristics

The climate of the study area is classified as subhumid with generally southerly prevailing winds (Patterson, 1984). Historical data from selected climate stations in southern Oklahoma (Climate Division 8) have been quality assured and summarized monthly as part of the Global Historical Climatology Network (National Centers for Environmental Information, 2021). These data were used to calculate and graph annual and monthly temperature and precipitation statistics for the study area (fig. 4, table 2). A locally weighted scatterplot smoothing line (Cleveland, 1979) was used to delineate periods of below- and above-mean annual temperature and precipitation (fig. 4). The mean annual temperature during 1895–2020 was 62.2 degrees Fahrenheit (°F) (fig. 4, table 2), and the mean annual temperature differed by about 2 °F from northwest to southeast across the study area (Oklahoma Climatological Survey, 2021a, b). The mean annual precipitation during 1895–2020 was 38.0 inches per year (in/yr) (fig. 4, table 2), and the mean annual precipitation increased by about 7.5 in/yr from northwest to southeast across the study area (Oklahoma Climatological Survey, 2021a, b). The mean annual precipitation for the 1980–2017 study period was 40.1 in/yr (fig. 4, table 2). The 1981–2002 period was an unprecedented wet period in the available record; above-mean precipitation was recorded in 18 of these 22 years (82 percent of the years), contributing to a higher mean precipitation during the study period compared to the mean for the period of record (fig. 4A). May is typically the wettest month, and January is typically the driest month (fig. 5A).

Figure 4. Graph shows mean precipitation and temperature for the period of record
                        and mean temperature and precipitation fluctuations.
Figure 4.

A, Precipitation and B, temperature, southern Oklahoma (Climate Division 8), 1895–2020 (National Centers for Environmental Information, 2021).

Table 2.    

Mean annual precipitation and mean annual temperature for selected periods, southern Oklahoma, 1895–2020.

[Reach and station locations shown on figures 1 and 2. --, data not available or not applicable]

Region or location Period Number of years Mean annual precipitation, in inches per year Mean annual temperature, in degrees Fahrenheit
Southern Oklahoma (Climate Division 8; National Centers for Environmental Information, 2021) 1895–2020 126 38.0 62.2
1895–1936 42 36.3 62.2
1937–78 42 37.2 62.2
1979–2020 42 40.5 62.3
1980–2017 38 40.1 62.4
Reach 3 of the Washita River aquifer (Chickasha, Okla., station USC00341750; National Centers for Environmental Information, 2022) 1980–2017 (missing daily values [1.3 percent of record] were assumed to equal the mean of available daily values) 38 36.2 --
Reach 4 of the Washita River aquifer (Pauls Valley, Okla., stations USC00346926 and USC00346931; National Centers for Environmental Information, 2022) 1980–2017 (missing daily values [3.3 percent of record] were assumed to equal the mean of available daily values) 38 38.6 --
Table 2.    Mean annual precipitation and mean annual temperature for selected periods, southern Oklahoma, 1895–2020.
Figure 5. Graph depicting monthly precipitation and temperatures during the period
                        of record and study period record.
Figure 5.

A, Mean monthly precipitation and B, mean monthly temperature, southern Oklahoma (Climate Division 8), 1895–2020 and 1980–2017 (National Centers for Environmental Information, 2021).

Multiyear to decadal droughts are common for the study area (fig. 4). The 1929–41 (Dust Bowl; Egan, 2006) and 1952–56 drought periods were among the most severe in Oklahoma in the 20th century; two shorter, less severe drought periods also occurred in the late 20th century during 1961–72 and 1976–81 (Tortorelli, 2008; Shivers and Andrews, 2013). The 21st century began with the 2002–6 and 2010–14 drought periods (Tortorelli, 2008; Shivers and Andrews, 2013). The most severe droughts on record developed from extended periods of below-mean precipitation paired with above-mean temperature.

Streamflow Characteristics and Trends

Daily streamflow data, with various and sometimes interrupted periods of record, were recorded at selected USGS streamgages in the study area (figs. 1 and 2) and summarized for the 1980–2017 study period (table 3). Streamflow measured at streamgages is the sum of runoff and base flow originating upstream; base flow is the component of streamflow that is supplied by the discharge of groundwater to streams (Barlow and Leake, 2012). For this report, streamflow-hydrograph data obtained from the NWIS database (USGS, 2020) were separated into runoff and base-flow components by using the standard Base-Flow Index code (Wahl and Wahl, 1995) in the USGS Groundwater Toolbox (Barlow and others, 2015). The Base-Flow Index code uses the minimum streamflow in a moving n-day window as a basis for hydrograph separation; for consistency, a 6-day window was used for all streamgages in this report. The 6-day window was selected by testing multiple n-day windows for the study period at USGS streamgage 07328100 Washita River at Alex, Okla. (hereinafter referred to as the “Alex gage”) and plotting the resulting mean base-flow percentage against n; a slope change was evident at 6 days. Base flow, as computed using this method, accounts for about 35–48 percent of annual streamflow for the period of record at streamgages on the Washita River main stem (table 3). This percentage, known as the base-flow index (BFI) (Wahl and Wahl, 1995; Barlow and others, 2015), generally increased over the period of record at streamgages on the Washita River (fig. 6AD). When summarized annually using the Kendall (1938) tau test with the Theil-Sen (Sen, 1968) slope estimator and an alpha value of 0.0500 set as the significance level, the upward trends in the BFI were found to be statistically significant (probability value [p-value] = 0.0054) at USGS streamgage 07326500 Washita River at Anadarko, Okla. (hereinafter referred to as the “Anadarko gage”), statistically significant (p-value = 0.0248) at the Alex gage, statistically significant (p-value = 0.0000) at USGS streamgage 07328500 Washita River near Pauls Valley, Okla. (hereinafter referred to as the “Pauls Valley gage”), and statistically significant (p-value = 0.0083) at USGS streamgage 07331000 Washita River near Dickson, Okla. (hereinafter referred to as the “Dickson gage”). Upward trends in base flow (fig. 6AD) were found to be not statistically significant (p-value = 0.0772) at the Anadarko gage, statistically significant (p-value = 0.0111) at the Alex gage, statistically significant (p-value = 0.0163) at the Pauls Valley gage, and statistically significant (p-value = 0.0177) at the Dickson gage. The overall pattern of upward trends in the BFI and base flow may have developed following the construction of more than 570 NRCS floodwater-retarding structures (dams; fig. 2) on tributaries of the Washita River within the study area. NRCS dam construction in this area began in 1948, peaked at 65 dams per year in 1963, greatly diminished to about 1 dam per year in the 1980s, and concluded in the early 1990s (U.S. Army Corps of Engineers, 2021). NRCS dams temporarily impound runoff and release it over several days with the purpose of reducing peak flows on the river main stem. The flows that are slowly and steadily released by these dams are indistinguishable from groundwater-derived base flows when they reach the main stem and are mostly classified as base flows by the Base-Flow Index code. Though most upward trends in the BFI and base flow were statistically significant for the streamgage periods of record, no statistically significant trends in base flow or BFI were detected for the 1980–2017 study period at streamgages on the Washita River main stem. One interpretation of the upward trends in the BFI and base flow for the period of record is that the construction of the floodwater-retarding structures caused a step change in apparent base flows in the Washita River.

Table 3.    

Annual streamflow, base-flow, and base-flow index values for selected U.S. Geological Survey streamgages in and near the Washita River aquifer study area, southern Oklahoma, for the various periods of record, 1903–2020, and for the study period, 1980–2017.

[Data from the U.S. Geological Survey (USGS) National Water Information System (USGS, 2020). Values computed by using the Base-Flow Index code (Wahl and Wahl, 1995) in the USGS Groundwater Toolbox (Barlow and others, 2015). Gage locations shown on figures 1 and 2. Gray shading indicates the data for the 1980–2017 study period. BFI, base-flow index; --, data not available or not applicable; POR, period of record]

Year USGS streamgage 07326500 Washita River at Anadarko, Okla. USGS streamgage 07328100 Washita River at Alex, Okla. USGS streamgage 07328500 Washita River near Pauls Valley, Okla. USGS streamgage 07331000 Washita River near Dickson, Okla.
Mean streamflow, in cubic feet per second Mean base flow Mean streamflow, in cubic feet per second Mean base flow Mean streamflow, in cubic feet per second Mean base flow Mean streamflow, in cubic feet per second Mean base flow
In cubic feet per second In percent (BFI) In cubic feet per second In percent (BFI) In cubic feet per second In percent (BFI) In cubic feet per second In percent (BFI)
1903 614.7 227.0 37 -- -- -- -- -- -- -- -- --
1904 257.7 128.9 50 -- -- -- -- -- -- -- -- --
1905 417.3 173.0 41 -- -- -- -- -- -- -- -- --
1906 385.5 207.6 54 -- -- -- -- -- -- -- -- --
1907 1,074.9 398.1 37 -- -- -- -- -- -- -- -- --
1908 -- -- -- -- -- -- -- -- -- -- -- --
1909 -- -- -- -- -- -- -- -- -- -- -- --
1910 -- -- -- -- -- -- -- -- -- -- -- --
1911 -- -- -- -- -- -- -- -- -- -- -- --
1912 -- -- -- -- -- -- -- -- -- -- -- --
1913 -- -- -- -- -- -- -- -- -- -- -- --
1914 -- -- -- -- -- -- -- -- -- -- -- --
1915 -- -- -- -- -- -- -- -- -- -- -- --
1916 -- -- -- -- -- -- -- -- -- -- -- --
1917 -- -- -- -- -- -- -- -- -- -- -- --
1918 -- -- -- -- -- -- -- -- -- -- -- --
1919 -- -- -- -- -- -- -- -- -- -- -- --
1920 -- -- -- -- -- -- -- -- -- -- -- --
1921 -- -- -- -- -- -- -- -- -- -- -- --
1922 -- -- -- -- -- -- -- -- -- -- -- --
1923 -- -- -- -- -- -- -- -- -- -- -- --
1924 -- -- -- -- -- -- -- -- -- -- -- --
1925 -- -- -- -- -- -- -- -- -- -- -- --
1926 -- -- -- -- -- -- -- -- -- -- -- --
1927 -- -- -- -- -- -- -- -- -- -- -- --
1928 -- -- -- -- -- -- -- -- -- -- -- --
1929 -- -- -- -- -- -- -- -- -- 1,545.4 561.1 36
1930 -- -- -- -- -- -- -- -- -- 1,556.6 581.8 37
1931 -- -- -- -- -- -- -- -- -- 1,045.9 375.5 36
1932 -- -- -- -- -- -- -- -- -- 1,731.7 540.4 31
1933 -- -- -- -- -- -- -- -- -- 1,447.6 368.0 25
1934 -- -- -- -- -- -- -- -- -- 836.2 332.3 40
1935 -- -- -- -- -- -- -- -- -- 2,249.9 731.7 33
1936 446.5 154.7 35 -- -- -- -- -- -- 1,031.2 425.4 41
1937 -- -- -- -- -- -- -- -- -- 671.3 201.5 30
1938 -- -- -- -- -- -- 765.2 371.3 49 1,700.9 449.9 26
1939 -- -- -- -- -- -- 308.5 98.4 32 372.1 130.0 35
1940 -- -- -- -- -- -- 305.4 86.8 28 804.8 176.1 22
1941 -- -- -- -- -- -- 1,717.3 593.6 35 3,405.6 1,043.6 31
1942 -- -- -- -- -- -- 1,280.7 622.1 49 2,805.8 1,011.5 36
1943 -- -- -- -- -- -- 852.8 398.8 47 1,856.7 484.5 26
1944 -- -- -- -- -- -- 681.4 282.3 41 1,069.8 431.8 40
1945 -- -- -- -- -- -- 1,569.1 460.1 29 4,066.8 846.7 21
1946 -- -- -- -- -- -- 849.3 302.4 36 2,111.0 580.7 28
1947 -- -- -- -- -- -- 1,232.7 381.0 31 2,178.1 769.4 35
1948 -- -- -- -- -- -- 575.9 207.8 36 1,004.1 326.9 33
1949 -- -- -- -- -- -- 1,197.5 407.7 34 1,777.4 635.5 36
1950 -- -- -- -- -- -- 961.9 354.2 37 1,977.6 521.6 26
1951 -- -- -- -- -- -- 1,001.0 319.3 32 1,717.9 532.4 31
1952 -- -- -- -- -- -- 370.4 115.7 31 573.0 165.5 29
1953 -- -- -- -- -- -- 427.3 77.4 18 910.9 149.8 16
1954 -- -- -- -- -- -- 362.8 132.2 36 892.1 268.9 30
1955 -- -- -- -- -- -- 714.5 151.4 21 1,091.5 210.0 19
1956 -- -- -- -- -- -- 170.7 63.5 37 248.5 86.9 35
1957 -- -- -- -- -- -- 1,663.7 667.4 40 3,616.3 977.8 27
1958 -- -- -- -- -- -- 390.7 185.6 48 776.7 324.6 42
1959 -- -- -- -- -- -- 998.0 265.1 27 1,145.5 308.3 27
1960 -- -- -- -- -- -- 790.1 437.4 55 1,407.8 595.0 42
1961 -- -- -- -- -- -- 813.9 363.8 45 1,171.4 468.9 40
1962 -- -- -- -- -- -- 794.9 274.6 35 1,239.0 341.3 28
1963 -- -- -- -- -- -- 240.6 114.0 47 361.6 153.3 42
1964 211.1 66.6 32 -- -- -- 365.9 101.8 28 653.6 157.7 24
1965 438.2 131.5 30 489.4 168.3 34 509.5 197.0 39 690.7 274.4 40
1966 137.0 95.5 70 195.6 113.7 58 217.7 121.3 56 334.2 137.7 41
1967 82.4 33.5 41 148.3 49.3 33 181.0 57.4 32 437.1 88.0 20
1968 290.7 90.0 31 358.1 132.3 37 541.9 175.4 32 1,366.8 348.3 25
1969 353.4 120.1 34 472.6 143.8 30 619.8 232.1 37 1,361.2 426.3 31
1970 92.8 52.3 56 162.6 72.7 45 276.9 105.8 38 1,024.5 234.1 23
1971 110.7 35.3 32 185.0 55.7 30 231.9 72.4 31 481.9 146.8 30
1972 84.9 36.2 43 137.8 61.3 45 225.1 88.7 39 528.6 174.3 33
1973 431.4 160.7 37 672.7 279.8 42 947.6 367.9 39 2,890.8 722.7 25
1974 298.1 125.3 42 451.8 210.9 47 599.3 266.1 44 1,803.1 537.5 30
1975 564.5 292.1 52 842.0 441.2 52 1,220.8 567.0 46 2,493.7 1,028.9 41
1976 181.7 119.7 66 275.4 183.8 67 320.0 228.5 71 618.8 325.3 53
1977 373.1 127.0 34 483.4 179.0 37 649.8 195.2 30 1,128.9 305.8 27
1978 228.2 95.4 42 378.3 131.4 35 505.0 186.5 37 890.9 188.4 21
1979 234.0 86.6 37 410.8 142.8 35 553.0 193.7 35 1,078.0 296.3 27
1980 318.0 132.3 42 408.9 194.9 48 555.2 258.9 47 863.5 355.5 41
1981 161.1 72.4 45 250.7 99.9 40 279.0 101.6 36 978.6 228.9 23
1982 685.9 308.9 45 882.3 377.8 43 1,045.4 600.8 57 2,368.5 1,178.4 50
1983 625.7 195.3 31 887.1 318.0 36 1,206.3 411.5 34 1,865.4 555.7 30
1984 141.3 118.2 84 301.7 168.7 56 439.9 255.8 58 930.7 418.3 45
1985 264.9 122.2 46 748.1 220.7 30 1,560.0 524.4 34 3,086.7 1,089.2 35
1986 913.4 420.3 46 -- -- -- 1,714.0 855.5 50 2,394.7 1,144.1 48
1987 1,172.7 777.1 66 -- -- -- 2,782.3 1,400.6 50 4,839.1 1,927.8 40
1988 648.4 403.4 62 -- -- -- 1,114.5 591.0 53 2,089.0 1,209.8 58
1989 943.4 461.1 49 1,172.5 494.8 42 1,294.4 513.8 40 2,309.2 891.7 39
1990 685.9 359.3 52 1,154.8 559.5 48 1,633.8 759.5 46 4,532.1 1,198.0 26
1991 456.3 178.8 39 811.1 271.3 33 1,061.1 475.2 45 2,784.0 794.0 29
1992 615.3 303.0 49 1,239.6 599.6 48 1,601.9 915.4 57 3,038.0 1,459.9 48
1993 1,092.6 444.8 41 1,631.3 1,087.1 67 2,184.7 1,289.2 59 4,428.9 1,823.0 41
1994 312.7 169.2 54 567.5 351.1 62 726.1 418.8 58 2,385.9 696.7 29
1995 1,170.0 364.8 31 1,675.4 763.9 46 2,231.8 951.2 43 4,036.5 1,487.6 37
1996 771.9 433.7 56 1,074.9 629.7 59 1,471.4 887.3 60 2,360.0 1,331.0 56
1997 1,452.3 990.6 68 1,631.1 1,094.6 67 2,062.7 1,307.3 63 3,095.1 1,806.9 58
1998 820.9 525.9 64 1,243.0 726.5 58 1,599.4 836.6 52 2,533.9 1,182.2 47
1999 493.0 278.4 56 696.4 333.3 48 874.3 462.9 53 1,398.0 717.6 51
2000 593.2 235.2 40 899.5 353.9 39 1,133.3 340.8 30 1,802.2 608.0 34
2001 517.0 328.9 64 770.1 459.5 60 1,003.6 564.0 56 2,289.0 1,031.2 45
2002 212.5 117.7 55 355.2 171.7 48 465.6 239.5 51 1,245.7 332.1 27
2003 209.3 145.5 69 293.0 191.5 65 416.5 275.0 66 794.2 398.2 50
2004 405.0 124.0 31 491.3 176.4 36 636.0 228.1 36 1,622.5 489.2 30
2005 384.4 188.9 49 413.3 238.7 58 490.2 291.5 59 1,276.4 437.7 34
2006 146.6 85.6 58 191.5 103.9 54 209.2 103.3 49 427.0 134.3 31
2007 1,486.6 491.6 33 2,122.6 646.2 30 3,083.3 1,251.6 41 6,410.2 2,164.4 34
2008 628.3 291.1 46 1,008.5 516.1 51 1,321.3 709.6 54 2,041.8 753.8 37
2009 391.6 218.5 56 596.0 291.8 49 794.5 467.4 59 2,011.5 575.2 29
2010 315.7 157.8 50 471.5 274.1 58 667.3 366.2 55 1,282.3 628.0 49
2011 121.5 65.0 54 188.0 89.3 48 229.1 109.1 48 403.8 147.9 37
2012 111.6 51.8 46 190.4 86.5 45 219.4 100.9 46 578.2 200.3 35
2013 97.9 43.0 44 236.6 79.0 33 295.9 125.3 42 539.9 178.7 33
2014 73.2 37.4 51 116.2 60.5 52 141.1 88.3 63 259.6 112.5 43
2015 809.5 328.8 41 1,322.1 504.3 38 1,578.7 616.5 39 3,444.0 1,179.4 34
2016 470.1 240.9 51 753.6 411.5 55 962.9 472.3 49 1,997.8 852.5 43
2017 563.5 274.8 49 772.8 358.0 46 871.8 496.5 57 1,588.8 647.7 41
2018 472.6 225.9 48 710.3 315.1 44 853.7 379.7 44 2,135.5 627.1 29
2019 1,319.6 776.6 59 1,898.3 1,198.6 63 2,103.7 1,288.4 61 3,494.2 1,684.7 48
2020 346.4 222.9 64 592.3 327.7 55 -- -- -- 2,300.3 669.7 29
Mean, 1980–2017 560.1 275.9 50 787.7 380.1 48 1,104.2 543.8 50 2,166.7 851.8 39
Mean, 1980–2017, in thousands of acre-feet per year 405.8 199.9 -- 570.6 275.4 -- 799.9 393.9 -- 1,569.7 617.1 --
Mean, POR 483.1 234.7 48 687.4 330.4 47 901.4 405.1 44 1,787.1 635.3 35
Mean, POR, in thousands of acre-feet per year 350.0 170.0 -- 498.0 239.4 -- 653.0 293.5 -- 1,294.7 460.3 --
Table 3.    Annual streamflow, base-flow, and base-flow index values for selected U.S. Geological Survey streamgages in and near the Washita River aquifer study area, southern Oklahoma, for the various periods of record, 1903–2020, and for the study period, 1980–2017.
Figure 6. Graphs show a general increase in base flow at four main-stem gages in the
                        study area during the study period.
Figure 6.

Annual base-flow and annual base-flow index values for U.S. Geological Survey streamgages A, 07326500 Washita River at Anadarko, Okla.; B, 07328100 Washita River at Alex, Okla.; C, 07328500 Washita River near Pauls Valley, Okla.; and D, 07331000 Washita River near Dickson, Okla., on the main stem of the Washita River in the Washita River aquifer study area, southern Oklahoma, for the various periods of record, 1903–2020 (U.S. Geological Survey, 2020). Lowess, locally weighted scatterplot smoothing.

Surface-water use in the Washita River aquifer study area consists of permitted diversions along the main stem of the Washita River and its associated tributaries (OWRB, 2019d). At the end of the study period (2017), 60 active, nontemporary permitted diversions on the Washita River main stem allowed for a maximum total withdrawal of 13,957 acre-ft/yr. All of those permitted diversions were for irrigation with the exception of three permitted diversions that were used for public supply, and all but 15 of those permitted diversions were in reach 3 of the Washita River aquifer. At the end of the study period (2017), 29 active, nontemporary permitted diversions on Washita River tributaries allowed for a maximum total withdrawal of 9,749 acre-ft/yr.

Groundwater Use

The OWRB permits and regulates groundwater use of more than 5 acre-ft/yr for domestic and agricultural purposes and groundwater use for irrigating more than 3 acres of land for growing gardens, orchards, or lawns in Oklahoma (Oklahoma Statutes §82–1020.1[2] [Oklahoma State Legislature, 2021a; OWRB, 2014]; Oklahoma Statutes §82–1020.3 [Oklahoma State Legislature, 2021c]). Groundwater-use data since 1980 are self-reported annually to the OWRB by permitted users (Christopher Neel, OWRB, written commun., 2019). In 2017, about 61 long-term temporary (regular) groundwater-use permits and about 66 prior-right groundwater-use permits were active for the Washita River aquifer (OWRB, 2019d). Each permit may be served by multiple wells that share the allocated groundwater use. Most groundwater-use permits for the Washita River aquifer were allocated for irrigation (38.0 percent) and public supply (61.4 percent), with some allocated for other uses (0.6 percent), including agricultural (nonirrigation); commercial; industrial; mining; power; recreation, fish, and wildlife; or uncategorized uses (figs. 7 and 8, table 4) (OWRB, 2019d; Rogers and others, 2023). Groundwater use for domestic supply (self-supplied directly to a home by a private well) was assumed to be a negligible part of the total groundwater use because the population of the land overlying the aquifer was predominantly rural or concentrated in cities supplied water by municipal or rural water districts. For the purposes of this report, all groundwater use was assumed to be consumptive use (that is, none of the groundwater that is withdrawn returns to the aquifer as recharge).

Figure 7. Maps showing that the abundance of groundwater wells occur in reach 4 when
                        compared to reach 3 of the Washita River aquifer.
Figure 7.

Land areas and wells permitted for groundwater use from A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, 2017.

Figure 8. Bar graphs and pie charts showing that groundwater use is predominantly
                        public supply and irrigation.
Figure 8.

Annual reported and total annual permitted groundwater use by use type from A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, 1967–2017.

Table 4.    

Annual reported groundwater use from reach 3 and reach 4 of the Washita River aquifer, southern Oklahoma, 1967–2017.

[Permit-level reported groundwater-use data from the Oklahoma Water Resources Board (OWRB) were aggregated by groundwater-use type in this table (Rogers and others, 2023) owing to restrictions of proprietary interest and permit-holder anonymity; table excludes groundwater use of less than 5 acre-feet per year for domestic and agricultural purposes and groundwater use for irrigation of fewer than 3 acres of land for growing of gardens, orchards, or lawns (Oklahoma Statutes §82–1020.3; Oklahoma State Legislature, 2021c). All values are in units of acre-feet per year. Totals may not equal sum of components because of independent rounding. Gray shading indicates the data for the 1980–2017 study period]

Reach 3 of the Washita River aquifer Reach 4 of the Washita River aquifer
Year Groundwater-use type Total groundwater use Total permitted and reported use in the previous 5 years Groundwater-use type Total groundwater use Total permitted and reported use in the previous 5 years
Irrigation Public supply Other Irrigation Public supply Other
1967 263.3 1,358.6 0.0 1,621.9 3,845.0 1,190.3 968.8 102.8 2,261.9 6,715.0
1968 398.0 0.0 18.7 416.7 3,934.0 475.6 449.8 41.7 967.1 6,757.0
1969 198.7 0.0 25.5 224.2 3,934.0 1,305.2 388.1 33.8 1,727.1 5,628.0
1970 246.7 0.0 25.0 271.7 4,262.0 1,290.7 543.7 82.4 1,916.8 5,636.0
1971 368.7 2,819.2 33.1 3,221.0 4,299.0 1,112.5 646.1 85.5 1,844.1 5,636.0
1972 480.8 0.0 4.9 485.7 4,299.0 1,097.4 807.9 7.6 1,912.9 6,478.0
1973 139.0 0.0 41.2 180.2 4,299.0 944.5 1,258.2 117.0 2,319.7 6,395.0
1974 207.0 0.0 37.7 244.7 4,210.0 899.5 1,104.8 76.2 2,080.5 6,449.0
1975 335.0 0.0 9.8 344.8 4,450.0 534.9 595.9 0.0 1,130.8 6,436.0
1976 314.0 0.0 0.0 314.0 4,282.0 819.3 1,157.7 0.0 1,977.0 6,544.0
1977 533.7 0.0 0.0 533.7 4,282.0 1,676.7 1,127.1 0.0 2,803.8 6,155.0
1978 603.8 0.0 0.0 603.8 4,282.0 2,071.7 1,193.6 0.0 3,265.3 6,788.0
1979 505.9 0.0 0.0 505.9 4,282.0 633.6 1,180.6 0.9 1,815.1 6,910.7
1980 1,111.0 0.0 0.9 1,111.9 3,999.0 990.4 1,355.4 0.5 2,346.3 6,865.7
1981 227.6 0.0 1.2 228.8 3,999.0 236.9 85.4 0.3 322.6 6,865.7
1982 237.0 0.0 1.2 238.2 3,999.0 456.2 715.5 0.0 1,171.7 7,025.7
1983 343.1 29.7 0.0 372.8 4,210.5 373.1 69.5 1.2 443.8 7,004.7
1984 365.8 37.3 0.0 403.1 4,173.5 810.9 684.1 1.3 1,496.3 6,701.7
1985 236.0 29.7 0.0 265.7 4,173.5 234.0 789.0 0.3 1,023.3 5,911.7
1986 433.6 27.0 0.0 460.6 4,173.5 332.3 728.5 0.0 1,060.8 5,875.3
1987 186.7 29.2 0.0 215.9 4,173.5 157.8 599.2 1.7 758.7 6,269.3
1988 207.5 88.5 0.0 296.0 1,639.5 131.6 607.2 0.0 738.8 6,109.3
1989 154.6 69.1 0.0 223.7 1,639.5 299.5 542.0 0.3 841.8 5,551.3
1990 618.3 0.0 0.0 618.3 1,639.5 58.5 624.8 0.2 683.5 5,491.3
1991 114.2 0.0 0.0 114.2 1,639.5 104.3 621.9 17.2 743.4 5,719.3
1992 20.0 82.9 0.0 102.9 1,595.5 78.7 473.4 0.1 552.2 5,518.9
1993 56.3 3,138.5 8.6 3,203.4 1,595.5 416.5 601.4 0.3 1,018.2 5,124.9
1994 100.8 8.7 2.2 111.7 1,610.5 210.2 633.4 0.1 843.7 4,924.9
1995 61.5 4.5 4.6 70.6 1,533.5 323.2 715.9 0.0 1,039.1 4,964.9
1996 32.5 11.4 15.7 59.6 1,533.5 213.3 130.9 0.0 344.2 4,964.9
1997 57.5 7.7 4.7 69.9 1,533.5 161.8 178.0 0.3 340.1 4,787.9
1998 33.3 20.0 2.1 55.4 1,633.5 173.7 213.3 19.2 406.2 4,839.8
1999 72.0 11.0 0.3 83.3 1,633.5 35.7 760.8 13.3 809.8 4,793.8
2000 96.0 7.7 0.0 103.7 1,550.5 170.1 585.6 1.5 757.2 5,371.5
2001 52.6 3.8 0.2 56.6 1,550.5 102.2 335.4 8.5 446.1 5,641.5
2002 10.0 2.7 0.1 12.8 1,550.5 44.4 636.3 2.8 683.5 5,641.5
2003 43.5 5.7 0.4 49.6 1,550.5 998.0 691.0 24.0 1,713.0 5,441.5
2004 35.0 12.9 0.8 48.7 1,520.5 263.1 579.6 0.7 843.4 5,354.5
2005 10.0 25.5 1.2 36.7 1,520.5 43.1 724.8 1.9 769.8 5,354.5
2006 32.5 25.2 0.8 58.5 1,520.5 119.7 333.9 24.3 477.9 5,356.0
2007 10.0 99.4 0.6 110.0 1,520.5 66.5 712.7 1.7 780.9 5,356.0
2008 10.0 130.0 1.3 141.3 1,520.5 389.8 913.9 1.4 1,305.1 5,252.0
2009 10.0 104.6 1.3 115.9 1,520.5 479.8 853.7 1.4 1,334.9 6,740.0
2010 10.0 25.0 1.2 36.2 1,520.5 355.5 1,036.9 1.8 1,394.2 6,740.0
2011 335.5 27.8 1.5 364.8 2,560.5 458.8 1,246.6 1.6 1,707.0 6,413.0
2012 353.1 108.6 1.1 462.8 2,560.5 273.1 653.1 1.6 927.8 6,352.0
2013 215.0 106.3 1.3 322.6 2,560.5 244.6 691.5 1.8 937.9 6,352.0
2014 332.5 107.8 1.8 442.1 3,160.5 345.9 274.5 49.9 670.3 6,352.0
2015 172.0 100.8 2.6 275.4 3,160.5 53.4 618.8 2.7 674.9 6,443.6
2016 158.9 99.0 2.1 260.0 3,160.5 249.8 870.0 2.4 1,122.2 6,443.6
2017 89.0 482.2 0.0 571.2 3,160.5 444.5 382.8 3.6 830.9 5,843.6
Mean annual, 1967–2017 220.4 181.3 5.0 406.7 2,822.7 489.3 680.3 14.5 1,184.0 5,966.4
Mean annual, 1980–2017 174.9 133.4 1.6 309.9 2,349.9 286.9 612.4 5.0 904.2 5,835.8
Table 4.    Annual reported groundwater use from reach 3 and reach 4 of the Washita River aquifer, southern Oklahoma, 1967–2017.

Mean annual groundwater use for the 1980–2017 study period was 1,214.1 acre-ft/yr, which includes 745.8 acre-ft/yr (61.4 percent) for public supply and 461.8 acre-ft/yr (38.0 percent) for irrigation (table 4) (Rogers and others, 2023). Mean annual groundwater use for the period of record (1967–2017) was 1,590.7 acre-ft/yr, which includes 861.6 acre-ft/yr (54.2 percent) for irrigation and 709.7 acre-ft/yr (44.6 percent) for public supply (table 4) (Rogers and others, 2023). The decrease in groundwater use, particularly in irrigation groundwater use after 1980, is likely a result of the OWRB implementing changes in methods of groundwater-use reporting. Since 1980, the OWRB has required irrigation permit holders to estimate and annually report their groundwater use in terms of number of applications and number of inches of water applied (Christopher Neel, OWRB, written commun., 2019). Prior to 1980, however, inches of water applied during irrigation were not required to be reported. As a result, the OWRB adopted rules to estimate the number of inches of water applied for pre-1980 data on the basis of the number of water applications. This estimate results in what appears to be a step-change decrease in irrigation groundwater use after 1980 (fig. 8). In 2017, regular permits accounted for about 58 percent and prior-right permits accounted for about 42 percent of reported groundwater use from the Washita River aquifer (OWRB, 2019d).

Mean annual groundwater use for the 1980–2017 study period was 309.9 acre-ft/yr in reach 3 and 904.2 acre-ft/yr in reach 4 of the Washita River aquifer (table 4) (Rogers and others, 2023). Mean annual groundwater use for the period of record 1967–2017 was 406.7 acre-ft/yr in reach 3 and 1,184.0 acre-ft/yr in reach 4 of the Washita River aquifer (table 4) (Rogers and others, 2023). Groundwater-use permits for public supply were predominantly for wells near Lindsay, Okla., or Pauls Valley, Okla., in reach 4 of the Washita River aquifer, whereas groundwater-use permits for irrigation were for wells scattered throughout the Washita River aquifer (fig. 7). Reach 4 accounted for 74.5 percent of the groundwater use for the 1980–2017 study period and 74.4 percent of the groundwater use for the period of record 1967–2017.

Historical groundwater yields from wells completed in the Washita River aquifer generally decreased from northwest to southeast as the aquifer materials thin; groundwater yields of about 350 gallons per minute (gal/min) were reported near Lindsay and groundwater yields of about 200 gal/min were reported near Davis in the 1950s (Leonard and others, 1958). Driller-estimated groundwater-yield values, which have been mostly recorded since the 1980s, averaged about 70 gal/min; about 50 percent of the values recorded were greater than 25 gal/min (the median value), and about 9 percent of the values recorded were greater than 200 gal/min (OWRB, 2019a). Historical groundwater-yield values in the Washita River aquifer obtained from the NWIS database (USGS, 2020) averaged about 240 gal/min. The median groundwater-yield value was 175 gal/min, and a groundwater-yield value of 1,200 gal/min was reported at one well northwest of Pauls Valley (USGS station 344625097181001). These groundwater-yield statistics likely were biased high because high-yield irrigation and public-supply wells were the most common type of well to have measured discharges (USGS, 2020).

Hydrogeology of the Washita River Aquifer and Surrounding Units

The alluvium and terrace deposits of the Washita River aquifer mostly lie atop Permian-age bedrock units of the northwest-southeast trending Anadarko Basin syncline axis where it approaches the Arbuckle Uplift in southern Oklahoma (fig. 9). The Permian-age bedrock units are generally composed of shale, siltstone, and fine-grained sandstone that serve as boundaries in relation to the alluvium and terrace deposits of the Washita River aquifer (fig. 10). The Arbuckle Uplift deformed and exposed pre-Permian Paleozoic-age sedimentary and igneous units that form the Arbuckle Mountains in the southeastern part of the study area (fig. 9) (Reeds, 1927; Ye and others, 1996). Where they cross the Arbuckle Mountains, the alluvium and terrace deposits of the Washita River aquifer narrow and thin (Kent and others, 1984) and are no longer classified as a major aquifer by the OWRB. In the discussion herein, the geologic units of the study area are presented in reverse chronological order, which is the order in which the units are crossed by the Washita River (upstream to downstream).

Alluvium and Terrace Deposits

Alluvium and terrace deposits of Quaternary age are the principal water-bearing units of the Washita River aquifer (fig. 9). These deposits are made up of alluvial gravel, sand, silt, and clay (figs. 9 and 10; Leonard and others, 1958). Most terrace deposits are hydraulically connected to the valley-bottom alluvium, but some are hydraulically disconnected from the other terrace deposits and the valley-bottom alluvium and thus are unable to store water (Leonard and others, 1958). The disconnected terrace deposits are most evident north of Anadarko and northwest of Pauls Valley (fig. 9). The alluvium and terrace deposits are typically about 2–3 mi wide, but where the Washita River Valley constricts, those deposits are less than 1 mi wide (fig. 11). The Washita River Valley is narrowest where it approaches the Arbuckle Uplift (figs. 9 and 11E). In that area, the alluvium and terrace deposits are only about 0.25 mi wide, which the OWRB regards as the southern boundary of the Washita River aquifer (fig. 9; OWRB, 2019a).

Figure 9. Map showing the Quaternary, Permian, and pre-Permian units of the Washita
                        River aquifer.
Figure 9.

Surficial geologic units and major structural features in A, reach 3 and B, reach 4 of the Washita River aquifer study area, southern Oklahoma.

Figure 10. Stratigraphic chart shows surficial geologic and hydrogeologic units of
                        Washita River aquifer study area, southern Oklahoma.
Figure 10.

Stratigraphic descriptions of surficial geologic and hydrogeologic units of the Washita River aquifer study area, southern Oklahoma.

Figure 11. Cross sections at five locations show altitude of the alluvium, bedrock,
                        and potentiometric surface, Washita River aquifer.
Figure 11.

Hydrogeologic cross sections of the Washita River aquifer, southern Oklahoma, at various locations A–E (fig. 9) in the study area.

Bedrock Units

Bedrock units of Permian age underlie all of reach 3 (fig. 9A) and most of reach 4 (fig. 9B) of the Washita River aquifer. In the northwesternmost part of the aquifer, west of Anadarko, the Whitehorse Group directly underlies the geologic units that contain the Washita River aquifer (fig. 9). The Whitehorse Group consists of the Rush Springs Formation and the Marlow Formation (fig. 10) but is shown in figure 9 as an undivided unit. The Rush Springs Formation is described as orange-brown, cross-bedded fine-grained sandstone, with some dolomite and gypsum beds (Ellis, 2018a; Neel and others, 2018). The thickness of this geologic unit ranges between 186 and 560 ft in the study area (Ellis, 2018a). Stratigraphically below the Rush Springs Formation is the Marlow Formation, which consists of orange-brown, fine-grained sandstone and siltstone interbedded with mudstone, gypsum, and dolomite (Kent and others, 1984; Ellis, 2018a). The Marlow Formation directly underlies the Washita River aquifer near Anadarko and along Sugar Creek (Kent and others, 1984).

Beneath the White Horse Group is the Permian-age El Reno Group, which consists of the Dog Creek Shale, the Chickasha Formation, and the Duncan Sandstone (fig. 10) but is shown in figure 9 as an undivided unit (Heran and others, 2003). The Dog Creek Shale underlies the Washita River aquifer in the western part of Grady County and consists of red, brown, and green gypsiferous shales interbedded with siltstone, sandstone, and dolomite (Kent and others, 1984; Heran and others, 2003; Ellis, 2018a). The Dog Creek Shale has an estimated thickness of 200 ft and is described as having very low permeability (Mogg and others, 1960; Ellis, 2018a), which serves to limit downward flow in the Washita River aquifer. The Chickasha Formation and Duncan Sandstone make up the stratigraphically lower part of the El Reno Group. The Chickasha Formation is characterized as a red-brown mudstone conglomerate with some shale and siltstone and is estimated to be about 200 ft thick (Kent and others, 1984; Ellis and others, 2017). The Duncan Sandstone has a similar thickness to the Chickasha Formation and the Dog Creek Shale and is characterized as a light gray to reddish cross-bedded fine-grained sandstone interbedded with shale (Heran and others, 2003; Ellis and others, 2017). The Duncan Sandstone directly underlies the Washita River aquifer as reach 3 transitions to reach 4 (Heran and others, 2003). The units of the El Reno Group contain the El Reno minor aquifer in the study area.

The Permian-age Hennessey Group directly underlies the Washita River aquifer in the northern part of Garvin County (Kent and others, 1984; fig. 9). The Hennessey Group is made up of the Bison Formation, Purcell Sandstone, and Fairmont Shale (fig. 10) but is shown in figure 9 as an undivided unit. The Bison Formation, estimated to be 50–90 ft thick in the study area, makes up the upper part of the Hennessey Group and is characterized as gray to red-brown shale containing calcium carbonate (Heran and others, 2003). The Purcell Sandstone is a red-brown, fine- to coarse-grained sandstone with some shale and mudstone conglomerate (Heran and others, 2003). The Purcell Sandstone is about 90–150 ft thick in the study area and lies directly above the Fairmont Shale. The Fairmont Shale is similar in color to the Purcell Sandstone, characterized as red-brown and blocky. The Fairmont Shale thickness is 40–80 ft in the study area.

The Permian-age Sumner Group, which includes the Garber Sandstone and Wellington Formation, underlies the Washita River aquifer northwest of Pauls Valley (figs. 9 and 10). The Garber Sandstone is characterized as red-brown fine-grained sandstone with a thickness of about 100–150 ft in the study area (Heran and others, 2003). The Wellington Formation is composed of red-brown shale with several 20–30-ft bituminous sandstones at the base and a thickness of about 100–200 ft decreasing southeastward in the study area (Heran and others, 2003). The Garber Sandstone and Wellington Formation contain the Garber-Wellington aquifer (also known as the Central Oklahoma aquifer, Mashburn and others, 2013) in the northeastern part of the study area, where the Garber Sandstone and Wellington Formation are thicker and more permeable. However, the Garber Sandstone and Wellington Formation are not considered to be a major aquifer near the Washita River aquifer (Hart, 1974; Patterson, 1984).

Underlying the Sumner Group is the Pontotoc Group, a mostly Permian-age red-brown to grayish shale with arkosic sandstone and limestone conglomerates (fig. 10). The Pontotoc Group directly underlies the Washita River aquifer in reach 4 where it approaches the Arbuckle Mountains (fig. 9) and has a relative thickness of 300–500 ft in the study area (Heran and others, 2003). The Pontotoc Group is not sufficiently permeable to yield large quantities of water in most locations and is not considered to be a major water-bearing unit in the study area; however, well yields of 60 gal/min have been reported in a small area in southwestern Garvin County (Hart, 1974).

Pre-Permian Paleozoic-age units make up the Arbuckle Mountains and underlie the southeasternmost part of the Washita River aquifer (fig. 9). These units are composed of gray-brown limestone with shale and some bituminous sandstone and limestone conglomerate (fig. 10) and were greatly faulted and folded by the Arbuckle Uplift (Reeds, 1927; Ye and others, 1996). Some of the pre-Permian Paleozoic-age units form the Arbuckle-Simpson aquifer in southern Oklahoma. It is uncertain, however, whether that aquifer is in direct hydraulic connection with the Washita River aquifer as defined in this report. Igneous units (Wichita Granite and Raggedy Mountain Gabbro Groups [fig. 9], consisting of diallage and olivine gabbro, anorthosite, perthite leucogranites, and diorite composing a layered intrusion [fig. 10]) are present within the study area southwest of Davis; however, these igneous units are not in direct hydraulic connection with the Washita River aquifer.

Groundwater Quality

Groundwater-quality data were collected in the study area during July 28–August 25, 2014, by the OWRB (OWRB, 2017). Groundwater was sampled from 14 wells in the alluvium and terrace deposits of the Washita River aquifer study area (fig. 9; site information is available in Rogers and others [2023]) and analyzed for selected constituents including physical properties (pH, temperature, and specific conductance), major ions, nutrients, and trace metals. Total dissolved solids (TDS) concentrations in the groundwater samples ranged from 215 to 1,070 milligrams per liter (mg/L) with a median concentration of 523 mg/L and a mean concentration of 594 mg/L. The U.S. Environmental Protection Agency (2017) has established a secondary drinking water standard of 500 mg/L for TDS. The State of Oklahoma, however, acknowledges a beneficial domestic use for general use (class II) groundwater with TDS concentrations of less than 3,000 mg/L and limited use (class III) groundwater with TDS concentrations of 3,000–5,000 mg/L (OWRB, 2015). Specific conductance ranged from 316 to 1,879 microsiemens per centimeter at 25 degrees Celsius (µS/cm at 25 °C) with a median value of 955 µS/cm at 25 °C and a mean of 1,051 µS/cm at 25 °C. The hardness ranged from 142 to 1,470 mg/L; the median hardness value was 493 mg/L, and the mean hardness value was 553 mg/L. All groundwater samples were classified as hard water (hardness greater than 180 mg/L) except for two samples that were collected from terrace deposits of northern Garvin County, which had hardness values of 142 and 170 mg/L.

In addition to the groundwater-quality data collected in 2014 from 14 wells by the OWRB (OWRB, 2017), groundwater-quality data from 3 wells sampled during 1947–53 were obtained from the NWIS database (USGS, 2020) (fig. 9). Groundwater from those wells was analyzed for selected constituents including physical properties (pH, temperature, and specific conductance), major ions, and nutrients. TDS concentrations in the groundwater samples ranged from 432 to 7,960 mg/L with a median concentration of 716 mg/L and a mean concentration of 1,471 mg/L. Specific conductance ranged from 195 to 6,990 µS/cm at 25 °C with a median value of 1,120 µS/cm at 25 °C and a mean value of 2,115 µS/cm at 25 °C. The hardness ranged from 160 to 900 mg/L; the median value was 622 mg/L, and the mean value was 219 mg/L. All groundwater samples were classified as hard water except for two samples, which had hardness values of 160 and 180 mg/L, from western Grady County.

Major-ion groundwater-quality data were converted to units of milliequivalents per liter and excluded from further analysis if the difference in the cation-anion balance for the sample was greater than 10 percent; one sample was excluded from the OWRB (2017) data, and two samples were excluded from the USGS (2020) data. The data from the remaining 13 OWRB (2017) and 12 USGS (2020) samples were plotted on a Piper (1944) diagram for visualization of groundwater types and trends (fig. 12). The northern parts of the Washita River aquifer are adjacent to Permian-age bedrock transitioning in the southern parts to the Pontotoc Group, Permian and pre-Permian Paleozoic-age red-brown to grayish shale with arkosic sandstone and limestone conglomerates near the Arbuckle Mountains (figs. 9 and 10). The groundwater samples show a downgradient transition in anions from a sulfate-rich water to a bicarbonate-rich water (fig. 12). Bicarbonate dominates the anions, particularly in the downgradient samples, resulting in the calcium-magnesium bicarbonate water type for a majority of the samples. A calcium-magnesium sulfate-chloride water type was identified in seven of the samples. These samples occurred mostly in the northwestern part of the aquifer (fig. 9A), where interactions with the Permian-age bedrock contribute the sulfate anion through dissolution of gypsum beds (fig. 10).

Figure 12. Piper diagram showing groundwater quality and how it changes from upstream
                        to downstream in the Washita River aquifer.
Figure 12.

Groundwater-quality samples of water produced from the Washita River aquifer, 1947–53 and 2014, southern Oklahoma.

Hydrogeologic Framework of the Washita River Aquifer

A hydrogeologic framework is a three-dimensional representation of an aquifer and how that aquifer interfaces with surrounding geologic units at a scale that captures the regional controls on groundwater flow (Smith and others, 2021). The hydrogeologic framework for the alluvium and terrace deposits of the Washita River aquifer includes updated definitions of the three-dimensional aquifer extent and potentiometric surface, as well as a description of the hydraulic and textural properties of aquifer materials. The hydrogeologic framework was used in the construction of the numerical groundwater-flow model of the Washita River aquifer (Rogers and others, 2023) described in this report.

Aquifer Extent and Thickness

The geographic extent of the Washita River aquifer (fig. 1) was updated from the OWRB (2019b) by using information obtained from relatively new 1:100,000-scale geologic maps (Miller and Stanley, 2004; Chang and Stanley, 2010; Stanley and Chang, 2012). An older 1:250,000-scale geologic map (Heran and others, 2003) was used for a part of the aquifer near Chickasha, Okla., for which finer scale geologic maps were unavailable. Compared to the relatively coarse scale of the older 1:250,000-scale geologic map of the study area, the finer 1:100,000-scale geologic maps showed a narrower extent of the alluvium and terrace deposits that contain the aquifer, so the updated aquifer area presented in this report is smaller than previously documented.

Where present, the top of the Washita River aquifer was defined as the land-surface altitude obtained from a 10-m (horizontal resolution) digital elevation model (DEM) (USGS, 2015) with filled depressions. The altitude of the base of the Washita River aquifer, which was equivalent to the top of the bedrock, was contoured at a 50-ft interval from bedrock depths obtained from drillers’ lithologic logs, well-completion reports, and test-hole data (OWRB, 2019a; USGS, 2020) and direct-push test holes (fig. 13; USGS, 2020). For each of these data sources, the altitude of the base of the aquifer was calculated by subtracting the measured bedrock depth from the land-surface altitude. For consistency, the land-surface altitude was obtained from the 10-m DEM, even when the data source provided a land-surface altitude.

Figure 13. Map showing the base of the Washita River aquifer decreasing in altitude
                        from northwest to southeast.
Figure 13.

The altitude of the base in A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma.

Selected wells and direct-push test holes in the Washita River aquifer were assumed to fully penetrate the aquifer (USGS, 2020; fig. 13), and their depths were assumed to be equal to the bedrock depth. Bedrock depths from these wells and test holes were given the highest priority because they had each been visited and checked by the authors of this report. Other wells and test holes included reports with drillers’ lithologic logs (OWRB, 2019a; fig. 13) that were searched for terms representing consolidated Permian-age bedrock units (such as “redbed” and “bedrock”). The altitude associated with the first occurrence of these terms in the logs was used as the altitude of the aquifer base. For logs that indicated that the well completely penetrated the Washita River aquifer (for example, the well log did not include a bedrock-related term), the altitude associated with the hole depth listed on the lithologic log was assumed to be the maximum (highest) possible altitude of the aquifer base at that location.

The mean thickness of alluvium and terrace deposits in the Washita River aquifer was about 54 ft, estimated by using the land-surface altitude and the aquifer-base altitude defined in this report. The mean thickness of alluvium and terrace deposits in reach 3 was about 56 ft and in reach 4 was about 53 ft.

Potentiometric Surface and Saturated Thickness

Potentiometric-surface maps show the altitude at which the water level would have stood in tightly cased wells at a specified time; the potentiometric surface is usually contoured or spatially interpolated from synoptic water-table-altitude measurements in many wells across an aquifer extent. Potentiometric-surface maps are used to indicate the general directions of groundwater flow in an aquifer. Groundwater generally flows perpendicular to potentiometric contours in the direction of decreasing contour altitudes (Freeze and Cherry, 1979).

In this report, the Washita River aquifer water table was assumed to be equivalent to the potentiometric surface because the aquifer is unconfined (Kent and others, 1984). The potentiometric surface for the Washita River aquifer was mapped by using groundwater-level measurements from selected wells and subtracting the measured depth to water from the land-surface altitude derived from a light detection and ranging (lidar) DEM (USGS, 2019). The selected wells were domestic and irrigation wells that were unused at the time of measurement. Methods described in Cunningham and Schalk (2011) were used to obtain groundwater-level measurements. A January 2020 potentiometric surface (fig. 14) for the Washita River aquifer was contoured primarily from synoptic groundwater-level altitudes measured during January 21–23, 2020, in 66 wells and supplemented with selected water-level altitudes measured mostly during February 2018 and December 1986–January 1987 in 24 wells (USGS, 2020; sites not listed in table 1; instead, site information is available in Rogers and others [2023]). Water-level altitudes measured in 2018 and 1986–87 were used to add control in areas where there was not a measurement in 2020. The mean saturated thickness of the Washita River aquifer was about 47 ft. The mean saturated thickness of reach 3 was about 50 ft and of reach 4 was about 45 ft (Rogers and others, 2023). Local flow in the Washita River aquifer was generally from topographically high areas toward the Washita River and major tributaries with regional flow in the aquifer generally from northwest to southeast (fig. 14).

Figure 14. Map shows potentiometric contours and accompanying groundwater flow arrows
                        as well as locations of water level measurements.
Figure 14.

Potentiometric-surface contours and general direction of groundwater flow in A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, 2020.

Hydraulic and Storage Properties of Aquifer Materials

The distribution and variability of textural and hydraulic properties of aquifer materials, especially the horizontal hydraulic conductivity, were assumed to be the primary controls on groundwater flow in the Washita River aquifer. Multiple methods were used to estimate the range and central tendency of horizontal hydraulic conductivity values in the aquifer. These methods included in-place estimation in test holes, a multiwell aquifer test, and summary of data in drillers’ lithologic logs.

Horizontal Hydraulic Conductivity Estimated in Test Holes

A direct-push Geoprobe hydraulic profiling tool (HPT; Geoprobe Systems, 2015) was used at six Geoprobe test-hole sites (table 5, fig. 1) during June 24–July 15, 2021, to obtain in-place estimates of horizontal hydraulic conductivity by using the McCall (2010) method. This method can be used in the saturated zone to estimate horizontal hydraulic conductivity values between 0.1 and 150 feet per day (ft/d). During direct-push drilling, the HPT injects water into the aquifer material while logging electrical conductivity, injection pressure (corrected for hydrostatic gradient in the saturated zone), and injection flow rate at 0.05-ft depth intervals (Geoprobe Systems, 2015). The Direct Image Viewer (version 3.0) software package (Geoprobe Systems, 2020) was used to calculate discrete horizontal hydraulic conductivity values with depth by using the ratio of the injection flow rate and the hydrostatic-corrected injection pressure. The HPT test holes were drilled to a depth of refusal, which was assumed to be the bedrock contact. The 5,595 discrete (every 0.05 ft) HPT horizontal hydraulic conductivity measurements from all six test holes ranged from 0.1 to 89.0 ft/d (fig. 15) with a mean of 12.2 ft/d (table 5; USGS, 2021). These horizontal hydraulic conductivity data may be underestimated, though; small algal mats from the injection-water tank became entrained in the injection line during logging, which may have partially obstructed flow through the injection line, possibly increasing the measured injection pressure and making the aquifer appear less permeable. Progressive obstruction of the injection line may explain why the first two drilled test holes (TH11 and TH12, table 5, fig. 15A) had greater mean horizontal hydraulic conductivity measurements (22.1 and 26.7 ft/d, respectively) than did the others. Most of the HPT test holes showed generally decreasing electrical conductivity and generally increasing horizontal hydraulic conductivity with depth, indicating that the aquifer materials tend to fine upward from fine to medium sand to silt and clay (fig. 9.2 in Boggs [1995]).

Table 5.    

Hydraulic properties calculated from Geoprobe hydraulic profiling tool logs at direct-push Geoprobe test holes in the Washita River aquifer, southern Oklahoma, 2020–21.

[Data from U.S. Geological Survey (USGS; 2020, 2021). M/D/Y, month/day/year; NAVD 88, North American Vertical Datum of 1988; HPT, hydraulic profiling tool; Kh, horizontal hydraulic conductivity; --, not available or not applicable]

Map name (fig. 1) Field identifier USGS station identifier USGS site name Latitude, in decimal degrees Longitude, in decimal degrees Reach Drilling date (M/D/Y) Land-surface altitude, in feet above NAVD 88 Well or hole depth, in feet HPT logging date (M/D/Y) HPT depth, in feet below land surface HPT depth to water, in feet below land surface HPT mean Kh, in feet per day
TH01 WYNNEWOOD3 343537097085801 02N-01E-36 DCC 1 34.5936 −97.1494 4 10/5/2020 826 40.6 7/13/2021 33.2 3.4 5.7
TH02 PAULSVALLEY1 344011097114301 02N-01E-03 CCB 1 34.6697 −97.1953 4 10/5/2020 864 48.6 -- -- -- --
TH03 PAULSVALLEY2 344053097101001 03N-01E-35 CDD 1 34.6814 −97.1696 4 10/2/2020 851 46.6 -- -- -- --
TH04 LINDSAY3 344842097283401 04N-03W-23 AAA 1 34.8118 −97.4762 4 9/11/2020 958 46.9 -- -- -- --
TH05 LINDSAY1 344916097311401 04N-03W-16 BDA 1 34.8212 −97.5205 4 10/2/2020 961 63.2 7/14/2021 63.0 15.6 6.8
TH06 LINDSAY4 345122097294001 05N-03W-34 DDD 1 34.8562 −97.4946 4 10/2/2020 956 74.1 -- -- -- --
TH07 BRADLEY1 345402097434901 05N-05W-16 DCB 1 34.9005 −97.7302 4 10/6/2020 1,019 49.2 7/14/2021 48.5 4.3 7.1
TH08 BRADLEY3 345422097421801 05N-05W-14 CBB 1 34.9061 −97.7051 4 10/9/2020 1,013 90.2 -- -- -- --
TH09 ALEX3 345749097484101 06N-06W-26 CBB 1 34.9637 −97.8114 3 10/1/2020 1,038 59.2 -- -- -- --
TH10 ALEX2 345817097470602 06N-06W-24 DCC 2 34.9715 −97.7849 3 10/1/2020 1,038 43.4 -- -- -- --
TH11 ALEX1 345829097484101 06N-06W-23 CCB 1 34.9748 −97.8115 3 10/9/2020 1,043 57.5 6/28/2021 58.9 19.2 22.1
TH12 CHICKASHA2 350513098021301 07N-08W-10 CCD 1 35.0870 −98.0369 3 10/1/2020 1,115 63.6 6/24/2021 63.6 2.8 26.7
TH13 ANADARKO2 350513098115501 07N-10W-12 DDD 1 35.0869 −98.1987 3 10/8/2020 1,164 69.2 7/15/2021 68.5 7.4 4.8
TH14 CHICKASHA1 350606098022501 07N-08W-03 CCC 1 35.1017 −98.0403 3 10/9/2020 1,123 49.3 -- -- -- --
TH15 ANADARKO1 350610098114701 07N-09W-06 CCC 1 35.1028 −98.1965 3 10/8/2020 1,165 98.0 -- -- -- --
TH16 CHICKASHA3 350700098022501 07N-08W-04 AAA 1 35.1167 −98.0403 3 10/1/2020 1,114 35.2 -- -- -- --
TH17 ANADARKO3 350701098121901 07N-10W-01 ABB 1 35.1170 −98.2054 3 10/8/2020 1,172 92.7 -- -- -- --
TH18 CHICKASHA4 350753098022001 08N-08W-28 DDD 1 35.1314 −98.0389 3 10/9/2020 1,117 47.3 -- -- -- --
Mean 55.9 8.8 12.2
Table 5.    Hydraulic properties calculated from Geoprobe hydraulic profiling tool logs at direct-push Geoprobe test holes in the Washita River aquifer, southern Oklahoma, 2020–21.
Figure 15. Line graphs show six test holes and electrical conductivity, maximum measurable
                           pressure, and hydraulic conductivity in each.
Figure 15.

Logs of electrical conductivity, maximum hydraulic profiling tool pressures, and horizontal hydraulic conductivity (U.S. Geological Survey, 2021) from selected Geoprobe test holes (U.S. Geological Survey, 2020) in A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, June 24–July 15, 2021.

Horizontal Hydraulic Conductivity Estimated From a Multiwell Aquifer Test

Multiwell aquifer tests are the best method for estimating mean hydraulic and storage properties of aquifer materials over a small area (Stallman, 1971). Data for two unpublished multiwell aquifer tests were found in the USGS paper archives in Oklahoma City, Okla. A multiwell aquifer test in the terrace north of Chickasha was of short duration (24 hours) and failed to induce adequate drawdown in the observation wells (an unknown distance from the pumping well) at a pumping rate of about 60 gal/min (Davis, 1955); this test could not be used for this report. Another multiwell aquifer test, conducted just east of Pauls Valley on March 21–22, 1950, also was of short duration (24 hours) but induced a drawdown of about 7–8 ft in observation wells about 15–30 ft from a production well discharging about 90 gal/min (Rogers and others, 2023). Using the Cooper and Jacob (1946) method, this test yielded a transmissivity value of 922 feet squared per day (ft2/d) (6,900 gallons per day per foot). A horizontal hydraulic conductivity value of about 13 ft/d was calculated by dividing the transmissivity value by the 70 ft of pre-pumping saturated thickness at the aquifer-test site. A specific yield value of about 0.04 (dimensionless) also was calculated from the Pauls Valley multiwell aquifer test. The Pauls Valley multiwell aquifer test was not able to drain the presumed coarser and more porous material at the base of the aquifer, so the calculated transmissivity, horizontal hydraulic conductivity, and specific yield values are likely underestimated for the whole saturated thickness at the aquifer-test site. One short duration multiwell aquifer test is not sufficient to attribute hydraulic and storage properties for an entire aquifer. However, the Pauls Valley multiwell aquifer test is useful as a piece of evidence in characterizing groundwater flow and storage in the aquifer. The data and field notes relevant to the Pauls Valley multiwell aquifer test were scanned and included in the data release that accompanies this report (Rogers and others, 2023).

Horizontal Hydraulic Conductivity Estimated From Lithologic Logs

Drillers’ lithologic logs that penetrated at least 20 ft of the aquifer also were used to estimate the horizontal hydraulic conductivity of aquifer materials. Textural terms in each lithologic log (OWRB, 2019a) were standardized, categorized, and converted to percentage-coarse-material values by using the methods of Mashburn and others (2013). Lithologic logs included terms such as “gravel,” “sand,” “silt,” and “clay” to describe the texture of unconsolidated alluvium and terrace deposits of the Washita River aquifer; however, terms used for lithologic descriptions varied among drillers. To simplify and standardize the lithologic logs, lithologic descriptions of unconsolidated deposits were reclassified into five lithologic categories (clay, silt, fine to medium sand, coarse sand, and gravel) that were assumed to have ranges of percentage coarse material (0–20, 20–40, 40–60, 60–80, and 80–100 percent coarse material, respectively). The midpoint of each range (10, 30, 50, 70, or 90 percent coarse material, respectively) was then assigned to each lithologic depth interval. The percentage-coarse-material value for each lithologic log was computed as the thickness-weighted mean of percentage-coarse-material values assigned to the unconsolidated lithologic categories in the log. The theoretical maximum percentage-coarse-material value for any lithologic log was 90 percent (all gravel), and the theoretical minimum percentage-coarse-material value for any lithologic log was 10 percent (all clay). A total of 260 lithologic logs (OWRB, 2019a) were used for the percentage-coarse-material analysis, and at least 161 of those logs fully penetrated the aquifer. Layers of predominantly gravel were noted in about 20 percent of the lithologic logs analyzed, and clay layers were noted in almost 100 percent of the logs analyzed.

Ellis and others (2017) developed the following equation to characterize the relation between horizontal hydraulic conductivity and the percentage-coarse-material value:

Kh
= (1.25 ×
Ps
) – 12.4
(1)
where

Kh

is the horizontal hydraulic conductivity, in feet per day; and

Ps

is the percentage-coarse-material value.

Equation 1 was used by Ellis and others (2017) to relate core-derived percentage-coarse-material values to horizontal hydraulic conductivity for alluvium and terrace deposits in the Canadian River aquifer of west-central Oklahoma; equation 1 is used in this report to estimate horizontal hydraulic conductivity values for lithologic logs in the Washita River aquifer. Lithologic-log horizontal hydraulic conductivity values of less than 0.1 ft/d were reassigned to that minimum value. The lithologic-log horizontal hydraulic conductivity values ranged from 0.1 to 90.1 ft/d with a mean of 32.9 ft/d.

Textural and Hydraulic Properties From Other Reports

Hemann (1985) measured horizontal hydraulic conductivity and specific yield values for core samples from reach 3 of the Washita River aquifer east of Anadarko by using a laboratory permeameter. The horizontal hydraulic conductivity measurements ranged from 0.1 to 184 ft/d and averaged about 23 ft/d. The specific yield measurements ranged from 0.03 to 0.22 and averaged about 0.10.

A hydrologic investigation and calibrated groundwater-flow model by Ellis (2018a, b) included simulations of groundwater flow in the alluvium and terrace deposits of the Washita River upgradient from and adjacent to the Washita River aquifer as described in this report. The Ellis (2018b) groundwater-flow model of the Rush Springs aquifer used calibrated horizontal hydraulic conductivity values that ranged from 60 to 90 ft/d for the most downgradient alluvium and terrace deposits, which are coincident with the upgradient extent of the study area for this report. The horizontal hydraulic conductivity values were comparable for the Washita River aquifer and Canadian River aquifer (fig. 43B in Ellis [2018a]). However, the Ellis (2018b) groundwater-flow model used calibrated specific yield values of 0.06–0.12 for the most downgradient Washita River alluvium and terrace deposits; those values were considerably less than the calibrated specific yield values (0.15–0.29) used for the most downgradient Canadian River alluvium and terrace deposits (fig. 46B in Ellis [2018a]). Ellis (2018a) related the difference in calibrated specific yield values to the textural difference in the source rocks eroded by (and contributing sediment to) the two streams; the Washita River alluvium and terrace deposits mostly originated from finer grained Permian-age sedimentary units in present-day Oklahoma, whereas the Canadian River alluvium and terrace deposits mostly originated from coarser grained Tertiary-age sedimentary units in present-day Oklahoma, Texas, and New Mexico (Hendricks, 1937).

The vertical anisotropy (ratio of horizontal to vertical hydraulic conductivity) and specific storage values have not been measured in the Washita River aquifer and were assumed to be comparable to those used in simulations of water availability in the nearby Salt Fork Red River aquifer (Smith and others, 2021), which used vertical anisotropy and specific storage values of 10.0 (dimensionless) and 1×10−5 inverse foot, respectively; these values were each within the range suggested by Domenico and Schwartz (1998) for unconsolidated aquifer materials like those of the Washita River aquifer.

Conceptual Groundwater-Flow Model

A conceptual groundwater-flow model is a simplified description (or diagram) of the major inflow and outflow sources (hydrologic boundaries) of a groundwater-flow system and includes a water-budget accounting of the estimated mean flows from those sources for a specified period. A conceptual groundwater-flow model (hereinafter referred to as the “conceptual model”) for the Washita River aquifer was developed to constrain the construction and calibration of a numerical groundwater-flow model (hereinafter referred to as the “numerical model”) for the aquifer. The conceptual-model water budget (table 6, fig. 16A) was used to estimate mean annual inflows to, and outflows from, the Washita River aquifer for the 1980–2017 study period and included a subaccounting of mean annual inflows and outflows for reach 3 (upgradient from the Alex gage) and reach 4 (downgradient from the Alex gage) (fig. 1). The conceptual-model (and numerical-model) aquifer area (modeled aquifer area, fig. 1) totaled 178,938 acres (table 6).

Table 6.    

Conceptual-model water budget of estimated mean annual inflows and outflows for the Washita River aquifer, southern Oklahoma, 1980–2017.

[All water-budget component values are in units of acre-feet per year. Totals may not equal sum of components because of independent rounding. Net streambed seepage, net lateral groundwater flow, and net change in groundwater storage represent the net effect of aquifer inflows and outflows. --, not quantified]

Reach 3 Reach 4 Total Percentage of water budget Notes
Modeled area in cells 4,948 7,231 12,179 -- Reach 3 and reach 4 are upgradient and downgradient, respectively, from the vicinity of U.S. Geological Survey streamgage 07328100 Washita River at Alex, Okla.
in acres 72,697 106,241 178,938 --
in percent 40.6 59.4 100.0 --
Inflow
Water-budget component Recharge 49,071 71,712 120,783 85.5 8.1 inches per year (table 7), or 20.2 percent of mean annual precipitation
Net lateral groundwater flow 6,763 13,797 20,560 14.5 Unknown; calculated as balance of water budget
Total inflow 55,835 85,509 141,344 100.0
Outflow
Water-budget component Net streambed seepage 49,041 75,130 124,171 87.8 Estimated from base-flow data at streamgages (table 3)
Saturated-zone evapotranspiration 6,484 9,476 15,960 11.3 Estimated as 1.33 feet per year1 times 12,000 acres2
Net change in groundwater storage -- -- -- -- Assumed to be a negligible part of water budget
Well withdrawals 310 903 1,213 0.9 Calculated from Oklahoma Water Resources Board reported groundwater-use data (table 4) (Rogers and others, 2023)
Total outflow 55,835 85,509 141,344 100.0
Table 6.    Conceptual-model water budget of estimated mean annual inflows and outflows for the Washita River aquifer, southern Oklahoma, 1980–2017.
1

Reduced slightly from an estimate of 1.46 feet per year by Ellis (2018a).

2

About 12,000 acres along the Washita River stream corridor were classified as wetland by the National Wetlands Inventory (U.S. Fish and Wildlife Service, 2014).

Figure 16. Graphs show the differences in water-budget components between the conceptual
                     model and the calibrated numerical model.
Figure 16.

Estimated mean annual inflows and outflows by water-budget component for the A, conceptual model and B, calibrated numerical model of the Washita River aquifer, southern Oklahoma, 1980–2017.

Estimated mean annual flows to hydrologic boundaries of the conceptual model have varying levels of uncertainty. Where possible, those estimated flows were based on field measurements from the study area and study period. In cases where field measurements were unavailable or too difficult or expensive to obtain, estimated flows of the conceptual model were assumed to be analogous to those of published conceptual models from similar aquifers in Oklahoma (Ryter and Correll, 2016; Ellis and others, 2017, 2020; Smith and others, 2017, 2021; Ellis, 2018a). The “notes” section of table 6 summarizes data sources and assumptions used to construct the conceptual-model water budget for the Washita River aquifer.

Hydrologic Boundaries

Hydrologic boundaries in the conceptual model represent actual sources (inflows) and sinks (outflows) of water to and from the aquifer. Boundaries that act as both inflows and outflows may be referred to as “net inflows” or “net outflows” depending on which flow component dominates.

Recharge

Recharge is the predominant inflow to the Washita River aquifer. Recharge is defined in this report as the amount of precipitation that infiltrates from the land surface through the unsaturated zone and reaches the water table over a given time. This definition of recharge includes irrigation return flows to groundwater. Other processes such as stream seepage or lateral groundwater flow from adjacent hydrogeologic units are not considered recharge and are accounted for separately in the conceptual-model water budget. Recharge rates are controlled by many factors including precipitation rate, land-surface gradient, soil and sediment permeability, evapotranspiration rates, and vegetation cover type. Although recharge rates are difficult to measure because of high spatial and temporal variability, methods that include environmental tracers, physical measurements, streamflow-hydrograph techniques, and computer codes can be used to estimate recharge rates. For this report, a groundwater-hydrograph-based water-table fluctuation (WTF) method (Healy and Cook, 2002) was used to estimate localized recharge rates for 2017–19, and a code-based water-balance-estimation technique (Westenbroek and others, 2010) was used to estimate spatially distributed recharge rates for the 1980–2017 study period.

Water-Table Fluctuation Method

The WTF method was the primary method used to estimate recharge to the Washita River aquifer. The WTF method assumes that rises in groundwater levels in unconfined aquifers that occurred during a relatively short period (hours to a few days) can be attributed to recharge arriving at the saturation zone following a period of precipitation. This method is most appropriately applied to groundwater wells with shallow water tables and hydrographs that display sharp rises in groundwater levels after precipitation (Healy and Cook, 2002). The WTF method cannot account for a steady rate of recharge or recharge from sources other than precipitation. Annual recharge (R), in inches per year, was estimated by using the following equation:

R
=
Sy
× Σ(
Δh
/
Δt
)
(2)
where

Sy

is the specific yield (dimensionless);

Δh

is the rise in groundwater-level altitude, in inches; and

Δt

is the change in time, in years.

Water-level hydrographs from three USGS continuous water-level recorder wells (W01, W02, and W05, table 1; USGS, 2020) were used to estimate annual recharge during 2017–19. W01 was located near Lindsay, W02 was located near Chickasha, and W05 was located near Pauls Valley (fig. 1). The water-level hydrographs from other continuous water-level recorder wells in the study area were not analyzed because they were affected by nearby pumping (in the case of W03, fig. 17A) or because they did not display sharp water-level rises (in the case of W04 and W06, fig. 17B). Annual precipitation data for the 2016–19 period of analysis were obtained from the Chickasha (CHIC) and Pauls Valley (PAUL) climate stations (fig. 1; Oklahoma Mesonet, 2019). Using the daily precipitation record from the nearest climate station and a specific yield of 0.10 from Hemann (1985), the annual recharge estimates for wells W02 and W05 in 2017 were 5.3 inches (in.; about 13.5 percent of the station annual precipitation) and 7.4 in. (about 18.5 percent of the station annual precipitation), respectively (table 7). The annual recharge estimates for wells W02 and W05 in 2018 were 9.0 in. (about 23.7 percent of the station annual precipitation) and 11.5 in. (about 22.8 percent of the station annual precipitation), respectively. In a year-long period from April 1, 2018, through March 31, 2019, the annual recharge estimate at well W01 was 10.3 in. (about 22.6 percent of the station annual precipitation). When those recharge estimates are normalized by the annual mean precipitation of 40.1 in/yr for the 1980–2017 study period and averaged, the resulting station-averaged mean annual recharge for 2017–19 in reach 3 and 4 of the Washita River aquifer was 8.1 in/yr (about 20.2 percent of mean annual precipitation for the 1980–2017 period). Multiplied by the 178,938-acre modeled area and unit converted, the WTF-calculated mean annual recharge for 2017–19 was 120,783 acre-ft/yr; this value was used as the conceptual-model recharge during 1980–2017 (table 6).

Figure 17. Line graph showing depth to water responses to daily precipitation in the
                              Washita River aquifer, August 2016–June 2019.
Figure 17.

Daily precipitation and depth to water in U.S. Geological Survey (USGS) continuous water-level recorder wells completed in A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, August 2016–June 2019.

Table 7.    

Summary of recharge amounts estimated using the water-table fluctuation method for the Washita River aquifer, southern Oklahoma, 2017–19.

[Continuous water-level recorder wells W03 and W04 were not suitable for analysis with the water-table fluctuation method (Healy and Cook, 2002). Dates for the periods are in month-day-year (MM-DD-YYYY) format; --, data not available or not applicable]

U.S. Geological Survey continuous water-level recorder well (fig. 1, table 1)
W02 W01 W05
Annual mean precipitation 1980–2017, in inches per year, southern Oklahoma, Climate Division 8 (National Centers for Environmental Information, 2021) 40.1 40.1 40.1
Climate station (fig. 2, table 1) CHIC PAUL PAUL
Estimated specific yield 0.10 0.10 0.10
Period 1: 01-01-2017 to 12-31-2017
Station annual precipitation, in inches 39.4 -- 40.1
Sum of water-level rises, in feet 4.4 -- 6.2
Recharge, in inches per year 5.3 -- 7.4
Recharge, percent of annual precipitation 13.5 -- 18.5
Recharge, in inches per year, normalized to mean annual precipitation during 1980–2017 5.4 -- 7.4
Period 2: 01-01-2018 to 12-31-2018
Station annual precipitation, in inches 37.9 -- 50.5
Sum of water-level rises, in feet 7.5 -- 9.6
Recharge, in inches per year 9.0 -- 11.5
Recharge, percent of annual precipitation 23.7 -- 22.8
Recharge, in inches per year, normalized to mean annual precipitation during 1980–2017 9.5 -- 9.2
Period 3: 04-01-2018 to 03-31-2019
Station annual precipitation, in inches -- 45.5 --
Sum of water-level rises, in feet -- 8.6 --
Recharge, in inches per year -- 10.3 --
Recharge, percent of annual precipitation -- 22.6 --
Recharge, in inches per year, normalized to mean annual precipitation during 1980–2017 -- 9.0 --
Mean annual recharge during 2017–19, in inches per year, normalized to mean annual precipitation 1980–2017 8.1
Table 7.    Summary of recharge amounts estimated using the water-table fluctuation method for the Washita River aquifer, southern Oklahoma, 2017–19.
Soil-Water-Balance Code

The Soil-Water-Balance code (SWB; Westenbroek and others, 2010) was used to estimate the amount and spatial distribution of daily groundwater recharge to the Washita River aquifer for each month of the 1980–2017 study period. SWB uses a modified Thornthwaite and Mather (1957) soil-water-balance method on a gridded data structure to compute the daily amount of infiltration (minus losses) that exceeds the storage capacity of the plant root zone. The soil-water-balance equation (modified from Westenbroek and others [2010]) has the form of the following equation:

R
= (
P
+
S
+
Ri
) – (
Int
+
Ro
+
Pet
) –
ΔSm
(3)
where

R

is recharge, in inches per day;

P

is precipitation, in inches per day;

S

is snowmelt, in inches per day;

Ri

is surface runoff inflow, in inches per day;

Int

is plant interception, in inches per day;

Ro

is surface runoff outflow, in inches per day;

Pet

is potential evapotranspiration, in inches per day; and

ΔSm

is the change in soil moisture, in inches per day.

SWB requires climate and landscape characteristic data inputs including precipitation, temperature, soil-water storage capacity, hydrologic soil group, surface-water flow direction, and land-cover type (Westenbroek and others, 2010). Each of these inputs was assigned to a user-specified grid that was 600 by 500 cells of 800 by 800 ft each. The landscape inputs were assumed to remain constant during the study period, but climate data inputs varied daily. Daily climate data including precipitation and minimum and maximum air temperature were obtained for 14 selected climate stations (fig. 2, table 1; National Centers for Environmental Information, 2022) in and near the study area and gridded using inverse-distance-weighted (Esri, 2021b) interpolation. Accumulated snowmelt was added on the basis of the daily minimum and maximum air temperatures. Potential evapotranspiration was calculated by using the Hargreaves and Samani (1985) method for a reference latitude range of 34.5–35.4 degrees. Surface runoff was routed downslope by using a flow-direction grid derived from a 10-m DEM (USGS, 2015). Depressions in the DEM were filled by using the ArcGIS Fill tool (Esri, 2021a) to ensure correct routing of surface runoff and to eliminate areas of internal drainage that can result in unrealistically high amounts of recharge. Soil properties (soil-water storage capacity and hydrologic soil group) were derived from the Soil Survey Geographic database (NRCS, 2022), which is an inventory of generalized soil characteristics. Land-cover types (Fry and others, 2011; Multi-Resolution Land Characteristics Consortium, 2011; fig. 3) were used in conjunction with hydrologic soil groups to partition daily precipitation into plant interception (Int) and surface runoff (Ri and Ro) components and assign plant root-zone depths. The root-zone depths for grass/pasture and cropland (the dominant land-cover types for land overlying the aquifer; fig. 3) varied with soil texture but ranged from about 0.8 to 1.5 ft (40 percent of the values used by Westenbroek and others [2010] for permeable glacial deposits in Wisconsin). The soil-water storage capacity, analogous to specific yield in the saturated zone, multiplied by the root-zone depth determined the maximum volume of water available in the root zone. Changes in soil moisture (ΔSm) exceeding the soil-water storage capacity were assumed to be recharge (R) to the saturated zone. Larger root-zone depths resulted in increased evapotranspiration of water from the root zone and decreased recharge; smaller root-zone depths resulted in decreased evapotranspiration of water from the root zone and increased recharge. Recharge from irrigation was not simulated by SWB but was assumed to be negligible given the relatively small amount of irrigation groundwater use in the study area. The SWB model for the Washita River aquifer study area is included in Rogers and others (2023).

The SWB-estimated mean annual recharge for the 1980–2017 study period was about 4.9 in/yr, or about 12.2 percent of the mean annual precipitation of 40.1 in/yr during 1980–2017 (fig. 18A). The minimum and maximum SWB-estimated annual recharge values were about 1.1 in/yr in 2003 and 11.6 in/yr in 2015, respectively. Recharge efficiency (mean monthly recharge as a percentage of mean monthly precipitation) was greatest from November to February, when evapotranspiration was at a minimum; recharge in these months was greater than 20 percent of precipitation (fig. 18B). Recharge efficiency was lowest in July and August, when evapotranspiration was at a maximum; recharge in these months was less than 6 percent of precipitation. Spatially, mean annual recharge for the study period was greatest in areas of alluvium near the Washita River and tributaries (fig. 19). The SWB-estimated mean annual recharge generally increased from northwest to southeast in the Washita River aquifer; the mean annual recharge for reach 4 was about 14 percent greater than mean annual recharge for reach 3 of the Washita River aquifer. Increased recharge from reach 3 to reach 4 is correlated with increased precipitation from northwest to southeast in the study area (Oklahoma Climatological Survey, 2021a, b).

Figure 18. Graphs show annual, monthly recharge estimated with Soil-Water-Balance
                              code, calibrated with model, and annual precipitation.
Figure 18.

A, Annual precipitation with annual recharge estimated by using the Soil-Water-Balance code (Westenbroek and others, 2010) and B, mean monthly precipitation with mean monthly recharge and evapotranspiration estimated by using the Soil-Water-Balance code for the Washita River aquifer, southern Oklahoma, 1980–2017.

Figure 19. Map showing the spatially distributed mean annual recharge over the Washita
                              River aquifer.
Figure 19.

Mean annual recharge estimated by using the Soil-Water-Balance code (Westenbroek and others, 2010) in A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, 1980–2017.

Mean annual recharge rates estimated from SWB were compared to published estimates of mean annual recharge for the Washita River aquifer. Kent and others (1984) estimated a mean annual recharge of 3.72 in/yr, or about 12 percent of the mean annual precipitation of 32.2 in/yr during a 20-year study period (July 1, 1973–July 1, 1993). That estimate was applied across the entire Washita River aquifer in Oklahoma and included reaches that are outside the scope of this report. Schipper (1983) estimated a mean annual recharge of 3.17 in/yr, or about 13 percent of the mean annual precipitation (24.90 in/yr), during the same 20-year study period for a groundwater model of reach 1 of the Washita River aquifer (in western Oklahoma). Ellis and others (2020) estimated a mean annual recharge for reach 1 of the Washita River aquifer (in western Oklahoma) of 3.15 in/yr, or about 14 percent of the mean annual precipitation of 23 in/yr, during a 36-year study period (1980–2015). The mean annual recharge rate of 4.9 in/yr that was estimated by using SWB for the study area in this report is thus within the range of published estimates of recharge for the Washita River aquifer, particularly when calculated as a percentage of the applicable mean annual precipitation.

The SWB-estimated mean annual recharge was, however, about 40 percent less than the WTF-estimated mean annual recharge used for the conceptual model. SWB is a numerical model, and SWB-estimated recharge must, therefore, be checked against and calibrated to field-based measurements or estimates of recharge like those estimated using WTF. In other reports (Smith and others, 2017, 2021; Ellis, 2018a; Ellis and others, 2020), SWB-estimated recharge was calibrated by decreasing root-zone depths until the SWB-estimated mean annual recharge approximated the conceptual-model recharge. That method was also applied in this report to a point but was halted when the mean root-zone depths began to decrease to less than 1.2 ft. Additional adjustment of SWB-estimated recharge was accomplished during calibration of the numerical model as detailed in the “Numerical Groundwater-Flow Model” section of this report.

Saturated-Zone Evapotranspiration

Evapotranspiration is the process by which water is transferred to the atmosphere directly through evaporation and indirectly through plant transpiration; most evapotranspiration occurs at the land surface where precipitation pools as surface water or where it infiltrates the soil unsaturated zone and becomes available to plant root zones (Lubczynski, 2009). The surface-water and unsaturated-zone components of evapotranspiration were not considered to be a part of the conceptual model for the Washita River aquifer because they occur before infiltrating precipitation has reached the saturated zone to become groundwater recharge. An additional component of evapotranspiration, however, occurs in areas of the aquifer where the saturated zone intersects the plant root zone, most commonly in lower lying or wetland areas along streams (Lubczynski, 2009); this component of evapotranspiration (hereinafter referred to as “saturated-zone evapotranspiration”) was an important part of the conceptual-model water budget.

Rates of saturated-zone evapotranspiration are difficult to estimate over a large area such as the study area but are expected to be roughly proportional to (1) the area where the saturated zone intersects the plant root zone, (2) the mean depth to groundwater in that area during the growing season, and (3) the mean rate of transpiration associated with the assemblage of plants in that area. The area where the saturated zone intersects the plant root zone probably is small compared to the entire modeled area and is mostly confined to the typically perennial Washita River stream corridor. About 12,000 acres along the Washita River stream corridor were classified as wetland (land area with frequently saturated or flooded soils [Cowardin and others, 1979]) by the National Wetlands Inventory (U.S. Fish and Wildlife Service, 2014). The saturated-zone component of evapotranspiration was assumed to be active during the growing season (April–October [National Agricultural Statistics Service, 2020; Oklahoma Climatological Survey, 2020]), greatest annually in wet and hot years, and greatest monthly in early summer (Scholl and others, 2005) when precipitation and temperature typically exceed their mean annual values (fig. 5).

By using the assumptions described in the previous paragraph, groundwater outflow by saturated-zone evapotranspiration could be estimated from daily water-level fluctuation data at wells with shallow depths to water according to methods of White (1932). Wells with continuously measured groundwater-level data were not available during the study period in the study area, but gage-height data from the Anadarko, Alex, and Pauls Valley gages (fig. 1) indicated daily declines in stream stage during daylight hours in summer low-flow conditions. These daily declines in stream stage (with rebounds at night) indicated that saturated-zone evapotranspiration was an active process in the Washita River aquifer, but the declines were too small to be accurately measured from the gage-height data. For this reason, the White (1932) methods were not used. Instead, an annual saturated-zone evapotranspiration rate of at least 1.33 feet per year (table 6, reduced slightly from an estimate of 1.46 feet per year by Ellis [2018a]) was assumed to be appropriate for the Washita River aquifer. If about 12,000 acres of wetland area had similar cover and depths to water, this assumed rate would correspond to an annual saturated-zone evapotranspiration outflow of 15,960 acre-ft/yr (table 6). This estimated annual saturated-zone evapotranspiration outflow was allocated to reach 3 and reach 4 in proportion to modeled area instead of wetland area for simplicity and consistency with other calculated components of the conceptual model.

Streambed Seepage

Base flow can be measured directly in streams when the runoff component of streamflow is at or near zero cubic feet per second (Garner and Bills, 2012). When base-flow measurements are collected at multiple locations over a short span of time, they are called synoptic base-flow (seepage-run) measurements. These measurements can be used to calculate net streambed seepage and classify stream reaches as gaining (exhibiting an increase in base flow between the upstream and downstream endpoints of the reach) or losing (exhibiting a decrease in base flow between the upstream and downstream endpoints of the reach). Streambed seepage rates were calculated as the difference in measured base flows (adjusted for tributary inflows) divided by the stream reach length between measurement locations.

Synoptic base-flow measurements were collected by using the methods of Rantz and others (1982) during a period of minimal runoff on January 23, 2018. These measurements, along with data from active USGS streamgages, served to delineate tributary inflow and base-flow conditions across the aquifer in January 2018. Tributaries along the Washita River alluvium contributed base flows ranging from 0.21 to 32.3 cubic feet per second (ft3/s; fig. 20). Tributaries with unmeasured base flows, many of which were visually observed as having minimal or no flows, were assumed to contribute no base flows. The base-flow measurements on the Washita River main stem indicated that there were gaining and losing reaches of the stream. The greatest base-flow gain (12.6 cubic feet per second per mile) occurred over a 5-mi reach south of Davis, and the greatest base-flow loss (–4.8 cubic feet per second per mile) occurred over a 9-mi reach near Lindsay.

Figure 20. Map showing the location of gaining and losing reaches of the Washita River
                           aquifer on January 23, 2018.
Figure 20.

Synoptic streamflow (base-flow) measurements and estimated streambed seepage (base-flow gain or loss) in gaining and losing reaches of A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma, January 23, 2018.

Net streambed seepage terms for the conceptual model were assumed to be net outflows from the aquifer and were estimated from mean annual base flows computed by using the Base-Flow Index code (Wahl and Wahl, 1995) in the USGS Groundwater Toolbox (Barlow and others, 2015) on selected streamgages in the study area. Net streambed seepage in reach 3 for the 1980–2017 study period (49,041 acre-ft/yr, table 6) was roughly estimated as the mean annual base flow at the Alex gage minus the mean annual base flows at the Anadarko gage and selected tributary streamgages: USGS streamgage 07327000 Sugar Creek near Gracemont, Okla.; USGS streamgage 07327550 Little Washita River east of Ninnekah, Okla.; and USGS streamgage 07328070 Winter Creek near Alex, Okla. (fig. 2, tables 1, 3; tributary base flows are summarized in Rogers and others [2023]). Net streambed seepage in reach 4 during 1980–2017 (75,130 acre-ft/yr, table 6) was estimated as the mean annual base flow at the downgradient end of the aquifer (estimated using data from the Dickson gage) minus the mean annual base flows at the Alex gage and tributary streamgages: USGS streamgage 07328180 North Criner Creek near Criner, Okla.; USGS streamgage 07329000 Rush Creek at Purdy, Okla.; USGS streamgage 07329700 Wildhorse Creek near Hoover, Okla.; and USGS streamgage 07329780 Honey Creek below Turner Falls near Davis, Okla. (fig. 2, tables 1, 3; tributary base flows are summarized in Rogers and others [2023]). When available, streamflow data from the 1980–2017 study period were used to estimate mean annual base flows for the conceptual model. Streamflow data from outside the study period were used to estimate mean annual base flows when other data were unavailable. The total net streambed seepage was 124,171 acre-ft/yr and accounted for 87.8 percent of total outflows in the conceptual-model water budget (table 6).

Well Withdrawals

Well withdrawals were assumed to equal the mean annual reported groundwater use for the 1980–2017 study period, or about 1,214 acre-ft/yr (table 4) (Rogers and others, 2023). About 74 percent of the well withdrawals for that period were from permitted wells in reach 4. Well withdrawals were greatest in dry and hot years because more water was required in those years to grow healthy crops. The altitude of the water table generally declines during dry and hot years (especially during extended droughts) and rises during wet and cool years (fig. 4). The degree to which the water table fluctuates annually at a location is related in part to the volume of nearby well withdrawals and the distribution (or concentration) of recharge near that location. No wells with annual water-level measurements made during 1980–2017 were available for the Washita River aquifer, so estimating net storage change in the aquifer was not possible. The net storage change of the Washita River aquifer was assumed to be a negligible component of the conceptual-model water budget.

Lateral Groundwater Flows

No data were available to estimate lateral groundwater flows (including vertical groundwater flows) across hydrogeologic units and political boundaries (administrative sections or reaches) of the Washita River aquifer. The lateral groundwater flows were from two assumed sources: (1) flow across the alluvium or intra-alluvial flow at the most upgradient and downgradient boundaries of the aquifer and (2) the exchange of groundwater between bedrock units and alluvium. Lateral groundwater flows from tributary alluvium and terrace deposits (those not considered to be part of the Washita River aquifer) were assumed to be a negligible part of the conceptual-model water budget. Net lateral flow was calculated as the difference between the aquifer inflows (recharge) and summed aquifer outflows (saturated-zone evapotranspiration, streambed seepage, and well withdrawals) and was used to balance the water budget. The net lateral flows for reach 3 and reach 4 of the Washita River aquifer accounted for 14.5 percent of the total inflows of the conceptual-model water budget (table 6). By comparison, net lateral flows accounted for 44 percent of the total inflows in the conceptual-model water budget for reach 1 of the Washita River aquifer (Ellis and others, 2020) and 65 percent of the total inflows in the conceptual-model water budget for the Canadian River alluvial aquifer (Ellis and others, 2017). Both of those alluvial aquifers, unlike the Washita River aquifer, were in contact with major bedrock aquifers. In the current study area, however, net lateral flows may be disproportionately greater in reach 4 than in reach 3 of the Washita River aquifer because the bedrock geologic units underlying reach 4 include the Garber Sandstone and Purcell Sandstone (members of the Sumner Group and Hennessey Group, respectively), which are among the coarsest lithologies in the study area (figs. 9 and 10).

Conceptual-Model Water Budget

The conceptual-model water budget (table 6) summarized mean water flows exchanged between each hydrologic boundary and the Washita River aquifer for the 1980–2017 study period. Components of the water budget were estimated from analyses of available data; however, where data were not available, assumptions were made by using published analogs as described in the “Hydrologic Boundaries” section. Recharge accounts for 85.5 percent of the conceptual-model inflows to the Washita River aquifer, and net streambed seepage accounts for 87.8 percent of the outflows from the Washita River aquifer. Net lateral groundwater flow (14.5 percent of inflows) and saturated-zone evapotranspiration (11.3 percent of outflows) were the only other components estimated to be greater than 5 percent of inflows or outflows. Well withdrawals accounted for less than 1 percent of conceptual-model water-budget outflows.

Numerical Groundwater-Flow Model

A finite-difference numerical model of the Washita River aquifer was constructed by using the USGS modular groundwater-flow model MODFLOW-2005 (Harbaugh, 2005) with the Newton formulation solver (MODFLOW-NWT, version 1.2.0; Niswonger and others, 2011) for improved solution of problems involving drying and rewetting. In the modular design of MODFLOW, each hydrologic boundary of the conceptual model, such as streambed seepage, recharge, or well withdrawals, is included as a boundary-condition package that, when activated, adds new inflow and outflow terms to the groundwater-flow equation being solved. Data inputs for each package are specified in machine-readable text files. Model space is discretized into cells, and the cell size is the finest resolution at which spatially varying properties (such as land-surface altitude or horizontal hydraulic conductivity) may be represented and varied. Model time is discretized into time steps within stress periods. The stress-period length is the finest resolution at which temporally varying inflows and outflows may be represented and varied, and the time-step length is the finest length of time for which model outputs may be written (Harbaugh, 2005). Selected numerical-model input values were adjusted to calibrate the model to available water-table-altitude and base-flow observations. The calibrated numerical groundwater-flow model inputs, outputs, metadata, directions for use, and ancillary data were published in a USGS data release (Rogers and others, 2023).

Spatial and Temporal Discretization

The model domain (fig. 21) of the Washita River aquifer was spatially discretized into 500 rows; 600 columns; 12,179 active cells measuring 800 by 800 ft each; and a single convertible layer (Harbaugh, 2005) based on the hydrogeologic framework described in this report. The cell size was chosen to minimize model-processing time while still representing the variability of properties being simulated. The chosen cell size also ensured that the narrowest part of the aquifer was represented by no fewer than three cells; model instability and water-level volatility often occur in narrow parts of the modeled area where groundwater flow is focused into fewer cells. The single convertible layer represented the undifferentiated Quaternary-age alluvium and terrace deposits with variable thickness determined from the hydrogeologic framework; the underlying bedrock was not represented as a layer. The altitude of the top of the aquifer was multiplied by 1.01 in the numerical model to prevent confined aquifer conditions that occur as a side effect of model discretization when the simulated water-table altitude exceeds the altitude of the top of the aquifer. This multiplier, which adds about 10 ft to the altitude of the top of the aquifer, has no effect on most parts of the model where the simulated water-table altitude is below the altitude of the top of the aquifer.

Figure 21. Maps show streamflow routing, evapotranspiration and recharge area, wells,
                        general-head boundaries of active modeled area.
Figure 21.

Model domain, active modeled area, hydrologic boundaries, and parameter zones for the numerical groundwater-flow model of A, reach 3 and B, reach 4 of the Washita River aquifer, southern Oklahoma.

The active modeled area (fig. 21) was initially derived from the Washita River aquifer extent (modified from OWRB [2019b]) as defined in the hydrogeologic framework of the Washita River aquifer. The active modeled area was expanded or contracted in some areas to ensure that each active cell was in connection with at least one other active cell. The active modeled area was further modified by inactivating some model cells in narrow areas along selected Washita River tributaries to improve model stability.

The numerical model was temporally discretized into 456 monthly transient stress periods (each with two time steps to improve model stability) representing the 1980–2017 study period (the numerical-modeling period [hereinafter referred to as the “modeling period”] coincides with the 1980–2017 study period). An initial steady-state stress period represented mean annual inflows to, and outflows from, the aquifer. The steady-state solution was used as the initial condition for subsequent transient stress periods, as well as some groundwater-availability scenarios. The numerical model was constructed by using length and time units of feet and days, respectively.

Simulation of Hydrologic Boundaries and Hydraulic Parameters

Hydrologic boundaries in the numerical model (fig. 21) define where and how water may enter or leave the model and include specified-flux and head-dependent boundaries (Harbaugh and others, 2000). Specified-flux boundaries were used to simulate recharge and well withdrawals. Head-dependent boundaries were used to simulate streambed seepage, saturated-zone evapotranspiration, and lateral groundwater flow exchanged with adjacent (upstream and downstream) alluvium and terrace deposits of the Washita River aquifer and with the surrounding bedrock units. When available, hydrologic data (along with data-based assumptions and analogs) were used to estimate or constrain precalibration model inputs for each hydrologic boundary.

Recharge and Distribution With the Soil-Water-Balance Code

Recharge to the Washita River aquifer was simulated using the Recharge package (Harbaugh and others, 2000). Recharge was spatially and temporally distributed for each month of the 1980–2017 study period by using outputs from SWB (Westenbroek and others, 2010). The SWB-output monthly recharge grids were converted to model units of feet per day and used as precalibration numerical-model inputs. These inputs were then scaled with multipliers during the numerical-model calibration to improve agreement between the simulated recharge rate and the conceptual-model recharge determined by the WTF method. The initial steady-state recharge-rate multiplier (rch001, table 8) of 1.66 was allowed to vary within about 3.0 percent (between 1.61 and 1.71); the initial transient recharge-rate multipliers (rch002–rch457, table 8) of 1.66 for each month in the period 1980–2017 were allowed to vary between 1.40 and 2.00 (Rogers and others, 2023). A narrow range was required for the steady-state recharge-rate multiplier to keep the numerical-model recharge closely aligned with the conceptual-model recharge.

Table 8.    

Calibration parameters for the numerical groundwater-flow model of the Washita River aquifer, southern Oklahoma, 1980–2017.

[A more detailed and longer version of this table is available in the associated data release (Rogers and others, 2023). GHB, general-head boundary; NAVD 88, North American Vertical Datum of 1988]

Parameter group Number of parameters Parameter names Parameter descriptions and units
evt 1 evtavg Steady-state saturated-zone evapotranspiration-rate multiplier (dimensionless)
12 evtjan–evtdec Transient (monthly) saturated-zone evapotranspiration-rate multiplier (dimensionless)
1 evtextd Saturated-zone evapotranspiration extinction (root-zone) depth, in feet
1 evttop Saturated-zone evapotranspiration reference altitude multiplier (dimensionless)
ghb 1 ghcup GHB conductance, in feet squared per day, boundary with upgradient alluvium and terrace deposits
3 ghc3a–ghc3c Reach 3 GHB conductance, in feet squared per day, boundary with bedrock units
5 ghc4a–ghc4e Reach 4 GHB conductance, in feet squared per day, boundary with bedrock units
1 ghcdn GHB conductance, in feet squared per day, boundary with downgradient alluvium and terrace deposits
1 ghdup GHB altitude in feet above NAVD 88, boundary with upgradient alluvium and terrace deposits
3 ghd3a–ghd3c Reach 3 GHB depth, in feet, boundary with bedrock units
5 ghd4a–ghd4e Reach 4 GHB depth, in feet, boundary with bedrock units
1 ghddn GHB altitude in feet above NAVD 88, boundary with downgradient alluvium and terrace deposits
hyd 146 hkpp001–hkpp146 Horizontal hydraulic conductivity, in feet per day, pilot points
rch 1 rch001 Steady-state recharge-rate multiplier (dimensionless)
456 rch002–rch457 Transient (monthly) recharge-rate multipliers (dimensionless)
sfr 28 sbk01–sbk28 Streambed hydraulic conductivities, in feet per day
28 sbb01–sbb28 Streambed thicknesses, in feet
28 sbd01–sbd28 Streambed depths, in feet
28 sbw01–sbw28 Streambed widths, in feet
sto 1 sy3 Reach 3 specific yield (dimensionless)
1 sy4 Reach 4 specific yield (dimensionless)
1 ss Specific storage, in inverse feet
Total 753
Table 8.    Calibration parameters for the numerical groundwater-flow model of the Washita River aquifer, southern Oklahoma, 1980–2017.

Lateral Groundwater Flows

Lateral groundwater flows between the general-head boundaries (GHBs) of the Washita River aquifer and the adjoining alluvium and terrace deposits and bedrock units were simulated by using the General-Head Boundary package (Harbaugh and others, 2000). The flow to or from a GHB is the product of the GHB conductance and the difference between the water-table altitude in the model cell and the specified GHB for that model cell. The GHB altitudes along the northwestern and southeastern model boundaries (where the aquifer is adjoined by the continuation of Washita River aquifer alluvium and terrace deposits) were estimated at about 1,200 and 740 ft, respectively, on the basis of the 2020 potentiometric surface of the aquifer (fig. 14). The altitudes specified for the GHBs remained constant throughout the simulation period. Lateral groundwater flows across the northwestern numerical-model boundary (“up,” fig. 21A) were likely net aquifer inflows, and lateral groundwater flows across the southeastern numerical-model boundary (“dn,” fig. 21B) were likely net aquifer outflows. For simplicity, GHBs simulating lateral groundwater flow exchange with the bedrock units (3a–3c and 4a–4e, fig. 21) were only assigned to model cells with side faces contacting those bedrock units; model cells that only contacted bedrock units through their bottom faces were not assigned a GHB. The GHB altitudes of those model cells were initially assigned an altitude 10 ft below land surface (assigned from a 10-m DEM [USGS, 2015]) and were allowed to vary in the range of 3–15 ft below land surface (ghd3a–ghd3c and ghd4a–ghd4e, table 8; Rogers and others, 2023). GHB conductance values were adjusted by zones during the numerical-model calibration and were allowed to vary in the range of 1,000–8,000 ft2/d for the contacts with adjoining alluvium and terrace deposits (ghcup and ghcdn, table 8; Rogers and others, 2023). GHB conductance values were allowed to vary in the range of 0.1–6,000 ft2/d for the bedrock contacts with the lower Hennessey Group and the Sumner Group, which included some of the coarsest bedrock lithologies (figs. 9 and 10) in the study area, and 0.1–2,000 ft2/d for the bedrock contacts with other units, which included generally finer bedrock lithologies (ghc3a–ghc3c and ghc4a–ghc4e, table 8; Rogers and others, 2023).

Streams

Selected named streams listed in the National Hydrography Dataset Plus (Horizon Systems Corporation, 2015) were simulated by using the Streamflow-Routing package, version 2 (SFR2) (Niswonger and Prudic, 2005). Only base flow was simulated in the SFR2 streams. Inflows for SFR2 stream cells included base flows routed from upstream segments, specified inflows (base flows) from tributary streams, and streambed seepage from the aquifer; outflows for SFR2 stream cells included base flows routed downstream and streambed seepage to the aquifer. SFR2 computes streambed seepage between the aquifer and the stream according to Darcy’s Law (Darcy, 1856; Bennett, 1976); the flow exchanged between the aquifer and the stream is the product of the streambed conductance and the difference between the water-table altitude and the stream stage, where the streambed conductance is the product of the hydraulic conductivity of the streambed sediments and the area of the stream channel divided by the streambed thickness. In SFR2, simulated base flows are calculated in each model cell and are routed downstream by segments, or group of cells with uniform hydraulic properties. Accumulated base flows are passed to the next downstream segment until flows are routed out of the active modeled area.

All SFR2 stream segments (1–28, fig. 21) were initially assigned a streambed thickness (sbb01–sbb28, table 8) of 2.0 ft and a streambed hydraulic conductivity (sbk01–sbk28, table 8) of 10.0 ft/d (based on Ryter and Correll [2016]; Ellis and others [2017]; Smith and others [2017, 2021]), but these properties were adjusted during calibration within ranges of 1.0–7.0 ft and 2.0–12.0 ft/d, respectively (Rogers and others, 2023). The streambed widths of stream segments (sbw01–sbw28, table 8) were estimated from aerial photographs (Microsoft Corporation, 2018) and ranged from about 40 to 600 ft (Rogers and others, 2023), gradually increasing downstream. The channel widths of some streams were increased to as much as 2.5 times the estimated values during calibration; increased stream widths resulted in decreased simulated stream stages and corresponding increases in the flow from the aquifer to the stream, allowing parts of the model to drain more efficiently. The streambed depth below land surface (sbd01–sbd28, table 8) was initially set to about 4.5–7.0 ft and was adjusted during calibration within the range of 1.0–8.0 ft (Rogers and others, 2023). The streambed depth was subtracted at model run time from the land-surface altitude as represented by a 10-m DEM (USGS, 2015). This streambed incision is necessary to compensate for lost spatial resolution of the stream channel (altitude averaging) caused by the large DEM cell size. When land-surface altitudes of actual features are represented by a DEM cell, the lowest altitudes (often stream channels) and highest altitudes (summits) are generalized to a single value—the mean land-surface altitude in the cell. The difference between the actual altitude of those features and the cell-averaged altitude depends on the local altitude relief and the cell size.

Monthly base-flow inflows at the boundary of the active modeled area were estimated for the Washita River (WASH), Sugar Creek (SUGA), Rock Hollow Creek (ROCK), Line Creek (LINE), West Bitter Creek (WBIT), East Bitter Creek (EBIT), Little Washita River (LITW), Roaring Creek (ROAR), Bear Creek (BEAR), Hybarger Creek (HYBA), Criner Creek (CRIN), Rush Creek (RUSH), Dry Sandy Creek (DSAN), Wildhorse Creek (WILD), and Honey Creek (HONE) (fig. 21). Streamgage records were used to estimate base-flow inflows where available (table 3). Where streamgage records were unavailable, the January 2018 synoptic streamflow (base-flow) measurements (fig. 20) were used to estimate some base-flow inflows. Base-flow inflows for some streams were estimated from streamgage records and base-flow measurements on similar streams nearby (USGS, 2020). Rock Creek, a relatively large tributary near the downgradient end of the model, was not simulated because it was dominated by periodic reservoir releases from Lake of the Arbuckles (fig. 21).

Saturated-Zone Evapotranspiration

Saturated-zone evapotranspiration was simulated by using the Evapotranspiration package (Harbaugh and others, 2000) and was expected to occur near streams and riparian wooded areas along the Washita River and perennial tributaries. Maximum rates of saturated-zone evapotranspiration were assumed to not exceed the difference between potential and actual evapotranspiration as computed by SWB; this assumption prevented the summed components of evapotranspiration from exceeding the potential evapotranspiration. Arrays representing the potential minus actual evapotranspiration for each monthly transient stress period (evtjan–evtdec, table 8) and the steady-state stress period (evtavg, table 8) were initially scaled by a factor of 0.6. During calibration, the scaling factors were allowed to vary independently in the range of 0.5–1.0 (Rogers and others, 2023). The evapotranspiration extinction (root-zone) depth (evtextd, table 8), or the depth below land surface at which the saturated zone becomes inaccessible to plants, was initially set at 1.0 ft in the active modeled area, which was consistent with the mean plant root-zone depth specified for grass/pasture and cropland (the dominant land-cover types, fig. 3) in SWB. During calibration, the evapotranspiration extinction depth was allowed to vary in the range of 0.8–1.2 ft. The saturated-zone evapotranspiration reference altitude (evttop, table 8) was initially set to the land-surface altitude of each cell but was adjusted during calibration by a multiplier that was allowed to vary in the range of 0.998–1.010. The narrow range on the reference altitude multiplier was necessary for model stability.

Well Withdrawals

Well withdrawals for the 1980–2017 study period were simulated by using the Well package (Harbaugh and others, 2000). Annual reported groundwater use for each permit was evenly distributed among all well locations recorded for that permit (OWRB, 2019d). In the case of permits for which there were no recorded wells, a single well was placed in the center of the land parcel recorded for the permit. Annual reported groundwater use was temporally split into model stress periods by using the mean monthly water demand distribution (fig. 22; OWRB, 2012a) from Oklahoma Comprehensive Water Plan water-management planning basins 14–16 (fig. 2; OWRB, 2019c). The monthly water demand for irrigation was greatest in the summer months; July–September accounted for about 76 percent of irrigation groundwater use (fig. 22). The monthly water demand for public supply also was greatest in the summer months; however, July–September only accounted for about 34 percent of public-supply groundwater use (fig. 22).

Figure 22. Graph showing irrigation drastically outpacing public supply for monthly
                           water demand during the months of May–August.
Figure 22.

Monthly water demand by groundwater-use type for Oklahoma Comprehensive Water Plan water-management planning basins 14–16, averaged, Washita River aquifer, southern Oklahoma.

Storage and Hydraulic Properties

The Upstream Weighting package of MODFLOW-2005 (Niswonger and others, 2011) was used to represent storage and hydraulic properties of the Washita River aquifer. Initial storage and hydraulic property values were similar to those used in numerical models for other alluvial aquifers in Oklahoma (Ellis and others, 2017, 2020; Smith and others, 2017) and identical to those used in the numerical model for the Salt Fork Red River aquifer (Smith and others, 2021). The storage and hydraulic property values were adjusted during model calibration but were held spatially uniform and temporally constant through all stress periods in all simulations. The specific yield (sy3 and sy4, table 8) and specific storage (ss, table 8) of the Washita River aquifer were initially set as 0.10 (dimensionless) and 1×10–5 inverse foot, respectively. Horizontal hydraulic conductivity (hkpp001–hkpp146, table 8) was represented by using 146 pilot points (1–146, fig. 23; Doherty, 2010) placed at 12,800-ft intervals (1 at every 16th column and row) within a 10,560-ft buffer of the aquifer boundary. (Pilot points are regularly spaced control points used to support an interpolation of cell-based horizontal hydraulic conductivity values across an active modeled area.) The pilot points were assigned an initial horizontal hydraulic conductivity value of 20 ft/d and interpolated before each model run to create the horizontal hydraulic conductivity array read by the numerical model (Doherty, 2010). During calibration, pilot-point horizontal hydraulic conductivity values were allowed to vary independently between 10 and 100 ft/d (Rogers and others, 2023).

Figure 23. Map shows that the highest hydraulic conductivity mainly occurs in the
                           northwestern part of the Washita River aquifer.
Figure 23.

Horizontal hydraulic conductivity pilot points and calibrated horizontal hydraulic conductivity in A, reach 3 and B, reach 4 for the numerical groundwater-flow model of the Washita River aquifer, southern Oklahoma.

Solver Settings and Budget Percentage Discrepancies

Most of the settings for the Newton solver were unchanged from suggested input values (Winston, 2018), which were converted to model units of feet and days. The head tolerance was increased to 1.5×10–3 ft, and the flux tolerance was increased to 5.0×104 cubic feet per day. These settings improved model stability while keeping the budget percentage discrepancies in the solution under 1.0 percent for all 457 stress periods; the largest budget percentage discrepancy (by absolute value) was about 0.16 percent in stress period 411 (Rogers and others, 2023).

Calibration

Model calibration is the process of systematically changing selected model input values (parameters) within limits to improve the fit between model-simulated data and observed data (calibration targets). The preferred calibration results (1) minimize the differences (residuals) between simulated and observed data and (2) conform to the predetermined conceptual model. The calibration process for the numerical model included both manual and automated adjustment of parameters. The manual calibration approach was primarily focused on aligning the numerical-model water budget, especially the recharge and streambed seepage components, to the conceptual-model water budget. The manual calibration also involved selected structural elements of the model, such as boundary elevations, that are not easily adjusted by automated calibration methods. The automated calibration approach was focused solely on minimizing residuals and used the PEST++ iterative ensemble smoother (White, 2018) to reduce run times associated with the calibration of highly parameterized models.

Calibration Targets

The suite of calibrated parameter values was evaluated on the basis of the minimization of an objective function that quantifies the difference between simulated and observed values included in the model (White, 2018). The objective function was calculated as the sum of squared weighted residuals for calibration targets in five observation groups: water-table-altitude observations, base-flow observations at the Alex gage, base-flow observations at the Pauls Valley gage, base-flow observations (estimated) at the downgradient end of the model, and the conceptual-model recharge rate (table 9). Observation weights were required to prevent calibration targets of any one observation group from dominating the objective function and, therefore, becoming the sole focus of parameter adjustments during automated calibration (Hill and Tiedeman, 2007; Doherty, 2010). The observation weights (1.0, 5.0 × 10–7, 1.7 × 10–7, 1.8 × 10–7, 2.5 × 10–2 for water-table-altitude observations, base-flow observations at the Alex gage, base-flow observations at the Pauls Valley gage, base-flow observations [estimated] at the downgradient end of the model, and conceptual-model recharge rate, respectively) were uniform for all observations from the transient simulation within each observation group; the base-flow observations from the steady-state simulation were given increased weights (200 times greater) relative to the base-flow observations from the transient simulation. The observation weights were chosen so that, at the beginning of automated calibration (that is, precalibration), the observation groups accounted for 24.0, 10.5, 55.3, 10.2, and 0.0 percent of the objective function, respectively (table 9). This weighting scheme reflects an assumption that base-flow observations (totaling 76.0 percent of the objective function) are more important than water-table-altitude observations for the purpose of determining the parameter values that are most influential in estimating groundwater availability (the primary subject of this report) in an unconfined alluvial aquifer. Base flows are generally more sensitive to changes in specific yield than they are to changes in water-table altitudes. Water-table altitudes are generally more sensitive to changes in horizontal hydraulic conductivity than to changes in specific yield.

Table 9.    

Components of the objective function for the automated calibration of the numerical groundwater-flow model of the Washita River aquifer, southern Oklahoma, 1980–2017.

[Objective function is calculated as the sum of squared weighted residuals. SS, steady-state simulation; TR, transient simulation; Min, minimum; Max, maximum; RMSE, root-mean-square error; OWRB, Oklahoma Water Resources Board; USGS, U.S. Geological Survey; E, base-10 exponent (for example, 1.7E-07 means 1.7 × 10−7)]

Observation group Source1 Number of observations Weight multiplier Precalibration residuals Calibrated residuals
SS
TR
SS
TR
Min Mean Max Interquartile range RMSE Objective function component
(percent in parentheses)
Min Mean Max Interquartile range RMSE Objective function component (percent in parentheses)
Water-table altitude, in feet (headobs) Water-table-altitude observations (OWRB, 2019a; USGS, 2020) 63
130
2.0E+00
2.0E+00
−22.5 −6.5 17.2 9.1 18.6 66,641
(24.0)
−19.5 −4.1 10.5 6.8 13.8 36,522
(24.9)
Base flow, in cubic feet per second (gageobs1) Alex gage (USGS, 2020) 1
432
1.0E-04
5.0E-07
−1,238.8 14.0 970.2 95.6 8.2 28,993
(10.5)
−1,228.1 −7.7 981.8 93.9 7.9 27,109
(18.5)
Base flow, in cubic feet per second (gageobs2) Pauls Valley gage (USGS, 2020) 1
456
3.4E-05
1.7E-07
−1,207.7 105.4 2,549.3 213.2 18.3 153,413
(55.3)
−1,224.2 45.6 2,489.8 213.8 10.8 53,763
(36.6)
Base flow, in cubic feet per second (gageobs3) Downgradient end of model; adjusted from Dickson gage (USGS, 2020) 1
456
3.6E-05
1.8E-07
−1,457.2 8.8 4,952.7 234.1 7.9 28,231
(10.2)
−1,433.3 −8.8 4,932.0 234.1 8.0 29,173
(19.9)
Conceptual-model recharge, in acre-feet per year (budget) Conceptual model (table 6) 1
1
2.0E-02
2.0E-02
353.7 353.8 353.8 0.0 7.1 100.1
(0.0)
−829.8 −416.5 −3.2 413.3 11.7 275.4
(0.2)
Total 1,542 Precalibration objective function: 277,378 Calibrated objective function: 146,842
Table 9.    Components of the objective function for the automated calibration of the numerical groundwater-flow model of the Washita River aquifer, southern Oklahoma, 1980–2017.
1

Alex gage refers to USGS streamgage 07328100 Washita River at Alex, Okla.; Pauls Valley gage refers to USGS streamgage 07328500 Washita River near Pauls Valley, Okla.; and Dickson gage refers to USGS streamgage 07331000 Washita River near Dickson, Okla. (fig. 1, table 1).

Water-Table-Altitude Observations

The Head Observation (HOB) package (Hill and others, 2000) was used to compare simulated water-table altitudes to observed water-table altitudes in the Washita River aquifer for the 1980–2017 modeling period. Observed water-table altitudes were calculated from depth-to-water measurements and the land-surface altitude obtained from a lidar DEM (USGS, 2019). For consistency, the land-surface altitude was obtained from the lidar DEM, even when the data source provided a land-surface altitude. Sites with water-table altitude observations were screened for location errors, and observations were filtered such that only one observation per cell per stress peri