Oklahoma Watershed Response Model Marc Petrequin GEGN 582 – Integrated Surface Hydrology Professor Reed Maxwell, Department of Hydrology, Colorado School of Mines December 9th, 2011 Abstract This project centered on the development of a linear reservoir response model for an unspecified watershed in Oklahoma. I was tasked with analyzing a rainfall hyetograph and corresponding watershed outflow from a past storm event, determining a range of appropriate parameters from this analysis, using these parameter to develop a series of different linear watershed models, calibrating these models to match the given hyetograph of the past rainfall event, and applying these models to predict the outflow of a future storm event. A significant portion of this project required assumptions inferred from the brief watershed description and the observed characteristics of the outflo ur e of the atershed’s respo se. My primary sources of methodology were the assigned texts for this course: Physical Hydrology, Second Edition by S. Lawrence Dingman and Hydrology: An Introduction by Wilfried Brutsaert, as well as the report Hydrographs by Single Linear Reservoir Model by Pedersen et al of the U.S. Army Corps of Engineers Hydrologic Engineering Center.            Assumptions The average daily precipitation (W) in the watershed for both datasets is roughly evenly distributed across the entire area. The barometric pressure in Little Park, Oklahoma, a state with an average elevation of 1300 feet above sea level, remained at a relatively constant 98.0 kPa for the given time periods. The land use in the watershed primarily constitutes pasture and rangeland cover with a small amount (roughly 20%) constituting row crops (straight row, good hydrologic conditions). The amount of agriculture and grazing lands in the watershed and its initial soil water content are significant enough to classify the surface as well-watered, and equilibrium conditions remain in effect for for purposes of determining the potential evaporation content (αPT = 1.26). For purposes of determining the SCS Curve Numbers, the soil in the watershed falls into hydrologic soil group B due to its relatively high hydraulic conductivity and storage capacity, and because the outflow at the end of the given time interval is slightly below that at the beginning (before the peak inflow), there was a rainfall event within the previous 5 days before the recording that places the soil wetness under Condition II ("average"). The tree cover in the watershed is evenly distributed. The watershed has a constant incoming base flow, which can be estimated using the straightline method in Dingman. This base flow is roughly the same for both time periods and does not contribute to basin storage. The initial water content of the watershed is roughly the same for both storm events. All increases to the recorded outflow after the beginning of the effective water input occur from storm events within the watershed. All reservoirs and streams in the watershed have the same response time. Groundwater seepage into watershed reservoirs and streams is negligible. O l the al ulated effe ti e pre ipitatio o tri utes to the storage of the atershed. Part 1: Quantitative Description The first phase of this project involved a rudimentary quantitative analysis of given data regarding a past storm event with the incorporation of concepts covered earlier in the semester. Computing the cumulative outflow required unit conversion to m3/d and summing the total amount. Deriving an estimate of the total effective rainfall (Weff), conversely, required a broad range of assumptions based on the given information. Using the SCS method illustrated in Page 445 of Dingman, I deri ed a hart of the atershed’s land cover and soil conditions based on the provided description (regarding grassland, agriculture, and forest), a d fro there deri ed a esti ate of the atershed’s Maximum Retention Capacity (Vmax). Combining this with the given total observed rainfall gave me an effective rainfall rate, which I multiplied by the entire watershed area to derive the cumulative rainfall amount. Because this method incorporates all storage in the watershed, using it made it easier to accommodate for additional watershed components that would diminish the effective rainfall and runoff amount (i.e. groundwater percolation into perched aquitards). Determining the total change in watershed storage Δ“ , as defined by Equation 4.55 in Brutseart, required a rough estimate of the daily evapotranspiration rate (ET) and runoff (R) (the amount of outflow in excess of base flow). Due to a lack of information on the local meteorological conditions, such as the average wind velocity and pressure gradient, I chose to calculate this amount based on the estimated potential evapotranspiration (PET) and from that a predicted actual evapotranspiration, as demonstrated in Chapter 7 of Dingman. The radiation-based equilibrium potential method employed by Priestley and Taylor (1972) seemed to fit best with the given data, which provided the average net longwave and shortwave radiation (K+L) for each day. In contrast, the Hamon method based on day length required more thorough use of trigonomic functions than EXCEL™ was capable of, and the Malmstrom method based on the mean monthly air temperature and vapor pressure seemed too broad an estimate for a time period of only 60 days. Wishing to maximize accuracy wherever possible, I derived for each day the pressure-temperature slope (Δ), the water density (ρw), the latent heat of vaporization (λv), and the psychometric constant (γ) using empirical methods based on the given temperature (T), and, assuming a constant alpha factor (αPT), the potential evapotranspiration rate. To convert this potential evapotranspiration rate to actual, a lack of reliable data on the average annual precipitation and soil water content led me to rely on the complementary (advection-aridity) approach used in Equation 7-70 of Dingman. Already relying on the assumption proposed by Brutsaert and Stricker that the PET is under equilibrium conditions (as employed in the Priestly and Taylor method), I used the net incoming radiation as the primary source of the energy flux and assumed such energy represented the difference between potential and actual. Naturally, this led to certain days having a negative evapotranspiration reading, primarily during or following rainfall, when a larger portion of the incoming solar energy would be absorbed by the increased surface water content or reflected by the increased cloud cover. To reduce the inaccuracy in the total evapotranspiration that ould e aused a assi e de poi t he applying such sparse negative amounts to the entire watershed, I coded the spreadsheet to simply return a 0 for any negative ET values, essentially assuming perfect equilibrium. Finally, determining the net change in storage also required determining the total runoff from the area, as opposed to the total outflow. This outflow over the recorded time appeared to taper off and steadily recede at a steady rate about 15 days after the end of the initial storm, to the point it dropped below the initial outflow in spite of some periods of light precipitation following the storm. This suggested that 1) the hydrograph was initially receding from a storm event before the recorded time period, 2) the periods of light precipitation following the initial storm event were too small to contribute to storm runoff, and 3) the watershed was undergoing a net loss of water shortage. Keeping these assumptions in mind, I used the straight-line method used in Page 375 of Dingman to compose an average constant base flow for the watershed and subtracted this from the outflow to determine the total runoff for the period. The second phase of Part 1 required interpreting the specific time instants on the hydrograph to determine a series of characteristic intervals, specifically the response lag time (TLR), the lag-to-peak time (TLP), the centroid lag time (TLC), and the concentration time (TC). The first two intervals were derived from observed instants in the recorded hydrographs of the precipitation and outflow of the watershed; for the centroid lag time, I relied on Equation 9-3 in Dingman to compute the centroids of effective water input and hydrograph response. Given the broad time intervals, I estimated the time of peak discharge (tpk) occurring roughly halfway between the 12th and 13th day of the recorded storm event. To determine the concentration time, Dingman provided four different methods over the course of the text; after comparing the results of each, I decided simply basing it on the time between the estimated end of response and end of effective water input intervals produced the most reasonable result. Cum. Rainfall Amount (Weff) (m3): Cu . Baseflo Cu . I flo Cu . Outflo ΣQ ΣQi 2179659.407 3 ): 366000 3 ): ΣQout Σ * As 2545659.407 (Weff + Base Flow) 3 ): Cu . Ru off A ou t ΣR 6498000 3 ): Cu . E apotra spiratio ΣET 6132000 (Outflow - Baseflow) 3 ): 1020035.737 Net Cha ge i “torage Δ“ ): Outflow/Inflow Ratio (Qout/Qin): Duration of Water Input (Tw) (d): Response Lag Time (TLR) (d): 3 -4972376.33 (P - R - ET, Brutsaert 4.55) 2.552580279 6 (tw0 - twe) 2 (tq0 - tw0) Rise Time (Tr) (d): 5 Rise Time (Tr) (d): (tpk - tq0) 3.121660607 (Dingman 9-48) Lag to Peak Time (TLP) (d): 7 (tpk - tq0) Centroid Lag Time (TLC) (d): 9.030922297 (tqc - twc) Centroid Lag to Peak Time (TLPC) (d): 0.830552952 (tpk - twc) Concentration Time (Tc) (d): 0.202767678 (Dingman Table 9-9, Watt and Chow) Concentration Time (Tc) (d): 1.38702343 Concentration Time (Tc) (d): 41.58739718 (Dingman 9-6) Concentration Time (Tc) (d): 15 (= 1.67 * TLPC, Dingman 441) (tqe - twe) Table 1: Results of Past Storm Analysis Equations: Maximum Retention Capacity (Dingman 9-47): Effective Rainfall (Dingman 9-46): Pressure-Temperature Slope (Dingman 7-6): Latent Heat of Vaporization (Dingman 7-8): Psychometric Constant (Dingman 7-13): Water Density (Dingman B-3): Potential Evapotranspiration (Dingman 7-65): Actual Evapotranspiration (Dingman 7-70): Centroid of Effective Water Input (Dingman 9-3): Change in Water Storage (Brutsaert 4.55): Part 2: Estimation of Future Storm Parameters The second phase of this project involved deriving estimates of the cumulative runoff amount (ΣR), the net change in storage, and the time of peak outflow for a future storm event, for which my instructor provided the daily precipitation, net radiation, and temperature. Deriving these estimates required use of the same assumptions as Part 1, and for consistency the same spreadsheet format. In this iteration, I estimated the cumulative outflow (ΣQout by multiplying the sum of cumulative effective rainfall and estimated base flow by the ratio of outflow to inflow (Qin) from Part 1, and derived the total runoff from the difference between the total outflow and the total base flow (again assuming ase flo re ai ed o sta t for oth the atershed’s i put a d output . The high Outflow:Inflow ratio used in both assessments implied the watershed was undergoing a massive net decrease in water storage, suggesting I had severely underestimated this base flow, as my first attempt at modeling the watershed response would prove. Finally, the beginning and end of the effective water input for this future storm event were observed on days 1 and 7, respectively; I estimated the beginning of the hydrograph rise and the time of peak discharge by adding the derived response lag time and time of rise (TR), respectively, from Part 1. Flow 3 Cu . Rai fall A ou t ΣWeff ): 3 Cu . Baseflo ΣQ ): 3 Cu . I flo ΣQi ): Outflow/Inflow Ratio (Qout/Qin): 31677338.78 138000 31815338.78 2.552580279 Cu . Outflo 81211206.34 (=Qin * (Qout/Qin)) 81073206.34 (Outflow - Baseflow) ΣQ 3 ): Cum. Ru off A ou t ΣR Cu . E apotra spiratio Net Cha ge i “torage Δ“ 3 ): ΣET 3 3 ): ): 0 -49395867.56 (Weff - ΣR - ΣET Instants Day Beginning of Hydro. Rise (tq0): 3 (Estimated, tw0 + TLR) Time of Peak Discharge (tpk): 6.121660607 (Estimated, tq0 + TR) Table 2: Estimates of Future Storm Part 3: Linear Watershed Modeling / Calibration As per the instructions, I began the composition of a linear reservoir model in the Java-based differential equation program Berkley-Madonna with the most rudimentary design laid out i oth Di g a ’s te t a d Pederse et al. i the U.S. ACOE handbook. This design accounted only for the initial outflow and the response time of the given watershed, the latter of which I had initially based on the centroid lag time calculated in Part 1. Needless to say, this initial model greatly exaggerated the expected outflow of the watershed, showing a peak discharge more than twice that of the recorded peak while having a residual outflow towards the end of the recorded period that dropped well below that of the recorded outflow. However, when basing the response time instead on the centroid lag time, as suggested on Page 401 in Dingman, the calculated peak outflow dropped by an order of magnitude, signifying my first customization of parameters for the project. Figure 1: Primary Basic Linear Model M e t step uilt upo this odel first addi g a Weff fu tio ased o the ratio of total precipitation to the effective precipitation previously calculated from the SCS method in Part 1 and m3 / d applied this to the outflow of the first reservoir. Essentially this first reservoir would come to represent the total catchment area of the basin, while the second reservoir that the Weff flow fed into would represent the total precipitation that could contribute to surface runoff (i.e. not absorbed by the crops grown in the aforementioned agricultural regions or evaporated directly from the tree cover). To further close the gap between the recorded and expected outflow and o tri ute to the odel’s le el of detail, I also added a second flow component based on the assumed base flow of the watershed, feeding directly to this second reservoir based on the given response time. Finally, for the sake of realism I added an evapotranspiration outflow from this new reservoir based on a rounded graph of the given ET data, based on the assumption that the most significant amount of evapotranspiration from tree and plant cover would emanate from this early stage, before feeding into surface and groundwater flow. This e asi odel ould ser e as the pri ar te plate for later ad a ed Figure 2: Secondary Basic Linear Model odels. 1800000 1600000 1400000 1200000 1000000 800000 600000 400000 200000 0 Recorded Outflow Calc. Outflow ET Effective Precip. 0 10 20 30 40 50 60 70 Time (d) Total Outflow, Calc. (m3/d): 29761257.68 Total Infiltration (m3/d): Outflow RSME: RSME % of Total Outflow: 85270.32165 (Inflow - Outflow - ET) 527987.9232 1.774077994 Figure 3: Secondary Linear Model Hydrograph The process of building more advanced linear reservoir models upon this basic template primarily consisted of extensive trial and error when adding different components to get a response curve that roughly mimicked that of the recorded outflow, followed by tweaking of various parameters to close the gap between the two. The relatively short response time and small ratio of total precipitation to outflow of the watershed suggested the presence of ulterior modes of outflow and secondary reservoirs with differing outflow velocities. Going under the assumption that the first of such modes would occur from the interception of precipitation from tree cover and subsequent evaporation, I first added such a fractional flow (of the same global response time of the watershed) from the first reservoir based on rough assumptions of tree cover regression modeling proposed in Chapter 7 of Dingman, ensuring only residual flow from the treetop would contribute to the effective precipitation ratio. Ultimately I reasoned that this first reservoir would primarily represent precipitation suspended by tree and crop cover, the fractional outflow of which would represent water directly evaporated from the plant surface (not considered in the e apotra spiratio flo . O e riti al ala i g ethod for this model involved giving this first reservoir a relatively high initial water content, essentially ensuring the watershed would still have an effective precipitation for the first few days of this study as the watershed responded to a previous storm event before the recorded time interval. I then added another reservoir after the first outflow with a higher initial water content and much slower outflow rate, representing two areas of starkly contrasting topography. The relati el large ra ge of Ma i g’s roughness coefficients in Table 9-6 of Dingman attributed to the various land types described in the basin (trees, row crops, grassland) supports the rationality of this assumption. Figure 4: First Advanced Model For the second advanced model, I significantly cut the evaporative outflow from tree cover, increased the estimated base flow, and added groundwater (GW) outflow components, as the watershed description of heterogeneous loamy soils with some tree and agriculture cover suggested the former aspect would have a stronger impact. Reasoning that GW infiltration would occur at different rates over different areas of the basin, and likely be highest following areas of greatest concentration (reservoirs), I added two GW flows of differing percolation rates on each of the secondary reservoirs, relying on the assumption that stream gain from groundwater compared to surface water runoff was negligible. The relatively high assumed hydraulic conductivity and field capacity of the soil and the relatively short response time of the watershed, arguably brought on by high initial water content from the assumed previous storm event, supports this hypothesis. Figure 5: Second Advanced Model After inputting the precipitation and evapotranspiration graphs of the recorded storm event into both models, I copied the readings of the output variables into a spreadsheet, plotted the calibrated outflow curve on top of the recorded outflow curves, and computed the respective cumulative outflow and infiltration (the total amount of inflow from base flow and effective precipitation, minus the amount of outflow and evapotranspiration). For the sake of accuracy in the second model, I also calculated the total infiltration based on the total calculated groundwater outflow of both effective reservoirs in the watershed and was pleased with the narrow range of the results. Most notably, by accounting for groundwater flow, the second model showed a net increase in storage, as opposed to the net decrease predicted in the first model. Finally, I used the Root Mean Squared Error (RMSE) method to quantify the accuracy of each calculated outflow curve to the recorded outflow curve; as anticipated, the second model showed a lower Root Mean Square Error (RSME) value, supporting the accuracy of its design. 1000000 m3/d 800000 600000 Recorded Outflow 400000 Calculated Outflow ET 200000 Effective Precipitation 0 0 10 20 30 40 50 60 70 Time (d) Total Outflow, Calc. (m3/d): 3 Total Infiltration (m /d): Outflow RSME: RSME % of Total Outflow: 7620421.104 -367131.4736 (Inflow - Outflow - ET) 48549.28741 0.637094548 Figure 6: Hydrograph of Past Storm, First Advanced Model 1400000 1200000 m3/d 1000000 800000 Recorded Outflow 600000 Calc. Outflow 400000 ET 200000 Effective Precipitation 0 0 10 20 30 40 50 60 70 Time (d) Total Outflow (m3/d): 7621039.416 Total Infiltration (m3/d): 9894895.462 (Total GW Flow) 3 Total Infiltration (m /d): 9820553.113 (Inflow - Outflow - ET) Outflow RSME: 39768.59027 RSME % of Total Outflow: 0.52182633 Figure 7: Hydrograph of Past Storm, Second Advanced Model If pressed to rely on one of the two models to predict the outflow of the watershed in response to a storm event, I would have to choose the latter model, incorporating canopy interception, multiple reservoirs with different surface streamflow rates, and multiple regions of varying groundwater infiltration. Not only were its results statistically closer to the observed outflow, but its prototype design seemed more realistic in a watershed where the sandy and loamy soils and scattered tree cover would ensure groundwater infiltration would play a much larger role than canopy interception in changes to the overall water storage. Additionally, due to the relatively large elevation gradient over the span of the watershed, one could interpret both the periods of highly different surface flow velocity and high rates of groundwater flow through porous media assumed in the model to be perfectly reasonable. Moreover, the response curve of the second model showed considerably less sensitivity towards changes in the initial water content of the three reservoirs in the watershed, further improving its accuracy in cases when initial watershed conditions were unknown. Equations: Linear Reservoir Model Outflow (Dingman 9B2-4): Root Mean Square Error: Part 4: Prediction of Future Runoff m3/d The final task of this project was to apply the two advanced linear models in the prediction of a future rainfall event (with a given rate of precipitation). A primary issue with calibrating such models for future events was a lack of data regarding the initial water content of the watershed, whether in the phreatic zone or the given reservoirs; for the sake of this experiment, I chose to assume roughly the same initial conditions applied. Using the new precipitation graph, I also calibrated both models with the new Weff:W ratio and predicted actual evapotranspiration, which, due to a lower overall rate of solar radiation input, turned out to be negligible for this time period. 40000000 35000000 30000000 25000000 20000000 15000000 10000000 5000000 0 Precipitation Effective Precipitation Outflow 0 5 10 15 20 25 Time (d) Cum. Runoff Amount (m3): 37879486.25 (Outflow - Baseflow) Net Cha ge i “torage Δ“ ): 40107916.65 (P - R - ET, Brutsaert 4.55) 7 Time of Peak Discharge (d) (tpk): 3 Figure 8: Prediction of Future Storm, 1st Advanced Model 40000000 35000000 m3/d 30000000 25000000 Precipitation 20000000 15000000 Effective Precipitation 10000000 Calc. Outflow 5000000 0 0 5 10 15 20 25 Time (d) Cum. Runoff Amount (m3): 22134947.23 (Outflow - Baseflow) Net Cha ge i “torage Δ“ ): 33850353.34 (P - R - ET, Brutsaert 4.55) 6 Time of Peak Discharge (d) (tpk): (Observed) 3 Figure 9: Prediction of Future Storm, 2nd Advanced Model Again, I would be compelled to favor the second model over the first for this prediction; it accommodated for a more thorough range of variables defining the watershed, and showed considerably less reliance on unknown factors, such as initial water content and tree cover. It also depicted a more realistic estimation of predicted watershed elements, such as the portion of precipitation intercepted by tree cover and the different rates of overland flow between reservoirs. Additionally, each of the given components of the second model could be easily customized to accommodate for any future statistical data gathered on the watershed (i.e. soil porosity, groundwater infiltration), whereas the first model would require more extensive alterations to its more limited range of parameters, further stretching the bounds of realism in its representation of the watershed. As such, this project illustrated to me the fine line every earth scientist and engineer must walk between accuracy and realism when making any kind of prediction relating to watershed hydrology. Sources: Brutsaert, Wilfried. Hydrology: An Introduction. Cambridge, UK: Cambridge University Press, 2010. Print. Dingman, S. Lawrence. Physical Hydrology. 2nd. Long Grove, IL: Waveland Press, Inc., 2008. Print. Pedersen, J.T.; Peters, J.C.; Helweg, O.J. Hydrographs by Single Linear Reservoir Model. Journal of Hydraulics Division, ASCE. 106.HY5 (1980): 837-852. Print.