UBC Theses and Dissertations

UBC Theses Logo

UBC Theses and Dissertations

A study of soil temperature and seedling survival in a forest clearcut Stathers, Robert John 1983

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

Item Metadata

Download

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

Full Text

A STUDY OF SOIL TEMPERATURE AND SEEDLING SURVIVAL IN A FOREST CLEARCUT By ROBERT JOHN STATHERS B.Sc.(Agr.)» The University of B r i t i s h Columbia, 1981  A THESIS SUBMITTED IN PARTIAL FULFILLMENT OF THE REQUIREMENTS FOR THE DEGREE OF MASTER OF SCIENCE ln THE FACULTY OF GRADUATE STUDIES Department of S o i l Science  We accept this thesis as conforming to the required standard  THE UNIVERSITY OF BRITISH COLUMBIA /  August 1983 (g) Robert John Stathers, 1 9 8 3  E-6  In p r e s e n t i n g  this  thesis i n partial  f u l f i l m e n t of the  r e q u i r e m e n t s f o r an a d v a n c e d d e g r e e a t t h e of B r i t i s h Columbia, I agree that it  freely  the L i b r a r y  a v a i l a b l e f o r r e f e r e n c e and s t u d y .  agree that p e r m i s s i o n f o r e x t e n s i v e for  University  s c h o l a r l y p u r p o s e s may  for  financial  shall  of  So'A  Sg-^ice.  The U n i v e r s i t y o f B r i t i s h 2075 W e s b r o o k P l a c e V a n c o u v e r , Canada V6T 1W5 Date  (2/79)  OCT. l l ,  R 63  thesis  Columbia  my  It is thesis  n o t be a l l o w e d w i t h o u t my  permission.  Department  further  be g r a n t e d by t h e h e a d o f  copying or p u b l i c a t i o n of this  gain  I  make  copying of t h i s  d e p a r t m e n t o r by h i s o r h e r r e p r e s e n t a t i v e s . understood that  shall  written  ABSTRACT  This thesis reports the results of a study of s o i l temperature and seedling survival on a steep south facing, high elevation forest clearcut on Vancouver Island, B.C., having a history of regeneration f a i l u r e s . During the f i r s t year of the study (1981) there were harsh weather conditions during the summer, while during the second year (1982) summer weather conditions were much milder.  Lethally high surface s o i l temperatures were  found to be the major cause of seedling mortality on this s i t e .  Of 2160  Douglas-fir, western hemlock and P a c i f i c s i l v e r f i r seedlings planted i n A p r i l 1981, 87, 33 and 17% respectively had survived by the end of the second growing season.  Shade cards s i g n i f i c a n t l y increased survival f o r  each of these species.  Of an additional 960 Douglas-fir, western hemlock,  and P a c i f i c s i l v e r f i r seedlings planted i n A p r i l 1982, 92, 74, and 88% respectively had survived by October 1982.  Shade cards, t r i c k l e  irriga-  tion, and shade cards and t r i c k l e i r r i g a t i o n combined only s i g n i f i c a n t l y increased survival i n western hemlock i n the seedlings planted i n 1982 because of the high o v e r a l l survival rate resulting from the mild summer weather  conditions.  A physically based s o i l temperature model was developed to ( i ) quantify the effects of the factors which influence the s o i l thermal regime and ( i i ) to model the time course of s o i l p r o f i l e temperatures from e a s i l y measured or estimated meteorological data and s i t e s p e c i f i c characteristics.  The model, based upon a f i n i t e difference solution of the heat flow  equation, energy and radiation balance, and aerodynamic theory, calculates surface and subsurface temperatures i n s o i l s having thermal properties that  Iff  vary with depth.  A s l i g h t l y modified version of the atmospheric  stability  correction method of Paulson (1970) was required to accurately predict surface temperatures under unstable atmospheric conditions. The model was  tested with data presented i n two other energy balance  studies reported i n the l i t e r a t u r e and found to estimate s o i l very accurately.  temperatures  In addition, good agreement was found between modelled  and measured s o i l p r o f i l e temperatures for 6 and 16 day periods i n the forest clearcut mentioned above. Energy balance theory was used to show that maximum surface s o i l temperatures i n a forest clearcut are not l i k e l y day.  Wind speed, aerodynamic  to exceed 85°C on a clear hot  roughness of the surface, and s o i l thermal  properties were also shown to strongly influence s o i l  temperatures.  TABLE OF CONTENTS Page ABSTRACT  JL  TABLE OF CONTENTS  J& JZD  LIST OF TABLES LIST OF FIGURES  IZfflt  LIST OF SYMBOLS  SJL  ACKNOWLEDGEMENTS  SSS  INTRODUCTION CHAPTER 1 .  ± EFFECT OF MICROCLIMATE MODIFICATION ON THE SURVIVAL OF 3 SPECIES OF CONIFER SEEDLINGS 4  1.1  Introduction  1.2  Experimental  Site  1.3  Experimental  Design and Procedure  1.4  Results and Discussion  CHAPTER 2 .  £  A.  Survival of Seedlings Planted i n 1981  JO  B.  Survival of Seedlings Planted i n 1982  13  C.  The Effect of Weather on Seedling Survival  73  D.  The E f f e c t s of Shade and I r r i g a t i o n on S o i l Temperature  1.5  $  Conclusions  2.1 25  FACTORS AFFECTING SURFACE SOIL TEMPERATURES IN A FOREST CLEARCUT  2.1  Introduction  2.2  Theory A.  The Energy and Radiation Balances of a Forest Clearcut  29  Z9  is: Page  2.3  B.  The S o i l Heat Flux  C.  Aerodynamic Equations Describing the Sensible Heat Flux, Latent Heat Flux and Momentum Flux i n the Turbulent Boundary Layer  Surface Energy Balance and Clearcut Microclimate ... 4/ a) b) c) d)  B.  C. D.  The Radiation Balance The S o i l Thermal Regime The S o i l Water Regime Boundary Layer Resistance to Heat, Water Vapour and Momentum Transfer  Comparison of Atmospheric S t a b i l i t y Correction Methods  4g  Examination of Micrometeorological Factors 53 §"J  Conclusions  MODELLING SOIL TEMPERATURE IN A FOREST CLEARCUT Introduction  3.2  Theory  3.3  Procedure  3.5  45 45  3.1  3.4  4/ ^£ 4.3  Relationship Between Surface Energy Balance Components and Surface Temperature  Affecting Surface Temperature  CHAPTER 3.  32-  Results and Discussion A.  2.4  31  50 67  A.  Assumptions Used i n the Model  64  B.  The Solution Method  65  Results and Discussion A.  Validation of the Model  57  B.  Sensitivity  "77  Conclusions  Analysis of the Model  a-7  Page CONCLUDING REMARKS  <\<\  REFERENCES  J 01  APPENDIX A. APPENDIX B.  Derivation of the 1-Dimensional F i n i t e Difference Solution of the S o i l Heat Conduction Equation  jQfa  Modelling Clear Sky Solar Irradiance on Variable Slopes and Aspects  //£.  APPENDIX C.  S o i l Thermal Properties  JJQ  APPENDIX D.  A Listing  of the Computer Program Used to Calculate  S o i l Temperatures and Surface Energy Fluxes  723  VTT  LIST OF TABLES Page CHAPTER 1. TABLE 1.1  Monthly percentage survival of seedlings which were planted i n A p r i l 1981  77  TABLE 1.2  Percent survival i n October of 1981 and 1982 of seedlings planted i n A p r i l 1981. Numbers followed by the same l e t t e r are not s i g n i f i c a n t l y different (t = 0.05) according to Duncan's New Multiple Range Test  11  TABLE 1.3 TABLE 1.4  Monthly percentage survival of seedlings which were planted i n A p r i l 1982 , Percent survival i n October 1982 of seedlings planted i n A p r i l 1981. Numbers followed by the same l e t t e r are not s i g n i f i c a n t l y different (t = 0.05) according to Duncan's New Multiple Range Test  H  15  CHAPTER 2. TABLE 2.1  Thermal Properties of Natural Materials ..  44  TABLE 2.2  Aerodynamic Properties of Natural  46  Surfaces  3zm  LIST OF FIGURES Page CHAPTER 1. FIGURE 1.1  Daily values of p r e c i p i t a t i o n , solar irradiance, average a i r temperature, and average r e l a t i v e humidity during the 1981 growing season  FIGURE 1.2  P r o f i l e s of gravimetric s o i l water content at selected times during the summer of 1981  FIGURE 1.3  Daily values of p r e c i p i t a t i o n , solar irradiance, average a i r temperature, and average r e l a t i v e humidity during the 1982 growing season  FIGURE 1.4  FIGURE 1.5  FIGURE 1.6  The diurnal course of s o i l temperature at several depths, surface s o i l temperature (measured with a fine (0.075 mm diameter) wire copper constantan thermo-couple) and a i r temperature at the 1 m height on August 7, 1982 The diurnal course of surface s o i l temperature adjacent to shaded, i r r i g a t e d , shaded and i r r i g a t e d , and control treatments on August 7, 1982. The course of a i r temperature at the 1 m height i s also shown The diurnal course of s o i l temperature i n the seedling root zone adjacent to shaded, i r r i g a t e d , shaded and i r r i g a t e d , and control treatments on August 7, 1982. Numbers above the p r o f i l e refer to the time of day (P.S.T.)  J7  20  23  Z4  26  CHAPTER 2. FIGURE 2.1  FIGURE 2.2  Integrated diabatic p r o f i l e correction factors for heat and momentum versus 6 or R^ (indices of atmospheric s t a b i l i t y ) f o r unstable, neutral, and stable atmospheric conditions Surface s o i l temperature versus solar irradiance calculated f o r a range of values of H + LE + G from energy balance theory (equation 2.7) assuming a = 0.18, e = 9.7, c = 0.85, and T = 25°C  40  D  s  a c  a  41  IX  FIGURE 2.3  Surface s o i l temperature versus wind speed calculated using energy balance and aerodynamic theory and the s t a b i l i t y correction methods of Paulson (1970), Thorn et al.(1975), Holbo (1971), Fuchs et al,(1968), and Hammel et al.(1982). Values were calculated assuming R = 1000 W m , G + LE = 150 W m , T = 25°C, ct = 0.18, e = 0.97, e = 0.85 and z = 0.01 m. The methods of Fuchs et a l and Hammel et al„give i d e n t i c a l results (see t e x t ) . . -2  s  -2  a  s  a  Q  FIGURE 2.4  Calculated versus measured surface s o i l temperatures on August 9, 1981. Calculated values were obtained from measured R , T , u, G and LE/R = 0.05 using the s t a b i l i t y correction methods of Paulson (1970), Thorn et al,(1975), Fuchs et al.(1968) and Hammel et al (1982). Surface s o i l temperatures were measured with an infra-red thermometer and a fine wire thermocouple s  a  Q  g  e  FIGURE 2.5  Surface s o i l temperature versus wind speed at at the 1 m height for various solar irradiances. Values were calculated assuming G + LE = 150 W m , T = 25°C, z = 0.01 m, a = 0.18, e = 0.97, and e = 0.85. -2  Q  a  Q  s  a  FIGURE 2.6  Surface s o i l temperature versus wind speed at a i r temperatures of 5, 25, and 45°C. Values were calculated assuming R = 1000 W m" , G + LE = 150 W m , z = 0.01 m, a = 0.18, e = 0.97, and e = 0.85. 2  s  -2  0  s  a  FIGURE 2.7  Surface s o i l temperatures versus wind speed at aerodynamic roughness lengths of 0.001, 0.01, and 0.1 m. Values were calculated assuming R = 1000 W m , G + LE = 150 W m , T = 25°C, a = 0.18, e = 0.97, and e = 0.85 s  -2  -2  Q  s  a  a  CHAPTER 3 . FIGURE 3.1  Flow chart showing the solution method used to calculate surface and subsurface s o i l temperatures.  FIGURE 3.2  Comparison of modelled and measured (Stearns 1967) temperatures i n a dry Peruvian desert s o i l having uniform thermal properties with depth on July 11, 1964  FIGURE 3.3  Comparison o f m o d e l l e d and measured (Novak 1981) t e m p e r a t u r e s a n d s u r f a c e e n e r g y f l u x e s on May 1 8 , 1979 i n a w e t s i l t l o a m s o i l a t A g a s s i z , B.C  FIGURE 3.4  C o m p a r i s o n o f m o d e l l e d and measured (Novak 1981) t e m p e r a t u r e s a n d s u r f a c e e n e r g y f l u x e s on J u l y 2 1 , 1979 i n a s i l t l o a m s o i l a t A g a s s i z , B.C. w h i c h h a d d r i e d near the s u r f a c e  FIGURE 3.5  C o m p a r i s o n o f m o d e l l e d ( 0 ) a n d m e a s u r e d (+) s o i l t e m p e r a t u r e s a t 4 d e p t h s i n t h e Cameron V a l l e y forest c l e a r - c u t experimental s i t e during the 6 h o t t e s t days i n 1981  FIGURE 3.6  C o m p a r i s o n o f m o d e l l e d ( 0 ) a n d m e a s u r e d (+) s o i l t e m p e r a t u r e s a t 3 d e p t h s i n t h e Cameron V a l l e y f o r e s t c l e a r - c u t e x p e r i m e n t a l s i t e d u r i n g 16 d a y s of v a r i a b l e weather c o n d i t i o n s i n J u l y 1981  FIGURE 3.7a  The e f f e c t o f t h e a e r o d y n a m i c r o u g h n e s s ( z ) on m o d e l l e d s u r f a c e s o i l t e m p e r a t u r e . Values of z u s e d w e r e 5 mm ( x ) , 10 mm ( + ) , a n d 50 mm (0) w h i l e a l l o t h e r v a r i a b l e s were h e l d c o n s t a n t . Q  Q  FIGURE 3.7b  A s f o r F i g u r e 3.7a b u t f o r t e m p e r a t u r e s a t t h e 0.1 m d e p t h  FIGURE 3.7c  As f o r F i g u r e 3.7a b u t f o r t e m p e r a t u r e s a t t h e 0.3 m d e p t h  FIGURE 3.7d  As f o r F i g u r e 3.7a b u t f o r t e m p e r a t u r e s a t t h e 0.5 m d e p t h  FIGURE 3.8a  The e f f e c t o f t h e l a t e n t h e a t f l u x d e n s i t y ( L E ) on m o d e l l e d s u r f a c e s o i l t e m p e r a t u r e . Values of L E u s e d w e r e 0.0 ( + ) , 0.05 R ( x ) , a n d 0.2 R (0) w h i l e a l l o t h e r v a r i a b l e s were h e l d c o n s t a n t . s  s  FIGURE 3.8b  As f o r F i g u r e 3.8a b u t f o r t e m p e r a t u r e s a t t h e 0.1 m d e p t h  FIGURE 3.8c  As f o r F i g u r e 3.8a b u t f o r t e m p e r a t u r e s a t t h e 0.3 m d e p t h  FIGURE 3.8d  As f o r F i g u r e 3.8a b u t f o r t e m p e r a t u r e s a t t h e 0.5 m d e p t h  Page FIGURE 3.9a  The effects of s o i l thermal properties on modeled surface s o i l temperature. Values of k and C i n the 0-0.02, 0.02-0.05 and 0.05-0.8 m s o i l layers were 1.0, 1.2 and 1.4 MJ m ° C and 0.19 0.29 and 0.67 W m °C respectively (X), values of k are the same as those given above but values of C are twice those given above (+) and values of k and C are twice those given above (0). A l l other variables were held constant. s  s  -3  _1  _ 1  _ 1  s  s  s  FIGURE 3.9b FIGURE 3.9c FIGURE 3.9d  s  As f o r Figure 3.9a but f o r temperatures at the 0.1 m depth  tyQ  As for Figure 3.9a but for temperatures at the 0.3 m depth  9/  As f o r Figure 3.9a but f o r temperatures at the 0.5 m depth  FIGURE 3.10a  FIGURE 3.10b  The effect of the temperature at the lower surface of the modelled s o i l slab (0.8 m depth) [lower boundary condition] on modelled surface s o i l temperature. Values of 10°C (+), 12°C (x), and 14°C (0) were used while a l l other variables were held constant  FIGURE 3.10d  92  93  As for Figure 3.10a but f o r temperatures at the 0.1 m depth  FIGURE 3.10c  Qfl  9^  As for Figure 3.10a but for temperatures at the 0.3 m depth  95  As for Figure 3.10a but f o r temperatures at the 0.5 m depth  96  ixn  L I S T O F SYMBOLS  C  volumetric  C  c l o u d cover  C  s  F G  volumetric  heat c a p a c i t y  (J m  °C  )  dimensionless heat c a p a c i t y o f s o i l  atmospheric s t a b i l i t y  (J m  correction factor  °C  )  dimensionless  s o i l heat f l u x d e n s i t y  (W m  )  H  s e n s i b l e heat f l u x d e n s i t y  (W m  )  H'  s e n s i b l e heat f l u x d e n s i t y under  (W m  )  Q  -  neutral  conditions Kg  eddy d i f f u s i v i t y  f o r s e n s i b l e heat  (m  s  )  KM  eddy d i f f u s i v i t y  f o r momentum (eddy v i s c o s i t y )  (m  s  )  Ky  eddy d i f f u s i v i t y  f o r water vapour  (m  s  )  L  Monin-Obukhov s t a b i l i t y  LE  l a t e n t heat f l u x d e n s i t y  LE'  l a t e n t heat f l u x d e n s i t y under  length  ( p c T u* /(kgH)) 3  p  a  neutral  (m) (W m  )  (W m  )  conditions Ri  gradient  R i c h a r d s o n number  R^  longwave i r r a d i a n c e  dimensionless —2 (W m ) _ o  R  n  net r a d i a t i o n f l u x d e n s i t y  (W m  )  R  s  solar irradiance  (W m  )  temperature  (°C) (°C, K)  T T  a  a i r temperature  T  a  mean a i r temperature between z and z  T  s  soil  s u r f a c e temperature  Q  (K) (°C, K)  -XZ7J  c  r a t i o of LE to R  dimensionless  s  C p  a  s p e c i f i c heat of a i r  (kj kg *  °C )  Cp  S  s p e c i f i c heat of s o i l  (kJ k g  °C )  e e  -  vapour p r e s s u r e of water i n a i r  - 1  - 1  (kPa)  vapour p r e s s u r e of water i n a i r a t s o i l  Q  - 1  surfce  (kPa)  e*(T )  s a t u r a t i o n vapour p r e s s u r e of water  (kPa)  g  gravitational acceleration  (m s  k  von Karman's c o n s t a n t  dimensionless  k  thermal c o n d u c t i v i t y  (W m  -1  "C"" )  thermal c o n d u c t i v i t y of s o i l  (W m  -1  °C )  u  wind v e l o c i t y  (m  u*  friction  (m s  z  depth below or h e i g h t above s o i l s u r f a c e  s  k  s  velocity  )  1  - 1  s *) -  - 1  ),  (both p o s i t i v e )  (m)  aerodynamic roughness l e n g t h  (m)  a  shortwave r e f l e c t i v i t y  dimensionless  Y  psychrometric  e  v i s c o u s d i s s i p a t i o n of k i n e t i c energy  dimensionless dimensionless  z  Q  (albedo)  constant  (Pa  _ 1  e  a  atmospheric  e  a c  c l e a r sky e m i s s i v i t y  dimensionless  G  S  soil  dimensionless  £ p  emissivity  °C )  surface emissivity  index of atmospheric a  d e n s i t y of a i r  stability  (z/L)  dimensionless (kg m  q  )  Ps  density of s o i l  (kg  a  Stephan-Boltzmann constant  (W m  T  momentum flux  (N  density  m" ) 3  -2  lC ) k  m ) -2  diabatic influence function f o r sensible •H  heat transfer  dimensionless  diabatic influence function for momentum •' IM  transfer  dimensionless  diabatic influence function for water <t>v  vapour transfer  dimensionless  ACKNOWLEDGEMENTS  I wish to express my s i n c e r e g r a t i t u d e and a p p r e c i a t i o n  to Dr. T.A.  B l a c k f o r h i s i n v a l u a b l e guidance and s u p e r v i s i o n d u r i n g my M.Sc. program. I a l s o wish to thank M a c M i l l a n B l o e d e l L t d . , Woodlands S e r v i c e s D i v i s i o n for  funding  this  research.  My a p p r e c i a t i o n goes t o the other members of my committee: Novak ( A s s i s t a n t P r o f e s s o r ,  Department o f S o i l S c i e n c e ) , Dr. T.M. B a l l a r d  ( P r o f e s s o r , Department of S o i l S c i e n c e ) , Silviculturist,  Dr. M.D.  MacMillan Bloedel  and Mr. B.G. Dunsworth  (Research  Ltd.).  I extend my thanks t o the biometeorology group a t U.B.C.:  Dr. D.L.  S p i t t l e h o u s e , Mr. F. K e l l i h e r , Mr. D.T. P r i c e , Mr. D.G. G i l e s , and i n p a r t i c u l a r , Mr. N.J. L i v i n g s t o n , f o r t h e i r v a l u a b l e  discussions  and a s s i s t -  ance . I a l s o wish to thank the Department of S o i l  Science  for their  accept-  ance of me i n t o the M.Sc. program, and Ms. D. Shunamon f o r her e x c e l l e n t word p r o c e s s i n g  of t h i s t h e s i s .  F i n a l l y , my wholehearted and s i n c e r e thanks goes t o my f a m i l y , and  daughter f o r t h e i r encouragement and p r e c i o u s  friendship.  wife  INTRODUCTION  S o i l temperature i s one of the major factors a f f e c t i n g conifer seedl i n g regeneration i n B r i t i s h Columbia.  In the north central i n t e r i o r  region, low s o i l temperatures i n h i b i t growth (Dobbs and McMinn 1977).  In  the southwest coastal region, l e t h a l l y high surface temperatures contribute to mortality (Arnott 1975). The problems associated with seedling regeneration on s i t e s prone to adverse microclimatic conditions have been recognized (Baker 1929; Various  Maguire 1955)  for a long time  and have prompted considerable  techniques such as shading (Hobbs et a l . 1980;  research.  Minore 1970), mulch-  ing  (Hermann 1964), i r r i g a t i n g (Livingston and Stathers 1982), and modify-  ing  the s o i l thermal properties (Hermann and Chilcote 1965)  to ameliorate  have been used  harsh environmental conditions.  The objectives of this study were to investigate p r a c t i c a l methods of improving seedling s u r v i v a l on a high-elevation south facing clearcut on Vancouver Island, B.C.  The experimental s i t e , located i n the Cameron  Valley near Port Alberni, t y p i c a l l y experiences harsh environmental conditions during the summer and was  replanted with d i f f e r e n t species at least 3  times before an adequate stocking l e v e l was  obtained.  Over the 2 year course of this study, 3 species of conifer seedlings were tested under a number of microclimate  modification treatments.  These  tests prompted a more detailed analysis of the factors influencing seedling microclimate was  and as a result a physically based computer simulation model  developed.  The results of these endeavours form the basis of the  chapters i n this thesis.  three  In Chapter 1, the results of seedling survival studies i n i t i a t e d i n 1981 and 1982 are reported.  The effects of weather and microclimate modi-  f i c a t i o n on survival are also discussed. In Chapter 2, the environmental factors influencing surface s o i l temperature are outlined.  Relationships between the energy and r a d i a t i o n  balances of the surface and surface s o i l temperature are developed.  Aero-  dynamic methods of c a l c u l a t i n g the sensible and latent heat fluxes from the surface are also discussed. In Chapter 3, a s o i l temperature model capable of predicting surface temperatures i s developed and tested.  The model i s used to ( i ) simulate  s o i l temperatures over extended time periods using easily measured climate data and ( i i ) simulate the effects of changing certain environmental v a r i ables on s o i l temperatures.  -3-  CHAPTER 1  EFFECT OF MICROCLIMATE MODIFICATION ON THE SURVIVAL OF THREE SPECIES OF CONIFER SEEDLINGS-  -4-  EFFECT OF MICROCLIMATE MODIFICATION ON  THE  SURVIVAL OF THREE SPECIES OF CONIFER SEEDLINGS  1.1  INTRODUCTION Plantation establishment  can be very d i f f i c u l t on the high elevation  south facing clearcuts on Vancouver Island.  The accumulated winter snow-  pack melts late i n the spring, leaving very l i t t l e time for newly planted seedlings to acclimatize to intense solar irradiance during the early summer.  The shallow, dense, coarse textured s o i l s commonly found at these  sites have a low water holding capacity and the seedling root zone dries rapidly a f t e r the snow melts.  When drought conditions p e r s i s t for more  than a few weeks during the summer, c r i t i c a l l y high surface s o i l temperatures and low s o i l water potentials can occur.  These conditions  are  thought to be the major cause of seedling mortality at high elevations i n this region (Arnott 1975). Many operational techniques have been applied to improve seedling s u r v i v a l where moisture and temperature stress occur.  The correct tree  species, provenance, and stock type selection have reduced seedling  mortal-  i t y ; however, under the extreme environmental conditions of summer drought, additional protective measures are necessary. microclimate  may  be the appropriate,  Modifying  the seedling  c o s t - e f f e c t i v e solution to the regen-  eration problems on these s i t e s . The objectives of this study were to assess the effects of moisture and temperature stress reduction on the s u r v i v a l of three seedling  species  on a high elevation south facing clearcut on Vancouver Island, B.C.  The  methods used to reduce stress included natural shade, a r t i f i c i a l shade,  -5-  i r r i g a t i o n , and a r t i f i c i a l shade and i r r i g a t i o n combined.  A r t i f i c i a l shade  from shade cards has been shown to increase survival on s t r e s s f u l sites i n Oregon (Hobbs et al,1980; Minore 1970) scribed operational treatment.  and thus was of interest as a pre-  Natural self shading by i n c l i n i n g the stem  to the south has also been e f f e c t i v e i n seedling nursery beds (Korstian and Fetherolf 1921)  and was  t i o n a l procedure. ative procedure  considered an inexpensive, easily applied opera-  I r r i g a t i o n has long been recognized as the most amelior-  for stress reduction (Baker 1929;  Maguire 1955); however,  i t i s usually not a cost e f f e c t i v e alternative i n remote areas.  Irrigation  treatments were applied i n this study to determine seedling performance under conditions of minimal water stress, and also to determine whether water or thermal stress caused seedling mortality.  1.2  EXPERIMENTAL SITE The experimental s i t e was  located at an elevation of 1150 m on a 30°  south facing slope i n the Cameron Valley near Mt. Arrowsmith, Vancouver Island, B r i t i s h Columbia (49° 13'N,  124° 36'W).  The area i s i n the Mari-  time Coastal Western Hemlock biogeoclimatic subzone (Klinka et al,1970). The s i t e was very representative of the surrounding area. The s o i l at this s i t e i s an orthic Humo-ferric  Podzol.  It has devel-  oped over compact g l a c i a l t i l l and contains a large proportion of coarse fragments («60% > 2mm  i n diameter by dry weight).  It has a sandy loam  _ o  texture.  The bulk density i s 1.1-1.3 Mg m  near the surface and 1.6-1.8  —3 Mg m  below the 0.5 m depth.  These c h a r a c t e r i s t i c s give i t a very low  available water storage capacity.  A thin (<0.05 m) discontinuous layer of  organic matter i s present on the surface.  ~6~  The s i t e was  harvested  i n 1974.  The stand was  comprised predominantly  of Douglas-fir (Pseudotsuga menziesii (Mirb.) Franko), western hemlock (Tsuga heterophylla (Raf.) Sarg.) and yellow cedar (Chamaecyparis nootkatensis (D. Don)  Spach).  A sparse canopy of fireweed  folium) and s a l a l (Gaultheria shallon) now  (Epilobium  covers the s i t e .  angusti-  Douglas-fir  seedlings were planted at least 3 times prior to this study before adequate stocking l e v e l was t a l l when this study was  1.3  These seedlings were less than 1 m  initiated.  EXPERIMENTAL DESIGN AND PROCEDURE In late A p r i l 1981,  site.  achieved.  an  High-elevation  2160  seedlings were planted on the experimental  provenances of Douglas-fir, western hemlock, and  P a c i f i c s i l v e r f i r were selected for the t r i a l . lings for each species were planted.  1-0  The same number of seed-  container grown styroplugs were  used exclusively because of the operational advantages of producing and planting this stock type.  Seedlings were maintained i n a dormant state i n  cold storage for 4 months prior to planting. Nine combinations of tree species and microclimate ments were tested.  modification t r e a t -  An a r t i f i c i a l shade treatment, a natural  self-shade  treatment, and an unshaded control treatment were tested for each of the 3 species.  The a r t i f i c i a l shade treatment was  imposed by placing a 0.2  x  0.3 m waxed cardboard shade card (manufactured by International Forest Supplies, Eugene, Oregon) on the southwest side, 0.05 sloping approximately 10° back over the seedling.  m from the stem, and  This orientation was  necessary to e f f e c t i v e l y shade the root c o l l a r region because of the high solar zenith angle during the middle of the day. treatment was  The natural  self-shade  imposed by i n c l i n i n g the seedling stem to the southwest so  -7-  that the foliage shaded the root c o l l a r region at midday.  Unshaded seed-  lings were planted v e r t i c a l l y and the root c o l l a r region was exposed to full  sunlight. Microsite selection and planting techniques were similar for a l l  species and treatments.  Seedlings were planted away from logging debris  and i n s c a r i f i e d mineral s o i l .  Good hydraulic contact was established  between the root system and the s o i l by c a r e f u l l y hand-compacting the disturbed s o i l near the roots. The seedlings were l a i d out on the s i t e i n a randomized complete block design.  Each of the 9 species-treatment combinations was replicated i n 24  blocks.  In every block, a treatment consisted of 10 seedlings, planted  1.5 m apart i n the downslope d i r e c t i o n .  Consequently, there were 240 seed-  lings per treatment on the s i t e . In site.  1982, 4 additional blocks of seedlings were planted on the same The same 3 species, provenance, stock type, and planting technique  were used; however 4 treatments were tested.  Treatments included shade  cards to reduce heat stress as i n 1981, t r i c k l e i r r i g a t i o n to reduce plant water stress, a combination of shade and i r r i g a t i o n to minimize water and heat stress, and an unshaded unwatered control. A t r i c k l e i r r i g a t i o n system (Drip-Eze Systems, E l Cajon, C a l i f o r n i a , U.S.A.) was used to deliver water through drip emitters to i n d i v i d u a l seedlings (Livingston 1984).  Water was applied when tensiometers on the s i t e  indicated that the s o i l water potential was less than -0.05 MPa.  Arti-  f i c i a l shade was imposed using the same method described i n the 1981 t r i a l . Analysis of variance was used to determine s i g n i f i c a n t differences in survival between blocks and treatments.  Duncan's New Multiple Range  -8-  Test (t = 0.05) was used to rank each treatment. Seedling survival was assessed monthly from May through October i n 1981 and 1982 i n order to relate mortality to climatic conditions.  A  weather station was maintained on the s i t e from May u n t i l mid-September i n both years since summer conditions were suspected to be the most important factor influencing seedling s u r v i v a l . temperature  Solar irradiance, net radiation, a i r  (1 m height), r e l a t i v e humidity (1 m height), wind speed (1 m  height), p r e c i p i t a t i o n , s o i l temperature  (0.5 m, 0.3 m, 0.05 m, 0.015 m  depths, and the surface) and s o i l heat flux (0.03 m) were recorded on battery-powered Campbell S c i e n t i f i c CR21 data loggers (Campbell S c i e n t i f i c , Logan, Utah, U.S.A.).  The data loggers were programmed to record hourly  averages and store the data on audio cassette tapes.  Data loggers were  housed i n p l a s t i c (Tupperware) containers and kept dry with s i l i c a g e l desiccant.  The p l a s t i c containers were kept i n watertight metal m i l i t a r y  ammunition boxes which were locked and chained to an underground pad.  concrete  The boxes were located i n a covered 0.6 m deep hole to protect the  data loggers from the harsh climatic conditions at the surface and make them less vulnerable to vandalism. The pyranometer (Model Li-200SZ-05, Li-Cor Inc., Lincoln,  Nebraska,  U.S.A.) was calibrated with a laboratory standard Kipp and Zonen pyranometer at the beginning and end of each season. S-l,  The net radiometer (Type  Swissteco Pty. Ltd., Melbourne, Australia) c a l i b r a t i o n was determined  by comparison with similar model laboratory standard net radiometer. Net radiation was measured during July and August of both years. Relative humidity was measured with a sulphonated polystyrene sensor (Model 201, Phys-Chemical Research Corp., New York, New York, U.S.A.).  -9-  It was found to have a ±3% accuracy over the range from about 10-97%. R a i n f a l l was measured with a tipping bucket raingauge (Model RG2501, Sierra Misco, Inc., Berkeley, C a l i f o r n i a , U.S.A.) which had a resolution of 1 mm. Wind speed was measured at the 1.0 m. height with a cup anemometer (Model 014A.Met-One, Sunneyvale, C a l i f o r n i a , U.S.A.) which had a s t a l l i n g speed of 0.5 m/s.  Air and deep s o i l temperatures were measured with thermistor  (Fenwall UUT-51J1) thermometers which were accurate to ±0.5°C between 0 and 45°C.  Temperatures at and near the s o i l surface were measured with fine  (0.075 mm diameter) wire copper constantan thermocouples which were r e f e r enced to a deep s o i l thermistor thermometer. also measured with an infrared thermometer  Surface temperatures were  (Instatherm Model 14-2200-4  close focus, Barnes Engineering Co., Stanford, Conn., U.S.A.).  S o i l heat  flux at the 0.03 m. depth was measured with 3 heat flux plates connected i n series (made by Dr. M.D. Novak, Dept. of S o i l Science, U.B.C. following Fuchs and Tanner (1968)). Measurements of the s o i l water content and potential of the seedling root zone were collected on a weekly basis throughout the summer.  A neu-  tron probe (Model 503 hydroprobe, Campbell P a c i f i c Nuclear, Pacheco, C a l i fornia, U.S.A.) was used to monitor the s o i l water content over a range of depths from 0.15 to 0.75 m.  Gravimetric sampling was used to determine  seedling root zone water contents and to c a l i b r a t e the neutron probe.  Soil  water potentials were measured on an HR33T dew point hygrometer with 40 s o i l psychrometers (Model PT51-10, Wescor Inc., Logan, Utah, U.S.A.) placed  i n the seedling root zone at the 0.15 m depth.  The psychrometers were  calibrated i n the laboratory using a range of KC1 solutions.  Care was tak-  en to ensure good hydraulic contact between the s o i l and the horizontally  -10-  oriented psychrometers.  Readings were taken i n the morning before  large  s o i l temperature gradients occurred, using the dew point mode and correcting the cooling c o e f f i c i e n t  (II ) for temperature V  variation.  Seedling shoot water potentials were measured throughout the summer with a Scholander pressure bomb.  Seedling leaf conductances were also  measured i n s i t u using a transient state chamber porometer developed for this study and described i n Livingston et al,(1983).  1.4  RESULTS AND DISCUSSION  A.  Survival of Seedlings Planted In 1981 Survival of seedlings planted i n A p r i l 1981 i s summarized on a monthly  basis i n Table 1.1 for each of the 9 species and treatment combinations. After the f i r s t growing season a s i g n i f i c a n t difference i n survival was found between species (Table 1.2).  Untreated Douglas-fir (91%)  survived  s i g n i f i c a n t l y better than untreated western hemlock (38%) which i n turn survived s i g n i f i c a n t l y better than untreated  P a c i f i c s i l v e r f i r (23%).  Shade cards s i g n i f i c a n t l y increased survival  i n western hemlock by 35%  and P a c i f i c s i l v e r f i r by 23%, however, survival was not s i g n i f i c a n t l y increased i n Douglas-fir (8%) because of i t s o v e r a l l high survival The inclined  self-shade treatment had no s i g n i f i c a n t  rate.  effect i n the f i r s t  year. After the second growing season, the survival  of Douglas-fir (76%) was  s t i l l s i g n i f i c a n t l y better than western hemlock (15%) and P a c i f i c s i l v e r f i r (8%).  Treatment effects  also became more noticeable.  Shade cards  s i g n i f i c a n t l y increased survival by 21% i n Douglas-fir, 44% i n western hemlock, and 27% i n P a c i f i c s i l v e r f i r . The inclined  self-shade  treatment  TABLE 1.1  Monthly percentage survival of seedlings which were planted In A p r i l 1981.  WESTERN HEMLOCK  DOUGLAS-FIR Date  PACIFIC SILVER FIR Shade Cards  Inclined Shade  100  100  100  100  100  100  96  97  98  95  95  94  95  96  94  95  95  91  94  92  93  96  93  92  90  92  88  87  88  88  94  96  88  84  84  85  81  77  79  79  95  94  96  78  52  43  58  53  34  29  38  98  94  91  94  73  45  38  52  46  25  23  31  May 82  98  91  88  92  63  25  25  38  38  13  11  21  June 82  97  89  79  88  63  25  21  35  38  13  9  20  July 82  97  89  78  87  62  25  19  34  38  10  9  19  Aug 82  97  89  78  87  62  25  19  34  35  10  8  17  Oct 82  97  89  76  87  59  25  15  33  35  10  8  17  Shade Cards  Inclined Shade  100  100  100  100  98  99  99  99  99  97  98  98  July 81  99  95  94  Aug 81  99  95  Sept 81  98  Oct 81  Inclined Shade  Apr 81  100  100  May 81  99  June 81  Unshaded Control  Unshaded Control  Unshaded Control  Total  Total  Total  Shade Cards  -/2_  TABLE 1.2  Percent survival i n October of 1981 and 1982 of seedlings planted i n A p r i l 1981. Numbers followed by the same l e t t e r are not s i g n i f i c a n t l y different (t = 0.05) according to Duncan's New Multiple Range Test.  October 1981  TREATMENT  Douglas-Fir  Western Hemlock  CONTROL UNTREATED  91 (a)  38 (c)  23 (d)  INCLINED STEM SELF SHADED  94 (a)  45 (c)  25 (d)  SHADE CARDS ARTIFICIAL SHADE  98 (a)  73 (b)  46 (c)  Pacific Silver Fir  October 1982 >  TREATMENT  Douglas-Fir  Western Hemlock  Pacific Silver Fir  CONTROL UNTREATED  76 (b)  15 ( f )  8 (f)  INCLINED STEM SELF SHADED  89 (a)  25 (e)  10 ( f )  SHADE CARDS ARTIFICIAL SHADE  97 (a)  59 (c)  35 (d)  -73"  s i g n i f i c a n t l y increased s u r v i v a l by 13% i n Douglas-fir and 10% i n western hemlock; however, s u r v i v a l was  increased by only 2% i n P a c i f i c s i l v e r f i r .  By the second year seedlings i n the i n c l i n e d treatment had responded t r o p i c a l l y , so there was  geo-  no practical, difference between this treatment and  the control.  This suggests that the s i g n i f i c a n t e f f e c t of the i n c l i n e d  treatment was  fortuitous rather than due  B.  Survival of Seedlings Planted i n  to shading.  1982  Survival of seedlings planted i n A p r i l 1982 basis i n Table 1.3  i s summarized on a monthly  for each of the 12 species and treatment combinations.  At the end of the growing season, trends s i m i l a r to those observed i n were seen i n untreated lock.  and shade card treated Douglas-fir and western hem-  In contrast, P a c i f i c s i l v e r f i r s u r v i v a l was Survival of untreated  f i r seedlings was survival by 5%, by 4%,  26%,  90%,  dramatically improved.  Douglas-fir, western hemlock, and P a c i f i c  42%,  and 69% respectively.  53%, and 21% respectively.  and 24% respectively.  cards Increased  s u r v i v a l by 0%,  western hemlock (Table 1.4).  silver  Shade cards increased  I r r i g a t i o n increased s u r v i v a l  The combination of i r r i g a t i o n and shade  48%,  and 29% respectively.  The only s i g n i f i c a n t difference between treatments was  found i n  Other treatment differences were i n s i g n i f i -  cant because of the higher o v e r a l l s u r v i v a l rate and the high deviation between blocks where survival was  C.  1981  standard  lower.  The E f f e c t of Weather on Seedling Survival Weather was  found to be a major factor influencing seedling s u r v i v a l  at the experimental s i t e .  The highest rates of seedling mortality  occurred  TABLE 1.3  Monthly percentage s u r v i v a l of s e e d l i n g s which were p l a n t e d l n A p r i l 1982.  WESTERN  D O U G L A S - F I R  P A C 1: F I C  HEMLOCK  S I L V E R  I' I R  Date  Unshaded Control  Shade Cards  Irrigation  Shade & Irrigation  Unshaded Control  Shade Cards  Irrigation  Shade & Irrigation  Unshaded Control  Shade Cards  Irrigation  Shade & Irrigation  Hay 82  100  100  100  100  100  100  100  100  100  100  100  100  June 82  96  99  96  93  96  95  93  90  98  91  93  98  J u l y 82  87  98  96  90  75  95  80  90  94  91  93  98  Aug 82  83  95  95  90  42  95  68  90  79  91  93  98  Sept 82  83  95  94  90  42  95  68  90  74  90  93  98  Oct  90  95  94  90  42  95  68  90  69  90  93  98  82  -J5-  TABLE 1.4  Percent survival i n October of 1982 of seedlings which were planted i n A p r i l 1982. Numbers followed by the same l e t t e r are not s i g n i f i c a n t l y different (t = 0.05) according to Duncan's New Multiple Range Test.  October 1982  Douglas-Fir  Western Hemlock  CONTROL UNTREATED  90 (a)  42.(b)  69 (a)  SHADE CARDS  95 (a)  95 (a)  90 (a)  IRRIGATION  94 (a)  68 (a)  93 (a)  SHADE CARDS AND IRRIGATION  90 (a)  90 (a)  98 (a)  TREATMENT  Pacific Silver Fir  '16-  during extended periods of hot dry weather when the s o i l p r o f i l e water content was  low.  The 1981  growing season was  characterized by cool cloudy wet weather  from the time the seedlings were planted u n t i l mid-July. u n t i l September less than 10 mm weather occurred. during this period.  of r a i n f e l l and exceptionally hot clear  Mid-day maximum a i r temperatures ranged from 25-37°C A summary of p r e c i p i t a t i o n , solar radiation, a i r temp-  erature, and r e l a t i v e humidity Figure 1.1.  From mid-July  during the 1981  growing season i s shown i n  The course of s o i l p r o f i l e water contents i n 1981  showed very  l i t t l e change i n root zone water storage during the early part of the summer because of the high r a i n f a l l .  From mid-July  drawdown occurred i n the top 0.4 m (Figure 1.2).  to July 27 a rapid  This depletion continued  at a slower rate during the hot weather i n August u n t i l minimum s o i l water contents of <5% by dry s o i l weight were measured i n the top 0.1 m. i n storage were observed to a depth of 0.75  m.  S o i l water potentials i n  the seedling root zone (0-0.15 m) were found to be as low as -1.8 during this period. t i a l was  Changes  MPa  Very l i t t l e measurable difference i n s o i l water poten-  found between the untreated, shaded, and inclined  Daytime seedling water potentials as low as -3.4  MPa  treatments.  were also measured  during this period (Livingston and Stathers 1981). The majority of seedling mortality occurred during the hot weather and dry s o i l conditions i n August.  The seedling assessment after this period  showed a 41% decline i n untreated western hemlock survival and a 50% decline i n untreated P a c i f i c s i l v e r f i r survival (Table 1.1). withstood only a few  these harsh conditions very well and i t s survival was percent.  Douglas-fir reduced by  -11'  PRECIPITATION (MM)  J  ll, , I.I 60'  eo  RELATIVE  SO  HUMIDITY  (%)  TOTAL SOLAR RADIATION  (MT/ M* )  2010  40 AVERAGE  AIR TEMPERATURE  C O  30 20 10  « ; T* 19 JUNE  FIGURE 1.1  24  29  4  9  1  4 JULY  19  24  29  4  9 14 AUGUST  19  D a l l y values of p r e c i p i t a t i o n , solar i r r a d i a n c e , average a l r temperature, and average r e l a t i v e humidity during the 1981 growing season.  -18-  "I  [  1  L....._""l  AU6..6.  Q. UJ O  0.4  r  .1 rl  0.2  JE  T  1  AUG 14  r  jJ~  -  J  [  JUL 27  j APR 24 (field  r  capacity )  LJ _J  fe 0.6 DC CL O CO  0.8 _L  FIGURE 1.2  _L  5  10  15  GRAVIMETRIC  SOIL WATER CONTENT  20 (%)  P r o f i l e s of gravimetric s o i l water content at selected times during the summer of 1981.  25  '79'  During the winter of 1981-82, seedling mortality did not change dramatically i n the Douglas-fir treatments. Western hemlock and P a c i f i c s i l v e r f i r survival were reduced by about 10%. The summer of 1982 was not as intensely hot and dry as i n 1981.  A  summary of the weather during the 1982 growing season i s given i n Figure 1.3.  Cool overcast conditions occurred through most of the summer with the  exception of a week i n June and two weeks i n mid-August.  In June, the s o i l  water content was s t i l l high enough to buffer the effects of the warm clear weather.  Less than 30 mm of rain f e l l during June, July and August, so the  clear hot weather i n August caused more seedling stress because of the dry soil profile.  Even though considerably more r a i n f e l l i n 1981 than 1982  s o i l water contents were similar i n August of both years because of the r a i n f a l l d i s t r i b u t i o n , the low water holding capacity of the s o i l , and the different d i s t r i b u t i o n of potential evapotranspiration rates i n each summer.  Measured root zone water potentials were also similar i n August of  both years. The duration of s t r e s s f u l conditions was much shorter i n 1982 than i n 1981. the  Even though s o i l water contents were similar, the overcast skies i n  summer of 1982 resulted i n very few days when surface s o i l temperatures  reached l e t h a l l y high l e v e l s .  The higher seedling survival rate i n 1982  r e f l e c t s these milder conditions.  The fact that s o i l water contents and  potentials i n the seedling root zone were very similar i n August of each year, while at the same time the intensity and duration of high surface temperatures was higher i n 1981 when s i g n i f i c a n t l y more mortality was observed, suggests that thermal stress was the major cause of seedling mortality for shade tolerant western hemlock and P a c i f i c  silver  -20-  FIGURE 1.3  Daily values of p r e c i p i t a t i o n , solar irradiance, average a i r temperature, and average r e l a t i v e humidity during the 1982 growing season.  f i r seedlings. The results of the 1982 experiment dence of t h i s .  (Table 1.4) provide further e v i -  The n u l l hypothesis of the 1982 experiment was that there  would be no difference i n survival between treatments.  There was  evidence  to disprove the n u l l hypothesis for western hemlock, i . e . a s i g n i f i c a n t improvement i n survival was obtained by shading, i r r i g a t i n g , and combined shading and i r r i g a t i n g .  If water stress had been the sole cause of mortal-  i t y , then the unwatered treatments should have shown lower survival; however, the shade treated seedlings survived as well as the seedlings i n the 2 treatments which received water.  If thermal stress had been the sole  cause of mortality, then unshaded treatments should have shown lower survival.  The unshaded-unwatered treatment did show s i g n i f i c a n t l y lower  s u r v i v a l , but the unshaded-watered treatment did not; however, surface s o i l temperatures did not reach l e t h a l levels i n the i r r i g a t e d treatment of  the wetter s o i l .  because  This suggests that thermal stress was the major cause  of seedling mortality.  D.  The E f f e c t s of Shade and I r r i g a t i o n on S o i l  Temperatures  I r r i g a t i o n , shade, and combined i r r i g a t i o n and shade treatments were shown to s i g n i f i c a n t l y increase seedling survival at the experimental site.  Shade cards were also shown to be more e f f e c t i v e as the severity of  drought increased. In this section the effects of each of these treatments on s o i l temperatures i s i l l u s t r a t e d using measurements collected on August 7,  1982.  This p a r t i c u l a r day was representative of summer drought conditions at the experimental s i t e .  The sky was clear, the a i r was warm, and the s o i l pro-  -22-  f i l e was very dry.  Similar conditions had been observed for the previous  week. The diurnal course of a i r temperature at the 1.0 m height, s o i l surface temperature (measured with a fine wire thermocouple), and s o i l temperatures at several depths below the untreated control treatment are shown in Figure 1.4.  Although a i r temperature only reached 26.5°C, the surface  temperature remained above 50°C f o r more than 7 hours. temperature of 64°C was recorded at 1500 P.S.T. maximum temperature was 32.5°C.  A maximum surface  At the 0.05 m depth, the  At the 0.15 m depth, the minimum and maxi-  mum s o i l temperatures were 15 and 20°C.  At the 0.3 and 0.5 m depths, the  temperatures remained constant at approximately 13 and 15°C respectively. These results show a strong contrast between temperatures at the surface and those below 0.15 m.  Deep s o i l temperatures were low for the o p t i -  mal growth of Douglas-fir seedlings (Brix 1967; Heninger and White 1974); however, at the same, time surface temperatures were l e t h a l l y high (Maquire 1955; Silen 1960; Baker 1929; Cleary 1978). Operational procedures which would increase the heat flow from the surface into the p r o f i l e would be b e n e f i c i a l on this s i t e .  This would  reduce the daytime surface maximum temperature, increase the nighttime surface minimum temperature, and increase the mean s o i l temperature at depth.  Factors affecting the s o i l thermal regime are discussed i n more  d e t a i l i n Chapter 2. The diurnal course of a i r temperature and surface temperature for the untreated control, shade, i r r i g a t i o n , and shade and i r r i g a t i o n treatments are shown i n Figure 1.5.  The combined  combined  shade and i r r i g a t i o n  treatment gave the greatest surface temperature reduction; about 10°C less  -23-  - — i  1  1  AIR  TEMPERATURE  "I  ( 1.0 m height)  SURFACE  TIME FIGURE 1.4  ( hr* P.S.T. )  The diurnal course of s o i l temperature at several depths, surface s o i l temperature (measured with a fine (0.075 mm diameter) wire copper constantan thermo-couple) and a i r temperature at the 1 m height on August 7, 1982.  r  FIGURE 1.5  1  1  1  TIME  P.S.T.  1  1  >  r  (hrs)  The diurnal course of surface s o i l temperature adjacent to shaded, i r r i g a t e d , shaded and i r r i g a t e d , and control treatments on August 7, 1982. The course of a i r temperature at the 1 m height i s also shown.  -25-  than a i r temperature and up to 50°C less than the untreated control.  The  i r r i g a t i o n treatment reduced midday surface temperatures by up to 35°C and the shade card treatment by up to 40°C.  The dramatic r i s e i n temperature  after 1500 P.S.T. i n the shaded treatments was a consequence of the shade card's shadow moving to the east of the seedling.  This rapid  temperature  r i s e even at this late hour i n the day indicates the importance of the correct shade card placement for e f f e c t i v e shading. The corresponding seedling root zone temperatures for each of these treatments are shown i n Figure 1.6.  Very l i t t l e temperature difference was  measured between treatments at depths below 0.05 m, even though very large differences were found at the surface.  This i s explained by the fact that  the treatments affected a r e l a t i v e l y small proportion of the s o i l surface, and therefore the treatment effect was overwhelmed at depth. A l l of the higher temperature amelioration treatments were e f f e c t i v e i n keeping surface temperatures below the c r i t i c a l l i m i t of 54°C (Cleary 1978;  Silen 1960) during periods of hot weather when the s o i l was dry.  This i s reflected i n the seedling survival s t a t i s t i c s ; p a r t i c u l a r l y the 1981 t r i a l where s i g n i f i c a n t differences i n survival were seen between untreated and shaded seedlings.  1.5  CONCLUSIONS Very low survival of planted 1-0 containerized seedling stock was  observed on a steep south facing high elevation clearcut on Vancouver Island.  At the end of the summer with harsh weather conditions, survival  of untreated Douglas-fir seedlings was found to be s i g n i f i c a n t l y higher than that of untreated western hemlock and P a c i f i c s i l v e r f i r seedlings.  -Z6-  1—1  1  1  \  1  1  T~7]~~ \  N  / / / V  SURFACE  / '  E o  ;  CL UJ Q  1  1 —  *S  TREATMENT  0  - 7 - 9  i  1  1  "  TEMPERATURE  H O U R S A B O V E 5 0 °C 1  i  i  1713  " 1  NO SHADE - WATER  *~~yr7 ^ r  i*r if  \ / /  SURFACE  // '* *  HOURS ABOVE 50°C - 0 HOURS ABOVE 30°C - 3 i I  1— V.  1  \ / 1 1 1*  O CC CL  '  I  |  I  l  i  1  1  ^  SX  TEMPERATURE  1  S H A D E  1  '  " ATER W  M  SURFACE  /I 1  CO —1  • 1  •  _ i  1  1  SHADE - NO WATER  /£'  SURFACE  I II  TEMPERATURE  HOURS ABOVE 5 0 ° C - 0 HOURS A B O V E 30°C 6  \ 1/ ** **  •  TEMPERATURE  HOURS ABOVE 50°C - 0 HOURS ABOVE 30°C - 1 1  \ 5 -  // * **  1  W  t  1  J-  1  SOIL TEMPERATURE FIGURE 1.6  1  HOURS ABOVE 30°C .  1  20  1  < i *»#/ if •i i /  1 -  1  (°C)  The d i u r n a l course of s o i l temperature i n the seedling root zone adjacent to shaded, i r r i g a t e d , shaded and i r r i gated, and c o n t r o l treatments on August 7, 1982. Numbers above the p r o f i l e r e f e r to the time of day ( P . S . T . ) .  -27-  T h i s was a l s o t r u e f o r these s e e d l i n g s a t the end of the second growing season. ing  For s e e d l i n g s  season, u n t r e a t e d  significantly  planted  D o u g l a s - f i r and P a c i f i c  b e t t e r than the u n t r e a t e d  the harsh summer, shade cards hemlock and P a c i f i c  of the m i l d e r  significantly  survival marginally.  The periods  majority  During  s u r v i v a l of western the second growing  increased  the s u r v i v a l  self-shade  t r i c k l e i r r i g a t i o n and  significantly  of seedling m o r t a l i t y occurred  increased  west-  d u r i n g hot droughty  i n the summer when r o o t zone s o i l was very  collar  survived  survival.  concluded t h a t thermal s t r e s s due t o very ling  grow-  F o r s e e d l i n g s p l a n t e d a t the  growing season, shade c a r d s ,  the combination o f these two treatments only ern hemlock  increased  During the harsh summer, an i n c l i n e d  treatment only i n c r e a s e d beginning  s i l v e r f i r seedlings  s i l v e r f i r s e e d l i n g s ; however, a f t e r  species.  of the m i l d e r  western hemlock s e e d l i n g s .  significantly  season f o l l o w i n g p l a n t i n g , shade cards of a l l t h r e e  a t the beginning  milder  high  dry.  I t was t e n t a t i v e l y  s u r f a c e s o i l and r o o t  temperatures was the major cause of m o r t a l i t y .  seed-  -£8~  CHAPTER 2  FACTORS AFFECTING SURFACE SOIL TEMPERATURE IN A FOREST CLEARC UT  FACTORS AFFECTING SURFACE SOIL TEMPERATURE IN A FOREST CLEARCUT  2.1  INTRODUCTION In Chapter 1 c o n i f e r o u s  s e e d l i n g m o r t a l i t y on a h i g h e l e v a t i o n ,  f a c i n g f o r e s t c l e a r c u t on Vancouver I s l a n d , B.C. with high The the h i g h ity  and  surface  soil  temperatures d u r i n g  periods  o b j e c t i v e of t h i s chapter i s to d i s c u s s surface  soil  determine the ameliorative weather and  the  soil  temperature w i l l  treatments t h a t can  surface  improve s e e d l i n g  s u r v i v a l under harsh  s e n s i b l e and  d e s c r i b i n g the combined r a d i a t i o n  are o u t l i n e d .  The  aerodynamic approach to  l a t e n t heat f l u x e s i s a l s o d i s c u s s e d .  This  theory  f o l l o w i n g chapter where i t i s used i n c o n j u n c t i o n  with a  from e a s i l y measured m e t e o r o l o g i c a l  subsurface s o i l  THEORY  A.  The Energy and Radiation Balance of a Forest Clearcut energy balance of a s o i l R n  H,  LE and  G  Q  are  n  surface  - H - LE - G  the net  temperatures  data.  2.2  where R ,  temperature i n  a i d i n the development of  s o i l heat flow model to e s t i m a t e s u r f a c e and  The  mortal-  conditions.  energy balances of a s u r f a c e  i s a p p l i e d i n the  f a c t o r s which cause  A b e t t e r u n d e r s t a n d i n g of the p h y s i c a l p r o c e s s e s which  To meet t h i s o b j e c t i v e , the p h y s i c s  estimating  associated  temperatures t h a t u l t i m a t e l y l e a d to s e e d l i n g  surface  soil  found to be  of summer drought.  to c a l c u l a t e the e f f e c t of these f a c t o r s on  a forest clearcut.  and  was  south  can Q  be w r i t t e n as =  follows:  0  radiation, sensible, latent,  (1) and  -30-  s o i l heat f l u x d e n s i t i e s r e s p e c t i v e l y a t the s u r f a c e . flux density  can be expressed as the sum of the net s o l a r and long wave  irradiances  as f o l l o w s : R  where R  S  (albedo),  N  =  R g ( l - a) + e ( R S  L  - aT! )  (2)  k s  i s the s o l a r i r r a d i a n c e , a i s the s o l a r r e f l e c t i o n c o e f f i c i e n t e  s  a i s the  i s the long wave e m i s s i v i t y o f the s u r f a c e ,  Stephan-Boltzmann c o n s t a n t , T RL  The net r a d i a t i o n  s  i s the long wave i r r a d i a n c e  i s the s o i l  surface  temperature i n K, and  from the sky. R L can be e s t i m a t e d as  follows: RL  where T  a  =  £a  a  T  a  (3)  k  i s the a i r temperature i n K a t s c r e e n h e i g h t and e  e f f e c t i v e e m i s s i v i t y of the s k y . The e m i s s i v i t y wavelength; however, i t can be t r e a t e d value f o r f a i r l y is  strongly  turn i s well (e ) a c  dependent upon the vapour d e n s i t y c o r r e l a t e d w i t h a i r temperature.  i s c l o s e l y approximated as f o l l o w s  where T  a  The  a c  =  i s a f u n c t i o n of  Under c l e a r s k i e s , e  a  of the atmosphere which i n The c l e a r sky e m i s s i v i t y  (Idso and J a c k s o n 1969):  1 - 0.261 exp (-7.77 x 10~^ T  2 a  )  (4)  i s i n °C.  atmospheric e m i s s i v i t y n o r m a l l y ranges from 0.75 to 0.85 under  clear skies.  Under cloudy s k i e s  the atmospheric e m i s s i v i t y  cause the e m i s s i v i t y of c l o u d s i s n e a r l y on  i s the  as a constant equal to some average  l a r g e wavebands (Campbell 1977).  e  a  1.0.  i s h i g h e r be-  The atmospheric  emissivity  cloudy days can be estimated by adding the energy emitted by the c l e a r  portions  of the sky t o the energy e m i t t e d by the c l o u d s as f o l l o w s  -31-  (Campbell 1977): e  a  "  e  ac  +  C  c  1  ~ a c " 4AT/T )  (5)  e  a  where C i s the f r a c t i o n of sky covered by cloud, AT i s the temperature difference between screen height and the cloud base a i r temperature and T  a  i s i n K.  From inspection of (5) i t can be seen that estimates of e  a  under cloudy conditions are quite i n s e n s i t i v e (±10%) to errors i n the e s t i mation of C and AT. e  Note when C = 0 ( i e . clear skies), (5) reduces to e = a  . -  ac  By combining equations (1) and (2) as follows: R (l-ct) + e s  s  (R - oTg^) - H - LE - G L  D  = 0  (6)  i t becomes evident that the absorbed solar and long wave radiation are d i s sipated at the s o i l surface as either emitted long wave radiation, sensible heat, latent heat, or s o i l heat. Upon rearrangement of (6) the s o i l surface temperature can be expressed  as follows: T  B.  s  =  [ R ( l - a) + e  s  s  R  L  - H - LE - G ] / e Q  s  a)  (7)  0 , 2 5  The S o i l Heat Flux The surface s o i l heat flux density can be found from Fourier's Law of  heat conduction as follows: G where k the  s  D  =  -k  s  (3T/3z)| = o  (8)  z  i s the apparent thermal conductivity of the s o i l (which includes  combined  effects of pure conduction and vapour d i s t i l l a t i o n  temperature gradients), z i s s o i l depth, and 3T/8z|  z  =  due to  o i s the s o i l  -32-  temperature  gradient at the surface.  The negative sign indicates that the  flux i s directed into the s o i l ( G > 0) when 3T/3z < 0 ( i . e .  temperature  0  i s decreasing with depth).  S o i l heat flow i s quite accurately described by  (8) because radiative and convective heat transfer are n e g l i g i b l e i n the soil.  Heat transfer due to vapour movement i s s i m i l a r l y i n s i g n i f i c a n t  because the thermal and isothermal vapour fluxes tend to cancel during the day and are very small at night (Kimball et a l . 1976). A quantitative method of calculating G  Q  from a numerical solution of  the heat flow equation w i l l be presented i n the following chapter; however, to a f i r s t approximation, daytime average G  Q  i s reasonably  by 15% of daytime average R , while during the nighttime G s  approximated  C.  approximately Q  i s well  by Rn (Novak 1981).  Aerodynamic Equations Describing the Sensible Heat Flux, Latent Heat Flux, and Momentum Flux i n the Turbulent Boundary Layer Sensible and latent heat transfer occur predominantly by convection i n  the turbulent boundary layer.  Therefore, they are closely related to the  momentum flux to the surface.  Equations describing the sensible heat,  latent heat, and momentum flux densities are respectively:  where K^, Ky and  H  =  -p  LE  =  -(p  T  =  p  a  c  a  K  H  (3T/3z)  c /y) K  a  K  p a  p a  M  v  (9)  (3e/3z)  (3u/3z)  (10) (ID  are the eddy d i f f u s i v i t i e s for heat, vapour and  momentum respectively, pa i s the density of a i r , Cp i s the s p e c i f i c heat a  of a i r , y i s the psychrometric constant, e i s the vapour pressure of water in a i r , u i s the windspeed, and z i s the height above the ground.  The sign  -33"  convention used above i s determined from the d i r e c t i o n of transfer. the day  During  the gradients are negative (lapse) and H and LE are directed away  from the surface (positive f l u x e s ) .  At night the gradients are p o s i t i v e  (inversion) and H and LE are directed toward the surface The momentum f l u x , also known as the shearing toward the  (negative f l u x e s ) .  stress, i s always directed  surface.  Using aerodynamic theory, the eddy d i f f u s i v i t i e s i n (9), (10) and  (11)  are defined as follows: K  H  =  k u* z/>  K  v  =  k u^ z/(|)  V  (13)  =  k u^ z/cJ)  (14)  %  (12)  H  M  where k i s von Karman's constant, u^  (=  (T/p )^" ) 5  a  i s the f r i c t i o n  v e l o c i t y (the tangential v e l o c i t y of an eddy i n the turbulent boundary layer (Thorn 1975)), and <j>, <[>y and cj) are the diabatic influence H  M  functions for heat, vapour and momentum respectively.  The diabatic i n f l u -  ence functions account for the effects of atmospheric s t a b i l i t y on turbulent transfer of H, LE, and  the  T.  The atmosphere i s i n a state of neutral equilibrium when the v e r t i c a l temperature gradient i s exactly equal to the rate of cooling (due to expansion) of a parcel of a i r r i s i n g a d i a b a t i c a l l y (-0.01°C m ). -1  During the  day, when the s o i l surface i s warmed by solar radiation, the a i r temperature near the surface decreases rapidly with height.  The atmosphere i s  then i n a state of unstable equilibrium because parcels of a i r which are displaced upward continue to ascend since they are less dense than the surrounding a i r .  At night when the s o i l surface temperature i s lower than  the a i r temperature because of longwave radiative cooling, buoyancy forces  -34-  oppose any o r i g i n a l a i r displacement and the atmosphere in a state of stable equilibrium.  i s considered to be  In unstable conditions turbulence i s  enhanced by buoyancy and i n stable conditions i t i s suppressed.  Under  neutral conditions the diabatic influence functions equal 1, under unstable conditions they are less than 1, and under stable conditions they are greater than 1 (Webb 1965). The state of the atmosphere  can be described quantitatively by the  r e l a t i o n between the production of turbulent k i n e t i c energy by buoyancy and mechanical forces.  This i s described by the turbulent k i n e t i c energy  balance equation as follows: — p  "  e  a  — 3z  +  p  a  g H c T p a  u a  ^  where g i s the g r a v i t a t i o n a l acceleration, e i s the viscous d i s s i p a t i o n of k i n e t i c energy and T  a  i s i n degrees K.  The f i r s t term i n this equation  describes the mechanical production of turbulent k i n e t i c energy and the second term describes the buoyant production of turbulent k i n e t i c energy. The r a t i o of the buoyant to mechanical turbulent k i n e t i c energy production terms i s used as the basis of two different indices of atmospheric s t a b i l i t y ; the Richardson number, R i , (Richardson 1920) and the MoninObukhov scaling or s t a b i l i t y length, L, (Monin and Obukhov 1954). The gradient Richardson number, derived from (15), can be written i n d i f f e r e n t i a l form as follows: Ri  =  (g/T ) (3T/3z)  /  a  (3u/3z)  2  (16)  and i s calculated i n f i n i t e difference form as follows: Ri where T  a  =  ( g / T ) [ ( T - T ) / u ] [ln ( z / z ) ] z 2  a  s  a  0  (17)  i s the mean a i r temperature between the measurement heights z  -35-  and z . Q  The surface roughness length, z , i s the distance above the Q  ground where the wind speed extrapolates to zero when plotted against logarithm of height.  the  It can also be thought of as the size of the smallest  eddy i n the turbulent boundary layer (Thorn 1975). The Monin-Obukhov scaling length, also derived from (15), i s c a l c u l a t ed as follows: L  =  (p  a  c  p a  T  u* )  /  3  a  (k g H).  (18)  L i s defined i n terms of the sensible heat flux and i s therefore independent of height whereas Ri i s defined i n terms of the a i r temperature gradient and varies with height.  A more commonly used form of L i s the s t a b i l i -  ty index, £ (= z/L), which has been found to be related to the Richardson number as follows (Dyer and Hicks 1970): Ri  =  c (W^M )-  Measurements by Dyer and Hicks suggest that (<|>jl/(j)M  =  1.00  ± 0.06)  <$>y[  ) i s nearly unity  over a very wide range of atmospheric s t a b i l -  i t y , which implies that Ri « £. al,(1975).  (19)  2  This has also been suggested by Thorn et  Paulson (1970) and Businger (1972, 1975)  for unstable conditions (£ or Ri < 0).  suggest that Ri = c,  For stable conditions (£ or Ri > 0)  they suggest that the relationship between Ri and £ i s approximated as follows: C  -  Ri / (1 - 4.7  Ri).  (20)  Barker and Baxter (1975) used a similar method of c a l c u l a t i n g £ from R i . The indices of s t a b i l i t y are used to calculate the diabatic influence functions as outlined below.  It should be noted that these are semi-empir-  i c a l relationships which have been derived from the analysis of many large  - 3 6 -  data sets as well as dimensional analysis.  Physically quantitative  relationships describing turbulence do not yet e x i s t . For unstable atmospheric conditions the relationship between the diabatic influence functions and £ i s described as follows (Businger (1966); Pandolfo (1966); Dyer and Hicks (1970)): *H  *  *V  =  <f>M 2  =  (1-16C)"  :  0-5  Z < 0.  (21)  For stable atmospheric conditions, Businger et a l . (1971) found that: *H  ™  *V  =  *M  =  (°-  7 4  +  The a p p l i c a b i l i t y of (21) i s questionable tions where free convection  4  '  7  5) : C > 0.  (22)  under extremely unstable condi-  i s the dominant mode of heat transfer (Dyer and  Hicks 1970). Many d i f f e r e n t methods of incorporating diabatic influence  functions  i n the aerodynamic equations have appeared i n the l i t e r a t u r e i n the past few decades.  There are b a s i c a l l y three d i f f e r e n t approaches.  Two of these  are integrated p r o f i l e methods and the t h i r d i s based upon the r a t i o s of the eddy d i f f u s i v i t i e s . The e a r l i e s t of the integrated p r o f i l e methods was introduced by Lettau (1962) and was used by Fuchs et a l . (1968) to determine the sensible heat flux density.  Assuming s i m i l a r i t y between momentum and sensible heat  transfer (KJJ = K^) Fuchs et a l . calculated the sensible heat flux density as follows: H  =  P  a  Cp  a  k  2  u (T - T ) s  a  /  l n ( z / z ) - lfa) Q  2  (23)  where ty^, the integrated diabatic p r o f i l e correction factor for momentum transfer i s defined as:  -37-  *M Lettau numerically  z (4»M " 1) (z + z ) o  /  =  _ 1  d  0  integrated (24) •  using:  = (1 - 18 R i ) " 0  M  (24)  z-  (25)  2 5  as suggested by Panofsky, Blackadar, and McVehil (1960. The second integrated diabatic p r o f i l e approach, developed by Paulson (1970), d i f f e r s from Lettau's method i n that separate diabatic influence By substituting (12) into (9)  functions are used for heat and momentum. and integrating with respect to z from z  to z, the sensible heat flux  Q  density i s expressed as follows: H  =  Pa p a * ( s " a ) c  k  u  T  ln(z/z ) - ^ )  /  T  H  Q  (26)  where T|»  =  h  /  (1 - $ ) r  1  (27)  dr..  H  o s i m i l a r l y by combining (11) and (14) to  The f r i c t i o n v e l o c i t y i s obtained give: = and integrating from z  Q  u  (6u/6z)  k(z - z ) / <|> Q  M  (28)  to z to give: =  4  ku  /  In ( z / z ) - I } J 0  m  (29)  where C  Tpn  Integrated  =  J o  (1 - <|) ) V, M  d?  , (3 ) 0  solutions of (27) and (30) i n terms of x, were obtained for  unstable atmospheric conditions (using 21) by Paulson (1970) as follows:  f  H  =  2 In [(1+x ) 12] 2  : c < 0  (31)  -38-  ifa  =  2 £n [(1+x) / 2] + In [(1+x ) / 2] - 2 t a n x + tr/2 2  _1  : c. < 0 (32)  where  x = (1 - 1 6 ? ) ° ' .  (33)  25  For  stable atmospheric conditions Businger (1975) suggested that the  diabatic p r o f i l e correction factors for heat and momentum were i d e n t i c a l and could be calculated as follows: =  =  4  -  7  5 : 5 > 0.  (34)  There are several variations on each of the two integrated diabatic p r o f i l e methods quoted i n the l i t e r a t u r e .  Some of the more notable recent  work i n this f i e l d has been done by Berkowicz and Prahm (1982), Hammel et a l . (1982), Businger (1972, 1975), and Barker and Baxter (1982).  Berkowic  and Prahm used (31), (32), and (33) to estimate H from routine meteorologi cal data.  Hammel et a l . assumed  ipjj = ipjj and  model s o i l temperatures i n a wheat f i e l d .  used (31) and (33) to  It should be noted that both  Businger and Hammel et a l . have l e f t out the m u l t i p l i c a t i o n factor of 2 i n (31).  This was probably a typographical error since integration of (27)  shows that Paulson's work i s correct.  Berkowicz and Prahm also confirm  this. The third basic approach used to calculate sensible and latent heat fluxes was developed independently by Holbo (1971) and Thorn et a l . (1975). It was used by Leckie et a l . (1981) to estimate sensible heat flux from remotely sensed surface temperature.  This non-integrated p r o f i l e method i  based upon the s t a b i l i t y dependent m u l t i p l i c a t i v e correction factor, "F", which i s incorporated i n the basic sensible heat flux equation as follows: H  =  p  a  c  p a  k  2  u (T - T ) F s  a  /  l n (z/z ) . 2  Q  (35)  -39-  The correction factor F i s defined as follows: F  =  (ch C^M)  =  -1  H  (R  - G)  n  /  Q  (H'+LE')  (36)  where H' and LE' are the sensible and latent heat flux densities calculated assuming neutral atmospheric conditions. Using the same diabatic influence functions from equations (21) Ri * ?) and  (with  (22), Thorn (1975) has shown that:  and  F  =  (1-16  Ri) -  F  =  (1 - 5 R i )  0  7 5  2  :  Ri < 0  (37)  :  Ri > 0  (38)  whereas Holbo (1971) found that for a bare pumice s o i l i n Oregon, U.S.A.: F  =  (1-34  Ri)  0 - 5 5  :  Ri < 0  (39)  which Thorn et a l . (1975) claim d i f f e r s from their F values by less than 10% through the moderately unstable range (-0.6  < Ri < 0).  The integrated diabatic correction factors used i n a l l of the integrated methods are shown i n Figure 2.1.  The diabatic influence functions  used by Thorn et al.(1975) i n their non-integrated those used by Paulson (1970). function = Fuchs et a l . ' s ^  method are i d e n t i c a l to  It i s of interest to note that  Lettau's  function « Businger's and Hammel et a l . ' s  i^H function = Paulson's ipg function * 2. The latent heat flux i s defined s i m i l a r l y i n a l l 3 methods since ij'H = <()Y and Kg = Ky.  Using the integrated diabatic p r o f i l e method  of Paulson (1970) the latent heat flux density i s calculated as follows: LE = p where e  Q  a  c  p a  k u^ ( e  Q  - e) / [ l n ( z / z ) - i|>] y 0  v  - e i s the vapour pressure difference between the s o i l  and the atmosphere.  (40) surface  This equation i s easily applied when the s o i l i s moist  -40-  -2.5  -2.0 INDEX OF  -1.5 STABILITY  1.0  -0.5  ( Rl or £ )  FIGURE 2 . 1 Integrated diabatic p r o f i l e correction factors for heat and momentum versus 6 or Ri (indices of atmospheric s t a b i l i t y ) for unstable, neutral, and stable atmospheric conditions.  since e  Q  = e * ( T ) , i . e . the saturation vapour pressure at surface s  temperature, but i t becomes very d i f f i c u l t  to estimate e  Q  as the s o i l  surface dries.  2.3  RESULTS AND DISCUSSION  A.  Surface Energy Balance and Clearcut Microclimate The energy balance of a forest clearcut i s a function of the net solar  and longwave irradiances at the surface, the thermal and moisture c h a r a c t e r i s t i c s of the s o i l p r o f i l e , the atmospheric gradients of temperature, vapour pressure and wind, and the surface resistance to heat and vapour exchange.  a)  The Radiation Balance The net solar and longwave irradiances have the greatest single  influence upon the microclimate of a forest clearcut.  Solar radia-  tion i s highly variable i n i n t e n s i t y throughout the year, with tude, and with slope and aspect.  lati-  South-facing slopes receive s i g n i f i -  cantly more radiation than other aspects i n the northern hemisphere, p a r t i c u l a r l y when the solar elevation angle i s low i n the winter. Increasing slope angle increases the irradiance on south-facing slopes and decreases  i t on north-facing slopes.  The t o t a l daily irradiance  on east and west facing slopes i s similar, but the irradiance on the west facing slope lags the east-facing slope by a few hours.  Solar  radiation i s attentuated and backscattered by clouds and particulate matter i n the atmosphere.  In B r i t i s h Columbia there are t y p i c a l l y  -4Z-  fewer hours of sunshine i n June than July because June usually has more cloudy days (Hay 1979). The solar r e f l e c t i o n c o e f f i c i e n t (albedo) of the surface  influ-  ences the net solar irradiance and can quite easily be used to manipulate the energy balance.  Typical albedos of s o i l s range from 0.06 to  about 0.6 depending on the angle of the sun, the colour of the surface and the s o i l water content. The net longwave irradiance depends upon the temperatures and e m i s s i v i t i e s of the sky and the surface. usually range from 0.95 to 0.98,  Surface emissivities of s o i l  increasing with s o i l water content.  Atmospheric e m i s s i v i t i e s range from about 0.75 to 1.0 depending upon sky  cloud cover and the temperature and humidity of the lower atmos-  phere .  b)  The S o i l Thermal Regime The s o i l thermal regime i s dependent upon (1) the surface heat  flux density; G  Q  (the quantity of heat passing through a unit area  of surface s o i l per unit time), (2) the s o i l thermal conductivity; k  s  (the s o i l heat flux density per unit temperature gradient), and  (3) the volumetric heat capacity of the s o i l ; C  s  (the change i n heat  content of a unit s o i l volume required to cause a 1°C temperature change). The s o i l thermal conductivity and volumetric heat capacity have a very s i g n i f i c a n t influence upon s o i l p r o f i l e temperatures.  Both of  these thermal properties are functions of the s o i l bulk density and the  r e l a t i v e proportions of the a i r , organic, mineral, and water  -43-  fractions of the s o i l .  Table 2.1 shows the values of the thermal  properties of each of these materials.  Since the thermal conductivity  and volumetric heat capacity of water are so much higher than the other s o i l constituents, p a r t i c u l a r l y a i r , the s o i l water content has a large influence upon the s o i l thermal regime.  The values of both  thermal properties for organic s o i l are also much lower than for mineral s o i l and as a result the thermal response of a layered forest s o i l can be s i g n i f i c a n t l y modified by removing  the upper organic  horizon.  c)  The S o i l Water Regime S o i l water i n d i r e c t l y influences the surface energy balance i n 3  separate ways.  As mentioned  above, the values of the thermal conduct-  i v i t y and volumetric heat capacity of water are much higher than those of a i r . As the s o i l drains, pore water i s replaced by a i r and consequently the thermal conductivity and volumetric heat capacity of the s o i l decrease. The latent heat flux density i s also dependent upon the surface s o i l water content.  During hot dry periods when there i s a large  evaporative demand, the surface s o i l dries, forming a mulch which r e s t r i c t s vapour movement.  The water content of the s o i l p r o f i l e and  the s o i l hydraulic properties determine the supply of water to the surface. The t h i r d effect occurs when p r e c i p i t a t i o n having a temperature different from that of the s o i l percolates through the s o i l p r o f i l e . This i s an example of heat transport resulting from mass flow.  TABLE 2.1  Thermal properties of natural materials (after DeVries 1963)  Remarks  Water  4°C  still  Air  10°C  still  turbulent  Mg  m  -3  k  C  Density Material  MJ m  -3  C  _ 1  W m . C -1  0.57  1.00  4.18  0.0012  0.0012  0.025  0.0012  0.0012  «125  Sandy S o i l  dry  1.6  1.28  0.3  40% pore space  saturated  2.0  2.96  2.2  Clay S o i l  dry  1.6  1.42  0.25  saturated  2.0  3.10  1.58  dry  0.3  0.58  0.06  saturated  1.1  4.02  0.50  40% pore space Peat S o i l 80% pore space  - 1  -45-  d)  Boundary Layer Resistance to Heat, Water Vapour and Momentum Transfer Heat, water vapour, and momentum transfer are strongly influenced  by the resistance of the laminar and turbulent boundary layers of the clearcut surface.  Transfer i n the laminar boundary layer occurs by  molecular d i f f u s i o n .  Beyond the laminar boundary layer turbulence  causes convective transfer to become the dominant transfer process. The resistance of the turbulent boundary layer depends upon the aerodynamic roughness of the surface and the wind speed.  Extended  f l a t surfaces which are aerodynamically smooth develop r e l a t i v e l y thick boundary layers.  As the height of the roughness elements  increases eddies induced by turbulence near the surface break down the boundary layer.  The aerodynamic roughness lengths of various natural  surfaces are given i n Table 2.2. As the wind velocity increases, turbulent mixing increases heat and vapour transfer.  The a i r tempera-  ture, water vapour pressure and wind speed d i r e c t l y above the surface depend to a certain extent on the nature of the l o c a l surface; however, they are also a function of macroclimate and upstream weather processes.  B.  Relationship Between Surface Energy Balance Components and Surface Temperature Equation (7) can be used to study the relationship between the surface energy balance and surface temperature.  This i s i l l u s t r a t e d  i n Figure 2.2 where s o i l surface temperatures are shown for a range of  -46"  TABLE 2.2  Surface  Aerodynamic properties of natural surfaces (after Oke 1978).  Remarks  Water  S t i l l - open sea  Ice  smooth  Roughness length (m) 0.1 - 10.0 x 10 -5 _  0.1 x 10~  h  Snow  0.5 - 10.0 x 10"  Sanddesert  0.0003  Soils  0.001 - 0.01  Grass  0.02 - 0.1 m t a l l 0.25 - 1.0 m t a l l  0.003 - 0.01 0.04 - 0.10  Agricultural Crops  0.04 - 0.2  Orchards  0.5 - 1.0  Forests  1.0 - 6.0  Sources:  Sutton (1953),  Szeicz et a l (1969),  Kraus (1972).  -47-  H + L E + G  (Wm-*)  SOLAR IRRADIANCE (Wm* ) 1  FIGURE 2.2  S u r f a c e s o i l temperature ed f o r a range o f v a l u e s balance theory ( e q u a t i o n 0.97, e - 0.85, and T a  c  a  versus s o l a r i r r a d i a n c e c a l c u l a t o f H + LE + G from energy 2.7) assuming o = 0.18, e = 25°C. Q  8  solar irradiances from 0-1100 W m  when the combined sensible, latent and  s o i l heat flux densities increase from 0 to 600 W m  . A constant a i r  temperature of 25°C, albedo of 0.18, atmospheric emissivity of 0.85 and surface emissivity of 0.97 were used for these calculations.  The maximum  surface temperature attainable under these conditions i s represented by the upper l i n e , where only radiative heat loss occurs.  These calculations show  that the maximum surface temperatures attainable for solar irradiances of 200, 400, 600, 800, and 1100 W m ively.  -2  are 40, 60, 80, 95, and 110°C respect-  Even under conditions of an extremely dry, thick LFH surface and  calm atmospheric conditions surface temperatures would f a l l short of these values.  As the combined sensible, latent, and s o i l heat flux densities  increase by 100 W m  the surface temperature drops by approximately 10°C  under these conditions. be less than 300 W m  The combined t o t a l of H, LE, and G  Q  would rarely  during a mid-summer day, hence surface temperatures  should not exceed about 85°C.  C.  Comparison of Atmospheric S t a b i l i t y Correction Methods As indicated e a r l i e r , there are three methods of s t a b i l i t y correcting  aerodynamic  estimates of the sensible heat flux density.  The following  comparison of these methods was undertaken because the most recent, scient i f i c a l l y rigorous approach suggested by Paulson (1970) requires a complex and time-consuming  i t e r a t i v e procedure to estimate separate diabatic cor-  rection factors for heat and momentum.  Because the other integrated pro-  f i l e method, f i r s t derived by Lettau (1962) to describe momentum transfer and later used by Fuchs et a l . (1968) to describe heat transfer, assumes that the diabatic correction factors for heat and momentum are i d e n t i c a l i t i s much easier to implement.  The method developed by Thorn et a l . (1975)  -49-  and Holbo (1971) i s appealing because of i t s s i m p l i c i t y and ease of a p p l i cation; however, because i t i s a non-integrated  p r o f i l e method that was  developed for use at a single point some distance from the surface where gradients are smaller, i t s application near a surface subject to very unstable conditions was questioned. The results of a comparison of these methods applied under i d e n t i c a l environmental conditions are shown i n Figure 2.3. following variables remained constant: LE = 150 W m , -2  0.01  a = 0.18, e  s  = 0.97, e  T a  a  In each case, the  = 25°C, R  = 1000 W m ,  G +  -2  s  = 0.85, z = 1.0 m, and z = Q  m while the wind speed was allowed to drop from 10 to 0.5 m s  - 1  so  that performance under very unstable atmospheric conditions could be determined. The upper s o l i d l i n e i n Figure 2.3 shows the calculated surface temperature when s t a b i l i t y correction was not used ( i . e . (7) solved i t e r a t i v e l y using (26) with \|>JJ = 0) .  It i s evident from the large difference between  the corrected and uncorrected methods that s t a b i l i t y correction i s important, i . e . that free convection i s a very important under unstable conditions.  heat transfer process  The integrated p r o f i l e methods of Paulson,  Fuchs et a l , and Hammel et a l , a l l predict surface temperatures within at least a few degrees of each other over the whole range of wind speeds.  The  method outlined by Paulson was solved using (7), (26), (29), (31), (32), and (33).  Hammel et a l . ' s method s i m i l a r l y used (7), (26), (29) and (33);  however i t assumes that ipj^ = i|>g and that the without  the m u l t i p l i c a t i o n factor of 2.  function equals (31)  Fuchs et al.'s method u t i l i z e s (7),  (23), (24), and (25) which gives a \JJ function very similar to that used m  by Hammel et a l . (as seen i n Figure 2.1).  The non-integrated methods of  -50-  I  •  0  I  -  i  I  2  i  1  4  1  6  1  8  1  1—  10  WIND SPEED ( m s-' ) FIGURE 2.3  Surface s o i l temperature versus wind speed c a l c u l a t e d using energy balance and aerodynamic theory and the s t a b i l i t y c o r r e c t i o n methods of Paulson (1970), Thom et a l (1975), Holbo (1971), Fuchs et a l (1968), and Hammel et a l (1982). Values were c a l c u l a t e d assuming Rg = 1000 W m , G + LE - 150 W m " , T - 25°C, o - 0.18, e = 0.97, e «• 0.85 and z = 0.01 m. The methods of Fuchs et a l and Hammel et a l give i d e n t i c a l r e s u l t s (see t e x t ) . - 2  2  a  s  a  Q  -5/-  Thom et a l . and Holbo were calculated using (7), (35), (37), (38), and (39). The non-integrated p r o f i l e method used by Thom and Holbo breaks down under very unstable conditions near the surface.  At low wind speeds sur-  face temperatures were underestimated by up to 20°C.  This occurred because  the "F" factor increased too quickly for negative Richardson numbers, hence the sensible heat flux was overestimated and as a result the surface temperature was seriously underestimated.  This method, though well suited to  applictions away from the surface, does not accurately predict fluxes near the surfce where gradients are highly non-linear under unstable conditions. A further test of each method involved the comparison of measured and modelled surface temperatures at the experimental s i t e described i n Chapter 1.  Measurements of surface temperature were taken throughout August 9,  1981 with an infrared thermometer and a fine wire thermocouple.  Surface  temperatures were modelled using measured a i r temperature, solar i r r a d i ance, wind speed, and s o i l heat flux density.  The latent heat flux density  was estimated to be similar to that measured by Ballard et a l . (1977) on another dry south facing clearcut i n southwestern B.C. A 1:1 plot showing the correspondence between modelled and measured surface temperatures using each of the methods outlined above i s shown i n Figure 2.4.  The integrated p r o f i l e methods of Paulson, Fuchs et a l . , and  Hammel et a l . s a t i s f a c t o r i l y predicted surface temperatures throughout the range of atmospheric s t a b i l i t y conditions encountered on this day. The non-integrated diabatic p r o f i l e methods of Thom et a l . and Holbo successf u l l y predicted temperatures under stable atmospheric conditions, but as i n the previous test they underestimated surface temperatures by up to 15°C under unstable atmospheric conditions.  -SZ-  1  -70 o  LU  tr  1  1  PAULSON METHOD TH0M/H0LB0 METHODS FUCHS ET A L . METHOD  A o —a  HAMEL ET AL. METHOD ALL METHODS IDENTICAL-  •- • •  I<  rr _ _ LU 50  . 6  A  oo  LU  LU  o  2  rr 30  hi  00  LINE  o LU  _l -J LU Q  O  10  10  30  50  MEASURED SURFACE TEMPERATURE FIGURE 2.4  70 (°C )  C a l c u l a t e d v e r s u s measured s u r f a c e s o i l temperatures on August 9, 1981. C a l c u l a t e d v a l u e s were o b t a i n e d from measured Rg, T , u, G and LE/R = 0.05 u s i n g the s t a b i l i t y c o r r e c t i o n methods o f P a u l s o n (1970), Thorn e t a l (1975), Fuchs e t a l (1968) and Hammel e t a l (1982). S u r f a c e s o i l temperatures were measured w i t h an i n f r a - r e d thermometer and a f i n e wire thermocouple. a  Q  8  -S3'  Further micrometeorological studies are required to determine whether Paulson's method can be improved upon for sloping clearcuts.  The non-  integrated p r o f i l e method of Thorn et a l . and Holbo i s not recommended for this application.  D.  Examination of Micrometeorological Factors A f f e c t i n g Surface Temperature In this section, the effects of changes i n solar irradiance (energy  input), windspeed, aerodynamic roughness and a i r temperature on surface temperature are compared.  The method of c a l c u l a t i o n makes use of the  Paulson integrated diabatic p r o f i l e correction procedure tested and recommended i n the previous section.  As indicated e a r l i e r , this involves the  interative solution of (7), (26), (29), (31), (32) and (33).  The results  of the calculations are shown i n Figures 2.5, 2.6 and 2.7 where the effects of solar irradiance, a i r temperature, and surface roughness on surface temperature are shown for wind speeds ranging from 0.5 to 10 m s * . -  upper, middle and middle curves of Figures 2.5, 2.6 and 2.7 are the same as the Paulson curve of Figure 2.3.  The  respectively  Figure 2.5 shows that  higher wind speeds are required to reduce surface temperatures as the solar irradiance increases.  Figure 2.6 shows that the surface temperature  approaches a i r temperature as the wind speed increases.  This i s the r e s u l t  of decreasing turbulent boundary resistance with increasing wind speed.  To  a good approximation, i f a l l other variables are constant, a r i s e i n a i r temperature results i n the same r i s e i n surface temperature.  Figure 2.4  shows that higher wind speeds are required to reduce surface temperatures as the aerodynamic roughness of the surface decreases. The roughness of the  -54-  80 h  o  I.  60-  LU  WIND SPEED FIGURE 2.5  (m/s)  Surface s o i l temperature versus wind speed at at the 1 m height for various s o l a r i r r a d i a n c e s . Values were calcu l a t e d assuming G + LE = 150 W uT , T - 25°C, z 0.01 m, a - 0.18, e - 0.97, and e - 0.85. 0  a  8  a  Q  -55-  80h  a  i,  60  AIR TEMPERATURE = 45 °C  40  AIR TEMPERATURE = 25 ° C  20  AIR TEMPERATURE = 5 ° C  UJ  cc < CC  £ ui Ul  o cc  CO  WIND SPEED FIGURE 2.6  10  8  6 (m/s)  S u r f a c e s o i l temperature v e r s u s wind speed a t a i r tempe r a t u r e s o f 5, 25, and 45°C. V a l u e s were c a l c u l a t e d assuming R - 1000 W m" , G + LE » 150 W m~ , z 0.01 m, a » 0.18, e = 0.97, and e - 0.85. 2  2  s  Q  8  a  -S6'  WIND SPEED FIGURE 2.7  (m/s)  S u r f a c e s o i l temperatures versus wind speed a t a e r o dynamic roughness l e n g t h s o f 0.001, 0.01, and_0.1 m. V a l u e s were c a l c u l a t e d assuming Rg » 1000 W m , G LE - 150 W n T , T - 25°C, a = 0.18, e = 0.97, and e » 0.85. -  0  2  a  a  8  -57-  clearcut i n this study was estimated to be between 0.01 and 0.07 m.  At low  wind speeds, surface temperature rises sharply i n response to decreased roughness.  2.4  CONCLUSIONS By combining energy and radiation balance theory, an equation was  derived that enabled surface temperatures to be related to the major energy balance components i n a forest clearcut. mum  This expression showed that maxi-  surface temperatures occur when a l l absorbed radiation i s dissipated  r a d i a t i v e l y ( i . e . H + LE + G H, LE, or G  0  =0).  An increase of 10 Wm  i n the sum  was found to cause a surface temperature reduction of  Q  of approximately 1°C. A method was developed for calculating surface temperatures which combined  the above expression with aerodynamic theory of sensible heat  transfer through the turbulent boundary layer.  The aerodynamic theory used  to calculate H was shown to require integrated diabatic p r o f i l e correction functions to correct for the strongly unstable atmospheric conditions i n clearcuts.  While 3 different approaches tested gave satisfactory r e s u l t s ,  the integrated diabatic p r o f i l e method of Paulson (1970) i s recommended for estimating fluxes from near surface gradients.  The non-integrated diabatic  p r o f i l e methods of Thorn et a l . (1975) and Holbo (1971) were shown to overestimate the sensible heat flux density under very unstable atmospheric conditions. It should be noted that the analysis i n this section i s somewhat s i m p l i s t i c i n section D since a change i n one variable (such as R , s  H, LE, or G ) Q  i s generally accompanied  T, a  by changes i n one or more other  -56-  variables.  These interactive effects w i l l be accounted for i n the model  described i n Chapter 3.  -59-  CHAPTER 3  MODELLING SOIL TEMPERATURE IN A FOREST CLEARCUT  -60-  MODELLING SOIL TEMPERATURE IN A FOREST CLEARCUT  3.1  INTRODUCTION The importance of s o i l temperature on the survival and growth of coni-  ferous seedlings was discussed i n Chapter 1.  In Chapter 2 energy balance,  radiation balance and aerodynamic theory were used to describe the factors affecting surface s o i l temperatures.  In this chapter the above theory i s  used, i n conjunction with a numerical solution of the heat flow equation, to model the time course of s o i l p r o f i l e temperatures.  A physically based  model of the s o i l thermal regime has an important application i n determining whether high surface s o i l temperatures would occur under certain environmental conditions and to assess the effects of microclimate modification on seedling root zone temperatures.  The a v a i l a b i l i t y of a validated  soil  temperature model would also greatly reduce the need for expensive long term programs of s o i l temperature measurement. In the past two decades s o i l temperature modelling has received considerable attention i n the l i t e r a t u r e .  Fourier analysis was used by  Stearns (1967), Ballard (1972), and MacKinnon (1975) to describe s o i l temperatures i n desert, alpine, and a g r i c u l t u r a l s o i l s respectively.  A  finite  difference solution of the heat flow equation was used by Hanks et a l . (1971) and Wierenga and de Witt (1970) to describe s o i l temperatures i n a g r i c u l t u r a l s o i l s ; however, as i n the previous method i t was necessary to specify the time course of surface s o i l temperature.  Gupta et a l . (1981)  used the method described by Hanks et a l . and an empirical relationship between a i r and surface temperature to model temperatures of a s o i l covered with various thicknesses of corn mulch with reasonable success. More  -61-  recently, Hammel et a l . (1982) developed a physically based method of e s t i mating the surface s o i l temperature which required the f i n i t e solution of the heat and water flow equations for s o i l , and  difference aerodynamic  theory (see Chapter 2); however, their midday modelled and measured surface temperatures differed by as much as 20°C. The objectives of this chapter are (1) to develop a physically based model for calculating the time course of surface and subsurface temperatures from easily measured or estimated meteorological and s o i l data; (2) to validate the model with measured data from two other energy balance studies reported i n the l i t e r a t u r e as well as s o i l temperatures measured i n the  forest clearcut described i n Chapter 1; and (3) to use a s e n s i t i v i t y  analysis approach to show how various environmental factors influence the s o i l thermal regime. 3.2  THEORY The p a r t i a l d i f f e r e n t i a l equation describing the flow of heat i n s o i l  can be written as follows: - 3G/3z  =  p  s  c  p s  3T/3t  ' (1)  where G i s the s o i l heat flux density, T i s the s o i l temperature, z i s the s o i l depth, t i s the time, p Cp  S  s  i s the s o i l density ( s o l i d plus water), and  i s the s p e c i f i c heat of s o i l .  The equation i s derived from the heat  balance of a unit volume of s o i l where the change i n heat storage ( p Cp  S  s  3T/3t) i s equated with the divergence of the heat flux density  (-3G/3z). The s o i l heat flux density i s described by Fourier's Law as follows: G where k  s  =  -k  s  3T/3z  i s the apparent thermal conductivity of the s o i l .  (2) The negative  -61'  s i g n i m p l i e s t h a t the heat f l u x i s p o s i t i v e  (downward) when the  temperature  g r a d i e n t i s n e g a t i v e i f z i n c r e a s e s w i t h depth. Combining e q u a t i o n s (1) and C where C  s  (= p  Philip  s  s  (2) g i v e s :  3T/3t  =  3(k  s  3T/3z)/3z  C p ) i s the v o l u m e t r i c heat c a p a c i t y of the S  the e f f e c t s of vapour movement on temperature;  was  not used because  l a r l y when the s o i l  relatively i s dry.  little  however, h i s approach  heat flows by t h i s p r o c e s s , p a r t i c u -  K i m b a l l et a l . (1972) suggest that heat  to vapour movement i s i n s i g n i f i c a n t because  vapour  soil.  (1957) suggested a more g e n e r a l heat f l o w e q u a t i o n i n c o r p o r a t -  ing  due  (3)  transfer  the thermal and i s o t h e r m a l  f l u x e s r o u g h l y c a n c e l d u r i n g the day and are very s m a l l at n i g h t .  Non-conductive apparent  heat t r a n s f e r p r o c e s s e s can a l s o be accounted f o r by u s i n g the  thermal c o n d u c t i v i t y of the s o i l  r a t h e r than the t r u e  thermal  conductivity. The  implicit  (Richtmeyer 1957)  f i n i t e difference  s o l u t i o n of the heat d i f f u s i o n e q u a t i o n  i n v o l v e s the d i s c r e t i z a t i o n of the ( z , t ) s o i l  profile  plane i n t o a r e c t a n g u l a r g r i d of p o i n t s j = 1,2,3... L along the z a x i s n = 1,2,3... a l o n g the t a x i s .  and  The d i s t a n c e between the v e r t i c a l nodes i s  Az and between the time steps i s A t .  The  temperature  at any g i v e n g r i d  point i s T j . n  In  f i n i t e d i f f e r e n c e form, (3) becomes:  n n-1 T - T  V/  1_  1  Az  At  -  k,  V/  2  n-1 n n-1 1 /T +T -T 2Az \ J+l i+1 J  (  n-1 n n-1 n - T T + T - T J J  n -T 1  (4)  -63-  for a l l i n t e r i o r nodes j = 2 to L - l . When the upper and lower boundary conditions T  n  and T , n  are s p e c i f i e d a system of L equations i n L  unknowns i s established.  These equations form a tridiagonal matrix which  can be solved quite e f f i c i e n t l y using Gaussian elimination. The derivation and method of solving (4) are outlined i n more d e t a i l i n Appendix A.  The advantages of using a f i n i t e difference approximation  rather than an a n a l y t i c a l solution of the heat flow equation are that ( i ) the s o i l thermal properties can vary with depth and i n time and  ( i i ) the  upper and lower boundary conditions can be a r b i t r a r y rather than periodic functions of time.  The major d i f f i c u l t y with this type of model involves  the s p e c i f i c a t i o n of the upper and lower boundary conditions. boundary condition i s quite easily simulated  by modelling  The lower  the s o i l p r o f i l e  to a depth where the temperature does not change appreciably with time. Over a period of a few weeks this depth i s about 1 m.  On an annual basis a  p r o f i l e depth of a few meters would be required. The surface boundary condition i s much more d i f f i c u l t to estimate because i t depends upon many f a c t o r s .  Using the radiation and energy  balance theory outlined i n Chapter 2 the s o i l surface temperature can  be  calculated as follows: T  s  =  ( R ( l - a) + e s  s  R  L  - H - LE - G ) Q  /  e  s  a  0. 25  (5)  By ( i ) estimating H by the method described i n Chapter 2, ( i i ) approximating LE using a simple empirical method (described l a t e r ) , ( i i i ) calculating G (4), and  Q  and  using (2) i n f i n i t e difference form after solving  (iv) knowing the solar irradiance, a i r temperature, windspeed,  albedo, and atmospheric and surface e m i s s i v i t i e s , i t i s possible to  -64-  i t e r a t i v e l y calculate the surface temperature.  3.3  PROCEDURE The s o i l temperature simulations were run on an 8 b i t microcomputer  (Campbell S c i e n t i f i c Inc., Logan, Utah, U.S.A., Model C2000 with a 4 MHz Z80 CPU) having only 64 K bytes of random access memory.  The model was  written i n the basic programming language (Microsoft Mbasic) and run i n a compiled  form to increase the operating e f f i c i e n c y .  The microcomputer con-  tained just enough memory to simulate a s o i l p r o f i l e depth of 0.8 m. It took between 10 to 15 minutes of computer time to simulate one actual day using 24 one hour time steps.  The model could be easily transferred to a  mainframe computer where the execution time would increase by several orders of magnitude and memory limitations would not determine the depth of the s o i l p r o f i l e .  A.  The program source code i s l i s t e d i n Appendix D.  Assumptions Used i n the Model Three assumptions were made to simplify the solution of the surface  s o i l temperature.  The f i r s t was that the latent heat flux density could be  modelled empirically using: LE = c R  s  where c i s an experimentally determined c o e f f i c i e n t . supported (1975).  (6) This r e l a t i o n s h i p i s  by the measurements of LE i n a forest clearcut by Ballard et a l . During a period of summer drought the latent heat flux density was  approximately  i n phase with the solar irradiance during the day and was  n e g l i g i b l e at night. s o i l conditions.  Their data suggested that c « 0.05 under dry surface  Under moist surface s o i l conditions, the p o t e n t i a l  -65-  evaporation  formula of Jensen and Haise (1963) (LE = 0.025°C  3°C)R ) could be used to give a more accurate S  The empirical modelling present  (T + a  estimate of LE.  of LE i s a shortcoming i n the model at the  time; however, the model i s set up i n such a way that this problem  can be overcome.  A numerical s o i l water flow model running  could be used to calculate the vapour pressure face (Hammel et a l . 1982).  concurrently  of the a i r at the s o i l  sur-  This would allow an aerodynamic estimate of the  latent heat flux density to be made using the diabatic p r o f i l e correction factors calculated for sensible heat transfer and an additional i t e r a t i o n routine. The  second assumption was that the semi-empirical  diabatic influence  functions described i n Chapter 2 could be used to correct the aerodynamic estimate of the sensible heat flux density under non-neutral conditions on a sloping forest clearcut. semi-empirical, sets obtained  Although the diabatic influence functions are  relationships, which have been derived from many large data  over generally horizontal land surfaces and dimensional  analysis, they quite accurately correct for the effects of buoyancy on a sloping clearcut as shown i n Chapter 2. A third assumption was that equations (4) and (5) i n Chapter 2 could adequately estimate the atmospheric emissivity under cloudy reasonable estimates of e  a  skies.  Very  were obtained using this method and even a  r e l a t i v e l y large margin of error (eg. ± 0.05) would not seriously a f f e c t the surface energy balance.  B.  The Solution Method The solution involves the sequential estimation of s o i l temperatures  -66-  at time (t = 2) based upon the temperature d i s t r i b u t i o n i n the s o i l p r o f i l e at  time (t = 1) and estimates of the upper and lower boundary temperatures  at  time (t = 2 ) . The model i s i n i t i a t e d by representing the s o i l p r o f i l e as a series of  nodes i n a two-dimensional array; one dimension for depth and the other for time.  Each of the nodes at time (t = 1) i s given an i n i t i a l  temperature  (uniform through the p r o f i l e ) , thermal conductivity and volumetric heat capacity. are  Values of the thermal conductivity and volumetric heat capacity  also specified at time (t = 2).  Solar irradiance, a i r temperature,  wind speed, and % sky cloud cover are then read from a data f i l e for time (t = 2).  The latent heat flux density and atmospheric emissivity are semi-  empirically determined from these variables.  Rough estimates of the  surface temperature and s o i l heat flux density are also made. Using an i t e r a t i v e procedure, convergence between the aerodynamic and energy balance estimates of the sensible heat flux density i s found at a unique surface temperature.  This temperature i s used as the upper boundary  condition i n the s o i l temperature model at time (t = 2).  S o i l temperatures  between the upper and lower boundaries and the near surface s o i l heat flux density (using the f i n i t e difference form of (2)) are then calculated. the  If  estimated, ( i . e . the value used i n (5)), and calculated s o i l heat flux  densities are not equal, a new estimate of the s o i l heat flux density i s determined and the process i s repeated.  When the estimated and calculated  s o i l heat flux densities converge, the correct s o i l temperature d i s t r i b u tion i s attained.  Pertinent data i s then stored or printed, the s o i l temperature  distri-  bution at the next time step i s set equal to that at the present time step, new meteorological data i s read from a f i l e , and the process i s repeated.  The flow chart i n Figure 3.1 shows this solution method.  3.4  RESULTS AND DISCUSSION  A.  V a l i d a t i o n of the Model Extensive testing was undertaken to determine the v a l i d i t y of the  model.  The f i n i t e difference s o i l heat flow component was tested f i r s t by  comparison with the a n a l y t i c a l solution for heat flow through a i m slab having uniform thermal properties. condition, very good agreement (T =  Using a sinusoidal upper  thick boundary  ±0.1°C) was found between the two  solution methods over the range of thermal properties and surface  tempera-  ture amplitudes observed i n s o i l s . The f i n i t e difference s o i l heat flow component of the model was evaluated using Fourier analyzed s o i l temperatures reported by Stearns (1967).  These data were collected on July 11, 1964 i n a dry sandy Peruvian  desert s o i l having approximately constant thermal properties through the soil profile.  In this test, the temperature measured at the 1 mm depth by  Stearns was used as the upper boundary  i n the f i n i t e difference solution.  Modelled and measured temperatures at the 0.02, are  compared i n Figure 3.2.  0.1, and 0.2 m s o i l depths  Excellent agreement was found through the  whole diurnal period at each of the depths considering (1) measurement errors of ±0.5°C, (2) errors i n the estimated thermal properties and i n the s o i l homogeneity assumption, (3) errors i n the temperature probe location (±1 cm), and (4) errors associated with the Fourier transcription of Steam's data (±0.5°C).  ~6B-  INITIALIZE MODEL T ( r . l ) , ks(z,t), C ( i , t ) B  1 INPUT Ra, T and 2 u at t a  I CALCULATE p , c , L E , and c ESTIMATE G ( e e t ) , T a  p a  0  a  B  CALCULATE R , e <TT , Rn, and H (energy balance) <>  L  B  8  Estimate New G(est)  CALC T(z,2) and G  Set  0  T ( r , l ) - T(t,2)  END  FIGURE 3.1  Flow chart showing the s o l u t i o n method used to c a l c u l a t e surface and subsurface s o i l temperatures.  -69-  i  FIGURE 3.2  8 TIME  12 (hrs)  Comparison of modelled and measured (Stearns 1967) temperatures i n a dry Peruvian desert s o i l having uniform thermal properties w i t h depth on J u l y 11, 1964.  -70'  These two tests suggest that the s o i l heat flow component of the model i s very accurate.  It was assumed that the f i n i t e difference model worked  equally well for a s o i l p r o f i l e having thermal properties that varied with depth. The complete model, using the energy balance and aerodynamic theory to estimate the surface s o i l temperature, was tested with measured energy fluxes and s o i l temperatures reported by Novak (1981) for a large a g r i c u l t u r a l f i e l d at Agassiz, B.C.  Tests were conducted on two diurnal periods;  one on May 18, 1979, a clear day when the s o i l water content was near f i e l d capacity, and the other on July 21, 1979, another clear day after the s o i l had dried considerably near the surface.  Novak (personal  supplied hourly values of the solar irradiance,°air  communication)  temperature, wind  speed, and latent heat flux density i n addition to the s o i l thermal propert i e s , albedo, surface emissivity, and temperature at the bottom of the modelled s o i l p r o f i l e (0.8 m).  The surface roughness length was estimated  to be 0.005 m using values quoted i n the l i t e r a t u r e for a smooth a g r i c u l tural f i e l d .  The i n i t i a l p r o f i l e temperature d i s t r i b u t i o n was obtained by  running the model for one day, starting from a temperature p r o f i l e set equal to the lower boundary condition.  This procedure gave an i n i t i a l  temperature d i s t r i b u t i o n very similar to that measured by Novak. The comparison of modelled and measured energy fluxes and s o i l temperatures i s shown i n Figure 3.3 for the May 18 test.  In general, excellent  agreement was found between modelled and measured values except for a few hours near midday.  During this time the modelled sensible heat flux was  overestimated by up to 30 W m  and as a result the modelled surface  temperature was underestimated by about 2.5°C.  Modelled and measured  -71'  1000  m i l  • *  «. IS  M  «.  —  0  •  -200  •  i  I  J  U  L  I  J  L.  T — i — r  55  45  - i — i — i — - T — i — r  r  T  r-  SOIL DEPTH MODELLED MEASURED 0.0 m O 0.025 m —— •• • 0.1 m — — 0.2 m A 0.5 m 1.0 m  •  •  —  A?  ^J.I--U-i.'.'-'.^L  12 TIME  FIGURE 3.3  18  24  ( HRS. RS.T. )  Comparison of modelled and measured (Novak 1981) temperatures and surface energy fluxes on May 18, 1979 i n a wet s i l t loam s o i l at A g a s s i z , B . C .  -7Z'  temperatures at the 0.1, 0.2 and 0.5 m depths were i n near perfect agreement.  There was a discrepancy of about 3°C between the midday modelled and  measured temperatures at the 0.025 m depth.  This was p a r t i a l l y due to the  error i n the modelled surface temperature but i t also could have been due to an overestimate of the volumetric heat capacity i n this layer. The results of the July 21 test are shown i n Figure 3.4.  Excellent  agreement between modelled and measured energy fluxes was obtained throughout the entire day.  Modelled and measured surface temperature were also i n  close agreement during the daytime; however, modelled values were up to 3°C less than measured values i n the early morning.  Subsoil temperatures at  the 0.1, 0.2, and 0.5 m depths were a l l modelled to within less than 1°C through the whole day.  The temperature at the 0.025 m depth was underesti-  mated during the daytime and overestimated during the nighttime, suggesting as i n the May 18 test, that the volumetric heat capacity just below the surface was i n c o r r e c t l y estimated. The results shown i n Figures 3.3 and 3.4 suggest that the model i s capable of accurately simulating the time course of energy fluxes, surface and subsurface temperatures from very few data inputs.  The r e l a t i v e l y  small differences between modelled and measured s o i l temperatures are probably less than the s p a t i a l v a r i a b i l i t y that could be expected within a few meters of the temperature sensors. A l l of the variables used as model inputs except the latent heat flux density are routine climate station measurements.  The objective of the  next test was to determine i f the simple empirical LE model (6) could be used to model LE with s u f f i c i e n t accuracy to model the surface temperature.  Theory outlined i n Chapter 2 suggests that a change i n the latent  -73-  -200 55  o'  1  45  A A  /  /  \o<  A *  /  0  V  35 ec  /A  //  < BC  z ii!  25  L  //  _  o  f as.  9%  /-  ° pTTr _ _ _ _  CO  •  '  o  v N  • o*.-  SOIL DEPTH MODELLED MEASURED O ao m Ojoes m A  15  0 . 1m 0.2 m 0.9 m 1 . 0m  •  •  6  t  I  1  12  o 1  1  L  18  24  TIME (HR& R&T.) FIGURE  3.4  Comparison of modelled and measured (Novak 1981) temperatures and surface energy fluxes on J u l y 2 1 , 1979 In a s l i t loam s o i l at A g a s s i z , B . C . which had d r i e d near the surface.  -74-  heat flux density of 10 W m  -  change.  would cause roughly a 1°C surface temperature  It should be possible then, to use an empirical LE model to calcu-  late the surface temperature to within at least a few degrees. The test of the simple LE relationship for calculating s o i l  tempera-  tures from standard climate data was conducted on the clearcut s i t e described i n Chapter 1 for a period of the s i x hottest clear sunny days i n August 1981.  The s o i l p r o f i l e had dried considerably from the surface to a  depth of well over 0.3 m.  To account for this effect on the s o i l thermal  properties i n the model, the s o i l was divided into 4 layers.  The thermal  conductivity and volumetric heat capacity i n the 0-0.02, 0.02-0.05, 0.050.25 and 0.25-0.8 m layers were 0.3, 0.3, 0.5, 0.7 W m 2.0, 2.2 MJ m  - 3 o  C  - 1  _ l o  C  _ 1  and 1.2, 1.6,  respectively.  The temperature measured at the 0.8 m depth was 15.0°C and did not change s i g n i f i c a n t l y during the six days. boundary condition. the  surface.  A roughness length of 0.05 m was used to characterize  Measured  solar irradiance, a i r temperature, and wind speed  were used as inputs to the model. assumed to be 0.05 i n (6). 0.95  This was used as the lower  Since the s o i l surface was dry, c was  The surface emissivity and albedo used were  and 0.16 respectively. The results of this simulation are shown i n Figure 3.5.  measured temperatures at the surface and the 0.15,  Modelled and  0.3, and 0.5 m s o i l  depths generally showed agreement to within less than 2°C.  The increase i n  mean temperature with time was also followed well by the model. led  The model-  amplitude at the 0.15 m depth was underestimated by approximately 2°C;  however, the modelled mean daily temperature at this depth showed good agreement with the measured values.  The modelled amplitude at the 0.3 m  -75-  80 60 40  SOIL SURFACE r *  20 0 40  0.15 m DEPTH  u 30 e  "St. «"  20  LU  or z> »< or LU LU O  r  T  CO  10 0 40  1  -i  1  1  r  0.3 m DEPTH  30 20 10 0 40  T  1  1  1  0.5 m DEPTH  30 20 10 24  FIGURE 3.5  120 72 96 ( hrs RS.T.) 6 7 8 9 DATE (AUGUST, 1981) 48 TIME  144 10  Comparison of modelled (0) and measured (+) s o i l temperatures at 4 depths i n the Cameron Valley forest clear-cut experimental s i t e during the 6 hottest days i n 1981.  76-  depth were s l i g h t l y greater than the measured values, but again the mean daily temperature was well modelled.  The discrepancy between modelled and  measured amplitudes at these depths suggests that the estimated s o i l thermal properties were s l i g h t l y i n error. These results show that a simple function describing the latent heat flux density can be incorporated i n the model so that easily measured climate station data can be used to simulate s o i l temperatures under dry s o i l conditions. Another objective was to test the model for a longer simulation period and under variable weather conditions.  A sixteen day i n t e r v a l beginning on  July 2, 1981 was chosen for this test.  The weather for the week preceding  July 2 was warm and dry. the  From July 2 to the July 5 i t was hot and dry.  In  early morning of July 6 a f r o n t a l storm approached and cloudy cool  weather occurred u n t i l the July 11. this period.  Less than 10 mm of rain f e l l during  From July 12 u n t i l the July 17 a gradual warming trend  occurred and no rain f e l l .  As i n the last test, four layers were used to  represent the change i n thermal properties with s o i l depth.  The thermal  conductivity and volumetric heat capacity i n the modelled 0-0.02, 0.020.05, 0.05-0.25, and 0.25-0.8 m layers were 0.19, — 3  and 0.78, at  0.95,  1.4, 1.6 MJ m  0.23,  0.76,  0.96 W m C _1  _1  — 1  C  respectively.  The temperature measured  the 0.8 m depth was 8.0°C and the surface roughness, albedo, and surface  emissivity remained the same as i n the previous test.  Measured solar  irradiance, a i r temperature, and wind speed were used as model inputs.  The  latent heat flux density was estimated to be 5 percent of the solar i r r a d i ance, as i n the previous test.  The cloudy sky atmospheric emissivity  (equations (4) and (5), Chapter 2) function was used throughout the  -77-  i n t e r v a l and i t was assumed that the sky was completely cloudy from July 5 u n t i l July 11. The results of the test are shown i n Figure 3.6.  Modelled and measur-  ed s o i l temperatures at the surface, 0.1, and 0.2 m s o i l depths showed very good agreement.  Modelled surface temperature showed an accurate response  to the change i n the weather.  Modelled subsoil temperatures also r e f l e c t e d  the pattern of heat storage, heat loss, then heat storage during the consecutive clear hot, cloudy cool, then clear hot periods.  Modelled and  measured temperatures d i f f e r e d at the 0.1 and 0.2 m depths for the f i r s t two days.  This reflected inaccuracies i n the i n i t i a l temperature  t i o n i n the p r o f i l e rather than a deficiency i n the model.  distribu-  The results  shown i n Figure 3.6 further suggest that standard climate data can be used to model s o i l surface and subsurface temperatures with very good accuracy. On days when the latent heat flux density i s higher than i n the two cases modelled here either an increase i n the c values used i n (6) or the method of Jensen and Haise (1963) could be used to more accurately simulate the course of LE.  More research i s necessary to determine the effect of  s o i l moisture content on the relationship between LE and R  B.  s  i n clearcuts.  S e n s i t i v i t y Analysis of the Model The effects of the aerodynamic  roughness of the surface, the latent  heat flux density, the s o i l thermal properties, and the temperature at the bottom of the s o i l p r o f i l e on modelled s o i l p r o f i l e temperatures were determined using a s e n s i t i v i t y analysis.  Changes i n the s o i l thermal  regime due to microclimate modification and the necessary accuracy of e s t i mated s i t e dependent variables can be determined from these t e s t s .  -76-  60  -i  r  SURFACE  40 — 20 O  < *  o  W 30 Q_  _l  I  1  L.  -i  1  1  r  u  _i  i  0.2 m depth  6  TIME FIGURE 3.6  I  _1  8  (days)  10  12  i  14  1  1-  16  JULY 1981  Comparison of modelled (0) and measured (+) s o i l temperatures at 3 depths i n the Cameron Valley forest clear-cut experimental s i t e during 16 days of variable weather conditions i n July 1981.  -19-  A s i x day test period was used so that both cumulative and effects could be observed. standard  tests (and represented  s  The standard  conditions used i n a l l of the  by the (x) symbol i n Figures 3.7-3.10) were as f o l -  measured solar irradiance, a i r temperature, and wind speed from  August 5-10, e  The test was conducted by selecting a set of  conditions and allowing the variable of interest to change while  holding a l l others constant.  lows:  immediate  = 0.97  1981, c = 0.05,  z = 0.005 m, T(0.8 m) = 12.0°C, a = 0.16, Q  and the volumetric heat capacity and thermal conductivity of  the 0.0-0.02, 0.02-0.05, 0.05-0.8 m s o i l layers were estimated to be 1.0, 1.2, 1.4 MJ m ^" , and 0.19, 0.29, 0.67 W m -3  1  - 1 0  ^  respectively.  1  The effect of the aerodynamic roughness length on s o i l p r o f i l e temperatures i s shown i n Figures 3.7a-d.  Increasing z from 5 to 10 mm Q  reduced midday surface s o i l temperatures by about 5°C.  Increasing z  from 5 to 50 mm reduced midday surface temperatures by about 15°C.  Q  Very  l i t t l e effect was noticed at the surface at night because z influences Q  sensible and latent heat transfer and these fluxes are small at night. The aerodynamic roughness had a large influence on s o i l warming at depth. Decreasing the roughness length caused less of the available energy at the surface to be partitioned as sensible heat, and as a result more s o i l warming  occurred.  This effect was p a r t i c u l a r l y noticeable at the 0.3 and 0.5 m  depths as shown i n Figures 3.7c and  3.7d.  Hammel et a l . (1982) also noted the significance of surface roughness and suggested that i t was the most important factor involved i n the effect of t i l l a g e on s o i l temperature.  Management practices which increase the  aerodynamic roughness of the surface would be b e n e f i c i a l on forest c l e a r cuts prone to high surface s o i l temperatures.  -80-  The effect of the latent heat flux density upon s o i l p r o f i l e temperatures i s shown i n Figures 3.8a-d. c = 0.0, 0.05 rates.  and 0.2  S o i l temperatures were calculated using  i n (6) to give a wide range of s o i l evaporation  Surface s o i l temperatures decreased as the latent heat flux  increased.  Midday surface s o i l temperatures were 10°C lower when LE was _ o  increased from 0 to about 200 W m  under these conditions.  This appears  to disagree with the observation i n Chapter 2 that a 100 W m~ H + LE + G causes a 10°C increase i n T . s  decrease i n  Both observations are correct.  The simulation i n this chapter shows that while LE i s decreased by 200 W m , -2  H + G  has increased by about 100 W m~ . 2  Q  Novak (1981) observ-  ed this effect i n his energy balance-soil temperature study of a bare s o i l . Increasing the latent heat flux density also had the effect of reducing the t o t a l amount of available energy at the surface for s o i l heat transfer so that subsurface s o i l temperatures did not increase as r a p i d l y . The effect on s o i l p r o f i l e temperatures of doubling the volumetric heat capacity or doubling both the volumetric heat capacity and the thermal conductivity are shown i n Figures 3.9a-d.  Roughly, doubling C  and k  s  s  would be the result of increasing the volumetric water content from about 10 to 20% i n a coarse textured mineral s o i l .  These changes i n thermal  properties had s i g n i f i c a n t effects upon s o i l temperatures (up to 10°C). Figure 3.9 shows that increased moisture content at the surface decreases T  s  mainly because of i t s effect of increasing k  assuming LE i s constant.  s  rather than C  s  Increasing the thermal conductivity increased  both the rate of s o i l warming and the amplitude of the temperature wave at depth.  While increasing the volumetric heat capacity s l i g h t l y reduced  surface s o i l temperatures i t , markedly decreased temperature at depth.  At  T  o 60  M  W M  UJ CC  o o  •  o  o  oo • o o o  o  H 401  cc UJ  OL UJ  Ml  20 "8 o 0  o  -oj*  SOIL SURFACE  CO  7 DATE (days) FIGURE 3.7a  8 AUGUST 1981  10  9  The e f f e c t of the aerodynamic roughness ( z ) on modelled surface s o i l temperature. Values of z used were 5 mm (x) , 10 mm (+), and 50 mm (0) while a l l other variables were held constant. Q  Q  0  0.1 m DEPTH  60  O  LU  or  H 40 or LU  •"**•.«  «:::«  o.  50  o°°°o * « 0  ?  °°Sssgss°  O to  DATE  FIGURE 3.7b  7 (days)  8 AUGUST 1981  As f o r F i g u r e 3.7a but f o r temperatures a t the 0.1 m depth.  10  TEMPERATURE  > 01  o 0Q  c 1-1 a  OJ  » a* c rt  rti O H  § (D » rt C M  (0  CD  » rt  n> o  S  a. n  T3  -£0-  TEMPERATURE  5  6  7 DATE  FIGURE 3.8a  (days)  8 AUGUST 1981  The e f f e c t of the latent heat flux density (LE) modelled surface s o i l temperature. Values of LE were 0.0 (+), 0.05 Rg ( x ) , and 0.2 Rj, (0) while other variables were held constant.  7 DATE FIGURE 3.8b  (days)  8  9  AUGUST 1981  As f o r F i g u r e 3.8a but f o r temperatures a t the 0.1 m depth.  0  60  © UJ  cc 3  £ 40  cc UJ Q.  UJ  gjSSoOoooooogggggJgggoooO gggSoSoSSSooSsgsaigScS^  20  )88«  88i  o  CO  6  7 DATE  FIGURE 3.8c  (days)  8  9  AUGUST 1981  As f o r F i g u r e 3.8b but f o r temperatures a t the 0.3 m depth.  DATE FIGURE 3.8d  (days)  AUGUST 1981  As f o r F i g u r e 3.8b but f o r temperatures a t the 0.5 m depth.  0  • •  60  o 111  CC  H 40  I  cr  oa  UJ CL  £20 o  » «o **»*  **«»** *M  0oo  2 :8 !  s  SOIL SURFACE  7 DATE (days) FIGURE 3.9a  8 AUGUST 1981  The e f f e c t s of s o i l thermal properties on modeled surface s o i l temperature. Values of k and C i n the 0-0.02, 0.02-0.05 and 0.05-0.8 m s o i l layers were 1.0, 1.2 and 1.4 raj m ° C and 0.19 0.29 and 0.67 W m °C~ r e s p e c t i v e l y ( X ) , values of k are the same as those given above but values of C are twice those given above (+) and values of k and C are twice those given above ( 0 ) . A l l other v a r i a b l e s were held constant. g  - 3  s  - 1  - 1  l  8  s  s  g  o 60 O  UJ  cc 3  < 40 CC Ul CL 2 UJ  *****  8  8  ***** gooogg  20  o  CO  7 DATE ( d a y s )  FIGURE 3.9b  8 AUGUST 1981  9  As f o r F i g u r e 3.9a but f o r temperatures a t the 0.1 m depth.  10  0.3 m DEPTH  O60 O  LU  cr I N O  cr  «  O  CO  6  FIGURE 3.9c  DATE  7 (days)  8 AUGUST 1981  9  As f o r F i g u r e 3.9a but f o r temperatures a t the 0.3 m depth.  0.5 m DEPTH 0 O  60  UJ  cc ? 40  cc UJ  o.  I  20  o  CO  7 DATE ( d a y s ) FIGURE 3.9d  8 AUGUST 1981  9  As f o r F i g u r e 3.9a but f o r temperatures a t the 0.5 m depth.  10  o 60 o UJ CC 3 I<  • 40  cc Ul Q_ 2  UJ  20 SOIL SURFACE  o  i  CO  7 DATE ( d a y s ) VTCURE 10a FIGURE 3 3.10a  8 AUGUST 1981  9  The The effect of^the temperature ^ ^ ^ at the lower surface of d e p t h )  [ l o w e r  b  o  u  n  d  condition] on modelled surface s o i l temperature. Values of 10°C (+), 12°C (x), and 14°C (0) were used while a l l other variables were held constant.  10  o 60 O  tr I<  40  Q  tr  r  *  ' t o  0  «9  LU Q.  '  8 8  9 9  W  2  LU _ _  {_ 20  o  CO  6  7 DATE  FIGURE 3.10b  (days)  8 AUGUST 1981  As f o r F i g u r e 3.10a but f o r temperatures a t the 0.1 m depth,  *o  NO  DATE  FIGURE 3.10c  (days)  AUGUST 1981  As f o r F i g u r e 3.10a but f o r temperatures at the 0.3 m depth.  DATE  FIGURE 3.10d  (days)  AUGUST 1981  As f o r F i g u r e 3.10a but f o r temperatures a t the  -97-  depth there i s a compensating effect of increasing both k since increasing k increasing C  s  s  s  and  C, s  w i l l tend to increase the mean temperature and  w i l l tend to decrease i t .  Thermal properties obviously have a large influence upon s o i l thermal regimes.  On sites which are prone to high surface temperatures, management  practices, such as scalping off the surface organic layer or mixing the surface organic layer with the mineral s o i l below which tend to increase C  s  or k  s  w i l l reduce l e t h a l l y high daytime surface temperatures.  The f i n a l s e n s i t i v i t y test involved the effect of the lower boundary s o i l temperature upon temperatures at other depths i n the s o i l p r o f i l e . The test was conducted because of the d i f f i c u l t y i n obtaining measured deep s o i l temperatures, p a r t i c u l a r l y on remote sites where estimates may be i n error by a few degrees.  The test was conducted by using 10, 12, and 14°C  as the lower boundary temperatures and noting the effects at other depths. The results of the test are shown i n Figures 3.10a-d.  They show that the  modelled lower boundary temperature has no effect upon the surface ture and very l i t t l e effect upon root zone s o i l temperatures.  tempera-  Provided  that the lower boundary i s deep enough i n the s o i l , good estimates of s o i l temperature can be made even i f the lower boundary temperature i s s l i g h t l y in error.  3.5  CONCLUSIONS A f i n i t e difference numerical s o i l temperature model was developed to  calculate surface and subsurface temperatures i n s o i l s having thermal properties which vary with depth.  The model accurately predicted s o i l temp-  eratures i n a forest clearcut over extended periods of time using e a s i l y  '98'  measured climate station data.  Excellent agreement was also found using  data from two other energy balance studies reported i n the l i t e r a t u r e . A s e n s i t i v i t y analysis showed that the s o i l thermal properties, the latent heat flux density, and the aerodynamic  roughness of the surface and  the lower boundary s o i l temperature could be estimated with s u f f i c i e n t accuracy from relationships i n the l i t e r a t u r e when measurements of these variables were not available from the s i t e .  -99-  CONCLUDING REMARKS  Shade tolerant western hemlock and P a c i f i c s i l v e r f i r were not well adapted to the severe heat and water stress which occurred during the summer on a high elevation south facing slope.  Douglas-fir, a moderately  shade intolerant species, survived much better under these conditions. Modifying  the seedling microclimate  combination of these two  with shade cards, i r r i g a t i o n , or a  treatments s i g n i f i c a n t l y increased  survival.  Presumably, these treatments w i l l be more e f f e c t i v e during hot dry summers when there i s inadequate water i n the s o i l to supply plant requirements and suppress surface s o i l and seedling root c o l l a r temperatures.  The cost  effectiveness of microclimate modification should be c a r e f u l l y examined for different types of s i t e s .  The use of good stock, careful planting, and  attention to micro-site conditions may  also substantially increase s u r v i v a l  on stress prone s i t e s . Using the theory developed i n Chapter 2 i t i s possible to speculate about alternative microclimate modification techniques.  High surface  soil  and root c o l l a r temperatures can be decreased by either reducing the amount of solar radiation absorbed by the surface or d i s s i p a t i n g the absorbed energy more e f f i c i e n t l y as sensible, latent, or s o i l heat. Blocking the solar irradiance i s probably reducing the energy input. modifying  the most e f f e c t i v e means of  This can be accomplished to a lesser degree by  the solar r e f l e c t i v i t y of the surface.  Evaporative  cooling of  the surface can be substantially increased by i r r i g a t i n g dry s o i l ; however, this i s not a very p r a c t i c a l operational procedure.  The sensible heat flux  density i s strongly influenced by the aerodynamic roughness of the surface  wo-  and the wind speed. the  Designing cut blocks to maintain a i r flow and keeping  surface as rough as possible w i l l help reduce high surface  tures.  tempera-  The seedling thermal regime can also be improved by planting i n  mineral s o i l rather than s u r f i c i a l  organic material.  Scarifying increases  both the thermal conductivity and volumetric heat capacity of the s o i l and as a result surface temperature fluctuatons are reduced and s o i l tures at depth respond more rapidly to conditions at the surface.  temperaThe  procedures for increasing the sensible and s o i l heat flux densities are probably the most inexpensive and e a s i l y applied methods of improving the s o i l thermal regime; however, these measures may not be e f f e c t i v e enough on some s i t e s . The computer model described i n Chapter 3 can be used to predict the e f f e c t s of changing various s i t e c h a r a c t e r i s t i c s and weather conditions on s o i l temperatures.  Standard climate data was used to successfully model  s o i l temperatures for extended periods of time i n the clearcut.  With  minimal additional work the model could be used to model s o i l temperatures at remote s i t e s based upon existing climate and s o i l s data. The model could also be used to provide a s o i l temperature growing degree day index and could provide the means for a p r o b a b l i s t i c analysis of planting success.  Other p r a c t i c a l applications w i l l become apparent as the  model receives more use.  -707-  REFERENCES  Arnott, J.T. 1975. F i e l d performance of container-grown and bare root trees i n coastal B r i t i s h Columbia. Can. J . For. Res. 5: 186-194. Barker, E.H. and T.L. Baxter. 1975. A note on the computation of surface layer fluxes for use i n numerical modelling. J . Appl. Meteorol. 14: 620-622. Baker, F.S. 1929. Effect of excessively high temperatures on coniferous reproduction. J . For. 27: 949-975. Ballard, T.M. 1972. Subalpine s o i l temperature regimes i n southwestern B r i t i s h Columbia. A r c t i c and Alpine Research 4: 139-146. Ballard, T.M., T.A. Black and K.G. McNaughton. 1977. Summer energy balance and temperatures i n a forest clearcut i n south-western B r i t i s h Columbia. In: Energy, Water, and the Physical Environment of the S o i l . Report of the 6th B.C. S o i l Science Workshop, Richmond, B.C. A p r i l 20-21, 1977. Bercowicz, R. and L.P. Prahm. 1982. Sensible heat flux estimated from routine meteorological data by the resistance method. J . Appl. Meteorol. 21: 1845-1864. Brix, H. 1967. An analysis of dry matter production of Douglas-fir seedlings i n r e l a t i o n to temperature and l i g h t i n t e n s i t y . Can. Jour. Bot. 45: 2063-2072. Businger, J.A. 1975. Aerodynamics of Vegetated Surfaces. In: Heat and Mass Transfer i n the Biosphere, Part 1 - Transfer Processes i n the Plant Environment, (Devries, D.A. and N.H. Afgan, eds.), Scripta Book Co., John Wiley and Sons, New York. pp. 5-28. Businger, J.A., J.C. Wyngard, Y. Izume and E.F. Bradley. 1971. Fulx-prof i l e relationships i n the atmospheric surface layer. J . Atmos. S c i . 28: 181-189. Businger, J.A. 1966. Transfer of momentum and heat i n the planetary boundary layer. Proceedings of the Symposium on A r c t i c Heat Budget and Atmospheric C i r c u l a t i o n (The Rand Corporation), pp. 305-332. Campbell, G.S. 1977. An Introduction to Environmental Biophysics, Springer-Verlag, New York. Carslaw, H.S. and J.C. Jaeger. 1959. Conduction of Heat i n Solids, 2nd Ed., Oxford University Press, London. Cleary, B.D., R.D. Graves and R.K. Hermann. 1978. Regenerating Oregon's Forests. Oregon State University Extension Service, C o r v a l l i s , Oregon.  -702-  Cochran, P.H. 1969. Thermal properties and surface temperatures of seedbeds. P a c i f i c Northwest Forest and Range Expt. St., USDA Forest Service, Portland, Oregon. De V r i e s , D.A. 1963. Thermal Properties of S o i l s . In: Physics of Plant Environment, (Van Wijk, W.R., ed.), North Holland Publishing Co., Amsterdam, The Netherlands. Dobbs, R.C. and R.G. McMinn. 1977. Effects of scalping on s o i l temperature and growth of white spruce seedlings. In: Energy, Water and the Physical Environment of the S o i l . Report of the 6th B.C. S o i l Science Workshop, Richmond, B.C., A p r i l 20-21, 1977. Dyer, A.J. and B.B. Hicks. 1970. Flux-gradient relationships i n the constant flux layer. Quart. J.R. Met. Soc. 96: 715-721. Dyer, A.J. 1967. The turbulent transfer of heat and water vapour i n an unstable atmosphere. Quart. J . Roy. Meteor. Soc. 93: 501-508. Freeze, R.A. and J.A. Cherry. Englewood C l i f f s , New Jersey,  1979. Groundwater, Prentice-Hall, Inc., p. 544-545.  Fuchs, M. and C B . Tanner. 1968. C a l i b r a t i o n and f i e l d test of s o i l heat flux plates. S o i l S c i . Soc. Amer. Proc. 32: 326-328. Fuchs, M., C B . Tanner, G.W. T h u r t e l l and T.A. Black. 1968. Evaporation from Drying Surfaces Using the Combination Method. Agronomy Journal 61: 22-26. Geiger, R. 1965. Cambridge, Mass.  The Climate Near the Ground.  Harvard Univ. Press,  Gupta, S.C, J.K. Radke and W.E. Larson. 1981. Predicting temperatures of bare and residue covered s o i l s with and without a corn crop. S o i l S c i . Soc. Am. J. 45: 405-412. Hammel, J.E., R.I. Papendick and G.S. Campbell. 1981. Fallow T i l l a g e Effects on Evaporation and Seedzone Water Content i n a Dry Summer Climate. S o i l S c i . Soc. Am. J . 45: 1016-1022. Hanks, R.J., D.D. Austin and W.T. Ondrechen. 1971. S o i l temperature e s t i mation by a numerical method. S o i l S c i . Soc. Amer. Proc. 35: 665-667. Heninger, R.L. and D.P. White. 1974. Tree seedling growth at d i f f e r e n t s o i l temperatures. For. S c i . 20: 363-367. Hermann, R.K. 1964. J . For. 62: 98-101.  Paper mulch for reforestation i n Southwest Oregon.  Hermann, R.K. and W.W. Chilcote. 1965. Effect of seedbeds on germination and survival of Douglas-fir. For. Res. Lab. Oregon State University, C o r v a l l i s , Oregon, Research Paper 4, 28 p.  -703-  Hobbs, S.D., R.H. Byars, D.C. Henneman and CR. Frost. 1980. F i r s t year performance of 1-0 containerized Douglas-fir seedlings on droughty sites i n southwest Oregon. Forest Research Laboratory. Oregon State University, C o r v a l l i s , Oregon. Holbo, H.R. 1971. Development of a s t a b i l i t y correction for estimating Corrective Transfer. Technical Report No. 71-7, Department of Atmospheric Services, Oregon State University, C o r v a l l i s , Oregon. Idso,' S.B. and R.D. Jackson. 1969. Thermal radiation from the atmosphere. J . Geophys. Res. 74: 5397-5403. Jackson, R.D. and S.A. Taylor. 1965. Heat Transfer. In: Methods of S o i l Analysis, Part 1 - Physical and Mineralogical Properties Including S t a t i s t i c s of Measurement and Sampling. C.A. Black, ed. American Society of Agronomy, Inc., Madison, Wisconsin, pp. 349-360. Jensen, M.E. and H.R. Haise. 1963. Estimating evapotranspiration from solar r a d i a t i o n . J . I r r i g . Drain. Div., Am. Soc. Engr. 89: 15-41. Kimball, B.A., R.D. Jackson, R.J. Reginato, F.S. Nakayama and S.B. Idso. 1976. Comparison of f i e l d measured and calculated s o i l heat fluxes. S o i l S c i . Soc. Amer. J . 40: 18-25. Klinka, K., F.C Nuszdorfer and L. Skoda. 1979. Biogeoclimatic units of central and southern Vancouver Island. Prov. of B r i t i s h Columbia, Min. of Forests. 120 pp (1 map). Korstian, C F . and N.J. Fetherolf. 1921. Control of stem g i r d l e of spruce transplants caused by excessive heat. Phytopathology. 2: 485-490. Lavender, D.P. and W.S. Overton. 1972. Thermoperiods and s o i l temperatures as they affect growth and dormancy of Douglas-fir seedlings of d i f ferent geographic o r i g i n . For. Res. Lab. Sch. For., Oregon State Univers i t y , C o r v a l l i s , Oregon, Research Paper 13, 26 p. Leckie, D.C, T.A. Black and P.A. Murtha. 1981. Development and testing of a method of estimating sensible heat flux from natural surfaces using remotely sensed surface temperatures. Can. Jour, of Chem. Eng. 59: 189-194. Lettau, H.H. 1962. Notes on t h e o r e t i c a l models of p r o f i l e structure i n the diabatic influence layer. F i n a l Rep. Contract DA-36-039-SC-80282. USEPG, Ft. Huachuca, Arizona, and Dept. of Meteor., University of Wisconsin, Madison, Wisconson. pp. 195-226. L i s t , R.J. 1971. Smithsonian Meteorological Tables, 6th Ed., Washington, D.C, Smithsonian I n s t i t u t i o n Press. Livingston, N.J. and R.J. Stathers. 1982. Annual Report for Seedling Establishment Project #316-1. MacMillan Bloedel Ltd., Inhouse Report, 130 pp.  -/Of-  Livingston, N.J. and R.J. Stathers. 1983. Annual Report for Seedling Establishment Project #316-1. MacMillan Bloedel Ltd., Inhouse Report, 35 pp. Livingston, N.J., T.A. Black, D. Beames and B.G. Dunsworth. An instrument for measuring the average stomatal conductance of conifer seedlings. Can. J . For. Res. (submitted). Lumley, J.L. and H.A. Panofsky. 1964. The Structure of Atmospheric Turbulence. Interscience Monographs and Texts i n Physics and Astronomy, V o l . XII. Wiley, New York. MacKinnon, J.C. 1975. Modelling the thermal regime of a Nova Scotian s o i l under oats. Can. Agric. Eng. 18: 41-45. Maguire, W.P. 1955. Radiation, surface temperature, and seedling s u r v i val. For. S c i . 1(4): 277-285. Minore, D. 1970. Shade benefits Douglas-fir i n Southwestern Oregon cutover area. U.S.D.A. Forest Service Research Note. 2 p. P a c i f i c Northwest Forest and Range Exp. St., C o r v a l l i s , Oregon. Monin, A.S. and A.M. Obukhov. 1954. Dimensionless c h a r a c t e r i s t i c s of turbulence i n the surface layer. Akad. Nauk. SSSR. Geofiz. Inst., Tr. 163-187. Monteith, J.L. London. Munn, R.E.  1973.  1966.  P r i n c i p l e s of Environmental Physics.  Descriptive Micrometeorology.  24:  Edward Arnold,  Academic Press, New  York.  Novak, M.D. 1981. The moisture and thermal regimes of a bare s o i l i n the Lower Fraser Valley during spring. Ph.D. Thesis, University of B r i t i s h Columbia, Vancouver, B.C. Oke, T.R.  1978.  Boundary Layer Climates, Methuen and Co. Ltd., New  York.  Pandolfo, J.P. 1966. Wind and temperature p r o f i l e s for constant-flux boundry layers i n lapse conditions with a variable eddy conductivity to eddy v i s c o s i t y r a t i o . J . Atmos. S c i . 23: 495-502. Panofsky, H.A. measurements.  1963. Determination of stress from wind and temperature Quart. J.R. Meteor. Soc. 89: 85-94.  Panofsky, H., A.K. Blackadar and G.E. McVehil. 1960. P r o f i l e . Quart J . Roy. Meteor. Soc. 86: 390-395.  The Diabatic Wind  Paulson, C A . 1970. The mathematical representation of wind speed and temperature p r o f i l e s i n the unstable atmospheric surface layer. J . of Appl. Meteorol. 9: 857-861.  '105-  P h i l i p , J.R. 1957. Evaporation, and moisture and heat f i e l d s i n the P r i e s t l e y , C.H.B. 1959. Turbulent Transfer i n the Lower Atmosphere, Univ. Chicago Press, Chicago. Richardson, L.F. 1920. The supply of energy from and to atmospheric eddies. Proc. Roy. S o c , London. A97: 354-373. Richtmyer, R.D. 1957. Difference Methods f o r I n i t i a l Value Problems. Wiley-Interscience, New York. S i l e n , R.R. 1960. Lethal Surface temperatures and their interpretation for Douglas-fir. Ph.D. Thesis, Oregon State University, C o r v a l l i s , Oregon. Smith, F.H. and R.R. S i l e n . 1963. Anatomy of Heat Damaged Douglas-fir Seedlings. For. S c i . 9: 15-32. Stearns, C R . 1966. Micrometeorological Studies i n the Coastal Desert of Southern Peru, Ph.D. Thesis, Univ. of Wisconsin, Madison, Wisconsin. Steel, R.G. and J.H. T o r r i e . 1960. Principles and Procedures of S t a t i s t i c s With Special Reference to the B i o l o g i c a l Sciences. McGraw H i l l , New York. pp. 99-109. Thom, A.S. 1975. Momentum, Mass, and Heat Exchange of Plant Communities. In: "Vegetation and the Atmosphere", V o l . I, (Monteith, J.L., ed.), Academic Press, London, pp. 57-109. Thom, A.S. 1972. Momentum, mass and heat exchange of vegetation. J.R. Meteor. Soc. 98: 124  Quart.  Thom, A.S., J.B. Stewart, H.R. Olver and J.H.C. Gash. 1975. Comparison of aerodynamic and energy budget estimates of fluxes over a pine f o r e s t . Quart. J.R. Meteor. Soc. 101: 93-10.5. Van Wijk, W.R. and W.J. Derekson. S i n i s o i d a l Temperature Variation i n a Layered S o i l . In: Physics of Plant Environment, (Van Wijk, ed.), North Wierenga, P.J. and C.T. deWit. 1970. Simulation of heat transfer i n s o i l s . S o i l S c i . Soc. Amer. Proc. 34: 845-848.  -106-  APPENDIX  DERIVATION  OF T H E ONE DIMENSIONAL OF  THE SOIL  HEAT  A  FINITE  DIFFERENCE  FLOW E Q U A T I O N  SOLUTION  -107-  APPENDIX A Derivation of the One Dimensional F i n i t e Difference Solution Of the S o i l Heat Flow Equation  The p a r t i a l d i f f e r e n t i a l equation describing the d i f f u s i o n of heat i n a s o i l p r o f i l e can be written as follows: -3G 8z  =  p c s  (1)  3T_ 9t  p s  where T i s the temperature, t i s the time, z i s the depth, p s o i l bulk density, C p heat flux density.  S  s  i s the  i s the s p e c i f i c heat of s o i l , and G i s the s o i l  The equation i s derived from the heat balance of a unit  volume of s o i l where the change i n heat storage ( p c s  p s  6T/6t) i s equat-  ed with the divergence of the s o i l heat flux density (-6G/Sz). The s o i l heat flux density i s given by Fourier's G where  k  s  =  -k  s  Law:  <5T/6z  (2)  i s the s o i l thermal conductivity and 6T/6z i s the s o i l tempera-  ture gradient.  The negative sign implies a positive heat flux into the  s o i l when z i s measured p o s i t i v e l y downward from zero at the surface. Combining equations (1) and (2) gives:  C  where C  s  (= p c s  p s  s  9T 9t  =  9_ 3z  k  s  9T 3z  (3)  ) i s the volumetric heat capacity of the s o i l .  The i m p l i c i t f i n i t e difference solution of (3) (Richtmeyer  1957)  involves the d i s c r e t i z a t i o n of the (z,t) plane i n the s o i l p r o f i l e into a rectangular grid of points j = 1,2,3...L along the z axis and n = 1,2,3... along the t axis.  The distance between the v e r t i c a l nodes i s At.  The  solution at any given grid point i s T j . n  The p a r t i a l d i f f e r e n t i a l s i n equation (3) can be written i n difference form as: n- 1/2 /  n- / 1  1  3T/3t  =  C  S  s  n n-1 (T - T )/At j 3  2  j  and n- / 1  3  Tz  [k 3T/3z.] j  ^  s  [ ^  j  n- / [k ._ 1  s  +  i  /  2  2  ("/a.)  , n- / (3T/3z)_  2  1  1 / 2  j + 1 / 2  l  2  1/2  ]  where  3T 3z 3+ / 1  1_ Az  2  and n- 1/2 / 1  8T 3z  n-1 n T + T  n-1 n T + T  f n-1 n T + T  n-1 n T + T ( 3-1 3-1  j+1  j  1_ Az  j  j+1  i  "  5  The completed f i n i t e difference form of equation (3) becomes: n n-1 T - T 3 3 At  = !_  <  k  Az n- 1/2 / k.  'j+V  2  n-1  1  -  n-1 n n1 /T + T - T 2Az V j+1 j+1 j  1/T  n n-1 + T - T  n -T  C o l l e c t i n g similar terms and rearranging n- / 1  n- / 1  At  2  +  -T  S  T  j  1 +  /  :  - T  c  2Az'  +T  n- / 1  At  2  s,a/2 'j+V  V/  n-1 T  =  At  2  + C  2  (9)  Az'  j+l  2  2Az  2Az'  n - / At ..1,2  2  2Az'  1  2  1  ;  2  1  n - / At s . 1,2 1  n- /  n - / At n-1/2 / s . 1/2 + C J-V' 1  At  2  2Az'  2Az'  n-1  gives:  +  t  At  n-1 T 2Az  Z  2  Equation (9) can also be expressed as: -AT  n  +  j j+l  B  T  n  -  C  j J  n  T  J j "  1  =  D  3  (10)  where n-1/ /2 = k " At/2Az 1 ,2 D 1  A  j  1 /2  n-V C  j  B  D  j  =  A  (ID  2  j-V  =  k s j +  =  j  A  3  .1/2 n-1 / nT + C s  j j+l V j s  l/ +  2  A t  C  /  J  .2  (12)  2 A Z  +  C  s  3  (13)  (U)  -770-  When equations for the upper and lower boundary conditions ( j = 1 and j = L) are incorporated, the f i n i t e difference approximations i n equation (10) constitute a system of L l i n e a r equations i n L unknowns.The c o e f f i cients A,B,C  and D vary with the node ( j = 1, j = 2 to L -1, j = L ) , with  the time step (n = 1, n = 2, n > 2), and with changes i n the heat capacity and thermal conductivity of the p r o f i l e .  The variables on which A,B,C  and  D depend are the s o i l thermal properties, the node size i n the g r i d , and the temperatures at the upper and lower boundary  conditions.  The values of T j are calculated from the recurrence r e l a t i o n : n  T E  where  n j j D  and  j  =  E  =  A +  j J C  T  n  B  /  j  +  J+1  F  j  J-1  F  (15)  J  -  C  /  B  j j  E -  (16)  J-1 C  J  E  (17) J-1  The E and F c o e f f i c i e n t s are calculated from nodes j = 1 to j = L, and the T's are back-calculated from j = L to j = 1. The c o e f f i c i e n t s A-F at node j = 1 and A-F at node j = L are defined at each new time step as follows: A (1)  =  B (1)  0 =1  C (1)  =  0  D (1)  =  T ( t , 1)  E (1)  =  A (1) / B (1)  F (1)  =  D (1) / B (1)  ( l g )  A (L)  =  0  B (L)  =  1  C (L)  =  0  D (L)  =  T ( t , L)  E  (L)  =  A (L) / B (L)  F (L)  =  D (L) / B (L)  -7J2-  APPENDIX B  MODELLING CLEAR SKY SOLAR IRRADIANCE ON VARIABLE SLOPES AND ASPECTS  -773-  APPENDIX B. Modelling Clear Sky Solar Irradiance on Variable Slopes and Aspects  1.  INTRODUCTION Potential solar irradiance i s the maximum short wave irradiance poss  ble for any given observation s i t e and time.  It only occurs under clear  skies and varies with l a t i t u d e , a l t i t u d e , slope, aspect, and time of the day and year.  It can be calculated from the geometrical relationship  between the sun and earth and an empirical atmospheric radiation attenuation function as follows: ^t slope where S  t  s  i pe> 0  ^b slope>  a n c  ^b slope  =  ^ $d slope  a  r  e  t  $d slope  +  n  e  t  o  t  a  l  solar i r r a d i -  ance, direct and diffuse components respectively. For a horizontal surface, the direct component i s calculated as follows ( L i s t 1971): s  b horiz  =  s  p  c  o  s  (2)  z  where z i s the zenith angle of the sun. Sp, the direct solar beam irradiance can be empirically determined as follows ( L i s t 1967): S where Sp  D  p  =  a* S  (3)  p o  i s the e x t r a t e r r e s t r i a l direct solar beam irradiance (solar  constant = 1360 W m  ), a i s an atmospheric transmission c o e f f i c i e n t , and  i s the o p t i c a l airmass number. The atmospheric transmission c o e f f i c i e n t varies from 0.9 f o r a very clear atmosphere to 0.6 for a very hazy or smoggy atmosphere. value of (a) f o r clear days would be 0.84 (Campbell 1977).  A typical  The o p t i c a l  airmass number, the ratio of the slant path length of solar radiation to the zenith path length i n the atmosphere, for elevation angles greater than 10 degrees i s given as follows: m  =  (P/P ) / (cos z  (4)  Q  where (P/P ) , the atmospheric pressure at the observation s i t e divided by Q  sea l e v e l atmospheric pressure, corrects m for a l t i t u d e e f f e c t s . The zenith angle of the sun (z) can be found from: cos z  =  s i n c|>  sin e  +  cos §  cos e  cos 15 (t - t ) Q  (5)  where tj> i s the l a t i t u d e , e i s the solar declination, t i s the time of day, and t  0  i s the time of solar noon.  The horizontal direct short wave irradiance can be corrected for slope and aspect as follows (Hay 1979): Sb slope  where i ,  Sb  =  froriz cos l cos z  the angle of incidence of the sun's rays on the slope, i s found  from: cos i  =  cos s  cos z  +  sin s  sin z  cos (az-b)  where s i s the i n c l i n a t i o n angle of the slope, b i s the azimuth of the slope (aspect) and az, the azimuth of the sun, i s found from: cos az  =  s i n § cos z - s i n e cos <j> s i n z  ^  To a good approximation the clear sky diffuse radiation can be found from half the difference between the irradiance on a horizontal surface  -7/5-  below and above the atmosphere as follows ( L i s t 1971): S  d  =  0.5 S  (1 - a ) cos z m  p o  Sjj i s almost a constant f r a c t i o n of S  p o  (8)  over most of the range of  zenith angles (Campbell 1977), hence i t does not have to be corrected for slope and aspect, i . e . S^ slope " ^d* Reflected radiation might also increase the t o t a l solar irradiance by a few percent depending upon the sky view factor of the observation s i t e .  -us-  1200  900 E ui  < 600 Q <  tr cr <  o to 300 SOUTH 50* HORIZONTAL EAST 30* WEST SO* NORTH SO*  /  solar irradiance  200 JULIAN DAY  a t s o l a r noon on h o r i z o n t a l ,  and West 3 0 ° , and North 30° s l o p e s .  \ \  J  100  Calculated  \ *  s  L  300  south 3 0 ° , Ea  -117'  L i s t of Symbols  _o  S  t  s  iope  t o t a l solar irradiance on a sloping surface  =  W m  Sb slope  =  d i r e c t solar irradiance on a sloping surface  W m  b horiz  =  d i r e c t solar irradiance on a horizontal surface  W m  s  _  2  _  2  _2 S,j ^d slope  =  diffuse solar irradiance  W m  diffuse solar irradiance on a sloping surface  =  W m  _  2  _2 Sp  =  direct solar beam irradiance  Spo  =  solar constant (1360 ± 40 W m~ )  i  =  angle of incidence of sun s rays on a slope  degrees  z  =  zenith angle of sun  degrees  s  =  i n c l i n a t i o n of slope from horizontal  degrees  b  =  azimuth of slope (aspect)  degrees  az  =  azimuth of sun  degrees  p  =  l a t i t u d e of observation s i t e  degrees  e  =  declination of sun  degrees  P  =  atmospheric pressure at observation s i t e  Pa  =  atmospheric pressure at sea l e v e l  Pa  m  =  o p t i c a l airmass number  dimensionless  a  =  atmospheric transmission c o e f f i c i e n t  dimensionless  t  =  time of day  h  =  solar noon  h  P  t  Q  n  W m 2  1  W  m  -2  APPENDIX C  SOIL THERMAL PROPERTIES  APPENDIX C S o i l Thermal Properties  The  s o i l thermal properties used i n the simulations presented i n this  thesis were chosen on the basis of s o i l texture and s o i l water content. The volumetric  heat capacity was  calculated from the following equa-  t i o n (De Vries 1963): c  s  =  V^m  + o o + V^w x  c  +  x  a a  ^  c  where the subscripts m, o, w and a represent mineral s o i l , organic matter, water, and a i r , respectively, x i s the volume f r a c t i o n and C i s the volumetric heat capacity.  Volumetric  heat capacities of mineral s o i l s having  d i f f e r i n g bulk densities and water contents are shown i n Figure The  1.  thermal d i f f u s i v i t y of a mineral s o i l sample taken from the  Cameron Valley clearcut site was measured using the method described S t e r l i n g and Taylor (1965). tents of 0.0  and 0.2  by  Measurements taken at volumetric water con-  agreed very well with the d i f f u s i v i t i e s of sand quoted  by van Wijk and de Vries (1963) (Figure 2). The  thermal conductivity was  calculated from the volumetric  heat  capacity and the thermal d i f f u s i v i t y as follows: k  s  =  where k  s  and C  i s the volumetric heat capacity.  s  i s the thermal conductivity, K  K C S  s  S  (2)  i s the thermal d i f f u s i v i t y , The calculated values of thermal  conductivity are plotted against s o i l water content i n Figure 3.  -\io-  T  1  -i  r  S o i l v o l u m e t r i c heat c a p a c i t y v e r s u s s o i l v o l u m e t r i c water content mineral  soil.  for a  1.0  1  1  I  1  1  1  1  A  A  —  A  » 0.8 E  (0 I  O  —  o >-  0.6  u. t 0.4 o  O  o  A  o  O  O  0  O  O  —  i  < —  —  o  —  o  cc LU X  O  O  o —  CO Z>  A  0.2  1 . 1 1 1 02 0.1 SOIL VOLUMETRIC WATER CONTENT (mVm ) 1  1  1  3  S o i l thermal d i f f u s i v i t y versus volumetric water content of a sandy s o i l having a bulk density of 1.0 ( Q ) ,  1.3 ( O )» and 1.6 ( A ) Mg  m~ . 3  0  0.1  0.2  03  0.4  SOIL VOLUMETRIC WATER CONTENT (mVm5) S o i l thermal conductivity versus volumetric water content of a sandy -3 having a bulk density of 1.0 ( Q > »  1  ,  3  ( O )»  a  n  d  1.6 ( A ) Mg m  c>  APPENDIX D  COMPUTER PROGRAM USED TO CALCULATE SOIL TEMPERATURES AND SURFACE ENERGY FLUXES  SET UP MODEL INITIAL CONDITIONS AND CONSTANTS****************** OPEN"I",#1,"SOILTEST.DAT" • TO READ IN CLIMATE DATA OPEN"0",#2,"LONG.DAT" • TO OUTPUT SOILTEMP DATA DIM A(161),B(161),C(161),D(161),E(161),F(161),COND(161),RHOCP(161 TEMP(2,161),DEPTH(161) INPUT"MEAN DAILY SURFACE SOIL TEMPERATURE "jTBAR 270 INPUT"SURFACE ROUGHNESS LENGTH 20 (M) 0.001-lMM "jZO 280 INPUT"LE • X * RS INPUT X "; ALF 290 5.67E-08 STEPHAN BOLTZMANN CONSTANT SIGMA 310 .95 SURFACE EMISSIVITY ESURF 320 .16 SURFACE ALBEDO ALBDO 330 9.8 GRAV ACCEL GRAV 340 .41 VON KARMAN CONSTANT K 350 1 MEASUREMENT HEIGHT Zl 360 60 MINUTES BETWEEN DATA SAMPLES DT 380 161 NUMBER OF NODES IN THE PROFILE NODES 390 .005 80 CM PROFILE WITH .5 CM NODES DZ 400 2 TIME VARIABLE T 410 • 3.14159 3.14159 ' PI PI 420 E T Uf UP INITIAL AND BOUNDARY CONDITIONS OF THE SOIL PROFILE 500 REM SSET 510 REM THERMAL PROPERTIES OF THE PROFILE PRINT"FROM THE SURFACE DOWN TO THE 2 CM DEPTH " 515 INPUT"ENTER THE VOLUMETRIC HEAT CAPACITY 1.6 (MJ/M3 K) ";HC 518 INPUT"ENTER THE THERMAL CONDUCTIVITY 1.32 (W/M K) "jTC 520 FOR Z-NODES TO (NODES-4) STEP -1 522 COND(Z) - TC 'J/S M C 524 RHOCP(Z)- HC *lE+06 'J/M3 C 526 530 NEXT Z 535 PRINT"FROM THE 2 CM SOIL DEPTH DOWN TO THE 5 CM DEPTH " •;HC 538 INPUT"ENTER THE VOLUMETRIC HEAT CAPACITY 1.6 (MJ/M3 K) 540 INPUT"ENTER THE THERMAL CONDUCTIVITY 1.32 (W/M K) ";TC 542 FOR Z-(NODES-5) TO (NODES-10) STEP -1 544 COND(Z) - TC 'J/S M C 546 RHOCP(Z)- HC *lE+06 «J/M3 C 548 NEXT Z 571 PRINT"FROM THE 5 CM SOIL DEPTH DOWN TO THE 25 CM DEPTH jHC 575 INPUT"ENTER THE VOLUMETRIC HEAT CAPACITY 1.6 (MJ/M3 K) 576 INPUT"ENTER THE THERMAL CONDUCTIVITY 1.32 (W/M K) "jTC 580 FOR Z-(NODES-ll) TO (NODES-50) STEP -1 590 COND(Z) TC *J/S M C 600 RHOCP(Z)" HC *lE+06 'J/M3 C 610 NEXT Z ;HC 615 PRINT"FROM THE 25 CM SOIL DEPTH DOWN TO 80 CM " 620 INPUT"ENTER THE VOLUMETRIC HEAT CAPACITY 1.6 (MJ/M3 K) 630 INPUT"ENTER THE THERMAL CONDUCTIVITY 1.32 (W/M K) "»TC 780 FOR Z- (NODES-51) TO 0 STEP -1 790 COND(Z) TC 'J/S M C 800 RHOCP(Z)- HC*lE+06 'J/M3 C 810 NEXT Z 820 DD • SQR(COND(15)/RHOCP(15)*864001/PI) 830 840 PRINT'DAMPING DEPTH (CM) « "jDD Z • -DZ  200 230 240 250  REM  -IZ5-  850 POR I • NODES TO 1 STEP -1 860 Z - Z • DZ 870 DEPTH(I) - Z 880 TEMP(1,I)-TBAR • +10*EXP(-DEPTH(I)/DD)*SIN(2*PI - DEPTH(I)/DD) 885 TEMP(2,1)-TBAR 890 REM IF I<200 THEN PRINT USING"##.#### "jI,DEPTH(I)jTEMP(1,1),COND(I),RHC 900 IF I > NODES THEN LPRINT USING"If.f ";TEMP(l,I)j 910 NEXT I 911 LPRINT 912 PRINT 914 LPRINT" " j 915 LPRINT USING"\ \";"TIME">"TA"j"U"j"RI"> 917 LPRINT OSING"\ \";"RS";"RN";"G";"H";"LE"j 920 LPRINT USING"\ \";"T0";"T5";"T10">"T15";"T20"j"T25";"T30";"T35";"T40" 925 LPRINT 1000 REM*********************MAIN PROGRAM STARTS HERE**********"""** 1010 WHILE NOT EOF(l) 1020 INPUT!1,DAY,TIME,T0,T1,T10,T15,T30,T50,TA,RS,U 1030 PRINT" T« "jTIME" SOLAR - "RS» 1050 PRINT"TAIR - "TA" WIND • "U 1060 GOSUB 2500 'CALCULATE VARIABLES THAT DEPEND ON RS,TA,U 1120 WHILE ABS(GO-GA) > 1 1140 IF GA>G0 THEN GA=GA-(ABS(GA-G0))/2 'ESTIMATE A NEW 1145 IF G0>GA THEN GA«=GA+ (ABS (GA-G0)) /2 'VALUE OF GA 1160 GOSUB 3000 'ITERATE FOR H AND TS 1180 GOSUB 4000 'ITERATE FOR GO 1182 PRINT "RN,LE,GA,H1,TS"; 1183 PRINT USING"!!!!.!! ";RN,LE,GA,Hl,TS 1185 PRINT "RN,LE,G0,H1,TS"j 1187 PRINT USING"####.#t RN,LE,GO,HI,TEMP (T-l,NODES) 1192 PRINT "RI ,PHIH,PHIM,USTAR"j 1193 PRINT USING"!!I.lttt ";RI,PHIH,PHIM,USTAR 1194 REM PRINT TWOGASAGO,LASTGA,GA 1195 IF TWOGASAGO - GA THEN GA - (GA+LASTGA)/2 1196 TWOGASAGO - LASTGA 1197 LASTGA - GA 1199 REM PRINT TWOGASAGO,LASTGA,GA 1200 PRINT 1205 WEND 1210 PRINT 1220 REM CHANGE FROM ONE HOUR TO NEXT WHEN RN-H+LE+G**** 1230 FOR Z - 1 TO NODES 1240 TEMP (1,Z) - TEMP(2,2) 1250 NEXT Z 1 3 0 0 REM PRINT DATA TO FILE AND PRINTER 1320 LPRINT USING"!!.t ">TIME,TA,U,RI} 1330 LPRINT USING"!!!!.! "jRS,RN,H1,GA,LE; 1360 REM FOR Z - NODES TO 143 STEP -1  -IZ6-  LPRINT USING"it.it ";TEMP(T,Z); 1370 REM NEXT Z 1380 REM LPRINT:LPRINT" 1385 REM PRINTi 2,USING"ft"»TIME; 1387 FOR Z « NODES TO 61 STEP -10 1390 LPRINT USING" ##.I";TEMP(T,Z); 1400 PRINT42,USING" it"jDEPTH(Z)*100j 1405 PRINTi2,USING" f(.#"jTEMP(T,Z); 1407 NEXT Z 1410 LPRINT:PRINTi2,"" 1500 2400 WEND CLOSE : END 2440 2500 REM SUBROUTINE TO SET VARIABLES****************************** - 1.00055-(TA*(-.0001516)) 'DENSITY OF AIR RHOA 2520 « 1300.05 - (TA*4.07) 'SPECIFIC HEAT OF AIR CPA 2540 RHOCPA • RHOA*CPA 'HEAT CAPACITY OF AIR 2560 EAIR 2580 - 1-(.261*EXP((TA*2)*-.000777)) 'SKY EMISSIVITY RLD 2620 • EAIR*SIGMA*(273.2+TA)"4 'DOWNWARD LONGWAVE FLUX GA 2640 GO • 0.15*RS-40 'ESTIMATED GO 2650 2660 REM LE - GA-10 LE 2680 - .015*(TA+3)*RS HI 2740 H2 » ALF*RS 'NOTE ALF RANGES FROM 0 TO 1 2750 TS 2760 - 100 'RESET HI RETURN » 20 2880 'RESET H2 3000 • TA -5 "RESET TS 3010 RETURN 'TO MAIN PROGRAM 3020 REM SENSIBLE HEAT TRANSFER SUBROUTINE a************************** 3030 WHILE H2 < HI 3040 GOSUB 3500> 'ESTIMATE WHILE HI TS - H2 3050 IFGOSUB HI H2 >3500 H2 THEN 'ESTIMATETS H2+ 1 3060 WEND IF HI < H2 THEN TS = TS - .25 3070 WEND 3080 WHILE H2 < HI 3100 GOSUB 3500 'ESTIMATE H2 3110 IF HI > H2 THEN TS - TS + .05 3120 WEND 3140 WHILE H2 > HI 3160 3180 GOSUB 3500 'ESTIMATE H2 3200 IF HKH2 THEN TS • TS - .01 3210 WEND 3400 RETURN ' TO MAIN PROGRAM  -111-  3500 REM SUBROUTINE TO ESTIMATE H2************************************ 3520 TK - (TS+TA)/2 +273.15 3540 RLD - ESURF"EAIR*SIGMA*(273.15+TA)"4 3560 RLU • ESURF*SIGMA*(273.15+TS)~4 3580 HI •> RS*(1-ALBDO)+RLD-RLU-LE-GA 3600 RN - RS*(1-ALBDO)+RLD-RLU 3610 RI - -(GRAV/TK)*(TS-TA)*Z1/U~2 • REM TO MULT LOG(Z/Z0) TO Z 3612 IF (16*RI)>11 THEN RI - .0625 3615 X - <l-16*RIp.25 3620 PHIH - 2* LOG((l+X~2)/2) 3622 IF RI > 0 THEN PHIH - 4.7* ( RI/(1-4.7*RI)) 'PHIH • 4.7 ZETA 3625 PHIM - 2*LOG((l+X)/2) + LOG((1+X~2)/2) - 2* ATN(X) + PI/2 3627 IF RI > 0! THEN PHIM - PHIH 3630 USTAR - (K*U) / ( LOG(Z1/Z0) - PHIM) 3635 HTRANS - K*USTAR/(LOG(Z1/Z0) - PHIH) 3700 H2 - RHOCPA*HTRANS*(TS-TA) 3730 RETURN * TO SUBROUTINE ITERATION 4000 REM SUBROUTINE SOILTEMP CALCULATES GO AND PROFILES OF TEMPERATURE 4020 TEMP(T.NODES)- TS 4060 DZZ - 2* (DZ~2) 4080 DTT - 3600 'SECONDS PER HOUR 4100 A(l) * 0 4120 B(l) - 1 4140 C(l) - 0 4160 D(l) « TEMP(T,1) 4180 E(l) - A(l)/B(l) 4200 F ( l ) - D(l)/B(l) 4220 A(NODES)- 0 4240 B (NODES)- 1 4260 C(NODES)- 0 4280 D(NODES)- TEMP(T.NODES) 4300 E(NODES)- A(NODES)/B(NODES) 4320 F(NODES)- D(NODES)/B(NODES) 4340 FOR Z - 2 TO NODES-1 4360 C(Z) - COND(Z-l) *DTT / DZZ 4380 A(Z) • COND(Z) *DTT / DZZ 4400 B(Z) - RHOCP(Z) + A(Z) + C(Z) 4420 M - A(Z) * TEMP(T-1,Z+1) 4440 N » (RHOCP(Z) - A(Z) - C(Z)) * TEMP(T-1,Z) 4460 O - C(Z) * TEMP(T-1,Z-1) 4480 D(Z) - M + N + O 4500 NEXT Z 4520 FOR Z - 2 TO NODES-1 4540 BCE - B(Z) - ( C(Z)*E(Z-1) ) 4560 E(Z)- A(Z) / BCE 4580 P(Z)- (D(Z) + (C(Z) » F(Z-l))) / BCE 4600 NEXT Z 4620 FOR Z - NODES-1 TO 2 STEP -1 4640 TEMP(T.Z) •( E(Z) • TEMP(T,Z+1) ) + F(Z) 4660 NEXT Z 4740 GO - COND(NODES)* (TEMP(T.NODES) - TEMP(T,NODES-1)) / DZ 4750 GI - COND(NODES-1)*(TEMP(T,NODES-1) - TEMP (T.NODES-2)) / DZ 4760 REM PRINT"G0 - "jGO" GAIR » "»GA" TS - "»TS" T(T,TOP) • "% TEMP(T,1) 4800 RETURN 'TO MAIN PROGRAM  

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-0095699/manifest

Comment

Related Items