Open Collections

UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

A computer analysis of the flow of water and nutrients in agricultural soils as affected by subsurface… Richard, Paul François 1988

Your browser doesn't seem to have a PDF viewer, please download the PDF to view this item.

Item Metadata

Download

Media
831-UBC_1988_A1 R52.pdf [ 7.71MB ]
Metadata
JSON: 831-1.0098263.json
JSON-LD: 831-1.0098263-ld.json
RDF/XML (Pretty): 831-1.0098263-rdf.xml
RDF/JSON: 831-1.0098263-rdf.json
Turtle: 831-1.0098263-turtle.txt
N-Triples: 831-1.0098263-rdf-ntriples.txt
Original Record: 831-1.0098263-source.json
Full Text
831-1.0098263-fulltext.txt
Citation
831-1.0098263.ris

Full Text

A COMPUTER ANALYSIS OF THE FLOW OF WATER AND NUTRIENTS IN AGRICULTURAL SOILS AS AFFECTED BY SUBSURFACE DRAINAGE by PAUL FRANCOIS RICHARD B.Sc. (Agr. Eng.), M c G i l l U n i v e r s i t y , 1978 M.Sc, M c G i l l U n i v e r s i t y , 1981 A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF DOCTOR OF PHILOSOPHY in THE FACULTY OF GRADUATE STUDIES Resource Management Science We accept t h i s t h e s i s as conforming to the r e q u i r e d standard THE UNIVERSITY OF BRITISH COLUMBIA August 1988 ® Paul Frangois Richard, 1988 In presenting this thesis in partial fulfilment of the requirements for an advanced degree at the University of British Columbia, I agree that the Library shall make it freely available for reference and study. I further agree that permission for extensive copying of this thesis for scholarly purposes may be granted by the head of my department or by his or her representatives. It is understood that copying or publication of this thesis for financial gain shall not be allowed without my written permission. Department The University of British Columbia 1956 Main Mall Vancouver, Canada V6T 1Y3 DE-6 (3 /81 ) ABSTRACT A computer model was developed in order to determine the effects of drainage practices on nutrient losses from level agricultural soils. The model performs a daily simulation of the vertical flow of water, nitrogen, phosphorus, and heat, and of the growth of crops. A water f.low submodel calculates the depth of the water table based on daily predictions of evaporation, transpiration, flow to drains and ditches, and deep percolation. An original saturated-unsaturated flow algorithm is used to determine moisture infiltration, redistribution, and upward flow in the soil matrix, as well as bypassing flow in the soil macropores and horizontal flux between the soil matrix and the macropores, and surface runoff. Nutrient movement occurs by mass flow. Heat flow, nutrient biochemical transformations, and crop growth are determined by using well established relations. Field tests were carried out for a period of two years on an experimental site in the Lower Fraser Valley of British Columbia. The water table depth was measured on a continuous basis. Grab samples of drainwater and observation wells were obtained periodically and analyzed for nitrogen (N0 3 -N, NH4-N, and TKN) and phosphorus (P0 4 -P and TP). The field results show a decrease in the concentration of all nutrients over the sampling period, and provide evidence that denitrif ication and bypassing flow are important mechanisms affecting the nutrient balance of this soi l . These results were used to calibrate the model. An excellent - i i -f i t o f the o b s e r v e d water t a b l e p r o f i l e and an adequate f i t o f t h e o b s e r v e d d r a i n c o n c e n t r a t i o n o f n i t r a t e were o b t a i n e d . The s i m u l a t i o n r e v e a l e d t h a t b y p a s s i n g f l o w i s a v e r y i m p o r t a n t t r a n s f e r mechanism i n t h i s s o i l and must be i n c l u d e d i n o r d e r t o o b t a i n a s a t i s f a c t o r y f i t o f the e x p e r i m e n t a l d a t a . A s e n s i t i v i t y a n a l y s i s o f the model showed t h a t the p a t t e r n s o f m o i s t u r e f l o w have a p r e d o m i n a n t i n f l u e n c e on t h e r a t e o f n u t r i e n t l e a c h i n g . In p a r t i c u l a r , i t was found t h a t the n u t r i e n t c o n c e n t r a t i o n i n d r a i n water i s a s t r o n g f u n c t i o n o f t h e h y d r a u l i c c o n d u c t i v i t y o f t h e s o i l m a t r i x and of t h e h o r i z o n t a l d i s t a n c e between the s o i l macropores, which c o n t r o l the r a t i o o f m o i s t u r e f l o w i n t h e s o i l m a t r i x t o t h e macropore f l o w and the l a t e r a l d i f f u s i o n o f n u t r i e n t s between the s o i l m a t r i x and the macropores. The e f f e c t s o f f o u r d i f f e r e n t d r a i n a g e d e s i g n s on n u t r i e n t l o s s e s were s i m u l a t e d o v e r a p e r i o d o f two y e a r s f o r t h r e e d i f f e r e n t s o i l s and two d i f f e r e n t n u t r i e n t d i s t r i b u t i o n s i n t h e s o i l . I t was found t h a t t h e r e i s a l a r g e d i f f e r e n c e between t h e amount o f n u t r i e n t s l e a c h e d from d r a i n a g e systems u s i n g d i f f e r e n t d r a i n a g e c o e f f i c i e n t s . There was a l s o a l a r g e d i f f e r e n c e i n the r e s p o n s e o f two d r a i n a g e d e s i g n s b a s e d on t h e same d r a i n a g e c o e f f i c i e n t b u t u s i n g d i f f e r e n t d e p t h and s p a c i n g o f d r a i n s . T r a n s i e n t e f f e c t s , as d e t e r m i n e d by t h e i n i t i a l v e r t i c a l d i s t r i b u t i o n o f the n u t r i e n t s , were seen t o remain dominant o v e r the two y e a r d u r a t i o n o f the s i m u l a t i o n . The model was found t o be u s e f u l i n e x p l a i n i n g the a p p a r e n t c o n t r a d i c t i o n s found i n t h e l i t e r a t u r e a s s e s s i n g the e f f e c t s o f s u b s u r f a c e d r a i n a g e on n u t r i e n t l o s s e s . The r e s u l t s from t h e model show t h e s e e f f e c t s t o be s t r o n g l y s i t e and c o n d i t i o n s p e c i f i c . F u r t h e r m o r e , t h e model shows t h a t s o i l s and d r a i n a g e d e s i g n s t h a t produce s i m i l a r volumes o f d r a i n f l o w may e x h i b i t v e r y d i f f e r e n t l e a c h i n g r e s p o n s e s , and t h a t d r a i n a g e d e s i g n s e q u i v a l e n t f r o m a h y d r a u l i c s t a n d p o i n t can be v e r y d i s s i m i l a r i n t h e i r p o t e n t i a l f o r l e a c h i n g n u t r i e n t s . The model p r o v i d e s a t o o l which can be used t o d e t e r m i n e t h e a p p r o p r i a t e n e s s o f d i f f e r e n t d r a i n a g e d e s i g n s i n s o i l s where m i n i m i z i n g n u t r i e n t l o s s e s i s c r i t i c a l . - i v -TABLE OF CONTENTS Page ABSTRACT i i TABLE OF CONTENTS v LIST OF TABLES v i i LIST OF FIGURES v i i i LIST OF SYMBOLS x ACKNOWLEDGMENTS x i i i INTRODUCTION 1 CHAPTER 1 - AGRICULTURAL DRAINAGE AND WATER QUALITY: A REVIEW 1.1 - Introduction 4 1.2 - Mechanisms 4 1.3 - Field findings 10 1.4 - Computer simulations 17 1.5 - Conclusions 19 1.6 - References 22 CHAPTER 2 - MODELLING WATER AND NUTRIENT FLUXES IN DRAINED AGRICULTURAL SOILS 2.1 - Introduction 2.2 - Modelling approach 2.3 - Climate 2.4 - Water balance 2.5 - Nutrient balance 2.6 - Crop growth 2.7 - Discussion 2.8 - References CHAPTER 3 - A FIELD STUDY OF WATER AND NUTRIENT MOVEMENT IN AN AGRICULTURAL SOIL UNDER DRAINED AND UNDRAINED CONDITIONS 3.1 - Introduction 79 3.2 - Site and methods 81 3.3 - Results 85 3.4 - Discussion 94 3.5 - Conclusions 101 3.6 - References 102 29 30 36 40 52 62 67 74 - v -CHAPTER 4 - A COMPARISON OF SIMULATED AND OBSERVED WATER AND NUTRIENT MOVEMENT IN DRAINED AGRICULTURAL SOILS 4.1 - Introduction 105 4.2 - Methods 107 4.3 - Results 113 4.4 - Discussion 128 4.5 - Conclusions 138 4.6 - References 141 CHAPTER 5 - A SIMULATED COMPARISON OF THE EFFECTS OF DIFFERENT DRAINAGE PRACTICES ON NUTRIENT LOSSES 5.1 - Introduction 143 5.2 - Methods 145 5.3 - Results 150 5.4 - Discussion 158 5.5 - Conclusions 173 5.6 - References 175 CONCLUDING REMARKS 178 APPENDIX I - An example of moisture flow and 180 redistribution simulation APPENDIX II - A comparison of the dynamic 184 equilibrium water movement model with a numerical solution of Richard's equation APPENDIX III - Diffusions between the macropore 192 and micropore zones APPENDIX IV - A comparison of the model with 194 an analytical solution - vi -LIST OF TABLES Page TABLE 1-1 Composition of drain discharge and surface 11 runoff TABLE 1-2 Nutrient losses from different drainage 14 practices TABLE 2-1 List of the major assumptions used in 32 the model TABLE 3-1 Some physical and chemical characteristics 82 of the soil at the experimental site TABLE 3-2 Tillage and fertilization of test plots 84 TABLE 3-3 Average and range of sample concentrations 92 TABLE 3-4 Comparison between nutrient content of 93 drainwater and well water in drained and undrained fields, using the Wilcoxon matched pairs test TABLE 4-1 Summary of simulation parameters 112 TABLE 4-2 Effects of parameter variation on simulation 118 results, expressed as a ratio for 1984-85 TABLE 4-3 Effects of parameter variation on simulation 119 results, expressed as a ratio for 1984 only TABLE 4-4 Effects of using a simple bypassing flow 127 model or a simple homogeneous flow model on the simulation results, expressed as a ratio. TABLE 5-1 Summary of the conditions simulated 149 TABLE 5-2 Effect of drainage design and soil type on 153 water flux TABLE 5-3 Average losses from surface flow for the 154 drainage designs TABLE 5-4 Relative drainage losses using distribution 156 C, expressed as a ratio of the average TABLE 5-5 Relative drainage losses using distribution 157 T, expressed as a ratio of the average TABLE 5-6 Effect of drainage design and transfer 159 mechanisms on biochemical transformations TABLE 5-7 Average overall nutrient losses from 169 the different designs - vii -LIST OF FIGURES Page FIGURE 2-1 Schematic diagram of the model layout 33 FIGURE 2-2 Flowchart of the computer program 35 FIGURE 2-3 Flowchart of the subprogram CLIMATE 37 FIGURE 2-4 Flowchart of the subprogram WATER 41 FIGURE 2-5 Typical hydraulic conductivity 47 versus moisture content curve FIGURE 2-6 Diagram of the dynamic equilibrium 49 moisture content concept FIGURE 2-7 Flowchart of the subprogram NUTRIENT 54 FIGURE 2-8 Diagram of the nutrient balance 55 FIGURE 2-9 Flowchart of the subprogram CROP 63 FIGURE 3-1 Experimental field layout 83 FIGURE 3-2 Precipitation (PREC) and water table 86 depth below the surface versus time for drained and undrained plots FIGURE 3-3 Nitrate nitrogen concentration versus time 88 FIGURE 3-4 Ammonium nitrogen concentration versus time 89 FIGURE 3-5 Phosphate phosphorus concentration versus 90 time FIGURE 3-6 Estimated nitrate nitrogen drain loss versus 100 time FIGURE 4-1 Moisture characteristics of the soil used in 110 the simulation FIGURE 4-2 Simulated and observed water table profile of 114 the drained plot FIGURE 4-3 Simulated and observed water table profile of 115 the undrained plot FIGURE 4-4 Simulated and observed nitrate concentration 116 of drainage water - v i i i -FIGURE 4-5 The effects of varying the ped conductivity 120 CS on the simulated water table profile FIGURE 4-6 The effects of varying the macropore porosity 121 8PR on the simulated water table profile FIGURE 4-7 The effects of varying the equivalent ped 122 spacing PEDS on the simulated water table profile FIGURE 4-8 Simulating water table movement using a 125 simple bypassing flow model or a simple homogeneous flow model FIGURE 4-9 Simulating the nitrate content of drainage 126 water using a simple bypassing flow model or a simple homogeneous flow model FIGURE 5-1 Simulated water table profiles resulting 151 from design A (subdrains at 1.1m depth and 14m spacing), and design B (subdrains at 1.1m depth and 20m spacing). FIGURE 5-2 Simulated water table profiles resulting 152 from design C (subdrains at 0.7m depth and 10m spacing) design 0 (no subdrains, 1.1m deep ditches at 50m spacing) FIGURE A-l Example of moisture flow and redistribution 181 in the profile FIGURE A-2 Comparison of the model with a numerical 187 solution of an infiltration event FIGURE A-3 Comparison of the model with a numerical 188 solution of redistribution following infiItration FIGURE A-4 Comparison of infiltration and redistri- 190 bution for a soil with a larger hydraulic conductivity FIGURE A-5 Comparison of the model output with an 196 analytical solution - ix -LIST OF SYMBOLS A: nutrient mass (tonne/ha) AB: absorbtivity AK: adsorption coefficient (ha«m/tonne) ALOSS: nutrient transfer mass (tonne/ha*day) ALFA: permeability exponent ATU: accumulated thermal units (day-C) AVOL: volatilized mass (tonne/ha*day) BETA: temperature coefficient BEVAP: Priestley-Taylor coefficient BR: storm erosivity (tonne/ha«day) C: unsaturated hydraulic conductivity (m/day) CA: above drain conductivity (m/day) CB: below drain conductivity (m/day) CC: volumetric heat capacity (J/m3 C) CEQ: dissolved nutrient mass (tonne/ha»m) CM: adsorption capacity (tonne/ha*m) CMA: adsorbed nutrient mass (tonne/ha»m) CMF: soil cover factor CS: saturated ped hydraulic conductivity (m/day) CSAV: average hydraulic conductivity (m/day) DAY: day count DEC: declination (rad) DEQ: equivalent depth (m) DGRO: day count since growth onset DNIT: denitrification rate (tonne/ha-day) DPT: layer depth (m) DRF: drain flow (m/day) DRS: drain lateral spacing (m) DSEED: onset of growth (day) DT: time increment (day) DV: maximum volatilization depth (m) DWT: depth to water table (m) DZ: layer thickness (m) ED: dryness coefficient EI: storm energy (tonne/ha»day) EM: wetness coefficient ENR: nutrient enrichment ratio EPF: erosion control practices factor ES: process specific coefficient ET: temperature coefficient F,u: carbon to nitrogen conversion factor FIELDL: slope length (m) FLOODF: root decay factor (day-C)"1 FMIN: soil mineral fraction (m/m) GAM: slope factor exponent GRAD: pressure gradient (m/m) GROW: dry matter increment (tonne/ha) H: water table height above drains (m) HC: thermal conductivity (J/m.C.day) - x -HO: water level in collector ditch (m) HSS: day length (rad) K: decay coefficient (day"1) L: latitude (rad) LAI: leaf area index MC02: carbon mineralization rate (tonne/ha.day) PART: init ial growth partition coefficient PAV: root available phosphorus (tonne/ha) PEDS: lateral equivalent ped spacing (m) PET: potential evapotranspiration (m/day) PINF: potential infiltration (m/day) PLANTN: plant nutrient concentration PNMIN: minimum plant nutrient concentration PNOPT: optimum plant nutrient concentration POTRAN: potential transpiration (m/day) PP: growth partition coefficient Q: moisture transfer between layers (m/day) QSD: lateral infiltration (m/day) RAD: solar radiation evaporative power (m/day) RADIN: net solar radiation (m/day) RAIN: precipitation (m/day) RFF: surface runoff (m/day) ROOT: root mass (tonne/ha) RTDI: root depth increment (m) RTDMX: maximum root depth (m) S: moisture suction (m) SEF: soil erosion factor SLF: slope and length factor SLOF: root decay rate (tonne/ha.day) SLP: field slope SNM: snowmelt rate (m/day) SOILL: soil loss rate (tonne/ha.day) SOLIOL: adsorbed nutrients loss rate (tonne/ha.day) STRES: crop stress factor SUN: sunshine duration (day) SV: evaporation factor TAV: average temperature (C) TK: absolute temperature (K) TMIN: minimum temperature (C) TMAX: maximum temperature (C) TOTS: mass of adsorbed nutrient (tonne/ha) TPCON: growth conversion factor (tonne/ha«m) TRIG: solar geometry factor TRN: transpiration (m/day) TS: soil temperature (C) TSM: total soil mass (tonne/ha) TU: thermal units (day-C) TUREP: accumulated thermal units at onset of reproductive growth (day-C) 9: volumetric moisture content (mobile water) 8AV: availability factor 8DR: moisture content following unsaturated drainage - xi -9EQ: moisture content at static equilibrium 9FL: moisture content following infiltration 9PR: moisture content of intrapedal zone 6QD: moisture content of dynamic equilibrium 9S: saturated moisture content 9UP: moisture content following upward flux - xii -ACKNOWLEDGMENTS The author wishes to express his sincere gratitude towards the people involved in this study. Dr. L. Smith, Dr. M. Novak, Dr. J . de Vries, Dr. L. Lavkulich, and Dr. N. Nagpal provided guidance throughout th is work, along with Dr. S .T . Chieng, whose encouragement and gentle prodding are much appreciated. The participation and discussion of fellow students L. Elder and E. Charter was very helpful in all phases of this work, and their contribution is gratefully acknowledged. This study would not have been feasible without the generous financial support of the Natural Sciences and Engineering Research Council. The deciphering and typing skil ls of A. Jones were also very helpful. I also wish to express my deepest gratitude to my friends and family, whose encouragement and support was much needed. - x i i i -1 INTRODUCTION Subsurface drainage is becoming increasingly widespread as a practice to remove excess water and improve crop growing conditions and f ie ld accessibility in humid areas. The effects of this practice on the hydraulic regimes of the soils are well documented. However, the overall effects of subsurface drainage on nutrient losses from field soils are s t i l l poorly understood. This makes it diff icult to assess the impact of this drainage practice on the environment, where nutrient losses may cause eutrophication problems, and on soil fert i l i ty status. This situation arises because several factors compound the overall results consequent to the use of this form of drainage. The mechanisms of nutrient losses are altered by subsurface drainage, because the proportions of surface flow to percolation tend to change as a result. Percolation losses typically consist of soluble forms of nutrients, while surface losses are composed of both soluble, adsorbed, and particulate forms. As a result, not only the magnitude of nutrient loss, but also the relative proportion of the different species, are affected by drainage practices. This situation is complicated by the fact that the leaching mechanisms are dissimilar from soil to so i l , as hydraulic conductivity, soil adsorption characteristics, and the presence of cracks and macropores vary and may even be affected by the drainage practice i t se l f . Climatic conditions also make site to site comparisons d i f f icu l t . Finally, changes in cropping practices 2 usually follow drainage improvements, with the result that changes in fer t i l izat ion, soil cover, and other variables can mask the effects that are specific to drainage. As a consequence, it is necessary to find a method where the mechanisms of nutrient losses specific to drainage can be examined in isolation. A computer simulation model of nutrients and water flow in the soils can be such a method, provided that the model is sufficiently comprehensive so as not to exclude any relevant transfer mechanism, and provided that these mechanisms can be simulated in a realistic fashion. Such a model, then, makes it possible to determine the relative importance of these mechanisms, the conditions under which they are important, and if there are designs of drainage systems that can minimize these losses. Consequently, a computer model was developed in order to address the above questions. The validity of this model was assessed using experimental results from a field site located in the Lower Fraser Valley of British Columbia. The model was used to compare various combinations of drainage designs and leaching mechanisms on nutrient losses in the context of this f ie ld . This thesis describes the result of this process, outlined in five chapters presented in paper format. Each chapter addresses a specific aspect of the problem, as outlined below. The introductory chapter presents a review of the literature relevant to this problem in order to set the context of this work. Its purpose is also to examine whether a consensus can be found amongst the various researchers as to the effects of drainage, and 3 whether there currently exists simulation models that can be adapted to the present task. The second chapter describes the computer model that was developed as part of this study. Its salient features and their relevance to the problem at hand are described and compared with other models when appropriate. The third chapter is a report of the experimental procedure used to calibrate and test the model. An experimental field with both drained and undrained plots was monitored for two years. Water table depth was read on a continuous basis and grab samples of drain and well water were collected and analyzed for nitrogen and phosphorus. In the fourth chapter, field and model results are compared. The quality of the f i t after calibration and the importance of the transfer mechanisms are discussed with the help of a sensitivity analysis. The final chapter presents the results obtained from comparing various drainage designs using soils characterized by different transfer mechanisms, including flow in macropores (bypassing flow). These results are compared with the literature to determine their validity. Their limitations and implications for drainage design are discussed. 4 CHAPTER 1 Agricultural Drainage and Water Quality: A Review 1.1 Introduction In humid areas, surface or subsurface drainage is often required to remove excess water from agricultural lands. At present there is l i t t le agreement with respect to the effect of the various types of drainage on soil fert i l i ty and water pollution by nutrients. For instance, Walter et a l . (1979), in their review of conservation practices, concluded that the only well-known effect of subsurface drainage was to increase nitrate losses while Switzer-Howse (1983), reviewing practices that could improve water quality in the Canadian Great Lakes basin, stated that the effect of subsurface drainage on nutrient movement and water quality is s t i l l unknown. The purpose of this chapter is to review the literature pertaining to drainage and nutrient movement in order to identify and discuss f ie ld studies where the effects of drainage on nutrients can be ascertained. It also discusses mechanisms via which drainage alters nutrient movement, and computer models that incorporate these mechanisms. 1.2 Mechanisms Agricultural drainage may consist of surface drainage, in which fields are profiled to favor runoff of ponded water, and subsurface drainage, which lowers the water table by means of deep 5 ditches or subdrains. Often land use is intensified to justify the investment in drainage. This in tens i f ica t ion is usually accompanied by increased nutrient pollution. For instance, Coote et a l . (1982) showed that subsurface drainage is positively correlated with erosion and nitrate losses in several Ontario watersheds. In this case, the correlation between subsurface drainage and row crops masks the specific effects due to drainage alone. Consequently, specific mechanisms must be inspected. Drainage affects nutrient movement and transformations in the so i l , particularly in the root zone, via four distinct mechanisms related to water movement. First, surface runoff can be altered by drainage; surface drainage promotes runoff, while subsurface drainage can reduce it under specific circumstances. Secondly, subsurface drainage provides an outlet for percolating water and may reduce the residence time of this water in the soi l . A third mechanism, distinct from the second, resides in the fact that subsurface drainage provides an outlet for bypassing flow if cracks and channels extend down from topsoil to drain depth. Because water flow involves solute flow and erosion from surface flow, any alteration of water flow patterns will affect nutrient movements. Finally, subsurface drainage reduces the occurrence of flooding and high water tables and promotes soil aeration. This fourth mechanism favors an oxidized environment in the so i l , which results in increased mineralization, nitr if ication, and phosphate sorption, and in decreased d e n i t r i f i c a t i o n , as well as increased transpiration and nutrient uptake for most farm crops. 6 Surface drainage shapes the land so as to eliminate surface storage and promotes surface flow via shallow ditches. Since the erosion rate is considered to increase with slope and, for individual storms, with runoff volume and flow rate (Williams et a l . , 1984), surface drainage can be expected to result in higher nutrient losses. Subsurface drainage can reduce surface runoff by increasing the soil storage capacity for infiltration water. It is well-known (Henninger et a l . , 1976; Edwards and Amerman, 1984) that naturally well-drained soils do not produce as much surface runoff as poorly drained soils. Subsurface drainage of poorly drained soils can produce the same effect. This has been shown to be the case by simulation (Whiteley and Ghate, 1979). In their review, Irwin and Whiteley (1983) found that there is compelling evidence for this mechanism from field studies. This was also seen by Istok and Kling (1983), who also found that subsurface drainage was effective in reducing erosion. Baden and Eggelsmann (1968) observed a reduction of surface runoff from t i le drained organic soi ls; they found that the storage volume was further augmented by increased transpiration. Subsurface drainage can also increase i n f i l t r a t i o n by improving the hydraulic properties of the soi l . The increase in infiltration rate following mole drainage has been ascribed to the soil cracks caused by moling (Robinson and Beven, 1983; Rennes et a l . , 1976). Newman and Robinson (1983) found that subsurface drainage reduced the incidence of surface poaching, while Hundal et 7 a l . (1976) and Yadav and Leyton (1960) found that drainage resulted in an increase of the soil hydraulic conductivity. Breckenridge et a l . (1985) demonstrated that subsurface drainage can increase infiltration rate by venting air trapped behind the wetting front. Since subsurface drainage increases percolation, it also increases leaching; there is ample evidence for this in the literature (for instance, in the reviews of Baker and Johnson, 1976, or Hore and Broughton, 1976). Thus the concentration of soluble species such as nitrates in drain effluent tend to reflect the concentration of this species adsorbed above the drains. In contrast, strongly adsorbed species, such as orthophosphates, generally have a drainwater concentration at equilibrium with the adsorbed concentration of the soil at drain depth (Sharpley et a l . , 1977; Baker et a l . , 1978). It is very difficult to determine the effect of subdrains on deep percolation and groundwater pollution. Gilliam et a l . (1979) found from a water balance that there is less deep percolation where subdrains are used. This would indicate that subsurface drainage may lessen deep percolation leaching to groundwater. Devitt et a l . (1976) and Cooke and Williams (1970) found lower nitrate levels below drains than above. MacGregor et a l . (1974) found similar results in t i le drained clay loam, in contrast to deep nitrate penetration in an untiled clay loam. Bypassing flow is now recognized as an important mechanism for the transfer of water and nutrients in soils (see the review of White, 1985). Because of cracks, root channels, and other 8 macropores, infiltrating water can rapidly percolate down without displacing the water held in the soil matrix, or without saturating this matrix. If such bypassing channels extend down to drain l eve l , subdrains can provide a ready outlet for such flow. Armstrong and Arrowsmith (1986) proposed that drainwater consists of the independent additions of bypassing flow and matrix flow; thus drainwater nutrient concentration should ref lect the composition of both sources. Several studies found that the orthophosphate concentration of drainwater tends to increase with flow rate, following fertilization and rain (Sharpley and Syers, 1979) or waste water application (Karlen et a l . , 1976; Hergert et a l . , 1981). Conversely, in similar situations (Sharpley and Syers, 1981; Hergert et a l . , 1981) it was found that nitrate levels decrease with increasing drain flow rates, indicating a dilution effect. This arises because matrix water is generally richer in nitrates and poorer in phosphates than bypassing water, whether the latter originates from irrigation or rainfal l . This situation is not confined to soils with pronounced shrinkage cracks; Coles and Trudgil l (1985) found evidence of bypassing flow in weakly structured soi ls , while Kanchanasut and Scotter (1982) demonstrated the contribution of root channels to bypassing flow. Bypassing flow also varies seasonally (Reid and Parkinson, 1984), being more prevalent during summer because of shrinkage cracks. Disruption of the soil structure by tillage would tend to increase matrix flow at the expense of bypassing flow, which may contribute to the variation in time of leaching patterns. Accordingly Bergstrom 9 (1987) and Roberts et a l . (1986) found that tillage was a more important factor than ferti l ization in causing high nitrate levels in drainwater. Lowering the water table, by subsurface drainage or other means, promotes soil aeration, which affects biochemical reactions in the so i l . High water tables slow the oxidation of soil organic matter, which accumulates in poorly drained soils. Draining such soils accelerates mineralization, which releases nutrients to soil water, particularly in the case of organic soils (Miller, 1979). While proceeding at a slower rate, oxidation of organic matter under submerged conditions depletes soil oxygen and increases denitrification. Thus poorly drained soils leach few nitrates, both because of reduced flow rate and denitrification (Gambrell et a l . , 1975). By contrast, aerated conditions in well drained soils favor nitrate formation from nitr i f ication and oxidation. For instance, Gast et a l . (1974) found that nitrate concentration increases as the lateral distance to subdrains decreases. Meek et a l . (1970) found that denitrification was complete below the water table in soil columns; in a subsequent paper (Willardson et a l . , 1972) they showed that submerging drain outlets was effective in reducing nitrate loadings. Flooding and draining cycles can also increase phosphorus leaching, i f flooding conditions sufficiently reduce the redox potential of the so i l . This arises from the release of phosphates adsorbed or complexed by iron or calcium compounds solubilized by the reducing conditions (Ponnamperuma, 1972). This was found to be the case when drained 10 fields are flooded for rice cultivation (Johnston et a l . , 1965), to control subsidence (Reddy, 1983), or for experimental purposes (Hill and Shawney, 1981). Drainage also indirectly influences nutrient transformations by promoting vigorous root growth, nutrient uptake, and nitrogen fixation in legumes (Dowding, 1981). Drainwater can also contain appreciable amounts of sediment, and sediment-adsorbed phosphates (Sharpley and Syers, 1979; Bottcher et a l . , 1981; Culley et a l . , 1983). It is not clear at present whether such sediments originate from the vicinity of the drains or are carried in bypassing flow, as suggested by Bottcher et a l . (1981). 1.3 Field Findings Several researchers have analyzed nutrient concentration from surface runoff and drain water. Table 1-1 summarizes these results. It can be seen that, in general, drainwater is more concentrated in nitrates, but less concentrated in other nutrients, than surface runoff. It is clear from these results that the composition of drainwater is quite different from that of surface runoff; thus, any alteration of the volume ratio of surface runoff and subdrain flow would markedly affect nutrient loading to surface waters. A few studies have been designed such that direct comparisons between different drainage regimes are possible. Schwab et a l . (1980) reported on several years of monitoring a silty clay site in Ohio which is drained by deep (lm) drains, shallow (0.5m) drains or 11 Table 1-1 Composition of Drain Discharge and Surface Runoff Site Description Parameters, mg/1 Type TN N03-N TP P04-•P SS Ref Florida, spodic sand » shallow tillage S 0.32 0.19 Rogers et a l . deep t i l lage 1 D 5.40 0.40 1976 S 0.71 0.25 D 1.89 0.13 Georgia, loamy sand, S 0.34 Jackson et a l . corn D 9.40 1973 Georgia, loamy sand, S 0.47 Hubbard and varied crops D 8.75 Sheridan, 1983 Iowa, s i l t loams, S 2.4 1.8 0.49 12350 Baker et a l . varied crops2 D 12.2 12.1 0.12 — 1978 Iowa, s i l t loams, S 11.0 1.01 0.013 1060 Hanway and corn3 D 18.0 0.028 0.004 — Laflen, 1974 Iowa, loess, S 34.3 2.0 0.82 0.29 20750 Burwel1 et a l . corn and forage D 13.7 13.5 0.06 0.06 — 1977 Ireland, gleyed S 1.7 0.7 Burke et a l . , clay pasture D 2.8 0.08 1974 Louisiana, clay S 0.67 0.17 146 Bengtson et a l . loam, corn D 0.96 0.015 — 1984 New York, s i l t loam, S 4.1 3.6 0.19 Klausner et a l . varied crops2 D 19.9 19.6 0.01 1974; Zwerman et a l . , 1972 New York, s i l t loam, S 18 2.5 . Hergert et a l . , manured corn4 D 44 0.78 1981 New Zealand, s i l t S 3.4 0.3 2.81 1.27 1780 Sharpley and loam, forage D 8.7 6.6 0.27 0.13 270 Syers, 1979, 19; Ohio, silty clay S 18.1 5.0 1236 Schwab et a l . loam, tillage D 15.7 0.6 880 1973 tillage no-ti l l S 14.5 1.8 277 D 8.6 0.5 177 Ontario, clay soi ls, S 0.19 0.077 71 Culley et a l . , various crops D 0.23 0.12 90 1983 Ontario, clay loams, S 5.7 1.72 Phillips et a l . corn 2 , 4 D 14.6 0.06 1981 Vermont, s i l t loams, S various crops D 1.1 5.7 0.70 0.02 Benoit, 1973 12 Notes: TN: t o t a l n i t r o g e n TP: t o t a l phosphorus SS: suspended s o l i d s S: s u r f a c e r u n o f f D: d r a i n d i s c h a r g e ^ e e p t i l l a g e impeded s u b s u r f a c e d r a i n a g e 2TN i n c l u d e s NH4-N and N03-N o n l y 3TN i n c l u d e s NH4-N, N02-N, and N03-N o n l y ; SS i n c l u d e s d i s s o l v e d s o l i d s . S p e c i f i c f l o w e v e n t s d r a i n w a t e r from t h e s e p l o t s . Sediment l o s s e s i n d r a i n w a t e r from t h e s e p l o t s was s i m i l a r t o sediment l o s s e s from the s u r f a c e d r a i n e d p l o t s d u r i n g major s t o r m s , a l t h o u g h t h e r e was no e v i d e n c e t h a t t h i s was a r e s u l t o f b y p a s s i n g f l o w . 13 surface drains. As can be seen in Table 1-2, the surface drained plots had the least nitrate losses and the most phosphorus and sediment losses. Phosphorus losses were higher in the deeper drained f ield because of the larger sediment concentration in drainwater from these plots. Sediment losses in drainwater from these plots were similar to sediment losses from the surface drained plots during major storms, although there was no evidence that this was a result of bypassing flow. Gilliam and Skaggs (1986), in North Carolina, found that subdrainage increases nitrate losses but decreases phosphorus losses. The comparison was between fields with surface drainage, deep ditches 100 m apart, or deep ditches and subsurface tubing. By contrast, Bengtson et a l . (1984) found that all losses were reduced by subsurface drainage. They attributed this reduction to the lowering of the water table in fal l and winter, when the most severe soil erosion takes place on bare, waterlogged soils for their Louisiana conditions. Drainage markedly increased both nitrogen and phosphorus losses in a Florida sandy citrus grove (Rogers et a l . , 1976). All plots under comparison were subdrained, but different ti l lage treatments caused some partial clogging of some plots. More surface runoff was seen in the clogged plots, but this source of nutrient loss was much less than the leaching losses in the well drained plots. The opposite situation was found by Campbell et a l . (1985), for a sandy soil growing potatoes also in Florida. They compared nutrient losses from a field with shallow subdrains to a 14 Table 1-2 Nutrient Losses from Different Drainage Practices Losses, kg/ha-yr. Description N03-N TN TP SS Ref Surface drained Shallow drains Deep drains 5.4 11.2 14.9 _ 1.3 0.8 1.2 887 218 395 Schwab et a l . , 1980 Surface drained Deep ditches Sub drains 3.7 15.7 32.4 13.6 20.0 42.1 0.53 0.33 0.21 - Gilliam and Skaggs, 1986 Surface drained Sub drains -7.0 6.1 3.1 2.1 1088 911 Bengtson et a l . , 1984 Furrow drained Sub drains 4.53 2.75 4.81* 3.12* 1.10+ 0.43 - Campbellt et a l . , 1985 Undrained Mole drains 5.2 15.4 - 2.2 0.6 - Burke et a l . , 1974 Undrained Subdrains Sub & Mole drains 27.4 21.9 18.2 27.7* 22.4* 19.0* -— Armstrong1* et a l . , 1984 * N03-N and NH4-N only + P04-P only t 13 weeks only x 10 months only 15 field with water furrows. All losses were found to be reduced by subsurface drainage. The lower nitrate concentrations in drainflow were unexpected; it was presumed that the high water tables of the ridge and furrow system increased the likelihood of nitrate leaching from the raised beds to the furrows. Unfortunately, the effect of denitrification was not assessed. In a detailed study in England, Harris et a l . (1984) found that a t i l e and mole drainage system significantly increased nitrate loadings, which reached an average of 51.7 kg/ha-yr. This is in contrast to an undrained f ie ld , where nitrate losses were negligible. Ammonium and phosphorus losses were small in both cases. This situation resulted from the large drainflow volume, on a site where surface runoff is normally small. By contrast, Armstrong (1984) found that subdrainage, and particularly mole and t i le drainage, reduced nitrate loadings because of the large decrease in surface runoff. Nitrate concentration remained unusually high in surface runoff from the undrained plot from this grassland field in England. Similarly, Burke et a l . (1974) found that mole drainage increased nitrate losses but reduced phosphorus losses from a gley field in Ireland. Seuna and Kauppi (1980) compared a t i le drained watershed to a watershed drained by open ditches and found that nitrogen losses significantly increased, while phosphate losses were also higher. This was attributed to the fact that drainage affected soil frost patterns in their Finn site. Finally a Czech study (Ferda and Novak, 1976) found that subsurface drainage of peat soi ls , as 16 compared to undrained, decreased phosphate concentrations from 0.18 mg/1 to 0.03 mg/1. Nitrate concentrations in both cases were negligible, possibly due to the organic matter of the so i l . Thus, it can be seen from the summary of Table 1-2 that it is difficult to generalize about the role of subsurface drainage on nutrient pollution. Most cases report an increase in nitrate loading, typical of leaching, and a decrease in phosphate loading, characteristically originating from surface losses. However, a high concentration of phosphorus in drainwater can arise, due to low adsorption, cracks, or high levels of native phosphates. For instance, Thornley and Bos (1985) found that drainwater quality from Ontario pastures was comparable to raw sewage; such cases present a clear example of drainage as an environmentally damaging practice. In contrast, where drainage is effective in reducing surface runoff, its role is beneficial, as seen above. However, drainage may not be effective in reducing runoff from very intense storms, which cause the most erosion. Baker et a l . (1978) found that a single localized storm on a pasture plot caused more runoff and nutrient losses than what was collected over a year from nearby row crop plots. Consequently, Bottcher et a l . (1981) suggested that surface runoff should be controlled where subsurface drainage provides an outlet for excess water. They showed that a field with raised borders to keep surface runoff in substantially reduced sediment and nutrient losses as compared to conventional drainage. 17 1.4 Computer Simulations A large number of computer models have been developed to simulate water and nutrient flow in soils. (See, for instance, the reviews of Tanji [1982], Haan [1982], or DeJong [1981]). However, few of these models were developed to consider the long term simulation of f ie ld conditions; fewer s t i l l are flexible and comprehensive enough to allow for an analysis of the effects of drainage practices on nutrient loading. Several models, such as CNS (Tubbs and Haith, 1981), CREAMS (Knisel, 1980), or particularly EPIC (Williams et a l . , 1984) concentrate on nutrient loadings from erosion while providing a generalised description of soil mechanisms affecting nutrient movement and transformations. Subsurface drainage can be accounted for in these models in the calculation of percolation, in a simplified manner. These models have not been used at present to simulate specifically the effects of subsurface drainage. Some models have been designed specifically to simulate mass transfer of solutes to drain t i les. Van Ommen (1985) developed a solution to show that solute concentration in drain lines can be represented using a concept of water flow cascading from layer to layer, akin to the models used to simulate hydrographic response of watersheds. Bypassing flow is included but several biological mechanisms, as well as surface runoff, are not. This contrasts with the model of Johnsson et a l . (1987), which emphasizes biochemical processes. Their model simulates water and nitrogen flow to drains in a layered soi l . Surface runoff is included in 18 the water balance calculations but is assumed to have a negligible nutrient concentration. Bypassing flow is not simulated. However, neither of these models was used to investigate the effects of drainage on nutrient flow. Duffy et a l . (1975) presented a model to simulate nitrogen movement in t i l e - d r a i n e d s o i l s . This model includes a comprehensive simulation of biochemical processes, but surface runoff and bypassing flow are not considered. Again, an analysis of the effects of varying drainage parameters on the response of the model was not attempted. Skaggs and coworkers modified a water balance model for drainage design (Skaggs, 1978) specifically to compare the effects of various drainage systems on nutrient loss. In this model, excess water is removed by a combination of surface and subsurface drainage. Biochemical processes simulation was not attempted, and soil nutrient was assigned an average concentration. Bypassing flow was not included. Under these conditions, it was found (Deal et a l . , 1986) that a combination of good subsurface drainage and poor surface drainage resulted in increased nitrate and TKN loadings when compared with the reverse combination. The reverse was true for phosphorus loadings, except for organic soils. They also showed that controlling drainage by means of weirs can reduce nitrate loads. In another study (Skaggs et a l . , 1982) it was found that subsurface drainage can substantially reduce surface erosion, even in relatively flat lands. Soil loss from surface runoff was reduced by increasing drain depth and reducing drain spacing. 19 The model of Kanwar and others (Kanwar and Johnson, 1984; Kanwar et a l . , 1983) was developed specif ical ly to simulate nitrogen dynamics of soils with subsurface drainage. This model includes a comprehensive set of water flow mechanisms (including surface runoff and deep percolation, but not bypassing flow) and of nutrient biochemical transformations. Again, it was found that drain depth and spacing had a pronounced effect on nitrate losses in drainwater. Unfortunately, the effects of subsurface drainage on surface runoff were not reported. Whitmore and Addiscott (1986) presented a model where bypassing flow and diffusion between mobile and immobile water are the main transfer mechanisms. Their model also featured a comprehensive simulation of nitrogen biochemical transformations. In this simulation, it was found that the presence or absence of mole drains made only a relatively small difference in the amount of nitrate leached from the soils simulated. 1.5 Conclusions Although there is a large body of published work on the movement of nutrients in soils under varying drainage conditions, l i t t le of this work has focussed specifically on the effects of drainage on nutrient movement. This arises from the difficulty of isolating the effects of drainage from other factors. Differences in crops, soi ls, and climate from site to site mask any effects that could be attributable to drainage. Even in a similar soil and location, a comparison between 20 drained and undrained site may be di f f icul t because drainage affects crop yield and timing of operations. When such comparisons are possib le , however, the evidence from the f i e ld seems contradictory. In some cases, subsurface drainage increases the overall loss of nutrients, both nitrogen and phosphorus, while the opposite can also be seen. Nevertheless, it seems that, in general, subsurface drainage will increase nitrate losses while reducing losses of soil and soil erosion, phosphorus and other species of nitrogen. To what extent this is verified depends on the type and prevalence of the mechanisms through which these nutrients are lost. In particular, the proportions of surface runoff, bypassing flow, and soil throughflow as nutrient carrying mechanisms seem to be of fundamental importance. As it appears diff icult to generalize from field evidence, computer simulation offers a promising avenue to investigate the specific effects of drainage independently from other variables. The predictive ability of models complex enough to encompass all physical, chemical and biological phenomena either is low, or requires such intense calibration as to be very restricted in application. This is particularly true of nitrogen models, as shown in the comparisons performed by de Willigen and Neeteson (1985). Nevertheless, even relatively simple models can be used for realistic comparisons of drainage scenarios. Such is the case, for instance, in the work of Skaggs and coworkers, whose computer-generated comparisons were found to agree with field observations (Gilliam and Skaggs, 1986). 21 At present there is no computer model that includes a comprehensive simulation of each nutrient transfer mechanism outlined in the previous sections. However, none of these mechanisms can be a priori excluded if an overall assessment of the effects of drainage is to be conducted using a computer simulation. Therefore, there is a need for this type of comprehensive computer simulation. The present study addresses this need. 22 1.6 References - Armstrong, A . C , 1984. The hydrology and water quality of a drained clay catchment, in T.P. Bart and D.G. Walling, eds., Catchment experiments in fluvial hydrology. Geo books, Norwich, 153-168. - Armstrong, A.C. and R. Arrowsmith, 1986. Field evidence for a bi-porous soil water regime in clay soils. Agric. Water Manag., 11, 117-125. - Baden, W. and R. Eggelsmann, 1968. The hydrologic budget of highbogs in the Atlantic region. Proc. 3rd. Int. Peat Congress, Quebec, 206-11. - Baker, J . L . , H.P. Johnson, M.A. Boncherding, and V.R. Payne, 1978. Nutrient and pesticide movement from field to stream, in R.C. Loehr, D.A. Haith, M.F. Walter and C S . Martin, eds. Best management practices for agriculture and silviculture. Ann Arbor Science, 213-246. - Baker, J . L . and H.P. Johnson, 1976. Impact of subsurface drainage on water quality. 3rd Int. Drainage Symp., ASAE, St. Joseph, 91-98. - Bengtson, R.L., C E . Carter, H.F. Morris, and J.G. Kowalczuk, 1984. Reducing water pollution with subsurface drainage. Trans. ASAE, 27, 80-83. - Benoit, G.R., 1973. Effect of agricultural management of wet sloping soil on nitrate and phosphorus in surface and subsurface water. Water Resour. Res. 9, 1296-1303. - Bergstrom, L., 1987. Nitrate leaching from annual and perennial crops in t i le drained plots and lysimeters. J . Environ. Qual., 16, 11-18. - Bottcher, A.B. , E.J. Monke, and L.F. Huggins, 1981. Nutrient and sediment loadings from a subsurface drainage system. Trans. ASAE, 24, 1221-1226. - Breckenridge, R.P., A.R. Jarrett, and J.R. Hoover, 1985. Runoff reduction by venting with shallow subsurface drainage. Trans. ASAE, 28, 476-479. - Burke, W., J . Mulqueen and P. Butler, 1974. Fertilizer losses in drainage water from a surface water gley. Ir. J . Agric. Res. 13, 203-214. - Burwell, R.E. , G.E. Schuman, H.G. Heinemann, and R.G. Spomer, 1977. Nitrogen and phosphorus movement from agricultural 23 watershed. J . Soil Water Cons. 32, 226-230. - Campbell, K.L. , J .S . Rogers and D.R. Hensol, 1985. Drainage water quality from potato production. Trans. ASAE, 28, 1798-1801. - Coles, N. and S. Trudgi l l , 1985. The movement of nitrate f e r t i l i z e r from the so i l surface to drainage waters by preferential flow in weakly structured soils. A g r i c , Ecosystems and Environ., 13, 241-259. - Cooke, G.W. and R.J.B. Williams, 1970. Losses of nitrogen and phosphorus from agricultural land. Water Treat. Exam. 19, 253-276. - Coote, D.R., E.M. MacDonald, W.T. Dickinson, R.C. Ostry and R. Frank, 1982. Agriculture and water quality in the Canadian Great Lakes basin. J . Environ. Qual. 11, 473-481. - Culley, J . L . B . , E.F. Bolton and V. Bernyk, 1983. Suspended solids and phosphorus loads from a clay so i l . J . Environ. Qual. 12, 493-503. - Deal, S .C . , J.W. Gilliam, R.W. Skaggs and K.D. Konyha, 1986. Prediction of nitrogen and phosphorus losses as related to agricultural drainage system design. Agric. Ecosyst. Environ. 18, 37-51. - De Jong, R., 1981. Soil water models: a review. Agriculture Canada Research Branch, LRRI no. 123. - Devitt, D., J . Letey, L.J . Lund, and J.W. Blair, 1976. Nitrate-nitrogen movement through soil as affected by soil profile characteristics. J . Environ. Qual. 5, 282-288. - deWilligen, P. and J . J . Neeteson, 1985. Comparison of six simulation models for the nitrogen cycle in the so i l . Fertilizer Res., 8, 157-171. - Dowding, P. 1981. The effect of art i f icial drainage on nitrogen cycle processes, in F.E. Clark and T. Rosswall, eds. Terrestrial nitrogen cycles. Ecol. Bull. (Stockholm) 33, 615-625. - Duffy, J . , C. Chung, C. Boast and M. Franklin, 1975. A simulation model of biophysiochemical transformations of nitrogen in tile-drained Corn Belt so i l . J . Environ. Qual., 4, 477-486. - Edwards, W.M., and CR. Amerman, 1984. Subsoil chracteristics influence hydrologic response to no-tillage. Trans. ASAE, 27, 1055-1058. 24 Ferda, J . and M. Novak, 1976. The effect of ameliorative measures on the changes of the qual i ty of surface and groundwaters in peat so i ls . Proc. 5th Int. Peat Congress, Poznan, 118-127. Gambre l l , R . P . , J.W. G i l l i a m and S . B . Weed, 1975 . Denitrification in subsoils of the North Carolina coastal plain as affected by soil drainage. J . Environ. Qual. 11, 311-316. Gast, R.G., W.W. Nelson, and J.M. MacGregor, 1974. Nitrate and chloride accumulation and distribution in fertil ized tile-drained soils. J . Environ. Qual. 3, 209-213. Gilliam, J.W. and R.W. Skaggs, 1986. Controlled agricultural drainage to maintain water quality. J . Irrig. Drain. Eng., 112, 254-263. Gil l iam, J.W., R.W. Skaggs, and S.B. Weed, 1979. Drainage control to diminish nitrate loss from agricultural f ields. J . Environ. Qual. 8, 137-142. Haan, C.T., ed., 1982. Hydrologic modeling of small watersheds. ASAE monograph no. 5, ASAE, St. Joseph. Hanway, J . J . and J.M. Laflen, 1974. Plant nutrient losses from ti le outlet terraces. J . Environ. Qual. 3, 351-356. Harris, G.L. , M.J. Goss, R.J. Dowdell, K.R. Howse and P. Morgan, 1984. A study of mole drainage with simplified cultivation for autumn-sown crops on a clay so i l . J . Agric. Sci . Camb., 102, 561-581. Henninger, D.L., G.W. Petersen and E.T. Engman, 1976. Surface soil moisture relationships to surface runoff. Soil Sc. Soc. Am. J . , 40, 773-776. Hergert, G.W., D.R. Bouldin, S.D. Klausner and P.J. Zwerman, 1981. Phosphorus concentration - water flow interactions in t i le effluent from manured land. J . Environ. Qual. 10, 338-344. Hi l l , D.E. and B.L. Shawney, 1981. Removal of phosphorus from wastewater by soil under aerobic and anaerobic conditions. J . Environ. Qual. 10, 401-405. Hore, F.R. and R.S. Broughton, 1976. Overview of environmental effects of agricultural drainage in Canada. Environ. Asp. Irr. Dr., ASCE Symp., Ottawa, 138-154. Hubbard, R.K., and J.M. Sheridan, 1983. Water and nitrate-nitrogen losses from a small, upland, coastal plain watershed. J . Environ. Qual. 12, 291-295. 25 Hundal, S.S. , G.O. Schwab and G.S. Taylor, 1976. Drainage system effects on physical properties of a Lakebed clay so i l . Soil Sc. Soc. Am. J . , 40, 300-305. Irwin, R.W. and H.R. Whiteley, 1983. Effects of land drainage on stream flow. Can. Water Resources J . , 8, 88-103. Istok, J.D. and G.F. Kling, 1983. Effect of subsurface drainage on runoff and sediment yield from an agricultural watershed in western Oregon. J . Hydrology, 65, 274-291. Jackson, W.A., L.E. Asmussen, E.W. Hauser, and A.W. White, 1973. Nitrate in surface and subsurface flow from a small agricultural watershed. J . Environ. Qual. 2, 480-482. Johnsson, H., L. Bergstrom, P.E. Jansson and K. Paustian, 1987. Simulated nitrogen dynamics and losses in a layered agricultural so i l . Agric. Ecosys. Environ., 18, 333-356. Johnston, W.R., F. Ittahadieh, R.M. Daum and A.F. Pillsbury, 1965. Nitrogen and phosphorus in t i le drainage effluent. Soil Sc. Soc. Am. Proc. 29, 287-290. Kanchanasut, P. and D.R. Scotter, 1982. Leaching patterns in soil under pasture and crop. Aust. J . Soil Res. 20, 193-202. Kanwar, R.S., and H.P. Johnson, 1984. Sensitivity analysis of a drainage simulation model. Water Resources Bul l . , 20, 339-342. Kanwar, R.S., H.P. Johnson and J .L . Baker, 1983. Comparison of simulated and measured nitrate losses in t i le effluent. Trans. ASAE, 26, 1451-1457. Karlen, D.L., M.L. Vitosh, and R.J. Kunze, 1976. Irrigation of corn with simulated municipal sewage effluent. J . Environ. Qual., 5, 269-273. Klausner, S.D., P.J. Zwerman and D.F. E l l i s , 1974. Surface runoff losses of soluble nitrogen and phosphorus under two systems of soil management. J . Environ. Qual. 3 42-46. Knisel, W.D., ed . , 1980. CREAMS: a f ie ld-scale model for chemicals, runoff and erosion from agricultural management systems. USDA Cons. Res. Rep., no. 26. MacGregor, J .M. , G.R. Blake and S.D. Evans, 1974. Mineral nitrogen movement into subsoils following continued annual fertilization for corn. Soil Sc. Soc. Am. Proc. 38, 110-112. Meek, B.D., L.B. Grass, L.S. Willardson and A.J . Mackenzie, 1970. 26 Nitrate transformation in a column with a controlled water table. Soil Sc. Soc. Am. Proc. 34, 235-239. Miller, M.H., 1979. Contribution of nitrogen and phosphorus to subsurface drainage water from intensively cropped mineral and organic soils in Ontario. J . Environ. Qual. 8, 42-47. Newman, M.D., and M. Robinson, 1983. Effects of agricultural drainage on upland stream flow. J . Environ. Manag. 17, 333-348. Phill ips, P.A., J .L.B. Culley, F.R. Hore and N.K. Patni, 1981. Pollution potentials and corn yields from selected rates and timing of liquid manure applications. Trans. ASAE, 24 139-144. Ponnamperuna, F.N., 1972. The chemistry of submerged soils. Adv. Agron. 24, 29-96. Reddy, K.R., 1983. Soluble phosphorus release from organic soils. Agric. Ecosyst. and Environ. 9, 373-382. Reid, I. and R.J. Parkinson, 1984. Seasonal changes in so i l -water redistribution processes affecting drain flow. In J . Bouma and P.A.C. Raats, eds., Water and solute movement in heavy clay soils. ILRI pub. 37, Wageningen, 156-159. Rennes, A . , R.W. Tillman, J.K. Syers and D.G. Bowler, 1976. Effect of mole drainage on surface runoff from a soil under pasture. N.Z. J . Agric. Res. 20, 45-49. Roberts, G. , J.A. Hudson and J.R. Blackie, 1986. Effect of upland pasture improvement on nutrient release in flows from a natural lysimeter and a field drain. Agric. Water. Manag., 11, 231-245. Robinson, M. and K.J. Beven, 1983. The effect of mole drainage on the hydrological response of a swelling clay s o i l . J . Hydrology, 64, 205-228. Rogers, J . S . , C H . Stewart, D.V. Calvert and R.S. Mansell, 1976. Water quality from a subsurface drained citrus grove. 3rd Nat. Drainage Symp. ASAE pub. 1-77, 99-103. Schwab, G.O., N.R. Fausey and D.G. Kopcak, 1980. Sediment and chemical content of agricultural drainage water. Trans. ASAE, 23, 1446-1449. Schwab, G.O., E.O. McLean, A . C Waldron, R.K. White and D.W. Michener, 1973. Quality of drainage water from a heavy-textured soi l . Trans. ASAE, 16, 1104-1107. 27 Seuna, P. and L. Kauppi, 1980. Influence of sub-drainage on water quality and quantity in a cultivated area in Finland. The influence of man on the hydrological regime. Symp. Proc , IAHS pub. 130, 177-187. Sharpley, A.N. , R.W. Tillman and J.K. Syers, 1977. Use of laboratory extraction data to predict losses of dissolved inorganic phosphate in surface runoff and t i le drainage. J . Environ. Qual. 6, 33-36. Sharpley, A.N. and J.K. Syers, 1979. Phosphorus inputs into a stream draining an agricultural watershed. Water Air Soil Poll . 16, 417-428. Sharpley, A.N. and J.K. Syers, 1981. Amounts and relative significance of runoff types in the transport of nitrogen into a stream draining an agricultural watershed. Water Air Soil Poll . 15, 299-308. Skaggs, R.W., A. Nassehzadeh-Tabrizi and G.R. Foster, 1982. Subsurface drainage effects on erosion. J . Soil Water Cons. 37, 167-172. Skaggs, R.W., 1978. A water management model for shallow water table s o i l s . Tech. Rep. 134. Water Res. Inst., U. North Carolina, Raleigh. Switzer-Howse, K.D., 1983. Agricultural management practices for improved water quality in the Canadian Great Lakes basin. Agr. Can., Res. Branch publ. 1983-24-E. Tanji, K.K., 1982. Modeling of the soil nitrogen cycle. In F.J. Stevenson, ed . , Nitrogen in Agricultural S o i l s , Agronomy monograph 22, ASA-CSSA-SSSA, Madison, 721-772. Thornley, S. and A.W. Bos, 1985. Effects of livestock wastes and agricultural drainage on water quality: an Ontario case study. J . Soil Water Conservation, 40, 173-175. Tubbs, L . J . , and D.A. Haith, 1981. Simulation model for agricultural non-point-source pollution. Water Pollution Control Fed. J . , 53, 1425-1433. Van Ommen, H . C , 1985. Systems approach to an unsaturated-saturated groundwater quality model, including adsorption, decomposition and bypass. Agric. Water Manag., 10, 193-203. Walter, M.F., T.S. Steenhuis and D.A. Haith, 1979. Nonpoint source pol lut ion control by soi l and water conservation practices. Trans. ASAE, 22, 834-840. 28 - White, R.E., 1985. The influence of macropores on the transport of dissolved and suspended matter through so i l , in B.A. Stewart, ed., Advances in soil science, vol. 3, Springer Verlag, New York, 95-120. - Whiteley, H.R. and S.R. Ghate, 1979. Sources and amounts of overland runoff from rain on three small watersheds. Can. Agric. Eng. 21, 1-7. - Whitmore, A.P. and T.M. Addiscott, 1986. Computer simulation of winter leaching losses of nitrate from soils cropped with winter wheat. Soil use manag., 2, 26-30. - Willardson, L.S., B.O. Meek, L.B. Grass, G.L. Dickey and J.W. Bailey, 1972. Nitrate reduction with submerged drains. Trans. ASAE, 15, 84-90. -Wil l iams, J .R . , C A . Jones and P.T. Dyke, 1984. A modeling approach to determining the relationship between erosion and soil productivity. Trans. ASAE, 27, 129-144. - Yadav, J.S.P. and L. Leyton, 1960. Effect of drainage on certain physical properties of a heavy clay soi l . J . Soil S c i . , 11, 305-312. - Zwerman, P . J . , T. Greweling, S.D. Klausner and D.J. Lathwell, 1972. Nitrogen and phosphorus content of water from t i le drains at two levels of management and ferti l ization. Soil Sc. Soc. Am. Proc. 36, 134-136. 29 CHAPTER 2 Modelling the effects of drainage on water flow and  nutrient transformations 2.1 Introduction Field drainage is often required in humid agricultural areas to remove excess water and improve crop production. The effects of surface and subsurface drainage on water movement are well documented, but the implications of drainage for nutrient movement and non point source pollution are poorly understood (Switzer-Howse, 1983). A review of f ie ld studies (see previous chapter) where different drainage methods are compared indicates that no universal conclusion can be derived from field studies. Subsurface drainage effects range from an overall increase in nitrogen and phosphorus loadings (Deal et a l . , 1986) to an overall reduction in both (Bengston et a l . , 1986), depending on the specific situations compared. Since these results are so s i te -spec i f ic , i t is necessary to examine the various nutrient transfer mechanisms involved and their relative contribution to nutrient loading if these results are to be explained and extrapolated to other sites. Computer simulation enables one to examine the several nutrient transfer mechanisms (surface runoff, bypassing flow, drainflow, biochemical transformations, etc.) that control nutrient losses from s o i l s . This makes it possible to ascertain the relative contribution of these mechanisms under different drainage regimes. Several models (see previous chapter) simulate nutrient 30 movement in farm soils where subsurface drainage is practiced. These models address specific aspects of nutrient transport and consequently omit or simplify some transfer mechanisms, such as surface runoff or bypassing flow. The purpose of this chapter is to describe an original computer model developed with the specific objective of assessing the overall effect of various drainage regimes on nutrient loading. 2.2 Modelling Approach The overall purpose for which this model is designed is the comparison of expected nutrient losses resulting from different drainage practices. Accordingly, the modelling strategy adopted consists of using a relatively detailed and responsive subprogram to simulate water and nutrient flow, while simulating processes less sensitive to changes in drainage practices in a simplified manner. This approach may reduce the power of the model to accurately predict responses in the f ie ld; however, this relatively low predicting ability is a common feature of nutrient behaviour models, even when highly sophisticated (de Willi gen and Neeteson, 1985). In contrast, remarkably accurate comparisons can be achieved by simpler models, such as is seen in the work of Deal et al. (1986). This is because in a comparison between computer runs, errors common to both runs tend to cancel out. In order to establish meaningful comparisons, simulation must be carried out over a relatively long duration of operations. This requirement imposed a minimum time increment of one day, which 31 restricts the choice of approaches to the water flow simulation. Numerical solutions of Richards' equation remain accurate only when a small time increment is used, because of the strong nonlinearity of the flow equation. (This is not the case for the heat transfer equation, for which a finite difference approach was retained). Furthermore, numerical formulations that include hysteresis (e.g., Dane and Wierenga, 1975) or bypassing flow (e .g . , Sine and Agneessens, 1980) become rapidly cumbersome. Accordingly in this model water flow is simulated using Darcy's law with simplifying assumptions regarding the gradients of pressure potentials. (Table 2-1 l ists these and other assumptions used in the model). Darcy's law governs the water flow in a double one-dimensional parallel domains system, as shown in the diagram of Figure 2-1. This system consists of two zones, the micropores (the matrix or soil ped) zone and the macropores (or channels) zone. Vertical, parallel flows occur in both zones, each of which are discretized in horizontal layers. Each of these layers is assumed to be homogeneous with respect to its state variables (such as moisture and nutrient content, temperature, etc.) and i ts parameters. These may vary from layer to layer, and the choice of number and thickness of layers is left to the operator. Horizontal flow is assumed to take place between the macropore and micropore zone in the form of lateral infiltration or seepage between the peds and the bypassing channels, providing the linkage between the two parallel domains. As can be seen in figure 2-1, the assumption of vertical flow arises from the observation that 32 Table 2-1 List of Assumptions - The flow domain is one-dimensional and vertical. Flow occurs in two dist inct , parallel zones of micropores and macropores. Horizontal flow can occur between the two zones. - The flow parameters are invariant with time. - The flow domain is discretized in horizontal layers. Each layer in each zone is homogeneous and can be characterized by a single value for each water content and nutrient content. - There is continuity in the macropore channels. - Darcy's law governs the flow of water, except in the unsaturated portion of the macropore zone, where water flow is assigned a fixed velocity. - Nutrient flow occurs by advection. Dispersion and diffusion are negligible, and numerical dispersion does not significantly alter the results. - There is no nutrient adsorption in the macropore zone. - Biochemical reactions are described by simple f i rst -order (mineralization, nitrification) or zero-order (denitrification) reactions with variable coefficients. Figure 2-1: Schematic Dlagraa of the Model Layout 34 the flow streamlines above the drains are mostly vertical. Beiow the drains the streamlines are curvilinear and characterized by different travel times, as pointed out by Jury (1978). In order to account for this situation while retaining a one-dimensional formulation, the layer immediately below the drain is given a thickness large enough to include most of the streamlines. In this layer water percolating from above is transferred laterally from the micropore zone to the macropore zone, then percolating from the macropores into the drains. The vertical movement of nutrients is assumed to occur solely as a result of mass flow (advection), while diffusion between the two zones is included with advection in the calculation of horizontal flow over the full depth of the profile. Adsorption and biochemical transformations of phosphate, nitrate, and ammonium, as well as organic forms of nitrogen, phosphorus and carbon are simulated in each layer using empirical relationships from the. literature. Likewise, crop growth and development are simulated in a simplistic, empirical fashion. However, care has been taken to select relationships that include the effects of soil water content, 1n order to make these simulations responsive to different drainage treatments. The organization of these different components of the model are shown in the flow chart of figure 2-2. Since the purpose of the model is to compare the effects of different drainage procedures on nutrient flow, the emphasis is placed on a detailed simulation of the water flow. Consequently the main thrust of the program consists of modelling nutrient advection as determined by evaporation, transpiration, runoff, and i n i t i a l i z e read parameters z n time loop I / read weather 1 c a l l CLIMATE c a l l WATER c a l l NUTRIENT c a l l CROPS 1 p r i n t t o t a l s ( 5 t ° P ) p r i n t d a i l y values Figure 2-2: Flowchart of the Computer Program 36 flow to subsurface drains as well as unsaturated, saturated and bypassing flows simulated as distinct mechanisms. 2.3 Climate As seen in the flowchart of f ig . 2-3, the model reads daily values of minimum and maximum temperature (TMIN and TMAX), precipitation (PPT), and sunshine hours (SUN). These data are used to calculate potential evapotranspiration (PET) and soi l temperatures (TS). The total solar radiation (RADIN) reaching the canopy on a cloudless day can be calculated, following Penning de Vries and Van Laar (1982), as RADIN = TRIG • (1 + 0.034 • cos(DAY))/138.62 2-1 where TRIG = cos(L)«cos(DEC)«sin(HSS) + HSS«sin(L)•sin (DEC) 2-2 DEC = -0.41 • cos(DAY) 2-3 and HSS = arccos(-tan(L) • tan(DEC)) 2-4 where DAY is the Julian day of the year from January 1, L is the latitude and DEC is the declination. The net radiation RAD is found by RAD = RADIN • AB • (0.2 + 0.1048 • SUN/HSS) 2-5 where AB is the absorptivity, which is averaged over the relative area covered by the crop leaves, the residue cover on the so i l , and the soil i tself , taking soil surface wetness into account. This calculate potential evapotranspiration T calculate potential transpiration soil layers calculate actual transpiration 1 calculate root water uptake calculate actual evaporation r soil layers update surface conditions t ^ I return J calculate soil temperatures Figure 2-3: Flowchart of the Subprogram CLIMATE 38 equation assumes that net radiation on a fully overcast day is one f i f th of that of a cloudless day. The Priestley-Taylor equation (Priestley and Taylor, 1972) is used to calculate potential evapotranspiration from net radiation, as PET = RAD • (1 + BEVAP) • SV/(SV + 0.66) 2-6 where SV = (5304/TK2) • EXP (21.255 - 5304/TK) 2-7 TK is the average air absolute temperature, and BEVAP is a constant reflecting the geographical location. Novak (1981) found the Priestley-Taylor equation reliable for the Lower Fraser Valley. The temperatures of each soil layer are calculated using an implicit finite difference formula of the heat conduction equation 1n the vertical dimension. The heat capacity CC and thermal conductivity HC of each soil layer are updated daily as a function of moisture content, given by CC = (2 • FMIN + 4.184 • 9) • (1.0E+6) 2-8 and HC = 17 + 344 • 6 / (6S + 3 • 0) 2-9 where 6 and 6S are the actual and saturated moisture content, and FMIN the volumetric mineral fraction of each soil layer. These functions are adapted from data presented by Hlllel (1980a). The thermal properties of the surface layer are modified to account for crop residue or snow cover. Advective transfer of heat is assumed to be negligible, since the water flux through the profile seldom 39 exceeds 1 mm/h (Kimball and Jackson, 1975). The soil profile is assumed to be insulated at the bottom (if the profile is sufficiently deep), while the top boundary condition is given by using the temperature of the air film at the soil surface. This film temperature is equal to the average daily air temperature, unless the actual evapotranspiration is less than potential. In that case, the film temperature is increased proportionally to the difference, up to a maximum of five degrees, to account for radiative heating. When soil temperatures reach below freezing, sensible heat is converted to latent heat until soil water is fully frozen in the layer under consideration; the reverse occurs during thawing. Soil layer temperatures are set at zero while freezing or thawing is taking place. Once a soil layer is frozen, its hydraulic conductivity is set to zero. It is assumed that this oversimplification of soil freezing phenomena is adequate for the purposes of the model. Snow is allowed to accumulate while air temperatures remain below freezing, during which all precipitation is accounted as snowfall. Above freezing, potential snowmelt (SNM) expressed as a water depth is SNM = TAV • (0.0045 + 0.013 • RAIN) + 0.00025 2-10 according to the equation presented by Kattelmann et al. (1985). where TAV is the daily average temperature. 40 2.4 Water Balance As seen in figure 2.1, water movement occurs in two zones, the micropore and macropore domains. Water can also be held in a third zone labelled the immobile zone, which represents the volume of water that is unavailable for evapotranspiration or percolation because of the high matrix potential under which it is held. This immobile zone completes the total volume of soil water and contributes to nutrient dilution but is not involved in the flow of water. The computation of water movement consists of two phases. In the f i rst one the fluxes into and out of the domain, namely, evaporation, transpiration, rainfall and flows to drains, ditches, and deep seepage, are calculated. This provides a set of init ial and boundary conditions from which water redistribution within the domain is calculated. This redistribution procedure, which includes in f i l t ra t ion , unsaturated flow (upward or downward), bypassing flow, saturated flow below the water table and adjustment of the water table depth, is calculated using a simplified dynamic equilibrium approach (outlined further) as opposed to a numerical formulation of Richard's equation. This sequence of operation is shown in figure 2.4. Evaporation and transpiration are calculated separately as each affect nutrient flow in a different manner. The ratio of potential transpiration (POTRAN) to PET is calculated as dependent 0 calculate rainfall, snowmelt 41 soil layers calculate drainflow, seepage predict water table depth calculate unsaturated static equilibrium potential upward flux update moisture content soil layers calculate surface runoff calculate water table depth ^ return ^ update static equilibrium calculate ped flux update moisture content calculate bypassing flux Figure 2-4: Flowchart of the Subprogram WATER 42 on the leaf area index (LAI), according to POTRAN = PET • 5 • LAI / (2.8 + 4 • LAI) 2-11 which is a f i t of data presented by De Smedt et a l . (1981). POTRAN is set to PET when the LAI grows larger than 2.8 (after Stewart, 1981). Actual transpiration (TRN) depends on available soil water, root density, and POTRAN. Stewart et a l . (1985) presented graphical results of a detailed simulation of moisture extraction by roots. In their graphs the ratio of actual to potential transpiration can be found as a function of root density, soil moisture, and" PET. The relationships expressed by these graphs were approximated with two empirical equations, 2-12 and 2-13. (Since these equations are simple numerical approximations, they are not dimensionally homogeneous and must be used with the units specified in the nomenclature). The root uptake (or transpiration, TRN) is calculated as TRN = POTRAN • 9AV/(50»PET) 2-12 for any value of 0AV inferior to 50'PET, where 6AV is the ratio of actual over potential plant available water. This quantity is calculated with 9AV = (9/9S) • 16 • ROOT/DZ 2-13 1 + 15 • ROOT/DZ where ROOT is the root mass in a given layer. 9AV is calculated for each layer in which there is uptake of water by roots and a 43 depth-weighted average of these values is used in equation 2-12. The depth of the soil over which 9AV is averaged is determined by either the root depth or the water table depth, whichever is smallest. This assumes that root metabolism is inhibited by flooding below the water table, which has been observed in a number of f ield situations (Williamson and Kriz, 1970; Cannell, 1977). Transpired moisture is removed from each layer in proportion of the value of 9AV calculated in each. Potential evaporation (POTEV) is calculated as the difference between PET and TRN. Ponded water is evaporated without restrictions until POTEV is reached, but the actual evaporation rate (EVP) from the soil is assumed to be proportional to the moisture content of the surface soil (down to a depth of ten centimeters), in the intrapedal zone. If the water table is high, capillary flow upwards will replenish some of the water lost to evaporation and transpiration. When used with a deep water table, this model gives results comparable to the well-known model of Ritchie (1972). Water movement below the water table is modelled as saturated flow to drainpipes, ditches, and deep seepage. Broughton and Foroud (1978) showed that Hooghoudt's equation for steady-state drainage (Luthin, 1978) gives satisfactory results when used in a transient simulation with daily time increments. This equation is DRF = 4 • CA • (H - HD)2 + 8 • CB • DEQ • (H - HD) DRS2 2-14 44 As seen in figure 2-1, H is the water table height above drain-depth, HO is the height of water above the drains at the outlets, DEQ is the effective flow depth below the drains, and DRS is the spacing between drains. DRF is the flow rate out of the drains (or d i tches) . CA and CB are the average saturated hydraulic conductivities of the flow regions above and below the drains, respectively. This average is calculated using a depth weighted arithmetic mean since, as soil stratification and heterogeneity between the layers increase, the horizontal component of flow tends to dominate. The total hydraulic conductivity is calculated as the sum of the saturated conductivity of the interpedal zone and the conductivity of the intrapedal (micropores) zone, as these flows occur in parallel. The same equation is used, with the appropriate spacing and conductivity parameters, to calculate flow to ditches. Oitchflow and drainflow are computed separately and summed when ditches are relatively closely spaced and run perpendicular to the drains. Lagace et a l . (1982) showed that this summation does not result in too large an overestimation of flow; a correction can be made to drain and ditch spacing to account for this overestimation. Equation 2-14 is also used to calculate subirrigation, when HD is higher than H, with the appropriate sign inversion. However, when the water table reaches deeper than drain or ditch level and subirrigation is practiced, the depth of lateral flow is then calculated as the difference between the drain or ditch depth and the equivalent depth DEQ. Deep seepage is assigned a constant 45 value by the operator. Once the i n i t i a l and boundary c o n d i t i o n s of evapotranspiration, saturated drainage (or subirrigation), and rainfall (or snowmelt) are known, the water is redistributed within the different soil layers. Redistribution is calculated using Darcy's law rather than empirical equations, using simplifying assumptions to reduce computational requirements. The expected final moisture content in each layer is calculated for different conditions, from which the dynamic equilibrium moisture content is selected. A detailed example of this procedure, outlined below, can be found in Appendix I. The top layer is located above the soil surface, where ponded water and rainfall intercepted by the canopy is stored. The water content at saturation of this layer is modified during the simulation to account for crop growth, residues, and tillage since this value represents the maximum volume available for ponding. All other layers are below the soil surface and are characterized by a moisture content 0 (less or equal to saturation, 0S), that represents water held in the interpedal and intrapedal zones. The height of each layer above the water table determines an equilibrium moisture content (9EQ) for each layer. 9EQ is the sum of the water held in the macropores (interpedal water), and of the water in the intrapedal zone, which is characterized by a conventional moisture content vs. suction curve. Because of hysteresis in this curve, 9EQ takes on a different value for wetting or drying (9EQW or 9EQD). 46 The hydraulic conductivity is calculated using a discontinuous function, as shown in figure 2-5, with separate values for the macropores and micropores zones. The intrapedal conductivity C is calculated, after Gardner (1960), as C = CS • (9/9PR)A L F A 2-15 where CS is the intrapedal saturated conductivity, ALFA is an empirical coefficient, and 6PR is the saturation value of the intrapedal zone. When 9 is larger than 9PR, C becomes equal to CS, while the conductivity of the interpedal zone, relevant only in this situation, is given a constant value, CK. By adjusting the values of 9PR, CS and CK, the relative importance of bypassing flow and throughflow of water in the soil can be altered. Moisture redistribution, in any layer, takes place either in upward or downward flow, depending on the prevailing conditions. Before downward flow is calculated (if any), the model computes upward flow whenever a layer has a moisture content inferior to 9EQW. This flux is assured to take place under a constant gradient and results, after one time increment, in the new water content 9UP. The*smallest of the values of 9EQW and 9UP is assigned to the layer. Downward flux is then computed and the resulting dynamic equilibrium value, 9QD, takes on either the value of the init ial moisture content, of the static equilibrium 9EQD, or of QDR or 0FL, which represent values obtained from continuous drainage or continuous inf i l t rat ion, respectively, as detailed below. The largest of these values is used for 9QD since it is the f irst 47 HYDRAULIC CONDUCTIVITY Figure 2-5: Typical Hydraulic Conductivity Versus Moisture Content Curve 48 limiting value encountered as drainage proceeds. If this value of 6QD is larger than the matrix saturation value, OPR, the resulting excess water is drained by macropore flow. This sequence of calculations is diagrammed in figure 2-6. When there is negligible flow of moisture into a layer from above, this layer will drain to a value of 9DR as This equation can be obtained by integrating Darcy's law for vertical flow under a unit gradient for a duration DT, using equation 2-15 to express unsaturated conductivity. Nielsen et a l . (1973) found the assumption of unit gradient satisfactory to model unsaturated drainage. Moisture content may also reach a value 8FL which is at equilibrium with unsaturated flux into the layer from above under a unit gradient, following Rubin (1966), as where Q is the depth of water flowing during the interval DT. The assumption of unit gradient will fai l when the layer is at static equilibrium with the water table. In this case, no water will drain from the layer below the value of 6EQ. Otherwise, the 1 1-ALFA 2-16 1 2-17 a 0/ 6FL 9QD dz T 9EQD WTD V — SDR —1 r— 6PR es - I n i t i a l c o n d i t i o n s S t a t i c e q u i l i b r i u m - S t eady d r a i n a g e F l u x e q u i l i b r i u m - Dynamic e q u i l i b r i u m 6DR < 6EQD < 9FL .*. 9QD • GPL M o i s t u r e c o n t e n t Figure 2-6: Diagram of the Dynamic Equilibrium Moisture Content Concept 50 dynamic equilibrium value 8QD becomes the larger of 0FL or 9DR. The layer will drain to this value, provided drainage is not restricted by the layer below i t . The maximum flux PINF that can infiltrate a given layer (in intrapedal flow) is calculated as where GRAD is the average gradient (over time) under which infiltration is taking place. GRAD is set at unity unless the intrapedal zone becomes saturated. In that case, GRAD is calculated as the average obtained from the Green and Ampt model (Hil lel , 1980b). This value of GRAD is obtained by dividing the infiltration volume predicted by the Green and Ampt model after a time DT by the value obtained with a unit gradient, ie, CS • DT. However, since the average gradient cannot be determined for vertical in f i l t ra t ion without i terat ing, this model uses the gradient obtained for horizontal in f i l t ra t ion , and adds this gradient to the vertical unit gradient, as which approximates the vertical infiltration solution. S is the suction ahead of the wetting front. However, when the conductivity of the layer in which water is infiltrating is much larger than that of the layer above i t , as for instance below a soil crust or a plowpan, GRAD is set constant as GRAD = (DZ + S)/DZ 2-20 PINF = CS • GRAD • DT 2-18 GRAD = 1 + 2-19 51 Water is cascaded downwards from layer to layer. The resulting distribution may result in some layers with a water content superior to 6PR, and with water ponded in the top layer. This excess water is cascaded down in the interpedal zone, using CK as the maximum flux that can infiltrate into the macropore channels of a layer. Lateral infiltration from interpedal to intrapedal zone QSD is calculated, again using Green and Ampt, as QSD = J(Z • CS • S • (8PR-6) • DT)' • (DZ/PEDS) 2-21 where PEDS is the average spacing between macropore channels (assumed to be f lat , wide cracks). Following this procedure, any layer with a moisture content higher than 9S is adjusted back to 9S, and the excess moisture cascaded upwards to the next layer above. In effect, this procedure corrects excessive downward flow predictions. Excess water in the f irst layer is removed and disposed of as surface runoff. Upward flow is assumed to take place only in the intrapedal zone, causing the moisture content to increase to a value of 9UP as 9UP = 9 + CSAV • S - DWT • DT DWT DZ 2-22 where CSAV is the average unsaturated conductivity between two layers, DWT is the height above the water table, and S is the average suction of the layer during the flow duration DT. This is a simple expression of Darcy's law at steady-state, using the 52 wetting limb of the moisture retention curve to calculate S. Since upward flow usually takes place at a relatively uniform rate, a steady-state expression gives satisfactory results, as shown by Ragab and Amer (1986). To avoid excessive flow values the upper limit that a layer can reach from upward flow is set equal to water content at static equilibrium, for wetting. The amount of water cascaded down to the saturated zone is substracted from the amount removed by capillary upward flow and saturated drainage. The difference between the amounts determines the updated depth of the water table. Water-table depth is a determining factor in the calculation of moisture redistribution, because of the static equilibrium values 6EQ. In order to avoid feedback oscillations, of the calculated depth, the values of 8EQ are calculated according to a predicted water table depth rather than to the actual depth at the previous time increment. The validity of the dynamic equilibrium approach was tested by comparing two simple i n f i l t r a t i o n -redistribution events with the solution of Richard's equation as approximated by a finite difference procedure. It was found that the differences between the two solutions were minimal. Details of the validation procedure are provided in Appendix II. 2.5 Nutrient Balance The model simulates the movement and transformations of nitrogen, phosphorus and carbon according to their distribution in nine distinct types of nutrients. The amounts of each type of 53 nutrient 1n each layer are referred to as pools. A separate pool 1s assigned to nitrate, ammonium and orthophosphate. Organic matter is accounted for by the six other pools, which keep track of moderately stable organic matter (carbon, nitrogen and phosphorus), and of labile organic matter, which represents microbial biomass and rapidly decaying wastes. Fig. 2-7 shows the flowchart of the nutrient submodel, while a diagram of the processes is shown in figure 2-8. Nutrient movement is simulated as an advection process, whereby the nutrients are carried by water flow. Vertical diffusion or disperson are not explicitly simulated; however, some dispersion occurs as a consequence of bypassing flow and numerical dispersion as show by Van Hoorn (1981). Vertical nutrient movement occurs simultaneously 1n the interpedal and intrapedal zone, while lateral movement between the zones occurs either from lateral Infiltration of water or from diffusion. Diffusion between the intrapedal and interpedal zones is simulated using a mixing algorithm similar to that used for surface runoff, as described below. The dissolved nutrient content of surface runoff and ponding is calculated as that resulting from mixing uniformly runoff water and soil water held within one centimeter of the surface. This mixing is assumed to result from diffusion of solutes between water at the surface of the matrix and water held in i t . This common procedure has been reviewed by Ahuja (1986) and the theoretical basis for i t can be found in Appendix III. This procedure similarly alters the nutrient concentration of bypassing flow 54 soi l layers assimilation carbon oxidation mineralisation nitri f icat ion denitrification volatil ization - / s o i l layers \_ \ nutrients / micropore mass flux macropore mass flux lateral mass flux surface erosion nutrients surface losses Figure 2-7: Flowchart of the Subprogram NUTRIENT 55 F e r t i l i z a t i o n and p l a n t decay V F e r t i l i z a t i o n Min N u t r i e n t p o o l s M i n L e a c h i n g losses N03, P04 Legend OCS: s t a b l e carbon p o o l ONS: s t a b l e n i t r o g e n p o o l OPS: s t a b l e phosphorus p o o l M i n : m i n e r a l i z a t i o n V o l : v o l a t i l i z a t i o n OCL: l a b i l e ca rbon p o o l ONL: l a b i l e n i t r o g e n p o o l O P L : l a b i l e phosphorus p o o l N i t : n i t r i f i c a t i o n D e n : d e n i t r i f i c a t i o n Figure 2-8: Diagram of the Nutrient Balance 56 infiltrating the interpedal zone from the surface. The degree of mixing between intrapedal water and interpedal flow is the ratio of the mixing thickness, again taken as one centimeter, tp the average ped spacing (PEOS). Nutrient adsorption in the interpedal zone is assumed to be negligible, while it is calculated in the intrapedal zone according to the Langmuir isotherm (Hi l lel , 1980a) as CEQ/CMA = CEQ/CM + 1/(AK • CM) 2-23 where CEQ and CMA are the nutrient mass, dissolved and adsorbed, respectively, per unit volume, CM is the adsorption capacity at saturation and AK is the adsorption coefficient. This procedure greatly simplifies the computational requirements of this submodel, but imposes some limitations, which are addressed in section 2.7. The composition of drain water (or of water draining to a ditch) reflects the relative contributions of intra and interpedal flows. Drainflow is taken from the water of the soil layers located between drain depth and equivalent depth (DEQ), from both the interpedal and intrapedal zones, in proportion to the relative magnitude of bypassing flow and throughflow. This procedure ensures adequate mixing of the water reaching the drain tube from lateral as well as vertical flow and is numerically similar to that discussed by Jury (1978). Nitrate and orthophosphate are also removed by root uptake. Like water uptake, nutrient uptake is assumed to place only in the root zone above the water table. Nitrate absorption is calculated 57 as a mass flow carried by transpiration water. Liao and Bartholomew (1974) verified that ignoring diffusion of nitrate to the roots does not result in significant errors under normal conditions. However, this simplification is not valid for phosphorus, because of the low concentration of this strongly adsorbed species in soi l water. Instead, the approach of absorption as controlled by crop needs used by Williams et a l . (1984) is followed here. It assumes that phosphate will diffuse to the root mass according to root uptake, provided that enough phosphate is available in the soil layer. This mass of phosphate PAV available to the roots is calculated as PAV = CEQ • 6 • DZ • 6AV 2-24 where 9AV is an availability factor dependent on root density and soil moisture given by equation 2-13. Surface runoff can result in a loss of dissolved nutrients, as discussed previously, but also in soil erosion and the consequent loss of adsorbed nutrients even in level areas. Soil loss (SOILL) is predicted using the universal soil loss equation as modified by Onstad and Forster (1975) for single storms, calculated as SOILL = BR • SEF • CMF • EPF • SLF 2-25 where BR is the storm erosivity, and SEF, CMF, and EPF are the soil erosion factor, the soil cover factor, and the erosion control practices factor respectively, calculated according to the soil conservation service guidelines (Wischmeir and Smith, 1978). This 58 equation is commonly used, even for relatively flat topography (e.g., Williams et a l . , 1984; Deal et a l . , 1986). SLF is the field slope and length factor, calculated following Williams et a l . (1984) as SLF = (0.045 • FIELDL)GAM • (65.41 • SLP2 + 4.56 • SLP + 0.015) 2-26 with GAM = 0.6 • (1 - EXP(-35.835 • SLP) 2-27 where FIELDL is the slope length and SLP the slope. The storm erosivity is given by BR = EI + 1558 • RFF • (RAIN/DT)0-33 2-28 with EI = (2176.2 + 351.26 • LN (RAIN/DT)) • RAIN2/DT 2-29 This equation is adapted from the model EPIC (Williams et a l . , 1984) and takes into account both the erosivity of rainfall and of runoff. RAIN/DT is the average rainfall rate and the coefficients reflect the assumption that the peak half-hour rainfall rate is twice that of average, and absolute peak is three times that value. Erosion arising from snowmelt without precipitation is modelled by setting EI to zero and using the snowmelt value as RAIN in eq. 2-28. Once the rate of soil loss is determined, the loss of adsorbed phosphorus and nitrogen from the various pools can be calculated as S0LIDL = TOTS • (S0ILL/TSM) • ENR 2-30 59 where SOLIDL is loss of nutrient from a total TOTS adsorbed on the soil surface (including surface residues) down to a depth of one centimeter, representing a mass of soil TSM. ENR is the enrichment ratio which is assumed to be constant and equal to 2, a value commonly encountered in field situations (Dean, 1983). Soil loss may also occur as a result of drainflow. In the model, sediment loss in drains is assumed to be proportional to drainflow, while nutrient loss in drain sediment is calculated using eq. 2-30, with the appropriate values for TSM and TOTS. Nutrients are transferred from one pool to another according to the rate of biochemical reactions taking place in each soil layer. These reactions are usually microbially mediated and thus strongly affected by environmental factors. Following Johnsson et a l . (1987), this model assumes that a single coefficient can be used for all reactions, as ALOSS = A • K • EM • ED • ET • ES 2-31 where ALOSS is the daily loss of a nutrient from pool A, K is a decay constant specific to a given reaction and EM, ED, and ET are environmental coefficients calculated as 0.8 < 9/9S < 1 2-32 o < e/es < 0.3 2-33 2-34 where EM accounts for excessive soil moisture, ED for lack of EM = 2 • 9/9S - 1.4 ED = 3.33 • e/es E T = B E T A T S - 2 0 60 moisture, and ET for soil temperature, with BETA a temperature coefficient. The decay constants K can be derived from f irst order decay constant published elsewhere (e.g. Overcash et a l . , 1980) when it is recognized that ALOSS = A • K = A • (1 - EXP (-K0 • DT)) 2-35 where K0 is a f irst order decay constant and ALOSS is the loss from pool A during the interval DT. Although the use of a temperature coefficient ET is an unusual formulation, it can be shown that, for small values of ET and K, it approximates the exponential increase of reactivity with temperature of Arrhenius' equation. Equations 2-31 to 2-34 are adapted from Johnsson et a l . (1987). ES is a coefficient used only in the mineralization and nitr i f icat ion calculations. The b iochemical t ransformat ions modelled include volatil ization, immobilization, mineralization, nitrification and denitrification. Mineralization from the stable organic pools into carbon dioxide, ammonium, and orthophosphate can be calculated directly using equation 2-31, with the appropriate coefficient, and setting ES equal to one. Assimilation of mineral nutrients into the stable organic pools is neglected, and these pools are replenished solely by addition of organic waste and from crop and root decay. Mineralization from the labile organic pools is calculated in a similar way, except that immobilization occurs simultaneously. This model assumes that the soil microbial fraction of the labile 61 biomass has a constant ratio of carbon to nitrogen to phosphorus. This microbial biomass fraction can be determined by finding the limiting constituent, whether carbon, nitrogen, or phosphorus. Mineralization of organic carbon is then calculated by equation 2-31, where ES becomes the ratio of microbial biomass to total biomass in the labile pools. In effect, this assumes that biomass decay is proportional to microbial biomass. Mineralization of nitrogen and phosphorus is calculated in a similar fashion, except if one of these is limiting. In this case, the mineralized fraction is assumed to be recycled internally such that no net mineralization occurs; instead, assimilation from the mineral pools takes place, again at a rate proportional to ES, the ratio of microbial to total unstable organic matter. However, the mineralization coefficient is not necessarily the same as the immobilization coefficient. Immobilization of nitrogen is taken f irst from the ammonium pool, then from the nitrate pool when ammonium is depleted. Nitrification also proceeds according to equation 2-31, except that it is assumed to be completely inhibited when the soil is saturated. ES is then calculated as ES = 20 • (1 - 9/9S) 0.95 < 9/9S < 1 2-36 Thus nitrification is calculated as a single step reaction and the production of nitrite is not simulated. Denitrification is simulated as a one step process, using a modification of the equation of Johnsson et a l . (1987). It assumes 62 that denitrification is proportional to the carbon mineralization rate and related to the presence of oxygen, which can be accounted for as DNIT = A • FN • MC02 • (9/9S - 0.7) 2 0.7 < 9/9S < 1 2-37 where A is the nitrate pool, FN is a conversion factor, MC02 is the daily rate of carbon mineralization, and DNIT is the daily denitrification. Volatilization of ammonium (AVOL), as ammonia, is assumed to occur down to a depth DV as AVOL = K • (1 - DPT/DV) • CEQ • 9 • DZ 2-38 where DPT is the layer depth and CEQ • 9 • DZ represents the total amount of ammonium or ammonia in solution in the depth DV. The volatilization constant K is varied according to the acidity of the soi l , but variations due to weather are not considered. 2.6 Crop Growth Fig. 2-9 shows the architecture of the subprogram CROP. Plant growth is initiated at a given day DSEED which corresponds to emergence, and is assumed to be proportional to transpired water TRN, a situation commonly observed (Hanks, 1983). This is expressed as GROW = TRN • TPCON • STRES 2-39 Figure 2-9: Flowchart of the Subprogram CROP 64 where GROW is daily growth (in dry matter increments), TPCON is a conversion coefficient, and STRES accounts for stress due to lack of nutrients. Other causes of growth reduction are either ignored or assumed to be incorporated in TRN. Nutrient stress is calculated as STRES = (PLANTN - PNMIN) / (PNOPT - PNMIN) PNMIN < PLANTN < PNOPT 2-40 where PLANTN is the plant concentration of nitrogen or phosphorus, PNMIN is the concentration below which growth is completely stopped, and PNOPT is the concentration at which no stunting effects are observed. This equation is a simplification of the procedure followed by Williams, et a l . (1984) Plant growth is calculated in two phases. During the init ial phase, vegetative growth only is simulated, and the dry matter increment is used to develop roots and leaves. During the second phase dry matter is assumed to be stored only in reproductive organs and all vegetative growth is stopped. The transition from the f irst phase to the second occurs after a given number of heat units TUREP has been accumulated since DSEED. Similarly, growth stops completely ( i .e . , the plant dies or enters dormancy) after a limiting number of days or heat units has been accumulated since growth, or after freezing temperatures are encountered. Stewart (1981) found that the leaf area index LAI of most crops grown in Canada could be represented by the relation 65 LAI = (2.77 X 10" 5) • DGRO 3.115 2-41 which is used here. DGRO represents the number of time increments counted from DSEED. After LAI reaches a value of 2.8 m2/m2 canopy coverage is complete, and although leaves continue to grow, the ratio of potential evaporation to transpiration (eq. 2-11) becomes „n-ii'T^artition of growth between root and shoot, PP, varies during the vegetative growth period and is calculated as PP = PART • (1 - ATU/TUREP) 2-42 where PART is the partition coefficient at the onset of growth and ATU is the accumulated heat units since DSEED. This equation is a simplification of the partition equations given by Penning de Vries and Van Laar (1982). Borg and Grimes (1986) showed that root depth can be described with a simple sinusoidal function. This model uses the derivative of their relationship to calculate the potential root depth increment, RTDI, as RTDI = 1.515 • RTDMX • (TU/TUREP) • where TU is the number of heat units during the time increment, and RTDMX is the maximum rooting depth reached after the accumulated heat units, ATU, total a value of TUREP. Roots will grow by this amount unless restricted by the presence of a shallow water table, which stops root growth. Daily increments of weight are added to cos(3.03 • ATU/TUREP - 1.47) 2-43 66 the roots during the vegetative period and are distributed evenly through the root zone (or down to the depth of the water table). Roots submerged by a rising water table are assumed to decay at a rate given by SLOF = ROOT • FLOODF • TS 2-44 where SLOF is the daily amount of sloughing, ROOT is the root mass in a soil layer, FLOODF is a decay coefficient that varies with crop sensitivity to flooding and TS is the soil temperature. This equation is original and gives the operator the option of altering the roots distribution as a response to flooding if there is documentation to derive the value of FLOODF for a particular crop. When the crop dies, it is partitioned between a harvested part and a residue part. The residue part itself is partitioned between the pools of stable and unstable organic matter of the soil and soil surface. Using a different partition coefficient for carbon, nitrogen and phosphorus makes 1t possible to mimic the effects of nutrient translocation in the end results while dispensing with the task of simulating such translocation for each time increment. This end partitioning method is also used to simulate the effects of a nitrogen fixing crop. In this case, there is no nitrogen stress on crop growth, but the crop nevertheless absorbs nitrate by mass flow from the transpiration stream. The nitrogen concentration of the crop at death is fixed at a minimum, the balance assumed to have originated from fixation. Nitrogen is returned to the soil by using the appropriate partitioning 67 coefficients at crop death. 2.7 Discussion In any attempt at simulating a situation as complex as that presented here, it is to be expected that the comprehensiveness of the model must bear a price in terms of lack of accuracy. This is partially because a number of mechanisms simulated are poorly known or quantified, such as the biological dynamics of plant or microorganism growth. However, even the more familiar phenomena of soil water flow or advection-dispersion of solutes must include simplications, if only because of restrictions in the number of spatial dimensions, or of the requirements of discretization of the time and space domains. The problem, therefore, consists of demonstrating that the errors and inaccuracies caused by the simplifications do not affect the validity of the conclusions derived from the model. This is easier to achieve when developing a model for comparative, as opposed to predictive, purposes. Therefore it was felt that the relatively crude computation of biological phenomena is adequate, because the sensitivity of nutrient movement to these (including n i t r i f i c a t i o n and denitrification) is low, as shown by Kanwar and Johnson (1984). Conversely, water flow must be simulated with a certain degree of accuracy, and flow mechanisms cannot be excluded without the relevance of the conclusions of the model to field situations becoming questionable. For the same reason, simulations of a relatively long duration must be carried out, which imposed on this 68 model the water balance approach as opposed to the more rigorous, but time intensive numerical solution to the water flow equations. Water balance models have problems of their own, as can be seen below. A number of water balance models (Chieng, 1975; Skaggs, 1978; Makkink and Van Heemst, 1975; Madramootoo and Broughton, 1987) determine the position of the water table and use an empirical or physical relationship to determine the rate of flow between the saturated and unsaturated zone. This type of model cannot be used to simulate nutrient flow without modifications, because of the lack of discretization of the domain. Several models use a flow domain discretized into horizontal layer and calculate water infiltration into, and drainage out of, these layers according to empirical functions (Baier et a l . , 1979; Wilkie, 1981; Williams et a l . , 1984; Thomas and Beasley, 1986). Such functions provide realistic simulations but must be properly calibrated as they rely on such concepts as field capacity which become calibration parameters (Carneiro Da Silva and De Jong, 1986). Furthermore, a delay must often be assumed to take place between infiltration and drainage; Wilkie (1981) found that the common assumption of instantaneous drainage can lead to an overprediction of leaching losses. This delaying function, in effect, accounts for the fact that both infiltration and drainage are continuous, although they are computed as a discrete step. As an alternative, some models (e.g. Bulley et a l . , 1980) use Darcy's law to calculate drainage under a unit gradient. However, this 69 assumption fails for soils with shallow water tables. Furthermore, in all cases a separate function must be provided to calculate upward flux from the water table. While this flux may be small compared to drainage, Makkink and Van Heemst (1975) showed its overriding importance in simulating summer water table depth fluctuations. (They found that a small change in the upward flux calculation can lower or raise the water table by several meters over a season.) The present model consists of a combination of physical and empirical approaches. Like an empirical model, it assigns an equilibrium moisture content to each layer. However, this equilibrium content is determined using physical laws of water flow, according to the existing flux situation. This makes it possible to dispense with delay functions and preset field capacity values. Because the flux draining out of a layer is also compared with the inf i l t rat ion capacity of the layer below i t , this formulation allows for a realistic simulation of layered soi ls. Finally, this approach also makes a physically-based simulation of bypassing flow possible in a flexible manner. Thus, for instance, the model calculates equilibrium moisture conditions during an infiltration event rather than calculating inf i l t rat ion to use as a set of in i t ia l conditions prior to drainage. The assumptions used in calculating infiltration are similar to those of Green and Ampt (in Hi 1 le i , 1980b); the Green and Ampt model can easily be adapted to varying conditions (Shirmohammadi and Skaggs, 1985) and may be more accurate than the 70 SCS curve number method even in soils for which a sharp suction front does not exist (Rawls and Brakensiek, 1986). Bypassing flow is easily included with such an approach by a calculation of fluxes in parallel. This model combines the downward flow and horizontal in f i l t ra t ion approach of Beven and Germann (1981) with the aggregate diffusion concept of Addiscott (1981) and the parallel flow idea of Armstrong and Arrowsmith (1986). Keng and Lin (1982) showed that the hydraulic conductivity of common Canadian soils behaved as a macropore-micropore network, which further reinforces the need to include bypassing flow as a transfer mechanism. Another advantage of using this physically based dynamic equilibrium concept lies in the fact that the model produces results which converge towards analytical solutions at steady state, for any set of initial conditions. For instance, a steady low-intensity rain in an initially dry profile will result in a moisture distribution similar to that given by the solution of Rubin (1966). Likewise, a steady state solution of static equilibrium will eventually be achieved when no moisture is added to or removed from the prof i le . Furthermore, the use of equilibrium values imposes a set of limits to the permissible values of the moisture content in a layer. This prevents overpredictions of the flux in the profile which could otherwise result in an unstable, oscillating solution. This feature makes it possible to calculate separately upward flux, downward. flux, and bypassing flux while avoiding the danger of creating unbounded oscillations. 71 Nutrient flow is calculated in this model as carried by water advection. While this approach offers the advantages of simplicity, it often fai ls to predict leaching accurately, in particular because of the stochastic nature and spatial variation of both water and nutrient flow. Accordingly, Jury (1982) proposed that a stat ist ical description of the flow mechanism, called Transfer Function Model, may present a better representation of the field situation. Addiscott et a l . (1986) tested this model but concluded that the improved accuracy may not warrant the calibration effort. Furthermore, White (1985a) found that this model could not adequately describe bypassing flow in dry soi ls. Accordingly, the simpler approach was retained in this model. It was assumed that the Langmuir isotherm could be used to model the adsorption of both phosphate, ammonium and organic species. This isotherm can predict a value for saturation and is computationally simple. Although it has obvious limitations (for instance, when modeling ammonium movement in alkaline soils) , this procedure improves the versatility of the model when compared to a simple designation of some nutrients as completely immobile (eg. Bulley, et a l . , 1980). Adsorption along the macropores is assumed to be driven by diffusion into the micropore matrix, and there is no provision for adsorption sites within the macropores. This approach was retained for its simplicity and remains adequate when diffusion within the matrix is rapid. However, when the matrix consists of large structural units between cracks, there may be an important horizontal concentration gradient within the matrix. As 72 discussed by White (1985b), this situation will clearly affect the extent by which solutes may be added or removed from the matrix by bypassing flow. In this case, the assumption of homogeneity within a layer is inadequate. This limitation can be partially remedied by adjusting the parameter PEDS but a more realistic approach would require a more complex model which would also be very difficult to validate. The effects of flooding on crop growth and development are often simulated using empirical relationships or indexes such as the SDI (Battacharya and Broughton, 1979). Since an accurate simulation of crop yield is not the object of this model, the effects of high water tables on crops are calculated only indirectly. The model assumes that roots do not grow, nor are active beyond the water table; any reduction in yield arises subsequently from a reduction in water or nutrients uptake, including that caused by denitrif ication. Although largely untested, this approach is closer to the physical reality. A relatively similar approach, which also includes the effects of drainage on crop emergence, resulted in a successful simulation of potato yields (van Wijk and Feddes, 1986). A commonplace test for verifying a model consists of comparing its output against an analytical solution when modeling a simple situation. To do this on this type of comprehensive model is difficult because several of its components do not lend themselves to an analytical solution (eg, crop growth). Furthermore, since analytical solutions form an integral part of the model (in the 73 water flow subsection), a comparison of the model output with these solutions only ensures that these have been transcribed properly in the model. However, the overall response of the model can be tested against the analytical solution developed by Van Ommen (1985), which describes the concentration of a solute in field drains as a function of time. The solution assumes a steady state water flow through the profile and includes saturated, unsaturated, and bypassing flows. When the model was modified to reproduce these conditions, using a non-adsorbed, conservative solute, the model output was nearly identical to the analytical solution. Details of the comparison are given in Appendix IV. In summary, this chapter presents a description of an original model designed for long term simulation of water and nutrient movements in humid agricultural soils. Its main purpose is to compare the effects of different drainage practices on nutrient losses under various conditions, including soils where bypassing flow is an important mechanism. The model features a comprehensive simulation of water and nutrient flow and transformation mechanisms, including crop growth simulation. In contrast with other models of water and nutrient flow, the emphasis of this model is on a detailed simulation of water flows, for which an original and flexible approach has been developed. 74 -9 REFERENCES Addiscott, T.M., P.J. Heys and A.P. Whitmore, 1986. Application of simple leaching models in heterogeneous soils. Geoderma, 38, 185-194. Addiscott, T.M., 1981. Leaching of nitrates in structured soils in M.J. Frissel and J.A. van Veen, eds. Simulation of nitrogen behaviour of soil-plant systems. Centre Agric. Pub. Docum., Wageningen. Ahuja, L.R., 1986. Characterization and modeling of chemical transfer to runoff. In B.A. Stewart, ed., Adv. Soil S c i . , 4. Springer-Verlag, New York, 148-188. Armstrong, A.C. and R. Arrowsmith, 1986. Field evidence for a bi-porous soil water regime in clay soils. Agric. Water Manag., 11, 117-125. Baier, W., J.A. Dyer and W.R. Sharp, 1979. The versatile soil moisture budget. Agrometeorol. Sect. LRR1 tech. bul l . 87, Research Branch, Agriculture Canada. Battacharya, A.K., and R.S. Broughton, 1979. Spacing drains to maximize revenue increases. Trans ASAE, 22, 351-356. Bengtson, R.L., C E . Carter, H.F. Morris, and J.G. Kowalczuk, 1984. Reducing water pollution with subsurface drainage. Trans. ASAE, 27, 80-83. Beven, K. and P. Germann, 1981. Water flow in soil macropores. J . Soil S c i . , 32, 15-29. Borg, H., and D.W. Grimes, 1986. Depth development of roots with time: an empirical description. Trans. ASAE, 29, 194-197. Broughton, R.S. and N. Foroud, 1978. A model to predict water table depths for flat lands. Can. Agric. Eng. 20, 81-86. Bulley, N.R., R.A. Bertrand and B. Cappelaere, 1980. Manure nitrogen management: a model. In Livestock waste: a renewable resource. Proc. 4th int. symp. livestock wastes. ASAE, St. Joseph. Cannell, R.Q., 1977. Soil aeration and compaction in relation to root growth and soil management. In T.H. Coaker, ed., Applied Biology, vol. II, Academic Press, New York, 1-86. Carneiro Da Silva, C , and E. DeJong, 1986. Comparison of two computer models for predicting soil water in a tropical monsoon climate. Agr. For. Meteor., 36, 249-262. 75 Chieng, S.T., 1975. The effects of subsurface drain depths and drainage rates on water table levels. Unpub. M.Sc. thesis, McGill Univ., Montreal. Dane, J .H . , and P.J. Wierenga, 1975. Effect of hysteresis on the prediction of infi ltration, redistribution and drainage of water in a layered so i l . J . Hydrol. 25, 229-242. Deal, S . C , J.W. Gilliam, R.W. Skaggs, and K.D. Konyha, 1986. Prediction of nitrogen and phosphorus losses as related to agricultural drainage system design. Agric. Ecosyst. Environ. 18, 37-51. Dean, J .D . , 1983. Potency factors and loading functions for predicting agricultural nonpoint source pollution. In F.W. Schaller and G.W. Barley, eds., Agricultural management and water quality. Iowa State U. press, Ames, 155-177. F. De Smedt, S. Al Khafaf and P.J. Wierenga, 1981. Simulation of water flow in plants and soils. In Dubois, D.M., ed., Progress in Ecological Engineering and Management by Mathematical Modelling. Cebedoc ed., Liege, 849-863. deWilligen, P. and J . J . Neeteson, 1985. Comparison of six simulation models for the nitrogen cycle in the so i l . Fertilizer Res., 8, 157-171. Gardner, W.R., 1960. Soil water relations in arid and semi-arid conditions. UNESCO 15, 37-61. Hanks, R.J. , 1983. Yield and water use relationships. In H.M. Taylor, W.R. Jordan and T.R. Sinclair , eds., Limitations to efficient water use in crop production. ASA, CSSA, SSSA, Madison, 393-412. H i l l e l , D., 1980a. Fundamentals of soil physics. Academic Press, New York. H i l l e l , D., 1980b. Applications of soil physics. Academic Press, New York. Johnsson, H., L. Bergstrom, P.E. Jansson and K. Paustian, 1987. Simulated nitrogen dynamics and losses in a layered agricultural soi l . Agric. Ecosys. Environ., 18, 333-356. Jury, W.A., 1982. Simulation of solute transport using a transfer function model. Water Resource Res. 18, 363-368. Jury, W.A., 1978. Estimating the influence of soil residence time on effluent water quality. In P.F. Pratt, ed., Management 76 of nitrogen in irrigated agriculture, Sacramento, 265-290. Kanwar, R.S. and H.P. Johnson, 1984. Sensitivity analysis of a drainage simulation model. Water Resources Bul l . , 20, 339-342. Kattelmann, R.C., N.H. Berg and M.K. Pack, 1985. Estimating regional snow water equivalent with a simple simulation model. Water Resource Bul l . , 21, 273-280. Keng, J.C.W., and C.S. Lin, 1982. A two line approximation of hydraulic conductivity for structured soils. Can. Agric. Eng., 24, 77-80. Kimball , B.A. and R.D. Jackson, 1975. Soil heat flux determination: a null-alignment method. Agric. Meteorol. 15, 1-9. Lagace, R., R.W. Skaggs and J .E . Parsons, 1982. Predicting water table drawdown for two-dimensional drainage. Proc. 4th Nat. Drainage Symp. ASAE, St. Joseph, 6-15. Liao, C.F.M. and W.V. Bartholomew, 1974. Relation between nitrate adsorption and water transpiration by corn. Soil Sc. Soc. Am. Proc. 38, 472-477. Luthin, J .N . , 1978. Drainage Engineering. Robert E. Krieger Pub., Huntingdon. Makkink, G.F., and M.D.J, van Heemst, 1975. Simulation of the water balance of arable land and pastures. Centre for Agr. Pub. Docum., Wageningen. Madramootoo, C.A. and R.S. Broughton, 1987. A computer simulation model of surface and subsurface flows from agricultural areas. Can. Water Resources J . , 12, 30-43. Nielsen, D.R., J.W. Biggar and K.T. Erb, 1973. Spatial variability of field measured soil water properties. Hilgardia, 42, 215-259. Novak, M.D., 1981. The moisture and thermal regions of a bare soil in the Lower Fraser Valley during spring. Unpub. Ph.D. thesis, U. British Columbia. Onstad, C.A. and G.R. Foster, 1975. Erosion modelling on a watershed. Trans. ASAE 18, 288-292. Overcash, M.R., K.R. Reddy and R. Khaleel, 1980. Chemical and biological processes influencing the modelling of nonpoint source water quality from land areas receiving wastes. Proc. Hydrol. Transport Modeling. ASAE pub. 4-80, St. Joseph, 1-11. 77 Penning de Vr ies , F.W.T. and H.H. van Laar, eds. , 1982. Simulation of plant growth and crop production. Center Agr. Pub. Docum. Wageningen. Priestley, L.H.B. and R.J. Taylor, 1972. On the assessment of surface heat flux and evaporation using large scale parameters. Mon. Weather Rev. 100, 81-92. Ragab, R.A. and F. Amer, 1986. Estimating water table contribution to the water supply of maize. Agric. Water Manag. 11, 221-230. Rawls, W.J. and O.L. Brakensiek, 1986. Comparison between Green-Ampt and curve number runoff predictions. Trans. ASAE, 29, 1597-1602. Ritchie, J .T . , 1972. Model for predicting evaporation from a row crop with incomplete cover. Water Resource Res., 8, 1204-1212. Rubin, J . , 1966. Theory of rainfall uptake by soils ini t ia l ly drier than their f ield capacity and its applications. Water Resource Res., 2, 739-749. Shirmohammadi, A . , and R.W. Skaggs, 1985. Predict ing infiltration for shallow water table soils with different surface covers. Trans. ASAE 28, 1829-1837. Sine, L. and J.P. Agneessens, 1980. Modelisation de la migration d'elements dans les sols. Pedologie, 30, 341-350. Skaggs, R.W., 1978. A water management model for shallow water table soi ls . Tech. Rep. 134, Water Resources Res. Inst., U. North Carolina, Raleigh. Stewart, D.W., L.W. Dwyer and R.L. Desjardins, 1985. The effect of available soil water and root density on actual and potential transpiration relationships. Can. Agric. Eng., 27, 7-11. Stewart, R.B., 1981. Modeling methodology for assessing crop production potentials in Canada. Agr. Can. Res. Branch, Tech. Bull. 96. Switzer-Howse, K.O., 1983. Agricultural management practices for improved water quality in the Canadian Great Lakes basin. Agr. Can., Res. Branch, pub. 1983-24-E. Thomas, D.L. and D.B. Beasley, 1986. A physically-based forest hydrology model. Trans. ASAE, 29, 962-981. Van Hoorn, J.W., 1981. Salt movement, leaching efficiency, and leaching requirement. Agric. Water Manag., 4, 409-428. 78 - Van Ommen, H . C , 1985. Systems approach to an unsaturated-saturated groundwater quality model, including adsorption, decomposition and bypass. Agric. Water Manag., 10, 193-203. - Van Wijk, A.L.M., and R.A. Feddes, 1986. Simulating effects of soil type and drainage on arable crop yield. Proc. Agric. Wat. Manag., Balkena, 97-112. - White, R.E., 1985a. The analysis of solute breakthrough curves to predict water redistribution during unsteady flow through undisturbed structured clay soi l . J . Hydrol. 79, 21-35. - White, R.E., 1985b. The influence of macropores on the transport of dissolved and suspended matter through so i l , in B.A. Stewart, ed., Advances in soil science, Springer Verlag, 95-120. - Wilkie, K.I., 1981. Nitrogen ferti l izer management with the aid of a soil nitrogen model. Unpub. Ph.D. Thesis, Techn. U. Nova Scotia, Halifax. -Wil l iams, J .R. , C A . Jones and P.T. Dyke, 1984. A modeling approach to determining the relationship between erosion and soil productivity. Trans. ASAE, 27, 129-144. - Williamson, R.E. and G.J. Kriz, 1970. Response of agricultural crops to f looding, depth of water table and soil gaseous composition. Trans. ASAE, 13, 216-220. - Wischmeir, W.H. and D.D. Smith, 1978. Predicting rainfall erosion losses: a guide to conservation planning. USDA Agric. Handbook 537. 79 CHAPTER 3 A field study of water and nutrient  movement in an agricultural soil under  drained and undrained conditions 3.1 Introduction In humid areas, drainage is commonly used to remove excess water and provide suitable conditions for crop growth and field trafficabi1ity. Subsurface drainage achieves this by promoting aeration of the root zone, lowering the water table, increasing percolation and improving soil strength. For this reason, it also tends to increase the leaching of nutrients, particularly nitrate, the formation of which requires aerated conditions and which is readily soluble and poorly adsorbed by the so i l . This effect has long been recognized (Walter et a l . , 1979), but the specific mechanisms that influence the loss of nitrates and other nutrients are d i f f icul t to isolate. The introduction of subsurface drainage in the soil-water-crop system interacts with numerous variables; the different conditions of weather, soil type, and crop management preclude any systematic comparisons between data from different work. The nitrate concentration of drainwater may vary from 2.8 mg/L to above 100 mg/L (Burke et a l . , 1974; Devitt et a l . , 1976), and that of orthophosphate from a trace to above 1.0 mg/L (Culley et a l . , 1983). Similarly, the response time of drainage concentration to a surface applied input varies greatly and may result from the action of different transport mechanisms; for instance, Iqbal and Warkentin (1981) found no change in drain 80 water following a surface application of manure, while Hergert et a l . (1981) found a rapid and transient response under the same conditions. Finally, comparative studies show that the effects of subsurface drainage on nutrient losses can range from an overall reduction (Bengtson et a l . , 1984) to an overall increase (Seuna and Kauppi, 1980). These conflicting results may be attributed to the fact that numerous transport mechanisms are affected by subsurface drainage. Drainage may affect the aeration status of the soi l , altering its nutrient composition. It may also change the hydrologic balance of the soi l , reducing surface runoff and increasing percolation, with a consequent increase of leaching and decrease in erosion. Finally, it may also provide an outlet for macropore flow, which bypasses the soil matrix. The existence of this mechanism has been suggested by Bottcher et a l . , (1981), who found unexpectedly high concentrations of a strongly adsorbed pesticide in drain water. Furthermore, subsurface drainage may favourize the formation of bracks, worm holes, and aggregation (Hundal et a l . , 1976), compounding this effect. This chapter presents the results of a field study undertaken to estimate the effects of subsurface drainage on nutrient losses for the conditions of the Lower Fraser Valley in British Columbia. The objectives were to determine the importance of nutrient losses in drainage water and to investigate the mechanisms by which subsurface drainage affects the movement of nutrients in this context. 81 3.2 Site and Methods The experimental site was located in the Lower Fraser Valley of British Columbia, an area characterized by flat lowlands and a humid climate. The soil at the test site is classified as a Humic Luvic Gleysol of the Ladner series, developed on fine textured (si l ty clay loam) deltaic deposits (Driehuyzen, 1983). Some physical and chemical characteristics of the soil are shown in Table 3-1. The test area, depicted in Fig. 3-1, consisted of three drained plots (A, B, C) and one undrained plot (D) of approximately 5000 m2 of surface area each. The drain lines consisted of 100 mm plastic corrugated tubing spaced 14 m apart and lying at a depth of 1.10 m below soil surface. Hog fuel (wood chips) was used as a f i l ter material in backfilling. A submersible pump was installed in the ditch to provide year-round drainage. The topography is flat and surface runoff is virtually non-existent. (No surface runoff was observed during the experimental period). Two of the drained plots (B and C) were subirrigated in the summer during periods of water deficit by introducing the water into the ditch and the drain tubes. Silage corn, potatoes, strawberries and hay (rye grass, timothy and clover) were grown on each plot. Table 3-2 details the fer t i l izat ion and timing schedule. The delay in cultivating plot D is due to lower trafficabi 1 ity on the undrained plot. The plots were established in 1982, prior to which no cultivation had taken place for about forty years (since the Second World War) and the area was overgrown with wild grasses. Details about the site, its management and its characteristics can be found 82 Table 3-1 Some Physical and Chemical Characteristics of the Soil at the Experimental Site (from Driehuyzen, 1983). Depth Sand Silt Clay Texture** B.D.+ pH* Carbon Total N m -% by wt.- kg/m3 — % by wt. — 0-0.3 3.5 73.3 '23.2 sil - 4.4 5.2 0.42 0.3-0.5 11.2 65.4 23.4 sil 1400 5.6 1.8 0.14 0.5-0.7 21.7 58.8 19.5 sil 1400 5.8 0.4 0.04 0.9-1.1 33.1 50.0 16.9 si 1-1 1370 4.5 0.4 0.03 ** Canadian system of soil classification + Bulk density 1:1 in water A PUMP LOCATION •4 SAMPLING OF DRAINS '00 metres X SAMPLING OF WELLS • ' Scale 1:2000 CO Figure 3-1: Experimental Field Layout 84 Table 3-2 Tillage and Fertilization of Test Plots PLOT CROP ACTIVITY DATES A,B,C Forage 35 kg N/ha 29/3/84, 14/8/84, 23/7/85 5/6/84* 5/6/85, D Forage 35 kg N/ha 5/6/84, 5/8/84, 23/7/85 11/7/84,* 11/5/85, A,B,C Forage 16 kg P/ha 5/6/84, 13/6/85 D Forage 16 kg P/ha 5/6/84, 13/6/85 A,B,C Corn Tillage 30/5/84, 9/5/85 0 Corn Ti1lage 21/6/84, 14/5/85 A,B,C Corn 69 kg 24 kg N/ha, 26 kg P/ha N/ha, 10 kg P/ha 31/5/84 13/5/85 D Corn 69 kg 24 kg N/ha, 15 kg P/ha N/ha, 10 kg P/ha 22/6/84 27/5/85 A.B.C Potatoes 69 kg 29 kg N/ha, 22 kg P/ha N/ha, 12 kg P/ha 31/5/84 5/6/85 0 Potatoes 69 kg 29 kg N/ha, 20 kg P/ha N/ha, 12 kg P/ha 22/5/84 5/6/85 A,B,C Potatoes Tillage 30/5/84, 2/5/85 D Potatoes Tillage 21/6/84, 24/5/85 *35 K N/ha applied at each date shown. 85 in Driehuyzen (1983), Chieng (1987), and Chieng et a l . , (1987). Piezometers and observation wells were installed in each plot to record water table height and allow sampling of groundwater. No barrier was installed between the plots. Sample collection took place from April 1984 to May 1986, taken twice a month on average. A vacuum hand pump was used to extract samples (500mls) from a depth of approximately 1.2 m. Grab samples of drain water were collected from drain outlets in plots A and B, but not from plot C as the outlets are partially submerged. The unfiltered samples were then frozen without delay and kept frozen until analysis. Concentrations of nitrate, ammonium and orthophosphate were measured with a Technicon autoanalyzer using reduction by hydrazine su l fa te , alkal ine phenol and ammonium molybdate methods, respectively (Technicon Industrial Systems, 1969, 1973a, 1973b). Acid digestion was performed on the samples using a micro-Kjeldahl technique and the digested samples were analysed for their content of total Kjeldahl nitrogen (excluding nitrate) and total phosphorus on the autoanalyser using alkaline phenol and ammonium molybdate methods (Technicon Industrial Systems, 1971, 1972). Drain flow rates were measured with a graduated cylinder and a stopwatch. The water table depth was recorded on a continuous basis. 3.3 Results Figure 3-2 shows the water table profile of plot A (drained) and plot D (undrained). As expected, the water table of undrained plot was markedly higher during the wet season. This is of WATER TABLE PROFILE oo Figure 3-2: Precipitation (PREC) and Water Table Depth Below the Surface Versus Time for Drained and Undrained Plots 87 particular importance in the spring, during which the water table of plot D remained high and delayed tillage and planting operations. Surface ponding in plot D was also common in winter, while it was infrequent in the other plots. Drain flow rates were measured during sample collection. There was a good correlation between flow rate and water table depth, as expected. However, flow rates from individual outlets differed. These differences most likely originated from spatial variability of the soil hydraulic conductivity; similar variations have been reported by other workers (Guitjens et al . , 1984). Further, there was a consistent difference between plot A and plot B. Both flow rates and water table depth measurements indicated that plot A responded faster than plot B. A similar water table depth produced a higher drain flow rate in plot A. Similarly, the water table of plot A rose more rapidly during rainfall and dropped faster following rainfall than that of plot B. It was suspected that the cause of this difference was the summer subirrigation in plots B and C, but more field trials are needed to confirm this. Nitrate, ammonium, and orthophosphate concentrations of the unfiltered samples are shown in Fig. 3-3, 3-4, and 3-5. In each figure the nutrient concentrations of drainwater, groundwater from the drained plots (averages of plots A and B), and groundwater from the undrained plots are depicted. It can be seen that nutrient concentrations were higher, in general, towards the beginning of the data collection period. This trend was observed in all cases. The arithmetic average and range of these values and of total NITRATE NITROGEN CONCENTRATION o - ORAIN WATER (DRAINED PLOTS) o - GROUNDWATER (DRAINED PLOTS) A - GROUNDWATER (UNDRAINED PLOTS) Figure 3-3: Nitrate Nitrogen Concentration 1n Unfiltered Samples Versus Time AMMONIUM NITROGEN CONCENTRATION • - DRAIN WATER (DRAINED PLOTS) o - GROUNDWATER (ORAINED PLOTS) A - GROUNDWATER (UNDRfllNEO PLOTS) i i r i f i i i i i i i i i i t i i i i i i i i i i t i i J F M f l r l J J A S O N D J F M R M J J A S O N D J F M A r l J 1984 TIME OF THE YEAR 1986 Figure 3-4: Ammonium Nitrogen Concentration 1n Unfiltered Samples Versus Time PHOSPHATE PHOSPHORUS CONCENTRATION J F M f l M J J f l S O N D J F M f i M 1984 J J f l S O N D J F M' fl' M' j ' TIME OF THE YEAR 1986 KO o Figure 3-5: Phosphate Phosphorus Concentration 1n Unfiltered Samples Versus Time 91 Kjeldahl nitrogen and total phosphorus are given in Table 3-3. Differences in the nutrient levels from the drain water, and the groundwater from both plots were assessed using a Wilcoxon signed-rank test for paired data. The results are shown in Table 3-4. Figure 3-3 shows that nitrate-nitrogen concentrations in the undrained plot were generally low, especially during spring and summer and particularly so in the last year. In contrast, N03-N concentration in drained plots during the course of the sampling period was relatively high, occasionally exceeding 10 mg/L N03-N in the early stages. The results also show that ground water from the drained plot contained higher nitrate concentrations than drainwater, except during the last winter of sampling, where the reverse was true. From Table 3-4 and Figure 3-3, it may be seen that groundwater from the drained plots and drainwater were significantly more concentrated in nitrate than groundwater from the undrained plot. The low level of significant difference between drainwater and groundwater from the drained plots masks the fact that groundwater was significantly less concentrated in N03-N in the 85-86 winter, and more concentrated otherwise. Ammonium nitrogen concentrations were relatively low (less than 2 mg/L) in both drainwater and groundwater from drained and undrained plots (Figure 3-4). Drainwater was found to contain significantly less NH4-N than groundwater from either drained or undrained plots. However, no significant difference in the NH4-N concentration in groundwater from drained and undrained plots were found. There were no significant differences in orthophosphate 92 Table 3-3 Average and range of sample concentrations, in mg/1. Nutrient Source Average Minimum Maximum N03* Drainwater 5.481 0.645 9.880 Drained plot 6.595 0.029 15.900 Undrained plot 2.236 0.000 10.610 NH4* Drainwater 0.321 0.070 1.650 Drained plot 0.481 0.140 1.760 Undrained plot 0.513 0.150 1.950 P04* Drainwater 0.069 0.000 1.330 Drained plot 0.082 0.000 0.920 Undrained plot 0.091 0.000 1.380 TKN** Drainwater 1.464 0.778 1.875 Drained plot 2.648 1.375 4.776 Undrained plot 2.487 1.806 3.122 TP** Drainwater 0.625 0.082 3.870 Drained plot 0.366 0.069 2.296 Undrained plot 0.453 0.194 1.593 * + Average of 52 samples Average of 14 samples, taken from Oct. 85 to May 86 Table 3-4 Comparison Between Nutrient Content of Drainwater and Well Water in Drained and Undrained Fields, Using the Wilcoxon Matched Pairs Test Nutrient Source * * Comparison N03-N NH4-N P04-P TKN TP DD vs GD p* < .0135 p < .001 N.S.* p < .001 N.S. DD vs GU p < .001 p < .001 N.S. p < .001 N.S. GD vs GU p < .001 N.S. N.S. N.S. N.S. * * DD = Drainwater of drained plots GD = Groundwater of drained plots GU = Groundwater of undrained plots * N.S. = No significant difference at 0.05 level + p < x = significantly different at level x 94 levels from the three sources, and orthophosphate levels were generally low. Higher levels of Ortho-P in ground and drain waters were occasionally encountered during the spring of 1984. Table 3-3 shows the average values and ranges of total Kjeldahl nitrogen and total phosphorus measured in the winter of 1985-86. The concentrations from all sources were relatively high, indicating a considerable contribution from organic forms of nitrogen and phosphorus. As with ammonium, TKN concentrations were significantly smaller in drainwater than in groundwater. There was no significant difference between all three sources for total phosphorus. A few samples (10% of the total) had concentrations higher than 1 mg/L. 3.4 Discussion One clear trend revealed by the data is the lowering of nutrient concentrations from the beginning of the data collection period to the end. This trend may arise from the transient effects of drainage, t i l lage, and cultivation. Draining and cultivating a field formerly untended or pastured causes rapid mineralization of the accumulated organic matter, which results in high nutrient levels in soil water for a few years following cultivation. This temporary effect has been observed in a variety of settings (Bergstrom, 1987; Roberts et a l . , 1986; Cameron and Wild, 1984; Devitt et a l . , 1976). Thus the nitrate concentration of drain water shows an overall decrease from the beginning to the end of the experimental period. 95 There is an increase in the f a l l , for both years, which may be caused by an init ial flushing of the nitrate accumulating in later summer and early fa l l from n i t r i f i ca t ion , as pointed out by Kowalenko (1987). However, the nitrate concentration of drainwater remains noticeable through the winter, which shows that the nitrate is not completely leached out at the end of the season. This suggests that there is some amount of preferential flow, or bypassing flow, in this s o i l , as more than one pore volume percolates through the soil during the winter. (One pore volume, to drain depth, was estimated to be equal to 0.5m, assuming all the soil water to be displaced. The winter precipitation, as seen in Figure 3-2, was obtained from a nearby weather station (Vancouver International Airport) and totalled, from September to Apri l , 823mm for 84-85 and 928mm for 85-86). During the period from April 84 to June 85, the nitrate concentration in the well water in the drained f i e ld was consistently higher than that of the drainwater, while the reverse was found during the following winter. The nitrate concentration of well water in plot D was always lower than in drainwater. Denitrification may be responsible for the low concentrations of plot D. Heavy mottling of the soil and a strong sulfurous odour emanating from the observation well of Plot D in the spring attest to a reducing and denitrifying environment. Gambrell et a l . (1975) and Devitt et a l . (1976) showed that denitrif ication is the main mechanism by which nitrate levels are lower in soils with impeded drainage. In these studies the comparison was made between 96 different soil types, where the heavier texture of some soils impeded drainage. Results from the present work would suggest that this situation may also arise as a sole effect of a subsurface drainage. Results from the wells in the drained plots are more puzzling. If the soil is homogeneous and the soil water is displaced in piston flow fashion, one would expect to see l i t t le difference between drain water and well water. This is clearly not the case, and this situation, similar to what has been experienced elsewhere (Karlen et a l . , 1976), can best be explained by the presence of bypassing flow. Observations of the soil surface and profile reveal large, deep cracks that do not completely close in winter. Further, there are numerous root channels, including many with concretions. Bypassing flow effectively divides soil water into two zones, one with rapidly flowing water held in a network of macropores, and the other where peds or aggregates hold the water relatively immobile. Drain water is usually composed of water from both zones, and its composition reflects the relative rate of percolation of each zone, and the dif fusion between zones (Armstrong and Arrowsmith, 1986). Thus a given amount of water percolating rapidly through large pores will leach fewer nutrients held in aggregates than the same amount applied slowly. This dependence of nutrient movement on flow rates was observed by Nagpal and Wiens (1984). In contrast to drains, observation wells sample a relatively small volume of soil and may or may not be connected to the network 97 of macropores. Differences between drainwater and well water nutrient concentrations are therefore to be expected. If well water consists mostly of ped water, that i s , water that is held relatively immobile in the micropores, one would expect well water to be relatively rich in nutrients. Also because of the relatively low permeability of micropore zone, there is l i t t le dilution from rainwater in ped water or well water as opposed to drainwater. Thus in this case, one may expect well water to show a higher concentration of ammonium and organic nitrogen, and of total orthophosphate. Nitrate levels in well water may or may not be higher than in drainwater, and will depend on the aeration status of well water. This seems to be the case in this experiment, at least for nitrogen. Ammonium and Kjeldahl nitrogen were consistently and significantly higher in well water than in drainwater over the sample collection period. Nitrate levels were also significantly higher in well water from the beginning of data collection to the spring of 1985; however, during the following winter, these levels were consistently lower than in drainwater. These systematic differences showed that drainwater and well water do not sample the same flow components. One possible hypothesis to explain the change of response of the wells after the summer of 1985 is that there may have been a change in the source of water sampled by these wells. The spring and summer of 1985 was much drier than the previous one (107mm of rain from May to August, compared to 210mm in 1984) and this may have caused a movement of fine soil 98 particles, or some other mechanism through which the wells may have become partially sealed or isolated from the macropore network. This would have altered the aeration status of the well water. The presence of a reducing environment in the wells of the drained plots was indicated by the smell of hydrogen sulfide from the well, which was first encountered in the winter of 1985-86. This implies that water from wells may be an unreliable indicator of the overall status of the soil water, since it samples a small volume of soil and may not adequately reflect the potentially different composition of macropore and micropore water. However, this may be a problem only if macropore and micropore water have a dramatically different residence time- in the soil, which is less likely to occur in a soil with impeded drainage, as in plot D. Nevertheless, the validity of conclusions derived from comparing well water from the undrained plot to drain water remains open to question; however, it is clear that the water reaching the drains does not originate from a reducing environment, while reduced conditions prevail in plot D, particularly during the spring, when a combination of high water tables (and flooding) and warm temperatures is found. There was no significant difference in orthophosphate (Figure 3-5) or total phosphorus (not shown) levels between drained and undrained plots. However, there were several instances, particularly in the spring of 1984, when drainwater or well water had unexpectedly high concentration of phosphorus. These are higher than what should be found after percolation through the mass 99 of the soil, since the soil has a large adsorption capacity for phosphate. Such concentrations may be further indications of bypassing flow in this soil. The importance of bypassing flow as a nutrient transfer mechanism has several implications. Fertilizer applied before a rain storm may be washed down to the drains through the macropore network. This may be the explanation for some of the high phosphorus values encountered. Conversely, since intrapedal flow and nutrient diffusion to macropores are fairly slow mechanisms, bypassing flow may result in less leaching, or leaching at a slower rate, than predicted by the theory of mass flow through homogeneous medium, if the nutrients are held in the micropores. Figure 3-6 shows the estimated daily nitrate loss from leaching to drainwater. This figure was constructed by calculating the daily drain flow as a function of the water table depth, using Hooghoudt's formula, and interpolating nitrate concentration linearly between data from sample collection. As expected, the figure reflects the decreasing nitrate levels over time. It also shows that nitrate leaching is largest in the fall, at the onset of the rainy season, but also that leaching is important throughout the winter, which indicates that there is residual nitrate in the soil in the spring. Thus the assumption, common in this area, of complete nitrate removal by winter leaching may be erroneous and lead to fertilization practices at above optimum levels. NITRRTE NITROGEN LOSS Figure 3-6: Estimated Nitrate Nitrogen Drain Loss Versus Time 101 3.5 Conclusions Short of performing tracer studies, it is difficult to provide definite evidence for specific mechanisms affecting nutrient movements and transformations at the present time. However, routine sampling of groundwater and drainwater as reported in this chapter can provide circumstantial evidence of the effects of some of these mechanisms. The conclusions that can be derived from the experimental data are as follows: • The transient situation following the disturbances caused by reclamation, drain installation and tillage results in a large amount of nutrient leaching. Transient effects remain evident three years after the initial installation. • There is indirect evidence that bypassing flow is an important mechanism which governs the loss of nutrients to drain water in this f ie ld. • Nitrate levels in drain water remain relatively constant over winter after an init ial flush in the f a l l . This indicates that there is residual nitrate in the soil in spring, as a possible result of both transient effects and bypassing flow. • There is indirect evidence that subsurface drainage may alter the aeration and redox status of the so i l , increasing the nitrate concentration of the water percolating to the drains. 102 3.6 References - Armstrong, A.C. and R. Arrowsmith, 1986. Field evidence for a bi-porous soil water regime in clay soils. Agric. Water Mgmt., 11, 117-125. - Bengtson, R.L., C E . Carter, H.F. Morris, and J.G. Kowalczuk, 1984. Reducing water pollution with subsurface drainage. Trans. ASAE, 27, 80-83. - Bergstrom, L., 1987. Nitrate leaching from annual and perennial crops in t i le drained plots and lysimeters. J . Env. Qual., 16, 11-18. - Bottcher, A.B. , E.J. Monke, and L.F. Huggins, 1981. Nutrient and sediment loadings from a subsurface drainage system. Trans. ASAE, 24, 1221-1226. - Burke, W., J . Mulqueen and P. Butler, 1974. Fertilizer losses in drainage water from a surface water gley. Ir. J . Agric. Res., 13, 203-214. - Cameron, K.C and A. Wild, 1984. Potential aquifer pollution from nitrate leaching following the plowing of temporary grassland. J . Env. Qual., 13, 274-278. - Chieng, S.T. , J . Keng and M.G. Driehuyzen, 1987. Effects of subsurface drainage and subirrigation on the yields of four crops. Can. Agric. Eng. 29, 21-26. - Chieng, S.T., 1987. Subirrigation and soil hydraulic properties. CSAE paper 87-308, CSAE, Ottawa, Ont. - Culley, J . L . B . , E.F. Bolton and V. Bernyk, 1983. Suspended solids and phosphorus loads from a clay soi l . J . Env. Qual., 12, 493-503. - Devitt, D., L. Letey, L.J. Lund and J.W. Blair, 1976. Nitrate-n i t rogen through s o i l as a f f e c t e d by s o i l p r o f i l e characteristics. J . Env. Qual., 5, 283-288. - Driehuyzen, M.G., 1983. Boundary Bay water control project, Unpub. report, B.C. Min. Agric. Fish. , Cloverdale, B.C. - G a m b r e l l , R . P . , J.W. G i l l i a m and S . B . Weed, 1975. Denitrification in subsoils of the North Carolina coastal plain as affected by soil drainage. J . Env. Qual., 4, 311-316. - Guitjens, J . C , P.S. Tsui and D.F. Thran, 1984. Quantity and quality variations in subsurface drainage. Trans. ASAE, 27, 425-428. 103 Hergert, G.W., D.R. Bouldin, S.D. Klausner, and P.J. Zwerman, 1981. Phosphorus concentration-water flow interactions in t i le effluent from manured land. J . Environ. Qual., 10, 338-344. Hundal, S.S. , G.O. Schwab and G.S. Taylor, 1976. Drainage system effects on physical properties of a lakebed clay so i l . Soil Sc. Soc. Am. J . , 40, 300-306. Iqbal, M.M., and B.P. Warkentin, 1981. Nitrogen and phosphorus content of flow from drains fertilized with cow manure slurry. Plant and Soi l , 63, 517-521. Karlen, D.L., M.L. Vitosh, and R.J. Kunze, 1976. Irrigation of corn with simulated municipal sewage effluent. J . Environ. Qual., 5, 269-273. Kowalenko, C.G., 1987. The dynamics of inorganic nitrogen in a Fraser Valley soil with and without spring or fa l l ammonium nitrate applications. Can. J . Soil Sc. 67, 367-382. Nagpal, N.K. and J.H. Wiens, 1984. Shawnigan Lake watershed study: Investigations of soil water quality below septic tank drain field and nutrient loading to the lake. Unpub. Report., Min. Env. and Parks, Victoria, B.C. Roberts, G. , J.A. Hudson and J.R. Blackie, 1986. Effect of upland pasture improvement on nutrient release in flows from a natural lysimeter and a field drain. Agric. Water Mgmt., 11, 231-245. Seuna, P. and L. Kauppi, 1980. Influence of sub-drainage on water quality and quantity in a cultivated area in Finland. In: The influence of man on the hydrological regime, Symp. Proc , IAHS pub. 130, 177-187. Technicon Industrial Systems, 1969. Nitrate and nitr i te in water. Indust. meth. 33-69W. Technicon I.S., Tarrytown, N.Y. Technicon Industrial Systems, 1971. Total nitrogen (Kjeldahl). Indust. Meth. 103-70A. Technicon I.S., Tarrytown, N.Y. Technicon Industiral Systems, 1971. Total phosphorus. Indus. Meth. 116-71W. Technicon I.S., Tarrytown, N.Y. Technicon Industrial Systems, 1973a. Ammonia in water and seawater. Indust. meth. 154-71W. Technicon I.S., Tarrytown, N.Y. 104 - Technicon Industrial Systems, 1973b. Orthophosphate in water and seawater. Indust. meth. 155-71W. Technicon I.S., Tarrytown, N.Y. -Walter, M.F., T.S. Steenhuis and D.A. Haith, 1979. Non-point source pol lut ion control by soi l and water conservation practices. Trans. ASAE, 28, 834-840. 105 CHAPTER 4 A comparison of simulated and measured water and  nutrient movement in drained and undrained  agricultural soils 4.1 Introduction Pollution of water from fertilizers and land degradation are two growing concerns associated with agriculture. These two problems arise from situations that are most affected by the movement of water in so i l s . Soil water flow controls such mechanisms as surface runoff, water erosion, or nutrient leaching, which may result in eutrophication, loss of land productivity, and pollution of surface and ground waters. Subsurface drainage alters the patterns of soil water flow, but its effects on nutrient movement and nutrient pollution remain largely unknown (Switzer-Howse, 1983). Subsurface drainage promotes soil aeration, increasing oxidation of organic nutrients and nitrification. It also tends to increase the rate of water flow (and leaching) through the s o i l , while possibly reducing surface runoff. It also provides an outlet for bypassing flow occurring in soil cracks or other macropores. Since it provides better conditions for crop root development, nutrient uptake may increase, counteracting the effects of leaching. Conversely, drainage is often associated with rowcrop cultivations, which are heavily fertilized and offer less protection against erosion. The interplay of these mechanisms causes field results to appear contradictory, as they depend on site specific conditions. A survey of the literature (see Chapter One) reveals that nitrate 106 and phosphate levels in drain water may reach values above 100 mg/L and 1 mg/L, respectively (Devitt et a l . , 1976; Culley et a l . , 1983), while remaining below detection levels in other situations (Burke et a l . , 1974). When compared with surface drained plots, plots equipped with subsurface drainage are found to exhibit a variety of responses, ranging from a reduction in nitrogen and phosphorus losses (Bengtson et a l . , 1984) to an increase in nitrogen losses (Harris et a l . , 1984). Simulating water and nutrient flows in the soil makes it possible to identify the conditions under which subsurface drainage increases or decreases nutrient losses and pollution. Such an approach was used by Deal et a l . (1986). Their simulation showed that, for most of the soils they modeled, nitrogen losses should increase but phosphorus losses should decrease as a result of subsurface drainage. Unfortunately, their model did not include bypassing flow as a transfer mechanism, nor were advection, crop growth, or chemical and biological nutrient transformations. Consequently, the applicability of the model and its aptitude at simulating different situations is somewhat limited. This chapter presents the results of the simulation of the water and nutrient movements in an experimental field where the responses of drained and undrained plots are monitored. The computer model consists of a simulation of a comprehensive set of mechanisms that affect nutrient and water status in humid, temperate farm soils. 107 4.2 Methods The results from the computer model are to be compared with measured values obtained from an experimental site located at Boundary Bay in the South West Lower Fraser Valley of British Columbia. The climate is moderate and humid, with relatively dry summers. The research plots were established in 1982 on an abandoned site after clearing and installing subsurface drainage on some of the plots. During the data collection period of 1984-1985, corn, grass, potatoes and strawberries were grown on the s i l t clay loam gleysol plots. Water table levels were recorded continuously and grab samples of drain water and of water from observation wells were taken per iodica l ly . More details of the experimental conditions are discussed elsewhere (see Chapter Three). A summary of the computer model (described in Chapter Two) is presented here to clarify the calibration procedure. The computer model consists of a water flow submodel, a nutrient submodel, and a crop growth submodel. The crop growth part consists of a simplified simulation of leaf cover growth and root extension as a function of transpiration. Results from this submodel are used to calculate evapotranspiration and nutrient uptake. The nutrient submodel simulates the rate of oxidation of organic matter, mineralization of organic nitrogen and phosphorus, immobilization of nitrogen and phosphorus by microbial biomass, nitrif ication, denitrification, and ammonia volatilization. These processes are calculated according to f i r s t order kinetics, modulated by soil water content and temperature. Erosion is 108 calculated with a modified version of the universal soil loss equation. In order to calculate the transfer of heat, nutrients, and water, the soil is divided into discrete horizontal layers. The nutrients move between layers by water advection. Adsorption is calculated using a Langmuir isotherm. The water flow submodel computes water movement between the layers, which is then used in the nutrient submodel. The water flow submodel uses daily weather data to simulate potential and actual evapotranspiration, and snowmelt. Drain flow and seepage are calculated based on the depth of the water table. The model uses an original infiltration-redistribution routine that simulates bypassing flow and ped flow as parallel mechanisms. Rainfall and snowmelt are distributed between ped infiltration and bypassing flow down soil macropores to the water table. The concept of flux equilibrium outlined in Chapter Two is used to calculate redistribution in the ped as well as upward capillary flow, both of which are strongly dependent on the water table depth and on the moisture retention curve. The permeability of the soil macropores is set at a constant value, CK, while the permeability of the micropore domain, C, is calculated as C = fJS • (9/9PR)A L F A 4-1 where CS is the permeability when the micropores are water saturated, which happens when the water content of a layer, 9, reaches a value of 9PR. The water holding capacity of the 109 macropores is given by the difference between total saturation (6S) and 6PR. Since they are very difficult to measure in the f ield, CS, 6PR, and the exponent ALFA are treated here as fitted parameters. Lateral transfer between macropores and micropores arises from horizontal infiltration and diffusional mixing between the two domains. Horizontal infiltration is calculated using the Green and Ampt model while diffusion is assumed to take place over a finite horizontal distance from the cracks. The extent of this lateral transfer is controlled by the spacing index PEDS, which represents the effective distance between the macropore channels. PEDS must also be estimated as a fitted parameter. (As it has no effect on the hydrologic functions this parameter cannot be considered to represent a measurable physical distance between actual cracks in a well defined geometry.) Since bypassing flow and ped flow are parallel, the sum of CS and CK represents the saturated hydraulic conductivity of the soi l . This value was estimated from data obtained for the experimental site by Chao (1987) and Dryiehuyzen (1983). Drainflow rate measurements were used to estimate the hydraulic conductivity below drain depth. The moisture retention curves for this soil were determined by Gibb (1987) and Chao (1987); figure 4-1 shows the values used for this simulation. The simulation was carried out for the period of January 1984 to December 1985, using the weather data of nearby Vancouver Airport meteorological station. The weather during 1984 was 110 MOISTURE RETENTION CURVE .10 0.15 0.45 0.50 MOISTURE CONTENT (M/M) Figure 4-1: Moisture Characteristics of the Soil Used 1n the Simulation I l l considerably wetter than in 1985 (1280 mm of precipitation versus 800 mm, respectively), and slightly cooler. Field management data (fertilization, operations timing) were taken from field records. Crop growth was simulated using parameters for corn calculated from data found in Stoskopf (1981). Nutrient kinetic coefficients were obtained from Overcash et a l . (1980) and Johnsson et a l . (1987). The values of the main model parameters can be found in Table 4-1, along with the values selected for the fitted parameters. The model was calibrated by adjusting the values of CS, ALFA, 8PR, and PEDS to obtain an adequate f i t of the water table profile for the period of January 1st to October 1st, 1984. The quality of the f i t at that point was determined by visual inspection. After calibration, the f i t obtained for the subsequent period was inspected to determine the extent to which the water flow submodel can be considered validated. After this process of calibration and partial validation, the quality of the f i t was assessed using the Wilcoxon matched pairs rank test and the Spearman correlation coefficient (Siege!, 1956). The nutrient content of the soil at the beginning of the simulation was unavailable and consequently the init ial concentrations of nitrate, ammonium and phosphate were treated as fitted parameters. These nutrients were assured to be distributed uniformly through the profile down to drain depth. These values were adjusted in order to obtain a reasonable f i t of the observed nutrient concentration of drainwater. For both water flow and nutrient flow, no attempt was made at an optimisation of the f i t because of the excessive requirements in computer time. 112 Table 4-1 Summary of Simulation Parameters Fixed Parameters Value Drain depth 1.1 m Drain spacing 14 m Equivalent depth to impermeable layer 2.1 m Overall hydraulic conductivity (0-0.7 m) 0.5 m/d Overall hydraulic conductivity (below .7 m) 0.03 m/d Depression storage 0.02 m Effective rooting depth 1.1 m Mineralization rate from stable organic matter 0.00001 per day Mineralization rate from unstable organic matter 0.0001 per day Nitrification rate 0.15 per day Thickness of layers (cms) 1,9,10,10,20,20 20,20,100,1000 Fitted Parameters i Ped conductivity CS 0.001 m/d Conductivity exponent ALFA 10.0 Macropore volume 9PR 0.03 m/m Equivalent ped spacing PEDS 5.0 m Initial soil nitrate (above drains) 0.014 tonne/ha Initial soil ammonium (above drains) 0.043 tonne/ha Initial soil phosphate (above drains) 0.028 tonne/ha 113 4.3 Results Predicted and actual water table depths are shown in Fig. 4-2 for the drained plot, and in Fig. 4-3 for the undrained plot. In general there is a very good agreement between predicted and actual values. Predicted and actual water table heights (daily average) were compared for statistical significance. For the drained case, the Spearman correlation coefficient (rho) is equal to 0.917 and there are no significant differences (p > 0.05) between the simulated and observed values using the Wilcoxon matched pairs test. Although the simulated and observed values of the water table height for the undrained case were found to be significantly different (p < 0.05, using the Wilcoxon test), there was a very good correlation between the two, the Spearman's rho being equal to 0.924. Fig. 4-4 shows the predicted nitrate concentration of drainage water along with actual concentrations from grab samples. The Wilcoxon matched pairs test reveals no significant differences between actual and predicted values. The Spearman correlation coefficient is low (rho = .584), while significant (p < 0.001). The predicted ammonium levels in drain water are in general too low by one order of magnitude, and are consequently significantly different at the .001 level. Further, no significant correlation can be found. Predicted levels of phosphorus are on average higher than actual, but there is no significant difference between the two, using the Wilcoxon test at the 0.001 level. Again, there is no significant correlation between predicted and actual values. WATER TABLE PROFILE M f l M J J R S O N D 1984 J F M f l M J J f l S O N TIME OP THE YEAR 1985 Figure 4-2: Simulated and Observed Water Table Profile of the Drained Plot WATER TABLE PROFILE J F M fl M j J flSONDJFMRHJJflSOND 198* TIME OP THE YEAR 1985 Figure 4-3: Simulated and Observed Water Table Profile of the Undrained Plot NITRRTE CONCENTRATION OBSERVED VALUES SIMULATED VALUES NO DRAIN FLOW Figure 4-4: Simulated and Observed Nitrate Concentration of Drainage Water 117 In order to compare the importance of the various water flow mechanisms on nutrient losses, a sensitivity analysis was carried out on the fitted parameters CS, 9PR, PEDS, and ALFA. These parameters were increased or reduced in value, which resulted in different values of water movement (evaporation, transpiration, drainflow) and nutrient losses by drainflow over the two year simulation. Table 4-2 shows the results expressed as a ratio of the values obtained using the original parameters (given in Table 4-1). The results obtained after only one year (1984) of simulation are presented in Table 4-3. The contrast between some of the values, particularly nutrient losses, il lustrates the difference in response from one year to the next. Figures 4-5, 4-6, and 4-7 show the profile of the time response of the water table when CS, 9PR, and PEDS are varied. The pertinence of the parallel flow approach used in this model can be evaluated by comparing the performance of the model with the results obtained using two other models. In the f irst of these two the modeling approach used the assumption of homogeneous flow through the so i l . In this case there is no bypassing flow and the soil behaves as a homogeneous medium with a high hydraulic conductivity, equal to the value given by the addition of CS and CK in the previous model. This approach is similar to that used in the models of Duffy et a l . (1975) or Kanwar et a l . (1983). In the second case the approach used is the opposite to the former, in that vertical flow is assumed to only take place in channels and macropores, and thus this model uses only bypassing flow. 118 Table 4-2 Effects of Parameter Variation on Simulation Results, Expressed as a Ratio, for 1984-85 Variable Parameter ALFA PEDS (m) CS (m/d) 6PR 7.5 13.5 900 0.2 .005 .002 .06 .015 Evaporation 1.03 0.98 1.00 1.00 1.04 0.93 1.00 1.00 Transpiration 1.03 0.99 1.00 1.01 1.04 0.95 1.01 1.00 Drainflow 0.98 1.01 1.00 0.99 0.96 1.06 1.01 0.99 Drain losses of NH4 0.89 1.13 0.82 1.58 1.46 0.84 1.05 1.04 N03 0.97 0.99 0.87 1.30 1.35 0.72 1.18 0.93 P0 4 0.98 1.05 0.98 1.16 1.19 0.99 1.48 0.91 119 Table 4-3 Effects of Parameter Variation on Simulation Results, Expressed as a Ratio, for 1984 only Variable Parameters ALFA PEDS (m) CS (m/d) 9PR .75 13.5 900 0.2 .005 .0002 .06 .015 Evaporating 1.02 0.99 1 .00 1.00 1 .04 0. .96 0. 98 Transpiration 1.02 1.00 1 .00 1.00 1 .00 1. .00 1. 00 Drainflow 0.99 1.00 1 .00 0.99 0 .98 1. .01 1. 03 Drain losses of NH4 0.94 1.08 0 .81 1.76 1 .67' 0. .73 1. 15 N03 0.97 1.02 0 .83 1.50 1 .54 0. .69 1. 24 P0 4 1.00 1.00 0 .91 1.36 1 .44 0, .87 1. 50 1.01 1.00 0.98 0.96 0.90 0.86 WATER TABLE PROFILE CONDUCTIVITY CS - 5 .0 HM/D CONDUCTIVITY CS - 0 .2 MM/0 Figure 4-5: The effects of varying the Matrix conductivity CS on the simulated water table profile WATER TABLE PROFILE Figure 4-6: The Effects of Varying the Macropore Porosity OPR on the Simulated Water Table Profile WATER TABLE PROFILE Figure 4 -7 : The Effects of Varying the Equivalent Ped Spacing PEDS on the Simulated Water Table Profile 123 Nutrients move Into or out of the matrix through horizontal diffusion only. This approach is akin to that used by Addiscott (1981). The two models used in this comparison (referred to here as the homogeneous flow model and bypassing flow only model) were developed by altering the present model. In the f irst case, the bypassing flow component was removed. This was accomplished by setting the porosity and conductivity of the macropore domain to zero. The matrix, which now becomes the sole flow domain, is assigned a saturated conductivity CS equal to the overall saturated conductivity of the soi l . The unsaturated conductivity parameter ALFA is then recalibrated, while the moisture retention curve is left unaltered. The saturated hydraulic conductivity and the moisture retention curve remain unchanged since they are measured parameters, leaving only ALFA to be evaluated as a f i t ted parameter. (Clearly, however, a soil dominated by bypassing flow and a soil where homogeneous flow is prevalent would have a very di f ferent s t ructure , conduct ivi ty , and moisture retention characteristics. This procedure, then, may be seen as testing whether a computer model based on the assumption of homogeneous flow can produce a sat isfactory f i t of the data, while acknowledging that such an assumption may be inappropriate and contradict the physical reality of the soil modelled). In the second case, a situation where water can only move in the macropores 1s created simply by setting CS, the matrix conductivity, to zero. In this case, the shape of the moisture retention curve becomes irrelevant, as there is no moisture flow regardless of the gradient. (Again, this model represents a departure from the physical reality of the system.) Fig. 4-8 shows 124 the simulated water table profile resulting from these two models, while the predicted nitrate concentration in drain flow can be seen in f ig . 4-9. Table 4-4 summarizes the results. Since it was established by Kanwar and Johnson (1984) and Caussade et a l . (1984) that changes to the water regime cause changes in nutrient movement of far greater magnitude than changes in biochemical kinetic coefficients do, no sensitivity analysis on these coefficients was performed. Surface runoff was predicted to occur only once (on January 3rd, 1984) in all simulations of the drained f ie ld. Since no surface runoff was actually observed during the simulation period, this occurrence is probably due to a failure to accurately estimate the soil moisture content at the beginning of the simulation runs. Therefore the simulated results for surface runoff may not be meaningful and are not discussed further. Because of the simplicity of its design, the model was found to require only an average of five CPU seconds for two year simulation run on batch mode of the University's Amdahl 470 V/8 computer. WATER TABLE PROFILE HOMOGENEOUS PLOW MODEL BYPASSING FLOW ONLY MODEL in ri' J F M R M J J f l S O N D J F M f l M J J R S O N D 1984 TIME OF THE YERR 1985 Figure 4-8: Simulating Hater Table Movement Using a Simple Bypassing Flow Model or a Simple Homogeneous Flow Model NITRATE CONCENTRATION J F M R M J J R S O N D J F M f i M J J R S O N D 1984 TIME OF THE YERR 1985 Figure 4-9: Simulating the Nitrate Content of Ordlnage Water Using a Simple Bypassing Flow Model or a Simple Homogeneous Flow Model Table 4-4 Effects of Using a Simple Bypassing Flow Model or a Simple Homogeneous Flow Model on the Simulation Results, Expressed as a Ratio V a r i a b l e Model Homogenous flow Bypassing flow only Overall 1984 only Overall 1984 only Evaporation 1.12 1.12 0.65 0.77 Transpiration 1.07 1.00 0.92 1.00 Drainflow 0.92 0.96 1.19 1.06 Drain loss of: NH4 1.56 1.85 0.72 0.61 NO3 1.35 1.61 0.77 0.75 PO4 1.06 1.32 0.57 0.48 128 4.4 Discussion A visual inspection of Figs. 4-2 and 4-3 reveals that the model is generally very successful at predicting the depth of the water table for the whole simulation period. The sequence of troughs and peaks is well reproduced throughout, and there is good agreement between. the predicted and actual water table depth. However, there are some significant discrepancies which illustrate the limitations of the model. The weather data was obtained from the Vancouver International Airport station located approximately 16 kilometers from the experimental site. The geographic variability of precipitation, in particular, may be responsible for such discrepancies as the difference in the water table peaks following the storms of May 1st, 1984 and April 25th, 1985, where the predicted rise is much sharper than actually observed. However, while it is important for specific events, this source of errors is unlikely to alter the overall performance of the model very significantly. The model assumes all precipitation to fall as snow when the average air temperature is below freezing; further, there is no account of snowmelt by heat conduction from the warm soi l . This simplification may have caused an overprediction of water table depth in the early December of both years, which is more pronounced in the undrained plot. The model predicts deeper water tables between storms than observed, particularly during the winter of 84-85, for both drained and undrained plots. This may be due to an overestimation of 129 drainage and percolation, or of surface evaporation. This overprediction may also have resulted from an erroneous estimation of the maximum depth of ponding; Madramootoo et a l . (1987) discussed the effects of this variable on soil hydrologic response. During the growing season, the water table drops largely as a result of transpiration and upward capillary flow, as shown by Makkink and Van Heemst (1975). This, in turn, depends on the growth rate, depth, and distribution of the roots. Root growth is simulated in a simplified manner, which may account for the observed differences in the shape of the water table curve. The model makes no adjustment in root growth pattern for water availability or soil temperatures. Consequently, the model cannot simulate a deeper root growth in the drier summer of 1985; although the model predicts a deeper water table in 1985, this is solely due to capillary flux and underestimates the depth reached. Hillel and Talpaz (1976) showed the importance of adjusting the simulation of root growth to water availabil i ty. However, in designing the present model it was felt that the expected gain in accuracy did not warrant the added complexity and calibration effort required to simulate root growth in a more realistic fashion. The model underestimates the rise of the water table in autumn, in both years, and for both the drained and undrained cases. This can arise from an overestimation of the total evapotranspiration accumulated over the growing season. However, this may also result from overestimating the rate at which the micropores become saturated. The model assumes that equilibrium is 130 attained instantaneously ( i .e . , within a time increment of one day) between water in the micropores and static water f i l l i n g the macropores. While lateral infiltration from water flowing in the macropores is calculated using the approach of Green and Ampt, lateral infiltration from static water in the macropores is allowed to fully saturate the micropores within a day, because of the increased gradient. This observation may not be adequate in the context of this f ie ld. The peaks resulting from the storms of early September 1984 and 1985 exhibit a very sharp rate of rise and drawdown. While the sharp rise is possible if the soil at that depth is nearly saturated, the drawdown rate cannot be explained solely by upward capillary flux and deep percolation. Lateral infiltration from the cracks network into the peds probably contributes largely to such a rapid drawdown. This mechanism could also explain the very rapid rate of rise of the water table in the f a l l . Thus i t appears that the differences between actual and predicted values result mainly from a failure to obtain very accurate model parameters and weather data. Nevertheless, the f i t of the model to actual values remains very good, as indicated by the results of the Spearman test. Further, the response of the simulated water table indicates that both bypassing flow and bulk flow must be included as water movement mechanisms in order to obtain a successful f i t . This is made evident by comparing the simulated profiles shown in Figs. 4-2 and 4-8. The three models (namely, the current model, the model 131 including solely bypassing flow, and the homogeneous flow model) show l i t t le difference in the prediction of water table depths over winter. This is to be expected since there is l i t t l e or no difference in the values of the overall saturated conductivity and, effectively, of the drainable porosity used in all three models. However, there are major differences in the summer responses of the models. A model that uses bypassing flow as its sole mode of vertical water movement has no provision for capillary flow to the root zone. In this simulation, deep percolation is assumed to be n i l , and therefore the water table does not drop below the drain level. This, of course, conflicts with the observed profile. It is not possible to adjust the model parameters to obtain a satisfactory f i t over the entire simulation period. Thus, it is clear that bulk flow cannot be excluded from a successful simulation of water movement in this type of soil characterized by high water tables. However, it is equally clear from observing late summer and autumn profiles that bypassing flow cannot be dismissed as a water movement mechanism for the present situation. The September storms can be seen to cause sharp rises of the actual water table. This -.sudden response cannot be duplicated when bypassing flow is neglected. This is due to the fact that there is l i t t le moisture in the unsaturated zone during late summer, as a result of evapotranspiration. Because of high moisture tensions, any precipitation is held by the upper layers of the soil and does not percolate to the water table. Consequently, the water table only 132 exhibits the progressive drop due to capi l lary flow towards the root zone, unti l the prof i le has been nearly replenished by autumn rains. In order to reproduce the observed response of the water t a b l e , either the water holding capacity of the s o i l , or the evapotranspiration must be reduced, or the upward flux must be increased. However, the latter results in unreal ist ical ly deep water tables, while the former two need to be reduced to such an extent as to become completely unrea l is t ic . This shows that under the present experimental conditions, bypassing flow is an important mechanism, at least in summer and autumn. While the model produces an excellent f i t of the observed water table p r o f i l e , the same cannot be said for the observed nutrient concentration of drain water. Predicted values of ammonium and orthophosphate do not f i t the data (according to the results of the Spearman and Wilcoxon tes ts ) , while the f i t of nitrate concentrations, as seen in F ig . 4-4, is only passable. Part of the d i f f i cu l ty l ies in the quality of the f i e ld data i t s e l f . Because of the cost of analysis, samples could not be collected at very frequent intervals; furthermore, such problems as occasional pump fa i lure and ditch flooding prevented some sampling. Further, as they are grab samples, they indicate an instantaneous concentration, as compared to the daily average concentration given by the model. A few values were obtained while drainflow was a mere t r i c k l e , in late spring. The phosphorus data contains a few values best described as outl iers (see Chapter Three). F ina l ly , nu t r i en t concentrat ion from dra in to drain was d i f f e r e n t , 133 indicating the importance of spatial variability. Predicted ammonium and orthophosphate concentrations in soil water are markedly affected by the choice of the adsorption parameters. Adjusting these parameters provides a decent f i t for orthophosphate (except for the outliners mentioned above). However, no parameter adjustment was found to provide a satisfactory f i t of ammonium concentrations, which are much higher than the model prediction, unless the initial concentration of soil ammonium was raised to an unrealistic level. Partial inhibition of nitrification in this fairly acidic soil (pH of 4.4 in topsoil), loss of topsoil particles through cracks, or contribution from microfauna that may be present in the drains are possible explanations for this lack of f i t , but none of these mechanisms were included in the simulation. The same difficulty may have been experienced by Johnsson et a l . (1987), as they found it necessary to maintain the ratio of ammonium to nitrate at an arbitrary level in order to obtain an adequate simulation of ammonium concentrations. The simulation of nitrate concentrations provides an adequate f i t of the data in terms of magnitude and overall trend. The quality of the f i t is comparable to that reported from the models of Duffy et a l . (1975), Kanwar et a l . (1983), and Johnsson et a l . (1987). However, this comparison does not imply that the present model is successfully validated; i t merely highlights the di f f icul t ies involved in simulating soil nitrogen behaviour. Furthermore, Kanwar and Johnson (1984) found their model to be most 134 sensitive to the init ial soil nitrogen content, which they treated (as does the present model) as a fitted parameter. An inspection of Figs. 4-4 and 4-9 makes it possible to contrast the effects of the different leaching mechanisms. Fig. 4-9 represents the concentrations resulting from bulk flow, which effectively acts as piston flow (with numerical dispersion). The smoothness of this curve contrasts with those of Kanwar et a l . (1983) and Duffy et a l . (1975) because these workers used a variable weighting factor to account for macropore flow. It can be seen from Fig . 4-9 that piston flow is characterized by an increased concentration in autumn followed by a progressive drop. Although there is a phase shift, these trends are found in the actual data. However, piston flow fails to predict the magnitude and variability of the data collected in Winter 84-85 (December 15 to May 15). By contrast, the inclusion of bypassing flow, as seen in Fig. 4-4, results in a more accurate prediction of both the magnitude and the variability of the data during the autumn and winter of 84-85. However, in both springs of 84 and 85, the bypassing flow model f a i l s to predict a smooth decline in concentration but rather exhibits the opposite trend. It is possible that the importance of bypassing flow varies seasonally; this has been observed in clay soils by Robinson and Beven (1983). A plot of concentrations resulting from using solely bypassing flow (as seen in Fig. 4-9) results in a more jagged profile, and in a poorer f i t . Thus, again, it seems that neither mechanism alone is sufficient to explain the observed values; further, there is an 135 indication that the relative importance of these mechanisms may not be constant. The sensitivity analysis, as seen in Figs. 4-4 to 4-6 and Tables 4-2 to 4-4, i l lustrates the relative importance of the dif ferent mechanisms on water flow and nutrient leaching. Intrapedal conductivity CS, conductivity exponent ALFA (not shown), and ped spacing PEDS have very l i t t le effect on the water table profile in winter. A decrease in the macropore volume 9PR causes the water table peaks to be somewhat more jagged in winter. However, differences are evident in summer. Since CS and ALFA determine the unsaturated conductivity, they control the rate of upward capillary flux, and thus the depth of the water table in summer. Macropore volume has a similar effect because it controls the volume of water available for capillary flow, per unit depth. The equivalent ped spacing PEDS controls the extent of lateral infiltration from macropores into peds. Accordingly, when PEDS is small , lateral i n f i l t r a t i o n prevents bypassing water from percolating to the water table. The resulting effect on the water table is similar to that of excluding bypassing flow, although the water distribution in the unsaturated zone is different. (It should be pointed out that the parameter PEDS can be varied independently of the bypass hydraulic conductivity CK in this model. However, it is expected that these two parameters be inversely correlated in an actual s o i l , because an increasing frequency of bypassing channels should increase the hydraulic conductivity of the profile i f the dimensions of the channels 136 remain relatively constant.) In general, the effect of varying these parameters on evaporation, transpiration, and drainflow is minimal. As expected, during the drier 1985 the lower unsaturated conductivity (from a low CS or a high ALFA), caused a reduction in evapotranspiration because of a reduction in upward flow. This is even more pronounced when CS is set to zero (Table 4-4). This lower evapotranspiration implies that less water is needed to replenish soil moisture, which results in an increased drainflow. This is also seen as a result of decreased lateral infiltration, though this effect is much less pronounced. However, the variation of these parameters had a very pronounced effect on nutrient losses. Decreasing bulk flow markedly decreases leaching losses, as can be seen in Tables 4-2, 4-3, and 4-4. This effect is more pronounced when CS is reduced than when ALFA is increased, because the latter only affects flow in the unsaturated mode, which is a small fraction of the total flow. Losses of nitrate are more sensitive to a variation in bulk flow conductivity than phosphate, probably because phosphate is strongly adsorbed and leaches slowly. In contrast, ammonium appears to be more sensitive; this reflects the fact that the rate of ammonification increases with an increase in bulk flow, probably because of improved aeration. The sensitivity of leaching losses to a variation in bulk flow is more pronounced at the end of one year (Table 4-3) than after two years (Table 4-2); this arises because the soil nutrients are partially depleted at the end of the 137 f i rst year. Because of the dynamic nature of the simulation, i t is difficult to eliminate such transient effects from the comparison. PEDS controls not only lateral inf i l trat ion but also the degree of mixing between macropore and micropore water. This effect is strong over both years for all nutrients. Superficially, the effects of decreasing PEDS appears similar to that of increasing bulk flow, but there are important differences. Bulk flow causes a vertical, downward movement of the solutes (except in summer) whereas lateral mixing causes mostly a reduction in nutrients without much altering their vertical distr ibution. Lateral mixing and infi l tration is controlled by the macropore volume 6PR, though in a less important fashion than by PEDS. In effect, a larger 6PR results in a longer average residence time of macropore water, which consequently increases the time available for mixing. Consequently, this results in more leaching losses, as there is more diffusion of nutrients to the macropore zone. In contrast to the other cases, phosphate is more sensitive to a change in 9PR than the nitrogen species. This is probably because a larger macropore volume provides both a longer duration for diffusion and a larger dilution volume, thereby affecting more importantly a strongly adsorbed nutrient. Timing effects are also important. In both cases of lower conductivity and large ped spacing, drains start flowing earlier in autumn because less water is needed to recharge the soil below drain level . These early flows cause nutrient losses to be temporarily higher than when PEDS or ALFA are small or CS is large. 138 Although nutrient losses from the latter cases eventually surpass the former, these losses would be even larger if this timing factor was of less importance. Evidently, the validity of the above discussion depends on the soundness of the assumptions on which the model is based. Nevertheless, these results unquestionably show the importance of bypassing flow, bulk flow, and lateral mixing as d ist inct mechanisms. These mechanisms can combine to produce a relatively equivalent hydrological balance (in terms of evapotranspiration and drainflow), but a very different nutrient balance. These findings have significant implications for the management of nutrients in drained soi ls, as management practices, such as ti l lage, may affect the relative importance of these mechanisms. The results of the model, in general, can be seen to justify the choice of the water table as the output variable for calibration. The response of the water table is very sensitive to the combination of water flow mechanisms selected during the drier months. Such a response provides a convenient criterion by which the importance of bypassing flow, in par t icu la r , can be ascertained. In the absence of such a criterion, it would be difficult to distinguish the effects of the different water flow mechanisms based on nutrient flux alone. 4.5 Conclusions In general, it can be seen that the model is successful at simulating the behaviour of the flow data obtained from the 139 experimental f ie ld . The model also shows sufficient f lexibi l i ty to investigate the relative effects of different water and nutrient movement mechanisms, even if the ability of the model to predict nutrient flow remains largely untested. Specific conclusions that can be derived from the calibration and sensitivity analysis are as follows: • The model provides a very good f i t of the observed water table prof i le and an adequate f i t of the observed nitrate concentration in drain water. The observed ammonium and phosphate concentration of drain water are not successfully simulated. • Both bypassing flow and bulk flow are important moisture transfer mechanisms of the experimental site. It is not possible to successfully simulate the behaviour of the experimental variables while excluding either one of these mechanisms. • Combinations of moisture transfer mechanisms that produce equivalent values of yearly hydrological balance, in terms of drainage and evapotranspiration, can result in very different values of yearly leaching losses of a l l nutrients under consideration. • In soi ls where lateral mixing between macropores and micropores is small, leaching is increased by an increase of the ratio of bulk flow to bypassing flow, everything else being equal. Leaching is also increased, to a lesser degree, when the ratio of macropore to micropore volume increases. In soils where bypassing flow is dominant, leaching increases with an increase of lateral mixing between micropores and macropores, everything else being the 140 same. In summary, the relative importance of different moisture and nutrient transfer mechanisms can be seen to have a significant influence on leaching losses through drainage. Consequently, the effects of different drainage practices on leaching losses can be expected to vary accordingly. This also implies that drainage losses may be altered by tillage and fertilization practices, if these affect the soil ratio of bulk flow to bypassing flow, or the distribution of nutrients between the micropore and the macropore domains. 141 4 . 6 References - Bengtson, R .L . , C E . Carter, H.F. Morris, and J.G. Kowalczuk, 1984. Reducing water pollution with subsurface drainage. Trans. ASAE, 27, 80-83. - Burke, W., J . Mulqueen and P. Butler, 1974. Fertilizer losses in drainage water from a surface water gley. Ir. J . Agric. Res., 13, 203-214. - Chao, E . C . Y . , 1987. Water table depth simulation for flat agricultural land under subsurface drainage and subirrigation practices. Unpub. M.Sc. Thesis, Univ. British Columbia. - Caussade, B., D. Pierre and M. Prat, 1984. Modele d'etude et de gestion de la teneur en azote des eaux de surface dans un bassin versant sous culture. J . Hydrol., 73, 105-128. - Culley, J . L . B . , E.F. Bolton and V. Bernyk, 1983. Suspended solids and phosphorus loads from a clay soil . J . Environ. Qual. 12, 493-503. - Deal, S . C , J.W. Gilliam, R.W. Skaggs, and K.D. Konyha, 1986. Prediction of nitrogen and phosphorus losses as related to agricultural drainage system design. Agric. Ecosyst. Environ. 18, 37-51. - Devitt, D., L. Letey, L . J . Lund and J.W. Blair, 1976. Nitrate-n i trogen through s o i l as af fected by s o i l p r o f i l e characteristics. J . Env. Qual., 5, 283-288. -Duf fy , J . , C. Chung, C. Boast and M. Franklin, 1975. A simulation model of biophysiochemical transformations of nitrogen in tile-drained Corn Belt soil . J . Environ. Qual., 4, 477-486. - Driehuyzen, M.G., 1983. Boundary Bay water control project. Unpublished report, B.C. Ministry of Agric . and F i s h . , Cloverdale, B.C. - Gibb, A . J . , 1986. The effect of water table lowering rate on the drainable porosity of agricultural soils. Unpub. B.Sc. Thesis, Univ. British Columbia, Vancouver. - Harris, G .L . , M.J. Goss, R.J. Dowdell, K.R. Howse and P. Morgan, 1984. A study of mole drainage with simplified cultivation for autumn-sown crops on a clay soil . J . Agric. Sci. Camb., 102, 561-581. 142 H111el, D., and H. Talpaz, 1976. Simulation of root growth and its effect on the pattern of soil water uptake by a non uniform root system. Soil Sc i . , 121, 307-312. Johnsson, H., L. Bergstrom, P.E. Jansson and K. Faustian, 1987. Simulated nitrogen dynamics and losses in a layered agricultural so i l . Agric. Ecosys. Environ., 18, 333-356. Kanwar, R.S., and H.P. Johnson, 1984. Sensitivity analysis of a drainage simulation model. Water Resources Bul l . , 20, 339-342. Kanwar, R.S., H.P. Johnson, and J .L . Baker, 1983. Comparison of simulated and measured nitrate losses in t i le effluent. Trans. ASAE, 26, 1451-1457. Madramootoo, C.A., R.S. Broughton, S.O. Prasher and F. Peyrow, 1987. Development and testing of an agricultural land drainage simulation model. Computers and Electronics in Agric. 2, 15-30. Makkink, G.F. and M.D.J, van Heemst, 1975. Simulation of the water balance of arable land and pastures. Centre for Agr. Pub. Docum., Wageningen. Overcash, M.R., K.R. Reddy and R. Khalel, 1980. Chemical and biological processes influencing the modelling of non point source water quality from land areas receiving wastes. Proc. Hydrol. Transport Modeling. ASAE pub. 4-80, St. Joseph, 1-11. Robinson, M., and K.J. Beven, 1983. The effect of mole drainage on the hydrological response of a clay soi l . J . Hydrol., 64, 205-223. Siegel, S. , 1956. Non-parametric statistics for the behavioural sciences. McGraw Hi l l , New York. Stoskopf, N.C., 1981. Understanding crop production. Reston Pub. Co., Reston. Switzer-Howse, K.D., 1983. Agricultural management practices for improved water quality in the Canadian Great Lakes basin. Agr. Can., Res. Branch publ. 1983-24-E. 143 CHAPTER 5 A simulated comparison of the effects of different drainage  practices on nutrient losses 5.1 Introduction In humid areas, drainage of agricultural land is often required in order to improve crop growing conditions and trafficability. Surface drainage, a combination of land forming and shallow ditches, promotes removal of excess water by surface runoff. Subsurface drainage relies on underdrains to lower the water table and provide an outlet for percolation of excess water. Subsurface drainage also tends to reduce surface runoff, because of the increased storage volume available in a soil with a lower water table, provided that the storm intensity does not exceed the infi ltrabil ity of the soi l . Surface drainage is associated with erosion and nutrient losses via surface flow, while subsurface drainage tends to promote leaching losses of nutrients. Accordingly, it is difficult to determine the overall impact of subsurface drainage on nutrient losses, as stated by Switzer-Howse (1983). Furthermore, experimental evidence shows that this impact varies largely with site and experimental conditions, and can range from an overall increase to an overall decrease in loss of nutrients (as seen in Chapter One). Consequently, it is dif f icult to draw general conclusions about the impact of subsurface drainage on the environment based solely on field evidence. Deal et a l . (1986) showed that the computer simulation of water and nutrient flow in 144 drained soils provides a method to derive generalizations regarding the impact of drainage, because it makes paired comparisons possible over a wide range of situations. Percolation of water and nutrients in field soils is complex and leaching patterns are often best represented by a statistical distribution which reflects the effects of spatial variability and dispersivity. For structured soils White (1985) showed that this distribution is characterized by the effects of two distinct flow domains: that of the flow within the micropores of the peds and that of macropore flow. The relative importance of these flows in paral lel and of the linkage between them has a determining influence on the amount of nutrients that is leached to drainflow (see Chapter Four). Subsurface drainage is usually designed on the basis of a desired drainage coefficient, which represents the amount of water that can be removed daily by a drainage system under the design conditions (e.g., Van Der Gulik et a l . , 1986). Different drain depths and spacings can be combined to give a similar drainage coefficient. Kanwar and Johnson (1984) and Skaggs and Gil 11am (1981) showed that nitrate losses are increased by a narrower spacing and deeper depth of drains. Skaggs et a l . (1982) also showed that surface runoff and erosion can be decreased by either increasing drain depth or decreasing drain spacing. Unfortunately, the effects of drain depth and spacing on the overall nutrient budget were not investigated. Also, 1t is not possible to compare, on the basis of their results, depth and spacing combinations that 145 result in the same drainage coefficient. In both works, macropore flow, bypassing the peds, was not explicitly simulated; consequently, whether their conclusions apply to strongly structured soils remains unknown. Also, the effects of drainage design on mechanisms that affect nutrient losses only indirectly, such as denitrification, were not assessed by these workers. This chapter presents the results of a simulation of water and nutrient flow mechanisms designed to address these questions. 5.2 Method A computer model was developed to provide a comprehensive simulation of the main mechanisms that affect moisture and nutrient transfer in flat agricultural soils with humid climate. Crop growth, nutrient flows and bioreactions, and water flow are the main components of the model. The crop submodel determines the ratio of evaporation to transpiration, root growth, and nutrient uptake. Such bioreactions as mineralization, immobilization, nitrif ication and denitrification are simulated in the nutrient submodel. Nutrient movement is calculated as occurring by mass flow of solutes (after adsorption) from bulk flow, bypassing flow, lateral flow, or surface runoff. In the water submodel, water is redistributed within the soil peds according to the da i ly inputs of p r e c i p i t a t i o n , evapotranspiration, and drainage from percolation. The resulting flux can be upward or downward. The maximum micropore, or ped, 146 flux is restricted by the ped conductivity, denoted by CS at saturation. If the precipitation rate exceeds this maximum flux the excess water is routed down the macropores, again provided that it does not exceed the conductivity of the macropores, CK. When the water table has reached the soil surface, or when CK is exceeded, the excess moisture is allowed to pond on the soil surface. Surface runoff occurs when the ponding capacity is exceeded. The solute concentration of surface runoff and ponded water is assured to reach equilibrium, by diffusion, with the water held in the pores of the upper centimeter of the soi l . Nutrient losses by erosion are simulated using a modified version of the universal soil loss equation (Wischmeir and Smith, 1978) and the appropriate enrichment ratio. Water percolating down the macropores is also assumed to mix with micropore water. The extent of this mixing is determined by the ratio of the mixing depth, taken as one centimetre, to the equivalent spacing between the peds, PEDS. If the micropores are unsaturated lateral inf i l t rat ion is also simulated. More details on the computer model can be found elsewhere (see Chapter Two). The model was calibrated against experimental results obtained from a drained silty clay loam site located at Boundary Bay in the Lower Fraser Valley of British Columbia for 1984 and 1985. It was found that results from the model gave a very good f i t of the water table profile and an adequate f i t of the nitrate concentration in drain water, while underestimating the ammonium concentration 147 actually observed (see Chapter Four). The calibration showed that by-passing flow is an important mechanism for this soi l . Four different drainage designs are compared in this study. Design A consists of the drainage installation as it currently exists in the experimental f ield, with the drains lying at a depth of 1.1 metres and a spacing of 14 metres. In design B, the drains are 20 metres apart, also at a depth of 1.1 metres, which results in a drainage coefficient for design B of nearly half that of design A. Design C uses the same drainage coefficient as design A, but both spacing and depth are reduced to 10 metres and 0.7 metre respectively. Design D relies on surface drainage and uses a maximum ponding depth of one centimetre, which is half that of the other three designs. Subdrains are not used in design D, but 1.1 metre deep ditches, spaced 50 metres apart, provide an outlet for a small amount of percolation. In all four designs deep percolation is set at zero, in order to isolate the effects of drainage from those of deep percolation. The flow parameters assigned to this soil are such that bypassing flow is a dominant mechanism that largely controls nutrient leaching. In order to determine whether this characteristic affects the response of the model to different drainage designs, the simulation of designs A, B, and C was repeated using two different soils. These two hypothetical soils are similar in all respects to the experimental soil except with respect to the relative importance of bypassing flow, bulk flow, and lateral mixing. In soil B the value of the ped conductivity CS 148 is increased from 0.001 m/d to 0.005 m/d, which increases the ratio of bulk flow to bypassing flow. In soil C, lateral mixing between macropores and micropores is increased by reducing the equivalent ped spacing PEDS from 5.0 metres to 0.2 metre. A sensitivity analysis showed that both these changes increase leaching; in the case of soil B, the increase in bulk flow is such as to make this soil nearly equivalent to one where bypassing flow is negligible (see Chapter Four). Since the simulation is carried out for only two years, transient effects cannot be dismissed and therefore the init ial nutrient distribution may affect the response of the model to the different drainage designs. Two initial distributions of nutrients were used in these simulations. In the f i rs t distribution (referred to here as distribution C) mineral nutrients (nitrate, ammonium, phosphate) are uniformly distributed through the profile. The concentrations used are those which gave an adequate f i t of the observed f ie ld results, as shown in Chapter Four. Organic nutrients are distributed with a concentration that tapers with depth, reproducing the observed organic matter distribution. The second distribution, called T, is similar to distribution C except that it is truncated at a depth of 0.7 metre, below which the concentration of all nutrients is set to zero. The simulation was carried out for two years using the 1984 and 1985 weather data from Vancouver International Airport weather station. There were 1280 mm of precipitation in 1984 as compared to 800 mm in 1985. Table 5-1 gives a summary of the simulated 149 Table 5-1 Summary of the Conditions Simulated Drainage Design Drain Spacing Drain Depth A 14m 1.1m B 20m 1.1m C 10m 0.7m D* 50m 1.1m Soil Type Ped Conductivity (CS) Ped Spacing (PEDS) A 1 mm/d 5.0m B 5 mm/d 5.0m C 1 mm/d 0.2m •design D uses ditches as opposed to pipe drains 150 conditions. 5.3 Results The effects of the four drainage designs on the water table profile in soil A can be seen in figures 5-1 and 5-2. The frequency of occurrence and duration of high water tables is reduced in designs C and A when compared with design B, and likewise when design B is compared with design D. High water tables drop faster in design C than in design A, but are slightly more frequent (e.g., on June 28, 84 and March 24, 85). The differences in the water table profile between soils A, B, and C are negligible during the wet season for all designs. Table 5-2 shows the water balance of the different soils and drainage combinations. Evapotranspiration is affected l i t t le by the different designs, although higher values are seen in design C. The major difference between the designs is in the ratio of surface runoff to percolation, the average values of which are 0.010, 0.078, 0.012, and 0.325 for designs A, B, C, and D, respectively. Percolation here includes both drainflow and seepage to ditches. No runoff occurred in all three underdrained designs during the relatively dry year of 1985. There is l i t t le difference in the hydrological response of the three soils for all designs. The nutrient losses arising from surface flow are summarized in Table 5-3. Again, there is l i t t le difference between the three soils and, obviously, none between the two distributions. Nitrogen losses include organic species, nitrate, and ammonium in both WATER TABLE PROFILE GN fl (NORMAL SPACING) 1?H . J WIDE" sPfiC 1 KG I. J F M f l M J J R S O N D J F M f l M J J f l S O N D 1984 TIME OF THE YERR 1985 Figure 5-1: Simulated Water Table Profiles Resulting From Design A (Subdrains at 1.1M Depth and 14M Spacing) and Design B (Subdrains at 1.1M Depth and 20M Spacing) WATER TABLE PROFILE DESIGN C (NARROW AND SHRLLOW DRAINS) Figure 5-2: Simulated Water Table Profile Resulting From Design C (Subdrains at 0.7M Depth and 10M Spacing) and Design D (No Subdrains, 1.1M Deep Ditches at 50M Spacing) 153 Table 5-2 Drainage and Soil Type Effect on Water Fluxes Drainage Soil Evap Tran Runf Drainf Seepf metres A A .399 .608 .011 1.051 .019 B A .403 .605 .078 .970 .028 C A .403 .649 .012 1.005 .008 A B .414 - .629 .017 1.010 .019 B B .420 .630 .082 .926 .028 C B .426 .650 .011 .981 .008 A C .402 .613 .019 1.035 .018 B C .407 .611 .088 .951 .027 C C .404 .650 .012 1.003 .008 D A .432 .613 .254 .000 .781 154 Table 5-3 Average Losses from Surface Flow for the Drainage Designs Type Design A B C D Runoff (m) .016 .083 .012 .254 Soil loss (tonne/ha) .045 .237 .038 .914 Carbon loss (kg/ha) 11.97 62.72 10.02 241.01 Nitrogen loss (kg/ha) .334 1.717 .278 6.535 Phosphorus loss (kg/ha) .189 .947 .157 3.431 155 adsorbed and solute forms. Phosphorus losses include organic and mineral forms, both adsorbed and solute. Losses of soluble forms are roughly proportional to runoff volume, while adsorbed or particulate losses increase faster with large runoff volumes. Tables 5-4 and 5-5 show the nutrient losses in drain water from designs A, B, and C for the three soils and the two distributions, expressed as a proportion of the average loss. In general, soils B and C leach more nutrients than soil A, as expected. However, the opposite is also seen, for instance in the case of designs A and B, where more ammonium and phosphorus is leached by soil A than soils B and C when the distribution T is used. The importance of the in i t ia l nutrient distribution can be seen when comparing the effects of the designs. When using distribution C, design B is found to leach the largest amount of nitrate and ammonium, while the leaching losses of phosphorus are relatively similar. In contrast, when using distribution T, design B leaches the lowest amount of both nitrate, ammonium and phosphate. Also, a larger amount of ammonium and phosphate is leached from design C than from design A when distribution T is used. The effects of drainage design on leaching losses are relatively independent of the effects of soil type, when the uniform distribution C is used. By contrast, there is a strong interaction between the effects of drainage design and the effects of soil type, for ammonium and phosphorus, when the nutrients are layered according to the truncated distribution T. Using this distribution, selecting design B over design A causes a far larger 156 Table 5-4 Relative Drain Losses Using Distribution C, as a Ratio of the Average Drainage Type A B C Soil Avg. N03 .91 1.23 1.18 1.11 A NH4 .77 1.12 1.21 1.03 P04 .92 1.09 1.06 1.02 N03 1.00 1.39 1.31 1.23 B NH4 .95 1.26 1.35 1.19 P04 .90 .99 .98 .96 N03 .51 .76 .70 .66 C NH4 .49 .91 .93 .78 P04 .94 1.06 1.07 1.02 N03 .81 1.13 1.06 1.00 Avg NH4 .74 1.10 1.16 1.00 P04 .92 1.05 1.04 1.00 Average losses are: 70.771 kg/ha N03-N 2.549 kg/ha NH4-N 4.105 kg/ha P04-P 157 Table 5-5 Relative Drain Losses Using Distribution T, as a Ratio of the Average Drainage Type Soil A B C Avg N03 .74 1.32 1.15 1.07 A NH4 .92 .74 .79 .82 PO4 1.36 .76 .47 .86 N03 .60 1.11 .95 .89 B NH4 .77 .26 .51 .51 P04 1.12 .29 .20 .54 N03 .73 1.30 1.11 1.05 C NH4 1.24 1.85 1.91 1.67 P04 2.02 1.38 1.41 1.60 N03 .69 1.24 1.07 1.00 Avg NH4 .98 .95 1.07 1.00 P04 1.50 .81 .69 1.00 Average losses are: 18.550 kg/ha N03-N .337 kg/ha NH4-N 1.318 kg/ha P04-P 158 reduction of losses of ammonium and phosphate in soils B and C than in soil A. Similarly, opting for design C instead of design B results in a more important increase of ammonium and phosphate loss in soils B and C than in soil A. Table 5-6 shows the effects of soil type and drainage design on mineralization, n i t r i f i ca t ion , and deni t r i f ica t ion , for distribution C. There is l i t t l e difference between the two distributions on these terms. Mineralization and nitrification are greatest in design C, followed by designs A, B, and D. Denitrification, conversely, is largest in design C, followed by designs D, B, and A. The effect of soil type is also noticeable, there being more nitr i f icat ion and mineralization, and less denitrification in soil B than soil C, and more in soil C than soil A. Denitrification is largest in soil C, followed by soils A and B. Seepage losses to ditches for design D are small, equal to 0.030, 1.946, and 0.096 kg/ha for ammonium, ni trate, and orthophosphate, respectively. 5.4 Discussion There are few unexpected results found in the simulated hydraulic response of the various combinations of soils and drainage designs. The frequency of occurrence and duration of high water tables decreases when the drainage design is based on a larger drainage coefficient. The total percolation volume also reflects the drainage coefficient. Evapotranspiration is affected 159 Table 5-6 Effect of Drainage Design and Transfer Mechanisms on Biochemical Transformations Soil Drainage Transformation (kg/ha) N Miner. P Miner. Nitrif. Denitr. A A 34.81 15.80 123.31 29.49 A B 34.12 15.51 121.92 31.62 A C 38.49 17.48 126.67 35.90 B A 41.42 18.83 128.54 27.08 B B 40.65 18.49 127.02 28.67 B C 42.62 19.36 129.68 30.72 C A 38.40 17.40 125.86 30.38 C B 37.38 16.96 124.09 32.24 C C 41.17 18.64 128.42 35.51 A D 33.50 15.23 116.22 32.60 160 l i t t l e by the drainage design because drainage is of l i t t le consequence in this soil during the growing season, when there is a high evapotranspiration rate. Since the overall water flux from the profile is nearly similar, the balance consists of surface runoff, which increases when a lower drainage coefficient is selected. These findings agree well with those found in the literature (e.g., the review of Irwin and Whiteley, 1983). Surface runoff is clearly reduced after using a high drainage coefficient in design. Runoff occurred during 44 days in design D, as compared to 9 days for design B and only on one day for designs A and C. Eleven of the runoff events of design D had a flow rate larger than one centimetre per day. All the runoff events were caused by the precipitation rate exceeding the percolation rate, bringing the water table to the surface. There was no instance of the precipitation exceeding the infiltrability of the so i l . This situation arises from the type of weather used in this simulation, characterized by storms of long duration and gentle intensity. However, it must be pointed out that the simulations were run with the assumption that daily precipitation is distributed evenly over the twenty-four hour period. This assumption is imposed by the choice of a daily time increment in the model, and may appear unrealistic in that it does not reflect the hourly variability of precipitation, where short intense downpours may easily exceed temporarily the inf i l t rabi l i ty of the soil and cause ponding. However, the depression storage in the field prevents runoff from occurring before a minimum depth of ponding is reached and this 161 prevents sudden surges of runoff. This mechanism, the gentle nature of the climate, and the fact that runoff was not actually experienced in the calibration field (design A) combine to make the above assumption workable. Nevertheless, it must be remembered that this assumption plays an important part in generating the differences in surface runoff found between the different designs. Although the hydraulic behaviour of designs A and C are generally similar, there are a few noticeable differences. The drainage coefficient is calculated on the basis of a water table depth of 25 centimetres. Design C, because of its shallow drains, exhibits a faster response than design A when the water table is above that depth, as can be seen by inspection of the steady state drainage equation (e.g. , from Van Der Gulik et a l . , 1986). Consequently ponded water resulting from the water table reaching the surface is drained faster in this design, which may explain the smaller volume of surface runoff predicted in soils B and C. However, since drainage ceases at a shallower depth, there is less storage volume available for rain water, with the consequent higher water table rise seen, for instance, on June 29, 84. Because of this factor, it is difficult to generalize on the effects that this design would have on surface runoff, if the simulation duration was extended. Another point which arises from the comparison is the fact that the predicted transpiration is higher in design C (by 6%, on average). This may be a consequence of the fact that drainage ceases at a shallower depth, thereby conserving moisture for the 162 crops. This fact, and the shorter duration of high water tables, may indicate that design C is better than design A from a point of view of crop yields and soil husbandry; however, the model is poorly equipped to investigate this aspect. The difference in the water table profiles of designs C and A during the summer is an artifact of the model (caused by the way upward flow is calculated) but this does not affect the present conclusions. As discussed previously (see Chapter Four), a decrease in the ratio of bypassing flow to bulk flow, which characterizes soils B and C, causes a small decrease in drain flow and an increase in evapotranspiration. This arises because more moisture is found in the unsaturated zone of these two soils than in soil A; this is independent of drainage design. For the same reason, soils B and C produce a larger volume of runoff. The accuracy of the predicted surface losses is unverified because these losses are directly related to the predicted amount of surface runoff, which is subject to the provisos mentioned above. Furthermore, as there was no runoff from the experimental f ie ld , the erosion parameters could only be estimated (with the exception of soi l erosivi ty , measured by Peterson [1986]). Nevertheless, the erosion and nutrient losses, predicted using well-accepted theory (adapted from Williams et a l . , 1984), give results comparable to those found elsewhere (for instance, see Table 1-2 in Chapter One). Erosion losses are comparatively larger in design D, because of the large number of runoff events with a high flow rate. Losses of adsorbed and particulate nutrients are 163 also high accordingly. Losses of organic carbon, nitrogen, and phosphorus follow approximately the proportions of composition of the topsoil organic matter. Ammonium and phosphate losses remain high even when erosion losses are low, because these nutrients desorb from the surface layer when in contact with the relatively large ponding and runoff volume. This does not occur with nitrate, because this nutrient is leached from the surface before it can mix with ponded water. Solute losses remain high even from a low runoff volume because of the initial flushing effect experienced at the onset of a runoff event. The responses from all three soils are similar and reflect the volume of runoff. Like surface losses, drainage losses of nutrients are comparable to those reported in the literature. However, because of the high initial nutrient levels in the soil and the low levels of fert i l izat ion, transient effects prevail in this simulation. Also, the much wetter climate of 1984 combines with these transient effects to produce much higher nutrient losses during that year. It is this dynamic and transient nature that produces some of the unexpected results that arise from the comparison of the different drainage designs. Skaggs and Gilliam (1981) and Kanwar and Johnson (1984) found that nitrate losses increased as a result of selecting narrower drain spacings. The latter authors also found such losses to increase when using deeper drains. Their sensitivity analysis showed these increases to be larger than the corresponding increase in drainflow volume. These results appear to contradict the 164 findings of the present study, as reported in Tables 5-4 and 5-5, and make it necessary to examine the leaching mechanisms in detail in order to explain such discrepancies. Changing the spacing of the drains alters the geometry of the flow patterns underneath the drains. Although the overall flow rate diminishes as a result of the wider spacing, the deeper layers of the profile contribute proportionately more water than at a narrower spacing because of the adjustment of the flow lines (unless these are constrained by a shallow impermeable layer). Consequently, if the whole profile is initially assigned a nitrate concentration constant through the profile, there is a larger mass of nitrate that can potentially leach to drainwater if the spacing of these drains is increased. This becomes particularly manifest when the soil nitrate is not fully replenished by fertilization and transient effects are dominant, as is the case in Table 5-4, calculated using the uniform distribution C. However, when the in i t ia l mass of nitrate that can leach to drains is constant regardless of the drainage geometry, as in the case of distribution T, this effect disappears and nitrate losses are reduced when a wider spacing is selected. This reduction arises from the concomitant reduction in drain flow, but also from changes in nitrification and denitrification. This transient effect is not found in Skaggs and Gi l l iam (1981) because the nitrate concentration of groundwater was set at a constant value throughout the simulation. Similarly, this effect was not seen in the sensitivity analysis of Kanwar and Johnson (1984) because the 165 nitrate mass available for drain leaching was not affected by the drainage geometry. Likewise the reduction in leaching losses consequent to a reduction in drain depth results from a reduction of drain flow in the findings of these workers. In the present study, drain flow volume is kept relatively constant by reducing both depth and spacing. When the uniform distribution C is used, the leaching losses from design C are smaller than A, again because the contributing volume is smaller, although it is depleted at a faster rate. This effect disappears when using distribution T; in this case, the nitrate losses of designs A and C are similar. Ammonium and phosphate losses follow the same pattern as nitrate, except for the effect of adsorption. When the uniform distribution C is used, the effects of soil type or drainage design are similar for these nutrients, except that adsorption tends to buffer the transient effects, especially for phosphate. By contrast, adsorption magnifies the difference between the designs using the layered distribution T. In this case, there is a larger volume of soil that adsorbs and retards the flux of these nutrients in designs A and B than in design C. The type of leaching mechanism has a marked influence on the response of the soil to the type of drainage design selected. When bypassing flow is the dominant mechanism, as in soil A, the differences between designs A, B, and C are much less pronounced than when the micropore water is more mobile, whether vertically in soil B or laterally in soil C. This interaction is particularly 166 evident in the case of distribution T for adsorbed nutrients. When bypassing flow is dominant, there is relatively l i t t le interaction between the nutrient-rich macropore water and the nutrient poor micropore water from the lower soil layers. Thus in this case the difference in leaching response is mostly related to the difference in drain flowrates from the various designs. On the other hand, the much larger micropore volume that can effectively dilute and adsorb nutrients in designs A and B has an important effect when micropore water is mobile, as in soils B and C. This effect disappears in distribution C because there is no adsorption taking place, as the adsorption sites are already saturated. In a soil where either micropore flow or bypassing flow is inexistant, this effect would be magnified. Thus water flow mechanisms appear to be as important as drain flow rates in determining leaching losses. This situation was also observed by Johnsson et a l . (1987), when investigating the effects of soil layering. The improved aeration status of the designs using a larger drainage coefficient is , in all likelihood, responsible for the more aggressive mineralization and nitrification witnessed in these designs. However, both mineralization, n i t r i f i ca t ion , and denitrif ication are most pronounced in design C; this is found in al l three so i ls , and would appear to be consistent with the observation that more moisture is conserved in this design that in design C, in the unsaturated zone. Thus this design creates more favourable conditions for biochemical reactions, but also increased the occurrence of wet conditions, which enhances denitrification. 167 As expected, denitrification is much larger in design D than design A, despite the fact that there are fewer nitrates present in the profile of design D because of the reduced nitrification. Design B is intermediate between A and D. Because of its low ped conductivity and lateral infiltration, soil A is drier in summer than the other two. Consequently there is , comparatively, a reduced biological activity in this so i l , and mineralization and n i t r i f ica t ion proceed at a slower rate. Moisture content during the dry months is also a plausible explanation for the difference between soils B and C in that respect. The reasons for the differences in the denitrification rate of the three so i ls are more d i f f i c u l t to ident i fy . Denitrification is a function of the denitrification rate (which depends on the moisture content and the carbon mineralization rate) and of the concentration of nitrate, which depends on leaching losses and nitrification rates. These, in turn, also depend on the drainage design, which makes it difficult to isolate the effects of the soil type on denitrification. It is unlikely that the biochemical reaction rates were much influenced by the differences 1n soil temperatures associated with the different designs. At planting time, these differences, as calculated by the model, were less than half a degree centigrade at 15 centimetres of depth. Interestingly, while design B produces slightly cooler than A or C, design D produces the warmest temperatures. This is probably due to the decreased albedo of the wetter so i l , which counteracts the effects of the reduced thermal 168 diffusivity. Because of the lack of verifiability of the surface losses predictions, it is difficult to compare the overall losses from the four designs. Table 5-7 shows the overall losses for soil A and distribution T for the various designs. It can be seen that the total losses of phosphorus and Kjeldahl nitrogen are reduced as a result of intensive drainage, while the losses of mineral nutrients are much higher for those designs. The lack of verification of the surface nutrient loss predictions makes it difficult to conclude in favor of one design or another in terms of total loss of nutrients. However, despite the possible lack of accuracy, it appears clear that subsurface \ drainage causes an increase in the loss of mineral nutrients, under the conditions simulated here. This response is comparable to the field results obtained by Harris et a l . (1984). This similarity is encouraging because the soil type, topography, and climate of their experimental site is comparable to those simulated in this study. By contrast, Schwab et a l . (1980), Bengtson et a l . (1984), and Gilliam and Skaggs (1987) used field sites exposed to violent summer storms associated with severe surface losses. In the work of Bengtson et a l . (1984), in particular, a single, very large storm was responsible for the bulk of the total losses, which were caused mostly by surface flow. Since these were smaller in the drained field of their Louisiana fieldsite, subsurface drainage was seen to have a positive effect for all nutrients. On the other hand, the higher losses from an 169 Table 5-7 Average Overall Nutrient Losses from Drainage Selection Source A Design B C D Runoff (kg/ha) NO3-N .000 .000 .000 .000 NH4-N .014 .034 .012 .058 P04-P .024 .121 .028 .241 TKN .250 1.569 .304 6.535 TP .142 .897 .172 3.431 Percolation (kg/ha) NO3-N 13.730 11.074 13.590 1.946 NH4-N .309 .258 .419 .030 P04-P 1.791 1.474 2.656 .096 TKN .372 .320 .513 .030 TP 1.802 1.480 2.673 .096 Total (kg/ha) NO3-N 13.730 11.074 13.590 1.946 NH4-N .323 .292 .431 .088 P04-P 1.815 1.595 2.684 .337 TKN .622 1.889 .817 6.565 TP 1.944 2.377 2.845 3.527 170 undrained f ie ld than from a drained one in an English study (Armstrong, 1984) cannot be explained on the basis of the climate. In their study surface runoff had an unusually high concentration of nitrates, which was responsible for the results; this situation cannot be duplicated by this model. Culley et a l . (1983) found larger losses of phosphorus from shallow than from deep drains. This finding is related to the adsorption characteristics of their soil and is in agreement with the results from this model. In contrast, Schwab et a l . (1980) found increased losses of nitrate and phosphate from deeper drains. Their findings represent long term averages, which eliminates the possibility of transient effects; instead, this loss increase is simply the result of an increase of drain flow volumes, which is also in agreement with the present findings. There are several implications for drainage design that result from the model. It is clear, under the conditions simulated, that the design of a drainage system according to a higher drainage coefficient results in more nutrient losses to surface waters. This is aggravated if shallow drains with a narrow spacing are selected. In this case, leaching will be worsened if the nutrients are concentrated near the soil surface as opposed to being uniformly distributed, which is more likely to be the case in field soils. Conversely, these effects will be less pronounced, and the drainage design less cr i t ical , if the main flow mechanism in the soil is bypassing flow in soil cracks. For instance, this implies that mole drains, despite their shallowness and narrow spacing, may 171 not lead to very large nutrient losses because mole drainage Induces crack formation and bypassing flow (Leeds-Harrison et a l . , 1982). This reasoning appears to be supported by the findings of Armstrong (1984). These conclusions apply to the transient situation simulated by the model. Under steady state conditions, when the soil nutrient content remains relatively constant from year to year, there is no particular reason, at least in the light of the present findings and others, to suppose that a shallow and narrow drain placement would cause higher losses of nutrients if these are distributed (and replenished) uniformly through the soil profile. In fact, higher denitrification levels, as observed in design C, may produce the opposite. However, even if long term equilibrium is achieved, much of the leaching would remain of a transient nature; thus a large storm following fertilization would cause more leaching in a design using shallow and narrow drains than the opposite, i f no surface runoff was produced. Since the timing of the leaching losses is often critical in terms of environmental impact, the faster response of this type of design may be considered undesirable. Further work of simulation over a longer duration would be required to confirm this assertion. It may be observed that the transient nature of the simulation biases the comparison of designs A and D in favour of the latter. In this simulation, it is assumed that the subsoil has a nutrient content that is initially negligible, which would be characteristic of a soil of low fert i l i ty status. Consequently, seepage losses to 172 ditches from design D remain low while nutrients are accumulating in the subsoil. At steady state, the differences between these designs may be far less pronounced than under the present transient conditions. Nevertheless, the large leaching losses often encountered when subsurface drains are f i rs t installed (for instance, as reported by Bergstrom, 1987) may possibly be of more environmental significance than steady state losses. It is appropriate at this point to examine the plausibility of the results discussed above in the light of the ability of the model to reproduce field results. While the model can provide an very good prediction of the water table depth and an adequate f i t of the nitrate data, it has failed to reproduce the observed concentrations of ammonium and phosphate in drain water. This failure may be the result of a poor calibration, a poor quality of the field data, or of the omission or poor understanding of some mechanism affecting ammonium or phosphate behaviour. Consequently, caution must be exercised when extending the model comparisons to field situations. In a strict sense, the model predictions remain applicable only to field sites which are characterized by nutrient movement mechanisms rigorously similar to those included in the model. This close correspondence between such a model and the physical reality is difficult to establish; accordingly, this model has a poor predictive ability, a trait not untypical of this type of model (de Willigen and Neeteson, 1985). Consequently, the conclusions must remain limited in their applicability. However, as the model consists of a synthesis of recognized mechanisms and 173 has the ability of explaining field results from the literature, it is felt that comparisons between simulations of this model s t i l l provide a meaningful insight to the workings of the mechanisms examined and their implications for the design of drainage systems. 5-5 Conclusions The results from the simulation model demonstrate that the choice and design of drainage system can have a determining influence on the magnitude of field leaching losses under transient conditions. Specific findings, which apply only to the conditions simulated and assumptions used here, can be summarized as follows: • Leaching losses arising from drainage are larger when a higher drainage coefficient is selected. Overall losses, including surface losses, may also be higher, although the evidence remains tentative in this case. • Different designs based on the same drainage coefficient produce different leaching patterns. The magnitude and direction of this difference is strongly dependent on nutrient distribution and soil flow mechanisms. • Under transient conditions, the initial distribution of the s o i l nutrients af fects markedly the leaching patterns. Specifically, when nutrients are uniformly distributed through the soil profile, drains with narrow spacing produce less drain losses of nutrients, while the reverse is found when the nutrients are init ial ly concentrated near the soil surface. • Soils where bypassing flow is important, and lateral mixing 174 of micropore water with bypassing flow is negligible, are less sensitive to changes in drainage design in terms of leaching losses of adsorbed nutrients than soils which have more mobile micropore water. • Nutr ient t ransformat ion rates ( m i n e r a l i z a t i o n , nitrification, and denitrification) are sensitive both to drainage design and to mechanisms of flow in the soi l . 175 0 5.6 References - Addiscott, T.M., 1981. Leaching of nitrates in structured soils in M.J. Frissel and J.A. van Veen, eds. Simulation of nitrogen behaviour of soil-plant systems. Centre Agric. Pub. Docum., Wageningen. -Armstrong, A .C . , 1984. The hydrology and water quality of a drained day catchment, in T.P. Bart and D.G. Walling, eds., Catchment experiments in fluvial hydrology. Geo books, Norwich, 153-168. - Bengtson, R.L., C E . Carter, H.F. Morris, and J.G. Kowalczuk, 1984. Reducing water pollution with subsurface drainage. Trans. ASAE, 27, 80-83. - Bergstrom, L., 1987. Nitrate leaching from annual and perennial crops in t i le drained plots and lysimeters. J . Environ. Qual., 16, 11-18. - Burke, W., J . Mulqueen and P. Butler, 1974. Fertilizer losses in drainage water from a surface water gley. Ir. J . Agric. Res. 13, 203-214. - Culley, J . L . B . , E.F. Bolton and V. Bernyk, 1983. Suspended solids and phosphorus loads from a clay soi l . J . Environ. Qual. 12, 493-503. - Deal, S.C. , J.W. Gilliam, R.W. Skaggs, and K.D. Konyha, 1986. Prediction of nitrogen and phosphorus losses as related to agricultural drainage system design. Agric. Ecosyst. Environ. 18, 37-51. - Devitt, D., J . Letey, L.J. Lund and J.W. Blair, 1976. Nitrate-nitrogen movement through soil as affected by soil profile characteristics. J . Environ. Qual. 5, 282-288. - de Willigen, P., and J . J . Neeteson, 1985. Comparison of six simulation models for the nitrogen cycle in the soi l . Fertilizer Res., 8, 157-171. - Duffy, J . , C. Chung, C. Boast and M. Franklin, 1975. A simulation model of biophysiochemical transformations of nitrogen in tile-drained Corn Belt soi l . J . Environ. Qual., 4, 477-486. - Gilliam, J.W. and R.W. Skaggs, 1986. Controlled agricultural drainage to maintain water quality. J . Irrig. Drain. Eng., 112, 254-263. 176 Harris, G.L., M.J. Goss, R.J. Dowdell, R.R. Howse and P. Morgan, 1984. A study of mole drainage with simplified cultivation for autumn-sown crops on a clay soi l . J . Agric. Sci. Camb., 102, 561-581. Irwin, R.W. and H.R. Whiteley, 1983. Effects of land drainage on stream flow. Can. Water Resour. J . , 8, 85-103. Johnsson, H., L. Bergstrom, P.E. Jansson and K. Paustian, 1987. Simulated nitrogen dynamics and losses in a layered agricultural so i l . Agric. Ecosys. Environ., 18, 333-356. Kanwar, R.S., and H.P. Johnson, 1984. Sensitivity analysis of a drainage simulation model. Water Resources Bul l . , 20, 339-342. Kanwar, R.S., H.P. Johnson and J.L. Baker, 1983. Comparison of simulated and measured nitrate losses in ti le effluent. Trans. ASAE, 26, 1451-1457. Leeds-Harrison, P.W., G. Spoor and R.J. Goodwin, 1982. Water flow to mole drains. J . Agric. Eng. Res. 27, 81-91. Makkink, G.F., and M.D.J, van Heemst, 1975. Simulation of the water balance of arable land and pastures. Centre for Agr. Pub. Docum., Wageningen. Overcash, M.R., R.R. Reddy and R. Khaleel, 1980. Chemical and biological processes influencing the modelling of nonpoint source water quality from land areas receiving wastes. Proc. Hydrol. Transport Modeling. ASAE pub. 4-80, St. Joseph, 1-11. Petersen, A.P., 1986. Effect of rainfall erosivity and field slope on erosion of a Fraser Valley soi l . Unpub. B.Sc. Thesis, Univ. British Columbia. Robinson, M. and K.J. Beven, 1983. The effect of mole drainage on the hydrological response of a swelling clay s o i l . J . Hydrology, 64, 205-228. Schwab, G.O., N.R. Fausey and D.G. Kopcak, 1980. Sediment and chemical content of agricultural drainage water. Trans. ASAE, 23, 1446-1449. Skaggs, R.W. and J.W. Gilliam, 1981. Effect of drainage system design and operation on nitrate transport. Trans. ASAE, 24, 929-934. Switzer-Howse, K.D., 1983. Agricultural management practices for improved water quality in the Canadian Great Lakes basin. Agr. Can., Res. Branch publ. 1983-24-E. 177 - Van der Gulik, T.W., S.T. Chieng and F.J. Smith, 1986. B.C. Agricultural Drainage Manual. B.C. Min. Agric. Fish. , Agric. Eng. Branch, Abbotsford. - Williams, J .R. , C.A. Jones and P.T. Dyke, 1984. A modeling approach to determining the relationship between erosion and soil productivity. Trans. ASAE, 27, 129-144. - White, R.G., 1985. The analysis of solute breakthrough curves to predict water redistribution during unsteady flow through undisturbed structured clay soils. J . Hydrol., 79, 21-35. - Wischmeir, W.M. and D.D. Smith, 1978. Predicting rainfall erosion losses: a guide to conservation planning. USDA Agric. Handbook 537. 178 Concluding Remarks The conclusions derived from the results of a model can only be as valid as the premises used to develop the model. The present work cannot escape this fundamental fact and the conclusions and findings discussed here are limited in their application to the limits of validity of the assumptions underlying the model. It has been assumed that the mechanisms simulated here present a realistic picture of what can actually happen in field soi ls, even if the parameters that govern the relative magnitude and extent of these mechanisms cannot be evaluated with any degree of accuracy. This assumption fails if there are important mechanisms that have been omitted, or if the spatial distribution of the soil properties is such that a one-dimensional simulation cannot reproduce field phenomena in meaningful fashion. The quality of the f i t between the model predictions and the results from the experimental field offer some degree of assurance that this assumption is justified in the present context and, furthermore, that the model parameters present some degree of veracity. However, independently of this fact, the main objective of this work was to demonstrate that drainage practices can have an important effect on nutrient losses. This simulation shows that there can be a case where this assertion is verified, whether or not the simulated case corresponds closely to the experimental conditions. More appropriately, the type of drainage design has been shown 179 by this simulation to play an important role in determining the extent of the losses of nutrients. What conditions determine the importance of drainage design has also been partially investigated and the simulation has shown conclusively that the mechanisms that determine soil moisture movement, and the vertical distribution of soil nutrients, are important factors, at least under transient conditions. Several more factors need to be investigated in order to properly circumscribe the question of the impact of drainage practices on nutrient losses. Simulation for a longer duration would be required to test the validity of the present conclusions under steady state conditions (if these can be reached). Tillage practices, ferti l ization timing and method of application, crop type and cover, to name a few, are all factors that can be expected to interact with the response of leaching to drainage practices. The effect of soil type (as opposed to simply transfer mechanism type) is undoubtedly significant, as is the effect of the type of climate. To address these factors would require a simulation effort that was not possible within the scope of this work, while also requiring improvements of the limitations of the model. Despite these limitations, the present work shows that the simulation methodology used here can be adapted successfully to this type of problem and that the model developed here provides a useful and flexible tool when used within its range of application. 180 Appendix I An example of moisture flow and redistribution simulation The sequence of moisture content profiles depicted in f ig . A-l i l lustrates a typical example of moisture flow in the s o i l , including inf i l t rat ion, redistribution, and bypass. The soil chosen for the simulation has a saturated moisture content of 0.24 (6S) while its micropore matrix is saturated at a moisture content of 0.22 (6PR). The saturated conductivity of the matrix is 3mm/day from the surface down to a depth of 30 cm and 1 mm/day below. This layering is arbitrary and was chosen for the purpose of illustrating the generation of bypassing flow. The parameter PEDS was assigned a very large value so that vertical flow can be clearly outlined. All other parameters are similar to those found in table 4.1. Init ial ly, the soil is assumed to have a uniform moisture content of 0.16 above drain depth (1.1m deep) while being saturated below. Evapotranspiration is set to zero and the simulation is carried out under a constant rainfall rate of 2 mm/day. In response to this rainfall rate and to the initial conditions, the water infiltrates from the surface, bringing the top layers to flux equilibrium (6FL = 0.21) while upward flow in the lower unsaturated layers progressively lowers the water table. This can be seen for the first seven days. Eventually the wetting front reaches the lower permeability horizon. At this point water infiltrates into this lower horizon MOISTURE PROFILES 181 Figure A-l Example of moisture flow and redistribution 1n the profile (Two viewing angles are provided 1n order to h igh l igh t the upward c a p i l l a r y flow and the Infiltration front, as well as the water table depth, denoted by V.) 182 MOISTURE P R O F I L E S Figure A-l Example of moisture flow and redistribution 1n the profile (Two viewing angles are provided In order to h igh l igh t the upward c a p i l l a r y flow and the Infiltration front, as well as the water table depth, denoted by V.) 183 at a much reduced rate, while it starts accumulating above i t . This quickly saturates the matrix above this horizon, causing seepage into the macropore channels. This seepage water rapidly percolates down the profile in the channels, bypassing the lower permeability matrix and causing a rise in the water table. Meanwhile, the infiltration front is progressing through the lower permeability matrix, along with continuing upward flux from the saturated zone. This sequence happens in days 8 to 24. The infiltration front slowly saturates the matrix of the lower horizon while the water table rises above the drain level. As the water table level rises steady state is eventually achieved as the drain flow rate, determined by the water table level, reaches a value equal to the rainfall rate. At this point the flow is fully conducted in unsaturated mode in the upper horizon and in parallel flow mode (equal amounts of matrix flow and bypassing flow) in the lower horizon. 184 Appendix II A comparison of the dynamic equilibrium water movement model with a numerical solution of Richard's equation The infiltration-redistribution segment of the model was tested by comparing it with a finite difference formulation of which is Richard's equation expressed in terms of pressure potential, V. This equation is solved for where R is the infiltration rate at the top of the profile and Z is the depth of the prof i le . If the rainfall rate exceeds the inf i l t rabi l i ty of the so i l , the top boundary condition can be replaced by K ¥ = 0 @z = 0 The profile of the moisture content vs. pressure potential curve is assumed to be represented by 185 0 = 0S • B* which results in 96 = C = 0S • ln(B) • B* ay The conductivity IC is assumed to follow equation 2-15, which is expressed in terms of pressure potential, as K = Ks • B 0 1* The solution to Richard's equation was approximated using a f inite difference method. The solution is generated using the Thomas algorithm to solve the matrix of flow parameters calculated from the Crank-Nicholson formulation. For each time increment the values of K and C are estimated using a simple iterative predictor-corrector procedure. Details of this method of solution can be found, for instance, in Remson et a l . (1971). The infiltration-redistribution segment of the computer model was isolated so that the matrix flow algorithm can be examined separately from the effects of evapotranspiration, drain flow, or bypassing flow. Thus only infiltration, redistribution, and upward flow are simulated, according to the dynamic equilibrium concept outlined in Chapter 2. Both the model and the f in i te difference solution were determined for a 2.5m deep profile, segmented in 10 layers of 0.25 m of thickness. The profile was assured to be homogeneous, with 186 the following properties: Saturated moisture content, 0S = 0.24 Saturated conductivity, Ks = 0.010 m/day Moisture curve parameter 6 = 1.28 Conductivity exponent = 10 These parameters are similar to those used in the calibrated model (in chapters 4 and 5). The value of the coefficient 8 is such that it produces a moisture tension curve equivalent to that of figure 4-1. The value of the matrix conductivity Ks was increased by an order of magnitude from the original in order to provide a faster infiltration rate. The simulation was carried out for 20 days, with a rainfall rate of 5 mm/day for the f irst 10 days and no rain afterwards (R = 0 for boundary condition). The dynamic equilibrium model uses a one day time increment, while the finite difference solution is calculated every 0.05 day. The results of the comparison can be seen in f ig . A-2 and A-3. Both solutions give very similar results during the infiltration event, although the dynamic equilibrium solution predicts a sharper infi l tration front. Likewise, during the redistribution period (fig A-3) the finite difference solution determines a very smooth redistr ibution pro f i le , while the dynamic equilibrium model solution shows a routing of water through the layers with a progressively attenuated bulge. (In f ig . A-3, the initial moisture content curve refers to the moisture content on day 10 as calculated by the simulation model.) Both solution methods have inherent drawbacks. The dynamic M O I S T U R E C O N T E N T P R O F I L E S Figure A-2 Comparison of the model with a numerical solution of an infiltration event M O I S T U R E C O N T E N T P R O F I L E S Q. O m 0.1 0.2 ° FIN. OIF. SOLUTION, T - 15 DAYS * FIN. OIF. SOLUTION, T - 20 DAYS SIMULATION MOOD., f - 1? DAYS ^IfflL^ONJlUULLj T - IS DAYS " MOISTURE CONTENT — i — 0.3 0.4 CO CO Figure A-3 Comparison of the model with a numerical solution of redistribution following Infiltration 189 equilibrium model produces an inflexion point in the moisture curve at 2.25 m of depth. This is due to an underprediction of the depth of the water table, which causes a small amount of upward flow of water. This water is then removed from the layer where the water table is located, causing this small drop in moisture content. Because of the rounding error accumulated over the large number of calculations and of the nonlinearity of the moisture-pressure relationship, the finite difference solution does not conserve water. After 20 days there is 9.7 mm of water unaccounted for (lost from the profile). However, this amounts to only 1.8% of the total, which does not markedly impair the validity of the solution. It can be seen that the dynamic equilibrium model overpredicts the rate of redistribution. This arises because this model assumes that a unity gradient governs the unsaturated redistribution; however, it can easily be seen that this is not the case after 20 days. (A unity gradient implies that the moisture content remains constant with depth). Nevertheless the magnitude of the error remains small and is likely to be negligible if bypassing flow is an important solute transfer mechanism, as is the case in this study. This overprediction may become more significant in a soil with a high matrix conductivity, however. To test this possibility the comparison was repeated using a hydraulic conductivity and a rainfall rate increased to 100 mm/day and 50 mm/day, respectively, and restricting the rainfall duration to one day. The results from the model and the finite difference approximation are shown in fig. M O I S T U R E C O N T E N T P R O F I L E S —-INIJJflL MOISTURE: CONTENT ° rpg. crTr. SOLUTION, T - J r w o F I N . O I F . SOLUTION, T - 2 DAYS SIMULATION MODEL, f - 1 DRY MOISTURE CONTENT 0.3 0.1 Snarj2" ° , f 1 n f 1 l t r a t 1 ° " and redistribution for a soil with a larger hydraulic conductivity 191 A-4. Both solutions give similar results following infi ltration. The redistribution profiles after two days are also quite similar, with the exception that the model predicts an increased transmission rate to the water table, with its subsequent r ise. Again, this probably arises from the assumption of a unity gradient in the flux model. Thus it can be concluded that the dynamic flux equilibrium model, using a time increment of one day, can provide a satisfactory approximation of the solution to Richard's equation, at least for the condition tested. Reference Remson, I., G.M. Hornberger and F .J . Molz, 1971. Numerical methods in subsurface hydrology. Wiley-Interscience, New York. 192 Appendix III Diffusion between the macropore and micropore zones In this model, a simple mixing depth approach is used to model diffusion. The purpose of this appendix 1s to show the analytical rationale for such an approach. Water held in the ped is assumed to mix with macropore water in proportion to a ratio of widths, where W 1s the mixing width and L is the total width of the ped. This 1s tantamount to where C a v g , an average dimensionless concentration, represents the ratio of C i and C f , the initial and final ped concentrations, and C 0 , the Initial macropore concentration. This equation can be compared to the analytical solution for diffusion at the boundary of two semi-infinite domains, in one dimension, when one domain has a unit dimensionless concentration and the other a concentration of zero. For a constant diffusion coefficient, the solution gives (Luikov, 1968, p. 402): 1 C(x,t) = — (1 - erf X) where X = x / (2 /Dt ) where x is the distance from the boundary, t the time, and D a diffusion coefficient. If we assume that the ped width, I, is suff iciently wide to justify the assumption of semi-infinite domain, then the average concentration at a time t and distance x is 193 4 C(x,t)dx ^avg J dx o Using a series expansion of the error function, this averaging can be easily carried out to yield 1 1 / X X3 X5 X7 C a v g " 2 " / I T \ 2~ 4-3" 6»5«2!~ 8-7»3! When evaluated at L, and given that C a v g = W/L, one gets for W L I / L2 L4 L6 W = — - - — + 2 /n \ 4 /Dt 96(Dt)K 5 1920(Dt)2,5 and thus it can be seen that there is a value of the mixing width W that can be used to represent the average concentration resulting from diffusion for any set of values of L, D and t. Reference Luikov, A.V., 1968. Analytical heat diffusion theory. Academic Press, New York. 194 Appendix IV A comparison of the model with an analytical solution Van Ommen (1985) has published an analytical solution to determine the solute concentration of drainwater for an idealized situation. The flow domain he used is a one-dimensional set of layers similar to that shown in f ig. 2-1, except that there is no macropore domain. Instead, vert ical flow is arb i t rar i ly partitioned between a matrix flow and a bypassing flow, where the matrix flow occurs through the layered domain above the drains, while the bypassing flow short-circuits the layers above the drains and is added directly to the larger layer below the drains. One consequence of this simplification is that there is no flow or solute movement between macropores and micropores. Water flow is assumed to be at steady-state, while solutes flow in response to a unit pulse at the soil surface of arbitrary duration. Solute flow arises from advection and mixing caused by the discretization of the flow domain. The solute concentration at the drain outlet, c(t), can be calculated as: c(t) = l( l-f) • (X + u-Y) + f • X] • exp (-ut) where X = exp(ut) - 1 and Y = l_n - vm • exp((u-v) • t) • Z + n 195 with tn-m • (-l)J Z = I j=o (m-j)I • (u-v)J + 1 and m = n - k u = Q/(9 • 0Z)b v = Q • (1 - f)/(9 • DZ)a In this equation n refers to the number of layers above the drain level, less one; f refers to the ratio of bypassing flow to the total flow; Q denotes the flow rate through the so i l , while t is the time; f inally, 9 and DZ represent the moisture content and the thickness of the layers above (subscripted a) and below (subscripted b) the drain l ine. This solution is obtained by assuming that the solutes are cascaded from one layer to the next, and that each layer has a homogeneous solute concentration. The drain solute concentration predicted by this equation can be seen in f ig . A-5 for values of f of 0.0 (no bypass flow), 1.0 (bypass flow only), and 0.5 (equal proportions). The rain has a constant rate of 2 mm/day and has a unite solute concentration for the initial hundred days. The drain line is at a depth of 1.1 m, above which four layers of equal thickness are assigned a moisture content of 0.19. The layer below drains has a thickness of lm and a moisture content of 0.2. The simulation model was modified in order to reproduce the situation analyzed by Van Ommen (1985). The sections dealing with daily weather, evapotranspiration, plant growth, and biochemistry were removed, leaving the water flow and solute flow subsections. The value of the parameter PEDS was increased in order to minimize BREAKTHROUGH CURVES AT DRAIN • ANALYTICAL, F-0.0 SIMULATION, CS-2.0 o ANALYTICAL, F-0.5 SIMULATION. CS-1.0 A ANALYTICAL, F-1.0 SIMULATION, CS-0.0 INPUT PULSE i l j i i i 1 0.0 100.0 200.0 300.0 400 DAYS Figure A-5 Comparison of the model output with an analyt ical solution 197 the flux between the macropore and micropore zones. The matrix conductivity CS was assigned the value of 2 mm/day, 1 mm/day, or 0 mm/day in order to reproduce the effect of the ratio f in the analytical solution. The geometry of the layering and the water content was identical to that used in the analytical solution. The results of the simulation runs can be compared in f ig . A-5. As seen in this figure, the results from the model are nearly identical to the analytical solution when there is no bypassing flow. This is to be expected as the situation simulated is the same in both cases. By contrast, some difference between the two solutions can be seen when bypassing flow is included. This arises because the model simulates bypassing flow along discrete layers, which causes numerical dispersion and retardation similar to that occurring in the matrix zone. Accordingly, the model does not reproduce the situation described by the analytical solution, but proposes a different way of representing bypassing flow, which may arguably be considered more realistic. If the model is further simplified so that bypassing flow is directly added to the lower layer, as in Van Ommen's solution, its results and the analytical solution are, again, nearly identical. Reference Van Ommen, H.C., 1985. Systems approach to an unsaturated-saturated groundwater quality model, including adsorption, decomposition and bypass. Agric. Water Mana., 10, 113-203. 

Cite

Citation Scheme:

        

Citations by CSL (citeproc-js)

Usage Statistics

Share

Embed

Customize your widget with the following options, then copy and paste the code below into the HTML of your page to embed this item in your website.
                        
                            <div id="ubcOpenCollectionsWidgetDisplay">
                            <script id="ubcOpenCollectionsWidget"
                            src="{[{embed.src}]}"
                            data-item="{[{embed.item}]}"
                            data-collection="{[{embed.collection}]}"
                            data-metadata="{[{embed.showMetadata}]}"
                            data-width="{[{embed.width}]}"
                            async >
                            </script>
                            </div>
                        
                    
IIIF logo Our image viewer uses the IIIF 2.0 standard. To load this item in other compatible viewers, use this url:
http://iiif.library.ubc.ca/presentation/dsp.831.1-0098263/manifest

Comment

Related Items