<%BANNER%>

Geologic, Vegetative and Climatic Controls on Coupled Hydrologic Processes in a Complex River Basin

MISSING IMAGE

Material Information

Title:
Geologic, Vegetative and Climatic Controls on Coupled Hydrologic Processes in a Complex River Basin Lessons Learned From a Fully Integrated Hydrologic Model
Physical Description:
1 online resource (206 p.)
Language:
english
Creator:
Srivastava, Vibhava
Publisher:
University of Florida
Place of Publication:
Gainesville, Fla.
Publication Date:

Thesis/Dissertation Information

Degree:
Doctorate ( Ph.D.)
Degree Grantor:
University of Florida
Degree Disciplines:
Agricultural and Biological Engineering
Committee Chair:
Graham, Wendy Dimbero
Committee Members:
Martin, Jonathan Bowman
Munoz-Carpena, Rafael
Cohen, Matthew J
Annable, Michael D
Maxwell, Reed M

Subjects

Subjects / Keywords:
integrated -- model
Agricultural and Biological Engineering -- Dissertations, Academic -- UF
Genre:
Agricultural and Biological Engineering thesis, Ph.D.
bibliography   ( marcgt )
theses   ( marcgt )
government publication (state, provincial, terriorial, dependent)   ( marcgt )
born-digital   ( sobekcm )
Electronic Thesis or Dissertation

Notes

Abstract:
This dissertation documents the first implementation of an integrated transient 3D surface water-groundwater-land surface process model, ParFlow. CLM, to evaluate the interacting geologic, climatic and vegetative controls on water budget components and streamflow generation processes over the Santa Fe River Basin in North Central Florida.  Model predictions indicate that evapotranspiration (ET) is the most important water balance component in the basin comprising 77% of rainfall.  Geologic conditions and vegetative properties were found to exert primary control on the spatial variability of streamflow generation processes in the basin through their influence on the balance between rainfall, ET, runoff and infiltration processes.  Climatic variability was found to provide primary control on the temporal variability of streamflow generation processes. Model predictions indicate that in the upper basin more than 95% of streamflow is generated by recent near-stream rainfall. In contrast, in the lower basin the majority of streamflow is contributed by the Upper Floridan Aquifer, with the fraction of subsurface flow averaging approximately 77% at the outlet of basin. A global sensitivity analysis of the model revealed that the permeability of the Intermediate Aquifer System is the most influential factor driving hydrologic responsethroughout the SFRB. Particle tracking experiments predicted that the median age of streamflow in the upper basin ranges from approximately 1 day at the peak of storm hydrographs to approximately 7 days at the end of stormflow recession,with travel time distributions that vary over time but are generally well-fit with log-normal distributions at the peak of the storm hydrograph. The median age of subsurface contributions to streamflow in the lower portion of the basin was predicted to be approximately 17 years, and the travel time distribution for the subsurface contribution is well-fit by a gamma distribution showing fractal properties that do not vary significantly over time.  The fraction of new stormwater versus old groundwater in the streamflow in the unconfined region, and thus the shape of the total streamflow travel time distribution, varies as a function antecedent conditions, storm magnitude, time during the storm, and assumptions regarding the contrast in hydraulic conductivity between high permeability zones and the porous matrix.
General Note:
In the series University of Florida Digital Collections.
General Note:
Includes vita.
Bibliography:
Includes bibliographical references.
Source of Description:
Description based on online resource; title from PDF title page.
Source of Description:
This bibliographic record is available under the Creative Commons CC0 public domain dedication. The University of Florida Libraries, as creator of this bibliographic record, has waived all rights to it worldwide under copyright law, including all related and neighboring rights, to the extent allowed by law.
Statement of Responsibility:
by Vibhava Srivastava.
Thesis:
Thesis (Ph.D.)--University of Florida, 2013.
Local:
Adviser: Graham, Wendy Dimbero.

Record Information

Source Institution:
UFRGP
Rights Management:
Applicable rights reserved.
Classification:
lcc - LD1780 2013
System ID:
UFE0045666:00001

MISSING IMAGE

Material Information

Title:
Geologic, Vegetative and Climatic Controls on Coupled Hydrologic Processes in a Complex River Basin Lessons Learned From a Fully Integrated Hydrologic Model
Physical Description:
1 online resource (206 p.)
Language:
english
Creator:
Srivastava, Vibhava
Publisher:
University of Florida
Place of Publication:
Gainesville, Fla.
Publication Date:

Thesis/Dissertation Information

Degree:
Doctorate ( Ph.D.)
Degree Grantor:
University of Florida
Degree Disciplines:
Agricultural and Biological Engineering
Committee Chair:
Graham, Wendy Dimbero
Committee Members:
Martin, Jonathan Bowman
Munoz-Carpena, Rafael
Cohen, Matthew J
Annable, Michael D
Maxwell, Reed M

Subjects

Subjects / Keywords:
integrated -- model
Agricultural and Biological Engineering -- Dissertations, Academic -- UF
Genre:
Agricultural and Biological Engineering thesis, Ph.D.
bibliography   ( marcgt )
theses   ( marcgt )
government publication (state, provincial, terriorial, dependent)   ( marcgt )
born-digital   ( sobekcm )
Electronic Thesis or Dissertation

Notes

Abstract:
This dissertation documents the first implementation of an integrated transient 3D surface water-groundwater-land surface process model, ParFlow. CLM, to evaluate the interacting geologic, climatic and vegetative controls on water budget components and streamflow generation processes over the Santa Fe River Basin in North Central Florida.  Model predictions indicate that evapotranspiration (ET) is the most important water balance component in the basin comprising 77% of rainfall.  Geologic conditions and vegetative properties were found to exert primary control on the spatial variability of streamflow generation processes in the basin through their influence on the balance between rainfall, ET, runoff and infiltration processes.  Climatic variability was found to provide primary control on the temporal variability of streamflow generation processes. Model predictions indicate that in the upper basin more than 95% of streamflow is generated by recent near-stream rainfall. In contrast, in the lower basin the majority of streamflow is contributed by the Upper Floridan Aquifer, with the fraction of subsurface flow averaging approximately 77% at the outlet of basin. A global sensitivity analysis of the model revealed that the permeability of the Intermediate Aquifer System is the most influential factor driving hydrologic responsethroughout the SFRB. Particle tracking experiments predicted that the median age of streamflow in the upper basin ranges from approximately 1 day at the peak of storm hydrographs to approximately 7 days at the end of stormflow recession,with travel time distributions that vary over time but are generally well-fit with log-normal distributions at the peak of the storm hydrograph. The median age of subsurface contributions to streamflow in the lower portion of the basin was predicted to be approximately 17 years, and the travel time distribution for the subsurface contribution is well-fit by a gamma distribution showing fractal properties that do not vary significantly over time.  The fraction of new stormwater versus old groundwater in the streamflow in the unconfined region, and thus the shape of the total streamflow travel time distribution, varies as a function antecedent conditions, storm magnitude, time during the storm, and assumptions regarding the contrast in hydraulic conductivity between high permeability zones and the porous matrix.
General Note:
In the series University of Florida Digital Collections.
General Note:
Includes vita.
Bibliography:
Includes bibliographical references.
Source of Description:
Description based on online resource; title from PDF title page.
Source of Description:
This bibliographic record is available under the Creative Commons CC0 public domain dedication. The University of Florida Libraries, as creator of this bibliographic record, has waived all rights to it worldwide under copyright law, including all related and neighboring rights, to the extent allowed by law.
Statement of Responsibility:
by Vibhava Srivastava.
Thesis:
Thesis (Ph.D.)--University of Florida, 2013.
Local:
Adviser: Graham, Wendy Dimbero.

Record Information

Source Institution:
UFRGP
Rights Management:
Applicable rights reserved.
Classification:
lcc - LD1780 2013
System ID:
UFE0045666:00001


This item has the following downloads:


Full Text

PAGE 1

1 GEOLOGIC, VEGETATIVE AND CLIMATIC CONTROLS ON COUPLED HYDROLOGIC PROCESSES IN A COMPLEX RIVER BASIN: LESSONS LEARNED FROM A FULLY INTEGRATED HYDROLOGIC MODEL By VIBHAVA S R IVASTAVA A DISSERTATION PRESENTED TO THE GRADUATE S CHOOL OF THE UNIVERSITY OF FLORIDA IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY UNIVERSITY OF FLORIDA 2013

PAGE 2

2 2013 Vibhava Srivastava

PAGE 3

3 To my family

PAGE 4

4 ACKNOWLEDGEMENTS Fi rst and foremost I would like to thank GOD for everything that led to the fulfillment of all that a student can desire passion for research, a knowledgeable nice and perfect advisor, supportive family and friends. It has been a long journey during which I learnt valuable life lessons on perseverance, passion, and people. I am lucky to her enough for her patience, advice, and time amongst many other things. She taught me th ree most valuable lessons that I will live with for the rest of my life. First, work on things that you feel passionate about. Second, never give up. Third, you need love and support of all the people in your personal and professional life to succeed. Dr. Graham Thanks for everything! I am grateful to the members of my graduate committee Dr. Rafa el Munoz Carpena, Dr. Jonathan Martin, Dr. Matt hew Cohen, Dr. Michael Annable, and Dr. Reed Maxwell for their valuable time and interest in completing this rese arch. Their timely advice and input helped to improve the quality of this research. Special thanks to Dr. Carpena for inspiring new ideas and never letting me settle for less. a nd Sarita, and my sister Smita. If it was not for their encouragement and support I would have never made it to the United States. Their love and blessings kept me going when things got tough during my stay in the USA. I am lucky to have a friend in my wi fe, Richa, who has been there for me throughout my research. She inspires me with her sincerity and passion towards her work and with her perseverance. I am grateful to her for being there, listening to me when I got frustrated, cheering me up in bad times praying for me when I struggled,

PAGE 5

5 and for unconditionally loving me. My love and hugs to my son Parth whose smiles gave me endurance to take on all the challenges of the past three years. I am grateful to my Water Institute Family, Kathleen, Lisette, S yeWoon, Rob, Mary, Jason, and Wes for all the moral support and love that they have extended to me. Big thanks to Kathleen for working tirelessly to collect and manage all the data needed for my research. During my student life, I met several awesome peopl e. Many have since become good friends and are an integral part of my life. I extend my gratitude to Sumit Sen and Madhukar Varshney who have been like elder brothers to me. Thanks for your patience and advice during our long phone conversations. Thanks to Shailey for your words of encouragement and for being there for me and my family when we needed you. Thanks to Purushottum for standing by me during the most difficult days of my life. Thanks to Rob deRooij for accompanying me on countless coffee breaks a nd long discussions on personal and professional issues. At last I would like to thanks everyone who I have not mentioned above and who have been helpful to me during my time at the University of Florida.

PAGE 6

6 TABLE OF CONTENTS Page ACKNOWLEDGEMENTS ................................ ................................ ............................... 4 LIST OF TABLES ................................ ................................ ................................ ............ 8 LIST OF FIGURES ................................ ................................ ................................ .......... 9 ABSTRACT ................................ ................................ ................................ ................... 13 CHAPTER 1 INTRODUCTION ................................ ................................ ................................ .... 15 2 GEOLOGIC AND CLIMATI C CONTROLS ON STREAM FLOW GENERATION PROCESSES IN A COMPL EX EOGENETIC KARST B ASIN ................................ 21 2.1 Background ................................ ................................ ................................ ....... 21 2 .2 Santa Fe River Basin Physical System and Conceptual Model ..................... 24 2.3 Methods ................................ ................................ ................................ ............ 26 2.3.1 River Source Water Mixing Dynamics: End M ember Mixing Analyses .. 26 2.3.2 Integrated 3D Model ParFlow.CLM ................................ ..................... 29 2.3.2.1 Model setup, parameterization and initialization ....................... 32 2.3.2.2 Estimating surface water and groundwater contributions to streamflow using ParFlow.CLM ................................ ................ 37 2.3.2.3 Benchmarking ParFlow.C LM ................................ .................... 37 2.4 Results and Discussion ................................ ................................ ..................... 38 2.4.1 Streamflow ................................ ................................ ............................. 38 2.4.2 Gr oundwater Elevation ................................ ................................ .......... 39 2.4.3 Surface Versus Subsurface Contributions to Streamflow ...................... 42 2.4.4 Analysis of Water Budget Components for the SFRB ............................ 44 2.4.5 Geological Control Over Spatial Pattern of Monthly ET Losses and Subsurface Recharge ................................ ................................ ............ 46 2.4.6 Vegetation a nd Geological Control Over Water Table Depth and ET Dynamics ................................ ................................ ............................... 47 2.5 Summary and Conclusions ................................ ................................ ............... 48 3 GLOBAL SENSITIVITY A NALYSIS OF AN INTEGRATED HYDROL OGIC MODEL INSIGHT ON GEOLOGIC AND VEGETATIVE CONTR OLS OVER THE HYDROLOGIC FUNCT IONING OF A LARGE CO MPLEX BASIN ................. 68 3.1 Background ................................ ................................ ................................ ....... 68 3.2 Methods ................................ ................................ ................................ ............ 73 3.2.1 Hydro Geologic Characteristics of the Santa Fe River Basin ................ 73 3.2.2 Integrated 3 D Model ParFlow.CLM ................................ ..................... 75 3.2.3 Baseline Model ................................ ................................ ...................... 78

PAGE 7

7 3.2.4 Morris Screening Method ................................ ................................ ....... 80 3.2.5 GSA Experiment Design ................................ ................................ ........ 81 3.2.6 Measures of Goodness of Model Fit ................................ ...................... 83 3.3 Results and Discussion ................................ ................................ ..................... 84 3.3.1 Baseline Model Performance ................................ ................................ 84 3.3.2 Global Sensitivity Analysis ................................ ................................ ..... 86 3.3.2.1 Evapotranspiration ................................ ................................ .... 87 3.3.2.2 Total streamflow ................................ ................................ ....... 88 3.3.2.3 Peak streamflow ................................ ................................ ....... 89 3.3.2.4 Storm hydrograph recession ................................ ..................... 90 3.3.2.5 Mean and coefficient of variation of groundwater elevation ...... 91 3.3.2.6 Surface water groundwater interactions ................................ .. 92 3.3.3 Goodness of Fit of the GSA Model Runs ................................ ............... 93 3.4 Conclusions ................................ ................................ ................................ ...... 97 4 TOPOGRAPHIC, GEOLOGIC AND CLIMATIC CONTROLS OVER RIVER WATER SOURCES, MIXING DYNAMICS AND TRAVEL TIME DISTRIBUTIONS IN A LARGE COMPLEX BASIN ................................ ............... 126 4.1 Background ................................ ................................ ................................ ..... 126 4.2 Santa Fe River Basin ................................ ................................ ...................... 132 4.3 Methods ................................ ................................ ................................ .......... 133 4.3.1 Integrated Hydrologic Modeling ................................ ........................... 133 4.3.2 Surface Subsurface Particle Tracking ................................ ................. 135 4.3.3 River Water Source Mix ing Dynamics and Travel Time Distribution .... 136 4.4 Results and Discussion ................................ ................................ ................... 138 4.4.1 Streamflow Age Distribution Under Time Varyin g Hydrologic Conditions ................................ ................................ ............................ 138 4.4.2 Streamflow Origin and Flowpaths Over Contrasting Hydrogeologic Regions Under Different Hydrologic Conditions ................................ .. 141 4.4.3 Geologic Controls Over Water Travel Time and Flow Paths in an Unconfined Aquifer ................................ ................................ .............. 144 4.4.4 On Time Invariant Assumptions in TTD Modeling ............................... 146 4.5 Summary and Conclusions ................................ ................................ ............. 148 5 CONCLUSIONS, CON T RIBUTIONS AND RECOMMENDATIONS FOR FUTURE REASEARCH ................................ ................................ ........................ 183 APPENDIX: CHAPTER 4 ADDITIONAL FIGURES ................................ ..................... 190 LIST OF REFERENCES ................................ ................................ ............................. 195 BIOGRAPHICAL SKETCH ................................ ................................ .......................... 206

PAGE 8

8 LIST OF TABLES Table Page 2 1 Summary of ParFlow input variables and their relevant source of information. .. 52 2 2 Subsurface flow contributions as percent of total streamflow at river locations downstream of Cody Escarpment (CE). ................................ ............................. 5 2 2 3 Annual total of major water balance components for the SFRB. ........................ 53 3 1 Physical characteristics of the Santa Fe River Basin. ................................ ....... 102 3 2 Input data used in ParFlow.CLM for the SFRB. ................................ ................ 103 3 3 ParFlow parameters included in GSA by Morris method. ................................ 104 3 4 CLM parameters included in GSA by Morris method. ................................ ....... 105 3 5 Summary of relevant information on intial and boundary condition used in baseline ParFlow.CLM. ................................ ................................ .................... 106 3 6 ParFlow.CLM output variables assessed in the GSA. ................................ ...... 107 3 7 Morris mean and standard deviation based rank of individual parameters for mean daily ET Losses. ................................ ................................ ..................... 108 3 8 Morris mean and sigma based rank of individual parameters for total and peak streamflow. ................................ ................................ .............................. 109 3 9 Morris mean and sigma based rank of individual parameters for streamflow recession characteristics. ................................ ................................ ................. 110 3 10 Morris mean and standard deviation based rank of individual parameters for mean daily groundwater elevation. ................................ ................................ ... 112 3 11 Mo rris mean and standard deviation based rank of individual parameters for coefficient of variation daily groundwater elevation. ................................ ......... 113 3 12 Morris mean and standard deviation based rank of individ ual parameters for total baseflow contribution received by three unconfined river locations. ......... 114 4 1 1500. ................................ ................................ ................................ ................ 151 4 2 Parameter values for best fit function for overall CDFs at station 1500 and ................................ ..................... 151

PAGE 9

9 LIST OF FIGURES Figure Page 2 1 Map of the model domain.. ................................ ................................ ................. 54 2 2 3D mesh used to define the domain. ................................ ................................ .. 55 2 3 Spatial distribution of land surface properties over the SFRB. ........................... 56 2 4 Comparison of observed and model predicted daily streamflow measurements at various streamgage locations a cross the domain. ................. 57 2 5 Comparison of observed and model simulated groundwater surface elevation at various locations across the domain. ................................ .............................. 58 2 6 Potentiometric surface of UFAS with respect to mean sea level. ....................... 59 2 7 Time series of daily total subsurface flow contributions received by all the upstream river reaches from a given gage location. ................................ ........... 60 2 8 Comparison of EMMA based versus ParFlow.CLM simulated ratio of subsurface flow contributions to total streamflow at three UR river locations. .... 61 2 9 Comparison of measured and ParFlow.CLM predicted streamflow for calendar years 2002 2003 (left column) and 2007 2008 (right column). ............ 63 2 10 Monthly averaged water budget components for the study domain over entire simulation time period ................................ ................................ ......................... 64 2 11 Spatial distribution of monthly averaged Evapotranspiration (ET; mm) across the study d omain over entire simulation period. ................................ ................. 65 2 12 Spatial distribution of average monthly subsurface recharge (mm) across the study domain over the entire simulation period. ................................ ................. 66 2 13 Controls on average daily ET (mm) and water table depth (m) dynamics in the SFRB.. ................................ ................................ ................................ .......... 67 3 1 Map of the model domain. ................................ ................................ ................ 115 3 2 3D mesh used to define the domain. ................................ ................................ 116 3 3 Comparison of observed and model predicted daily streamflow measurements at various streamgage locations across the domain. ............... 117 3 4 Comparison of observed and model simulated groundwater elevation at various locations across the domain. ................................ ................................ 118

PAGE 10

10 3 5 effects of each parameter. ................................ ................................ ................ 119 3 6 entary effects of each parameter. ................................ ................................ ................ 120 3 7 effects of each parameter. ................................ ................................ ................ 120 3 8 Flow chart of separation of behavioral models from non behavioral models based on streamflow simulation by 340 GSA model runs. ............................... 121 3 9 Comparison of observed and model simulated total streamflow at USGS 1500 (top) and USGS 2500 (bottom). ................................ ............................... 122 3 10 Co mparison of observed and model simulated total streamflow at wells E2, E3, E8, and E10. ................................ ................................ .............................. 123 3 11 Range and median values for most sensitive parameter identified during the Morris Analysis. ................................ ................................ ................................ 124 3 12 Influence of most sensitive parameters i dentified during Morris global sensitivity analysis on median values of output variables of interest. ............... 125 4 1 Map of the model domain.. ................................ ................................ ............... 152 4 2 3D mesh used to define the domain. ................................ ................................ 153 4 3 Time series of ParFlow.CLM simulated daily flow at station 1500 and station 2500. ................................ ................................ ................................ ................ 154 4 4 Hydrograph points where particle tracking experiments were conducted during year 2003 and year 2004 storm events at station 1500. ........................ 156 4 5 Hydrograph points where particle tracking experiments were conducted during year 2003 and year 2004 storm events. ................................ ................ 158 4 6 Water particle age distribution at station 1500 for various points on storm ................................ ................................ ........................... 160 4 7 Water particle age distribution at station 2500 for various points on storm ................................ .............................. 162 4 8 W ater particle median age versus flow for base case model at various times during storm hydrographs ................................ ................................ ................. 164 4 9 Spatial extent from which water particles arrived at station 1500 during storm event ................................ ................................ ................................ ........... 165

PAGE 11

11 4 10 Spatial extent from which water particles arrived at station 1500 during storm ................................ ................................ ................................ ........... 166 4 11 Spatial extent ................................ ................................ ................................ 167 4 12 2500 during storm event ................................ ................................ ............. 168 4 13 ................................ ................................ ...................... 169 4 14 Spat ................................ ................................ ................................ 170 4 15 2500 during s ................................ ................................ ............. 171 4 16 ................................ ................................ ...................... 172 4 17 Water particle age distribution at station 2500 at various points on storm e ................................ ................................ ........ 173 4 18 Water particle age distribution at station 2500 at various points du ring ................................ ................................ .. 174 4 19 Temporal dynamics of various sources of water at 2500 during various time within a storm hydrograph for base case and its variants. ................................ 175 4 20 Water particle median age versus flow for base case model variant at various ................................ ................................ .... 175 4 21 Spatial extent of water particle source location during point e4 of storm event ................................ ................................ ................................ ..................... 176 4 22 Comparison of measured and ParFlow.CLM simulated total streamflow and subsurface contributions. ................................ ................................ .................. 177 4 23 Water particle age distribution at station 2500 ................................ .................. 178 4 24 Water particle age distribution at station 1500. ................................ ................. 179 4 25 Water particle age distribution at station 2500 during various time in storm ................................ ................................ .......... 180 A 1 Water particle median age versus flow for base case model variants. ............. 190

PAGE 12

12 A 2 Temporal dynamics of various sources of water at 2500 during various times within a storm hydrograph for base case and its variants. ................................ 191 A 3 Water particle age distribution at station 2500 at various points on storm ................................ ................................ ................................ ........... 192 A 4 Spatial extent of water particle s 193 A 5 ................. 194

PAGE 13

13 Abstract of Dissertation Pr esented to the Graduate School of the University of Florida in Partial Fulfillment of the Requirements for the Degree of Doctor of Philosophy GEOLOGIC, VEGETATIVE AND CLIMATIC CONTROLS ON COUPLED HYDROLOGIC PROCESSES IN A COMPLEX RIVER BAS IN: LESSONS LEARNED FROM A FULLY INTEGRATED HYDROLOGIC MODEL By Vibhava Srivastava August 2013 Chair: Wendy D. Graham Major: Agricultural and Biological Engineering This dissertation documents the first implementation of a n i ntegrated transient 3D surf ace water groundwater land surface process model, ParFlow.CLM, to evaluate the interacting geologic, climatic and vegetative controls on water budget components and streamflow generation processes over the Santa Fe River Basin in North Central Florida. Mo del predictions indicate that evapotranspiration (ET) is the most important water balance component in the basi n comprising 77% of rainfall. Geologic conditions and vegetative properties w ere found to exert primary control on the spatial variability of st reamf low generation processes in the basin through their influence on the balance between rainfall, ET, runoff and infiltration processes. Climatic variability was found to provide primary control on the temporal variability of streamflow generation proce sses Model predictions indicate that i n the upper basin more than 95% of streamflow is generated by recent near stream rainfall In contrast, in the lower basin the majority of streamflow is contributed by the Upper Floridan Aquifer with the fraction of subsurface flow averaging approximately 77% at the outlet of basin. A global sensitivity analysis of

PAGE 14

14 the model revealed that the permeability of the Intermediate Aquifer System is the most influential factor driving hydrologic response throughout the SFRB Particle tracking experiments predicted that the median age of streamflow in the upper basin ranges from approximately 1 day at the peak of storm hydrographs to approximately 7 days at the end of stormflow recession, with travel time distributions that vary over time but are generally well fit with log normal distributions at the peak of the storm hydrograph. The median age of subsurface contributions to streamflow in the lower portion of the basin was predicted to be approximately 17 years, and the tra vel time distribution for the subsurface contribution is well fit by a gamma distribution showing fractal properties that do not vary significantly over time. The fraction of new stormwater versus old groundwater in the streamflow in the unconfined region and thus the shape of the total streamflow travel time distribution, varies as a function antecedent conditions, storm magnitude, time during the storm, and assumptions regarding the contrast in hydraulic conductivity between high permeability zones and the porous matrix.

PAGE 15

15 CHAPTER 1 INTRODUCTION Increased pressure on water resources, due to expanding urban areas and a growing population, has raised the concern about the sustainability of water resources and the proper functioning of ecosystems dependent on th ese resources. This has increased the demand for robust predictive tools that can help water managers and policy makers with current and future water management decisions in large basins which are comprised of a complex mix of natural and urban environment s. All hydrological systems consists of three major flow domains namely, the atmosphere, the surface and the subsurface which are interconnected with continuous exchange of water between them (Furman, 2008). However, due to disparity in temporal and spati al scales associated with the flow processes within each domain the three domains are commonly modeled as separate entities. As such modeling is focused on one flow domain and the remaining two domains are simplified and applied as boundary conditions, res ulting in incomplete representation of actual physical processes by the models (Furman, 2008). Hydrological processes such as runoff, subsurface flow (both saturated and unsaturated), and evapotranspiration losses co vary in space and time (Huntington a nd Niswonger, 2012). Therefore, it is imperative to use hydrological models which facilitate integration of all hydrological processes when modeling objectives are geared towards exploring their behavior and interactions at various space and time scales un der varying physical settings as well as possible future land use, water use and climatic conditions. (e.g. Huntington and Niswonger, 2012; Bolger et al., 2011; Goderniaux et al., 2009).

PAGE 16

16 Although the first blueprint for an integrated hydrologic model was o utlined in the year 1969 (Freeze and Harlan, 1969), it was not until recently that significant advances were made towards the development of sophisticated, robust hydrological models that exploit advanced parallel computational capabilities and represent t he surface subsurface near land surface hydrologic processes in a seamlessly integrated and physically plausible manner (e.g. ParFlow (Ashby and Falgout, 1996; Kollet and M axwell, 2008a; Kollet and Maxwell, 2006), InHm (VanderKwaak, 1999), MODHMS (Panday and Huyakorn, 2004), HydroGeoSphere (Therrien et al., 2010) ). Past applications of these fully integrated models have been mostly limited to studying controls over rainfall runoff generation for hillslopes or small catchments (area ~ 10 0 10 1 km 2 e.g. Jones et al., 2006). A few applications of these fully integrated models have involved simulation of hydrologic behavior of medium to large scale watersheds (area ~ 10 2 km 2 e.g Kollet and Maxwell, 2008; Bolger et al., 2011; Ferguson and Maxwell, 2010). Modeling studies explicitly simulating evapotranspiration or land surface processes in the modeling framework have rarely been reported (e.g. Kollet and Maxwell, 2008; Li et al., 2008; Goderniaux et al., 2009). The application of fully integrated models to understand surface groundwater interactions and controls over river water source mixing dynamics in large river basins is still in its infancy. Although there have been stud ies on integrated model applications to understand surface groundwater dynamics in idealized test cases (e.g. Frei et al., 2009), applications over real large basins are rare (e.g. Huntington and Niswonger, 2012; Werner et al., 2006). Huntington and Niswon ger (2012) used the GSFLOW model to understand the linkage between snowmelt timing and surface water groundwater

PAGE 17

17 interactions in a 54 km 2 watershed at eastern Sierra Nevada. Werner et al. (2006 ) used MODHMS to i dentify river reaches receiving groundwater influx and spatio temporal characteristics of groundwater surface flow dynamics in the ~420 km 2 Pioneer Valley, Australia. There remains a need for studies which can test the capabilities of fully integrated mode ls at even larger (area ~10 3 km 2 ) and more complex river basins covering multiple geologic formulations with spatio temporally varying connectivity between surface water and groundwater flow systems. The Santa Fe River Basin (SFRB) in North Central Flo rida provides an ideal test bed to study controls over coupled surface subsurface processes and more specifically surface groundwater interactions in large river basins. This basin offers a complex mix of hydro geomorphic conditions which makes it an ideal setting to study and characterize important river flow producing mechanisms and their interactions under varying geological, and physiographic settings. Depending upon the combination of geology, physiography, and land cover, the sources of river water an d their mixing dynamics change dramatically along the entire length of the Santa Fe River. In the upper two thirds of the basin, hydrologic processes are dominated by surface runoff and surficial stores (wetlands and lakes) which is reflected in high varia bility in the flow regime, and low dissolved mineral contents (Ritorto et al., 2009; Bailly Comte et al., 2011; Bailly Comte et al., 2010; Mo ore et al., 2009) In the lower one third of the basin direct mixing between surface and ground water occurs and there are virtually no stream networks feeding the river. The SFRB is a well instrumented basin with various state and federal government agen cies collecting hydrologic and water quality data at daily to monthly temporal scales and it has also been the focus of numerous field scale

PAGE 18

18 studies. Thus there is a wealth of information available which can be used to develop and validate hydrologic model s. This dissertation documents the implementation of a fully integrated 3D surface water groundwater land surface process modeling platform, ParFlow.CLM, to evaluate the interacting geologic, climatic and vegetative controls on streamflow generation pro cesses in a complex eogenetic karst basin in North Central Florida. The model predictions were benchmarked against observations and used to understand hydrologic functioning of the basin, river water sources origins, travels times and mixing dynamics throu ghout the basin. The Parflow.CLM model developed for the SFRB significantly improves over the existing steady state groundwater modeling platforms available for the region. Findings of the study are presented in the next three chapters, each of which meets requirements for independent publication in the peer reviewed scientific literature including background, study area and data, methodology, results, and summary. The key research questions addressed in this dissertation are: 1. What are the major water budge t components in the basin and how do they vary in space and time? How do surface and subsurface water contributions to streamflow vary along the river? 2. What Parflow.CLM parameters are ET, groundwater level, streamflow, and subsurface contributions to stre amflow most sensitive to? How does parametric sensitivity, and interaction among sensitive parameters, vary across the basin? 3. How do streamflow sources, flow paths and travel times vary along the river and betwe en base flow and storm events? What processes and parameters exert primary control on travel time distributions? Chapter 2 documents the development and application of a three dimensional, fully integrated transient modeling platform, ParFlow.CLM, for the SFRB. The model was used to evaluate the int eracting geologic, climatic and vegetative controls on streamflow generation processes and spatiotemporal water budget components over

PAGE 19

19 the basin. Results of this study underscore the usefulness of combining end member mixing model results with integrated h ydrologic models to identify sources of model error while simulating surface water groundwater interactions. This study provides the first implementation of a fully integrated transient modeling platform in the region, which is important in light of the important current issues of declines in aquifer and streamflows and nutrient enrichment of surface and groundwater throughout the study area. Chapter 3 documents a global sensitivity analysis (GSA) experiment designed to investigate the sensitivity of pred ictions from ParFlow.CLM to uncertainties in a wide range of parameters governing the coupled surface, subsurface and land atmosphere processes. This application highlights the efficiency of using the Morris method of GSA to understand interactions among c oupled non linear surface groundwater land atmosphere processes, to identify variables most relevant for hydrologic process of interest, to understand interactions among these variables, and to identify those parameters which contribute highest uncertainty to model predictions The result of the Morris GSA was successful in identifying spatial and temporal variability among the dominant hydrological processes and sensitive parameters for this large complex basin. T he findings demonstrate important nonlinear interactions among geologic, soil and vegetation properties on land atmosphere, surface and subsurface processes across large scales. Moreover, based on our findings it can be concluded that any future sensitivity and uncertainty analysis on the model sho uld account for spatial variability in the parameter values within the contrasting hydro geologic regions of the basin. Chapter 4 documents results from the application of Slim Fast, a fully coupled surface subsurface particle tracking scheme in conjunctio n with ParFlow.CLM, to the

PAGE 20

20 SFRB to identify dominant sources of river water, their travel time, flow paths and mixing dynamics under varying hydrologic conditions. The results of particle tracking experiments identify first order controls over river water sources origins, flowpaths, and travel times and how these controls, and associated source water dynamics vary with space and time. The final chapter of this dissertation presents a summary of the findings of this study along with some general conclusions Based on the findings of this study recommendations are made for future research needed to build on the research presented in this dissertation.

PAGE 21

21 CHAPTER 2 GEOLOGIC AND CLIMATI C CONTROLS ON STREAM FLOW GENERATION PROCESSES IN A COMPL EX EOGENETIC KARST B ASIN 2.1 Backgr ound Streamflow at any given location and time is comprised of surface and subsurface contributions from various sources. The ability to identify the factors controlling these contributions is key to successfully understand the stores, fluxes, flowpaths an d travel times of water and solutes through hydrologic systems. Hydrological processes are highly non linear and interactive making it difficult to predict the emergent behavior of the system under alternative external stressors. Representation of these c oupled processes in an integrated physical model provides an efficient way to explore their behavior and interactions at various space and time scales under varying physical settings as well as possible future land use, water use and climatic conditions (VanderKwaak and Loague, 2001 ; Loague et al. 2005 ; Maxwell and Mi ller, 2005 ; Li et al. 2008 ; Jones et al. 2008 ; Sudicky et al. 2008) For example, Bolger et al. ( 2011) developed a fully integrated surface subsurface HydroGeosphere model to study historic hydrologic conditions in the San Joaquin Valley in California. M odel simulations helped understand evapotranspiration (ET) losses and associated root zone processes in the vadose zone as well as surface water groundwater interactions along the rivers and within wetlands areas during historic and pre development condi tions. Goderniaux et al. ( 2009) developed a fully integrated surface subsurface HydroGeosphere model for the Geer Basin in eastern Belgium to study impacts of climate change on groundwater resources. The integrated flow simulations indicated potential for significant reductions in groundwater levels and surface water flow rates in the basin by the year 2080. Sudicky et al. ( 2008) developed a fully

PAGE 22

22 integrated InHM model for the Laurel Creek watershed in southern Ontario, Canada. The model was used to understand the transport of surficial contaminants through the fully integrated surface subsurface environment of the basin for a variety of precipitation patterns. In addition to ph ysical modeling, data driven approaches such as end member mixing analyses (EMMA; Hooper et al. 1990 ; Christophersen and Hooper, 1992 ) which use natural and/or isotopic tracers offer alt ernative ways to explore the hydrological functioning of natural systems under varying hydrogeological settings. For example, Banks et al. ( 2011) analyzed time series of hydrologic data such as river flow, groundwater elevation, rainfall; water chemistry data such as temperature, electric conductivity, pH, total dissolved solids; stable isotope data along with tracer based techniques to understand spatiotemporal connectivity between surface and groundwater in Rocky Rive r catchment in South Australia. They found that the connection between surface and groundwater varied along the length of the river and was affected by hydrogeological and hydroclimatic conditions within the basin. Capel l et al. ( 2011) used multivariate tracer analysis to study dominant runoff generation processes in two contrasting geological regions within the North Esk catchment, north east Scotland. Results showed a contrast in major sources of streamflow water with near surface processes dominant in uplands and groundwater fed baseflow dominant in lower parts of the catchment. Although effective for exploration and to gain system insight, the data driven nature of these techniques limits their effectiveness as predic tive tools. In recent years there have been an increasing number of studies that have used data driven approaches within a physical modeling framework as additional diagnostic

PAGE 23

23 tests for testing the representation of hydrologic process within simplistic mo dels for small basins. For example, Seibert and McDonnell ( 2002) reservoir volumes to improve the interna l consistency and performance of a conceptual three box model for the ~3.8ha Maimai basin, in New Zealand. They concluded that the use of the soft data in addition to the hard data including runoff and groundwater elevation significantly improved their mod uncertainty McMillan et al. ( 2011) demonstrated the use of available rainfall, flow and soil moisture data for parameter estimation as well as identification of more scientifi cally defensible model structures for the ~0.25 km 2 Satellite Right subcatchment of the Mahurangi river basin, New Zealand Werner et al. ( 2006) used field observations with an integrated model, MODHMS, to identi fy river reaches receiving groundwater influx in the ~420 km 2 Pioneer Valley, Australia. However, they did not investigate the spatio temporal characteristics of the groundwater surface water flow dynamics in the basin In this study we use a fully integr ated 3D surface water groundwater land surface model ParFlow .CLM (Maxwell and Miller, 2005 ; Kollet and Maxwell 2008a ; Kollet and M axwell, 2006) to evaluate the interacting geologic, climatic and vegetative controls on streamflow generation processes in a complex 3700 km 2 eogenetic karst basin in North Central Florida. The period of study January 1, 2000 through December 31, 2008 was selected to include diverse climatologica l conditions spanning several extreme events. For instance year 2004 was an extremely wet period during which four hurricanes and a tropical storm crossed the state of Florida causing flooding over much of the study area. Conversely, years 2006 2007 were one of the driest periods ever recorded for the

PAGE 24

24 basin and were followed by relatively wet conditions as a result of a tropical storm which occurred in September 2008. In addition to traditional model evaluation criterion, such as comparing field observati ons to model simulated streamflow and groundwater groundwater interactions over space and time using EMMA that we performed using observed specific conductivity (SC) differences among surface and s ubsurface water sources throughout the domain. After benchmarking ParFlow against field observations, we used the model to address following questions about integrated hydrologic processes in this complex regional river basin; (1) how do spatiotemporal patterns of surface water: groundwater streamflow fractions vary throughout the basin; (2) how do geology, land cover and climatic variability interact to control surface water groundwater interactions in the basin, and (3) how do ge ology, land cover and climatic variability control spatiotemporal patterns of major hydrologic water budget components such as ET and groundwater recharge in the basin 2.2 Santa Fe River Basin P hysical S ystem and C onceptual M odel The 3700 km 2 Santa Fe Rive r Basin (SFRB) consists of two distinct yet linked hydrogeologic units ( Figure 2 1): the upper confined region (CR) and the lower unconfined region (UR) which are separated by a topographic break known as the Cody Escarpment (CE). The climate of the SFRB i s warm and humid with mean annual precipitation of 1356 mm, receiving most of its rain from mid may through mid October. The rainy season is followed by a cool, dry period from mid October through mid February and a moderately cool, wet period from mid Feb ruary through mid May characterized by periodic extra tropical fronts. The basin experiences rare/episodic

PAGE 25

25 extreme events in form of hurricanes, and seasonal and multi annual events caused by the El Nino Southern Oscillations (ENSO) phenomenon. In the upp er two thirds of the SFRB the Surficial Aquifer System (SAS) is a perched aquifer separated from the regional Upper Floridan Aquifer System (UFAS) by the presence of a low hydraulic conductivity Intermediate Aquifer System (IAS; Water resources associates, 2007 ) composed of Hawthorne clays ( Figure 2 2 A ). Major soils covering the CR are classified as poorly or very poorly drained with very low vertical hydraulic conductivity resulting in low infiltrati on (Arthur et al. 2005) Stream networks, lakes and swamp systems are well developed in the region (Upchurch, 2007) The abundance of swamps, wetlands, forests, and poorly dra ined soil result in high water table conditions even during extended dry periods. The CE is formed by marine, fluvial and karst related erosion of the Hawthorn clay which make s up the IAS (Upchurch 2007) The SAS is generally absent in this region and the IAS, if present, is in the form of local clay rich strata (Upchurch 2007) The soil types are a mix of poor or well drained depending on location. The CE is characterized by the presence of numerous sinkholes, sinking streams, siphons, springs and other karst features which have the potential to rapidly recharge the UFA S (Upchurch 2007) During baseflow conditions the upper Santa Fe River (SFR) is captured in its entirety by the SFR Sink (hereafter referred to as the sink) in the CE region. The river emerges from the SFR Rise (hereafter referred to as the rise), located approximately 6 km downstream of the sink, as a first order magn itude spring. Hydrological behavior of the river in the transition zone is highly dynamic with no flow entering the sink during extreme drought conditions and with a portion of the river

PAGE 26

26 occasionally bypassing the sink as overland flow to the lower SFR dur ing large flood events. In lower third of the SFRB both the SAS and the IAS are absent and the UFAS system becomes an unconfined aquifer, with a thin overlying sand layer and direct exchange with surface water in the river channels. Major soils covering th e UR region are classified as well or excessively well drained with high hydraulic conductivity resulting in high infiltration rates and virtually no surface runoff (Arthur et al. 2005) Thus rainfall in excess of ET rapidly moves through the vadose zone to the UFA, which feeds the river via series of springs and diffuse groundwater discharge. The UFA in this region is known to have major conduits creating internal subsurface drainage networks connected to the main ri ver ( Figure 2 2B ), some of which extend beyond the topographic boundary of the basin (Upchurch 2007 ; Upchurch et al. 2008 ; Meyer et al. 2008) In contrast to teleogenetic karst systems that have low intergranular porosity, because they have undergone deep burial and alterations, eogenetic karst systems that exist in Florida have not undergone deep burial and have retained their higher matrix porosity (Vacher and Mylroie 2002 ; Florea and Vacher 2007) As such primary matrix porosity in the SFRB plays an important role in water storage, flow through the subsurface environment, and exchange with the river system 2.3 Metho ds 2.3.1 River S ource W ater M ixing D ynamics: End Member Mixing Analyses During low flow conditions continuously flowing rivers typically receive baseflow contributions from the subsurface environment which helps sustain their flow. Because of their origin and a ge, the hydrologic modeling literature often refers to these baseflow (e.g. Desmarais and Rojstaczer

PAGE 27

27 2002) River source water characteristics change during sto rm events when rain water, (e.g. Desmarais and Rojstaczer 2002) typically dominates the river hydrographs. In this study our goal is to understand factors controlling contributions from the UFAS. Many hydrologic studies have used SC of water as a natural tracer in hydrog raph separation analysis (Nakamura, 1971 ; Pilgrim et al. 1979 ; Matsubayashi et al. 1993 ; Cox et al. 2007) Hoop er et al. ( 1990) and Christophersen and Hooper ( 1992) detail the mathematical basis and assumptions of EMMA for source water identification. For brevity, in this paper we provide a brief description of the me thod, equations used and information pertinent to application of the method in the SFRB. Rain water in the SFRB is known to have low ionic concentration and therefore, low SC (Ritorto et al. 2009) In contrast water in the subsurface environment of UFAS, because of continuous dissolution of surrounding matrix, has higher ionic concentration and therefore high SC (Langston et al. 2012) Given that SC values representative of ibutions in total streamflow (Eq. (2 3) ) can be approximated by solving binary mass balance Eqs. (2 1 ) and (2 2) (Desmarais and Rojstaczer 2002) (2 1) (2 2)

PAGE 28

28 (2 3) groundwater present in streamflow at a given location and time. SC SW SC GW SC river water, respectively. SC measurements at station 1500 in the CR ( Figure 2 1) during high stream flow conditions caused by hurricanes in October 2004 were used to represent SC SW values (0.07 m s/cm) in Eqs. (2 2), and (2 3) This value was chosen to represent SC value for R and UR) were flooded with rain water received from four hurricanes that passed through the study area over very short time. The SC measurements in the SFR during the peak of October 2004 storm event ranged from 0.07 m s/cm (at stations 1000, 1500) to 0.1 m s/cm (at station 2800). To account for known effects of local variations in geology on (Gulley et al. 2012) we used C for each UR river locations. This approach of using independent end member estimates for distinct locations has been previously used successfully to account for soil and aquifer heterogeneity (e.g. Soulsby and Dunn, 2 003) For each UR river location, SC values corresponding to lowest five percentile flow values were averaged, to represent SC GW The SC end member values thus calculated for SC GW at stations 1975, 2500, 2800, respectively, were 0.47, 0.39, and 0.37 m s/cm All the above mentioned calculations used discrete (monthly) SC measurements obtained from the Suwannee River Water Management District

PAGE 29

29 ( SRWMD ) to calculate surface and groundwater fractions at stations 1975, 2500 and 2800 over the entire duration of thi s study 2.3.2 Integrated 3D M odel ParFlow .CLM We used ParFlow to simulate the 3D integrated surface water groundwater system in the SFRB. Full details on model physics and how the surface and subsurface components are integrated can be found elsewhere ( Ashby and Falgout 1996 ; Jones and Woodward 2001 ; Kollet and Maxwell 2006) In this section we present a brief summary of equations used by the model. To simulate flow through v ariably saturated subsurface environment, ParFlow solves the 3D Richards equation using a cell centered finite difference scheme in space and an implicit Euler scheme in time. Overland flow is simulated using 2D Kinematic wave equation and Manning s equati on. Surface and subsurface equations are integrated using a free surface overland flow boundary condition ( Eqs. (2 4 ), (2 5), and (2 6) ) ( 2 4) ( 2 5) ( 2 6) p is the subsurface pressure head [L]; S w p ) is the degree of saturation [ p as defined by van Genuchten model (Van Genuchten 1980) ; S s is the specific storage [L 1 ]; is the porosity [ ]; k(x) is saturated hydraulic conductivity [LT 1 ]; k r is the relative permeability as a function of subsurface pressure

PAGE 30

30 head[ ]; z is depth below the land surface [L]; q s is a general source/sink term [T 1 ]; q r represent s rainfall and evapo rat ive fluxes [LT 1 ]; || A,B || indicates the greater of A and B; S f,x and S f,y are the friction slopes in x and y direction [ ]; n is the Manning s coefficient [TL 1/3 ]; an v is the depth averaged velocity vector of surface runoff [LT 1 ] ; o in Eqs. (2 5) and (2 6) is the pressure head or water ponding depth at the land surface [L] ParFlow .CLM incorporates a modified version of the Common Land Model (CLM ; Dai et al. 2003 ) into ParFlow to simulate near land surface energy and water balance. CLM simulates evaporation losses from the ground surface and vegetation canopy; transpiration losses from plants; snow accumulation and melt processes; and latent, sensible, and ground heat fluxes (Kollet and Maxwell 2008a) Near land surface atmospheric fluxes are calculated as a function of atmospheric variables that CLM requires as input (precipitation, air temperature, pressure, wind speed, specific humidity, and solar radi ation) and soil moisture calculated by ParFlow Details of equations used to simulate various water and energy fluxes in CLM and coupling of the two models are provided elsewhere (Maxwell and Miller 2005 ; Kollet et al. 2009; Kollet and Maxwell, 2008a) ; however, equations relevant to work presented here are summarized in the following. At any given land surface location, the available net radiation (given by sum of both short and long wave radiations) is dissipated through the sum of latent, sensible, and soil heat fluxes: ( 2 7) w here R n is net radiation at the land surface (W/m 2 ); LE is latent heat (W/m 2 ), H is sensible heat flux (W/m 2 ), and G is t he ground heat flux (W/m 2 ). The energy balance

PAGE 31

31 at the land surface is coupled to water balance in the subsurface environment through the soil moisture content at or close to ground surface Each of the terms on the right hand side of Eq. (2 7 ) are independently calculated by CLM usi ng atmospheric forcing, soil, vegetation characteristics and soil moisture provided by ParFlow LE is calculated as: ( 2 8) ( 2 9) where, L v is lat ent heat of evaporation (J/kg) and E is the sum of evaporative fluxes from the foliage ( E c (kg/m 2 s) if vegetation is present) and ground E g (kg/m 2 s) Evaporative fluxes from the ground are cal culated as: (2 10) Where a is the intrinsic density of air (kg/m 3 ); q g is the air specific humidity at the ground surface (kg/kg); q a is the air specific humidity at reference height z q obtained from atmospheric forcing (kg/kg); r d is the aerodynamics resistance at z q (s/m) Evap orative fluxes from the foliage are estimated as the sum of evaporative flux from wet foliage E w (kg/m 2 s) and total transpiration flux E tr (kg/m 2 s ) from the root zone cells: ( 2 11 ) ( 2 12 ) ( 2 1 3 ) Where f root,j* is effective root fraction in root zone layer j

PAGE 32

32 (2 14) (2 15) (2 16 ) Where f is the vegetation fraction ( ); L SAI is the stem plus leaf area index ( is the step function (one for positive and zero for zero and negative arguments); is the potential evaporation from wet foliage (kg/m 2 s) ; is the we tted fraction of the canopy ( ); L d is the dry fraction of foliage surface ( ); r b is the conductance of heat and vapor flux from leaves (s/m); and r s is the stomatal resistance (s/m). i fc wp respectively, are soil moisture conte nt in layer i, moisture content at field capacity, and permanent wilting point; f root,j is root fraction parameter which varies exponentially with the depth from land surface z ( ) in root zone layer j ; a and b are coefficients dependent upon vegetation ty pe. In ParFlow.CLM the term q s in subsurface mass balance equation (Eq. (4)) is expressed as (2 17) Where q g (s 1 ) is the flux of water infiltrating at the land surface after accounting for precipitation, canopy throughfall, and/or surface runoff 2.3.2.1 Model s etup, p arameterization and i nitialization In this study ParFlow .CLM wa s developed for a domain that extends beyond the SFRB topographic boundaries in order to include major springsheds that are thought to provide subsurface f low to the river in the UR of the basin (Meyer et al. 2008 ; Upchurch et al. 2008) The bottom of the domain was set at 70 m below mean sea level to

PAGE 33

33 include the top of the active subsurface zone believed to contribute flow to the SFR. A to discretize the domain resulting in 56x50x150 cells, in x, y, and z dimensions, respectively. Thus the domain included a total of 420,000 rectilinear elements covering a total land surface area of 6300 km 2 ( Figure 2 1). The land surface elevation in the domain ranges from 4 m to 74 m above mean sea level, resulting in a maximum domain depth of 144m. A 10 meter (m) resol ution digital elevation model (DEM) dataset obtained from the United States Geological Survey (USGS) was resampled to a resolution of 1500 m using bilinear interpolation approximation in the Spatial Analyst Module in ArcInfo ( Figure 2 3 A ). The main channe ls, as represented by the National Hydrography Dataset (NHD), were incised into the DEM by setting the bottom elevation of the channels equal to the elevations used in a HECRAS model previously developed for the region ( Water resources associates 2007) and personal communication Clay Coarsey, ( SRWMD, 2012 ). Manning s surface roughness coefficients (n) needed for ParFlow were taken from the n values assigned to the channels in the existing HECRAS model ( Water resources associates 2007) and personal communication Clay Coarsey, ( SRWMD, 2012 ), with the same values assigned to the channels and the floodplains by region (see Table 2 1). In the CR the n value approximately corresponds to that of channel with heavy brush and timber (Chow et al. 1988) In the UR the value corresponds to something between a clean straight and clean winding channel (Chow et al. 1988)

PAGE 34

3 4 Information regarding the spatial extent and depths of the SAS, IAS and UFAS in the region were taken from maps available from the Floridan Aquifer Vulnerability Assessment (FAVA) conducted by the Florida Geological Su rvey (FGS ; Arthur et al. 2005 ). These maps were available at resolutions ranging from 30m to 390m and were all resampled to resolution of 1500 m using bilinear interpolation approximation in the Spatial Analyst Modu le in ArcInfo. Where conflicts arose between the USGS DEM used to define the land surface and the top elevations, horizontal extents and thicknesses of the SAS and IAS obtained from FAVA, it was assumed that SAS occupied at least the top 8 meters below the USGS defined land surface over the region where the IAS was mapped. The resulting geological representation is shown in Figure 2 2 A The UFAS is known to have networks of conduits in the UR and CE regions of the domain. A digital map of the approximate lo cations of these conduits was obtained from (Meyer et al. 2008) who compiled their best estimates of conduit locations based on results of dye tracer studies, groundwater chemistry and potentiometric surface measureme nts in the domain ( Figure 2 2 B ). All the conduits were assumed to be 10m in height and to occur between depths of 5 m to 15 m below mean sea level, based on previous field investigations (Langston et al. 2012) Since ParFlow does not contain algorithms for turbulent conduit flow, conduits were represented as high hydraulic conductivity zones in this study. This assumption should be reasonable given the goal of the study to estimate controls on regional variabi lity of streamflow generation processes and water budget components over a large river basin and the resulting, relatively coarse 1500 m horizontal discretization of the model. Finer discretization and more accurate representation of turbulent conduit flow processes would be required to

PAGE 35

35 investigate smaller scale interactions (i.e on the order of meters to 10s of meters) between the river, conduits and surrounding karst matrix. The subsurface environment was assumed to be isotropic and homogeneous within ea ch aquifer type (SAS, IAS, UFAS, conduits). The effective aquifer properties for each type were estimated based on literature values from previous modeling, observational and experimental studies in the region and are summarized in Table 2 1. A 30 m reso lution digital map depicting land cover in the domain for was obtained from the Florida Fish and Wildlife Conservation Commission for the year 2003 ( Figure 2 3 B ). These data were resampled to 1500 m resolution using the nearest neighbor interpolation appr oximation in the Spatial Analyst Module in ArcInfo. The domain was initially described by 22 different land cover categories. Land cover categories were merged/reclassified into nine distinct land cover types (i.e. evergreen needleleaf forests (31%), grass lands (14%), open shrublands (13%), permanent wetland (13%), mixed forests (11%), urban areas (8%), croplands (4%), barren or sparsly vegetated (3%), and water bodies (2%)) to correspond to the vegetation types included in ParFlow .CLM default vegetation pr operty database, which follow the standard of the International Geosphere Biosphere Program (IGBP Running et al. 1994 ). Spatially distributed climatic forcing data obtained from the North American Land Data Assimil ation System (NLDAS ; Cosgrove et al. 2003 ) were used to drive the model. The NLDAS dataset is available at an hourly temporal resolution and at a spatial resolution of 1 2 km x 1 2 km. The forcing data was spatially resampled to a resolution of 1500 m by using a bilinear interpolation approximation.

PAGE 36

36 To initialize ParFlow requires estimation of the initial pressure distribution throughout the domain. The May 2002 UFAS potentiometric surface map published by the SRWMD was used as the initial condition in the UFAS. This is closest date to the beginning of the study period for which a potentiometric surface map was available. For the SAS and IAS no groundwater elevation maps were available. Therefore, to initialize the m odel the water table in the SAS was set at a depth equal to 10% of total SAS thickness (i.e. SAS was 90% saturated everywhere in the CR), and the IAS was assumed to be in hydrostatic equilibrium with the UFAS. Constant pressure boundary conditions equal t o the initial conditions were applied along all lateral domain boundaries. A no flux boundary condition was applied to the domain bottom and an overland flow boundary condition was applied at the land surface To minimize the influence of initial condition assumptions the model was spun up by repeatedly simulating the entire 9 year study period (of 2000 through 2008 ) using the spatially distributed NLDAS forcing data until a quasi dynamic equilibrium was reached, i.e. until the difference between the beginn ing and ending subsurface water storage (over the 9 year period run) dropped to less than <0.05 % of average subsurface storage and the simulated streamflow and groundwater elevations time series at multiple locations across the domain showed minimal vis ual difference between two consecutive 9 year runs. These criteria were met by the end of the second 9 year simulation period, thus all the model results presented in this paper are the model results obtained from round 3 of model spinup

PAGE 37

37 2.3.2.2 Estimating s urfac e w ater and g roundwater c ontributions to s treamflow u sing ParFlow .CLM The groundwater fraction of flow at a given river location on a given day includes the subsurface contribution that the location receives from its entire upstream contributing area on th at day. A portion of this contribution is received at upstream channel locations which then travels downstream as channel flow and a portion is received as local subsurface contribution. In this study we approximated instantaneous daily subsurface flow at a particular river location as the sum of subsurface flow contributions received over the entire upstream river reaches on that day. The subsurface flow to/from each channel cell was calculated using the saturated hydraulic conductivity and hydraulic head gradients provided by ParFlow at each river cell in the Richards equation on daily basis. This methodology neglects travel time within the stream channel, which may have some impact on instantaneous comparisons with the EMMA, but should not affect event sc ale, monthly, seasonal or annual subsurface contribution estimate s 2.3.2.3 Benchmarking ParFlow .CLM The adequacy of model simulated regional surface and subsurface flow processes within the domain was tested using traditional benchmarks including comparison of m easured and model simulated daily streamflow and groundwater surface elevations, at multiple locations across the study domain ( Figure 2 1). In addition comparison between ParFlow .CLM simulated and EMMA estimated groundwater flow contributions at multiple river locations was used to evaluate the accuracy of model simulated groundwater surface water interactions under varying hydrologic conditions. Daily streamflow data were obtained from the USGS and

PAGE 38

38 SRMWD, daily groundwater elevation data were obtained fr om the SRWMD, and monthly SC measurements were obtained from the SRWMD. The extreme climate variability (e.g. hurricanes in 2004 versus historic droughts in 2007) during the study period allows testing of the robustness of the model under highly variable hydrologic conditions. It should be noted that results presented in this paper were obtained without formal model calibration. Rather than optimally f itting a model structure and parameters to a calibration dataset we aim ed to develop a model that incor porates physically plausible hydrological characteristics of the domain based on best known existing information, and to use this model to 1) test our conceptual model of the hydrologic functioning of the basin, 2) gain improved insights about the interac tions of coupled hydrologic processes in the basin, and 3) to determine where additional field experiments may be needed to refine both the conceptual and numerical models of the basin 2.4 Results and D iscussion 2.4.1 Streamflow Measured and simulated flow hydrogr aphs at 5 locations along the SFR and 1 location on New River (station 1000 in Figure 2 1), a major tributary to SFR, were compared to evaluate ParFlow .CLM flow predictions ( Figure 2 4). Stations 1000 and 1500 are located within the CR of the domain, Stati on 1898 is located along the CE, and Stations 1975, 2500 and 2800 are located within the UR of the domain ( Figure 2 1). Model simulated streamflow at all CR and CE locations showed good agreement with measured data during both low and high flow conditions ( Figure 2 4 stations 1000, 1500, 1898). Timing and peaks of hydrographs were accurately simulated by the model

PAGE 39

39 at all locations with very few exceptions. ParFlow .CLM underpredicted flow during a major event in middle of year 2003 2004 at station 1000; how ever comparison of the hydrographs at other confined locations during same time period indicates that the poor performance at the given location might be attributed to inadequate representation of rainfall around station 1000. In general model predictions at all UR stations showed reasonable agreement with measured data during high flow conditions ( Figure 2 4 stations 1975, 2500, 2800). At station 1975, the model consistently simulated higher than observed streamflow during peaks of storm events. This is l which cannot simulate small scale short duration interactions between conduits and the porous matrix that are known to occur in the region between station 1898 and 1975. Nevertheless, the model adequately simulated the transition from high to low flow conditions at 1975. At stations 2500 and 2800 in the UR the model adequately simulated the timing and peaks of storm events but consistently underpredicted streamflow during the end of the falling limbs of s torm hydrographs. Thus the model as currently configured, is missing a transient source of groundwater to the streams in the UR. This could be due to the constant fixed pressure head boundary conditions applied at the lateral boundaries of the domain whi ch may result in some loss of groundwater from the domain during high rainfall periods, or inaccuracies associated with the extent and thickness of the IAS and the location of the CE boundary which could result in less water infiltrating into the UFAS in the upper portions of the UR 2.4.2 Groundwater E levation ParFlow .CLM simulated groundwater elevation in the UFAS were evaluated at six CR locations (wells E1 through E 6 in Figure 2 1) and three UR locations (wells E 7

PAGE 40

40 through E 9 in Figure 2 1). With exception s of wells E3, E7 and E 8 the model simulated groundwater surface elevations were in good agreement with field observations (Figure 2 5) Both measured and model simulated groundwater elevation at all CR wells show smooth annual variations corresponding to slow rise and fall of regional groundwater surface in response to wet and dry hydrologic periods. In contrast, both measured and simulated groundwater elevation at all unconfined well locations indicate direct response to individual rain events as well as smoother annual variations corresponding to the slow rise and fall of regional groundwater surface in response to wet and dry hydrologic periods. Well E3 and E 8 are close to the CE and therefore highly sensitive to the degree of confinement and extent of karstification along the CE, which are uncertain. Due to their locations these wells are expected to show greater response (as compared to well confined well locations) to large storm events. However the piezometric head at these locations is always high er than observed and showed even more response to storm events than measured data. We hypothesize that the difference between the model simulated and measured elevations at these location is due to uncertainty in extent of confinement and conduit locations in this area and perhaps the coarse resolution of the model. Well E 7 is located in a physiographic region known as Wacassassa Flats (Figure 2 1) It is assumed that there are remnants of the IAS in this region, which causes the observed local piezometric surface high in the region. However these IAS remnants are not mapped in the FGS aquifer datasets used in this study, in fact there is no published information on thickness of IAS in this region Therefore we did not include them in the

PAGE 41

41 conceptual model underlying ParFlow As a result the model was unable to capture the observed significant rise in the potentiometric surface in Wacassassa Flats after storm events and showed marginal rise in groundwater elevation in response to rise in regional groundwater surface rise. These results, as well as those for wells E3 and E 8 underscore the need for more accurate, spatially explicit hydrogeologic characterization of the region. In addition to model evaluation of groundwater elevations at point locations we eval uated the model simulated regional UFAS response to wet and dry hydrologic conditions resulting from extreme climate events. Here we present the piezometric surface of the UFAS during wet conditions caused by hurricanes in year 2004 ( Figure 2 6 A, B, and C ) and during extremely dry conditions in years 2006 and 2007 ( Figure 2 6 D, E, and F ). In general the piezometric surface shows regional flow from east to west, with the river draining the UFAS in the UR The transient effects of rainfall on UFAS groundwat er elevations are visible during peaks of hurricanes in year 2004 as well as tropical storm Fay in year 2008 ( Figure 2 6 B and F ). The effect is most prominent around high hydraulic conductivity zones, representative of conduits that breach the IAS in the CR, in the northwest portion of the domain. It is less prominent around high hydraulic conductivity zones representative of conduits in the UR in the south central part of the domain, and is not visible in the well confined western portion of the basin. January 2005, a few months after the end of the 2004 hurricane season, shows smoothing of the regional piezometric surface, but the piezometric head around the high hydraulic conductivity zones is still significantly elevated compared to the pre hurricane condition. These observations suggest that high hydraulic conductivity zones in the

PAGE 42

42 UFAS provide hotspots for rapid transport of water and surficial contaminants into the regional groundwater system during large storm events, particularly where they breach the IAS. Comparison of Figures 2 6 D and E shows that as the drought conditions extended during 2006 2007 the river drained water from regions further away, extending into the CR 2.4.3 Surface V ersus S ubsurface C ontributions to S treamflow Comparison of measure d versus predicted ratios of total daily subsurface flow to total streamflow was used as an additional diagnostic test of the adequacy of the conceptual model and parameterization underlying ParFlow .CLM over a range of hydrologic conditions. As shown in Fi gure 2 7, river reaches upstream of station 1898 receive no significant contributions from the subsurface, indicating rainfall and near stream overland flow to be the dominant source of streamflow in the CR of the domain. At station 1898 the river shows a net loss of stream flow to the subsurface during storm events as streamflow generated from surface processes in the CR is lost to the UFAS through high conductivity breaches in the IAS beneath the river channel. Thereafter, each downstream station in the U R (1975 through 2800) shows a consistent gain in subsurface flow contribution. Streamflow at all gages in the UR shows an immediate decrease in subsurface contributions immediately after the onset of largest storm events, as the surface pulse of storm flo w passes and hydraulic gradients from the local UFAS to the stream are reduced. This is followed by an extended period of elevated subsurface contributions from the UFAS to the river after the storm. Table 2 2 summarizes the net annual gain/loss of subsu rface flow at various river locations. Note that for station 1898 at the CE, values represent the total gain/loss in all river cells upstream of that gage location, whereas for the remaining UR stations the values

PAGE 43

43 represent total gain/loss of subsurface fl ow downstream of station 1898. The negative values of subsurface flow contributions for station 1898 indicate that in all years except 2006 the river shows a net loss of surface water generated in the CR to the subsurface environment. Figure 2 8 compare s the ratio of EMMA estimated (using SC measurements and observed streamflow) and ParFlow .CLM predicted groundwater to total streamflow ratios at the three UR river locations. This figure indicates that the model accurately simulates both the timing and ma gnitude of the surface to groundwater mixing ratio groundwater contributions from the UFAS d uring low flow conditions. Subsurface contributions decline during storm events due to both the high volume of streamflow coming off the CR as a result of the storms and the reduction of the hydraulic gradient towards the river as the surface water pulse passes. Both of these mechanisms are accurately represented by ParFlow .CLM, for example see prominent drop in subsurface flow contributions during 2004 hurricanes in Figure 2 7 and corresponding drop in groundwater to total streamflow ratio in Figure 2 8. Note that at station 1975, directly downstream of the CE, groundwater to total streamflow ratios are significantly underestimated by ParFlow during low flow periods. However during these periods total streamflow is extremely low (See Figures 2 4 and 2 9) making the groundwater to surface water ratio highly sensitive to small errors. Figure 2 9 shows observed and model simulated streamflow, along with total subsurface flow contributions as estimated by EMMA and simulated by ParFlow at

PAGE 44

44 three UR locations f or two different periods: 2002 2003 (a wet year 2003 following an equally wet year 2002; see total annual rainfall in Table 2 3) and 2007 2008 (a wet year 2008 following two dry years 2006 2007). These results show that ParFlow more accurately simulates both total streamflow response to storm events, and groundwater contributions to streamflow, after dry antecedent conditions (2008) than after wet antecedent conditions (2003), particularly for stations 2500 and 2800. Comparison with EMMA results indicate that ParFlow underestimates groundwater contribution to streamflow, especially when wet years follow wet years. This suggests an underestimate of interannual storage in the UFAS system, with unrealistic losses either occurring across the southwestern doma in boundar y or by ET (rather than recharge to the UFAS) in the transition region between the CR and UR where the extent thickness and hydraulic conductivity of the IAS is highly uncertain 2.4.4 Analysis of W ater B udget C omponents for the S FRB Use of fully inte grated models makes it possible to elucidate significant spatiotemporal patterns in coupled near land surface energy and hydrologic processes and resulting water balance components. In this section we present annual ( Table 2 3) and monthly averaged water b alance components (rainfall, ET, runoff, surface storage (S t s ) and subsurface storage (S t ss ), Figure 2 10) integrated over the study domain and compare them to commonly reported values for the study area wherever possible. Average annual rainfall received by the domain during the study period (January 2000 through December 2008) was 1253 mm as compared to the 1356 mm long term average annual rainfall reported for the study area (Schneider et al. 2008) The average annual ET loss simulated by ParFlow .CLM was 982 mm (78% of applied rainfall) as compared to long term average annual ET losses of 1041 mm (77% of

PAGE 45

45 rainfall as estimated by (Schneider et al. 2008) ). Thus ParFlow.CL M accurately simulates average ET fluxes from the domain, the most important water balance component after rainfall. Monthly averaged values for all water balance components for the domain are shown in Figure 2 10. A negative value in the figure indicates water loss from the domain and a positive value indicates water gain in the domain. Averaged monthly rainfall follow s a bimodal pattern with highest rainfall occurring during months of March and June, and 57% of the rainfall occurring during the months of June through September. Average monthly ET losses were highest during summer months (peaking during the month of July) and lowest during winter months of December January. In general the rainfall and ET patterns during the simulation period are consiste nt with long term trend reported for the basin (e.g. Tripathi, 2006) R unoff S t s S t ss and groundwater recharge all show a strong dependence on seasonal and interannual variability in rainfall and ET patterns over the basin ( Figure 2 10; Table 2 3). In general S t s and runoff from the domain outlet followed the rainfall pattern with wet months (e.g. months of August and September) and wet years (e.g. 2004 and 2005) showing an overall increase in S t s and higher runoff losses in response to higher rainfall over the domain. Patterns in S t ss indicate higher recharge during months of June through August in response to higher rainfall rates, despite higher ET losses. Net discharge from the subsurface occurs during periods o f rainfall deficit in April May and October and November (i.e. when ET > rainfall). Annual water balance components, summarized in Table 2 3, show a strong influence of interannual climate variability, with loss of S t s and S t ss from the domain

PAGE 46

46 during dry years (2000, 2001, 2006, 2007) and gain of S t s and S t ss in the domain during wet years (2002, 2004, 2005, and 2008). It is interesting to note the dependence of water balance dynamics on antecedent conditions. Year 2005, a wet year following a wet year, s howed a smaller gain in S t ss and more runoff and ET losses compared to year 2004 which received similar amount of rainfall. Thus the high antecedent S t ss conditions at the beginning of 2005 resulted in higher runoff and ET losses. Similarly year 2003 (an average rainfall year following an average rainfall year) showed more runoff and ET losses as compared to year 2002 (an average year following a dry year) resulting in very small net gain of storage in the domain. In contrast years such as 2002 and 2008 w hich followed dry conditions showed a net gain in S t ss and smaller runoff and ET losses compared to 2003 which received similar amount of rainfall but followed an average rainfall year 2.4.5 Geological C ontrol O ver S patial P attern of M onthly ET L osses and S ubs urface R echarge Spatial distribution of average monthly ET losses from the domain and average monthly groundwater recharge (net change in S t ss ) are shown in Figures 2 11 and 2 12, respectively. High ET losses and low groundwater recharge occur in the CR wh ere the presence of a confining unit and relatively low hydraulic conductivity result in a high water table and poorly drained conditions. In contrast the absence of the confining layer over the UFAS and the presence of high hydraulic conductivity zones (i e. conduits) in the UR create lower water tables and a well drained subsurface environment. Thus the contrast in the geology between the two regions strongly controls the balance between ET and groundwater recharge, which in turn controls streamflow gene ration

PAGE 47

47 mechanisms, i.e. dominance of surface flow in the CR and dominance subsurface flow in the UR. Monthly variations in climate modulate the water balance components, but do not alter the sharp contrast between the two geologic regimes. The months of J une through September receive more rainfall than ET losses and showed a net subsurface recharge ( Figures 2 10 and 2 12). Consistently higher recharge in the UR than in the CR during these summer months occurs due to drainage of subsurface storage by the hi ghly permeable UFAS to river. During April May and October November ET losses equal or exceed rainfall inputs and thus groundwater discharge to the river results in a loss of storage from the subsurface 2.4.6 Vegetation and G eological C ontrol O ver W ater T able D epth and ET D ynamics In addition to geologic controls on energy and hydrologic feedbacks between land surface and atmosphere (through control on water table depth), vegetation also play s an important role (Kollet and Maxwell, 2008a) Figure 2 13 A shows the relationship between averaged daily water table depths (WTD) from land surface and averaged daily ET loss for each land surface cell in the model domain. Blue and red circles, respectively, are cells in the CR an d UR of the basin. As discussed in section 4.5 the shallower WTD in the CR results in higher ET losses as compared to the UR. Similar to results presented by (Kollet and Maxwell 2008a) we observe a strong dependenc y between ET losses and WTD up to an approximate critical WTD of 5 m. Beyond this critical depth ET losses become independent of WTD and are limited to withdrawal from moisture available in the unsaturated subsurface environment. Figure 2 13 B shows that l andcover exerts a secondary control on ET and WTD, with grass

PAGE 48

48 lands showing the lowest ET losses and forests and wetlands showing the highest ET losses. These findings illustrate a strong link between potential changes in anthropogenic activities (via chan ges in land cover) to changes in the most significant component in the water balance for the basin i.e. ET losses 2.5 Summary and C onclusions ParFlow .CLM an integrated 3D surface water groundwater land surface hydrologic model, was used to evaluate the inte racting geologic climatic and vegetative controls on streamflow generation processes and spatiotemporal water budget components over a large, complex eogenetic karst basin in North Central Florida. Comparison to field observations throughout the basin sho wed that the model successfully simulated the spatial and temporal dynamics of streamflow, ET fluxes and groundwater elevations throughout th e domain across extreme hydrologic conditions and without any formal model calibration Results of this study ind icate that geologic heterogeneity exerts primary control on streamflow generation processes in this basin through its influence on water table depth and ET fluxes. In the upper basin, where the regional karst aquifer is overlain by a thick confining layer, a persistent high surficial water table results in high ET and low groundwater recharge. Low topographic relief and relatively low surficial aquifer hydraulic conductivity result in both low surface and low subsurface connectivity of the larger watershed with the streams in the upper basin despite the high water table. Thus streamflow is episodic, with more than 95% generated by recent near stream rainfall that temporarily raises the water table above the land surface near stream channels, and virtually n o inter event subsurface baseflow. In the lower basin, where the karst aquifer is unconfined, deeper water tables result in lower ET and higher groundwater

PAGE 49

49 recharge. There is little surface connectivity between the watershed and streams in the lower basin and thus virtually all surface contributions to streamflow originate in the upper confined basin. However the high karst aquifer hydraulic conductivity results in extensive subsurface connectivity with the streams in the lower basin. As a result the majori ty of the streamflow in the lower basin is subsurface flow originating as diffuse infiltration through the epikarst, with the fraction of subsurface flow increasing from 52% in the upper portion of the UR to 77% at the outlet of basin. Climatic variabili ty was found to provide a secondary control on streamflow generation processes, resulting in significant seasonal and interannual variability in both the timing and sources of streamflow. In the upper confined basin stream flow occurs primarily as a resul t of extra tropical cold fronts in the spring and tropical storms in the fall when ET fluxes are relatively low. Persistent summer streamflow only occurs in the confined basin when consecutive years of above average rainfall result in water storage in the region in excess of high summer evapotranspirative demand. In the lower basin the fraction of surface flow (i.e. flow contribution from the upper basin) ranges from 12% in late spring (May) to 55% in early spring and fall (March and September), and ranges from an annual average of 52% during periods of consecutive wet years to 23% during periods of consecutive dry years. The ratio of surface water groundwater contributions t o streamflow estimated from EMMA were used to test the surface water ground wa ter dynamics simulated by ParFlow.CLM. Results indicated that the model accurately simulated both the timing and magnitude of the surface to groundwater mixing ratio response to storms at all UR stations. The model accurately simulated the dominant fractio n of

PAGE 50

50 during peaks of all storm events and the dominant fraction contributions from the UFAS during low flow conditions. However comparison of total surface water groundwater contributions to streamflow estimated f rom EMMA to those simulated by ParFlow .CLM revealed that ParFlow consistently underestimates groundwater contribution to streamflow after wet antecedent hydrologic conditions. We hypothesize that is due to underprediction of interannual storage in the UFAS system, caused by errors in the representation of the exten t, thickness and hydraulic conductivity of the IAS in both the Cody Escarpment and Wacassassa Flats regions of the domain These results underscore the usefulness of combining integrated hydrolog ic modeling with insights from data driven EMMA to develop a quantitative, predictive understanding of surface water groundwater flow dynamics within a complex eogenetic karst basin, and to indicate where additional field measurements may be necessary to refine geologic heterogeneity and to improve numerical model predictions. To the best of our knowledge this is the first time EMMA based hydrograph separation results have been used to evaluate the transient and spatially distributed river source water mi xing dynamics simulated by a fully integrated hydrologic model for such a large and complex basin. In the next phase of this research we will use global sensitivity/uncertainty analyses (Saltelli et al. 2005) to conf irm key model processes and parameters governing model prediction errors identified in this study. Data assimilation techniques such as Generalized Likelihood Uncertainty Analysis (Beven 1993) or Ensemble Kalman Filt ering (Graham and McLaughlin 1991) will then be used

PAGE 51

51 with new high frequency groundwater level, streamflow, and specific conductivity data taken at key locations in the basin to improve model parameters and predic tions. The fully integrated ParFlow.CLM model presented here represents a baseline scenario that provides a strong basis to quantitatively predict the impacts of major changes in landuse (e.g. wetlands and forests cleared for agricultural or urban develo pment), water use and climate on stores, fluxes, flowpaths and travel times of water in the SFRB. These predictions are essential to inform holistic land and water resource planning in the region in order to provide reliable water supply for human uses as well environmental flows that are protective of aquatic ecosystems. Furthermore accurate predictions of surface versus groundwater contributions to streamflow will allow prediction of transport and transformations of ecologically relevant solutes such as carbon, nitrogen and phosphorus in the spring and river systems in the region

PAGE 52

52 Table 2 1. Summary of ParFlow input variables and their relevant source of information. Table 2 2. Subsurface flow contributions as percent of total streamflow at river loc ations downstream of Cody Es carp ment (CE).

PAGE 53

53 Table 2 3. Annual total of major water balance components for the SFRB

PAGE 54

54 Figure 2 1 Map of the model domain Also shown are the extent s of the Santa Fe River Basin, locations of the river channels monito ring wells USGS stream gage s and major hydrogeological regions in basin

PAGE 55

55 Figure 2 2. 3D mesh used to define the domain A) the extent of the three major hydrogeological regions based on Floridan Aquifer Vulnerability Assessment (FAVA) dataset within the mesh along with their indicator variables Surficial Aquifer System blue region Intermediate Aquifer System green region and Upper Floridan Aquifer System red region. B) Locations of high hydraulic conductivity zones (blue cells) and ma in channels (red cells). ( a ) 150 m A B

PAGE 56

56 A B Figure 2 3. Spatial distribution of land surface properties over the SFRB. A) E levation B) M ajo r landuse classes

PAGE 57

57 Figure 2 4. C omparison of observed and model predicted daily streamflow measurements at various stream gage locations across the domain

PAGE 58

58 Figure 2 5 Comparison of observed and model simulated groundwater surface elevation at various locations across the domain Wells 1 through 6 are in the Upper Confined Region ( CR ) and wells 7 through 9 are in the Lower Unconfined Region ( UR )

PAGE 59

59 A B C D E F Figure 2 6 Potentiometric surface of UFAS with respect to mean sea level. A ) pre 2004 hurricane B) peak during 2004 hurricanes C) post 2004 hurricanes D) duri ng dry period in 2006 E) during end of extended dry period of 2006 2007 F) peak during tropical storm Fay. Blank space on northeast corner of domain indicates that the UFAS is located at greater depth than the model domain

PAGE 60

60 Figure 2 7 Time series of daily total subsurface flow contributions received by all the upstream river reaches from a given gage location. Note that time series for stations downstream to usgs1898 show only the daily subsurface flow contributions received downstream of usgs1898.

PAGE 61

61 Figure 2 8 Comparison of EMMA based versus ParFlow .CLM simulated ratio of subsurface flow contributions to total streamflow at three UR river locations.

PAGE 62

62 Figure 2 8 Continued.

PAGE 63

63 Figure 2 9. C omparison of measured and ParFlow .CLM predic te d streamflow for calendar years 2002 2003 (left column) and 2007 2008 (right column). Also shown are EMMA estimated (yellow circles) and ParFlow .CLM predicted (light blue line) daily total subsurface flow contributions at three UR river locations

PAGE 64

64 Fig ure 2 10. Monthly averaged water budget components for the study domain over entire simulation time period

PAGE 65

65 Figure 2 11. Spatial distribution of monthly averaged Evapotranspiration ( ET; mm) across the study domain over entire simulation period.

PAGE 66

66 Figure 2 12. Spatial distribution of average monthly subsurface recharge (mm) across the study domain over the entire simulation period.

PAGE 67

67 A B Figure 2 13. Controls on average daily ET (mm) and water table depth (m) dynamics in the SFRB. A) G eolog y B) L anduse.

PAGE 68

68 CHAPTER 3 GLOBAL SENSITIVITY A NALYSIS OF AN INTEGR ATED HYDROLOGIC MODE L INSIGHT ON GEOLOGIC AND VEGETATIVE CONTR OLS OVER THE HYDROLO GIC FUNCTIO NING OF A LARGE COMP LEX BASIN 3.1 Background Increasing human impact on water resources has led to numerous studies on hydrologic functioning of large basin s and more importantly on hydrologic system response to external stressor s such as groundwater extracti on, natural or synthetic contaminants, changes in climate and landuse conditions (Ferguson and Maxwell 2010 Goderniaux et al. 2009 Bolger et al 2011 ; Du et al. 2012) Establishing an accurate cause and effect relationship between any external stressor and corresponding hydrologic response in a basin requires an accurate accounting of dominant hydrologic processes as well as their spatiotempor al variations and interactions within the basin H ydrologic models are commonly used to test and improve our understating about dominant hydrologic processes within small to very large complex basins (Jones et al. 2008 ; Kollet and Maxwell 2008a ; VanderKwaak and Loague 2001) Recent advances in hydrologic modeling have established significant feedback between surface water, groundwater and land surface processes (Ferguson and Maxwell 2010 ; Goderniaux et al. 2009 ; Kollet and Maxwell 2008a ; Kollet and Maxwell 2006) and thereby the need for integrated representatio n of various surface and subsurface flow processes that coexist within hydrologic systems. This integration is achieved by using highly complex and fully coupled hydrologic models (e.g. HydroGeoSphere, lnhm, ParFlow CLM ) which have high computational requi rements. Moreover, the large number of inputs and complex mathematical algorithms within these fully integrated

PAGE 69

69 models makes it difficult to i dentify parameters controlling dominant hydrologic processes simulated by the model and their interactions over sp ace and time S ensitivity analysis (SA) is performed in hydrologic modeling studies to (1) gain insight on the working of complex models and identify dominant hydrologic processes within the study domain (e.g. Nossent and Bauwens 2012) ; (2) identify parameters that have most influence over dominant hydrologic processes identified in step 1 (e.g. Francone et al. 2012) ; (3) help prioritiz e factors to be included (by process of fac tor fixing) in the parameter estimation or model optimization process so that computational resources are not exhausted on insensitive parameters (e.g. Linhoss et al. 2012) ; and (4) identify parameters that interact with other parameters in the coupled non linear system (e.g. Linhoss et al. 2012) S tudies introducing new SA technique s or applying existing techniques to range of hydrological problems ar e extensively documented (e.g. Linhoss et al. 2012 ; Muoz Carpena et al. 2007 ; Saltelli et al. 2004 ; Van Griensven et al. 2006) The need for cost effective SA and calibration of complex distributed models is also well documented (e.g. Foglia et al. 2009) Despite the critical role that SA plays in the process of hydrologic model development and application, formal SA is often av oided in studies involving large scale, fully integrated hydrological models primarily because of their e xtremely high computation cost Instead, typically manual calibration and occasionally SA are performed on these models. F or example, Huntington and Niswonger ( 2012) manually calibrated the GSFLOW model for a 54 km 2 domain covering three snowmelt dominated watersheds in the E astern Sierra Nevada. The model was then used to study the relationship between snowmelt timing and hydrologic

PAGE 70

70 processes such as streamflow and groundwater recharge, and evapotranspiration losses. The modelers also evaluated the effect of variations in hydraulic conductivity on spatial distribution of groundwa ter head values. Jones et al. ( 2008) manually calibrated a steady state InHM model developed for the 75 km 2 Laurel Creek Watershed in Southern Ontario, Canada. Discrete groundwater heads at 50 observation wells were us ed to fit model simulated groundwater heads to field observations. However, no SA was performed on the model. Rihani et al. ( 2010) developed a P arFlow. CLM application for an idealized hillslope (5000 m x 100m x 80 310m). They executed 14 different model setup s to study the effect of land cover, atmospheric conditions, and subsurface formation and terrain properties over major land surface and subsurface processes. Goderniaux et al. (2009) calibrated HydrogeoSphere m odel developed for 465 km 2 Geer basin, Belgium. No SA was performed on the model which was used to study the impact of climate change scenario o n groundwater reserves. SA methods can be classified into (1) local or OAT in which the value of one of the par ameter s is changed within the allowable parameter range (while keeping other parameters fixed) and the resulting change in model output is observed; (2) GSA methods ( e g method of SOBOL ( Sobol, 1993 ) ; extended Fouri er A mplitude Sensitivity Test (FAST; Saltelli 1999 ) ) in which all the parameter values are changed simultaneously within the allowable parameter range and the resu lting change in model output is observed. Both the methods have their merits and shortcomings. The OAT method is easy to implement however the sensitivity measure is for a single location within the entire parameter space and therefore does not account for parameter interactions which may lead to misle ading results in highly nonlinear systems with

PAGE 71

71 interacting parameters. The GSA methods randomly or systematically s ample the entire parameter space and generate a large number of model input parameter combinations to overcome the abovementioned issues with OAT methods However the large number of model runs needed in GSA methods can make them computationally prohibitive when evaluating fully integrated model for large basin. Although providing qualitative results, global screening methods (e.g. Morris met hod Morris 1991 or Latin Hypercube One factor At a Time (LH OAT) techniques Van Griensven et al. 2006 ), are a more feasible option than quantitative GSA methods wh en deali ng with computationally expensive models because of the considerably small er number of model executions required to effectively identify the most sensitive parameters in the model (Yang 2011) For example, Nossent and Bauwens (2012) conducted an extensive LH OAT evaluation on a SWAT model developed for the 580 km 2 Klein Nete catchment. They used LH OAT to better understand SWAT model p arameters effecting flow and water quality processes. Francone et al. (2012) presented the first application of the Morris method on the land surface model UTOPIA to investigate the influence of model input parameters such as leaf area index, maximum vegetation cover, and rooting depth on selected energy and hydrological budget components f or a vineyard located in Cocconato, N orthern Italy. Many studies have also demonstrated the use of global screening method s in a two step GSA and uncertainty analysis (UA) of hydrological and water quality models (Muoz Carpena et al. 2007 ; Yang et al. 2012 ; Linhoss et al. 2012) T he t wo steps used in the se stud ies involved (1) use of a screening technique such as Morris method (e.g. Muoz Carpena et al. 2007 ; Yang et al. 2012) to identify the most important parameters and reduce the

PAGE 72

72 dimensionality of the parameter space by the process of factor fixing; and (2) further investigation of the shortlisted paramete rs using a quantitative and more exhaustive GSA approach such as a variance based technique (e.g. Muoz Carpena et al. 2007 used extended the FAST ). In this paper we apply the GSA screening technique developed by Morris (Morris, 1991 ) to the integrated hydrologic model ParFlow .CLM (Ashby and Falgout 1996 ; Kollet and Maxwell 2006 ; Maxwell and Miller 2005) previously developed f or the Santa Fe R iver Basin (SFRB) N orth Central Florida USA (Srivastava et al., 2013 a ) to investigate the sensitivity of predictions to uncertainties in a wide range of parameters governing coupled surface, subsurface and land atmosphere processes in the basin. The Morris method was selected over other GSA techniques because of the extremely high computation requirements of the ParFlow.CLM application built for the SFRB (discussed later in this paper). The previously developed application of ParFlow.CLM for the SFRB (hereafter referred to as the baseline model) accurately reproduced regional and seasonal patterns of evapotranspiration, streamflow and groundwater elevation, as well as surface groundwater mixing ratios in the UR of the basin. However, it consistently underpredicted streamflow recession and baseflow at unconfined river locations. The specific objectives of this study were to (1) investigate the sensitivity of hydrologic processes throughout the basin to vegetation, land surface, soil and ge ologic parameter uncertainty during both flood and base flow conditions ; (2) identify the most sensitive and interactive parameters governing s urface and groundwater flow particularly for streamflow recession and groundwater surface water interactions in the UR of the

PAGE 73

73 basin; (3) identify ranges of sensitive parameter values for the subset of models that performed satisfactorily throughout the basin over the study period; and (4) make recommendations for the design of more computationally intensive spatiall y distributed parameter estimation studies and field measurement campaigns to improve model predictions 3.2 Methods 3.2.1 Hydro G eologic C haracteristics of the Santa Fe River Basin The SFRB is a mixed use basin with a drainage area of about 3700 km 2 covering nine counties in N orth C entral Florida. The basin spans three distinct yet linked hydrogeologic regions (Figu re 3 1; Schneider et al. 2008 ) known as the confined region (CR), semi confined region (SCR) and the unconfine d region (UR). The three regions are name d based on the degree of confinement of the Upper Floridan Aquifer System (UFAS) a regional aquifer system comprised mostly of limestone, that underlies the entire study area ( see red region Figure 3 2 A ). The degre e of confinement exerts important control over both the regional surface and groundwater system flow dynamics. Table 3 1 summarizes the physical attributes of the basin relevant to this study. In the CR the UFAS is overlain with a thick clay rich confini ng unit (the Intermediate Aquifer System, IAS, depicted as the green region in Figure 3 2 A ), and the Surficial Aquifer System (SAS) which consists of unconfined, saturated sands (Schneider et al., 2008; depicted as the blue region in Figure 3 2 A ). Poorly d rained conditions in the CR result in shallow water table depths and well developed network of surface channels and wetlands. Annual average recharge from the SAS to the confined UFAS is estimated to be less than 0.3 m/year (Schneider et al., 2008).

PAGE 74

74 In th e UR, where the confining unit on top of the UFAS is completely eroded, limestone rocks are overlain with thin layer of sand which results in direct rapid recharge of the UFAS by infiltrating rainwater (Arthur, 2005). There are virtually no surface stream s feeding the main river in this region, however the aquifer is known to have major conduits which creat e an internal subsurface drainage network connected to the river (Meyer et al. 2008 ; Upchur ch 2007) Unlike teleogenetic karst systems with low intergranular porosity the eogentic karst features in Florida have retained high matrix porosity (Vacher and Mylroie 2002 ; Florea and Vacher 2007) As such primary matrix porosity in the SFRB plays an important role in water storage, flow through the subsurface environment, and exchange with the river system The region that marks the transition from confined to unconfined conditions (known as the Cody Escarpment, CE) is characterized by the presence of numerous sinkholes, sinking streams, siphons, springs and other karst features which have the potential to rapidly recharge the UFA S (Upchurch, 2007) During baseflow conditions the upper Sant a Fe River (SFR) is captured in its entirety by the SFR Sink in the CE region. The river emerges from the SFR Rise, located approximately 6 km downstream of the sink, as a first order magnitude spring. Hydrological behavior of the river in the transition z one is highly dynamic with no flow entering the sink during extreme drought conditions and with a portion of the river occasionally bypassing the sink as overland flow to the lower SFR during large flood events. Note that in some regions of the SFRB the U FAS is mapped as semi confined (SCR) because of a leaky or discontinuous IAS (Schneider et al., 2008; see area enclosed by green lines in Figure 3 1). The recharge

PAGE 75

75 to the UFAS is highly variable in these regions and is focused mainly around sinkholes and s wallets 3.2.2 Integrated 3D M odel ParFlow.CLM In this study we used the fully integrated surface water groundwater land surface process model ParFlow.CLM (Maxw ell and Miller 2005 ; Kollet and Maxwell 2006 ; Ashby and Falgout, 1996 ; Jones and Woodward 2001) to simulate hydrologic conditions within the SFRB. Here we present a brief summary of equations used by the model that are relevant to this study. Full deta ils on model physics and how the surface and subsurface components are integrated can be found in the references cited above Variably saturated groundwater flow is simulated by ParFlow, which solves t he 3D Richards equation using a cell centered finite d ifference scheme in space and an implicit Euler scheme in time (Eq 3 1) Surface flow simulations and routing are performed by solving the 2D Kinematic wave equation (Eq 3 2) and Manning s equation (Eq 3 3) Surface and subsurface equations are integrated using a free surface overland flow boundary condition. (3 1 ) (3 2 ) (3 3 ) p is the subsurface pressure head [L]; S w p ) is the degree of saturatio n [ p as defined by van Genuchten model (Van Genuchten 1980) ; S s is the specific storage [L 1 ]; is the porosity [ ]; k(x) is saturated hydraulic

PAGE 76

76 conductivity [LT 1 ]; k r is the relative permeability as a function of subsurface pressure head[ ]; z is depth below the land surface [L]; q s is a general source/sink term [T 1 ]; q r represents rainfall and evapo rat ive fluxes [LT 1 ]; || A,B || indicates the greater of A an d B; S f,x and S f,y are the friction slopes in x and y direction [ ]; n is the Manning s coefficient [TL 1/3 ]; an v is the depth averaged velocity vector of surface runoff [LT 1 ] ; o in Eqs. (3 2) and (3 3) is the pressure head or water ponding depth at the land surface [L] The coupled land surface energy and water balance is simulated by a modified version of the Common Land Model (CLM ( Dai et al. 2003 ) ) that is incorporated in ParFlow. CLM and it requires hourly atmospheric variables such as precipitation, air temperature, pressure, wi nd speed, specific humidity, and solar radiation along with soil moisture estimates (which are provided by ParFl ow) to simulate near land surface atmospheric fluxes such as evapotranspiration loss (Kollet and Maxwell 2008a) Details of equations used in CLM to simulate all the water and energy fluxes are provided elsewhere (Dai et al. 2003 ; Kollet et al. 2009 ; Maxwell and Miller 2005) E quations relevant to work presented here are summarized in the following. Energy balance equation for any given po int in a basin is given by (3 4 ) w here R n is net radiation at the land surface (W/m 2 ); LE is latent heat (W/m 2 ), H is sensible heat flux (W/m 2 ), and G is the ground heat flux (W/m 2 ). The water balance in the subsurface environ ment and the energy balance at the land surface are coupled together via the soil moisture content at or close to ground surface Each of the terms on the right hand side of Eq. (3 4) are independently calculated by CLM usi ng

PAGE 77

77 atmospheric forcing soil, vegetation characteristics and soil moisture provided by ParFlow LE is calculated as: (3 5 ) (3 6 ) where, L v is lat ent heat of evaporation (J/kg) and E is the sum of evaporative flux es from the foliage ( E c (kg/m 2 s) if vegetation is present) and ground E g (kg/m 2 s) Evaporative fluxes from the ground are calculated as: (3 7) Where a is the intrinsic density of air (kg/m 3 ); q g is the air specific humidity at the ground surface (kg/kg); q a is the air specific humidity at reference height z q obtained from atmospheric forcing (kg/kg); r d is the aerodynamics resistance at z q (s/m) Evap orative fluxes from the foliage are estimated as the sum of evaporative flux from wet foliage E w (kg/m 2 s) and total transpiration flux E tr (kg/m 2 s ) from the root zone cells: ( 3 8 ) ( 3 9 ) ( 3 1 0 ) Where f root,j* is effective root fraction in root zone layer j (3 11) (3 12)

PAGE 78

78 (3 13 ) Where f is the vegetation fraction ( ); L SAI is the stem plus leaf area index ( is the step function (one for positive and zero for zero and negative arguments); is the potential evaporation from wet foliage (kg/m 2 s) ; is the we tted fraction of the canopy ( ); L d is the dry fraction of foliage surface ( ); r b is the conductance of heat and vapor flux from leaves (s/m); and r s is the stomatal resistance (s/m). i fc wp respectively, are soil moisture conte nt in layer i, moisture content at field capacity, and permanent wilting point; f root,j is root fraction parameter which varies exponentially with the depth from land surface z ( ) in root zone layer j ; a and b are coefficients dependent upon vegetation ty pe. In ParFlow.CLM the term q s in subsurface mass balance equation (Eq. (1)) is expressed as (3 14) Where q g (s 1 ) is the flux of water infiltrating at the land surface after accounting for precipitation, canopy throughfall, and/or surface runoff 3.2.3 Baseline M odel The baseline scenario for the GSA was established using the ParFlow.CLM model for a n area that covers and extends beyond the topographic boundaries of the SFRB to include major springsheds that are known to contribute subsurface flow to the river in the UR of the SFRB (Fig ure 3 1; Srivastava et al., 2013 a ). From here on any reference to SFRB implies the extended boundary study domain shown in Figure 3 1. The land surface elevation in the model domain ranges from about 4 m to 74 m above mean sea level, resulting in a maxim um domain depth of 144m. The bottom of the domain was set at 70 m below mean sea level to include the top of the active

PAGE 79

79 subsurface zone believed to contribute flow to the SFR. A lateral discretization ( resulting in 56x50x144 cells, in x, y, and z dimensions, respectively. Thus the domain included a total of 403,200 rectilinear elements covering a total land sur face area of 6300 km 2 (Figure 3 1). Table 3 2 summarize s the spatial input data used to develop ParFlow. CLM for the SFRB Conduit locations were taken from Meyer et al (2008) and were assumed to be 10m in height and to occur between depths of 5 m to 15 m below mean sea level based on previous field investigations (Langston et al. 2012) Since ParFlow does not contain algorithms for turbulent conduit flow, conduits were represented as high hydraulic conductivi ty zones in this study ( Ando et al. 2003 ; Tsang et al. 1996 ) This assumption should be reasonable given the goal of the study is to assess the influence of these high hydraul ic conductivity zones on seasonal hydrologic budgets and large scale streamflow and groundwater interactions in the basin. The subsurface environment was assumed to be isotropic and homogeneous within each aquifer type (SAS, IAS, UFAS, conduits). The model was parameterized using literature values or values estimated based on previous groundwater flow modeling in the basin (see references in Tables 3 3 and 3 4). In the absence of study area specific information, default values of CLM related vegetation para meters were used in the model. Table 3 boundary conditions. To minimize the influence of initial condition assumptions the baseline model was spun up by repeatedly simulating 9 years (ye ars 2000 through 2008) of hydrologi cal conditions within the basin. S patially distributed hourly North American Land Data

PAGE 80

80 Assimilation System ( NLDAS ) forcing data was used to drive the model until a quasi dynamic equilibrium was reached, i.e. until the dif ference between the beginning and ending subsurface water storage (over the 9 year period run) dropped to less than <0.05 % of average subsurface storage and the simulated streamflow and groundwater elevations time series at multiple locations across the d omain showed minimal visual difference betwe en two consecutive 9 year runs. These criteria were met by the end of the second 9 year simulation period, thus all the baseline model results presented in this paper are the model results obtained from round 3 o f model spin up 3.2.4 Morris S creening M ethod Morris (Morris 1991) dimensional p level grid (Campolongo et al. 2007) ; k being the number of independent parameters included in the analysis (i.e. the model requires an input parameter vector i selected to span the specified value range for each parameter B etween consecutive model executions in the GSA only one parameter value is changed and the corresponding change in mo del output is observed. This individual change in model output (y) in response to a change in the i th parameter value in is known as the elementary effect (EE) of the i th input on output y and is calculated as: ( 3 15) is a value in {1/(p 1/(p 1)}. For this study p = 8. A f inite sample of EE i is obtained by randomly sampling different end

PAGE 81

81 EE i samples are cal culated A larg e value of indicates an important influence of the corresponding input on model output and a large value indicates that the corresponding input exert s a non linear effect on model output or is involved in interactions with other model pa rameters (Saltelli et al. 2004) Campolongo et al. (2007) proposed an improvement in the estimate of index of Morris by averaging the absolute values of the elementary effects (*) to eliminate the effect of op posite signs in non monotonic models. Also, the new index * was shown to be an acceptable substitute for variance based total index (Campolongo et al. 2007) In this paper we have presented all the results using * which we refer to in the following as the Morris mean 3.2.5 GSA E xperiment D esign A total of 19 ParFlow and 14 CLM related parameters were considered in the GSA. The number of simulations required in Morris Analysis is given by: N = r(k+1) ( 3 16) Where r is the sampling size for search trajectory and r = 10 used in this study is known to produce satisfactory results (Carpena et al ., 2010) and k is the number of parameter s (33 in this study). This resulted in a total of 340 model runs (see Tab les 3 3 and 3 4). All ParFlow parameters for each geological region represented in the baseline model were included in the GSA. However due to the large computational requirement to execute a single ParFlow run (discussed in the end of this section ) and the fact that CLM parameters varie s with the 9 vegetation types occurring in the region a re duced set of CLM vegetation parameters was obtained by (1) fixing the values of insensitive parameters (to their respective default values) based on results from a OAT SA performed on all CLM parameters for ParFlow .CLM simulations on a single

PAGE 82

82 homogeneous soil column; and (2) limiting the analysis to only the vegetation parameters corresponding to major land cover types found in the basin (i.e. evergreen needleleaf forests (31%), grasslands (14%), open shrublands (13%), upland mixed forests (11%) isolated wetland and floodplain mixed forests wetlands and water ( 13 %)). All the CLM input parameters related to energy budget calculations were fixed to their default val ues. The Morris method based sampling scheme in Simlab 2.2 software was used to generate 340 different combinations of parameter values using ranges summarized in Tables 3 3 and 3 4. The r ange s of ParFlow parameter values used were based on the published literature or estimated from previous modeling studies in the SFRB. In absence of reliable basin specific information, CLM parameters were generated using a 50% variation in their default values. The pressure head distribution from the last time step of the baseline model was used to initialize each of the 340 models which were run for a four year simulation period, 1/1/2000 through 12/31/2003 T o minimize the transient effect of change in parameter value s from the baseline model values, we analyzed the model simulated outputs for only the last 15 months out of the total simulation period This period includes a significant storm event during 2/1/2003 5/30/2003 well defined pre storm low baseflow conditions during 10/1/2002 1/31/2003 and a variable po st storm period (6/1/2003 12/31/2003). This period excludes the e xtreme hurricane, tropical storm, and drought events that occurred during years 2004 through 2008 in the baseline model to eliminate any possible bias in the analysis produced by unusual hy drologic conditions during these extreme events. The output variables assessed in the sensitivity analysis

PAGE 83

83 are summarized in Table 3 6 Streamf low was assessed at 6 locations, and UFAS groundwater elevation s were assessed at 10, where observations were ava ilable (Figure 3 1). Five out of the ten output variables (peak flow, flow ratios 1 3, and maximum streamflow loss to the subsurface) were studied for the storm event that occurred during 10/1/2002 5/31/2003 (Figure 3 3). Groundwater flow contributions to streamflow were assessed only at three unconfined river locations. ET estimates were assessed for the two major contrasting geological regions 3.2.6 Measures of Goodness of Model Fit In addition to calculating the elementary effect of each parameter for eac h selected model output, we calculated goodness of fit measures for streamflow and groundwater level predictions at each of the observation locations for each of the 340 model runs. We used the Nash Sutcliffe model efficiency coefficient (NSC, Eq 3 1 7 ; Nash and Sutcliffe 1970 ) for streamflow and log streamflow to evaluate goodness of fit for high streamflows and low streamflows, respectively. For groundwater elevations we used the coefficient of determination (R 2 Eq 3 1 8 ; Legates and McCabe 1999 ) and percent bias (PBIAS, Eq. 3 1 9 ; Gupta et al. 1999 ) to evaluate goodness of fit. (3 1 7 ) (3 1 8 )

PAGE 84

84 ( 3 1 9 ) where , ,and ,respectively, are the i th observation, the i th simulated value, mean of observed data, and the m ean of simulated data; n is the Beven and Binley 1992 ; Beven and Freer 2001 ) for both str eamflow and log streamflow. Threshold values of 0.7 for R 2 and 0.35 for PBIAS were The model execution time depended on the difficultly level posed by combination of parameter val ues and for this study ranged from 63 hours to 201 hours for each model run (about 260 models finished within 100 hrs) 3.3 R e sults and D iscussion 3.3.1 Baseline M odel P erformance Detailed analyses of baseline model simulated streamflow, surface groundwater interact ions and groundwater elevation at multiple locations, along with analysis of major water budget components acr oss the basin, were presented by Srivastava et al. (2013 a ) Here we briefly present comparison of field observations and model simulated streamflo w at two locations (1500 and 2500), and groundwater elevations at four locations (E2, E3, E8, E9) to highlight the processes that are, and are not, well predicted by the baseline model. Understanding the sensitivity of these processes to uncertain paramet er values is an important first step to improving the predictive capability of the model.

PAGE 85

85 In general model performance at all the stations within a given hydro geological region (i.e. CR and UR) was comparable. Stations 1000, 1500 and 1898 are located in the CR of the basin where streamflow is rainfall dominated, with zero or insignificant measured flow in the absence of rain (see mean and coefficient of variation (CV) values for CR streamflow stations in Table 3 1 and the time series of streamflow at 1500 in Figure 3 3). In contrast Stations 1975, 2500 and 2800 located in the UR of the basin, receive substantial streamflow contributions from groundwater throughout the year (see mean and CV values for UR streamflow stations in Table 3 1 and the time series of streamflow at 2500 in Figure 3 3). For brevity we focus our streamflow analyses on stations 1500 and 2500, as representative of behaviors in the CR and UR, respectively. In general, m odel simulated streamflow at USGS 1500 in the CR showed good agreement with measured data during both low and high flow conditions ( Figure 3 3 ). Timing and peaks of hydrographs were accurately simulated by the model during all storm events and insignificant groundwater contributions resulted in zero streamflow in absence of rain. At station 2500 in the UR the model simulated streamflow showed reasonable agreement with measured data during high flow conditions T he model adequately simulated the timing and peaks of storm events but consistently underpredicted flow during stor m recession periods as well as during low flow conditions. These observations indicate that the baseline model is missing a near stream transient s ource of groundwater in the UR that is released after storm events, as well as underestimating the steady gro undwater baseflow. Srivastava et al (2013 a ) hypothesized that these problems may be attributed to inaccuracies in the extent and

PAGE 86

86 thickness of the IAS particularly in areas where the CR transitions to the UR, which could result in insufficient water infi ltrating into the UFAS in the these regions. Figure 3 4 compares model simulated and measured groundwater elevation at four UFAS locations; E1 and E3 in the UR, E8 in the SCR and E10 in the UR. These wells span the various behaviors observed by Srivastav a et al (2013 a ) throughout the SFRB. Wells E2 and E10 simulate both the mean groundwater level and the transient groundwater response accurately. Well E3 over estimates both the mean groundwater level and the magnitude of the transient response, and Well E8 underestimates both the mean groundwater level and the magnitude of the transient response. Srivastava et al. (2013 a ) attributed the d ifference s between the model simulated and measured elevations at these location to uncertainty in the extent of and hy draulic conductivity of the IAS a s well as uncertainty in conduit locations and properties in the vicinity of the affected wells. The goal of this GSA is test the Srivastava et al. (2013 a ) hypotheses and to provide insight on the relative influences of al l ParFlow.CLM parameters on evapotranspiration, streamflow predictions, groundwater levels, and groundwater surface water interactions throughout the basin 3.3.2 Global S ensitivity A nalysis In the following sections Morris analys e s are presented both in t abu lar and graphical format s. A table is provided for each output of interest which summarizes the ranking of all parameters based on their Morris mean *, and standard deviation ,values Parameters with the highest values (i.e. most sensitive or most in teractive) are assigned a rank of one. Non influential parameters (i.e. parameters with zero values for * and are represented by a blank cell ) In addition scatter plots of the actual values

PAGE 87

87 gnitudes of the sensitivities and interactions across rank 3.3.2.1 Evapotranspiration Evapotranspiration accounts for 75 80% of the water budget in the SFRB ( Srivastava et al. 2013 a ; Schneider et al 2008) and is thus an extremely important component of the w ater budget. Table 3 7 presents the parameter rankings based on weighted mean daily evapotranspiration occurring in the CR and UR of the basin over the study period. In the CR the hydraulic conductivity of the intermediate aquifer sy stem (ias_k) was found to be the most sensitive and most interactive parameter influencing ET, through its control on the surficial water table depth in the CR. Vegetation properties of mixed wetland and floodplain forests (lai_mf) and evergreen needlelea f forests (lai_enf) were the next most sensitive parameters for ET, with lai_mf also showing a high interactive effect. The location of mixed wetland and floodplain forests in higher water table regions makes them more influential in determining average d aily ET over the CR than the evergreen needleleaf forests, even though the evergreen needle forests occupy a greater land area (Table 3 1). These findings are consistent with previous results that show a strong dependency between ET losses and water table depth for shallow water table regions (Srivastava et al. 2013 a; Kollet and Maxwell 2008a) Table 3 7 indicates that surficial aquifer hydraulic conductivity (sas_k) is also sensitive an d interactive, but Figure 3 7 shows that ias_k and lai_mf are by far the most important parameters affecting ET in the CR. In contrast, Table 3 7 shows that wilting point (wp) is the most sensitive and interactive parameter affecting ET in the UR, followe d by the hydraulic conductivity of the upper Floridan aquifer system (fas_k), and the van Genuchten n parameter (vg_n).

PAGE 88

88 ET in the UR also shows some sensitivity to maximum LAI of mixed forest (lai_mf) and open shrubland (lai_os), but Figure 3 7 shows tha t wp is by far the most sensitive and interactive parameter. Thus in both the CR and the UR, ET is more sensitive to soil/geologic properties than vegetative properties. However in the energy limited CR, where water tables are shallow, ET is most sensitiv ity to properties of the confining layer that control the height of the surficial water table. In the moisture limited UR, where water tables are deep, ET is more sensitive to surficial soil/geologic parameters controlling unsaturated moisture retention 3.3.2.2 Total s treamflow Table s 3 8 and 3 9 summarize the parameter rankings based on Morris mean and s igma for stream flow characteristics at each of the streamgage locations. Analyzing multiple characteristics of daily streamflow time series at multiple locatio ns within the basin is useful to identify the spatial and temporal variability in dominant hydrological processes affecting streamflow across the basin. In general total streamflow at all locations across the domain w as significantly influenced by the sam e parameters affecting ET losses with some differences in the order of importance (Table 3 8 and Figures 3 5 A and 3 5 B ). This high influence of ET on total streamflow i s expected in the basin where ET is the biggest water budget component besides rainfall (Table 3 1) and show s the importance of an integrated modeling approach which accounts for important bi directional feedbacks between land surface subsurface components of hydrologic cycle. ET exerts significant control over available surface and subsur face storage for rain fall, which in turn affects the amount of water contributing to streamflow via surface, shallow or deep subsurface flow paths.

PAGE 89

89 In the CR total streamflow showed most sensitivity to lai_mf, followed by ias_k, with ias_k showing by far the most interactive effects (Figure 3 5 A ). Streamflow also showed some sensitivity to lai_enf, but significantly less than sensitivity to lai_mf. Again this is likely due to the proximity of the wetland and floodplain mixed forests to the river channels in the CR and is in agreement with previous studies in the basin which hypothesize that that stream flow in the CR is generated primarily from direct rain or surface runoff received from near stream locations (Upchurc h 2007) In the UR total streamflow is most sensitive to fas_k followed by lai_mf and wp, indicating the primary control of the Floridan aquifer hydraulic conductivity properties on groundwater contributions to the stream and the secondary control of E T on total water available in the surface subsurface system. In general the influence of land cover based ET parameters on total streamflow is more pronounced in the CR where land cover related parameters ( such as lai _mf, lai _enf, z0_enf, and z0_mf ) showed higher rankings relative to their UR rankings. In contrast in the UR the influence of soil based ET parameters (wp and vg_n) gained relative importance in comparison to rankings in the CR. Again t his difference is attributed to the contrast in soil/ geolo gical characteristics which control the water available to plants for ET in the two regions 3.3.2.3 Peak s treamflow Peak streamflow rates at all confined and unconfined locations showed highest CR_mannings ) follow ed by ias_k (Figure 3 5 C ), with the exception of USGS2800 where the order of sensitivity of these parameters was reversed (Table 3 8). This indicates that streamflow peaks at unconfined locations generally originate from rainfall in the CR and are contro lled by the rainfall runoff response in that region. The strong influence of IAS aquifer is again due

PAGE 90

90 to its influence on water table location in the CR, and confirms that saturation excess, rather than Hortonian streamflow generation processes dominate i n this basin ET related parameters (lai_mf, lai_en, wp, and z0_enf) showed influence over the peak flow rate at all the locations (Table 3 8 ) by virtue o f their control over intra event surface and subsurface storage 3.3.2.4 Storm h ydrograph r ecession Flow r atios 1 3 (defined in Table 3 6) reflect the rate of peak flow recession for the February May 2003 storm event. In general, a slower decay rate implies a system with heterogeneous sources, travel paths and travel times for event water reaching the stream, often associated with a higher influence of subsurface contributions during recession of the storm hydrograph (Panagopoulos and Lambrakis 2006) Investigating the sensitivity of flow rati os 1 3 to land surfa ce, vegetative and geologic parameters provides insight into the surface and groundwater flow process es that play important role s over the duration of storm flow recession. In general, for all confined and unconfined locations, ratio 1 showed highest influe nce CR_mannings (with highest sensitivity and interactions ) emphasizing the dominance of CR streamflow generation and transport processes during early hydrograph recession even at unconfined locations (Table 3 9 and Figure 3 5 D ) For stations in the UR, Manning s coefficient in the UR channel (UR_mannings) ranked second closely followed by fas_k and and then ias_k. The sensitivity of ratio 1 to fas_k in the UR impl ies that even during the early period of hydrograph recession groundwater flow parameters are influential For ratio s 2 and 3 the relative ranking of the parameters change s, but the same four parameters retain the top four ranks. Hydraulic conductivity of the Floridan aquifer (fas_k) becomes

PAGE 91

91 the most sensitive parameter for ratio 3, confirming an increased influence of Floridan groundwater flow on late hydrograph recession as baseflow conditions are approached In the CR, stations 1000 and 1500 showed large scatter in ranks from ratio1 to ratio 3 with CR_mannings maintain ing the highest direct influence during entire hydrograph recession period. At station 1000 ET related parameter gained importance for ratio 3 (e.g. i ncreased ranking of lai_mf, vg _n, sai_mf, lai_g, and wp) indicating higher influence of ET losses on streamflow during later par ts of the recession period. Station 1500 showed increased influence of shallow subsurface flow contribution (in comparison to station 1000) with increased sensitivity to SAS properties ( r ank 2 based on sensitivity and r ank 1 based on interactions for sas_k ). Insignificant effect of upper Floridan aquifer system on streamflow throughout the region is evident from consistent low ranking of UFAS related parameters 3.3.2.5 Mean and c oefficient of v ariation of g roundwater e levation Table 3 10 and Figures 3 6 A and 3 6 B summarize the parameter rankings based on Morris mean and s igma for average daily potentiometric surface elevation at multiple Floridan well locations across the domain. Spatial patterns in potentiometric surface elevation across the upper Floridan aquife r control the direction of regional groundwater flow in the study domain and affect streamflow in the UR of the domain. In general, a ll CR wells (E1 E7) show high potentiometric head sensitivity to ias_k and fas_k, and some sensitivity to fas _ss. Wells E3 and E4 that are located in the vicinity of conduits also show high sensitivity to conduit hydraulic conductivity (c_k). The high influence of ias_k on all CR wells demonstrates the strong control the IAS exerts over recharge to the UFAS. Wells E9 and E1 0 in the UR showed dominant sensitivity to fas_k, with some sensitivity to wp which has a major influence on ET and thus recharge in the UR. Wells

PAGE 92

92 close to the boundary of UR (wells E2 E6) also showed high rank for wp indicating the influence of UR recharg e on the potentiometric surface in these locations Well E8 in the SCR showed highest sensitivity to the hydraulic conductivity of the semi confined region (SCR), followed by fas_k and wp. Table 3 11 and Figures 3 6C and 3 6D summarize the parameter ranki ngs based on Morris mean and s igma for the coefficient of variation (CV) of potentiometric surface at multiple well locations across the domain. In general, fas_ss, fas_k and ias_k were the most sensitive and interactive parameters for potentiometic surfac e CV at all wells, with fas_k being the most sensitive parameter in the UCR and fas_ss being the most sensitive in the CR. Wells E2 and E3 in the CR also showed some sensitivity to conduit head is more sensitive to the storage properties of the aquifer resulting from compressibility of the water and porous matrix (specific storage) than the storage properties of the aquifer that result from raising the water table (porosity), even in the UR At well E8 in the SCR geologic properties in that local region (sc_ss, sc_k, and sc_n) were the most influential 3.3.2.6 Surface w ater g roundwater i nteractions Exchange of water between river channels and the underlying UFAS was investigated to determine t he most influential parameters control ling groundwater exchange with the unconfined river locations. It has been previously established that the CR stream locations d oes not receive any groundwater flow contribution from the UFAS and no significant groundw ater contribution from the SAS (Upchurch 2007) Table 3 12 and Figures 3 5 E and 3 5 F summarize the rank of parameters most influential to (1) total baseflow contributions for the 15 month study period, and (2) maxi mum surface water lost to groundwater during the peak of the February May 2003 storm event. For

PAGE 93

93 contributions that all unconfined river cells upstream of a given location receive each day. For instance, total baseflow received by station 1975 was calculated as the total groundwater received by the river between stations 1898 and 1975 and total baseflow received by station 2500 was calculated as the total groundwater contrib utions received by all channel cells between 1898 and 2500. Total baseflow contribution received by all UR locations w as most significantly influenced by fas _k followed by c _k showing the importance of groundwater flow from the UFAS including high hydraul ic conductivity zones. The E T related parameter wp again shows high influence on the total groundwater contributions to streamflow at all three locations (rank 3 for stations 2500 and 2800 and rank 4 for station 1975). Maximum loss of surface water to the UFAS was significantly influenced by fas ss UR_mannings, and CR_mannings (Ranked 1 3) at all locations Highest sensitivity to fas ss shows the importance of subsurface storage in the near river environment. Higher specific storage allows more water to b e transferred to, and stored in, the aquifer without increasing the potentiometric head, and thus without reducing the head gradient between the river and the aquifer. S ensitivity to the two roughness coefficients shows the importance of peak river levels and thus again higher gradients between the river and aquifer 3.3.3 Goodness of F it of the GSA M odel R uns Hydrologic representation of large scale basins in integrated models is by necessity a simplification of reality in that homogeneous, isotropic conditio ns are commonly assumed within the major hydro geologic units in the model domain, an assumption which is violated in all natural systems (McDonnell et al. 2007) Although

PAGE 94

94 inconsequential during SA (because SA is o nly designed to identify the parameters that a specific output variable are most sensitive to) these simplifying assumptions often make it difficult to reproduce observed spatially variable hydrologic responses throughout a basin. In addition there is ofte variety of different combinations of sensitive parameter values may achieve equally reliable model predictions (Beven and Binley 1992) Furthermore combinations of parameters that reliably reproduce hydrologic response in one region of the basin may not produce reliable predictions elsewhere, particularly for large basins like the SFRB which span contrasting hydro geologic conditions. In these cases allowing spatially variably parameter values within major hydro geologic regions may be required to improve model performance. GSA is a useful screening tool to both determine the sensitive parameters that should be included in more comprehensive spatially variable parameter estimati on techniques, and to help refine the range of parameters values that more reliably reproduce observed behavior in various regions of the basin. In this study we examined the goodness of fit of the 340 GSA model runs at the 6 stream gage locations and 11 well locations where observations were available based these performance criteria, were identified for each observation location. These behavioral models were then scr eened to identify those parameter sets that were behavioral at all observation locations within each subregion and then to identify the parameter sets that were behavioral at all locations in the SFRB (Figure 3 8) These results show that 162 models were f ound to be behavioral at all confined streamflow

PAGE 95

95 behavioral at station 1975, and only 8 were behavioral at all groundwater locations. The short term, small scale inter actions between conduit and aquifer matrix that are known to occur between stations 1898 and 1975 are not simulated by the model. In contrast to field observations where river water is lost to surrounding aquifer and recovered back at downstream locations, our model route the upstream contribution (with some interactions with surrounding aquifer) directly to the downstream river locations (i.e. 1975) via river channel. This results in model consistently simulating higher than observed streamflow during peak s of storm events and therefore making it difficult to obtain behavioral models at stations 1975 and downstream unconfined locations at the same time. As such when all streamflow stations are considered together there are no parameter sets that are behavio ral at all six locations However, if station 1975 is removed from the analysis then 28 out of 340 models were found to be behavioral for streamflow in remaining stations in the confined and unconfined regions of the basin, and of these 6 were also behavio ral at the 10 well locations. The small number of models that were identified as behavioral throughout the basin underscores the need to use a parameter estimation technique that can account for spatial variability in the parameter values within hydro geol ogic zones to improve the predictive performance of the model. Figure 3 9 shows a comparison between measured and model simulated streamflow for the baseline model, the best fit model for the location, and the set of models found to be behavioral for the majority of the basin, for stations 1500 and 2500. At station 2500, total groundwater contributions received by channels upstream of station 2500 are also shown (dotted lines) for each of the selected models, along with

PAGE 96

96 discrete estimates of total gro undwater contributions (yellow circles) obtained in a previous study using EMMA with legacy specific conductivity data (Srivastava et al., 2013). For station 1500 both the best fit model and the set of 6 models that are behavioral throughout the domain r eproduce streamflow quite well, although all models slightly underestimate the peak streamflow. At station 2500 the best fit model matches the pre storm base flow conditions, the peak streamflow and the base flow recession quite well, however when compar ed with the EMMA results it underestimates the groundwater contribution to streamflow during the early storm recession. Thus, even with locally best fit parameters, the model is unable to account for the near stream transient s ource of groundwater in the UR that is released after storm events. Model performance deteriorates slightly at station 2500 for the set of 6 models that are behavioral everywhere in the domain indicating that achieving good model results everywhere the domain, with parameters that a re homogeneous in a limited number of hydrogeologic zones, is difficult. It is interesting to note that the 6 behavioral models show a more significant drop in groundwater contributions to streamflow during the storm event than the best fit model or the b aseline model. Figure 3 10 shows a comparison between measured and model simulated groundwater elevation at wells E2, E3, E8, and E10, for the baseline and 6 behavioral model runs. These results indicate that the more behavioral models do not improve the mean groundwater levels significantly, but do improve the behavior of the groundwater coefficient of variation, particularly at E3. Additional insights can be gained by looking into the values of the most sensitive parameters for the behavioral models. Fo r example at station 1500, 95% (199 out of 209) of models with ias_k values less then ~10 6 m/hr were found to be

PAGE 97

97 behavioral where as only 44% (57 out of 131) of models with ias_k values greater than ~10 6 m/hr were found to be behavioral. Using a more st showed that 77% (160 out of 209) of models with ias_k values less then ~10 6 m/hr were found to be behavioral whereas only 18% (23 out of 131) of the models with ias_k values greater ~10 6 m/hr were found to be behavioral. Thus a low hydraulic conductivity IAS is required to produce the shallow water table needed to produce sufficient streamflow in the CR. Further investigation reveals that for models with ias_k values less than ~10 6 m/hr all the models that were non behaviora l are the models that have high lai_mf, lai_enf, z0_enf values indicating too much evaporation from the domain to produce sufficient streamflow. Similar analysis for station 2500 indicated that all models with extremely low values of ias_k (~10 9 m/hr) wer e non behaviora models (total of 141 look at NSC values for flow and log flow values indicated that an ias_k value of ~10 9 m/hr produced reasonable NSC values for flo w simulations but performed badly for NSC of logflow values indicating that ias_k values needs to be relaxed in order to improve baseflow predictions in the UR. Based on the results obtained in this study an ias_k value in between ~ 10 6 10 8 m/hr should give behavioral models for both the confined and unconfined locations 3.4 Conclusions The Morris method of screening was successfully applied and provided a reduced set of most sensitive and interactive parameters to be further investigated with more quant itative parameter estimation and uncertainty analysis techniques such as spatially distributed Ensemble Kalman Filtering (Graham and McLaughlin 1989 ; Graham and McLaughlin 1991) In general seven parameters (out of 33 parameters

PAGE 98

98 included in the analysis), namely the hydraulic conductivity of the IAS, the hydraulic conductivity and specific storage of the UFAS, maximum LAI for mixed forests, wilting confined and unconfined regions were identified to be most sensitive towards ET, streamflow and groundwater level predictions. Of these seven parameters all but maximum LAI for mixed forests showed interactive effects, with the hydraulic conductivity of the IAS being the most interactive parameter. The sensitivity and interaction of the hydraulic conductivity of the IAS underscores the importance of better mapping of the lateral and vertical extent of IAS in confined, transition and Wacasassa Flat regions The significance of ET parameters such as lai_mf, and wp underscores the importance of ET and highlights the significance of feedback between land atmosphere, surface and subsurface processes throughout the basin. Vegetative properties of land covers, p articularly those close to the river showed significant influence over several aspects of surface and groundwater flow. Vegetation properties were found to be more influential in the CR where ET is more energy limited than water limited due to the shallow surficial water table In UR w here water tables are and ET is water limited soil moisture retention properties such as wilting point and parameters of the van Genuchten equation were more influential. Studying the range and median values of most sensiti ve parameters identified during Morris analysis (Figure 3 11), for base case, all the model runs in GSA, and the behavioral models for the entire basin, provided some useful information for future sensitivity and uncertainty analysis of the model. Behavior al models were found to have high parameter values ias_k, fas_k, CR_mannings, UR_mannings, fas_ss as compared

PAGE 99

99 to the base case. In contrast wp and lai_mf were found to have lower values in the behavioral models (Figure 3 11). Further investigation on how t variables responded to variations in most sensitive parameter values further enhanced our understanding about the most sensitive parameters and their effects on output variables (Figure 3 12). Figure 3 12 illustrates the variations in so me of the model simulated output analysis. Evapotranspiration losses from the CR were found to decrease with an increase in ias_k values; this behavior might be attributed to the increase in water table depth with increased ias_k (Figures 3 12 A ). There was a dual trend in total and peak streamflow response at station 1500 and 2500 (Figure 3 12 B ), with both peak and total streamflow show an increase with increased ias_k values upto ~10 6 m/hr and thereafter there was an observed fall in peak and total streamflow values at both the stations. It was interesting to note that the water table was not significantly affected by an increase in ias_k values upto ~10 6 m/hr; after which i t showed significant increase with an increase in ias_k values. In light of the above observations we can attribute an initial increase in flow characteristics at both the locations to more water reaching the channels from IAS cells that were closer to the channels. When the ias_k values were increased above ~10 6 m/hr more water might have recharged the underlying UFAS increasing the water table depths from the land surface; greater water table depths would mean that more rain water will be lost to the sub surface storage and less water will arrive at the stream.

PAGE 100

100 Both the surface roughness coefficients were found to be sensitive towards peak flow and flow ratios in the UR. However, as clearly demonstrated in Figures 3 12 C and D, an increase in the value of these coefficients will improve the storm recessions at the cost of storm peaks. Therefore, caution must be taken to ensure that the storm recessions simulated by the model are not for the wrong reason (e.g. increasing surface roughness to hold water in t he channels for longer time to get slow recession). fas_k on the other hand showed an increase in total streamflow simulated at station 2500 along with an increase in flow ratio (i.e. slow er storm recession; Figure 3 12E ). An increase in fas_ss was found to significantly affect the total groundwater contributions received at station 2500 (Figure 3 12 F ). As such simultaneously varying both fas_k and fas_ss parameter values can provide a subsurface storage release mechanism similar to that observed in betwee n conduit aquifer; where fas_ss can store water lost from the river channels during storm peaks and fas_k can facilitate quick return of the stored water once the storm has passed. Parameters wp and lai_mf, respectively, showed a direct relationship with ET losses within the UR and CR of the basin. Increased value of wp meant that less water can be lost via vegetation whereas increase in lai_mf showed an increase in ET losses. To summarize above discussion, ias_k and lai_mf, respectively, are needed to be increased and decreased in order to achieve less ET losses and more recharge of UFAS in the CR of the basin. Whereas, fas_k and fas_ss should be varied simultaneously in order to obtain reliable storm recession simulations by the model. Although qualitativ e, the result of the Morris GSA was successful in identifying spatial and temporal variability among the dominant hydrological processes and

PAGE 101

101 sensitive parameters for this large complex basin. Moreover, the findings demonstrate important nonlinear interacti ons among geologic, soil and vegetation properties on land atmosphere, surface and subsurface processes across large scales. The results presented in this paper helped in minimizing the number of parameters that require more accurate estimations. Moreover, based on our findings it can be concluded that any future sensitivity and uncertainty analysis on the model should account for spatial variability in the parameter values within the contrasting hydro geologic regions of the basin

PAGE 102

102 Table 3 1. Physical cha racteristics of the Santa Fe River Basin Topography a Mean elevation (m) 36.2 Minimum elevation (m) 3.7 Maximum elevation (m) 74.1 Mean slope (%) 0.2 Maximum slope (%) 1.0 Major land use b Evergreen needleleaf (%) 31 Open sh rubland (%) 13 Grassland (%) 14 Mixed forest uplands (%) 11 Mixed forest isolated wetlands and floodplains (%) 13 Dominant soil type c Confined region sand to clayey sand Unconfined region sandy soil Hydro logic uni ts d Confined region SAS,IAS,UFAS Unconfined region UFAS Water budget components mean annual rainfall (mm) e 1252 mean annual evapotranspiration (mm) f 1041 Streamflow USGS 1000 (3.4,3.5 ) g Streamflow USGS 1500 (7.6,2.5 ) Stream flow USGS 1898 (8.3,2.7 ) Streamflow USGS 1975 (6.5,1.8 ) Streamflow USGS 2500 (30.1,0.8 ) Streamflow USGS 2800 (44.9,0.6 ) a 1500 m resolution elevation dataset used in this study (USGS, 2010) b FWCC land use data for year 2003 c,d Schneider et al., 2008 e NLDAS rain data 2000 2008 f Schneider et al., 2008 g (mean, coefficient of variation) daily USGS streamflow records for water years 2000 2008 SAS = Surficial Aquifer System IAS = Intermediate Aquifer System UFAS = Upper Floridan Aquifer System

PAGE 103

103 Table 3 2. Input data used in ParFlow .CLM for the SFRB

PAGE 104

104 Table 3 3. ParFlow parameters included in GSA by Morris method

PAGE 105

105 Table 3 4. CLM parameter s included in GSA by Morris method.

PAGE 106

106 Table 3 5. Summary of relevant information on intial and boundary condition used in baseline ParFlow.CLM Aquifer Condition applied Source Initial condition within major aquifers in domain SAS 90% saturated everywhere IAS Hydrostatic equilibrium with the UF AS UFAS May 2002 Potentiometric surface map SRWMD Boundary condition over domain boundaries East boundary Constant pressure equal to the initial condition West boundary Constant pressure equal to the initial condition North boundary Const ant pressure equal to the initial condition South boundary Constant pressure equal to the initial condition Domain bottom No flux Domain top Overland flow SRWMD = Suwannee River Water Management District

PAGE 107

107 Table 3 6. ParFlow .CLM output vari ables assessed in the GSA

PAGE 108

108 Table 3 7. Morris mean and standard deviation based rank of individual parameters for mean daily ET Losses Parameter Confined Unconfined Confined Unconfined c_k 28 17 28 10 c_n 23 23 c ss 28 30 28 30 CR_mannings 23 30 22 31 fc 11 10 15 7 ias_k 1 11 1 5 ias_n 20 19 9 13 ias_ss 20 23 14 23 lai_enf 3 9 10 14 lai_g 11 13 23 15 lai_mf 2 5 2 6 lai_os 7 6 13 8 sai_enf 9 19 12 20 sai_g 15 12 19 16 sai_mf 10 23 17 23 sai_os 16 19 17 20 sas_k 4 13 3 11 sas_n 19 19 sas_ss 20 19 sc_k 18 19 23 20 sc_n 24 23 25 23 sc_ss 28 23 28 23 fas_k 16 2 5 2 fas_n 17 18 fas_ ss 24 13 16 11 UR_mannings 13 19 vg_a 14 4 4 4 vg_n 11 3 7 3 wp 7 1 6 1 z0_enf 5 8 11 17 z0_g 24 23 25 23 z0_mf 6 7 8 9 z0_os 24 23 25 23

PAGE 109

109 Table 3 8. Morris mean and sigma based rank of individual parameters for total and peak streamflow

PAGE 110

110 Table 3 9. Morris mean and sigma based rank of individual parameters for streamflow recession characterist ics

PAGE 111

111 Table 3 9. Continued.

PAGE 112

112 Table 3 10. Morris mean and standard deviation based rank of individual parameters for mean daily groundwater elevation

PAGE 113

113 Table 3 11. Morris mean and standard deviation based rank of individual parameters for coefficient of variation daily groundwater elevation

PAGE 114

114 Table 3 12. Morris mean and standard deviation based rank of individual parameters for total baseflow contribution received by three unconfined river location s

PAGE 115

115 Figure 3 1. Map of the model domain Also show n are the extent s of the Santa Fe River Basin, locations of the river channel s, monitoring wells, USGS stream gages and major hydrogeological regions in basin.

PAGE 116

116 Figure 3 2 3D mesh used to define the domain A) the extent of the three major hydrogeologi cal regions based on Floridan Aquifer Vulnerability Assessment ( FAVA ) dataset within the mesh along with their indicator variables Surficial Aquifer System blue region Intermediate Aquifer System green region and Upper Floridan Aquifer System red region. B) Locations of high hydraulic conductivity zones (blue cells) and main channels (red cells) A B

PAGE 117

117 Figure 3 3 C omparison of observed a nd model predicted daily streamflow measurements at various stream gage locations across the domain

PAGE 118

118 Figure 3 4 Comparison of observed and model simulated groundwater elevation at various locations across the domain Wells E2 and E3 are in the Upper Confined Region ( CR ) and well s E8 and E10 are in the Semi Confined Region (SCR) and the Lower Unconfined Region ( UR )

PAGE 119

119 A B C D E F Figure 3 5. Scatter plots effects of each parameter. Results are presented for the simulation of A) T otal flow at USGS 1500 B) T otal flow at USGS 2500 C) P eak flow at USGS 1500 and USGS 2500 D) F low ratio 1 3 for USG S 2500. E) T otal groundwater contribution received by USGS 2500 F) E xchange of groundwater during peak of 2003 storm at USGS 2500. Labels are provided for the most significant parameters.

PAGE 120

120 A B C D Figure 3 distribution function of the elementary effects of each parameter. Results are presented for the simulation of A) M ean groundwater elevation at well E2 B) M ean groundwater elevation at well E10 C) C oefficient of variation (COV) at well E2 D) COV at well E10. Labels are provided for the most significant parameters. A B Figure 3 effects of each parameter. Results are presented for the simulation of e vapotranspiration losses from the A) C onfined region B) U nconfined region of the basin. Labels are provided for the most significant parameters.

PAGE 121

121 Figure 3 8. Flow chart of separation of behavioral models from non behavioral models based on streamflow simulation by 340 GSA model runs. Also shown is the number of behavioral models (dotted brown box) for groundwater elevation at all the ten well locations. 0.5 Groundwater elevation R 2

PAGE 122

122 Figure 3 9. C omparison of observed and model simulated total streamflow at USGS 1500 (top) and USGS 2500 (bottom). Model simulated flow are shown for best behavioral models consideri ng flow and groundwater elevation predictions at all confined and unconfined locations (green line) and baseline model (black line). Also shown are model simulated total groundwater contributions (dotted lines) and measured flow hydrograph separation resul ts using specific conductivity measurements and end member mixing modeling for USGS 2500.

PAGE 123

123 Figure 3 10. C omparison of observed and model simulated total streamflow at wells E2, E3, E8, and E10. Also shown are groundwater elevation simulated by six models found behavioral towards streamflow and groundwater elevation for everywhere in the basin (green line) and results for baseline model (black line).

PAGE 124

124 Figure 3 11. Range and median values for most sensitive parameter identified du ring the Morris Analysis. Red, blue and green bars, respectively, represent values for entire sample generate for Morris Analysis (all 340 models), models behavioral for flow everywhere in the basin (28 models), and models behavioral for both flow and grou ndwater elevation everywhere in the basin (6 models). Also shown is the value of parameter used in baseline model (black bar).

PAGE 125

125 A B C D E F G H Figure 3 12. Influence of most sensitive parameters identified during Morris global sensitivity analysis on median values of output variables of interest Influence of A) ias_k on total ET losses and mean water table depth in the confined region B) ias_k on total and peak streamflow at stations 1500 and 2500 C) CR 1500 and 2500 and flow ratios at station 2500 D) UR _mannings on peak flow and flow ratios at station 2500 E) fas_k on total streamflow, groundwater flow contributions and flow ratios at station 2500 F) fas_ss on total groundwater exchange at station 250 0 during 2003 storm G) lai_mf on total ET losses from confined region H) wp on ET losses from unconfined region.

PAGE 126

126 CHAPTER 4 TOP OGRAPHIC, GEOLOGI C AND CLIMATIC CONTROLS OVER RIVER WATER SOURCES, MIXING DYNAMICS AND TRAVEL TIME DISTRIBUTIONS IN A LARGE COMPLEX BASI N 4.1 Background Growing concern over degrading water quality in springs, rivers and aquifers has led to many studies focused on determining sources, travel paths, and travel times of pollutants to these water bodies (Darracq et al. 2010 ; Botter et al. 2008 ; Van der Velde et al. 2010) Strong link ages between water pollutants and specific anthropogenic activities have been established (Katz 2004 ; Katz et al. 2001 ; Katz and Griffin 2008 ; Mueller Warrant et al. 2012) and various water quality management measures have been proposed and implemented to protect and rest ore natural water bodies throughout the world over the past few decades (Duriancik et al. 2008 ; Maresch et al. 2008 ; Osmond et al. 2012) A recent review of U. S Nonpoint Source (NPS) watershed projects in last four decades indicated that little or no water quality improvement has been achieved in most of the watershed scale NPS studies across the country (Meals et al. 20 10) Uncertainty in groundwater travel paths and travel time distributions was found to be one of the important underlying cause s of failure in many case studies discussed in the paper. The review highlighted the need for better characterization of hydra ulic residence time, transport rate and flowpaths and more specifically groundwater travel time, amongst several other design and management factors involved in watershed scale NPS projects. The mean transit time of water in a basin is a single integrated measure that has been frequently used as a hydrological descriptor for various storage and flow characteristics of a basin (Basu et al. 2012 ; McGuire and McDonnell 2006) Moreover,

PAGE 127

127 the tr ansit time distribution (TTD) has been used to provide a robust stochastic description of hydrologic functioning of basins in a single curve (Botter et al. 2010) Previous efforts to identify first order control s on transit times in contrasting geographical regions indicated that TTD of water within a basin is influenced by varying climatic and landscape characteristics. For example, McGuire et al. ( 2005) found a link between various topographic indices and transit time in steep landscapes of Western Cascades, USA. Soulsby et al. ( 2006) found that soil hydrologic characteristics exert significant control on mean transit time (MTT) in a gla ciated catchment in Scotland. Although insightful most of these studies on transit times focused on small experimental catchments ( <10 km 2 ; e.g. Dunn et al. 2010 ; Roa Ga rca and Weiler 2010 ; Birkel et al. 2012) with few er studies focusing on meso scale (>100 km 2 ) and macro scale (>1000 km 2 ) basins (Darracq et al. 2010 ; Speed et al. 2010 ; Tetzlaff et al. 2011) Small scale experimental catchments ( e.g. Kirchner et al. 2000) have an important place in conceptual and numerical hydrologic model development as they provide testbeds for emerging theo ries and improve our understanding about hydrologic functioning of natural system s. H owever, most of the water resources management decisions are made in large scale basins where we posses limited pr edictive understanding (Dunn et al. 2008) A f ew studies have investigated how tracer dynamics and subsequently MTT evolve in large scale basin s that span contrasting geologic regions (Ogrin c et al. 2008 ; Tetzlaff et al. 2011 ; Koeniger et al. 2009) These studies focused on the dynamics of landscape controls and increased influence of groundwater inputs on MTT in the downstream direction where the lowland landscape feature become more imp ortant

PAGE 128

128 than upland landscape features. For example Tetzlaff et al. ( 2011) found an overall increase in MTT estimates as more groundwater mixe d with water from catchment headwaters. Koeniger et al. ( 2009) found that water ages in lowlands of 46,000 km 2 Weser catchment in Germany range d from 14 50 years and increase d in the downstream direction as the groundwater inputs gain ed importance over headwater contributions. One of the most discussed and critical assumption s in TTD modeling studies has been the time invariant characteristic of the functions used to represent TTD (Botter et al. 2010 ; Botter 2012 ; 1982 ; Van der Velde et al. 2 010 ; McMillan et al. 2012 ; Hrachowitz et al. 2010b ; Rinaldo et al. 2011) Kirchner et al. ( 2000) used spectral analysis to show that a time invariant Gamma distribution with fit both short (sub daily) to long (>1 year) term tracer travel times in 1 3.5 km 2 headwater catchments at Plynlimon Wales. Several other catchment studies reported satisfactory results while applying different time invariant functions (e.g. exponential, piston flow a mongst others) for TTD ( see benchmark review by McGuire and McDonnell 2006 for a comprehensive list ). Beven ( 2010) suggested the use of a hypothesis testing framework while modeling TTD, instead of using a pre assumed form of TTD. Recent work by (Botter 2012 ; Botter et al. 2010 ; Van de r Velde et al. 2010 ; Hrachowitz et al. 2010b ; Hrachowitz et al. 2010a) presented an argument that TTD are in general time variant and reflect variability in climate and hydrogeologic conditions within basins. Hra chowitz et al. ( 2010a) used long term hydrometeorological and tracer data available for contrasting catchments in the Scottish Highlands (6.8 9.6 km 2 in

PAGE 129

129 area) to study inter annual and intra Gamma distribut ion commonly used in TTD studies R esults showed a and rainfall intensity above catchment specific thresholds ; a relationship with hydrologic characteristic of the catchments. It was further suggested that time variant TTDs could help in more realistic representation of contaminant transport through catchments. Botter ( 2012) showed that water residence time pdfs (defined as pdfs of ages of the water particles stored in a catchment at a given time) reflect rainfall variability with mean resident time dynamics dependent on rain event frequency, whereas travel time pdfs (defined as pdfs of time water particles spent in a basin) are streamflow dependent. A numerical particle simulation study on a 6.6 km 2 lowland catchment in the Netherlands indicated high variability in TTD s mostly due to the influence of rainfall a nd evapotranspiration dynamics (Van der Velde et al. 2010) This study showed tha t the numerical TTDs were highly dynamic and irregular with spikes reflecting rainfall and evapotranspiration variability. W ith few exceptions both small and large scale TTD modeling studies have used predetermined time invariant functions to represent TTD ; the exceptions being studies which relied on generalized analytical models (Botter 2012 ; Botter et al. 2010 ; Rinaldo et al. 2011) or have applied time variant TT D in very small basin s ( <10 km 2 ; e.g., McGuire et al. 2007 ; McMillan et al. 2012 ; Hrachowitz et al. 2010a ; Van der Velde et al. 2010) Moreov er, previous applications of time variant TTD in modeling have mostly used simplistic lumped parameter models (e.g. McGuire et al. 2007 ; McMillan et al. 2012) W hile numerous studies have pr ovid ed useful insights on time variant TTDs and their effect on solute (or simply water) transport, relatively little advancement has been

PAGE 130

130 made in study ing hydrologic controls over time variant TTD in large scale basins (> 10 3 km 2 ) spanning contrasting hyd rogeologic regime s. Therefore modeling efforts that examine time variance in TTD on mixing dynamics of major water sources within large scale basins under varying hydrologic conditions are needed U nderstanding and predicting TTD of water in large basin s is important for water quality management because TTD s provide insights on how water is stored released and travels through a basin (McGuire and McDonnell 2006) which in turn directly controls contaminant transpo rt and delivery to receiving water bodies For instance, a high proportion of long travel times indicates longer time for biogeochemical reactions in the subsurface env ironment (McGuire and McDonnell 2006) and impl ies prolonged release of more persistent pollutants via diffuse ground water contributions to streams (Butscher and Huggenberger 2009) In contrast, systems with a high proportion of short travel times are more vulnerable to highly reactive contamin ants because the faster travel times in these systems do not allow sufficient time for pollution mitigating processes such as adsorption degradation and decay (Butscher and Huggenberger 2009) Studies of lar ge basins have reported dominance of fast surface processes in uplands regions and increased dominance of slow groundwater processes in lowland regions (e.g. Tetzlaff et a l. 2011 ; Ogrinc et al. 2008 ; Koeniger et al. 2009) These large basins are prone to dual vulnerability as they may receive both rapidly transported reactive contaminant from the upland regions and persistent contaminant s from regional groundwater sour ces. Numerical p article tracking schemes used in conjunction with integrated physical hydrologic models offer a relatively inexpensive alternative to extensive field sampling

PAGE 131

131 based studies to understand factors controlling TTDs in catchments. P article tr acking experiments typically introduc e solute particles (or water particles) into flow fields produced by a groundwater model and subsequently record the with time R esults obtained include the particle flow path s and transit time s th rough the hydrologic system. If introduced in sufficiently large number, particles TTDs can be obtained without defining the functional form of the TTD a priori (e.g. Kollet and Maxwell 2008b ; de Rooij et al. 2012) Recent advancement s in particle tracking schemes ha ve made it possible to track particles flow paths through fully integrated surface water groundwater atmosphere (de Rooij et al. 2012) A s such all the relevant physical processes such as three dimensional saturated or unsaturated flow, two dimensional overland flow, evapotranspirative fluxes from the land surface and shallow subsurface environment and infiltration into the subsurface can b e accounted for while evaluating hydrologic and geologic factors controlling TTDs in a basin This information may then be useful for predicting response of water bodies to changes in land and water management practices. In this paper we present result s obtained by applying Slim Fast, a novel particle tracking scheme (de Rooij et al. 2012) which tracks solute or water particles through a fully integrated surface water groundwater land surface atmosphere syst em, to a large complex basin in North Central Florida. Slim Fast works in conjunction with ParFlow.CLM, a fully integrated three dimensional variably saturated groundwater two dimensional surface water one dimensional land surface model that may be driven by spatiotemporally variable atmospheric forcing data. Particle tracking in transient flow fields was used to estimate time varying TTDs for water particles arriving at multiple

PAGE 132

132 river locations within the Santa Fe River Basin (SFRB) during wet and dry hy drologic conditions in order to provide insights on factors controlling sources, travel paths, and travel time distributions for streamflow throughout the basin. In particular particle tracking experiments were designed to investigate (1) how river wate r sources evolve, interact and integrate in large basins spanning contrasting hydro geological regimes; (2) how the spatiotemporal dynamics of dominant water sources, travel paths and TTD change under wet and dry hydrologic conditions; and (3) what topogra phic, geologic and climatic characteristics act as first order controls on TTDs and how these controlling factors vary with space and time in large basins 4.2 Santa Fe River Basin The SFRB is a mixed use basin spanning an area of about 3700 km 2 in North Cent ral Florida. A complete description of the basin and its hydrogeologic characteristics are presented elsewhere (Srivastava et al., 2013 a ). Briefly, the SFRB consists of three distinct yet interconnected hydrogeological units (Figure 4 1; Schneider et al., 2008); (1) a confined region (CR) where the regional karst aquifer (Upper Floridan Aquifer System (UFAS)) is confined by the clay dominated Intermediate Aquifer System (IAS) which is overlain by an unconfined sand dominated surficial aquifer system (SAS); (2) an unconfined region (UR) where due to erosion the IAS confining unit is completely missing and the UFAS is overlain by a thin sand layer ; and (3) the semi confined region (SCR) located where the CR transitions to the UR and the UFAS is variably conf ined by leaky or discontinuous clay rich beds. The CR is characterized by flat topography and poorly drained soil resulting in shallow water table conditions and a well developed stream network. In contrast, the UR has well drained soil which along with a major conduit network within the UFAS creates

PAGE 133

133 internal subsurface drainage toward the main river. Eogenetic karst features found in Florida have retained their high matrix porosity in comparison to the more commonly known teleogenetic karst systems with l ow intergranular porosity (Vacher and Mylroie 2002 ; Florea and Vacher 2007) Therefore, in addition to the conduit network, primary matrix porosity in the UFAS plays an important role in groundwat er storage and flow towards the river 4.3 Methods In this study, the TTDs for water particles were estimated following a two step experiment which included; (1) development and testing of a 3D integrated surface subsurface land process model ParFlow .CLM for t he SFRB (Srivastava et al., 2013 a ); (2) use of three dimensional subsurface pressure fields and two dimensional surface pressure fields obtained from the integrated model in S lim Fast to simulate the flow path and transit times of tagged water particles in the SFRB. Det ails regarding the methodology for developing the integrated model and its simulation results (step 1) are descri bed in Srivastava et al., (2013 a ). Details on the numerical technique used to develop Slim Fast are summarized in (de Rooij et al. 2012) For brevity here we provide details pertinent only to the work presented here 4.3.1 I ntegrated H ydrologic M odeling A n integrated three dimensional surface water groundwater land surface model w as previously deve loped for the SFRB (Srivastava et al., 2013 a ). The numerical platform used to develop the integrated model included : (1) ParFlow (Ashby and Falgout 19 96 ; Kollet and Maxwell 2006 ; Kollet and Maxwell 2008a ; Maxwell and Miller 2005) which simulates three dimensional variably saturated subsurface flow using Richards equation along with two dimensional overland flow using the Kinematic

PAGE 134

134 wave equation alon and (2) a modified version of Common Land Model (CLM; Dai et al. 2003 ; Maxwell and Miller 2005 ; Kollet et al. 2009) which calculates near land surfac e atmospheric fluxes, such as evapotranspiration losses, snow accumulation and melt processes, latent, sensible and ground heat fluxes, as a function of atmospheric variables such as precipitation, air temperature, pressure, wind speed, specific humidity, solar radiation and soil moisture in shallow subsurface and land surface. ParFlow simulates the moisture content in the shallow subsurface environment and pass es it to CLM which uses it along with various atmospheric variables (listed above) to calculat e energy and water balance in near land surface environment. CLM then returns net infiltration (or exfiltration in absence of rain) to ParFlow for use in calculat ing various surface and subsurface water balances and fluxes This process is continued iterat ively until the end of the simulation As described above, and in Srivastava et al. (2013 a ), SFRB is a complex and large eogenetic karst basin spanning three distinct yet interconnected hydrogeological regions with varying hydrologic behavior in terms of surface groundwater interactions. The ParFlow .CLM application for the basin was developed and tested for the period January 1, 2000 through December 31, 2008 which included diverse climatological conditions including several extreme events. For instance ye ar 2004 was an extremely wet period during which four major hurricanes and a major tropical storm crossed the state of Florida causing flooding over much of the study area. Conversely, years 2006 2007 were one of the driest periods ever recorded for the basin and were followed by

PAGE 135

135 relatively wet conditions as a result of a tropical storm which occurred in September 2008. The reliability of model predictions was established by validation of model predictions against long term field observations at multiple locations. In addition to commonly used criteria of streamflow and groundwater elevation predictions, high resolution (daily) as well as discrete (monthly) specific conductivity measurements based hydrograph separation results were used to evaluate surfac e and groundwater flow dynamics at multiple locations within the basin. Overall, the model predictions were found to provide a reasonable representation of regional hydrologic conditions over the study period. To study the effect of variability in subsurf ace properties on source water locations, travel and TTDs characteristics, the baseline model (described above) was rerun to obtain following three additional models; (1) baseline model with no high hydraulic conductivity model with higher hydraulic conductivity m/hr); and (3) four order of magnitude higher hydraulic conductivity in IAS (10 5 m/hr as compared to baseline model value10 9 m/hr). Pressu re fields from the first two models were used to study subsurface characteristics control over TTD in the UR whereas the last model were used to study the same for the CR 4.3.2 S urface S ubsurface P article T racking Particle tracking schemes accept pressure fiel ds from flow models to reconstruct flux and flow velocity fields within the domain which can be used to trace flow path s of small imaginary particles introduced in the flow field (Anderson and Woessner 1992) A m odified version of the Lagrangian particle tracking scheme Slim Fast (Maxwell et al.

PAGE 136

136 2003 ; Maxwell and Kastenberg 1999) was used in this study to simulate TTDs of water particles in the SFRB. Previous version s of Slim Fast used the classical advective particle tracking method of Pollock (Pollock 1988) and allowed tracking of particles only within three dimensional variably saturated subsurface flow f ields T he modified particle tracking scheme used here allow s the movement of particles through a fully coupled surface and subsurface flow system using flow fields generated by ParFlow (de Rooij et al. 2012) Pathl ines in both the subsurface and overland domain s are obtained using (Pollock 1988) which is a simple linear interpolation scheme that computes velocity components at a given particle location and subsequently move s the particle (within a given cell) to a new location (by multiplying the velocity components by a finite time step. Transfer of particles from the subsurface to the surface flow system or the atmosphere occurs if their trajectory inter sects the land surface. Transfer of particles from the overland to subsurface domain or the atmosphere is achieved by using pathline specific transfer probabilities determined based on the mass balance of pathline (de Rooij et al. 2012) The purely advective particle tracking scheme used in this study do es not account for local dispersion or molecular diffusion 4.3.3 River W ater S ource M ixing D ynamics and T ravel T ime D is tribution Particle tracking schemes, when used in conjunction with integrated hydrologic models, provide maps of the spatial extent of varying water sources, travel paths and travel times under transient hydrologic conditions. For the experiments conducte d in this study, 100,000 particles were placed in the river cell s for stations 1500 and 2500 in the CR and UR of the SFRB respectively (Figure 4 1) ParFlow .CLM generated pressure fields from the baseline ParFlow.CLM model ( Srivastava et al. 2013 a ) and i ts three

PAGE 137

137 variants (presented in section 3.1) were used to calculate daily cell velocities The cell velocities were reversed in time to transport the water particles back to their source location. The nine year daily c ell velocity sequence was repeated as necessary allow backtracking for 1000 year s. This approach is similar to that adopted by (Kollet and Maxwell 2008b) who repeatedly used one year of pressure field from a prev ious ParFlow simulation to obtain the age di stribution of groundwater particles in Little Washita Basin Slim Fast keeps track of and writes output files with information on particle locations as they move through the model domain with passage of time. W hen executed in backtrack mode the particles start from the end point (i.e. at river stations 1500 or 2500 in this study) and move backward in time towards their location of origin, either one of the domain boundaries or the atmosphere (which implies that the particle arrived as rain water within the domain). The location where a particle stops moving through the flow field indicates its entry point into the domain, and the time it takes to reach this entry point (from the starting river location) is its transit time in the domain. If large number of particles are introduced in the domain, transit times for all the particles can be used to construct the TTD of water arriving at a given river location at a given time. Releasing large numbers of particles at various locations and at multiple times allows investigation of the spatial variability as well as the effect of varying hydrologic conditions on the shape of TTDs. To study the effect of antecedent hydrologic conditions and storm magnitude on source wate r location and TTDs for stations 1500 and 2500 we conducted particle tracking experiments at multiple points along the streamflow hydrographs from three storms with varying antecedent hydrologic

PAGE 138

138 conditions and an extended drought period (storms a, b, and e; droughts d1 through d3 in Figures 4 3 4 5). Results obtained for stations 1500 and 2500 were used to deduce and represent hydrologic behavior in the confined and unconfined hydrogeologic regions within the SFRB, respectively 4.4 Results and D iscussion 4.4.1 S treamflow A ge D istribution U nder T ime V arying H ydrologic C onditions In the CR of the watershed streamflow is episodic with virtually no i nter event subsurface baseflow. Thus there was no streamflow during the drought period d1 d3 at station 1500 (Figure 4 3) so CR particle tracking experiments were only performed for subsurface storage w as relatively high (Figure 4 series of moderate intensity rain events that were distributed fairly uniformly over a period of 2 months (i.e. approximately from Feb 15 2003 (a1) through April 5 2003 (a7) in Figure 4 4). and were of shorter duration (Figure 4 4). Hydrologic conditions during the period when these storms occurred are summarized in Table 4 1. Cumulative distribution functions (CDFs) for water travel times to station 1500 for the three selected storms are shown in Figure 4 6. The m edian age of water arriving at station 1500 was less than 7 days throughout the rising and falling limbs of the hydrographs for each of the storms, and 95% of the water arriving at station 1500 was less than 15 days old. The short median travel times for all TTDs at station 1500 indicates the dominance of rapid surface flow paths in this region. The variation among the TTDs over the hydrograph and by storm in dicates the influence of rainfall spatial

PAGE 139

139 pattern, rainfall intensity and antecedent conditions on TTDs for streamflow in the CR of the basin. Slightly longer travel times predominate at the beginning and end of each storm event, and there is more variabi lity in the shapes of these TTDs across storms than those corresponding to peak streamflows. These results indicate that antecedent conditions do not have a strong control on TTDs for these large storms, however they do influence the total streamflow gene rated by the storms. For example despite approximately similar total rainfall, similar TTDs shapes, and roughly equivalent median water age (<7 days), less total streamflow and a lower streamflow to rainfall volume 1). These observations indicate that streamflow in the CR is triggered only when sufficient rainfall occurs to raise the water table above a particular threshold that creates direct floodplain channel connectivity The pattern of flood plain connectivity controls the TTD which is relatively time invariant, particularly at peak flow; however the volume of rainfall above the trigger threshold varies with storm magnitude and antecedent conditions and controls the tot al streamflow volume. At station 2500 in the UR TTDs showed a dramatically different shape from those for station 1500 (Figure 4 7). TTDs for the three storm events showed similar structure than ~15days, or ~0.04 years old) increasing rapidly with the rising limb, reaching its peak value during peak of the storm hydrographs, and dropping more slowly with the recession limb. The fraction of new water ranges from about 12% of total flow duri ng baseflow and drought conditions, up to about 90% of total flow during the largest storm peak (storm b). A break in the slope of the TTDs for all storms occurred at ~ 15 days, and indicated a

PAGE 140

140 relatively small fraction of moderately old (i.e. >15 days to 20 years) near stream subsurface contributions for all storm events (~ 6% during storm b to ~20% during storm e) to and the drought conditions (~44% during drought d1 d3). A second break the slope of the TTDs occurs at about 20 years reflecting a va riable fraction very old subsurface contribution to streamflow, ranging from ~ 4% at the peak of the largest water produces long tails in the TTDs at station 2500 indicativ e of fractal behavior similar to that observed in experimental catchments using environmental isotopes and geochemical tracers (e.g. Kirchner et al. 2000) TTDs during drought periods and in on for storm events showed remarkable resemblance to each other (Figure 4 occurring directly on t he stream channel. The median age of water age of stream water at station 2500 varies with position on the streamflow hydrograph for each storm (Figure 4 8). During baseflow conditions (before and after the storm event) the median age of water arriving a t station 2500 was similar to that during droughts d1 d3 (~ 16 years; Figure 4 8). Over the storm event the median age drops to <15 days at peak flows, with a more rapid drop and faster rebound in median age occurring for the wetter antecedent conditions (i.e. storms a and b). Storm e, with a drier antecedent condition and smaller streamflow to rainfall ratio, showed a slower drop in median age over the rising limb of the hydrograph, and a more rapid increase in median age over the receding limb of the hy drograph indicating a greater contribution of subsurface flow for this event. This

PAGE 141

14 1 pattern of shift towards younger streamflow with increased antecedent moisture and therefore increased connectivity in the basin is consistent with findings in previous stud ies (e.g. McGuire et al. 2007) 4.4.2 Streamflow O rigin and F lowpaths O ver C ontrasting H ydrogeologic R egions U nder D ifferent H ydrologic C onditions Although insightful, TTDs alone do not provide information regarding the spatial origins of various sources of water within a basin. Particle tracking schemes have an advantage over methods which rely on only experimental tracer data to estimate TTDs in that they not only provide information on particles travel times but also p rovide information on source locations and flowpaths. The temporal variations in spatial extent of water sources help to identify dominant topographic, geologic, climatic and hydrologic processes that control observed hydrologic behavior at any given time and location in the basin. Figures 4 9 and 4 10 illustrate the spatial extent of the contributing source area for streamflow at station 1500, as well as the time taken to travel from the point of origin ( as rainfall) to the streamflow station during diff 3). As discussed previously these storms vary in terms of antecedent hydrologic conditions (Table 4 1; Figure 4 3). For the storm b, with higher rainfall volume and wetter antecede nt conditions, a more rapid increase in contributing area occurs during the rising limb of the hydrograph. In contrast during storm e rainfall must initially fill the soil moisture deficit and raise the water table before generating significant runoff. H owever at peak flows, and at the end of the receding limb of the hydrographs, the contributing areas for the two storms look quite similar in spatial extent. In both cases, rainfall ultimately raises the water table to the land surface in the

PAGE 142

142 vicinity of t he stream, resulting in increased overland connectivity of river with the immediate floodplain resulting in high flow conditions. At the last point on the streamflow recession curve the contributing area is reduced from the peak, but is still substantially larger than the contributing area for the equivalent flow on the rising limb of the hydrograph, especially for storm e with the drier antecedent conditions. At this end of the storm, flow at station 1500 is sourced mainly by channel water collected durin g early part of storm events and surface contributions from low lying topographic depressions temporarily connected to the channel in the headwaters of the basin (e.g. Lake Santa Fe and the Santa Fe Swamp which constitute the headwaters of the Lower Santa Fe River). Based on the spatial extent, age distribution and the original source of particles identified b y particle tracking (Figure 4 9G and 4 10H ) recent rainfall falling directly on the floodplain was confirmed to be the source of streamflow for stati on 1500. In contrast to station 1500 in the upper CR of the SFRB, streamflow at station 2500 in the UR was found to be a mixture of water of widely varying ages and sources. Based on shape and characteristics of TTDs water arriving at 2500 was classified a s 1000 years). Figures 4 11 4 points along 14 4 Also shown in the plots is the age distribution within the water particles belonging to the three age categories. Comparison of spatial extent and location of water particles within CR of the basin that travels to the stream following fast surface flow paths, or as rain

PAGE 143

143 falling directly on the unconfined stream channel, and that its spatial evolution with time follows that of station 1500 (Figure 4 9 and 4 10). In the early part of the storm hydrographs (Figure 4 11 A B and 4 14 A B C ) the unconfined river is disconnected from the lower UR of the watershed, thus rain dir ectly over the channels is the only on UR of the basin then infiltrated through the vadose zone to recharge the UFAS, with very small portion arriving from the CR via stream channels (Figures 4 12, 4 15). Within the UFAS this water either moved directly to the nearby stream channel or to nearby high hydraulic conductivity zones which then rapidly delivered it to the stream. The geometry of the high hydraulic conductivity regions and the contrast in hydraulic conductivity between the slow porous matrix and the high hydraulic conductivity regions explain the wide range of age of the moderately old water (~15 days to 7300 days ( 20 years)), with water age increasing dramatically with distance from the high hydraulic conductivity zones (see Figure 4 2 for high hydraulic conductivity zone locations for ing to the color scale in Figures 4 12 and 4 15 the red particles are older than 20 years. Thus these particles are still moving backward in the domain toward their point of entry 20 years prior to the storm event. Blue to orange particles exited the dom ain at the location indicated by their position and the time (less than 20 years) indicated by their color. (Figures 4 13 and 4 16). This rainfall moved through the vado se zone to recharge the

PAGE 144

144 UFAS in regions farther from the high hydraulic conductivity zones, leading to longer travel paths and longer travel times to the high hydraulic conductivity zones and to the river. It should be noted again that in Figures 4 13 and 4 16 blue to orange colors indicate that the particles have arrived at their entry point into the domain (i.e. they are no longer moving). A very small fraction of the particles (red) are located with the confined UFAS, still traveling toward the domain b oundaries 1000 years prior to the study period 4.4.3 Geologic C ontrols O ver W ater T ravel T ime and F low P aths in a n U nconfined A quifer To investigate the control of hydraulic conductivity contrasts between the porous matrix and high hydraulic conductivity regi ons on source water locations, travel paths and TTDs, ParFlow Slim Fast hydraulic conductivity regions were assigned the same hydraulic conductivity as the porous matrix (9 m/ hydraulic conductivity hydraulic conductivity regions were assigned a hydraulic conductivity of 600 m/hr (versus 270 m/hour for the base case). Significant variations in the shape of TTDs and spatial pattern of part icle source locations were observed at station 2500 between these cases for storm events as well as during droughts (Figures 4 7, 4 17 and 4 18). Similar behavior was observed for all storm events, so for brevity only storm e and the drought d results wil l be discussed here. Analogous plots for storm b can be found in the Appendix Figure A 3. Comparison of TTDs at 2500 for the three cases shows that the fractions of new water, moderately old water and very old water vary differently throughout the storm hydrograph based on the hydraulic conductivity contrast between the high hydraulic

PAGE 145

145 conductivity zones and the porous matrix (Figures 4 17 and 4 19). At baseflow conditions (point e1 before and point e7 after the storm) the no conduit case shows the lowest fraction of new water and moderately old water, and the highest fraction of very old water compared to the other cases (Figure 4 19). In contrast during high flows (points e3 through e5) the no conduit case shows the highest fraction of new water, the lo west fraction of moderately old water, and about the same fraction of very old water as the other cases (Figure 4 19). The inflection points of the hydrograph (points e2 and e6) show transitional behavior between these two regimes. Figure 4 20 shows that t he influence of hydraulic conductivity contrast on median age of the streamflow is most pronounced on the falling limb of the hydrograph (i.e. point e6). Points at the beginning and end of the streamflow hydrograph (e1, e7) have a median age of 16 and 12 y ears (base case and conduit cases) to 26 years for the no conduit case. At high streamflows (points e3 e5) the median water age is <15 days for all cases. However at the inflection point of the falling limb (point e6) the median age ranges from <15 days for the no conduit case, to < 60 days for the base case, to >1 year for the high hydraulic conductivity case. The increase in conduit hydraulic conductivity also increases the spatial extent of the influence of high hydraulic conductivity zones by draining more of the surrounding matrix and thus facilitating more 21). This hy draulic conductivity 19) and by the increased spatial distribution of moderately old water blue particles in the UR for storm e (Fig ure 4 21), storm b (Appendix ) and the drought condition (Appendix ).

PAGE 146

146 fractions of water contributed to streamflow during storm event e reveals useful insights on the influence of conduit hydraulic conductivity on surface groundwater interactions within the UR of the basin (Figure 4 22). Srivastava et al., (2013a) compared ParFlow predicted versus EMMA estimated subsurface contributions to streamflow during storm events for the base case model and concluded that the model was missing an important transient source of subsurface water during the storm recession period. Figure 4 22 compares EMMA versus ParFlow predictions of subsurface contributions to stream flow for the three conduit hydraulic conductivity cases during 2008 which includes storm e. This figure shows that increasing the hydraulic conductivity of the conduits i ncreases the subsurface contribution to stream flow during the recession of the streamflow curve and thus improves the fit of model predictions of groundwater contributions to EMMA results. This indicates that transient source of old water missing from the base case streamflow recession curves may be released by high hydraulic conductivity conduits that are not well represented in the base case. Comparison of TTDs during storm e at 1500 and 2500 for the high IAS hydraulic conductivity test case showed that the TTDs for the base case and its high hydraulic conductivity variant are similar. However, the TTDs at 2500 showed that an increase in IAS hydraulic conductivity decreased the fraction of new water originating in the CR and increase in moderately old and very old water fractions at 2500 (Figures 4 23) 4.4.4 On T ime I n variant A ssumption s in TTD M odeling Most of the studies of TTDs rely on inverse modeling techniques where rainfall and tracer time series are used to estimate the TTD of water (or some solute of i nterest) through an experimental or modeled catchment (Kirchner et al. 2001, McGuire and

PAGE 147

147 McDonnell 2006) These methods require a pre defined functional form of the TTD which is often ass umed to be invariant in time, an assumption which has been debated in the literature (Botter et al. 2010 ; Botter 2012 ; Van der Velde et al. 2010) Recent work by (Van der Velde et al. 2010 ; Rinaldo et al. 2011 ; Botter et al. 2010b ; Botter 2012) showed that the TTD are not smooth and may vary with precipitation patters and antecedent conditions. Hrachowitz et al. ( 2010a) showed that the beta parameter of the gamma function fit to TTDs could be related to precipitation intensity and that the TTDs can vary int er annually. McMillan et al. ( 2012) found that in the 0.9 km 2 Loch Ard Burn catchment a time invariant TTD was not able to reproduce the expected tracer response at shorter time scales. Longer time scale response s, on the other hand, were found to be reasonably well represented by a steady state TTD. All these studies focused on small catchments (<10 km 2 ) and have not been tested at large scale. Figure 4 24 shows the TTD observed during various times within sto for station 1500 along with the best fit log normal TTD for each curve. As evident from Figure 4 24 the lognormal function was able to capture general shape and characteristics of the TTDs as long its parameters were allowed to vary over th e storm. Furthermore this figure and T able 4 2 shows that the lognormal TTDs fit to the high flow portions of the hydrograph, when the maximum contributing area of the flood plain is connected to the river, have quite similar parameters. TTDs for station two time periods, the first consisting of particles less than 15 days old and the second consisting of part icles older than 15 days; then (2) rescaling each CDF by the total

PAGE 148

148 number of particles that arrived each time period. Similar to station 1500, a lognormal were allowe d to change throughout the hydrographs (Figure 4 all storms and drought periods were well approximated by the gamma function with parameters that indicate fractal behavior and are fairly constant over time (Figure 4 25 and Table 4 and dry conditions is in agreement with previous findings from small scale s tudies (e.g. McMillan et al. 2012) 4.5 Summary and C onclusions Slim Fast, a fully coupled surface subsurface particle tracking scheme, was used in conjunction with ParFlow .CLM, an integrated 3D surface water ground water land surface hydrologic model, to study water source areas, travel paths and transit time s within the Santa Fe River Basin. In addition to improving our understanding of the hydrologic functioning of this large complex basin, TTDs and spatiotemporal dynamics of major flow components estimated using Slim Fast provided a valuable additional diagnostic tool for model evaluation and improvement. R esults indicate that in the CR of the SFRB channel floodplain connectivity, which depends on topography and pr operties of the confining layer, exerts strong and largely time invariant control on the shape of the TTD. However antecedent storm conditions and storm volume exert a transient and non linear (threshold) control on total volume of streamflow generated in this region In the UR of the SFRB the total streamflow TTDs were contribution that originates as rainfall and travels overland in the CR

PAGE 149

149 rainfall in the UR that recharges groundwater then move s to the stream via both localized high hydraulic conductivity flow paths (conduits) and the porous limestone matrix. The new water TTD is controlled by the same factors which control tota l streamflow TTDs in the CR whereas the shape of the old water controlled by the regional groundwater flow system and the contrast in hydraulic conductivity between high hydraulic conductivity regions and the porous matrix. The relative weight s of the new water TTD and old water TTD to the total streamflow TTD in the UR are a function of antecedent conditions and storm magnitude via their control on streamflow volumes generated in the CR These findings underscore the dual vulnerability of the Santa Fe River to the rapid transport of reactive contaminants (e.g. pesticides) generated in the CR during storm events as well as long term persistent transport of more stable contaminants (i.e. nitrates) that have been leached into the unconfined F loridan aquifer system groundwater system. T he media n age of streamflow in the CR of the basin ranged from approximately 1 day at the peak of storm hydrographs to approximately 7 days at the end of stormflow recession, with travel time distributions that varied over the hydrograph and between storms, but were generally well fit by log normal distributions with time varying parameters These results are similar to previous theoretical studies (e.g. Van der Velde et al. 2010) that at short time scale s TTD s are transient in nature and not well represented by a single function. The median age of subsurface contributions to streamflow in the lower portion of the basin wa s approximately 17 years, in close agreement to me an water ages estimated in springs that feed the river in this region using both chlorofluorocarbons (CCl 3 F, CCl 2 F 2 and C 2 Cl 3 F 3 ) and environmental isotopes

PAGE 150

150 ( 18 O/ 16 O, D/H, 13 C/ 12 C, 15 N/ 14 N) (Katz 2004 ; Katz et al. 2001 ; Katz and Griffin 2008) The travel time distribution of old subsurface contributions to streamflow in the unconfined area was well fit by a gamma distribution showing fractal properties that did not significa ntly vary over time. T he coupled ParFlow .CLM Slim Fast modeling system was used to explore the effect of geologic heterogeneity on surface water groundwater interactions in the SFRB. These experiments showed that increasing the hydraulic conductivity of the confining layer decreased total streamflow in both the CR and unconfined regions of the basin by lowering water tables in the CR. However changes in the hydraulic conductivity of the confining layer did not significantly affect the shape of the TTDs or the median age of streamflow in the CR for rainfall events large enough to induce direct flood plain stream connectivity. These experiments also showed that while increasing conduit hydraulic conductivity had a relatively minor effect on total streamfl ow in the UR (Srivastava et al 2013b) increasing conduit hydraulic conductivity significantly affected the spatial pattern of the sources of subsurface contributions to streamflow. Furthermore increasing conduit hydraulic conductivity increase d groundwa ter contributions to the recession limb of the streamflow hydrograph, producing behavior more similar to that observed from EMMA (Srivastava et al 2013a). Th ese results indicate that more accurate representation of conduit matri x interactions via finer d iscre tization and/or explicit representation of discrete conduits will likely improve model representation of streamflow and subsurface flow dynamics in the UR of the basin

PAGE 151

151 Table 4 ation 1500. Table 4 2. Parameter values for best fit function for overall CDFs at station 1500 and

PAGE 152

152 Figure 4 1 Map of the model domain Also shown are the extent s of the Santa Fe River Basin, locations of the river channel s, monitoring wells, USGS stream gages and major hydrogeological regions in basin.

PAGE 153

153 Figure 4 2 3D mesh used to define the domain A) the extent of the three major hydrogeological regions based on Floridan Aquifer Vulnerability Assessme nt ( FAVA ) dataset within the mesh along with their indicator variables Surficial Aquifer System blue region Intermediate Aquifer System green region and Upper Floridan Aquifer System red region. B) Locations of high hydraulic conductivity zones (b lue cells) and main channels (red cells) A B

PAGE 154

154 Figure 4 3 Time series of ParFlow.CLM simulated daily flow at station 1 500 and station 2500 Storm events a, b, and e are marked along with drought d. Also shown is mean daily rainfall received by the domain (inverted bars in top plot).

PAGE 155

155 Figure 4 3 Continued.

PAGE 156

156 Figure 4 4 Hydrograph points where particle tracking experiments were conducted during year 2003 and year 2004 storm events at station 1500. These storms 4 3 Numbers next to each point show date and flow value corresponding to each point. Also shown is mean daily rainfall received by the domain (inverted bars).

PAGE 157

157 Figure 4 4. Continued

PAGE 158

158 Figure 4 5 Hydrograph points where particle tracking experim ents were conducted during year 2003 and year 2004 storm events. These storms are referred to 4 3 Numbers next to each point show date and flow value corresponding to each point. Also shown is mean daily rainfall received by the domain (inverted bars).

PAGE 159

159 Figure 4 5 Continued.

PAGE 160

160 Figure 4 6 Water particle age distribution at station 1500 for various points on storm events and Note that the water age is reported in days at station 1500.

PAGE 161

161 Figure 4 6 Continued.

PAGE 162

162 Figure 4 7 Water particle age distribution at station 2500 for various points on storm for base case model. new moderately old very old

PAGE 163

163 Figure 4 7 Continued.

PAGE 164

164 A B C Figure 4 8 Water particle median age versus flow for base case model at vario us time s during storm hydrographs for A) S B) S C) S with black circles present data for the rising limb and green lines with red circles present data for the recession period; value corresponding to storm peak is repre sented by orange circle. Also shown are the median water ages for the three drought periods (dotted circle on upper left corner of each plot).

PAGE 165

165 A B C D E F G Water age (days) Figure 4 9 Spatial extent from which water particles arrived at station 1500 during b A F corresponds to points b 1 b6 in f igure 4 4 The color scheme shows the age (days) of all the particles arriving at station 1500 b (G ) Source of water arriving at station 1500 during peak of b b3 in f igure 4 4 ) with blue and red dots representing rainfall or surficial aquifer system.

PAGE 166

166 A B C D E F G H Water age (days) Figure 4 10 Spatial extent from which water particles arrived at station 1500 during A G corresponds to points e1 e7 in f igure 4 4 The color scheme shows the age (days) of all the particles arriving at station 1500 ( H ) Source of water arriving at station 1500 during peak of e e4 in f igure 4 4 ) with blue and red dots representing rainfall or surficial aquifer system.

PAGE 167

167 A B C D E F G Water age (days) Figure 4 11 es arrived at station 2500 b F corresponds to points b 1 b6 in f igure 4 5 The color scheme shows the age (days) of all the particles arriving at station 2500 b ( G ) Source of water arriving at station 2500 during p eak of b b3 in f igure 4 5 ) with blue and red dots representing rainfall or surficial aquifer system.

PAGE 168

168 A B C D E F Water age (days) Figure 4 12 b A F corresponds to points b 1 b6 in f igure 4 5 The color scheme shows the age (days) of all the particles arriving at station b

PAGE 169

169 A B C D E F Water age (days) G Water age (days) Figure 4 13 b A F corresponds to points b 1 b6 in f igure 4 5 The color scheme shows the age (days) of all the particles arriving at station 2500 b G) P oint b6 result ( F above) with new scale shown below the figure which shows full range of age distribution of particles at the end of 1000 years

PAGE 170

170 A B C D E F G H Water age (days) Figure 4 14 A G c orresponds to points e1 e7 in f igure 4 5 The color scheme sh ows the age (days) of all the particles arriving a t station 2500 ) Source of water arriving at s tation 2500 during point e4 in f igure 4 5 with blue and red dots representing rainfall or surficial aquifer system.

PAGE 171

171 A B C D E F G Water age (days) Figure 4 1 5 n 2500 dur G c orresponds to points e1 e7 in f igure 4 5 The color scheme shows the age (days) of all the particles arriving at station

PAGE 172

172 A B C D E F G H Water age (days) Figure 4 16 A G corresponds to points e1 e7 in f igure 4 5 The co lor scheme shows the age (days) of all the particles arriving at station 2500 ( H) P oint e4 result ( G above) with new scale shown below the figure which shows full range of age distribution of particles at the end of 1000 years

PAGE 173

173 A B Figu re 4 17 Water particle age dis tribution at station 2500 at various points on storm for base case variants. A) N o conduit B) H igh er hydraulic conductivity conduit

PAGE 174

174 A B Figure 4 18 Water particle age dis tribution at station 2500 at various p oints during for base case variants. A) N o conduit B) H igh er hydraulic conductivity conduit

PAGE 175

175 Figure 4 19 Temporal dynamics of various sources of water at 2500 during various time within a storm hydro graph for base case and its variants Variants are hydraulic conductivity A B Figure 4 20 Water particle median age versus flow for base case model variant at igher hy draulic conductivity conduit B) N o high hydraulic conductivity conduit Blue line with black circles present data for the rising limb and green lines with red circles present data for the recession period; value corresponding to storm peak is represented by orange circle. Also shown are the median water ages for the three drought periods (dotted circle on upper left corner of each plot).

PAGE 176

176 A B C Figure 4 2 1 Spatial extent of water particle source location during poi nt e4 of storm A) Basecase. B) B ase case with high er hydraulic conductivity conduits C) B ase case with no conduit. The color scheme shows the age

PAGE 177

177 Figure 4 22 Compariso n of measured and ParFlow.CLM simulated total streamflow and subsurface contributions Also shown are end member mixing analysis based

PAGE 178

178 Figure 4 23 Water particle age dis tribution at station 2500 Results shown for various for high intermediate aquifer hydraulic conductivity variant of base case.

PAGE 179

179 A B Figure 4 2 4 Water particle age dis tribution at station 1500 A) During B) Du ring storm event (broken line of various colors). Numbers in bracket are best fit parameter values for each distribution.

PAGE 180

180 A B Figure 4 25 Water particle age dis tribution at station 2500 during various time in storm A) New water. B) Old water. Also shown are best fit lognormal distribution (variable colors dotted lines on top plot) and Gamma distribution (broken red line in bottom plot). Numbers in bracket are param eter of various distribution fitted to CDFs.

PAGE 181

181 C D Figure 4 25 Continued.

PAGE 182

182 E Figure 4 25 Continued.

PAGE 183

183 CHAPTER 5 CONCLUSIONS, CON T RIBUTIONS AND RECOMMENDATIONS FOR FUTURE REASEARCH Understanding the hydrologic functioning of a river basin at any scale is a co mplex problem. Hydrologic models are commonly used to gain insights on the hydrologic functioning of river basin s to make water management decisions, and to design future field experiments to improve conceptual models and parameter estimation to enhance f uture modeling of the system. This study demonstrated the potential of a fully integrated modeling platform to simulate and understand complex surface water ground water near land surface hydrologic processes and their interactions in a large river bas in. The modeling platform provided a robust tool to better understand how various surface and subsurface water sources coexist and interact in space and time within the study area. McDonnell et al., (2010), Kirchner (2006) and many others in the hydrologi c modeling community have argued that there is a need to not only validate the model simulated outputs but also to validate that those outputs are produced for the right reasons. This requires for example, validating that the model is not only simulating t he total streamflow correctly but it is simulating the proportion of surface and groundwater water, and the age of water forming the streamflow correctly. It has been long recognized that the traditional model evaluation criteria (i.e. comparison of stream flow and groundwater elevation) do not contain enough information to fully validate whether the model simulated outputs are occurring for the right reasons. Tracer information and additional water chemistry measurements can provide extra insights to valida te models; however, field experiments and data collection can be expensive which limits their use. As such there is a need for modeling studies introducing and demonstrating the

PAGE 184

184 application of additional criteria that are inexpensive and facilitate compreh ensive testing of hydrologic models and their underlying conceptual models. A significant contribution of this study is that it demonstrated the potential use of groundwater flow fractions estimations obtained using high resolution river water specific co nductivity measurements as an additional model validation criterion. Recent developments in sensor technology allowed relatively inexpensive and high frequency measurement of specific conductivity of water at multiple river locations within the SFRB. This new validation criteria helped testing the adequacy of surface and groundwater flow interactions as simulated by fully integrated model over a large river basin. Based on the broader objectives defined in the beginning of this dissertation, the research and its findings can be partitioned into three major research components (RC s ). Major findings and conclusions on each of these major RCs are summarized below. RC1: What are the major water budget components in the basin and how do they vary in space and t ime? How do surface and subsurface water contributions to streamflow vary along the river? ET is the most important water balance component in the basin comprising 7 7 % of rainfall. The remaining water balance components namely, runoff, surface storage, and subsurface storage are governed by seasonal and interannual differences between rainfall and ET patterns over the basin. In general recharge to groundwater occurs primarily in June through September when rainfall significantly exceeds ET. Peaks in av erage streamflow occur both in March, when winter frontal rainfall is relatively high and ET losses are low, and in September at the end of the summer rainy season. ET varies spatially throughout the basin with water table depth and land cover. On averag e 85 % of the ET occurs in the confined region where water tables are shallow (averaging 0.84 m below land surface) and thus ET is primarily energy limited. In the unconfined region water tables are deeper (averaging 12.1 m below land surface) and thus gen erally moisture limited. In both the confined and

PAGE 185

185 the unconfined regions grass lands show the lowest ET losses and forests and wetlands show the highest ET losses. Streamflow in the confined region of the basin is episodic, with more than 95% generated by recent near stream rainfall that temporarily raises the water table above the land surface near stream channels. Virtually no inter event subsurface baseflow from the surficial aquifer occurs in the confined region. In the unconfined region of the basin there is little surface connectivity between the watershed and stream resulting in 52% of the streamflow in the upper portion being generated by subsurface contributions from F loridan aqu ifer that originated as diffuse infiltration through the epikarst. Surface contributions to streamflow are episodic and originate in the upper confined region of the basin. Geologic conditions were found to exert primary control on the spatial variability of streamflow generation processes in this basin through their influence on the balance between rainfall, evapotranspiration, runoff and infiltration processes. Climatic variability was found to provide primary control on the temporal variability of stream flow generation processes, resulting in significant seasonal and interannual variability in both the timing and sources of streamflow. RC2: What Parflow.CLM processes and parameters are ET, groundwater level, streamflow, and subsurface contributions to str eamflow most sensitive to? How does parametric sensitivity, and interaction among sensitive parameters, vary across the basin? The global sensitivity analysis of ParFlow.CLM identified seven parameters (out of 33 parameters included in the analysis) to be most sensitive towards ET, groundwater level, streamflow, and subsurface contributions to streamflow. These parameters include the hydraulic conductivity of the IAS, the hydraulic conductivity and specific storage of the UFAS, maximum LAI for mixed forest s, regions. Of the seven sensitive parameters identified during GSA, all but maximum LAI for mixed forests showed high interactive effects, with the hydraulic conductivity of the I AS being the most interactive parameter. The sensitivity and interaction of the hydraulic conductivity of the IAS underscores the importance of better mapping of the lateral and vertical extent of IAS in confined, transition and Wacasassa Flat regions. The hydrologic response of the confined region of the basin showed highest sensitivity to the hydraulic conductivity of the Intermediate Aquifer System confining unit (IAS). Changes in hydraulic conductivity of the IAS changed the

PAGE 186

186 surficial aquifer water tabl e depth in the confined region which in turn influenced average ET, total streamflow and peak streamflow in this region. After the hydraulic conductivity of the IAS, the maximum LAI of mixed forests in the confined region was the second most influential pa rameter controlling both average ET and total streamflow in the confined region. This illustrates the importance of vegetation properties in determining average hydrologic behavior in this region. Mannings coefficient in the confined region was the second most influential parameter controlling peak streamflow in this region. Total streamflow and subsurface contributions to streamflow in the unconfined region were most sensitive to the hydraulic conductivity of the Upper Floridan Aquifer System (UFAS). Tot al streamflow also showed significant sensitivity to maximum LAI of mixed forest and wilting point in the unconfined region. For subsurface contributions to streamflow in the unconfined region hydraulic conductivity of the conduit system was the second mo st sensitive parameter after hydraulic conductivity of the UFAS. Peak streamflow in the unconfined region was most sensitive to the hydraulic conductivity of the IAS and the mannings coefficient of the confined region, underscoring the dominance of the con fined region characteristics in determining streamflow response to major rainfall events throughout the entire basin. In the unconfined region ET was most sensitive to wilting point, followed by the UFAS hydraulic conductivity and parameters of the van Ge nuchten function, indicating that due to the generally lower water tables ET is moisture limited in this region. Mean groundwater levels and the magnitude of groundwater level fluctuation in the UFAS were gener ally most sensitive to UFAS hydraulic conduct ivity, UFAS specific storage and IAS hydraulic conductivity throughout the basin. The significance of ET parameters in controlling hydrologic behavior throughout the watershed underscores the importance of ET and highlights the significance of feedback b etween land atmosphere, surface and subsurface processes throughout the basin. This finding also highlights the significance of adequate representation of ET processes within the modeling framework; something which was consistently missing in previous regi onal groundwater models for the region. Analysis of hydrologic response (e.g. streamflow, groundwater elevation) of the model over 340 GSA model runs indicated that parameter values that produce reliable hydrologic response in one of region of the basin ma y not reproduce reliable predictions elsewhere. As such future sensitivity/uncertainty analysis studies should incorporate spatial variability in parameter values within major hydro geologic regions of the basin.

PAGE 187

187 RC3: How do streamflow sources, flow p aths and travel times vary along the river and between base flow and storm events? What processes and parameters exert primary control on travel time distributions? Streamflow in the confined region is episodic and travel time distributions are temporally variable. However Parflow Slim simulations indicate that all rainfall within the previous 30 days and having a median age of less than 7 days. Channel floodplain connectivi ty in the confined region, which depends on rainfall spatial patterns, antecedent hydrologic conditions and topography, controls the source area for streamflow, total volume of streamflow, and the shape of the travel time distributions. In general simulat ions indicate that even during the extreme 2004 hurricane season most confined region streamflow originates as rainfall within 15 00 m of the river. Travel time distributions in the unconfined region show a wide distribution of ages with two distinct comp Aquifer System contribution originating primarily from rainfall falling in the unconfined region, infiltrating throu gh the epikarst, and making its way through the porous matrix of the UFAS and higher permeability conduits to the river. the shape of the confined region travel time distrib ution and is well fit by a log normal travel time function with parameters that vary slightly between storms and with hydrograph position. invariant and well fit by a gamma functi on with parameters that indicate fractal behavior. The shape of the travel time distribution is sensitive to conduit geometry, the contrast between conduit and UFAS hydraulic conductivity, and the geometry of the model predicted extent of the unconfined c ontributing area. unconfined region varies with hydrologic condition, ranging from 100% old water during baseflow conditions to 10 30% old water at the peak flow major stor m hydrographs. The relative weight of new and old water travel time distributions is sensitive to IAS hydraulic conductivity which governs the total streamflow generated in the confined region. During baseflow conditions the median age of streamflow in th e unconfined region is predicted to be approximately 15 17 years, in close agreement to age previously determined by Chlorofluorocarbon and isotopic measurements taken at springs within the basin. During high streamflow events the median streamflow age dro ps rapidly to approximately 15 17 days at the peak of the hydrograph.

PAGE 188

188 The fully integrated ParFlow.CLM model developed during this study represents a baseline scenario that provides a strong basis to quantitatively predict the impacts of major changes in l anduse (e.g. wetlands and forests cleared for agricultural or urban development), groundwater extraction and climate on stores, fluxes, flowpaths and travel times of water in the SFRB. These predictions are essential to inform holistic land and water reso urce planning in the region in order to provide reliable water supply for human uses as well as environmental flows that are protective of aquatic ecosystems. Furthermore accurate predictions of surface versus groundwater contributions to streamflow will a llow prediction of transport and transformations of ecologically relevant solutes such as carbon, nitrogen and phosphorus in the spring and river systems in the region. The ParFlow.CLM model for the SFRB presented in this dissertation was developed using model parameters obtained from various state and federal agencies. However the Global Sensitivity Analysis conducted as a part of this research provides a foundation to design quantit ative spatially distributed optimal parameter estimation techniques (such as the Ensemble Kalman Filter) to reduce the uncertainty of model predictions using high frequency streamflow, groundwater level and end member mixing model estimates. In particul ar seven parameters namely the hydraulic conductivity of the IAS, the hydraulic conductivity and specific storage of the UFAS, maximum LAI for mixed regions were identified to be most sensitive towards ET, streamflow and groundwater

PAGE 189

189 level predictions and require more accurate estimation for future models developed for the basin Boundary conditions applied at the four lateral boundaries of the basin in the current model were assumed to be time invariant constant head conditions which allow Although these lateral boundaries were located far from the river, and outside the surface watershed, in o rder to minimize their effects on surface groundwater interactions, future modeling efforts should evaluate the effects of extending the lateral boundaries (especially the south boundary of the basin), and/or the adoption of time varying boundary condi tions, on model predictions. Similarly the effect of model grid resolution on simulated hydrologic behavior should be further investigated. The version of ParFlow used in this study required a uniform grid This, together with the large size of the st udy domain, required a relatively coarse grid resolution to maintain computational efficiency. The ParFlow model now has the capability to discretize domains using variable grid sizes. Thus the ParFlow.CLM predictions obtained in this study could be compa red with a model that incorporates higher resolution grid cells in hydrologically dynamic regions near the land surface, stream networks and high conductivity zones Finally, the ParFlow model developed here assumed that overland and streamflow would be w ell predicted using the kinematic wave equation and that conduits in the system could be represented by high conductivity porous media networks. Future work should test the sensitivity of model predictions to alternative overland flow physics and alternat ive conduit flow physics and geometry.

PAGE 190

190 APPENDIX CHAPTER 4 AD D ITIONAL FIGURES A B C D Figure A 1 Water particle median age versus flow for base case model variant s A) N o high hydraulic conductivity conduit hydraulic conduc tivity igh er hydraulic conductivity conduits hydraulic conductivity Blue line with black circles present data for the rising limb and green lines with red circles present data for the recession period; value corresponding to storm peak is represented by orange circle. Also shown are the median water ages for the three drought periods (dotted circle on upper left corner of each plot)

PAGE 191

191 Figure A 2. Temporal dynamics of vario us sources of water at 2500 during various time s within a storm hydro graph for base case and its variants Variants are hydraulic conductivity

PAGE 192

192 A B Figure A 3 Water particle age dis tribution at station 2500 at various points on storm A) N o conduit base case B) H igh er hydraulic conductivity conduit base case

PAGE 193

193 A B C Water age (days) Figure A 4 Spatial extent of water particle source location during peak of storm event b A) B asecase B) B ase case with high er hydraulic conductivity conduits C) B ase case with no conduit. The color scheme shows the age (days) of all the b

PAGE 194

194 A B C Water age (days) Figure A 5 . A) B asecase B) B ase case with high er hydraulic conductivity conduits. C ) B ase case with no conduit. The color scheme shows the age (days) of all the

PAGE 195

195 LIST OF REFERENCES Anderson M.P., Woessner W.W. 1992. Applied groundwater modeling: simulation of flow and advective transport. Academic Press: California; 38 1. Ando K., Kostner A., Neuman S.P. 2003. Stochastic continuum modeling of flow and transport in a crystalline rock mass: Fanay Augeres, France, revisited. Hydrogeology Journal 11 : 521 35. Arthur J.D., Baker A.E., Cichon J.R., Wood H.A.R., Rudin A. 2005. Florida, Aquifer aquifer systems. Florida Geological Survey ; 67 : 1 149. Ashby S.F., Falgout R.D. 1996. A parallel multigrid preconditioned conjugate gradient algorithm for gro undwater flow simulations. Nuclear Science and Engineering 124 : 145 59. Bailly Comte V., Martin J.B., Screaton E. 2011. Time variant cross correlation to assess residence time of water and impli cation for hydraulics of a sink rise karst system. Water Reso urces Research 47 : W05547. DOI: doi:10.1029/2010WR009613. Bailly Comte V., Martin J.B., Jourde H., Screaton E.J., Pistre S., Langston A. 2010. Water exchange and pressure transfer between conduits and matrix and their influence on hydrodynamics of two kars t aquifers with sinking streams. Journal of Hydrology 386 : 55 66. Banks E.W., Simmons C.T., Love A.J., Shand P. 2011. Assessing spatial and temporal connectivity between surface water and groundwater in a regional catchment: Implications for regional scal e water quantity and quality. Journal of Hydrology 404 : 30 49. Basu N.B., Jindal P., Schilling K.E., Wolter C.F., Takle E.S. 2012. Evaluation of analytical and numerical approaches for the estimation of groundwater travel time distribution. Journal of Hyd rology 475 : 65 73. Beven K.J. 2010. Preferential flows and travel time distributions: defining adequate hypothesis tests for hydrological process models. Hydrological Processes 24 : 1537 47. Beven K., Binley A. 1992. The future of distributed models: mode l calibration and uncertainty prediction. Hydrological Processes 6 : 279 98. Beven K., Freer J. 2001. Equifinality, data assimilation, and uncertainty estimation in mechanistic modelling of complex environmental systems using the GLUE methodology. Journal of hydrology 249 : 11 29.

PAGE 196

196 Beven K. 1993. Prophecy, reality and uncertainty in distributed hydrological modelling. Advances in Water Resources 16 : 41 51. DOI: 10.1016/0309 1708(93)90028 E. Birkel C., Soulsby C., Tetzlaff D. Dunn S., Spezia L. 2012. High f requency storm event isotope sampling revea ls time variant transit time distributions and influence of diurnal cycles. Hydrological Processes 26 : 308 16. Bolger B.L., Park Y., Unger A.J.A., Sudicky E.A. 2011 Simulating the pre development hydrologic cond itions in the San Joaquin Valley, California. Journal of Hydrology 411 : 322 30. Botter G., Peratoner F., Putti M., Zuliani A., Zonta R., Rinaldo A., Marani M. 2008. Observation and modeling of catchment scale solute transport in the hydrologic response: A tracer study. Water Resources Research 44 : W05409. Botter G. 2012. Catchment mixing processes and travel time distributions. Water Resources Research 48 : W05545. DOI: 10.1029/2011WR011160. Botter G., Bertuzzo E., Rinaldo A. 2010. Transport in the hydro logic response: Travel time distributions, soil moisture dynamics, and the old water paradox. Water Resources Research 46 : W03514. DOI: 10.1029/2009WR008371. Butscher C., Huggenberger P. 2009. Modeling the temporal variability of karst groundwater vulnera bility, with implications for climate change. Environmental science & technology 43 : 1665 9. Campolongo F., Cariboni J., Saltelli A. 2007. An effective screening design for sensitivity analysis of large models. Environmental Modelling & Software 22 : 1509 18. Capell R., Tetzlaff D., Malcolm I.A., Hartley A.J., Soulsby C. 2011. Using hydrochemical tracers to conceptualise hydrological function in a larger scale catchment draining contrasting geologic provinces. Journal of Hydrology 408 : 164 77. Chow V.T., Maidment D.R., Mays L.W. 1988. Applied hydrology McGraw Hill: New York; 57 6 Christophersen N., Hooper R.P. 1992. Multivariate analysis of stream water chemical data: The use of principal components analysis for the end member mixing problem. Water Resou rces Research 28 : 99 107. Cosgrove B.A., Lohmann D., Mitchell K.E., Houser P.R., Wood E.F., Schaake J.C., Robock A., Marshall C., Sheffield J., Duan Q. 2003. Real time and retrospective forcing in the North American Land Data Assimilation System (NLDAS) p roject. Journal of Geophysical Research 108 : 8842.

PAGE 197

197 Cox M.H., Su G.W., Constantz J. 2007. Heat, chloride, and specific conductance as ground water tracers near streams. Ground Water 45 : 187 95. Dai Y., Zeng X., Dickinson R.E., Baker I., Bonan G.B., Bosilo vich M.G., Denning A.S., Dirmeyer P.A., Houser P.R., Niu G. 2003. The common land model. Bulletin of the American Meteorological Society 84 : 1013 23. Darracq A., Destouni G., Persson K., Prieto C., Jarsj J. 2010. Quantification of advective solute travel times and mass transport through hydrological catchments. Environmental fluid mechanics 10 : 103 20. de Rooij R., Graham W., Maxwell R.M. 2012. A particle tracking scheme for simulating pathlines in coupled surface subsurface flows. Advances in Water Reso urces 52 : 7 18. Desmarais K., Rojstaczer S. 2002. Inferring source waters from measurements of carbonate spring response to storms. Journal of Hydrology 260 : 118 34. Du J., Qian L., Rui H., Zuo T., Zheng D., Xu Y., Xu C. 2012. Assessing the effects of ur banization on annual runoff and flood events using an integrated hydrological modelling system for Qinhuai River basin, China. Journal of Hydrology 464 465 : 127 139. Dunn S., Birkel C., Tetzlaff D., Soulsby C. 2010. Transit time distributions of a concept ual model: their characteristics and sensitivities. Hydrological Processes 24 : 1719 29. Dunn S., Freer J., Weiler M., Kirkby M., Seibert J., Quinn P., Lischeid G., Tetzlaff D., Soulsby C. 2008. Conceptualization in catchment modelling: simply learning? Hy drological Processes 22 : 2389 93. Duriancik L.F., Bucks D., Dobrowolski J.P., Drewes T., Eckles S.D., Jolley L., Kellogg R.L., Lund D., Makuch J.R., O'Neill M.P. 2008. The first five years of the Conservation Effects Assessment Project. Journal of Soil an d Water Conservation 63 : 185A 97A. Fenicia F., McDonnell J.J., Savenije H.H. 2008a. Learning from model improvement: On the contribution of complementary data to process understanding. Water Resources Research 44 : W06419. Fenicia F., Savenije H.H., Matge n P., Pfister L. 2008b. Understanding catchment behavior through stepwise model concept improvement. Water Resources Research 44 : W01402.

PAGE 198

198 Ferguson I.M., Maxwell R.M. 2010. Role of groundwater in watershed response and land surface feedbacks under climate change. Water Resources Research 46 : W00F02. Florea L.J., Vacher H.L. 2007. Eogenetic karst hydrology: Insights from the 2004 hurricanes, peninsular Florida. Ground Water 45 : 439 46. Foglia L., Hill M., Mehl S., Burlando P. 2009. Sensitivity analysis, ca libration, and testing of a distributed hydrological model using error based weighting and one objective function. Water Resources Research 45 : W06427. DOI: 10.1029/2008WR007255. Francone C., Cassardo C., Richiardone R., Confalonieri R. 2012. Sensitivity Analysis and Investigation of the Behaviour of the UTOPIA Land Surface Process Model: A Case Study for Vineyards in Northern Italy. Boundary Layer Meteorology : 1 12. Freeze R. A., Harlan R. L. 1969. Blueprint for a physically based disgitally simulated hy drologic response model. Journal of Hydrology 9: 237 258. Frhlich H.L., Breuer L., Frede H., Huisman J.A., Vach K.B. 2008. Water source characterization through spatiotemporal patterns of major, minor and trace element stream concentrations in a compl ex, mesoscale German catchment. Hydrological Processes 22 : 2028 43. Furman A. 2008. Modeling coupled surface subsurface flow processes: a review. Vadose Zone J 7: 741 756. Goderniaux P., Brouyre S., Fowler H.J., Blenkinsop S., Therrien R., Orban P., Da ssargues A. 2009. Large scale surface subsurface hydrological model to assess climate change impacts on groundwater reserves. Journal of Hydrology 373 : 122 38. Graham W.D., McLaughlin D.B. 1991. A stochastic model of solute transport in groundwater: Appli cation to the Borden, Ontario, Tracer Test. Water Resources Research 27 : 1345 59. Graham W., McLaughlin D. 1989. Stochastic analysis of nonstationary subsurface solute transport: 2. Conditional moments Water Resources Research 25 : 2331 55. Gulley J.D. Martin J.B., Spellman P., Moore P.J., Screaton E.J. 2012. Influence of partial confinement and Holocene river formation on groundwater flow and dissolution in the Florida carbonate platform. Hydrological Processes : DOI:10.1002/hyp.9601. DOI: 10.1002/hyp .9601.

PAGE 199

199 Gupta H.V., Sorooshian S., Yapo P.O. 1999. Status of automatic calibration for hydrologic models: Comparison with multilevel expert calibration. Journal of Hydrologic Engineering 4 : 135 43. Hooper R.P., Christophersen N., Peters N.E. 1990. Modellin g streamwater chemistry as a mixture of soilwater end members An application to the Panola Mountain catchment, Georgia, USA. Journal of Hydrology 116 : 321 43. Hrachowitz M., Soulsby C., Tetzlaff D., Malcolm I., Schoups G. 2010a. Gamma distribution models for transit time estimation in catchments: Physical interpretation of parameters and implications for time variant transit time assessment. Water Resources Research 46 : W10536. Hrachowitz M., Soulsby C., Tetzlaff D., Speed M. 2010b. Catchment transit time s and landscape controls does scale matter? Hydrological Processes 24 : 117 25. Huntington J.L., Niswonger R.G. 2012. Role of surface water and groundwater interactions on projected summertime streamflow in snow dominated regions: An integrated modeling ap proach. Water Resources Research 48 : W11524. Jones J.E., Woodward C.S. 2001. Newton Krylov multigrid solvers for large scale, highly heterogeneous, variably saturated flow problems. Advances in Water Resources 24 : 763 74. Jones J.P., Sudicky E.A., McLare n R.G. 2008. Application of a fully integrated surface subsurface flow model at the watershed scale: A case study. Water Resources Research 44 : W03407. DOI: 10.1029/2006WR005603. Katz B. 2004. Sources of nitrate contamination and age of water in large kar stic springs of Florida. Environmental Geology 46 : 689 706. Katz B.G., Bhlke J.K., Hornsby H.D. 2001. Timescales for nitrate contamination of spring waters, northern Florida, USA. Chemical Geology 179 : 167 86. Katz B.G., Griffin D.W. 2008. Using chemica l and microbiological indicators to track the impacts from the land application of treated municipal wastewater and other sources on groundwater quality in a karstic springs basin. Environmental Geology 55 : 801 21. Kirchner J.W. 2006. Getting the right an swers for the right reasons: Linking measurements, analyses, and models to advance the science of hydrology. Water Resources Research 42 : W03S04. Kirchner J.W., Feng X., Neal C. 2001. Catchment scale advection and dispersion as a mechanism for fractal sca ling in stream tracer concentrations. Journal of Hydrology 254 : 82 101.

PAGE 200

200 Kirchner J.W., Feng X., Neal C. 2000. Fractal stream chemistry and its implications for contaminant transport in catchments. Nature 403 : 524 7. Kirchner J.W., Feng X., Neal C., Robso n A.J. 2004. The fine structure of water quality dynamics: The (high frequency) wave of the future. Hydrological Processes 18 : 1353 9. Koeniger P., Leibundgut C., Stichler W. 2009. Spatial and temporal characterisation of stable isotopes in river water as indicators of groundwater contribution and Isotopes in environmental and health studies 45 : 289 302. Kollet S., Cvijanovic I., Schttemeyer D., Maxwell R., Moene A., Bayer P. 2009. The influence of rain sensible heat, sub surface heat convection and the lower temperature boundary condition on the energy balance at the land surface. Vadose Zone J 8 : 846 57. Kollet S.J., Maxwell R.M. 2008a. Capturing the influence of groundwater dynamics on land surface processes using an in tegrated, distributed watershed model. Water Resources Research 44 : W02402. Kollet S.J., Maxwell R.M. 2008b. Demonstrating fractal scaling of baseflow residence time distributions using a fully coupled groundwater and land surface model. Geophysical Resea rch Letters 35 : L07402. DOI: 10.1029/2008GL033215. Kollet S.J., Maxwell R.M. 2006. Integrated surface groundwater flow modeling: A free surface overland flow boundary condition in a parallel groundwater flow model. Advances in Water Resources 29 : 945 58. Langston A., Screaton E., Martin J., Bailly Comte V. 2012. Interactions of diffuse and focused allogenic recharge in an eogenetic karst aquifer (Florida, USA). Hydrogeology Journal 20 : 767 81. ess of hydrologic and hydroclimatic model validation. Water Resources Research 35 : 233 41. Li Q., Unger A.J.A., Sudicky E.A., Kassenaar D., Wexler E.J., Shikaze S. 2008. Simulating the multi seasonal response of a large scale watershed wi th a 3D physically based hydrologic model. Journal of Hydrology 357 : 317 36. Linhoss A., Muoz Carpena R., Kiker G., Hughes D. 2012. Hydrologic Modeling, Uncertainty, and Sensitivity in the Okavango Basin: Insights for Scenario Assessment: Case Study. Jo urnal of Hydrologic Engineering DOI: 10.1061/(ASCE)HE.1943 5584.0000755.

PAGE 201

201 Loague K., Heppner C.S., Abrams R.H., Carr A.E., VanderKwaak J.E., Ebel B.A. 2005. Further testing of the Integrated Hydrology Model (InHM): event based simulations for a small rang eland catchment located near Chickasha, Oklahoma. Hydrological Processes 19 : 1373 98. with the aid of environmental tracers: 1. Models and their applicability. Journal of hydrology 57 : 207 31. Maresch W., Walbridge M.R., Kugler D. 2008. Enhancing conservation on agricultural landscapes: A new direction for the Conservation Effects Assessment Project. Journal of Soil and Water Conservation 63 : 198A 203A. Matsubayashi U., Velasquez G.T., Takagi F. 1993. Hydrograph separation and flow analysis by specific electrical conductance of water. Journal of Hydrology 152 : 179 99. Maxwell R.M., Miller N.L. 2005. Development of a coupled land surface and groundwater model. Journal of Hydrometeorology 6 : 233 47. Maxwell R.M., Welty C., Tompson A.F. 2003. Streamline based simulation of virus transport resulting from long term artificial recharge in a heterogeneous aquifer. Advances in Water Resources 26 : 1075 96. Maxwell R., Kastenberg W. 1999. Stochastic environmental risk analysis: an integrated methodology for predicting cancer risk from contaminated groundwater. Stochastic Environmental Research and Risk Assessment 13 : 27 47. McDonnell J., Sivapalan M., Vach K., Dunn S., Grant G., Haggerty R., Hinz C., Hooper R., Kirchner J., Roderick M. 2007. Moving beyond heterogeneity and process complexity: A new vision for watershed hydrology. Water Resources Research 43 : W07301. DOI: 10.1029/2006WR005467. McGuire K.J., McDonnell J.J. 2006. A review and evaluation of catchment transit time modeling. Journal of Hydrology 330 : 543 63. McGuire K., McDonnell J., Weiler M., Kendall C., McGlynn B., Welker J., Seibert J. 2005. The role of topography on catchment scale water residence time. Water Resources Research 41 : W05002. DOI: 10.1029/2004WR003657. McGuire K., Weiler M., McDonnell J. 2007. Integrating tracer experiments with modeling to assess runoff processes and water transit times. Advances in Water Resou rces 30 : 824 37.

PAGE 202

202 McMillan H.K., Clark M.P., Bowden W.B., Duncan M., Woods R.A. 2011. Hydrological field data from a modeller's perspective: Part 1. Diagnostic tests for model structure. Hydrological Processes 25 : 511 22. McMillan H., Tetzlaff D., Clark M ., Soulsby C. 2012. Do time variable tracers aid the evaluation of hydrological model structure? A multimodel approach Water Resources Research 48 : W05501. DOI: 10.1029/2011WR011688. Meals D.W., Dressing S.A., Davenport T.E. 2010. Lag time in water quality response to best managem ent practices: A review. Journal of environmental quality 39 : 85 96. Meyer B.A., Kincaid T.R., Hazlett T.J. 2008. Modeling Karstic Controls on Watershed Scale Groundwater Flow in the Floridan Aquifer of North Florida. In Sinkholes and the Engineering and Environmental Impacts of Karst (GSP 183) Yuhr LB, Alexander EC, Beck BF (eds). ASCE : VA; 351 361. Moore P.J., Martin J.B., Screaton E.J. 2009. Geochemical and statistical evidence of recharge, mixing, and controls on spring discharge in an eogenetic karst aquifer. Journal of Hydrology 376 : 443 55. Morris M.D. 1991. Factorial sampling plans for preliminary computational experiments. Technometrics 33 : 161 74. Mueller Warrant G., Griffith S., Whittaker G., Banowetz G., Pfender W., Garcia T., Giannico G. 201 2. Impact of land use patterns and agricultural practices on water quality in the Calapooia River Basin of western Oregon. Journal of Soil and Water Conservation 67 : 183 201. Muoz Carpena R., Zajac Z., Kuo Y. 2007. Global sensitivity and uncertainty anal yses of the water quality model VFSMOD W. Trans.ASABE 50 : 1719 32. Nakamura R. 1971. Runoff analysis by electrical conductance of water. Journal of Hydrology 14 : 197 212. Nash J., Sutcliffe J. 1970. River flow forecasting through conceptual models part I A discussion of principles. Journal of hydrology 10 : 282 90. Nossent J., Bauwens W. 2012. Multi variable sensitivity and identifiability analysis for a complex environmental model in view of integrated water quantity and water quality modeling. Water sci ence and technology 65 : 539 49. 18 Journal of Hydrology 359 : 303 12.

PAGE 203

203 Osmond D., Meals D., Hoag D., Arabi M., Luloff A. Jennings G., McFarland M., Spooner J., Sharpley A., Line D. 2012. Improving conservation practices programming to protect water quality in agricultural watersheds: Lessons learned from the National Institute of Food and Agriculture Conservation Effects A ssessment Project. Journal of Soil and Water Conservation 67 : 122A 7A. Panagopoulos G., Lambrakis N. 2006. The contribution of time series analysis to the study of the hydrodynamic characteristics of the karst systems: Application on two typical karst aqu ifers of Greece (Trifilia, Almyros Crete). Journal of Hydrology 329 : 368 76. Pilgrim D.H., Huff D.D., Steele T.D. 1979. Use of specific conductance and contact time relations for separating flow components in storm runoff. Water Resources Research 15 : 32 9 39. Pollock D.W. 1988. Semianalytical computation of path lines for finite difference models. Ground Water 26 : 743 50. Rihani J.F., Maxwell R.M., Chow F.K. 2010. Coupling groundwater and land surface processes: Idealized simulations to identify effects of terrain and subsurface heterogeneity on land surface energy fluxes. Water Resources Research 46 : W12523. DOI: 10.1029/2010WR009111. Rinaldo A., Beven K.J., Bertuzzo E., Nicotina L., Davies J., Fiori A., Russo D., Botter G. 2011. Catchment travel time distributions and water flow in soils. Water Resources Research 47 : W07537. DOI: 10.1029/2011WR010478. Ritorto M., Screaton E., Martin J., Moore P. 2009. Relative importance and chemical effects of diffuse and focused recharge in an eogenetic karst aquifer: an example from the unconfined upper Floridan aquifer, USA. Hydrogeology J ournal 17 : 1687 98. Roa Garca M., Weiler M. 2010. Integrated response and transit time distributions of watersheds by combining hydrograph separation and long term transit time modeling. Hydrology and Earth System Sciences 14 : 1537 49. Running S.W., Lov eland T.R., Pierce L.L. 1994. A vegetation classification logic based on remote sensing for use in global biogeochemical models. Ambio 23 : 77 81. Saltelli A., Ratto M., Tarantola S., Campolongo F. 2005. Sensitivity analysis for chemical models. Chemical R eviews Columbus 105 : 2811 28. Saltelli A. 1999. Sensitivity analysis: Could better methods be used? Journal of Geophysical Research: Atmospheres (1984 2012) 104 : 3789 93.

PAGE 204

204 Saltelli A., Tarantola S., Campolongo F., Ratto M. 2004. Sensitivity analysis in pr actice: a guide to assessing scientific models Wiley; 219. Schneider J.W., Upchurch S.B., Chen J., Cain C. 2008. Simulation of groundwater flow in north central Florida and south central Georgia. Suwannee River Water Management District. : Live Oak; 92. S eibert J., McDonnell J.J. 2002. On the dialog between experimentalist and modeler in catchment hydrology: Use of soft data for multicriteria model calibration. Water Resources Research 38 : 1241. DOI: 10.1029/2001WR000978. Sobol I.M. 1993. Sensitivity esti mates for nonlinear mathematical models. Mathematical Modelling and Computational Experiments 1 : 407 14. Soulsby C., Dunn S.M. 2003. Towards integrating tracer studies in conceptual rainfall runoff models: recent insights from a sub arctic catchment in th e Cairngorm Mountains, Scotland. Hydrological Processes 17 : 403 16. Soulsby C., Tetzlaff D., Dunn S., Waldron S. 2006. Scaling up and out in runoff process understanding: insights from nested experimental catchment studies. Hydrological Processes 20 : 246 1 5. Speed M., Tetzlaff D., Soulsby C., Hrachowitz M., Waldron S. 2010. Isotopic and geochemical tracers reveal similarities in transit times in contrasting mesoscale catchments. Hydrological Processes 24 : 1211 24. Srivastava V., Graham W., Maxwell R. 201 3a. Geologic and climatic controls on streamflow generation processes in a complex eogenetic karst basin. Under preparation. Srivastava V., Graham W., Carpena R., Maxwell R. 2013b. Global sensitivity analysis of an integrated hydrologic model insights on g eologic and vegetative controls over the hydrologic functioning of a large complex basin. Under preparation Sudicky E.A., Jones J.P., Park Y.J., Brookfield A.E., Colautti D. 2008. Simulating complex flow and transport dynamics in an integrated surface sub surface modeling framework. Geosciences Journal 12 : 107 22. Tetzlaff D., Soulsby C., Hrachowitz M., Speed M. 2011. Relative influence of upland and lowland headwaters on the isotope hydrology and transit times of larger catchments. Journal of Hydrology 40 0 : 438 47. Tripathi N. 2006 Spatial water balance of the Suwannee River Basin, Florida: an application of El Nio Southern Oscillation Climatology University of Florida: Gainesville; 174.

PAGE 205

205 Tsang Y., Tsang C., Hale F., Dverstorp B. 1996. Tracer transport in a stochastic continuum model of fractured media. Water Resources Research 32 : 3077 92. Upchurch S.B. 2007. An introduction to the Cody Escarpment, North Central Florida. Suwannee River Water Manage ment District : Live Oak, Florida ; 16. Upchurch S., Chen J ., Cain C. 2008. Springsheds of the Santa Fe River Basin. Alachua county Department of Enviremental Protection, Gainesville, Fl. Springsheds of the Santa Fe River Basin. Alachua County (Florida) Departm ent of Environmental Protection: Gainesville, FL ; 93 Vacher H., Mylroie J.E. 2002. Eogenetic karst from the perspective of an equivalent porous medium. Carbonates and Evaporites 17 : 182 96. Van der Velde Y., De Rooij G., Rozemeijer J., Van Geer F., Broers H. 2010. Nitrate response of a lowland catchment: O n the relation between stream concentration and travel time distribution dynamics. Water Resources Research 46 : W11534. DOI: 10.1029/2010WR009105. Van Genuchten M.T. 1980. A closed form equation for predicting the hydraulic conductivity of unsaturated soi ls. Soil Science Society of America Journal 44 : 892 8. Van Griensven A., Meixner T., Grunwald S., Bishop T., Diluzio M., Srinivasan R. 2006. A global sensitivity analysis tool for the parameters of multi variable catchment models Journal of Hydrology 324 : 10 23. VanderKwaak J.E., Loague K. 2001. Hydrologic Response simulations for the R 5 catchment with a comprehensive physics based model. Water Resources Research 37 : 999 1013. Water resources associates 2007. MFL establishment for the upper Santa Fe River. Live Oak; 294. Werner A.D., Gallagher M.R., Weeks S.W. 2006. Regional scale, fully coupled modelling of stream aquifer interaction in a tropical catchment. Journal of Hydrology 328 : 497 510. Yang J. 2011. Convergence and uncertainty analyses in Monte Carlo based sensitivity analysis. Environmental Modelling & Software 26 : 444 57. Yang J., Liu Y., Yang W., Chen Y. 2012. Multi objective sensitivity analysis of a fully distributed hydrologic model WetS pa. Water Resources Management : 1 20.

PAGE 206

206 BIOGRAPHICAL SKETCH Vibhava Srivastava was born and grew up in the city of Allahabad, Utt ar Pradesh, India. He earned a Bachelor of Technology in a gricultural e ngineering from the Allahabad Agricultural Institute D eemed University in 2003. In January 2004 he joined the University of Arkansas from where he received a Master of Science in b iological and a gricultural e ngineering in August 2006. He joined the University of Florida in August 2006 and received a Doctor of Philosophy in August 2013. He is currently working as Post Doctorate researcher at the Water Institute, University of Florida.